All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Public Types | Public Member Functions | Private Attributes | List of all members
sbnd::CRTTrackMatchingAna Class Reference
Inheritance diagram for sbnd::CRTTrackMatchingAna:

Classes

struct  Config
 

Public Types

using Parameters = art::EDAnalyzer::Table< Config >
 

Public Member Functions

 CRTTrackMatchingAna (Parameters const &config)
 
virtual void beginJob () override
 
virtual void analyze (const art::Event &event) override
 
virtual void endJob () override
 
double DistToCrtHit (TVector3 trackPos, sbn::crt::CRTHit crtHit)
 

Private Attributes

art::InputTag fSimModuleLabel
 name of detsim producer More...
 
art::InputTag fCRTTrackLabel
 name of CRT producer More...
 
art::InputTag fTPCTrackLabel
 name of CRT producer More...
 
bool fVerbose
 print information about what's going on More...
 
CRTTrackMatchAlg trackAlg
 
CRTGeoAlg fCrtGeo
 
TPCGeoAlg fTpcGeo
 
CRTBackTracker fCrtBackTrack
 
CRTEventDisplay evd
 
TH1D * hAngle
 
TH1D * hMatchAngle
 
TH1D * hNoMatchAngle
 
TH1D * hDCA
 
TH1D * hMatchDCA
 
TH1D * hNoMatchDCA
 
TH2D * hMatchAngleDCA
 
TH2D * hNoMatchAngleDCA
 
TH1D * hEffAngleTotal
 
TH1D * hEffAngleReco
 
TH1D * hPurityAngleTotal
 
TH1D * hPurityAngleReco
 
TH1D * hEffDCATotal
 
TH1D * hEffDCAReco
 
TH1D * hPurityDCATotal
 
TH1D * hPurityDCAReco
 

Detailed Description

Definition at line 61 of file CRTTrackMatchingAna_module.cc.

Member Typedef Documentation

using sbnd::CRTTrackMatchingAna::Parameters = art::EDAnalyzer::Table<Config>

Definition at line 104 of file CRTTrackMatchingAna_module.cc.

Constructor & Destructor Documentation

sbnd::CRTTrackMatchingAna::CRTTrackMatchingAna ( Parameters const &  config)
explicit

Definition at line 163 of file CRTTrackMatchingAna_module.cc.

164  : EDAnalyzer(config)
165  , fSimModuleLabel (config().SimModuleLabel())
166  , fCRTTrackLabel (config().CRTTrackLabel())
167  , fTPCTrackLabel (config().TPCTrackLabel())
168  , fVerbose (config().Verbose())
169  , trackAlg (config().CRTTrackAlg())
170  , fCrtBackTrack (config().CrtBackTrack())
171  , evd (config().Evd())
172  {
173 
174  } //CRTTrackMatchingAna()
Verbose
BEGIN_PROLOG Verbose
Definition: simtestpulse_ICARUS.fcl:8
sbnd::CRTTrackMatchingAna::fSimModuleLabel
art::InputTag fSimModuleLabel
name of detsim producer
Definition: CRTTrackMatchingAna_module.cc:124
sbnd::CRTTrackMatchingAna::fCRTTrackLabel
art::InputTag fCRTTrackLabel
name of CRT producer
Definition: CRTTrackMatchingAna_module.cc:125
TPCTrackLabel
BEGIN_PROLOG don t mess with this TPCTrackLabel
Definition: crtt0matchingalg_icarus.fcl:6
sbnd::CRTTrackMatchingAna::fVerbose
bool fVerbose
print information about what&#39;s going on
Definition: CRTTrackMatchingAna_module.cc:127
sbnd::CRTTrackMatchingAna::fTPCTrackLabel
art::InputTag fTPCTrackLabel
name of CRT producer
Definition: CRTTrackMatchingAna_module.cc:126

Member Function Documentation

void sbnd::CRTTrackMatchingAna::analyze ( const art::Event &  event)
overridevirtual

Definition at line 210 of file CRTTrackMatchingAna_module.cc.

211  {
212 
213  // Fetch basic event info
214  if(fVerbose){
215  std::cout<<"============================================"<<std::endl
216  <<"Run = "<<event.run()<<", SubRun = "<<event.subRun()<<", Event = "<<event.id().event()<<std::endl
217  <<"============================================"<<std::endl;
218  }
219 
220  //----------------------------------------------------------------------------------------------------------
221  // GETTING PRODUCTS
222  //----------------------------------------------------------------------------------------------------------
223 
224  // Get g4 particles
225  auto particleHandle = event.getValidHandle<std::vector<simb::MCParticle>>(fSimModuleLabel);
226 
227  // Get CRT hits from the event
228  art::Handle< std::vector<sbn::crt::CRTTrack>> crtTrackHandle;
229  std::vector<art::Ptr<sbn::crt::CRTTrack> > crtTrackList;
230  if (event.getByLabel(fCRTTrackLabel, crtTrackHandle))
231  art::fill_ptr_vector(crtTrackList, crtTrackHandle);
232 
233  // Get reconstructed tracks from the event
234  auto tpcTrackHandle = event.getValidHandle<std::vector<recob::Track>>(fTPCTrackLabel);
235  art::FindManyP<recob::Hit> findManyHits(tpcTrackHandle, event, fTPCTrackLabel);
236 
237  fCrtBackTrack.Initialize(event);
238 
239  //----------------------------------------------------------------------------------------------------------
240  // TRUTH MATCHING
241  //----------------------------------------------------------------------------------------------------------
242 
243  std::map<int, simb::MCParticle> particles;
244  // Loop over the true particles
245  for (auto const& particle: (*particleHandle)){
246 
247  // Make map with ID
248  int partID = particle.TrackId();
249  particles[partID] = particle;
250 
251  }
252 
253  std::vector<sbn::crt::CRTTrack> crtTracks;
254  std::map<int, std::vector<sbn::crt::CRTTrack>> crtTrackMap;
255  int track_i = 0;
256  double minTrackTime = 99999;
257  double maxTrackTime = -99999;
258 
259  for(auto const& track : (*crtTrackHandle)){
260  crtTracks.push_back(track);
261  int trueID = fCrtBackTrack.TrueIdFromTrackId(event, track_i);
262  track_i++;
263  if(trueID != -99999) crtTrackMap[trueID].push_back(track);
264 
265  double trackTime = (double)(int)track.ts1_ns * 1e-3;
266  if(trackTime < minTrackTime) minTrackTime = trackTime;
267  if(trackTime > maxTrackTime) maxTrackTime = trackTime;
268  }
269 
270  //----------------------------------------------------------------------------------------------------------
271  // ANGLE AND DISTANCE OF CLOSEST APPROACH ANALYSIS
272  //----------------------------------------------------------------------------------------------------------
273 
274  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(event);
275  auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(event, clockData);
276 
277  // Loop over reconstructed tracks
278  for (auto const& tpcTrack : (*tpcTrackHandle)){
279  // Get the associated hits
280  std::vector<art::Ptr<recob::Hit>> hits = findManyHits.at(tpcTrack.ID());
281  int trackTrueID = RecoUtils::TrueParticleIDFromTotalRecoHits(clockData, hits, false);
282 
283  if(particles.find(trackTrueID) == particles.end()) continue;
284  // Only consider primary muons
285  if(!(std::abs(particles[trackTrueID].PdgCode()) == 13 && particles[trackTrueID].Mother() == 0)) continue;
286 
287  // Only consider particles withing time window of reco track
288  double trueTime = particles[trackTrueID].T() * 1e-3;
289  if(trueTime < minTrackTime || trueTime > maxTrackTime) continue;
290 
291  //----------------------------------------------------------------------------------------------------------
292  // SINGLE ANGLE CUT ANALYSIS
293  //----------------------------------------------------------------------------------------------------------
294  // Find the closest track by angle
295  std::pair<sbn::crt::CRTTrack, double> closestAngle = trackAlg.ClosestCRTTrackByAngle(detProp,tpcTrack, crtTracks, event);
296  if(closestAngle.second != -99999){
297  hAngle->Fill(closestAngle.second);
298  }
299  // Is crt track matched to that tpc track
300  if(closestAngle.second != -99999){
301  int crtTrueID = fCrtBackTrack.TrueIdFromTotalEnergy(event, closestAngle.first);
302  if(crtTrueID == trackTrueID && crtTrueID != -99999){
303  hMatchAngle->Fill(closestAngle.second);
304  }
305  else{
306  hNoMatchAngle->Fill(closestAngle.second);
307  }
308  }
309 
310  int nbinsAngle = hEffAngleTotal->GetNbinsX();
311  for(int i = 0; i < nbinsAngle; i++){
312  double angleCut = hEffAngleTotal->GetBinCenter(i);
313 
314  // Fill total efficiency histogram with each cut if track matches any hits
315  if(crtTrackMap.find(trackTrueID) != crtTrackMap.end()){
316  hEffAngleTotal->Fill(angleCut);
317 
318  // If closest hit is below limit and track matches any hits then fill efficiency
319  if(closestAngle.second < angleCut && closestAngle.second != -99999){
320  hEffAngleReco->Fill(angleCut);
321  }
322  }
323 
324  // Fill total purity histogram with each cut if closest hit is below limit
325  if(closestAngle.second < angleCut && closestAngle.second != -99999){
326  hPurityAngleTotal->Fill(angleCut);
327 
328  // If closest hit is below limit and matched time is correct then fill purity
329  double trackTime = closestAngle.first.ts1_ns * 1e-3;
330  if(particles.find(trackTrueID) != particles.end()){
331  if(std::abs(trackTime - trueTime) < 2.){
332  hPurityAngleReco->Fill(angleCut);
333  }
334  }
335  }
336  }
337 
338  //----------------------------------------------------------------------------------------------------------
339  // SINGLE DCA CUT ANALYSIS
340  //----------------------------------------------------------------------------------------------------------
341  // Find the closest track by average distance of closest approach
342  std::pair<sbn::crt::CRTTrack, double> closestDCA = trackAlg.ClosestCRTTrackByDCA(detProp,tpcTrack, crtTracks, event);
343  if(closestDCA.second != -99999){
344  hDCA->Fill(closestDCA.second);
345  }
346  if(closestDCA.second != -99999){
347  int crtTrueID = fCrtBackTrack.TrueIdFromTotalEnergy(event, closestDCA.first);
348  if(crtTrueID == trackTrueID && crtTrueID != -99999){
349  hMatchDCA->Fill(closestDCA.second);
350  }
351  else{
352  hNoMatchDCA->Fill(closestDCA.second);
353  }
354  }
355 
356  int nbinsDCA = hEffDCATotal->GetNbinsX();
357  for(int i = 0; i < nbinsDCA; i++){
358  double DCAcut = hEffDCATotal->GetBinCenter(i);
359 
360  // Fill total efficiency histogram with each cut if track matches any hits
361  if(crtTrackMap.find(trackTrueID) != crtTrackMap.end()){
362  hEffDCATotal->Fill(DCAcut);
363 
364  // If closest hit is below limit and track matches any hits then fill efficiency
365  if(closestDCA.second < DCAcut && closestDCA.second != -99999){
366  hEffDCAReco->Fill(DCAcut);
367  }
368  }
369 
370  // Fill total purity histogram with each cut if closest hit is below limit
371  if(closestDCA.second < DCAcut && closestDCA.second != -99999){
372  hPurityDCATotal->Fill(DCAcut);
373 
374  // If closest hit is below limit and matched time is correct then fill purity
375  double trackTime = closestDCA.first.ts1_ns * 1e-3;
376  if(particles.find(trackTrueID) != particles.end()){
377  if(std::abs(trackTime - trueTime) < 2.){
378  hPurityDCAReco->Fill(DCAcut);
379  }
380  }
381  }
382  }
383 
384 
385  //----------------------------------------------------------------------------------------------------------
386  // JOINT DCA AND ANGLE CUT ANALYSIS
387  //----------------------------------------------------------------------------------------------------------
388  std::vector<sbn::crt::CRTTrack> possTracks = trackAlg.AllPossibleCRTTracks(detProp,tpcTrack, crtTracks, event);
389  for(auto const& possTrack : possTracks){
390  int crtTrueID = fCrtBackTrack.TrueIdFromTotalEnergy(event, possTrack);
391  double angle = trackAlg.AngleBetweenTracks(tpcTrack, possTrack);
392  double DCA = trackAlg.AveDCABetweenTracks(detProp, tpcTrack, possTrack, event);
393  if(crtTrueID == trackTrueID && crtTrueID != -99999){
394  hMatchAngleDCA->Fill(angle, DCA);
395  }
396  else{
397  hNoMatchAngleDCA->Fill(angle, DCA);
398  }
399  }
400 
401  }
402 
403  /*
404  evd.SetDrawCrtTracks(true);
405  evd.SetDrawCrtData(true);
406  evd.SetDrawTpcTracks(true);
407  evd.SetDrawTrueTracks(true);
408  evd.SetDrawTpc(true);
409  evd.SetTrueId(114);
410  evd.Draw(event);
411 */
412 
413  } // CRTTrackMatchingAna::analyze()
sbnd::CRTBackTracker::TrueIdFromTrackId
int TrueIdFromTrackId(const art::Event &event, int track_i)
Definition: sbndcode/sbndcode/CRT/CRTUtils/CRTBackTracker.cc:448
sbnd::CRTTrackMatchingAna::fSimModuleLabel
art::InputTag fSimModuleLabel
name of detsim producer
Definition: CRTTrackMatchingAna_module.cc:124
track
process_name use argoneut_mc_hitfinder track
Definition: standard_reco.fcl:36
sbnd::CRTTrackMatchAlg::AngleBetweenTracks
double AngleBetweenTracks(recob::Track tpcTrack, sbn::crt::CRTTrack crtTrack)
Definition: CRTTrackMatchAlg.cc:347
sbnd::CRTTrackMatchAlg::AveDCABetweenTracks
double AveDCABetweenTracks(recob::Track tpcTrack, sbn::crt::CRTTrack crtTrack, double shift)
Definition: CRTTrackMatchAlg.cc:364
sbnd::CRTBackTracker::TrueIdFromTotalEnergy
int TrueIdFromTotalEnergy(const art::Event &event, const sbnd::crt::CRTData &data)
Definition: sbndcode/sbndcode/CRT/CRTUtils/CRTBackTracker.cc:281
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
sbnd::CRTTrackMatchAlg::ClosestCRTTrackByDCA
std::pair< sbn::crt::CRTTrack, double > ClosestCRTTrackByDCA(detinfo::DetectorPropertiesData const &detProp, recob::Track tpcTrack, std::vector< sbn::crt::CRTTrack > crtTracks, const art::Event &event, double minAngle=0.)
Definition: CRTTrackMatchAlg.cc:259
sbnd::CRTTrackMatchingAna::fCRTTrackLabel
art::InputTag fCRTTrackLabel
name of CRT producer
Definition: CRTTrackMatchingAna_module.cc:125
RecoUtils::TrueParticleIDFromTotalRecoHits
int TrueParticleIDFromTotalRecoHits(detinfo::DetectorClocksData const &clockData, const std::vector< art::Ptr< recob::Hit > > &hits, bool rollup_unsaved_ids=1)
Definition: icaruscode/icaruscode/CRT/CRTUtils/RecoUtils.cc:83
sbnd::CRTTrackMatchAlg::AllPossibleCRTTracks
std::vector< sbn::crt::CRTTrack > AllPossibleCRTTracks(detinfo::DetectorPropertiesData const &detProp, recob::Track tpcTrack, std::vector< sbn::crt::CRTTrack > crtTracks, const art::Event &event)
Definition: CRTTrackMatchAlg.cc:166
sbnd::CRTTrackMatchingAna::fVerbose
bool fVerbose
print information about what&#39;s going on
Definition: CRTTrackMatchingAna_module.cc:127
sbnd::CRTTrackMatchAlg::ClosestCRTTrackByAngle
std::pair< sbn::crt::CRTTrack, double > ClosestCRTTrackByAngle(detinfo::DetectorPropertiesData const &detProp, recob::Track tpcTrack, std::vector< sbn::crt::CRTTrack > crtTracks, const art::Event &event, double minDCA=0.)
Definition: CRTTrackMatchAlg.cc:217
e
do i e
Definition: updateVersionList.sh:68
sbnd::CRTTrackMatchingAna::fTPCTrackLabel
art::InputTag fTPCTrackLabel
name of CRT producer
Definition: CRTTrackMatchingAna_module.cc:126
angle
finds tracks best matching by angle
Definition: trackfinderalgorithms.fcl:64
cout
BEGIN_PROLOG could also be cout
Definition: messageservice.fcl:10
detsim::detProp
auto const detProp
Definition: SimWire_module.cc:157
void sbnd::CRTTrackMatchingAna::beginJob ( )
overridevirtual

Definition at line 177 of file CRTTrackMatchingAna_module.cc.

178  {
179 
180  // Access tfileservice to handle creating and writing histograms
181  art::ServiceHandle<art::TFileService> tfs;
182 
183  hAngle = tfs->make<TH1D>("Angle", "", 50, 0, 2);
184  hMatchAngle = tfs->make<TH1D>("MatchAngle", "", 50, 0, 2);
185  hNoMatchAngle = tfs->make<TH1D>("NoMatchAngle", "", 50, 0, 2);
186 
187  hDCA = tfs->make<TH1D>("DCA", "", 50, 0, 150);
188  hMatchDCA = tfs->make<TH1D>("MatchDCA", "", 50, 0, 150);
189  hNoMatchDCA = tfs->make<TH1D>("NoMatchDCA", "", 50, 0, 150);
190 
191  hMatchAngleDCA = tfs->make<TH2D>("MatchAngleDCA", "", 20, 0, 2, 20, 0, 150);
192  hNoMatchAngleDCA = tfs->make<TH2D>("NoMatchAngleDCA", "", 20, 0, 2, 20, 0, 150);
193 
194  hEffAngleTotal = tfs->make<TH1D>("EffAngleTotal", "", 20, 0, 1);
195  hEffAngleReco = tfs->make<TH1D>("EffAngleReco", "", 20, 0, 1);
196  hPurityAngleTotal = tfs->make<TH1D>("PurityAngleTotal", "", 20, 0, 1);
197  hPurityAngleReco = tfs->make<TH1D>("PurityAngleReco", "", 20, 0, 1);
198 
199  hEffDCATotal = tfs->make<TH1D>("EffDCATotal", "", 20, 0, 100);
200  hEffDCAReco = tfs->make<TH1D>("EffDCAReco", "", 20, 0, 100);
201  hPurityDCATotal = tfs->make<TH1D>("PurityDCATotal", "", 20, 0, 100);
202  hPurityDCAReco = tfs->make<TH1D>("PurityDCAReco", "", 20, 0, 100);
203 
204  // Initial output
205  if(fVerbose) std::cout<<"----------------- CRT T0 Matching Ana Module -------------------"<<std::endl;
206 
207  } // CRTTrackMatchingAna::beginJob()
sbnd::CRTTrackMatchingAna::fVerbose
bool fVerbose
print information about what&#39;s going on
Definition: CRTTrackMatchingAna_module.cc:127
tfs
art::ServiceHandle< art::TFileService > tfs
Definition: NuVertexChargeTree_module.cc:254
cout
BEGIN_PROLOG could also be cout
Definition: messageservice.fcl:10
double sbnd::CRTTrackMatchingAna::DistToCrtHit ( TVector3  trackPos,
sbn::crt::CRTHit  crtHit 
)
void sbnd::CRTTrackMatchingAna::endJob ( )
overridevirtual

Definition at line 416 of file CRTTrackMatchingAna_module.cc.

416  {
417 
418 
419  } // CRTTrackMatchingAna::endJob()

Member Data Documentation

CRTEventDisplay sbnd::CRTTrackMatchingAna::evd
private

Definition at line 135 of file CRTTrackMatchingAna_module.cc.

CRTBackTracker sbnd::CRTTrackMatchingAna::fCrtBackTrack
private

Definition at line 134 of file CRTTrackMatchingAna_module.cc.

CRTGeoAlg sbnd::CRTTrackMatchingAna::fCrtGeo
private

Definition at line 131 of file CRTTrackMatchingAna_module.cc.

art::InputTag sbnd::CRTTrackMatchingAna::fCRTTrackLabel
private

name of CRT producer

Definition at line 125 of file CRTTrackMatchingAna_module.cc.

art::InputTag sbnd::CRTTrackMatchingAna::fSimModuleLabel
private

name of detsim producer

Definition at line 124 of file CRTTrackMatchingAna_module.cc.

TPCGeoAlg sbnd::CRTTrackMatchingAna::fTpcGeo
private

Definition at line 132 of file CRTTrackMatchingAna_module.cc.

art::InputTag sbnd::CRTTrackMatchingAna::fTPCTrackLabel
private

name of CRT producer

Definition at line 126 of file CRTTrackMatchingAna_module.cc.

bool sbnd::CRTTrackMatchingAna::fVerbose
private

print information about what's going on

Definition at line 127 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hAngle
private

Definition at line 138 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hDCA
private

Definition at line 142 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hEffAngleReco
private

Definition at line 150 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hEffAngleTotal
private

Definition at line 149 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hEffDCAReco
private

Definition at line 155 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hEffDCATotal
private

Definition at line 154 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hMatchAngle
private

Definition at line 139 of file CRTTrackMatchingAna_module.cc.

TH2D* sbnd::CRTTrackMatchingAna::hMatchAngleDCA
private

Definition at line 146 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hMatchDCA
private

Definition at line 143 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hNoMatchAngle
private

Definition at line 140 of file CRTTrackMatchingAna_module.cc.

TH2D* sbnd::CRTTrackMatchingAna::hNoMatchAngleDCA
private

Definition at line 147 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hNoMatchDCA
private

Definition at line 144 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hPurityAngleReco
private

Definition at line 152 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hPurityAngleTotal
private

Definition at line 151 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hPurityDCAReco
private

Definition at line 157 of file CRTTrackMatchingAna_module.cc.

TH1D* sbnd::CRTTrackMatchingAna::hPurityDCATotal
private

Definition at line 156 of file CRTTrackMatchingAna_module.cc.

CRTTrackMatchAlg sbnd::CRTTrackMatchingAna::trackAlg
private

Definition at line 129 of file CRTTrackMatchingAna_module.cc.

The documentation for this class was generated from the following file:
Generated by   doxygen 1.8.5