From 1a6da4d8366f839d352bf56df20ed1c8cd121dfd Mon Sep 17 00:00:00 2001 From: David Nathan Brown Date: Fri, 17 Nov 2023 21:20:04 -0800 Subject: [PATCH 1/4] Add entries to the MCSim branch for particles that contribute hits to the track, up to the configurable limit (default 1). add nhits and nactive leafs to the MCSim branch. Restore the gen and proc leafs --- inc/InfoMCStructHelper.hh | 2 +- inc/SimInfo.hh | 12 ++++-- src/InfoMCStructHelper.cc | 71 +++++++++++++++++++---------------- src/TrkAnaTreeMaker_module.cc | 3 +- 4 files changed, 49 insertions(+), 39 deletions(-) diff --git a/inc/InfoMCStructHelper.hh b/inc/InfoMCStructHelper.hh index e0afdc0..ac0f608 100644 --- a/inc/InfoMCStructHelper.hh +++ b/inc/InfoMCStructHelper.hh @@ -60,7 +60,7 @@ namespace mu2e { void fillTrkInfoMC(const KalSeed& kseed, const KalSeedMC& kseedmc, TrkInfoMC& trkinfomc); void fillTrkInfoMCDigis(const KalSeed& kseed, const KalSeedMC& kseedmc, TrkInfoMC& trkinfomc); void fillHitInfoMC(const KalSeedMC& kseedmc, TrkStrawHitInfoMC& tshinfomc, const TrkStrawHitMC& tshmc); - void fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector& siminfos, int n_generations); + void fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector& siminfos, int n_generations, int n_match); void fillVDInfo(KalSeed const& kseed, const KalSeedMC& kseedmc, std::vector& vdinfos); void fillHitInfoMCs(const KalSeedMC& kseedmc, std::vector& tshinfomcs); void fillCaloClusterInfoMC(CaloClusterMC const& ccmc, CaloClusterInfoMC& ccimc); diff --git a/inc/SimInfo.hh b/inc/SimInfo.hh index 2ded05a..474af5d 100644 --- a/inc/SimInfo.hh +++ b/inc/SimInfo.hh @@ -13,13 +13,17 @@ namespace mu2e // general info about the SimParticle which was simulated struct SimInfo { bool valid = false; // true/false as to whether the data is valid + Int_t nhits = 0; // number of matched hits in the track generated by this particle + Int_t nactive = 0; // number of active matched hits in the track generated by this particle + Int_t rank = -1; // rank of this particle, by # of associated hits Int_t pdg = -1; // true PDG code + Int_t proc = -1; // process code Int_t gen = -1; // generator code - Float_t time = -1.0; // time of this step - XYZVectorF mom = XYZVectorF(); // particle momentum at the start of this step - XYZVectorF pos = XYZVectorF(); // particle position at the start of this step + Float_t time = -1.0; // Origin time of the SimParticle + XYZVectorF mom = XYZVectorF(); // origin momentumof the SimParticle + XYZVectorF pos = XYZVectorF(); // origin position of the SimParticle MCRelationship prirel = MCRelationship::none; // relationship to the event primary particles - MCRelationship trkrel = MCRelationship::none; // relationship to the particle that created the track + MCRelationship trkrel = MCRelationship::none; // relationship to the particle that created hits in the track void reset() { *this = SimInfo(); } }; } diff --git a/src/InfoMCStructHelper.cc b/src/InfoMCStructHelper.cc index 63f1c00..d02c1fa 100644 --- a/src/InfoMCStructHelper.cc +++ b/src/InfoMCStructHelper.cc @@ -125,39 +125,45 @@ namespace mu2e { tshinfomc.doca = -1*dperp; } - void InfoMCStructHelper::fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector& siminfos, int n_generations) { - auto trkprimaryptr = kseedmc.simParticle().simParticle(_spcH); - auto trkprimary = trkprimaryptr->originParticle(); - - auto current_sim_particle_ptr = trkprimaryptr; - auto current_sim_particle = trkprimary; - if (n_generations == -1) { // means do all generations - n_generations = std::numeric_limits::max(); - } + void InfoMCStructHelper::fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector& siminfos, int n_generations, int n_match) { + for(int imatch = 0 ; imatch < std::min(n_match,static_cast(kseedmc.simParticles().size())); ++imatch) { + auto trkprimaryptr = kseedmc.simParticle(imatch).simParticle(_spcH); + auto trkprimary = trkprimaryptr->originParticle(); + auto current_sim_particle_ptr = trkprimaryptr; + auto current_sim_particle = trkprimary; + if (n_generations == -1) { // means do all generations + n_generations = std::numeric_limits::max(); + } - for (int i_generation = 0; i_generation < n_generations; ++i_generation) { - SimInfo sim_info; - fillSimInfo(current_sim_particle, sim_info); - sim_info.trkrel = MCRelationship(current_sim_particle_ptr, trkprimaryptr); + for (int i_generation = 0; i_generation < n_generations; ++i_generation) { + SimInfo sim_info; + fillSimInfo(current_sim_particle, sim_info); + sim_info.trkrel = MCRelationship(current_sim_particle_ptr, trkprimaryptr); + sim_info.rank = imatch; + if(i_generation == 0){ + sim_info.nhits = kseedmc.simParticle(imatch)._nhits; + sim_info.nactive = kseedmc.simParticle(imatch)._nactive; + } - auto bestprimarysp = primary.primarySimParticles().front(); - MCRelationship bestrel; - for(auto const& spp : primary.primarySimParticles()){ - MCRelationship mcrel(current_sim_particle_ptr, spp); - if(mcrel > bestrel){ - bestrel = mcrel; - bestprimarysp = spp; + auto bestprimarysp = primary.primarySimParticles().front(); + MCRelationship bestrel; + for(auto const& spp : primary.primarySimParticles()){ + MCRelationship mcrel(current_sim_particle_ptr, spp); + if(mcrel > bestrel){ + bestrel = mcrel; + bestprimarysp = spp; + } } - } - sim_info.prirel = bestrel; + sim_info.prirel = bestrel; - siminfos.push_back(sim_info); - if (current_sim_particle.parent().isNonnull()) { - current_sim_particle_ptr = current_sim_particle.parent(); - current_sim_particle = current_sim_particle_ptr->originParticle(); - } - else { - break; // this particle doesn't have a parent + siminfos.push_back(sim_info); + if (current_sim_particle.parent().isNonnull()) { + current_sim_particle_ptr = current_sim_particle.parent(); + current_sim_particle = current_sim_particle_ptr->originParticle(); + } + else { + break; // this particle doesn't have a parent + } } } @@ -175,12 +181,12 @@ namespace mu2e { } } if (!already_added) { + auto trkprimaryptr = kseedmc.simParticle().simParticle(_spcH); sim_info.trkrel = MCRelationship(spp, trkprimaryptr); sim_info.prirel = MCRelationship(spp, spp); siminfos.push_back(sim_info); } } - } @@ -194,12 +200,11 @@ namespace mu2e { } void InfoMCStructHelper::fillSimInfo(const SimParticle& sp, SimInfo& siminfo) { - GeomHandle det; - siminfo.valid = true; + if(sp.genParticle().isNonnull())siminfo.gen = sp.genParticle()->generatorId().id(); + siminfo.proc = sp.creationCode(); siminfo.pdg = sp.pdgId(); - siminfo.gen = sp.creationCode(); siminfo.time = sp.startGlobalTime(); siminfo.mom = XYZVectorF(sp.startMomentum()); siminfo.pos = XYZVectorF(det->toDetector(sp.startPosition())); diff --git a/src/TrkAnaTreeMaker_module.cc b/src/TrkAnaTreeMaker_module.cc index 58a5df3..e572864 100644 --- a/src/TrkAnaTreeMaker_module.cc +++ b/src/TrkAnaTreeMaker_module.cc @@ -113,6 +113,7 @@ namespace mu2e { fhicl::Atom filltrkpid{Name("fillTrkPID"), Comment("Switch to turn on filling of the full TrkPIDInfo for this set of tracks"), false}; fhicl::Atom required{Name("required"), Comment("True/false if you require this type of track in the event"), false}; fhicl::Atom genealogyDepth{Name("genealogyDepth"), Comment("The depth of the genealogy information you want to keep"), 1}; + fhicl::Atom matchDepth{Name("matchDepth"), Comment("The depth into the MC true particle matching you want to keep"), 1}; }; struct BranchConfig { @@ -854,7 +855,7 @@ namespace mu2e { _infoMCStructHelper.fillTrkInfoMC(kseed, kseedmc, _allMCTIs.at(i_branch)); auto& mcvdis = _allMCVDInfos.at(i_branch); _infoMCStructHelper.fillVDInfo(kseed, kseedmc, mcvdis); - _infoMCStructHelper.fillAllSimInfos(kseedmc, primary, _allMCSimTIs.at(i_branch), branchConfig.options().genealogyDepth()); + _infoMCStructHelper.fillAllSimInfos(kseedmc, primary, _allMCSimTIs.at(i_branch), branchConfig.options().genealogyDepth(), branchConfig.options().matchDepth()); if(_conf.diag() > 1 || (_conf.fillhits() && branchConfig.options().fillhits())){ _infoMCStructHelper.fillHitInfoMCs(kseedmc, _allTSHIMCs.at(i_branch)); From 0d63d948381382c1470acfebcc41507f53efd3ea Mon Sep 17 00:00:00 2001 From: David Nathan Brown Date: Sun, 19 Nov 2023 00:21:27 -0600 Subject: [PATCH 2/4] Add DeP trigger config. Allow -1 as a match depth --- fcl/prolog_trigger.fcl | 33 ++++++++++++++++++++++++++------- inc/SimInfo.hh | 4 ++-- src/InfoMCStructHelper.cc | 19 ++++++++++++------- 3 files changed, 40 insertions(+), 16 deletions(-) diff --git a/fcl/prolog_trigger.fcl b/fcl/prolog_trigger.fcl index cf1ae1d..4d247d7 100644 --- a/fcl/prolog_trigger.fcl +++ b/fcl/prolog_trigger.fcl @@ -6,13 +6,6 @@ BEGIN_PROLOG TrkAnaTrigger : { TrkAnaTT : { @table::TrkAnaTreeMaker - candidate : { - @table::DeM - options : { - @table::AllOpt - fillHits : true - } - } supplements : [] ExtraMCStepCollectionTags : [] diagLevel : 2 @@ -37,5 +30,31 @@ TrkAnaTrigger : { } } } +TrkAnaTrigger : { + @table::TrkAnaTrigger + TrkAnaDeMTT : { + @table::TrkAnaTrigger.TrkAnaTT + candidate : { + @table::DeM + options : { + @table::AllOpt + fillHits : true + } + } + } + + TrkAnaDePTT : { + @table::TrkAnaTrigger.TrkAnaTT + candidate : { + @table::DeP + options : { + @table::AllOpt + fillHits : true + } + } + } + + +} END_PROLOG diff --git a/inc/SimInfo.hh b/inc/SimInfo.hh index 474af5d..3124be7 100644 --- a/inc/SimInfo.hh +++ b/inc/SimInfo.hh @@ -13,8 +13,8 @@ namespace mu2e // general info about the SimParticle which was simulated struct SimInfo { bool valid = false; // true/false as to whether the data is valid - Int_t nhits = 0; // number of matched hits in the track generated by this particle - Int_t nactive = 0; // number of active matched hits in the track generated by this particle + Int_t nhits = -1; // number of matched hits in the track generated by this particle + Int_t nactive = -1; // number of active matched hits in the track generated by this particle Int_t rank = -1; // rank of this particle, by # of associated hits Int_t pdg = -1; // true PDG code Int_t proc = -1; // process code diff --git a/src/InfoMCStructHelper.cc b/src/InfoMCStructHelper.cc index d02c1fa..e76d3a4 100644 --- a/src/InfoMCStructHelper.cc +++ b/src/InfoMCStructHelper.cc @@ -126,24 +126,24 @@ namespace mu2e { } void InfoMCStructHelper::fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector& siminfos, int n_generations, int n_match) { + // interpret -1 as no llimit + if (n_generations == -1) { + n_generations = std::numeric_limits::max(); + } + if (n_match == -1 ) { + n_match = std::numeric_limits::max(); + } for(int imatch = 0 ; imatch < std::min(n_match,static_cast(kseedmc.simParticles().size())); ++imatch) { auto trkprimaryptr = kseedmc.simParticle(imatch).simParticle(_spcH); auto trkprimary = trkprimaryptr->originParticle(); auto current_sim_particle_ptr = trkprimaryptr; auto current_sim_particle = trkprimary; - if (n_generations == -1) { // means do all generations - n_generations = std::numeric_limits::max(); - } for (int i_generation = 0; i_generation < n_generations; ++i_generation) { SimInfo sim_info; fillSimInfo(current_sim_particle, sim_info); sim_info.trkrel = MCRelationship(current_sim_particle_ptr, trkprimaryptr); sim_info.rank = imatch; - if(i_generation == 0){ - sim_info.nhits = kseedmc.simParticle(imatch)._nhits; - sim_info.nactive = kseedmc.simParticle(imatch)._nactive; - } auto bestprimarysp = primary.primarySimParticles().front(); MCRelationship bestrel; @@ -155,6 +155,11 @@ namespace mu2e { } } sim_info.prirel = bestrel; + // only count hits for direct contributors + if(i_generation == 0){ + sim_info.nhits = kseedmc.simParticle(imatch)._nhits; + sim_info.nactive = kseedmc.simParticle(imatch)._nactive; + } siminfos.push_back(sim_info); if (current_sim_particle.parent().isNonnull()) { From cf3b1b30df0c4f12f50f98ec21be463bbb2f1291 Mon Sep 17 00:00:00 2001 From: David Nathan Brown Date: Mon, 20 Nov 2023 10:05:42 -0600 Subject: [PATCH 3/4] Change branch name to 'trk' --- fcl/prolog_trigger.fcl | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/fcl/prolog_trigger.fcl b/fcl/prolog_trigger.fcl index 4d247d7..721e1e1 100644 --- a/fcl/prolog_trigger.fcl +++ b/fcl/prolog_trigger.fcl @@ -35,7 +35,8 @@ TrkAnaTrigger : { TrkAnaDeMTT : { @table::TrkAnaTrigger.TrkAnaTT candidate : { - @table::DeM + branch: "trk" + suffix : "DeM" options : { @table::AllOpt fillHits : true @@ -46,7 +47,8 @@ TrkAnaTrigger : { TrkAnaDePTT : { @table::TrkAnaTrigger.TrkAnaTT candidate : { - @table::DeP + branch: "trk" + suffix : "DeP" options : { @table::AllOpt fillHits : true From 38c6bfaa6f014c3af2b1e5219c383cf36483c9c3 Mon Sep 17 00:00:00 2001 From: David Nathan Brown Date: Tue, 21 Nov 2023 16:05:47 -0600 Subject: [PATCH 4/4] Add index leaf to allow navigation to sub-ranked genealogy sim info. Allow matching multiple sim particles by default --- fcl/prolog.fcl | 4 +++- inc/SimInfo.hh | 1 + src/InfoMCStructHelper.cc | 4 +++- 3 files changed, 7 insertions(+), 2 deletions(-) diff --git a/fcl/prolog.fcl b/fcl/prolog.fcl index 9d072e8..cd1ca37 100644 --- a/fcl/prolog.fcl +++ b/fcl/prolog.fcl @@ -97,13 +97,15 @@ dioLLWeight : { } genCountLogger: { module_type: GenEventCountReader } -AllOpt : { fillMC : true +AllOpt : { + fillMC : true trkqual : "TrkQual" fillTrkQual : true trkpid : "TrkPID" fillTrkPID : true fillHits : true genealogyDepth : -1 + matchDepth : -1 } DeM : { input : "KK" diff --git a/inc/SimInfo.hh b/inc/SimInfo.hh index 3124be7..0615e42 100644 --- a/inc/SimInfo.hh +++ b/inc/SimInfo.hh @@ -20,6 +20,7 @@ namespace mu2e Int_t proc = -1; // process code Int_t gen = -1; // generator code Float_t time = -1.0; // Origin time of the SimParticle + Int_t index; // index into the SimInfo vector XYZVectorF mom = XYZVectorF(); // origin momentumof the SimParticle XYZVectorF pos = XYZVectorF(); // origin position of the SimParticle MCRelationship prirel = MCRelationship::none; // relationship to the event primary particles diff --git a/src/InfoMCStructHelper.cc b/src/InfoMCStructHelper.cc index e76d3a4..8159644 100644 --- a/src/InfoMCStructHelper.cc +++ b/src/InfoMCStructHelper.cc @@ -160,7 +160,8 @@ namespace mu2e { sim_info.nhits = kseedmc.simParticle(imatch)._nhits; sim_info.nactive = kseedmc.simParticle(imatch)._nactive; } - + // record the index this object will have + sim_info.index = siminfos.size(); siminfos.push_back(sim_info); if (current_sim_particle.parent().isNonnull()) { current_sim_particle_ptr = current_sim_particle.parent(); @@ -189,6 +190,7 @@ namespace mu2e { auto trkprimaryptr = kseedmc.simParticle().simParticle(_spcH); sim_info.trkrel = MCRelationship(spp, trkprimaryptr); sim_info.prirel = MCRelationship(spp, spp); + sim_info.index = siminfos.size(); siminfos.push_back(sim_info); } }