truepion_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: truepion
3 // File: truepion_module.cc
4 //
5 // Extract protoDUNE useful information, do a first pre-selection and save output to a flat tree
6 //
7 // Some parts are copied from the beam module example
8 //
9 // Georgios Christodoulou - georgios.christodoulou at cern.ch
10 // Heng-Ye Liao modified for his protonanalysis - liao@phys.ksu.edu
11 //Ajib Paudel Modified for the pion truth cross-section studies
12 ////////////////////////////////////////////////////////////////////////
13 
14 #include "art/Framework/Core/EDAnalyzer.h"
15 #include "art/Framework/Core/ModuleMacros.h"
16 #include "art/Framework/Principal/Event.h"
17 #include "art/Framework/Principal/Handle.h"
18 #include "art/Framework/Principal/Run.h"
19 #include "art/Framework/Principal/SubRun.h"
20 #include "art_root_io/TFileService.h"
21 #include "canvas/Utilities/InputTag.h"
22 #include "canvas/Persistency/Common/FindManyP.h"
23 #include "fhiclcpp/ParameterSet.h"
24 #include "messagefacility/MessageLogger/MessageLogger.h"
25 
26 #include "larcore/Geometry/Geometry.h"
27 #include "larcoreobj/SimpleTypesAndConstants/geo_types.h"
28 #include "lardata/DetectorInfoServices/DetectorClocksService.h"
29 #include "larsim/MCCheater/BackTrackerService.h"
30 #include "larsim/MCCheater/ParticleInventoryService.h"
31 #include "lardataobj/RecoBase/Hit.h"
32 #include "lardataobj/RecoBase/Track.h"
33 #include "lardataobj/RecoBase/Shower.h"
34 #include "lardataobj/RecoBase/PFParticle.h"
35 #include "nusimdata/SimulationBase/MCParticle.h"
36 #include "nusimdata/SimulationBase/MCTruth.h"
37 #include "lardataobj/AnalysisBase/CosmicTag.h"
38 #include "lardataobj/AnalysisBase/T0.h"
39 #include "larreco/RecoAlg/TrackMomentumCalculator.h"
40 #include "lardataobj/AnalysisBase/Calorimetry.h"
41 
42 #include "dune/DuneObj/ProtoDUNEBeamEvent.h"
43 //#include "dune/EventGenerator/ProtoDUNEbeamDataProducts/ProtoDUNEbeamsim.h"
44 //#include "dune/EventGenerator/ProtoDUNEbeamDataProducts/ProtoDUNEBeamInstrument.h"
45 
46 #include "protoduneana/Utilities/ProtoDUNETrackUtils.h"
47 #include "protoduneana/Utilities/ProtoDUNEShowerUtils.h"
48 #include "protoduneana/Utilities/ProtoDUNETruthUtils.h"
49 #include "protoduneana/Utilities/ProtoDUNEPFParticleUtils.h"
50 #include "dune/Protodune/singlephase/DataUtils/ProtoDUNEDataUtils.h"
51 
52 #include "larsim/MCCheater/ParticleInventoryService.h"
53 
54 // ROOT includes
55 #include "TTree.h"
56 #include "TFile.h"
57 #include "TString.h"
58 #include "TTimeStamp.h"
59 #include "TH2.h"
60 // C++ Includes
61 #include <stdio.h>
62 #include <stdlib.h>
63 #include <string>
64 #include <sstream>
65 #include <cmath>
66 #include <algorithm>
67 #include <iostream>
68 #include <vector>
69 #include "TComplex.h"
70 #include <fstream>
71 #include "TPaveStats.h"
72 #include <iostream>
73 #include <string>
74 #include "math.h"
75 #include <iterator>
76 #include <map>
77 
78 // Maximum number of beam particles to save
79 const int NMAXDAUGTHERS = 30;
80 double m_pi=0.1395;//GeV/c^2
81 double prim_energy=0.0;
82 //double geant_energy=0;
83 namespace protoana {
84  class truepion;
85 }
86 
87 
88 class protoana::truepion : public art::EDAnalyzer {
89 public:
90 
91  explicit truepion(fhicl::ParameterSet const & p);
92 
93  truepion(truepion const &) = delete;
94  truepion(truepion &&) = delete;
95  truepion & operator = (truepion const &) = delete;
96  truepion & operator = (truepion &&) = delete;
97 
98  virtual void beginJob() override;
99  virtual void endJob() override;
100 
101  // Required functions.
102  void analyze(art::Event const & evt) override;
103 
104 private:
105 
106  // Helper utility functions
107  protoana::ProtoDUNEDataUtils dataUtil;
108  protoana::ProtoDUNEPFParticleUtils pfpUtil;
109  protoana::ProtoDUNETrackUtils trackUtil;
110  protoana::ProtoDUNETruthUtils truthUtil;
111 
112  // Initialise tree variables
113  void Initialise();
114 
115  // Fill cosmics tree
116  void FillCosmicsTree(art::Event const & evt, std::string pfParticleTag);
117 
118  // fcl parameters
119  const art::InputTag fBeamModuleLabel;
120  std::string fCalorimetryTag;
121  std::string fTrackerTag;
122  std::string fShowerTag;
123  std::string fPFParticleTag;
124  std::string fGeneratorTag;
125  bool fVerbose;
126 
127  geo::GeometryCore const * fGeometry;
128 
129  TTree *fPandoraBeam;
130  //TTree *fPandoraCosmics;
131 
132  // Tree variables
133  int fRun;
134  int fSubRun;
135  int fevent;
136  double fTimeStamp;
137 
138 
139  double fbeamtrackMomentum;
140  double fbeamtrackP[3]; //Px/Py/Pz
141  double fbeamtrackEnergy;
142  double fbeamtrackPos[3];
143  double fbeamtrackDir[3];
144  double fbeamtrackTime;
145  int fbeamtrackPdg;
146  int fbeamtrackID;
147 
148  //int NumberBeamTrajectoryPoints;
149  std::vector<float> beamtrk_x;
150  std::vector<float> beamtrk_y;
151  std::vector<float> beamtrk_z;
152  std::vector<float> beamtrk_Px;
153  std::vector<float> beamtrk_Py;
154  std::vector<float> beamtrk_Pz;
155  std::vector<float> beamtrk_Eng;
156 
157 
158 
159  // Reconstructed tracks/showers information
160  double fvertex[3];
161  double fsecvertex[3];
162 
163  int fisprimarytrack;
164  int fisprimaryshower;
165  double fprimaryBDTScore;
166  int fprimaryNHits;
167  double fprimaryTheta;
168  double fprimaryPhi;
169  double fprimaryLength;
170  double fprimaryMomentum;
171  double fprimaryEndMomentum;
172  double fprimaryEndPosition[3];
173  double fprimaryStartPosition[3];
174  double fprimaryEndDirection[3];
175  double fprimaryStartDirection[3];
176  double fprimaryOpeningAngle;
177  int fprimaryShowerBestPlane;
178  double fprimaryShowerEnergy;
179  double fprimaryShowerMIPEnergy;
180  double fprimaryShowerdEdx;
181  double fprimaryMomentumByRangeProton;
182  double fprimaryMomentumByRangeMuon;
183  double fprimaryKineticEnergy[3];
184  double fprimaryRange[3];
185  int fprimaryID;
186 
187 
188 
189 
190  //*********************truth information*******************************************//
191  int fprimary_truth_TrackId;
192  int fprimary_truth_Pdg;
193  int ftruthpdg;
194  double fprimary_truth_StartPosition[4];
195  double fprimary_truth_EndPosition[4];
196  std::string fprimary_truth_EndProcess;
197  std::string truth_last_process;
198  double fprimary_truth_P;
199  double fprimary_truth_Momentum[4];
200  double fprimary_truth_EndMomentum[4];
201  double fprimary_truth_Pt;
202  double fprimary_truth_Mass;
203  double fprimary_truth_Theta;
204  double fprimary_truth_Phi;
205  int fprimary_truth_Process;
206  int fprimary_truth_Isbeammatched;
207  int fprimary_truth_NDaughters;
208  double fprimary_truth_tracklength;
209 
210 
211  //interaction point details
212  std::vector<double> interactionX;
213  std::vector<double> interactionY;
214  std::vector<double> interactionZ;
215  std::vector<std::string> interactionProcesslist;
216  std::vector<double> interactionAngles;
217  std::vector<double> Zintersection;
218  std::vector<double> timeintersection;
219 
220  //calo info
221  std::vector< std::vector<double> > primtrk_dqdx;
222  std::vector< std::vector<double> > primtrk_resrange;
223  std::vector< std::vector<double> > primtrk_dedx;
224  std::vector<double> primtrk_range;
225  std::vector< std::vector<double> > primtrk_hitx;
226  std::vector< std::vector<double> > primtrk_hity;
227  std::vector< std::vector<double> > primtrk_hitz;
228  std::vector< std::vector<double> > primtrk_pitch;
229  std::vector<std::vector<double> > primtrk_truth_Z;
230  std::vector<std::vector<double> > primtrk_truth_Eng;
231  std::vector<std::vector<double> > primtrk_truth_trkide;
232 
233 
234 
235  // Daughters from primary
236  int fNDAUGHTERS;
237  int fisdaughtertrack[NMAXDAUGTHERS];
238  int fisdaughtershower[NMAXDAUGTHERS];
239  int fdaughterNHits[NMAXDAUGTHERS];
240  double fdaughterTheta[NMAXDAUGTHERS];
241  double fdaughterPhi[NMAXDAUGTHERS];
242  double fdaughterLength[NMAXDAUGTHERS];
243  double fdaughterMomentum[NMAXDAUGTHERS];
244  double fdaughterEndMomentum[NMAXDAUGTHERS];
245  double fdaughterEndPosition[NMAXDAUGTHERS][3];
246  double fdaughterStartPosition[NMAXDAUGTHERS][3];
247  double fdaughterEndDirection[NMAXDAUGTHERS][3];
248  double fdaughterStartDirection[NMAXDAUGTHERS][3];
249  double fdaughterOpeningAngle[NMAXDAUGTHERS];
250  double fdaughterShowerEnergy[NMAXDAUGTHERS];
251  double fdaughterShowerMIPEnergy[NMAXDAUGTHERS];
252  double fdaughterShowerdEdx[NMAXDAUGTHERS];
253  int fdaughterShowerBestPlane[NMAXDAUGTHERS];
254  double fdaughterMomentumByRangeProton[NMAXDAUGTHERS];
255  double fdaughterMomentumByRangeMuon[NMAXDAUGTHERS];
256  double fdaughterKineticEnergy[NMAXDAUGTHERS][3];
257  double fdaughterRange[NMAXDAUGTHERS][3];
258  int fdaughterID[NMAXDAUGTHERS];
259  //double fdaughterT0[NMAXDAUGTHERS];
260 
261  int fdaughter_truth_TrackId[NMAXDAUGTHERS];
262  int fdaughter_truth_Pdg[NMAXDAUGTHERS];
263  double fdaughter_truth_StartPosition[NMAXDAUGTHERS][4];
264  double fdaughter_truth_EndPosition[NMAXDAUGTHERS][4];
265  double fdaughter_truth_P[NMAXDAUGTHERS];
266  double fdaughter_truth_Momentum[NMAXDAUGTHERS][4];
267  double fdaughter_truth_EndMomentum[NMAXDAUGTHERS][4];
268  double fdaughter_truth_Pt[NMAXDAUGTHERS];
269  double fdaughter_truth_Mass[NMAXDAUGTHERS];
270  double fdaughter_truth_Theta[NMAXDAUGTHERS];
271  double fdaughter_truth_Phi[NMAXDAUGTHERS];
272  int fdaughter_truth_Process[NMAXDAUGTHERS];
273 
274  double minX = -360.0;
275  double maxX = 360.0;
276  double minY =0.0;
277  double maxY = 600.0;
278  double minZ = 0.0;
279  double maxZ = 695.0;
280 
281 };
282 
283 
284 protoana::truepion::truepion(fhicl::ParameterSet const & p)
285  :
286  EDAnalyzer(p),
287  dataUtil(p.get<fhicl::ParameterSet>("DataUtils")),
288  fBeamModuleLabel(p.get< art::InputTag >("BeamModuleLabel")),
289  fCalorimetryTag(p.get<std::string>("CalorimetryTag")),
290  fTrackerTag(p.get<std::string>("TrackerTag")),
291  fShowerTag(p.get<std::string>("ShowerTag")),
292  fPFParticleTag(p.get<std::string>("PFParticleTag")),
293  fGeneratorTag(p.get<std::string>("GeneratorTag")),
294  fVerbose(p.get<bool>("Verbose"))
295 {
296 
297 }
298 
299 void protoana::truepion::beginJob(){
300 
301  art::ServiceHandle<art::TFileService> tfs;
302 
303  fPandoraBeam = tfs->make<TTree>("PandoraBeam", "Beam events reconstructed with Pandora");
304  fPandoraBeam->Branch("run", &fRun, "run/I");
305  fPandoraBeam->Branch("subrun", &fSubRun, "subrun/I");
306  fPandoraBeam->Branch("event", &fevent, "event/I");
307  fPandoraBeam->Branch("timestamp", &fTimeStamp, "timestamp/D");
308  fPandoraBeam->Branch("beamtrackMomentum", &fbeamtrackMomentum, "beamtrackMomentum/D");
309  fPandoraBeam->Branch("beamtrackP", &fbeamtrackP, "beamtrackP[3]/D");
310  fPandoraBeam->Branch("beamtrackEnergy", &fbeamtrackEnergy, "beamtrackEnergy/D");
311  fPandoraBeam->Branch("beamtrackPos", &fbeamtrackPos, "beamtrackPos[3]/D");
312  fPandoraBeam->Branch("beamtrackDir", &fbeamtrackDir, "beamtrackDir[3]/D");
313  fPandoraBeam->Branch("beamtrackTime", &fbeamtrackTime, "beamtrackTime/D");
314  fPandoraBeam->Branch("beamtrackPdg", &fbeamtrackPdg, "beamtrackPdg/I");
315  fPandoraBeam->Branch("beamtrackID", &fbeamtrackID, "beamtrackID/I");
316 
317  //fPandoraBeam->Branch("NumberBeamTrajectoryPoints", &NumberBeamTrajectoryPoints, "NumberBeamTrajectoryPoints/I");
318  fPandoraBeam->Branch("beamtrk_x",&beamtrk_x);
319  fPandoraBeam->Branch("beamtrk_y",&beamtrk_y);
320  fPandoraBeam->Branch("beamtrk_z",&beamtrk_z);
321  fPandoraBeam->Branch("beamtrk_Px",&beamtrk_Px);
322  fPandoraBeam->Branch("beamtrk_Py",&beamtrk_Py);
323  fPandoraBeam->Branch("beamtrk_Pz",&beamtrk_Pz);
324  fPandoraBeam->Branch("beamtrk_Eng",&beamtrk_Eng);
325 
326 
327  fPandoraBeam->Branch("vertex", &fvertex, "vertex[3]/D");
328  fPandoraBeam->Branch("secvertex", &fsecvertex, "secvertex[3]/D");
329  fPandoraBeam->Branch("isprimarytrack", &fisprimarytrack, "isprimarytrack/I");
330  fPandoraBeam->Branch("isprimaryshower", &fisprimaryshower, "isprimaryshower/I");
331  fPandoraBeam->Branch("primaryBDTScore", &fprimaryBDTScore, "primaryBDTScore/D");
332  fPandoraBeam->Branch("primaryNHits", &fprimaryNHits, "fprimaryNHits/I");
333  fPandoraBeam->Branch("primaryTheta", &fprimaryTheta, "primaryTheta/D");
334  fPandoraBeam->Branch("primaryPhi", &fprimaryPhi, "primaryPhi/D");
335  fPandoraBeam->Branch("primaryLength", &fprimaryLength, "primaryLength/D");
336  fPandoraBeam->Branch("primaryMomentum", &fprimaryMomentum, "primaryMomentum/D");
337  fPandoraBeam->Branch("primaryEndMomentum", &fprimaryEndMomentum, "primaryEndMomentum/D");
338  fPandoraBeam->Branch("primaryEndPosition", &fprimaryEndPosition, "primaryEndPosition[3]/D");
339  fPandoraBeam->Branch("primaryStartPosition", &fprimaryStartPosition, "primaryStartPosition[3]/D");
340  fPandoraBeam->Branch("primaryEndDirection", &fprimaryEndDirection, "primaryEndDirection[3]/D");
341  fPandoraBeam->Branch("primaryStartDirection", &fprimaryStartDirection, "primaryStartDirection[3]/D");
342  fPandoraBeam->Branch("primaryOpeningAngle", &fprimaryOpeningAngle, "primaryOpeningAngle/D");
343  fPandoraBeam->Branch("primaryID", &fprimaryID, "primaryID/I");
344  fPandoraBeam->Branch("primaryShowerBestPlane", &fprimaryShowerBestPlane, "primaryShowerBestPlane/I");
345  fPandoraBeam->Branch("primaryShowerEnergy", &fprimaryShowerEnergy, "primaryShowerEnergy/I");
346  fPandoraBeam->Branch("primaryShowerMIPEnergy", &fprimaryShowerMIPEnergy, "primaryShowerMIPEnergy/I");
347  fPandoraBeam->Branch("primaryShowerdEdx", &fprimaryShowerdEdx, "primaryShowerdEdx/I");
348  fPandoraBeam->Branch("primaryMomentumByRangeProton", &fprimaryMomentumByRangeProton, "primaryMomentumByRangeProton/D");
349  fPandoraBeam->Branch("primaryMomentumByRangeMuon", &fprimaryMomentumByRangeMuon, "primaryMomentumByRangeMuon/D");
350  fPandoraBeam->Branch("primaryKineticEnergy", &fprimaryKineticEnergy, "primaryKineticEnergy[3]/D");
351  fPandoraBeam->Branch("primaryRange", &fprimaryRange, "primaryRange[3]/D");
352 
353 
354  fPandoraBeam->Branch("primary_truth_TrackId", &fprimary_truth_TrackId, "primary_truth_TrackId/I");
355  fPandoraBeam->Branch("primary_truth_Pdg", &fprimary_truth_Pdg, "primary_truth_Pdg/I");
356  fPandoraBeam->Branch("truthpdg", &ftruthpdg, "truthpdg/I");
357 
358 
359  fPandoraBeam->Branch("primary_truth_StartPosition", &fprimary_truth_StartPosition, "primary_truth_StartPosition[4]/D");
360  fPandoraBeam->Branch("primary_truth_EndPosition", &fprimary_truth_EndPosition, "primary_truth_EndPosition[4]/D");
361  fPandoraBeam->Branch("primary_truth_Momentum", &fprimary_truth_Momentum, "primary_truth_Momentum[4]/D");
362  fPandoraBeam->Branch("primary_truth_EndMomentum", &fprimary_truth_EndMomentum, "primary_truth_EndMomentum[4]/D");
363  fPandoraBeam->Branch("primary_truth_P", &fprimary_truth_P, "primary_truth_P/D");
364  fPandoraBeam->Branch("primary_truth_Pt", &fprimary_truth_Pt, "primary_truth_Pt/D");
365  fPandoraBeam->Branch("primary_truth_Mass", &fprimary_truth_Mass, "primary_truth_Mass/D");
366  fPandoraBeam->Branch("primary_truth_Theta", &fprimary_truth_Theta, "primary_truth_Theta/D");
367  fPandoraBeam->Branch("primary_truth_Phi", &fprimary_truth_Phi, "primary_truth_Phi/D");
368  fPandoraBeam->Branch("primary_truth_Process", &fprimary_truth_Process, "primary_truth_Process/I");
369  fPandoraBeam->Branch("primary_truth_Isbeammatched", &fprimary_truth_Isbeammatched, "primary_truth_Isbeammatched/I");
370  fPandoraBeam->Branch("primary_truth_NDaughters", &fprimary_truth_NDaughters, "primary_truth_NDaughters/I");
371  fPandoraBeam->Branch("primary_truth_EndProcess", &fprimary_truth_EndProcess);
372  fPandoraBeam->Branch("truth_last_process", &truth_last_process);
373  fPandoraBeam->Branch("primary_truth_tracklength", &fprimary_truth_tracklength, "primary_truth_tracklength/D");
374 
375  fPandoraBeam->Branch("interactionX",&interactionX);
376  fPandoraBeam->Branch("interactionY",&interactionY);
377  fPandoraBeam->Branch("interactionZ",&interactionZ);
378  fPandoraBeam->Branch("interactionProcesslist",&interactionProcesslist);
379  fPandoraBeam->Branch("interactionAngles",&interactionAngles);
380  fPandoraBeam->Branch("Zintersection",&Zintersection);
381  fPandoraBeam->Branch("timeintersection",&timeintersection);
382 
383 
384  fPandoraBeam->Branch("NDAUGHTERS", &fNDAUGHTERS, "NDAUGHTERS/I");
385  fPandoraBeam->Branch("isdaughtertrack", &fisdaughtertrack, "isdaughtertrack[NDAUGHTERS]/I");
386  fPandoraBeam->Branch("isdaughtershower", &fisdaughtershower, "isdaughtershower[NDAUGHTERS]/I");
387  fPandoraBeam->Branch("daughterNHits", &fdaughterNHits, "daughterNHits[NDAUGHTERS]/I");
388  fPandoraBeam->Branch("daughterTheta", &fdaughterTheta, "daughterTheta[NDAUGHTERS]/D");
389  fPandoraBeam->Branch("daughterPhi", &fdaughterPhi, "daughterPhi[NDAUGHTERS]/D");
390  fPandoraBeam->Branch("daughterLength", &fdaughterLength, "daughterLength[NDAUGHTERS]/D");
391  fPandoraBeam->Branch("daughterMomentum", &fdaughterMomentum, "daughterMomentum[NDAUGHTERS]/D");
392  fPandoraBeam->Branch("daughterEndMomentum", &fdaughterEndMomentum, "daughterEndMomentum[NDAUGHTERS]/D");
393  fPandoraBeam->Branch("daughterEndPosition", &fdaughterEndPosition, "daughterEndPosition[NDAUGHTERS][3]/D");
394  fPandoraBeam->Branch("daughterStartPosition", &fdaughterStartPosition, "daughterStartPosition[NDAUGHTERS][3]/D");
395  fPandoraBeam->Branch("daughterStartDirection", &fdaughterStartDirection, "daughterStartDirection[NDAUGHTERS][3]/D");
396  fPandoraBeam->Branch("daughterEndDirection", &fdaughterEndDirection, "daughterEndDirection[NDAUGHTERS][3]/D");
397  fPandoraBeam->Branch("daughterOpeningAngle", &fdaughterOpeningAngle, "daughterOpeningAngle[NDAUGHTERS]/D");
398  fPandoraBeam->Branch("daughterShowerBestPlane", &fdaughterShowerBestPlane, "daughterShowerBestPlane[NDAUGHTERS]/D");
399  fPandoraBeam->Branch("daughterShowerEnergy", &fdaughterShowerEnergy, "daughterShowerEnergy[NDAUGHTERS]/D");
400  fPandoraBeam->Branch("daughterShowerMIPEnergy", &fdaughterShowerMIPEnergy, "daughterShowerMIPEnergy[NDAUGHTERS]/D");
401  fPandoraBeam->Branch("daughterShowerdEdx", &fdaughterShowerdEdx, "daughterShowerdEdx[NDAUGHTERS]/D");
402  fPandoraBeam->Branch("daughterMomentumByRangeProton", &fdaughterMomentumByRangeProton,"daughterMomentumByRangeProton[NDAUGHTERS]/D");
403  fPandoraBeam->Branch("daughterMomentumByRangeMuon", &fdaughterMomentumByRangeMuon, "daughterMomentumByRangeMuon[NDAUGHTERS]/D");
404  fPandoraBeam->Branch("daughterKineticEnergy", &fdaughterKineticEnergy, "daughterKineticEnergy[NDAUGHTERS][3]/D");
405  fPandoraBeam->Branch("daughterRange", &fdaughterRange, "daughterRange[NDAUGHTERS][3]/D");
406  fPandoraBeam->Branch("daughterID", &fdaughterID, "daughterID[NDAUGHTERS]/I");
407  // fPandoraBeam->Branch("daughterT0", &fdaughterT0, "daughterT0[NDAUGHTERS]/D");
408 
409  fPandoraBeam->Branch("daughter_truth_TrackId", &fdaughter_truth_TrackId, "daughter_truth_TrackId[NDAUGHTERS]/I");
410  fPandoraBeam->Branch("daughter_truth_Pdg", &fdaughter_truth_Pdg, "daughter_truth_Pdg[NDAUGHTERS]/I");
411  fPandoraBeam->Branch("daughter_truth_StartPosition", &fdaughter_truth_StartPosition, "daughter_truth_StartPosition[NDAUGHTERS][4]/D");
412  fPandoraBeam->Branch("daughter_truth_EndPosition", &fdaughter_truth_EndPosition, "daughter_truth_EndPosition[NDAUGHTERS][4]/D");
413  fPandoraBeam->Branch("daughter_truth_Momentum", &fdaughter_truth_Momentum, "daughter_truth_Momentum[NDAUGHTERS][4]/D");
414  fPandoraBeam->Branch("daughter_truth_EndMomentum", &fdaughter_truth_EndMomentum, "daughter_truth_EndMomentum[NDAUGHTERS][4]/D");
415  fPandoraBeam->Branch("daughter_truth_P", &fdaughter_truth_P, "daughter_truth_P[NDAUGHTERS]/D");
416  fPandoraBeam->Branch("daughter_truth_Pt", &fdaughter_truth_Pt, "daughter_truth_Pt[NDAUGHTERS]/D");
417  fPandoraBeam->Branch("daughter_truth_Mass", &fdaughter_truth_Mass, "daughter_truth_Mass[NDAUGHTERS]/D");
418  fPandoraBeam->Branch("daughter_truth_Theta", &fdaughter_truth_Theta, "daughter_truth_Theta[NDAUGHTERS]/D");
419  fPandoraBeam->Branch("daughter_truth_Phi", &fdaughter_truth_Phi, "daughter_truth_Phi[NDAUGHTERS]/D");
420  fPandoraBeam->Branch("daughter_truth_Process", &fdaughter_truth_Process, "daughter_truth_Process[NDAUGHTERS]/I");
421 
422 
423  // hKE_truth_reco = tfs->make<TH2D>("hKE_truth_reco","KE truth vs reco;true KE in MeV;reconstructed KE in MeV" , 3000, -100, 1400 ,3000, -100, 1400);
424  // hEndZ_truth_reco = tfs->make<TH2D>("hEndZ_truth_reco","EndZ truth vs reco;true EndZ (cm);reconstructed EndZ(cm)" , 695, 0, 695 ,695,0,695);
425 
426  fPandoraBeam->Branch("primtrk_dqdx",&primtrk_dqdx);
427  fPandoraBeam->Branch("primtrk_dedx",&primtrk_dedx);
428  fPandoraBeam->Branch("primtrk_resrange",&primtrk_resrange);
429  fPandoraBeam->Branch("primtrk_range",&primtrk_range);
430  fPandoraBeam->Branch("primtrk_hitx",&primtrk_hitx);
431  fPandoraBeam->Branch("primtrk_hity",&primtrk_hity);
432  fPandoraBeam->Branch("primtrk_hitz",&primtrk_hitz);
433  fPandoraBeam->Branch("primtrk_pitch",&primtrk_pitch);
434  ///////////////////////////////////////////////////
435  fPandoraBeam->Branch("primtrk_truth_Z",&primtrk_truth_Z);//primary track true Z positions
436  fPandoraBeam->Branch("primtrk_truth_Eng",&primtrk_truth_Eng);//primary track true Energy deposited for each Z position
437  fPandoraBeam->Branch("primtrk_truth_trkide",&primtrk_truth_trkide);//primary track true Energy deposited for each Z position
438 }//begin job
439 
440 void protoana::truepion::analyze(art::Event const & evt){
441 
442  // Initialise tree parameters
443  Initialise();
444 
445  art::ServiceHandle<cheat::ParticleInventoryService> pi_serv;
446  fGeometry = &*(art::ServiceHandle<geo::Geometry>());
447  // const sim::ParticleList& plist=pi_serv->ParticleList();
448 
449  int beamid=-9999;
450  int truthid=-999;
451 
452  fRun = evt.run();
453  fSubRun = evt.subRun();
454  fevent = evt.id().event();
455  art::Timestamp ts = evt.time();
456  if (ts.timeHigh() == 0){
457  TTimeStamp ts2(ts.timeLow());
458  fTimeStamp = ts2.AsDouble();
459  }
460  else{
461  TTimeStamp ts2(ts.timeHigh(), ts.timeLow());
462  fTimeStamp = ts2.AsDouble();
463  }
464 
465  bool beamTriggerEvent = false;
466  // If this event is MC then we can check what the true beam particle is
467  auto mcTruths = evt.getValidHandle<std::vector<simb::MCTruth>>(fGeneratorTag);
468  if(!evt.isRealData()){
469  // mcTruths is basically a pointer to an std::vector of simb::MCTruth objects. There should only be one
470  // of these, so we pass the first element into the function to get the good particle
471  const simb::MCParticle* geantGoodParticle = truthUtil.GetGeantGoodParticle((*mcTruths)[0],evt);
472  //std::cout<<"geantGoodParticle "<<geantGoodParticle.size()<<std::endl;
473  if(geantGoodParticle != 0x0){
474  std::cout<<"geant good particle loop "<<std::endl;
475  std::cout << "Found GEANT particle corresponding to the good particle with pdg = " << geantGoodParticle->PdgCode()
476  << " , track id = " << geantGoodParticle->TrackId()
477  << " , Vx/Vy/Vz = " << geantGoodParticle->Vx() << "/"<< geantGoodParticle->Vy() << "/" << geantGoodParticle->Vz()
478 
479  << std::endl;
480 
481  std::vector<double> tmp_primtrk_truth_Z;
482  std::vector<double> tmp_primtrk_truth_Eng;
483  std::vector<double> tmp_primtrk_truth_trkide;
484 
485  beamTriggerEvent = true;
486  fbeamtrackPos[0] = geantGoodParticle->Vx();
487  fbeamtrackPos[1] = geantGoodParticle->Vy();
488  fbeamtrackPos[2] = geantGoodParticle->Vz();
489  fbeamtrackMomentum = geantGoodParticle->P();
490  fbeamtrackP[0] = geantGoodParticle->Px();
491  fbeamtrackP[1] = geantGoodParticle->Py();
492  fbeamtrackP[2] = geantGoodParticle->Pz();
493  // fbeamtrackEnergy = geantGoodParticle->E();
494  fbeamtrackPdg = geantGoodParticle->PdgCode();
495  fbeamtrackTime = geantGoodParticle->T();
496  fbeamtrackID = geantGoodParticle->TrackId();
497  fprimary_truth_TrackId = geantGoodParticle->TrackId();
498  fprimary_truth_Pdg = geantGoodParticle->PdgCode();
499  beamid = geantGoodParticle->TrackId();
500  fprimary_truth_StartPosition[3] = geantGoodParticle->T();
501  fprimary_truth_EndPosition[3] = geantGoodParticle->EndT();
502  fprimary_truth_P = geantGoodParticle->P();
503  fprimary_truth_Momentum[3] = geantGoodParticle->E();
504  fprimary_truth_Pt = geantGoodParticle->Pt();
505  fprimary_truth_Mass = geantGoodParticle->Mass();
506  fprimary_truth_EndMomentum[0] = geantGoodParticle->EndPx();
507  fprimary_truth_EndMomentum[1] = geantGoodParticle->EndPy();
508  fprimary_truth_EndMomentum[2] = geantGoodParticle->EndPz();
509  fprimary_truth_EndMomentum[3] = geantGoodParticle->EndE();
510  fprimary_truth_Theta = geantGoodParticle->Momentum().Theta();
511  fprimary_truth_Phi = geantGoodParticle->Momentum().Phi();
512  fprimary_truth_NDaughters = geantGoodParticle->NumberDaughters();
513 
514  ////////////////////////////////////////
515  truth_last_process=geantGoodParticle->EndProcess();
516  ////here is another new parameter added
517 
518  fprimary_truth_Process = int(geantGoodParticle->Trajectory().ProcessToKey(geantGoodParticle->Process()));
519  prim_energy=0;
520  for(size_t i_s=0; i_s < geantGoodParticle->NumberTrajectoryPoints(); i_s++){ //loop over beam tracks
521  // if(geantGoodParticle->Position(i_s).Z()>0) break;
522  beamtrk_x.push_back(geantGoodParticle->Position(i_s).X());
523  beamtrk_y.push_back(geantGoodParticle->Position(i_s).Y());
524  beamtrk_z.push_back(geantGoodParticle->Position(i_s).Z());
525 
526  beamtrk_Px.push_back(geantGoodParticle->Momentum(i_s).X());
527  beamtrk_Py.push_back(geantGoodParticle->Momentum(i_s).Y());
528  beamtrk_Pz.push_back(geantGoodParticle->Momentum(i_s).Z());
529 
530  beamtrk_Eng.push_back(geantGoodParticle->Momentum(i_s).E()-geantGoodParticle->Mass());
531  if(geantGoodParticle->Position(i_s).Z()<0) prim_energy=1000*(geantGoodParticle->Momentum(i_s).E()-geantGoodParticle->Mass());
532  if(geantGoodParticle->Position(i_s).Z()<0) fbeamtrackEnergy=prim_energy; //correct energy at the beginning of track
533  } //loop over beam trks
534 
535  //new section
536  art::ServiceHandle<cheat::BackTrackerService> bt_serv;
537 
538 
539  art::ServiceHandle<geo::Geometry> geom;
540  simb::MCTrajectory truetraj=geantGoodParticle->Trajectory();
541  auto thisTrajectoryProcessMap1 = truetraj.TrajectoryProcesses();
542  if (thisTrajectoryProcessMap1.size()){
543  for(auto const& couple: thisTrajectoryProcessMap1){
544  // int_label=truetraj.KeyToProcess(couple.second);
545  fprimary_truth_EndPosition[0]=((truetraj.at(couple.first)).first).X();
546  fprimary_truth_EndPosition[1]=((truetraj.at(couple.first)).first).Y();
547  fprimary_truth_EndPosition[2]=((truetraj.at(couple.first)).first).Z();
548  fprimary_truth_EndProcess=truetraj.KeyToProcess(couple.second);
549  fprimary_truth_Momentum[0]=((truetraj.at(couple.first)).second).X();
550  fprimary_truth_Momentum[1]= ((truetraj.at(couple.first)).second).Y();
551  fprimary_truth_Momentum[2]=((truetraj.at(couple.first)).second).Z();
552  break;
553  }
554  }
555 
556  ////saving the complete information of all the interactions
557  std::cout<<"interaction map size "<<thisTrajectoryProcessMap1.size()<<std::endl;
558  if (thisTrajectoryProcessMap1.size()){
559  for(auto const& couple1: thisTrajectoryProcessMap1){
560 
561  if ((truetraj.KeyToProcess(couple1.second)).find("CoulombScat")!= std::string::npos) continue;
562  // Let's check if the interaction is in the the TPC
563  auto interactionPos4D = (truetraj.at(couple1.first)).first ;
564  if (interactionPos4D.Z() < minZ || interactionPos4D.Z() > maxZ ) continue;
565  else if (interactionPos4D.X() < minX || interactionPos4D.X() > maxX ) continue;
566  else if (interactionPos4D.Y() < minY || interactionPos4D.Y() > maxY ) continue;
567  interactionX.push_back(((truetraj.at(couple1.first)).first).X());
568  interactionY.push_back(((truetraj.at(couple1.first)).first).Y());
569  interactionZ.push_back(((truetraj.at(couple1.first)).first).Z());
570  interactionProcesslist.push_back(truetraj.KeyToProcess(couple1.second));
571  std::cout<<"number of interactions "<<thisTrajectoryProcessMap1.size()<<std::endl;
572  std::cout<<"int X, Y, Z and process "<<((truetraj.at(couple1.first)).first).X()<<" "<<((truetraj.at(couple1.first)).first).Y()<<" "<<((truetraj.at(couple1.first)).first).Z()<<" "<<truetraj.KeyToProcess(couple1.second)<<std::endl;
573  ///get the interaction angle here
574  double interactionAngle = 999999.; // This needs to be changed
575  //--------------------- Int Angle ---------------------------
576  // Try to retreive the interaction angle
577  auto prevInteractionPos4D = (truetraj.at(couple1.first-1)).first ;
578  auto prevInteractionPos3D = prevInteractionPos4D.Vect() ;
579  auto interactionPos3D = interactionPos4D.Vect() ;
580  auto distanceBtwPoint = interactionPos3D - prevInteractionPos3D;
581  //Let's try to see if the next point exists
582  if (truetraj.size() > couple1.first + 1) {
583  // The particle doesn't die. No need to check for anything else.
584  auto nextInteractionPos4D = (truetraj.at(couple1.first+1)).first ;
585  auto nextInteractionPos3D = nextInteractionPos4D.Vect() ;
586  auto distanceBtwPointNext = nextInteractionPos3D - interactionPos3D;
587  interactionAngles.push_back(TMath::ACos(distanceBtwPointNext.Dot(distanceBtwPoint)/(distanceBtwPointNext.Mag()*distanceBtwPoint.Mag() ) ));
588  continue;
589  }
590  interactionAngles.push_back(interactionAngle);
591  }
592  }
593 
594  geo::View_t view = geom->View(2);
595  auto simIDE_prim=bt_serv->TrackIdToSimIDEs_Ps(geantGoodParticle->TrackId(),view);
596  std::map<double, sim::IDE> orderedSimIDE;
597  for (auto& ide : simIDE_prim) orderedSimIDE[ide->z]= *ide;
598  auto inTPCPoint = truetraj.begin();
599  auto Momentum0 = inTPCPoint->second;
600  auto old_iter = orderedSimIDE.begin();
601  double tlen=0.0;
602  double xi=0.0;double yi=0.0;double zi=0.0;
603  int count=0;
604 
605 
606  for ( auto iter= orderedSimIDE.begin(); iter!= orderedSimIDE.end(); iter++,old_iter++){
607  auto currentIde = iter->second;
608  if(currentIde.z<minZ) continue;
609  else if (currentIde.x < minX || currentIde.x > maxX ) continue;
610  else if (currentIde.y < minY || currentIde.y > maxY ) continue;
611  tmp_primtrk_truth_Z.push_back(currentIde.z);
612  tmp_primtrk_truth_Eng.push_back(currentIde.energy);
613  tmp_primtrk_truth_trkide.push_back(currentIde.trackID);
614  if(count==0){
615  fprimary_truth_StartPosition[0] = currentIde.x;
616  fprimary_truth_StartPosition[1] = currentIde.y;
617  fprimary_truth_StartPosition[2] = currentIde.z;
618  }
619  if(currentIde.trackID>=0){
620  if(count>0){
621  tlen=tlen+TMath::Sqrt(std::pow(currentIde.x-xi,2)+std::pow(currentIde.y-yi,2)+std::pow(currentIde.z-zi,2));
622  }//if count
623  xi=currentIde.x;yi=currentIde.y;zi=currentIde.z;
624  count++;
625  }//trackid>0 loop
626  }// iter loop
627  fprimary_truth_tracklength=tlen;
628  primtrk_truth_Z.push_back(tmp_primtrk_truth_Z);
629  primtrk_truth_Eng.push_back(tmp_primtrk_truth_Eng);
630  primtrk_truth_trkide.push_back(tmp_primtrk_truth_trkide);
631 
632  tmp_primtrk_truth_Z.clear();
633  tmp_primtrk_truth_Eng.clear();
634  tmp_primtrk_truth_trkide.clear();
635  //new section
636 
637  }// geantGoodParticle
638  }//is real data loop
639 
640  /*
641  // Now we want to access the output from Pandora. This comes in the form of particle flow objects (recob::PFParticle).
642  // The primary PFParticles are those we want to consider and these PFParticles then have a hierarchy of daughters that
643  // describe the whole interaction of a given primary particle
644  //
645  // / daughter track
646  // /
647  // primary track /
648  // ---------------- ---- daughter track
649  // \
650  // /\-
651  // /\\-- daughter shower
652  //
653  // The above primary PFParticle will have links to three daughter particles, two track-like and one shower-like
654  */
655 
656  // Track momentum algorithm calculates momentum based on track range
657  trkf::TrackMomentumCalculator trmom;
658  //trmom.SetMinLength(100);
659 
660  // Get all of the PFParticles, by default from the "pandora" product
661  auto recoParticles = evt.getValidHandle<std::vector<recob::PFParticle>>(fPFParticleTag);
662  std::cout << "All primary pfParticles = " << pfpUtil.GetNumberPrimaryPFParticle(evt,fPFParticleTag) << std::endl;
663 
664  // We'd like to find the beam particle. Pandora tries to do this for us, so let's use the PFParticle utility
665  // to look for it. Pandora reconstructs slices containing one (or sometimes more) primary PFParticles. These
666  // are tagged as either beam or cosmic for ProtoDUNE. This function automatically considers only those
667  // PFParticles considered as primary
668  // std::vector<recob::PFParticle*> pfParticles = pfpUtil.GetPFParticlesFromBeamSlice(evt,fPFParticleTag);
669 
670 
671 
672  auto pfParticles = pfpUtil.GetPFParticlesFromBeamSlice(evt,fPFParticleTag);
673  //cluster information
674 
675  std::vector<art::Ptr<recob::Track> > tracklist;
676  auto trackListHandle = evt.getHandle< std::vector<recob::Track> >("pandoraTrack");
677  if(trackListHandle) art::fill_ptr_vector(tracklist, trackListHandle);
678 
679  std::vector<art::Ptr<recob::PFParticle> > pfplist;
680  auto PFPListHandle = evt.getHandle< std::vector<recob::PFParticle> >("pandora");
681  if(PFPListHandle) art::fill_ptr_vector(pfplist, PFPListHandle);
682 
683  // get information about the hits
684  std::vector<art::Ptr<recob::Cluster>> clusterlist;
685  auto clusterListHandle = evt.getHandle< std::vector<recob::Cluster> >("pandora");
686  if(clusterListHandle)
687  art::fill_ptr_vector(clusterlist, clusterListHandle);
688 
689  art::FindManyP<recob::Cluster> fmcp(PFPListHandle,evt,"pandora");
690  art::FindManyP<recob::Track> pftrack(PFPListHandle,evt,"pandoraTrack");
691 
692  std::cout<<"number of pfp_particles "<<pfplist.size()<<std::endl;
693  std::cout<<" size of pfParticles "<<pfParticles.size()<<std::endl;
694 
695  /* for(size_t p1=0;p1<pfplist.size();p1++){
696  std::vector<art::Ptr<recob::Track>> trk=pftrack.at(p1);
697  if(trk.size()) std::cout<<" trk key "<<trk[0].key()<<std::endl;
698 
699  std::vector<art::Ptr<recob::Cluster>> allClusters=fmcp.at(p1);
700  std::cout<<"cluster size for each particle "<<allClusters.size()<<std::endl;
701  for(size_t c1=0;c1<allClusters.size();c1++){
702  std::cout<<" cluster ID "<<allClusters[c1]->ID();
703  std::cout<<" plane number "<<allClusters[c1]->Plane().Plane;
704  std::cout<<" TPC number "<<allClusters[c1]->Plane().TPC;
705  std::cout<<" start wire "<<allClusters[c1]->StartWire();
706  std::cout<<" end wire "<<allClusters[c1]->EndWire();
707  std::cout<<" start tick "<<allClusters[c1]->StartTick();
708  std::cout<<" end tick "<<allClusters[c1]->EndTick();
709  }
710  }*/
711 
712 
713 
714 
715  // We can now look at these particles
716  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService>()->DataFor(evt);
717  for(const recob::PFParticle* particle : pfParticles){
718  // Pandora's BDT beam-cosmic score
719  fprimaryBDTScore = (double)pfpUtil.GetBeamCosmicScore(*particle,evt,fPFParticleTag);
720  // NHits associated with this pfParticle
721  fprimaryNHits = (pfpUtil.GetPFParticleHits(*particle,evt,fPFParticleTag)).size();
722  // of this particle might be more helpful. These return null pointers if not track-like / shower-like
723  const recob::Track* thisTrack = pfpUtil.GetPFParticleTrack(*particle,evt,fPFParticleTag,fTrackerTag);
724  const recob::Shower* thisShower = pfpUtil.GetPFParticleShower(*particle,evt,fPFParticleTag,fShowerTag);
725  /////new line added here
726  // std::vector<art::Ptr<recob::Cluster>> allClusters=fmcp.at(particle);
727  // std::cout<<allClusters.size();
728  /* for(size_t c1=0;c1<allClusters.size();c1++){
729 
730 
731  }*/
732 
733 
734  if(thisTrack != 0x0){
735  // Get the true mc particle
736  const simb::MCParticle* mcparticle = truthUtil.GetMCParticleFromRecoTrack(clockData, *thisTrack, evt, fTrackerTag);
737  if(mcparticle!=0x0){
738  std::cout<<"ftruth pdg "<<mcparticle->PdgCode()<<std::endl;
739  ftruthpdg=mcparticle->PdgCode();
740 
741  truthid=mcparticle->TrackId();
742  fprimary_truth_Isbeammatched=0;
743  if(beamid==truthid) fprimary_truth_Isbeammatched=1;
744  }//mcparticle loop
745 
746 
747 
748 
749 
750  fisprimarytrack = 1;
751  fisprimaryshower = 0;
752  fprimaryID = thisTrack->ID();
753  std::cout<<"this Track track ID "<<thisTrack->ID()<<std::endl;
754  fprimaryTheta = thisTrack->Theta();
755  fprimaryPhi = thisTrack->Phi();
756  fprimaryLength = thisTrack->Length();
757  fprimaryMomentum = thisTrack->StartMomentum();
758  fprimaryEndMomentum = thisTrack->EndMomentum();
759  fprimaryEndPosition[0] = thisTrack->End().X();
760  fprimaryEndPosition[1] = thisTrack->End().Y();
761  fprimaryEndPosition[2] = thisTrack->End().Z();
762  fprimaryStartPosition[0] = thisTrack->Start().X();
763  fprimaryStartPosition[1] = thisTrack->Start().Y();
764  fprimaryStartPosition[2] = thisTrack->Start().Z();
765  fprimaryEndDirection[0] = thisTrack->Trajectory().EndDirection().X();
766  fprimaryEndDirection[1] = thisTrack->Trajectory().EndDirection().Y();
767  fprimaryEndDirection[2] = thisTrack->Trajectory().EndDirection().Z();
768  fprimaryStartDirection[0] = thisTrack->Trajectory().StartDirection().X();
769  fprimaryStartDirection[1] = thisTrack->Trajectory().StartDirection().Y();
770  fprimaryStartDirection[2] = thisTrack->Trajectory().StartDirection().Z();
771  fprimaryMomentumByRangeMuon = trmom.GetTrackMomentum(thisTrack->Length(),13);
772  fprimaryMomentumByRangeProton = trmom.GetTrackMomentum(thisTrack->Length(),2212);
773  std::vector<anab::Calorimetry> calovector = trackUtil.GetRecoTrackCalorimetry(*thisTrack, evt, fTrackerTag, fCalorimetryTag);
774  if(calovector.size() != 3)
775  std::cerr << "WARNING::Calorimetry vector size for primary is = " << calovector.size() << std::endl;
776  std::vector<double> tmp_primtrk_dqdx;
777  std::vector<double> tmp_primtrk_resrange;
778  std::vector<double> tmp_primtrk_dedx;
779  std::vector<double> tmp_primtrk_hitx;
780  std::vector<double> tmp_primtrk_hity;
781  std::vector<double> tmp_primtrk_hitz;
782  std::vector<double> tmp_primtrk_pitch;
783  for (auto & calo : calovector) {
784  if (calo.PlaneID().Plane == 2){ //only collection plane
785  primtrk_range.push_back(calo.Range());
786  for (size_t ihit = 0; ihit < calo.dQdx().size(); ++ihit){ //loop over hits
787  tmp_primtrk_dqdx.push_back(calo.dQdx()[ihit]);
788  tmp_primtrk_resrange.push_back(calo.ResidualRange()[ihit]);
789  tmp_primtrk_dedx.push_back(calo.dEdx()[ihit]);
790  tmp_primtrk_pitch.push_back(calo.TrkPitchVec()[ihit]);
791 
792  const auto &primtrk_pos=(calo.XYZ())[ihit];
793  tmp_primtrk_hitx.push_back(primtrk_pos.X());
794  tmp_primtrk_hity.push_back(primtrk_pos.Y());
795  tmp_primtrk_hitz.push_back(primtrk_pos.Z());
796  } //loop over hits
797  } //only collection plane
798  }//calovector
799 
800  if(tmp_primtrk_dqdx.size()!=0){
801  if(tmp_primtrk_hitz[0]>tmp_primtrk_hitz[tmp_primtrk_hitz.size()-1]){
802  std::reverse(tmp_primtrk_hitz.begin(),tmp_primtrk_hitz.end());
803  std::reverse(tmp_primtrk_pitch.begin(),tmp_primtrk_pitch.end());
804  std::reverse(tmp_primtrk_dedx.begin(),tmp_primtrk_dedx.end());
805  std::reverse(tmp_primtrk_dqdx.begin(),tmp_primtrk_dqdx.end());
806  std::reverse(tmp_primtrk_resrange.begin(),tmp_primtrk_resrange.end());
807  }
808  primtrk_dqdx.push_back(tmp_primtrk_dqdx);
809  primtrk_resrange.push_back(tmp_primtrk_resrange);
810  primtrk_dedx.push_back(tmp_primtrk_dedx);
811  primtrk_hitx.push_back(tmp_primtrk_hitx);
812  primtrk_hity.push_back(tmp_primtrk_hity);
813  primtrk_hitz.push_back(tmp_primtrk_hitz);
814  primtrk_pitch.push_back(tmp_primtrk_pitch);
815  }
816  tmp_primtrk_dqdx.clear();
817  tmp_primtrk_resrange.clear();
818  tmp_primtrk_dedx.clear();
819  tmp_primtrk_hitx.clear();
820  tmp_primtrk_hity.clear();
821  tmp_primtrk_hitz.clear();
822  tmp_primtrk_pitch.clear();
823 
824  for(size_t k = 0; k < calovector.size() && k<3; k++){
825  fprimaryKineticEnergy[k] = calovector[k].KineticEnergy();
826  fprimaryRange[k] = calovector[k].Range();
827  //const std::vector< double > & dedxvec = calovector[k].dEdx();
828  }
829 
830  //**************Section for dE/dx correction**************************//
831  std::vector<float> Stw, Endw, Stt, Endt, Stwires, Endwires, Stticks, Endticks, TPCb, TPCcl;
832  Stw.clear(); Endw.clear(); Stt.clear(); Endt.clear(); Stwires.clear(); Endwires.clear(); Stticks.clear(); Endticks.clear(); TPCb.clear(); TPCcl.clear();
833  float den;
834  float numw, numt,wire_no,ticks_no;
835  for(size_t p1=0;p1<pfplist.size();p1++){
836  std::vector<art::Ptr<recob::Track>> trk=pftrack.at(p1);
837  std::vector<art::Ptr<recob::Cluster>> allClusters=fmcp.at(p1);
838  for(size_t c1=0;c1<allClusters.size();c1++){
839  if(allClusters[c1]->Plane().Plane!=2) continue;
840  if(trk.size() && int(trk[0].key())==fprimaryID){
841  Stw.push_back(allClusters[c1]->StartWire());
842  Endw.push_back(allClusters[c1]->EndWire());
843  Stt.push_back(allClusters[c1]->StartTick());
844  Endt.push_back(allClusters[c1]->EndTick());
845  TPCb.push_back(allClusters[c1]->Plane().TPC);
846  }
847  else{
848  Stwires.push_back(allClusters[c1]->StartWire());
849  Endwires.push_back(allClusters[c1]->EndWire());
850  Stticks.push_back(allClusters[c1]->StartTick());
851  Endticks.push_back(allClusters[c1]->EndTick());
852  TPCcl.push_back(allClusters[c1]->Plane().TPC);
853  }
854  }
855  }
856  for(size_t clt=0;clt<Stw.size();clt++){
857  for(size_t cl1=0;cl1<Stwires.size();cl1++){
858  if(TPCcl[cl1]!=TPCb[clt]) continue;
859  std::cout<<"tpc are equal "<<std::endl;
860  den=(Stw[clt]-Endw[clt])*(Stticks[cl1]-Endticks[cl1])-(Stt[clt]-Endt[clt])*(Stwires[cl1]-Endwires[cl1]);
861  if(den==0) continue;
862  numw=(Stw[clt]*Endt[clt]-Stt[clt]*Endw[clt])*(Stwires[cl1]-Endwires[cl1])-(Stw[clt]-Endw[clt])*(Stwires[cl1]*Endticks[cl1]-Stticks[cl1]*Endwires[cl1]);
863  numt=(Stw[clt]*Endt[clt]-Stt[clt]*Endw[clt])*(Stticks[cl1]-Endticks[cl1])-(Stt[clt]-Endt[clt])*(Stwires[cl1]*Endticks[cl1]-Stticks[cl1]*Endwires[cl1]);
864  wire_no=numw/den;
865  ticks_no=numt/den;
866  // std::cout<<"wireno and ticks not solution "<<wire_no<<" "<<ticks_no<<std::endl;
867  if(((Stw[clt]<wire_no && Endw[clt]>wire_no)||(Stw[clt]>wire_no && Endw[clt]<wire_no))&&((Stt[clt]<ticks_no && Endt[clt]>ticks_no)||(Stt[clt]>ticks_no && Endt[clt]<ticks_no)) && ((Stwires[cl1]<wire_no && Endwires[cl1]>wire_no)||(Stwires[cl1]>wire_no && Endwires[cl1]<wire_no)) && ((Stticks[cl1]<ticks_no && Endticks[cl1]>ticks_no)||(Stticks[cl1]>ticks_no && Endticks[cl1]<ticks_no)))
868  {
869  std::cout<<"intersection wire and ticks are "<<std::round(wire_no)<<" "<<ticks_no<<" Stw Endw StT EndT "<<Stwires[cl1]<<" "<<Endwires[cl1]<<" "<<Stticks[cl1]<<" "<<Endticks[cl1]<<std::endl;
870  double xyzStart[3];
871  double xyzEnd[3];
872  unsigned int wireno=std::round(wire_no);
873  geo::WireID wireid(0,TPCb[clt],2,wireno);
874  fGeometry->WireEndPoints(0,TPCb[clt],2,wireno, xyzStart, xyzEnd);
875  std::cout<<"Z position of intersection = "<<xyzStart[2]<<" "<<xyzEnd[2]<<" "<<wireno<<std::endl;
876  Zintersection.push_back(xyzStart[2]);
877  timeintersection.push_back(ticks_no);
878  }
879 
880  }
881  }
882 
883 
884 
885  ////*****************section for dE/dx correction*******************************////
886 
887 
888  }//this track
889 
890  else if(thisShower != 0x0){
891  fisprimarytrack = 0;
892  fisprimaryshower = 1;
893 
894  fprimaryID = thisShower->ID();
895  fprimaryLength = thisShower->Length();
896  fprimaryShowerBestPlane = thisShower->best_plane();
897  fprimaryOpeningAngle = thisShower->OpenAngle();
898  fprimaryStartPosition[0] = thisShower->ShowerStart().X();
899  fprimaryStartPosition[1] = thisShower->ShowerStart().Y();
900  fprimaryStartPosition[2] = thisShower->ShowerStart().Z();
901  fprimaryStartDirection[0] = thisShower->Direction().X();
902  fprimaryStartDirection[1] = thisShower->Direction().Y();
903  fprimaryStartDirection[2] = thisShower->Direction().Z();
904  if( (thisShower->Energy()).size() > 0 )
905  fprimaryShowerEnergy = thisShower->Energy()[0]; // thisShower->best_plane()
906  if( (thisShower->MIPEnergy()).size() > 0 )
907  fprimaryShowerMIPEnergy = thisShower->MIPEnergy()[0];
908  if( (thisShower->dEdx()).size() > 0 )
909  fprimaryShowerdEdx = thisShower->dEdx()[0];
910  }
911  else{
912  std::cout << "INFO::Primary pfParticle is not track or shower. Skip!" << std::endl;
913  continue;
914  }
915 
916  // Find the particle vertex. We need the tracker tag here because we need to do a bit of
917  // additional work if the PFParticle is track-like to find the vertex.
918  const TVector3 vtx = pfpUtil.GetPFParticleVertex(*particle,evt,fPFParticleTag,fTrackerTag);
919  fvertex[0] = vtx.X(); fvertex[1] = vtx.Y(); fvertex[2] = vtx.Z();
920 
921  // Now we can look for the interaction point of the particle if one exists, i.e where the particle
922  // scatters off an argon nucleus. Shower-like objects won't have an interaction point, so we can
923  // check this by making sure we get a sensible position
924  const TVector3 interactionVtx = pfpUtil.GetPFParticleSecondaryVertex(*particle,evt,fPFParticleTag,fTrackerTag);
925  fsecvertex[0] = interactionVtx.X(); fsecvertex[1] = interactionVtx.Y(); fsecvertex[2] = interactionVtx.Z();
926 
927  // Maximum number of daugthers to be processed
928  if(particle->NumDaughters() > NMAXDAUGTHERS)
929  std::cout << "INFO::Number of daughters is " << particle->NumDaughters() << ". Only the first NMAXDAUGTHERS are processed." << std::endl;
930 
931  // Let's get the daughter PFParticles... we can do this simply without the utility
932  for(const int daughterID : particle->Daughters()){
933  // Daughter ID is the element of the original recoParticle vector
934  const recob::PFParticle *daughterParticle = &(recoParticles->at(daughterID));
935  std::cout << "Daughter " << daughterID << " has " << daughterParticle->NumDaughters() << " daughters" << std::endl;
936 
937  const recob::Track* daughterTrack = pfpUtil.GetPFParticleTrack(*daughterParticle,evt,fPFParticleTag,fTrackerTag);
938  const recob::Shower* daughterShower = pfpUtil.GetPFParticleShower(*daughterParticle,evt,fPFParticleTag,fShowerTag);
939 
940  if(daughterTrack != 0x0){
941  fisdaughtertrack[fNDAUGHTERS] = 1;
942  fisdaughtershower[fNDAUGHTERS] = 0;
943  fdaughterTheta[fNDAUGHTERS] = daughterTrack->Theta();
944  fdaughterPhi[fNDAUGHTERS] = daughterTrack->Phi();
945  fdaughterLength[fNDAUGHTERS] = daughterTrack->Length();
946  fdaughterMomentum[fNDAUGHTERS] = daughterTrack->StartMomentum();
947  fdaughterEndMomentum[fNDAUGHTERS] = daughterTrack->EndMomentum();
948  fdaughterStartPosition[fNDAUGHTERS][0] = daughterTrack->Trajectory().Start().X();
949  fdaughterStartPosition[fNDAUGHTERS][1] = daughterTrack->Trajectory().Start().Y();
950  fdaughterStartPosition[fNDAUGHTERS][2] = daughterTrack->Trajectory().Start().Z();
951  fdaughterEndPosition[fNDAUGHTERS][0] = daughterTrack->Trajectory().End().X();
952  fdaughterEndPosition[fNDAUGHTERS][1] = daughterTrack->Trajectory().End().Y();
953  fdaughterEndPosition[fNDAUGHTERS][2] = daughterTrack->Trajectory().End().Z();
954  fdaughterStartDirection[fNDAUGHTERS][0] = daughterTrack->Trajectory().StartDirection().X();
955  fdaughterStartDirection[fNDAUGHTERS][1] = daughterTrack->Trajectory().StartDirection().Y();
956  fdaughterStartDirection[fNDAUGHTERS][2] = daughterTrack->Trajectory().StartDirection().Z();
957  fdaughterEndDirection[fNDAUGHTERS][0] = daughterTrack->Trajectory().EndDirection().X();
958  fdaughterEndDirection[fNDAUGHTERS][1] = daughterTrack->Trajectory().EndDirection().Y();
959  fdaughterEndDirection[fNDAUGHTERS][2] = daughterTrack->Trajectory().EndDirection().Z();
960 
961  fdaughterMomentumByRangeMuon[fNDAUGHTERS] = trmom.GetTrackMomentum(daughterTrack->Length(),13);
962  fdaughterMomentumByRangeProton[fNDAUGHTERS] = trmom.GetTrackMomentum(daughterTrack->Length(),2212);
963 
964  std::vector<anab::Calorimetry> daughtercalovector = trackUtil.GetRecoTrackCalorimetry(*daughterTrack, evt, fTrackerTag, fCalorimetryTag);
965  if(daughtercalovector.size() != 3)
966  std::cerr << "WARNING::Calorimetry vector size for daughter is = " << daughtercalovector.size() << std::endl;
967 
968  for(size_t k = 0; k < daughtercalovector.size() && k<3; k++){
969  fdaughterKineticEnergy[fNDAUGHTERS][k] = daughtercalovector[k].KineticEnergy();
970  fdaughterRange[fNDAUGHTERS][k] = daughtercalovector[k].Range();
971  }
972 
973  // Get the true mc particle
974  const simb::MCParticle* mcdaughterparticle = truthUtil.GetMCParticleFromRecoTrack(clockData, *daughterTrack, evt, fTrackerTag);
975  if(mcdaughterparticle != 0x0){
976  fdaughter_truth_TrackId[fNDAUGHTERS] = mcdaughterparticle->TrackId();
977  fdaughter_truth_Pdg[fNDAUGHTERS] = mcdaughterparticle->PdgCode();
978  fdaughter_truth_StartPosition[fNDAUGHTERS][0] = mcdaughterparticle->Vx();
979  fdaughter_truth_StartPosition[fNDAUGHTERS][1] = mcdaughterparticle->Vy();
980  fdaughter_truth_StartPosition[fNDAUGHTERS][2] = mcdaughterparticle->Vz();
981  fdaughter_truth_StartPosition[fNDAUGHTERS][3] = mcdaughterparticle->T();
982  fdaughter_truth_EndPosition[fNDAUGHTERS][0] = mcdaughterparticle->EndX();
983  fdaughter_truth_EndPosition[fNDAUGHTERS][1] = mcdaughterparticle->EndY();
984  fdaughter_truth_EndPosition[fNDAUGHTERS][2] = mcdaughterparticle->EndZ();
985  fdaughter_truth_EndPosition[fNDAUGHTERS][3] = mcdaughterparticle->EndT();
986  fdaughter_truth_P[fNDAUGHTERS] = mcdaughterparticle->P();
987  fdaughter_truth_Momentum[fNDAUGHTERS][0] = mcdaughterparticle->Px();
988  fdaughter_truth_Momentum[fNDAUGHTERS][1] = mcdaughterparticle->Py();
989  fdaughter_truth_Momentum[fNDAUGHTERS][2] = mcdaughterparticle->Pz();
990  fdaughter_truth_Momentum[fNDAUGHTERS][3] = mcdaughterparticle->E();
991  fdaughter_truth_Pt[fNDAUGHTERS] = mcdaughterparticle->Pt();
992  fdaughter_truth_Mass[fNDAUGHTERS] = mcdaughterparticle->Mass();
993  fdaughter_truth_EndMomentum[fNDAUGHTERS][0] = mcdaughterparticle->EndPx();
994  fdaughter_truth_EndMomentum[fNDAUGHTERS][1] = mcdaughterparticle->EndPy();
995  fdaughter_truth_EndMomentum[fNDAUGHTERS][2] = mcdaughterparticle->EndPz();
996  fdaughter_truth_EndMomentum[fNDAUGHTERS][3] = mcdaughterparticle->EndE();
997  fdaughter_truth_Theta[fNDAUGHTERS] = mcdaughterparticle->Momentum().Theta();
998  fdaughter_truth_Phi[fNDAUGHTERS] = mcdaughterparticle->Momentum().Phi();
999  fdaughter_truth_Process[fNDAUGHTERS] = int(mcdaughterparticle->Trajectory().ProcessToKey(mcdaughterparticle->Process()));
1000  std::cout << "Daughter Process = " << (mcdaughterparticle->Process()).c_str()
1001  << " , mother = " << mcdaughterparticle->Mother()
1002  << std::endl;
1003  }
1004  }
1005  else if(daughterShower != 0x0){
1006  fisdaughtertrack[fNDAUGHTERS] = 0;
1007  fisdaughtershower[fNDAUGHTERS] = 1;
1008  fdaughterLength[fNDAUGHTERS] = daughterShower->Length();
1009  fdaughterShowerBestPlane[fNDAUGHTERS] = daughterShower->best_plane();
1010  fdaughterOpeningAngle[fNDAUGHTERS] = daughterShower->OpenAngle();
1011  fdaughterStartPosition[fNDAUGHTERS][0] = daughterShower->ShowerStart().X();
1012  fdaughterStartPosition[fNDAUGHTERS][1] = daughterShower->ShowerStart().Y();
1013  fdaughterStartPosition[fNDAUGHTERS][2] = daughterShower->ShowerStart().Z();
1014  fdaughterStartDirection[fNDAUGHTERS][0] = daughterShower->Direction().X();
1015  fdaughterStartDirection[fNDAUGHTERS][1] = daughterShower->Direction().Y();
1016  fdaughterStartDirection[fNDAUGHTERS][2] = daughterShower->Direction().Z();
1017  if( (daughterShower->Energy()).size() > 0 )
1018  fdaughterShowerEnergy[fNDAUGHTERS] = daughterShower->Energy()[0]; // thisShower->best_plane()
1019  if( (daughterShower->MIPEnergy()).size() > 0 )
1020  fdaughterShowerMIPEnergy[fNDAUGHTERS] = daughterShower->MIPEnergy()[0];
1021  if( (daughterShower->dEdx()).size() > 0 )
1022  fdaughterShowerdEdx[fNDAUGHTERS] = daughterShower->dEdx()[0];
1023  }
1024  else{
1025  std::cout << "INFO::Daughter pfParticle is not track or shower. Skip!" << std::endl;
1026  continue;
1027  }
1028 
1029  fdaughterID[fNDAUGHTERS] = daughterID;
1030  // NHits associated with this pfParticle
1031  fdaughterNHits[fNDAUGHTERS] = (pfpUtil.GetPFParticleHits(*daughterParticle,evt,fPFParticleTag)).size();
1032  // T0
1033  // std::vector<anab::T0> pfdaughterT0vec = pfpUtil.GetPFParticleT0(*daughterParticle,evt,fPFParticleTag);
1034  // if(!pfT0vec.empty())
1035  // fdaughterT0[fNDAUGHTERS] = pfdaughterT0vec[0].Time();
1036 
1037  fNDAUGHTERS++;
1038 
1039  // Only process NMAXDAUGTHERS
1040  if(fNDAUGHTERS > NMAXDAUGTHERS) break;
1041 
1042  }
1043 
1044  break;
1045  }//particle loop pfparticle
1046 
1047  // Fill trees
1048  if(beamTriggerEvent)
1049  fPandoraBeam->Fill();
1050 
1051  //fPandoraCosmics->Fill();
1052 }//analyzer
1053 
1054 void protoana::truepion::endJob(){
1055 
1056 }
1057 
1058 void protoana::truepion::FillCosmicsTree(art::Event const & evt, std::string pfParticleTag){
1059 
1060  // To fill
1061 
1062 }
1063 
1064 
1065 
1066 void protoana::truepion::Initialise(){
1067 
1068  fRun = -999;
1069  fSubRun = -999;
1070  fevent = -999;
1071  fTimeStamp = -999.0;
1072 
1073  for(int k=0; k < 3; k++){
1074  fvertex[k] = -999.0;
1075  fsecvertex[k] = -999.0;
1076  fprimaryEndPosition[k] = -999.0;
1077  fprimaryStartPosition[k] = -999.0;
1078  fprimaryEndDirection[k] = -999.0;
1079  fprimaryStartDirection[k] = -999.0;
1080  fprimaryKineticEnergy[k] = -999.0;
1081  fprimaryRange[k] = -999.0;
1082  }
1083  fbeamtrackMomentum = -999.0;
1084  fbeamtrackEnergy = 999.0;
1085  fbeamtrackPdg = -9999;
1086  fbeamtrackTime = -999.0;
1087  fbeamtrackID = -999;
1088  for(int l=0; l < 3; l++){
1089  fbeamtrackP[l] = -999.0;
1090  fbeamtrackPos[l] = -999.0;
1091  fbeamtrackDir[l] = -999.0;
1092  }
1093 
1094  //NumberBeamTrajectoryPoints=0;
1095  beamtrk_x.clear();
1096  beamtrk_y.clear();
1097  beamtrk_z.clear();
1098  beamtrk_Px.clear();
1099  beamtrk_Py.clear();
1100  beamtrk_Pz.clear();
1101  beamtrk_Eng.clear();
1102 
1103 
1104 
1105  fisprimarytrack = 0;
1106  fisprimaryshower = 0;
1107  fNDAUGHTERS = 0;
1108 
1109  fprimaryBDTScore = -999.0;
1110  fprimaryNHits = -999;
1111  fprimaryTheta = -999.0;
1112  fprimaryPhi = -999.0;
1113  fprimaryLength = -999.0;
1114  fprimaryMomentum = -999.0;
1115  fprimaryEndMomentum = -999.0;
1116  fprimaryOpeningAngle = -999.0;
1117  fprimaryShowerBestPlane = -999;
1118  fprimaryShowerEnergy = -999;
1119  fprimaryShowerMIPEnergy = -999;
1120  fprimaryShowerdEdx = -999;
1121  fprimaryID = -999;
1122  fprimaryMomentumByRangeProton = -999.0;
1123  fprimaryMomentumByRangeMuon = -999.0;
1124 
1125  fprimary_truth_TrackId = -999;
1126  fprimary_truth_Pdg = -999;
1127  ftruthpdg=-999;
1128  fprimary_truth_P = -999.0;
1129  fprimary_truth_Pt = -999.0;
1130  fprimary_truth_Mass = -999.0;
1131  fprimary_truth_Theta = -999.0;
1132  fprimary_truth_Phi = -999.0;
1133  fprimary_truth_Process = -999;
1134  fprimary_truth_Isbeammatched = -999;
1135  fprimary_truth_NDaughters = -999;
1136  fprimary_truth_tracklength=-9999;
1137  fprimary_truth_EndProcess="";
1138  truth_last_process="";
1139  for(int l=0; l < 4; l++){
1140  fprimary_truth_StartPosition[l] = -999.0;
1141  fprimary_truth_EndPosition[l] = -999.0;
1142  fprimary_truth_Momentum[l] = -999.0;
1143  fprimary_truth_EndMomentum[l] = -999.0;
1144  }
1145 
1146  for(int k=0; k < NMAXDAUGTHERS; k++){
1147  fisdaughtertrack[k] = -999;
1148  fisdaughtershower[k] = -999;
1149  fdaughterNHits[k] = -999;
1150  fdaughterTheta[k] = -999.0;
1151  fdaughterPhi[k] = -999.0;
1152  fdaughterLength[k] = -999.0;
1153  fdaughterMomentum[k] = -999.0;
1154  fdaughterEndMomentum[k] = -999.0;
1155  for(int l=0; l < 3; l++){
1156  fdaughterEndPosition[k][l] = -999.0;
1157  fdaughterStartPosition[k][l] = -999.0;
1158  fdaughterEndDirection[k][l] = -999.0;
1159  fdaughterStartDirection[k][l] = -999.0;
1160  fdaughterKineticEnergy[k][l] = -999.0;
1161  fdaughterRange[k][l] = -999.0;
1162  }
1163  fdaughterOpeningAngle[k] = -999.0;
1164  fdaughterShowerBestPlane[k] = -999;
1165  fdaughterShowerEnergy[k] = -999;
1166  fdaughterShowerMIPEnergy[k] = -999;
1167  fdaughterShowerdEdx[k] = -999;
1168  fdaughterMomentumByRangeProton[k] = -999.0;
1169  fdaughterMomentumByRangeMuon[k] = -999.0;
1170  fdaughterID[k] = -999;
1171  // fdaughterT0[k] = -999;
1172 
1173  fdaughter_truth_TrackId[k] = -999;
1174  fdaughter_truth_Pdg[k] = -999;
1175  fdaughter_truth_P[k] = -999.0;
1176  fdaughter_truth_Pt[k] = -999.0;
1177  fdaughter_truth_Mass[k] = -999.0;
1178  fdaughter_truth_Theta[k] = -999.0;
1179  fdaughter_truth_Phi[k] = -999.0;
1180  fdaughter_truth_Process[k] = -999;
1181  for(int l=0; l < 4; l++){
1182  fdaughter_truth_StartPosition[k][l] = -999.0;
1183  fdaughter_truth_EndPosition[k][l] = -999.0;
1184  fdaughter_truth_Momentum[k][l] = -999.0;
1185  fdaughter_truth_EndMomentum[k][l] = -999.0;
1186  }
1187  }
1188 
1189  primtrk_dqdx.clear();
1190  primtrk_resrange.clear();
1191  primtrk_dedx.clear();
1192  primtrk_range.clear();
1193  primtrk_hitx.clear();
1194  primtrk_hity.clear();
1195  primtrk_hitz.clear();
1196  primtrk_pitch.clear();
1197  primtrk_truth_Z.clear();
1198  primtrk_truth_Eng.clear();
1199  primtrk_truth_trkide.clear();
1200  interactionX.clear();
1201  interactionY.clear();
1202  interactionZ.clear();
1203  interactionProcesslist.clear();
1204  interactionAngles.clear();
1205  Zintersection.clear();
1206  timeintersection.clear();
1207 }
1208 
1209 DEFINE_ART_MODULE(protoana::truepion)
simb::MCParticle::E
double E(const int i=0) const
Definition: MCParticle.h:233
protoana::truepion::beamtrk_Py
std::vector< float > beamtrk_Py
Definition: truepion_module.cc:153
protoana::truepion::fdaughterMomentumByRangeProton
double fdaughterMomentumByRangeProton[NMAXDAUGTHERS]
Definition: truepion_module.cc:254
recob::Shower::best_plane
int best_plane() const
Definition: Shower.h:200
protoana::truepion::fdaughterEndMomentum
double fdaughterEndMomentum[NMAXDAUGTHERS]
Definition: truepion_module.cc:244
protoana::truepion::beginJob
virtual void beginJob() override
Definition: truepion_module.cc:299
simb::MCParticle::NumberTrajectoryPoints
unsigned int NumberTrajectoryPoints() const
Definition: MCParticle.h:218
protoana::truepion::fdaughter_truth_Theta
double fdaughter_truth_Theta[NMAXDAUGTHERS]
Definition: truepion_module.cc:270
recob::Shower::ShowerStart
const TVector3 & ShowerStart() const
Definition: Shower.h:192
simb::MCParticle::Position
const TLorentzVector & Position(const int i=0) const
Definition: MCParticle.h:219
recob::Track::EndMomentum
double EndMomentum() const
Definition: Track.h:144
simb::MCParticle::PdgCode
int PdgCode() const
Definition: MCParticle.h:212
protoana::truepion::fdaughter_truth_EndMomentum
double fdaughter_truth_EndMomentum[NMAXDAUGTHERS][4]
Definition: truepion_module.cc:267
art::Timestamp::timeLow
constexpr std::uint32_t timeLow() const
Definition: Timestamp.h:29
protoana::truepion::fprimary_truth_EndMomentum
double fprimary_truth_EndMomentum[4]
Definition: truepion_module.cc:200
protoana::truepion::primtrk_truth_Eng
std::vector< std::vector< double > > primtrk_truth_Eng
Definition: truepion_module.cc:230
simb::MCParticle::Py
double Py(const int i=0) const
Definition: MCParticle.h:231
protoana::truepion::primtrk_hitx
std::vector< std::vector< double > > primtrk_hitx
Definition: truepion_module.cc:225
recob::PFParticle::NumDaughters
int NumDaughters() const
Returns the number of daughter particles flowing from this one.
Definition: PFParticle.h:89
protoana::ProtoDUNETruthUtils::GetGeantGoodParticle
const simb::MCParticle * GetGeantGoodParticle(const simb::MCTruth &genTruth, const art::Event &evt) const
Definition: ProtoDUNETruthUtils.cxx:860
simb::MCParticle::EndE
double EndE() const
Definition: MCParticle.h:244
simb::MCParticle::EndZ
double EndZ() const
Definition: MCParticle.h:228
protoana::truepion::primtrk_pitch
std::vector< std::vector< double > > primtrk_pitch
Definition: truepion_module.cc:228
recob::Shower::Length
double Length() const
Definition: Shower.h:201
protoana::truepion::primtrk_truth_trkide
std::vector< std::vector< double > > primtrk_truth_trkide
Definition: truepion_module.cc:231
protoana::ProtoDUNEPFParticleUtils::GetPFParticleShower
const recob::Shower * GetPFParticleShower(const recob::PFParticle &particle, art::Event const &evt, const std::string particleLabel, const std::string showerLabel) const
Get the shower associated to this particle. Returns a null pointer if not found.
Definition: ProtoDUNEPFParticleUtils.cxx:496
art::DataViewImpl::getHandle
Handle< PROD > getHandle(SelectorBase const &) const
Definition: DataViewImpl.h:382
geo::View_t
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
string
std::string string
Definition: nybbler.cc:12
protoana::truepion::fPFParticleTag
std::string fPFParticleTag
Definition: truepion_module.cc:123
protoana::truepion::fdaughter_truth_Pt
double fdaughter_truth_Pt[NMAXDAUGTHERS]
Definition: truepion_module.cc:268
simb::MCParticle::Trajectory
const simb::MCTrajectory & Trajectory() const
Definition: MCParticle.h:253
simb::MCParticle::Mother
int Mother() const
Definition: MCParticle.h:213
recob::Track::Trajectory
const recob::TrackTrajectory & Trajectory() const
Access to the stored recob::TrackTrajectory.
Definition: Track.h:98
protoana::truepion::truepion
truepion(fhicl::ParameterSet const &p)
Definition: truepion_module.cc:284
protoana::ProtoDUNETrackUtils::GetRecoTrackCalorimetry
std::vector< anab::Calorimetry > GetRecoTrackCalorimetry(const recob::Track &track, art::Event const &evt, const std::string trackModule, const std::string caloModule) const
Get the Calorimetry(s) from a given reco track.
Definition: ProtoDUNETrackUtils.cxx:72
simb::MCTrajectory::at
const value_type & at(const size_type i) const
Definition: MCTrajectory.h:175
protoana::truepion::primtrk_range
std::vector< double > primtrk_range
Definition: truepion_module.cc:224
art::Timestamp::timeHigh
constexpr std::uint32_t timeHigh() const
Definition: Timestamp.h:34
simb::MCParticle::Mass
double Mass() const
Definition: MCParticle.h:239
cet::pow
constexpr T pow(T x)
Definition: pow.h:72
simb::MCTrajectory::KeyToProcess
std::string KeyToProcess(unsigned char const &key) const
Definition: MCTrajectory.cxx:117
protoana::ProtoDUNEPFParticleUtils::GetPFParticleHits
const std::vector< const recob::Hit * > GetPFParticleHits(const recob::PFParticle &particle, art::Event const &evt, const std::string particleLabel) const
Get the hits associated to the PFParticle.
Definition: ProtoDUNEPFParticleUtils.cxx:564
recob::Shower::Energy
const std::vector< double > & Energy() const
Definition: Shower.h:195
simb::MCParticle::Px
double Px(const int i=0) const
Definition: MCParticle.h:230
protoana::ProtoDUNEPFParticleUtils::GetPFParticleVertex
const TVector3 GetPFParticleVertex(const recob::PFParticle &particle, art::Event const &evt, const std::string particleLabel, const std::string trackLabel) const
Function to find the interaction vertex of a primary PFParticle.
Definition: ProtoDUNEPFParticleUtils.cxx:352
protoana::truepion::beamtrk_Pz
std::vector< float > beamtrk_Pz
Definition: truepion_module.cc:154
protoana::truepion::endJob
virtual void endJob() override
Definition: truepion_module.cc:1054
protoana::truepion::fdaughter_truth_StartPosition
double fdaughter_truth_StartPosition[NMAXDAUGTHERS][4]
Definition: truepion_module.cc:263
protoana::truepion::truth_last_process
std::string truth_last_process
Definition: truepion_module.cc:197
recob::tracking::Plane::Plane
Plane(const Point_t &planePos, const Vector_t &planeDir)
Constructor from reference position on the plane and direction orthogonal to the plane.
Definition: TrackingPlane.h:61
protoana::truepion::fprimaryStartDirection
double fprimaryStartDirection[3]
Definition: truepion_module.cc:175
cheat::BackTrackerService::TrackIdToSimIDEs_Ps
std::vector< const sim::IDE * > TrackIdToSimIDEs_Ps(int const &id) const
Definition: BackTrackerService_service.cc:134
std
STL namespace.
protoana::truepion::interactionZ
std::vector< double > interactionZ
Definition: truepion_module.cc:214
protoana::truepion::fdaughterOpeningAngle
double fdaughterOpeningAngle[NMAXDAUGTHERS]
Definition: truepion_module.cc:249
art::EDAnalyzer::EDAnalyzer
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:25
simb::MCParticle::Process
std::string Process() const
Definition: MCParticle.h:215
protoana::truepion::fisdaughtertrack
int fisdaughtertrack[NMAXDAUGTHERS]
Definition: truepion_module.cc:237
MCParticle.h
Particle class.
simb::MCTrajectory::ProcessToKey
unsigned char ProcessToKey(std::string const &process) const
Definition: MCTrajectory.cxx:97
simb::MCParticle::EndY
double EndY() const
Definition: MCParticle.h:227
protoana::truepion::beamtrk_Px
std::vector< float > beamtrk_Px
Definition: truepion_module.cc:152
protoana::truepion::fdaughterRange
double fdaughterRange[NMAXDAUGTHERS][3]
Definition: truepion_module.cc:257
protoana::truepion::fdaughter_truth_P
double fdaughter_truth_P[NMAXDAUGTHERS]
Definition: truepion_module.cc:265
protoana::truepion::fdaughterShowerMIPEnergy
double fdaughterShowerMIPEnergy[NMAXDAUGTHERS]
Definition: truepion_module.cc:251
l
static QStrList * l
Definition: config.cpp:1044
simb::MCParticle::NumberDaughters
int NumberDaughters() const
Definition: MCParticle.h:217
NMAXDAUGTHERS
const int NMAXDAUGTHERS
Definition: truepion_module.cc:79
Geometry.h
art framework interface to geometry description
simb::MCParticle::TrackId
int TrackId() const
Definition: MCParticle.h:210
protoana::truepion::fdaughterStartDirection
double fdaughterStartDirection[NMAXDAUGTHERS][3]
Definition: truepion_module.cc:248
protoana::truepion::beamtrk_Eng
std::vector< float > beamtrk_Eng
Definition: truepion_module.cc:155
protoana::truepion::fprimaryKineticEnergy
double fprimaryKineticEnergy[3]
Definition: truepion_module.cc:183
protoana::truepion::fdaughterKineticEnergy
double fdaughterKineticEnergy[NMAXDAUGTHERS][3]
Definition: truepion_module.cc:256
protoana::truepion::fdaughterShowerdEdx
double fdaughterShowerdEdx[NMAXDAUGTHERS]
Definition: truepion_module.cc:252
protoana::truepion::fprimary_truth_StartPosition
double fprimary_truth_StartPosition[4]
Definition: truepion_module.cc:194
protoana::truepion::fdaughterLength
double fdaughterLength[NMAXDAUGTHERS]
Definition: truepion_module.cc:242
util::size
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
Definition: StdUtils.h:92
protoana::truepion::interactionX
std::vector< double > interactionX
Definition: truepion_module.cc:212
art::DataViewImpl::isRealData
bool isRealData() const
Definition: DataViewImpl.cc:113
protoana::truepion::beamtrk_z
std::vector< float > beamtrk_z
Definition: truepion_module.cc:151
simb::MCParticle::Pt
double Pt(const int i=0) const
Definition: MCParticle.h:236
protoana::truepion::fdaughterID
int fdaughterID[NMAXDAUGTHERS]
Definition: truepion_module.cc:258
recob::Track::Phi
double Phi() const
Definition: Track.h:178
protoana::truepion::primtrk_hity
std::vector< std::vector< double > > primtrk_hity
Definition: truepion_module.cc:226
recob::Track::Length
double Length(size_t p=0) const
Access to various track properties.
Definition: Track.h:167
protoana::truepion::fGeometry
geo::GeometryCore const * fGeometry
Definition: truepion_module.cc:127
art::DataViewImpl::time
Timestamp time() const
Definition: DataViewImpl.cc:106
recob::Shower::dEdx
const std::vector< double > & dEdx() const
Definition: Shower.h:203
DEFINE_ART_MODULE
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:67
recob::Shower::MIPEnergy
const std::vector< double > & MIPEnergy() const
Definition: Shower.h:198
simb::MCParticle::EndProcess
std::string EndProcess() const
Definition: MCParticle.h:216
protoana::truepion::fprimary_truth_EndPosition
double fprimary_truth_EndPosition[4]
Definition: truepion_module.cc:195
recob::Track::Start
Point_t const & Start() const
Access to track position at different points.
Definition: Track.h:123
protoana::truepion::primtrk_dedx
std::vector< std::vector< double > > primtrk_dedx
Definition: truepion_module.cc:223
protoana::truepion::operator=
truepion & operator=(truepion const &)=delete
protoana::truepion::Zintersection
std::vector< double > Zintersection
Definition: truepion_module.cc:217
wirecell.dfp.graph.key
def key(type, name=None)
Definition: graph.py:13
simb::MCParticle::EndPz
double EndPz() const
Definition: MCParticle.h:243
protoana::truepion::fprimary_truth_Momentum
double fprimary_truth_Momentum[4]
Definition: truepion_module.cc:199
protoana::truepion::fCalorimetryTag
std::string fCalorimetryTag
Definition: truepion_module.cc:120
recob::Track::Theta
double Theta() const
Access to spherical or geographical angles at vertex or at any point.
Definition: Track.h:176
protoana::truepion::fdaughter_truth_TrackId
int fdaughter_truth_TrackId[NMAXDAUGTHERS]
Definition: truepion_module.cc:261
simb::MCParticle::P
double P(const int i=0) const
Definition: MCParticle.h:234
protoana::ProtoDUNEPFParticleUtils::GetPFParticleTrack
const recob::Track * GetPFParticleTrack(const recob::PFParticle &particle, art::Event const &evt, const std::string particleLabel, const std::string trackLabel) const
Get the track associated to this particle. Returns a null pointer if not found.
Definition: ProtoDUNEPFParticleUtils.cxx:476
protoana::truepion::fdaughter_truth_Phi
double fdaughter_truth_Phi[NMAXDAUGTHERS]
Definition: truepion_module.cc:271
protoana::truepion::fdaughter_truth_Mass
double fdaughter_truth_Mass[NMAXDAUGTHERS]
Definition: truepion_module.cc:269
recob::Shower::OpenAngle
double OpenAngle() const
Definition: Shower.h:202
simb::MCParticle::T
double T(const int i=0) const
Definition: MCParticle.h:224
art::DataViewImpl::getValidHandle
ValidHandle< PROD > getValidHandle(InputTag const &tag) const
Definition: DataViewImpl.h:441
test.p
p
Definition: test.py:223
protoana::truepion::fdaughterEndPosition
double fdaughterEndPosition[NMAXDAUGTHERS][3]
Definition: truepion_module.cc:245
protoana::truepion::primtrk_truth_Z
std::vector< std::vector< double > > primtrk_truth_Z
Definition: truepion_module.cc:229
protoana::truepion::primtrk_dqdx
std::vector< std::vector< double > > primtrk_dqdx
Definition: truepion_module.cc:221
art::DataViewImpl::subRun
SubRunNumber_t subRun() const
Definition: DataViewImpl.cc:78
protoana::truepion::fdaughterPhi
double fdaughterPhi[NMAXDAUGTHERS]
Definition: truepion_module.cc:241
protoana::ProtoDUNETruthUtils::GetMCParticleFromRecoTrack
const simb::MCParticle * GetMCParticleFromRecoTrack(detinfo::DetectorClocksData const &clockData, const recob::Track &track, art::Event const &evt, std::string trackModule) const
Definition: ProtoDUNETruthUtils.cxx:609
protoana::truepion::fdaughterMomentumByRangeMuon
double fdaughterMomentumByRangeMuon[NMAXDAUGTHERS]
Definition: truepion_module.cc:255
protoana::ProtoDUNEPFParticleUtils::GetPFParticlesFromBeamSlice
const std::vector< const recob::PFParticle * > GetPFParticlesFromBeamSlice(art::Event const &evt, const std::string particleLabel) const
Return the pointers for the PFParticles in the beam slice. Returns an empty vector is no beam slice w...
Definition: ProtoDUNEPFParticleUtils.cxx:175
protoana::truepion::fdaughterTheta
double fdaughterTheta[NMAXDAUGTHERS]
Definition: truepion_module.cc:240
recob::Shower::Direction
const TVector3 & Direction() const
Definition: Shower.h:189
protoana::truepion::fisdaughtershower
int fisdaughtershower[NMAXDAUGTHERS]
Definition: truepion_module.cc:238
art::DataViewImpl::run
RunNumber_t run() const
Definition: DataViewImpl.cc:71
simb::MCParticle::EndT
double EndT() const
Definition: MCParticle.h:229
geo::GeometryCore
Description of geometry of one entire detector.
Definition: GeometryCore.h:1468
protoana::truepion::fGeneratorTag
std::string fGeneratorTag
Definition: truepion_module.cc:124
protoana::truepion::fprimaryEndPosition
double fprimaryEndPosition[3]
Definition: truepion_module.cc:172
simb::MCTrajectory::TrajectoryProcesses
ProcessMap const & TrajectoryProcesses() const
Definition: MCTrajectory.h:188
geo::GeometryCore::View
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
Definition: GeometryCore.cxx:601
geo_types.h
Definition of data types for geometry description.
protoana::truepion::pfpUtil
protoana::ProtoDUNEPFParticleUtils pfpUtil
Definition: truepion_module.cc:108
protoana::truepion::fdaughterShowerBestPlane
int fdaughterShowerBestPlane[NMAXDAUGTHERS]
Definition: truepion_module.cc:253
protoana::truepion::fprimaryMomentumByRangeMuon
double fprimaryMomentumByRangeMuon
Definition: truepion_module.cc:182
protoana::truepion::timeintersection
std::vector< double > timeintersection
Definition: truepion_module.cc:218
protoana::truepion::beamtrk_x
std::vector< float > beamtrk_x
Definition: truepion_module.cc:149
keras_to_tensorflow.int
int
Definition: keras_to_tensorflow.py:69
protoana::truepion::fdaughter_truth_Pdg
int fdaughter_truth_Pdg[NMAXDAUGTHERS]
Definition: truepion_module.cc:262
protoana::truepion::fdaughterMomentum
double fdaughterMomentum[NMAXDAUGTHERS]
Definition: truepion_module.cc:243
recob::TrackTrajectory::EndDirection
Vector_t EndDirection() const
Returns the direction of the trajectory at the last point.
Definition: TrackTrajectory.h:321
protoana::truepion::fdaughterStartPosition
double fdaughterStartPosition[NMAXDAUGTHERS][3]
Definition: truepion_module.cc:246
protoana::truepion::fdaughter_truth_EndPosition
double fdaughter_truth_EndPosition[NMAXDAUGTHERS][4]
Definition: truepion_module.cc:264
simb::MCParticle::Vx
double Vx(const int i=0) const
Definition: MCParticle.h:221
recob::Track::ID
int ID() const
Definition: Track.h:198
Hit.h
Declaration of signal hit object.
recob::Track::StartMomentum
double StartMomentum() const
Definition: Track.h:143
recob::PFParticle
Hierarchical representation of particle flow.
Definition: PFParticle.h:44
protoana::truepion::FillCosmicsTree
void FillCosmicsTree(art::Event const &evt, std::string pfParticleTag)
Definition: truepion_module.cc:1058
m_pi
double m_pi
Definition: truepion_module.cc:80
protoana::truepion::interactionAngles
std::vector< double > interactionAngles
Definition: truepion_module.cc:216
simb::MCTrajectory::size
size_type size() const
Definition: MCTrajectory.h:166
protoana::truepion::fdaughter_truth_Momentum
double fdaughter_truth_Momentum[NMAXDAUGTHERS][4]
Definition: truepion_module.cc:266
simb::MCParticle::EndPy
double EndPy() const
Definition: MCParticle.h:242
reverse
static unsigned int reverse(QString &chars, unsigned char *level, unsigned int a, unsigned int b)
Definition: qstring.cpp:11649
protoana::truepion::primtrk_hitz
std::vector< std::vector< double > > primtrk_hitz
Definition: truepion_module.cc:227
recob::TrackTrajectory::End
Point_t const & End() const
Returns the position of the last valid point of the trajectory [cm].
Definition: TrackTrajectory.h:257
simb::MCParticle::Momentum
const TLorentzVector & Momentum(const int i=0) const
Definition: MCParticle.h:220
simb::MCParticle::Pz
double Pz(const int i=0) const
Definition: MCParticle.h:232
Track.h
Provides recob::Track data product.
simb::MCParticle::Vz
double Vz(const int i=0) const
Definition: MCParticle.h:223
protoana::truepion::fprimaryStartPosition
double fprimaryStartPosition[3]
Definition: truepion_module.cc:173
protoana::truepion::fprimary_truth_EndProcess
std::string fprimary_truth_EndProcess
Definition: truepion_module.cc:196
geo::GeometryCore::WireEndPoints
void WireEndPoints(geo::WireID const &wireid, double *xyzStart, double *xyzEnd) const
Fills two arrays with the coordinates of the wire end points.
Definition: GeometryCore.cxx:1243
protoana::truepion::fdaughterShowerEnergy
double fdaughterShowerEnergy[NMAXDAUGTHERS]
Definition: truepion_module.cc:250
protoana::ProtoDUNEPFParticleUtils::GetPFParticleSecondaryVertex
const TVector3 GetPFParticleSecondaryVertex(const recob::PFParticle &particle, art::Event const &evt, const std::string particleLabel, const std::string trackLabel) const
Function to find the secondary interaction vertex of a primary PFParticle.
Definition: ProtoDUNEPFParticleUtils.cxx:407
art::EventID::event
EventNumber_t event() const
Definition: EventID.h:116
recob::Track::End
Point_t const & End() const
Definition: Track.h:125
reco_momentum_tuples.count
int count
Definition: reco_momentum_tuples.py:202
protoana::truepion::analyze
void analyze(art::Event const &evt) override
Definition: truepion_module.cc:440
simb::MCParticle::EndPx
double EndPx() const
Definition: MCParticle.h:241
protoana::truepion::trackUtil
protoana::ProtoDUNETrackUtils trackUtil
Definition: truepion_module.cc:109
protoana::truepion::dataUtil
protoana::ProtoDUNEDataUtils dataUtil
Definition: truepion_module.cc:107
tca::evt
TCEvent evt
Definition: DataStructs.cxx:7
protoana::truepion::interactionProcesslist
std::vector< std::string > interactionProcesslist
Definition: truepion_module.cc:215
art::get
auto const & get(AssnsNode< L, R, D > const &r)
Definition: AssnsNode.h:115
art::fill_ptr_vector
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:297
simb::MCParticle::EndX
double EndX() const
Definition: MCParticle.h:226
protoana::ProtoDUNEPFParticleUtils::GetNumberPrimaryPFParticle
unsigned int GetNumberPrimaryPFParticle(art::Event const &evt, const std::string particleLabel) const
Get the number of primary PFParticles.
Definition: ProtoDUNEPFParticleUtils.cxx:24
recob::TrackTrajectory::StartDirection
Vector_t StartDirection() const
Returns the direction of the trajectory at the first point.
Definition: TrackTrajectory.h:317
util::quantities::second
second_as<> second
Type of time stored in seconds, in double precision.
Definition: spacetime.h:85
trkf::Plane
recob::tracking::Plane Plane
Definition: TrackState.h:17
recob::TrackTrajectory::Start
Point_t const & Start() const
Returns the position of the first valid point of the trajectory [cm].
Definition: TrackTrajectory.h:253
protoana::truepion::fdaughterNHits
int fdaughterNHits[NMAXDAUGTHERS]
Definition: truepion_module.cc:239
prim_energy
double prim_energy
Definition: truepion_module.cc:81
protoana::truepion::interactionY
std::vector< double > interactionY
Definition: truepion_module.cc:213
bool
int bool
Definition: qglobal.h:345
protoana::truepion::fdaughter_truth_Process
int fdaughter_truth_Process[NMAXDAUGTHERS]
Definition: truepion_module.cc:272
recob::Shower::ID
int ID() const
Definition: Shower.h:187
protoana::truepion::truthUtil
protoana::ProtoDUNETruthUtils truthUtil
Definition: truepion_module.cc:110
protoana::truepion::primtrk_resrange
std::vector< std::vector< double > > primtrk_resrange
Definition: truepion_module.cc:222
protoana::truepion::fBeamModuleLabel
const art::InputTag fBeamModuleLabel
Definition: truepion_module.cc:119
protoana::ProtoDUNEPFParticleUtils::GetBeamCosmicScore
float GetBeamCosmicScore(const recob::PFParticle &particle, art::Event const &evt, const std::string particleLabel) const
Access the BDT output used to decide if a slice is beam-like or cosmic-like.
Definition: ProtoDUNEPFParticleUtils.cxx:182
protoana::truepion::fprimaryEndDirection
double fprimaryEndDirection[3]
Definition: truepion_module.cc:174
art::Event::id
EventID id() const
Definition: Event.cc:34
simb::MCParticle::Vy
double Vy(const int i=0) const
Definition: MCParticle.h:222
recob::Track
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a "fitted" track:
Definition: Track.h:49
calo
calorimetry
Definition: CaloChecker_module.cc:43
trkf::TrackMomentumCalculator::GetTrackMomentum
double GetTrackMomentum(double trkrange, int pdg) const
Definition: TrackMomentumCalculator.cxx:144
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
protoana::truepion::beamtrk_y
std::vector< float > beamtrk_y
Definition: truepion_module.cc:150
protoana::truepion::fdaughterEndDirection
double fdaughterEndDirection[NMAXDAUGTHERS][3]
Definition: truepion_module.cc:247
protoana::truepion::fprimaryMomentumByRangeProton
double fprimaryMomentumByRangeProton
Definition: truepion_module.cc:181
Generated by   doxygen 1.8.11