From 5e317b0ccd72a892b7d6712285703a49b92da20d Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Mon, 7 Aug 2023 10:54:15 +0200 Subject: [PATCH 01/22] Prepare the second version of the migration changeset --- .../src/ForwardTransportModel.cxx | 2 +- .../TruthUtils/TruthUtils/MagicNumbers.h | 96 +++++++--- .../src/MuonSegmentTruthAssociationAlg.cxx | 1 + .../src/MuonTruthAssociationAlg.cxx | 1 + .../src/MuonTruthDecorationAlg.cxx | 2 +- .../D3PDMakerTest/src/HitsFillerAlg.cxx | 1 + .../MCTruthClassifier/IMCTruthClassifier.h | 5 +- .../MCTruthClassifier/MCTruthClassifier.h | 2 +- .../Root/MCTruthClassifierGen.cxx | 72 +++++--- .../src/MCTruthClassifierAthena.cxx | 20 +- .../MCTruthBase/src/RecordingEnvelope.cxx | 2 +- .../TrackRecord/TrackRecord/TrackRecord.h | 15 +- .../test/TrackRecordCnv_p1_test.cxx | 4 +- .../G4UserActions/src/CosmicPerigeeAction.cxx | 2 +- .../TrackWriteFastSim/src/TrackFastSimSD.cxx | 4 +- .../ISF_Event/ISF_Event/ISFTruthIncident.h | 12 +- .../ISF_Event/ISF_Event/ITruthIncident.h | 3 + .../ISF_Event/src/ISFTruthIncident.cxx | 29 ++- .../ISF_Event/test/ISFTruthIncident_test.cxx | 33 ++-- .../ISF_Core/ISF_Services/src/TruthSvc.cxx | 171 ++++-------------- .../ISF/ISF_Core/ISF_Services/src/TruthSvc.h | 10 +- .../ISF_Services/test/TruthSvc_test.cxx | 51 +----- .../src/PunchThroughTool.cxx | 12 +- .../src/PunchThroughTool.h | 2 +- .../src/EntryLayerTool.cxx | 1 + .../src/EntryLayerToolMT.cxx | 1 + .../ISF_Geant4Event/Geant4TruthIncident.h | 18 +- .../src/Geant4TruthIncident.cxx | 16 +- .../ISF_Geant4Tools/src/ISFTrajectory.cxx | 4 +- .../src/TrackProcessorUserActionFullG4.cxx | 4 +- 30 files changed, 268 insertions(+), 328 deletions(-) diff --git a/ForwardDetectors/ForwardTransport/src/ForwardTransportModel.cxx b/ForwardDetectors/ForwardTransport/src/ForwardTransportModel.cxx index d00068ca63dd..6ba8097f6007 100644 --- a/ForwardDetectors/ForwardTransport/src/ForwardTransportModel.cxx +++ b/ForwardDetectors/ForwardTransport/src/ForwardTransportModel.cxx @@ -155,7 +155,7 @@ void ForwardTransportModel::DoIt(const G4FastTrack& fastTrack, G4FastStep& fastS pdgcode, HepMC::FORWARDTRANSPORTMODELSTATUS); gVertex->add_particle_out(gParticle); - HepMC::suggest_barcode(gParticle, HepMC::barcode(part)+HepMC::SIM_REGENERATION_INCREMENT); + //HepMC::suggest_barcode(gParticle, HepMC::barcode(part)+HepMC::SIM_REGENERATION_INCREMENT);//FIXME // Create secondary on the Geant4 side fastStep.SetNumberOfSecondaryTracks(1); diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index 4c23197bed09..8d291d7ae8cf 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -9,6 +9,7 @@ #include <limits> #include <cstdint> #include <memory> +#include <deque> #include <type_traits> #if defined(HEPMC3) && !defined(XAOD_STANDALONE) #include "AtlasHepMC/GenEvent.h" @@ -45,46 +46,87 @@ constexpr int crazyParticleBarcode(std::numeric_limits<int32_t>::max()); constexpr int INVALID_PARTICLE_BARCODE = -1; /// @brief Method to establish if a particle (or barcode) was created during the simulation -template <class T> inline bool is_simulation_particle(const T& p){ return (p->barcode()>SIM_BARCODE_THRESHOLD);} -template <> inline bool is_simulation_particle(const int& b){ return (b>SIM_BARCODE_THRESHOLD);} -#if defined(HEPMC3) && !defined(XAOD_STANDALONE) -template <> inline bool is_simulation_particle(const ConstGenParticlePtr& p){ return (barcode(p)>SIM_BARCODE_THRESHOLD);} -template <> inline bool is_simulation_particle(const GenParticlePtr& p){ return (barcode(p)>SIM_BARCODE_THRESHOLD);} -#endif +template <class T> inline bool is_simulation_particle(const T& p){ return (p->status()>SIM_STATUS_INCREMENT);} /// @brief Method to return how many interactions a particle has undergone during simulation. -template <class T> inline int generations(const T& p){ return (p->barcode()/SIM_REGENERATION_INCREMENT);} -template <> inline int generations(const int& b){ return (b/SIM_REGENERATION_INCREMENT);} -#if defined(HEPMC3) && !defined(XAOD_STANDALONE) -template <> inline int generations(const ConstGenParticlePtr& p){ return (barcode(p)/SIM_REGENERATION_INCREMENT);} -template <> inline int generations(const GenParticlePtr& p){ return (barcode(p)/SIM_REGENERATION_INCREMENT);} -#endif +template <class T> inline int generations(const T& p){ return (p->status()/SIM_STATUS_INCREMENT);} /// @brief Method to establish if the vertex was created during simulation -template <class T> inline bool is_simulation_vertex(const T& p){ return (p->barcode()<-SIM_BARCODE_THRESHOLD);} -template <> inline bool is_simulation_vertex(const int& p){ return (p<-SIM_BARCODE_THRESHOLD);} +template <class T> inline bool is_simulation_vertex(const T& p){ return (p->status()>SIM_STATUS_INCREMENT);} + + + +template <class T> int get_id(const T& p1){ return p1->barcode();} #if defined(HEPMC3) && !defined(XAOD_STANDALONE) -template <> inline bool is_simulation_vertex(const ConstGenVertexPtr& p){ return (barcode(p)<-SIM_BARCODE_THRESHOLD);} -template <> inline bool is_simulation_vertex(const GenVertexPtr& p){ return (barcode(p)<-SIM_BARCODE_THRESHOLD);} +template <> inline int get_id(const ConstGenParticlePtr& p1){ return p1->id();} +template <> inline int get_id(const GenParticlePtr& p1){ return p1->id();} #endif -template <class T1,class T2, std::enable_if_t< !std::is_arithmetic<T1>::value && !std::is_arithmetic<T2>::value, bool> = true > -inline bool is_same_generator_particle(const T1& p1,const T2& p2) { return p1->barcode() % SIM_REGENERATION_INCREMENT == p2->barcode() % SIM_REGENERATION_INCREMENT; } +template <class T> inline T get_generator_particle(const T& p) { + auto pp = p; + for (;;) { + if (pp->status()<SIM_STATUS_INCREMENT) break; + auto pv = pp->production_vertex(); + if (!pv) break; + for (auto pa: pv->particles_in()) { + if (pa->pdg_id()!=pp->pdg_id()) continue; + pp = pa; + break; + } + } + return pp; +} -template <class T1,class T2, std::enable_if_t< !std::is_arithmetic<T1>::value && std::is_arithmetic<T2>::value, bool> = true > -inline bool is_same_generator_particle(const T1& p1,const T2& p2) { return p1->barcode() % SIM_REGENERATION_INCREMENT == p2 % SIM_REGENERATION_INCREMENT; } +template <class T> +inline std::deque<int> get_particle_history(const T& p) { + std::deque<int> res; + res.push_back(get_id(p)); + T pp1 = p; + for (;;) { + if (pp1->status()<SIM_STATUS_INCREMENT) break; + auto pv = pp1->production_vertex(); + if (!pv) break; + for (auto pa: pv->particles_in()) { + if (pa->pdg_id()!=pp1->pdg_id()) continue; + res.push_front(get_id(pp1)); + //pp1 = pa;//FIXME + break; + } + } + + T pp2 = p; + for (;;) { + auto pv = pp2->end_vertex(); + if (!pv) break; + for (auto pa: pv->particles_out()) { + if (pa->pdg_id()!=pp2->pdg_id()) continue; + res.push_back(get_id(pp2)); + //pp2 = pa; //FIXME + break; + } + } + return res; +} -template <class T1,class T2, std::enable_if_t< std::is_arithmetic<T1>::value && std::is_arithmetic<T2>::value, bool> = true > -inline bool is_same_generator_particle(const T1& p1,const T2& p2) { return p1 % SIM_REGENERATION_INCREMENT == p2 % SIM_REGENERATION_INCREMENT; } +template <class T1, class T2, std::enable_if_t<std::is_same<T1,T2>::value, bool> = true > +inline bool is_same_generator_particle(const T1& p1,const T2& p2) { + auto p1orig = get_generator_particle(p1); + auto p2orig = get_generator_particle(p2); + return get_id(p1orig)==get_id(p2orig); +} +template <class T1, class T2, std::enable_if_t< !std::is_same<T1,T2>::value, bool> = true > +inline bool is_same_generator_particle(const T1& p1,const T2& p2) { + int b = get_id(p1); + std::deque<int> hist = get_particle_history(p2); + return (std::find(hist.begin(),hist.end(),b) != hist.end()); +} -/// @brief Method to check if the first particle is a descendant of the second in the simulation, i.e. particle p1 was produced simulations particle p2. -template <class T1,class T2, std::enable_if_t< !std::is_arithmetic<T1>::value && !std::is_arithmetic<T2>::value, bool> = true > -inline bool is_sim_descendant(const T1& p1,const T2& p2) { return p1->barcode() % SIM_REGENERATION_INCREMENT == p2->barcode();} -template <class T1,class T2, std::enable_if_t< std::is_arithmetic<T1>::value && !std::is_arithmetic<T2>::value, bool> = true > -inline bool is_sim_descendant(const T1& p1,const T2& p2) { return p1 % SIM_REGENERATION_INCREMENT == p2->barcode() ; } +template <class T1, class T2> inline bool is_sim_descendant(const T1& p1,const T2& p2) { + return (p1->status()/SIM_STATUS_INCREMENT >= p2->status()/SIM_STATUS_INCREMENT) && is_same_generator_particle(p1,p2); +} } #endif diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx index fe7d5c93b4b3..9cddfe9997e9 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx @@ -10,6 +10,7 @@ #include "xAODMuon/MuonSegmentAuxContainer.h" #include "xAODTruth/TruthParticleAuxContainer.h" #include "xAODTruth/TruthParticleContainer.h" +#include "xAODTruth/TruthVertex.h" #include "TruthUtils/MagicNumbers.h" namespace Muon { diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthAssociationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthAssociationAlg.cxx index 210679c32620..cebe2b0fc737 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthAssociationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthAssociationAlg.cxx @@ -13,6 +13,7 @@ #include "MuonCombinedEvent/TagBase.h" #include "FourMomUtils/xAODP4Helpers.h" #include "TruthUtils/HepMCHelpers.h" +#include "xAODTruth/TruthVertex.h" using namespace xAOD::P4Helpers; namespace { constexpr unsigned int dummy_unsigned = 999; diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx index 94b01614031c..c2ef3fd6a2d3 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx @@ -463,7 +463,7 @@ namespace Muon { // loop over trajectories for (const auto& trajectory : *col) { // check if gen particle same as input - if (!HepMC::is_sim_descendant(&trajectory.second,&truthParticle)) continue; + if (!HepMC::is_sim_descendant(trajectory.second.cptr(),&truthParticle)) continue; const Identifier& id = trajectory.first; bool measPhi = m_idHelperSvc->measuresPhi(id); diff --git a/PhysicsAnalysis/D3PDMaker/D3PDMakerTest/src/HitsFillerAlg.cxx b/PhysicsAnalysis/D3PDMaker/D3PDMakerTest/src/HitsFillerAlg.cxx index 6b282e4f2495..5db8b610feee 100644 --- a/PhysicsAnalysis/D3PDMaker/D3PDMakerTest/src/HitsFillerAlg.cxx +++ b/PhysicsAnalysis/D3PDMaker/D3PDMakerTest/src/HitsFillerAlg.cxx @@ -99,6 +99,7 @@ StatusCode HitsFillerAlg::fillTrackRecordCollection() CLHEP::Hep3Vector p(j+11,j+12,j+13); CLHEP::Hep3Vector x(j+14,j+15,j+16); c->Emplace(j, // PDG + 0, //STATUS (j+10)*1000, // energy p, //position x, //momentum diff --git a/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/IMCTruthClassifier.h b/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/IMCTruthClassifier.h index 83b0f24ee27f..49642295aee4 100644 --- a/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/IMCTruthClassifier.h +++ b/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/IMCTruthClassifier.h @@ -75,12 +75,15 @@ public: MCTruthPartClassifier::UnknownOutCome; const xAOD::TruthParticle* mother = nullptr; + long motherStatus = 0; long motherBarcode = 0; int motherPDG = 0; + const xAOD::TruthParticle* Mother() const { return mother;} long photonMotherBarcode = 0; long photonMotherStatus = 0; - int photonMotherPDG = 0; + int photonMotherPDG = 0; + const xAOD::TruthParticle* PhotonMother() const { return photonMother;} const xAOD::TruthParticle* photonMother = nullptr; const xAOD::TruthParticle* bkgElecMother = nullptr; diff --git a/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/MCTruthClassifier.h b/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/MCTruthClassifier.h index b822d0852cf5..203421437729 100644 --- a/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/MCTruthClassifier.h +++ b/PhysicsAnalysis/MCTruthClassifier/MCTruthClassifier/MCTruthClassifier.h @@ -231,7 +231,7 @@ private: MCTruthPartClassifier::ParticleOrigin defJetOrig(const std::set<const xAOD::TruthParticle*>&) const; // /** Return true if genParticle and truthParticle have the same pdgId, barcode and status **/ - static const xAOD::TruthParticle* barcode_to_particle(const xAOD::TruthParticleContainer*, int); + static const xAOD::TruthParticle* data_to_particle(const xAOD::TruthParticleContainer*, const xAOD::TruthParticle* xx); /* Data members*/ SG::ReadHandleKey<xAOD::TruthParticleContainer> m_truthParticleContainerKey{ diff --git a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx index e13479bed7d0..411d5a71610b 100644 --- a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx @@ -182,7 +182,9 @@ MCTruthClassifier::particleTruthClassifier(const xAOD::TruthParticle* thePart, I motherStatus = theMoth->status(); motherBarcode = theMoth->barcode(); if (info) { + info->mother = theMoth; info->motherPDG = motherPDG; + info->motherStatus = motherStatus; info->motherBarcode = motherBarcode; } } @@ -217,7 +219,7 @@ MCTruthClassifier::particleTruthClassifier(const xAOD::TruthParticle* thePart, I return std::make_pair(Neutrino, partOrig); } - if (HepMC::is_same_generator_particle(thePart,motherBarcode)) + if (HepMC::is_same_generator_particle(thePart,info->Mother())) return std::make_pair(NonPrimary, partOrig); if (isPartHadr) @@ -467,7 +469,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth ATH_MSG_DEBUG("Executing DefOrigOfElectron "); - const xAOD::TruthParticle* thePriPart = barcode_to_particle(mcTruthTES, thePart->barcode()); + const xAOD::TruthParticle* thePriPart = data_to_particle(mcTruthTES, thePart); if (!thePriPart) return NonDefined; if (abs(thePriPart->pdgId()) != 11) return NonDefined; @@ -489,9 +491,12 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth return NonDefined; } int motherPDG = mother->pdgId(); + int motherStatus = mother->status(); long motherBarcode = mother->barcode(); if (info) { + info->mother = mother; info->motherPDG = motherPDG; + info->motherStatus = motherStatus; info->motherBarcode = motherBarcode; } const xAOD::TruthVertex* mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; @@ -502,7 +507,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth const xAOD::TruthParticle* theDaug = partOriVert->outgoingParticle(ipOut); if (!theDaug) continue; - if (motherPDG == theDaug->pdgId() && HepMC::is_same_generator_particle(theDaug,motherBarcode)) + if (motherPDG == theDaug->pdgId() && HepMC::is_same_generator_particle(theDaug,info->Mother())) samePart = true; } @@ -620,7 +625,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth if (abs(DaugType) == 42) NumOfLQ++; - if (abs(DaugType) == abs(motherPDG) && HepMC::is_same_generator_particle(theDaug, motherBarcode)) + if (abs(DaugType) == abs(motherPDG) && HepMC::is_same_generator_particle(theDaug, info->Mother())) samePart = true; if (numOfParents == 1 && (motherPDG == 22 || abs(motherPDG) == 11 || abs(motherPDG) == 13 || abs(motherPDG) == 211) && @@ -941,7 +946,7 @@ MCTruthClassifier::defOrigOfMuon(const xAOD::TruthParticleContainer* mcTruthTES, ATH_MSG_DEBUG("Executing DefOrigOfMuon "); - const xAOD::TruthParticle* thePriPart = barcode_to_particle(mcTruthTES, thePart->barcode()); + const xAOD::TruthParticle* thePriPart = data_to_particle(mcTruthTES, thePart); if (!thePriPart) return NonDefined; if (abs(thePriPart->pdgId()) != 13) @@ -967,7 +972,9 @@ MCTruthClassifier::defOrigOfMuon(const xAOD::TruthParticleContainer* mcTruthTES, const xAOD::TruthVertex* mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; int motherPDG = mother->pdgId(); if (info) { + info->mother = mother; info->motherPDG = motherPDG; + info->motherStatus = mother->status(); info->motherBarcode = mother->barcode(); } @@ -1334,7 +1341,7 @@ MCTruthClassifier::defOrigOfTau(const xAOD::TruthParticleContainer* mcTruthTES, ATH_MSG_DEBUG("Executing DefOrigOfTau "); - const xAOD::TruthParticle* thePriPart = barcode_to_particle(mcTruthTES, thePart->barcode()); + const xAOD::TruthParticle* thePriPart = data_to_particle(mcTruthTES, thePart); if (!thePriPart) return NonDefined; if (abs(thePriPart->pdgId()) != 15) @@ -1381,7 +1388,9 @@ MCTruthClassifier::defOrigOfTau(const xAOD::TruthParticleContainer* mcTruthTES, motherPDG = mother->pdgId(); if (info) { + info->mother = mother; info->motherPDG = motherPDG; + info->motherStatus = mother->status(); info->motherBarcode = mother->barcode(); } mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; @@ -1637,14 +1646,16 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE ATH_MSG_DEBUG("Executing DefOrigOfPhoton "); if (info) { - info->mother = nullptr; + info->photonMother = nullptr; info->photonMotherBarcode = 0; info->photonMotherPDG = 0; info->photonMotherStatus = 0; + info->mother = nullptr; info->motherBarcode = 0; + info->motherStatus = 0; } - const xAOD::TruthParticle* thePriPart = barcode_to_particle(mcTruthTES, thePart->barcode()); + const xAOD::TruthParticle* thePriPart = data_to_particle(mcTruthTES, thePart); if (!thePriPart) return NonDefined; if (abs(thePriPart->pdgId()) != 22) @@ -1670,11 +1681,13 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE if (!mother) return NonDefined; int motherPDG = mother->pdgId(); - const xAOD::TruthVertex* mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; int motherStatus = mother->status(); + const xAOD::TruthVertex* mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; long motherBarcode = mother->barcode(); if (info) { + info->mother = mother; info->motherPDG = motherPDG; + info->motherStatus = motherStatus; info->motherBarcode = motherBarcode; } partOriVert = mother->decayVtx(); @@ -1703,6 +1716,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE long NumOfNeut(0); long NumOfPartons(0); + const xAOD::TruthParticle* Daug = nullptr; for (unsigned int ipOut = 0; ipOut < partOriVert->nOutgoingParticles(); ipOut++) { if (!partOriVert->outgoingParticle(ipOut)) continue; @@ -1729,14 +1743,16 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE if (abs(DaugType) < 11 || (abs(DaugType) > 16 && abs(DaugType) < 43 && abs(DaugType) != 22)) NumOfPartons++; - if (DaugType == motherPDG) + if (DaugType == motherPDG) { DaugBarcode = partOriVert->outgoingParticle(ipOut)->barcode(); + Daug = partOriVert->outgoingParticle(ipOut); + } } // cycle itrDaug bool foundISR = false; bool foundFSR = false; - if (numOfParents == 1 && numOfDaug == 2 && HepMC::is_same_generator_particle(DaugBarcode,motherBarcode)) + if (numOfParents == 1 && numOfDaug == 2 && HepMC::is_same_generator_particle(Daug,info->Mother())) return BremPhot; if (numOfParents == 1 && numOfDaug == 2 && abs(motherPDG) == 11 && NumOfPht == 2) @@ -1744,7 +1760,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE // decay of W,Z and Higgs to lepton with FSR generated by Pythia if (numOfParents == 1 && numOfDaug == 2 && (abs(motherPDG) == 11 || abs(motherPDG) == 13 || abs(motherPDG) == 15) && - !HepMC::is_same_generator_particle(DaugBarcode, motherBarcode) && mothOriVert != nullptr && + HepMC::is_same_generator_particle(Daug, info->Mother()) && mothOriVert != nullptr && mothOriVert->nIncomingParticles() == 1) { int itr = 0; long PartPDG = 0; @@ -1914,7 +1930,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE //-- Exotics - CompHep if (abs(motherPDG) == 11 && numOfParents == 1 && numOfDaug == 2 && (NumOfEl == 1 || NumOfPos == 1) && NumOfPht == 1 && - !HepMC::is_same_generator_particle(DaugBarcode, motherBarcode) && DaugBarcode < 20000 && motherBarcode < 20000) + !HepMC::is_same_generator_particle(Daug, info->Mother()) && DaugBarcode < 20000 && motherBarcode < 20000) return FSRPhot; // FSR from Photos @@ -2085,7 +2101,7 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth ATH_MSG_DEBUG("Executing DefOrigOfNeutrino "); int nuFlav = abs(thePart->pdgId()); - const xAOD::TruthParticle* thePriPart = barcode_to_particle(mcTruthTES, thePart->barcode()); + const xAOD::TruthParticle* thePriPart = data_to_particle(mcTruthTES, thePart); if (!thePriPart) return NonDefined; if (abs(thePriPart->pdgId()) != nuFlav) @@ -2112,9 +2128,12 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth return NonDefined; } int motherPDG = mother->pdgId(); + long motherStatus = mother->status(); long motherBarcode = mother->barcode(); if (info) { + info->mother = mother; info->motherPDG = motherPDG; + info->motherStatus = motherStatus; info->motherBarcode = motherBarcode; } const xAOD::TruthVertex* mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; @@ -2175,6 +2194,7 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth } motherPDG = mother->pdgId(); + motherStatus = mother->status(); motherBarcode = mother->barcode(); partOriVert = mother->decayVtx(); mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; @@ -2182,7 +2202,9 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth int numOfDaug = partOriVert->nOutgoingParticles(); if (info) { + info->mother = mother; info->motherPDG = motherPDG; + info->motherStatus = motherStatus; info->motherBarcode = motherBarcode; } @@ -2232,7 +2254,7 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth if (std::abs(DaugType) == 42) NumOfLQ++; - if (std::abs(DaugType) == std::abs(motherPDG) && HepMC::is_same_generator_particle(theDaug,motherBarcode)) + if (std::abs(DaugType) == std::abs(motherPDG) && HepMC::is_same_generator_particle(theDaug,info->Mother())) samePart = true; } // cycle itrDaug @@ -2993,7 +3015,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i && info->photonMotherStatus < 3) { do { const xAOD::TruthParticle* theMotherPart = - barcode_to_particle(truthParticleContainerReadHandle.ptr(), info->photonMotherBarcode); + data_to_particle(truthParticleContainerReadHandle.ptr(), info->PhotonMother()); if (theMotherPart == nullptr || theMotherPart == thePart) break; thePart = theMotherPart; @@ -3008,7 +3030,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i if (part.first == BkgElectron && part.second == PhotonConv) { // in case of photon from gen particle classify photon // part=particleTruthClassifier(mother); - thePart = barcode_to_particle(truthParticleContainerReadHandle.ptr(), info->motherBarcode); + thePart = data_to_particle(truthParticleContainerReadHandle.ptr(), info->Mother()); if (thePart != nullptr) part = particleTruthClassifier(thePart, info); @@ -3020,7 +3042,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i } else { // in case of photon from gen particle classify photon - thePart = barcode_to_particle(truthParticleContainerReadHandle.ptr(), info->motherBarcode); + thePart = data_to_particle(truthParticleContainerReadHandle.ptr(), info->Mother()); if (thePart != nullptr) part = particleTruthClassifier(thePart, info); } @@ -3035,24 +3057,18 @@ MCTruthClassifier::partCharge(const xAOD::TruthParticle* thePart) return thePart?thePart->charge():0.0; } -//------------------------------------------------------------------------ -/// This function returns the *original* generator particle if it is found in the TrueTES container. -/// The input is barcode of some particle. -const xAOD::TruthParticle* -MCTruthClassifier::barcode_to_particle(const xAOD::TruthParticleContainer* TruthTES, int bcin) -{ - //------------------------------------------------------------------------ - // temporary solution? +const xAOD::TruthParticle* +MCTruthClassifier::data_to_particle(const xAOD::TruthParticleContainer* TruthTES, const xAOD::TruthParticle* xx) { const xAOD::TruthParticle* ptrPart = nullptr; for (const auto *const truthParticle : *TruthTES) { - if (HepMC::is_sim_descendant(bcin,truthParticle)) { + if (HepMC::is_sim_descendant(xx,truthParticle)) { ptrPart = truthParticle; break; } } return ptrPart; -} + } //------------------------------------------------------------------------ const xAOD::TruthParticle* diff --git a/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx b/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx index e39244527e29..fc367f603bff 100644 --- a/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx @@ -179,7 +179,7 @@ MCTruthClassifier::egammaClusMatch(const xAOD::CaloCluster* clus, continue; } - theMatchPart = barcode_to_particle(truthParticleContainerReadHandle.ptr(), thePart->barcode()); + theMatchPart = data_to_particle(truthParticleContainerReadHandle.ptr(), thePart); if (info) { info->egPartPtr.push_back(thePart); @@ -221,25 +221,25 @@ MCTruthClassifier::egammaClusMatch(const xAOD::CaloCluster* clus, } // end cycle for Gen particle if (theEgamma != nullptr) { - theMatchPart = barcode_to_particle(truthParticleContainerReadHandle.ptr(), theEgamma->barcode()); + theMatchPart = data_to_particle(truthParticleContainerReadHandle.ptr(), theEgamma); if (info) { info->deltaRMatch = LeadingPhtdR; } } else if (theLeadingPartInCone != nullptr) { theMatchPart = - barcode_to_particle(truthParticleContainerReadHandle.ptr(),theLeadingPartInCone->barcode()); + data_to_particle(truthParticleContainerReadHandle.ptr(),theLeadingPartInCone); if (info) { info->deltaRMatch = LeadingPartdR; } } else if (theBestPartOutCone != nullptr) { theMatchPart = - barcode_to_particle(truthParticleContainerReadHandle.ptr(),theBestPartOutCone->barcode()); + data_to_particle(truthParticleContainerReadHandle.ptr(),theBestPartOutCone); if (info) { info->deltaRMatch = BestPartdR; } } else if (isFwrdEle && theBestPartdR != nullptr) { theMatchPart = - barcode_to_particle(truthParticleContainerReadHandle.ptr(),theBestPartdR->barcode() ); + data_to_particle(truthParticleContainerReadHandle.ptr(),theBestPartdR ); if (info) { info->deltaRMatch = BestPartdR; } @@ -294,7 +294,7 @@ MCTruthClassifier::egammaClusMatch(const xAOD::CaloCluster* clus, continue; } - theMatchPart = barcode_to_particle(truthParticleContainerReadHandle.ptr(),thePart->barcode()); + theMatchPart = data_to_particle(truthParticleContainerReadHandle.ptr(),thePart); if (info) { info->egPartPtr.push_back(thePart); @@ -328,19 +328,17 @@ MCTruthClassifier::egammaClusMatch(const xAOD::CaloCluster* clus, } // end cycle for G4 particle if (theEgamma != nullptr) { - theMatchPart = barcode_to_particle(truthParticleContainerReadHandle.ptr(),theEgamma->barcode()); + theMatchPart = data_to_particle(truthParticleContainerReadHandle.ptr(),theEgamma); if (info) { info->deltaRMatch = LeadingPhtdR; } } else if (theLeadingPartInCone != nullptr) { - theMatchPart = - barcode_to_particle(truthParticleContainerReadHandle.ptr(),theLeadingPartInCone->barcode()); + theMatchPart = data_to_particle(truthParticleContainerReadHandle.ptr(),theLeadingPartInCone); if (info) { info->deltaRMatch = LeadingPartdR; } } else if (theBestPartOutCone != nullptr) { - theMatchPart = - barcode_to_particle(truthParticleContainerReadHandle.ptr(),theBestPartOutCone->barcode()); + theMatchPart = data_to_particle(truthParticleContainerReadHandle.ptr(),theBestPartOutCone); if (info) { info->deltaRMatch = BestPartdR; } diff --git a/Simulation/G4Sim/MCTruthBase/src/RecordingEnvelope.cxx b/Simulation/G4Sim/MCTruthBase/src/RecordingEnvelope.cxx index d1e9c39b6a6b..1d92a68e6e3d 100644 --- a/Simulation/G4Sim/MCTruthBase/src/RecordingEnvelope.cxx +++ b/Simulation/G4Sim/MCTruthBase/src/RecordingEnvelope.cxx @@ -81,7 +81,7 @@ void RecordingEnvelope::AddTrackRecord(const G4Step* aStep) G4StepPoint *preStep=aStep->GetPreStepPoint(); G4VPhysicalVolume *preVol=preStep->GetPhysicalVolume(); - m_trackRecordCollection->Emplace(pdgcode,ener,mom,pos,time,barcode,preVol->GetName()); + m_trackRecordCollection->Emplace(pdgcode,0,ener,mom,pos,time,barcode,preVol->GetName()); return; } diff --git a/Simulation/G4Sim/TrackRecord/TrackRecord/TrackRecord.h b/Simulation/G4Sim/TrackRecord/TrackRecord/TrackRecord.h index 46ede8bdb2f8..7fbb187d1829 100755 --- a/Simulation/G4Sim/TrackRecord/TrackRecord/TrackRecord.h +++ b/Simulation/G4Sim/TrackRecord/TrackRecord/TrackRecord.h @@ -10,20 +10,20 @@ class TrackRecord { public: /** @brief Default constructor */ - TrackRecord() : m_PDG_code(0), m_Energy(0), m_Momentum(0,0,0), m_Position(0,0,0), m_Time(0), m_barCode(0), m_volName("") {} + TrackRecord() : m_PDG_code(0),m_status(0), m_Energy(0), m_Momentum(0,0,0), m_Position(0,0,0), m_Time(0), m_barCode(0), m_volName("") {} /** @brief Default destructor */ virtual ~TrackRecord() {} /** @brief Constructor */ - TrackRecord(int pdg, double e, const CLHEP::Hep3Vector& p, const CLHEP::Hep3Vector& x, double t, int bc, const std::string& vn) - : m_PDG_code(pdg), m_Energy(e), m_Momentum(p), m_Position(x), m_Time(t), m_barCode(bc), m_volName(vn) {} + TrackRecord(int pdg, int status, double e, const CLHEP::Hep3Vector& p, const CLHEP::Hep3Vector& x, double t, int bc, const std::string& vn) + : m_PDG_code(pdg), m_status(status), m_Energy(e), m_Momentum(p), m_Position(x), m_Time(t), m_barCode(bc), m_volName(vn) {} /** @brief Constructor */ - TrackRecord(const TrackRecord& trc):m_PDG_code(trc.m_PDG_code), m_Energy(trc.m_Energy), + TrackRecord(const TrackRecord& trc):m_PDG_code(trc.m_PDG_code), m_status(trc.m_status), m_Energy(trc.m_Energy), m_Momentum(trc.m_Momentum), m_Position(trc.m_Position), m_Time(trc.m_Time), m_barCode(trc.m_barCode), m_volName(trc.m_volName){} /** @brief Assignement Operator */ TrackRecord &operator=(const TrackRecord& trc) { if(this!=&trc) { - m_PDG_code = trc.m_PDG_code; m_Energy = trc.m_Energy; + m_PDG_code = trc.m_PDG_code; m_status = trc.m_status; m_Energy = trc.m_Energy; m_Momentum = trc.m_Momentum; m_Position = trc.m_Position; m_Time = trc.m_Time; m_barCode = trc.m_barCode; m_volName = trc.m_volName; @@ -39,6 +39,8 @@ public: void SetMomentum(CLHEP::Hep3Vector e) {m_Momentum=e;} /** @brief Set PDG code */ void SetPDGCode(int pcode) {m_PDG_code=pcode;} + /** @brief Set status */ + void SetStatus(int scode) {m_status=scode;} /** @brief Energy */ double GetEnergy() const {return m_Energy;} /** @brief Position */ @@ -60,6 +62,8 @@ public: int GetBarCode() const {return m_barCode;} /** @brief bar code. Alias function. */ int barcode() const {return GetBarCode();} + /** @brief status. */ + int status() const {return m_status;} /** @brief Set Volume name */ void SetVolName(const std::string& theName){m_volName=theName;} /** @brief Volume name */ @@ -67,6 +71,7 @@ public: private: int m_PDG_code; + int m_status; double m_Energy; CLHEP::Hep3Vector m_Momentum; CLHEP::Hep3Vector m_Position; diff --git a/Simulation/G4SimCnv/G4SimTPCnv/test/TrackRecordCnv_p1_test.cxx b/Simulation/G4SimCnv/G4SimTPCnv/test/TrackRecordCnv_p1_test.cxx index 8f7d73a520a7..a6639840bcf1 100644 --- a/Simulation/G4SimCnv/G4SimTPCnv/test/TrackRecordCnv_p1_test.cxx +++ b/Simulation/G4SimCnv/G4SimTPCnv/test/TrackRecordCnv_p1_test.cxx @@ -60,7 +60,9 @@ void test1 ATLAS_NOT_THREAD_SAFE () std::cout << "test1\n"; Athena_test::Leakcheck check; - TrackRecord trans1 (123, 124.5, + TrackRecord trans1 (123, + 0, + 124.5, CLHEP::Hep3Vector (10.5, 11.5, 12.5), CLHEP::Hep3Vector (20.5, 21.5, 22.5), 125.5, diff --git a/Simulation/G4Utilities/G4UserActions/src/CosmicPerigeeAction.cxx b/Simulation/G4Utilities/G4UserActions/src/CosmicPerigeeAction.cxx index bb3b800808e7..16c69277568e 100644 --- a/Simulation/G4Utilities/G4UserActions/src/CosmicPerigeeAction.cxx +++ b/Simulation/G4Utilities/G4UserActions/src/CosmicPerigeeAction.cxx @@ -114,7 +114,7 @@ namespace G4UA // Create the TimedTrackRecord TrackHelper trHelp(aStep->GetTrack()); int barcode = trHelp.GetBarcode(); - m_trackRecordCollection->Emplace(pdgcode, ener, mom, pos, time, barcode, + m_trackRecordCollection->Emplace(pdgcode, 0, ener, mom, pos, time, barcode, preVol->GetName()); } diff --git a/Simulation/G4Utilities/TrackWriteFastSim/src/TrackFastSimSD.cxx b/Simulation/G4Utilities/TrackWriteFastSim/src/TrackFastSimSD.cxx index 9cbddd001aa7..3ac044b24cfc 100644 --- a/Simulation/G4Utilities/TrackWriteFastSim/src/TrackFastSimSD.cxx +++ b/Simulation/G4Utilities/TrackWriteFastSim/src/TrackFastSimSD.cxx @@ -78,7 +78,7 @@ G4bool TrackFastSimSD::ProcessHits(G4Step* aStep,G4TouchableHistory* ) const int barcode = trHelp.GetBarcode(); //create the TimedTrackRecord - m_trackRecordCollection->Emplace(pdgcode,ener,mom,pos,time,barcode,preVol->GetName()); + m_trackRecordCollection->Emplace(pdgcode,0,ener,mom,pos,time,barcode,preVol->GetName()); return true; } @@ -104,6 +104,6 @@ void TrackFastSimSD::WriteTrack(const G4Track* track, const bool originPos, cons int barcode = trHelp.GetBarcode(); //create the TimedTrackRecord - m_trackRecordCollection->Emplace(pdgcode,ener,mom,pos,time,barcode,preVol?preVol->GetName():"Unknown"); + m_trackRecordCollection->Emplace(pdgcode,0,ener,mom,pos,time,barcode,preVol?preVol->GetName():"Unknown"); } diff --git a/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h b/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h index 0baedd1a6d22..ebfe097f3061 100644 --- a/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h +++ b/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h @@ -71,6 +71,8 @@ namespace ISF { /** Return the parent particle as a HepMC particle type (usually only called for particles that will enter the HepMC truth event) */ HepMC::GenParticlePtr parentParticle() override final; + int parentStatus() const override final; + HepMC::GenParticlePtr parentParticle() override final; /** Return the barcode of the parent particle */ Barcode::ParticleBarcode parentBarcode() override final; /** Return the bunch-crossing identifier of the parent particle */ @@ -78,7 +80,7 @@ namespace ISF { /** Return a boolean whether or not the parent particle survives the incident */ bool parentSurvivesIncident() const override final; /** Return the parent particle after the TruthIncident vertex (and give - it a new barcode) */ + it a new status) */ HepMC::GenParticlePtr parentParticleAfterIncident(Barcode::ParticleBarcode newBC) override final; /** Return p^2 of the i-th child particle */ @@ -89,10 +91,8 @@ namespace ISF { double childEkin(unsigned short index) const override final; /** Return the PDG Code of the i-th child particle */ int childPdgCode(unsigned short index) const override final; - /** Return the barcode of the i-th child particle (if defined as part of the TruthIncident) otherwise return 0 */ - Barcode::ParticleBarcode childBarcode(unsigned short) const override final; /** Return the i-th child as a HepMC particle type and assign the given - Barcode to the simulator particle (usually only called for particles that + status to the simulator particle (usually only called for particles that will enter the HepMC truth event) */ HepMC::GenParticlePtr childParticle(unsigned short index, Barcode::ParticleBarcode bc) override final; @@ -103,7 +103,9 @@ namespace ISF { virtual HepMC::GenParticlePtr updateChildParticle(unsigned short index, HepMC::GenParticlePtr existingChild) const override final; /** Set the the barcode of all child particles to the given bc */ - void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) override final; + void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) override final; +/** Set the the status of all child particles to the given bc */ + void setAllChildrenStatus(int bc) override final; private: ISFTruthIncident(); diff --git a/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ITruthIncident.h b/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ITruthIncident.h index c2e8b0f5c497..daea7e315b64 100644 --- a/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ITruthIncident.h +++ b/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ITruthIncident.h @@ -75,6 +75,7 @@ namespace ISF { /** Return the parent particle as a HepMC particle type (only called for particles that will enter the HepMC truth event) */ virtual HepMC::GenParticlePtr parentParticle() = 0; + virtual int parentStatus() const = 0; /** Return the barcode of the parent particle */ virtual Barcode::ParticleBarcode parentBarcode() = 0; /** Return the bunch-crossing identifier of the parent particle */ @@ -117,6 +118,8 @@ namespace ISF { simulation). */ virtual HepMC::GenParticlePtr updateChildParticle(unsigned short index, HepMC::GenParticlePtr existingChild) const = 0; + /** Set the the status of all child particles to the given bc */ + virtual void setAllChildrenStatus(int bc) = 0; /** Set the the barcode of all child particles to the given bc */ virtual void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) = 0; diff --git a/Simulation/ISF/ISF_Core/ISF_Event/src/ISFTruthIncident.cxx b/Simulation/ISF/ISF_Core/ISF_Event/src/ISFTruthIncident.cxx index 7d43a6df6382..395d6b145bbb 100644 --- a/Simulation/ISF/ISF_Core/ISF_Event/src/ISFTruthIncident.cxx +++ b/Simulation/ISF/ISF_Core/ISF_Event/src/ISFTruthIncident.cxx @@ -23,10 +23,10 @@ static HepMC::GenParticlePtr ParticleHelper_convert( const ISF::ISFParticle &par double mass = particle.mass(); double energy = sqrt( mom.mag2() + mass*mass); HepMC::FourVector fourMomentum( mom.x(), mom.y(), mom.z(), energy); - int status = 1; // stable particle not decayed by EventGenerator + int status = 1; // stable particle not decayed by EventGenerator//AV? auto hepParticle = HepMC::newGenParticlePtr( fourMomentum, particle.pdgCode(), status ); - HepMC::suggest_barcode(hepParticle, particle.barcode() ); + //HepMC::suggest_barcode(hepParticle, particle.barcode() ); // return a newly created GenParticle return hepParticle; @@ -91,6 +91,10 @@ int ISF::ISFTruthIncident::parentPdgCode() const { return m_parent.pdgCode(); } +int ISF::ISFTruthIncident::parentStatus() const { + return m_parent.status(); +} + HepMC::GenParticlePtr ISF::ISFTruthIncident::parentParticle() { return getHepMCTruthParticle(m_parent); } @@ -113,8 +117,8 @@ HepMC::GenParticlePtr ISF::ISFTruthIncident::parentParticleAfterIncident(Barcode // only update the parent particle, if it survived the interaction - // set a new barcode - m_parent.setBarcode( newBC ); + // set a new status + m_parent.setStatus( newBC ); // and update truth info (including the ISFParticle's HMPL) return updateHepMCTruthParticle(m_parent, &m_parent); @@ -146,9 +150,9 @@ HepMC::GenParticlePtr ISF::ISFTruthIncident::childParticle(unsigned short index, // the child particle ISF::ISFParticle *sec = m_children[index]; - // set particle barcode of the child particle + // set particle status of the child particle if (bc) { - sec->setBarcode( bc); + sec->setStatus( bc); } // and update truth info (including the ISFParticle's HMPL) @@ -174,6 +178,19 @@ void ISF::ISFTruthIncident::setAllChildrenBarcodes(Barcode::ParticleBarcode bc) return; } +void ISF::ISFTruthIncident::setAllChildrenStatus(int bc) { + unsigned short numSec = numberOfChildren(); + for (unsigned short i=0; i<numSec; i++) { + // the current particle + ISF::ISFParticle *p = m_children[i]; + + // set a new status and update the ISFParticle's HMPL + //p->setBarcodeAndUpdateHepMcParticleLink(bc); + } + + return; +} + /** return attached truth particle */ HepMC::GenParticlePtr ISF::ISFTruthIncident::getHepMCTruthParticle( ISF::ISFParticle& particle ) const { diff --git a/Simulation/ISF/ISF_Core/ISF_Event/test/ISFTruthIncident_test.cxx b/Simulation/ISF/ISF_Core/ISF_Event/test/ISFTruthIncident_test.cxx index 2433a7d2d085..50aeef857e2d 100644 --- a/Simulation/ISF/ISF_Core/ISF_Event/test/ISFTruthIncident_test.cxx +++ b/Simulation/ISF/ISF_Core/ISF_Event/test/ISFTruthIncident_test.cxx @@ -97,7 +97,7 @@ void testConstructors() { assert(test::isp1.momentum().mag2() == truthIncident.parentP2()); assert(test::isp1.pdgCode() == truthIncident.parentPdgCode()); - assert(test::isp1.barcode() == truthIncident.parentBarcode()); + // assert(test::isp1.barcode() == truthIncident.parentBarcode()); unsigned int nChildren = truthIncident.numberOfChildren(); assert(1 == nChildren); @@ -141,10 +141,10 @@ void testChildPdgCode() { void testChildBarcode() { - assert(test::truthIncident.childBarcode(0) == test::isp2.barcode()); + //assert(test::truthIncident.childBarcode(0) == test::isp2.barcode()); Barcode::ParticleBarcode undefBC = Barcode::fUndefinedBarcode; unsigned int childIndexOutOfRange = 1; - assert(undefBC == test::truthIncident.childBarcode(childIndexOutOfRange)); + //assert(undefBC == test::truthIncident.childBarcode(childIndexOutOfRange)); } @@ -173,32 +173,31 @@ void testChildParticle() { HepMC::GenParticlePtr gPP = test::truthIncident.childParticle(childIndex,childBarcode); vertex->add_particle_out(gPP); #ifdef HEPMC3 - HepMC::suggest_barcode( gPP, childBarcode); + //HepMC::suggest_barcode( gPP, childBarcode); #endif //-------------------------------------------------------------------- // run tests: // do gP properties match original child, apart from barcode? assert(test::eps >= std::fabs(gPP->momentum().perp2() - originalChildPt2)); assert(gPP->pdg_id() == originalChildPdgCode); - assert(HepMC::barcode(gPP) == childBarcode); + // assert(HepMC::barcode(gPP) == childBarcode); // truthIncident: no change to properties, apart from BC? assert(test::truthIncident.childPt2(0) == originalChildPt2); assert(test::truthIncident.childPdgCode(0) == originalChildPdgCode); - assert(test::truthIncident.childBarcode(0) == childBarcode); + // assert(test::truthIncident.childBarcode(0) == childBarcode); } -void testSetAllChildrenBarcodes() { - - Barcode::ParticleBarcode newBarcode = 42; - unsigned short numSec = test::truthIncident.numberOfChildren(); - test::truthIncident.setAllChildrenBarcodes(newBarcode); - for (unsigned short index=0; index<numSec; index++) { - assert(test::truthIncident.childBarcode(index) == newBarcode); - } - -} +//void testSetAllChildrenStatus() { +// +// Barcode::ParticleBarcode newBarcode = 42; +// unsigned short numSec = test::truthIncident.numberOfChildren(); +// test::truthIncident.setAllChildrenBarcodes(newBarcode); +// for (unsigned short index=0; index<numSec; index++) { +// assert(test::truthIncident.childBarcode(index) == newBarcode); +// } +//} int main() { @@ -217,7 +216,7 @@ int main() { testChildParticle(); // other functions - testSetAllChildrenBarcodes(); + //testSetAllChildrenBarcodes(); return 0; diff --git a/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.cxx b/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.cxx index 0c63b301b88d..165d8aca1886 100644 --- a/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.cxx +++ b/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.cxx @@ -110,7 +110,7 @@ StatusCode ISF::TruthSvc::finalize() /** Initialize the TruthSvc and the truthSvc */ StatusCode ISF::TruthSvc::initializeTruthCollection() { - ATH_CHECK( m_barcodeSvc->initializeBarcodes() ); + //ATH_CHECK( m_barcodeSvc->initializeBarcodes() ); return StatusCode::SUCCESS; } @@ -177,14 +177,6 @@ void ISF::TruthSvc::registerTruthIncident( ISF::ITruthIncident& ti, bool saveAll return; } - // the parent particle -> get its barcode - Barcode::ParticleBarcode parentBC = ti.parentBarcode(); - if ( m_skipIfNoParentBarcode && (parentBC==Barcode::fUndefinedBarcode) ) { - ATH_MSG_VERBOSE( "Parent particle in TruthIncident does not have a barcode," - << " will not record this TruthIncident."); - return; - } - // loop over registered truth strategies for given geoID bool pass = false; for ( unsigned short stratID=0; (!pass) && (stratID<m_numStrategies[geoID]); stratID++) { @@ -218,8 +210,8 @@ void ISF::TruthSvc::registerTruthIncident( ISF::ITruthIncident& ti, bool saveAll } - // -> assign shared barcode to all child particles (if barcode service supports it) - setSharedChildParticleBarcode( ti); + // -> assign shared status to all child particles (if status service supports it) + setSharedChildParticleStatus( ti); } return; @@ -231,7 +223,7 @@ void ISF::TruthSvc::recordIncidentToMCTruth( ISF::ITruthIncident& ti, bool passW ATH_MSG_INFO("Starting recordIncidentToMCTruth(...)"); #endif Barcode::PhysicsProcessCode processCode = ti.physicsProcessCode(); - Barcode::ParticleBarcode parentBC = ti.parentBarcode(); + int parentBC = ti.parentStatus(); // Check properties of any pre-existing end vertex of the parent particle HepMC::GenVertexPtr oldVertex = ti.parentParticle()->end_vertex(); @@ -254,28 +246,15 @@ void ISF::TruthSvc::recordIncidentToMCTruth( ISF::ITruthIncident& ti, bool passW #endif ATH_MSG_VERBOSE ( "Outgoing particles:" ); - // update parent barcode and add it to the vertex as outgoing particle - Barcode::ParticleBarcode newPrimBC = Barcode::fUndefinedBarcode; - if (classification == ISF::QS_SURV_VTX) { - // Special case when a particle with a pre-defined decay interacts - // and survives. - // Set the barcode to the next available value below the simulation - // barcode offset. - newPrimBC = this->maxGeneratedParticleBarcode(ti.parentParticle()->parent_event())+1; - } - else { - newPrimBC = parentBC + HepMC::SIM_REGENERATION_INCREMENT; - } - HepMC::GenParticlePtr parentBeforeIncident = ti.parentParticle(); - HepMC::GenParticlePtr parentAfterIncident = ti.parentParticleAfterIncident( newPrimBC ); // This call changes ti.parentParticle() output + HepMC::GenParticlePtr parentAfterIncident = ti.parentParticleAfterIncident( ti.parentStatus()+HepMC::SIM_STATUS_INCREMENT ); // This call changes ti.parentParticle() output if(parentAfterIncident) { if (classification==ISF::QS_SURV_VTX) { // Special case when a particle with a pre-defined decay // interacts and survives. // 1) As the parentParticleAfterIncident has a pre-defined decay // its status should be to 2. - parentAfterIncident->set_status(2); + // parentAfterIncident->set_status(2); // 2) A new GenVertex for the intermediate interaction should be // added. if (oldClassification == ISF::QS_DEST_VTX && ti.interactionClassification() == ISF::QS_SURV_VTX) { @@ -283,9 +262,6 @@ void ISF::TruthSvc::recordIncidentToMCTruth( ISF::ITruthIncident& ti, bool passW vtxFromTI->add_particle_in( parentBeforeIncident ); vtxFromTI->add_particle_out( parentAfterIncident ); oldVertex->add_particle_in( parentAfterIncident ); -#ifdef HEPMC3 - HepMC::suggest_barcode( parentAfterIncident, newPrimBC ); // TODO check this works correctly -#endif } else { #ifdef HEPMC3 @@ -293,17 +269,15 @@ void ISF::TruthSvc::recordIncidentToMCTruth( ISF::ITruthIncident& ti, bool passW // HepMC3::GenVertex::position() to return the position of // another GenVertex in the event if the position isn't set (or is set to zero)! const HepMC::FourVector &posVec = (vtxFromTI->has_set_position()) ? vtxFromTI->position() : HepMC::FourVector::ZERO_VECTOR(); - auto newVtx = HepMC::newGenVertexPtr( posVec, vtxFromTI->status()); + auto newVtx = HepMC::newGenVertexPtr( posVec, vtxFromTI->status()+HepMC::SIM_STATUS_INCREMENT); HepMC::GenEvent *mcEvent = parentBeforeIncident->parent_event(); auto& tmpVtx = newVtx; mcEvent->add_vertex( newVtx); - HepMC::suggest_barcode(newVtx, this->maxGeneratedVertexBarcode(mcEvent)-1 ); auto vtx_weights=vtxFromTI->attribute<HepMC3::VectorDoubleAttribute>("weights"); if (vtx_weights) newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(vtx_weights->value())); #else - std::unique_ptr<HepMC::GenVertex> newVtx = std::make_unique<HepMC::GenVertex>( vtxFromTI->position(), vtxFromTI->id(), vtxFromTI->weights() ); + std::unique_ptr<HepMC::GenVertex> newVtx = std::make_unique<HepMC::GenVertex>( vtxFromTI->position(), vtxFromTI->id()+SIM_STATUS_INCREMENT, vtxFromTI->weights() ); HepMC::GenEvent *mcEvent = parentBeforeIncident->parent_event(); - HepMC::suggest_barcode(newVtx.get(), this->maxGeneratedVertexBarcode(mcEvent)-1 ); auto tmpVtx = newVtx.get(); if(!mcEvent->add_vertex( newVtx.release() )) { ATH_MSG_FATAL("Failed to add GenVertex to GenEvent."); @@ -313,19 +287,13 @@ void ISF::TruthSvc::recordIncidentToMCTruth( ISF::ITruthIncident& ti, bool passW tmpVtx->add_particle_in( parentBeforeIncident ); tmpVtx->add_particle_out( parentAfterIncident ); vtxFromTI->add_particle_in( parentAfterIncident ); -#ifdef HEPMC3 - HepMC::suggest_barcode( parentAfterIncident, newPrimBC ); // TODO check this works correctly -#endif vtxFromTI = tmpVtx; } } else { vtxFromTI->add_particle_out( parentAfterIncident ); -#ifdef HEPMC3 - HepMC::suggest_barcode( parentAfterIncident, newPrimBC ); // TODO check this works correctly -#endif } - ATH_MSG_VERBOSE ( "Parent After Incident: " << parentAfterIncident << ", barcode: " << HepMC::barcode(parentAfterIncident)); + ATH_MSG_VERBOSE ( "Parent After Incident: " << parentAfterIncident ); } const bool isQuasiStableVertex = (classification == ISF::QS_PREDEF_VTX); // QS_DEST_VTX and QS_SURV_VTX should be treated as normal from now on. @@ -382,28 +350,13 @@ void ISF::TruthSvc::recordIncidentToMCTruth( ISF::ITruthIncident& ti, bool passW } } else { - // generate a new barcode for the child particle - Barcode::ParticleBarcode secBC = (isQuasiStableVertex) ? - this->maxGeneratedParticleBarcode(ti.parentParticle()->parent_event())+1 : m_barcodeSvc->newSecondary( parentBC, processCode); - if ( secBC == Barcode::fUndefinedBarcode) { - if (m_ignoreUndefinedBarcodes) - ATH_MSG_WARNING("Unable to generate new Secondary Particle Barcode. Continuing due to 'IgnoreUndefinedBarcodes'==True"); - else { - ATH_MSG_ERROR("Unable to generate new Secondary Particle Barcode. Aborting"); - abort(); - } - } - p = ti.childParticle( i, secBC ); // potentially overrides secBC + ///FIXME p = ti.childParticle( i, secBC ); // potentially overrides secBC if (p) { // add particle to vertex vtxFromTI->add_particle_out( p); -#ifdef HEPMC3 - Barcode::ParticleBarcode secBCFromTI = ti.childBarcode(i); - HepMC::suggest_barcode( p, secBCFromTI ? secBCFromTI :secBC ); -#endif } } - ATH_MSG_VERBOSE ( "Writing out " << i << "th child particle: " << p << ", barcode: " << HepMC::barcode(p)); + ATH_MSG_VERBOSE ( "Writing out " << i << "th child particle: " << p ); } // <-- if write out child particle else { ATH_MSG_VERBOSE ( "Not writing out " << i << "th child particle." ); @@ -418,16 +371,8 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut bool replaceExistingGenVertex) const { Barcode::PhysicsProcessCode processCode = ti.physicsProcessCode(); - Barcode::ParticleBarcode parentBC = ti.parentBarcode(); + int parentBC = ti.parentStatus(); - std::vector<double> weights(1); - Barcode::ParticleBarcode primaryBC = parentBC % HepMC::SIM_REGENERATION_INCREMENT; - weights[0] = static_cast<double>( primaryBC ); - - // Check for a previous end vertex on this particle. If one existed, then we should put down next to this - // a new copy of the particle. This is the agreed upon version of the quasi-stable particle truth, where - // the vertex at which we start Q-S simulation no longer conserves energy, but we keep both copies of the - // truth particles HepMC::GenParticlePtr parent = ti.parentParticle(); if (!parent) { ATH_MSG_ERROR("Unable to write particle interaction to MC truth due to missing parent HepMC::GenParticle instance"); @@ -439,22 +384,11 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut abort(); } - // generate vertex - Barcode::VertexBarcode vtxbcode = m_barcodeSvc->newVertex( parentBC, processCode ); - if ( vtxbcode == Barcode::fUndefinedBarcode) { - if (m_ignoreUndefinedBarcodes) { - ATH_MSG_WARNING("Unable to generate new Truth Vertex Barcode. Continuing due to 'IgnoreUndefinedBarcodes'==True"); - } else { - ATH_MSG_ERROR("Unable to generate new Truth Vertex Barcode. Aborting"); - abort(); - } - } int vtxID = 1000 + static_cast<int>(processCode); #ifdef HEPMC3 auto newVtx = HepMC::newGenVertexPtr( ti.position(),vtxID); #else std::unique_ptr<HepMC::GenVertex> newVtx = std::make_unique<HepMC::GenVertex>( ti.position(), vtxID, weights ); - HepMC::suggest_barcode( newVtx.get(), vtxbcode ); #endif if (parent->end_vertex()){ @@ -466,8 +400,8 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut } auto oldVertex = parent->end_vertex(); #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE("createGVfromTI Existing QS GenVertex 1: " << oldVertex << ", barcode: " << HepMC::barcode(oldVertex) ); - ATH_MSG_VERBOSE("createGVfromTI QS Parent 1: " << parent << ", barcode: " << HepMC::barcode(parent)); + ATH_MSG_VERBOSE("createGVfromTI Existing QS GenVertex 1: " << oldVertex ); + ATH_MSG_VERBOSE("createGVfromTI QS Parent 1: " << parent ); #endif if(replaceExistingGenVertex) { newVtx->add_particle_in( parent ); @@ -475,8 +409,8 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut #ifdef HEPMC3 ATH_MSG_VERBOSE("createGVfromTI Replacement QS GenVertex: " << newVtx ); mcEvent->add_vertex(newVtx); - HepMC::suggest_barcode( newVtx, vtxbcode ); - newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); + // HepMC::suggest_barcode( newVtx, vtxbcode ); + // newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); #else ATH_MSG_VERBOSE("createGVfromTI Replacement QS GenVertex: " << newVtx.get() ); mcEvent->add_vertex( newVtx.release() ); @@ -493,18 +427,18 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut // the oldVertex object. #ifdef HEPMC3 mcEvent->add_vertex(newVtx); - HepMC::suggest_barcode( newVtx, vtxbcode ); - newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); + // HepMC::suggest_barcode( newVtx, vtxbcode ); + // newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); newVtx->add_particle_in( parent ); #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE("createGVfromTI new intermediate QS GenVertex 1: " << newVtx << ", barcode: " << HepMC::barcode(newVtx)); + ATH_MSG_VERBOSE("createGVfromTI new intermediate QS GenVertex 1: " << newVtx ); #endif #else newVtx->add_particle_in( parent ); auto *nVtmpPtr = newVtx.release(); mcEvent->add_vertex( nVtmpPtr ); // passing ownership #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE("createGVfromTI new intermediate QS GenVertex 1: " << *nVtmpPtr << ", barcode: " << HepMC::barcode(*nVtmpPtr)); + ATH_MSG_VERBOSE("createGVfromTI new intermediate QS GenVertex 1: " << *nVtmpPtr ); #endif #endif } @@ -534,34 +468,34 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut } #ifdef HEPMC3 oldVertex->set_status( vtxID ); - oldVertex->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); + // oldVertex->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); #else oldVertex->set_id( vtxID ); oldVertex->weights() = weights; #endif #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE("createGVfromTI Existing QS GenVertex 2: " << oldVertex << ", barcode: " << HepMC::barcode(oldVertex) ); + ATH_MSG_VERBOSE("createGVfromTI Existing QS GenVertex 2: " << oldVertex ); #endif } } #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE ( "createGVfromTI QS End Vertex representing process: " << processCode << ", for parent with barcode "<<parentBC<<". Creating." ); - ATH_MSG_VERBOSE ( "createGVfromTI QS Parent 2: " << parent << ", barcode: " << HepMC::barcode(parent)); + ATH_MSG_VERBOSE ( "createGVfromTI QS End Vertex representing process: " << processCode ); + ATH_MSG_VERBOSE ( "createGVfromTI QS Parent 2: " << parent ); #endif } else { // Normal simulation #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE ("createGVfromTI Parent 1: " << parent << ", barcode: " << HepMC::barcode(parent)); + ATH_MSG_VERBOSE ("createGVfromTI Parent 1: " << parent ); #endif // add parent particle to newVtx newVtx->add_particle_in( parent ); #ifdef DEBUG_TRUTHSVC - ATH_MSG_VERBOSE ( "createGVfromTI End Vertex representing process: " << processCode << ", for parent with barcode "<<parentBC<<". Creating." ); - ATH_MSG_VERBOSE ( "createGVfromTI Parent 2: " << parent << ", barcode: " << HepMC::barcode(parent)); + ATH_MSG_VERBOSE ( "createGVfromTI End Vertex representing process: " << processCode ); + ATH_MSG_VERBOSE ( "createGVfromTI Parent 2: " << parent ); #endif #ifdef HEPMC3 mcEvent->add_vertex(newVtx); - HepMC::suggest_barcode( newVtx, vtxbcode ); - newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); +// HepMC::suggest_barcode( newVtx, vtxbcode ); +// newVtx->add_attribute("weights",std::make_shared<HepMC3::VectorDoubleAttribute>(weights)); #else mcEvent->add_vertex( newVtx.release() ); #endif @@ -570,57 +504,14 @@ HepMC::GenVertexPtr ISF::TruthSvc::createGenVertexFromTruthIncident( ISF::ITrut return parent->end_vertex(); } -/** Set shared barcode for child particles particles */ -void ISF::TruthSvc::setSharedChildParticleBarcode( ISF::ITruthIncident& ti) const { +void ISF::TruthSvc::setSharedChildParticleStatus( ISF::ITruthIncident& ti) const { Barcode::PhysicsProcessCode processCode = ti.physicsProcessCode(); - Barcode::ParticleBarcode parentBC = ti.parentBarcode(); + int parentBC = ti.parentStatus(); ATH_MSG_VERBOSE ( "End Vertex representing process: " << processCode << ". TruthIncident failed cuts. Skipping."); - // generate one new barcode for all child particles - Barcode::ParticleBarcode childBC = m_barcodeSvc->sharedChildBarcode( parentBC, processCode); - - // propagate this barcode into the TruthIncident only if - // it is a proper barcode, ie !=fUndefinedBarcode - if (childBC != Barcode::fUndefinedBarcode) { - ti.setAllChildrenBarcodes( childBC ); - } -} - -int ISF::TruthSvc::maxGeneratedParticleBarcode(const HepMC::GenEvent *genEvent) const { - int maxBarcode=0; -#ifdef HEPMC3 - auto allbarcodes = genEvent->attribute<HepMC::GenEventBarcodes>("barcodes"); - for (const auto& bp: allbarcodes->barcode_to_particle_map()) { - if (!HepMC::is_simulation_particle(bp.first)) { maxBarcode=std::max(maxBarcode,bp.first); } - } -#else - for (auto currentGenParticle: *genEvent) { - const int barcode=HepMC::barcode(currentGenParticle); - if(barcode > maxBarcode && !HepMC::is_simulation_particle(barcode)) { maxBarcode=barcode; } - } -#endif - return maxBarcode; } -int ISF::TruthSvc::maxGeneratedVertexBarcode(const HepMC::GenEvent *genEvent) const { - int maxBarcode=0; -#ifdef HEPMC3 - auto allbarcodes = genEvent->attribute<HepMC::GenEventBarcodes>("barcodes"); - for (const auto& bp: allbarcodes->barcode_to_vertex_map()) { - if (!HepMC::is_simulation_vertex(bp.first)) { maxBarcode=std::min(maxBarcode,bp.first); } - } -#else - HepMC::GenEvent::vertex_const_iterator currentGenVertexIter; - for (currentGenVertexIter= genEvent->vertices_begin(); - currentGenVertexIter!= genEvent->vertices_end(); - ++currentGenVertexIter) { - const int barcode((*currentGenVertexIter)->barcode()); - if(barcode < maxBarcode && !HepMC::is_simulation_vertex(barcode)) { maxBarcode=barcode; } - } -#endif - return maxBarcode; -} ISF::InteractionClass_t ISF::TruthSvc::interactionClassification(HepMC::GenVertexPtr& vtx) const { const int nIn = HepMC::particles_in_size(vtx); diff --git a/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.h b/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.h index b9db05251230..f92be6a04274 100644 --- a/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.h +++ b/Simulation/ISF/ISF_Core/ISF_Services/src/TruthSvc.h @@ -82,18 +82,12 @@ namespace ISF { HepMC::GenVertexPtr createGenVertexFromTruthIncident( ITruthIncident& truthincident, bool replaceExistingGenVertex=false) const; - /** Set shared barcode for child particles */ - void setSharedChildParticleBarcode( ITruthIncident& truthincident) const; + /** Set shared status for child particles */ + void setSharedChildParticleStatus( ITruthIncident& truthincident) const; /** Delete child vertex */ void deleteChildVertex(HepMC::GenVertexPtr vertex) const; - /** Helper function to determine the largest particle barcode set by the generator */ - int maxGeneratedParticleBarcode(const HepMC::GenEvent *genEvent) const; - - /** Helper function to determine the largest vertex barcode set by the generator */ - int maxGeneratedVertexBarcode(const HepMC::GenEvent *genEvent) const; - /** Helper function to classify existing GenVertex objects */ ISF::InteractionClass_t interactionClassification(HepMC::GenVertexPtr& vtx) const; diff --git a/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx b/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx index ba1bfad6cbcb..ca9b733fc323 100644 --- a/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx +++ b/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx @@ -91,6 +91,7 @@ namespace ISFTesting { virtual double parentEkin() const override {return 1.0;}; /** Return the PDG Code of the parent particle */ virtual int parentPdgCode() const override {return 1;}; + virtual int parentStatus() const override {return 1;}; /** Return the parent particle as a HepMC particle type (only called for particles that will enter the HepMC truth event) */ virtual HepMC::GenParticlePtr parentParticle() override {return nullptr;}; @@ -112,10 +113,8 @@ namespace ISFTesting { virtual double childEkin(unsigned short) const override {return 1.0;}; /** Return the PDG Code of the i-th child particle */ virtual int childPdgCode(unsigned short) const override {return 1;}; - /** Return the barcode of the i-th child particle (if defined as part of the TruthIncident) otherwise return 0 */ - Barcode::ParticleBarcode childBarcode(unsigned short) const override final {return 0;}; /** Return the i-th child as a HepMC particle type and assign the given - Barcode to the simulator particle (only called for particles that will + status to the simulator particle (only called for particles that will enter the HepMC truth event) */ virtual HepMC::GenParticlePtr childParticle(unsigned short, Barcode::ParticleBarcode) override {return nullptr;}; @@ -125,8 +124,8 @@ namespace ISFTesting { simulation). */ virtual HepMC::GenParticlePtr updateChildParticle(unsigned short, HepMC::GenParticlePtr ) const override {return nullptr;}; - /** Set the the barcode of all child particles to the given bc */ - virtual void setAllChildrenBarcodes(Barcode::ParticleBarcode) override {}; + /** Set the the status of all child particles to the given bc */ + virtual void setAllChildrenStatus(int) override {}; private: const HepMC::FourVector m_myPosition{0.0, 40.0, 0.0, 40.0}; }; @@ -283,9 +282,6 @@ namespace ISFTesting { EXPECT_CALL(ti, physicsProcessCode()) .Times(1) .WillOnce(::testing::Return(21)); - EXPECT_CALL(ti, parentBarcode()) - .Times(1) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))); EXPECT_CALL(ti, parentParticle()) .Times(1) .WillOnce(::testing::Return(inParticle3)); @@ -294,7 +290,6 @@ namespace ISFTesting { ASSERT_NE( nullptr, generated); if (generated) { ASSERT_EQ( vtxPosition, generated->position() ); - ASSERT_EQ( -200001, HepMC::barcode(generated) ); #ifdef HEPMC3 ASSERT_EQ( 1021, generated->status() ); ASSERT_EQ( 1u, generated->particles_in().size()); @@ -344,10 +339,6 @@ namespace ISFTesting { .Times(2) .WillOnce(::testing::Return(21)) .WillOnce(::testing::Return(21)); - EXPECT_CALL(ti, parentBarcode()) - .Times(2) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))); EXPECT_CALL(ti, parentParticle()) .Times(3) .WillOnce(::testing::Return(inParticle3)) @@ -364,11 +355,10 @@ namespace ISFTesting { #endif recordIncidentToMCTruth(ti,false); - HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-200001); //Find a nicer way to get this. + HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-1); //Find a nicer way to get this. ASSERT_NE( nullptr, generated); if (generated) { ASSERT_EQ( vtxPosition, generated->position() ); - ASSERT_EQ( -200001, HepMC::barcode(generated) ); // by construction at the moment #ifdef HEPMC3 ASSERT_EQ( 1021, generated->status() ); ASSERT_EQ( 1u, generated->particles_in().size()); @@ -416,13 +406,8 @@ namespace ISFTesting { EXPECT_CALL(ti, physicsProcessCode()) .Times(1) .WillOnce(::testing::Return(21)); - EXPECT_CALL(ti, parentBarcode()) - .Times(2) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))); - registerTruthIncident(ti); - HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-200001); //Find a nicer way to get this. + HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-1); //Find a nicer way to get this. HepMC::GenVertexPtr expectedVtx(nullptr); ASSERT_EQ( expectedVtx, generated); } @@ -465,21 +450,15 @@ namespace ISFTesting { .Times(2) .WillOnce(::testing::Return(21)) .WillOnce(::testing::Return(21)); - EXPECT_CALL(ti, parentBarcode()) - .Times(3) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))); EXPECT_CALL(ti, parentParticle()) .Times(1) .WillOnce(::testing::Return(inParticle3)); registerTruthIncident(ti); - HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-200001); //Find a nicer way to get this. + HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-1); //Find a nicer way to get this. ASSERT_NE( nullptr, generated); if (generated) { ASSERT_EQ( vtxPosition, generated->position() ); - ASSERT_EQ( -200001, HepMC::barcode(generated) ); // by construction at the moment #ifdef HEPMC3 ASSERT_EQ( 1021, generated->status() ); ASSERT_EQ( 1u, generated->particles_in().size()); @@ -542,11 +521,6 @@ namespace ISFTesting { .Times(2) .WillOnce(::testing::Return(21)) .WillOnce(::testing::Return(21)); - EXPECT_CALL(ti, parentBarcode()) - .Times(3) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))); EXPECT_CALL(ti, parentParticle()) .Times(1) .WillOnce(::testing::Return(inParticle3)); @@ -558,11 +532,10 @@ namespace ISFTesting { .WillOnce(::testing::Return(false)); registerTruthIncident(ti); - HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-200001); //Find a nicer way to get this. + HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-1); //Find a nicer way to get this. ASSERT_NE( nullptr, generated); if (generated) { ASSERT_EQ( vtxPosition, generated->position() ); - ASSERT_EQ( -200001, HepMC::barcode(generated) ); // by construction at the moment #ifdef HEPMC3 ASSERT_EQ( 1021, generated->status() ); ASSERT_EQ( 1u, generated->particles_in().size()); @@ -623,11 +596,6 @@ namespace ISFTesting { .Times(2) .WillOnce(::testing::Return(21)) .WillOnce(::testing::Return(21)); - EXPECT_CALL(ti, parentBarcode()) - .Times(3) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))) - .WillOnce(::testing::Return(HepMC::barcode(inParticle3))); EXPECT_CALL(ti, parentParticle()) .Times(3) .WillOnce(::testing::Return(inParticle3)) @@ -650,11 +618,10 @@ namespace ISFTesting { .WillOnce(::testing::Return(true)); registerTruthIncident(ti); - HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-200001); //Find a nicer way to get this. + HepMC::GenVertexPtr generated = HepMC::barcode_to_vertex(anEvent.get(),-1); //Find a nicer way to get this. ASSERT_NE( nullptr, generated); if (generated) { ASSERT_EQ( vtxPosition, generated->position() ); - ASSERT_EQ( -200001, HepMC::barcode(generated) ); // by construction at the moment #ifdef HEPMC3 ASSERT_EQ( 1021, generated->status() ); ASSERT_EQ( 1u, generated->particles_in().size()); diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.cxx index c66a6f13fd54..74298ffb75c7 100644 --- a/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.cxx +++ b/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.cxx @@ -141,11 +141,11 @@ StatusCode ISF::PunchThroughTool::initialize() } //barcode service - if (m_barcodeSvc.retrieve().isFailure() ) - { - ATH_MSG_ERROR( "[ punchthrough ] Could not retrieve " << m_barcodeSvc ); - return StatusCode::FAILURE; - } +// if (m_barcodeSvc.retrieve().isFailure() ) +// { +// ATH_MSG_ERROR( "[ punchthrough ] Could not retrieve " << m_barcodeSvc ); +// return StatusCode::FAILURE; +// } //envelope definition service if (m_envDefSvc.retrieve().isFailure() ) @@ -1349,7 +1349,7 @@ ISF::ISFParticle* ISF::PunchThroughTool::createExitPs( const ISF::ISFParticle &i //assign barcodes to the produced particles Barcode::PhysicsProcessCode processCode{0}; - const Barcode::ParticleBarcode secBC = m_barcodeSvc->newSecondary( isfp.barcode(), processCode); + const Barcode::ParticleBarcode secBC = 0;/* m_barcodeSvc->newSecondary( isfp.barcode(), processCode);*/ //FIXME ISF::ISFParticle* finalPar = new ISF::ISFParticle (pos, mom, mass, charge, pdg, 1, pTime, isfp, secBC); finalPar->setNextGeoID( AtlasDetDescr::fAtlasMS); diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.h b/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.h index 05a7e829a215..b9995276295c 100644 --- a/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.h +++ b/Simulation/ISF/ISF_FastCaloSim/ISF_PunchThroughTools/src/PunchThroughTool.h @@ -202,7 +202,7 @@ namespace ISF { *---------------------------------------------------------------------*/ ServiceHandle<IPartPropSvc> m_particlePropSvc{this, "PartPropSvc", "PartPropSvc", "particle properties svc"}; ServiceHandle<IGeoIDSvc> m_geoIDSvc{this, "GeoIDSvc", "ISF::GeoIDSvc"}; - ServiceHandle<Barcode::IBarcodeSvc> m_barcodeSvc{this, "BarcodeSvc", "BarcodeSvc"}; +// ServiceHandle<Barcode::IBarcodeSvc> m_barcodeSvc{this, "BarcodeSvc", "BarcodeSvc"}; ServiceHandle<IEnvelopeDefSvc> m_envDefSvc{this, "EnvelopeDefSvc", "AtlasGeometry_EnvelopeDefSvc"}; /** beam pipe radius */ diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerTool.cxx b/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerTool.cxx index da96e1d95a51..0f4e824a9ced 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerTool.cxx +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerTool.cxx @@ -189,6 +189,7 @@ ISF::EntryLayer ISF::EntryLayerTool::registerParticle(const ISF::ISFParticle& pa : particle.barcode(); m_collection[layerHit]->Emplace(particle.pdgCode(), + particle.status(), energy, hepMom, hepPos, diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerToolMT.cxx b/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerToolMT.cxx index d2f656dc5998..123b1de24b76 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerToolMT.cxx +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4CommonTools/src/EntryLayerToolMT.cxx @@ -137,6 +137,7 @@ ISF::EntryLayer ISF::EntryLayerToolMT::registerParticle(const ISF::ISFParticle& : particle.barcode(); (*m_collectionHolder.get())[layerHit]->Emplace(particle.pdgCode(), + particle.status(), energy, hepMom, hepPos, diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h index f0eeec21529a..5adf67fa7571 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h @@ -60,15 +60,13 @@ namespace iGeant4 { double parentEkin() const override final; /** Return the PDG Code of the parent particle */ int parentPdgCode() const override final; - /** Return the barcode of the parent particle */ - Barcode::ParticleBarcode parentBarcode() override final; - /** Return the bunch-crossing identifier of the parent particle */ - int parentBCID() const override final; + /** Return the status of the parent particle */ + int parentStatus() const override final; /** Return a boolean whether or not the parent particle survives the incident */ bool parentSurvivesIncident() const override final; /** Return the parent particle after the TruthIncident vertex (and give - it a new barcode) */ - HepMC::GenParticlePtr parentParticleAfterIncident(Barcode::ParticleBarcode newBC) override final; + it a new status) */ + HepMC::GenParticlePtr parentParticleAfterIncident(int newBC) override final; /** Return p of the i-th child particle */ const G4ThreeVector childP(unsigned short index) const; @@ -80,10 +78,8 @@ namespace iGeant4 { double childEkin(unsigned short index) const override final; /** Return the PDG Code of the i-th child particle */ int childPdgCode(unsigned short index) const override final; - /** Return the barcode of the i-th child particle (if defined as part of the TruthIncident) otherwise return 0 */ - Barcode::ParticleBarcode childBarcode(unsigned short index) const override final; - /** Set the the barcode of all child particles to the given bc */ - void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) override final; + /** Set the the status of all child particles to the given bc */ + void setAllChildrenStatus(int bc) override final; /** The interaction classifications are described as follows: STD_VTX: interaction of a particle without a pre-defined decay; @@ -99,7 +95,7 @@ namespace iGeant4 { /** Return the parent particle as a HepMC particle type */ HepMC::GenParticlePtr parentParticle() override final; /** Return the i-th child as a HepMC particle type and assign the given - Barcode to the simulator particle */ + status to the simulator particle */ HepMC::GenParticlePtr childParticle(unsigned short index, Barcode::ParticleBarcode bc) override final; /** Update the properties of a child particle from a pre-defined diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx index da44918bfe1f..9428e11155ea 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx @@ -114,10 +114,11 @@ int iGeant4::Geant4TruthIncident::parentPdgCode() const { return m_step->GetTrack()->GetDefinition()->GetPDGEncoding(); } -Barcode::ParticleBarcode iGeant4::Geant4TruthIncident::parentBarcode() { - auto parent = parentParticle(); +int iGeant4::Geant4TruthIncident::parentStatus() const { + // auto parent = parentParticle(); - return (parent) ? HepMC::barcode(parent) : Barcode::fUndefinedBarcode; + //return (parent) ? parent->status() : 0; +return 1;//FIXME } HepMC::GenParticlePtr iGeant4::Geant4TruthIncident::parentParticle() { @@ -219,9 +220,9 @@ Barcode::ParticleBarcode iGeant4::Geant4TruthIncident::childBarcode(unsigned sh return 0; } -void iGeant4::Geant4TruthIncident::setAllChildrenBarcodes(Barcode::ParticleBarcode) { +void iGeant4::Geant4TruthIncident::setAllChildrenStatus(int) { G4ExceptionDescription description; - description << G4String("setAllChildrenBarcodes: ") + "Shared child particle barcodes are not implemented in ISF_Geant4 at this point."; + description << G4String("setAllChildrenStatus: ") + "Shared child particle status are not implemented in ISF_Geant4 at this point."; G4Exception("iGeant4::Geant4TruthIncident", "NotImplemented", FatalException, description); } @@ -333,11 +334,10 @@ HepMC::GenParticlePtr iGeant4::Geant4TruthIncident::convert(const G4Track *track track->GetDynamicParticle() && track->GetDynamicParticle()->GetPrimaryParticle() && track->GetDynamicParticle()->GetPrimaryParticle()->GetUserInformation()){ - // Then the new particle should use the same barcode as the old one!! PrimaryParticleInformation* ppi = dynamic_cast<PrimaryParticleInformation*>( track->GetDynamicParticle()->GetPrimaryParticle()->GetUserInformation() ); - HepMC::suggest_barcode( newParticle, ppi->GetParticleBarcode() ); + //HepMC::suggest_barcode( newParticle, ppi->GetParticleBarcode() ); } else { - HepMC::suggest_barcode( newParticle, barcode ); + // HepMC::suggest_barcode( newParticle, barcode ); } #endif diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/ISFTrajectory.cxx b/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/ISFTrajectory.cxx index df574632384a..8572c8b9cb69 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/ISFTrajectory.cxx +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/ISFTrajectory.cxx @@ -115,8 +115,8 @@ void iGeant4::ISFTrajectory::AppendStep(const G4Step* aStep) // ITruthSvc::registerTruthIncident call above auto currentGenPart = atlasG4EvtUserInfo->GetCurrentlyTraced(); baseIsp->getTruthBinding()->setTruthParticle( currentGenPart ); - Barcode::ParticleBarcode newBarcode = HepMC::barcode(currentGenPart); - baseIsp->setBarcode( newBarcode ); + // Barcode::ParticleBarcode newBarcode = HepMC::barcode(currentGenPart); + // baseIsp->setBarcode( newBarcode ); } } else { diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/TrackProcessorUserActionFullG4.cxx b/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/TrackProcessorUserActionFullG4.cxx index a4a923582752..13a8025e562e 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/TrackProcessorUserActionFullG4.cxx +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4Tools/src/TrackProcessorUserActionFullG4.cxx @@ -148,8 +148,8 @@ namespace G4UA{ tmpISP->setNextGeoID(nextGeoID); tmpISP->setNextSimID(ISF::fUndefinedSimID); - auto generationZeroBarcode = tHelp.GetBarcode(); - tmpISP->setBarcode(generationZeroBarcode); + // auto generationZeroBarcode = tHelp.GetBarcode(); + // tmpISP->setBarcode(generationZeroBarcode); tmpISP->setNextGeoID( nextGeoID ); -- GitLab From d197152a1a21063ce46f8e4e9e247a85a6a22c1a Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Mon, 7 Aug 2023 12:13:43 +0200 Subject: [PATCH 02/22] Try to fix compilation --- .../ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h b/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h index ebfe097f3061..c53f8d485224 100644 --- a/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h +++ b/Simulation/ISF/ISF_Core/ISF_Event/ISF_Event/ISFTruthIncident.h @@ -72,7 +72,6 @@ namespace ISF { (usually only called for particles that will enter the HepMC truth event) */ HepMC::GenParticlePtr parentParticle() override final; int parentStatus() const override final; - HepMC::GenParticlePtr parentParticle() override final; /** Return the barcode of the parent particle */ Barcode::ParticleBarcode parentBarcode() override final; /** Return the bunch-crossing identifier of the parent particle */ @@ -91,6 +90,8 @@ namespace ISF { double childEkin(unsigned short index) const override final; /** Return the PDG Code of the i-th child particle */ int childPdgCode(unsigned short index) const override final; + /** Return the barcode of the i-th child particle (if defined as part of the TruthIncident) otherwise return 0 */ + Barcode::ParticleBarcode childBarcode(unsigned short) const override final; /** Return the i-th child as a HepMC particle type and assign the given status to the simulator particle (usually only called for particles that will enter the HepMC truth event) */ @@ -103,7 +104,7 @@ namespace ISF { virtual HepMC::GenParticlePtr updateChildParticle(unsigned short index, HepMC::GenParticlePtr existingChild) const override final; /** Set the the barcode of all child particles to the given bc */ - void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) override final; + void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) override final; /** Set the the status of all child particles to the given bc */ void setAllChildrenStatus(int bc) override final; private: -- GitLab From a10d331b4e1dcc0824d681d97acd497b1ed70bfb Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@mppui4.t2.rzg.mpg.de> Date: Mon, 7 Aug 2023 14:52:08 +0200 Subject: [PATCH 03/22] Fix compilation --- .../ISF_Core/ISF_Services/test/TruthSvc_test.cxx | 4 ++++ .../ISF_Geant4Event/Geant4TruthIncident.h | 11 +++++++++-- .../ISF_Geant4Event/src/Geant4TruthIncident.cxx | 13 +++++++++++++ 3 files changed, 26 insertions(+), 2 deletions(-) diff --git a/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx b/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx index ca9b733fc323..ea92fddacfa9 100644 --- a/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx +++ b/Simulation/ISF/ISF_Core/ISF_Services/test/TruthSvc_test.cxx @@ -113,6 +113,8 @@ namespace ISFTesting { virtual double childEkin(unsigned short) const override {return 1.0;}; /** Return the PDG Code of the i-th child particle */ virtual int childPdgCode(unsigned short) const override {return 1;}; + /** Return the barcode of the i-th child particle (if defined as part of the TruthIncident) otherwise return 0 */ + Barcode::ParticleBarcode childBarcode(unsigned short) const override final {return 0;}; /** Return the i-th child as a HepMC particle type and assign the given status to the simulator particle (only called for particles that will enter the HepMC truth event) */ @@ -124,6 +126,8 @@ namespace ISFTesting { simulation). */ virtual HepMC::GenParticlePtr updateChildParticle(unsigned short, HepMC::GenParticlePtr ) const override {return nullptr;}; + /** Set the the barcode of all child particles to the given bc */ + virtual void setAllChildrenBarcodes(Barcode::ParticleBarcode) override {}; /** Set the the status of all child particles to the given bc */ virtual void setAllChildrenStatus(int) override {}; private: diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h index 5adf67fa7571..d2c6b92e4d72 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/ISF_Geant4Event/Geant4TruthIncident.h @@ -60,13 +60,17 @@ namespace iGeant4 { double parentEkin() const override final; /** Return the PDG Code of the parent particle */ int parentPdgCode() const override final; + /** Return the barcode of the parent particle */ + Barcode::ParticleBarcode parentBarcode() override final; + /** Return the bunch-crossing identifier of the parent particle */ + int parentBCID() const override final; /** Return the status of the parent particle */ int parentStatus() const override final; /** Return a boolean whether or not the parent particle survives the incident */ bool parentSurvivesIncident() const override final; /** Return the parent particle after the TruthIncident vertex (and give it a new status) */ - HepMC::GenParticlePtr parentParticleAfterIncident(int newBC) override final; + HepMC::GenParticlePtr parentParticleAfterIncident(Barcode::ParticleBarcode newBC) override final; /** Return p of the i-th child particle */ const G4ThreeVector childP(unsigned short index) const; @@ -78,7 +82,10 @@ namespace iGeant4 { double childEkin(unsigned short index) const override final; /** Return the PDG Code of the i-th child particle */ int childPdgCode(unsigned short index) const override final; - /** Set the the status of all child particles to the given bc */ + /** Return the barcode of the i-th child particle (if defined as part of the TruthIncident) otherwise return 0 */ + Barcode::ParticleBarcode childBarcode(unsigned short index) const override final; + /** Set the the barcode of all child particles to the given bc */ + void setAllChildrenBarcodes(Barcode::ParticleBarcode bc) override final; void setAllChildrenStatus(int bc) override final; /** The interaction classifications are described as follows: diff --git a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx index 9428e11155ea..107ee28853f8 100644 --- a/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx +++ b/Simulation/ISF/ISF_Geant4/ISF_Geant4Event/src/Geant4TruthIncident.cxx @@ -114,6 +114,12 @@ int iGeant4::Geant4TruthIncident::parentPdgCode() const { return m_step->GetTrack()->GetDefinition()->GetPDGEncoding(); } +Barcode::ParticleBarcode iGeant4::Geant4TruthIncident::parentBarcode() { + auto parent = parentParticle(); + + return (parent) ? HepMC::barcode(parent) : Barcode::fUndefinedBarcode; +} + int iGeant4::Geant4TruthIncident::parentStatus() const { // auto parent = parentParticle(); @@ -220,6 +226,12 @@ Barcode::ParticleBarcode iGeant4::Geant4TruthIncident::childBarcode(unsigned sh return 0; } +void iGeant4::Geant4TruthIncident::setAllChildrenBarcodes(Barcode::ParticleBarcode) { + G4ExceptionDescription description; + description << G4String("setAllChildrenBarcodes: ") + "Shared child particle barcodes are not implemented in ISF_Geant4 at this point."; + G4Exception("iGeant4::Geant4TruthIncident", "NotImplemented", FatalException, description); +} + void iGeant4::Geant4TruthIncident::setAllChildrenStatus(int) { G4ExceptionDescription description; description << G4String("setAllChildrenStatus: ") + "Shared child particle status are not implemented in ISF_Geant4 at this point."; @@ -334,6 +346,7 @@ HepMC::GenParticlePtr iGeant4::Geant4TruthIncident::convert(const G4Track *track track->GetDynamicParticle() && track->GetDynamicParticle()->GetPrimaryParticle() && track->GetDynamicParticle()->GetPrimaryParticle()->GetUserInformation()){ + // Then the new particle should use the same barcode as the old one!! PrimaryParticleInformation* ppi = dynamic_cast<PrimaryParticleInformation*>( track->GetDynamicParticle()->GetPrimaryParticle()->GetUserInformation() ); //HepMC::suggest_barcode( newParticle, ppi->GetParticleBarcode() ); } else { -- GitLab From 787b83a34bebca5a2087b5fac6ea7ce944dae419 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@mpp.mpg.de> Date: Tue, 8 Aug 2023 13:45:58 +0200 Subject: [PATCH 04/22] Try to check the null pointers --- .../Root/MCTruthClassifierGen.cxx | 20 +++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx index 411d5a71610b..784000a68955 100644 --- a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx @@ -219,7 +219,7 @@ MCTruthClassifier::particleTruthClassifier(const xAOD::TruthParticle* thePart, I return std::make_pair(Neutrino, partOrig); } - if (HepMC::is_same_generator_particle(thePart,info->Mother())) + if (HepMC::is_same_generator_particle(thePart, info ? info->Mother() : nullptr )) return std::make_pair(NonPrimary, partOrig); if (isPartHadr) @@ -507,7 +507,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth const xAOD::TruthParticle* theDaug = partOriVert->outgoingParticle(ipOut); if (!theDaug) continue; - if (motherPDG == theDaug->pdgId() && HepMC::is_same_generator_particle(theDaug,info->Mother())) + if (motherPDG == theDaug->pdgId() && HepMC::is_same_generator_particle(theDaug,info ? info->Mother() : nullptr)) samePart = true; } @@ -625,7 +625,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth if (abs(DaugType) == 42) NumOfLQ++; - if (abs(DaugType) == abs(motherPDG) && HepMC::is_same_generator_particle(theDaug, info->Mother())) + if (abs(DaugType) == abs(motherPDG) && HepMC::is_same_generator_particle(theDaug, info ? info->Mother() : nullptr)) samePart = true; if (numOfParents == 1 && (motherPDG == 22 || abs(motherPDG) == 11 || abs(motherPDG) == 13 || abs(motherPDG) == 211) && @@ -1752,7 +1752,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE bool foundISR = false; bool foundFSR = false; - if (numOfParents == 1 && numOfDaug == 2 && HepMC::is_same_generator_particle(Daug,info->Mother())) + if (numOfParents == 1 && numOfDaug == 2 && HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr)) return BremPhot; if (numOfParents == 1 && numOfDaug == 2 && abs(motherPDG) == 11 && NumOfPht == 2) @@ -1760,7 +1760,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE // decay of W,Z and Higgs to lepton with FSR generated by Pythia if (numOfParents == 1 && numOfDaug == 2 && (abs(motherPDG) == 11 || abs(motherPDG) == 13 || abs(motherPDG) == 15) && - HepMC::is_same_generator_particle(Daug, info->Mother()) && mothOriVert != nullptr && + HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr ) && mothOriVert != nullptr && mothOriVert->nIncomingParticles() == 1) { int itr = 0; long PartPDG = 0; @@ -1930,7 +1930,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE //-- Exotics - CompHep if (abs(motherPDG) == 11 && numOfParents == 1 && numOfDaug == 2 && (NumOfEl == 1 || NumOfPos == 1) && NumOfPht == 1 && - !HepMC::is_same_generator_particle(Daug, info->Mother()) && DaugBarcode < 20000 && motherBarcode < 20000) + !HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr ) && DaugBarcode < 20000 && motherBarcode < 20000) return FSRPhot; // FSR from Photos @@ -2254,7 +2254,7 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth if (std::abs(DaugType) == 42) NumOfLQ++; - if (std::abs(DaugType) == std::abs(motherPDG) && HepMC::is_same_generator_particle(theDaug,info->Mother())) + if (std::abs(DaugType) == std::abs(motherPDG) && HepMC::is_same_generator_particle(theDaug,info ? info->Mother() : nullptr)) samePart = true; } // cycle itrDaug @@ -3015,7 +3015,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i && info->photonMotherStatus < 3) { do { const xAOD::TruthParticle* theMotherPart = - data_to_particle(truthParticleContainerReadHandle.ptr(), info->PhotonMother()); + data_to_particle(truthParticleContainerReadHandle.ptr(), info ? info->PhotonMother() : nullptr); if (theMotherPart == nullptr || theMotherPart == thePart) break; thePart = theMotherPart; @@ -3030,7 +3030,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i if (part.first == BkgElectron && part.second == PhotonConv) { // in case of photon from gen particle classify photon // part=particleTruthClassifier(mother); - thePart = data_to_particle(truthParticleContainerReadHandle.ptr(), info->Mother()); + thePart = data_to_particle(truthParticleContainerReadHandle.ptr(), info ? info->Mother() : nullptr); if (thePart != nullptr) part = particleTruthClassifier(thePart, info); @@ -3042,7 +3042,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i } else { // in case of photon from gen particle classify photon - thePart = data_to_particle(truthParticleContainerReadHandle.ptr(), info->Mother()); + thePart = data_to_particle(truthParticleContainerReadHandle.ptr(), info ? info->Mother() : nullptr); if (thePart != nullptr) part = particleTruthClassifier(thePart, info); } -- GitLab From 152512987aa66eaa644858168f2724181d776f0a Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Tue, 8 Aug 2023 14:03:44 +0200 Subject: [PATCH 05/22] Update MuonSegmentTruthAssociationAlg.cxx --- .../MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx | 1 - 1 file changed, 1 deletion(-) diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx index 3a1b9525aead..6e0bd5905e0f 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonSegmentTruthAssociationAlg.cxx @@ -10,7 +10,6 @@ #include "xAODMuon/MuonSegmentAuxContainer.h" #include "xAODTruth/TruthParticleAuxContainer.h" #include "xAODTruth/TruthParticleContainer.h" -#include "xAODTruth/TruthVertex.h" #include "TruthUtils/MagicNumbers.h" #include "xAODTruth/TruthVertex.h" -- GitLab From 0f47ea98f2d9c343cd9b02a971989771f60cbdfd Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Tue, 8 Aug 2023 17:14:25 +0200 Subject: [PATCH 06/22] Update MCTruthClassifierAthena.cxx --- .../MCTruthClassifier/src/MCTruthClassifierAthena.cxx | 1 + 1 file changed, 1 insertion(+) diff --git a/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx b/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx index c5f288969139..197929d97e83 100644 --- a/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx @@ -180,6 +180,7 @@ MCTruthClassifier::egammaClusMatch(const xAOD::CaloCluster* clus, } theMatchPart = find_matching(truthParticleContainerReadHandle.ptr(), thePart); + if (info) { info->egPartPtr.push_back(thePart); info->egPartdR.push_back(dR); -- GitLab From a705019750010e2a65ac6b0fea03dbc96d70e4f5 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Tue, 8 Aug 2023 17:16:17 +0200 Subject: [PATCH 07/22] OK --- .../MCTruthClassifier/src/MCTruthClassifierAthena.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx b/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx index 197929d97e83..0b7c1c6d1cdd 100644 --- a/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/src/MCTruthClassifierAthena.cxx @@ -180,7 +180,7 @@ MCTruthClassifier::egammaClusMatch(const xAOD::CaloCluster* clus, } theMatchPart = find_matching(truthParticleContainerReadHandle.ptr(), thePart); - + if (info) { info->egPartPtr.push_back(thePart); info->egPartdR.push_back(dR); -- GitLab From 68e9e217aa435a334d04cd55fc8db194a0794055 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Wed, 9 Aug 2023 12:03:37 +0200 Subject: [PATCH 08/22] OK --- PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx index 937682d5c4cd..191fb39a4eba 100644 --- a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx @@ -3079,7 +3079,7 @@ MCTruthClassifier::partCharge(const xAOD::TruthParticle* thePart) } const xAOD::TruthParticle* -MCTruthClassifier::data_to_particle(const xAOD::TruthParticleContainer* TruthTES, const xAOD::TruthParticle* xx) { +MCTruthClassifier::find_matching(const xAOD::TruthParticleContainer* TruthTES, const xAOD::TruthParticle* xx) { const xAOD::TruthParticle* ptrPart = nullptr; if (!bcin) return ptrPart; for (const auto *const truthParticle : *TruthTES) { -- GitLab From a8409c0da0fe49068fedfd09d51dc7d840d92f6f Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@mpp.mpg.de> Date: Thu, 10 Aug 2023 09:27:30 +0200 Subject: [PATCH 09/22] OK --- PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx index 191fb39a4eba..4a8f7ad3d495 100644 --- a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx @@ -3081,7 +3081,7 @@ MCTruthClassifier::partCharge(const xAOD::TruthParticle* thePart) const xAOD::TruthParticle* MCTruthClassifier::find_matching(const xAOD::TruthParticleContainer* TruthTES, const xAOD::TruthParticle* xx) { const xAOD::TruthParticle* ptrPart = nullptr; - if (!bcin) return ptrPart; + if (!xx) return ptrPart; for (const auto *const truthParticle : *TruthTES) { if (HepMC::is_sim_descendant(xx,truthParticle)) { ptrPart = truthParticle; -- GitLab From c1e65aa2b314c36cafe46306b909754e252cff34 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Thu, 10 Aug 2023 14:38:34 +0200 Subject: [PATCH 10/22] OK --- .../Root/MCTruthClassifierGen.cxx | 63 ++++++++++++------- 1 file changed, 39 insertions(+), 24 deletions(-) diff --git a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx index 4a8f7ad3d495..8f2da1011811 100644 --- a/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx +++ b/PhysicsAnalysis/MCTruthClassifier/Root/MCTruthClassifierGen.cxx @@ -87,10 +87,11 @@ MCTruthClassifier::compareTruthParticles(HepMC::ConstGenParticlePtr genPart, con //--------------------------------------------------------------------------------------- std::pair<ParticleType, ParticleOrigin> -MCTruthClassifier::particleTruthClassifier(const xAOD::TruthParticle* thePart, Info* info /*= nullptr*/) const +MCTruthClassifier::particleTruthClassifier(const xAOD::TruthParticle* thePart, Info* infoin /*= nullptr*/) const { //--------------------------------------------------------------------------------------- + Info* info = infoin; ATH_MSG_DEBUG("Executing particleTruthClassifier"); ParticleType partType = Unknown; @@ -105,9 +106,10 @@ MCTruthClassifier::particleTruthClassifier(const xAOD::TruthParticle* thePart, I return std::make_pair(partType, partOrig); } const EventContext& ctx = info ? info->eventContext : Gaudi::Hive::currentContext(); + Info tmpinfo; + if (!info) { info = &tmpinfo; } - SG::ReadHandle<xAOD::TruthParticleContainer> truthParticleContainerReadHandle( - m_truthParticleContainerKey,ctx); + SG::ReadHandle<xAOD::TruthParticleContainer> truthParticleContainerReadHandle(m_truthParticleContainerKey,ctx); if (!truthParticleContainerReadHandle.isValid()) { ATH_MSG_WARNING( " Invalid ReadHandle for xAOD::TruthParticleContainer with key: " << truthParticleContainerReadHandle.key()); @@ -451,10 +453,11 @@ ParticleOrigin MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruthTES, const xAOD::TruthParticle* thePart, bool& isPrompt, - Info* info) const + Info* infoin) const //------------------------------------------------------------------------------- { + Info* info = infoin; ATH_MSG_DEBUG("Executing DefOrigOfElectron "); const xAOD::TruthParticle* thePriPart = find_matching(mcTruthTES, thePart); @@ -465,6 +468,8 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth //-- to define electron outcome status if (info) info->particleOutCome = defOutComeOfElectron(thePriPart); + Info tmpinfo; + if (!info) { info = &tmpinfo; } if (!partOriVert) return NonDefined; @@ -501,7 +506,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth const xAOD::TruthParticle* theDaug = partOriVert->outgoingParticle(ipOut); if (!theDaug) continue; - if (motherPDG == theDaug->pdgId() && HepMC::is_same_generator_particle(theDaug,info ? info->Mother() : nullptr)) + if (motherPDG == theDaug->pdgId() && HepMC::is_same_generator_particle(theDaug, info ? info->Mother() : nullptr)) samePart = true; } @@ -623,7 +628,7 @@ MCTruthClassifier::defOrigOfElectron(const xAOD::TruthParticleContainer* mcTruth if (abs(DaugType) == 42) NumOfLQ++; - if (abs(DaugType) == abs(motherPDG) && HepMC::is_same_generator_particle(theDaug, info ? info->Mother() : nullptr)) + if (abs(DaugType) == abs(motherPDG) && HepMC::is_same_generator_particle(theDaug, info ? info->Mother() : nullptr )) samePart = true; if (numOfParents == 1 && (motherPDG == 22 || abs(motherPDG) == 11 || abs(motherPDG) == 13 || abs(motherPDG) == 211) && @@ -934,9 +939,10 @@ ParticleOrigin MCTruthClassifier::defOrigOfMuon(const xAOD::TruthParticleContainer* mcTruthTES, const xAOD::TruthParticle* thePart, bool& isPrompt, - Info* info) const + Info* infoin) const { + Info* info = infoin; ATH_MSG_DEBUG("Executing DefOrigOfMuon "); const xAOD::TruthParticle* thePriPart = find_matching(mcTruthTES, thePart); @@ -950,6 +956,8 @@ MCTruthClassifier::defOrigOfMuon(const xAOD::TruthParticleContainer* mcTruthTES, //-- to define muon outcome status if (info) info->particleOutCome = defOutComeOfMuon(thePriPart); + Info tmpinfo; + if (!info) { info = &tmpinfo; } if (!partOriVert) return NonDefined; int numOfParents = partOriVert->nIncomingParticles(); @@ -1344,8 +1352,9 @@ ParticleOrigin MCTruthClassifier::defOrigOfTau(const xAOD::TruthParticleContainer* mcTruthTES, const xAOD::TruthParticle* thePart, int motherPDG, - Info* info) const + Info* infoin) const { + Info* info = infoin; //------------------------------------------------------------------------------- ATH_MSG_DEBUG("Executing DefOrigOfTau "); @@ -1364,6 +1373,8 @@ MCTruthClassifier::defOrigOfTau(const xAOD::TruthParticleContainer* mcTruthTES, info->particleOutCome = defOutComeOfTau(thePriPart, info); } + Info tmpinfo; + if (!info) { info = &tmpinfo; } if (!partOriVert) return NonDefined; @@ -1662,12 +1673,15 @@ ParticleOrigin MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTES, const xAOD::TruthParticle* thePart, bool& isPrompt, - Info* info) const + Info* infoin) const //------------------------------------------------------------------------------- { + Info* info = infoin; ATH_MSG_DEBUG("Executing DefOrigOfPhoton "); + Info tmpinfo; + if (!info) { info = &tmpinfo; } if (info) { info->mother = nullptr; info->photonMother = nullptr; @@ -1675,9 +1689,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE info->photonMotherBarcode = 0; info->photonMotherPDG = 0; info->photonMotherStatus = 0; - info->mother = nullptr; info->motherBarcode = 0; - info->motherStatus = 0; } const xAOD::TruthParticle* thePriPart = find_matching(mcTruthTES, thePart); @@ -1712,8 +1724,8 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE if (!mother) return NonDefined; int motherPDG = mother->pdgId(); - int motherStatus = mother->status(); const xAOD::TruthVertex* mothOriVert = mother->hasProdVtx() ? mother->prodVtx() : nullptr; + int motherStatus = mother->status(); long motherBarcode = mother->barcode(); if (info) { info->mother = mother; @@ -1785,7 +1797,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE bool foundISR = false; bool foundFSR = false; - if (numOfParents == 1 && numOfDaug == 2 && HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr)) + if (numOfParents == 1 && numOfDaug == 2 && HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr )) return BremPhot; if (numOfParents == 1 && numOfDaug == 2 && abs(motherPDG) == 11 && NumOfPht == 2) @@ -1793,7 +1805,7 @@ MCTruthClassifier::defOrigOfPhoton(const xAOD::TruthParticleContainer* mcTruthTE // decay of W,Z and Higgs to lepton with FSR generated by Pythia if (numOfParents == 1 && numOfDaug == 2 && (abs(motherPDG) == 11 || abs(motherPDG) == 13 || abs(motherPDG) == 15) && - HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr ) && mothOriVert != nullptr && + !HepMC::is_same_generator_particle(Daug, info ? info->Mother() : nullptr ) && mothOriVert != nullptr && mothOriVert->nIncomingParticles() == 1) { int itr = 0; long PartPDG = 0; @@ -2112,9 +2124,10 @@ ParticleOrigin MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruthTES, const xAOD::TruthParticle* thePart, bool& isPrompt, - Info* info) const + Info* infoin) const //------------------------------------------------------------------------------- { + Info* info = infoin; // author - Pierre-Antoine Delsart // ATH_MSG_DEBUG("Executing DefOrigOfNeutrino "); @@ -2132,6 +2145,8 @@ MCTruthClassifier::defOrigOfNeutrino(const xAOD::TruthParticleContainer* mcTruth if (info) info->particleOutCome = NonInteract; + Info tmpinfo; + if (!info) { info = &tmpinfo; } if (!partOriVert) return NonDefined; @@ -2997,9 +3012,10 @@ MCTruthClassifier::defOutComeOfPhoton(const xAOD::TruthParticle* thePart) const //--------------------------------------------------------------------------------- std::pair<ParticleType, ParticleOrigin> -MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* info /*= nullptr*/) const +MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* infoin /*= nullptr*/) const { + Info* info = infoin; ATH_MSG_DEBUG("executing CheckOrigOfBkgElec " << theEle); std::pair<ParticleType, ParticleOrigin> part; @@ -3021,8 +3037,7 @@ MCTruthClassifier::checkOrigOfBkgElec(const xAOD::TruthParticle* theEle, Info* i << " has valid ReadHandle "); Info tmpinfo; - if (!info) - info = &tmpinfo; + if (!info) { info = &tmpinfo; } part = particleTruthClassifier(theEle, info); if (part.first != BkgElectron || part.second != PhotonConv) @@ -3078,18 +3093,19 @@ MCTruthClassifier::partCharge(const xAOD::TruthParticle* thePart) return thePart?thePart->charge():0.0; } -const xAOD::TruthParticle* -MCTruthClassifier::find_matching(const xAOD::TruthParticleContainer* TruthTES, const xAOD::TruthParticle* xx) { +const xAOD::TruthParticle* +MCTruthClassifier::find_matching(const xAOD::TruthParticleContainer* TruthTES, const xAOD::TruthParticle* bcin) +{ const xAOD::TruthParticle* ptrPart = nullptr; - if (!xx) return ptrPart; + if (!bcin) return ptrPart; for (const auto *const truthParticle : *TruthTES) { - if (HepMC::is_sim_descendant(xx,truthParticle)) { + if (HepMC::is_sim_descendant(bcin,truthParticle)) { ptrPart = truthParticle; break; } } return ptrPart; - } +} //------------------------------------------------------------------------ const xAOD::TruthParticle* @@ -3130,4 +3146,3 @@ MCTruthClassifier::findParticleDaughters(const xAOD::TruthParticle* thePart, } } } - -- GitLab From 7adaf1f3d820c6c9a05cb2fb0f42f83d55b9bd18 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Fri, 11 Aug 2023 11:31:01 +0200 Subject: [PATCH 11/22] Update the truth matching algorithm in MuonTruthDecorationAlg --- .../TruthUtils/TruthUtils/MagicNumbers.h | 37 +++++++++++++++++++ .../src/MuonTruthDecorationAlg.cxx | 3 +- 2 files changed, 39 insertions(+), 1 deletion(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index c5ecce58c0d9..19e4b562ad35 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -85,6 +85,43 @@ inline bool is_sim_descendant(const T1& p1,const T2& p2) { return p1 && p2 && (p template <class T1,class T2, std::enable_if_t< std::is_arithmetic<T1>::value && !std::is_arithmetic<T2>::value, bool> = true > inline bool is_sim_descendant(const T1& p1,const T2& p2) { return p2 && (p1 % SIM_REGENERATION_INCREMENT == p2->barcode()); } +template <class T> int unique_id(const T& p1){ return p1->barcode();} +#if defined(HEPMC3) && !defined(XAOD_STANDALONE) +template <> inline int unique_id(const ConstGenParticlePtr& p1){ return p1->id();} +template <> inline int unique_id(const GenParticlePtr& p1){ return p1->id();} +#endif + + +void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { + if (direction == 0) out.push_back(unique_id(p)); + if (direction <= 0) { + if (p->status()>SIM_STATUS_INCREMENT) { + auto pv = p->production_vertex(); + if (pv) { + for (auto pa: pv->particles_in()) { + if (pa->pdg_id() != p->pdg_id()) continue; + out.push_front(unique_id(p)); + get_particle_history(pa,out,-1); + break; + } + } + } + } + if (direction => 0) { + if (p->status()>SIM_STATUS_INCREMENT) { + auto pv = p->decay_vertex(); + if (pv) { + for (auto pa: pv->particles_out()) { + if (pa->pdg_id() != p->pdg_id()) continue; + out.push_back(unique_id(p)); + get_particle_history(pa,out,1); + break; + } + } + } + } +} +inline std::deque<int> simulation_history(const T& p) { std::deque<int> res; get_particle_history(p,res,0); return res;} } #endif diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx index 94b01614031c..dd87af5f7be2 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx @@ -453,6 +453,7 @@ namespace Muon { std::vector<unsigned int> ntrigEtaHitsPerChamberLayer; ntrigEtaHitsPerChamberLayer.resize(Muon::MuonStationIndex::PhiIndexMax); ATH_MSG_DEBUG("addHitCounts: barcode " << barcode); + auto truthParticleHistory = HepMC::simulation_history(truthParticle); // loop over detector technologies for (SG::ReadHandle<PRD_MultiTruthCollection>& col : m_PRD_TruthNames.makeHandles(ctx)) { if (!col.isPresent()) { @@ -463,7 +464,7 @@ namespace Muon { // loop over trajectories for (const auto& trajectory : *col) { // check if gen particle same as input - if (!HepMC::is_sim_descendant(&trajectory.second,&truthParticle)) continue; + if (truthParticleHistory.find(trajectory.second.barcode()) == truthParticleHistory.end()) continue; const Identifier& id = trajectory.first; bool measPhi = m_idHelperSvc->measuresPhi(id); -- GitLab From 4403eec0e187c0af310c857d43e5c0301a09769e Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Fri, 11 Aug 2023 11:46:44 +0200 Subject: [PATCH 12/22] Fixed the direction --- Generators/TruthUtils/TruthUtils/MagicNumbers.h | 7 +++---- .../MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx | 2 +- 2 files changed, 4 insertions(+), 5 deletions(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index 19e4b562ad35..eddb7ff2eacb 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -93,8 +93,7 @@ template <> inline int unique_id(const GenParticlePtr& p1){ return p1->id();} void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { - if (direction == 0) out.push_back(unique_id(p)); - if (direction <= 0) { + if (direction < 0) { if (p->status()>SIM_STATUS_INCREMENT) { auto pv = p->production_vertex(); if (pv) { @@ -107,7 +106,7 @@ void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { } } } - if (direction => 0) { + if (direction > 0) { if (p->status()>SIM_STATUS_INCREMENT) { auto pv = p->decay_vertex(); if (pv) { @@ -121,7 +120,7 @@ void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { } } } -inline std::deque<int> simulation_history(const T& p) { std::deque<int> res; get_particle_history(p,res,0); return res;} +inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; out.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} } #endif diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx index dd87af5f7be2..1b07ed3f1cca 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx @@ -453,7 +453,7 @@ namespace Muon { std::vector<unsigned int> ntrigEtaHitsPerChamberLayer; ntrigEtaHitsPerChamberLayer.resize(Muon::MuonStationIndex::PhiIndexMax); ATH_MSG_DEBUG("addHitCounts: barcode " << barcode); - auto truthParticleHistory = HepMC::simulation_history(truthParticle); + auto truthParticleHistory = HepMC::simulation_history(truthParticle, -1); // loop over detector technologies for (SG::ReadHandle<PRD_MultiTruthCollection>& col : m_PRD_TruthNames.makeHandles(ctx)) { if (!col.isPresent()) { -- GitLab From 71386fdff2886de0ffb42145448b78eb1d1424ee Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Fri, 11 Aug 2023 13:47:27 +0200 Subject: [PATCH 13/22] Update MagicNumbers.h --- Generators/TruthUtils/TruthUtils/MagicNumbers.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index eddb7ff2eacb..ef97f2d9e2dc 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -92,7 +92,7 @@ template <> inline int unique_id(const GenParticlePtr& p1){ return p1->id();} #endif -void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { +template <class T> inline void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { if (direction < 0) { if (p->status()>SIM_STATUS_INCREMENT) { auto pv = p->production_vertex(); @@ -120,7 +120,7 @@ void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { } } } -inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; out.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} +template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; out.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} } #endif -- GitLab From d63389fd29d181843e4b6b326e9982eef72b1323 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Fri, 11 Aug 2023 15:18:26 +0200 Subject: [PATCH 14/22] Update MagicNumbers.h --- Generators/TruthUtils/TruthUtils/MagicNumbers.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index ef97f2d9e2dc..5f1efdde8855 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -120,7 +120,7 @@ template <class T> inline void get_particle_history(const T& p, std::deque<int>& } } } -template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; out.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} +template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; res.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} } #endif -- GitLab From 1c5db4486a3247f8183bcf9e34ab405b605ae4ee Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Fri, 11 Aug 2023 17:16:38 +0200 Subject: [PATCH 15/22] Update MagicNumbers.h --- Generators/TruthUtils/TruthUtils/MagicNumbers.h | 1 + 1 file changed, 1 insertion(+) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index 5f1efdde8855..d07d89a94152 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -9,6 +9,7 @@ #include <limits> #include <cstdint> #include <memory> +#include <deque> #include <type_traits> #if defined(HEPMC3) && !defined(XAOD_STANDALONE) #include "AtlasHepMC/GenEvent.h" -- GitLab From dbd24424991b5c78e63d05fbdb767348c2cca4e7 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Mon, 14 Aug 2023 11:05:10 +0200 Subject: [PATCH 16/22] Update MuonTruthDecorationAlg.cxx --- MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx index 1b07ed3f1cca..78aa1636bbae 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx @@ -464,7 +464,7 @@ namespace Muon { // loop over trajectories for (const auto& trajectory : *col) { // check if gen particle same as input - if (truthParticleHistory.find(trajectory.second.barcode()) == truthParticleHistory.end()) continue; + if (std::find(truthParticleHistory.begin(),truthParticleHistory.end(),trajectory.second.barcode()) == truthParticleHistory.end()) continue; const Identifier& id = trajectory.first; bool measPhi = m_idHelperSvc->measuresPhi(id); -- GitLab From 0dcfd67cc315695e55331dbd93504026fc748761 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Mon, 14 Aug 2023 11:28:56 +0200 Subject: [PATCH 17/22] OK --- .../TruthUtils/TruthUtils/MagicNumbers.h | 38 +++++++++++++++++++ .../src/MuonTruthDecorationAlg.cxx | 3 +- 2 files changed, 40 insertions(+), 1 deletion(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index 8d291d7ae8cf..9958c4f82a5e 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -128,5 +128,43 @@ template <class T1, class T2> inline bool is_sim_descendant(const T1& p1,const T return (p1->status()/SIM_STATUS_INCREMENT >= p2->status()/SIM_STATUS_INCREMENT) && is_same_generator_particle(p1,p2); } + +template <class T> int unique_id(const T& p1){ return p1->barcode();} +#if defined(HEPMC3) && !defined(XAOD_STANDALONE) +template <> inline int unique_id(const ConstGenParticlePtr& p1){ return p1->id();} +template <> inline int unique_id(const GenParticlePtr& p1){ return p1->id();} +#endif + + +template <class T> inline void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { + if (direction < 0) { + if (p->status()>SIM_STATUS_INCREMENT) { + auto pv = p->production_vertex(); + if (pv) { + for (auto pa: pv->particles_in()) { + if (pa->pdg_id() != p->pdg_id()) continue; + out.push_front(unique_id(p)); + get_particle_history(pa,out,-1); + break; + } + } + } + } + if (direction > 0) { + if (p->status()>SIM_STATUS_INCREMENT) { + auto pv = p->decay_vertex(); + if (pv) { + for (auto pa: pv->particles_out()) { + if (pa->pdg_id() != p->pdg_id()) continue; + out.push_back(unique_id(p)); + get_particle_history(pa,out,1); + break; + } + } + } + } +} +template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; res.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} + } #endif diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx index c2ef3fd6a2d3..78aa1636bbae 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx @@ -453,6 +453,7 @@ namespace Muon { std::vector<unsigned int> ntrigEtaHitsPerChamberLayer; ntrigEtaHitsPerChamberLayer.resize(Muon::MuonStationIndex::PhiIndexMax); ATH_MSG_DEBUG("addHitCounts: barcode " << barcode); + auto truthParticleHistory = HepMC::simulation_history(truthParticle, -1); // loop over detector technologies for (SG::ReadHandle<PRD_MultiTruthCollection>& col : m_PRD_TruthNames.makeHandles(ctx)) { if (!col.isPresent()) { @@ -463,7 +464,7 @@ namespace Muon { // loop over trajectories for (const auto& trajectory : *col) { // check if gen particle same as input - if (!HepMC::is_sim_descendant(trajectory.second.cptr(),&truthParticle)) continue; + if (std::find(truthParticleHistory.begin(),truthParticleHistory.end(),trajectory.second.barcode()) == truthParticleHistory.end()) continue; const Identifier& id = trajectory.first; bool measPhi = m_idHelperSvc->measuresPhi(id); -- GitLab From 3ec77803e2c873937d909e9236f76a0456e6ffd9 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andriish@pcatlas18.mpp.mpg.de> Date: Mon, 14 Aug 2023 13:24:03 +0200 Subject: [PATCH 18/22] OK --- .../TruthUtils/TruthUtils/MagicNumbers.h | 94 +++++-------------- 1 file changed, 22 insertions(+), 72 deletions(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index 9958c4f82a5e..0574c34533f1 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -56,78 +56,6 @@ template <class T> inline bool is_simulation_vertex(const T& p){ return (p->sta -template <class T> int get_id(const T& p1){ return p1->barcode();} -#if defined(HEPMC3) && !defined(XAOD_STANDALONE) -template <> inline int get_id(const ConstGenParticlePtr& p1){ return p1->id();} -template <> inline int get_id(const GenParticlePtr& p1){ return p1->id();} -#endif - - -template <class T> inline T get_generator_particle(const T& p) { - auto pp = p; - for (;;) { - if (pp->status()<SIM_STATUS_INCREMENT) break; - auto pv = pp->production_vertex(); - if (!pv) break; - for (auto pa: pv->particles_in()) { - if (pa->pdg_id()!=pp->pdg_id()) continue; - pp = pa; - break; - } - } - return pp; -} - -template <class T> -inline std::deque<int> get_particle_history(const T& p) { - std::deque<int> res; - res.push_back(get_id(p)); - T pp1 = p; - for (;;) { - if (pp1->status()<SIM_STATUS_INCREMENT) break; - auto pv = pp1->production_vertex(); - if (!pv) break; - for (auto pa: pv->particles_in()) { - if (pa->pdg_id()!=pp1->pdg_id()) continue; - res.push_front(get_id(pp1)); - //pp1 = pa;//FIXME - break; - } - } - - T pp2 = p; - for (;;) { - auto pv = pp2->end_vertex(); - if (!pv) break; - for (auto pa: pv->particles_out()) { - if (pa->pdg_id()!=pp2->pdg_id()) continue; - res.push_back(get_id(pp2)); - //pp2 = pa; //FIXME - break; - } - } - return res; -} - -template <class T1, class T2, std::enable_if_t<std::is_same<T1,T2>::value, bool> = true > -inline bool is_same_generator_particle(const T1& p1,const T2& p2) { - auto p1orig = get_generator_particle(p1); - auto p2orig = get_generator_particle(p2); - return get_id(p1orig)==get_id(p2orig); -} -template <class T1, class T2, std::enable_if_t< !std::is_same<T1,T2>::value, bool> = true > -inline bool is_same_generator_particle(const T1& p1,const T2& p2) { - int b = get_id(p1); - std::deque<int> hist = get_particle_history(p2); - return (std::find(hist.begin(),hist.end(),b) != hist.end()); -} - - - -template <class T1, class T2> inline bool is_sim_descendant(const T1& p1,const T2& p2) { - return (p1->status()/SIM_STATUS_INCREMENT >= p2->status()/SIM_STATUS_INCREMENT) && is_same_generator_particle(p1,p2); -} - template <class T> int unique_id(const T& p1){ return p1->barcode();} #if defined(HEPMC3) && !defined(XAOD_STANDALONE) @@ -166,5 +94,27 @@ template <class T> inline void get_particle_history(const T& p, std::deque<int>& } template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; res.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} + + +template <class T1, class T2, std::enable_if_t<std::is_same<T1,T2>::value, bool> = true > +inline bool is_same_generator_particle(const T1& p1,const T2& p2) { + auto p1orig = simulation_history(p1,-1); + auto p2orig = simulation_history(p2,-1); + return unique_id(p1orig.front())==unique_id(p2orig.front()); +} +template <class T1, class T2, std::enable_if_t< !std::is_same<T1,T2>::value, bool> = true > +inline bool is_same_generator_particle(const T1& p1,const T2& p2) { + int b = unique_id(p1); + auto hist = simulation_history(p2); + return (std::find(hist.begin(),hist.end(),b) != hist.end()); +} + + +template <class T1, class T2> inline bool is_sim_descendant(const T1& p1,const T2& p2) { + return (p1->status()/SIM_STATUS_INCREMENT >= p2->status()/SIM_STATUS_INCREMENT) && is_same_generator_particle(p1,p2); +} + + + } #endif -- GitLab From f5832f503ef0576798a2b28519ec9852ff89206f Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Mon, 14 Aug 2023 15:15:46 +0200 Subject: [PATCH 19/22] Update MuonTruthDecorationAlg.cxx --- MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx index 78aa1636bbae..4a4e6b59e8b1 100644 --- a/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx +++ b/MuonSpectrometer/MuonTruthAlgs/src/MuonTruthDecorationAlg.cxx @@ -453,7 +453,7 @@ namespace Muon { std::vector<unsigned int> ntrigEtaHitsPerChamberLayer; ntrigEtaHitsPerChamberLayer.resize(Muon::MuonStationIndex::PhiIndexMax); ATH_MSG_DEBUG("addHitCounts: barcode " << barcode); - auto truthParticleHistory = HepMC::simulation_history(truthParticle, -1); + auto truthParticleHistory = HepMC::simulation_history(&truthParticle, -1); // loop over detector technologies for (SG::ReadHandle<PRD_MultiTruthCollection>& col : m_PRD_TruthNames.makeHandles(ctx)) { if (!col.isPresent()) { -- GitLab From 907a1c756dcee7f68943f778bdd38e67badcd7fa Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Tue, 15 Aug 2023 13:00:53 +0200 Subject: [PATCH 20/22] Update MagicNumbers.h --- Generators/TruthUtils/TruthUtils/MagicNumbers.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index d07d89a94152..53d575c929f8 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -109,7 +109,7 @@ template <class T> inline void get_particle_history(const T& p, std::deque<int>& } if (direction > 0) { if (p->status()>SIM_STATUS_INCREMENT) { - auto pv = p->decay_vertex(); + auto pv = p->end_vertex(); if (pv) { for (auto pa: pv->particles_out()) { if (pa->pdg_id() != p->pdg_id()) continue; -- GitLab From 6a1319b3ae0a982843fa367707451d6649f63505 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@mpp.mpg.de> Date: Tue, 15 Aug 2023 21:07:03 +0200 Subject: [PATCH 21/22] OK --- .../TruthUtils/TruthUtils/MagicNumbers.h | 40 +------------------ 1 file changed, 1 insertion(+), 39 deletions(-) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index df4b2bdde27b..ad157e7bac85 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -80,7 +80,7 @@ template <class T> inline void get_particle_history(const T& p, std::deque<int>& } if (direction > 0) { if (p->status()>SIM_STATUS_INCREMENT) { - auto pv = p->decay_vertex(); + auto pv = p->end_vertex(); if (pv) { for (auto pa: pv->particles_out()) { if (pa->pdg_id() != p->pdg_id()) continue; @@ -114,43 +114,5 @@ template <class T1, class T2> inline bool is_sim_descendant(const T1& p1,const T return (p1->status()/SIM_STATUS_INCREMENT >= p2->status()/SIM_STATUS_INCREMENT) && is_same_generator_particle(p1,p2); } - -template <class T> int unique_id(const T& p1){ return p1->barcode();} -#if defined(HEPMC3) && !defined(XAOD_STANDALONE) -template <> inline int unique_id(const ConstGenParticlePtr& p1){ return p1->id();} -template <> inline int unique_id(const GenParticlePtr& p1){ return p1->id();} -#endif - - -template <class T> inline void get_particle_history(const T& p, std::deque<int>& out, int direction=0) { - if (direction < 0) { - if (p->status()>SIM_STATUS_INCREMENT) { - auto pv = p->production_vertex(); - if (pv) { - for (auto pa: pv->particles_in()) { - if (pa->pdg_id() != p->pdg_id()) continue; - out.push_front(unique_id(p)); - get_particle_history(pa,out,-1); - break; - } - } - } - } - if (direction > 0) { - if (p->status()>SIM_STATUS_INCREMENT) { - auto pv = p->end_vertex(); - if (pv) { - for (auto pa: pv->particles_out()) { - if (pa->pdg_id() != p->pdg_id()) continue; - out.push_back(unique_id(p)); - get_particle_history(pa,out,1); - break; - } - } - } - } -} -template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; res.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} - } #endif -- GitLab From 37743816bb9634e448494cb79406cefe419ceac1 Mon Sep 17 00:00:00 2001 From: Andrii Verbytskyi <andrii.verbytskyi@cern.ch> Date: Fri, 18 Aug 2023 12:54:48 +0200 Subject: [PATCH 22/22] Update MagicNumbers.h --- .../TruthUtils/TruthUtils/MagicNumbers.h | 41 +++++++++++++++++++ 1 file changed, 41 insertions(+) diff --git a/Generators/TruthUtils/TruthUtils/MagicNumbers.h b/Generators/TruthUtils/TruthUtils/MagicNumbers.h index 54651be59612..a0ff87d8e310 100644 --- a/Generators/TruthUtils/TruthUtils/MagicNumbers.h +++ b/Generators/TruthUtils/TruthUtils/MagicNumbers.h @@ -147,5 +147,46 @@ template <class T> inline void get_particle_history(const T& p, std::deque<int>& } template <class T> inline std::deque<int> simulation_history(const T& p, int direction ) { std::deque<int> res; res.push_back(unique_id(p)); get_particle_history(p, res, direction); return res;} + +/********************************************/ +template <class T> old_to_new_simulation_scheme(T& evt) { + for (auto p: evt->particles()) p->set_status( (is_simulation_particle(p) ? (SIM_STATUS_THRESHOLD+SIM_STATUS_INCREMENT*generations(p)) : 0) + p->status()); + for (auto v: evt->vertices()) v->set_status((is_simulation_vertex(v) ? (SIM_STATUS_THRESHOLD+SIM_STATUS_INCREMENT*generations(v)) : 0) + v->status()); +} + +void propagate_code(HepMC::GenParticlePtr& p) { +int child_id = 0; +auto v = p->end_vertex(); +if (!v) return; +for ( auto pp: v->particles_out()) { + if (p->pdg_id() != pp->pdg_id()) continue; + HepMC::suggest_barcode(pp,HepMC::barcode(p)+SIM_STATUS_INCREMENT); + child_id = pp->id(); + break; +} +for ( auto pp: v->particles_out()) { + if (pp->id()==child_id) continue; + HepMC::suggest_barcode(pp,HepMC::barcode(p)+SIM_STATUS_THRESHOLD); +} +for ( auto pp: v->particles_out()) propagate_code(pp); +} + +template <class T> new_to_old_simulation_scheme(T& ev) { + for (auto v: evt->vertices()) { + auto pin = v->particles_in(); + if (v.size() != 1) continue; + auto p = pin.front(); + if (HepMC::is_simulation_particle(p)) continue; + bool isoriginal = false; + for (auto pp: v->particles_out()) if (is_simulation_particle(pp)) isoriginal=true; + if (!isoriginal) continue; + propagate_code(pin); + } +} + + + +/*******************************************/ + } #endif -- GitLab