ToFSingleCounterAnalysis_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: ToFSingleCounterAnalysis
3 // Plugin Type: producer (art v2_12_01)
4 // File: ToFSingleCounterAnalysis_module.cc
5 //
6 // Aidan Medcalf (University of Dallas) <amedcalf@udallas.edu>
7 // 10 July 2019
8 //
9 // Description: Single-counter ToF statistics
10 ////////////////////////////////////////////////////////////////////////
11 
12 #include "art/Framework/Core/EDAnalyzer.h"
13 #include "art/Framework/Core/ModuleMacros.h"
14 #include "art/Framework/Principal/Event.h"
15 #include "art/Framework/Principal/Handle.h"
16 #include "art/Framework/Principal/Run.h"
17 #include "art/Framework/Principal/SubRun.h"
18 #include "canvas/Utilities/InputTag.h"
19 #include "fhiclcpp/ParameterSet.h"
20 #include "messagefacility/MessageLogger/MessageLogger.h"
21 #include "canvas/Persistency/Common/Ptr.h"
22 #include "canvas/Persistency/Common/PtrVector.h"
23 #include "art/Framework/Services/Registry/ServiceHandle.h"
24 #include "art/Framework/Services/Optional/TFileService.h"
25 #include "art/Framework/Services/Optional/TFileDirectory.h"
26 
27 #include "BeamlineRecoBase/BeamlineDigit.h"
28 #include "Utilities/AssociationUtil.h"
29 #include "BeamlineUtils/BeamlineChannelMap.h"
30 #include "BeamlineReco/ToFClusterAlg.h"
31 
32 #include "TH1D.h"
33 #include "TH2D.h"
34 
35 #include <iostream>
36 //#include "ToFPositionFinder.h"
37 
38 namespace beamlinereco {
39  class ToFSingleCounterAnalysis : public art::EDAnalyzer {
40 
41  public:
42  ToFSingleCounterAnalysis(fhicl::ParameterSet const &p);
43 
44  void reconfigure(const fhicl::ParameterSet &p);
45  void beginJob() override;
46  void analyze(const art::Event &e);
47 
48  private:
49  double calculateTime(ToFPulseCluster c, std::string method) const;
50 
51  // FHiCl parameters
52  art::InputTag fToFDigitLabel;
53  double fSamplingInterval;
54  double fHitClusterWindow;
55  std::string fClusterTimeMethod;
56 
57  double fHistdtBound;
58  double fHistdABound;
59  double fHistdIBound;
60  unsigned int fHistdtNBins;
61  unsigned int fHistdANBins;
62  unsigned int fHistdINBins;
63  unsigned int fHistdAdtNBins;
64  unsigned int fHistdIdtNBins;
65  unsigned int fHistdAdINBins;
66 
67  // Histograms
68 
69  std::map<beamlinegeo::ChannelID, TH1D*> hAmplitude_AllHits_channel;
70  std::map<beamlinegeo::ChannelID, TH1D*> hIntegral_AllHits_channel;
71  std::map<beamlinegeo::ChannelID, TH1D*> hStartTime_AllHits_channel;
72  std::map<beamlinegeo::ChannelID, TH1D*> hWidth_AllHits_channel;
73  std::map<beamlinegeo::ChannelID, TH2D*> hAmplitude_vs_Integral_AllHits_channel;
74  std::map<beamlinegeo::ChannelID, TH2D*> hAmplitude_vs_Width_AllHits_channel;
75 
76  TH1I *hUSHitMultiplicity;
77  TH1I *hDSHitMultiplicity;
78  TH1I *hSiPMHitMultiplicity;
79 
80  TH1I *hUSFirstHitDetector;
81  TH1I *hDSFirstHitDetector;
82  TH1I *hSiPMFirstHitDetector;
83 
84  TH1D *hUSClusterTime;
85  TH1D *hDSClusterTime;
86  TH1D *hSiPMClusterTime;
87 
88  TH1D *hUSdt[4][4];
89  TH1D *hUSdA[4][4];
90  TH1D *hUSdI[4][4];
91  TH2D *hUSdAdt[4][4];
92  TH2D *hUSdIdt[4][4];
93  TH2D *hUSdAdI[4][4];
94 
95  TH1D *hDSdt[4][4];
96  TH1D *hDSdA[4][4];
97  TH1D *hDSdI[4][4];
98  TH2D *hDSdAdt[4][4];
99  TH2D *hDSdIdt[4][4];
100  TH2D *hDSdAdI[4][4];
101 
102  TH1D *hSiPMdt[4][4];
103  TH1D *hSiPMdA[4][4];
104  TH1D *hSiPMdI[4][4];
105  TH2D *hSiPMdAdt[4][4];
106  TH2D *hSiPMdIdt[4][4];
107  TH2D *hSiPMdAdI[4][4];
108 
109  // Services
110 
111  art::ServiceHandle<beamlineutil::BeamlineChannelMap> fChannelMap;
112  art::ServiceHandle<art::TFileService> tfs;
113 
114  };
115 }
116 
117 static inline unsigned int onlineChannel(ChannelID id) {
118  if (id.System == beamlinegeo::BeamlineSystem::ToF) {
119  switch (id.Detector) {
120  case beamlinegeo::ToFCounter::US: return 8 + id.Channel;
121  case beamlinegeo::ToFCounter::DS: return 12 + id.Channel;
122  case beamlinegeo::ToFCounter::DSSiPM: return 4 + id.Channel;
123  default: return 99999;
124  }
125  } else if (id.System == beamlinegeo::BeamlineSystem::Cherenkov) {
126  return id.Channel;
127  }
128  return 99999;
129 }
130 
131 static inline std::vector<art::Ptr<brb::BeamlineDigit>>::iterator getEarliest(std::vector<art::Ptr<brb::BeamlineDigit>> v) {
132  auto dit = std::min_element(v.begin(), v.end(),
133  [](art::Ptr<brb::BeamlineDigit> p1, art::Ptr<brb::BeamlineDigit> p2){
134  return p1->StartTimeInNanoSec() < p2->StartTimeInNanoSec();
135  });
136  return dit;
137 }
138 
139 // Calculate single time for a cluster
140 double beamlinereco::ToFSingleCounterAnalysis::calculateTime(ToFPulseCluster c, std::string method) const {
141  if (c.pulses.size() == 0) {
142  return std::numeric_limits<double>::quiet_NaN();
143  }
144 
145  double t = std::numeric_limits<double>::quiet_NaN();
146 
147  if (method == "Mean") {
148  t = 0;
149  for (auto p : c.pulses) { t += p->StartTimeInNanoSec(); }
150  t = t / (double(c.pulses.size()));
151 
152  } else if (method == "Weighted_MaxAmpNorm") {
153  t = 0;
154  auto dmaxamp = *std::max_element(c.pulses.begin(), c.pulses.end(),
155  [](art::Ptr<brb::BeamlineDigit> p1, art::Ptr<brb::BeamlineDigit> p2){
156  return p1->AmplitudeInADC() > p2->AmplitudeInADC();
157  });
158  double maxamp = dmaxamp->AmplitudeInADC();
159  double w = 0;
160  for (auto p : c.pulses) {
161  t += p->StartTimeInNanoSec() * (p->AmplitudeInADC() / maxamp);
162  w += p->AmplitudeInADC() / maxamp;
163  }
164  t = t / w;
165 
166  } else if (method == "First") {
167  auto d = *std::min_element(c.pulses.begin(), c.pulses.end(),
168  [](art::Ptr<brb::BeamlineDigit> p1, art::Ptr<brb::BeamlineDigit> p2){
169  return p1->StartTimeInNanoSec() < p2->StartTimeInNanoSec();
170  });
171  t = d->StartTimeInNanoSec();
172 
173  } else {
174  throw art::Exception(art::errors::Configuration)
175  << "Unknown time of flight cluster time calculation method: " << method << "." << std::endl;
176  }
177 
178  return t;
179 }
180 
181 beamlinereco::ToFSingleCounterAnalysis::ToFSingleCounterAnalysis(const fhicl::ParameterSet &p) : EDAnalyzer(p)
182 {
183  reconfigure(p);
184 }
185 
186 void beamlinereco::ToFSingleCounterAnalysis::reconfigure(const fhicl::ParameterSet &p) {
187  fToFDigitLabel = p.get<art::InputTag>("ToFDigitLabel");
188  fSamplingInterval = p.get<double>("TimeSamplingInterval");
189  fHitClusterWindow = p.get<double>("HitClusterWindow");
190  fClusterTimeMethod = p.get<std::string>("ClusterTimeMethod", "First");
191 
192  fHistdtBound = p.get<double>("HistdtBound");
193  fHistdABound = p.get<double>("HistdABound");
194  fHistdIBound = p.get<double>("HistdIBound");
195  fHistdtNBins = p.get<unsigned int>("HistdtNBins");
196  fHistdANBins = p.get<unsigned int>("HistdANBins");
197  fHistdINBins = p.get<unsigned int>("HistdINBins");
198  fHistdAdtNBins = p.get<unsigned int>("HistdAdtNBins");
199  fHistdIdtNBins = p.get<unsigned int>("HistdIdtNBins");
200  fHistdAdINBins = p.get<unsigned int>("HistdIdtNBins");
201 }
202 
203 void beamlinereco::ToFSingleCounterAnalysis::beginJob() {
204  double endtime = 1024.*fSamplingInterval;
205 
206  //std::cout << "Channels in map:" << std::endl;
207  for (auto channel : fChannelMap->ChannelIDs()) {
208  unsigned int online_channel = onlineChannel(channel);
209  //std::cout << "\t" << online_channel << " " << channel << std::endl;
210  art::TFileDirectory chdir = tfs->mkdir(Form("Channel%d",online_channel));
211  this->hAmplitude_AllHits_channel[channel] = chdir.make<TH1D>(Form("hAmplitude_AllHits_channel%d",online_channel),
212  Form("Amplitude Channel %d (all hits)",online_channel),
213  200, 0, 4000);
214  this->hIntegral_AllHits_channel[channel] = chdir.make<TH1D>(Form("hIntegral_AllHits_channel%d",online_channel),
215  Form("Charge Integral Channel %d (all hits)",online_channel),
216  200, 0, 100000);
217  this->hStartTime_AllHits_channel[channel] = chdir.make<TH1D>(Form("hStartTime_AllHits_channel%d",online_channel),
218  Form("Start Time Channel %d (all hits)",online_channel),
219  200, 0, endtime);
220  this->hWidth_AllHits_channel[channel] = chdir.make<TH1D>(Form("hWidth_AllHits_channel%d",online_channel),
221  Form("Width Channel %d (all hits)",online_channel),
222  200, 0, 20);
223  this->hAmplitude_vs_Integral_AllHits_channel[channel] = chdir.make<TH2D>(Form("hAmplitude_vs_Integral_AllHits_channel%d",online_channel),
224  Form("Amplitude vs Integral Channel %d (all hits)",online_channel),
225  200, 0, 100000, 200, 0, 4000);
226  this->hAmplitude_vs_Width_AllHits_channel[channel] = chdir.make<TH2D>(Form("hAmplitude_vs_Width_AllHits_channel%d",online_channel),
227  Form("Amplitude vs Width Channel %d (all hits)",online_channel),
228  200, 0, 20, 200, 0, 4000);
229  }
230 
231  this->hUSClusterTime = tfs->make<TH1D>("hUSClusterTime", "US Cluster Start Time", 200, 0, endtime);
232  this->hDSClusterTime = tfs->make<TH1D>("hDSClusterTime", "DS Cluster Start Time", 200, 0, endtime);
233  this->hSiPMClusterTime = tfs->make<TH1D>("hSiPMClusterTime", "SiPM Cluster Start Time", 200, 0, endtime);
234 
235  this->hUSHitMultiplicity = tfs->make<TH1I>("hUSHitMultiplicity", "US Hits per Event", 15, 0, 15);
236  this->hDSHitMultiplicity = tfs->make<TH1I>("hDSHitMultiplicity", "DS Hits per Event", 15, 0, 15);
237  this->hSiPMHitMultiplicity = tfs->make<TH1I>("hSiPMHitMultiplicity", "SiPM Hits per Event", 15, 0, 15);
238 
239  this->hUSFirstHitDetector = tfs->make<TH1I>("hUSFirstHitDetector", "US First Hit Detector", 4, 1, 5);
240  this->hDSFirstHitDetector = tfs->make<TH1I>("hDSFirstHitDetector", "DS First Hit Detector", 4, 1, 5);
241  this->hSiPMFirstHitDetector = tfs->make<TH1I>("hSiPMFirstHitDetector", "SiPM First Hit Detector", 4, 1, 5);
242 
243  art::TFileDirectory usdir = tfs->mkdir("US deltas");
244  art::TFileDirectory dsdir = tfs->mkdir("DS deltas");
245  art::TFileDirectory sipmdir = tfs->mkdir("SiPM deltas");
246 
247  for (int i=0; i<4; i++) {
248  for (int j=i+1; j<4; j++) {
249  art::TFileDirectory usijdir = tfs->mkdir(Form("US delta %d-%d", i+1, j+1));
250  art::TFileDirectory dsijdir = tfs->mkdir(Form("DS delta %d-%d", i+1, j+1));
251  art::TFileDirectory sipmijdir = tfs->mkdir(Form("SiPM delta %d-%d", i+1, j+1));
252 
253  hUSdt[i][j] = usijdir.make<TH1D>(Form("hUSdt%d%d",i+1,j+1), Form("ToF dt between channels %d and %d",i+8,j+8),
254  fHistdtNBins, -fHistdtBound, fHistdtBound);
255  hUSdA[i][j] = usijdir.make<TH1D>(Form("hUSdA%d%d",i+1,j+1), Form("ToF dA between channels %d and %d",i+8,j+8),
256  fHistdANBins, -fHistdABound, fHistdABound);
257  hUSdI[i][j] = usijdir.make<TH1D>(Form("hUSdI%d%d",i+1,j+1), Form("ToF dI between channels %d and %d",i+8,j+8),
258  fHistdINBins, -fHistdIBound, fHistdIBound);
259  hUSdAdt[i][j] = usijdir.make<TH2D>(Form("hUSdAdt%d%d",i+1,j+1), Form("ToF dA vs dt channels %d and %d",i+8,j+8),
260  fHistdAdtNBins, -fHistdtBound, fHistdtBound,
261  fHistdAdtNBins, -fHistdABound, fHistdABound);
262  hUSdIdt[i][j] = usijdir.make<TH2D>(Form("hUSdIdt%d%d",i+1,j+1), Form("ToF dI vs dt channels %d and %d",i+8,j+8),
263  fHistdIdtNBins, -fHistdtBound, fHistdtBound,
264  fHistdIdtNBins, -fHistdIBound, fHistdIBound);
265  hUSdAdI[i][j] = usijdir.make<TH2D>(Form("hUSdAdI%d%d",i+1,j+1), Form("ToF dA vs dI channels %d and %d",i+8,j+8),
266  fHistdAdINBins, -fHistdIBound, fHistdIBound,
267  fHistdAdINBins, -fHistdABound, fHistdABound);
268 
269  hDSdt[i][j] = dsijdir.make<TH1D>(Form("hDSdt%d%d",i+1,j+1), Form("ToF dt between channels %d and %d",i+12,j+12),
270  fHistdtNBins, -fHistdtBound, fHistdtBound);
271  hDSdA[i][j] = dsijdir.make<TH1D>(Form("hDSdA%d%d",i+1,j+1), Form("ToF dA between channels %d and %d",i+12,j+12),
272  fHistdANBins, -fHistdABound, fHistdABound);
273  hDSdI[i][j] = dsijdir.make<TH1D>(Form("hDSdI%d%d",i+1,j+1), Form("ToF dI between channels %d and %d",i+12,j+12),
274  fHistdINBins, -fHistdIBound, fHistdIBound);
275  hDSdAdt[i][j] = dsijdir.make<TH2D>(Form("hDSdAdt%d%d",i+1,j+1), Form("ToF dA vs dt channels %d and %d",i+12,j+12),
276  fHistdAdtNBins, -fHistdtBound, fHistdtBound,
277  fHistdAdtNBins, -fHistdABound, fHistdABound);
278  hDSdIdt[i][j] = dsijdir.make<TH2D>(Form("hDSdIdt%d%d",i+1,j+1), Form("ToF dI vs dt channels %d and %d",i+12,j+12),
279  fHistdIdtNBins, -fHistdtBound, fHistdtBound,
280  fHistdIdtNBins, -fHistdIBound, fHistdIBound);
281  hDSdAdI[i][j] = dsijdir.make<TH2D>(Form("hDSdAdI%d%d",i+1,j+1), Form("ToF dA vs dI channels %d and %d",i+12,j+12),
282  fHistdAdINBins, -fHistdIBound, fHistdIBound,
283  fHistdAdINBins, -fHistdABound, fHistdABound);
284 
285  hSiPMdt[i][j] = sipmijdir.make<TH1D>(Form("hSiPMdt%d%d",i+1,j+1), Form("ToF dt between channels %d and %d",i+4,j+4),
286  fHistdtNBins, -fHistdtBound, fHistdtBound);
287  hSiPMdA[i][j] = sipmijdir.make<TH1D>(Form("hSiPMdA%d%d",i+1,j+1), Form("ToF dA between channels %d and %d",i+4,j+4),
288  fHistdANBins, -fHistdABound, fHistdABound);
289  hSiPMdI[i][j] = sipmijdir.make<TH1D>(Form("hSiPMdI%d%d",i+1,j+1), Form("ToF dI between channels %d and %d",i+4,j+4),
290  fHistdINBins, -fHistdIBound, fHistdIBound);
291  hSiPMdAdt[i][j] = sipmijdir.make<TH2D>(Form("hSiPMdAdt%d%d",i+1,j+1), Form("ToF dA vs dt channels %d and %d",i+4,j+4),
292  fHistdAdtNBins, -fHistdtBound, fHistdtBound,
293  fHistdAdtNBins, -fHistdABound, fHistdABound);
294  hSiPMdIdt[i][j] = sipmijdir.make<TH2D>(Form("hSiPMdIdt%d%d",i+1,j+1), Form("ToF dI vs dt channels %d and %d",i+4,j+4),
295  fHistdIdtNBins, -fHistdtBound, fHistdtBound,
296  fHistdIdtNBins, -fHistdIBound, fHistdIBound);
297  hSiPMdAdI[i][j] = sipmijdir.make<TH2D>(Form("hSiPMdAdI%d%d",i+1,j+1), Form("ToF dA vs dI channels %d and %d",i+4,j+4),
298  fHistdAdINBins, -fHistdIBound, fHistdIBound,
299  fHistdAdINBins, -fHistdABound, fHistdABound);
300  }
301  }
302 }
303 
304 void beamlinereco::ToFSingleCounterAnalysis::analyze(const art::Event &e)
305 {
306  // Load ToF digits
307  art::Handle<std::vector<brb::BeamlineDigit>> tofDigitHandle;
308  std::vector<art::Ptr<brb::BeamlineDigit>> tofDigits;
309  if (e.getByLabel(fToFDigitLabel, tofDigitHandle))
310  art::fill_ptr_vector(tofDigits, tofDigitHandle);
311 
312  // process digits
313 
314  unsigned int USHits = 0;
315  unsigned int DSHits = 0;
316  unsigned int SiPMHits = 0;
317 
318  //std::cout << "------------------------------ run " << e.run() << " subrun " << e.subRun() << " event " << e.event() << std::endl;
319  for (auto d : tofDigits) {
320  //unsigned int ch = onlineChannel(d->ChannelID());
321  //std::cout << "Digit " << d << " Online Channel " << ch << " Channel " << d->ChannelID();
322  //std::cout << " Associated hStartTime_AllHits_channel: " << hStartTime_AllHits_channel[d->ChannelID()]->GetName() << std::endl << std::endl;
323 
324  switch(d->ChannelID().Detector) {
325  case beamlinegeo::ToFCounter::US: ++USHits; break;
326  case beamlinegeo::ToFCounter::DS: ++DSHits; break;
327  case beamlinegeo::ToFCounter::DSSiPM: ++SiPMHits; break;
328  default: break;
329  }
330 
331  //std::cout << "width: " << d->WidthInNanoSec() << std::endl;
332  hAmplitude_AllHits_channel[d->ChannelID()]->Fill(d->AmplitudeInADC());
333  hIntegral_AllHits_channel[d->ChannelID()]->Fill(abs(d->AreaInADCNanoSec()));
334  hStartTime_AllHits_channel[d->ChannelID()]->Fill(d->StartTimeInNanoSec());
335  hWidth_AllHits_channel[d->ChannelID()]->Fill(d->WidthInNanoSec());
336  hAmplitude_vs_Integral_AllHits_channel[d->ChannelID()]->Fill(abs(d->AreaInADCNanoSec()), d->AmplitudeInADC());
337  hAmplitude_vs_Width_AllHits_channel[d->ChannelID()]->Fill(d->WidthInNanoSec(), d->AmplitudeInADC());
338  }
339 
340  // Process clusters
341 
342  // Find clusters
343  ToFClusterAlg clusterAlg;
344  clusterAlg.FindClusters(tofDigits, fHitClusterWindow);
345 
346  std::vector<ToFPulseCluster> usclusters = clusterAlg.GetClusters(beamlinegeo::ToFCounter::US);
347  std::vector<ToFPulseCluster> dsclusters = clusterAlg.GetClusters(beamlinegeo::ToFCounter::DS);
348  std::vector<ToFPulseCluster> sipmclusters = clusterAlg.GetClusters(beamlinegeo::ToFCounter::DSSiPM);
349 
350  for (auto c : usclusters) {
351  //std::cout << "----------------------------" << std::endl;
352  //for (auto a : c.pulses) std::cout << a->Channel() << " " << a->StartTimeInNanoSec() << std::endl;
353  unsigned int earliest = onlineChannel((*getEarliest(c.pulses))->ChannelID()) - 8+1;
354  //std::cout << "Earliest: " << earliest + 8-1 << " (" << earliest << ")" << std::endl;
355  hUSFirstHitDetector->Fill(Form("%d",earliest), earliest);
356  hUSClusterTime->Fill(calculateTime(c, fClusterTimeMethod));
357  double times[4] = { -1, -1, -1, -1 };
358  double amps[4] = { -1, -1, -1, -1 };
359  double areas[4] = { -1, -1, -1, -1 };
360  for (auto d : c.pulses) {
361  times[onlineChannel(d->ChannelID()) - 8] = d->StartTimeInNanoSec();
362  amps[onlineChannel(d->ChannelID()) - 8] = d->AmplitudeInADC();
363  areas[onlineChannel(d->ChannelID()) - 8] = d->AreaInADCNanoSec();
364  }
365  for (int i=0; i<4; i++) {
366  for (int j=i+1; j<4; j++) {
367  if ((times[i] >= 0) && (times[j] >= 0)) {
368  double dA = amps[j] - amps[i];
369  double dI = areas[j] - areas[i];
370  double dt = times[j] - times[i];
371  hUSdt[i][j]->Fill(dt);
372  hUSdA[i][j]->Fill(dA);
373  hUSdI[i][j]->Fill(dI);
374  hUSdAdt[i][j]->Fill(dt, dA);
375  hUSdIdt[i][j]->Fill(dt, dI);
376  hUSdAdI[i][j]->Fill(dI, dA);
377  }
378  }
379  }
380  }
381 
382  for (auto c : dsclusters) {
383  unsigned int earliest = onlineChannel((*getEarliest(c.pulses))->ChannelID()) - 12+1;
384  hDSFirstHitDetector->Fill(Form("%d",earliest), earliest);
385  hDSClusterTime->Fill(calculateTime(c, fClusterTimeMethod));
386  double times[4] = { -1, -1, -1, -1 };
387  double amps[4] = { -1, -1, -1, -1 };
388  double areas[4] = { -1, -1, -1, -1 };
389  for (auto d : c.pulses) {
390  times[onlineChannel(d->ChannelID()) - 12] = d->StartTimeInNanoSec();
391  amps[onlineChannel(d->ChannelID()) - 12] = d->AmplitudeInADC();
392  areas[onlineChannel(d->ChannelID()) - 12] = d->AreaInADCNanoSec();
393  }
394  for (int i=0; i<4; i++) {
395  for (int j=i+1; j<4; j++) {
396  if ((times[i] >= 0) && (times[j] >= 0)) {
397  double dA = amps[j] - amps[i];
398  double dI = areas[j] - areas[i];
399  double dt = times[j] - times[i];
400  hDSdt[i][j]->Fill(dt);
401  hDSdA[i][j]->Fill(dA);
402  hDSdI[i][j]->Fill(dI);
403  hDSdAdt[i][j]->Fill(dt, dA);
404  hDSdIdt[i][j]->Fill(dt, dI);
405  hDSdAdI[i][j]->Fill(dI, dA);
406  }
407  }
408  }
409  }
410 
411  for (auto c : sipmclusters) {
412  unsigned int earliest = onlineChannel((*getEarliest(c.pulses))->ChannelID()) - 4+1;
413  hSiPMFirstHitDetector->Fill(Form("%d",earliest), earliest);
414  hSiPMClusterTime->Fill(calculateTime(c, fClusterTimeMethod));
415  double times[4] = { -1, -1, -1, -1 };
416  double amps[4] = { -1, -1, -1, -1 };
417  double areas[4] = { -1, -1, -1, -1 };
418  for (auto d : c.pulses) {
419  times[onlineChannel(d->ChannelID()) - 4] = d->StartTimeInNanoSec();
420  amps[onlineChannel(d->ChannelID()) - 4] = d->AmplitudeInADC();
421  areas[onlineChannel(d->ChannelID()) - 4] = d->AreaInADCNanoSec();
422  }
423  for (int i=0; i<4; i++) {
424  for (int j=i+1; j<4; j++) {
425  if ((times[i] >= 0) && (times[j] >= 0)) {
426  double dA = amps[j] - amps[i];
427  double dI = areas[j] - areas[i];
428  double dt = times[j] - times[i];
429  hSiPMdt[i][j]->Fill(dt);
430  hSiPMdA[i][j]->Fill(dA);
431  hSiPMdI[i][j]->Fill(dI);
432  hSiPMdAdt[i][j]->Fill(dt, dA);
433  hSiPMdIdt[i][j]->Fill(dt, dI);
434  hSiPMdAdI[i][j]->Fill(dI, dA);
435  }
436  }
437  }
438  }
439 
440  hUSHitMultiplicity->Fill(USHits);
441  hDSHitMultiplicity->Fill(DSHits);
442  hSiPMHitMultiplicity->Fill(SiPMHits);
443 }
444 
445 DEFINE_ART_MODULE(beamlinereco::ToFSingleCounterAnalysis)
plot_validation_datamc.p2
p2
Definition: plot_validation_datamc.py:199
beamlinereco::ToFSingleCounterAnalysis::calculateTime
double calculateTime(ToFPulseCluster c, std::string method) const
Definition: ToFSingleCounterAnalysis_module.cc:140
beamlinereco::ToFClusterAlg::GetClusters
std::vector< ToFPulseCluster > GetClusters(beamlinegeo::ToFCounter counter) const
Definition: ToFClusterAlg.h:75
beamlinereco::ToFSingleCounterAnalysis::hAmplitude_vs_Integral_AllHits_channel
std::map< beamlinegeo::ChannelID, TH2D * > hAmplitude_vs_Integral_AllHits_channel
Definition: ToFSingleCounterAnalysis_module.cc:73
beamlinereco::ToFSingleCounterAnalysis::hWidth_AllHits_channel
std::map< beamlinegeo::ChannelID, TH1D * > hWidth_AllHits_channel
Definition: ToFSingleCounterAnalysis_module.cc:72
beamlinereco::ToFSingleCounterAnalysis::tfs
art::ServiceHandle< art::TFileService > tfs
Definition: ToFSingleCounterAnalysis_module.cc:112
getEarliest
static std::vector< art::Ptr< brb::BeamlineDigit > >::iterator getEarliest(std::vector< art::Ptr< brb::BeamlineDigit >> v)
Definition: ToFSingleCounterAnalysis_module.cc:131
p
const char * p
Definition: xmltok.h:285
brb::BeamlineDigit::AmplitudeInADC
int AmplitudeInADC() const
Definition: BeamlineDigit.cxx:56
art::TFileDirectory::mkdir
TFileDirectory mkdir(std::string const &dir, std::string const &descr="")
beamlineutil::BeamlineChannelMap::ChannelIDs
std::vector< ChannelID > ChannelIDs() const
All the ChannelIDs in this channel mapping.
Definition: BeamlineChannelMap_service.cc:168
brb::BeamlineDigit::StartTimeInNanoSec
double StartTimeInNanoSec() const
Definition: BeamlineDigit.cxx:36
abs
void abs(TH1 *hist)
Definition: wrong_sign_uncertainty.C:716
beamlinereco::ToFClusterAlg::FindClusters
void FindClusters(std::vector< art::Ptr< brb::BeamlineDigit >> digits, double clusterWindow)
Definition: ToFClusterAlg.h:43
DEFINE_ART_MODULE
DEFINE_ART_MODULE(TestTMapFile)
dt
Definition: Cand.cxx:23
beamlinereco::ToFSingleCounterAnalysis::hIntegral_AllHits_channel
std::map< beamlinegeo::ChannelID, TH1D * > hIntegral_AllHits_channel
Definition: ToFSingleCounterAnalysis_module.cc:70
onlineChannel
static unsigned int onlineChannel(ChannelID id)
Definition: ToFSingleCounterAnalysis_module.cc:117
BeamlineDigit.h
Encapsulation of reconstructed digitizer &#39;hits&#39;. Used for ToF PMTs and SiPMs, and Cherenkov and Muon ...
beamlinereco::ToFSingleCounterAnalysis::hAmplitude_vs_Width_AllHits_channel
std::map< beamlinegeo::ChannelID, TH2D * > hAmplitude_vs_Width_AllHits_channel
Definition: ToFSingleCounterAnalysis_module.cc:74
fhicl::ParameterSet::get
T get(std::string const &key) const
Definition: ParameterSet.h:231
plot_validation_datamc.c
c
Definition: plot_validation_datamc.py:198
d
Float_t d
Definition: plot.C:236
calib::j
const double j
Definition: BetheBloch.cxx:29
beamlinereco::ToFSingleCounterAnalysis::hStartTime_AllHits_channel
std::map< beamlinegeo::ChannelID, TH1D * > hStartTime_AllHits_channel
Definition: ToFSingleCounterAnalysis_module.cc:71
art::EDAnalyzer::EDAnalyzer
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
plot_timing_data.times
dictionary times
Definition: plot_timing_data.py:11
Event.h
art::Exception
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
art::TFileDirectory::make
T * make(ARGS...args) const
Definition: TFileDirectory.h:82
art::DataViewImpl::getByLabel
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
beamlinereco::ToFSingleCounterAnalysis::fChannelMap
art::ServiceHandle< beamlineutil::BeamlineChannelMap > fChannelMap
Definition: ToFSingleCounterAnalysis_module.cc:111
beamlinereco::ToFPulseCluster::pulses
std::vector< art::Ptr< brb::BeamlineDigit > > pulses
Definition: ToFClusterAlg.h:16
beamlinereco::ToFSingleCounterAnalysis::hAmplitude_AllHits_channel
std::map< beamlinegeo::ChannelID, TH1D * > hAmplitude_AllHits_channel
Definition: ToFSingleCounterAnalysis_module.cc:69
registry_explorer.v
v
Definition: registry_explorer.py:22
art::fill_ptr_vector
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:464
e
Float_t e
Definition: plot.C:35
w
Float_t w
Definition: plot.C:20
art::Ptr
Definition: fwd.h:28
plot_validation_datamc.p1
p1
Definition: plot_validation_datamc.py:199
BeamlineChannelMap.h
Channel mapping service which may be used to interpret the channels which are read out by the various...
string
enum BeamMode string
Generated on Mon Mar 1 2021 14:14:20 for NOvA by   doxygen 1.8.11