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