LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
SpacePointAna_module.cc
Go to the documentation of this file.
1 //
2 // Name: SpacePointAna_module.cc
3 //
4 // Purpose: Module SpacePointAna.
5 //
6 //
7 // Configuration parameters.
8 //
9 // HitModuleLabel: Hit module label.
10 // UseClusterHits: If true, use clustered hits (otherwise all hits).
11 // ClusterModuleLabel: Cluster module label.
12 // UseMC: Use MC truth information.
13 // SpacePointAlgTime: SpacePointAlg configuration (loose time cuts).
14 // SpacePointAlgSep: SpacePointAlg configuration (loose separation cuts).
15 // SpacePointAlgDefault: SpacePointAlg configuration (default cuts).
16 //
17 // Created: 2-Aug-2011 H. Greenlee
18 //
19 
24 #include "art_root_io/TFileService.h"
27 
39 
40 #include "TH1F.h"
41 #include "TH2F.h"
42 
43 #include "range/v3/view.hpp"
44 
45 using lar::to_element;
46 using ranges::views::transform;
47 
48 namespace trkf {
49 
50  class SpacePointAna : public art::EDAnalyzer {
51  public:
52  explicit SpacePointAna(fhicl::ParameterSet const& pset);
53 
54  private:
55  void analyze(const art::Event& evt) override;
56 
57  void bookHistograms(bool mc);
58 
59  // Fcl Attributes.
60 
61  const SpacePointAlg fSptalgTime; // Algorithm object (increased time cut).
62  const SpacePointAlg fSptalgSep; // Algorithm object (increased sepataion cut).
63  const SpacePointAlg fSptalgDefault; // Algorithm object (default cuts).
64  std::string fHitModuleLabel;
66  std::string fClusterModuleLabel;
67  bool fUseMC;
68  double fMinX; // Minimum x.
69  double fMaxX; // Maximum x.
70  double fMinY; // Minimum y.
71  double fMaxY; // Maximum y.
72  double fMinZ; // Minimum z.
73  double fMaxZ; // Maximum z.
74 
75  // Histograms.
76 
77  bool fBooked; // Have histograms been booked yet?
78  TH1F* fHDTUE; // U-drift electrons time difference.
79  TH1F* fHDTVE; // V-drift electrons time difference.
80  TH1F* fHDTWE; // W-drift electrons time difference.
81  TH1F* fHDTUPull; // U-drift electrons time pull.
82  TH1F* fHDTVPull; // V-drift electrons time pull.
83  TH1F* fHDTWPull; // W-drift electrons time pull.
84  TH1F* fHDTUV; // U-V time difference.
85  TH1F* fHDTVW; // V-W time difference.
86  TH1F* fHDTWU; // W-U time difference.
87  TH2F* fHDTUVU; // U-V time difference vs. U.
88  TH2F* fHDTUVV; // U-V time difference vs. V.
89  TH2F* fHDTVWV; // V-W time difference vs. V.
90  TH2F* fHDTVWW; // V-W time difference vs. W.
91  TH2F* fHDTWUW; // W-U time difference vs. W.
92  TH2F* fHDTWUU; // W-U time difference vs. U.
93  TH1F* fHS; // Spatial separation.
94  TH1F* fHchisq; // Space point chisquare.
95  TH1F* fHx; // X position.
96  TH1F* fHy; // Y position.
97  TH1F* fHz; // Z position.
98  TH1F* fHAmpU; // U hit amplitude.
99  TH1F* fHAmpV; // V hit amplitude.
100  TH1F* fHAmpW; // W hit amplitude.
101  TH1F* fHAreaU; // U hit area.
102  TH1F* fHAreaV; // V hit area.
103  TH1F* fHAreaW; // W hit area.
104  TH1F* fHSumU; // U hit sum ADC.
105  TH1F* fHSumV; // V hit sum ADC.
106  TH1F* fHSumW; // W hit sum ADC.
107  TH1F* fHMCdx; // X residual (reco vs. mc truth).
108  TH1F* fHMCdy; // Y residual (reco vs. mc truth).
109  TH1F* fHMCdz; // Z residual (reco vs. mc truth).
110  TH1F* fHMCxpull; // X pull (reco vs. mc truth).
111  TH1F* fHMCypull; // Y pull (reco vs. mc truth).
112  TH1F* fHMCzpull; // Z pull (reco vs. mc truth).
113 
114  // Statistics.
115 
117  };
118 
120 
121  SpacePointAna::SpacePointAna(const fhicl::ParameterSet& pset)
122  //
123  // Purpose: Constructor.
124  //
125  // Arguments: pset - Module parameters.
126  //
127  : EDAnalyzer(pset)
128  , fSptalgTime(pset.get<fhicl::ParameterSet>("SpacePointAlgTime"))
129  , fSptalgSep(pset.get<fhicl::ParameterSet>("SpacePointAlgSep"))
130  , fSptalgDefault(pset.get<fhicl::ParameterSet>("SpacePointAlgDefault"))
131  , fHitModuleLabel(pset.get<std::string>("HitModuleLabel"))
132  , fUseClusterHits(pset.get<bool>("UseClusterHits"))
133  , fClusterModuleLabel(pset.get<std::string>("ClusterModuleLabel"))
134  , fUseMC(pset.get<bool>("UseMC"))
135  , fMinX(pset.get<double>("MinX", -1.e10))
136  , fMaxX(pset.get<double>("MaxX", 1.e10))
137  , fMinY(pset.get<double>("MinY", -1.e10))
138  , fMaxY(pset.get<double>("MaxY", 1.e10))
139  , fMinZ(pset.get<double>("MinZ", -1.e10))
140  , fMaxZ(pset.get<double>("MaxZ", 1.e10))
141  , fBooked(false)
142  , fHDTUE(0)
143  , fHDTVE(0)
144  , fHDTWE(0)
145  , fHDTUPull(0)
146  , fHDTVPull(0)
147  , fHDTWPull(0)
148  , fHDTUV(0)
149  , fHDTVW(0)
150  , fHDTWU(0)
151  , fHDTUVU(0)
152  , fHDTUVV(0)
153  , fHDTVWV(0)
154  , fHDTVWW(0)
155  , fHDTWUW(0)
156  , fHDTWUU(0)
157  , fHS(0)
158  , fHchisq(0)
159  , fHx(0)
160  , fHy(0)
161  , fHz(0)
162  , fHAmpU(0)
163  , fHAmpV(0)
164  , fHAmpW(0)
165  , fHAreaU(0)
166  , fHAreaV(0)
167  , fHAreaW(0)
168  , fHSumU(0)
169  , fHSumV(0)
170  , fHSumW(0)
171  , fHMCdx(0)
172  , fHMCdy(0)
173  , fHMCdz(0)
174  , fHMCxpull(0)
175  , fHMCypull(0)
176  , fHMCzpull(0)
177  , fNumEvent(0)
178  {
179  mf::LogInfo("SpacePointAna") << "SpacePointAna configured with the following parameters:\n"
180  << " HitModuleLabel = " << fHitModuleLabel << "\n"
181  << " UseClusterHits = " << fUseClusterHits << "\n"
182  << " ClusterModuleLabel = " << fClusterModuleLabel << "\n"
183  << " UseMC = " << fUseMC;
184  }
185 
187  {
188  if (!fBooked) {
189  fBooked = true;
190 
193  art::TFileDirectory dir = tfs->mkdir("sptana", "SpacePointAna histograms");
194 
195  unsigned int nwiresU = 0, nwiresV = 0, nwiresW = 0;
196 
197  // Figure out the number of wires in U, V, and W planes.
198 
199  // Loop over cryostats, tpcs, and planes.
200 
201  for (auto const& pgeom : geom->Iterate<geo::PlaneGeo>()) {
202  unsigned int const nwires = pgeom.Nwires();
203  switch (pgeom.View()) {
204  case geo::kU: nwiresU = nwires; break;
205  case geo::kV: nwiresV = nwires; break;
206  case geo::kZ: nwiresW = nwires; break;
207  default: {
208  }
209  }
210  }
211 
212  if (mc && fUseMC) {
213  fHDTUE = dir.make<TH1F>("MCDTUE", "U-Drift Electrons Time Difference", 100, -5., 5.);
214  fHDTVE = dir.make<TH1F>("MCDTVE", "V-Drift Electrons Time Difference", 100, -5., 5.);
215  fHDTWE = dir.make<TH1F>("MCDTWE", "W-Drift Electrons Time Difference", 100, -5., 5.);
216  fHDTUPull = dir.make<TH1F>("MCDTUPull", "U-Drift Electrons Time Pull", 100, -50., 50.);
217  fHDTVPull = dir.make<TH1F>("MCDTVPull", "V-Drift Electrons Time Pull", 100, -50., 50.);
218  fHDTWPull = dir.make<TH1F>("MCDTWPull", "W-Drift Electrons Time Pull", 100, -50., 50.);
219  }
220  if (!fSptalgTime.merge()) {
221  fHDTUV = dir.make<TH1F>("DTUV", "U-V time difference", 100, -20., 20.);
222  fHDTVW = dir.make<TH1F>("DTVW", "V-W time difference", 100, -20., 20.);
223  fHDTWU = dir.make<TH1F>("DTWU", "W-U time difference", 100, -20., 20.);
224  fHDTUVU = dir.make<TH2F>(
225  "DTUVU", "U-V time difference vs. U", 100, 0., double(nwiresU), 100, -20., 20.);
226  fHDTUVV = dir.make<TH2F>(
227  "DTUVV", "U-V time difference vs. V", 100, 0., double(nwiresV), 100, -20., 20.);
228  fHDTVWV = dir.make<TH2F>(
229  "DTVWV", "V-W time difference vs. V", 100, 0., double(nwiresV), 100, -20., 20.);
230  fHDTVWW = dir.make<TH2F>(
231  "DTVWW", "V-W time difference vs. W", 100, 0., double(nwiresW), 100, -20., 20.);
232  fHDTWUW = dir.make<TH2F>(
233  "DTWUW", "W-U time difference vs. W", 100, 0., double(nwiresW), 100, -20., 20.);
234  fHDTWUU = dir.make<TH2F>(
235  "DTWUU", "W-U time difference vs. U", 100, 0., double(nwiresU), 100, -20., 20.);
236  fHS = dir.make<TH1F>("DS", "Spatial Separation", 100, -2., 2.);
237  }
238  fHchisq = dir.make<TH1F>("chisq", "Chisquare", 100, 0., 20.);
239 
240  fHx = dir.make<TH1F>("xpos", "X Position", 100, 0., 2. * geom->DetHalfWidth());
241  fHy =
242  dir.make<TH1F>("ypos", "Y Position", 100, -geom->DetHalfHeight(), geom->DetHalfHeight());
243  fHz = dir.make<TH1F>("zpos", "Z Position", 100, 0., geom->DetLength());
244  fHAmpU = dir.make<TH1F>("ampU", "U Hit Amplitude", 50, 0., 50.);
245  fHAmpV = dir.make<TH1F>("ampV", "V Hit Amplitude", 50, 0., 50.);
246  fHAmpW = dir.make<TH1F>("ampW", "W Hit Amplitude", 50, 0., 50.);
247  fHAreaU = dir.make<TH1F>("areaU", "U Hit Area", 100, 0., 500.);
248  fHAreaV = dir.make<TH1F>("areaV", "V Hit Area", 100, 0., 500.);
249  fHAreaW = dir.make<TH1F>("areaW", "W Hit Area", 100, 0., 500.);
250  fHSumU = dir.make<TH1F>("sumU", "U Hit Sum ADC", 100, 0., 500.);
251  fHSumV = dir.make<TH1F>("sumV", "V Hit Sum ADC", 100, 0., 500.);
252  fHSumW = dir.make<TH1F>("sumW", "W Hit Sum ADC", 100, 0., 500.);
253  if (mc && fUseMC) {
254  fHMCdx = dir.make<TH1F>("MCdx", "X MC Residual", 100, -1., 1.);
255  fHMCdy = dir.make<TH1F>("MCdy", "Y MC Residual", 100, -1., 1.);
256  fHMCdz = dir.make<TH1F>("MCdz", "Z MC Residual", 100, -1., 1.);
257  fHMCxpull = dir.make<TH1F>("MCxpull", "X MC Pull", 100, -50., 50.);
258  fHMCypull = dir.make<TH1F>("MCypull", "Y MC Pull", 100, -50., 50.);
259  fHMCzpull = dir.make<TH1F>("MCzpull", "Z MC Pull", 100, -50., 50.);
260  }
261  }
262  }
263 
265  {
266  ++fNumEvent;
267 
268  // Make sure histograms are booked.
269 
270  bool mc = !evt.isRealData();
271  bookHistograms(mc);
272 
273  // Get Services.
274 
276 
277  // Get Hits.
278 
280 
281  if (fUseClusterHits) {
282 
283  // Get clusters.
284 
286  evt.getByLabel(fClusterModuleLabel, clusterh);
287 
288  // Get hits from all clusters.
290 
291  if (clusterh.isValid()) {
292  int nclus = clusterh->size();
293 
294  for (int i = 0; i < nclus; ++i) {
295  art::Ptr<recob::Cluster> pclus(clusterh, i);
296  std::vector<art::Ptr<recob::Hit>> clushits = fm.at(i);
297  int nhits = clushits.size();
298  hits.reserve(hits.size() + nhits);
299 
300  for (std::vector<art::Ptr<recob::Hit>>::const_iterator ihit = clushits.begin();
301  ihit != clushits.end();
302  ++ihit) {
303  hits.push_back(*ihit);
304  }
305  }
306  }
307  }
308  else {
309 
310  // Get unclustered hits.
311 
313  evt.getByLabel(fHitModuleLabel, hith);
314  if (hith.isValid()) {
315  int nhits = hith->size();
316  hits.reserve(nhits);
317 
318  for (int i = 0; i < nhits; ++i)
319  hits.push_back(art::Ptr<recob::Hit>(hith, i));
320  }
321  }
322 
323  // Fill histograms that don't depend on space points.
324 
325  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(evt);
326  auto const detProp =
328 
329  if (mc && fUseMC) {
330 
332 
333  // Loop over hits and fill hit-electron time difference histogram.
334 
335  for (auto const& hitPtr : hits) {
336  const recob::Hit& hit = *hitPtr;
337 
338  double tpeak = hit.PeakTime();
339  double terr = hit.SigmaPeakTime();
340 
341  bool tav_ok = true;
342  double tav = 0.;
343  try {
344  std::vector<double> hitxyz = bt_serv->HitToXYZ(clockData, hit);
345  tav = detProp.ConvertXToTicks(
346  hitxyz[0], hit.WireID().Plane, hit.WireID().TPC, hit.WireID().Cryostat);
347  }
348  catch (cet::exception& x) {
349  tav_ok = false;
350  }
351  if (tav_ok) {
352  if (hit.View() == geo::kU) {
353  fHDTUE->Fill(tpeak - tav);
354  fHDTUPull->Fill((tpeak - tav) / terr);
355  }
356  else if (hit.View() == geo::kV) {
357  fHDTVE->Fill(tpeak - tav);
358  fHDTVPull->Fill((tpeak - tav) / terr);
359  }
360  else if (hit.View() == geo::kZ) {
361  fHDTWE->Fill(tpeak - tav);
362  fHDTWPull->Fill((tpeak - tav) / terr);
363  }
364  else
365  throw cet::exception("SpacePointAna") << "Bad view = " << hit.View() << "\n";
366  }
367  }
368  }
369 
370  std::vector<recob::SpacePoint> spts1; // For time histograms.
371  std::vector<recob::SpacePoint> spts2; // For separation histogram.
372  std::vector<recob::SpacePoint> spts3; // Default cuts.
373 
374  // If nonzero time cut is specified, make space points using that
375  // time cut (for time histograms).
376 
377  if (!fSptalgTime.merge()) {
378  if (mc && fUseMC)
379  fSptalgTime.makeMCTruthSpacePoints(clockData, detProp, hits, spts1);
380  else
381  fSptalgTime.makeSpacePoints(clockData, detProp, hits, spts1);
382 
383  // Report number of space points.
384 
385  MF_LOG_DEBUG("SpacePointAna")
386  << "Found " << spts1.size() << " space points using special time cut.";
387  }
388 
389  // If nonzero separation cut is specified, make space points using that
390  // separation cut (for separation histogram).
391 
392  if (!fSptalgSep.merge()) {
393  if (mc && fUseMC)
394  fSptalgSep.makeMCTruthSpacePoints(clockData, detProp, hits, spts2);
395  else
396  fSptalgSep.makeSpacePoints(clockData, detProp, hits, spts2);
397 
398  // Report number of space points.
399 
400  MF_LOG_DEBUG("SpacePointAna")
401  << "Found " << spts2.size() << " space points using special seperation cut.";
402  }
403 
404  // Make space points using default cuts.
405 
406  if (mc && fUseMC)
407  fSptalgDefault.makeMCTruthSpacePoints(clockData, detProp, hits, spts3);
408  else
409  fSptalgDefault.makeSpacePoints(clockData, detProp, hits, spts3);
410 
411  // Report number of space points.
412 
413  MF_LOG_DEBUG("SpacePointAna") << "Found " << spts3.size()
414  << " space points using default cuts.";
415 
416  if (!fSptalgTime.merge()) {
417 
418  // Loop over space points and fill time histograms.
419 
420  for (auto const& spt : spts1) {
421  if (spt.XYZ()[0] < fMinX || spt.XYZ()[0] > fMaxX || spt.XYZ()[1] < fMinY ||
422  spt.XYZ()[1] > fMaxY || spt.XYZ()[2] < fMinZ || spt.XYZ()[2] > fMaxZ)
423  continue;
424 
425  // Get hits associated with this SpacePoint.
426 
428 
429  // Make a double loop over hits and fill hit time difference histograms.
430 
431  for (auto const& hit1 : spthits | transform(to_element)) {
432  geo::WireID hit1WireID = hit1.WireID();
433  unsigned int tpc1 = hit1WireID.TPC;
434  unsigned int plane1 = hit1WireID.Plane;
435  unsigned int wire1 = hit1WireID.Wire;
436 
437  geo::View_t view1 = hit1.View();
438  double t1 = fSptalgTime.correctedTime(detProp, hit1);
439 
440  for (auto const& hit2 : spthits | transform(to_element)) {
441  geo::WireID hit2WireID = hit2.WireID();
442  unsigned int tpc2 = hit2WireID.TPC;
443  unsigned int plane2 = hit2WireID.Plane;
444  unsigned int wire2 = hit2WireID.Wire;
445 
446  // Require same tpc, different view.
447 
448  if (tpc1 == tpc2 && plane1 != plane2) {
449 
450  geo::View_t view2 = hit2.View();
451  double t2 = fSptalgTime.correctedTime(detProp, hit2);
452 
453  if (view1 == geo::kU) {
454  if (view2 == geo::kV) {
455  fHDTUV->Fill(t1 - t2);
456  fHDTUVU->Fill(double(wire1), t1 - t2);
457  fHDTUVV->Fill(double(wire2), t1 - t2);
458  }
459  if (view2 == geo::kZ) {
460  fHDTWU->Fill(t2 - t1);
461  fHDTWUW->Fill(double(wire2), t2 - t1);
462  fHDTWUU->Fill(double(wire1), t2 - t1);
463  }
464  }
465  if (view1 == geo::kV) {
466  if (view2 == geo::kZ) {
467  fHDTVW->Fill(t1 - t2);
468  fHDTVWV->Fill(double(wire1), t1 - t2);
469  fHDTVWW->Fill(double(wire2), t1 - t2);
470  }
471  if (view2 == geo::kU) {
472  fHDTUV->Fill(t2 - t1);
473  fHDTUVU->Fill(double(wire2), t2 - t1);
474  fHDTUVV->Fill(double(wire1), t2 - t1);
475  }
476  }
477  if (view1 == geo::kZ) {
478  if (view2 == geo::kU) {
479  fHDTWU->Fill(t1 - t2);
480  fHDTWUW->Fill(double(wire1), t1 - t2);
481  fHDTWUU->Fill(double(wire2), t1 - t2);
482  }
483  if (view2 == geo::kV) {
484  fHDTVW->Fill(t2 - t1);
485  fHDTVWV->Fill(double(wire2), t2 - t1);
486  fHDTVWW->Fill(double(wire1), t2 - t1);
487  }
488  }
489  }
490  }
491  }
492  }
493 
494  // Loop over space points and fill seperation histograms.
495 
496  for (auto const& spt : spts2) {
497  if (spt.XYZ()[0] < fMinX || spt.XYZ()[0] > fMaxX || spt.XYZ()[1] < fMinY ||
498  spt.XYZ()[1] > fMaxY || spt.XYZ()[2] < fMinZ || spt.XYZ()[2] > fMaxZ)
499  continue;
500 
501  // Get hits associated with this SpacePoint.
502 
504 
505  // Fill separation histogram.
506 
507  double sep = fSptalgSep.separation(spthits);
508  fHS->Fill(sep);
509  }
510  }
511 
512  // Loop over default space points and fill histograms.
513 
514  for (auto const& spt : spts3) {
515  if (spt.XYZ()[0] < fMinX || spt.XYZ()[0] > fMaxX || spt.XYZ()[1] < fMinY ||
516  spt.XYZ()[1] > fMaxY || spt.XYZ()[2] < fMinZ || spt.XYZ()[2] > fMaxZ)
517  continue;
518 
519  fHchisq->Fill(spt.Chisq());
520  fHx->Fill(spt.XYZ()[0]);
521  fHy->Fill(spt.XYZ()[1]);
522  fHz->Fill(spt.XYZ()[2]);
523 
524  // Get hits associated with this SpacePoint.
525 
526  std::vector<art::Ptr<recob::Hit>> spthits;
528  for (auto const& ptr : av_spthits) {
529  spthits.push_back(ptr);
530  }
531 
532  // Fill single hit histograms.
533 
534  for (auto const& hitPtr : spthits) {
535  const recob::Hit& hit = *hitPtr;
536 
537  geo::View_t view = hit.View();
538 
539  if (view == geo::kU) {
540  fHAmpU->Fill(hit.PeakAmplitude());
541  fHAreaU->Fill(hit.Integral());
542  fHSumU->Fill(hit.SummedADC());
543  }
544  if (view == geo::kV) {
545  fHAmpV->Fill(hit.PeakAmplitude());
546  fHAreaV->Fill(hit.Integral());
547  fHSumV->Fill(hit.SummedADC());
548  }
549  if (view == geo::kZ) {
550  fHAmpW->Fill(hit.PeakAmplitude());
551  fHAreaW->Fill(hit.Integral());
552  fHSumW->Fill(hit.SummedADC());
553  }
554  }
555 
556  if (mc && fUseMC) {
557 
559 
560  try {
561  std::vector<double> mcxyz = bt_serv->SpacePointHitsToWeightedXYZ(clockData, spthits);
562  fHMCdx->Fill(spt.XYZ()[0] - mcxyz[0]);
563  fHMCdy->Fill(spt.XYZ()[1] - mcxyz[1]);
564  fHMCdz->Fill(spt.XYZ()[2] - mcxyz[2]);
565  if (spt.ErrXYZ()[0] > 0.)
566  fHMCxpull->Fill((spt.XYZ()[0] - mcxyz[0]) / std::sqrt(spt.ErrXYZ()[0]));
567  if (spt.ErrXYZ()[2] > 0.)
568  fHMCypull->Fill((spt.XYZ()[1] - mcxyz[1]) / std::sqrt(spt.ErrXYZ()[2]));
569  if (spt.ErrXYZ()[5] > 0.)
570  fHMCzpull->Fill((spt.XYZ()[2] - mcxyz[2]) / std::sqrt(spt.ErrXYZ()[5]));
571  }
572  catch (cet::exception&) {
573  }
574  }
575  }
576  }
577 }
Float_t x
Definition: compare.C:6
details::range_type< T > Iterate() const
Initializes the specified ID with the ID of the first cryostat.
Definition: GeometryCore.h:541
void reserve(size_type n)
Definition: PtrVector.h:337
constexpr to_element_t to_element
Definition: ToElement.h:25
TTree * t1
Definition: plottest35.C:26
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
Planes which measure V.
Definition: geo_types.h:136
Declaration of signal hit object.
Length_t DetHalfWidth(TPCID const &tpcid=tpc_zero) const
Returns the half width of the active volume of the specified TPC.
CryostatID_t Cryostat
Index of cryostat.
Definition: geo_types.h:211
STL namespace.
Planes which measure Z direction.
Definition: geo_types.h:138
float Integral() const
Integral under the calibrated signal waveform of the hit, in tick x ADC units.
Definition: Hit.h:244
geo::View_t View() const
View for the plane of the hit.
Definition: Hit.h:276
WireID_t Wire
Index of the wire within its plane.
Definition: geo_types.h:563
geo::WireID const & WireID() const
Initial tdc tick for hit.
Definition: Hit.h:280
bool isRealData() const
Definition: Event.cc:53
std::vector< double > HitToXYZ(detinfo::DetectorClocksData const &clockData, const recob::Hit &hit) const
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.cc:6
Length_t DetLength(TPCID const &tpcid=tpc_zero) const
Returns the length of the active volume of the specified TPC.
const art::PtrVector< recob::Hit > & getAssociatedHits(const recob::SpacePoint &spt) const
const SpacePointAlg fSptalgSep
float PeakAmplitude() const
The estimated amplitude of the hit at its peak, in ADC units.
Definition: Hit.h:232
void makeMCTruthSpacePoints(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, const art::PtrVector< recob::Hit > &hits, std::vector< recob::SpacePoint > &spts) const
double separation(const art::PtrVector< recob::Hit > &hits) const
bool isValid() const noexcept
Definition: Handle.h:203
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
Definition: DumpUtils.h:289
Planes which measure U.
Definition: geo_types.h:135
void hits()
Definition: readHits.C:15
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:435
parameter set interface
void makeSpacePoints(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, const art::PtrVector< recob::Hit > &hits, std::vector< recob::SpacePoint > &spts) const
Geometry information for a single wire plane.The plane is represented in the geometry by a solid whic...
Definition: PlaneGeo.h:78
void analyze(const art::Event &evt) override
bool merge() const noexcept
Definition: SpacePointAlg.h:87
PlaneID_t Plane
Index of the plane within its TPC.
Definition: geo_types.h:481
double correctedTime(detinfo::DetectorPropertiesData const &detProp, const recob::Hit &hit) const
Declaration of cluster object.
size_type size() const
Definition: PtrVector.h:302
Detector simulation of raw signals on wires.
TTree * t2
Definition: plottest35.C:36
float PeakTime() const
Time of the signal peak, in tick units.
Definition: Hit.h:220
decltype(auto) get(T &&obj)
ADL-aware version of std::to_string.
Definition: StdUtils.h:120
SpacePointAna(fhicl::ParameterSet const &pset)
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Encapsulate the construction of a single detector plane.
const SpacePointAlg fSptalgTime
TDirectory * dir
Definition: macro.C:5
std::vector< double > SpacePointHitsToWeightedXYZ(detinfo::DetectorClocksData const &clockData, std::vector< art::Ptr< recob::Hit >> const &hits) const
#define MF_LOG_DEBUG(id)
float SummedADC() const
The sum of calibrated ADC counts of the hit (0. by default)
Definition: Hit.h:240
const SpacePointAlg fSptalgDefault
constexpr WireID()=default
Default constructor: an invalid TPC ID.
float SigmaPeakTime() const
Uncertainty for the signal peak, in tick units.
Definition: Hit.h:224
Length_t DetHalfHeight(TPCID const &tpcid=tpc_zero) const
Returns the half height of the active volume of the specified TPC.
2D representation of charge deposited in the TDC/wire plane
Definition: Hit.h:46
TCEvent evt
Definition: DataStructs.cxx:8
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:399
Algorithm for generating space points from hits.
art framework interface to geometry description
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
Encapsulate the construction of a single detector plane.