include/Pythia8/StringFragmentation.h

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