PYTHIA  8.312
HadronLevel.h
1 // HadronLevel.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2024 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 main class for hadron-level generation.
7 // HadronLevel: handles administration of fragmentation and decay.
8 
9 #ifndef Pythia8_HadronLevel_H
10 #define Pythia8_HadronLevel_H
11 
12 #include "Pythia8/Basics.h"
13 #include "Pythia8/BoseEinstein.h"
14 #include "Pythia8/ColourTracing.h"
15 #include "Pythia8/DeuteronProduction.h"
16 #include "Pythia8/Event.h"
17 #include "Pythia8/FragmentationFlavZpT.h"
18 #include "Pythia8/FragmentationSystems.h"
19 #include "Pythia8/HadronWidths.h"
20 #include "Pythia8/HiddenValleyFragmentation.h"
21 #include "Pythia8/Info.h"
22 #include "Pythia8/JunctionSplitting.h"
23 #include "Pythia8/LowEnergyProcess.h"
24 #include "Pythia8/SigmaLowEnergy.h"
25 #include "Pythia8/MiniStringFragmentation.h"
26 #include "Pythia8/NucleonExcitations.h"
27 #include "Pythia8/ParticleData.h"
28 #include "Pythia8/ParticleDecays.h"
29 #include "Pythia8/PartonVertex.h"
30 #include "Pythia8/PhysicsBase.h"
31 #include "Pythia8/PythiaStdlib.h"
32 #include "Pythia8/RHadrons.h"
33 #include "Pythia8/Settings.h"
34 #include "Pythia8/StringFragmentation.h"
35 #include "Pythia8/TimeShower.h"
36 #include "Pythia8/UserHooks.h"
37 
38 namespace Pythia8 {
39 
40 //==========================================================================
41 
42 // The HadronLevel class contains the top-level routines to generate
43 // the transition from the partonic to the hadronic stage of an event.
44 
45 class HadronLevel : public PhysicsBase {
46 
47 public:
48 
49  // Constructor.
50  HadronLevel() = default;
51 
52  // Initialize HadronLevel classes as required.
53  bool init( TimeShowerPtr timesDecPtr, RHadrons* rHadronsPtrIn,
54  DecayHandlerPtr decayHandlePtr, vector<int> handledParticles,
55  StringIntPtr stringInteractionsPtrIn, PartonVertexPtr partonVertexPtrIn,
56  SigmaLowEnergy& sigmaLowEnergyIn,
57  NucleonExcitations& nucleonExcitationsIn);
58 
59  // Get pointer to StringFlav instance (needed by BeamParticle).
60  StringFlav* getStringFlavPtr() {return &flavSel;}
61 
62  // Generate the next event.
63  bool next(Event& event);
64 
65  // Try to decay the specified particle. Returns false if decay failed.
66  bool decay( int iDec, Event& event) { return
67  (event[iDec].isFinal() && event[iDec].canDecay() && event[iDec].mayDecay())
68  ? decays.decay( iDec, event) : true;}
69 
70  // Special routine to allow more decays if on/off switches changed.
71  bool moreDecays(Event& event);
72 
73  // Perform rescattering. Return true if new strings must be hadronized.
74  bool rescatter(Event& event);
75 
76  // Prepare and pick process for a low-energy hadron-hadron scattering.
78  int pickLowEnergyProcess(int idA, int idB, double eCM, double mA, double mB);
79 
80  // Special routine to do a low-energy hadron-hadron scattering.
81  bool doLowEnergyProcess(int i1, int i2, int procTypeIn, Event& event) {
82  if (!lowEnergyProcess.collide( i1, i2, procTypeIn, event)) {
83  loggerPtr->ERROR_MSG("low energy collision failed");
84  return false;
85  }
86  return true;
87  }
88 
89  // Tell if we did an early user-defined veto of the event.
90  bool hasVetoedHadronize() const {return doHadronizeVeto; }
91 
92 protected:
93 
94  virtual void onInitInfoPtr() override{
95  registerSubObject(flavSel);
96  registerSubObject(pTSel);
97  registerSubObject(zSel);
98  registerSubObject(stringFrag);
99  registerSubObject(ministringFrag);
100  registerSubObject(decays);
101  registerSubObject(lowEnergyProcess);
102  registerSubObject(boseEinstein);
103  registerSubObject(hiddenvalleyFrag);
104  registerSubObject(junctionSplitting);
105  registerSubObject(deuteronProd);
106  }
107 
108 private:
109 
110  // Constants: could only be changed in the code itself.
111  static const double MTINY;
112 
113  // Initialization data, read from Settings.
114  bool doHadronize{}, doDecay{}, doPartonVertex{}, doBoseEinstein{},
115  doDeuteronProd{}, allowRH{}, closePacking{}, doNonPertAll{};
116  double mStringMin{}, pNormJunction{}, widthSepBE{}, widthSepRescatter{};
117  vector<int> nonPertProc{};
118 
119  // Configuration of colour-singlet systems.
120  ColConfig colConfig;
121 
122  // Colour and mass information.
123  vector<int> iParton{}, iJunLegA{}, iJunLegB{}, iJunLegC{},
124  iAntiLegA{}, iAntiLegB{}, iAntiLegC{}, iGluLeg{};
125  vector<double> m2Pair{};
126 
127  // The generator class for normal string fragmentation.
128  StringFragmentation stringFrag;
129 
130  // The generator class for special low-mass string fragmentation.
131  MiniStringFragmentation ministringFrag;
132 
133  // Try ministring fragmentation also if normal fails.
134  bool tryMiniAfterFailedFrag{};
135 
136  // The generator class for normal decays.
137  ParticleDecays decays;
138 
139  // The generator class for Bose-Einstein effects.
140  BoseEinstein boseEinstein;
141 
142  // The generator class for deuteron production.
143  DeuteronProduction deuteronProd;
144 
145  // Classes for flavour, pT and z generation.
146  StringFlav flavSel;
147  StringPT pTSel;
148  StringZ zSel;
149 
150  // Class for colour tracing.
151  ColourTracing colTrace;
152 
153  // Junction splitting class.
154  JunctionSplitting junctionSplitting;
155 
156  // The RHadrons class is used to fragment off and decay R-hadrons.
157  RHadrons* rHadronsPtr;
158 
159  // Special class for Hidden-Valley hadronization. Not always used.
160  HiddenValleyFragmentation hiddenvalleyFrag;
161  bool useHiddenValley{};
162 
163  // Special case: colour-octet onium decays, to be done initially.
164  bool decayOctetOnia(Event& event);
165 
166  // Trace colour flow in the event to form colour singlet subsystems.
167  // Option to keep junctions, needed for rope hadronization.
168  bool findSinglets(Event& event, bool keepJunctions = false);
169 
170  // Class to displace hadron vertices from parton impact-parameter picture.
171  PartonVertexPtr partonVertexPtr;
172 
173  // Hadronic rescattering.
174  class PriorityNode;
175  bool doRescatter{}, scatterManyTimes{}, scatterQuickCheck{},
176  scatterNeighbours{}, delayRegeneration{};
177  double b2Max, tauRegeneration{};
178  void queueDecResc(Event& event, int iStart,
179  priority_queue<HadronLevel::PriorityNode>& queue);
180  int boostDir;
181  double boost;
182  bool doBoost;
183  bool useVelocityFrame;
184 
185  // User veto performed right after string hadronization.
186  bool doHadronizeVeto;
187 
188  // The generator class for low-energy hadron-hadron collisions.
189  LowEnergyProcess lowEnergyProcess;
190  int impactModel{};
191  double impactOpacity{};
192 
193  // Cross sections for low-energy processes.
194  SigmaLowEnergy* sigmaLowEnergyPtr = {};
195 
196  // Nucleon excitations data.
197  NucleonExcitations* nucleonExcitationsPtr = {};
198 
199  // Class for event geometry for Rope Hadronization. Production vertices.
200  StringRepPtr stringRepulsionPtr;
201  FragModPtr fragmentationModifierPtr;
202 
203  // Extract rapidity pairs.
204  vector< vector< pair<double,double> > > rapidityPairs(Event& event);
205 
206  // Calculate the rapidity for string ends, protected against too large y.
207  double yMax(Particle pIn, double mTiny) {
208  double temp = log( ( pIn.e() + abs(pIn.pz()) ) / max( mTiny, pIn.mT()) );
209  return (pIn.pz() > 0) ? temp : -temp; }
210  double yMax(Vec4 pIn, double mTiny) {
211  double mTemp = pIn.m2Calc() + pIn.pT2();
212  mTemp = (mTemp >= 0.) ? sqrt(mTemp) : -sqrt(-mTemp);
213  double temp = log( ( pIn.e() + abs(pIn.pz()) ) / max( mTiny, mTemp) );
214  return (pIn.pz() > 0) ? temp : -temp; }
215 
216 };
217 
218 //==========================================================================
219 
220 } // end namespace Pythia8
221 
222 #endif // Pythia8_HadronLevel_H
Definition: LowEnergyProcess.h:28
bool collide(int i1, int i2, int typeIn, Event &event, Vec4 vtx=Vec4(), Vec4 vtx1=Vec4(), Vec4 vtx2=Vec4())
Produce outgoing primary hadrons from collision of incoming pair.
Definition: LowEnergyProcess.cc:105
bool moreDecays(Event &event)
Special routine to allow more decays if on/off switches changed.
Definition: HadronLevel.cc:336
Definition: NucleonExcitations.h:23
Definition: PhysicsBase.h:27
void registerSubObject(PhysicsBase &pb)
Register a sub object that should have its information in sync with this.
Definition: PhysicsBase.cc:56
The ParticleDecays class contains the routines to decay a particle.
Definition: ParticleDecays.h:57
The Event class holds all info on the generated event.
Definition: Event.h:453
Definition: HadronLevel.cc:19
Definition: JunctionSplitting.h:30
The StringPT class is used to select select transverse momenta.
Definition: FragmentationFlavZpT.h:322
Gets cross sections for hadron-hadron collisions at low energies.
Definition: SigmaLowEnergy.h:22
Definition: StringFragmentation.h:119
StringFlav * getStringFlavPtr()
Get pointer to StringFlav instance (needed by BeamParticle).
Definition: HadronLevel.h:60
Definition: RHadrons.h:29
int pickLowEnergyProcess(int idA, int idB, double eCM, double mA, double mB)
Pick process for a low-energy hadron-hadron scattering.
Definition: HadronLevel.cc:380
HadronLevel()=default
Constructor.
The StringZ class is used to sample the fragmentation function f(z).
Definition: FragmentationFlavZpT.h:264
Definition: HadronLevel.h:45
bool next(Event &event)
Generate the next event.
Definition: HadronLevel.cc:172
Definition: HiddenValleyFragmentation.h:119
bool rescatter(Event &event)
Perform rescattering. Return true if new strings must be hadronized.
Definition: HadronLevel.cc:617
bool doLowEnergyProcess(int i1, int i2, int procTypeIn, Event &event)
Special routine to do a low-energy hadron-hadron scattering.
Definition: HadronLevel.h:81
bool hasVetoedHadronize() const
Tell if we did an early user-defined veto of the event.
Definition: HadronLevel.h:90
ColourTracing class. It is used to trace colours within the event record.
Definition: ColourTracing.h:22
bool initLowEnergyProcesses()
Prepare and pick process for a low-energy hadron-hadron scattering.
Definition: HadronLevel.cc:356
Definition: Event.h:32
The ColConfig class describes the colour configuration of the whole event.
Definition: FragmentationSystems.h:60
Definition: BoseEinstein.h:48
Logger * loggerPtr
Pointer to logger.
Definition: PhysicsBase.h:87
bool decay(int iDec, Event &event)
Perform a decay of a single particle.
Definition: ParticleDecays.cc:118
virtual void onInitInfoPtr() override
Definition: HadronLevel.h:94
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:84
The DeuteronProduction class.
Definition: DeuteronProduction.h:22
Definition: MiniStringFragmentation.h:30
bool init(TimeShowerPtr timesDecPtr, RHadrons *rHadronsPtrIn, DecayHandlerPtr decayHandlePtr, vector< int > handledParticles, StringIntPtr stringInteractionsPtrIn, PartonVertexPtr partonVertexPtrIn, SigmaLowEnergy &sigmaLowEnergyIn, NucleonExcitations &nucleonExcitationsIn)
Initialize HadronLevel classes as required.
Definition: HadronLevel.cc:56
bool decay(int iDec, Event &event)
Try to decay the specified particle. Returns false if decay failed.
Definition: HadronLevel.h:66
Definition: Basics.h:32