TriggerRateAna_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: TriggerRateAna
3 // Plugin Type: analyzer (art v2_12_01)
4 // File: TriggerRateAna_module.cc
5 //
6 // Generated at Mon Jun 10 10:50:11 2019 by Dung Phan using cetskelgen
7 // from cetlib version v3_06_00.
8 ////////////////////////////////////////////////////////////////////////
9 
20 #include "fhiclcpp/ParameterSet.h"
21 
22 // nova
27 #include "BeamlineRecoBase/ToF.h"
29 #include "RawData/RawBeamline.h"
31 
32 // ROOT
33 #include "TCanvas.h"
34 #include "TFile.h"
35 #include "TStyle.h"
36 #include "TROOT.h"
37 #include "TGraph.h"
38 #include "TH1D.h"
39 #include "TH2D.h"
40 #include "TLatex.h"
41 #include "TLegend.h"
42 #include "TMath.h"
43 #include "TTree.h"
44 
45 // stl
46 #include <iostream>
47 #include <vector>
48 
49 namespace novatb {
50 class TriggerRateAna;
51 }
52 
54  public:
55  explicit TriggerRateAna(fhicl::ParameterSet const& p);
56  void reconfigure(const fhicl::ParameterSet& p);
57  void analyze(art::Event const& e) override;
58  void beginJob() override;
59  void beginRun(art::Run const& r) override;
60  void endJob() override;
61 
62  void decimal2binary(unsigned int decimalNumber, unsigned int* binaryNumber);
63  std::vector<double> findCluster(unsigned int tofIndex);
64  // double calculate_tof_inns(double distance_incm, double mass_inMeV, double
65  // momentum_inMeV);
66 
67  private:
68  unsigned int _RunNumber;
69 
74 
75  std::vector<double> tof_p1_p4;
76  std::vector<double> momentum;
77  std::vector<bool> isElectron;
78 
79  unsigned int countTriggerPad4in4;
80  unsigned int countTriggerWC3in4;
81  unsigned int countTotalTriggers;
83  unsigned int countElectronEvents;
84  unsigned int totalEvents;
85 
87  unsigned int _SelectedRun;
88  double _MomentumFlip;
89  unsigned int _MagnetCurrent;
90 
91  std::vector<double> startTimeOfHitInChannel_4;
92  std::vector<double> startTimeOfHitInChannel_5;
93  std::vector<double> startTimeOfHitInChannel_6;
94  std::vector<double> startTimeOfHitInChannel_7;
95  std::vector<double> startTimeOfHitInChannel_8;
96  std::vector<double> startTimeOfHitInChannel_9;
97  std::vector<double> startTimeOfHitInChannel_10;
98  std::vector<double> startTimeOfHitInChannel_11;
99  std::vector<double> startTimeOfHitInChannel_12;
100  std::vector<double> startTimeOfHitInChannel_13;
101  std::vector<double> startTimeOfHitInChannel_14;
102  std::vector<double> startTimeOfHitInChannel_15;
103 
109 
110  TH1D* hTriggerBit;
111  std::string TriggerBitName[16] = {"WC1IGNOREME", "WC2IGNOREME", "WC3IGNOREME", "WCCOINC",
112  "PAD1" , "PAD2" , "PAD3" , "PAD4" ,
113  "READOUTGATE", "SWYDGATE" , "MWPCOFF" , "DIGBUSY",
114  "WC1" , "WC2" , "WC3" , "WC4"};
115 
116  TGraph* pidPlot_noE;
117  TGraph* pidPlot_E;
118 
119  // services
121 };
122 
124  : EDAnalyzer(p) {
125  reconfigure(p);
126 }
127 
129  _ToFInputTag = "tofcfdhitfinder";
130  _WCInputTag = "wctrackreco";
131  _CherenkovInputTag = "cherenkovreco";
132  _TriggerInputTag = "daq:Trigger";
133 
134  _CherenkovDelayAcceptance = p.get<double>("CherenkovDelayAcceptance");
135  _SelectedRun = p.get<unsigned int>("SelectedRun");
136 }
137 
139  hTrigPattern = new TH1D("hTrigPattern", ";Trigger Pattern;", 6, 0., 6.);
140  hTOF_MNC_0500A = new TH1D("hTOF_MNC_0500A",
141  "Magnet Current at 500A;TOF (ns);", 100, 0, 100.);
142  hTOF_MNC_1000A = new TH1D("hTOF_MNC_1000A",
143  "Magnet Current at 1000A;TOF (ns);", 100, 0, 100.);
144  hTOF_MNC_1500A = new TH1D("hTOF_MNC_1500A",
145  "Magnet Current at 1500A;TOF (ns);", 100, 0, 100.);
146  hTOF_MNC_2000A = new TH1D("hTOF_MNC_2000A",
147  "Magnet Current at 2000A;TOF (ns);", 100, 0, 100.);
148 
149  hTriggerBit = new TH1D("hTriggerBit",
150  ";;", 16, 0, 16);
151  hTriggerBit->GetXaxis()->CanExtend();
152  for (unsigned int i = 0; i < 16; i++) {
153  hTriggerBit->Fill(TriggerBitName[i].c_str(), 1);
154  }
155 
156 
158  countTriggerWC3in4 = 0;
159  countTotalTriggers = 0;
162  totalEvents = 0;
163  tof_p1_p4.clear();
164  momentum.clear();
165  isElectron.clear();
166 }
167 
169  _RunNumber = r.run();
170  if (_RunNumber == 2854 || _RunNumber == 2855 || _RunNumber == 2856 ||
171  _RunNumber == 2857 || _RunNumber == 2861 || _RunNumber == 2862) {
172  _MomentumFlip = 1920.;
173  _MagnetCurrent = 2000;
174  } else if (_RunNumber == 2853 || _RunNumber == 2858) {
175  _MomentumFlip = 1480.;
176  _MagnetCurrent = 1500;
177  } else if (_RunNumber == 2851 || _RunNumber == 2852 || _RunNumber == 2860) {
178  _MomentumFlip = 980.;
179  _MagnetCurrent = 1000;
180  } else {
181  _MomentumFlip = 450.;
182  _MagnetCurrent = 500;
183  }
184 }
185 
187  if (true) {
188  // if (_RunNumber == _SelectedRun) {
190  std::vector<art::Ptr<brb::BeamlineDigit> > ToFs;
191  if (e.getByLabel(_ToFInputTag, ToFHandle))
192  art::fill_ptr_vector(ToFs, ToFHandle);
193 
195  std::vector<art::Ptr<brb::WCTrack> > WCs;
196  if (e.getByLabel(_WCInputTag, WCHandle))
197  art::fill_ptr_vector(WCs, WCHandle);
198 
200  std::vector<art::Ptr<brb::BeamlineDigit> > Cherenkovs;
201  if (e.getByLabel(_CherenkovInputTag, CherenkovHandle))
202  art::fill_ptr_vector(Cherenkovs, CherenkovHandle);
203 
205  std::vector<art::Ptr<rawdata::RawBeamlineTrigger> > Triggers;
206  if (e.getByLabel(_TriggerInputTag, TriggerHandle))
207  art::fill_ptr_vector(Triggers, TriggerHandle);
208 
221  for (auto itr = ToFs.begin(); itr != ToFs.end(); itr++) {
222  unsigned int channel = fChannelMap->OnlineDigitChannel((*itr)->ChannelID());
223  switch (channel) {
224  case (4): {
225  startTimeOfHitInChannel_4.push_back((*itr)->StartTimeInNanoSec());
226  break;
227  }
228  case (5): {
229  startTimeOfHitInChannel_5.push_back((*itr)->StartTimeInNanoSec());
230  break;
231  }
232  case (6): {
233  startTimeOfHitInChannel_6.push_back((*itr)->StartTimeInNanoSec());
234  break;
235  }
236  case (7): {
237  startTimeOfHitInChannel_7.push_back((*itr)->StartTimeInNanoSec());
238  break;
239  }
240  case (8): {
241  startTimeOfHitInChannel_8.push_back((*itr)->StartTimeInNanoSec());
242  break;
243  }
244  case (9): {
245  startTimeOfHitInChannel_9.push_back((*itr)->StartTimeInNanoSec());
246  break;
247  }
248  case (10): {
249  startTimeOfHitInChannel_10.push_back((*itr)->StartTimeInNanoSec());
250  break;
251  }
252  case (11): {
253  startTimeOfHitInChannel_11.push_back((*itr)->StartTimeInNanoSec());
254  break;
255  }
256  case (12): {
257  startTimeOfHitInChannel_12.push_back((*itr)->StartTimeInNanoSec());
258  break;
259  }
260  case (13): {
261  startTimeOfHitInChannel_13.push_back((*itr)->StartTimeInNanoSec());
262  break;
263  }
264  case (14): {
265  startTimeOfHitInChannel_14.push_back((*itr)->StartTimeInNanoSec());
266  break;
267  }
268  case (15): {
269  startTimeOfHitInChannel_15.push_back((*itr)->StartTimeInNanoSec());
270  break;
271  }
272  default: {
273  break;
274  }
275  }
276  } // end ToFs elem loop
277 
278  std::vector<double> hitStartTOFUS = findCluster(0);
279  std::vector<double> hitStartTOFDS = findCluster(1);
280  bool passSelection = true;
281  unsigned int nHitsTOFUS = hitStartTOFUS.size();
282  unsigned int nHitsTOFDS = hitStartTOFDS.size();
283  if (!(nHitsTOFUS == 1 && nHitsTOFDS == 1)) passSelection = false;
284  if (nHitsTOFUS == 1 && nHitsTOFDS == 1) {
285  if (hitStartTOFUS.at(0) - 10.8 >= hitStartTOFDS.at(0))
286  passSelection = false;
287  if (hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8 <= 39)
288  passSelection = false;
289 
290  switch (_MagnetCurrent) {
291  case (500): {
292  if (hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8 >= 100)
293  passSelection = false;
294  break;
295  }
296  case (1000): {
297  if (hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8 >= 70)
298  passSelection = false;
299  break;
300  }
301  case (1500): {
302  if (hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8 >= 58)
303  passSelection = false;
304  break;
305  }
306  case (2000): {
307  if (hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8 >= 55)
308  passSelection = false;
309  break;
310  }
311  default: {
312  break;
313  }
314  }
315  }
316 
317  totalEvents++;
318  for (auto trigItr = Triggers.begin(); trigItr != Triggers.end();
319  trigItr++) {
320  unsigned int pattern = (*trigItr)->Pattern();
321 
322  unsigned int binaryNumber[32];
323  decimal2binary(pattern, &binaryNumber[0]);
324  // std::cout << pattern << ", ";
325  int sumPattern = 0;
326  for (unsigned int i = 0; i < 16; i++) {
327  if (binaryNumber[i] == 1) hTriggerBit->Fill(TriggerBitName[i].c_str(), 1);
328  sumPattern += binaryNumber[i];
329  // std::cout << binaryNumber[i] << ", ";
330  }
331  // std::cout << std::endl;
332  if (sumPattern == 2) countTotalTriggers++;
333  }
334 
335  // For event that pass selection, get momentum
336 
337  // if (passSelection) {
338  // countPassSelectionEvents++;
339 
340  // if (WCs.size() == 1) {
341  // for (auto itr = WCs.begin(); itr != WCs.end(); itr++) {
342  // double mom = (*itr)->Momentum();
343  // // double mom_diff = TMath::Abs(mom - _MomentumFlip);
344  // // if (_MagnetCurrent != 500) {
345  // // mom = mom < _MomentumFlip ? _MomentumFlip + mom_diff :
346  // _MomentumFlip - mom_diff;
347  // // }
348  // momentum.push_back(mom);
349  // double tof = hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8;
350  // tof_p1_p4.push_back(tof);
351 
352  // switch (_MagnetCurrent) {
353  // case (500): {
354  // hTOF_MNC_0500A->Fill(tof - 10.8);
355  // break;
356  // }
357  // case (1000): {
358  // hTOF_MNC_1000A->Fill(tof - 10.8);
359  // break;
360  // }
361  // case (1500): {
362  // hTOF_MNC_1500A->Fill(tof - 10.8);
363  // break;
364  // }
365  // case (2000): {
366  // hTOF_MNC_2000A->Fill(tof - 10.8);
367  // break;
368  // }
369  // default: {break;}
370  // }
371 
372  // bool foundElectron = false;
373  // for (auto CherenkovItr = Cherenkovs.begin(); CherenkovItr !=
374  // Cherenkovs.end(); CherenkovItr++) {
375  // double chevStart = (*CherenkovItr)->StartTimeInNanoSec();
376  // if (hitStartTOFUS.at(0) - 10.8 > chevStart) continue;
377  // if (hitStartTOFDS.at(0) > chevStart) continue;
378  // if (chevStart - hitStartTOFDS.at(0) > _CherenkovDelayAcceptance)
379  // continue; foundElectron = true; countElectronEvents++;
380  // }
381  // isElectron.push_back(foundElectron);
382  // } // end WCs elem loop
383  // } // make sure only one and only one track is associated with this
384  // trigger
385 
386  // for (auto trigItr = Triggers.begin(); trigItr != Triggers.end();
387  // trigItr++) {
388  // countTotalTriggers++;
389  // unsigned int pattern = (*trigItr)->Pattern();
390 
391  // unsigned int binaryNumber[32];
392  // decimal2binary(pattern, &binaryNumber[0]);
393  // if (!binaryNumber[0] && binaryNumber[1] && binaryNumber[2] &&
394  // binaryNumber[3]) {
395  // hTrigPattern->Fill(0.5);
396  // } else if (binaryNumber[0] && !binaryNumber[1] && binaryNumber[2] &&
397  // binaryNumber[3]) {
398  // hTrigPattern->Fill(1.5);
399  // } else if (binaryNumber[0] && binaryNumber[1] && !binaryNumber[2] &&
400  // binaryNumber[3]) {
401  // hTrigPattern->Fill(2.5);
402  // } else if (binaryNumber[0] && binaryNumber[1] && binaryNumber[2] &&
403  // !binaryNumber[3]) {
404  // hTrigPattern->Fill(3.5);
405  // } else if (binaryNumber[0] && binaryNumber[1] && binaryNumber[2] &&
406  // binaryNumber[3]) {
407  // hTrigPattern->Fill(4.5);
408  // } else if (binaryNumber[4] && binaryNumber[5] && binaryNumber[6] &&
409  // binaryNumber[7]) {
410  // hTrigPattern->Fill(5.5);
411  // }
412  // }
413  // }
414  }
415 }
416 
417 std::vector<double> novatb::TriggerRateAna::findCluster(unsigned int tofIndex) {
418  std::vector<double> hitStartTOF;
419 
420  if (tofIndex == 0) {
421  for (unsigned int iNS = 0; iNS < 400; iNS++) {
422  double coincWindowStart = (double)iNS;
423  double coincWindowEnd = (double)iNS + 5.;
424 
425  bool coincidence_channel_8 = false;
426  bool coincidence_channel_9 = false;
427  bool coincidence_channel_10 = false;
428  bool coincidence_channel_11 = false;
429  double coincidence_time_channel_8 = 0;
430  double coincidence_time_channel_9 = 0;
431  double coincidence_time_channel_10 = 0;
432  double coincidence_time_channel_11 = 0;
433 
434  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_8.size();
435  iElem++) {
436  if (coincWindowStart < startTimeOfHitInChannel_8.at(iElem) &&
437  startTimeOfHitInChannel_8.at(iElem) < coincWindowEnd) {
438  coincidence_channel_8 = true;
439  coincidence_time_channel_8 = startTimeOfHitInChannel_8.at(iElem);
440  }
441  }
442 
443  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_9.size();
444  iElem++) {
445  if (coincWindowStart < startTimeOfHitInChannel_9.at(iElem) &&
446  startTimeOfHitInChannel_9.at(iElem) < coincWindowEnd) {
447  coincidence_channel_9 = true;
448  coincidence_time_channel_9 = startTimeOfHitInChannel_9.at(iElem);
449  }
450  }
451 
452  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_10.size();
453  iElem++) {
454  if (coincWindowStart < startTimeOfHitInChannel_10.at(iElem) &&
455  startTimeOfHitInChannel_10.at(iElem) < coincWindowEnd) {
456  coincidence_channel_10 = true;
457  coincidence_time_channel_10 = startTimeOfHitInChannel_10.at(iElem);
458  }
459  }
460 
461  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_11.size();
462  iElem++) {
463  if (coincWindowStart < startTimeOfHitInChannel_11.at(iElem) &&
464  startTimeOfHitInChannel_11.at(iElem) < coincWindowEnd) {
465  coincidence_channel_11 = true;
466  coincidence_time_channel_11 = startTimeOfHitInChannel_11.at(iElem);
467  }
468  }
469 
470  // check coinciden 3 out of 4 or 4 out of 4
471  bool tofus_coinc_8_9_10 =
472  (coincidence_channel_8 && coincidence_channel_9 &&
473  coincidence_channel_10);
474  bool tofus_coinc_8_9_11 =
475  (coincidence_channel_8 && coincidence_channel_9 &&
476  coincidence_channel_11);
477  bool tofus_coinc_8_10_11 =
478  (coincidence_channel_8 && coincidence_channel_10 &&
479  coincidence_channel_11);
480  bool tofus_coinc_9_10_11 =
481  (coincidence_channel_9 && coincidence_channel_10 &&
482  coincidence_channel_11);
483  bool tofus_coinc = (tofus_coinc_8_9_10 || tofus_coinc_8_9_11 ||
484  tofus_coinc_8_10_11 || tofus_coinc_9_10_11);
485  if (tofus_coinc) {
486  double sum = coincidence_time_channel_8 + coincidence_time_channel_9 +
487  coincidence_time_channel_10 + coincidence_time_channel_11;
488  if (tofus_coinc_8_9_10 && tofus_coinc_8_9_11 && tofus_coinc_8_10_11 &&
489  tofus_coinc_9_10_11) {
490  hitStartTOF.push_back(sum / 4);
491  } else {
492  if (tofus_coinc_8_9_10) hitStartTOF.push_back(sum / 3);
493  if (tofus_coinc_8_9_11) hitStartTOF.push_back(sum / 3);
494  if (tofus_coinc_8_10_11) hitStartTOF.push_back(sum / 3);
495  if (tofus_coinc_9_10_11) hitStartTOF.push_back(sum / 3);
496  }
497  iNS = (unsigned int)coincWindowEnd;
498  } else
499  continue;
500  } // loop over 400ns step of 1ns
501  } // apply for US TOF
502 
503  else if (tofIndex == 1) {
504  for (unsigned int iNS = 0; iNS < 400; iNS++) {
505  double coincWindowStart = (double)iNS;
506  double coincWindowEnd = (double)iNS + 5.;
507 
508  bool coincidence_channel_12 = false;
509  bool coincidence_channel_13 = false;
510  bool coincidence_channel_14 = false;
511  bool coincidence_channel_15 = false;
512  double coincidence_time_channel_12 = 0;
513  double coincidence_time_channel_13 = 0;
514  double coincidence_time_channel_14 = 0;
515  double coincidence_time_channel_15 = 0;
516 
517  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_12.size();
518  iElem++) {
519  if (coincWindowStart < startTimeOfHitInChannel_12.at(iElem) &&
520  startTimeOfHitInChannel_12.at(iElem) < coincWindowEnd) {
521  coincidence_channel_12 = true;
522  coincidence_time_channel_12 = startTimeOfHitInChannel_12.at(iElem);
523  }
524  }
525 
526  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_13.size();
527  iElem++) {
528  if (coincWindowStart < startTimeOfHitInChannel_13.at(iElem) &&
529  startTimeOfHitInChannel_13.at(iElem) < coincWindowEnd) {
530  coincidence_channel_13 = true;
531  coincidence_time_channel_13 = startTimeOfHitInChannel_13.at(iElem);
532  }
533  }
534 
535  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_14.size();
536  iElem++) {
537  if (coincWindowStart < startTimeOfHitInChannel_14.at(iElem) &&
538  startTimeOfHitInChannel_14.at(iElem) < coincWindowEnd) {
539  coincidence_channel_14 = true;
540  coincidence_time_channel_14 = startTimeOfHitInChannel_14.at(iElem);
541  }
542  }
543 
544  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_15.size();
545  iElem++) {
546  if (coincWindowStart < startTimeOfHitInChannel_15.at(iElem) &&
547  startTimeOfHitInChannel_15.at(iElem) < coincWindowEnd) {
548  coincidence_channel_15 = true;
549  coincidence_time_channel_15 = startTimeOfHitInChannel_15.at(iElem);
550  }
551  }
552 
553  // check coinciden 3 out of 4 or 4 out of 4
554  bool tofus_coinc_12_13_14 =
555  (coincidence_channel_12 && coincidence_channel_13 &&
556  coincidence_channel_14);
557  bool tofus_coinc_12_13_15 =
558  (coincidence_channel_12 && coincidence_channel_13 &&
559  coincidence_channel_15);
560  bool tofus_coinc_12_14_15 =
561  (coincidence_channel_12 && coincidence_channel_14 &&
562  coincidence_channel_15);
563  bool tofus_coinc_13_14_15 =
564  (coincidence_channel_13 && coincidence_channel_14 &&
565  coincidence_channel_15);
566  bool tofus_coinc = (tofus_coinc_12_13_14 || tofus_coinc_12_13_15 ||
567  tofus_coinc_12_14_15 || tofus_coinc_13_14_15);
568  if (tofus_coinc) {
569  double sum = coincidence_time_channel_12 + coincidence_time_channel_13 +
570  coincidence_time_channel_14 + coincidence_time_channel_15;
571  if (tofus_coinc_12_13_14 && tofus_coinc_12_13_15 &&
572  tofus_coinc_12_14_15 && tofus_coinc_13_14_15) {
573  hitStartTOF.push_back(sum / 4);
574  } else {
575  if (tofus_coinc_12_13_14) hitStartTOF.push_back(sum / 3);
576  if (tofus_coinc_12_13_15) hitStartTOF.push_back(sum / 3);
577  if (tofus_coinc_12_14_15) hitStartTOF.push_back(sum / 3);
578  if (tofus_coinc_13_14_15) hitStartTOF.push_back(sum / 3);
579  }
580  iNS = (unsigned int)coincWindowEnd;
581  } else
582  continue;
583  } // loop over 400ns step of 1ns
584  } // apply for US TOF
585 
586  return hitStartTOF;
587 }
588 
590  std::cout << "Total Triggers: " << totalEvents << std::endl;
591  std::cout << "Weird Triggers: " << countTotalTriggers << std::endl;
592 
593  std::cout << "countPassSelectionEvents: " << countPassSelectionEvents
594  << std::endl;
595  std::cout << "countElectronEvents: " << countElectronEvents << std::endl;
596 
597  gStyle->SetOptStat(0);
598 
599  TCanvas* c = new TCanvas("c", "c", 2400, 1200);
600  hTriggerBit->Draw();
601  c->SaveAs("TriggerBit.pdf");
602 
603  // TCanvas* tof_c = new TCanvas("tof_c", "tof_c", 1200, 1200);
604  // tof_c->Divide(2,2);
605  // tof_c->cd(1);
606  // hTOF_MNC_0500A->Draw();
607  // tof_c->cd(2);
608  // hTOF_MNC_1000A->Draw();
609  // tof_c->cd(3);
610  // hTOF_MNC_1500A->Draw();
611  // tof_c->cd(4);
612  // hTOF_MNC_2000A->Draw();
613  // tof_c->SaveAs("TOF.pdf");
614 
615  // const unsigned int nTotalPoints = momentum.size();
616 
617  // unsigned int countElectrons = 0;
618  // for (unsigned int i = 0; i < isElectron.size(); i++) {
619  // if (isElectron.at(i)) countElectrons++;
620  // }
621  // const unsigned int nElectronPoints = countElectrons;
622  // const unsigned int nNotElectronPoints = nTotalPoints - countElectrons;
623 
624  // double* tof_array_noE = new double[nNotElectronPoints];
625  // double* momentum_array_noE = new double[nNotElectronPoints];
626  // double* tof_array_E = new double[nElectronPoints];
627  // double* momentum_array_E = new double[nElectronPoints];
628  // unsigned int index_E = 0;
629  // unsigned int index_noE = 0;
630  // for (unsigned int i = 0; i < nTotalPoints; i++) {
631  // if (isElectron.at(i)) {
632  // *(tof_array_E + index_E) = tof_p1_p4.at(i);
633  // *(momentum_array_E + index_E) = momentum.at(i);
634  // index_E++;
635  // } else {
636  // *(tof_array_noE + index_noE) = tof_p1_p4.at(i);
637  // *(momentum_array_noE + index_noE) = momentum.at(i);
638  // index_noE++;
639  // }
640  // }
641 
642  // // double xmin = 100;
643  // // double xmax = 3000;
644  // // double ymin = 40;
645  // // double ymax = 100;
646 
647  // // double* momentum_for_particle_plots = new double[1000];
648  // // double* tof_for_deuteron_plots = new double[1000];
649  // // double* tof_for_proton_plots = new double[1000];
650  // // double* tof_for_kaon_plots = new double[1000];
651  // // double* tof_for_muon_plots = new double[1000];
652  // // double* tof_for_pion_plots = new double[1000];
653  // // double* tof_for_electron_plots = new double[1000];
654 
655  // // double momentum_step = (xmax - xmin) / 1000;
656  // // double distance_p1_p4_incm = 1316.22;
657  // // double mass_deuteron_inMeV = 1875.6;
658  // // double mass_proton_inMeV = 938.272;
659  // // double mass_kaon_inMeV = 493.677;
660  // // double mass_pion_inMeV = 139.570;
661  // // double mass_muon_inMeV = 105.658;
662  // // double mass_electron_inMeV = 0.511;
663  // // for (unsigned int i = 0; i < 1000; i++) {
664  // // *(momentum_for_particle_plots + i) = xmin + (double)i * momentum_step;
665  // // *(tof_for_deuteron_plots + i) =
666  // // calculate_tof_inns(distance_p1_p4_incm, mass_deuteron_inMeV,
667  // // *(momentum_for_particle_plots + i));
668  // // *(tof_for_proton_plots + i) =
669  // // calculate_tof_inns(distance_p1_p4_incm, mass_proton_inMeV,
670  // // *(momentum_for_particle_plots + i));
671  // // *(tof_for_kaon_plots + i) =
672  // // calculate_tof_inns(distance_p1_p4_incm, mass_kaon_inMeV,
673  // // *(momentum_for_particle_plots + i));
674  // // *(tof_for_pion_plots + i) =
675  // // calculate_tof_inns(distance_p1_p4_incm, mass_pion_inMeV,
676  // // *(momentum_for_particle_plots + i));
677  // // *(tof_for_muon_plots + i) =
678  // // calculate_tof_inns(distance_p1_p4_incm, mass_muon_inMeV,
679  // // *(momentum_for_particle_plots + i));
680  // // *(tof_for_electron_plots + i) =
681  // // calculate_tof_inns(distance_p1_p4_incm, mass_electron_inMeV,
682  // // *(momentum_for_particle_plots + i));
683  // // }
684 
685  // // TCanvas* canvas = new TCanvas("canvas", "canvas", 2400, 2400);
686  // pidPlot_noE =
687  // new TGraph(nNotElectronPoints, momentum_array_noE, tof_array_noE);
688  // pidPlot_noE->SetTitle("");
689  // pidPlot_noE->SetMarkerColor(kBlack);
690  // pidPlot_noE->SetMarkerStyle(20);
691  // pidPlot_noE->SetMarkerSize(1);
692  // pidPlot_noE->Draw("AP");
693 
694  // // TGraph* deuteronPlot = new TGraph(1000, momentum_for_particle_plots,
695  // // tof_for_deuteron_plots); TGraph* protonPlot = new TGraph(1000,
696  // // momentum_for_particle_plots, tof_for_proton_plots); TGraph* kaonPlot =
697  // // new TGraph(1000, momentum_for_particle_plots, tof_for_kaon_plots); TGraph*
698  // // pionPlot = new TGraph(1000, momentum_for_particle_plots,
699  // // tof_for_pion_plots); TGraph* muonPlot = new TGraph(1000,
700  // // momentum_for_particle_plots, tof_for_muon_plots); TGraph* electronPlot =
701  // // new TGraph(1000, momentum_for_particle_plots, tof_for_electron_plots);
702 
703  // if (nElectronPoints > 0) {
704  // pidPlot_E = new TGraph(nElectronPoints, momentum_array_E, tof_array_E);
705  // pidPlot_E->SetTitle("");
706  // pidPlot_E->SetMarkerColor(kRed);
707  // pidPlot_E->SetMarkerStyle(20);
708  // pidPlot_E->SetMarkerSize(1);
709  // pidPlot_E->Draw("P SAME");
710  // }
711 
712  // TFile* graphFile = new TFile("graphFile_10pB.root", "RECREATE");
713  // graphFile->cd();
714  // pidPlot_E->SetName("pidPlot_E");
715  // pidPlot_E->Write();
716  // pidPlot_noE->SetName("pidPlot_noE");
717  // pidPlot_noE->Write();
718  // graphFile->Write();
719  // graphFile->Close();
720 
721  // // deuteronPlot->SetLineColor(kMagenta);
722  // // protonPlot->SetLineColor(kRed);
723  // // kaonPlot->SetLineColor(kGreen - 2);
724  // // pionPlot->SetLineColor(kBlue + 1);
725  // // muonPlot->SetLineColor(kOrange - 2);
726  // // electronPlot->SetLineColor(kViolet - 1);
727 
728  // // deuteronPlot->Draw("L SAME");
729  // // protonPlot->Draw("L SAME");
730  // // kaonPlot->Draw("L SAME");
731  // // pionPlot->Draw("L SAME");
732  // // muonPlot->Draw("L SAME");
733  // // electronPlot->Draw("L SAME");
734 
735  // // TLegend* legend = new TLegend(0.5, 0.5, 0.8, 0.87);
736  // // legend->SetBorderSize(0);
737  // // legend->AddEntry(pidPlot_noE, "data", "p");
738  // // if (nElectronPoints > 0) {
739  // // legend->AddEntry(pidPlot_E, "data e", "p");
740  // // }
741  // // legend->AddEntry(deuteronPlot, "D", "l");
742  // // legend->AddEntry(protonPlot, "p", "l");
743  // // legend->AddEntry(kaonPlot, "K", "l");
744  // // legend->AddEntry(pionPlot, "#pi", "l");
745  // // legend->AddEntry(muonPlot, "#mu", "l");
746  // // legend->AddEntry(electronPlot, "e", "l");
747  // // legend->Draw();
748 
749  // // pidPlot_noE->GetXaxis()->SetTitle("Momentum (MeV/c)");
750  // // pidPlot_noE->GetXaxis()->CenterTitle();
751  // // pidPlot_noE->GetXaxis()->SetLimits(xmin, xmax);
752  // // pidPlot_noE->GetYaxis()->SetTitle("TOF (ns)");
753  // // pidPlot_noE->GetYaxis()->CenterTitle();
754  // // pidPlot_noE->GetYaxis()->SetRangeUser(ymin, ymax);
755 
756  // // auto prelim = TLatex(.9, .95, "NO#nuA Preliminary");
757  // // prelim.SetTextColor(4);
758  // // prelim.SetNDC();
759  // // prelim.SetTextSize(1.5 / 30.);
760  // // prelim.SetTextAlign(32);
761  // // prelim.Draw();
762 
763  // // canvas->SaveAs("pid.pdf");
764 }
765 
766 void novatb::TriggerRateAna::decimal2binary(unsigned int decimalNumber,
767  unsigned int* binaryNumber) {
768  for (int index = 31; index >= 0; index--) {
769  unsigned int reference = TMath::Power(2, (int)index);
770  if (decimalNumber >= reference) {
771  *(binaryNumber + index) = 1;
772  decimalNumber = decimalNumber - reference;
773  } else {
774  *(binaryNumber + index) = 0;
775  }
776  }
777 }
778 
779 // double novatb::TriggerRateAna::calculate_tof_inns(double distance_incm,
780 // double mass_inMeV, double momentum_inMeV) {
781 // double gamma = TMath::Sqrt((TMath::Power(momentum_inMeV, 2) +
782 // TMath::Power(mass_inMeV, 2)) / TMath::Power(mass_inMeV, 2)); double beta =
783 // TMath::Sqrt(1 - TMath::Power(1 / gamma, 2)); double velocity_mps =
784 // 299792458 * beta; double tof_s = (distance_incm / 100) / velocity_mps;
785 // return tof_s * 1E9;
786 // }
787 
void beginRun(art::Run const &r) override
std::vector< double > startTimeOfHitInChannel_6
std::vector< double > startTimeOfHitInChannel_9
const char * p
Definition: xmltok.h:285
TriggerRateAna(fhicl::ParameterSet const &p)
std::vector< bool > isElectron
DEFINE_ART_MODULE(TestTMapFile)
RunNumber_t run() const
Definition: Run.h:47
unsigned int OnlineDigitChannel(ChannelID channel) const
Online digiziter channel number for this offline ChannelID.
std::vector< double > startTimeOfHitInChannel_14
Definition: Run.h:31
Encapsulation of reconstructed digitizer &#39;hits&#39;. Used for ToF PMTs and SiPMs, and Cherenkov and Muon ...
void analyze(art::Event const &e) override
T get(std::string const &key) const
Definition: ParameterSet.h:231
Encapsulation of reconstructed Time-of-Flight (ToF) information. Part of beamline reconstruction for ...
std::vector< double > startTimeOfHitInChannel_4
Encapsulation of reconstructed Wire Chamber track. Part of beamline reconstruction for NOvA test beam...
std::vector< double > startTimeOfHitInChannel_5
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
std::vector< double > findCluster(unsigned int tofIndex)
OStream cout
Definition: OStream.cxx:6
std::vector< double > startTimeOfHitInChannel_15
std::vector< double > startTimeOfHitInChannel_11
void reconfigure(const fhicl::ParameterSet &p)
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
std::vector< double > startTimeOfHitInChannel_13
std::vector< double > startTimeOfHitInChannel_8
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
std::vector< double > startTimeOfHitInChannel_7
std::vector< double > startTimeOfHitInChannel_10
art::ServiceHandle< beamlineutil::BeamlineChannelMap > fChannelMap
TRandom3 r(0)
Raw data definitions for beamline data used in NOvA test beam experiment.
std::vector< double > tof_p1_p4
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:464
Double_t sum
Definition: plot.C:31
Encapsulation of reconstructed PID information from detectors in the beamline (ToF, WCs, Cherenkov). Part of beamline reconstruction for NOvA test beam.
Float_t e
Definition: plot.C:35
std::vector< double > momentum
void decimal2binary(unsigned int decimalNumber, unsigned int *binaryNumber)
Channel mapping service which may be used to interpret the channels which are read out by the various...
std::vector< double > startTimeOfHitInChannel_12
Encapsulation of reconstructed track in the muon stack downstream of test beam detector. Part of beamline reconstruction for NOvA test beam.