PYTHIA  8.312
SimpleSpaceShower.h
1 // SimpleSpaceShower.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 // Header file for the original simple spacelike initial-state showers.
7 // SpaceDipoleEnd: data on a radiating dipole end in ISR.
8 // SimpleSpaceShower: handles the showering description.
9 
10 #ifndef Pythia8_SimpleSpaceShower_H
11 #define Pythia8_SimpleSpaceShower_H
12 
13 #include "Pythia8/SpaceShower.h"
14 #include "Pythia8/SimpleWeakShowerMEs.h"
15 
16 namespace Pythia8 {
17 
18 //==========================================================================
19 
20 // Data on radiating dipole ends, only used inside SimpleSpaceShower.
21 
23 
24 public:
25 
26  // Constructor.
27  SpaceDipoleEnd( int systemIn = 0, int sideIn = 0, int iRadiatorIn = 0,
28  int iRecoilerIn = 0, double pTmaxIn = 0., int colTypeIn = 0,
29  int chgTypeIn = 0, int weakTypeIn = 0, int MEtypeIn = 0,
30  bool normalRecoilIn = true, int weakPolIn = 0,
31  int iColPartnerIn = 0, int idColPartnerIn = 0) :
32  system(systemIn), side(sideIn), iRadiator(iRadiatorIn),
33  iRecoiler(iRecoilerIn), pTmax(pTmaxIn), colType(colTypeIn),
34  chgType(chgTypeIn), weakType(weakTypeIn), MEtype(MEtypeIn),
35  normalRecoil(normalRecoilIn), weakPol(weakPolIn),
36  iColPartner(iColPartnerIn), idColPartner(idColPartnerIn),
37  nBranch(0), idDaughter(), idMother(), idSister(), iFinPol(), x1(), x2(),
38  m2Dip(), pT2(), z(), xMo(), Q2(), mSister(), m2Sister(), pT2corr(),
39  pT2Old(0.), zOld(0.5), asymPol(), m2IF(), mColPartner(),
40  pAccept() { }
41 
42  // Store values for trial emission.
43  void store( int idDaughterIn, int idMotherIn, int idSisterIn,
44  double x1In, double x2In, double m2DipIn, double pT2In, double zIn,
45  double xMoIn, double Q2In, double mSisterIn, double m2SisterIn,
46  double pT2corrIn, int iColPartnerIn, double m2IFIn, double mColPartnerIn)
47  {idDaughter = idDaughterIn; idMother = idMotherIn;
48  idSister = idSisterIn; x1 = x1In; x2 = x2In; m2Dip = m2DipIn;
49  pT2 = pT2In; z = zIn; xMo = xMoIn; Q2 = Q2In; mSister = mSisterIn;
50  m2Sister = m2SisterIn; pT2corr = pT2corrIn; iColPartner = iColPartnerIn;
51  m2IF = m2IFIn; mColPartner = mColPartnerIn;}
52 
53  // Basic properties related to evolution and matrix element corrections.
54  int system, side, iRadiator, iRecoiler;
55  double pTmax;
56  int colType, chgType, weakType, MEtype;
57  bool normalRecoil;
58  int weakPol, iColPartner, idColPartner;
59 
60  // Properties specific to current trial emission.
61  int nBranch, idDaughter, idMother, idSister, iFinPol;
62  double x1, x2, m2Dip, pT2, z, xMo, Q2, mSister, m2Sister, pT2corr,
63  pT2Old, zOld, asymPol, m2IF, mColPartner;
64 
65  // Properties needed for the evaluation of parameter variations
66  double pAccept;
67 
68 } ;
69 
70 //==========================================================================
71 
72 // The SimpleSpaceShower class does spacelike showers.
73 
75 
76 public:
77 
78  // Constructor.
79  SimpleSpaceShower() : rescatterFail(), gamma2qqbar(), hasWeaklyRadiated(),
80  iSysSel(), pTmaxFudge(), doQCDshower(), doQEDshowerByQ(), doQEDshowerByL(),
81  useSamePTasMPI(), doWeakShower(), doMEcorrections(), doMEafterFirst(),
82  doPhiPolAsym(), doPhiPolAsymHard(), doPhiIntAsym(), doRapidityOrder(),
83  useFixedFacScale(), doSecondHard(), canVetoEmission(), hasUserHooks(),
84  alphaSuseCMW(), singleWeakEmission(), vetoWeakJets(), weakExternal(),
85  doRapidityOrderMPI(), doMPI(), doDipoleRecoil(), doPartonVertex(),
86  pTmaxMatch(), pTdampMatch(), alphaSorder(), alphaSnfmax(), alphaEMorder(),
87  nQuarkIn(), enhanceScreening(), weakMode(), pT0paramMode(), pTdampFudge(),
88  mc(), mb(), m2c(), m2b(), renormMultFac(), factorMultFac(),
89  fixedFacScale2(), alphaSvalue(), alphaS2pi(), Lambda3flav(), Lambda4flav(),
90  Lambda5flav(), Lambda3flav2(), Lambda4flav2(), Lambda5flav2(), pT0Ref(),
91  ecmRef(), ecmPow(), pTmin(), sCM(), eCM(), pT0(), pTminChgQ(), pTminChgL(),
92  pT20(), pT2min(), pT2minChgQ(), pT2minChgL(), pTweakCut(), pT2weakCut(),
93  pTmaxFudgeMPI(), strengthIntAsym(), weakEnhancement(), mZ(), gammaZ(),
94  thetaWRat(), mW(), gammaW(), weakMaxWt(), vetoWeakDeltaR2(), sideA(),
95  twoHard(), dopTlimit1(), dopTlimit2(), dopTdamp(), tChannel(),
96  doUncertaintiesNow(), iNow(), iRec(), idDaughter(), nRad(), idResFirst(),
97  idResSecond(), xDaughter(), x1Now(), x2Now(),
98  m2ColPair(), mColPartner(), m2ColPartner(), m2Dip(), m2Rec(), pT2damp(),
99  pTbegRef(), pdfScale2(), doTrialNow(), canEnhanceEmission(),
100  canEnhanceTrial(), canEnhanceET(), iDipNow(), iSysNow(), dipEndNow(),
101  iDipSel(), dipEndSel() { beamOffset = 0; pdfMode = 0; }
102 
103  // Destructor.
104  virtual ~SimpleSpaceShower() override {}
105 
106  // Initialize generation. Possibility to force re-initialization by hand.
107  virtual void init(BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn)
108  override;
109 
110  // Find whether to limit maximum scale of emissions, and whether to dampen.
111  virtual bool limitPTmax( Event& event, double Q2Fac = 0.,
112  double Q2Ren = 0.) override;
113 
114  // Prepare system for evolution; identify ME.
115  virtual void prepare( int iSys, Event& event, bool limitPTmaxIn = true)
116  override;
117 
118  // Update dipole list after each FSR emission.
119  virtual void update( int iSys, Event& event, bool hasWeakRad = false)
120  override;
121 
122  // Select next pT in downwards evolution.
123  virtual double pTnext( Event& event, double pTbegAll, double pTendAll,
124  int nRadIn = -1, bool doTrialIn = false) override;
125 
126  // ME corrections and kinematics that may give failure.
127  virtual bool branch( Event& event) override;
128 
129  // Print dipole list; for debug mainly.
130  virtual void list() const override;
131 
132  // Initialize data members for calculation of uncertainty bands.
133  virtual bool initUncertainties() override;
134 
135  // Initialize data members for application of enhancements.
136  virtual bool initEnhancements() override;
137 
138  // Flag for failure in branch(...) that will force a retry of parton level.
139  virtual bool doRestart() const override {return rescatterFail;}
140 
141  // Tell if latest scattering was a gamma->qqbar.
142  virtual bool wasGamma2qqbar() override { return gamma2qqbar; }
143 
144  // Tell whether ISR has done a weak emission.
145  virtual bool getHasWeaklyRadiated() override {return hasWeaklyRadiated;}
146 
147  // Tell which system was the last processed one.
148  virtual int system() const override {return iSysSel;}
149 
150  // Potential enhancement factor of pTmax scale for hardest emission.
151  virtual double enhancePTmax() const override {return pTmaxFudge;}
152 
153  // Functions to directly extract the probability of no emission between two
154  // scales. These functions are not used in the Pythia core code, but can be
155  // used by external programs to interface with the shower directly.
156  double noEmissionProbability( double pTbegAll, double pTendAll, double m2dip,
157  int id, int type, double s = -1., double x = -1.) override;
158  double pTnext( vector<SpaceDipoleEnd> dipEnds, Event event, double pTbegAll,
159  double pTendAll, double m2dip, int id, int type, double s = -1.,
160  double x = -1.);
161  int pdfMode;
162 
163 private:
164 
165  // Constants: could only be changed in the code itself.
166  static const int MAXLOOPTINYPDF;
167  static const double MCMIN, MBMIN, CTHRESHOLD, BTHRESHOLD, EVALPDFSTEP,
168  TINYPDF, TINYKERNELPDF, TINYPT2, HEAVYPT2EVOL, HEAVYXEVOL,
169  EXTRASPACEQ, LAMBDA3MARGIN, PT2MINWARN, LEPTONXMIN, LEPTONXMAX,
170  LEPTONPT2MIN, LEPTONFUDGE, WEAKPSWEIGHT, HEADROOMQ2Q, HEADROOMQ2G,
171  HEADROOMG2G, HEADROOMG2Q, HEADROOMHQG, REJECTFACTOR, PROBLIMIT;
172 
173  // Store properties to be returned by methods.
174  bool rescatterFail, gamma2qqbar, hasWeaklyRadiated;
175  int iSysSel;
176  double pTmaxFudge;
177 
178  // Initialization data, normally only set once.
179  bool doQCDshower, doQEDshowerByQ, doQEDshowerByL, useSamePTasMPI,
180  doWeakShower, doMEcorrections, doMEafterFirst, doPhiPolAsym,
181  doPhiPolAsymHard, doPhiIntAsym, doRapidityOrder, useFixedFacScale,
182  doSecondHard, canVetoEmission, hasUserHooks, alphaSuseCMW,
183  singleWeakEmission, vetoWeakJets, weakExternal, doRapidityOrderMPI,
184  doMPI, doDipoleRecoil, doPartonVertex;
185  int pdfModeSave;
186  int pTmaxMatch, pTdampMatch, alphaSorder, alphaSnfmax, alphaEMorder,
187  nQuarkIn, enhanceScreening, weakMode, pT0paramMode;
188  double pTdampFudge, mc, mb, m2c, m2b, renormMultFac, factorMultFac,
189  fixedFacScale2, alphaSvalue, alphaS2pi, Lambda3flav, Lambda4flav,
190  Lambda5flav, Lambda3flav2, Lambda4flav2, Lambda5flav2, pT0Ref,
191  ecmRef, ecmPow, pTmin, sCM, eCM, pT0, pTminChgQ, pTminChgL, pT20,
192  pT2min, pT2minChgQ, pT2minChgL, pTweakCut, pT2weakCut, pTmaxFudgeMPI,
193  strengthIntAsym, weakEnhancement, mZ, gammaZ, thetaWRat, mW, gammaW,
194  weakMaxWt, vetoWeakDeltaR2;
195 
196  // alphaStrong and alphaEM calculations.
197  AlphaStrong alphaS;
198  AlphaEM alphaEM;
199 
200  // Weak matrix elements used for corrections both of ISR and FSR.
201  SimpleWeakShowerMEs simpleWeakShowerMEs;
202 
203  // Some current values.
204  bool sideA, twoHard, dopTlimit1, dopTlimit2, dopTdamp, tChannel,
205  doUncertaintiesNow;
206  int iNow, iRec, idDaughter, nRad, idResFirst, idResSecond;
207  double xDaughter, x1Now, x2Now, m2ColPair, mColPartner, m2ColPartner,
208  m2Dip, m2Rec, pT2damp, pTbegRef, pdfScale2;
209 
210  // Bookkeeping of enhanced actual or trial emissions (see EPJC (2013) 73).
211  bool doTrialNow, canEnhanceEmission, canEnhanceTrial, canEnhanceET;
212  string splittingNameNow, splittingNameSel;
213  map< double, pair<string,double> > enhanceFactors;
214  void storeEnhanceFactor(double pT2, string name, double enhanceFactorIn)
215  { enhanceFactors.insert(make_pair(pT2,make_pair(name,enhanceFactorIn)));}
216 
217  // List of emissions in different sides in different systems:
218  vector<int> nRadA,nRadB;
219 
220  // All dipole ends
221  vector<SpaceDipoleEnd> dipEnd;
222 
223  // List of 2 -> 2 momenta for external weak setup.
224  vector<Vec4> weakMomenta;
225 
226  // Pointers to the current and hardest (so far) dipole ends.
227  int iDipNow, iSysNow;
228  SpaceDipoleEnd* dipEndNow;
229  int iDipSel;
230  SpaceDipoleEnd* dipEndSel;
231 
232  // Evolve a QCD dipole end.
233  void pT2nextQCD( double pT2begDip, double pT2endDip);
234 
235  // Evolve a QCD and QED dipole end near heavy quark threshold region.
236  void pT2nearThreshold( BeamParticle& beam, double m2Massive,
237  double m2Threshold, double xMaxAbs, double zMinAbs,
238  double zMaxMassive, int iColPartner);
239 
240  // Evolve a QED dipole end.
241  void pT2nextQED( double pT2begDip, double pT2endDip);
242 
243  // Evolve a Weak dipole end.
244  void pT2nextWeak( double pT2begDip, double pT2endDip);
245 
246  // Find class of ME correction.
247  int findMEtype( int iSys, Event& event, bool weakRadiation = false);
248 
249  // Provide maximum of expected ME weight; for preweighting of evolution.
250  double calcMEmax( int MEtype, int idMother, int idDaughterIn);
251 
252  // Provide actual ME weight for current branching.
253  double calcMEcorr(int MEtype, int idMother, int idDaughterIn, double M2,
254  double z, double Q2,double m2Sister);
255 
256  // Provide actual ME weight for t-channel weak emissions.
257  double calcMEcorrWeak(int MEtype, double m2, double z,
258  double pT2, Vec4 pMother, Vec4 pB, Vec4 pDaughter,
259  Vec4 pB0, Vec4 p1, Vec4 p2, Vec4 pSister);
260 
261  // Find coefficient of azimuthal asymmetry from gluon polarization.
262  void findAsymPol( Event& event, SpaceDipoleEnd* dip);
263 
264  // Find a possible colour partner in the case of dipole recoil.
265  int findColPartner(Event& event, int iSideA, int iSideB, int iSystem);
266 
267  // Calculate uncertainty-band weights for accepted/rejected trial branching.
268  void calcUncertainties(bool accept, double pAcceptIn, double pT20in,
269  double enhance, SpaceDipoleEnd* dip, Particle* motherPtr,
270  Particle* sisterPtr);
271 
272 };
273 
274 //==========================================================================
275 
276 } // end namespace Pythia8
277 
278 #endif // Pythia8_SimpleSpaceShower_H
double pAccept
Properties needed for the evaluation of parameter variations.
Definition: SimpleSpaceShower.h:66
Definition: StandardModel.h:23
The Event class holds all info on the generated event.
Definition: Event.h:453
void store(int idDaughterIn, int idMotherIn, int idSisterIn, double x1In, double x2In, double m2DipIn, double pT2In, double zIn, double xMoIn, double Q2In, double mSisterIn, double m2SisterIn, double pT2corrIn, int iColPartnerIn, double m2IFIn, double mColPartnerIn)
Store values for trial emission.
Definition: SimpleSpaceShower.h:43
Definition: BeamParticle.h:133
int nBranch
Properties specific to current trial emission.
Definition: SimpleSpaceShower.h:61
Data on radiating dipole ends, only used inside SimpleSpaceShower.
Definition: SimpleSpaceShower.h:22
virtual bool wasGamma2qqbar() override
Tell if latest scattering was a gamma->qqbar.
Definition: SimpleSpaceShower.h:142
Definition: StandardModel.h:106
virtual bool doRestart() const override
Flag for failure in branch(...) that will force a retry of parton level.
Definition: SimpleSpaceShower.h:139
virtual int system() const override
Tell which system was the last processed one.
Definition: SimpleSpaceShower.h:148
Definition: Event.h:32
Definition: SimpleWeakShowerMEs.h:26
virtual bool getHasWeaklyRadiated() override
Tell whether ISR has done a weak emission.
Definition: SimpleSpaceShower.h:145
double m2(const Vec4 &v1)
The squared invariant mass of one or more four-vectors.
Definition: Basics.cc:605
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
virtual double enhancePTmax() const override
Potential enhancement factor of pTmax scale for hardest emission.
Definition: SimpleSpaceShower.h:151
SpaceDipoleEnd(int systemIn=0, int sideIn=0, int iRadiatorIn=0, int iRecoilerIn=0, double pTmaxIn=0., int colTypeIn=0, int chgTypeIn=0, int weakTypeIn=0, int MEtypeIn=0, bool normalRecoilIn=true, int weakPolIn=0, int iColPartnerIn=0, int idColPartnerIn=0)
Constructor.
Definition: SimpleSpaceShower.h:27
virtual ~SimpleSpaceShower() override
Destructor.
Definition: SimpleSpaceShower.h:104
int system
Basic properties related to evolution and matrix element corrections.
Definition: SimpleSpaceShower.h:54
The SpaceShower class does spacelike showers.
Definition: SpaceShower.h:33
SimpleSpaceShower()
Constructor.
Definition: SimpleSpaceShower.h:79
Definition: Basics.h:32
The SimpleSpaceShower class does spacelike showers.
Definition: SimpleSpaceShower.h:74