PYTHIA  8.315
StringFragmentation.h
1 // StringFragmentation.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2025 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 
6 // This file contains the classes for string fragmentation.
7 // StringEnd: keeps track of the fragmentation step.
8 // StringFragmentation: is the top-level class.
9 
10 #ifndef Pythia8_StringFragmentation_H
11 #define Pythia8_StringFragmentation_H
12 
13 #include "Pythia8/FragmentationModel.h"
14 
15 namespace Pythia8 {
16 
17 //==========================================================================
18 
19 // The StringEnd class contains the information related to
20 // one of the current endpoints of the string system.
21 // Only to be used inside StringFragmentation, so no private members.
22 
23 class StringEnd {
24 
25 public:
26 
27  // Constructor.
28  StringEnd() : particleDataPtr(), flavSelPtr(), pTSelPtr(), zSelPtr(),
29  fromPos(), thermalModel(), mT2suppression(), iEnd(), iMax(), idHad(),
30  iPosOld(), iNegOld(), iPosNew(), iNegNew(), hadSoFar(), colOld(), colNew(),
31  pxOld(), pyOld(), pxNew(), pyNew(), pxHad(), pyHad(), mHad(), mT2Had(),
32  zHad(), GammaOld(), GammaNew(), xPosOld(), xPosNew(), xPosHad(), xNegOld(),
33  xNegNew(), xNegHad(), aLund(), bLund(), iPosOldPrev(), iNegOldPrev(),
34  colOldPrev(), pxOldPrev(), pyOldPrev(), GammaOldPrev(), xPosOldPrev(),
35  xNegOldPrev(), mVecRatio(1.), tinyEq(), pT2tiny() {}
36 
37  // Save pointers.
38  void init( ParticleData* particleDataPtrIn, StringFlav* flavSelPtrIn,
39  StringPT* pTSelPtrIn, StringZ* zSelPtrIn, Settings& settings) {
40  particleDataPtr = particleDataPtrIn; flavSelPtr = flavSelPtrIn;
42  pTSelPtr = pTSelPtrIn; zSelPtr = zSelPtrIn;
43  bLund = zSelPtr->bAreaLund(); aLund = zSelPtr->aAreaLund();
44  thermalModel = settings.flag("StringPT:thermalModel");
45  mT2suppression = settings.flag("StringPT:mT2suppression");
46  closePacking = settings.flag("ClosePacking:doClosePacking"); }
47 
48  // Set up initial endpoint values from input.
49  void setUp(bool fromPosIn, int iEndIn, int idOldIn, int iMaxIn,
50  double pxIn, double pyIn, double GammaIn, double xPosIn,
51  double xNegIn, int colIn, double mVecRatioIn);
52 
53  // Fragment off one hadron from the string system, in flavour and pT.
54  void newHadron(double kappaModifier, bool forbidPopcornNow = false,
55  double strangeJunc = 0., double probQQmod = 1.);
56 
57  // Fragment off one hadron from the string system, in momentum space,
58  // by taking steps either from positive or from negative end.
59  Vec4 kinematicsHadron(StringSystem& system, StringVertex& newVertex,
60  double zHadIn);
61 
62  // Generate momentum for some possible next hadron, based on mean values
63  // to get an estimate for rapidity and pT.
64  Vec4 kinematicsHadronTmp(StringSystem system, Vec4 pRem, double phi,
65  double mult);
66 
67  // Update string end information after a hadron has been removed.
68  void update();
69  void storePrev();
70  void updateToPrev();
71 
72  // Constants: could only be changed in the code itself.
73  static const double TINY, PT2SAME, MEANMMIN, MEANM, MEANPT;
74 
75  // Pointer to the particle data table.
77 
78  // Pointers to classes for flavour, pT and z generation.
80  StringPT* pTSelPtr;
81  StringZ* zSelPtr;
82 
83  // Local copy of flavSelPtr for modified flavour selection.
85 
86  // Data members.
87  bool fromPos, thermalModel, mT2suppression, closePacking;
88  int iEnd, iMax, idHad, iPosOld, iNegOld, iPosNew, iNegNew, hadSoFar,
89  colOld, colNew;
90  double pxOld, pyOld, pxNew, pyNew, pxHad, pyHad, mHad, mT2Had, zHad,
91  GammaOld, GammaNew, xPosOld, xPosNew, xPosHad, xNegOld, xNegNew,
92  xNegHad, aLund, bLund;
93  int iPosOldPrev, iNegOldPrev, colOldPrev;
94  double pxOldPrev, pyOldPrev, GammaOldPrev, xPosOldPrev, xNegOldPrev,
95  mVecRatio, tinyEq, pT2tiny;
96  FlavContainer flavOld, flavNew, flavOldPrev;
97  Vec4 pHad, pSoFar;
98 
99 };
100 
101 //==========================================================================
102 
103 // The StringFragmentation class contains the top-level routines
104 // to fragment a colour singlet partonic system.
105 
107 
108 public:
109 
110  // Constructor.
112  FragmentationModel(), flavRopePtr(), closePacking(),
113  setVertices(), constantTau(), smearOn(), traceColours(false),
114  hadronVertex(), stopMass(), stopNewFlav(), stopSmear(),
115  pNormJunction(), pMaxJunction(), eBothLeftJunction(),
116  eMaxLeftJunction(), eMinLeftJunction(), mJoin(), bLund(),
117  closePackingFluxRatio(1.), closePackingPT20(1.), pT20(),
118  xySmear(), maxSmear(), maxTau(), kappaVtx(), mc(), mb(),
119  hasJunction(), isClosed(), iPos(), iNeg(), nExtraJoin(),
120  w2Rem(), stopMassNow(), mVecRatio(1.), closedM2max(),
121  idDiquark(), legMin(), legMid() {}
122 
123  // Initialize and save pointers.
124  bool init(StringFlav* flavSelPtrIn = nullptr, StringPT* pTSelPtrIn = nullptr,
125  StringZ* zSelPtrIn = nullptr, FragModPtr fragModPtrIn = nullptr) override;
126 
127  // Do the fragmentation: driver routine.
128  bool fragment(int iSub, ColConfig& colConfig, Event& event,
129  bool isDiff = false, bool systemRecoil = true) override;
130 
131  // Local copy of flavSelPtr for modified flavour selection.
133 
134  // Find the boost matrix to the rest frame of a junction.
135  Vec4 junctionRestFrame(const Vec4& p0, const Vec4& p1, const Vec4& p2,
136  const bool angleCheck = true) const;
137 
138  // Set the vector mass ratio.
139  void setMVecRatio(double mVecRatioIn) {mVecRatio = mVecRatioIn;}
140 
141 private:
142 
143  // Constants: could only be changed in the code itself.
144  static const int NTRYFLAV, NTRYJOIN, NSTOPMASS,
145  NTRYJNMATCH, NTRYJRFEQ, NTRYSMEAR, MAXVETOFINTWO;
146  static const double FACSTOPMASS, CLOSEDM2MAX, CLOSEDM2FRAC, EXPMAX,
147  MATCHPOSNEG, M2MINJRF, EMINJRF, EEXTRAJNMATCH,
148  MDIQUARKMIN, CONVJRFEQ, CHECKPOS;
149 
150  // Pointer to flavour-composition-changing ropes.
151  FragModPtr flavRopePtr;
152 
153  // Initialization data, read from Settings.
154  bool closePacking, setVertices, constantTau, smearOn,
155  traceColours, hardRemn, doStrangeJunc;
156  int hadronVertex;
157  double stopMass, stopNewFlav, stopSmear, pNormJunction, pMaxJunction,
158  eBothLeftJunction, eMaxLeftJunction, eMinLeftJunction,
159  mJoin, bLund, closePackingFluxRatio, closePackingPT20,
160  qqSupPar, qqSupAnti, pT20, xySmear, maxSmear, maxTau,
161  kappaVtx, mc, mb, dampPopcorn, aRemn, bRemn, strangeJuncParm;
162 
163  // Data members.
164  bool hasJunction, isClosed;
165  int iPos, iNeg, nExtraJoin;
166  double w2Rem, stopMassNow, kappaModifier, probQQmod, mVecRatio, closedM2max;
167  Vec4 pSum, pRem, pJunctionHadrons;
168 
169  // UserHooks flags.
170  bool doChangeFragPar = false, doVetoFrag = false;
171 
172  // List of partons in string system.
173  vector<int> iParton, iPartonMinLeg, iPartonMidLeg, iPartonMax;
174 
175  // Vertex information from the fragmentation process.
176  vector<StringVertex> stringVertices, legMinVertices, legMidVertices;
177  StringVertex newVertex;
178 
179  // Variables used for JRF finding.
180  vector<Vec4> listJRF;
181  vector<double> weightJRF;
182  int iLeg[3], idLeg[3], legEnd[3];
183  double weightSum, pSumJRF, m2Leg[3];
184  Vec4 pLeg[3];
185  bool lastJRF, endpoint[3];
186 
187  // Boost from/to rest frame of a junction to original frame.
188  RotBstMatrix MfromJRF, MtoJRF;
189 
190  // Information on diquark created at the junction.
191  int idDiquark;
192 
193  // Fictitious opposing partons in JRF: string ends for vertex location.
194  Vec4 pMinEnd, pMidEnd;
195 
196  // Temporary event record for the produced particles.
197  Event hadrons;
198 
199  // Information on the system of string regions.
200  StringSystem system, systemMin, systemMid;
201 
202  // Information on the two current endpoints of the fragmenting
203  // system, with backup copies if there is a veto.
204  StringEnd posEnd, negEnd, posEndSave, negEndSave;
205 
206  // Find region where to put first string break for closed gluon loop.
207  vector<int> findFirstRegion(int iSub, const ColConfig& colConfig,
208  const Event& event) const;
209 
210  // Set flavours and momentum position for initial string endpoints.
211  void setStartEnds(int idPos, int idNeg, const StringSystem& systemNow,
212  int legNow = 3);
213 
214  // Check remaining energy-momentum whether it is OK to continue.
215  bool energyUsedUp(bool fromPos);
216 
217  // Produce the final two partons to complete the system.
218  bool finalTwo(bool fromPos, const Event& event, bool usedPosJun,
219  bool usedNegJun);
220 
221  // Final region information.
222  Vec4 pPosFinalReg, pNegFinalReg, eXFinalReg, eYFinalReg;
223 
224  // Set hadron production points in space-time picture.
225  bool setHadronVertices(Event& event);
226 
227  // Construct a special joining region for the final two hadrons.
228  StringRegion finalRegion();
229 
230  // Store the hadrons in the normal event record, ordered from one end.
231  void store(Event& event);
232 
233  // Fragment off two of the string legs in to a junction.
234  bool fragmentToJunction(Event& event,
235  vector< vector< pair<double,double> > >& rapPairs);
236 
237  // Initially considered legs from the junction.
238  int legMin, legMid;
239 
240  // Functions used in JRF iterative procedure.
241  bool collinearPair(Event& event);
242  bool perturbedJRF(Event& event);
243  int updateLegs(Event& event, Vec4 vJunIn, bool juncCoM = false);
244  double updateWeights(double pSmall, Vec4 vJunIn);
245  void nextParton(Event& event, int leg);
246 
247  // Join extra nearby partons when stuck.
248  int extraJoin(double facExtra, Event& event);
249 
250  // Get the number of nearby strings given the energies.
251  void kappaEffModifier(StringSystem& systemNow,
252  StringEnd end, bool fromPos, vector<int> partonList,
253  vector< vector< pair<double,double> > >& rapPairs,
254  double mRem, Event& event);
255 
256  double yMax(Particle pIn, double mTiny) {
257  double temp = log( ( pIn.e() + abs(pIn.pz()) ) / max( mTiny, pIn.mT()) );
258  return (pIn.pz() > 0) ? temp : -temp; }
259 
260 };
261 
262 //==========================================================================
263 
264 } // end namespace Pythia8
265 
266 #endif // Pythia8_StringFragmentation_H
Definition: FragmentationSystems.h:225
bool flag(string keyIn)
Give back current value, with check that key exists.
Definition: Settings.cc:1602
virtual double aAreaLund()
a and b fragmentation parameters needed in some operations.
Definition: FragmentationFlavZpT.h:303
bool fromPos
Data members.
Definition: StringFragmentation.h:87
The Event class holds all info on the generated event.
Definition: Event.h:408
Vec4 kinematicsHadron(StringSystem &system, StringVertex &newVertex, double zHadIn)
Definition: StringFragmentation.cc:132
void newHadron(double kappaModifier, bool forbidPopcornNow=false, double strangeJunc=0., double probQQmod=1.)
Fragment off one hadron from the string system, in flavour and pT.
Definition: StringFragmentation.cc:67
void setUp(bool fromPosIn, int iEndIn, int idOldIn, int iMaxIn, double pxIn, double pyIn, double GammaIn, double xPosIn, double xNegIn, int colIn, double mVecRatioIn)
Set up initial endpoint values from input.
Definition: StringFragmentation.cc:39
The StringPT class is used to select select transverse momenta.
Definition: FragmentationFlavZpT.h:333
static const double PT2SAME
Assume two (eX, eY) regions are related if pT2 differs by less.
Definition: StringFragmentation.h:73
StringFlav * flavSelPtr
Pointers to classes for flavour, pT and z generation.
Definition: StringFragmentation.h:79
ParticleData * particleDataPtr
Pointer to the particle data table.
Definition: StringFragmentation.h:76
Definition: StringFragmentation.h:106
Definition: FragmentationSystems.h:186
StringFragmentation()
Constructor.
Definition: StringFragmentation.h:111
Definition: FragmentationSystems.h:136
StringEnd()
Constructor.
Definition: StringFragmentation.h:28
The StringZ class is used to sample the fragmentation function f(z).
Definition: FragmentationFlavZpT.h:268
void updateToPrev()
Update string end information to previous string break.
Definition: StringFragmentation.cc:572
Definition: StringFragmentation.h:23
Vec4 kinematicsHadronTmp(StringSystem system, Vec4 pRem, double phi, double mult)
Definition: StringFragmentation.cc:321
void storePrev()
Store old string end information.
Definition: StringFragmentation.cc:554
Definition: Event.h:32
StringFlav flavSelNow
Local copy of flavSelPtr for modified flavour selection.
Definition: StringFragmentation.h:84
The ColConfig class describes the colour configuration of the whole event.
Definition: FragmentationSystems.h:60
void update()
Update string end information after a hadron has been removed.
Definition: StringFragmentation.cc:536
void setMVecRatio(double mVecRatioIn)
Set the vector mass ratio.
Definition: StringFragmentation.h:139
Definition: Basics.h:254
Definition: FragmentationFlavZpT.h:33
FragmentationModel is the base class for handling fragmentation algorithms.
Definition: FragmentationModel.h:28
double phi(const Vec4 &v1, const Vec4 &v2)
phi is azimuthal angle between v1 and v2 around z axis.
Definition: Basics.cc:704
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
The StringFlav class is used to select quark and hadron flavours.
Definition: FragmentationFlavZpT.h:76
StringFlav flavSelNow
Local copy of flavSelPtr for modified flavour selection.
Definition: StringFragmentation.h:132
static const double MEANMMIN
Definition: StringFragmentation.h:73
void init(ParticleData *particleDataPtrIn, StringFlav *flavSelPtrIn, StringPT *pTSelPtrIn, StringZ *zSelPtrIn, Settings &settings)
Save pointers.
Definition: StringFragmentation.h:38
This class holds a map of all ParticleDataEntries.
Definition: ParticleData.h:422
static const double TINY
Constants: could only be changed in the code itself.
Definition: StringFragmentation.h:73
Definition: Basics.h:32
Definition: Settings.h:196