All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PFParticleMonitoring_module.cc
Go to the documentation of this file.
1 /**
2  * @file larpandora/LArPandoraAnalysis/PFParticleMonitoring_module.cc
3  *
4  * @brief Analysis module for created particles
5  *
6  */
7 
8 #include "art/Framework/Core/EDAnalyzer.h"
9 #include "art/Framework/Core/ModuleMacros.h"
10 
11 #include "TTree.h"
12 
13 #include "larpandora/LArPandoraInterface/LArPandoraHelper.h"
14 
15 #include <string>
16 
17 //------------------------------------------------------------------------------------------------------------------------------------------
18 
19 namespace lar_pandora {
20 
21  /**
22  * @brief PFParticleMonitoring class
23  */
24  class PFParticleMonitoring : public art::EDAnalyzer {
25  public:
26  /**
27  * @brief Constructor
28  *
29  * @param pset
30  */
31  PFParticleMonitoring(fhicl::ParameterSet const& pset);
32 
33  /**
34  * @brief Destructor
35  */
36  virtual ~PFParticleMonitoring();
37 
38  void beginJob();
39  void endJob();
40  void analyze(const art::Event& evt);
41  void reconfigure(fhicl::ParameterSet const& pset);
42 
43  private:
44  typedef std::set<art::Ptr<recob::PFParticle>> PFParticleSet;
45  typedef std::set<art::Ptr<simb::MCParticle>> MCParticleSet;
46  typedef std::set<art::Ptr<simb::MCTruth>> MCTruthSet;
47 
48  /**
49  * @brief Build mapping from true neutrinos to hits
50  *
51  * @param truthToParticles the input mapping from true event to true particles
52  * @param trueParticlesToHits the input mapping from true particles to hits
53  * @param trueNeutrinosToHits the output mapping from trues event to hits
54  * @param trueHitsToNeutrinos the output mappign from hits to true events
55  */
56  void BuildTrueNeutrinoHitMaps(const MCTruthToMCParticles& truthToParticles,
57  const MCParticlesToHits& trueParticlesToHits,
58  MCTruthToHits& trueNeutrinosToHits,
59  HitsToMCTruth& trueHitsToNeutrinos) const;
60 
61  /**
62  * @brief Build mapping from reconstructed neutrinos to hits
63  *
64  * @param recoParticleMap the input mapping from reconstructed particle and particle ID
65  * @param recoParticlesToHits the input mapping from reconstructed particles to hits
66  * @param recoNeutrinosToHits the output mapping from reconstructed particles to hits
67  * @param recoHitsToNeutrinos the output mapping from reconstructed hits to particles
68  */
69  void BuildRecoNeutrinoHitMaps(const PFParticleMap& recoParticleMap,
70  const PFParticlesToHits& recoParticlesToHits,
71  PFParticlesToHits& recoNeutrinosToHits,
72  HitsToPFParticles& recoHitsToNeutrinos) const;
73 
74  /**
75  * @brief Perform matching between true and reconstructed neutrino events
76  *
77  * @param recoNeutrinosToHits the mapping from reconstructed neutrino events to hits
78  * @param trueHitsToNeutrinos the mapping from hits to true neutrino events
79  * @param matchedNeutrinos the output matches between reconstructed and true neutrinos
80  * @param matchedNeutrinoHits the output matches between reconstructed neutrinos and hits
81  */
82  void GetRecoToTrueMatches(const PFParticlesToHits& recoNeutrinosToHits,
83  const HitsToMCTruth& trueHitsToNeutrinos,
84  MCTruthToPFParticles& matchedNeutrinos,
85  MCTruthToHits& matchedNeutrinoHits) const;
86 
87  /**
88  * @brief Perform matching between true and reconstructed neutrino events
89  *
90  * @param recoNeutrinosToHits the mapping from reconstructed neutrino events to hits
91  * @param trueHitsToNeutrinos the mapping from hits to true neutrino events
92  * @param matchedNeutrinos the output matches between reconstructed and true neutrinos
93  * @param matchedNeutrinoHits the output matches between reconstructed neutrinos and hits
94  * @param recoVeto the veto list for reconstructed particles
95  * @param trueVeto the veto list for true particles
96  */
97  void GetRecoToTrueMatches(const PFParticlesToHits& recoNeutrinosToHits,
98  const HitsToMCTruth& trueHitsToNeutrinos,
99  MCTruthToPFParticles& matchedNeutrinos,
100  MCTruthToHits& matchedNeutrinoHits,
101  PFParticleSet& recoVeto,
102  MCTruthSet& trueVeto) const;
103 
104  /**
105  * @brief Perform matching between true and reconstructed particles
106  *
107  * @param recoParticlesToHits the mapping from reconstructed particles to hits
108  * @param trueHitsToParticles the mapping from hits to true particles
109  * @param matchedParticles the output matches between reconstructed and true particles
110  * @param matchedHits the output matches between reconstructed particles and hits
111  */
112  void GetRecoToTrueMatches(const PFParticlesToHits& recoParticlesToHits,
113  const HitsToMCParticles& trueHitsToParticles,
114  MCParticlesToPFParticles& matchedParticles,
115  MCParticlesToHits& matchedHits) const;
116 
117  /**
118  * @brief Perform matching between true and reconstructed particles
119  *
120  * @param recoParticlesToHits the mapping from reconstructed particles to hits
121  * @param trueHitsToParticles the mapping from hits to true particles
122  * @param matchedParticles the output matches between reconstructed and true particles
123  * @param matchedHits the output matches between reconstructed particles and hits
124  * @param recoVeto the veto list for reconstructed particles
125  * @param trueVeto the veto list for true particles
126  */
127  void GetRecoToTrueMatches(const PFParticlesToHits& recoParticlesToHits,
128  const HitsToMCParticles& trueHitsToParticles,
129  MCParticlesToPFParticles& matchedParticles,
130  MCParticlesToHits& matchedHits,
131  PFParticleSet& recoVeto,
132  MCParticleSet& trueVeto) const;
133 
134  /**
135  * @brief Count the number of reconstructed hits in a given wire plane
136  *
137  * @param view the wire plane ID
138  * @param hitVector the input vector of reconstructed hits
139  */
140  int CountHitsByType(const int view, const HitVector& hitVector) const;
141 
142  /**
143  * @brief Find the start and end points of the true particle in the active region of detector
144  *
145  * @param trueParticle the input true particle
146  * @param startT the true start point
147  * @param endT the true end point
148  */
149  void GetStartAndEndPoints(const art::Ptr<simb::MCParticle> trueParticle,
150  int& startT,
151  int& endT) const;
152 
153  /**
154  * @brief Find the length of the true particle trajectory through the active region of the detector
155  *
156  * @param trueParticle the input true particle
157  * @param startT the true start point
158  * @param endT the true end point
159  */
160  double GetLength(const art::Ptr<simb::MCParticle> trueParticle,
161  const int startT,
162  const int endT) const;
163 
164  TTree* m_pRecoTree; ///<
165 
166  int m_run; ///<
167  int m_event; ///<
168  int m_index; ///<
169 
170  int m_nMCParticles; ///<
171  int m_nNeutrinoPfos; ///<
172  int m_nPrimaryPfos; ///<
173  int m_nDaughterPfos; ///<
174 
175  int m_mcPdg; ///<
176  int m_mcNuPdg; ///<
177  int m_mcParentPdg; ///<
178  int m_mcPrimaryPdg; ///<
179  int m_mcIsNeutrino; ///<
180  int m_mcIsPrimary; ///<
181  int m_mcIsDecay; ///<
182  int m_mcIsCC; ///<
183 
184  int m_pfoPdg; ///<
185  int m_pfoNuPdg; ///<
186  int m_pfoParentPdg; ///<
187  int m_pfoPrimaryPdg; ///<
188  int m_pfoIsNeutrino; ///<
189  int m_pfoIsPrimary; ///<
190  int m_pfoIsStitched; ///<
191 
192  int m_pfoTrack; ///<
193  int m_pfoVertex; ///<
194  double m_pfoVtxX; ///<
195  double m_pfoVtxY; ///<
196  double m_pfoVtxZ; ///<
197  double m_pfoEndX; ///<
198  double m_pfoEndY; ///<
199  double m_pfoEndZ; ///<
200  double m_pfoDirX; ///<
201  double m_pfoDirY; ///<
202  double m_pfoDirZ; ///<
203  double m_pfoLength; ///<
204  double m_pfoStraightLength; ///<
205 
206  int m_mcVertex; ///<
207  double m_mcVtxX; ///<
208  double m_mcVtxY; ///<
209  double m_mcVtxZ; ///<
210  double m_mcEndX; ///<
211  double m_mcEndY; ///<
212  double m_mcEndZ; ///<
213  double m_mcDirX; ///<
214  double m_mcDirY; ///<
215  double m_mcDirZ; ///<
216  double m_mcEnergy; ///<
217  double m_mcLength; ///<
218  double m_mcStraightLength; ///<
219 
220  double m_completeness; ///<
221  double m_purity; ///<
222 
223  int m_nMCHits; ///<
224  int m_nPfoHits; ///<
225  int m_nMatchedHits; ///<
226 
227  int m_nMCHitsU; ///<
228  int m_nMCHitsV; ///<
229  int m_nMCHitsW; ///<
230 
231  int m_nPfoHitsU; ///<
232  int m_nPfoHitsV; ///<
233  int m_nPfoHitsW; ///<
234 
235  int m_nMatchedHitsU; ///<
236  int m_nMatchedHitsV; ///<
237  int m_nMatchedHitsW; ///<
238 
239  int
240  m_nTrueWithoutRecoHits; ///< True hits which don't belong to any reconstructed particle - "available"
241  int
242  m_nRecoWithoutTrueHits; ///< Reconstructed hits which don't belong to any true particle - "missing"
243 
244  double m_spacepointsMinX; ///<
245  double m_spacepointsMaxX; ///<
246 
247  std::string m_hitfinderLabel; ///<
248  std::string m_trackLabel; ///<
249  std::string m_particleLabel; ///<
250  std::string m_backtrackerLabel; ///<
251  std::string m_geantModuleLabel; ///<
252 
253  bool m_useDaughterPFParticles; ///<
254  bool m_useDaughterMCParticles; ///<
255  bool m_addDaughterPFParticles; ///<
256  bool m_addDaughterMCParticles; ///<
257 
258  bool m_recursiveMatching; ///<
259  bool m_printDebug; ///< switch for print statements (TODO: use message service!)
260  bool
261  m_disableRealDataCheck; ///< Whether to check if the input file contains real data before accessing MC information
262  };
263 
264  DEFINE_ART_MODULE(PFParticleMonitoring)
265 
266 } // namespace lar_pandora
267 
268 //------------------------------------------------------------------------------------------------------------------------------------------
269 // implementation follows
270 
271 #include "art/Framework/Principal/Event.h"
272 #include "art/Framework/Principal/Handle.h"
273 #include "art/Framework/Services/Registry/ServiceHandle.h"
274 #include "art_root_io/TFileDirectory.h"
275 #include "art_root_io/TFileService.h"
276 #include "fhiclcpp/ParameterSet.h"
277 #include "messagefacility/MessageLogger/MessageLogger.h"
278 
279 #include "larcore/Geometry/Geometry.h"
280 #include "larcorealg/Geometry/CryostatGeo.h"
281 #include "larcorealg/Geometry/PlaneGeo.h"
282 #include "larcorealg/Geometry/TPCGeo.h"
283 #include "lardata/DetectorInfoServices/DetectorPropertiesService.h"
284 #include "lardata/DetectorInfoServices/LArPropertiesService.h"
285 #include "lardata/Utilities/AssociationUtil.h"
286 #include "lardataobj/AnalysisBase/T0.h"
287 #include "lardataobj/RecoBase/Cluster.h"
288 #include "lardataobj/RecoBase/Hit.h"
289 #include "lardataobj/RecoBase/PFParticle.h"
290 #include "lardataobj/RecoBase/SpacePoint.h"
291 #include "lardataobj/RecoBase/Track.h"
292 #include "lardataobj/RecoBase/Vertex.h"
293 #include "nusimdata/SimulationBase/MCParticle.h"
294 #include "nusimdata/SimulationBase/MCTruth.h"
295 
296 #include <iostream>
297 
298 namespace lar_pandora {
299 
300  PFParticleMonitoring::PFParticleMonitoring(fhicl::ParameterSet const& pset)
301  : art::EDAnalyzer(pset)
302  {
303  this->reconfigure(pset);
304  }
305 
306  //------------------------------------------------------------------------------------------------------------------------------------------
307 
308  PFParticleMonitoring::~PFParticleMonitoring() {}
309 
310  //------------------------------------------------------------------------------------------------------------------------------------------
311 
312  void
313  PFParticleMonitoring::reconfigure(fhicl::ParameterSet const& pset)
314  {
315  m_trackLabel = pset.get<std::string>("TrackModule", "pandoraTracks");
316  m_particleLabel = pset.get<std::string>("PFParticleModule", "pandora");
317  m_hitfinderLabel = pset.get<std::string>("HitFinderModule", "gaushit");
318  m_backtrackerLabel = pset.get<std::string>("BackTrackerModule", "gaushitTruthMatch");
319  m_geantModuleLabel = pset.get<std::string>("GeantModule", "largeant");
320 
321  m_useDaughterPFParticles = pset.get<bool>("UseDaughterPFParticles", false);
322  m_useDaughterMCParticles = pset.get<bool>("UseDaughterMCParticles", true);
323  m_addDaughterPFParticles = pset.get<bool>("AddDaughterPFParticles", true);
324  m_addDaughterMCParticles = pset.get<bool>("AddDaughterMCParticles", true);
325 
326  m_recursiveMatching = pset.get<bool>("RecursiveMatching", false);
327  m_printDebug = pset.get<bool>("PrintDebug", false);
328  m_disableRealDataCheck = pset.get<bool>("DisableRealDataCheck", false);
329  }
330 
331  //------------------------------------------------------------------------------------------------------------------------------------------
332 
333  void
334  PFParticleMonitoring::beginJob()
335  {
336  mf::LogDebug("LArPandora") << " *** PFParticleMonitoring::beginJob() *** " << std::endl;
337 
338  //
339  art::ServiceHandle<art::TFileService const> tfs;
340 
341  m_pRecoTree = tfs->make<TTree>("pandora", "LAr Reco vs True");
342  m_pRecoTree->Branch("run", &m_run, "run/I");
343  m_pRecoTree->Branch("event", &m_event, "event/I");
344  m_pRecoTree->Branch("index", &m_index, "index/I");
345  m_pRecoTree->Branch("nMCParticles", &m_nMCParticles, "nMCParticles/I");
346  m_pRecoTree->Branch("nNeutrinoPfos", &m_nNeutrinoPfos, "nNeutrinoPfos/I");
347  m_pRecoTree->Branch("nPrimaryPfos", &m_nPrimaryPfos, "nPrimaryPfos/I");
348  m_pRecoTree->Branch("nDaughterPfos", &m_nDaughterPfos, "nDaughterPfos/I");
349  m_pRecoTree->Branch("mcPdg", &m_mcPdg, "mcPdg/I");
350  m_pRecoTree->Branch("mcNuPdg", &m_mcNuPdg, "mcNuPdg/I");
351  m_pRecoTree->Branch("mcParentPdg", &m_mcParentPdg, "mcParentPdg/I");
352  m_pRecoTree->Branch("mcPrimaryPdg", &m_mcPrimaryPdg, "mcPrimaryPdg/I");
353  m_pRecoTree->Branch("mcIsNeutrino", &m_mcIsNeutrino, "mcIsNeutrino/I");
354  m_pRecoTree->Branch("mcIsPrimary", &m_mcIsPrimary, "mcIsPrimary/I");
355  m_pRecoTree->Branch("mcIsDecay", &m_mcIsDecay, "mcIsDecay/I");
356  m_pRecoTree->Branch("mcIsCC", &m_mcIsCC, "mcIsCC/I");
357  m_pRecoTree->Branch("pfoPdg", &m_pfoPdg, "pfoPdg/I");
358  m_pRecoTree->Branch("pfoNuPdg", &m_pfoNuPdg, "pfoNuPdg/I");
359  m_pRecoTree->Branch("pfoParentPdg", &m_pfoParentPdg, "pfoParentPdg/I");
360  m_pRecoTree->Branch("pfoPrimaryPdg", &m_pfoPrimaryPdg, "pfoPrimaryPdg/I");
361  m_pRecoTree->Branch("pfoIsNeutrino", &m_pfoIsNeutrino, "pfoIsNeutrino/I");
362  m_pRecoTree->Branch("pfoIsPrimary", &m_pfoIsPrimary, "pfoIsPrimary/I");
363  m_pRecoTree->Branch("pfoIsStitched", &m_pfoIsStitched, "pfoIsStitched/I");
364  m_pRecoTree->Branch("pfoTrack", &m_pfoTrack, "pfoTrack/I");
365  m_pRecoTree->Branch("pfoVertex", &m_pfoVertex, "pfoVertex/I");
366  m_pRecoTree->Branch("pfoVtxX", &m_pfoVtxX, "pfoVtxX/D");
367  m_pRecoTree->Branch("pfoVtxY", &m_pfoVtxY, "pfoVtxY/D");
368  m_pRecoTree->Branch("pfoVtxZ", &m_pfoVtxZ, "pfoVtxZ/D");
369  m_pRecoTree->Branch("pfoEndX", &m_pfoEndX, "pfoEndX/D");
370  m_pRecoTree->Branch("pfoEndY", &m_pfoEndY, "pfoEndY/D");
371  m_pRecoTree->Branch("pfoEndZ", &m_pfoEndZ, "pfoEndZ/D");
372  m_pRecoTree->Branch("pfoDirX", &m_pfoDirX, "pfoDirX/D");
373  m_pRecoTree->Branch("pfoDirY", &m_pfoDirY, "pfoDirY/D");
374  m_pRecoTree->Branch("pfoDirZ", &m_pfoDirZ, "pfoDirZ/D");
375  m_pRecoTree->Branch("pfoLength", &m_pfoLength, "pfoLength/D");
376  m_pRecoTree->Branch("pfoStraightLength", &m_pfoStraightLength, "pfoStraightLength/D");
377  m_pRecoTree->Branch("mcVertex", &m_mcVertex, "mcVertex/I");
378  m_pRecoTree->Branch("mcVtxX", &m_mcVtxX, "mcVtxX/D");
379  m_pRecoTree->Branch("mcVtxY", &m_mcVtxY, "mcVtxY/D");
380  m_pRecoTree->Branch("mcVtxZ", &m_mcVtxZ, "mcVtxZ/D");
381  m_pRecoTree->Branch("mcEndX", &m_mcEndX, "mcEndX/D");
382  m_pRecoTree->Branch("mcEndY", &m_mcEndY, "mcEndY/D");
383  m_pRecoTree->Branch("mcEndZ", &m_mcEndZ, "mcEndZ/D");
384  m_pRecoTree->Branch("mcDirX", &m_mcDirX, "mcDirX/D");
385  m_pRecoTree->Branch("mcDirY", &m_mcDirY, "mcDirY/D");
386  m_pRecoTree->Branch("mcDirZ", &m_mcDirZ, "mcDirZ/D");
387  m_pRecoTree->Branch("mcEnergy", &m_mcEnergy, "mcEnergy/D");
388  m_pRecoTree->Branch("mcLength", &m_mcLength, "mcLength/D");
389  m_pRecoTree->Branch("mcStraightLength", &m_mcStraightLength, "mcStraightLength/D");
390  m_pRecoTree->Branch("completeness", &m_completeness, "completeness/D");
391  m_pRecoTree->Branch("purity", &m_purity, "purity/D");
392  m_pRecoTree->Branch("nMCHits", &m_nMCHits, "nMCHits/I");
393  m_pRecoTree->Branch("nPfoHits", &m_nPfoHits, "nPfoHits/I");
394  m_pRecoTree->Branch("nMatchedHits", &m_nMatchedHits, "nMatchedHits/I");
395  m_pRecoTree->Branch("nMCHitsU", &m_nMCHitsU, "nMCHitsU/I");
396  m_pRecoTree->Branch("nMCHitsV", &m_nMCHitsV, "nMCHitsV/I");
397  m_pRecoTree->Branch("nMCHitsW", &m_nMCHitsW, "nMCHitsW/I");
398  m_pRecoTree->Branch("nPfoHitsU", &m_nPfoHitsU, "nPfoHitsU/I");
399  m_pRecoTree->Branch("nPfoHitsV", &m_nPfoHitsV, "nPfoHitsV/I");
400  m_pRecoTree->Branch("nPfoHitsW", &m_nPfoHitsW, "nPfoHitsW/I");
401  m_pRecoTree->Branch("nMatchedHitsU", &m_nMatchedHitsU, "nMatchedHitsU/I");
402  m_pRecoTree->Branch("nMatchedHitsV", &m_nMatchedHitsV, "nMatchedHitsV/I");
403  m_pRecoTree->Branch("nMatchedHitsW", &m_nMatchedHitsW, "nMatchedHitsW/I");
404  m_pRecoTree->Branch("nTrueWithoutRecoHits", &m_nTrueWithoutRecoHits, "nTrueWithoutRecoHits/I");
405  m_pRecoTree->Branch("nRecoWithoutTrueHits", &m_nRecoWithoutTrueHits, "nRecoWithoutTrueHits/I");
406  m_pRecoTree->Branch("spacepointsMinX", &m_spacepointsMinX, "spacepointsMinX/D");
407  m_pRecoTree->Branch("spacepointsMaxX", &m_spacepointsMaxX, "spacepointsMaxX/D");
408  }
409 
410  //------------------------------------------------------------------------------------------------------------------------------------------
411 
412  void
413  PFParticleMonitoring::endJob()
414  {}
415 
416  //------------------------------------------------------------------------------------------------------------------------------------------
417 
418  void
419  PFParticleMonitoring::analyze(const art::Event& evt)
420  {
421  if (m_printDebug) std::cout << " *** PFParticleMonitoring::analyze(...) *** " << std::endl;
422 
423  m_run = evt.run();
424  m_event = evt.id().event();
425  m_index = 0;
426 
427  m_nMCParticles = 0;
428  m_nNeutrinoPfos = 0;
429  m_nPrimaryPfos = 0;
430  m_nDaughterPfos = 0;
431 
432  m_mcPdg = 0;
433  m_mcNuPdg = 0;
434  m_mcParentPdg = 0;
435  m_mcPrimaryPdg = 0;
436  m_mcIsNeutrino = 0;
437  m_mcIsPrimary = 0;
438  m_mcIsDecay = 0;
439  m_mcIsCC = 0;
440 
441  m_pfoPdg = 0;
442  m_pfoNuPdg = 0;
443  m_pfoParentPdg = 0;
444  m_pfoPrimaryPdg = 0;
445  m_pfoIsNeutrino = 0;
446  m_pfoIsPrimary = 0;
447  m_pfoIsStitched = 0;
448  m_pfoTrack = 0;
449  m_pfoVertex = 0;
450  m_pfoVtxX = 0.0;
451  m_pfoVtxY = 0.0;
452  m_pfoVtxZ = 0.0;
453  m_pfoEndX = 0.0;
454  m_pfoEndY = 0.0;
455  m_pfoEndZ = 0.0;
456  m_pfoDirX = 0.0;
457  m_pfoDirY = 0.0;
458  m_pfoDirZ = 0.0;
459  m_pfoLength = 0.0;
460  m_pfoStraightLength = 0.0;
461 
462  m_mcVertex = 0;
463  m_mcVtxX = 0.0;
464  m_mcVtxY = 0.0;
465  m_mcVtxZ = 0.0;
466  m_mcEndX = 0.0;
467  m_mcEndY = 0.0;
468  m_mcEndZ = 0.0;
469  m_mcDirX = 0.0;
470  m_mcDirY = 0.0;
471  m_mcDirZ = 0.0;
472  m_mcEnergy = 0.0;
473  m_mcLength = 0.0;
474  m_mcStraightLength = 0.0;
475 
476  m_completeness = 0.0;
477  m_purity = 0.0;
478 
479  m_nMCHits = 0;
480  m_nPfoHits = 0;
481  m_nMatchedHits = 0;
482  m_nMCHitsU = 0;
483  m_nMCHitsV = 0;
484  m_nMCHitsW = 0;
485  m_nPfoHitsU = 0;
486  m_nPfoHitsV = 0;
487  m_nPfoHitsW = 0;
488  m_nMatchedHitsU = 0;
489  m_nMatchedHitsV = 0;
490  m_nMatchedHitsW = 0;
491 
492  m_nTrueWithoutRecoHits = 0;
493  m_nRecoWithoutTrueHits = 0;
494 
495  m_spacepointsMinX = 0.0;
496  m_spacepointsMaxX = 0.0;
497 
498  if (m_printDebug) {
499  std::cout << " Run: " << m_run << std::endl;
500  std::cout << " Event: " << m_event << std::endl;
501  }
502 
503  // Collect Hits
504  // ============
505  HitVector hitVector;
506  LArPandoraHelper::CollectHits(evt, m_hitfinderLabel, hitVector);
507 
508  if (m_printDebug) std::cout << " Hits: " << hitVector.size() << std::endl;
509 
510  // Collect SpacePoints and SpacePoint <-> Hit Associations
511  // =======================================================
512  SpacePointVector spacePointVector;
513  SpacePointsToHits spacePointsToHits;
514  HitsToSpacePoints hitsToSpacePoints;
515  LArPandoraHelper::CollectSpacePoints(
516  evt, m_particleLabel, spacePointVector, spacePointsToHits, hitsToSpacePoints);
517 
518  if (m_printDebug) std::cout << " SpacePoints: " << spacePointVector.size() << std::endl;
519 
520  // Collect Tracks and PFParticle <-> Track Associations
521  // ====================================================
522  TrackVector recoTrackVector;
523  PFParticlesToTracks recoParticlesToTracks;
524  LArPandoraHelper::CollectTracks(evt, m_trackLabel, recoTrackVector, recoParticlesToTracks);
525 
526  if (m_printDebug) std::cout << " Tracks: " << recoTrackVector.size() << std::endl;
527 
528  // Collect TOs and PFParticle <-> T0 Associations
529  // ==============================================
530  T0Vector t0Vector;
531  PFParticlesToT0s particlesToT0s;
532  LArPandoraHelper::CollectT0s(evt, m_particleLabel, t0Vector, particlesToT0s);
533 
534  // Collect Vertices and PFParticle <-> Vertex Associations
535  // =======================================================
536  VertexVector recoVertexVector;
537  PFParticlesToVertices recoParticlesToVertices;
538  LArPandoraHelper::CollectVertices(
539  evt, m_particleLabel, recoVertexVector, recoParticlesToVertices);
540 
541  if (m_printDebug) std::cout << " Vertices: " << recoVertexVector.size() << std::endl;
542 
543  // Collect PFParticles and match Reco Particles to Hits
544  // ====================================================
545  PFParticleVector recoParticleVector;
546  PFParticleVector recoNeutrinoVector;
547  PFParticlesToHits recoParticlesToHits;
548  HitsToPFParticles recoHitsToParticles;
549 
550  LArPandoraHelper::CollectPFParticles(evt, m_particleLabel, recoParticleVector);
551  LArPandoraHelper::SelectNeutrinoPFParticles(recoParticleVector, recoNeutrinoVector);
552  LArPandoraHelper::BuildPFParticleHitMaps(
553  evt,
554  m_particleLabel,
555  recoParticlesToHits,
556  recoHitsToParticles,
557  (m_useDaughterPFParticles ? (m_addDaughterPFParticles ? LArPandoraHelper::kAddDaughters :
558  LArPandoraHelper::kUseDaughters) :
559  LArPandoraHelper::kIgnoreDaughters));
560 
561  if (m_printDebug) {
562  std::cout << " RecoNeutrinos: " << recoNeutrinoVector.size() << std::endl;
563  std::cout << " RecoParticles: " << recoParticleVector.size() << std::endl;
564  }
565 
566  // Collect MCParticles and match True Particles to Hits
567  // ====================================================
568  MCParticleVector trueParticleVector;
569  MCTruthToMCParticles truthToParticles;
570  MCParticlesToMCTruth particlesToTruth;
571  MCParticlesToHits trueParticlesToHits;
572  HitsToMCParticles trueHitsToParticles;
573 
574  if (m_disableRealDataCheck || !evt.isRealData()) {
575  LArPandoraHelper::CollectMCParticles(evt, m_geantModuleLabel, trueParticleVector);
576  LArPandoraHelper::CollectMCParticles(
577  evt, m_geantModuleLabel, truthToParticles, particlesToTruth);
578 
579  LArPandoraHelper::BuildMCParticleHitMaps(
580  evt,
581  m_geantModuleLabel,
582  hitVector,
583  trueParticlesToHits,
584  trueHitsToParticles,
585  (m_useDaughterMCParticles ? (m_addDaughterMCParticles ? LArPandoraHelper::kAddDaughters :
586  LArPandoraHelper::kUseDaughters) :
587  LArPandoraHelper::kIgnoreDaughters));
588 
589  if (trueHitsToParticles.empty()) {
590  if (m_backtrackerLabel.empty())
591  throw cet::exception("LArPandora") << " PFParticleMonitoring::analyze - no sim channels "
592  "found, backtracker module must be set in FHiCL "
593  << std::endl;
594 
595  LArPandoraHelper::BuildMCParticleHitMaps(
596  evt,
597  m_geantModuleLabel,
598  m_hitfinderLabel,
599  m_backtrackerLabel,
600  trueParticlesToHits,
601  trueHitsToParticles,
602  (m_useDaughterMCParticles ? (m_addDaughterMCParticles ? LArPandoraHelper::kAddDaughters :
603  LArPandoraHelper::kUseDaughters) :
604  LArPandoraHelper::kIgnoreDaughters));
605  }
606  }
607 
608  if (m_printDebug) {
609  std::cout << " TrueParticles: " << particlesToTruth.size() << std::endl;
610  std::cout << " TrueEvents: " << truthToParticles.size() << std::endl;
611  std::cout << " MatchedParticles: " << trueParticlesToHits.size() << std::endl;
612  }
613 
614  if (trueParticlesToHits.empty()) {
615  m_pRecoTree->Fill();
616  return;
617  }
618 
619  // Build Reco and True Particle Maps (for Parent/Daughter Navigation)
620  // =================================================================
621  MCParticleMap trueParticleMap;
622  PFParticleMap recoParticleMap;
623 
624  LArPandoraHelper::BuildMCParticleMap(trueParticleVector, trueParticleMap);
625  LArPandoraHelper::BuildPFParticleMap(recoParticleVector, recoParticleMap);
626 
627  m_nMCParticles = trueParticlesToHits.size();
628  m_nNeutrinoPfos = 0;
629  m_nPrimaryPfos = 0;
630  m_nDaughterPfos = 0;
631 
632  // Count reconstructed particles
633  for (PFParticleVector::const_iterator iter = recoParticleVector.begin(),
634  iterEnd = recoParticleVector.end();
635  iter != iterEnd;
636  ++iter) {
637  const art::Ptr<recob::PFParticle> recoParticle = *iter;
638 
639  if (LArPandoraHelper::IsNeutrino(recoParticle)) { m_nNeutrinoPfos++; }
640  else if (LArPandoraHelper::IsFinalState(recoParticleMap, recoParticle)) {
641  m_nPrimaryPfos++;
642  }
643  else {
644  m_nDaughterPfos++;
645  }
646  }
647 
648  // Match Reco Neutrinos to True Neutrinos
649  // ======================================
650  PFParticlesToHits recoNeutrinosToHits;
651  HitsToPFParticles recoHitsToNeutrinos;
652  HitsToMCTruth trueHitsToNeutrinos;
653  MCTruthToHits trueNeutrinosToHits;
654  this->BuildRecoNeutrinoHitMaps(
655  recoParticleMap, recoParticlesToHits, recoNeutrinosToHits, recoHitsToNeutrinos);
656  this->BuildTrueNeutrinoHitMaps(
657  truthToParticles, trueParticlesToHits, trueNeutrinosToHits, trueHitsToNeutrinos);
658 
659  MCTruthToPFParticles matchedNeutrinos;
660  MCTruthToHits matchedNeutrinoHits;
661  this->GetRecoToTrueMatches(
662  recoNeutrinosToHits, trueHitsToNeutrinos, matchedNeutrinos, matchedNeutrinoHits);
663 
664  for (MCTruthToHits::const_iterator iter = trueNeutrinosToHits.begin(),
665  iterEnd = trueNeutrinosToHits.end();
666  iter != iterEnd;
667  ++iter) {
668  const art::Ptr<simb::MCTruth> trueEvent = iter->first;
669  const HitVector& trueHitVector = iter->second;
670 
671  if (trueHitVector.empty()) continue;
672 
673  if (!trueEvent->NeutrinoSet()) continue;
674 
675  const simb::MCNeutrino trueNeutrino(trueEvent->GetNeutrino());
676  const simb::MCParticle trueParticle(trueNeutrino.Nu());
677 
678  m_mcIsCC = ((simb::kCC == trueNeutrino.CCNC()) ? 1 : 0);
679  m_mcPdg = trueParticle.PdgCode();
680  m_mcNuPdg = m_mcPdg;
681  m_mcParentPdg = 0;
682  m_mcPrimaryPdg = 0;
683  m_mcIsNeutrino = 1;
684  m_mcIsPrimary = 0;
685  m_mcIsDecay = 0;
686 
687  m_mcVertex = 1;
688  m_mcVtxX = trueParticle.Vx();
689  m_mcVtxY = trueParticle.Vy();
690  m_mcVtxZ = trueParticle.Vz();
691  m_mcEndX = m_mcVtxX;
692  m_mcEndY = m_mcVtxY;
693  m_mcEndZ = m_mcVtxZ;
694  m_mcDirX = trueParticle.Px() / trueParticle.P();
695  m_mcDirY = trueParticle.Py() / trueParticle.P();
696  m_mcDirZ = trueParticle.Pz() / trueParticle.P();
697  m_mcEnergy = trueParticle.E();
698  m_mcLength = 0.0;
699  m_mcStraightLength = 0.0;
700 
701  m_nMCHits = trueHitVector.size();
702  m_nMCHitsU = this->CountHitsByType(geo::kU, trueHitVector);
703  m_nMCHitsV = this->CountHitsByType(geo::kV, trueHitVector);
704  m_nMCHitsW = this->CountHitsByType(geo::kW, trueHitVector);
705 
706  m_pfoPdg = 0;
707  m_pfoNuPdg = 0;
708  m_pfoParentPdg = 0;
709  m_pfoPrimaryPdg = 0;
710  m_pfoIsNeutrino = 0;
711  m_pfoIsPrimary = 0;
712  m_pfoIsStitched = 0;
713  m_pfoTrack = 0;
714  m_pfoVertex = 0;
715  m_pfoVtxX = 0.0;
716  m_pfoVtxY = 0.0;
717  m_pfoVtxZ = 0.0;
718  m_pfoEndX = 0.0;
719  m_pfoEndY = 0.0;
720  m_pfoEndZ = 0.0;
721  m_pfoDirX = 0.0;
722  m_pfoDirY = 0.0;
723  m_pfoDirZ = 0.0;
724  m_pfoLength = 0.0;
725  m_pfoStraightLength = 0.0;
726 
727  m_nPfoHits = 0;
728  m_nPfoHitsU = 0;
729  m_nPfoHitsV = 0;
730  m_nPfoHitsW = 0;
731 
732  m_nMatchedHits = 0;
733  m_nMatchedHitsU = 0;
734  m_nMatchedHitsV = 0;
735  m_nMatchedHitsW = 0;
736 
737  m_nTrueWithoutRecoHits = 0;
738  m_nRecoWithoutTrueHits = 0;
739 
740  m_spacepointsMinX = 0.0;
741  m_spacepointsMaxX = 0.0;
742 
743  m_completeness = 0.0;
744  m_purity = 0.0;
745 
746  for (HitVector::const_iterator hIter1 = trueHitVector.begin(),
747  hIterEnd1 = trueHitVector.end();
748  hIter1 != hIterEnd1;
749  ++hIter1) {
750  if (recoHitsToNeutrinos.find(*hIter1) == recoHitsToNeutrinos.end())
751  ++m_nTrueWithoutRecoHits;
752  }
753 
754  MCTruthToPFParticles::const_iterator pIter1 = matchedNeutrinos.find(trueEvent);
755  if (matchedNeutrinos.end() != pIter1) {
756  const art::Ptr<recob::PFParticle> recoParticle = pIter1->second;
757  m_pfoPdg = recoParticle->PdgCode();
758  m_pfoNuPdg = m_pfoPdg;
759  m_pfoParentPdg = m_pfoPdg;
760  m_pfoPrimaryPdg = 0;
761  m_pfoIsNeutrino = 1;
762  m_pfoIsPrimary = 0;
763 
764  if (!LArPandoraHelper::IsNeutrino(recoParticle))
765  std::cout << " Warning: Found neutrino with an invalid PDG code " << std::endl;
766 
767  PFParticlesToHits::const_iterator pIter2 = recoNeutrinosToHits.find(recoParticle);
768  if (recoParticlesToHits.end() != pIter2) {
769  const HitVector& recoHitVector = pIter2->second;
770 
771  for (HitVector::const_iterator hIter2 = recoHitVector.begin(),
772  hIterEnd2 = recoHitVector.end();
773  hIter2 != hIterEnd2;
774  ++hIter2) {
775  if (trueHitsToNeutrinos.find(*hIter2) == trueHitsToNeutrinos.end())
776  ++m_nRecoWithoutTrueHits;
777  }
778 
779  MCTruthToHits::const_iterator pIter3 = matchedNeutrinoHits.find(trueEvent);
780  if (matchedNeutrinoHits.end() != pIter3) {
781  const HitVector& matchedHitVector = pIter3->second;
782 
783  m_nPfoHits = recoHitVector.size();
784  m_nPfoHitsU = this->CountHitsByType(geo::kU, recoHitVector);
785  m_nPfoHitsV = this->CountHitsByType(geo::kV, recoHitVector);
786  m_nPfoHitsW = this->CountHitsByType(geo::kW, recoHitVector);
787 
788  m_nMatchedHits = matchedHitVector.size();
789  m_nMatchedHitsU = this->CountHitsByType(geo::kU, matchedHitVector);
790  m_nMatchedHitsV = this->CountHitsByType(geo::kV, matchedHitVector);
791  m_nMatchedHitsW = this->CountHitsByType(geo::kW, matchedHitVector);
792 
793  PFParticlesToVertices::const_iterator pIter4 =
794  recoParticlesToVertices.find(recoParticle);
795  if (recoParticlesToVertices.end() != pIter4) {
796  const VertexVector& vertexVector = pIter4->second;
797  if (!vertexVector.empty()) {
798  if (vertexVector.size() != 1 && m_printDebug)
799  std::cout << " Warning: Found particle with more than one associated vertex "
800  << std::endl;
801 
802  const art::Ptr<recob::Vertex> recoVertex = *(vertexVector.begin());
803  double xyz[3] = {0.0, 0.0, 0.0};
804  recoVertex->XYZ(xyz);
805 
806  m_pfoVertex = 1;
807  m_pfoVtxX = xyz[0];
808  m_pfoVtxY = xyz[1];
809  m_pfoVtxZ = xyz[2];
810  }
811  }
812  }
813  }
814  }
815 
816  m_purity =
817  ((m_nPfoHits == 0) ? 0.0 :
818  static_cast<double>(m_nMatchedHits) / static_cast<double>(m_nPfoHits));
819  m_completeness =
820  ((m_nPfoHits == 0) ? 0.0 :
821  static_cast<double>(m_nMatchedHits) / static_cast<double>(m_nMCHits));
822 
823  if (m_printDebug)
824  std::cout << " MCNeutrino [" << m_index << "]"
825  << " trueNu=" << m_mcNuPdg << ", truePdg=" << m_mcPdg
826  << ", recoNu=" << m_pfoNuPdg << ", recoPdg=" << m_pfoPdg
827  << ", mcHits=" << m_nMCHits << ", pfoHits=" << m_nPfoHits
828  << ", matchedHits=" << m_nMatchedHits
829  << ", availableHits=" << m_nTrueWithoutRecoHits << std::endl;
830 
831  m_pRecoTree->Fill();
832  ++m_index; // Increment index number
833  }
834 
835  // Match Reco Particles to True Particles
836  // ======================================
837  MCParticlesToPFParticles matchedParticles;
838  MCParticlesToHits matchedParticleHits;
839  this->GetRecoToTrueMatches(
840  recoParticlesToHits, trueHitsToParticles, matchedParticles, matchedParticleHits);
841 
842  // Compare true and reconstructed particles
843  for (MCParticlesToHits::const_iterator iter = trueParticlesToHits.begin(),
844  iterEnd = trueParticlesToHits.end();
845  iter != iterEnd;
846  ++iter) {
847  const art::Ptr<simb::MCParticle> trueParticle = iter->first;
848  const HitVector& trueHitVector = iter->second;
849 
850  if (trueHitVector.empty()) continue;
851 
852  m_mcPdg = trueParticle->PdgCode();
853  m_mcNuPdg = 0;
854  m_mcParentPdg = 0;
855  m_mcPrimaryPdg = 0;
856  m_mcIsNeutrino = 0;
857  m_mcIsPrimary = 0;
858  m_mcIsDecay = 0;
859  m_mcIsCC = 0;
860 
861  m_pfoPdg = 0;
862  m_pfoNuPdg = 0;
863  m_pfoParentPdg = 0;
864  m_pfoPrimaryPdg = 0;
865  m_pfoIsNeutrino = 0;
866  m_pfoIsPrimary = 0;
867  m_pfoIsStitched = 0;
868  m_pfoTrack = 0;
869  m_pfoVertex = 0;
870  m_pfoVtxX = 0.0;
871  m_pfoVtxY = 0.0;
872  m_pfoVtxZ = 0.0;
873  m_pfoEndX = 0.0;
874  m_pfoEndY = 0.0;
875  m_pfoEndZ = 0.0;
876  m_pfoDirX = 0.0;
877  m_pfoDirY = 0.0;
878  m_pfoDirZ = 0.0;
879  m_pfoLength = 0.0;
880  m_pfoStraightLength = 0.0;
881 
882  m_mcVertex = 0;
883  m_mcVtxX = 0.0;
884  m_mcVtxY = 0.0;
885  m_mcVtxZ = 0.0;
886  m_mcEndX = 0.0;
887  m_mcEndY = 0.0;
888  m_mcEndZ = 0.0;
889  m_mcDirX = 0.0;
890  m_mcDirY = 0.0;
891  m_mcDirZ = 0.0;
892  m_mcEnergy = 0.0;
893  m_mcLength = 0.0;
894  m_mcStraightLength = 0.0;
895 
896  m_completeness = 0.0;
897  m_purity = 0.0;
898 
899  m_nMCHits = 0;
900  m_nMCHitsU = 0;
901  m_nMCHitsV = 0;
902  m_nMCHitsW = 0;
903 
904  m_nPfoHits = 0;
905  m_nPfoHitsU = 0;
906  m_nPfoHitsV = 0;
907  m_nPfoHitsW = 0;
908 
909  m_nMatchedHits = 0;
910  m_nMatchedHitsU = 0;
911  m_nMatchedHitsV = 0;
912  m_nMatchedHitsW = 0;
913 
914  m_nTrueWithoutRecoHits = 0;
915  m_nRecoWithoutTrueHits = 0;
916 
917  m_spacepointsMinX = 0.0;
918  m_spacepointsMaxX = 0.0;
919 
920  // Set true properties
921  try {
922  int startT(-1);
923  int endT(-1);
924  this->GetStartAndEndPoints(trueParticle, startT, endT);
925 
926  // vertex and end positions
927  m_mcVertex = 1;
928  m_mcVtxX = trueParticle->Vx(startT);
929  m_mcVtxY = trueParticle->Vy(startT);
930  m_mcVtxZ = trueParticle->Vz(startT);
931  m_mcEndX = trueParticle->Vx(endT);
932  m_mcEndY = trueParticle->Vy(endT);
933  m_mcEndZ = trueParticle->Vz(endT);
934 
935  const double dx(m_mcEndX - m_mcVtxX);
936  const double dy(m_mcEndY - m_mcVtxY);
937  const double dz(m_mcEndZ - m_mcVtxZ);
938 
939  m_mcStraightLength = std::sqrt(dx * dx + dy * dy + dz * dz);
940  m_mcLength = this->GetLength(trueParticle, startT, endT);
941 
942  // energy and momentum
943  const double Ptot(trueParticle->P(startT));
944 
945  if (Ptot > 0.0) {
946  m_mcDirX = trueParticle->Px(startT) / Ptot;
947  m_mcDirY = trueParticle->Py(startT) / Ptot;
948  m_mcDirZ = trueParticle->Pz(startT) / Ptot;
949  m_mcEnergy = trueParticle->E(startT);
950  }
951  }
952  catch (cet::exception& e) {
953  }
954 
955  // Get the true parent neutrino
956  MCParticlesToMCTruth::const_iterator nuIter = particlesToTruth.find(trueParticle);
957  if (particlesToTruth.end() == nuIter)
958  throw cet::exception("LArPandora") << " PFParticleMonitoring::analyze --- Found a true "
959  "particle without any ancestry information ";
960 
961  const art::Ptr<simb::MCTruth> trueEvent = nuIter->second;
962 
963  if (trueEvent->NeutrinoSet()) {
964  const simb::MCNeutrino neutrino(trueEvent->GetNeutrino());
965  m_mcNuPdg = neutrino.Nu().PdgCode();
966  m_mcIsCC = ((simb::kCC == neutrino.CCNC()) ? 1 : 0);
967  }
968 
969  // Get the true 'parent' and 'primary' particles
970  try {
971  const art::Ptr<simb::MCParticle> parentParticle(
972  LArPandoraHelper::GetParentMCParticle(trueParticleMap, trueParticle));
973  const art::Ptr<simb::MCParticle> primaryParticle(
974  LArPandoraHelper::GetFinalStateMCParticle(trueParticleMap, trueParticle));
975  m_mcParentPdg = ((parentParticle != trueParticle) ? parentParticle->PdgCode() : 0);
976  m_mcPrimaryPdg = primaryParticle->PdgCode();
977  m_mcIsPrimary = (primaryParticle == trueParticle);
978  m_mcIsDecay = ("Decay" == trueParticle->Process());
979  }
980  catch (cet::exception& e) {
981  }
982 
983  // Find min and max X positions of space points
984  bool foundSpacePoints(false);
985 
986  for (HitVector::const_iterator hIter1 = trueHitVector.begin(),
987  hIterEnd1 = trueHitVector.end();
988  hIter1 != hIterEnd1;
989  ++hIter1) {
990  const art::Ptr<recob::Hit> hit = *hIter1;
991 
992  HitsToSpacePoints::const_iterator hIter2 = hitsToSpacePoints.find(hit);
993  if (hitsToSpacePoints.end() == hIter2) continue;
994 
995  const art::Ptr<recob::SpacePoint> spacepoint = hIter2->second;
996  const double X(spacepoint->XYZ()[0]);
997 
998  if (!foundSpacePoints) {
999  m_spacepointsMinX = X;
1000  m_spacepointsMaxX = X;
1001  foundSpacePoints = true;
1002  }
1003  else {
1004  m_spacepointsMinX = std::min(m_spacepointsMinX, X);
1005  m_spacepointsMaxX = std::max(m_spacepointsMaxX, X);
1006  }
1007  }
1008 
1009  // Count number of available hits
1010  for (HitVector::const_iterator hIter1 = trueHitVector.begin(),
1011  hIterEnd1 = trueHitVector.end();
1012  hIter1 != hIterEnd1;
1013  ++hIter1) {
1014  if (recoHitsToParticles.find(*hIter1) == recoHitsToParticles.end())
1015  ++m_nTrueWithoutRecoHits;
1016  }
1017 
1018  // Match true and reconstructed hits
1019  m_nMCHits = trueHitVector.size();
1020  m_nMCHitsU = this->CountHitsByType(geo::kU, trueHitVector);
1021  m_nMCHitsV = this->CountHitsByType(geo::kV, trueHitVector);
1022  m_nMCHitsW = this->CountHitsByType(geo::kW, trueHitVector);
1023 
1024  MCParticlesToPFParticles::const_iterator pIter1 = matchedParticles.find(trueParticle);
1025  if (matchedParticles.end() != pIter1) {
1026  const art::Ptr<recob::PFParticle> recoParticle = pIter1->second;
1027  m_pfoPdg = recoParticle->PdgCode();
1028  m_pfoNuPdg = LArPandoraHelper::GetParentNeutrino(recoParticleMap, recoParticle);
1029  m_pfoIsPrimary = LArPandoraHelper::IsFinalState(recoParticleMap, recoParticle);
1030 
1031  const art::Ptr<recob::PFParticle> parentParticle =
1032  LArPandoraHelper::GetParentPFParticle(recoParticleMap, recoParticle);
1033  m_pfoParentPdg = parentParticle->PdgCode();
1034 
1035  const art::Ptr<recob::PFParticle> primaryParticle =
1036  LArPandoraHelper::GetFinalStatePFParticle(recoParticleMap, recoParticle);
1037  m_pfoPrimaryPdg = primaryParticle->PdgCode();
1038 
1039  PFParticlesToHits::const_iterator pIter2 = recoParticlesToHits.find(recoParticle);
1040  if (recoParticlesToHits.end() == pIter2)
1041  throw cet::exception("LArPandora")
1042  << " PFParticleMonitoring::analyze --- Found a reco particle without any hits ";
1043 
1044  const HitVector& recoHitVector = pIter2->second;
1045 
1046  for (HitVector::const_iterator hIter2 = recoHitVector.begin(),
1047  hIterEnd2 = recoHitVector.end();
1048  hIter2 != hIterEnd2;
1049  ++hIter2) {
1050  if (trueHitsToParticles.find(*hIter2) == trueHitsToParticles.end())
1051  ++m_nRecoWithoutTrueHits;
1052  }
1053 
1054  MCParticlesToHits::const_iterator pIter3 = matchedParticleHits.find(trueParticle);
1055  if (matchedParticleHits.end() == pIter3)
1056  throw cet::exception("LArPandora") << " PFParticleMonitoring::analyze --- Found a "
1057  "matched true particle without matched hits ";
1058 
1059  const HitVector& matchedHitVector = pIter3->second;
1060 
1061  m_nPfoHits = recoHitVector.size();
1062  m_nPfoHitsU = this->CountHitsByType(geo::kU, recoHitVector);
1063  m_nPfoHitsV = this->CountHitsByType(geo::kV, recoHitVector);
1064  m_nPfoHitsW = this->CountHitsByType(geo::kW, recoHitVector);
1065 
1066  m_nMatchedHits = matchedHitVector.size();
1067  m_nMatchedHitsU = this->CountHitsByType(geo::kU, matchedHitVector);
1068  m_nMatchedHitsV = this->CountHitsByType(geo::kV, matchedHitVector);
1069  m_nMatchedHitsW = this->CountHitsByType(geo::kW, matchedHitVector);
1070 
1071  PFParticlesToVertices::const_iterator pIter4 = recoParticlesToVertices.find(recoParticle);
1072  if (recoParticlesToVertices.end() != pIter4) {
1073  const VertexVector& vertexVector = pIter4->second;
1074  if (!vertexVector.empty()) {
1075  if (vertexVector.size() != 1 && m_printDebug)
1076  std::cout << " Warning: Found particle with more than one associated vertex "
1077  << std::endl;
1078 
1079  const art::Ptr<recob::Vertex> recoVertex = *(vertexVector.begin());
1080  double xyz[3] = {0.0, 0.0, 0.0};
1081  recoVertex->XYZ(xyz);
1082 
1083  m_pfoVertex = 1;
1084  m_pfoVtxX = xyz[0];
1085  m_pfoVtxY = xyz[1];
1086  m_pfoVtxZ = xyz[2];
1087  }
1088  }
1089 
1090  PFParticlesToTracks::const_iterator pIter5 = recoParticlesToTracks.find(recoParticle);
1091  if (recoParticlesToTracks.end() != pIter5) {
1092  const TrackVector& trackVector = pIter5->second;
1093  if (!trackVector.empty()) {
1094  if (trackVector.size() != 1 && m_printDebug)
1095  std::cout << " Warning: Found particle with more than one associated track "
1096  << std::endl;
1097 
1098  const art::Ptr<recob::Track> recoTrack = *(trackVector.begin());
1099  const auto& vtxPosition = recoTrack->Vertex();
1100  const auto& endPosition = recoTrack->End();
1101  const auto& vtxDirection = recoTrack->VertexDirection();
1102 
1103  m_pfoTrack = 1;
1104  m_pfoVtxX = vtxPosition.x();
1105  m_pfoVtxY = vtxPosition.y();
1106  m_pfoVtxZ = vtxPosition.z();
1107  m_pfoEndX = endPosition.x();
1108  m_pfoEndY = endPosition.y();
1109  m_pfoEndZ = endPosition.z();
1110  m_pfoDirX = vtxDirection.x();
1111  m_pfoDirY = vtxDirection.y();
1112  m_pfoDirZ = vtxDirection.z();
1113  m_pfoStraightLength = (endPosition - vtxPosition).R();
1114  m_pfoLength = recoTrack->Length();
1115  }
1116  }
1117 
1118  m_pfoIsStitched = (particlesToT0s.end() != particlesToT0s.find(recoParticle));
1119  }
1120 
1121  m_purity =
1122  ((m_nPfoHits == 0) ? 0.0 :
1123  static_cast<double>(m_nMatchedHits) / static_cast<double>(m_nPfoHits));
1124  m_completeness =
1125  ((m_nPfoHits == 0) ? 0.0 :
1126  static_cast<double>(m_nMatchedHits) / static_cast<double>(m_nMCHits));
1127 
1128  if (m_printDebug)
1129  std::cout << " MCParticle [" << m_index << "]"
1130  << " trueNu=" << m_mcNuPdg << ", truePdg=" << m_mcPdg
1131  << ", recoNu=" << m_pfoNuPdg << ", recoPdg=" << m_pfoPdg
1132  << ", mcHits=" << m_nMCHits << ", pfoHits=" << m_nPfoHits
1133  << ", matchedHits=" << m_nMatchedHits
1134  << ", availableHits=" << m_nTrueWithoutRecoHits << std::endl;
1135 
1136  m_pRecoTree->Fill();
1137  ++m_index; // Increment index number
1138  }
1139  }
1140 
1141  //------------------------------------------------------------------------------------------------------------------------------------------
1142 
1143  void
1144  PFParticleMonitoring::BuildTrueNeutrinoHitMaps(const MCTruthToMCParticles& truthToParticles,
1145  const MCParticlesToHits& trueParticlesToHits,
1146  MCTruthToHits& trueNeutrinosToHits,
1147  HitsToMCTruth& trueHitsToNeutrinos) const
1148  {
1149  for (MCTruthToMCParticles::const_iterator iter1 = truthToParticles.begin(),
1150  iterEnd1 = truthToParticles.end();
1151  iter1 != iterEnd1;
1152  ++iter1) {
1153  const art::Ptr<simb::MCTruth> trueNeutrino = iter1->first;
1154  const MCParticleVector& trueParticleVector = iter1->second;
1155 
1156  for (MCParticleVector::const_iterator iter2 = trueParticleVector.begin(),
1157  iterEnd2 = trueParticleVector.end();
1158  iter2 != iterEnd2;
1159  ++iter2) {
1160  const MCParticlesToHits::const_iterator iter3 = trueParticlesToHits.find(*iter2);
1161  if (trueParticlesToHits.end() == iter3) continue;
1162 
1163  const HitVector& hitVector = iter3->second;
1164 
1165  for (HitVector::const_iterator iter4 = hitVector.begin(), iterEnd4 = hitVector.end();
1166  iter4 != iterEnd4;
1167  ++iter4) {
1168  const art::Ptr<recob::Hit> hit = *iter4;
1169  trueHitsToNeutrinos[hit] = trueNeutrino;
1170  trueNeutrinosToHits[trueNeutrino].push_back(hit);
1171  }
1172  }
1173  }
1174  }
1175 
1176  //------------------------------------------------------------------------------------------------------------------------------------------
1177 
1178  void
1179  PFParticleMonitoring::BuildRecoNeutrinoHitMaps(const PFParticleMap& recoParticleMap,
1180  const PFParticlesToHits& recoParticlesToHits,
1181  PFParticlesToHits& recoNeutrinosToHits,
1182  HitsToPFParticles& recoHitsToNeutrinos) const
1183  {
1184  for (PFParticleMap::const_iterator iter1 = recoParticleMap.begin(),
1185  iterEnd1 = recoParticleMap.end();
1186  iter1 != iterEnd1;
1187  ++iter1) {
1188  const art::Ptr<recob::PFParticle> recoParticle = iter1->second;
1189  const art::Ptr<recob::PFParticle> recoNeutrino =
1190  LArPandoraHelper::GetParentPFParticle(recoParticleMap, recoParticle);
1191 
1192  if (!LArPandoraHelper::IsNeutrino(recoNeutrino)) continue;
1193 
1194  const PFParticlesToHits::const_iterator iter2 = recoParticlesToHits.find(recoParticle);
1195  if (recoParticlesToHits.end() == iter2) continue;
1196 
1197  const HitVector& hitVector = iter2->second;
1198 
1199  for (HitVector::const_iterator iter3 = hitVector.begin(), iterEnd3 = hitVector.end();
1200  iter3 != iterEnd3;
1201  ++iter3) {
1202  const art::Ptr<recob::Hit> hit = *iter3;
1203  recoHitsToNeutrinos[hit] = recoNeutrino;
1204  recoNeutrinosToHits[recoNeutrino].push_back(hit);
1205  }
1206  }
1207  }
1208 
1209  //------------------------------------------------------------------------------------------------------------------------------------------
1210 
1211  void
1212  PFParticleMonitoring::GetRecoToTrueMatches(const PFParticlesToHits& recoNeutrinosToHits,
1213  const HitsToMCTruth& trueHitsToNeutrinos,
1214  MCTruthToPFParticles& matchedNeutrinos,
1215  MCTruthToHits& matchedNeutrinoHits) const
1216  {
1217  PFParticleSet recoVeto;
1218  MCTruthSet trueVeto;
1219 
1220  this->GetRecoToTrueMatches(recoNeutrinosToHits,
1221  trueHitsToNeutrinos,
1222  matchedNeutrinos,
1223  matchedNeutrinoHits,
1224  recoVeto,
1225  trueVeto);
1226  }
1227 
1228  //------------------------------------------------------------------------------------------------------------------------------------------
1229 
1230  void
1231  PFParticleMonitoring::GetRecoToTrueMatches(const PFParticlesToHits& recoNeutrinosToHits,
1232  const HitsToMCTruth& trueHitsToNeutrinos,
1233  MCTruthToPFParticles& matchedNeutrinos,
1234  MCTruthToHits& matchedNeutrinoHits,
1235  PFParticleSet& vetoReco,
1236  MCTruthSet& vetoTrue) const
1237  {
1238  bool foundMatches(false);
1239 
1240  for (PFParticlesToHits::const_iterator iter1 = recoNeutrinosToHits.begin(),
1241  iterEnd1 = recoNeutrinosToHits.end();
1242  iter1 != iterEnd1;
1243  ++iter1) {
1244  const art::Ptr<recob::PFParticle> recoNeutrino = iter1->first;
1245  if (vetoReco.count(recoNeutrino) > 0) continue;
1246 
1247  const HitVector& hitVector = iter1->second;
1248 
1249  MCTruthToHits truthContributionMap;
1250 
1251  for (HitVector::const_iterator iter2 = hitVector.begin(), iterEnd2 = hitVector.end();
1252  iter2 != iterEnd2;
1253  ++iter2) {
1254  const art::Ptr<recob::Hit> hit = *iter2;
1255 
1256  HitsToMCTruth::const_iterator iter3 = trueHitsToNeutrinos.find(hit);
1257  if (trueHitsToNeutrinos.end() == iter3) continue;
1258 
1259  const art::Ptr<simb::MCTruth> trueNeutrino = iter3->second;
1260  if (vetoTrue.count(trueNeutrino) > 0) continue;
1261 
1262  truthContributionMap[trueNeutrino].push_back(hit);
1263  }
1264 
1265  MCTruthToHits::const_iterator mIter = truthContributionMap.end();
1266 
1267  for (MCTruthToHits::const_iterator iter4 = truthContributionMap.begin(),
1268  iterEnd4 = truthContributionMap.end();
1269  iter4 != iterEnd4;
1270  ++iter4) {
1271  if ((truthContributionMap.end() == mIter) ||
1272  (iter4->second.size() > mIter->second.size())) {
1273  mIter = iter4;
1274  }
1275  }
1276 
1277  if (truthContributionMap.end() != mIter) {
1278  const art::Ptr<simb::MCTruth> trueNeutrino = mIter->first;
1279 
1280  MCTruthToHits::const_iterator iter5 = matchedNeutrinoHits.find(trueNeutrino);
1281 
1282  if ((matchedNeutrinoHits.end() == iter5) || (mIter->second.size() > iter5->second.size())) {
1283  matchedNeutrinos[trueNeutrino] = recoNeutrino;
1284  matchedNeutrinoHits[trueNeutrino] = mIter->second;
1285  foundMatches = true;
1286  }
1287  }
1288  }
1289 
1290  if (!foundMatches) return;
1291 
1292  for (MCTruthToPFParticles::const_iterator pIter = matchedNeutrinos.begin(),
1293  pIterEnd = matchedNeutrinos.end();
1294  pIter != pIterEnd;
1295  ++pIter) {
1296  vetoTrue.insert(pIter->first);
1297  vetoReco.insert(pIter->second);
1298  }
1299 
1300  if (m_recursiveMatching)
1301  this->GetRecoToTrueMatches(recoNeutrinosToHits,
1302  trueHitsToNeutrinos,
1303  matchedNeutrinos,
1304  matchedNeutrinoHits,
1305  vetoReco,
1306  vetoTrue);
1307  }
1308 
1309  //------------------------------------------------------------------------------------------------------------------------------------------
1310 
1311  void
1312  PFParticleMonitoring::GetRecoToTrueMatches(const PFParticlesToHits& recoParticlesToHits,
1313  const HitsToMCParticles& trueHitsToParticles,
1314  MCParticlesToPFParticles& matchedParticles,
1315  MCParticlesToHits& matchedHits) const
1316  {
1317  PFParticleSet recoVeto;
1318  MCParticleSet trueVeto;
1319 
1320  this->GetRecoToTrueMatches(
1321  recoParticlesToHits, trueHitsToParticles, matchedParticles, matchedHits, recoVeto, trueVeto);
1322  }
1323 
1324  //------------------------------------------------------------------------------------------------------------------------------------------
1325 
1326  void
1327  PFParticleMonitoring::GetRecoToTrueMatches(const PFParticlesToHits& recoParticlesToHits,
1328  const HitsToMCParticles& trueHitsToParticles,
1329  MCParticlesToPFParticles& matchedParticles,
1330  MCParticlesToHits& matchedHits,
1331  PFParticleSet& vetoReco,
1332  MCParticleSet& vetoTrue) const
1333  {
1334  bool foundMatches(false);
1335 
1336  for (PFParticlesToHits::const_iterator iter1 = recoParticlesToHits.begin(),
1337  iterEnd1 = recoParticlesToHits.end();
1338  iter1 != iterEnd1;
1339  ++iter1) {
1340  const art::Ptr<recob::PFParticle> recoParticle = iter1->first;
1341  if (vetoReco.count(recoParticle) > 0) continue;
1342 
1343  const HitVector& hitVector = iter1->second;
1344 
1345  MCParticlesToHits truthContributionMap;
1346 
1347  for (HitVector::const_iterator iter2 = hitVector.begin(), iterEnd2 = hitVector.end();
1348  iter2 != iterEnd2;
1349  ++iter2) {
1350  const art::Ptr<recob::Hit> hit = *iter2;
1351 
1352  HitsToMCParticles::const_iterator iter3 = trueHitsToParticles.find(hit);
1353  if (trueHitsToParticles.end() == iter3) continue;
1354 
1355  const art::Ptr<simb::MCParticle> trueParticle = iter3->second;
1356  if (vetoTrue.count(trueParticle) > 0) continue;
1357 
1358  truthContributionMap[trueParticle].push_back(hit);
1359  }
1360 
1361  MCParticlesToHits::const_iterator mIter = truthContributionMap.end();
1362 
1363  for (MCParticlesToHits::const_iterator iter4 = truthContributionMap.begin(),
1364  iterEnd4 = truthContributionMap.end();
1365  iter4 != iterEnd4;
1366  ++iter4) {
1367  if ((truthContributionMap.end() == mIter) ||
1368  (iter4->second.size() > mIter->second.size())) {
1369  mIter = iter4;
1370  }
1371  }
1372 
1373  if (truthContributionMap.end() != mIter) {
1374  const art::Ptr<simb::MCParticle> trueParticle = mIter->first;
1375 
1376  MCParticlesToHits::const_iterator iter5 = matchedHits.find(trueParticle);
1377 
1378  if ((matchedHits.end() == iter5) || (mIter->second.size() > iter5->second.size())) {
1379  matchedParticles[trueParticle] = recoParticle;
1380  matchedHits[trueParticle] = mIter->second;
1381  foundMatches = true;
1382  }
1383  }
1384  }
1385 
1386  if (!foundMatches) return;
1387 
1388  for (MCParticlesToPFParticles::const_iterator pIter = matchedParticles.begin(),
1389  pIterEnd = matchedParticles.end();
1390  pIter != pIterEnd;
1391  ++pIter) {
1392  vetoTrue.insert(pIter->first);
1393  vetoReco.insert(pIter->second);
1394  }
1395 
1396  if (m_recursiveMatching)
1397  this->GetRecoToTrueMatches(recoParticlesToHits,
1398  trueHitsToParticles,
1399  matchedParticles,
1400  matchedHits,
1401  vetoReco,
1402  vetoTrue);
1403  }
1404 
1405  //------------------------------------------------------------------------------------------------------------------------------------------
1406 
1407  int
1408  PFParticleMonitoring::CountHitsByType(const int view, const HitVector& hitVector) const
1409  {
1410  int nHits(0);
1411 
1412  for (HitVector::const_iterator iter = hitVector.begin(), iterEnd = hitVector.end();
1413  iter != iterEnd;
1414  ++iter) {
1415  const art::Ptr<recob::Hit> hit = *iter;
1416  if (hit->View() == view) ++nHits;
1417  }
1418 
1419  return nHits;
1420  }
1421 
1422  //------------------------------------------------------------------------------------------------------------------------------------------
1423 
1424  void
1425  PFParticleMonitoring::GetStartAndEndPoints(const art::Ptr<simb::MCParticle> particle,
1426  int& startT,
1427  int& endT) const
1428  {
1429  art::ServiceHandle<geo::Geometry const> theGeometry;
1430 
1431  bool foundStartPosition(false);
1432 
1433  const int numTrajectoryPoints(static_cast<int>(particle->NumberTrajectoryPoints()));
1434 
1435  for (int nt = 0; nt < numTrajectoryPoints; ++nt) {
1436  try {
1437  double pos[3] = {particle->Vx(nt), particle->Vy(nt), particle->Vz(nt)};
1438  unsigned int which_tpc(std::numeric_limits<unsigned int>::max());
1439  unsigned int which_cstat(std::numeric_limits<unsigned int>::max());
1440  theGeometry->PositionToTPC(pos, which_tpc, which_cstat);
1441 
1442  // TODO: Apply fiducial cut due to readout window
1443 
1444  endT = nt;
1445  if (!foundStartPosition) {
1446  startT = endT;
1447  foundStartPosition = true;
1448  }
1449  }
1450  catch (cet::exception& e) {
1451  continue;
1452  }
1453  }
1454 
1455  if (!foundStartPosition) throw cet::exception("LArPandora");
1456  }
1457 
1458  //------------------------------------------------------------------------------------------------------------------------------------------
1459 
1460  double
1461  PFParticleMonitoring::GetLength(const art::Ptr<simb::MCParticle> particle,
1462  const int startT,
1463  const int endT) const
1464  {
1465  if (endT <= startT) return 0.0;
1466 
1467  double length(0.0);
1468 
1469  for (int nt = startT; nt < endT; ++nt) {
1470  const double dx(particle->Vx(nt + 1) - particle->Vx(nt));
1471  const double dy(particle->Vy(nt + 1) - particle->Vy(nt));
1472  const double dz(particle->Vz(nt + 1) - particle->Vz(nt));
1473  length += sqrt(dx * dx + dy * dy + dz * dz);
1474  }
1475 
1476  return length;
1477  }
1478 
1479 } //namespace lar_pandora
lar_pandora::LArPandoraHelper::BuildPFParticleHitMaps
static void BuildPFParticleHitMaps(const PFParticleVector &particleVector, const PFParticlesToSpacePoints &particlesToSpacePoints, const SpacePointsToHits &spacePointsToHits, PFParticlesToHits &particlesToHits, HitsToPFParticles &hitsToParticles, const DaughterMode daughterMode=kUseDaughters)
Build mapping between PFParticles and Hits using PFParticle/SpacePoint/Hit maps.
Definition: LArPandoraHelper.cxx:524
lar_pandora::PFParticleMonitoring::m_pfoIsStitched
int m_pfoIsStitched
Definition: PFParticleMonitoring_module.cc:190
lar_pandora::PFParticleMonitoring::m_pfoIsNeutrino
int m_pfoIsNeutrino
Definition: PFParticleMonitoring_module.cc:188
lar_pandora::PFParticleMonitoring::m_nNeutrinoPfos
int m_nNeutrinoPfos
Definition: PFParticleMonitoring_module.cc:171
lar_pandora::PFParticleMonitoring::analyze
void analyze(const art::Event &evt)
Definition: PFParticleMonitoring_module.cc:419
lar_pandora::HitsToMCParticles
std::map< art::Ptr< recob::Hit >, art::Ptr< simb::MCParticle > > HitsToMCParticles
Definition: LArPandoraHelper.h:99
lar_pandora::PFParticleMonitoring::CountHitsByType
int CountHitsByType(const int view, const HitVector &hitVector) const
Count the number of reconstructed hits in a given wire plane.
Definition: PFParticleMonitoring_module.cc:1408
lar_pandora::LArPandoraHelper::BuildPFParticleMap
static void BuildPFParticleMap(const PFParticleVector &particleVector, PFParticleMap &particleMap)
Build particle maps for reconstructed particles.
Definition: LArPandoraHelper.cxx:1201
lar_pandora::PFParticleMonitoring::m_mcParentPdg
int m_mcParentPdg
Definition: PFParticleMonitoring_module.cc:177
lar_pandora::LArPandoraHelper::GetParentNeutrino
static int GetParentNeutrino(const PFParticleMap &particleMap, const art::Ptr< recob::PFParticle > daughterParticle)
Return the parent neutrino PDG code (or zero for cosmics) for a given reconstructed particle...
Definition: LArPandoraHelper.cxx:1396
lar_pandora::PFParticleMonitoring::m_useDaughterPFParticles
bool m_useDaughterPFParticles
Definition: PFParticleMonitoring_module.cc:253
CryostatGeo.h
Encapsulate the construction of a single cyostat.
lar_pandora::PFParticleMonitoring::m_pfoTrack
int m_pfoTrack
Definition: PFParticleMonitoring_module.cc:192
lar_pandora::PFParticleMonitoring::m_mcVtxX
double m_mcVtxX
Definition: PFParticleMonitoring_module.cc:207
lar_pandora::PFParticleMonitoring::m_mcDirX
double m_mcDirX
Definition: PFParticleMonitoring_module.cc:213
lar_pandora::PFParticleMonitoring::GetRecoToTrueMatches
void GetRecoToTrueMatches(const PFParticlesToHits &recoNeutrinosToHits, const HitsToMCTruth &trueHitsToNeutrinos, MCTruthToPFParticles &matchedNeutrinos, MCTruthToHits &matchedNeutrinoHits) const
Perform matching between true and reconstructed neutrino events.
Definition: PFParticleMonitoring_module.cc:1212
lar_pandora::PFParticleMonitoring::m_pfoIsPrimary
int m_pfoIsPrimary
Definition: PFParticleMonitoring_module.cc:189
lar_pandora::PFParticleMonitoring::PFParticleSet
std::set< art::Ptr< recob::PFParticle > > PFParticleSet
Definition: PFParticleMonitoring_module.cc:44
lar_pandora::PFParticleMonitoring::m_pfoVtxZ
double m_pfoVtxZ
Definition: PFParticleMonitoring_module.cc:196
lar_pandora::LArPandoraHelper::IsNeutrino
static bool IsNeutrino(const art::Ptr< recob::PFParticle > particle)
Determine whether a particle has been reconstructed as a neutrino.
Definition: LArPandoraHelper.cxx:1444
lar_pandora::PFParticleMonitoring::m_backtrackerLabel
std::string m_backtrackerLabel
Definition: PFParticleMonitoring_module.cc:250
lar_pandora::MCTruthToHits
std::map< art::Ptr< simb::MCTruth >, HitVector > MCTruthToHits
Definition: LArPandoraHelper.h:91
geo::kV
Planes which measure V.
Definition: geo_types.h:130
lar_pandora::PFParticleMonitoring::m_purity
double m_purity
Definition: PFParticleMonitoring_module.cc:221
Hit.h
Declaration of signal hit object.
lar_pandora::MCParticlesToMCTruth
std::map< art::Ptr< simb::MCParticle >, art::Ptr< simb::MCTruth > > MCParticlesToMCTruth
Definition: LArPandoraHelper.h:93
lar_pandora::PFParticleMonitoring::m_pfoDirY
double m_pfoDirY
Definition: PFParticleMonitoring_module.cc:201
lar_pandora::PFParticleMonitoring::m_geantModuleLabel
std::string m_geantModuleLabel
Definition: PFParticleMonitoring_module.cc:251
lar_pandora::PFParticleMonitoring::m_mcVtxY
double m_mcVtxY
Definition: PFParticleMonitoring_module.cc:208
lar_pandora::LArPandoraHelper::GetFinalStatePFParticle
static art::Ptr< recob::PFParticle > GetFinalStatePFParticle(const PFParticleMap &particleMap, const art::Ptr< recob::PFParticle > daughterParticle)
Return the final-state parent particle by navigating up the chain of parent/daughter associations...
Definition: LArPandoraHelper.cxx:1248
lar_pandora::PFParticleMonitoring::m_pfoEndY
double m_pfoEndY
Definition: PFParticleMonitoring_module.cc:198
lar_pandora::PFParticleMonitoring::m_nPfoHitsU
int m_nPfoHitsU
Definition: PFParticleMonitoring_module.cc:231
lar_pandora::PFParticleMonitoring::m_nMatchedHitsU
int m_nMatchedHitsU
Definition: PFParticleMonitoring_module.cc:235
lar_pandora::PFParticleMonitoring::m_nMatchedHitsW
int m_nMatchedHitsW
Definition: PFParticleMonitoring_module.cc:237
lar_pandora::MCParticlesToPFParticles
std::map< art::Ptr< simb::MCParticle >, art::Ptr< recob::PFParticle > > MCParticlesToPFParticles
Definition: LArPandoraHelper.h:96
lar_pandora::PFParticleMonitoring::m_pfoParentPdg
int m_pfoParentPdg
Definition: PFParticleMonitoring_module.cc:186
lar_pandora::PFParticleMonitoring::m_mcLength
double m_mcLength
Definition: PFParticleMonitoring_module.cc:217
lar_pandora::PFParticleMonitoring::m_trackLabel
std::string m_trackLabel
Definition: PFParticleMonitoring_module.cc:248
lar_pandora::PFParticleMonitoring::m_addDaughterMCParticles
bool m_addDaughterMCParticles
Definition: PFParticleMonitoring_module.cc:256
lar_pandora::PFParticleMonitoring::m_pfoVtxY
double m_pfoVtxY
Definition: PFParticleMonitoring_module.cc:195
lar_pandora::PFParticleMonitoring::m_mcIsCC
int m_mcIsCC
Definition: PFParticleMonitoring_module.cc:182
hit
process_name hit
Definition: cheaterreco.fcl:51
X
then echo echo For and will not be changed by echo further linking echo echo B echo The symbol is in the uninitialized data multiple common symbols may appear with the echo same name If the symbol is defined the common echo symbols are treated as undefined references For more echo details on common see the discussion of warn common echo in *Note Linker see the discussion of warn common echo in *Note Linker such as a global int variable echo as opposed to a large global array echo echo I echo The symbol is an indirect reference to another symbol This echo is a GNU extension to the a out object file format which is echo rarely used echo echo N echo The symbol is a debugging symbol echo echo R echo The symbol is in a read only data section echo echo S echo The symbol is in an uninitialized data section for small echo objects echo echo T echo The symbol is in the the normal defined echo symbol is used with no error When a weak undefined symbol echo is linked and the symbol is not the value of the echo weak symbol becomes zero with no error echo echo W echo The symbol is a weak symbol that has not been specifically echo tagged as a weak object symbol When a weak defined symbol echo is linked with a normal defined the normal defined echo symbol is used with no error When a weak undefined symbol echo is linked and the symbol is not the value of the echo weak symbol becomes zero with no error echo echo echo The symbol is a stabs symbol in an a out object file In echo this the next values printed are the stabs other echo the stabs desc and the stab type Stabs symbols are echo used to hold debugging information For more echo see *Note or object file format specific echo echo For Mac OS X
Definition: find_global_symbol.sh:107
lar_pandora::PFParticleMonitoring::m_mcIsNeutrino
int m_mcIsNeutrino
Definition: PFParticleMonitoring_module.cc:179
lar_pandora::PFParticleMonitoring::m_nMCHitsU
int m_nMCHitsU
Definition: PFParticleMonitoring_module.cc:227
lar_pandora::PFParticleMonitoring::m_mcPdg
int m_mcPdg
Definition: PFParticleMonitoring_module.cc:175
lar_pandora::PFParticleMonitoring::m_pRecoTree
TTree * m_pRecoTree
Definition: PFParticleMonitoring_module.cc:164
lar_pandora::LArPandoraHelper::GetParentPFParticle
static art::Ptr< recob::PFParticle > GetParentPFParticle(const PFParticleMap &particleMap, const art::Ptr< recob::PFParticle > daughterParticle)
Return the top-level parent particle by navigating up the chain of parent/daughter associations...
Definition: LArPandoraHelper.cxx:1216
lar_pandora::PFParticleMonitoring::m_pfoPrimaryPdg
int m_pfoPrimaryPdg
Definition: PFParticleMonitoring_module.cc:187
lar_pandora::PFParticlesToTracks
std::map< art::Ptr< recob::PFParticle >, TrackVector > PFParticlesToTracks
Definition: LArPandoraHelper.h:77
lar_pandora::PFParticleMonitoring::m_mcDirZ
double m_mcDirZ
Definition: PFParticleMonitoring_module.cc:215
lar_pandora::PFParticleMonitoring::m_mcPrimaryPdg
int m_mcPrimaryPdg
Definition: PFParticleMonitoring_module.cc:178
lar_pandora::PFParticlesToVertices
std::map< art::Ptr< recob::PFParticle >, VertexVector > PFParticlesToVertices
Definition: LArPandoraHelper.h:81
lar_pandora::PFParticleMap
std::map< int, art::Ptr< recob::PFParticle > > PFParticleMap
Definition: LArPandoraHelper.h:105
lar_pandora::PFParticleMonitoring::m_pfoVtxX
double m_pfoVtxX
Definition: PFParticleMonitoring_module.cc:194
lar_pandora::PFParticleMonitoring::BuildTrueNeutrinoHitMaps
void BuildTrueNeutrinoHitMaps(const MCTruthToMCParticles &truthToParticles, const MCParticlesToHits &trueParticlesToHits, MCTruthToHits &trueNeutrinosToHits, HitsToMCTruth &trueHitsToNeutrinos) const
Build mapping from true neutrinos to hits.
Definition: PFParticleMonitoring_module.cc:1144
lar_pandora::PFParticleMonitoring::m_pfoEndX
double m_pfoEndX
Definition: PFParticleMonitoring_module.cc:197
lar_pandora::MCParticleMap
std::map< int, art::Ptr< simb::MCParticle > > MCParticleMap
Definition: LArPandoraHelper.h:109
lar_pandora::PFParticleMonitoring::m_useDaughterMCParticles
bool m_useDaughterMCParticles
Definition: PFParticleMonitoring_module.cc:254
lar_pandora::PFParticleMonitoring::m_pfoPdg
int m_pfoPdg
Definition: PFParticleMonitoring_module.cc:184
lar_pandora::PFParticleMonitoring::MCTruthSet
std::set< art::Ptr< simb::MCTruth > > MCTruthSet
Definition: PFParticleMonitoring_module.cc:46
lar_pandora::PFParticleMonitoring::m_hitfinderLabel
std::string m_hitfinderLabel
Definition: PFParticleMonitoring_module.cc:247
lar_pandora::PFParticleMonitoring::m_event
int m_event
Definition: PFParticleMonitoring_module.cc:167
lar_pandora::PFParticleMonitoring::m_nMatchedHits
int m_nMatchedHits
Definition: PFParticleMonitoring_module.cc:225
lar_pandora::LArPandoraHelper::CollectSpacePoints
static void CollectSpacePoints(const art::Event &evt, const std::string &label, SpacePointVector &spacePointVector, SpacePointsToHits &spacePointsToHits)
Collect the reconstructed SpacePoints and associated hits from the ART event record.
Definition: LArPandoraHelper.cxx:121
lar_pandora::PFParticleMonitoring::m_mcVertex
int m_mcVertex
Definition: PFParticleMonitoring_module.cc:206
lar_pandora::PFParticleMonitoring::m_nPrimaryPfos
int m_nPrimaryPfos
Definition: PFParticleMonitoring_module.cc:172
geo::kU
Planes which measure U.
Definition: geo_types.h:129
ana::Current::kCC
Charged-current interactions.
Definition: IPrediction.h:38
lar_pandora::PFParticleMonitoring::GetLength
double GetLength(const art::Ptr< simb::MCParticle > trueParticle, const int startT, const int endT) const
Find the length of the true particle trajectory through the active region of the detector.
Definition: PFParticleMonitoring_module.cc:1461
lar_pandora::PFParticleMonitoring::m_nTrueWithoutRecoHits
int m_nTrueWithoutRecoHits
True hits which don&#39;t belong to any reconstructed particle - &quot;available&quot;.
Definition: PFParticleMonitoring_module.cc:240
lar_pandora::PFParticlesToT0s
std::map< art::Ptr< recob::PFParticle >, T0Vector > PFParticlesToT0s
Definition: LArPandoraHelper.h:103
lar_pandora::HitsToMCTruth
std::map< art::Ptr< recob::Hit >, art::Ptr< simb::MCTruth > > HitsToMCTruth
Definition: LArPandoraHelper.h:100
lar_pandora::LArPandoraHelper::GetParentMCParticle
static art::Ptr< simb::MCParticle > GetParentMCParticle(const MCParticleMap &particleMap, const art::Ptr< simb::MCParticle > daughterParticle)
Return the top-level parent particle by navigating up the chain of parent/daughter associations...
Definition: LArPandoraHelper.cxx:1286
lar_pandora::PFParticleMonitoring::BuildRecoNeutrinoHitMaps
void BuildRecoNeutrinoHitMaps(const PFParticleMap &recoParticleMap, const PFParticlesToHits &recoParticlesToHits, PFParticlesToHits &recoNeutrinosToHits, HitsToPFParticles &recoHitsToNeutrinos) const
Build mapping from reconstructed neutrinos to hits.
Definition: PFParticleMonitoring_module.cc:1179
lar_pandora::PFParticleMonitoring::m_nPfoHitsV
int m_nPfoHitsV
Definition: PFParticleMonitoring_module.cc:232
lar_pandora::MCParticlesToHits
std::map< art::Ptr< simb::MCParticle >, HitVector > MCParticlesToHits
Definition: LArPandoraHelper.h:94
lar_pandora::PFParticleMonitoring::m_recursiveMatching
bool m_recursiveMatching
Definition: PFParticleMonitoring_module.cc:258
lar_pandora::PFParticleMonitoring::m_pfoLength
double m_pfoLength
Definition: PFParticleMonitoring_module.cc:203
lar_pandora::PFParticleVector
std::vector< art::Ptr< recob::PFParticle > > PFParticleVector
Definition: LArPandoraHelper.h:65
lar_pandora::MCTruthToMCParticles
std::map< art::Ptr< simb::MCTruth >, MCParticleVector > MCTruthToMCParticles
Definition: LArPandoraHelper.h:90
lar_pandora::MCTruthToPFParticles
std::map< art::Ptr< simb::MCTruth >, art::Ptr< recob::PFParticle > > MCTruthToPFParticles
Definition: LArPandoraHelper.h:92
lar_pandora::SpacePointVector
std::vector< art::Ptr< recob::SpacePoint > > SpacePointVector
Definition: LArPandoraHelper.h:59
lar_pandora::PFParticlesToHits
std::map< art::Ptr< recob::PFParticle >, HitVector > PFParticlesToHits
Definition: LArPandoraHelper.h:83
lar_pandora::PFParticleMonitoring::m_printDebug
bool m_printDebug
switch for print statements (TODO: use message service!)
Definition: PFParticleMonitoring_module.cc:259
lar_pandora::LArPandoraHelper::CollectVertices
static void CollectVertices(const art::Event &evt, const std::string &label, VertexVector &vertexVector, PFParticlesToVertices &particlesToVertices)
Collect the reconstructed PFParticles and associated Vertices from the ART event record.
Definition: LArPandoraHelper.cxx:492
lar_pandora::PFParticleMonitoring::m_nMCHits
int m_nMCHits
Definition: PFParticleMonitoring_module.cc:223
lar_pandora::PFParticleMonitoring::PFParticleMonitoring
PFParticleMonitoring(fhicl::ParameterSet const &pset)
Constructor.
Definition: PFParticleMonitoring_module.cc:300
lar_pandora::HitsToSpacePoints
std::map< art::Ptr< recob::Hit >, art::Ptr< recob::SpacePoint > > HitsToSpacePoints
Definition: LArPandoraHelper.h:97
lar_pandora::PFParticleMonitoring::m_mcEndX
double m_mcEndX
Definition: PFParticleMonitoring_module.cc:210
lar_pandora::MCParticleVector
std::vector< art::Ptr< simb::MCParticle > > MCParticleVector
Definition: LArPandoraHelper.h:67
lar_pandora::PFParticleMonitoring::m_spacepointsMaxX
double m_spacepointsMaxX
Definition: PFParticleMonitoring_module.cc:245
lar_pandora::PFParticleMonitoring::m_completeness
double m_completeness
Definition: PFParticleMonitoring_module.cc:220
lar_pandora::PFParticleMonitoring::m_spacepointsMinX
double m_spacepointsMinX
Definition: PFParticleMonitoring_module.cc:244
lar_pandora::TrackVector
std::vector< art::Ptr< recob::Track > > TrackVector
Definition: LArPandoraHelper.h:63
Cluster.h
Declaration of cluster object.
lar_pandora::LArPandoraHelper::CollectTracks
static void CollectTracks(const art::Event &evt, const std::string &label, TrackVector &trackVector, PFParticlesToTracks &particlesToTracks)
Collect the reconstructed PFParticles and associated Tracks from the ART event record.
Definition: LArPandoraHelper.cxx:330
lar_pandora::PFParticleMonitoring::m_mcEnergy
double m_mcEnergy
Definition: PFParticleMonitoring_module.cc:216
lar_pandora::LArPandoraHelper::CollectPFParticles
static void CollectPFParticles(const art::Event &evt, const std::string &label, PFParticleVector &particleVector)
Collect the reconstructed PFParticles from the ART event record.
Definition: LArPandoraHelper.cxx:96
lar_pandora::PFParticleMonitoring::GetStartAndEndPoints
void GetStartAndEndPoints(const art::Ptr< simb::MCParticle > trueParticle, int &startT, int &endT) const
Find the start and end points of the true particle in the active region of detector.
Definition: PFParticleMonitoring_module.cc:1425
Track.h
Provides recob::Track data product.
lar_pandora::LArPandoraHelper::CollectHits
static void CollectHits(const art::Event &evt, const std::string &label, HitVector &hitVector)
Collect the reconstructed Hits from the ART event record.
Definition: LArPandoraHelper.cxx:72
lar_pandora::PFParticleMonitoring::beginJob
void beginJob()
Definition: PFParticleMonitoring_module.cc:334
lar_pandora::PFParticleMonitoring::m_nMatchedHitsV
int m_nMatchedHitsV
Definition: PFParticleMonitoring_module.cc:236
lar_pandora::PFParticleMonitoring::m_pfoVertex
int m_pfoVertex
Definition: PFParticleMonitoring_module.cc:193
lar_pandora::LArPandoraHelper::SelectNeutrinoPFParticles
static void SelectNeutrinoPFParticles(const PFParticleVector &inputParticles, PFParticleVector &outputParticles)
Select reconstructed neutrino particles from a list of all reconstructed particles.
Definition: LArPandoraHelper.cxx:707
lar_pandora::PFParticleMonitoring::m_mcDirY
double m_mcDirY
Definition: PFParticleMonitoring_module.cc:214
lar_pandora::PFParticleMonitoring::m_pfoDirX
double m_pfoDirX
Definition: PFParticleMonitoring_module.cc:200
lar_pandora::PFParticleMonitoring::m_particleLabel
std::string m_particleLabel
Definition: PFParticleMonitoring_module.cc:249
lar_pandora::HitVector
std::vector< art::Ptr< recob::Hit > > HitVector
Definition: LArPandoraHelper.h:58
lar_pandora::PFParticleMonitoring::m_mcEndY
double m_mcEndY
Definition: PFParticleMonitoring_module.cc:211
lar_pandora::PFParticleMonitoring::m_pfoStraightLength
double m_pfoStraightLength
Definition: PFParticleMonitoring_module.cc:204
lar_pandora::PFParticleMonitoring::m_mcNuPdg
int m_mcNuPdg
Definition: PFParticleMonitoring_module.cc:176
PlaneGeo.h
Encapsulate the construction of a single detector plane.
lar_pandora::PFParticleMonitoring::m_pfoNuPdg
int m_pfoNuPdg
Definition: PFParticleMonitoring_module.cc:185
lar_pandora::HitsToPFParticles
std::map< art::Ptr< recob::Hit >, art::Ptr< recob::PFParticle > > HitsToPFParticles
Definition: LArPandoraHelper.h:98
lar_pandora::PFParticleMonitoring::m_mcIsDecay
int m_mcIsDecay
Definition: PFParticleMonitoring_module.cc:181
lar_pandora::PFParticleMonitoring::m_pfoDirZ
double m_pfoDirZ
Definition: PFParticleMonitoring_module.cc:202
lar_pandora::LArPandoraHelper::kUseDaughters
Definition: LArPandoraHelper.h:126
lar_pandora::PFParticleMonitoring::m_mcStraightLength
double m_mcStraightLength
Definition: PFParticleMonitoring_module.cc:218
lar_pandora::PFParticleMonitoring::m_nMCHitsV
int m_nMCHitsV
Definition: PFParticleMonitoring_module.cc:228
e
do i e
Definition: updateVersionList.sh:68
lar_pandora::VertexVector
std::vector< art::Ptr< recob::Vertex > > VertexVector
Definition: LArPandoraHelper.h:62
lar_pandora::PFParticleMonitoring::reconfigure
void reconfigure(fhicl::ParameterSet const &pset)
Definition: PFParticleMonitoring_module.cc:313
lar_pandora::PFParticleMonitoring::m_run
int m_run
Definition: PFParticleMonitoring_module.cc:166
lar_pandora::PFParticleMonitoring::endJob
void endJob()
Definition: PFParticleMonitoring_module.cc:413
lar_pandora::PFParticleMonitoring::MCParticleSet
std::set< art::Ptr< simb::MCParticle > > MCParticleSet
Definition: PFParticleMonitoring_module.cc:45
lar_pandora::PFParticleMonitoring::m_nPfoHits
int m_nPfoHits
Definition: PFParticleMonitoring_module.cc:224
geo::kW
Planes which measure W (third view for Bo, MicroBooNE, etc).
Definition: geo_types.h:131
tfs
art::ServiceHandle< art::TFileService > tfs
Definition: NuVertexChargeTree_module.cc:254
lar_pandora::PFParticleMonitoring::m_mcVtxZ
double m_mcVtxZ
Definition: PFParticleMonitoring_module.cc:209
lar_pandora::PFParticleMonitoring::m_nMCParticles
int m_nMCParticles
Definition: PFParticleMonitoring_module.cc:170
tca::evt
TCEvent evt
Definition: DataStructs.cxx:8
lar_pandora::PFParticleMonitoring::m_index
int m_index
Definition: PFParticleMonitoring_module.cc:168
lar_pandora::T0Vector
std::vector< art::Ptr< anab::T0 > > T0Vector
Definition: LArPandoraHelper.h:72
lar_pandora::PFParticleMonitoring
PFParticleMonitoring class.
Definition: PFParticleMonitoring_module.cc:24
lar_pandora::PFParticleMonitoring::m_nDaughterPfos
int m_nDaughterPfos
Definition: PFParticleMonitoring_module.cc:173
lar_pandora::LArPandoraHelper::CollectT0s
static void CollectT0s(const art::Event &evt, const std::string &label, T0Vector &t0Vector, PFParticlesToT0s &particlesToT0s)
Collect a vector of T0s from the ART event record.
Definition: LArPandoraHelper.cxx:776
lar_pandora::PFParticleMonitoring::m_nMCHitsW
int m_nMCHitsW
Definition: PFParticleMonitoring_module.cc:229
lar_pandora::PFParticleMonitoring::m_addDaughterPFParticles
bool m_addDaughterPFParticles
Definition: PFParticleMonitoring_module.cc:255
lar_pandora::LArPandoraHelper::BuildMCParticleHitMaps
static void BuildMCParticleHitMaps(const art::Event &evt, const HitVector &hitVector, const SimChannelVector &simChannelVector, HitsToTrackIDEs &hitsToTrackIDEs)
Collect the links from reconstructed hits to their true energy deposits.
Definition: LArPandoraHelper.cxx:917
lar_pandora::LArPandoraHelper::CollectMCParticles
static void CollectMCParticles(const art::Event &evt, const std::string &label, MCParticleVector &particleVector)
Collect a vector of MCParticle objects from the ART event record.
Definition: LArPandoraHelper.cxx:831
LArPandoraHelper.h
helper function for LArPandoraInterface producer module
lar_pandora::LArPandoraHelper::GetFinalStateMCParticle
static art::Ptr< simb::MCParticle > GetFinalStateMCParticle(const MCParticleMap &particleMap, const art::Ptr< simb::MCParticle > daughterParticle)
Return the final-state parent particle by navigating up the chain of parent/daughter associations...
Definition: LArPandoraHelper.cxx:1315
lar_pandora::PFParticleMonitoring::m_nRecoWithoutTrueHits
int m_nRecoWithoutTrueHits
Reconstructed hits which don&#39;t belong to any true particle - &quot;missing&quot;.
Definition: PFParticleMonitoring_module.cc:242
lar_pandora::LArPandoraHelper::kAddDaughters
Definition: LArPandoraHelper.h:127
Geometry.h
art framework interface to geometry description
cout
BEGIN_PROLOG could also be cout
Definition: messageservice.fcl:10
lar_pandora::LArPandoraHelper::IsFinalState
static bool IsFinalState(const PFParticleMap &particleMap, const art::Ptr< recob::PFParticle > daughterParticle)
Determine whether a particle has been reconstructed as a final-state particle.
Definition: LArPandoraHelper.cxx:1420
lar_pandora::PFParticleMonitoring::m_nPfoHitsW
int m_nPfoHitsW
Definition: PFParticleMonitoring_module.cc:233
lar_pandora::SpacePointsToHits
std::map< art::Ptr< recob::SpacePoint >, art::Ptr< recob::Hit > > SpacePointsToHits
Definition: LArPandoraHelper.h:89
lar_pandora::LArPandoraHelper::kIgnoreDaughters
Definition: LArPandoraHelper.h:125
lar_pandora::PFParticleMonitoring::m_disableRealDataCheck
bool m_disableRealDataCheck
Whether to check if the input file contains real data before accessing MC information.
Definition: PFParticleMonitoring_module.cc:261
TPCGeo.h
Encapsulate the construction of a single detector plane.
lar_pandora::PFParticleMonitoring::~PFParticleMonitoring
virtual ~PFParticleMonitoring()
Destructor.
Definition: PFParticleMonitoring_module.cc:308
lar_pandora::PFParticleMonitoring::m_mcEndZ
double m_mcEndZ
Definition: PFParticleMonitoring_module.cc:212
lar_pandora::PFParticleMonitoring::m_mcIsPrimary
int m_mcIsPrimary
Definition: PFParticleMonitoring_module.cc:180
lar_pandora::PFParticleMonitoring::m_pfoEndZ
double m_pfoEndZ
Definition: PFParticleMonitoring_module.cc:199
lar_pandora::LArPandoraHelper::BuildMCParticleMap
static void BuildMCParticleMap(const MCParticleVector &particleVector, MCParticleMap &particleMap)
Build particle maps for true particles.
Definition: LArPandoraHelper.cxx:1185
Generated by   doxygen 1.8.5