11 #ifndef Pythia8_ColourReconnection_H 12 #define Pythia8_ColourReconnection_H 14 #include "Pythia8/Basics.h" 15 #include "Pythia8/BeamParticle.h" 16 #include "Pythia8/Event.h" 17 #include "Pythia8/FragmentationFlavZpT.h" 18 #include "Pythia8/Info.h" 19 #include "Pythia8/ParticleData.h" 20 #include "Pythia8/StringFragmentation.h" 21 #include "Pythia8/PartonDistributions.h" 22 #include "Pythia8/PartonSystems.h" 23 #include "Pythia8/PythiaStdlib.h" 24 #include "Pythia8/Settings.h" 25 #include "Pythia8/StringLength.h" 26 #include "Pythia8/StringInteractions.h" 42 int colReconnectionIn = 0,
bool isJunIn =
false,
bool isAntiJunIn =
false,
43 bool isActiveIn =
true,
bool isRealIn =
false) :
col(colIn), iCol(iColIn),
44 iAcol(iAcolIn), colReconnection(colReconnectionIn), isJun(isJunIn),
45 isAntiJun(isAntiJunIn),isActive(isActiveIn), isReal(isRealIn)
46 {iColLeg = 0; iAcolLeg = 0; printed =
false; p1p2 = 0.;}
48 double mDip(
Event & event) {
49 if (isJun || isAntiJun)
return 1E9;
50 else return m(event[iCol].p(),event[iAcol].p());
54 int col, iCol, iAcol, iColLeg, iAcolLeg, colReconnection;
55 bool isJun, isAntiJun, isActive, isReal, printed;
56 weak_ptr<ColourDipole> leftDip{}, rightDip{};
57 vector<weak_ptr<ColourDipole> > colDips, acolDips;
67 inline bool operator<(
const ColourDipolePtr& d1,
const ColourDipolePtr& d2) {
68 return ( d1 && d2 ? d1->index < d2->index : !d1 );}
82 dipsOrig() {
for(
int i = 0;i < 3;++i) {
83 dips[i] = ju.dips[i]; dipsOrig[i] = ju.dipsOrig[i];}
86 this->Junction::operator=(ju);
87 for(
int i = 0;i < 3;++i) {
88 dips[i] = ju.dips[i]; dipsOrig[i] = ju.dipsOrig[i];
93 ColourDipolePtr getColDip(
int i) {
return dips[i];}
94 void setColDip(
int i, ColourDipolePtr dip) {dips[i] = dip;}
95 ColourDipolePtr dips[3];
96 ColourDipolePtr dipsOrig[3];
112 ColourDipolePtr dip3In = 0, ColourDipolePtr dip4In = 0,
int modeIn = 0,
113 double lambdaDiffIn = 0) {
114 dips.push_back(dip1In); dips.push_back(dip2In);
115 dips.push_back(dip3In); dips.push_back(dip4In);
116 mode = modeIn; lambdaDiff = lambdaDiffIn;
120 cout <<
"mode: " << mode <<
" " <<
"lambdaDiff: " << lambdaDiff << endl;
121 for (
int i = 0;i < int(dips.size()) && dips[i] != 0;++i) {
122 cout <<
" "; dips[i]->list(); }
125 vector<ColourDipolePtr> dips;
143 vector<vector<ColourDipolePtr> > dips;
144 vector<bool> colEndIncluded, acolEndIncluded;
145 vector<ColourDipolePtr> activeDips;
151 void listActiveDips();
168 singleReconOnly(), lowerLambdaOnly(), allowDiqJunCR(), nSys(),
169 nReconCols(), swap1(), swap2(), reconnectMode(), flipMode(),
170 timeDilationMode(), eCM(), sCM(), pT0(), pT20Rec(), pT0Ref(), ecmRef(),
171 ecmPow(), reconnectRange(), m0(), mPseudo(), m2Lambda(), fracGluon(),
172 dLambdaCut(), timeDilationPar(), timeDilationParGeV(), tfrag(), blowR(),
173 blowT(), rHadron(), kI(), dipMaxDist(), nColMove() {}
180 {beamAPtr = beamAPtrIn; beamBPtr = beamBPtrIn;}
183 bool next(
Event & event,
int oldSize);
188 static const double MINIMUMGAIN, MINIMUMGAINJUN, TINYP1P2;
189 static const int MAXRECONNECTIONS;
192 bool allowJunctions, sameNeighbourCol, singleReconOnly, lowerLambdaOnly;
194 int nSys, nReconCols, swap1, swap2, reconnectMode, flipMode,
196 double eCM, sCM, pT0, pT20Rec, pT0Ref, ecmRef, ecmPow, reconnectRange,
197 m0, mPseudo, m2Lambda, fracGluon, dLambdaCut, timeDilationPar,
198 timeDilationParGeV, tfrag, blowR, blowT, rHadron, kI;
202 vector<ColourDipolePtr> dipoles, usedDipoles;
208 void addDipole(
int colIn = 0,
int iColIn = 0,
int iAcolIn = 0,
209 int colReconnectionIn = 0,
bool isJunIn =
false,
bool isAntiJunIn =
false,
210 bool isActiveIn =
true,
bool isRealIn =
false) {
211 dipoles.push_back(make_shared<ColourDipole>(colIn, iColIn, iAcolIn,
212 colReconnectionIn, isJunIn, isAntiJunIn, isActiveIn, isRealIn));
213 dipoles.back()->index = ++dipoleIndex;
218 dipoles.push_back(make_shared<ColourDipole>(dipole));
219 dipoles.back()->index = ++dipoleIndex;
223 vector<ColourJunction> junctions;
224 vector<ColourParticle> particles;
225 vector<TrialReconnection> junTrials, dipTrials;
226 vector<vector<int> > iColJun;
227 map<int,double> formationTimes;
236 bool nextNew(
Event & event,
int oldSize);
239 void swapDipoles(ColourDipolePtr dip1, ColourDipolePtr dip2,
243 void setupDipoles(
Event& event,
int iFirst = 0);
246 void makePseudoParticle(ColourDipolePtr dip,
int status,
247 bool setupDone =
false);
251 bool getJunctionIndices(ColourDipolePtr dip,
int &iJun,
int &i0,
int &i1,
252 int &i2,
int &junLeg0,
int &junLeg1,
int &junLeg2);
255 void makeAllPseudoParticles(
Event & event,
int iFirst = 0);
258 void updateEvent(
Event& event,
int iFirst = 0);
260 double calculateStringLength( ColourDipolePtr dip,
261 vector<ColourDipolePtr> & dips);
264 double calculateStringLength(
int i,
int j);
268 double calculateJunctionLength(
int i,
int j,
int k);
273 double calculateDoubleJunctionLength(
int i,
int j,
int k,
int l);
278 bool findJunctionParticles(
int iJun, vector<int>& iParticles,
279 vector<bool> &usedJuns,
int &nJuns, vector<ColourDipolePtr> &dips);
282 void singleReconnection( ColourDipolePtr dip1, ColourDipolePtr dip2);
285 void singleJunction(ColourDipolePtr dip1, ColourDipolePtr dip2);
288 void singleJunction(ColourDipolePtr dip1, ColourDipolePtr dip2,
289 ColourDipolePtr dip3);
292 void listChain(ColourDipolePtr dip);
295 void listAllChains();
298 void listDipoles(
bool onlyActive =
false,
bool onlyReal =
false);
301 void listParticles();
304 void listJunctions();
310 void checkRealDipoles(
Event& event,
int iFirst);
313 double mDip(ColourDipolePtr dip);
318 Vec4 getVProd(ColourDipolePtr dip,
bool anti)
const;
322 Vec4 getVProd(
int iJun, ColourDipolePtr dip,
bool anti)
const;
326 bool checkDist(ColourDipolePtr dip1, ColourDipolePtr dip2);
330 bool findAntiNeighbour(ColourDipolePtr& dip);
334 bool findColNeighbour(ColourDipolePtr& dip);
337 bool checkJunctionTrials();
352 double getLambdaDiff(ColourDipolePtr dip1, ColourDipolePtr dip2,
353 ColourDipolePtr dip3, ColourDipolePtr dip4,
int mode);
356 double getLambdaDiff(ColourDipolePtr dip1, ColourDipolePtr dip2);
359 void updateDipoleTrials();
362 void updateJunctionTrials();
365 bool checkTimeDilation(ColourDipolePtr dip1 = 0, ColourDipolePtr dip2 = 0,
366 ColourDipolePtr dip3 = 0, ColourDipolePtr dip4 = 0);
369 bool checkTimeDilation(
Vec4 p1,
Vec4 p2,
double t1,
double t2);
372 Vec4 getDipoleMomentum(ColourDipolePtr dip);
375 void addJunctionIndices(
int iSinglePar, vector<int> &iPar,
376 vector<int> &usedJuncs);
379 void setupFormationTimes(
Event & event);
382 double getJunctionMass(
Event & event,
int col);
385 void addJunctionIndices(
Event & event,
int iSinglePar,
386 vector<int> &iPar, vector<int> &usedJuncs);
389 bool reconnectMPIs(
Event& event,
int oldSize);
394 vector<int> iReduceCol, iExpandCol;
398 vector<double> lambdaijMove;
401 double lambda12Move(
int i,
int j) {
402 int iAC = iReduceCol[i];
int jAC = iReduceCol[j];
403 return lambdaijMove[nColMove * min( iAC, jAC) + max( iAC, jAC)];
407 double lambda123Move(
int i,
int j,
int k) {
408 int iAC = iReduceCol[i];
int jAC = iReduceCol[j];
int kAC = iReduceCol[k];
409 return lambdaijMove[nColMove * min( iAC, jAC) + max( iAC, jAC)]
410 + lambdaijMove[nColMove * min( iAC, kAC) + max( iAC, kAC)]
411 - lambdaijMove[nColMove * min( jAC, kAC) + max( jAC, kAC)];
415 bool reconnectMove(
Event& event,
int oldSize);
418 bool reconnectTypeCommon(
Event& event,
int oldSize);
421 map<double,pair<int,int> > reconnectTypeI(
Event& event,
422 vector<vector<ColourDipole> > &dips,
Vec4 decays[2]);
426 map<double,pair<int,int> > reconnectTypeII(
Event& event,
427 vector<vector<ColourDipole> > &dips,
Vec4 decays[2]);
431 double determinant3(vector<vector< double> >& vec);
Definition: ColourReconnection.h:36
TrialReconnection class.
Definition: ColourReconnection.h:107
The Event class holds all info on the generated event.
Definition: Event.h:453
Definition: BeamParticle.h:133
void reassignBeamPtrs(BeamParticle *beamAPtrIn, BeamParticle *beamBPtrIn)
New beams possible for handling of hard diffraction.
Definition: ColourReconnection.h:179
Definition: StringFragmentation.h:119
int col
Members.
Definition: ColourReconnection.h:54
Definition: StringInteractions.h:78
ColourDipole(int colIn=0, int iColIn=0, int iAcolIn=0, int colReconnectionIn=0, bool isJunIn=false, bool isAntiJunIn=false, bool isActiveIn=true, bool isRealIn=false)
Constructor.
Definition: ColourReconnection.h:41
StringLength class. It is used to calculate the lambda measure.
Definition: StringLength.h:23
ColourParticle class.
Definition: ColourReconnection.h:137
Definition: ColourReconnection.h:75
The ColourReconnection class handles the colour reconnection.
Definition: ColourReconnection.h:162
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
double m(const Vec4 &v1)
Invariant mass and its square.
Definition: Basics.cc:576
void list()
Printing function, mainly intended for debugging.
Definition: ColourReconnection.cc:38
bool operator<(const ColourDipolePtr &d1, const ColourDipolePtr &d2)
Comparison operator by index for two dipole pointers.
Definition: ColourReconnection.h:67
ColourReconnection()
Constructor.
Definition: ColourReconnection.h:167