Merge pull request #116 from AndrewEdmonds11/v05-dev

v05 development

README.md

@@ -2,11 +2,11 @@
 
 ## Introduction
 
-TrkAna is a track-based ntuple for Mu2e analyses. Each entry in the ntuple corresponds to a single reconstructed track. 
+TrkAna is an event-based ntuple for Mu2e analyses. Each entry in the ntuple corresponds to a single Mu2e event. 
 
 ## Useful Links
 
 * [Mu2eWiki page](https://mu2ewiki.fnal.gov/wiki/TrkAna)
 * [spreadsheet explaining branch and leaf names](https://docs.google.com/spreadsheets/d/1IMZVgj74vxGTWqAFRtFEEcYChujJyc38VcopZh52jhc/edit)
 * [version history of TrkAna](https://mu2ewiki.fnal.gov/wiki/TrkAna#Version_History_.2F_TrkAna_Musings)
-* [tutorial](tutorial/README.md)
+* [tutorial](tutorial/README.md) (some but not all updated for v5)

fcl/TrkAnaReco.fcl

@@ -22,27 +22,19 @@ physics.TrkAnaTrigPath : [ PBIWeight ]
 physics.TrkAnaEndPath : [ @sequence::TrkAnaReco.EndSequence ]
 
 # Include more information (MC, full TrkQual and TrkPID branches)
-physics.analyzers.TrkAnaNeg.candidate.options : @local::AllOpt
-physics.analyzers.TrkAnaPos.candidate.options : @local::AllOpt
+# TODO: add these options back
+#physics.analyzers.TrkAnaNeg.candidate.options : @local::AllOpt
+#physics.analyzers.TrkAnaPos.candidate.options : @local::AllOpt
 
 # for hit level diagnostics, set diagLevel to 2
-physics.analyzers.TrkAnaNeg.diagLevel : 1
-physics.analyzers.TrkAnaNeg.FillMCInfo : true
-physics.analyzers.TrkAnaNeg.FillCRVHits : true
-physics.analyzers.TrkAnaNeg.FillTriggerInfo : false
-physics.analyzers.TrkAnaNeg.FitType : LoopHelix
-
-physics.analyzers.TrkAnaPos.diagLevel : 1
-physics.analyzers.TrkAnaPos.FillMCInfo : true
-physics.analyzers.TrkAnaPos.FillCRVHits : true
-physics.analyzers.TrkAnaPos.FillTriggerInfo : false
-physics.analyzers.TrkAnaPos.FitType : LoopHelix
-
-physics.analyzers.TrkAnaNeg.FillTrkQualInfo : false
-physics.analyzers.TrkAnaNeg.FillTrkPIDInfo : false
-
-physics.analyzers.TrkAnaPos.FillTrkQualInfo : false
-physics.analyzers.TrkAnaPos.FillTrkPIDInfo : false
+physics.analyzers.TrkAna.diagLevel : 1
+physics.analyzers.TrkAna.FillMCInfo : true
+physics.analyzers.TrkAna.FillCRVHits : true
+physics.analyzers.TrkAna.FillTriggerInfo : false
+physics.analyzers.TrkAna.FitType : LoopHelix
+
+physics.analyzers.TrkAna.FillTrkQualInfo : false
+physics.analyzers.TrkAna.FillTrkPIDInfo : false
 
 
 services.TFileService.fileName: "nts.owner.trkana-reco.version.sequencer.root"

fcl/TrkAnaRecoEnsemble-Data.fcl

@@ -1,7 +1,6 @@
 #include "TrkAna/fcl/TrkAnaReco.fcl"
 
-physics.analyzers.TrkAnaNeg.FillMCInfo : false
-physics.analyzers.TrkAnaPos.FillMCInfo : false
+physics.analyzers.TrkAna.FillMCInfo : false
 
 physics.TrkAnaEndPath : [ @sequence::TrkAnaReco.EndSequenceNoMC ]
 

fcl/TrkAnaReco_ceSimReco.fcl

@@ -2,11 +2,10 @@
 
 process_name : TrkAnaCeSimReco
 
-physics.analyzers.TrkAnaNeg.supplements : [  ]
-physics.analyzers.TrkAnaPos.supplements : [  ]
+physics.analyzers.TrkAna.supplements : [  ]
 
-physics.analyzers.TrkAnaNeg.FillTriggerInfo : false
+physics.analyzers.TrkAna.FillTriggerInfo : false
 physics.TrkAnaTrigPath : [ PBIWeight, TrkQualDeM, TrkPIDDeM ]
-physics.TrkAnaEndPath : [ TrkAnaNeg, genCountLogger ]
+physics.TrkAnaEndPath : [ TrkAna, genCountLogger ]
 
 services.TFileService.fileName: "nts.owner.trkana-ce-sim-reco.version.sequencer.root"

fcl/prolog.fcl

@@ -111,6 +111,7 @@ AllOpt : {
 DeM : { input : "KK"
   branch : "dem"
   suffix : "DeM"
+  options : { fillMC : true   genealogyDepth : -1 }
 }
 UeM : { input : "KK"
   branch : "uem"
@@ -127,6 +128,7 @@ UmuM : { input : "KK"
 DeP : { input : "KK"
   branch : "dep"
   suffix : "DeP"
+  options : { fillMC : true   genealogyDepth : -1 }
 }
 UeP : { input : "KK"
   branch : "uep"
@@ -148,6 +150,8 @@ Ext : { input : "KK"
 
 TrkAnaTreeMaker : {
   module_type : TrkAnaTreeMaker
+  branches : [ @local::DeM, @local::DeP, @local::UeM, @local::UeP,
+    @local::DmuM, @local::DmuP, @local::UmuM, @local::UmuP ]
   RecoCountTag : "SelectRecoMC"
   PBITag : "PBISim"
   PBTTag : "EWMProducer"
@@ -171,6 +175,8 @@ TrkAnaTreeMaker : {
   ProcessEmptyEvents : false
   FillCRVHits : true
   FillCRVPulses : false
+  FillCaloMC : true
+  FillHelixInfo : false
   PrimaryParticleTag : "compressRecoMCs"
   KalSeedMCAssns : "SelectRecoMC"
   CaloClusterMCTag : "compressRecoMCs"
@@ -179,6 +185,7 @@ TrkAnaTreeMaker : {
     MinGoodMomFraction : 0.9
     MaxVDDt : 5 # ns
   }
+  ExtraMCStepCollectionTags : [ "compressRecoMCs:protonabsorber", "compressRecoMCs:stoppingtarget" ]
 }
 
 TrkAnaReco : {
@@ -190,28 +197,18 @@ TrkAnaReco : {
   }
 
   analyzers : {
-    TrkAnaNeg : { @table::TrkAnaTreeMaker
-      candidate : @local::DeM
-      supplements : [ @local::UeM, @local::DmuM ]
-      ExtraMCStepCollectionTags : [ "compressRecoMCs:protonabsorber", "compressRecoMCs:stoppingtarget" ]
-    }
-    TrkAnaPos : { @table::TrkAnaTreeMaker
-      candidate : @local::DeP
-      supplements : [ @local::UeP, @local::DmuP ]
-      ExtraMCStepCollectionTags : [ "compressRecoMCs:protonabsorber", "compressRecoMCs:stoppingtarget" ]
-    }
+    TrkAna : { @table::TrkAnaTreeMaker  }
+
     TrkAnaExt : { @table::TrkAnaTreeMaker
-      candidate : @local::Ext
-      supplements : [ ]
-      ExtraMCStepCollectionTags : [ "compressRecoMCs:protonabsorber", "compressRecoMCs:stoppingtarget" ]
+      branches : [ @local::Ext ]
     }
 
     genCountLogger : @local::genCountLogger
   }
 
   TrigSequence : [ PBIWeight, @sequence::TrkQualProducersPath, @sequence::TrkPIDProducersPath ]
-  EndSequenceNoMC : [ TrkAnaNeg, TrkAnaPos ]
-  EndSequence : [ TrkAnaNeg, TrkAnaPos, genCountLogger ]
+  EndSequenceNoMC : [ TrkAna ]
+  EndSequence : [ TrkAna, genCountLogger ]
 
 }
 

inc/InfoMCStructHelper.hh

@@ -57,15 +57,15 @@ namespace mu2e {
       InfoMCStructHelper(const Config& conf);
       void updateEvent(const art::Event& event);
 
-      void fillTrkInfoMC(const KalSeed& kseed, const KalSeedMC& kseedmc, TrkInfoMC& trkinfomc);
+      void fillTrkInfoMC(const KalSeed& kseed, const KalSeedMC& kseedmc, std::vector<TrkInfoMC>& all_trkinfomcs);
       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<SimInfo>& siminfos, int n_generations, int n_match);
-      void fillVDInfo(KalSeed const& kseed, const KalSeedMC& kseedmc, std::vector<MCStepInfo>& vdinfos);
-      void fillHitInfoMCs(const KalSeedMC& kseedmc, std::vector<TrkStrawHitInfoMC>& tshinfomcs);
-      void fillCaloClusterInfoMC(CaloClusterMC const& ccmc, CaloClusterInfoMC& ccimc);
+      void fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector<std::vector<SimInfo>>& all_siminfos, int n_generations, int n_match);
+      void fillVDInfo(KalSeed const& kseed, const KalSeedMC& kseedmc, std::vector<std::vector<MCStepInfo>>& all_vdinfos);
+      void fillHitInfoMCs(const KalSeedMC& kseedmc, std::vector<std::vector<TrkStrawHitInfoMC>>& all_tshinfomcs);
+      void fillCaloClusterInfoMC(CaloClusterMC const& ccmc, std::vector<CaloClusterInfoMC>& ccimc);
       void fillExtraMCStepInfos(KalSeedMC const& kseedmc, StepPointMCCollection const& mcsteps,
-      MCStepInfos& mcsic, MCStepSummaryInfo& mcssi);
+                                std::vector<MCStepInfos>& mcsics, std::vector<MCStepSummaryInfo>& mcssis);
   };
 }
 

inc/InfoStructHelper.hh

@@ -48,17 +48,17 @@ namespace mu2e {
     void fillHitCount(const StrawHitFlagCollection& flags, HitCount& hitcount);
     void fillHitCount(RecoCount const& nrec, HitCount& hitcount);
 
-    void fillTrkInfo(const KalSeed& kseed,TrkInfo& trkinfo);
-    void fillTrkFitInfo(const KalSeed& kseed,std::vector<TrkFitInfo>& trkfitinfos);
+    void fillTrkInfo(const KalSeed& kseed,std::vector<TrkInfo>& trkinfo);
+    void fillTrkFitInfo(const KalSeed& kseed,std::vector<std::vector<TrkFitInfo>>& trkfitinfos);
     void fillTrkInfoHits(const KalSeed& kseed,TrkInfo& trkinfo);
     void fillTrkInfoStraws(const KalSeed& kseed,TrkInfo& trkinfo);
-    void fillLoopHelixInfo(const KalSeed& kseed, std::vector<LoopHelixInfo>& lhis);
-    void fillCentralHelixInfo(const KalSeed& kseed, std::vector<CentralHelixInfo>& chis);
-    void fillKinematicLineInfo(const KalSeed& kseed, std::vector<KinematicLineInfo>& klis);
+    void fillLoopHelixInfo(const KalSeed& kseed, std::vector<std::vector<LoopHelixInfo>>& lhis);
+    void fillCentralHelixInfo(const KalSeed& kseed, std::vector<std::vector<CentralHelixInfo>>& chis);
+    void fillKinematicLineInfo(const KalSeed& kseed, std::vector<std::vector<KinematicLineInfo>>& klis);
 
-    void fillHitInfo(const KalSeed& kseed, std::vector<TrkStrawHitInfo>& tshinfos );
-    void fillMatInfo(const KalSeed& kseed, std::vector<TrkStrawMatInfo>& tminfos );
-    void fillCaloHitInfo(const KalSeed& kseed, TrkCaloHitInfo& tchinfo );
+    void fillHitInfo(const KalSeed& kseed, std::vector<std::vector<TrkStrawHitInfo>>& tshinfos );
+    void fillMatInfo(const KalSeed& kseed, std::vector<std::vector<TrkStrawMatInfo>>& tminfos );
+    void fillCaloHitInfo(const KalSeed& kseed, std::vector<TrkCaloHitInfo>& tchinfo );
     void fillTrkQualInfo(const TrkQual& tqual, TrkQualInfo& trkqualInfo);
     void fillTrkPIDInfo(const TrkCaloHitPID& tchp, const KalSeed& kseed, TrkPIDInfo& trkpidInfo);
     void fillHelixInfo(art::Ptr<HelixSeed> const& hptr, HelixInfo& hinfo);

src/InfoMCStructHelper.cc

@@ -47,9 +47,10 @@ namespace mu2e {
     _onSpill = (ewMarker.spillType() == EventWindowMarker::SpillType::onspill);
   }
 
-  void InfoMCStructHelper::fillTrkInfoMC(const KalSeed& kseed, const KalSeedMC& kseedmc, TrkInfoMC& trkinfomc) {
+  void InfoMCStructHelper::fillTrkInfoMC(const KalSeed& kseed, const KalSeedMC& kseedmc, std::vector<TrkInfoMC>& all_trkinfomcs) {
     // use the primary match of the track
     // primary associated SimParticle
+    TrkInfoMC trkinfomc;
     if(kseedmc.simParticles().size() > 0){
       auto const& simp = kseedmc.simParticles().front();
       trkinfomc.valid = true;
@@ -58,6 +59,7 @@ namespace mu2e {
     }
 
     fillTrkInfoMCDigis(kseed, kseedmc, trkinfomc);
+    all_trkinfomcs.push_back(trkinfomc);
   }
 
   void InfoMCStructHelper::fillTrkInfoMCDigis(const KalSeed& kseed, const KalSeedMC& kseedmc, TrkInfoMC& trkinfomc) {
@@ -123,7 +125,9 @@ namespace mu2e {
     tshinfomc.doca = -1*dperp;
   }
 
-  void InfoMCStructHelper::fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector<SimInfo>& siminfos, int n_generations, int n_match) {
+  void InfoMCStructHelper::fillAllSimInfos(const KalSeedMC& kseedmc, const PrimaryParticle& primary, std::vector<std::vector<SimInfo>>& all_siminfos, int n_generations, int n_match) {
+    std::vector<SimInfo> siminfos;
+
     // interpret -1 as no llimit
     if (n_generations == -1) {
       n_generations = std::numeric_limits<int>::max();
@@ -192,6 +196,8 @@ namespace mu2e {
         siminfos.push_back(sim_info);
       }
     }
+
+    all_siminfos.push_back(siminfos);
   }
 
 
@@ -217,8 +223,8 @@ namespace mu2e {
     siminfo.endpos = XYZVectorF(det->toDetector(sp.endPosition()));
   }
 
-  void InfoMCStructHelper::fillVDInfo(const KalSeed& kseed, const KalSeedMC& kseedmc, std::vector<MCStepInfo>& vdinfos) {
-    vdinfos.clear();
+  void InfoMCStructHelper::fillVDInfo(const KalSeed& kseed, const KalSeedMC& kseedmc, std::vector<std::vector<MCStepInfo>>& all_vdinfos) {
+    std::vector<MCStepInfo> vdinfos;
     const auto& vdsteps = kseedmc._vdsteps;
     const auto& inters = kseed.intersections();
     double tmin = std::numeric_limits<float>::max();
@@ -268,19 +274,22 @@ namespace mu2e {
       vdinfos[imin].early = true;
       vdinfos[imax].late = true;
     }
+    all_vdinfos.push_back(vdinfos);
   }
 
-  void InfoMCStructHelper::fillHitInfoMCs(const KalSeedMC& kseedmc, std::vector<TrkStrawHitInfoMC>& tshinfomcs) {
-    tshinfomcs.clear();
+  void InfoMCStructHelper::fillHitInfoMCs(const KalSeedMC& kseedmc, std::vector<std::vector<TrkStrawHitInfoMC>>& all_tshinfomcs) {
+    std::vector<TrkStrawHitInfoMC> tshinfomcs;
 
     for(const auto& i_tshmc : kseedmc._tshmcs) {
       TrkStrawHitInfoMC tshinfomc;
       fillHitInfoMC(kseedmc, tshinfomc, i_tshmc);
       tshinfomcs.push_back(tshinfomc);
     }
+    all_tshinfomcs.push_back(tshinfomcs);
   }
 
-  void InfoMCStructHelper::fillCaloClusterInfoMC(CaloClusterMC const& ccmc, CaloClusterInfoMC& ccimc) {
+  void InfoMCStructHelper::fillCaloClusterInfoMC(CaloClusterMC const& ccmc, std::vector<CaloClusterInfoMC>& ccimcs) {
+    CaloClusterInfoMC ccimc;
     auto const& edeps = ccmc.energyDeposits();
     ccimc.nsim = edeps.size();
     ccimc.etot = ccmc.totalEnergyDep();
@@ -291,13 +300,14 @@ namespace mu2e {
       ccimc.tprimary = primary.time();
       ccimc.prel = primary.rel();
     }
+    ccimcs.push_back(ccimc);
   }
 
   void InfoMCStructHelper::fillExtraMCStepInfos(KalSeedMC const& kseedmc, StepPointMCCollection const& mcsteps,
-      MCStepInfos& mcsic, MCStepSummaryInfo& mcssi) {
+                                                std::vector<MCStepInfos>& mcsics, std::vector<MCStepSummaryInfo>& mcssis) {
+    MCStepInfos mcsic;
+    MCStepSummaryInfo mcssi;
     GeomHandle<DetectorSystem> det;
-    mcssi.reset();
-    mcsic.clear();
     MCStepInfo mcsi;
     // only count the extra steps associated with the primary MC truth match
     auto simp = kseedmc.simParticle().simParticle(_spcH);
@@ -338,5 +348,7 @@ namespace mu2e {
       mcsic.push_back(mcsi);
       mcssi.addStep(mcsi);
     }
+    mcsics.push_back(mcsic);
+    mcssis.push_back(mcssi);
   }
 }