From 442c864a89e0195547f4cc1bcf81ec91fd254a85 Mon Sep 17 00:00:00 2001 From: Thomas Madlener Date: Mon, 6 Dec 2021 15:56:26 +0100 Subject: [PATCH 1/4] Use auto where possible to avoid explicitly mentioning edm4hep types --- DDG4/edm4hep/Geant4Output2EDM4hep.cpp | 11 ++++------- 1 file changed, 4 insertions(+), 7 deletions(-) diff --git a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp index 504fb8e4f..0fbbe8d15 100644 --- a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp +++ b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp @@ -281,7 +281,6 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { size_t cnt = 0; map p_ids; vector p_part(pm.size(),0); - vector p_edm4hep(pm.size()); // First create the particles for(ParticleMap::const_iterator i=pm.begin(); i!=pm.end();++i, ++cnt) { int id = (*i).first; @@ -289,7 +288,7 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { PropertyMask mask(p->status); // std::cout << " ********** mcp status : 0x" << std::hex << p->status << ", mask.isSet(G4PARTICLE_GEN_STABLE) x" << std::dec << mask.isSet(G4PARTICLE_GEN_STABLE) <create(); mcp.setPDG(p->pdgID); float ps_fa[3] = {float(p->psx/CLHEP::GeV),float(p->psy/CLHEP::GeV),float(p->psz/CLHEP::GeV)}; @@ -337,17 +336,15 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { mcp.setSpin(p->spin); mcp.setColorFlow(p->colorFlow); - mcpc->push_back(mcp); p_ids[id] = cnt; p_part[cnt] = p; - p_edm4hep[cnt] = mcp; } // Now establish parent-daughter relationships for(size_t i=0, n=p_ids.size(); i::iterator k; const Geant4Particle* p = p_part[i]; - edm4hep::MCParticle q = p_edm4hep[i]; + auto q = (*mcpc)[i]; const Geant4Particle::Particles& dau = p->daughters; for(Geant4Particle::Particles::const_iterator j=dau.begin(); j!=dau.end(); ++j) { int idau = *j; @@ -356,7 +353,7 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { continue; } int iqdau = (*k).second; - edm4hep::MCParticle qdau = p_edm4hep[iqdau]; + auto qdau = (*mcpc)[iqdau]; qdau.addToParents(q); } const Geant4Particle::Particles& par = p->parents; @@ -367,7 +364,7 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { continue; } int iqpar = (*k).second; - edm4hep::MCParticle qpar = p_edm4hep[iqpar]; + auto qpar = (*mcpc)[iqpar]; q.addToParents(qpar); } } From 3f12c7d4b20889a13bc30cf6912a5664cac41ae0 Mon Sep 17 00:00:00 2001 From: Thomas Madlener Date: Tue, 7 Dec 2021 11:08:21 +0100 Subject: [PATCH 2/4] Use map to keep track of collections and avoid const_cast --- DDG4/edm4hep/Geant4Output2EDM4hep.cpp | 52 ++++++++++++++------------- 1 file changed, 27 insertions(+), 25 deletions(-) diff --git a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp index 0fbbe8d15..9fc61b294 100644 --- a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp +++ b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp @@ -43,7 +43,7 @@ #include #include #include - +#include /// Namespace for the AIDA detector description toolkit namespace dd4hep { @@ -88,6 +88,8 @@ namespace dd4hep { std::map< std::string, std::string > m_eventParametersString; bool m_FirstEvent = true ; + std::unordered_map m_collections; + /// create the podio collections for the particles and hits void createCollections(OutputContext& ctxt) ; /// Data conversion interface for MC particles to EDM4hep format @@ -273,9 +275,7 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { typedef Geant4ParticleMap::ParticleMap ParticleMap; const ParticleMap& pm = particles->particleMap; size_t nparts = pm.size(); - edm4hep::MCParticleCollection* mcpc = - const_cast( - &m_store->get("MCParticles")); + auto mcpc = static_cast(m_collections["MCParticles"]); if ( nparts > 0 ) { size_t cnt = 0; @@ -396,8 +396,8 @@ void Geant4Output2EDM4hep::saveEvent(OutputContext& ctxt) { // this does not compile as create() is we only get a const ref - need to review PODIO EventStore API // auto& evtHCol = m_store->get("EventHeader") ; - // auto evtHdr = evtHCol.create() ; - auto* evtHCol = const_cast(&m_store->get("EventHeader") ); +// auto evtHdr = evtHCol.create() ; + auto* evtHCol = static_cast(m_collections["EventHeader"]); auto evtHdr = evtHCol->create() ; evtHdr.setRunNumber(runNumber); @@ -435,15 +435,12 @@ void Geant4Output2EDM4hep::saveCollection(OutputContext& /*ctxt*/, G4VH Geant4ParticleMap* pm = context()->event().extension(false); - edm4hep::MCParticleCollection* mcpc = - const_cast( - &m_store->get("MCParticles")); + auto* mcpc = static_cast(m_collections["MCParticles"]); //------------------------------------------------------------------- if( typeid( Geant4Tracker::Hit ) == coll->type().type() ){ - edm4hep::SimTrackerHitCollection* sthc = - const_cast(&m_store->get(colName)); + auto* sthc = static_cast(m_collections[colName]); for(unsigned i=0 ; i < nhits ; ++i){ auto sth = sthc->create() ; @@ -481,15 +478,11 @@ void Geant4Output2EDM4hep::saveCollection(OutputContext& /*ctxt*/, G4VH Geant4Sensitive* sd = coll->sensitive(); int hit_creation_mode = sd->hitCreationMode(); - edm4hep::SimCalorimeterHitCollection* sCaloHitColl = - const_cast( - &m_store->get(colName)); + auto* sCaloHitColl = static_cast(m_collections[colName]); colName += "Contributions" ; - edm4hep::CaloHitContributionCollection* sCaloHitContColl = - const_cast( - &m_store->get(colName)); + auto* sCaloHitContColl = static_cast(m_collections[colName]); for(unsigned i=0 ; i < nhits ; ++i){ @@ -536,14 +529,17 @@ void Geant4Output2EDM4hep::saveCollection(OutputContext& /*ctxt*/, G4VH } } - void Geant4Output2EDM4hep::createCollections(OutputContext& ctxt){ - m_store->create("MCParticles"); + auto* mcColl = new edm4hep::MCParticleCollection(); + m_collections.emplace("MCParticles", mcColl); + m_store->registerCollection("MCParticles", mcColl); m_file->registerForWrite("MCParticles"); printout(DEBUG,"Geant4Output2EDM4hep","+++ created collection MCParticles" ); - m_store->create("EventHeader") ; + auto* evtHeader = new edm4hep::EventHeaderCollection(); + m_collections.emplace("EventHeader", evtHeader); + m_store->registerCollection("EventHeader", evtHeader); m_file->registerForWrite("EventHeader"); printout(DEBUG,"Geant4Output2EDM4hep","+++ created collection EventHeader" ); @@ -566,22 +562,28 @@ void Geant4Output2EDM4hep::createCollections(OutputContext& ctxt){ if( typeid( Geant4Tracker::Hit ) == coll->type().type() ){ - auto& sthc = m_store->create(colName); + auto* sthc = new edm4hep::SimTrackerHitCollection(); + m_collections.emplace(colName, sthc); + m_store->registerCollection(colName, sthc); m_file->registerForWrite(colName); - auto& sthc_md = m_store->getCollectionMetaData( sthc.getID() ); + auto& sthc_md = m_store->getCollectionMetaData( sthc->getID() ); sthc_md.setValue("CellIDEncodingString", sd_enc); printout(DEBUG,"Geant4Output2EDM4hep","+++ created collection %s",colName.c_str() ); } else if( typeid( Geant4Calorimeter::Hit ) == coll->type().type() ){ - auto& schc = m_store->create(colName); + auto* schc = new edm4hep::SimCalorimeterHitCollection(); + m_collections.emplace(colName, schc); + m_store->registerCollection(colName, schc); m_file->registerForWrite(colName); - auto& schc_md = m_store->getCollectionMetaData( schc.getID() ); + auto& schc_md = m_store->getCollectionMetaData( schc->getID() ); schc_md.setValue("CellIDEncodingString", sd_enc); printout(DEBUG,"Geant4Output2EDM4hep","+++ created collection %s",colName.c_str() ); colName += "Contributions" ; - m_store->create(colName); + auto* chContribColl = new edm4hep::CaloHitContributionCollection(); + m_collections.emplace(colName, chContribColl); + m_store->registerCollection(colName, chContribColl); m_file->registerForWrite(colName); printout(DEBUG,"Geant4Output2EDM4hep","+++ created collection %s",colName.c_str() ); From 81107453a154a5e87cd437d467086189717577ef Mon Sep 17 00:00:00 2001 From: Thomas Madlener Date: Tue, 7 Dec 2021 13:23:49 +0100 Subject: [PATCH 3/4] Switch to range based for loops where possible --- DDG4/edm4hep/Geant4Output2EDM4hep.cpp | 48 ++++++++++++++------------- 1 file changed, 25 insertions(+), 23 deletions(-) diff --git a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp index 9fc61b294..15c2bd838 100644 --- a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp +++ b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp @@ -274,17 +274,18 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { typedef detail::ReferenceBitMask PropertyMask; typedef Geant4ParticleMap::ParticleMap ParticleMap; const ParticleMap& pm = particles->particleMap; - size_t nparts = pm.size(); auto mcpc = static_cast(m_collections["MCParticles"]); - if ( nparts > 0 ) { + if ( pm.size() > 0 ) { size_t cnt = 0; + // Mapping of ids in the ParticleMap to indices in the MCParticle collection map p_ids; - vector p_part(pm.size(),0); + vector p_part; + p_part.reserve(pm.size()); // First create the particles - for(ParticleMap::const_iterator i=pm.begin(); i!=pm.end();++i, ++cnt) { - int id = (*i).first; - const Geant4ParticleHandle p = (*i).second; + for (const auto& iParticle : pm) { + int id = iParticle.first; + const Geant4ParticleHandle p = iParticle.second; PropertyMask mask(p->status); // std::cout << " ********** mcp status : 0x" << std::hex << p->status << ", mask.isSet(G4PARTICLE_GEN_STABLE) x" << std::dec << mask.isSet(G4PARTICLE_GEN_STABLE) <spin); mcp.setColorFlow(p->colorFlow); - p_ids[id] = cnt; - p_part[cnt] = p; + p_ids[id] = cnt++; + p_part.push_back(p); } // Now establish parent-daughter relationships - for(size_t i=0, n=p_ids.size(); i::iterator k; + for(size_t i=0; i < p_ids.size(); ++i) { const Geant4Particle* p = p_part[i]; auto q = (*mcpc)[i]; - const Geant4Particle::Particles& dau = p->daughters; - for(Geant4Particle::Particles::const_iterator j=dau.begin(); j!=dau.end(); ++j) { - int idau = *j; - if ( (k=p_ids.find(idau)) == p_ids.end() ) { // Error!!! + + for (const auto& idau : p->daughters) { + const auto k = p_ids.find(idau); + if (k == p_ids.end()) { printout(FATAL,"Geant4Conversion","+++ Particle %d: FAILED to find daughter with ID:%d",p->id,idau); continue; } @@ -356,16 +356,18 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { auto qdau = (*mcpc)[iqdau]; qdau.addToParents(q); } - const Geant4Particle::Particles& par = p->parents; - for(Geant4Particle::Particles::const_iterator j=par.begin(); j!=par.end(); ++j) { - int ipar = *j; // A parent ID iof -1 means NO parent, because a base of 0 is perfectly leagal! - if ( ipar>=0 && (k=p_ids.find(ipar)) == p_ids.end() ) { // Error!!! - printout(FATAL,"Geant4Conversion","+++ Particle %d: FAILED to find parent with ID:%d",p->id,ipar); - continue; + + for (const auto& ipar : p->parents) { + if (ipar >= 0) { // A parent ID of -1 means NO parent, because a base of 0 is perfectly legal + const auto k = p_ids.find(ipar); + if (k == p_ids.end()) { + printout(FATAL,"Geant4Conversion","+++ Particle %d: FAILED to find parent with ID:%d",p->id,ipar); + continue; + } + int iqpar = (*k).second; + auto qpar = (*mcpc)[iqpar]; + q.addToParents(qpar); } - int iqpar = (*k).second; - auto qpar = (*mcpc)[iqpar]; - q.addToParents(qpar); } } } From d9294130b30aa425d155b5cf6d12aecaea253d4e Mon Sep 17 00:00:00 2001 From: Thomas Madlener Date: Tue, 7 Dec 2021 13:26:14 +0100 Subject: [PATCH 4/4] Set daugther relations as well This doesn't happen automatically in EDM4hep but does in LCIO! --- DDG4/edm4hep/Geant4Output2EDM4hep.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp index 15c2bd838..554937cbc 100644 --- a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp +++ b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp @@ -355,6 +355,7 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { int iqdau = (*k).second; auto qdau = (*mcpc)[iqdau]; qdau.addToParents(q); + q.addToDaughters(qdau); } for (const auto& ipar : p->parents) { @@ -367,6 +368,7 @@ void Geant4Output2EDM4hep::saveParticles(Geant4ParticleMap* particles) { int iqpar = (*k).second; auto qpar = (*mcpc)[iqpar]; q.addToParents(qpar); + qpar.addToDaughters(q); } } }