include/Pythia8/StringFragmentation.h

0001 // StringFragmentation.h is a part of the PYTHIA event generator.
0002 // Copyright (C) 2025 Torbjorn Sjostrand.
0003 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
0004 // Please respect the MCnet Guidelines, see GUIDELINES for details.
0005
0006 // This file contains the classes for string fragmentation.
0007 // StringEnd: keeps track of the fragmentation step.
0008 // StringFragmentation: is the top-level class.
0009
0010 #ifndef Pythia8_StringFragmentation_H
0011 #define Pythia8_StringFragmentation_H
0012
0013 #include "Pythia8/FragmentationModel.h"
0014
0015 namespace Pythia8 {
0016
0017 //==========================================================================
0018
0019 // The StringEnd class contains the information related to
0020 // one of the current endpoints of the string system.
0021 // Only to be used inside StringFragmentation, so no private members.
0022
0023 class StringEnd {
0024
0025 public:
0026
0027   // Constructor.
0028   StringEnd() : particleDataPtr(), flavSelPtr(), pTSelPtr(), zSelPtr(),
0029     fromPos(), thermalModel(), mT2suppression(), iEnd(), iMax(), idHad(),
0030     iPosOld(), iNegOld(), iPosNew(), iNegNew(), hadSoFar(), colOld(), colNew(),
0031     pxOld(), pyOld(), pxNew(), pyNew(), pxHad(), pyHad(), mHad(), mT2Had(),
0032     zHad(), GammaOld(), GammaNew(), xPosOld(), xPosNew(), xPosHad(), xNegOld(),
0033     xNegNew(), xNegHad(), aLund(), bLund(), iPosOldPrev(), iNegOldPrev(),
0034     colOldPrev(), pxOldPrev(), pyOldPrev(), GammaOldPrev(), xPosOldPrev(),
0035     xNegOldPrev(), mVecRatio(1.), tinyEq(), pT2tiny() {}
0036
0037   // Save pointers.
0038   void init( ParticleData* particleDataPtrIn, StringFlav* flavSelPtrIn,
0039     StringPT* pTSelPtrIn, StringZ* zSelPtrIn, Settings& settings) {
0040     particleDataPtr = particleDataPtrIn; flavSelPtr = flavSelPtrIn;
0041     flavSelNow = *flavSelPtr;
0042     pTSelPtr = pTSelPtrIn; zSelPtr = zSelPtrIn;
0043     bLund = zSelPtr->bAreaLund(); aLund = zSelPtr->aAreaLund();
0044     thermalModel   = settings.flag("StringPT:thermalModel");
0045     mT2suppression = settings.flag("StringPT:mT2suppression");
0046     closePacking = settings.flag("ClosePacking:doClosePacking"); }
0047
0048   // Set up initial endpoint values from input.
0049   void setUp(bool fromPosIn, int iEndIn, int idOldIn, int iMaxIn,
0050     double pxIn, double pyIn, double GammaIn, double xPosIn,
0051     double xNegIn, int colIn, double mVecRatioIn);
0052
0053   // Fragment off one hadron from the string system, in flavour and pT.
0054   void newHadron(double kappaModifier, bool forbidPopcornNow = false,
0055     double strangeJunc = 0., double probQQmod = 1.);
0056
0057   // Fragment off one hadron from the string system, in momentum space,
0058   // by taking steps either from positive or from negative end.
0059   Vec4 kinematicsHadron(StringSystem& system, StringVertex& newVertex,
0060     double zHadIn);
0061
0062   // Generate momentum for some possible next hadron, based on mean values
0063   // to get an estimate for rapidity and pT.
0064   Vec4 kinematicsHadronTmp(StringSystem system, Vec4 pRem, double phi,
0065     double mult);
0066
0067   // Update string end information after a hadron has been removed.
0068   void update();
0069   void storePrev();
0070   void updateToPrev();
0071
0072   // Constants: could only be changed in the code itself.
0073   static const double TINY, PT2SAME, MEANMMIN, MEANM, MEANPT;
0074
0075   // Pointer to the particle data table.
0076   ParticleData* particleDataPtr;
0077
0078   // Pointers to classes for flavour, pT and z generation.
0079   StringFlav*   flavSelPtr;
0080   StringPT*     pTSelPtr;
0081   StringZ*      zSelPtr;
0082
0083   // Local copy of flavSelPtr for modified flavour selection.
0084   StringFlav    flavSelNow;
0085
0086   // Data members.
0087   bool   fromPos, thermalModel, mT2suppression, closePacking;
0088   int    iEnd, iMax, idHad, iPosOld, iNegOld, iPosNew, iNegNew, hadSoFar,
0089          colOld, colNew;
0090   double pxOld, pyOld, pxNew, pyNew, pxHad, pyHad, mHad, mT2Had, zHad,
0091          GammaOld, GammaNew, xPosOld, xPosNew, xPosHad, xNegOld, xNegNew,
0092          xNegHad, aLund, bLund;
0093   int    iPosOldPrev, iNegOldPrev, colOldPrev;
0094   double pxOldPrev, pyOldPrev, GammaOldPrev, xPosOldPrev, xNegOldPrev,
0095          mVecRatio, tinyEq, pT2tiny;
0096   FlavContainer flavOld, flavNew, flavOldPrev;
0097   Vec4   pHad, pSoFar;
0098
0099 };
0100
0101 //==========================================================================
0102
0103 // The StringFragmentation class contains the top-level routines
0104 // to fragment a colour singlet partonic system.
0105
0106 class StringFragmentation : public FragmentationModel {
0107
0108 public:
0109
0110   // Constructor.
0111   StringFragmentation() :
0112     FragmentationModel(), flavRopePtr(), closePacking(),
0113     setVertices(), constantTau(), smearOn(), traceColours(false),
0114     hadronVertex(), stopMass(), stopNewFlav(), stopSmear(),
0115     pNormJunction(), pMaxJunction(), eBothLeftJunction(),
0116     eMaxLeftJunction(), eMinLeftJunction(), mJoin(), bLund(),
0117     closePackingFluxRatio(1.), closePackingPT20(1.), pT20(),
0118     xySmear(), maxSmear(), maxTau(), kappaVtx(), mc(), mb(),
0119     hasJunction(), isClosed(), iPos(), iNeg(), nExtraJoin(),
0120     w2Rem(), stopMassNow(), mVecRatio(1.), closedM2max(),
0121     idDiquark(), legMin(), legMid() {}
0122
0123   // Initialize and save pointers.
0124   bool init(StringFlav* flavSelPtrIn = nullptr, StringPT* pTSelPtrIn = nullptr,
0125     StringZ* zSelPtrIn = nullptr, FragModPtr fragModPtrIn = nullptr) override;
0126
0127   // Do the fragmentation: driver routine.
0128   bool fragment(int iSub, ColConfig& colConfig, Event& event,
0129     bool isDiff = false, bool systemRecoil = true) override;
0130
0131   // Local copy of flavSelPtr for modified flavour selection.
0132   StringFlav flavSelNow;
0133
0134   // Find the boost matrix to the rest frame of a junction.
0135   Vec4 junctionRestFrame(const Vec4& p0, const Vec4& p1, const Vec4& p2,
0136     const bool angleCheck = true) const;
0137
0138   // Set the vector mass ratio.
0139   void setMVecRatio(double mVecRatioIn) {mVecRatio = mVecRatioIn;}
0140
0141 private:
0142
0143   // Constants: could only be changed in the code itself.
0144   static const int    NTRYFLAV, NTRYJOIN, NSTOPMASS,
0145                       NTRYJNMATCH, NTRYJRFEQ, NTRYSMEAR, MAXVETOFINTWO;
0146   static const double FACSTOPMASS, CLOSEDM2MAX, CLOSEDM2FRAC, EXPMAX,
0147                       MATCHPOSNEG, M2MINJRF, EMINJRF, EEXTRAJNMATCH,
0148                       MDIQUARKMIN, CONVJRFEQ, CHECKPOS;
0149
0150   // Pointer to flavour-composition-changing ropes.
0151   FragModPtr  flavRopePtr;
0152
0153   // Initialization data, read from Settings.
0154   bool   closePacking, setVertices, constantTau, smearOn,
0155          traceColours, hardRemn, doStrangeJunc;
0156   int    hadronVertex;
0157   double stopMass, stopNewFlav, stopSmear, pNormJunction, pMaxJunction,
0158          eBothLeftJunction, eMaxLeftJunction, eMinLeftJunction,
0159          mJoin, bLund, closePackingFluxRatio, closePackingPT20,
0160          qqSupPar, qqSupAnti, pT20, xySmear, maxSmear, maxTau,
0161          kappaVtx, mc, mb, dampPopcorn, aRemn, bRemn, strangeJuncParm;
0162
0163   // Data members.
0164   bool   hasJunction, isClosed;
0165   int    iPos, iNeg, nExtraJoin;
0166   double w2Rem, stopMassNow, kappaModifier, probQQmod, mVecRatio, closedM2max;
0167   Vec4   pSum, pRem, pJunctionHadrons;
0168
0169   // UserHooks flags.
0170   bool doChangeFragPar = false, doVetoFrag = false;
0171
0172   // List of partons in string system.
0173   vector<int> iParton, iPartonMinLeg, iPartonMidLeg, iPartonMax;
0174
0175   // Vertex information from the fragmentation process.
0176   vector<StringVertex> stringVertices, legMinVertices, legMidVertices;
0177   StringVertex newVertex;
0178
0179   // Variables used for JRF finding.
0180   vector<Vec4>   listJRF;
0181   vector<double> weightJRF;
0182   int    iLeg[3], idLeg[3], legEnd[3];
0183   double weightSum, pSumJRF, m2Leg[3];
0184   Vec4   pLeg[3];
0185   bool   lastJRF, endpoint[3];
0186
0187   // Boost from/to rest frame of a junction to original frame.
0188   RotBstMatrix MfromJRF, MtoJRF;
0189
0190   // Information on diquark created at the junction.
0191   int    idDiquark;
0192
0193   // Fictitious opposing partons in JRF: string ends for vertex location.
0194   Vec4 pMinEnd, pMidEnd;
0195
0196   // Temporary event record for the produced particles.
0197   Event hadrons;
0198
0199   // Information on the system of string regions.
0200   StringSystem system, systemMin, systemMid;
0201
0202   // Information on the two current endpoints of the fragmenting
0203   // system, with backup copies if there is a veto.
0204   StringEnd posEnd, negEnd, posEndSave, negEndSave;
0205
0206   // Find region where to put first string break for closed gluon loop.
0207   vector<int> findFirstRegion(int iSub, const ColConfig& colConfig,
0208     const Event& event) const;
0209
0210   // Set flavours and momentum position for initial string endpoints.
0211   void setStartEnds(int idPos, int idNeg, const StringSystem& systemNow,
0212     int legNow = 3);
0213
0214   // Check remaining energy-momentum whether it is OK to continue.
0215   bool energyUsedUp(bool fromPos);
0216
0217   // Produce the final two partons to complete the system.
0218   bool finalTwo(bool fromPos, const Event& event, bool usedPosJun,
0219     bool usedNegJun);
0220
0221   // Final region information.
0222   Vec4 pPosFinalReg, pNegFinalReg, eXFinalReg, eYFinalReg;
0223
0224   // Set hadron production points in space-time picture.
0225   bool setHadronVertices(Event& event);
0226
0227   // Construct a special joining region for the final two hadrons.
0228   StringRegion finalRegion();
0229
0230   // Store the hadrons in the normal event record, ordered from one end.
0231   void store(Event& event);
0232
0233   // Fragment off two of the string legs in to a junction.
0234   bool fragmentToJunction(Event& event,
0235     vector< vector< pair<double,double> > >& rapPairs);
0236
0237   // Initially considered legs from the junction.
0238   int legMin, legMid;
0239
0240   // Functions used in JRF iterative procedure.
0241   bool   collinearPair(Event& event);
0242   bool   perturbedJRF(Event& event);
0243   int    updateLegs(Event& event, Vec4 vJunIn, bool juncCoM = false);
0244   double updateWeights(double pSmall, Vec4 vJunIn);
0245   void   nextParton(Event& event, int leg);
0246
0247   // Join extra nearby partons when stuck.
0248   int extraJoin(double facExtra, Event& event);
0249
0250   // Get the number of nearby strings given the energies.
0251   void kappaEffModifier(StringSystem& systemNow,
0252     StringEnd end, bool fromPos, vector<int> partonList,
0253     vector< vector< pair<double,double> > >& rapPairs,
0254     double mRem, Event& event);
0255
0256   double yMax(Particle pIn, double mTiny) {
0257     double temp = log( ( pIn.e() + abs(pIn.pz()) ) / max( mTiny, pIn.mT()) );
0258     return (pIn.pz() > 0) ? temp : -temp; }
0259
0260 };
0261
0262 //==========================================================================
0263
0264 } // end namespace Pythia8
0265
0266 #endif // Pythia8_StringFragmentation_H