PYTHIA  8.312
VinciaQED.h
1 // VinciaQED.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2024 Peter Skands, 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 QED antenna-shower class and auxiliary
7 // classes. Main author is Rob Verheyen.
8 
9 #ifndef Pythia8_VinciaQED_H
10 #define Pythia8_VinciaQED_H
11 
12 // PYTHIA 8 headers.
13 #include "Pythia8/BeamParticle.h"
14 #include "Pythia8/Event.h"
15 #include "Pythia8/StandardModel.h"
16 #include "Pythia8/PartonSystems.h"
17 
18 // VINCIA headers.
19 #include "Pythia8/VinciaCommon.h"
20 #include "Pythia8/VinciaWeights.h"
21 
22 namespace Pythia8 {
23 
24 //==========================================================================
25 
26 // Class for QED emissions.
27 
29 
30 public:
31 
32  // Friends for internal private members.
33  friend class QEDemitSystem;
34 
35  // Constuctor.
36  QEDemitElemental() : rndmPtr(nullptr), partonSystemsPtr(nullptr), q2Sav(0.),
37  zetaSav(0.), phiSav(0.), sxjSav(0.), syjSav(0.), alpha(0.), c(0.),
38  hasTrial(false), x(0), y(0), idx(0), idy(0), mx2(0.), my2(0.),
39  ex(0.), ey(0.), m2Ant(0.), sAnt(0.), QQ(0.), isII(false), isIF(false),
40  isFF(false), isRF(false), isIA(true), isDip(false), shh(0.),
41  isInitPtr(false), isInit(false), verbose(1) {;}
42 
43  // Initialize the pointers.
44  void initPtr(Rndm* rndmPtrIn, PartonSystems* partonSystemsPtrIn);
45  // Initialize.
46  void init(Event &event, int xIn, int yIn, double shhIn,
47  double verboseIn);
48  // Initialize.
49  void init(Event &event, int xIn, vector<int> iRecoilIn, double shhIn,
50  double verboseIn);
51  // Generate a trial point.
52  double generateTrial(Event &event, double q2Start, double q2Low,
53  double alphaIn, double cIn);
54 
55 private:
56 
57  // Random pointer.
58  Rndm* rndmPtr{};
59 
60  // Parton system pointer.
61  PartonSystems* partonSystemsPtr{};
62 
63  // Trial variables.
64  double q2Sav, zetaSav, phiSav, sxjSav, syjSav, alpha, c;
65  bool hasTrial;
66 
67  // Particle indices.
68  int x, y;
69  // Recoiler indices.
70  vector<int> iRecoil;
71  // IDs.
72  int idx, idy;
73  // Particle masses.
74  double mx2, my2;
75  // Particle energies.
76  double ex, ey;
77  // Antenna invariant mass.
78  double m2Ant;
79  // Antenna dot product.
80  double sAnt;
81  // The negative of the product of charges.
82  double QQ;
83 
84  // Type switches.
85  bool isII, isIF, isFF, isRF, isIA, isDip;
86 
87  // Hadronic invariant mass.
88  double shh;
89 
90  // Initialization.
91  bool isInitPtr, isInit;
92  int verbose;
93 
94 };
95 
96 //==========================================================================
97 
98 // Base class for QED systems.
99 
100 class QEDsystem {
101 
102  public:
103 
104  // Constructor.
105  QEDsystem() : infoPtr(nullptr), partonSystemsPtr(nullptr),
106  particleDataPtr(nullptr), rndmPtr(nullptr), settingsPtr(nullptr),
107  loggerPtr(nullptr), vinComPtr(nullptr), isInitPtr(false), iSys(-1),
108  verbose(0), jNew(0), shat(0.) {;}
109 
110  // Destructor.
111  virtual ~QEDsystem() = default;
112 
113  // Initialize pointers.
114  void initPtr(Info* infoPtrIn, ParticleData* particleDataPtrIn,
115  PartonSystems* partonSystemsPtrIn, Rndm* rndmPtrIn,
116  Settings* settingsPtrIn, VinciaCommon* vinComPtrIn);
117 
118  // Initialise settings for current run.
119  virtual void init(BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn,
120  int verboseIn) = 0;
121  virtual void setVerbose(int verboseIn) { verbose = verboseIn; }
122  // Prepare a parton system for evolution.
123  virtual void prepare(int iSysIn, Event &event, double q2CutIn,
124  bool isBelowHadIn, vector<double> evolutionWindowsIn, AlphaEM alIn) = 0;
125  // Build parton system.
126  virtual void buildSystem(Event &event) = 0;
127  // Generate a trial scale.
128  virtual double q2Next(Event &event, double q2Start) = 0 ;
129  // Generate kinematics and check for veto of branching.
130  virtual bool acceptTrial(Event &event) = 0;
131  // Update the envent.
132  virtual void updateEvent(Event &event) = 0;
133  // Update the parton systems.
134  virtual void updatePartonSystems();
135  // Print information about the system.
136  virtual void print() = 0;
137 
138  // Methods to tell which type of brancher this is.
139  virtual bool isSplitting() {return false;};
140  virtual bool isInitial() {return false;};
141 
142  protected:
143 
144  // Pointers.
145  Info* infoPtr{};
146  PartonSystems* partonSystemsPtr{};
147  ParticleData* particleDataPtr{};
148  Rndm* rndmPtr{};
149  Settings* settingsPtr{};
150  Logger* loggerPtr{};
151  VinciaCommon* vinComPtr{};
152  bool isInitPtr;
153 
154  // Event system.
155  int iSys;
156  vector<Vec4> pNew;
157 
158  // Verbose setting.
159  int verbose;
160 
161  // Information for partonSystems.
162  int jNew;
163  map<int,int> iReplace;
164  double shat;
165 
166 };
167 
168 //==========================================================================
169 
170 // Class for a QED emission system.
171 
172 class QEDemitSystem : public QEDsystem {
173 
174 public:
175 
176  QEDemitSystem() : shh(-1.), cMat(0.), trialIsVec(false), beamAPtr(nullptr),
177  beamBPtr(nullptr), qedMode(-1), qedModeMPI(-1), useFullWkernel(false),
178  isBelowHad(false), emitBelowHad(false), q2Cut(-1.), isInit(false),
179  TINYPDF(-1.), kMapTypeFinal(0) {;}
180 
181  // Initialise settings for current run.
182  void init(BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn, int verboseIn)
183  override;
184  // Prepare a parton system for photon emission evolution
185  void prepare(int iSysIn, Event &event, double q2CutIn, bool isBelowHadIn,
186  vector<double> evolutionWindowsIn, AlphaEM alIn) override;
187  // Set up antenna pairing for incoherent mode.
188  void buildSystem(Event &event) override;
189  // Generate a trial scale.
190  double q2Next(Event &event, double q2Start) override;
191  // Generate kinematics and check veto.
192  bool acceptTrial(Event &event) override;
193  // Update the envent
194  void updateEvent(Event &event) override;
195  // Branching type.
196  bool isInitial() override {return eleTrial->isII || eleTrial->isIF;}
197  // Print the QED emit internal system.
198  void print() override;
199 
200  // Trial antenna function.
201  double aTrial(QEDemitElemental* ele, double sxj, double syj, double sxy);
202  // Physical antenna function.
203  double aPhys (QEDemitElemental* ele, double sxj, double syj, double sxy);
204  // Ratio between PDFs.
205  double pdfRatio(bool isA, double eOld, double eNew, int id, double Qt2);
206 
207  private:
208 
209  // Event system.
210  double shh;
211 
212  // Internal storage.
213  vector<vector<QEDemitElemental> > eleMat;
214  vector<int> iCoh;
215  double cMat;
216  vector<QEDemitElemental> eleVec;
217 
218  // AlphaEM.
219  AlphaEM al;
220 
221  // Evolution window.
222  vector<double> evolutionWindows;
223 
224  // Trial pointer.
225  QEDemitElemental* eleTrial{};
226  bool trialIsVec;
227 
228  // Pointers.
229  BeamParticle* beamAPtr{};
230  BeamParticle* beamBPtr{};
231 
232  // Settings.
233  int qedMode, qedModeMPI;
234  bool useFullWkernel, isBelowHad, emitBelowHad;
235  double q2Cut;
236 
237  // Initialization.
238  bool isInit;
239 
240  // PDF check.
241  double TINYPDF;
242 
243  // Recoil strategy.
244  int kMapTypeFinal;
245 
246  // Global recoil momenta.
247  Vec4 pRecSum;
248  vector<Vec4> pRec;
249  vector<int> iRec;
250 
251 };
252 
253 //==========================================================================
254 
255 // Class for trial QED splittings.
256 
258 
259 public:
260 
261  // Friends for internal private members.
262  friend class QEDsplitSystem;
263 
264  // Constructor.
265  QEDsplitElemental(Event &event, int iPhotIn, int iSpecIn):
266  iPhot(iPhotIn), iSpec(iSpecIn), ariWeight(0) {
267  m2Ant = max(VinciaConstants::PICO,
268  m2(event[iPhotIn], event[iSpecIn]));
269  sAnt = max(VinciaConstants::PICO,
270  2.*event[iPhotIn].p()*event[iSpecIn].p());
271  m2Spec = max(0., event[iSpecIn].m2());
272  }
273 
274  // Kallen function.
275  double getKallen() {return m2Ant/(m2Ant - m2Spec);}
276 
277 private:
278 
279  // Internal members.
280  int iPhot, iSpec;
281  double m2Spec, m2Ant, sAnt;
282  double ariWeight;
283 };
284 
285 //==========================================================================
286 
287 // Class for a QED splitting system.
288 
289 class QEDsplitSystem : public QEDsystem {
290 
291 public:
292 
293  // Constructor.
295  totIdWeight(-1.), hasTrial(false),
296  q2Trial(-1.), zTrial(-1.), phiTrial(-1.), idTrial(0), nQuark(-1),
297  nLepton(-1), q2Max(-1.), q2Cut(-1.), isBelowHad(false),
298  beamAPtr(nullptr), beamBPtr(nullptr), isInit(false), kMapTypeFinal(0) {;}
299 
300  // Initialize.
301  void init(BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn, int verboseIn)
302  override;
303  // Prepare list of final-state photons - with recoilers - for splittings.
304  void prepare(int iSysIn, Event &event, double q2CutIn, bool isBelowHadIn,
305  vector<double> evolutionWindowsIn, AlphaEM alIn) override;
306  // Build the splitting system.
307  void buildSystem(Event &event) override;
308  // Generate a scale for the system.
309  double q2Next(Event &event, double q2Start) override;
310  // Generate kinematics and check veto.
311  bool acceptTrial(Event &event) override;
312  // Update Event.
313  void updateEvent(Event &event) override;
314  // Branching type: isSplitting() = true.
315  bool isSplitting() override { return true;}
316  // Print the system.
317  void print() override;
318 
319 private:
320 
321  // AlphaEM.
322  AlphaEM al;
323 
324  // Evolution window.
325  vector<double> evolutionWindows;
326 
327  // Weights for splitting IDs.
328  vector<int> ids;
329  vector<double> idWeights;
330  double totIdWeight;
331 
332  // Antennae.
333  vector<QEDsplitElemental> eleVec;
334 
335  // Trial variables.
336  bool hasTrial;
337  double q2Trial, zTrial, phiTrial, idTrial;
338  QEDsplitElemental* eleTrial{};
339 
340  // Settings.
341  int nQuark, nLepton;
342  double q2Max, q2Cut;
343  bool isBelowHad;
344 
345  // Pointers.
346  BeamParticle* beamAPtr;
347  BeamParticle* beamBPtr;
348 
349  // Initialization.
350  bool isInit;
351 
352  // Recoil strategy.
353  int kMapTypeFinal;
354 
355 };
356 
357 //==========================================================================
358 
359 // Class for a QED conversion system.
360 
361 class QEDconvSystem : public QEDsystem {
362 
363 public:
364 
365  // Constructor.
366  QEDconvSystem() : totIdWeight(-1.), maxIdWeight(-1.), shh(-1.), s(-1.),
367  iA(-1), iB(-1), isAPhot(false), isBPhot(false), hasTrial(false),
368  iPhotTrial(-1), iSpecTrial(-1), q2Trial(-1.), zTrial(-1.), phiTrial(-1.),
369  idTrial(-1), nQuark(-1), q2Cut(-1.),
370  isBelowHad(false), beamAPtr(nullptr), beamBPtr(nullptr),isInit(false),
371  TINYPDF(-1.) {;}
372 
373  // Initialize.
374  void init(BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn, int verboseIn)
375  override;
376  // Prepare for backwards-evolution of photons.
377  void prepare(int iSysIn, Event &event, double q2CutIn, bool isBelowHadIn,
378  vector<double> evolutionWindowsIn, AlphaEM alIn) override;
379  // Build the system.
380  void buildSystem(Event &event) override;
381  // Generate a trial scale.
382  double q2Next(Event &event, double q2Start) override;
383  // Generate kinematics and check veto.
384  bool acceptTrial(Event &event) override;
385  // Update.
386  void updateEvent(Event &event) override;
387  // Branching type: isInitial() = true.
388  bool isInitial() override {return true;};
389  // Print.
390  void print() override;
391 
392 private:
393 
394  // Trial pdf ratios for conversion.
395  map<int,double> Rhat;
396 
397  // AlphaEM.
398  AlphaEM al;
399 
400  // Evolution window.
401  vector<double> evolutionWindows;
402 
403  // Weights for conversion IDs.
404  vector<int> ids;
405  vector<double> idWeights;
406  double totIdWeight, maxIdWeight;
407  double shh;
408 
409  // Antenna parameters.
410  double s;
411  int iA, iB;
412  bool isAPhot, isBPhot;
413 
414  // Trial variables.
415  bool hasTrial;
416  int iPhotTrial, iSpecTrial;
417  double q2Trial, zTrial, phiTrial, idTrial;
418 
419  // Settings.
420  int nQuark;
421  double q2Cut;
422  bool isBelowHad;
423 
424  // Pointers.
425  BeamParticle* beamAPtr;
426  BeamParticle* beamBPtr;
427 
428  // Initialization.
429  bool isInit;
430 
431  // PDF check.
432  double TINYPDF;
433 
434  // Global recoil momenta
435  vector<Vec4> pRec;
436  vector<int> iRec;
437 
438 };
439 
440 //==========================================================================
441 
442 // Base class for Vincia's QED and EW shower modules.
443 
445 
446 public:
447 
448  // Default constructor.
449  VinciaModule(): verbose(-1), isInitPtr(false), isInitSav(false) {;};
450 
451  // Default destructor.
452  virtual ~VinciaModule() = default;
453 
454  // Initialise pointers (called at construction time).
455  virtual void initPtr(Info* infoPtrIn, VinciaCommon* vinComPtrIn) = 0;
456 
457  //Some early initialisation (needed for EW shower).
458  virtual void load() {;}
459 
460  // Initialise settings for current run (called as part of Pythia::init()).
461  virtual void init(BeamParticle* beamAPtrIn = 0, BeamParticle* beamBPtrIn = 0)
462  = 0;
463  bool isInit() {return isInitSav;}
464 
465  // Select helicities for a system of particles.
466  virtual bool polarise(vector<Particle>&) {return false;}
467 
468  // Prepare to shower a system.
469  virtual bool prepare(int iSysIn, Event &event, bool isBelowHadIn) = 0;
470 
471  // Update shower system each time something has changed in event.
472  virtual void update(Event &event, int iSys) = 0;
473 
474  // Set verbosity level.
475  virtual void setVerbose(int verboseIn) {verbose = verboseIn;}
476 
477  // Clear everything (optionally for specific system).
478  virtual void clear(int iSys = -1) = 0;
479 
480  // Generate a trial scale.
481  virtual double q2Next(Event &event, double q2Start, double q2End) = 0;
482 
483  // Which system does the winner belong to?
484  virtual int sysWin() = 0;
485 
486  // Get information about the latest branching.
487  virtual bool lastIsSplitting() = 0;
488  virtual bool lastIsInitial() = 0;
489  virtual bool lastIsResonanceDecay() {return false;}
490 
491  // Generate kinematics and check veto.
492  virtual bool acceptTrial(Event &event) = 0;
493 
494  // Update event after branching accepted.
495  virtual void updateEvent(Event &event) = 0;
496 
497  // Update partonSystems after branching accepted.
498  virtual void updatePartonSystems(Event &event) = 0;
499 
500  // End scales.
501  virtual double q2minColoured() = 0;
502  virtual double q2min() = 0;
503 
504  // Get number of branchers / systems.
505  virtual unsigned int nBranchers() = 0;
506  virtual unsigned int nResDec() = 0;
507 
508  // Members.
509  BeamParticle* beamAPtr{};
510  BeamParticle* beamBPtr{};
511  Info* infoPtr{};
512  ParticleData* particleDataPtr{};
513  Logger* loggerPtr{};
514  PartonSystems* partonSystemsPtr{};
515  Rndm* rndmPtr{};
516  Settings* settingsPtr{};
517  VinciaCommon* vinComPtr{};
518 
519  protected:
520 
521  int verbose;
522  bool isInitPtr, isInitSav;
523 
524 };
525 
526 //==========================================================================
527 
528 // Class for performing QED showers.
529 
530 class VinciaQED : public VinciaModule {
531 
532 public:
533 
534  // Friends for internal private members.
535  friend class VinciaFSR;
536 
537  // Constructor.
538  VinciaQED() {;}
539 
540  // Initialise pointers (called at construction time).
541  void initPtr(Info* infoPtrIn, VinciaCommon* vinComPtrIn) override;
542  // Initialise settings for current run (called as part of Pythia::init()).
543  void init(BeamParticle* beamAPtrIn = 0, BeamParticle* beamBPtrIn = 0)
544  override;
545  // Prepare to shower a system.
546  bool prepare(int iSysIn, Event& event, bool isBelowHadIn) override;
547  // Update QED shower system(s) each time something has changed in event.
548  void update(Event& event, int iSys) override;
549  // Set or change verbosity level, and propagate to QED systems.
550  virtual void setVerbose(int verboseIn) override {
551  verbose = verboseIn;
552  emptyQEDemitSystem.setVerbose(verboseIn);
553  emptyQEDsplitSystem.setVerbose(verboseIn);
554  emptyQEDconvSystem.setVerbose(verboseIn);
555  }
556  // Clear everything, or specific system.
557  void clear(int iSys = -1) override {
558  if (iSys < 0) {emitSystems.clear(); splitSystems.clear();
559  convSystems.clear();}
560  else {emitSystems.erase(iSys); splitSystems.erase(iSys);
561  convSystems.erase(iSys);}
562  qedTrialSysPtr = nullptr;}
563 
564  // Generate a trial scale.
565  double q2Next(Event& event, double q2Start, double) override;
566  // Return the system window.
567  int sysWin() override {return iSysTrial;}
568  // Information about last branching.
569  bool lastIsSplitting() override {
570  if (qedTrialSysPtr != nullptr) return qedTrialSysPtr->isSplitting();
571  else return false;}
572  bool lastIsInitial() override {
573  if (qedTrialSysPtr != nullptr) return qedTrialSysPtr->isInitial();
574  else return false;}
575  // Generate kinematics and check veto.
576  bool acceptTrial(Event &event) override;
577  // Update event after branching accepted.
578  void updateEvent(Event &event) override;
579  // Update partonSystems after branching accepted.
580  void updatePartonSystems(Event &event) override;
581  // Return scales.
582  double q2minColoured() override {return q2minColouredSav;}
583  double q2min() override {return q2minSav;}
584 
585  // Getter for number of systems.
586  unsigned int nBranchers() override {
587  int sizeNow = emitSystems.size();
588  sizeNow = max(sizeNow, (int)splitSystems.size());
589  sizeNow = max(sizeNow, (int)convSystems.size());
590  return sizeNow;}
591 
592  // This module does not implement resonance decays.
593  unsigned int nResDec() override { return 0; }
594 
595 private:
596 
597  // Get Q2 for QED system.
598  template <class T>
599  void q2NextSystem(map<int, T>& QEDsystemList, Event& event, double q2Start);
600 
601  // Systems.
602  QEDemitSystem emptyQEDemitSystem;
603  QEDsplitSystem emptyQEDsplitSystem;
604  QEDconvSystem emptyQEDconvSystem;
605  map< int, QEDemitSystem> emitSystems;
606  map< int, QEDsplitSystem> splitSystems;
607  map< int, QEDconvSystem> convSystems;
608 
609  // Settings.
610  bool doQED, doEmission;
611  int nGammaToLepton, nGammaToQuark;
612  bool doConvertGamma, doConvertQuark;
613 
614  // Scales.
615  double q2minSav, q2minColouredSav;
616 
617  // Trial information
618  int iSysTrial;
619  double q2Trial;
620  QEDsystem* qedTrialSysPtr{};
621 
622  // AlphaEM.
623  AlphaEM al;
624 
625  // Evolution windows
626  vector<double> evolutionWindows;
627 
628 };
629 
630 //==========================================================================
631 
632 } // end namespace Pythia8
633 
634 #endif // Pythia8_VinciaQED_H
bool isInitial() override
Branching type.
Definition: VinciaQED.h:196
virtual bool isSplitting()
Methods to tell which type of brancher this is.
Definition: VinciaQED.h:139
int jNew
Information for partonSystems.
Definition: VinciaQED.h:162
Class for a QED emission system.
Definition: VinciaQED.h:172
bool isSplitting() override
Branching type: isSplitting() = true.
Definition: VinciaQED.h:315
Definition: Info.h:45
The Event class holds all info on the generated event.
Definition: Event.h:453
int iSys
Event system.
Definition: VinciaQED.h:155
bool isInitial() override
Branching type: isInitial() = true.
Definition: VinciaQED.h:388
Definition: BeamParticle.h:133
QEDsplitSystem()
Constructor.
Definition: VinciaQED.h:294
unsigned int nResDec() override
This module does not implement resonance decays.
Definition: VinciaQED.h:593
void init(Event &event, int xIn, int yIn, double shhIn, double verboseIn)
Initialize.
Definition: VinciaQED.cc:35
unsigned int nBranchers() override
Getter for number of systems.
Definition: VinciaQED.h:586
Base class for Vincia&#39;s QED and EW shower modules.
Definition: VinciaQED.h:444
virtual bool polarise(vector< Particle > &)
Select helicities for a system of particles.
Definition: VinciaQED.h:466
friend class QEDemitSystem
Friends for internal private members.
Definition: VinciaQED.h:33
double generateTrial(Event &event, double q2Start, double q2Low, double alphaIn, double cIn)
Generate a trial point.
Definition: VinciaQED.cc:132
virtual void setVerbose(int verboseIn) override
Set or change verbosity level, and propagate to QED systems.
Definition: VinciaQED.h:550
Class for a QED conversion system.
Definition: VinciaQED.h:361
Definition: StandardModel.h:106
Class for performing QED showers.
Definition: VinciaQED.h:530
Definition: Logger.h:23
QEDconvSystem()
Constructor.
Definition: VinciaQED.h:366
virtual void load()
Some early initialisation (needed for EW shower).
Definition: VinciaQED.h:458
VinciaModule()
Default constructor.
Definition: VinciaQED.h:449
VinciaQED()
Constructor.
Definition: VinciaQED.h:538
Definition: VinciaFSR.h:582
Definition: Basics.h:386
Definition: VinciaCommon.h:494
QEDsplitElemental(Event &event, int iPhotIn, int iSpecIn)
Constructor.
Definition: VinciaQED.h:265
Base class for QED systems.
Definition: VinciaQED.h:100
double q2minColoured() override
Return scales.
Definition: VinciaQED.h:582
QEDemitElemental()
Constuctor.
Definition: VinciaQED.h:36
bool lastIsSplitting() override
Information about last branching.
Definition: VinciaQED.h:569
Class for trial QED splittings.
Definition: VinciaQED.h:257
Class for QED emissions.
Definition: VinciaQED.h:28
double m2(const Vec4 &v1)
The squared invariant mass of one or more four-vectors.
Definition: Basics.cc:605
void clear(int iSys=-1) override
Clear everything, or specific system.
Definition: VinciaQED.h:557
The PartonSystems class describes the whole set of subcollisions.
Definition: PartonSystems.h:42
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
void initPtr(Rndm *rndmPtrIn, PartonSystems *partonSystemsPtrIn)
Initialize the pointers.
Definition: VinciaQED.cc:24
double getKallen()
Kallen function.
Definition: VinciaQED.h:275
This class holds a map of all ParticleDataEntries.
Definition: ParticleData.h:422
int verbose
Verbose setting.
Definition: VinciaQED.h:159
virtual void setVerbose(int verboseIn)
Set verbosity level.
Definition: VinciaQED.h:475
Definition: Basics.h:32
int sysWin() override
Return the system window.
Definition: VinciaQED.h:567
Definition: Settings.h:195
QEDsystem()
Constructor.
Definition: VinciaQED.h:105
Class for a QED splitting system.
Definition: VinciaQED.h:289