MCParticleShowerMatching_module.cc
Go to the documentation of this file.
1 /////////////////////////////////////////////////////////////////////////////
2 /// Class: MCParticleShowerMatching
3 /// Module Type: producer
4 /// File: MCParticleShowerMatching_module.cc
5 ///
6 /// Author: Wesley Ketchum
7 /// E-mail address: wketchum@fnal.gov
8 ///
9 /// This module uses existing hit<-->MCParticle assns to make assns of showers
10 /// to mcparticles. It also stores some idea of cleanliness of the matching.
11 /// Note: it's probably better to do fancier things with the hit<->particle
12 /// info, but this is a start of sorts.
13 ///
14 /// Input: recob::Showers, recob::Hit collection, recob::Hit<--->simb::MCparticle assns
15 /// Output: recob::Shower/simb::MCParticle assns, with BackShowererMatchingData.
16 ///
17 /////////////////////////////////////////////////////////////////////////////
18 
19 // Framework includes
24 #include "canvas/Persistency/Common/FindManyP.h"
27 #include "fhiclcpp/ParameterSet.h"
28 
29 #include <memory>
30 #include <iostream>
31 #include <iterator>
32 
33 // LArSoft
38 
39 namespace t0 {
40  class MCParticleShowerMatching;
41 }
42 
44 public:
46  // The destructor generated by the compiler is fine for classes
47  // without bare pointers or other resource use.
48 
49  // Plugins should not be copied or assigned.
54 
55  // Required functions.
56  void produce(art::Event & e) override;
57 
58 private:
59 
64 
65 };
66 
67 
69  : EDProducer{p}
70 {
71  fShowerModuleLabel = p.get<art::InputTag>("ShowerModuleLabel");
72  fShowerHitAssnLabel = p.get<art::InputTag>("ShowerHitAssnLabel",fShowerModuleLabel);
73  fHitModuleLabel = p.get<art::InputTag>("HitModuleLabel");
74  fHitParticleAssnLabel = p.get<art::InputTag>("HitParticleAssnLabel");
75 
76  produces< art::Assns<recob::Shower , simb::MCParticle, anab::BackTrackerMatchingData > > ();
77 }
78 
80 {
81  if (evt.isRealData()) return;
82 
83  //auto mcpartHandle = evt.getValidHandle< std::vector<simb::MCParticle> >("largeant");
84  std::unique_ptr< art::Assns<recob::Shower, simb::MCParticle, anab::BackTrackerMatchingData > > MCPartShowerassn( new art::Assns<recob::Shower, simb::MCParticle, anab::BackTrackerMatchingData >);
85 
86 
87  double maxe = -1;
88  double tote = 0;
89  // int trkid = -1;
90  //int maxtrkid = -1;
91  //double maxn = -1;
92  //double totn = 0;
93  //int maxntrkid = -1;
94 
96  std::unordered_map<int,double> trkide;
97 
98  art::Handle< std::vector<recob::Shower> > showerListHandle;
99  evt.getByLabel(fShowerModuleLabel,showerListHandle);
100 
101  art::Handle< std::vector<recob::Hit> > hitListHandle;
102  evt.getByLabel(fHitModuleLabel,hitListHandle);
103 
104  if(!showerListHandle.isValid()){
105  std::cerr << "Shower handle is not valid!" << std::endl;
106  return;
107  }
108 
109  if(!hitListHandle.isValid()){
110  std::cerr << "Hit handle is not valid!" << std::endl;
111  return;
112  }
113 
114  auto const& showerList(*showerListHandle);
115  art::FindManyP<recob::Hit> fmtht(showerListHandle, evt, fShowerHitAssnLabel);
116  //auto const& mcpartList(*mcpartHandle);
117 
118  for(size_t i_t=0; i_t<showerList.size(); ++i_t){
119  art::Ptr<recob::Shower> shwPtr(showerListHandle,i_t);
120  trkide.clear();
121  tote = 0; maxe=-1; art::Ptr<simb::MCParticle> maxp;
122 
123  std::vector< art::Ptr<recob::Hit> > allHits = fmtht.at(i_t);
124 
125  std::vector<anab::BackTrackerHitMatchingData const*> bthmd_vec;
126  std::vector< art::Ptr<simb::MCParticle> > matchedParticlePtrs;
127 
128  art::FindManyP<simb::MCParticle,anab::BackTrackerHitMatchingData>
129  particles_per_hit(hitListHandle, evt, fHitParticleAssnLabel);
130 
131  for(size_t i_h=0; i_h<allHits.size(); ++i_h){
132  bthmd_vec.clear(); matchedParticlePtrs.clear();
133  particles_per_hit.get(allHits[i_h].key(),matchedParticlePtrs,bthmd_vec);
134 
135  for(size_t i_p=0; i_p<matchedParticlePtrs.size(); ++i_p){
136  trkide[ matchedParticlePtrs[i_p]->TrackId() ] += bthmd_vec[i_p]->energy;
137  tote += bthmd_vec[i_p]->energy;
138  if( trkide[ matchedParticlePtrs[i_p]->TrackId() ] > maxe ){
139  maxe = trkide[ matchedParticlePtrs[i_p]->TrackId() ];
140  maxp = matchedParticlePtrs[i_p];
141  }
142  }//end loop over particles per hit
143 
144  }//end loop over hits
145 
146  btdata.cleanliness = maxe/tote;
147  if(maxe>0)
148  MCPartShowerassn->addSingle(shwPtr, maxp, btdata);
149 
150  }//end loop over showers
151 
152  evt.put(std::move(MCPartShowerassn));
153 } // Produce
154 
code to link reconstructed objects back to the MC truth information
int maxp
Definition: tracks.py:196
MCParticleShowerMatching & operator=(MCParticleShowerMatching const &)=delete
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.h:20
Particle class.
bool isValid() const noexcept
Definition: Handle.h:191
bool isRealData() const
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Definition: DataViewImpl.h:633
const double e
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:67
def key(type, name=None)
Definition: graph.py:13
def move(depos, offset)
Definition: depos.py:107
p
Definition: test.py:223
ProductID put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
Definition: DataViewImpl.h:686
Declaration of signal hit object.
MCParticleShowerMatching(fhicl::ParameterSet const &p)
TCEvent evt
Definition: DataStructs.cxx:7
QTextStream & endl(QTextStream &s)