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     doFlavBeforePT(), fromPos(), iEnd(), iMax(), idHad(), iPosOld(),
0030     iNegOld(), iPosNew(), iNegNew(), iPosNewTmp(), iNegNewTmp(),
0031     hadSoFar(), colOld(), colNew(),
0032     pxOld(), pyOld(), pxNew(), pyNew(), pxHad(), pyHad(), mHad(),
0033     mT2Had(), zHad(), GammaOld(), GammaNew(), xPosOld(), xPosNew(),
0034     xPosHad(), xNegOld(), xNegNew(), xNegHad(), aLund(), bLund(),
0035     iPosOldPrev(), iNegOldPrev(), colOldPrev(), pxOldPrev(),
0036     pyOldPrev(), GammaOldPrev(), xPosOldPrev(), xNegOldPrev(),
0037     mVecRatio(1.), tinyEq(), pT2tiny() {}
0038
0039   // Save pointers.
0040   void init( ParticleData* particleDataPtrIn, StringFlav* flavSelPtrIn,
0041     StringPT* pTSelPtrIn, StringZ* zSelPtrIn, Settings& settings,
0042     bool doFlavBeforePTin = true) { particleDataPtr = particleDataPtrIn;
0043     flavSelPtr = flavSelPtrIn; flavSelNow = *flavSelPtr; pTSelPtr = pTSelPtrIn;
0044     zSelPtr = zSelPtrIn; doFlavBeforePT = doFlavBeforePTin;
0045     bLund = zSelPtr->bAreaLund(); aLund = zSelPtr->aAreaLund();
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   doFlavBeforePT, fromPos, closePacking;
0088   int    iEnd, iMax, idHad, iPosOld, iNegOld, iPosNew, iNegNew, iPosNewTmp,
0089          iNegNewTmp, hadSoFar, 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(), doFlavBeforePT(true),
0113     closePacking(), setVertices(), constantTau(), smearOn(),
0114     traceColours(false), hadronVertex(), stopMass(), stopNewFlav(),
0115     stopSmear(), pNormJunction(), pMaxJunction(), eJunctionCutoff(),
0116     mJunctionCutoff(), eBothLeftJunction(),
0117     eMaxLeftJunction(), eMinLeftJunction(), mJoin(), bLund(),
0118     closePackingFluxRatio(1.), closePackingPT20(1.), pT20(),
0119     xySmear(), maxSmear(), maxTau(), kappaVtx(), mc(), mb(),
0120     hasJunction(), isClosed(), iPos(), iNeg(), nExtraJoin(),
0121     w2Rem(), stopMassNow(), mVecRatio(1.),
0122     idDiquark(), legMin(), legMid() {}
0123
0124   // Set order of flavour and pT selection in StringEnd.
0125   void setFlavBeforePT( bool doFlavBeforePTin) {
0126     doFlavBeforePT = doFlavBeforePTin; }
0127
0128   // Initialize and save pointers.
0129   bool init(StringFlav* flavSelPtrIn = nullptr, StringPT* pTSelPtrIn = nullptr,
0130     StringZ* zSelPtrIn = nullptr, FragModPtr fragModPtrIn = nullptr) override;
0131
0132   // Do the fragmentation: driver routine.
0133   bool fragment(int iSub, ColConfig& colConfig, Event& event,
0134     bool isDiff = false, bool systemRecoil = true) override;
0135
0136   // Local copy of flavSelPtr for modified flavour selection.
0137   StringFlav flavSelNow;
0138
0139   // Find the boost matrix to the rest frame of a junction.
0140   Vec4 junctionRestFrame(const Vec4& p0, const Vec4& p1, const Vec4& p2,
0141     const bool angleCheck = true) const;
0142
0143   // Set the vector mass ratio.
0144   void setMVecRatio(double mVecRatioIn) {mVecRatio = mVecRatioIn;}
0145
0146 private:
0147
0148   // Constants: could only be changed in the code itself.
0149   static const int    NTRYFLAV, NTRYJOIN, NSTOPMASS,
0150                       NTRYJNMATCH, NTRYJRFEQ, NTRYSMEAR, MAXVETOFINTWO;
0151   static const double FACSTOPMASS, CLOSEDM2MAX, CLOSEDM2FRAC, EXPMAX,
0152                       MATCHPOSNEG, M2MINJRF, EMINJRF, EEXTRAJNMATCH,
0153                       MDIQUARKMIN, CONVJRFEQ, CHECKPOS;
0154
0155   // Pointer to flavour-composition-changing ropes.
0156   FragModPtr  flavRopePtr;
0157
0158   // Initialization data, read from Settings.
0159   bool   doFlavBeforePT, closePacking, setVertices, constantTau, smearOn,
0160          traceColours, hardRemn, doStrangeJunc;
0161   int    hadronVertex;
0162   double stopMass, stopNewFlav, stopSmear, pNormJunction, pMaxJunction,
0163          eJunctionCutoff, mJunctionCutoff,
0164          eBothLeftJunction, eMaxLeftJunction, eMinLeftJunction,
0165          mJoin, bLund, closePackingFluxRatio, closePackingPT20,
0166          qqSupPar, qqSupAnti, pT20, xySmear, maxSmear, maxTau,
0167          kappaVtx, mc, mb, dampPopcorn, aRemn, bRemn, strangeJuncParm;
0168
0169   // Data members.
0170   bool   hasJunction, isClosed;
0171   int    iPos, iNeg, nExtraJoin;
0172   double w2Rem, stopMassNow, kappaModifier, probQQmod, mVecRatio;
0173   Vec4   pSum, pRem, pJunctionHadrons;
0174
0175   // UserHooks flags.
0176   bool doChangeFragPar = false, doVetoFrag = false;
0177
0178   // List of partons in string system.
0179   vector<int> iParton, iPartonMinLeg, iPartonMidLeg, iPartonMax;
0180
0181   // Vertex information from the fragmentation process.
0182   vector<StringVertex> stringVertices, legMinVertices, legMidVertices;
0183   StringVertex newVertex;
0184
0185   // Variables used for JRF finding.
0186   vector<Vec4>   listJRF;
0187   vector<double> weightJRF;
0188   int    iLeg[3], idLeg[3], legEnd[3];
0189   double weightSum, pSumJRF, m2Leg[3];
0190   Vec4   pLeg[3];
0191   bool   lastJRF, endpoint[3];
0192
0193   // Boost from/to rest frame of a junction to original frame.
0194   RotBstMatrix MfromJRF, MtoJRF;
0195
0196   // Information on diquark created at the junction.
0197   int    idDiquark;
0198
0199   // Fictitious opposing partons in JRF: string ends for vertex location.
0200   Vec4 pMinEnd, pMidEnd;
0201
0202   // Temporary event record for the produced particles.
0203   Event hadrons;
0204
0205   // Information on the system of string regions.
0206   StringSystem system, systemMin, systemMid;
0207
0208   // Information on the two current endpoints of the fragmenting
0209   // system, with backup copies if there is a veto.
0210   StringEnd posEnd, negEnd, posEndSave, negEndSave;
0211
0212   // Find region where to put first string break for closed gluon loop.
0213   vector<int> findFirstRegion(int iSub, const ColConfig& colConfig,
0214     const Event& event) const;
0215
0216   // Set flavours and momentum position for initial string endpoints.
0217   void setStartEnds(int idPos, int idNeg, const StringSystem& systemNow,
0218     int legNow = 3);
0219
0220   // Check remaining energy-momentum whether it is OK to continue.
0221   bool energyUsedUp(bool fromPos);
0222
0223   // Produce the final two partons to complete the system.
0224   bool finalTwo(bool fromPos, const Event& event, bool usedPosJun,
0225     bool usedNegJun);
0226
0227   // Final region information.
0228   Vec4 pPosFinalReg, pNegFinalReg, eXFinalReg, eYFinalReg;
0229
0230   // Set hadron production points in space-time picture.
0231   bool setHadronVertices(Event& event);
0232
0233   // Construct a special joining region for the final two hadrons.
0234   StringRegion finalRegion();
0235
0236   // Store the hadrons in the normal event record, ordered from one end.
0237   void store(Event& event);
0238
0239   // Fragment off two of the string legs in to a junction.
0240   bool fragmentToJunction(Event& event,
0241     vector< vector< pair<double,double> > >& rapPairs);
0242
0243   // Initially considered legs from the junction.
0244   int legMin, legMid;
0245
0246   // Functions used in JRF iterative procedure.
0247   bool   collinearPair(Event& event);
0248   bool   perturbedJRF(Event& event);
0249   int    updateLegs(Event& event, Vec4 vJunIn, bool juncCoM = false);
0250   double updateWeights(double pSmall, Vec4 vJunIn);
0251   void   nextParton(Event& event, int leg);
0252
0253   // Join extra nearby partons when stuck.
0254   int extraJoin(double facExtra, Event& event);
0255
0256   // Get the number of nearby strings given the energies.
0257   void kappaEffModifier(StringSystem& systemNow,
0258     StringEnd end, vector<int> partonList,
0259     vector< vector< pair<double,double> > >& rapPairs,
0260     Vec4 pRemNow, Event& event);
0261
0262   double yMax(Particle pIn, double mTiny) {
0263     double temp = log( ( pIn.e() + abs(pIn.pz()) ) / max( mTiny, pIn.mT()) );
0264     return (pIn.pz() > 0) ? temp : -temp; }
0265
0266   double yMax(Vec4 pIn, double mTiny) {
0267     double mTemp = pIn.m2Calc() + pIn.pT2();
0268     mTemp = (mTemp >= 0.) ? sqrt(mTemp) : -sqrt(-mTemp);
0269     double temp = log( ( pIn.e() + abs(pIn.pz()) ) / max( mTiny, mTemp) );
0270     return (pIn.pz() > 0) ? temp : -temp; }
0271
0272 };
0273
0274 //==========================================================================
0275
0276 } // end namespace Pythia8
0277
0278 #endif // Pythia8_StringFragmentation_H