diff --git a/DDG4/include/DDG4/Geant4PrimaryHandler.h b/DDG4/include/DDG4/Geant4PrimaryHandler.h index 1c43cc44bf586212842cf6ccedd16f939ac09a27..b4a098ab402fd54f02260558508d778258af7930 100644 --- a/DDG4/include/DDG4/Geant4PrimaryHandler.h +++ b/DDG4/include/DDG4/Geant4PrimaryHandler.h @@ -21,6 +21,23 @@ namespace dd4hep { /// Namespace for the Geant4 based simulation part of the AIDA detector description toolkit namespace sim { + /// Geant4PrimaryConfig to hold configuration for PrimaryHandlers + class Geant4PrimaryConfig { + public: + /// particles with these PDG IDs are not passed to geant for simulation + std::set<int> m_rejectPDGs = {1, 2, 3, 4, 5, 6, 21, 23, 24, 25}; + /// particles with these PDG IDs are not passed to geant for simulation if their properTime is zero + std::set<int> m_zeroTimePDGs = {11, 13, 15, 17}; + + std::string toString() const { + std::stringstream str; + str << "\nRejectPDGs: "; + for (int i: m_rejectPDGs) { str << i << ", "; } + str << "\nzeroTimePDGs: "; + for (int i: m_zeroTimePDGs) { str << i << ", "; } + return str.str(); + } + }; /// Geant4Action to convert the particle information to Geant4 /** Convert the primary interaction (object \tt{Geant4PrimaryInteraction} object @@ -41,8 +58,7 @@ namespace dd4hep { virtual void operator()(G4Event* event); public: - /// particles with these PDG IDs are not passed to geant for simulation - std::set<int> m_rejectPDGs={1,2,3,4,5,6,21,23,24}; + Geant4PrimaryConfig m_primaryConfig; }; } // End namespace sim diff --git a/DDG4/python/DDSim/Helper/Physics.py b/DDG4/python/DDSim/Helper/Physics.py index 150d713285c10846e512f2ea2212785066e89801..8f7863a0d082fccbf8719f3d385d70846790453c 100644 --- a/DDG4/python/DDSim/Helper/Physics.py +++ b/DDG4/python/DDSim/Helper/Physics.py @@ -14,6 +14,7 @@ class Physics( ConfigHelper ): self.decays = True self._pdgfile = None self._rejectPDGs = {1,2,3,4,5,6,21,23,24,25} + self._zeroTimePDGs = {11, 13, 15, 17} @property def rejectPDGs( self ): @@ -24,7 +25,18 @@ class Physics( ConfigHelper ): return self._rejectPDGs @rejectPDGs.setter def rejectPDGs( self, val ): - self._rejectPDGs = val + self._rejectPDGs = self.makeSet(val) + + @property + def zeroTimePDGs(self): + """Set of PDG IDs for particles that should not be passed to Geant4 if their properTime is 0. + + The properTime of 0 indicates a documentation to add FSR to a lepton for example. + """ + return self._zeroTimePDGs + @zeroTimePDGs.setter + def zeroTimePDGs(self, val): + self._zeroTimePDGs = self.makeSet(val) @property def rangecut( self ): diff --git a/DDG4/src/Geant4InputHandling.cpp b/DDG4/src/Geant4InputHandling.cpp index f3de6c4d6e3507719b5ff8d83970320740ae009c..24383ac3585bae7a0c06a64e0798dfe509f4bca9 100644 --- a/DDG4/src/Geant4InputHandling.cpp +++ b/DDG4/src/Geant4InputHandling.cpp @@ -359,7 +359,7 @@ static vector< pair<Geant4Particle*,G4PrimaryParticle*> > getRelevant(set<int>& visited, map<int,G4PrimaryParticle*>& prim, Geant4PrimaryInteraction::ParticleMap& pm, - const set<int>& rejectPDGs, + const Geant4PrimaryConfig& primaryConfig, const Geant4ParticleHandle p) { typedef vector< pair<Geant4Particle*,G4PrimaryParticle*> > Primaries; @@ -384,12 +384,11 @@ getRelevant(set<int>& visited, double proper_time_Precision = pow(10.,-DBL_DIG)*me*fmax(fabs(p->time),fabs(dp->time)); bool isProperTimeZero = (proper_time <= proper_time_Precision); - const std::set<int> leptonPDGs{11,13,15,17}; // -- remove original if --- bool rejectParticle = not p.definition() // completely unknown to geant4 - or (rejectPDGs.count(abs(p->pdgID)) != 0) // quarks, gluon, "strings", W, Z etc. + or (primaryConfig.m_rejectPDGs.count(abs(p->pdgID)) != 0) // quarks, gluon, "strings", W, Z etc. or (isProperTimeZero and p.definition()->GetPDGStable() ) // initial state electrons, etc. - or (isProperTimeZero and leptonPDGs.count(abs(p->pdgID)) != 0 ) ; // charged 'documentation' leptons, e.g. in lepton pairs w/ FSR + or (isProperTimeZero and primaryConfig.m_zeroTimePDGs.count(abs(p->pdgID)) != 0 ) ; // charged 'documentation' leptons, e.g. in lepton pairs w/ FSR if (not rejectParticle) { map<int,G4PrimaryParticle*>::iterator ip4 = prim.find(p->id); G4PrimaryParticle* p4 = (ip4 == prim.end()) ? 0 : (*ip4).second; @@ -400,7 +399,7 @@ getRelevant(set<int>& visited, Primaries daughters; for(Geant4Particle::Particles::const_iterator i=dau.begin(); i!=dau.end(); ++i) { if ( visited.find(*i) == visited.end() ) { - Primaries tmp = getRelevant(visited,prim,pm,rejectPDGs,pm[*i]); + Primaries tmp = getRelevant(visited,prim,pm,primaryConfig,pm[*i]); daughters.insert(daughters.end(), tmp.begin(),tmp.end()); } } @@ -412,7 +411,7 @@ getRelevant(set<int>& visited, else { for(Geant4Particle::Particles::const_iterator i=dau.begin(); i!=dau.end(); ++i) { if ( visited.find(*i) == visited.end() ) { - Primaries tmp = getRelevant(visited,prim,pm,rejectPDGs,pm[*i]); + Primaries tmp = getRelevant(visited,prim,pm,primaryConfig,pm[*i]); res.insert(res.end(), tmp.begin(),tmp.end()); } } @@ -436,13 +435,9 @@ int dd4hep::sim::generatePrimaries(const Geant4Action* caller, set<int> visited; auto const* primHandler = dynamic_cast<const Geant4PrimaryHandler*>(caller); - auto const& rejectPDGs = primHandler ? primHandler->m_rejectPDGs : std::set<int>(); + auto const& primaryConfig = primHandler ? primHandler->m_primaryConfig : Geant4PrimaryConfig(); - caller->debug("Rejecting PDGs: %s", [&rejectPDGs]{ - std::stringstream str; - for (int i: rejectPDGs) { str << i << ", "; } - return str.str(); - }().c_str()); + caller->debug("PrimaryConfiguration:%s", primaryConfig.toString().c_str()); if ( interaction->locked ) { caller->abortRun("Locked interactions may not be used to generate primaries!", @@ -466,7 +461,7 @@ int dd4hep::sim::generatePrimaries(const Geant4Action* caller, mask.set(G4PARTICLE_HAS_SECONDARIES); } if ( p->parents.size() == 0 ) { - Primaries relevant = getRelevant(visited,prim,pm,rejectPDGs,p); + Primaries relevant = getRelevant(visited,prim,pm,primaryConfig,p); for(Primaries::const_iterator j=relevant.begin(); j!= relevant.end(); ++j) { Geant4ParticleHandle r = (*j).first; G4PrimaryParticle* p4 = (*j).second; diff --git a/DDG4/src/Geant4PrimaryHandler.cpp b/DDG4/src/Geant4PrimaryHandler.cpp index 8c98a3e620e757508cfaa9e1b35f942fbb8e77e3..4bb741052c7e2bd83c2bfd9450f64aeddb19b0fb 100644 --- a/DDG4/src/Geant4PrimaryHandler.cpp +++ b/DDG4/src/Geant4PrimaryHandler.cpp @@ -24,7 +24,8 @@ Geant4PrimaryHandler::Geant4PrimaryHandler(Geant4Context* ctxt, const std::strin : Geant4GeneratorAction(ctxt,nam) { InstanceCount::increment(this); - declareProperty("RejectPDGs", m_rejectPDGs); + declareProperty("RejectPDGs", m_primaryConfig.m_rejectPDGs); + declareProperty("ZeroTimePDGs", m_primaryConfig.m_zeroTimePDGs); } /// Default destructor