PYTHIA  8.312
Info.h
1 // Info.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 that keep track of generic event info.
7 // Info: contains information on the generation process and errors.
8 // Info: user interface to make parts of Info publicly available.
9 
10 #ifndef Pythia8_Info_H
11 #define Pythia8_Info_H
12 
13 #include "Pythia8/Basics.h"
14 #include "Pythia8/LHEF3.h"
15 #include "Pythia8/Logger.h"
16 #include "Pythia8/PythiaStdlib.h"
17 #include "Pythia8/SharedPointers.h"
18 #include "Pythia8/Weights.h"
19 
20 namespace Pythia8 {
21 
22 // Forward declaration of various classes.
23 class Settings;
24 class ParticleData;
25 class Rndm;
26 class BeamSetup;
27 class CoupSM;
28 class CoupSUSY;
29 class BeamParticle;
30 class PartonSystems;
31 class SigmaTotal;
32 class SigmaCombined;
33 class HadronWidths;
34 
35 // Forward declaration of HIInfo class.
36 class HIInfo;
37 
38 //==========================================================================
39 
40 // The Info class contains a mixed bag of information on the event
41 // generation activity, especially on the current subprocess properties,
42 // and on the number of errors encountered. This is used by the
43 // generation machinery.
44 
45 class Info {
46 
47 public:
48 
49  // Constructors.
50  Info() = default;
51  Info(bool) : Info() {}
52 
53  // Destructor for clean-up.
54  ~Info() {
55  if (hasOwnEventAttributes && eventAttributes != nullptr)
56  delete eventAttributes;
57  }
58 
59  // Assignment operator gives a shallow copy; no objects pointed to
60  // are copied.
61  Info & operator=(const Info &) = default;
62 
63  // Pointers to other class objects, carried piggyback by Info.
64 
65  // Set pointers to other class objects.
66  void setPtrs(Settings* settingsPtrIn, ParticleData* particleDataPtrIn,
67  Logger* loggerPtrIn, Rndm* rndmPtrIn, BeamSetup* beamSetupIn,
68  CoupSM* coupSMPtrIn, CoupSUSY* coupSUSYPtrIn,
69  PartonSystems* partonSystemsPtrIn, SigmaTotal* sigmaTotPtrIn,
70  SigmaCombined* sigmaCmbPtrIn, HadronWidths* hadronWidthsPtrIn,
71  WeightContainer* weightContainerPtrIn) {
72  settingsPtr = settingsPtrIn; particleDataPtr = particleDataPtrIn;
73  loggerPtr = loggerPtrIn; rndmPtr = rndmPtrIn; beamSetupPtr = beamSetupIn;
74  coupSMPtr = coupSMPtrIn; coupSUSYPtr = coupSUSYPtrIn;
75  partonSystemsPtr = partonSystemsPtrIn; sigmaTotPtr = sigmaTotPtrIn;
76  sigmaCmbPtr = sigmaCmbPtrIn; hadronWidthsPtr = hadronWidthsPtrIn;
77  weightContainerPtr = weightContainerPtrIn; }
78 
79  // Pointer to the settings database.
81 
82  // Pointer to the particle data table.
84 
85  // Pointer to the logger.
87 
88  // Pointer to the random number generator.
90 
91  // Pointers to beam configuration.
93 
94  // Pointers to Standard Model and Beyond SM couplings.
96  CoupSUSY* coupSUSYPtr{};
97 
98  // Pointer to information on subcollision parton locations.
100 
101  // Pointers to the total/elastic/diffractive cross sections.
103  SigmaCombined* sigmaCmbPtr{};
104 
105  // Pointer to the hadron widths data table.
107 
108  // Pointer to the UserHooks object set for the run.
109  UserHooksPtr userHooksPtr{};
110 
111  // Pointer to information about a HeavyIons run and the current event.
112  // (Is NULL if HeavyIons object is inactive.)
114 
115  WeightContainer* weightContainerPtr{};
116 
117  // Listing of most available information on current event.
118  void list() const;
119 
120  // Beam particles (in rest frame). CM energy of event.
121  int idA() const {return idASave;}
122  int idB() const {return idBSave;}
123  double pzA() const {return pzASave;}
124  double pzB() const {return pzBSave;}
125  double eA() const {return eASave;}
126  double eB() const {return eBSave;}
127  double mA() const {return mASave;}
128  double mB() const {return mBSave;}
129  double eCM() const {return eCMSave;}
130  double s() const {return sSave;}
131 
132  // Warnings from initialization.
133  bool tooLowPTmin() const {return lowPTmin;}
134 
135  // Process name and code, and the number of final-state particles.
136  string name() const {return nameSave;}
137  int code() const {return codeSave;}
138  int nFinal() const {return nFinalSave;}
139 
140  // Are beam particles resolved, with pdf's? Are they diffractive?
141  bool isResolved() const {return isRes;}
142  bool isDiffractiveA() const {return isDiffA;}
143  bool isDiffractiveB() const {return isDiffB;}
144  bool isDiffractiveC() const {return isDiffC;}
145  bool isNonDiffractive() const {return isND;}
146  bool isElastic() const {return (codeSave == 102);}
147  // Retained for backwards compatibility.
148  bool isMinBias() const {return isND;}
149 
150  // Information for Les Houches Accord and reading files.
151  bool isLHA() const {return isLH;}
152  bool atEndOfFile() const {return atEOF;}
153 
154  // For nondiffractive and Les Houches Accord identify hardest subprocess.
155  bool hasSub(int i = 0) const {return hasSubSave[i];}
156  string nameSub(int i = 0) const {return nameSubSave[i];}
157  int codeSub(int i = 0) const {return codeSubSave[i];}
158  int nFinalSub(int i = 0) const {return nFinalSubSave[i];}
159 
160  // Incoming parton flavours and x values.
161  int id1(int i = 0) const {return id1Save[i];}
162  int id2(int i = 0) const {return id2Save[i];}
163  double x1(int i = 0) const {return x1Save[i];}
164  double x2(int i = 0) const {return x2Save[i];}
165  double y(int i = 0) const {return 0.5 * log( x1Save[i] / x2Save[i]);}
166  double tau(int i = 0) const {return x1Save[i] * x2Save[i];}
167 
168  // Hard process flavours, x values, parton densities, couplings, Q2 scales.
169  int id1pdf(int i = 0) const {return id1pdfSave[i];}
170  int id2pdf(int i = 0) const {return id2pdfSave[i];}
171  double x1pdf(int i = 0) const {return x1pdfSave[i];}
172  double x2pdf(int i = 0) const {return x2pdfSave[i];}
173  double pdf1(int i = 0) const {return pdf1Save[i];}
174  double pdf2(int i = 0) const {return pdf2Save[i];}
175  double QFac(int i = 0) const {return sqrtpos(Q2FacSave[i]);}
176  double Q2Fac(int i = 0) const {return Q2FacSave[i];}
177  bool isValence1() const {return isVal1;}
178  bool isValence2() const {return isVal2;}
179  double alphaS(int i = 0) const {return alphaSSave[i];}
180  double alphaEM(int i = 0) const {return alphaEMSave[i];}
181  double QRen(int i = 0) const {return sqrtpos(Q2RenSave[i]);}
182  double Q2Ren(int i = 0) const {return Q2RenSave[i];}
183  double scalup(int i = 0) const {return scalupSave[i];}
184 
185  // Kinematics of photons from lepton beams.
186  double xGammaA() const {return x1GammaSave;}
187  double xGammaB() const {return x2GammaSave;}
188  double Q2GammaA() const {return Q2Gamma1Save;}
189  double Q2GammaB() const {return Q2Gamma2Save;}
190  double eCMsub() const {return eCMsubSave;}
191  double thetaScatLepA() const {return thetaLepton1;}
192  double thetaScatLepB() const {return thetaLepton2;}
193  double sHatNew() const {return sHatNewSave;}
194  int photonMode() const {return gammaModeEvent;}
195 
196  // Information on VMD state inside a photon.
197  bool isVMDstateA() const {return isVMDstateAEvent;}
198  bool isVMDstateB() const {return isVMDstateBEvent;}
199  int idVMDA() const {return idVMDASave;}
200  int idVMDB() const {return idVMDBSave;}
201  double mVMDA() const {return mVMDASave;}
202  double mVMDB() const {return mVMDBSave;}
203  double scaleVMDA() const {return scaleVMDASave;}
204  double scaleVMDB() const {return scaleVMDBSave;}
205 
206  // Mandelstam variables (notation as if subcollision).
207  double mHat(int i = 0) const {return sqrt(sH[i]);}
208  double sHat(int i = 0) const {return sH[i];}
209  double tHat(int i = 0) const {return tH[i];}
210  double uHat(int i = 0) const {return uH[i];}
211  double pTHat(int i = 0) const {return pTH[i];}
212  double pT2Hat(int i = 0) const {return pTH[i] * pTH[i];}
213  double m3Hat(int i = 0) const {return m3H[i];}
214  double m4Hat(int i = 0) const {return m4H[i];}
215  double thetaHat(int i = 0) const {return thetaH[i];}
216  double phiHat(int i = 0) const {return phiH[i];}
217 
218  // Weight of current event; normally 1, but used for Les Houches events
219  // or when reweighting phase space selection. Conversion from mb to pb
220  // for LHA strategy +-4. Uncertainty variations can be accessed by
221  // providing an index >= 1 (0 = no variations). Also cumulative sum.
222  double weight(int i = 0) const;
223  double weightSum() const;
224  double lhaStrategy() const {return lhaStrategySave;}
225 
226  // Further access to uncertainty weights: number and labels.
227  int nWeights() const {
228  return weightContainerPtr->weightsShowerPtr->getWeightsSize() +
229  weightContainerPtr->weightsFragmentation.getWeightsSize() - 1;}
230  string weightLabel(int iWgt) const;
231  int nWeightGroups() const {return weightContainerPtr->
232  weightsShowerPtr->nWeightGroups() + weightContainerPtr->
233  weightsFragmentation.nWeightGroups();}
234  string getGroupName(int iGN) const;
235  double getGroupWeight(int iGW) const;
236 
237  // Number of times other steps have been carried out.
238  int nISR() const {return nISRSave;}
239  int nFSRinProc() const {return nFSRinProcSave;}
240  int nFSRinRes() const {return nFSRinResSave;}
241 
242  // Maximum pT scales for MPI, ISR and FSR (in hard process).
243  double pTmaxMPI() const {return pTmaxMPISave;}
244  double pTmaxISR() const {return pTmaxISRSave;}
245  double pTmaxFSR() const {return pTmaxFSRSave;}
246 
247  // Current evolution scale (for UserHooks).
248  double pTnow() const {return pTnowSave;}
249 
250  // Impact parameter picture, global information
251  double a0MPI() const {return a0MPISave;}
252 
253  // Impact parameter picture, as set by hardest interaction.
254  double bMPI() const {return (bIsSet) ? bMPISave : 1.;}
255  double enhanceMPI() const {return (bIsSet) ? enhanceMPISave : 1.;}
256  double enhanceMPIavg() const {return (bIsSet) ? enhanceMPIavgSave : 1.;}
257  double eMPI(int i) const {return (bIsSet) ? eMPISave[i] : 1.;}
258  double bMPIold() const {return (bIsSet) ? bMPIoldSave : 1.;}
259  double enhanceMPIold() const {return (bIsSet) ? enhanceMPIoldSave : 1.;}
260  double enhanceMPIoldavg() const {return (bIsSet)
261  ? enhanceMPIoldavgSave : 1.;}
262 
263  // Number of multiparton interactions, with code and pT for them.
264  int nMPI() const {return nMPISave;}
265  int codeMPI(int i) const {return codeMPISave[i];}
266  double pTMPI(int i) const {return pTMPISave[i];}
267  int iAMPI(int i) const {return iAMPISave[i];}
268  int iBMPI(int i) const {return iBMPISave[i];}
269 
270  // Cross section estimate, optionally process by process.
271  vector<int> codesHard();
272 
273  // Name of the specified process.
274  string nameProc(int i = 0) const {
275  if (i == 0) return "sum";
276  auto itr = procNameM.find(i);
277  if (itr != procNameM.end()) return itr->second;
278  loggerPtr->ERROR_MSG("process code not found", to_string(i));
279  return "unknown process";}
280 
281  // The number of phase-space points tried.
282  long nTried(int i = 0) const {
283  if (i == 0) return nTry;
284  auto itr = nTryM.find(i);
285  if (itr != nTryM.end()) return itr->second;
286  loggerPtr->ERROR_MSG("process code not found", to_string(i));
287  return 0;}
288 
289  // The number of selected hard processes.
290  long nSelected(int i = 0) const {
291  if (i == 0) return nSel;
292  auto itr = nSelM.find(i);
293  if (itr != nSelM.end()) return itr->second;
294  loggerPtr->ERROR_MSG("process code not found", to_string(i));
295  return 0;}
296 
297  // The number of accepted events.
298  long nAccepted(int i = 0) const {
299  if (i == 0) return nAcc;
300  auto itr = nAccM.find(i);
301  if (itr != nAccM.end()) return itr->second;
302  loggerPtr->ERROR_MSG("process code not found", to_string(i));
303  return 0;}
304 
305  // The estimated cross-section in units of mb.
306  double sigmaGen(int i = 0) const {
307  if (i == 0) return sigGen;
308  auto itr = sigGenM.find(i);
309  if (itr != sigGenM.end()) return itr->second;
310  loggerPtr->ERROR_MSG("process code not found", to_string(i));
311  return 0;}
312 
313  // The uncertainty on the estimated cross-section in units of mb.
314  double sigmaErr(int i = 0) const {
315  if (i == 0) return sigErr;
316  auto itr = sigErrM.find(i);
317  if (itr != sigErrM.end()) return itr->second;
318  //loggerPtr->ERROR_MSG("process code not found", to_string(code));
319  return 0;}
320 
321  // Counters for number of loops in various places.
322  int getCounter( int i) const {return counters[i];}
323 
324  // Set or increase the value stored in a counter.
325  void setCounter( int i, int value = 0) {counters[i] = value;}
326  void addCounter( int i, int value = 1) {counters[i] += value;}
327 
328  // Set initialization warning flag when too low pTmin in ISR/FSR/MPI.
329  void setTooLowPTmin(bool lowPTminIn) {lowPTmin = lowPTminIn;}
330 
331  // Set info on valence character of hard collision partons.
332  void setValence( bool isVal1In, bool isVal2In) {isVal1 = isVal1In;
333  isVal2 = isVal2In;}
334 
335  // Set and get some MPI/ISR/FSR properties needed for matching,
336  // i.e. mainly of internal relevance.
337  void hasHistory( bool hasHistoryIn) {hasHistorySave = hasHistoryIn;}
338  bool hasHistory() {return hasHistorySave;}
339  void zNowISR( double zNowIn) {zNowISRSave = zNowIn;}
340  double zNowISR() {return zNowISRSave;}
341  void pT2NowISR( double pT2NowIn) {pT2NowISRSave = pT2NowIn;}
342  double pT2NowISR() {return pT2NowISRSave;}
343 
344  // Return merging weight.
345  double mergingWeight(int i=0) const {
346  return weightContainerPtr->weightsMerging.getWeightsValue(i);}
347 
348  // Return the complete NLO weight.
349  double mergingWeightNLO(int i=0) const {
350  return weightContainerPtr->weightsMerging.getWeightsValue(i); }
351 
352  // Return an LHEF header
353  string header(const string &key) const {
354  if (headers.find(key) == headers.end()) return "";
355  else return headers.at(key);
356  }
357 
358  // Return a list of all header key names
359  vector<string> headerKeys() const;
360 
361  // Return the number of processes in the LHEF.
362  int nProcessesLHEF() const { return int(sigmaLHEFSave.size());}
363  // Return the cross section information read from LHEF.
364  double sigmaLHEF(int iProcess) const { return sigmaLHEFSave[iProcess];}
365 
366  // LHEF3 information: Public for easy access.
368 
369  // Save process information as read from init block of LHEF.
370  vector<double> sigmaLHEFSave;
371 
372  // Contents of the LHAinitrwgt tag
374 
375  // Contents of the LHAgenerator tags.
376  vector<LHAgenerator > *generators{};
377 
378  // A map of the LHAweightgroup tags, indexed by name.
379  map<string,LHAweightgroup > *weightgroups{};
380 
381  // A map of the LHAweight tags, indexed by name.
382  map<string,LHAweight > *init_weights{};
383 
384  // Store current-event Les Houches event tags.
386  map<string, string > *eventAttributes{};
387 
388  // The weights associated with this event, as given by the LHAwgt tags
389  map<string,double > *weights_detailed{};
390 
391  // The weights associated with this event, as given by the LHAweights tags
392  vector<double > *weights_compressed{};
393 
394  // Contents of the LHAscales tag.
396 
397  // Contents of the LHAweights tag (compressed format)
399 
400  // Contents of the LHArwgt tag (detailed format)
402 
403  // Vectorized version of LHArwgt tag, for easy and ordered access.
404  vector<double> weights_detailed_vector;
405 
406  // Value of the unit event weight read from a Les Houches event, necessary
407  // if additional weights in an otherwise unweighted input file are in
408  // relation to this number.
409  double eventWeightLHEF = {};
410 
411  // Set the LHEF3 objects read from the init and header blocks.
412  void setLHEF3InitInfo();
413  void setLHEF3InitInfo( int LHEFversionIn, LHAinitrwgt *initrwgtIn,
414  vector<LHAgenerator> *generatorsIn,
415  map<string,LHAweightgroup> *weightgroupsIn,
416  map<string,LHAweight> *init_weightsIn, string headerBlockIn );
417  // Set the LHEF3 objects read from the event block.
418  void setLHEF3EventInfo();
419  void setLHEF3EventInfo( map<string, string> *eventAttributesIn,
420  map<string, double > *weights_detailedIn,
421  vector<double > *weights_compressedIn,
422  LHAscales *scalesIn, LHAweights *weightsIn,
423  LHArwgt *rwgtIn, vector<double> weights_detailed_vecIn,
424  vector<string> weights_detailed_name_vecIn,
425  string eventCommentsIn, double eventWeightLHEFIn );
426 
427  // Retrieve events tag information.
428  string getEventAttribute(string key, bool doRemoveWhitespace = false) const;
429 
430  // Externally set event tag auxiliary information.
431  void setEventAttribute(string key, string value, bool doOverwrite = true) {
432  if (eventAttributes == nullptr) {
433  eventAttributes = new map<string,string>();
434  hasOwnEventAttributes = true;
435  }
436  map<string, string>::iterator it = eventAttributes->find(key);
437  if ( it != eventAttributes->end() && !doOverwrite ) return;
438  if ( it != eventAttributes->end() ) eventAttributes->erase(it);
439  eventAttributes->insert(make_pair(key,value));
440  }
441 
442  // Retrieve LHEF version
443  int LHEFversion() const { return LHEFversionSave; }
444 
445  // Retrieve initrwgt tag information.
446  unsigned int getInitrwgtSize() const;
447 
448  // Retrieve generator tag information.
449  unsigned int getGeneratorSize() const;
450  string getGeneratorValue(unsigned int n = 0) const;
451  string getGeneratorAttribute( unsigned int n, string key,
452  bool doRemoveWhitespace = false) const;
453 
454  // Retrieve rwgt tag information.
455  unsigned int getWeightsDetailedSize() const;
456  double getWeightsDetailedValue(string n) const;
457  string getWeightsDetailedAttribute(string n, string key,
458  bool doRemoveWhitespace = false) const;
459 
460  // Retrieve weights tag information.
461  unsigned int getWeightsCompressedSize() const;
462  double getWeightsCompressedValue(unsigned int n) const;
463  string getWeightsCompressedAttribute(string key,
464  bool doRemoveWhitespace = false) const;
465 
466  // Retrieve scales tag information.
467  string getScalesValue(bool doRemoveWhitespace = false) const;
468  double getScalesAttribute(string key) const;
469 
470  // Retrieve complete header block and event comments
471  // Retrieve scales tag information.
472  string getHeaderBlock() const { return headerBlock;}
473  string getEventComments() const { return eventComments;}
474 
475  // Set LHEF headers
476  void setHeader(const string &key, const string &val)
477  { headers[key] = val; }
478 
479  // Set abort in parton level.
480  void setAbortPartonLevel(bool abortIn) {abortPartonLevel = abortIn;}
481  bool getAbortPartonLevel() const {return abortPartonLevel;}
482 
483  // Get information on hard diffractive events.
484  bool hasUnresolvedBeams() const {return hasUnresBeams;}
485  bool hasPomPsystem() const {return hasPomPsys;}
486  bool isHardDiffractive() const {return isHardDiffA || isHardDiffB;}
487  bool isHardDiffractiveA() const {return isHardDiffA;}
488  bool isHardDiffractiveB() const {return isHardDiffB;}
489  double xPomeronA() const {return xPomA;}
490  double xPomeronB() const {return xPomB;}
491  double tPomeronA() const {return tPomA;}
492  double tPomeronB() const {return tPomB;}
493 
494  // History information needed to setup the weak shower for 2 -> n.
495  vector<int> getWeakModes() const {return weakModes;}
496  vector<pair<int,int> > getWeakDipoles() const {return weakDipoles;}
497  vector<Vec4> getWeakMomenta() const {return weakMomenta;}
498  vector<int> getWeak2to2lines() const {return weak2to2lines;}
499  void setWeakModes(vector<int> weakModesIn) {weakModes = weakModesIn;}
500  void setWeakDipoles(vector<pair<int,int> > weakDipolesIn)
501  {weakDipoles = weakDipolesIn;}
502  void setWeakMomenta(vector<Vec4> weakMomentaIn)
503  {weakMomenta = weakMomentaIn;}
504  void setWeak2to2lines(vector<int> weak2to2linesIn)
505  {weak2to2lines = weak2to2linesIn;}
506 
507  // Check if onia is included in showers (used for MPI process setup).
508  void setOniumShower(bool oniumShowerIn) {oniumShower = oniumShowerIn;}
509  bool getOniumShower() const {return oniumShower;}
510 
511  // From here on what used to be the private part of the class.
512 
513  // Allow conversion from mb to pb.
514  static const double CONVERTMB2PB;
515 
516  // Store common beam quantities.
517  int idASave{}, idBSave{};
518  double pzASave{}, eASave{}, mASave{}, pzBSave{}, eBSave{}, mBSave{},
519  eCMSave{}, sSave{};
520 
521  // Store initialization information.
522  bool lowPTmin;
523 
524  // Store common integrated cross section quantities.
525  long nTry{}, nSel{}, nAcc{};
526  double sigGen{}, sigErr{}, wtAccSum{};
527  map<int, string> procNameM;
528  map<int, long> nTryM, nSelM, nAccM;
529  map<int, double> sigGenM, sigErrM;
530  int lhaStrategySave{};
531 
532  // Store common MPI information.
533  double a0MPISave;
534 
535  // Store current-event quantities.
536  bool isRes{}, isDiffA{}, isDiffB{}, isDiffC{}, isND{}, isLH{},
537  hasSubSave[4], bIsSet{}, evolIsSet{}, atEOF{}, isVal1{}, isVal2{},
538  hasHistorySave{}, abortPartonLevel{}, isHardDiffA{}, isHardDiffB{},
539  hasUnresBeams{}, hasPomPsys{};
540  int codeSave{}, codeSubSave[4], nFinalSave{}, nFinalSubSave[4], nTotal{},
541  id1Save[4], id2Save[4], id1pdfSave[4], id2pdfSave[4], nMPISave{},
542  nISRSave{}, nFSRinProcSave{}, nFSRinResSave{};
543  double x1Save[4], x2Save[4], x1pdfSave[4], x2pdfSave[4], pdf1Save[4],
544  pdf2Save[4], Q2FacSave[4], alphaEMSave[4], alphaSSave[4],
545  Q2RenSave[4], scalupSave[4], sH[4], tH[4], uH[4], pTH[4], m3H[4],
546  m4H[4], thetaH[4], phiH[4], bMPISave{}, enhanceMPISave{},
547  enhanceMPIavgSave{}, bMPIoldSave{}, enhanceMPIoldSave{},
548  enhanceMPIoldavgSave{}, pTmaxMPISave{}, pTmaxISRSave{},
549  pTmaxFSRSave{}, pTnowSave{}, zNowISRSave{}, pT2NowISRSave{}, xPomA{},
550  xPomB{}, tPomA{}, tPomB{};
551  string nameSave{}, nameSubSave[4];
552  vector<int> codeMPISave, iAMPISave, iBMPISave;
553  vector<double> pTMPISave, eMPISave;
554 
555  // Variables related to photon kinematics.
556  bool isVMDstateAEvent{}, isVMDstateBEvent{};
557  int gammaModeEvent{}, idVMDASave{}, idVMDBSave{};
558  double x1GammaSave{}, x2GammaSave{}, Q2Gamma1Save{}, Q2Gamma2Save{},
559  eCMsubSave{}, thetaLepton1{}, thetaLepton2{}, sHatNewSave{},
560  mVMDASave{}, mVMDBSave{}, scaleVMDASave{}, scaleVMDBSave{};
561 
562  // Vector of various loop counters.
563  int counters[50];
564 
565  // Map for LHEF headers.
566  map<string, string> headers;
567 
568  // Strings for complete header block and event comments.
569  string headerBlock{}, eventComments{};
570 
571  // Set info on the two incoming beams: only from Pythia class.
572  void setBeamIDs( int idAin, int idBin) { idASave = idAin; idBSave = idBin;}
573  void setBeamA( int idAin, double pzAin, double eAin, double mAin) {
574  idASave = idAin; pzASave = pzAin; eASave = eAin; mASave = mAin;}
575  void setBeamB( int idBin, double pzBin, double eBin, double mBin) {
576  idBSave = idBin; pzBSave = pzBin; eBSave = eBin; mBSave = mBin;}
577  void setECM( double eCMin) {eCMSave = eCMin; sSave = eCMSave * eCMSave;}
578 
579  // Set info related to gamma+gamma subcollision.
580  void setX1Gamma( double x1GammaIn) { x1GammaSave = x1GammaIn; }
581  void setX2Gamma( double x2GammaIn) { x2GammaSave = x2GammaIn; }
582  void setQ2Gamma1( double Q2gammaIn) { Q2Gamma1Save = Q2gammaIn; }
583  void setQ2Gamma2( double Q2gammaIn) { Q2Gamma2Save = Q2gammaIn; }
584  void setTheta1( double theta1In) { thetaLepton1 = theta1In; }
585  void setTheta2( double theta2In) { thetaLepton2 = theta2In; }
586  void setECMsub( double eCMsubIn) { eCMsubSave = eCMsubIn; }
587  void setsHatNew( double sHatNewIn) { sHatNewSave = sHatNewIn; }
588  void setGammaMode( double gammaModeIn) { gammaModeEvent = gammaModeIn; }
589  // Set info on VMD state.
590  void setVMDstateA(bool isVMDAIn, int idAIn, double mAIn, double scaleAIn)
591  {isVMDstateAEvent = isVMDAIn; idVMDASave = idAIn; mVMDASave = mAIn;
592  scaleVMDASave = scaleAIn;}
593  void setVMDstateB(bool isVMDBIn, int idBIn, double mBIn, double scaleBIn)
594  {isVMDstateBEvent = isVMDBIn; idVMDBSave = idBIn; mVMDBSave = mBIn;
595  scaleVMDBSave = scaleBIn;}
596 
597  // Reset info for current event: only from Pythia class.
598  void clear() {
599  isRes = isDiffA = isDiffB = isDiffC = isND = isLH = bIsSet
600  = evolIsSet = atEOF = isVal1 = isVal2 = hasHistorySave
601  = isHardDiffA = isHardDiffB = hasUnresBeams = hasPomPsys = false;
602  codeSave = nFinalSave = nTotal = nMPISave = nISRSave = nFSRinProcSave
603  = nFSRinResSave = 0;
604  bMPISave = enhanceMPISave = enhanceMPIavgSave = bMPIoldSave
605  = enhanceMPIoldSave = enhanceMPIoldavgSave = 1.;
606  pTmaxMPISave = pTmaxISRSave = pTmaxFSRSave = pTnowSave = zNowISRSave
607  = pT2NowISRSave = 0.;
608  nameSave = " ";
609  for (int i = 0; i < 4; ++i) {
610  hasSubSave[i] = false;
611  codeSubSave[i] = nFinalSubSave[i] = id1pdfSave[i] = id2pdfSave[i]
612  = id1Save[i] = id2Save[i] = 0;
613  x1pdfSave[i] = x2pdfSave[i] = pdf1Save[i] = pdf2Save[i]
614  = Q2FacSave[i] = alphaEMSave[i] = alphaSSave[i] = Q2RenSave[i]
615  = scalupSave[i] = x1Save[i] = x2Save[i] = sH[i] = tH[i] = uH[i]
616  = pTH[i] = m3H[i] = m4H[i] = thetaH[i] = phiH[i] = 0.;
617  nameSubSave[i] = " ";
618  }
619  codeMPISave.resize(0); iAMPISave.resize(0); iBMPISave.resize(0);
620  pTMPISave.resize(0); eMPISave.resize(0); setHardDiff();
621  weightContainerPtr->clear();
622  }
623 
624  // Reset info arrays only for parton and hadron level.
625  int sizeMPIarrays() const { return codeMPISave.size(); }
626  void resizeMPIarrays(int newSize) { codeMPISave.resize(newSize);
627  iAMPISave.resize(newSize); iBMPISave.resize(newSize);
628  pTMPISave.resize(newSize); eMPISave.resize(newSize); }
629 
630  // Set info on the (sub)process: from ProcessLevel, ProcessContainer or
631  // MultipartonInteractions classes.
632  void setType( string nameIn, int codeIn, int nFinalIn,
633  bool isNonDiffIn = false, bool isResolvedIn = true,
634  bool isDiffractiveAin = false, bool isDiffractiveBin = false,
635  bool isDiffractiveCin = false, bool isLHAin = false) {
636  nameSave = nameIn; codeSave = codeIn; nFinalSave = nFinalIn;
637  isND = isNonDiffIn; isRes = isResolvedIn; isDiffA = isDiffractiveAin;
638  isDiffB = isDiffractiveBin; isDiffC = isDiffractiveCin; isLH = isLHAin;
639  nTotal = 2 + nFinalSave; bIsSet = false; hasSubSave[0] = false;
640  nameSubSave[0] = " "; codeSubSave[0] = 0; nFinalSubSave[0] = 0;
641  evolIsSet = false;}
642  void setSubType( int iDS, string nameSubIn, int codeSubIn,
643  int nFinalSubIn) { hasSubSave[iDS] = true; nameSubSave[iDS] = nameSubIn;
644  codeSubSave[iDS] = codeSubIn; nFinalSubSave[iDS] = nFinalSubIn;}
645  void setPDFalpha( int iDS, int id1pdfIn, int id2pdfIn, double x1pdfIn,
646  double x2pdfIn, double pdf1In, double pdf2In, double Q2FacIn,
647  double alphaEMIn, double alphaSIn, double Q2RenIn, double scalupIn) {
648  id1pdfSave[iDS] = id1pdfIn; id2pdfSave[iDS] = id2pdfIn;
649  x1pdfSave[iDS] = x1pdfIn; x2pdfSave[iDS] = x2pdfIn;
650  pdf1Save[iDS] = pdf1In; pdf2Save[iDS] = pdf2In; Q2FacSave[iDS] = Q2FacIn;
651  alphaEMSave[iDS] = alphaEMIn; alphaSSave[iDS] = alphaSIn;
652  Q2RenSave[iDS] = Q2RenIn; scalupSave[iDS] = scalupIn;}
653  void setScalup( int iDS, double scalupIn) {scalupSave[iDS] = scalupIn;}
654  void setKin( int iDS, int id1In, int id2In, double x1In, double x2In,
655  double sHatIn, double tHatIn, double uHatIn, double pTHatIn,
656  double m3HatIn, double m4HatIn, double thetaHatIn, double phiHatIn) {
657  id1Save[iDS] = id1In; id2Save[iDS] = id2In; x1Save[iDS] = x1In;
658  x2Save[iDS] = x2In; sH[iDS] = sHatIn; tH[iDS] = tHatIn;
659  uH[iDS] = uHatIn; pTH[iDS] = pTHatIn; m3H[iDS] = m3HatIn;
660  m4H[iDS] = m4HatIn; thetaH[iDS] = thetaHatIn; phiH[iDS] = phiHatIn;}
661  void setTypeMPI( int codeMPIIn, double pTMPIIn, int iAMPIIn = 0,
662  int iBMPIIn = 0, double eMPIIn = 1.) {codeMPISave.push_back(codeMPIIn);
663  pTMPISave.push_back(pTMPIIn); iAMPISave.push_back(iAMPIIn);
664  iBMPISave.push_back(iBMPIIn); eMPISave.push_back(eMPIIn); }
665 
666  // Set info on cross section: from ProcessLevel (reset in Pythia).
667  void sigmaReset() { nTry = nSel = nAcc = 0; sigGen = sigErr = wtAccSum = 0.;
668  procNameM.clear(); nTryM.clear(); nSelM.clear(); nAccM.clear();
669  sigGenM.clear(); sigErrM.clear();}
670  void setSigma( int i, string procNameIn, long nTryIn, long nSelIn,
671  long nAccIn, double sigGenIn, double sigErrIn, double wtAccSumIn) {
672  if (i == 0) {nTry = nTryIn; nSel = nSelIn; nAcc = nAccIn;
673  sigGen = sigGenIn; sigErr = sigErrIn; wtAccSum = wtAccSumIn;}
674  else {procNameM[i] = procNameIn; nTryM[i] = nTryIn; nSelM[i] = nSelIn;
675  nAccM[i] = nAccIn; sigGenM[i] = sigGenIn; sigErrM[i] = sigErrIn;} }
676  void addSigma( int i, long nTryIn, long nSelIn, long nAccIn, double sigGenIn,
677  double sigErrIn) { nTryM[i] += nTryIn; nSelM[i] += nSelIn;
678  nAccM[i] += nAccIn; sigGenM[i] += sigGenIn;
679  sigErrM[i] = sqrtpos(sigErrM[i]*sigErrM[i] + sigErrIn*sigErrIn); }
680 
681  // Set info on impact parameter: from PartonLevel.
682  void setImpact( double bMPIIn, double enhanceMPIIn, double enhanceMPIavgIn,
683  bool bIsSetIn = true, bool pushBack = false) {
684  if (pushBack) {bMPIoldSave = bMPISave; enhanceMPIoldSave = enhanceMPISave;
685  enhanceMPIoldavgSave = enhanceMPIavgSave;}
686  bMPISave = bMPIIn; enhanceMPISave = eMPISave[0] = enhanceMPIIn,
687  enhanceMPIavgSave = enhanceMPIavgIn, bIsSet = bIsSetIn;}
688 
689  // Set info on pTmax scales and number of evolution steps: from PartonLevel.
690  void setPartEvolved( int nMPIIn, int nISRIn) {
691  nMPISave = nMPIIn; nISRSave = nISRIn;}
692  void setEvolution( double pTmaxMPIIn, double pTmaxISRIn, double pTmaxFSRIn,
693  int nMPIIn, int nISRIn, int nFSRinProcIn, int nFSRinResIn) {
694  pTmaxMPISave = pTmaxMPIIn; pTmaxISRSave = pTmaxISRIn;
695  pTmaxFSRSave = pTmaxFSRIn; nMPISave = nMPIIn; nISRSave = nISRIn;
696  nFSRinProcSave = nFSRinProcIn; nFSRinResSave = nFSRinResIn;
697  evolIsSet = true;}
698 
699  // Set current pT evolution scale for MPI/ISR/FSR; from PartonLevel.
700  void setPTnow( double pTnowIn) {pTnowSave = pTnowIn;}
701 
702  // Set a0 from MultipartonInteractions.
703  void seta0MPI(double a0MPIin) {a0MPISave = a0MPIin;}
704 
705  // Set info whether reading of Les Houches Accord file at end.
706  void setEndOfFile( bool atEOFin) {atEOF = atEOFin;}
707 
708  // Set event weight, either for LHEF3 or for uncertainty bands. If weight
709  // has units (lhaStrategy 4 or -4), input in mb
710  void setWeight( double weightIn, int lhaStrategyIn) {
711  for (int i = 0; i < nWeights(); ++i)
712  weightContainerPtr->weightsShowerPtr->setValueByIndex(i,1.);
713  // Nominal weight in weightContainer saved in pb for lhaStrategy +-4
714  weightContainerPtr->setWeightNominal(
715  abs(lhaStrategyIn) == 4 ? CONVERTMB2PB*weightIn : weightIn);
716  lhaStrategySave = lhaStrategyIn;}
717  //
718  // Set info on resolved processes.
719  void setIsResolved(bool isResIn) {isRes = isResIn;}
720 
721  // Set info on hard diffraction.
722  void setHardDiff(bool hasUnresBeamsIn = false, bool hasPomPsysIn = false,
723  bool isHardDiffAIn = false, bool isHardDiffBIn = false,
724  double xPomAIn = 0., double xPomBIn = 0., double tPomAIn = 0.,
725  double tPomBIn = 0.) { hasUnresBeams = hasUnresBeamsIn;
726  hasPomPsys = hasPomPsysIn; isHardDiffA = isHardDiffAIn;
727  isHardDiffB = isHardDiffBIn;
728  xPomA = xPomAIn; xPomB = xPomBIn;
729  tPomA = tPomAIn; tPomB = tPomBIn;}
730 
731  // Move process information to an another diffractive system.
732  void reassignDiffSystem(int iDSold, int iDSnew);
733 
734  // Set information in hard diffractive events.
735  void setHasUnresolvedBeams(bool hasUnresBeamsIn)
736  {hasUnresBeams = hasUnresBeamsIn;}
737  void setHasPomPsystem(bool hasPomPsysIn) {hasPomPsys = hasPomPsysIn;}
738 
739  // Variables for weak and onia shower setup.
740  vector<int> weakModes, weak2to2lines;
741  vector<Vec4> weakMomenta;
742  vector<pair<int, int> > weakDipoles;
743  bool oniumShower{false};
744 
745  int numberOfWeights() const {
746  return weightContainerPtr->numberOfWeights(); }
747  double weightValueByIndex(int key=0) const {
748  return weightContainerPtr->weightValueByIndex(key); }
749  string weightNameByIndex(int key=0) const {
750  return weightContainerPtr->weightNameByIndex(key); }
751  vector<double> weightValueVector() const {
752  return weightContainerPtr->weightValueVector(); }
753  vector<string> weightNameVector() const {
754  return weightContainerPtr->weightNameVector(); }
755 
756 };
757 
758 //==========================================================================
759 
760 } // end namespace Pythia8
761 
762 #endif // Pythia8_Info_H
Info()=default
Constructors.
SigmaTotal * sigmaTotPtr
Pointers to the total/elastic/diffractive cross sections.
Definition: Info.h:102
void list() const
Listing of most available information on current event.
Definition: Info.cc:32
string name() const
Process name and code, and the number of final-state particles.
Definition: Info.h:136
void setPartEvolved(int nMPIIn, int nISRIn)
Set info on pTmax scales and number of evolution steps: from PartonLevel.
Definition: Info.h:690
BeamSetup * beamSetupPtr
Pointers to beam configuration.
Definition: Info.h:92
Rndm * rndmPtr
Pointer to the random number generator.
Definition: Info.h:89
bool hasUnresolvedBeams() const
Get information on hard diffractive events.
Definition: Info.h:484
UserHooksPtr userHooksPtr
Pointer to the UserHooks object set for the run.
Definition: Info.h:109
bool isVMDstateAEvent
Variables related to photon kinematics.
Definition: Info.h:556
string getScalesValue(bool doRemoveWhitespace=false) const
Retrieve scales tag information.
Definition: Info.cc:423
Definition: Info.h:45
double bMPI() const
Impact parameter picture, as set by hardest interaction.
Definition: Info.h:254
Definition: BeamSetup.h:33
bool hasOwnEventAttributes
Store current-event Les Houches event tags.
Definition: Info.h:385
void seta0MPI(double a0MPIin)
Set a0 from MultipartonInteractions.
Definition: Info.h:703
Definition: Weights.h:394
WeightsFragmentation weightsFragmentation
Fragmentation weights.
Definition: Weights.h:421
Logger * loggerPtr
Pointer to the logger.
Definition: Info.h:86
Definition: SigmaLowEnergy.h:135
vector< double > * weights_compressed
The weights associated with this event, as given by the LHAweights tags.
Definition: Info.h:392
double pTnow() const
Current evolution scale (for UserHooks).
Definition: Info.h:248
double pTmaxMPI() const
Maximum pT scales for MPI, ISR and FSR (in hard process).
Definition: Info.h:243
int counters[50]
Vector of various loop counters.
Definition: Info.h:563
Info & operator=(const Info &)=default
int idASave
Store common beam quantities.
Definition: Info.h:517
void setCounter(int i, int value=0)
Set or increase the value stored in a counter.
Definition: Info.h:325
CoupSM * coupSMPtr
Pointers to Standard Model and Beyond SM couplings.
Definition: Info.h:95
Definition: SigmaTotal.h:141
int nMPI() const
Number of multiparton interactions, with code and pT for them.
Definition: Info.h:264
vector< double > weights_detailed_vector
Vectorized version of LHArwgt tag, for easy and ordered access.
Definition: Info.h:404
bool isVMDstateA() const
Information on VMD state inside a photon.
Definition: Info.h:197
double mergingWeightNLO(int i=0) const
Return the complete NLO weight.
Definition: Info.h:349
map< string, string > headers
Map for LHEF headers.
Definition: Info.h:566
vector< int > codesHard()
Cross section estimate, optionally process by process.
Definition: Info.cc:212
double getWeightsValue(int iPos) const override
Modified weight getter to include first order weight.
Definition: Weights.h:257
double mHat(int i=0) const
Mandelstam variables (notation as if subcollision).
Definition: Info.h:207
LHAweights * weights
Contents of the LHAweights tag (compressed format)
Definition: Info.h:398
long nAccepted(int i=0) const
The number of accepted events.
Definition: Info.h:298
double sigmaErr(int i=0) const
The uncertainty on the estimated cross-section in units of mb.
Definition: Info.h:314
void reassignDiffSystem(int iDSold, int iDSnew)
Move process information to an another diffractive system.
Definition: Info.cc:453
Collect different scales relevant for an event.
Definition: LHEF3.h:281
map< string, double > * weights_detailed
The weights associated with this event, as given by the LHAwgt tags.
Definition: Info.h:389
The LHAweights struct represents the information in a weights tag.
Definition: LHEF3.h:245
Definition: Logger.h:23
string header(const string &key) const
Return an LHEF header.
Definition: Info.h:353
vector< LHAgenerator > * generators
Contents of the LHAgenerator tags.
Definition: Info.h:376
unsigned int getGeneratorSize() const
Retrieve generator tag information.
Definition: Info.cc:331
unsigned int getInitrwgtSize() const
Retrieve initrwgt tag information.
Definition: Info.cc:322
int id1pdf(int i=0) const
Hard process flavours, x values, parton densities, couplings, Q2 scales.
Definition: Info.h:169
bool hasSub(int i=0) const
For nondiffractive and Les Houches Accord identify hardest subprocess.
Definition: Info.h:155
void setWeight(double weightIn, int lhaStrategyIn)
Definition: Info.h:710
int LHEFversion() const
Retrieve LHEF version.
Definition: Info.h:443
The LHArwgt assigns a group-name to a set of LHAwgt objects.
Definition: LHEF3.h:473
void setOniumShower(bool oniumShowerIn)
Check if onia is included in showers (used for MPI process setup).
Definition: Info.h:508
long nSelected(int i=0) const
The number of selected hard processes.
Definition: Info.h:290
The LHAinitrwgt assigns a group-name to a set of LHAweightgroup objects.
Definition: LHEF3.h:511
void setEventAttribute(string key, string value, bool doOverwrite=true)
Externally set event tag auxiliary information.
Definition: Info.h:431
unsigned int getWeightsDetailedSize() const
Retrieve rwgt tag information.
Definition: Info.cc:362
int id1(int i=0) const
Incoming parton flavours and x values.
Definition: Info.h:161
LHAscales * scales
Contents of the LHAscales tag.
Definition: Info.h:395
Definition: Basics.h:386
bool isLHA() const
Information for Les Houches Accord and reading files.
Definition: Info.h:151
double a0MPI() const
Impact parameter picture, global information.
Definition: Info.h:251
void setVMDstateA(bool isVMDAIn, int idAIn, double mAIn, double scaleAIn)
Set info on VMD state.
Definition: Info.h:590
double mergingWeight(int i=0) const
Return merging weight.
Definition: Info.h:345
void setLHEF3EventInfo()
Set the LHEF3 objects read from the event block.
Definition: Info.cc:266
void setEndOfFile(bool atEOFin)
Set info whether reading of Les Houches Accord file at end.
Definition: Info.h:706
int nISR() const
Number of times other steps have been carried out.
Definition: Info.h:238
double a0MPISave
Store common MPI information.
Definition: Info.h:533
void setBeamIDs(int idAin, int idBin)
Set info on the two incoming beams: only from Pythia class.
Definition: Info.h:572
void setHasUnresolvedBeams(bool hasUnresBeamsIn)
Set information in hard diffractive events.
Definition: Info.h:735
long nTried(int i=0) const
The number of phase-space points tried.
Definition: Info.h:282
int numberOfWeights()
Functions to retrieve information stored in the subcategory members.
Definition: Weights.cc:871
HIInfo * hiInfo
Definition: Info.h:113
Definition: HadronWidths.h:22
void setWeightNominal(double weightNow)
The WeightContainer class.
Definition: Weights.cc:852
void clear()
Reset info for current event: only from Pythia class.
Definition: Info.h:598
double xGammaA() const
Kinematics of photons from lepton beams.
Definition: Info.h:186
string headerBlock
Strings for complete header block and event comments.
Definition: Info.h:569
void hasHistory(bool hasHistoryIn)
Definition: Info.h:337
map< string, LHAweight > * init_weights
A map of the LHAweight tags, indexed by name.
Definition: Info.h:382
map< string, LHAweightgroup > * weightgroups
A map of the LHAweightgroup tags, indexed by name.
Definition: Info.h:379
vector< string > weightNameVector()
Definition: Weights.cc:915
vector< int > getWeakModes() const
History information needed to setup the weak shower for 2 -> n.
Definition: Info.h:495
bool lowPTmin
Store initialization information.
Definition: Info.h:522
void setHardDiff(bool hasUnresBeamsIn=false, bool hasPomPsysIn=false, bool isHardDiffAIn=false, bool isHardDiffBIn=false, double xPomAIn=0., double xPomBIn=0., double tPomAIn=0., double tPomBIn=0.)
Set info on hard diffraction.
Definition: Info.h:722
vector< string > headerKeys() const
Return a list of all header key names.
Definition: Info.cc:224
~Info()
Destructor for clean-up.
Definition: Info.h:54
void setValueByIndex(int iPos, double val)
Functions to set values of weights.
Definition: Weights.h:95
int nProcessesLHEF() const
Return the number of processes in the LHEF.
Definition: Info.h:362
double sigmaGen(int i=0) const
The estimated cross-section in units of mb.
Definition: Info.h:306
Definition: StandardModel.h:135
void setX1Gamma(double x1GammaIn)
Set info related to gamma+gamma subcollision.
Definition: Info.h:580
unsigned int getWeightsCompressedSize() const
Retrieve weights tag information.
Definition: Info.cc:394
void clear()
Reset all members to default stage.
Definition: Weights.cc:939
bool isRes
Store current-event quantities.
Definition: Info.h:536
void setValence(bool isVal1In, bool isVal2In)
Set info on valence character of hard collision partons.
Definition: Info.h:332
long nTry
Store common integrated cross section quantities.
Definition: Info.h:525
int nWeights() const
Further access to uncertainty weights: number and labels.
Definition: Info.h:227
double sigmaLHEF(int iProcess) const
Return the cross section information read from LHEF.
Definition: Info.h:364
bool isResolved() const
Are beam particles resolved, with pdf&#39;s? Are they diffractive?
Definition: Info.h:141
vector< int > weakModes
Variables for weak and onia shower setup.
Definition: Info.h:740
ParticleData * particleDataPtr
Pointer to the particle data table.
Definition: Info.h:83
PartonSystems * partonSystemsPtr
Pointer to information on subcollision parton locations.
Definition: Info.h:99
void setIsResolved(bool isResIn)
Set info on resolved processes.
Definition: Info.h:719
string getHeaderBlock() const
Definition: Info.h:472
double eventWeightLHEF
Definition: Info.h:409
LHArwgt * rwgt
Contents of the LHArwgt tag (detailed format)
Definition: Info.h:401
HadronWidths * hadronWidthsPtr
Pointer to the hadron widths data table.
Definition: Info.h:106
void setImpact(double bMPIIn, double enhanceMPIIn, double enhanceMPIavgIn, bool bIsSetIn=true, bool pushBack=false)
Set info on impact parameter: from PartonLevel.
Definition: Info.h:682
static const double CONVERTMB2PB
From here on what used to be the private part of the class.
Definition: Info.h:514
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 setTooLowPTmin(bool lowPTminIn)
Set initialization warning flag when too low pTmin in ISR/FSR/MPI.
Definition: Info.h:329
LHAinitrwgt * initrwgt
Contents of the LHAinitrwgt tag.
Definition: Info.h:373
void setAbortPartonLevel(bool abortIn)
Set abort in parton level.
Definition: Info.h:480
vector< double > weightValueVector()
Definition: Weights.cc:888
bool isMinBias() const
Retained for backwards compatibility.
Definition: Info.h:148
WeightsMerging weightsMerging
Merging weights.
Definition: Weights.h:418
void sigmaReset()
Set info on cross section: from ProcessLevel (reset in Pythia).
Definition: Info.h:667
void setPTnow(double pTnowIn)
Set current pT evolution scale for MPI/ISR/FSR; from PartonLevel.
Definition: Info.h:700
void setPtrs(Settings *settingsPtrIn, ParticleData *particleDataPtrIn, Logger *loggerPtrIn, Rndm *rndmPtrIn, BeamSetup *beamSetupIn, CoupSM *coupSMPtrIn, CoupSUSY *coupSUSYPtrIn, PartonSystems *partonSystemsPtrIn, SigmaTotal *sigmaTotPtrIn, SigmaCombined *sigmaCmbPtrIn, HadronWidths *hadronWidthsPtrIn, WeightContainer *weightContainerPtrIn)
Pointers to other class objects, carried piggyback by Info.
Definition: Info.h:66
string getEventAttribute(string key, bool doRemoveWhitespace=false) const
Retrieve events tag information.
Definition: Info.cc:307
string nameProc(int i=0) const
Name of the specified process.
Definition: Info.h:274
int sizeMPIarrays() const
Reset info arrays only for parton and hadron level.
Definition: Info.h:625
bool tooLowPTmin() const
Warnings from initialization.
Definition: Info.h:133
vector< double > sigmaLHEFSave
Save process information as read from init block of LHEF.
Definition: Info.h:370
This class holds a map of all ParticleDataEntries.
Definition: ParticleData.h:422
void setLHEF3InitInfo()
Set the LHEF3 objects read from the init and header blocks.
Definition: Info.cc:235
int getCounter(int i) const
Counters for number of loops in various places.
Definition: Info.h:322
Settings * settingsPtr
Pointer to the settings database.
Definition: Info.h:80
WeightsShower * weightsShowerPtr
Definition: Weights.h:414
void setHeader(const string &key, const string &val)
Set LHEF headers.
Definition: Info.h:476
int LHEFversionSave
LHEF3 information: Public for easy access.
Definition: Info.h:367
Definition: HIInfo.h:27
void setType(string nameIn, int codeIn, int nFinalIn, bool isNonDiffIn=false, bool isResolvedIn=true, bool isDiffractiveAin=false, bool isDiffractiveBin=false, bool isDiffractiveCin=false, bool isLHAin=false)
Definition: Info.h:632
Definition: SusyCouplings.h:27
double weight(int i=0) const
Event weights and accumulated weight.
Definition: Info.cc:162
Definition: Settings.h:195
double sqrtpos(const double &x)
Avoid problem with negative square root argument (from roundoff).
Definition: PythiaStdlib.h:191
int idA() const
Beam particles (in rest frame). CM energy of event.
Definition: Info.h:121