PYTHIA  8.312
SigmaTotal.h
1 // SigmaTotal.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 classes for cross section parametrizations.
7 // SigmaTotAux : base class for the different parametrizations.
8 // SigmaTotal : top-level class, making use of the classes below.
9 // SigmaTotOwn : user-set values.
10 // SigmaSaSDL : Schuler and Sjostrand, based on Donnachie and Landshoff.
11 // SigmaMBR : Min Bias Rockefeller.
12 // SigmaABMST : Appleby, Barlow, Molson, Serluca and Toader.
13 // SigmaRPP : Review of Particle Physics 2014.
14 
15 #ifndef Pythia8_SigmaTotal_H
16 #define Pythia8_SigmaTotal_H
17 
18 #include "Pythia8/Info.h"
19 #include "Pythia8/ParticleData.h"
20 #include "Pythia8/PhysicsBase.h"
21 #include "Pythia8/PythiaStdlib.h"
22 #include "Pythia8/PythiaComplex.h"
23 #include "Pythia8/Settings.h"
24 
25 namespace Pythia8 {
26 
27 //==========================================================================
28 
29 // The SigmaTotAux is base class for the different parametrizations
30 // made accessible via SigmaTotal. Many variables are public, since
31 // they can only be accessed via the SigmaTotal methods anyway.
32 
33 class SigmaTotAux {
34 
35 public:
36 
37  // Constructor.
38  SigmaTotAux() : isExpEl(), hasCou(), sigTot(), rhoOwn(), sigEl(), bEl(),
39  sigTotCou(), sigElCou(), sigXB(), sigAX(), sigXX(), sigAXB(), idA(),
40  idB(), tryCoulomb(), chgSgn(), tAbsMin(), lambda(), phaseCst(),
41  particleDataPtr(), rndmPtr() {}
42 
43  // Destructor.
44  virtual ~SigmaTotAux() {};
45 
46  // Store pointers and initialize data members.
47  virtual void init(Info*) = 0;
48 
49  // Calculate total cros section only. Only implemented for SigmaSaSDL.
50  virtual bool calcTot( int, int, double) {return true;}
51 
52  // Calculate integrated total/elastic cross sections.
53  // Usage: calcTotEl( idAin, idBin, sIn, mAin, mBin).
54  virtual bool calcTotEl( int , int , double , double , double ) {return true;}
55 
56  // Store total and elastic cross section properties.
57  bool isExpEl, hasCou;
58  double sigTot, rhoOwn, sigEl, bEl, sigTotCou, sigElCou;
59 
60  // Differential elastic cross section, d(sigma_el) / dt.
61  virtual double dsigmaEl( double , bool = false, bool = false) {return 0.;}
62 
63  // Calculate integrated diffractive cross sections.
64  // Usage: calcDiff( idAin, idBin, sIn, mAin, mBin).
65  virtual bool calcDiff( int , int , double , double , double ) {return true;}
66 
67  // Store diffractive cross sections. Possibly also sigmaND.
68  double sigXB, sigAX, sigXX, sigAXB, sigNDtmp;
69 
70  // Differential single diffractive cross section,
71  // xi * d(sigma_SD) / (dxi dt).
72  virtual double dsigmaSD( double , double, bool = true, int = 0) {return 0.;}
73 
74  // Possibility to separate xi and t choices for diffraction.
75  virtual bool splitDiff() {return false;}
76 
77  // Differential double diffractive cross section,
78  // xi1 * xi2 * d(sigma_DD) / (dxi1 dxi2 dt).
79  virtual double dsigmaDD( double , double , double , int = 0) {return 0.;}
80 
81  // Differential central diffractive cross section,
82  // xi1 * xi2 * d(sigma_CD) / (dxi1 dxi2 dt1 dt2).
83  virtual double dsigmaCD( double , double , double , double, int ) {
84  return 0.;}
85 
86  // Minimal central diffractive mass.
87  virtual double mMinCD() {return 1.;}
88 
89  // Standard methods to find t range of a 2 -> 2 process
90  // and to check whether a given t value is in that range.
91  pair<double,double> tRange( double sIn, double s1In, double s2In,
92  double s3In, double s4In) {
93  double lambda12 = pow2( sIn - s1In - s2In) - 4. * s1In * s2In;
94  double lambda34 = pow2( sIn - s3In - s4In) - 4. * s3In * s4In;
95  if (lambda12 < 0. || lambda34 < 0.) return make_pair( 0., 0.);
96  double tLow = -0.5 * (sIn - (s1In + s2In + s3In + s4In) + (s1In - s2In)
97  * (s3In - s4In) / sIn + sqrtpos(lambda12 * lambda34) / sIn);
98  double tUpp = ( (s3In - s1In) * (s4In - s2In) + (s1In + s4In - s2In - s3In)
99  * (s1In * s4In - s2In * s3In) / sIn ) / tLow;
100  return make_pair( tLow, tUpp); }
101  bool tInRange( double tIn, double sIn, double s1In, double s2In,
102  double s3In, double s4In) {
103  pair<double, double> tRng = tRange( sIn, s1In, s2In, s3In, s4In);
104  return (tIn > tRng.first && tIn < tRng.second); }
105 
106  // Commonly used proton form factor.
107  double pFormFac(double tIn) {return (4. * SPROTON - 2.79 * tIn)
108  / ((4. * SPROTON - tIn) * pow2(1. - tIn / 0.71)); }
109 
110 protected:
111 
112  // Constants: could only be changed in the code itself.
113  static const int NPOINTS;
114  static const double ALPHAEM, CONVERTEL, MPROTON, SPROTON, MPION, SPION,
115  GAMMAEUL, TABSREF, TABSMAX, MINSLOPEEL;
116 
117  // Initialization data, normally only set once.
118  int idA, idB;
119 
120  // Add Coulomb corrections to the elastic cross section.
122  double chgSgn, tAbsMin, lambda, phaseCst;
123  virtual bool initCoulomb(Settings& settings,
124  ParticleData* particleDataPtrIn);
125  virtual bool addCoulomb();
126  virtual double dsigmaElCoulomb(double t);
127 
128  // Pointer to the particle data table.
130 
131  // Pointer to the random number generator.
133 
134 };
135 
136 //==========================================================================
137 
138 // The SigmaTotal class contains parametrizations of total, elastic and
139 // diffractive cross sections, and of the respective slope parameter.
140 
141 class SigmaTotal : public PhysicsBase {
142 
143 public:
144 
145  // Constructor.
146  SigmaTotal() : isCalc(false), ispp(), modeTotEl(), modeTotElNow(),
147  modeDiff(), modeDiffNow(), idAbsA(), idAbsB(), s(), sigND(),
148  sigTotElPtr(nullptr), sigDiffPtr(nullptr) {};
149 
150  // Destructor.
151  virtual ~SigmaTotal() { if (sigTotElPtr) delete sigTotElPtr;
152  if (sigDiffPtr) delete sigDiffPtr; }
153 
154  // Store pointers and initialize data members.
155  void init();
156 
157  // Calculate, or recalculate for new beams or new energy.
158  bool calc( int idA, int idB, double eCM);
159 
160  // Confirm that initialization worked.
161  bool hasSigmaTot() {return isCalc;}
162 
163  // Total and integrated elastic cross sections.
164  double sigmaTot() {return sigTotElPtr->sigTotCou;}
165  double rho() {return sigTotElPtr->rhoOwn;}
166  double sigmaEl() {return sigTotElPtr->sigElCou;}
167  bool bElIsExp() {return sigTotElPtr->isExpEl;}
168  double bSlopeEl() {return sigTotElPtr->bEl;}
169  bool hasCoulomb() {return sigTotElPtr->hasCou;}
170 
171  // Total and elastic cross section.
172  bool calcTotEl( int idAin, int idBin, double sIn, double mAin, double mBin) {
173  return sigTotElPtr->calcTotEl( idAin, idBin, sIn, mAin, mBin); }
174 
175  // Differential elastic cross section.
176  double dsigmaEl( double t, bool useCoulomb = false,
177  bool onlyPomerons = false) {
178  return sigTotElPtr->dsigmaEl( t, useCoulomb, onlyPomerons); }
179 
180  // Integrated diffractive cross sections.
181  double sigmaXB() const {return sigDiffPtr->sigXB;}
182  double sigmaAX() const {return sigDiffPtr->sigAX;}
183  double sigmaXX() const {return sigDiffPtr->sigXX;}
184  double sigmaAXB() const {return sigDiffPtr->sigAXB;}
185  double sigmaND() const {return sigND;}
186 
187  // Differential single diffractive cross section.
188  double dsigmaSD( double xi, double t, bool isXB = true, int step = 0) {
189  return sigDiffPtr->dsigmaSD( xi, t, isXB, step); }
190 
191  // Possibility to separate xi and t choices for diffraction.
192  virtual bool splitDiff() {return sigDiffPtr->splitDiff();}
193 
194  // Differential double diffractive cross section.
195  double dsigmaDD( double xi1, double xi2, double t, int step = 0) {
196  return sigDiffPtr->dsigmaDD( xi1, xi2, t, step); }
197 
198  // Differential central diffractive cross section.
199  double dsigmaCD( double xi1, double xi2, double t1, double t2, int step = 0)
200  { return sigDiffPtr->dsigmaCD( xi1, xi2, t1, t2, step); }
201 
202  // Minimal central diffractive mass.
203  double mMinCD() {return sigDiffPtr->mMinCD();}
204 
205  // Sample the VMD states for resolved photons.
206  void chooseVMDstates(int idA, int idB, double eCM, int processCode);
207 
208  // Standard methods to find t range of a 2 -> 2 process
209  // and to check whether a given t value is in that range.
210  pair<double,double> tRange( double sIn, double s1In, double s2In,
211  double s3In, double s4In) {
212  return sigDiffPtr->tRange( sIn, s1In, s2In, s3In, s4In); }
213  bool tInRange( double tIn, double sIn, double s1In, double s2In,
214  double s3In, double s4In) {
215  return sigDiffPtr->tInRange( tIn, sIn, s1In, s2In, s3In, s4In); }
216 
217 private:
218 
219  // Constants: could only be changed in the code itself.
220  static const double MMIN;
221 
222  // Initialization data, normally only set once.
223  bool isCalc, ispp;
224 
225  int modeTotEl, modeTotElNow, modeDiff, modeDiffNow, idAbsA, idAbsB,
226  idAOld, idBOld, modeTotElOld, modeDiffOld;
227  double s, sigND, eCMOld;
228 
229  // Pointer to class that handles total and elastic cross sections.
230  SigmaTotAux* sigTotElPtr;
231 
232  // Pointer to class that handles diffractive cross sections.
233  SigmaTotAux* sigDiffPtr;
234 
235 };
236 
237 //==========================================================================
238 
239 // The SigmaTotOwn class parametrizes total, elastic and diffractive
240 // cross sections by user settings.
241 
242 class SigmaTotOwn : public SigmaTotAux {
243 
244 public:
245 
246  // Constructor.
247  SigmaTotOwn() : dampenGap(), pomFlux(), s(), a0(), ap(), b0(), A1(), A2(),
248  A3(), a1(), a2(), a3(), bMinDD(), ygap(), ypow(), expPygap(), mMinCDnow(),
249  wtNow(), yNow(), yNow1(), yNow2(), b(), b1(), b2(), Q(), Q1(),
250  Q2() {};
251 
252  // Store pointers and initialize data members.
253  virtual void init(Info* infoPtrIn);
254 
255  // Calculate integrated total/elastic cross sections.
256  virtual bool calcTotEl( int idAin, int idBin, double , double , double);
257 
258  // Differential elastic cross section.
259  virtual double dsigmaEl( double t, bool useCoulomb = false, bool = false);
260 
261  // Calculate integrated diffractive cross sections.
262  virtual bool calcDiff( int , int , double sIn, double , double ) {
263  s = sIn; return true;}
264 
265  // Differential single diffractive cross section.
266  virtual double dsigmaSD( double xi, double t, bool = true, int = 0);
267 
268  // Differential double diffractive cross section.
269  virtual double dsigmaDD( double xi1, double xi2, double t, int = 0);
270 
271  // Differential central diffractive cross section.
272  virtual double dsigmaCD( double xi1, double xi2, double t1, double t2,
273  int = 0);
274 
275  // Minimal central diffractive mass.
276  virtual double mMinCD() {return mMinCDnow;}
277 
278 private:
279 
280  // Initialization data, normally only set once.
281  bool dampenGap;
282  int pomFlux;
283  double s, a0, ap, b0, A1, A2, A3, a1, a2, a3, bMinDD, ygap, ypow, expPygap,
284  mMinCDnow, wtNow, yNow, yNow1, yNow2, b, b1, b2, Q, Q1, Q2;
285 
286 };
287 
288 //==========================================================================
289 
290 // The SigmaSaSDL class parametrizes total, elastic and diffractive
291 // cross sections according to Schuler and Sjostrand, starting from
292 // Donnachie and Landshoff.
293 
294 class SigmaSaSDL : public SigmaTotAux {
295 
296 public:
297 
298  // Constructor.
299  SigmaSaSDL() : doDampen(), zeroAXB(), swapped(), sameSign(), idAbsA(),
300  idAbsB(), iProc(), iHadA(), iHadB(), iHadAtmp(), iHadBtmp(), iProcVP(),
301  iProcVV(), s(), mA(), mB(), bA(), bB(), maxXBOwn(), maxAXOwn(), maxXXOwn(),
302  maxAXBOwn(), epsSaS(), sigmaPomP(), mPomP(), pPomP(), sigAXB2TeV(),
303  mMin0(), cRes(), mRes0(), mMinCDnow(), alP2(), s0(), mMinXB(), mMinAX(),
304  mMinAXB(), mResXB(), mResAX(), sResXB(), sResAX(), wtNow(), mAtmp(),
305  mBtmp(), multVP(), multVV(), infoPtr() {};
306 
307  // Store pointers and initialize data members.
308  virtual void init(Info* infoPtrIn);
309 
310  // Fast method uniquely to return total cross section.
311  double sigmaTotal( int idAin, int idBin, double sIn, double mAin,
312  double mBin);
313 
314  // Calculate integrated total/elastic cross sections.
315  virtual bool calcTotEl( int idAin, int idBin, double sIn, double mAin,
316  double mBin);
317 
318  // Differential elastic cross section.
319  virtual double dsigmaEl( double t, bool useCoulomb = false, bool = false);
320 
321  // Calculate integrated diffractive cross sections.
322  virtual bool calcDiff( int idAin, int idBin, double sIn, double mAin,
323  double mBin);
324 
325  // Differential single diffractive cross section.
326  virtual double dsigmaSD( double xi, double t, bool isXB, int = 0);
327 
328  // Differential double diffractive cross section.
329  virtual double dsigmaDD( double xi1, double xi2, double t, int = 0);
330 
331  // Differential central diffractive cross section.
332  virtual double dsigmaCD( double xi1, double xi2, double t1, double t2,
333  int = 0);
334 
335  // Minimal central diffractive mass.
336  virtual double mMinCD() {return mMinCDnow;}
337 
338 private:
339 
340  // Constants: could only be changed in the code itself.
341  static const int IHADATABLE[], IHADBTABLE[], ISDTABLE[], IDDTABLE[], NVMD;
342  static const double EPSILON, ETA, X[], Y[], BETA0[], BHAD[], ALPHAPRIME,
343  CONVERTSD, CONVERTDD, VMDMASS[4], GAMMAFAC[4],
344  CSD[24][8], CDD[22][9], DIFFTHR, DIFFMULT;
345 
346  // Initialization data, normally only set once, and result of calculation.
347  bool doDampen, zeroAXB, swapped, sameSign;
348  int idAbsA, idAbsB, iProc, iHadA, iHadB, iHadAtmp[4],
349  iHadBtmp[4], iProcVP[4], iProcVV[4][4];
350  double s, mA, mB, bA, bB, maxXBOwn, maxAXOwn, maxXXOwn, maxAXBOwn, epsSaS,
351  sigmaPomP, mPomP, pPomP, sigAXB2TeV, mMin0, cRes, mRes0, mMinCDnow,
352  alP2, s0, mMinXB, mMinAX, mMinAXB, mResXB, mResAX, sResXB,
353  sResAX, wtNow, mAtmp[4], mBtmp[4], multVP[4], multVV[4][4];
354 
355  // Find combination of incoming beams.
356  bool findBeamComb( int idAin, int idBin, double mAin, double mBin);
357 
358  // Pointer to various information on the generation.
359  Info* infoPtr;
360 
361 };
362 
363 //==========================================================================
364 
365 // The SigmaMBR class parametrizes total and elastic cross sections
366 // according to the Minimum Bias Rockefeller (MBR) model.
367 
368 class SigmaMBR : public SigmaTotAux {
369 
370 public:
371 
372  // Constructor.
373  SigmaMBR() : s(), sigSD(), sigDD(), sigCD(), eps(), alph(), beta0gev(),
374  beta0mb(), sigma0mb(), sigma0gev(), m2min(), dyminSDflux(),
375  dyminDDflux(), dyminCDflux(), dyminSD(), dyminDD(), dyminCD(),
376  dyminSigSD(), dyminSigDD(), dyminSigCD(), a1(), a2(), b1(), b2(),
377  sdpmax(), ddpmax(), dpepmax() {};
378 
379  // Initialize data members.
380  virtual void init(Info* infoPtrIn);
381 
382  // Calculate integrated total/elastic cross sections.
383  virtual bool calcTotEl( int idAin, int idBin, double sIn, double , double );
384 
385  // Differential elastic cross section.
386  virtual double dsigmaEl(double t, bool useCoulomb = false, bool = false);
387 
388  // Calculate integrated diffractive cross sections.
389  virtual bool calcDiff( int , int , double sIn, double , double );
390 
391  // In MBR choice of xi and t values are separated.
392  virtual bool splitDiff() {return true;}
393 
394  // Differential single diffractive cross section.
395  virtual double dsigmaSD( double xi, double t, bool = true, int step = 0);
396 
397  // Differential double diffractive cross section.
398  virtual double dsigmaDD( double xi1, double xi2, double t, int step = 0);
399 
400  // Differential central diffractive cross section.
401  virtual double dsigmaCD( double xi1, double xi2, double t1, double t2,
402  int step = 0);
403 
404  // Minimal central diffractive mass.
405  virtual double mMinCD() {return sqrt(m2min);}
406 
407 private:
408 
409  // Integration of MBR cross sections and form factor approximation.
410  static const int NINTEG, NINTEG2;
411  static const double FFA1, FFA2, FFB1, FFB2;
412 
413  // Parameters of MBR model.
414  double s, sigSD, sigDD, sigCD, eps, alph, beta0gev,
415  beta0mb, sigma0mb, sigma0gev, m2min, dyminSDflux, dyminDDflux,
416  dyminCDflux, dyminSD, dyminDD, dyminCD, dyminSigSD, dyminSigDD,
417  dyminSigCD, a1, a2, b1, b2, sdpmax, ddpmax, dpepmax;
418 
419  // The error function erf(x) should normally be in your math library,
420  // but if not uncomment this simple parametrization by Sergei Winitzki.
421  //double erf(double x) { double x2 = x * x; double kx2 = 0.147 * x2;
422  // double tmp = sqrt(1. - exp(-x2 * (4./M_PI + kx2) / (1. + kx2)));
423  // return ((x >= 0.) ? tmp : -tmp); }
424 
425 };
426 
427 //==========================================================================
428 
429 // The SigmaABMST class parametrizes total and elastic cross sections
430 // according to Appleby, Barlow, Molson, Serluca and Toader (ABMST).
431 
432 class SigmaABMST : public SigmaTotAux {
433 
434 public:
435 
436  // Constructor.
437  SigmaABMST() : ispp(), dampenGap(), useBMin(), modeSD(), modeDD(), modeCD(),
438  s(), facEl(), m2minp(), m2minm(), alp0(), alpt(), s0(), c0(), ygap(),
439  ypow(), expPygap(), multSD(), powSD(), multDD(), powDD(), multCD(),
440  powCD(), mMinCDnow(), bMinSD(), bMinDD(), bMinCD() {};
441 
442  // Initialize data members.
443  virtual void init(Info* infoPtrIn);
444 
445  // Calculate integrated total/elastic cross sections.
446  virtual bool calcTotEl( int idAin, int idBin, double sIn, double , double );
447 
448  // Differential elastic cross section.
449  virtual double dsigmaEl( double t, bool useCoulomb = false,
450  bool onlyPomerons = false) {
451  return facEl * pow2(abs(amplitude( t, useCoulomb, onlyPomerons)));}
452 
453  // Calculate integrated diffractive cross sections.
454  virtual bool calcDiff( int idAin, int idBin, double sIn, double , double );
455 
456  // Differential single diffractive cross section.
457  virtual double dsigmaSD( double xi, double t, bool = true , int = 0);
458 
459  // Differential double diffractive cross section.
460  virtual double dsigmaDD( double xi1, double xi2, double t, int = 0);
461 
462  // Differential central diffractive cross section.
463  virtual double dsigmaCD( double xi1, double xi2, double t1, double t2,
464  int = 0);
465 
466  // Minimal central diffractive mass.
467  virtual double mMinCD() {return mMinCDnow;}
468 
469 private:
470 
471  // Constants: could only be changed in the code itself.
472  static const bool MCINTDD;
473  static const int NPOINTSTSD, NPOINTSTDD, NPOINTMCDD, NPOINTMCCD;
474  static const double EPSI[], ALPP[], NORM[], SLOPE[], FRACS[], TRIG[],
475  LAM2P, BAPPR[], LAM2FF, MRES[4], WRES[4], CRES[4],
476  AFAC[4], BFAC[4], CFAC[4], PPP[4], EPS[2], APR[2],
477  CPI[6], CNST[5], XIDIVSD, DXIRAWSD, DLNXIRAWSD,
478  XIDIVDD, DXIRAWDD, DLNXIRAWDD, BMCINTDD, BMCINTCD;
479 
480  // Initialization data, normally only set once.
481  bool ispp, dampenGap, useBMin;
482  int modeSD, modeDD, modeCD;
483  double s, facEl, m2minp, m2minm, alp0[2], alpt[3], s0, c0,
484  ygap, ypow, expPygap, multSD, powSD, multDD, powDD,
485  multCD, powCD, mMinCDnow, bMinSD, bMinDD, bMinCD;
486 
487  // The scattering amplitude, from which cross sections are derived.
488  complex amplitude( double t, bool useCoulomb = false,
489  bool onlyPomerons = false) ;
490 
491  // Auxiliary method for repetitive part of amplitude.
492  complex sModAlp( double sMod, double alpha) {
493  return exp(complex( 0., -0.5 * M_PI * alpha)) * pow( sMod, alpha); }
494 
495  // Core differential single diffractive cross section.
496  virtual double dsigmaSDcore(double xi, double t);
497 
498  // xi * d(sigma_SD) / (dxi dt) integrated over [t_min, t_max].
499  double dsigmaSDintT( double xi, double tMinIn, double tMaxIn);
500 
501  // d(sigma_SD) / (dxi dt) integrated over [xi_min, xi_max] * [t_min, t_max].
502  double dsigmaSDintXiT( double xiMinIn, double xiMaxIn, double tMinIn,
503  double tMaxIn );
504 
505  // d(sigma_DD) / (dxi1 dxi2 dt) integrated with Monte Carlo sampling.
506  double dsigmaDDintMC();
507 
508  // xi1 * xi2 * d(sigma_DD) / (dxi1 dxi2 dt) integrated over t.
509  double dsigmaDDintT( double xi1, double xi2, double tMinIn, double tMaxIn);
510 
511  // xi1 * d(sigma_DD) / (dxi1 dxi2 dt) integrated over xi2 and t.
512  double dsigmaDDintXi2T( double xi1, double xi2MinIn, double xi2MaxIn,
513  double tMinIn, double tMaxIn);
514 
515  // d(sigma_DD) / (dxi1 dxi2 dt) integrated over xi1, xi2 and t.
516  double dsigmaDDintXi1Xi2T( double xi1MinIn, double xi1MaxIn,
517  double xi2MinIn, double xi2MaxIn, double tMinIn, double tMaxIn);
518 
519  // d(sigma_CD) / (dxi1 dxi2 dt1 dt2) integrated with Monte Carlo sampling.
520  double dsigmaCDintMC();
521 
522 };
523 
524 //==========================================================================
525 
526 // The SigmaRPP class parametrizes total and elastic cross sections
527 // according to the fit in Review of Particle Physics 2014.
528 
529 class SigmaRPP : public SigmaTotAux {
530 
531 public:
532 
533  // Constructor.
534  SigmaRPP() : ispp(), s(), facEl() {};
535 
536  // Initialize data members.
537  virtual void init(Info* infoPtrIn) {
538  tryCoulomb = infoPtrIn->settingsPtr->flag("SigmaElastic:Coulomb");
539  tAbsMin = infoPtrIn->settingsPtr->parm("SigmaElastic:tAbsMin"); }
540 
541  // Calculate integrated total/elastic cross sections.
542  virtual bool calcTotEl( int idAin, int idBin, double sIn, double , double );
543 
544  // Differential elastic cross section.
545  virtual double dsigmaEl( double t, bool useCoulomb = false, bool = false) {
546  return facEl * pow2(abs(amplitude( t, useCoulomb))); }
547 
548 private:
549 
550  // Constants: could only be changed in the code itself.
551  static const double EPS1[], ALPP[], NORM[], BRPP[], KRPP[], LAM2FF;
552 
553  // Initialization data, normally only set once, and result of calculation.
554  bool ispp;
555  double s, facEl;
556 
557  // The scattering amplitude, from which cross sections are derived.
558  complex amplitude( double t, bool useCoulomb = false) ;
559 
560  // Auxiliary method for repetitive part of amplitude.
561  complex sModAlp( double sMod, double alpha) {
562  return exp(complex( 0., -0.5 * M_PI * alpha)) * pow( sMod, alpha); }
563 
564  // Auxiliary methods for Bessel J0 and J1 functions.
565  complex besJ0( complex x);
566  complex besJ1( complex x);
567 
568 };
569 
570 //==========================================================================
571 
572 } // end namespace Pythia8
573 
574 #endif // Pythia8_SigmaTotal_H
constexpr double pow2(const double &x)
Powers of small integers - for balance speed/code clarity.
Definition: PythiaStdlib.h:182
SigmaMBR()
Constructor.
Definition: SigmaTotal.h:373
std::complex< double > complex
Convenient typedef for double precision complex numbers.
Definition: PythiaComplex.h:17
static const double ALPHAEM
alpha_em(0).
Definition: SigmaTotal.h:114
bool isExpEl
Store total and elastic cross section properties.
Definition: SigmaTotal.h:57
bool flag(string keyIn)
Give back current value, with check that key exists.
Definition: Settings.cc:1478
SigmaTotal()
Constructor.
Definition: SigmaTotal.h:146
Definition: PhysicsBase.h:27
Definition: Info.h:45
virtual bool splitDiff()
Possibility to separate xi and t choices for diffraction.
Definition: SigmaTotal.h:75
bool tryCoulomb
Add Coulomb corrections to the elastic cross section.
Definition: SigmaTotal.h:121
pair< double, double > tRange(double sIn, double s1In, double s2In, double s3In, double s4In)
Definition: SigmaTotal.h:210
SigmaABMST()
Constructor.
Definition: SigmaTotal.h:437
double dsigmaDD(double xi1, double xi2, double t, int step=0)
Differential double diffractive cross section.
Definition: SigmaTotal.h:195
Definition: SigmaTotal.h:141
virtual double dsigmaElCoulomb(double t)
Coulomb contribution to the differential elastic cross sections.
Definition: SigmaTotal.cc:121
virtual bool addCoulomb()
Add Coulomb corrections to the elastic and total cross sections.
Definition: SigmaTotal.cc:67
virtual bool calcTotEl(int, int, double, double, double)
Definition: SigmaTotal.h:54
virtual bool initCoulomb(Settings &settings, ParticleData *particleDataPtrIn)
Initialize Coulomb correction parameters.
Definition: SigmaTotal.cc:46
pair< double, double > tRange(double sIn, double s1In, double s2In, double s3In, double s4In)
Definition: SigmaTotal.h:91
double dsigmaCD(double xi1, double xi2, double t1, double t2, int step=0)
Differential central diffractive cross section.
Definition: SigmaTotal.h:199
Definition: SigmaTotal.h:529
virtual void init(Info *)=0
Store pointers and initialize data members.
SigmaSaSDL()
Constructor.
Definition: SigmaTotal.h:299
double mMinCD()
Minimal central diffractive mass.
Definition: SigmaTotal.h:203
int idA
Initialization data, normally only set once.
Definition: SigmaTotal.h:118
Definition: SigmaTotal.h:242
Definition: SigmaTotal.h:294
virtual double mMinCD()
Minimal central diffractive mass.
Definition: SigmaTotal.h:336
virtual double mMinCD()
Minimal central diffractive mass.
Definition: SigmaTotal.h:405
virtual double dsigmaEl(double, bool=false, bool=false)
Differential elastic cross section, d(sigma_el) / dt.
Definition: SigmaTotal.h:61
Definition: Basics.h:386
virtual double dsigmaEl(double t, bool useCoulomb=false, bool onlyPomerons=false)
Differential elastic cross section.
Definition: SigmaTotal.h:449
double sigXB
Store diffractive cross sections. Possibly also sigmaND.
Definition: SigmaTotal.h:68
double dsigmaEl(double t, bool useCoulomb=false, bool onlyPomerons=false)
Differential elastic cross section.
Definition: SigmaTotal.h:176
Definition: SigmaTotal.h:33
double pFormFac(double tIn)
Commonly used proton form factor.
Definition: SigmaTotal.h:107
virtual bool splitDiff()
In MBR choice of xi and t values are separated.
Definition: SigmaTotal.h:392
virtual double dsigmaDD(double, double, double, int=0)
Definition: SigmaTotal.h:79
double sigmaXB() const
Integrated diffractive cross sections.
Definition: SigmaTotal.h:181
virtual bool calcTot(int, int, double)
Calculate total cros section only. Only implemented for SigmaSaSDL.
Definition: SigmaTotal.h:50
virtual bool splitDiff()
Possibility to separate xi and t choices for diffraction.
Definition: SigmaTotal.h:192
static const double TABSREF
Reference scale for nominal definition of t slope.
Definition: SigmaTotal.h:114
Definition: SigmaTotal.h:432
bool hasSigmaTot()
Confirm that initialization worked.
Definition: SigmaTotal.h:161
ParticleData * particleDataPtr
Pointer to the particle data table.
Definition: SigmaTotal.h:129
virtual bool calcDiff(int, int, double sIn, double, double)
Calculate integrated diffractive cross sections.
Definition: SigmaTotal.h:262
SigmaTotAux()
Constructor.
Definition: SigmaTotal.h:38
SigmaRPP()
Constructor.
Definition: SigmaTotal.h:534
virtual double dsigmaEl(double t, bool useCoulomb=false, bool=false)
Differential elastic cross section.
Definition: SigmaTotal.h:545
Rndm * rndmPtr
Pointer to the random number generator.
Definition: SigmaTotal.h:132
static const double MPROTON
Proton and pion masses, and their squares. Euler&#39;s constant.
Definition: SigmaTotal.h:114
virtual double mMinCD()
Minimal central diffractive mass.
Definition: SigmaTotal.h:467
virtual double mMinCD()
Minimal central diffractive mass.
Definition: SigmaTotal.h:87
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
Definition: SigmaTotal.h:368
double dsigmaSD(double xi, double t, bool isXB=true, int step=0)
Differential single diffractive cross section.
Definition: SigmaTotal.h:188
double sigmaTot()
Total and integrated elastic cross sections.
Definition: SigmaTotal.h:164
virtual ~SigmaTotAux()
Destructor.
Definition: SigmaTotal.h:44
virtual void init(Info *infoPtrIn)
Initialize data members.
Definition: SigmaTotal.h:537
virtual double dsigmaCD(double, double, double, double, int)
Definition: SigmaTotal.h:83
virtual bool calcDiff(int, int, double, double, double)
Definition: SigmaTotal.h:65
This class holds a map of all ParticleDataEntries.
Definition: ParticleData.h:422
static const int NPOINTS
Constants: could only be changed in the code itself.
Definition: SigmaTotal.h:113
Settings * settingsPtr
Pointer to the settings database.
Definition: Info.h:80
virtual ~SigmaTotal()
Destructor.
Definition: SigmaTotal.h:151
virtual double mMinCD()
Minimal central diffractive mass.
Definition: SigmaTotal.h:276
SigmaTotOwn()
Constructor.
Definition: SigmaTotal.h:247
static const double CONVERTEL
Conversion coefficients = 1/(16pi) * (mb <-> GeV^2).
Definition: SigmaTotal.h:114
virtual double dsigmaSD(double, double, bool=true, int=0)
Definition: SigmaTotal.h:72
Definition: Settings.h:195
bool calcTotEl(int idAin, int idBin, double sIn, double mAin, double mBin)
Total and elastic cross section.
Definition: SigmaTotal.h:172
double sqrtpos(const double &x)
Avoid problem with negative square root argument (from roundoff).
Definition: PythiaStdlib.h:191