|
27 | 27 | #include "CommonUtils/ConfigurationMacroHelper.h" |
28 | 28 |
|
29 | 29 | #include "TLorentzVector.h" // for TLorentzVector |
| 30 | +#include "TDatabasePDG.h" // for TDatabasePDG |
30 | 31 | #include "TParticle.h" // for TParticle |
31 | 32 | #include "TRefArray.h" // for TRefArray |
32 | 33 | #include "TVirtualMC.h" // for VMC |
|
36 | 37 | #include <cassert> |
37 | 38 | #include <cstddef> // for NULL |
38 | 39 | #include <cmath> |
| 40 | +#include <optional> |
| 41 | +#include <unordered_map> |
39 | 42 | #include "DetectorsBase/SimRandom.h" |
40 | 43 |
|
41 | 44 | using std::cout; |
42 | 45 | using std::endl; |
43 | 46 | using std::pair; |
44 | 47 | using namespace o2::data; |
45 | 48 |
|
| 49 | +namespace |
| 50 | +{ |
| 51 | +// a TParticle whose PDG record is looked up once per thread and PDG code, since every TDatabasePDG lookup takes ROOT's global lock |
| 52 | +TParticle makeParticle(int pdg, int status, int mother1, int mother2, int daughter1, int daughter2, |
| 53 | + double px, double py, double pz, double e, double vx, double vy, double vz, double t) |
| 54 | +{ |
| 55 | + // an empty optional marks a PDG code that is not in TDatabasePDG |
| 56 | + thread_local std::unordered_map<int, std::optional<TParticle>> prototypes; |
| 57 | + auto it = prototypes.find(pdg); |
| 58 | + if (it == prototypes.end()) { |
| 59 | + std::optional<TParticle> proto; |
| 60 | + if (TDatabasePDG::Instance()->GetParticle(pdg)) { |
| 61 | + proto.emplace(pdg, 0, -1, -1, -1, -1, 0., 0., 0., 0., 0., 0., 0., 0.); |
| 62 | + } |
| 63 | + it = prototypes.emplace(pdg, std::move(proto)).first; |
| 64 | + } |
| 65 | + if (!it->second) { // the mass then comes from the 4-momentum |
| 66 | + return TParticle(pdg, status, mother1, mother2, daughter1, daughter2, px, py, pz, e, vx, vy, vz, t); |
| 67 | + } |
| 68 | + TParticle p(*it->second); // the copy keeps the PDG record and the mass |
| 69 | + p.SetStatusCode(status); |
| 70 | + p.SetFirstMother(mother1); |
| 71 | + p.SetLastMother(mother2); |
| 72 | + p.SetFirstDaughter(daughter1); |
| 73 | + p.SetLastDaughter(daughter2); |
| 74 | + p.SetMomentum(px, py, pz, e); |
| 75 | + p.SetProductionVertex(vx, vy, vz, t); |
| 76 | + return p; |
| 77 | +} |
| 78 | +} // namespace |
| 79 | + |
46 | 80 | // small helper function to append to vector at arbitrary position |
47 | 81 | template <typename T, typename I> |
48 | 82 | void insertInVector(std::vector<T>& v, I index, T e) |
@@ -219,7 +253,7 @@ void Stack::PushTrack(Int_t toBeDone, Int_t parentId, Int_t pdgCode, Double_t px |
219 | 253 | // Int_t daughter1Id = -1; |
220 | 254 | // Int_t daughter2Id = -1; |
221 | 255 | Int_t iStatus = (proc == kPPrimary) ? is : trackId; |
222 | | - TParticle p(pdgCode, iStatus, parentId, secondparentId, daughter1Id, daughter2Id, px, py, pz, e, vx, vy, vz, time); |
| 256 | + TParticle p = makeParticle(pdgCode, iStatus, parentId, secondparentId, daughter1Id, daughter2Id, px, py, pz, e, vx, vy, vz, time); |
223 | 257 | p.SetPolarisation(polx, poly, polz); |
224 | 258 | p.SetWeight(weight); |
225 | 259 | p.SetUniqueID(proc); // using the unique ID to transfer process ID |
|
0 commit comments