StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Event.h
1 // Event.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2020 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 
6 // Header file for the Particle and Event classes.
7 // Particle: information on an instance of a particle.
8 // Junction: information on a junction between three colours.
9 // Event: list of particles in the current event.
10 
11 #ifndef Pythia8_Event_H
12 #define Pythia8_Event_H
13 
14 #include "Pythia8/Basics.h"
15 #include "Pythia8/ParticleData.h"
16 #include "Pythia8/PythiaStdlib.h"
17 
18 namespace Pythia8 {
19 
20 //==========================================================================
21 
22 // Forward references to ParticleDataEntry and ResonanceWidths classes.
23 class ParticleDataEntry;
24 class ResonanceWidths;
25 class Event;
26 
27 //==========================================================================
28 
29 // Particle class.
30 // This class holds info on a particle in general.
31 
32 class Particle {
33 
34 public:
35 
36  // Constructors.
37  Particle() : idSave(0), statusSave(0), mother1Save(0), mother2Save(0),
38  daughter1Save(0), daughter2Save(0), colSave(0), acolSave(0),
39  pSave(Vec4(0.,0.,0.,0.)), mSave(0.), scaleSave(0.), polSave(9.),
40  hasVertexSave(false), vProdSave(Vec4(0.,0.,0.,0.)), tauSave(0.),
41  pdePtr(0), evtPtr(0) { }
42  Particle(int idIn, int statusIn = 0, int mother1In = 0,
43  int mother2In = 0, int daughter1In = 0, int daughter2In = 0,
44  int colIn = 0, int acolIn = 0, double pxIn = 0., double pyIn = 0.,
45  double pzIn = 0., double eIn = 0., double mIn = 0.,
46  double scaleIn = 0., double polIn = 9.)
47  : idSave(idIn), statusSave(statusIn), mother1Save(mother1In),
48  mother2Save(mother2In), daughter1Save(daughter1In),
49  daughter2Save(daughter2In), colSave(colIn), acolSave(acolIn),
50  pSave(Vec4(pxIn, pyIn, pzIn, eIn)), mSave(mIn), scaleSave(scaleIn),
51  polSave(polIn), hasVertexSave(false), vProdSave(Vec4(0.,0.,0.,0.)),
52  tauSave(0.), pdePtr(0), evtPtr(0) { }
53  Particle(int idIn, int statusIn, int mother1In, int mother2In,
54  int daughter1In, int daughter2In, int colIn, int acolIn,
55  Vec4 pIn, double mIn = 0., double scaleIn = 0., double polIn = 9.)
56  : idSave(idIn), statusSave(statusIn), mother1Save(mother1In),
57  mother2Save(mother2In), daughter1Save(daughter1In),
58  daughter2Save(daughter2In), colSave(colIn), acolSave(acolIn),
59  pSave(pIn), mSave(mIn), scaleSave(scaleIn), polSave(polIn),
60  hasVertexSave(false), vProdSave(Vec4(0.,0.,0.,0.)), tauSave(0.),
61  pdePtr(0), evtPtr(0) { }
62  Particle(const Particle& pt) : idSave(pt.idSave),
63  statusSave(pt.statusSave), mother1Save(pt.mother1Save),
64  mother2Save(pt.mother2Save), daughter1Save(pt.daughter1Save),
65  daughter2Save(pt.daughter2Save), colSave(pt.colSave),
66  acolSave(pt.acolSave), pSave(pt.pSave), mSave(pt.mSave),
67  scaleSave(pt.scaleSave), polSave(pt.polSave),
68  hasVertexSave(pt.hasVertexSave), vProdSave(pt.vProdSave),
69  tauSave(pt.tauSave), pdePtr(pt.pdePtr), evtPtr(pt.evtPtr) { }
70  Particle& operator=(const Particle& pt) {if (this != &pt) {
71  idSave = pt.idSave; statusSave = pt.statusSave;
72  mother1Save = pt.mother1Save; mother2Save = pt.mother2Save;
73  daughter1Save = pt.daughter1Save; daughter2Save = pt.daughter2Save;
74  colSave = pt.colSave; acolSave = pt.acolSave; pSave = pt.pSave;
75  mSave = pt.mSave; scaleSave = pt.scaleSave; polSave = pt.polSave;
76  hasVertexSave = pt.hasVertexSave; vProdSave = pt.vProdSave;
77  tauSave = pt.tauSave; pdePtr = pt.pdePtr; evtPtr = pt.evtPtr; }
78  return *this; }
79 
80  // Destructor.
81  virtual ~Particle() {}
82 
83  // Member functions to set the Event and ParticleDataEntry pointers.
84  void setEvtPtr(Event* evtPtrIn) { evtPtr = evtPtrIn; setPDEPtr();}
85  void setPDEPtr(ParticleDataEntry* pdePtrIn = 0);
86 
87  // Member functions for input.
88  void id(int idIn) {idSave = idIn; setPDEPtr();}
89  void status(int statusIn) {statusSave = statusIn;}
90  void statusPos() {statusSave = abs(statusSave);}
91  void statusNeg() {statusSave = -abs(statusSave);}
92  void statusCode(int statusIn) {statusSave =
93  (statusSave > 0) ? abs(statusIn) : -abs(statusIn);}
94  void mother1(int mother1In) {mother1Save = mother1In;}
95  void mother2(int mother2In) {mother2Save = mother2In;}
96  void mothers(int mother1In = 0, int mother2In = 0)
97  {mother1Save = mother1In; mother2Save = mother2In;}
98  void daughter1(int daughter1In) {daughter1Save = daughter1In;}
99  void daughter2(int daughter2In) {daughter2Save = daughter2In;}
100  void daughters(int daughter1In = 0, int daughter2In = 0)
101  {daughter1Save = daughter1In; daughter2Save = daughter2In;}
102  void col(int colIn) {colSave = colIn;}
103  void acol(int acolIn) {acolSave = acolIn;}
104  void cols(int colIn = 0,int acolIn = 0) {colSave = colIn;
105  acolSave = acolIn;}
106  void p(Vec4 pIn) {pSave = pIn;}
107  void p(double pxIn, double pyIn, double pzIn, double eIn)
108  {pSave.p(pxIn, pyIn, pzIn, eIn);}
109  void px(double pxIn) {pSave.px(pxIn);}
110  void py(double pyIn) {pSave.py(pyIn);}
111  void pz(double pzIn) {pSave.pz(pzIn);}
112  void e(double eIn) {pSave.e(eIn);}
113  void m(double mIn) {mSave = mIn;}
114  void scale(double scaleIn) {scaleSave = scaleIn;}
115  void pol(double polIn) {polSave = polIn;}
116  void vProd(Vec4 vProdIn) {vProdSave = vProdIn; hasVertexSave = true;}
117  void vProd(double xProdIn, double yProdIn, double zProdIn, double tProdIn)
118  {vProdSave.p(xProdIn, yProdIn, zProdIn, tProdIn); hasVertexSave = true;}
119  void xProd(double xProdIn) {vProdSave.px(xProdIn); hasVertexSave = true;}
120  void yProd(double yProdIn) {vProdSave.py(yProdIn); hasVertexSave = true;}
121  void zProd(double zProdIn) {vProdSave.pz(zProdIn); hasVertexSave = true;}
122  void tProd(double tProdIn) {vProdSave.e(tProdIn); hasVertexSave = true;}
123  void vProdAdd(Vec4 vProdIn) {vProdSave += vProdIn; hasVertexSave = true;}
124  void tau(double tauIn) {tauSave = tauIn;}
125 
126  // Member functions for output.
127  int id() const {return idSave;}
128  int status() const {return statusSave;}
129  int mother1() const {return mother1Save;}
130  int mother2() const {return mother2Save;}
131  int daughter1() const {return daughter1Save;}
132  int daughter2() const {return daughter2Save;}
133  int col() const {return colSave;}
134  int acol() const {return acolSave;}
135  Vec4 p() const {return pSave;}
136  double px() const {return pSave.px();}
137  double py() const {return pSave.py();}
138  double pz() const {return pSave.pz();}
139  double e() const {return pSave.e();}
140  double m() const {return mSave;}
141  double scale() const {return scaleSave;}
142  double pol() const {return polSave;}
143  bool hasVertex() const {return hasVertexSave;}
144  Vec4 vProd() const {return vProdSave;}
145  double xProd() const {return vProdSave.px();}
146  double yProd() const {return vProdSave.py();}
147  double zProd() const {return vProdSave.pz();}
148  double tProd() const {return vProdSave.e();}
149  double tau() const {return tauSave;}
150 
151  // Member functions for output; derived int and bool quantities.
152  int idAbs() const {return abs(idSave);}
153  int statusAbs() const {return abs(statusSave);}
154  bool isFinal() const {return (statusSave > 0);}
155  int intPol() const;
156  bool isRescatteredIncoming() const {return
157  (statusSave == -34 || statusSave == -45 ||
158  statusSave == -46 || statusSave == -54);}
159 
160  // Member functions for output; derived double quantities.
161  double m2() const {return (mSave >= 0.) ? mSave*mSave
162  : -mSave*mSave;}
163  double mCalc() const {return pSave.mCalc();}
164  double m2Calc() const {return pSave.m2Calc();}
165  double eCalc() const {return sqrt(abs(m2() + pSave.pAbs2()));}
166  double pT() const {return pSave.pT();}
167  double pT2() const {return pSave.pT2();}
168  double mT() const {double temp = m2() + pSave.pT2();
169  return (temp >= 0.) ? sqrt(temp) : -sqrt(-temp);}
170  double mT2() const {return m2() + pSave.pT2();}
171  double pAbs() const {return pSave.pAbs();}
172  double pAbs2() const {return pSave.pAbs2();}
173  double eT() const {return pSave.eT();}
174  double eT2() const {return pSave.eT2();}
175  double theta() const {return pSave.theta();}
176  double phi() const {return pSave.phi();}
177  double thetaXZ() const {return pSave.thetaXZ();}
178  double pPos() const {return pSave.pPos();}
179  double pNeg() const {return pSave.pNeg();}
180  double y() const;
181  double eta() const;
182  double y(double mCut) const;
183  double y(double mCut, RotBstMatrix& M) const;
184  Vec4 vDec() const {return (tauSave > 0. && mSave > 0.)
185  ? vProdSave + tauSave * pSave / mSave : vProdSave;}
186  double xDec() const {return (tauSave > 0. && mSave > 0.)
187  ? vProdSave.px() + tauSave * pSave.px() / mSave : vProdSave.px();}
188  double yDec() const {return (tauSave > 0. && mSave > 0.)
189  ? vProdSave.py() + tauSave * pSave.py() / mSave : vProdSave.py();}
190  double zDec() const {return (tauSave > 0. && mSave > 0.)
191  ? vProdSave.pz() + tauSave * pSave.pz() / mSave : vProdSave.pz();}
192  double tDec() const {return (tauSave > 0. && mSave > 0.)
193  ? vProdSave.e() + tauSave * pSave.e() / mSave : vProdSave.e();}
194 
195  // Methods that can refer back to the event the particle belongs to.
196  virtual int index() const;
197  int iTopCopy() const;
198  int iBotCopy() const;
199  int iTopCopyId(bool simplify = false) const;
200  int iBotCopyId(bool simplify = false) const;
201  vector<int> motherList() const;
202  vector<int> daughterList() const;
203  vector<int> daughterListRecursive() const;
204  vector<int> sisterList(bool traceTopBot = false) const;
205  bool isAncestor(int iAncestor) const;
206  int statusHepMC() const;
207  bool isFinalPartonLevel() const;
208  bool undoDecay();
209 
210  // Further output, based on a pointer to a ParticleDataEntry object.
211  string name() const {
212  return (pdePtr != 0) ? pdePtr->name(idSave) : " ";}
213  string nameWithStatus(int maxLen = 20) const;
214  int spinType() const {
215  return (pdePtr != 0) ? pdePtr->spinType() : 0;}
216  int chargeType() const {
217  return (pdePtr != 0) ? pdePtr->chargeType(idSave) : 0;}
218  double charge() const {
219  return (pdePtr != 0) ? pdePtr->charge(idSave) : 0;}
220  bool isCharged() const {
221  return (pdePtr != 0) ? (pdePtr->chargeType(idSave) != 0) : false;}
222  bool isNeutral() const {
223  return (pdePtr != 0) ? (pdePtr->chargeType(idSave) == 0) : false;}
224  int colType() const {
225  return (pdePtr != 0) ? pdePtr->colType(idSave) : 0;}
226  double m0() const {
227  return (pdePtr != 0) ? pdePtr->m0() : 0.;}
228  double mWidth() const {
229  return (pdePtr != 0) ? pdePtr->mWidth() : 0.;}
230  double mMin() const {
231  return (pdePtr != 0) ? pdePtr->mMin() : 0.;}
232  double mMax() const {
233  return (pdePtr != 0) ? pdePtr->mMax() : 0.;}
234  double mSel() const {
235  return (pdePtr != 0) ? pdePtr->mSel() : 0.;}
236  double constituentMass() const {
237  return (pdePtr != 0) ? pdePtr->constituentMass() : 0.;}
238  double tau0() const {
239  return (pdePtr != 0) ? pdePtr->tau0() : 0.;}
240  bool mayDecay() const {
241  return (pdePtr != 0) ? pdePtr->mayDecay() : false;}
242  bool canDecay() const {
243  return (pdePtr != 0) ? pdePtr->canDecay() : false;}
244  bool doExternalDecay() const {
245  return (pdePtr != 0) ? pdePtr->doExternalDecay() : false;}
246  bool isResonance() const {
247  return (pdePtr != 0) ? pdePtr->isResonance() : false;}
248  bool isVisible() const {
249  return (pdePtr != 0) ? pdePtr->isVisible() : false;}
250  bool isLepton() const {
251  return (pdePtr != 0) ? pdePtr->isLepton() : false;}
252  bool isQuark() const {
253  return (pdePtr != 0) ? pdePtr->isQuark() : false;}
254  bool isGluon() const {
255  return (pdePtr != 0) ? pdePtr->isGluon() : false;}
256  bool isDiquark() const {
257  return (pdePtr != 0) ? pdePtr->isDiquark() : false;}
258  bool isParton() const {
259  return (pdePtr != 0) ? pdePtr->isParton() : false;}
260  bool isHadron() const {
261  return (pdePtr != 0) ? pdePtr->isHadron() : false;}
262  ParticleDataEntry& particleDataEntry() const {return *pdePtr;}
263 
264  // Member functions that perform operations.
265  void rescale3(double fac) {pSave.rescale3(fac);}
266  void rescale4(double fac) {pSave.rescale4(fac);}
267  void rescale5(double fac) {pSave.rescale4(fac); mSave *= fac;}
268  void rot(double thetaIn, double phiIn) {pSave.rot(thetaIn, phiIn);
269  if (hasVertexSave) vProdSave.rot(thetaIn, phiIn);}
270  void bst(double betaX, double betaY, double betaZ) {
271  pSave.bst(betaX, betaY, betaZ);
272  if (hasVertexSave) vProdSave.bst(betaX, betaY, betaZ);}
273  void bst(double betaX, double betaY, double betaZ, double gamma) {
274  pSave.bst(betaX, betaY, betaZ, gamma);
275  if (hasVertexSave) vProdSave.bst(betaX, betaY, betaZ, gamma);}
276  void bst(const Vec4& pBst) {pSave.bst(pBst);
277  if (hasVertexSave) vProdSave.bst(pBst);}
278  void bst(const Vec4& pBst, double mBst) {pSave.bst(pBst, mBst);
279  if (hasVertexSave) vProdSave.bst(pBst, mBst);}
280  void bstback(const Vec4& pBst) {pSave.bstback(pBst);
281  if (hasVertexSave) vProdSave.bstback(pBst);}
282  void bstback(const Vec4& pBst, double mBst) {pSave.bstback(pBst, mBst);
283  if (hasVertexSave) vProdSave.bstback(pBst, mBst);}
284  void rotbst(const RotBstMatrix& M, bool boostVertex = true) {pSave.rotbst(M);
285  if (hasVertexSave && boostVertex) vProdSave.rotbst(M);}
286  void offsetHistory( int minMother, int addMother, int minDaughter,
287  int addDaughter);
288  void offsetCol( int addCol);
289 
290 protected:
291 
292  // Constants: could only be changed in the code itself.
293  static const double TINY;
294 
295  // Properties of the current particle.
296  int idSave, statusSave, mother1Save, mother2Save, daughter1Save,
297  daughter2Save, colSave, acolSave;
298  Vec4 pSave;
299  double mSave, scaleSave, polSave;
300  bool hasVertexSave;
301  Vec4 vProdSave;
302  double tauSave;
303 
304  // Pointer to properties of the particle species.
305  // Should no be saved in a persistent copy of the event record.
306  // The //! below is ROOT notation that this member should not be saved.
307  // Event::restorePtrs() can be called to restore the missing information.
308  ParticleDataEntry* pdePtr;
309 
310  // Pointer to the whole event record to which the particle belongs (if any).
311  // As above it should not be saved.
312  Event* evtPtr;
313 
314 };
315 
316 // Invariant mass of a pair and its square.
317 // (Not part of class proper, but tightly linked.)
318 double m(const Particle&, const Particle&);
319 double m2(const Particle&, const Particle&);
320 
321 //==========================================================================
322 
323 // The junction class stores what kind of junction it is, the colour indices
324 // of the legs at the junction and as far out as legs have been traced,
325 // and the status codes assigned for fragmentation of each leg.
326 
327 class Junction {
328 
329 public:
330 
331  // Constructors.
332  Junction() : remainsSave(true), kindSave(0), colSave(), endColSave(),
333  statusSave() { }
334 
335  Junction( int kindIn, int col0In, int col1In, int col2In)
336  : remainsSave(true), kindSave(kindIn), colSave(), endColSave(),
337  statusSave() {
338  colSave[0] = col0In; colSave[1] = col1In; colSave[2] = col2In;
339  for (int j = 0; j < 3; ++j) {
340  endColSave[j] = colSave[j]; } }
341  Junction(const Junction& ju) : remainsSave(ju.remainsSave),
342  kindSave(ju.kindSave), colSave(), endColSave(), statusSave() {
343  for (int j = 0; j < 3; ++j) {
344  colSave[j] = ju.colSave[j]; endColSave[j] = ju.endColSave[j];
345  statusSave[j] = ju.statusSave[j]; } }
346  Junction& operator=(const Junction& ju) {if (this != &ju) {
347  remainsSave = ju.remainsSave; kindSave = ju.kindSave;
348  for (int j = 0; j < 3; ++j) { colSave[j] = ju.colSave[j];
349  endColSave[j] = ju.endColSave[j]; statusSave[j] = ju.statusSave[j]; } }
350  return *this; }
351 
352  // Set values.
353  void remains(bool remainsIn) {remainsSave = remainsIn;}
354  void col(int j, int colIn) {colSave[j] = colIn; endColSave[j] = colIn;}
355  void cols(int j, int colIn, int endColIn) {colSave[j] = colIn;
356  endColSave[j] = endColIn;}
357  void endCol(int j, int endColIn) {endColSave[j] = endColIn;}
358  void status(int j, int statusIn) {statusSave[j] = statusIn;}
359 
360  // Read out value.
361  bool remains() const {return remainsSave;}
362  int kind() const {return kindSave;}
363  int col(int j) const {return colSave[j];}
364  int endCol(int j) const {return endColSave[j];}
365  int status(int j) const {return statusSave[j];}
366 
367 private:
368 
369  // Kind, positions of the three ends and their status codes.
370  bool remainsSave;
371  int kindSave, colSave[3], endColSave[3], statusSave[3];
372 
373 };
374 
375 //==========================================================================
376 
377 // The Event class holds all info on the generated event.
378 
379 class Event {
380 
381 public:
382 
383  // Constructors.
384  Event(int capacity = 100) : startColTag(100), maxColTag(100),
385  savedSize(0), savedJunctionSize(0), savedPartonLevelSize(0),
386  scaleSave(0.), scaleSecondSave(0.),
387  headerList("----------------------------------------"),
388  particleDataPtr(0) { entry.reserve(capacity); }
389  Event& operator=(const Event& oldEvent);
390  Event(const Event& oldEvent) {*this = oldEvent;}
391 
392  // Initialize header for event listing, particle data table, and colour.
393  void init( string headerIn = "", ParticleData* particleDataPtrIn = 0,
394  int startColTagIn = 100) {
395  headerList.replace(0, headerIn.length() + 2, headerIn + " ");
396  particleDataPtr = particleDataPtrIn; startColTag = startColTagIn;}
397 
398  // Clear event record.
399  void clear() {entry.resize(0); maxColTag = startColTag;
400  savedPartonLevelSize = 0; scaleSave = 0.; scaleSecondSave = 0.;
401  clearJunctions();}
402  void free() {vector<Particle>().swap(entry); maxColTag = startColTag;
403  savedPartonLevelSize = 0; scaleSave = 0.; scaleSecondSave = 0.;
404  clearJunctions();}
405 
406  // Clear event record, and set first particle empty.
407  void reset() {clear(); append(90, -11, 0, 0, 0., 0., 0., 0., 0.);}
408 
409  // Overload index operator to access element of event record.
410  Particle& operator[](int i) {return entry.at(i);}
411  const Particle& operator[](int i) const {return entry.at(i);}
412 
413  // Implement standard references to elements in the particle array.
414  Particle& front() {return entry.front();}
415  Particle& at(int i) {return entry.at(i);}
416  Particle& back() {return entry.back();}
417 
418  // Implement iterators for the particle array.
419  vector<Pythia8::Particle>::iterator begin() { return entry.begin(); }
420  vector<Pythia8::Particle>::iterator end() { return entry.end(); }
421  vector<Pythia8::Particle>::const_iterator begin() const
422  { return entry.begin(); }
423  vector<Pythia8::Particle>::const_iterator end() const
424  { return entry.end(); }
425 
426  // Event record size.
427  int size() const {return entry.size();}
428 
429  // Put a new particle at the end of the event record; return index.
430  int append(Particle entryIn) {
431  entry.push_back(entryIn); setEvtPtr();
432  if (entryIn.col() > maxColTag) maxColTag = entryIn.col();
433  if (entryIn.acol() > maxColTag) maxColTag = entryIn.acol();
434  return entry.size() - 1;
435  }
436  int append(int id, int status, int mother1, int mother2, int daughter1,
437  int daughter2, int col, int acol, double px, double py, double pz,
438  double e, double m = 0., double scaleIn = 0., double polIn = 9.) {
439  entry.push_back( Particle(id, status, mother1, mother2, daughter1,
440  daughter2, col, acol, px, py, pz, e, m, scaleIn, polIn) ); setEvtPtr();
441  if (col > maxColTag) maxColTag = col;
442  if (acol > maxColTag) maxColTag = acol;
443  return entry.size() - 1;
444  }
445  int append(int id, int status, int mother1, int mother2, int daughter1,
446  int daughter2, int col, int acol, Vec4 p, double m = 0.,
447  double scaleIn = 0., double polIn = 9.) {
448  entry.push_back( Particle(id, status, mother1, mother2, daughter1,
449  daughter2, col, acol, p, m, scaleIn, polIn) ); setEvtPtr();
450  if (col > maxColTag) maxColTag = col;
451  if (acol > maxColTag) maxColTag = acol;
452  return entry.size() - 1;
453  }
454 
455  // Brief versions of append: no mothers and no daughters.
456  int append(int id, int status, int col, int acol, double px, double py,
457  double pz, double e, double m = 0., double scaleIn = 0.,
458  double polIn = 9.) { entry.push_back( Particle(id, status, 0, 0, 0, 0,
459  col, acol, px, py, pz, e, m, scaleIn, polIn) ); setEvtPtr();
460  if (col > maxColTag) maxColTag = col;
461  if (acol > maxColTag) maxColTag = acol;
462  return entry.size() - 1;
463  }
464  int append(int id, int status, int col, int acol, Vec4 p, double m = 0.,
465  double scaleIn = 0., double polIn = 9.) {entry.push_back( Particle(id,
466  status, 0, 0, 0, 0, col, acol, p, m, scaleIn, polIn) ); setEvtPtr();
467  if (col > maxColTag) maxColTag = col;
468  if (acol > maxColTag) maxColTag = acol;
469  return entry.size() - 1;
470  }
471 
472  // Set pointer to the event for a particle, by default latest one.
473  void setEvtPtr(int iSet = -1) {if (iSet < 0) iSet = entry.size() - 1;
474  entry[iSet].setEvtPtr( this);}
475 
476  // Add a copy of an existing particle at the end of the event record.
477  int copy(int iCopy, int newStatus = 0);
478 
479  // List the particles in an event.
480  void list(bool showScaleAndVertex = false,
481  bool showMothersAndDaughters = false, int precision = 3) const;
482 
483  // Remove last n entries.
484  void popBack(int nRemove = 1) { if (nRemove ==1) entry.pop_back();
485  else {int newSize = max( 0, size() - nRemove);
486  entry.resize(newSize);} }
487 
488  // Remove entries from iFirst to iLast, including endpoints, and fix history.
489  // (To the extent possible; history pointers in removed range are zeroed.)
490  void remove(int iFirst, int iLast, bool shiftHistory = true);
491 
492  // Restore all ParticleDataEntry* pointers in the Particle vector.
493  // Useful when a persistent copy of the event record is read back in.
494  void restorePtrs() { for (int i = 0; i < size(); ++i) setEvtPtr(i); }
495 
496  // Save or restore the size of the event record (throwing at the end).
497  void saveSize() {savedSize = entry.size();}
498  void restoreSize() {entry.resize(savedSize);}
499  int savedSizeValue() {return savedSize;}
500 
501  // Initialize and access colour tag information.
502  void initColTag(int colTag = 0) {maxColTag = max( colTag,startColTag);}
503  int lastColTag() const {return maxColTag;}
504  int nextColTag() {return ++maxColTag;}
505 
506  // Access scale for which event as a whole is defined.
507  void scale( double scaleIn) {scaleSave = scaleIn;}
508  double scale() const {return scaleSave;}
509 
510  // Need a second scale if two hard interactions in event.
511  void scaleSecond( double scaleSecondIn) {scaleSecondSave = scaleSecondIn;}
512  double scaleSecond() const {return scaleSecondSave;}
513 
514  // Find complete list of daughters.
515  // Note: temporarily retained for CMS compatibility. Do not use!
516  vector<int> daughterList(int i) const {return entry[i].daughterList();}
517 
518  // Return number of final-state particles, optionally charged only.
519  int nFinal(bool chargedOnly = false) const {
520  int nFin = 0;
521  for (int i = 0; i < size(); ++i)
522  if (entry[i].isFinal() && (!chargedOnly || entry[i].isCharged()))
523  ++nFin;
524  return nFin; }
525 
526  // Find separation in y, eta, phi or R between two particles.
527  double dyAbs(int i1, int i2) const {
528  return abs( entry[i1].y() - entry[i2].y() ); }
529  double detaAbs(int i1, int i2) const {
530  return abs( entry[i1].eta() - entry[i2].eta() ); }
531  double dphiAbs(int i1, int i2) const {
532  double dPhiTmp = abs( entry[i1].phi() - entry[i2].phi() );
533  if (dPhiTmp > M_PI) dPhiTmp = 2. * M_PI - dPhiTmp;
534  return dPhiTmp; }
535  double RRapPhi(int i1, int i2) const {
536  return sqrt( pow2(dyAbs(i1, i2)) + pow2(dphiAbs(i1, i2)) ); }
537  double REtaPhi(int i1, int i2) const {
538  return sqrt( pow2(detaAbs(i1, i2)) + pow2(dphiAbs(i1, i2)) ); }
539 
540  // Member functions for rotations and boosts of an event.
541  void rot(double theta, double phi)
542  {for (int i = 0; i < size(); ++i) entry[i].rot(theta, phi);}
543  void bst(double betaX, double betaY, double betaZ)
544  {for (int i = 0; i < size(); ++i) entry[i].bst(betaX, betaY, betaZ);}
545  void bst(double betaX, double betaY, double betaZ, double gamma)
546  {for (int i = 0; i < size(); ++i) entry[i].bst(betaX, betaY, betaZ,
547  gamma);}
548  void bst(const Vec4& vec)
549  {for (int i = 0; i < size(); ++i) entry[i].bst(vec);}
550  void rotbst(const RotBstMatrix& M, bool boostVertices = true)
551  {for (int i = 0; i < size(); ++i) entry[i].rotbst(M, boostVertices);}
552 
553  // Clear the list of junctions.
554  void clearJunctions() {junction.resize(0);}
555 
556  // Add a junction to the list, study it or extra input.
557  int appendJunction( int kind, int col0, int col1, int col2)
558  { junction.push_back( Junction( kind, col0, col1, col2) );
559  return junction.size() - 1;}
560  int appendJunction(Junction junctionIn) {junction.push_back(junctionIn);
561  return junction.size() - 1;}
562  int sizeJunction() const {return junction.size();}
563  bool remainsJunction(int i) const {return junction[i].remains();}
564  void remainsJunction(int i, bool remainsIn) {junction[i].remains(remainsIn);}
565  int kindJunction(int i) const {return junction[i].kind();}
566  int colJunction( int i, int j) const {return junction[i].col(j);}
567  void colJunction( int i, int j, int colIn) {junction[i].col(j, colIn);}
568  int endColJunction( int i, int j) const {return junction[i].endCol(j);}
569  void endColJunction( int i, int j, int endColIn)
570  {junction[i].endCol(j, endColIn);}
571  int statusJunction( int i, int j) const {return junction[i].status(j);}
572  void statusJunction( int i, int j, int statusIn)
573  {junction[i].status(j, statusIn);}
574  Junction& getJunction(int i) {return junction[i];}
575  const Junction& getJunction(int i) const {return junction[i];}
576  void eraseJunction(int i);
577 
578  // Save or restore the size of the junction list (throwing at the end).
579  void saveJunctionSize() {savedJunctionSize = junction.size();}
580  void restoreJunctionSize() {junction.resize(savedJunctionSize);}
581 
582  // List any junctions in the event; for debug mainly.
583  void listJunctions() const;
584 
585  // Save event record size at Parton Level, i.e. before hadronization.
586  void savePartonLevelSize() {savedPartonLevelSize = entry.size();}
587 
588  // Operator overloading allows to append one event to an existing one.
589  // Warning: particles should be OK, but some other information unreliable.
590  Event& operator+=(const Event& addEvent);
591 
592 private:
593 
594  // The Particle class needs to access particle data.
595  friend class Particle;
596 
597  // Constants: could only be changed in the code itself.
598  static const int IPERLINE;
599 
600  // Initialization data, normally only set once.
601  int startColTag;
602 
603  // The event: a vector containing all particles (entries).
604  // The explicit use of Pythia8:: qualifier patches a limitation in ROOT.
605  vector<Pythia8::Particle> entry;
606 
607  // The list of junctions.
608  // The explicit use of Pythia8:: qualifier patches a limitation in ROOT.
609  vector<Pythia8::Junction> junction;
610 
611  // The maximum colour tag of the event so far.
612  int maxColTag;
613 
614  // Saved entry and junction list sizes, for simple restoration.
615  int savedSize, savedJunctionSize, savedPartonLevelSize;
616 
617  // The scale of the event; linear quantity in GeV.
618  double scaleSave, scaleSecondSave;
619 
620  // Header specification in event listing (at most 40 characters wide).
621  string headerList;
622 
623  // Pointer to the particle data table.
624  // The //! below is ROOT notation that this member should not be saved.
625  ParticleData* particleDataPtr;
626 
627 };
628 
629 //==========================================================================
630 
631 } // end namespace Pythia8
632 
633 #endif // end Pythia8_Event_H
Event
Definition: AgUStep.h:26
Generated by   doxygen 1.8.5