All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
SeedAna_module.cc
Go to the documentation of this file.
1 //
2 // Name: SeedAna_module.cc
3 //
4 // Purpose: Module SeedAna.
5 //
6 // Configuration parameters.
7 //
8 // SeedModuleLabel: Seed module label.
9 // MCTrackModuleLabel: MCTrack module label.
10 // MinMCKE: Minimum MC particle kinetic energy.
11 // MatchColinearity: Minimum colinearity for mc-seed matching.
12 // MatchDisp: Maximum uv displacement for mc-seed matching.
13 //
14 // Created: 2-Aug-2011 H. Greenlee
15 //
16 
17 #include <cmath>
18 #include <iomanip>
19 #include <map>
20 #include <sstream>
21 
22 #include "art/Framework/Core/EDAnalyzer.h"
23 #include "art/Framework/Core/ModuleMacros.h"
24 #include "art/Framework/Principal/Event.h"
25 #include "art_root_io/TFileService.h"
26 #include "cetlib_except/exception.h"
27 #include "messagefacility/MessageLogger/MessageLogger.h"
28 
29 #include "larcore/Geometry/Geometry.h"
30 #include "lardata/DetectorInfoServices/DetectorPropertiesService.h"
31 #include "lardataobj/MCBase/MCTrack.h"
32 #include "lardataobj/RecoBase/Seed.h"
33 
34 #include "TH1F.h"
35 #include "TH2F.h"
36 #include "TMatrixD.h"
37 
38 namespace {
39 
40  // Calculate distance to boundary.
41  //----------------------------------------------------------------------------
42  double
43  bdist(const TVector3& pos, unsigned int /*tpc*/ = 0, unsigned int /*cstat*/ = 0)
44  {
45  // Get geometry.
46 
47  art::ServiceHandle<geo::Geometry const> geom;
48 
49  double d1 = pos.X(); // Distance to right side (wires).
50  double d2 = 2. * geom->DetHalfWidth() - pos.X(); // Distance to left side (cathode).
51  double d3 = pos.Y() + geom->DetHalfHeight(); // Distance to bottom.
52  double d4 = geom->DetHalfHeight() - pos.Y(); // Distance to top.
53  double d5 = pos.Z(); // Distance to front.
54  double d6 = geom->DetLength() - pos.Z(); // Distance to back.
55 
56  return std::min({d1, d2, d3, d4, d5, d6});
57  }
58 
59  // Find the closest matching mc trajectory point (sim::MCStep) for a given seed.
60  // Returned value is index of the trajectory point.
61  // Return -1 in case of no match.
62  int
63  mcmatch(detinfo::DetectorPropertiesData const& detProp,
64  const sim::MCTrack& mctrk,
65  const recob::Seed& seed)
66  {
67  // Get seed point.
68 
69  double pos[3];
70  double err[3];
71  seed.GetPoint(pos, err);
72 
73  // Calculate the x offset due to nonzero mc particle time.
74  double mctime = mctrk.Start().T(); // nsec
75  double mcdx = mctime * 1.e-3 * detProp.DriftVelocity(); // cm
76 
77  // Loop over trajectory points.
78 
79  int best_traj = -1;
80  double max_dist = 0.;
81  int ntraj = mctrk.size();
82  for (int itraj = 0; itraj < ntraj; ++itraj) {
83  const TLorentzVector& vec = mctrk[itraj].Position();
84  double dx = pos[0] - vec.X() - mcdx;
85  double dy = pos[1] - vec.Y();
86  double dz = pos[2] - vec.Z();
87  double dist = std::sqrt(dx * dx + dy * dy + dz * dz);
88  if (best_traj < 0 || dist < max_dist) {
89  best_traj = itraj;
90  max_dist = dist;
91  }
92  }
93  return best_traj;
94  }
95 
96  // Length of MCTrack.
97  // In this function, the extracted start and end momenta are converted to GeV
98  // (MCTrack stores momenta in Mev).
99  //----------------------------------------------------------------------------
100  double
101  length(detinfo::DetectorPropertiesData const& detProp,
102  const sim::MCTrack& mctrk,
103  double dx,
104  TVector3& start,
105  TVector3& end,
106  TVector3& startmom,
107  TVector3& endmom,
108  unsigned int /*tpc*/ = 0,
109  unsigned int /*cstat*/ = 0)
110  {
111  // Get services.
112 
113  art::ServiceHandle<geo::Geometry const> geom;
114 
115  // Get fiducial volume boundary.
116 
117  double xmin = 0.;
118  double xmax = 2. * geom->DetHalfWidth();
119  double ymin = -geom->DetHalfHeight();
120  double ymax = geom->DetHalfHeight();
121  double zmin = 0.;
122  double zmax = geom->DetLength();
123  double result = 0.;
124  TVector3 disp;
125  int n = mctrk.size();
126  bool first = true;
127 
128  for (int i = 0; i < n; ++i) {
129  TVector3 pos = mctrk[i].Position().Vect();
130 
131  // Make fiducial cuts. Require the particle to be within the physical volume of
132  // the tpc, and also require the apparent x position to be within the expanded
133  // readout frame.
134 
135  if (pos.X() >= xmin && pos.X() <= xmax && pos.Y() >= ymin && pos.Y() <= ymax &&
136  pos.Z() >= zmin && pos.Z() <= zmax) {
137  pos[0] += dx;
138  double ticks = detProp.ConvertXToTicks(pos[0], 0, 0, 0);
139  if (ticks >= 0. && ticks < detProp.ReadOutWindowSize()) {
140  if (first) {
141  start = pos;
142  startmom = 0.001 * mctrk[i].Momentum().Vect();
143  }
144  else {
145  disp -= pos;
146  result += disp.Mag();
147  }
148  first = false;
149  disp = pos;
150  end = pos;
151  endmom = 0.001 * mctrk[i].Momentum().Vect();
152  }
153  }
154  }
155 
156  return result;
157  }
158 
159  // Fill efficiency histogram assuming binomial errors.
160 
161  void
162  effcalc(const TH1* hnum, const TH1* hden, TH1* heff)
163  {
164  int nbins = hnum->GetNbinsX();
165  if (nbins != hden->GetNbinsX())
166  throw cet::exception("SeedAna") << "effcalc[" __FILE__ "]: incompatible histograms (I)\n";
167  if (nbins != heff->GetNbinsX())
168  throw cet::exception("SeedAna") << "effcalc[" __FILE__ "]: incompatible histograms (II)\n";
169 
170  // Loop over bins, including underflow and overflow.
171 
172  for (int ibin = 0; ibin <= nbins + 1; ++ibin) {
173  double num = hnum->GetBinContent(ibin);
174  double den = hden->GetBinContent(ibin);
175  if (den == 0.) {
176  heff->SetBinContent(ibin, 0.);
177  heff->SetBinError(ibin, 0.);
178  }
179  else {
180  double eff = num / den;
181  if (eff < 0.) eff = 0.;
182  if (eff > 1.) eff = 1.;
183  double err = std::sqrt(eff * (1. - eff) / den);
184  heff->SetBinContent(ibin, eff);
185  heff->SetBinError(ibin, err);
186  }
187  }
188 
189  heff->SetMinimum(0.);
190  heff->SetMaximum(1.05);
191  heff->SetMarkerStyle(20);
192  }
193 
194  // Fill multiplicity histogram.
195 
196  void
197  mulcalc(const TH1* hnum, const TH1* hden, TH1* hmul)
198  {
199  int nbins = hnum->GetNbinsX();
200  if (nbins != hden->GetNbinsX())
201  throw cet::exception("SeedAna") << "mulcalc[" __FILE__ "]: incompatible histograms (I)\n";
202  if (nbins != hmul->GetNbinsX())
203  throw cet::exception("SeedAna") << "mulcalc[" __FILE__ "]: incompatible histograms (II)\n";
204 
205  // Loop over bins, including underflow and overflow.
206 
207  for (int ibin = 0; ibin <= nbins + 1; ++ibin) {
208  double num = hnum->GetBinContent(ibin);
209  double den = hden->GetBinContent(ibin);
210  if (den == 0.) {
211  hmul->SetBinContent(ibin, 0.);
212  hmul->SetBinError(ibin, 0.);
213  }
214  else {
215  double mul = num / den;
216  if (mul < 0.) mul = 0.;
217  double err = std::sqrt((1. + mul) * mul / den);
218  hmul->SetBinContent(ibin, mul);
219  hmul->SetBinError(ibin, err);
220  }
221  }
222 
223  hmul->SetMinimum(0.);
224  hmul->SetMarkerStyle(20);
225  }
226 }
227 
228 namespace trkf {
229 
230  class SeedAna : public art::EDAnalyzer {
231  public:
232  // Embedded structs.
233 
234  // Struct for histograms that depend on seeds only.
235 
236  struct RecoHists {
237  RecoHists(const std::string& subdir);
238 
239  // Pure reco seed histograms.
240 
241  TH1F* fHx{nullptr}; // Seed x position.
242  TH1F* fHy{nullptr}; // Seed y position.
243  TH1F* fHz{nullptr}; // Seed z position.
244  TH1F* fHdist{nullptr}; // Seed distance to boundary.
245  TH1F* fHtheta{nullptr}; // Theta.
246  TH1F* fHphi{nullptr}; // Phi.
247  TH1F* fHtheta_xz{nullptr}; // Theta_xz.
248  TH1F* fHtheta_yz{nullptr}; // Theta_yz.
249  };
250 
251  // Struct for mc particles and mc-matched tracks.
252 
253  struct MCHists {
254  MCHists(const std::string& subdir);
255 
256  // Reco-MC matching.
257 
258  TH2F* fHduvcosth{nullptr}; // 2D mc vs. data matching, duv vs. cos(theta).
259  TH1F* fHcosth{nullptr}; // 1D direction matching, cos(theta).
260  TH1F* fHmcu{nullptr}; // 1D endpoint truth u.
261  TH1F* fHmcv{nullptr}; // 1D endpoint truth v.
262  TH1F* fHmcw{nullptr}; // 1D endpoint truth w.
263  TH1F* fHmcdudw{nullptr}; // Truth du/dw.
264  TH1F* fHmcdvdw{nullptr}; // Truth dv/dw.
265 
266  // Pure MC particle histograms (efficiency denominator).
267 
268  TH1F* fHmcstartx{nullptr}; // Starting x position.
269  TH1F* fHmcstarty{nullptr}; // Starting y position.
270  TH1F* fHmcstartz{nullptr}; // Starting z position.
271  TH1F* fHmcendx{nullptr}; // Ending x position.
272  TH1F* fHmcendy{nullptr}; // Ending y position.
273  TH1F* fHmcendz{nullptr}; // Ending z position.
274  TH1F* fHmctheta{nullptr}; // Theta.
275  TH1F* fHmcphi{nullptr}; // Phi.
276  TH1F* fHmctheta_xz{nullptr}; // Theta_xz.
277  TH1F* fHmctheta_yz{nullptr}; // Theta_yz.
278  TH1F* fHmcmom{nullptr}; // Momentum.
279  TH1F* fHmclen{nullptr}; // Length.
280 
281  // Matched seed histograms (multiplicity numerator).
282 
283  TH1F* fHmstartx{nullptr}; // Starting x position.
284  TH1F* fHmstarty{nullptr}; // Starting y position.
285  TH1F* fHmstartz{nullptr}; // Starting z position.
286  TH1F* fHmendx{nullptr}; // Ending x position.
287  TH1F* fHmendy{nullptr}; // Ending y position.
288  TH1F* fHmendz{nullptr}; // Ending z position.
289  TH1F* fHmtheta{nullptr}; // Theta.
290  TH1F* fHmphi{nullptr}; // Phi.
291  TH1F* fHmtheta_xz{nullptr}; // Theta_xz.
292  TH1F* fHmtheta_yz{nullptr}; // Theta_yz.
293  TH1F* fHmmom{nullptr}; // Momentum.
294  TH1F* fHmlen{nullptr}; // Length.
295 
296  // Matched seed histograms (efficiency numerator).
297 
298  TH1F* fHgstartx{nullptr}; // Starting x position.
299  TH1F* fHgstarty{nullptr}; // Starting y position.
300  TH1F* fHgstartz{nullptr}; // Starting z position.
301  TH1F* fHgendx{nullptr}; // Ending x position.
302  TH1F* fHgendy{nullptr}; // Ending y position.
303  TH1F* fHgendz{nullptr}; // Ending z position.
304  TH1F* fHgtheta{nullptr}; // Theta.
305  TH1F* fHgphi{nullptr}; // Phi.
306  TH1F* fHgtheta_xz{nullptr}; // Theta_xz.
307  TH1F* fHgtheta_yz{nullptr}; // Theta_yz.
308  TH1F* fHgmom{nullptr}; // Momentum.
309  TH1F* fHglen{nullptr}; // Length.
310 
311  // Multiplicity histograms.
312 
313  TH1F* fHmulstartx{nullptr}; // Starting x position.
314  TH1F* fHmulstarty{nullptr}; // Starting y position.
315  TH1F* fHmulstartz{nullptr}; // Starting z position.
316  TH1F* fHmulendx{nullptr}; // Ending x position.
317  TH1F* fHmulendy{nullptr}; // Ending y position.
318  TH1F* fHmulendz{nullptr}; // Ending z position.
319  TH1F* fHmultheta{nullptr}; // Theta.
320  TH1F* fHmulphi{nullptr}; // Phi.
321  TH1F* fHmultheta_xz{nullptr}; // Theta_xz.
322  TH1F* fHmultheta_yz{nullptr}; // Theta_yz.
323  TH1F* fHmulmom{nullptr}; // Momentum.
324  TH1F* fHmullen{nullptr}; // Length.
325 
326  // Efficiency histograms.
327 
328  TH1F* fHestartx{nullptr}; // Starting x position.
329  TH1F* fHestarty{nullptr}; // Starting y position.
330  TH1F* fHestartz{nullptr}; // Starting z position.
331  TH1F* fHeendx{nullptr}; // Ending x position.
332  TH1F* fHeendy{nullptr}; // Ending y position.
333  TH1F* fHeendz{nullptr}; // Ending z position.
334  TH1F* fHetheta{nullptr}; // Theta.
335  TH1F* fHephi{nullptr}; // Phi.
336  TH1F* fHetheta_xz{nullptr}; // Theta_xz.
337  TH1F* fHetheta_yz{nullptr}; // Theta_yz.
338  TH1F* fHemom{nullptr}; // Momentum.
339  TH1F* fHelen{nullptr}; // Length.
340  };
341 
342  explicit SeedAna(fhicl::ParameterSet const& pset);
343 
344  private:
345  void analyze(const art::Event& evt) override;
346  void endJob() override;
347 
348  // Fcl Attributes.
349 
350  std::string fSeedModuleLabel;
351  std::string fMCTrackModuleLabel;
352  int fDump; // Number of events to dump to debug message facility.
353  double fMinMCKE; // Minimum MC particle kinetic energy (GeV).
354  double fMinMCLen; // Minimum MC particle length in tpc (cm).
355  double fMatchColinearity; // Minimum matching colinearity.
356  double fMatchDisp; // Maximum matching displacement.
357  bool fIgnoreSign; // Ignore sign of mc particle if true.
358 
359  // Histograms.
360 
361  std::map<int, MCHists> fMCHistMap; // Indexed by pdg id.
362  std::map<int, RecoHists> fRecoHistMap; // Indexed by pdg id.
363 
364  // Statistics.
365 
366  int fNumEvent;
367  };
368 
369  DEFINE_ART_MODULE(SeedAna)
370 
371  // RecoHists methods.
372 
373  SeedAna::RecoHists::RecoHists(const std::string& subdir)
374  //
375  // Purpose: Initializing constructor.
376  //
377  {
378  // Get services.
379 
380  art::ServiceHandle<geo::Geometry const> geom;
381  art::ServiceHandle<art::TFileService> tfs;
382 
383  // Make histogram directory.
384 
385  art::TFileDirectory topdir = tfs->mkdir("seedana", "SeedAna histograms");
386  art::TFileDirectory dir = topdir.mkdir(subdir);
387 
388  // Book histograms.
389 
390  fHx =
391  dir.make<TH1F>("x", "X Position", 100, -2. * geom->DetHalfWidth(), 4. * geom->DetHalfWidth());
392  fHy = dir.make<TH1F>("y", "Y Position", 100, -geom->DetHalfHeight(), geom->DetHalfHeight());
393  fHz = dir.make<TH1F>("z", "Z Position", 100, 0., geom->DetLength());
394  fHdist =
395  dir.make<TH1F>("dist", "Position Distance to Boundary", 100, -10., geom->DetHalfWidth());
396  fHtheta = dir.make<TH1F>("theta", "Theta", 100, 0., 3.142);
397  fHphi = dir.make<TH1F>("phi", "Phi", 100, -3.142, 3.142);
398  fHtheta_xz = dir.make<TH1F>("theta_xz", "Theta_xz", 100, -3.142, 3.142);
399  fHtheta_yz = dir.make<TH1F>("theta_yz", "Theta_yz", 100, -3.142, 3.142);
400  }
401 
402  // MCHists methods.
403 
404  SeedAna::MCHists::MCHists(const std::string& subdir)
405  //
406  // Purpose: Initializing constructor.
407  //
408  {
409  // Get services.
410 
411  art::ServiceHandle<geo::Geometry const> geom;
412  art::ServiceHandle<art::TFileService> tfs;
413 
414  // Make histogram directory.
415 
416  art::TFileDirectory topdir = tfs->mkdir("seedana", "SeedAna histograms");
417  art::TFileDirectory dir = topdir.mkdir(subdir);
418 
419  // Book histograms.
420 
421  fHduvcosth =
422  dir.make<TH2F>("duvcosth", "Delta(uv) vs. Colinearity", 100, 0.95, 1., 100, 0., 1.);
423  fHcosth = dir.make<TH1F>("colin", "Colinearity", 100, 0.95, 1.);
424  fHmcu = dir.make<TH1F>("mcu", "MC Truth U", 100, -5., 5.);
425  fHmcv = dir.make<TH1F>("mcv", "MC Truth V", 100, -5., 5.);
426  fHmcw = dir.make<TH1F>("mcw", "MC Truth W", 100, -20., 20.);
427  fHmcdudw = dir.make<TH1F>("mcdudw", "MC Truth U Slope", 100, -0.2, 0.2);
428  fHmcdvdw = dir.make<TH1F>("mcdvdw", "MV Truth V Slope", 100, -0.2, 0.2);
429 
430  fHmcstartx = dir.make<TH1F>(
431  "mcxstart", "MC X Start Position", 10, -2. * geom->DetHalfWidth(), 4. * geom->DetHalfWidth());
432  fHmcstarty = dir.make<TH1F>(
433  "mcystart", "MC Y Start Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
434  fHmcstartz = dir.make<TH1F>("mczstart", "MC Z Start Position", 10, 0., geom->DetLength());
435  fHmcendx = dir.make<TH1F>(
436  "mcxend", "MC X End Position", 10, -2. * geom->DetHalfWidth(), 4. * geom->DetHalfWidth());
437  fHmcendy = dir.make<TH1F>(
438  "mcyend", "MC Y End Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
439  fHmcendz = dir.make<TH1F>("mczend", "MC Z End Position", 10, 0., geom->DetLength());
440  fHmctheta = dir.make<TH1F>("mctheta", "MC Theta", 20, 0., 3.142);
441  fHmcphi = dir.make<TH1F>("mcphi", "MC Phi", 10, -3.142, 3.142);
442  fHmctheta_xz = dir.make<TH1F>("mctheta_xz", "MC Theta_xz", 40, -3.142, 3.142);
443  fHmctheta_yz = dir.make<TH1F>("mctheta_yz", "MC Theta_yz", 40, -3.142, 3.142);
444  fHmcmom = dir.make<TH1F>("mcmom", "MC Momentum", 10, 0., 10.);
445  fHmclen = dir.make<TH1F>("mclen", "MC Particle Length", 10, 0., 1.1 * geom->DetLength());
446 
447  fHmstartx = dir.make<TH1F>("mxstart",
448  "Matched X Start Position",
449  10,
450  -2. * geom->DetHalfWidth(),
451  4. * geom->DetHalfWidth());
452  fHmstarty = dir.make<TH1F>(
453  "mystart", "Matched Y Start Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
454  fHmstartz = dir.make<TH1F>("mzstart", "Matched Z Start Position", 10, 0., geom->DetLength());
455  fHmendx = dir.make<TH1F>(
456  "mxend", "Matched X End Position", 10, -2. * geom->DetHalfWidth(), 4. * geom->DetHalfWidth());
457  fHmendy = dir.make<TH1F>(
458  "myend", "Matched Y End Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
459  fHmendz = dir.make<TH1F>("mzend", "Matched Z End Position", 10, 0., geom->DetLength());
460  fHmtheta = dir.make<TH1F>("mtheta", "Matched Theta", 20, 0., 3.142);
461  fHmphi = dir.make<TH1F>("mphi", "Matched Phi", 10, -3.142, 3.142);
462  fHmtheta_xz = dir.make<TH1F>("mtheta_xz", "Matched Theta_xz", 40, -3.142, 3.142);
463  fHmtheta_yz = dir.make<TH1F>("mtheta_yz", "Matched Theta_yz", 40, -3.142, 3.142);
464  fHmmom = dir.make<TH1F>("mmom", "Matched Momentum", 10, 0., 10.);
465  fHmlen = dir.make<TH1F>("mlen", "Matched Particle Length", 10, 0., 1.1 * geom->DetLength());
466 
467  fHgstartx = dir.make<TH1F>("gxstart",
468  "Good X Start Position",
469  10,
470  -2. * geom->DetHalfWidth(),
471  4. * geom->DetHalfWidth());
472  fHgstarty = dir.make<TH1F>(
473  "gystart", "Good Y Start Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
474  fHgstartz = dir.make<TH1F>("gzstart", "Good Z Start Position", 10, 0., geom->DetLength());
475  fHgendx = dir.make<TH1F>(
476  "gxend", "Good X End Position", 10, -2. * geom->DetHalfWidth(), 4. * geom->DetHalfWidth());
477  fHgendy = dir.make<TH1F>(
478  "gyend", "Good Y End Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
479  fHgendz = dir.make<TH1F>("gzend", "Good Z End Position", 10, 0., geom->DetLength());
480  fHgtheta = dir.make<TH1F>("gtheta", "Good Theta", 20, 0., 3.142);
481  fHgphi = dir.make<TH1F>("gphi", "Good Phi", 10, -3.142, 3.142);
482  fHgtheta_xz = dir.make<TH1F>("gtheta_xz", "Good Theta_xz", 40, -3.142, 3.142);
483  fHgtheta_yz = dir.make<TH1F>("gtheta_yz", "Good Theta_yz", 40, -3.142, 3.142);
484  fHgmom = dir.make<TH1F>("gmom", "Good Momentum", 10, 0., 10.);
485  fHglen = dir.make<TH1F>("glen", "Good Particle Length", 10, 0., 1.1 * geom->DetLength());
486 
487  fHmulstartx = dir.make<TH1F>("mulxstart",
488  "Multiplicity vs. X Start Position",
489  10,
490  -2. * geom->DetHalfWidth(),
491  4. * geom->DetHalfWidth());
492  fHmulstarty = dir.make<TH1F>("mulystart",
493  "Multiplicity vs. Y Start Position",
494  10,
495  -geom->DetHalfHeight(),
496  geom->DetHalfHeight());
497  fHmulstartz =
498  dir.make<TH1F>("mulzstart", "Multiplicity vs. Z Start Position", 10, 0., geom->DetLength());
499  fHmulendx = dir.make<TH1F>("mulxend",
500  "Multiplicity vs. X End Position",
501  10,
502  -2. * geom->DetHalfWidth(),
503  4. * geom->DetHalfWidth());
504  fHmulendy = dir.make<TH1F>("mulyend",
505  "Multiplicity vs. Y End Position",
506  10,
507  -geom->DetHalfHeight(),
508  geom->DetHalfHeight());
509  fHmulendz =
510  dir.make<TH1F>("mulzend", "Multiplicity vs. Z End Position", 10, 0., geom->DetLength());
511  fHmultheta = dir.make<TH1F>("multheta", "Multiplicity vs. Theta", 20, 0., 3.142);
512  fHmulphi = dir.make<TH1F>("mulphi", "Multiplicity vs. Phi", 10, -3.142, 3.142);
513  fHmultheta_xz = dir.make<TH1F>("multheta_xz", "Multiplicity vs. Theta_xz", 40, -3.142, 3.142);
514  fHmultheta_yz = dir.make<TH1F>("multheta_yz", "Multiplicity vs. Theta_yz", 40, -3.142, 3.142);
515  fHmulmom = dir.make<TH1F>("mulmom", "Multiplicity vs. Momentum", 10, 0., 10.);
516  fHmullen =
517  dir.make<TH1F>("mullen", "Multiplicity vs. Particle Length", 10, 0., 1.1 * geom->DetLength());
518 
519  fHestartx = dir.make<TH1F>("exstart",
520  "Efficiency vs. X Start Position",
521  10,
522  -2. * geom->DetHalfWidth(),
523  4. * geom->DetHalfWidth());
524  fHestarty = dir.make<TH1F>("eystart",
525  "Efficiency vs. Y Start Position",
526  10,
527  -geom->DetHalfHeight(),
528  geom->DetHalfHeight());
529  fHestartz =
530  dir.make<TH1F>("ezstart", "Efficiency vs. Z Start Position", 10, 0., geom->DetLength());
531  fHeendx = dir.make<TH1F>("exend",
532  "Efficiency vs. X End Position",
533  10,
534  -2. * geom->DetHalfWidth(),
535  4. * geom->DetHalfWidth());
536  fHeendy = dir.make<TH1F>(
537  "eyend", "Efficiency vs. Y End Position", 10, -geom->DetHalfHeight(), geom->DetHalfHeight());
538  fHeendz = dir.make<TH1F>("ezend", "Efficiency vs. Z End Position", 10, 0., geom->DetLength());
539  fHetheta = dir.make<TH1F>("etheta", "Efficiency vs. Theta", 20, 0., 3.142);
540  fHephi = dir.make<TH1F>("ephi", "Efficiency vs. Phi", 10, -3.142, 3.142);
541  fHetheta_xz = dir.make<TH1F>("etheta_xz", "Efficiency vs. Theta_xz", 40, -3.142, 3.142);
542  fHetheta_yz = dir.make<TH1F>("etheta_yz", "Efficiency vs. Theta_yz", 40, -3.142, 3.142);
543  fHemom = dir.make<TH1F>("emom", "Efficiency vs. Momentum", 10, 0., 10.);
544  fHelen =
545  dir.make<TH1F>("elen", "Efficiency vs. Particle Length", 10, 0., 1.1 * geom->DetLength());
546  }
547 
548  SeedAna::SeedAna(const fhicl::ParameterSet& pset)
549  //
550  // Purpose: Constructor.
551  //
552  // Arguments: pset - Module parameters.
553  //
554  : EDAnalyzer(pset)
555  , fSeedModuleLabel(pset.get<std::string>("SeedModuleLabel"))
556  , fMCTrackModuleLabel(pset.get<std::string>("MCTrackModuleLabel"))
557  , fDump(pset.get<int>("Dump"))
558  , fMinMCKE(pset.get<double>("MinMCKE"))
559  , fMinMCLen(pset.get<double>("MinMCLen"))
560  , fMatchColinearity(pset.get<double>("MatchColinearity"))
561  , fMatchDisp(pset.get<double>("MatchDisp"))
562  , fIgnoreSign(pset.get<bool>("IgnoreSign"))
563  , fNumEvent(0)
564  {
565  mf::LogInfo("SeedAna") << "SeedAna configured with the following parameters:\n"
566  << " SeedModuleLabel = " << fSeedModuleLabel << "\n"
567  << " MCTrackModuleLabel = " << fMCTrackModuleLabel << "\n"
568  << " Dump = " << fDump << "\n"
569  << " MinMCKE = " << fMinMCKE << "\n"
570  << " MinMCLen = " << fMinMCLen;
571  }
572 
573  void
574  SeedAna::analyze(const art::Event& evt)
575  //
576  // Purpose: Analyze method.
577  //
578  // Arguments: event - Art event.
579  //
580  {
581  ++fNumEvent;
582 
583  // Optional dump stream.
584 
585  std::unique_ptr<mf::LogInfo> pdump;
586  if (fDump > 0) {
587  --fDump;
588  pdump = std::unique_ptr<mf::LogInfo>(new mf::LogInfo("TrackAna"));
589  }
590 
591  // Make sure histograms are booked.
592 
593  bool mc = !evt.isRealData();
594 
595  // Get seed handle.
596 
597  art::Handle<std::vector<recob::Seed>> seedh;
598  evt.getByLabel(fSeedModuleLabel, seedh);
599 
600  // Seed->mc track matching map.
601 
602  std::map<const recob::Seed*, int> seedmap;
603 
604  if (mc) {
605 
606  // Get MCTracks.
607 
608  art::Handle<std::vector<sim::MCTrack>> mctrackh;
609  evt.getByLabel(fMCTrackModuleLabel, mctrackh);
610 
611  // Dump MCTracks.
612 
613  if (pdump) {
614  *pdump << "MC Tracks\n"
615  << " Id pdg x y z dx dy dz "
616  " p\n"
617  << "--------------------------------------------------------------------------------"
618  "-----------\n";
619  }
620 
621  // Loop over mc tracks, and fill histograms that depend only
622  // on mc particles.
623  auto const detProp =
624  art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(evt);
625 
626  for (std::vector<sim::MCTrack>::const_iterator imctrk = mctrackh->begin();
627  imctrk != mctrackh->end();
628  ++imctrk) {
629  const sim::MCTrack& mctrk = *imctrk;
630  int pdg = mctrk.PdgCode();
631  if (fIgnoreSign) pdg = std::abs(pdg);
632 
633  // Ignore everything except stable charged nonshowering particles.
634 
635  int apdg = std::abs(pdg);
636  if (apdg == 13 || // Muon
637  apdg == 211 || // Charged pion
638  apdg == 321 || // Charged kaon
639  apdg == 2212) { // (Anti)proton
640 
641  // Apply minimum energy cut.
642 
643  if (mctrk.Start().E() >= mctrk.Start().Momentum().Mag() + 1000. * fMinMCKE) {
644 
645  // Calculate the x offset due to nonzero mc particle time.
646 
647  double mctime = mctrk.Start().T(); // nsec
648  double mcdx = mctime * 1.e-3 * detProp.DriftVelocity(); // cm
649 
650  // Calculate the length of this mc particle inside the fiducial volume.
651 
652  TVector3 mcstart;
653  TVector3 mcend;
654  TVector3 mcstartmom;
655  TVector3 mcendmom;
656  double plen = length(detProp, mctrk, mcdx, mcstart, mcend, mcstartmom, mcendmom);
657 
658  // Apply minimum fiducial length cut. Always reject particles that have
659  // zero length in the tpc regardless of the configured cut.
660 
661  if (plen > 0. && plen > fMinMCLen) {
662 
663  // Dump MC particle information here.
664 
665  if (pdump) {
666  double pstart = mcstartmom.Mag();
667  double pend = mcendmom.Mag();
668  *pdump << "\nOffset" << std::setw(3) << mctrk.TrackID() << std::setw(6)
669  << mctrk.PdgCode() << " " << std::fixed << std::setprecision(2)
670  << std::setw(10) << mcdx << "\nStart " << std::setw(3) << mctrk.TrackID()
671  << std::setw(6) << mctrk.PdgCode() << " " << std::fixed
672  << std::setprecision(2) << std::setw(10) << mcstart[0] << std::setw(10)
673  << mcstart[1] << std::setw(10) << mcstart[2];
674  if (pstart > 0.) {
675  *pdump << " " << std::fixed << std::setprecision(3) << std::setw(10)
676  << mcstartmom[0] / pstart << std::setw(10) << mcstartmom[1] / pstart
677  << std::setw(10) << mcstartmom[2] / pstart;
678  }
679  else
680  *pdump << std::setw(32) << " ";
681  *pdump << std::setw(12) << std::fixed << std::setprecision(3) << pstart;
682  *pdump << "\nEnd " << std::setw(3) << mctrk.TrackID() << std::setw(6)
683  << mctrk.PdgCode() << " " << std::fixed << std::setprecision(2)
684  << std::setw(10) << mcend[0] << std::setw(10) << mcend[1] << std::setw(10)
685  << mcend[2];
686  if (pend > 0.01) {
687  *pdump << " " << std::fixed << std::setprecision(3) << std::setw(10)
688  << mcendmom[0] / pend << std::setw(10) << mcendmom[1] / pend
689  << std::setw(10) << mcendmom[2] / pend;
690  }
691  else
692  *pdump << std::setw(32) << " ";
693  *pdump << std::setw(12) << std::fixed << std::setprecision(3) << pend << "\n";
694  }
695 
696  // Fill histograms.
697 
698  if (fMCHistMap.count(pdg) == 0) {
699  std::ostringstream ostr;
700  ostr << "MC" << (fIgnoreSign ? "All" : (pdg > 0 ? "Pos" : "Neg")) << std::abs(pdg);
701  fMCHistMap.emplace(pdg, MCHists{ostr.str()});
702  }
703  const MCHists& mchists = fMCHistMap.at(pdg);
704 
705  double mctheta_xz = std::atan2(mcstartmom.X(), mcstartmom.Z());
706  double mctheta_yz = std::atan2(mcstartmom.Y(), mcstartmom.Z());
707 
708  mchists.fHmcstartx->Fill(mcstart.X());
709  mchists.fHmcstarty->Fill(mcstart.Y());
710  mchists.fHmcstartz->Fill(mcstart.Z());
711  mchists.fHmcendx->Fill(mcend.X());
712  mchists.fHmcendy->Fill(mcend.Y());
713  mchists.fHmcendz->Fill(mcend.Z());
714  mchists.fHmctheta->Fill(mcstartmom.Theta());
715  mchists.fHmcphi->Fill(mcstartmom.Phi());
716  mchists.fHmctheta_xz->Fill(mctheta_xz);
717  mchists.fHmctheta_yz->Fill(mctheta_yz);
718  mchists.fHmcmom->Fill(mcstartmom.Mag());
719  mchists.fHmclen->Fill(plen);
720 
721  // Loop over seeds and do matching.
722 
723  int nmatch = 0;
724  if (seedh.isValid()) {
725 
726  // Loop over seeds.
727 
728  int nseed = seedh->size();
729  for (int i = 0; i < nseed; ++i) {
730  art::Ptr<recob::Seed> pseed(seedh, i);
731  const recob::Seed& seed = *pseed;
732  if (seedmap.count(&seed) == 0) seedmap[&seed] = -1;
733 
734  // Get parameters of this seed.
735 
736  TVector3 pos;
737  TVector3 dir;
738  double err[3];
739  seed.GetPoint(&pos[0], err);
740  seed.GetDirection(&dir[0], err);
741 
742  // Calculate the global-to-local rotation matrix.
743  // Copied from Track.cxx.
744 
745  TMatrixD rot(3, 3);
746  double dirmag = dir.Mag();
747  double diryz = std::sqrt(dir.Y() * dir.Y() + dir.Z() * dir.Z());
748 
749  double sinth = dir.X() / dirmag;
750  double costh = diryz / dirmag;
751  double sinphi = 0.;
752  double cosphi = 1.;
753  if (diryz != 0) {
754  sinphi = -dir.Y() / diryz;
755  cosphi = dir.Z() / diryz;
756  }
757  rot(0, 0) = costh;
758  rot(1, 0) = 0.;
759  rot(2, 0) = sinth;
760  rot(0, 1) = sinth * sinphi;
761  rot(1, 1) = cosphi;
762  rot(2, 1) = -costh * sinphi;
763  rot(0, 2) = -sinth * cosphi;
764  rot(1, 2) = sinphi;
765  rot(2, 2) = costh * cosphi;
766 
767  // Get best matching mc trajectory point.
768 
769  int itraj = mcmatch(detProp, mctrk, seed);
770  if (itraj >= 0) {
771 
772  // Get mc relative position and direction at matching trajectory point.
773 
774  TVector3 mcpos = mctrk[itraj].Position().Vect() - pos;
775  TVector3 mcmom = mctrk[itraj].Momentum().Vect();
776  mcpos[0] += mcdx;
777 
778  // Rotate the momentum and position to the
779  // seed-local coordinate system.
780 
781  TVector3 mcmoml = rot * mcmom;
782  TVector3 mcposl = rot * mcpos;
783 
784  if (mcmoml.Z() < 0.) mcmoml = -mcmoml;
785  double costh = mcmoml.Z() / mcmoml.Mag();
786 
787  double u = mcposl.X();
788  double v = mcposl.Y();
789  double w = mcposl.Z();
790 
791  double pu = mcmoml.X();
792  double pv = mcmoml.Y();
793  double pw = mcmoml.Z();
794 
795  double dudw = pu / pw;
796  double dvdw = pv / pw;
797 
798  double u0 = u - w * dudw;
799  double v0 = v - w * dvdw;
800  double uv0 = std::sqrt(u0 * u0 + v0 * v0);
801 
802  // Fill matching histograms.
803 
804  mchists.fHduvcosth->Fill(costh, uv0);
805  if (std::abs(uv0) < fMatchDisp) {
806 
807  // Fill slope matching histograms.
808 
809  mchists.fHmcdudw->Fill(dudw);
810  mchists.fHmcdvdw->Fill(dvdw);
811  }
812  mchists.fHcosth->Fill(costh);
813  if (costh > fMatchColinearity) {
814 
815  // Fill displacement matching histograms.
816 
817  mchists.fHmcu->Fill(u0);
818  mchists.fHmcv->Fill(v0);
819  mchists.fHmcw->Fill(w);
820 
821  if (std::abs(uv0) < fMatchDisp) {
822 
823  // Now we have passed all matching cuts and we have a matching
824  // mc particle + seed pair.
825 
826  ++nmatch;
827  seedmap[&seed] = mctrk.TrackID();
828 
829  // Fill matched seed histograms (seed multiplicity).
830 
831  mchists.fHmstartx->Fill(mcstart.X());
832  mchists.fHmstarty->Fill(mcstart.Y());
833  mchists.fHmstartz->Fill(mcstart.Z());
834  mchists.fHmendx->Fill(mcend.X());
835  mchists.fHmendy->Fill(mcend.Y());
836  mchists.fHmendz->Fill(mcend.Z());
837  mchists.fHmtheta->Fill(mcstartmom.Theta());
838  mchists.fHmphi->Fill(mcstartmom.Phi());
839  mchists.fHmtheta_xz->Fill(mctheta_xz);
840  mchists.fHmtheta_yz->Fill(mctheta_yz);
841  mchists.fHmmom->Fill(mcstartmom.Mag());
842  mchists.fHmlen->Fill(plen);
843  }
844  }
845  }
846  }
847 
848  // If we found at least one matched seed, fill good
849  // particle histograms.
850 
851  if (nmatch > 0) {
852  mchists.fHgstartx->Fill(mcstart.X());
853  mchists.fHgstarty->Fill(mcstart.Y());
854  mchists.fHgstartz->Fill(mcstart.Z());
855  mchists.fHgendx->Fill(mcend.X());
856  mchists.fHgendy->Fill(mcend.Y());
857  mchists.fHgendz->Fill(mcend.Z());
858  mchists.fHgtheta->Fill(mcstartmom.Theta());
859  mchists.fHgphi->Fill(mcstartmom.Phi());
860  mchists.fHgtheta_xz->Fill(mctheta_xz);
861  mchists.fHgtheta_yz->Fill(mctheta_yz);
862  mchists.fHgmom->Fill(mcstartmom.Mag());
863  mchists.fHglen->Fill(plen);
864  }
865  }
866  }
867  }
868  }
869  }
870  }
871 
872  // Loop over seeds and fill reco-only seed histograms.
873 
874  if (seedh.isValid()) {
875 
876  // Loop over seeds.
877 
878  int nseed = seedh->size();
879 
880  if (nseed > 0 && pdump != 0) {
881  *pdump << "\nReconstructed Seeds\n"
882  << " MCid x y z dx dy dz "
883  " p\n"
884  << "--------------------------------------------------------------------------------"
885  "-----------\n";
886  }
887 
888  for (int i = 0; i < nseed; ++i) {
889  art::Ptr<recob::Seed> pseed(seedh, i);
890  const recob::Seed& seed = *pseed;
891 
892  // Fill histograms involving reco seeds only.
893 
894  TVector3 pos;
895  TVector3 dir;
896  double err[3];
897  seed.GetPoint(&pos[0], err);
898  seed.GetDirection(&dir[0], err);
899  double mdir = dir.Mag();
900  if (mdir != 0.) { dir *= (1. / mdir); }
901 
902  double dpos = bdist(pos);
903  double theta_xz = std::atan2(dir.X(), dir.Z());
904  double theta_yz = std::atan2(dir.Y(), dir.Z());
905 
906  // Dump seed information here.
907 
908  if (pdump) {
909  int mcid = seedmap[&seed];
910  *pdump << std::setw(15) << mcid << " " << std::fixed << std::setprecision(2)
911  << std::setw(10) << pos[0] << std::setw(10) << pos[1] << std::setw(10) << pos[2]
912  << " " << std::fixed << std::setprecision(3) << std::setw(10) << dir[0]
913  << std::setw(10) << dir[1] << std::setw(10) << dir[2] << "\n";
914  }
915 
916  // Fill histograms.
917 
918  if (fRecoHistMap.count(0) == 0) fRecoHistMap.emplace(0, RecoHists{"Reco"});
919  const RecoHists& rhists = fRecoHistMap.at(0);
920 
921  rhists.fHx->Fill(pos.X());
922  rhists.fHy->Fill(pos.Y());
923  rhists.fHz->Fill(pos.Z());
924  rhists.fHdist->Fill(dpos);
925  rhists.fHtheta->Fill(dir.Theta());
926  rhists.fHphi->Fill(dir.Phi());
927  rhists.fHtheta_xz->Fill(theta_xz);
928  rhists.fHtheta_yz->Fill(theta_yz);
929  }
930  }
931  }
932 
933  void
934  SeedAna::endJob()
935  //
936  // Purpose: End of job.
937  //
938  {
939  // Print summary.
940 
941  mf::LogInfo("SeedAna") << "SeedAna statistics:\n"
942  << " Number of events = " << fNumEvent;
943 
944  // Fill multiplicity histograms.
945 
946  for (std::map<int, MCHists>::const_iterator i = fMCHistMap.begin(); i != fMCHistMap.end();
947  ++i) {
948  const MCHists& mchists = i->second;
949  mulcalc(mchists.fHmstartx, mchists.fHmcstartx, mchists.fHmulstartx);
950  mulcalc(mchists.fHmstarty, mchists.fHmcstarty, mchists.fHmulstarty);
951  mulcalc(mchists.fHmstartz, mchists.fHmcstartz, mchists.fHmulstartz);
952  mulcalc(mchists.fHmendx, mchists.fHmcendx, mchists.fHmulendx);
953  mulcalc(mchists.fHmendy, mchists.fHmcendy, mchists.fHmulendy);
954  mulcalc(mchists.fHmendz, mchists.fHmcendz, mchists.fHmulendz);
955  mulcalc(mchists.fHmtheta, mchists.fHmctheta, mchists.fHmultheta);
956  mulcalc(mchists.fHmphi, mchists.fHmcphi, mchists.fHmulphi);
957  mulcalc(mchists.fHmtheta_xz, mchists.fHmctheta_xz, mchists.fHmultheta_xz);
958  mulcalc(mchists.fHmtheta_yz, mchists.fHmctheta_yz, mchists.fHmultheta_yz);
959  mulcalc(mchists.fHmmom, mchists.fHmcmom, mchists.fHmulmom);
960  mulcalc(mchists.fHmlen, mchists.fHmclen, mchists.fHmullen);
961  }
962  // Fill efficiency histograms.
963 
964  for (std::map<int, MCHists>::const_iterator i = fMCHistMap.begin(); i != fMCHistMap.end();
965  ++i) {
966  const MCHists& mchists = i->second;
967  effcalc(mchists.fHgstartx, mchists.fHmcstartx, mchists.fHestartx);
968  effcalc(mchists.fHgstarty, mchists.fHmcstarty, mchists.fHestarty);
969  effcalc(mchists.fHgstartz, mchists.fHmcstartz, mchists.fHestartz);
970  effcalc(mchists.fHgendx, mchists.fHmcendx, mchists.fHeendx);
971  effcalc(mchists.fHgendy, mchists.fHmcendy, mchists.fHeendy);
972  effcalc(mchists.fHgendz, mchists.fHmcendz, mchists.fHeendz);
973  effcalc(mchists.fHgtheta, mchists.fHmctheta, mchists.fHetheta);
974  effcalc(mchists.fHgphi, mchists.fHmcphi, mchists.fHephi);
975  effcalc(mchists.fHgtheta_xz, mchists.fHmctheta_xz, mchists.fHetheta_xz);
976  effcalc(mchists.fHgtheta_yz, mchists.fHmctheta_yz, mchists.fHetheta_yz);
977  effcalc(mchists.fHgmom, mchists.fHmcmom, mchists.fHemom);
978  effcalc(mchists.fHglen, mchists.fHmclen, mchists.fHelen);
979  }
980  }
981 }
pmtana::std
double std(const std::vector< short > &wf, const double ped_mean, size_t start, size_t nsample)
Definition: UtilFunc.cxx:42
pdg
var pdg
Definition: selectors.fcl:14
trkf::SeedAna::fMCTrackModuleLabel
std::string fMCTrackModuleLabel
Definition: SeedAna_module.cc:351
evd::err
EResult err(const char *call)
Definition: WebEVDServer.cxx:140
recob::Seed::GetPoint
void GetPoint(double *Pt, double *Err) const
Definition: Seed.cxx:108
get
get()
Definition: PMTconfigurationExtraction_module.cc:168
trkf::SeedAna::endJob
void endJob() override
Definition: SeedAna_module.cc:934
trkf::SeedAna::fSeedModuleLabel
std::string fSeedModuleLabel
Definition: SeedAna_module.cc:350
larg4::bool
static constexpr bool
Definition: OpFastScintillation.hh:588
sim::MCStep::T
double T() const
Definition: MCStep.h:45
util::quantities::ticks
tick ticks
Alias for common language habits.
Definition: electronics.h:78
trkf::SeedAna::MCHists::MCHists
MCHists(const std::string &subdir)
Definition: SeedAna_module.cc:404
trkf::SeedAna::analyze
void analyze(const art::Event &evt) override
Definition: SeedAna_module.cc:574
xmin
process_name pandoraGausCryo1 vertexChargeCryo1 vertexStubCryo1 xmin
Definition: run_vertexcharge_icarus.fcl:109
trkf::SeedAna::fMatchColinearity
double fMatchColinearity
Definition: SeedAna_module.cc:355
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
icarus::ns::util::end
auto end(FixedBins< T, C > const &) noexcept
Definition: FixedBins.h:585
detinfo::DetectorPropertiesData::ConvertXToTicks
double ConvertXToTicks(double X, int p, int t, int c) const
Definition: DetectorPropertiesData.cc:99
gen_crt_frags.pv
tuple pv
Definition: gen_crt_frags.py:1419
seed
unsigned int seed
Definition: NeutronOsc_module.cc:129
detinfo::DetectorPropertiesData::DriftVelocity
double DriftVelocity(double efield=0., double temperature=0.) const
cm/us
Definition: DetectorPropertiesData.cc:23
MCTrack.h
Class def header for mctrack data container.
trkf::SeedAna::fRecoHistMap
std::map< int, RecoHists > fRecoHistMap
Definition: SeedAna_module.cc:362
gen_crt_frags.w
tuple w
Definition: gen_crt_frags.py:1416
sim::MCTrack::PdgCode
int PdgCode() const
Definition: MCTrack.h:41
sim::MCStep::Momentum
const TLorentzVector & Momentum() const
Definition: MCStep.h:38
dropbox.dir
tuple dir
Definition: dropbox.py:28
trkf::SeedAna::SeedAna
SeedAna(fhicl::ParameterSet const &pset)
Definition: SeedAna_module.cc:548
trkf::SeedAna::fMCHistMap
std::map< int, MCHists > fMCHistMap
Definition: SeedAna_module.cc:361
larg4::dist
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Definition: OpFastScintillation.hh:563
sim::MCStep::E
double E() const
Definition: MCStep.h:49
sim::MCTrack::Start
const MCStep & Start() const
Definition: MCTrack.h:44
tfs
art::ServiceHandle< art::TFileService > tfs
Definition: NuVertexChargeTree_module.cc:254
sim::MCTrack::TrackID
unsigned int TrackID() const
Definition: MCTrack.h:42
tca::evt
TCEvent evt
Definition: DataStructs.cxx:8
recob::Seed::GetDirection
void GetDirection(double *Dir, double *Err) const
Definition: Seed.cxx:98
Geometry.h
art framework interface to geometry description
detsim::detProp
auto const detProp
Definition: SimWire_module.cc:157
trkf::SeedAna::RecoHists::RecoHists
RecoHists(const std::string &subdir)
Definition: SeedAna_module.cc:373
Generated by   doxygen 1.8.5