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  (void)passSelection; // quiet down the compiler for now
338  // if (passSelection) {
339  // countPassSelectionEvents++;
340 
341  // if (WCs.size() == 1) {
342  // for (auto itr = WCs.begin(); itr != WCs.end(); itr++) {
343  // double mom = (*itr)->Momentum();
344  // // double mom_diff = TMath::Abs(mom - _MomentumFlip);
345  // // if (_MagnetCurrent != 500) {
346  // // mom = mom < _MomentumFlip ? _MomentumFlip + mom_diff :
347  // _MomentumFlip - mom_diff;
348  // // }
349  // momentum.push_back(mom);
350  // double tof = hitStartTOFDS.at(0) - hitStartTOFUS.at(0) + 10.8;
351  // tof_p1_p4.push_back(tof);
352 
353  // switch (_MagnetCurrent) {
354  // case (500): {
355  // hTOF_MNC_0500A->Fill(tof - 10.8);
356  // break;
357  // }
358  // case (1000): {
359  // hTOF_MNC_1000A->Fill(tof - 10.8);
360  // break;
361  // }
362  // case (1500): {
363  // hTOF_MNC_1500A->Fill(tof - 10.8);
364  // break;
365  // }
366  // case (2000): {
367  // hTOF_MNC_2000A->Fill(tof - 10.8);
368  // break;
369  // }
370  // default: {break;}
371  // }
372 
373  // bool foundElectron = false;
374  // for (auto CherenkovItr = Cherenkovs.begin(); CherenkovItr !=
375  // Cherenkovs.end(); CherenkovItr++) {
376  // double chevStart = (*CherenkovItr)->StartTimeInNanoSec();
377  // if (hitStartTOFUS.at(0) - 10.8 > chevStart) continue;
378  // if (hitStartTOFDS.at(0) > chevStart) continue;
379  // if (chevStart - hitStartTOFDS.at(0) > _CherenkovDelayAcceptance)
380  // continue; foundElectron = true; countElectronEvents++;
381  // }
382  // isElectron.push_back(foundElectron);
383  // } // end WCs elem loop
384  // } // make sure only one and only one track is associated with this
385  // trigger
386 
387  // for (auto trigItr = Triggers.begin(); trigItr != Triggers.end();
388  // trigItr++) {
389  // countTotalTriggers++;
390  // unsigned int pattern = (*trigItr)->Pattern();
391 
392  // unsigned int binaryNumber[32];
393  // decimal2binary(pattern, &binaryNumber[0]);
394  // if (!binaryNumber[0] && binaryNumber[1] && binaryNumber[2] &&
395  // binaryNumber[3]) {
396  // hTrigPattern->Fill(0.5);
397  // } else if (binaryNumber[0] && !binaryNumber[1] && binaryNumber[2] &&
398  // binaryNumber[3]) {
399  // hTrigPattern->Fill(1.5);
400  // } else if (binaryNumber[0] && binaryNumber[1] && !binaryNumber[2] &&
401  // binaryNumber[3]) {
402  // hTrigPattern->Fill(2.5);
403  // } else if (binaryNumber[0] && binaryNumber[1] && binaryNumber[2] &&
404  // !binaryNumber[3]) {
405  // hTrigPattern->Fill(3.5);
406  // } else if (binaryNumber[0] && binaryNumber[1] && binaryNumber[2] &&
407  // binaryNumber[3]) {
408  // hTrigPattern->Fill(4.5);
409  // } else if (binaryNumber[4] && binaryNumber[5] && binaryNumber[6] &&
410  // binaryNumber[7]) {
411  // hTrigPattern->Fill(5.5);
412  // }
413  // }
414  // }
415  }
416 }
417 
418 std::vector<double> novatb::TriggerRateAna::findCluster(unsigned int tofIndex) {
419  std::vector<double> hitStartTOF;
420 
421  if (tofIndex == 0) {
422  for (unsigned int iNS = 0; iNS < 400; iNS++) {
423  double coincWindowStart = (double)iNS;
424  double coincWindowEnd = (double)iNS + 5.;
425 
426  bool coincidence_channel_8 = false;
427  bool coincidence_channel_9 = false;
428  bool coincidence_channel_10 = false;
429  bool coincidence_channel_11 = false;
430  double coincidence_time_channel_8 = 0;
431  double coincidence_time_channel_9 = 0;
432  double coincidence_time_channel_10 = 0;
433  double coincidence_time_channel_11 = 0;
434 
435  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_8.size();
436  iElem++) {
437  if (coincWindowStart < startTimeOfHitInChannel_8.at(iElem) &&
438  startTimeOfHitInChannel_8.at(iElem) < coincWindowEnd) {
439  coincidence_channel_8 = true;
440  coincidence_time_channel_8 = startTimeOfHitInChannel_8.at(iElem);
441  }
442  }
443 
444  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_9.size();
445  iElem++) {
446  if (coincWindowStart < startTimeOfHitInChannel_9.at(iElem) &&
447  startTimeOfHitInChannel_9.at(iElem) < coincWindowEnd) {
448  coincidence_channel_9 = true;
449  coincidence_time_channel_9 = startTimeOfHitInChannel_9.at(iElem);
450  }
451  }
452 
453  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_10.size();
454  iElem++) {
455  if (coincWindowStart < startTimeOfHitInChannel_10.at(iElem) &&
456  startTimeOfHitInChannel_10.at(iElem) < coincWindowEnd) {
457  coincidence_channel_10 = true;
458  coincidence_time_channel_10 = startTimeOfHitInChannel_10.at(iElem);
459  }
460  }
461 
462  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_11.size();
463  iElem++) {
464  if (coincWindowStart < startTimeOfHitInChannel_11.at(iElem) &&
465  startTimeOfHitInChannel_11.at(iElem) < coincWindowEnd) {
466  coincidence_channel_11 = true;
467  coincidence_time_channel_11 = startTimeOfHitInChannel_11.at(iElem);
468  }
469  }
470 
471  // check coinciden 3 out of 4 or 4 out of 4
472  bool tofus_coinc_8_9_10 =
473  (coincidence_channel_8 && coincidence_channel_9 &&
474  coincidence_channel_10);
475  bool tofus_coinc_8_9_11 =
476  (coincidence_channel_8 && coincidence_channel_9 &&
477  coincidence_channel_11);
478  bool tofus_coinc_8_10_11 =
479  (coincidence_channel_8 && coincidence_channel_10 &&
480  coincidence_channel_11);
481  bool tofus_coinc_9_10_11 =
482  (coincidence_channel_9 && coincidence_channel_10 &&
483  coincidence_channel_11);
484  bool tofus_coinc = (tofus_coinc_8_9_10 || tofus_coinc_8_9_11 ||
485  tofus_coinc_8_10_11 || tofus_coinc_9_10_11);
486  if (tofus_coinc) {
487  double sum = coincidence_time_channel_8 + coincidence_time_channel_9 +
488  coincidence_time_channel_10 + coincidence_time_channel_11;
489  if (tofus_coinc_8_9_10 && tofus_coinc_8_9_11 && tofus_coinc_8_10_11 &&
490  tofus_coinc_9_10_11) {
491  hitStartTOF.push_back(sum / 4);
492  } else {
493  if (tofus_coinc_8_9_10) hitStartTOF.push_back(sum / 3);
494  if (tofus_coinc_8_9_11) hitStartTOF.push_back(sum / 3);
495  if (tofus_coinc_8_10_11) hitStartTOF.push_back(sum / 3);
496  if (tofus_coinc_9_10_11) hitStartTOF.push_back(sum / 3);
497  }
498  iNS = (unsigned int)coincWindowEnd;
499  } else
500  continue;
501  } // loop over 400ns step of 1ns
502  } // apply for US TOF
503 
504  else if (tofIndex == 1) {
505  for (unsigned int iNS = 0; iNS < 400; iNS++) {
506  double coincWindowStart = (double)iNS;
507  double coincWindowEnd = (double)iNS + 5.;
508 
509  bool coincidence_channel_12 = false;
510  bool coincidence_channel_13 = false;
511  bool coincidence_channel_14 = false;
512  bool coincidence_channel_15 = false;
513  double coincidence_time_channel_12 = 0;
514  double coincidence_time_channel_13 = 0;
515  double coincidence_time_channel_14 = 0;
516  double coincidence_time_channel_15 = 0;
517 
518  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_12.size();
519  iElem++) {
520  if (coincWindowStart < startTimeOfHitInChannel_12.at(iElem) &&
521  startTimeOfHitInChannel_12.at(iElem) < coincWindowEnd) {
522  coincidence_channel_12 = true;
523  coincidence_time_channel_12 = startTimeOfHitInChannel_12.at(iElem);
524  }
525  }
526 
527  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_13.size();
528  iElem++) {
529  if (coincWindowStart < startTimeOfHitInChannel_13.at(iElem) &&
530  startTimeOfHitInChannel_13.at(iElem) < coincWindowEnd) {
531  coincidence_channel_13 = true;
532  coincidence_time_channel_13 = startTimeOfHitInChannel_13.at(iElem);
533  }
534  }
535 
536  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_14.size();
537  iElem++) {
538  if (coincWindowStart < startTimeOfHitInChannel_14.at(iElem) &&
539  startTimeOfHitInChannel_14.at(iElem) < coincWindowEnd) {
540  coincidence_channel_14 = true;
541  coincidence_time_channel_14 = startTimeOfHitInChannel_14.at(iElem);
542  }
543  }
544 
545  for (unsigned int iElem = 0; iElem < startTimeOfHitInChannel_15.size();
546  iElem++) {
547  if (coincWindowStart < startTimeOfHitInChannel_15.at(iElem) &&
548  startTimeOfHitInChannel_15.at(iElem) < coincWindowEnd) {
549  coincidence_channel_15 = true;
550  coincidence_time_channel_15 = startTimeOfHitInChannel_15.at(iElem);
551  }
552  }
553 
554  // check coinciden 3 out of 4 or 4 out of 4
555  bool tofus_coinc_12_13_14 =
556  (coincidence_channel_12 && coincidence_channel_13 &&
557  coincidence_channel_14);
558  bool tofus_coinc_12_13_15 =
559  (coincidence_channel_12 && coincidence_channel_13 &&
560  coincidence_channel_15);
561  bool tofus_coinc_12_14_15 =
562  (coincidence_channel_12 && coincidence_channel_14 &&
563  coincidence_channel_15);
564  bool tofus_coinc_13_14_15 =
565  (coincidence_channel_13 && coincidence_channel_14 &&
566  coincidence_channel_15);
567  bool tofus_coinc = (tofus_coinc_12_13_14 || tofus_coinc_12_13_15 ||
568  tofus_coinc_12_14_15 || tofus_coinc_13_14_15);
569  if (tofus_coinc) {
570  double sum = coincidence_time_channel_12 + coincidence_time_channel_13 +
571  coincidence_time_channel_14 + coincidence_time_channel_15;
572  if (tofus_coinc_12_13_14 && tofus_coinc_12_13_15 &&
573  tofus_coinc_12_14_15 && tofus_coinc_13_14_15) {
574  hitStartTOF.push_back(sum / 4);
575  } else {
576  if (tofus_coinc_12_13_14) hitStartTOF.push_back(sum / 3);
577  if (tofus_coinc_12_13_15) hitStartTOF.push_back(sum / 3);
578  if (tofus_coinc_12_14_15) hitStartTOF.push_back(sum / 3);
579  if (tofus_coinc_13_14_15) hitStartTOF.push_back(sum / 3);
580  }
581  iNS = (unsigned int)coincWindowEnd;
582  } else
583  continue;
584  } // loop over 400ns step of 1ns
585  } // apply for US TOF
586 
587  return hitStartTOF;
588 }
589 
591  std::cout << "Total Triggers: " << totalEvents << std::endl;
592  std::cout << "Weird Triggers: " << countTotalTriggers << std::endl;
593 
594  std::cout << "countPassSelectionEvents: " << countPassSelectionEvents
595  << std::endl;
596  std::cout << "countElectronEvents: " << countElectronEvents << std::endl;
597 
598  gStyle->SetOptStat(0);
599 
600  TCanvas* c = new TCanvas("c", "c", 2400, 1200);
601  hTriggerBit->Draw();
602  c->SaveAs("TriggerBit.pdf");
603 
604  // TCanvas* tof_c = new TCanvas("tof_c", "tof_c", 1200, 1200);
605  // tof_c->Divide(2,2);
606  // tof_c->cd(1);
607  // hTOF_MNC_0500A->Draw();
608  // tof_c->cd(2);
609  // hTOF_MNC_1000A->Draw();
610  // tof_c->cd(3);
611  // hTOF_MNC_1500A->Draw();
612  // tof_c->cd(4);
613  // hTOF_MNC_2000A->Draw();
614  // tof_c->SaveAs("TOF.pdf");
615 
616  // const unsigned int nTotalPoints = momentum.size();
617 
618  // unsigned int countElectrons = 0;
619  // for (unsigned int i = 0; i < isElectron.size(); i++) {
620  // if (isElectron.at(i)) countElectrons++;
621  // }
622  // const unsigned int nElectronPoints = countElectrons;
623  // const unsigned int nNotElectronPoints = nTotalPoints - countElectrons;
624 
625  // double* tof_array_noE = new double[nNotElectronPoints];
626  // double* momentum_array_noE = new double[nNotElectronPoints];
627  // double* tof_array_E = new double[nElectronPoints];
628  // double* momentum_array_E = new double[nElectronPoints];
629  // unsigned int index_E = 0;
630  // unsigned int index_noE = 0;
631  // for (unsigned int i = 0; i < nTotalPoints; i++) {
632  // if (isElectron.at(i)) {
633  // *(tof_array_E + index_E) = tof_p1_p4.at(i);
634  // *(momentum_array_E + index_E) = momentum.at(i);
635  // index_E++;
636  // } else {
637  // *(tof_array_noE + index_noE) = tof_p1_p4.at(i);
638  // *(momentum_array_noE + index_noE) = momentum.at(i);
639  // index_noE++;
640  // }
641  // }
642 
643  // // double xmin = 100;
644  // // double xmax = 3000;
645  // // double ymin = 40;
646  // // double ymax = 100;
647 
648  // // double* momentum_for_particle_plots = new double[1000];
649  // // double* tof_for_deuteron_plots = new double[1000];
650  // // double* tof_for_proton_plots = new double[1000];
651  // // double* tof_for_kaon_plots = new double[1000];
652  // // double* tof_for_muon_plots = new double[1000];
653  // // double* tof_for_pion_plots = new double[1000];
654  // // double* tof_for_electron_plots = new double[1000];
655 
656  // // double momentum_step = (xmax - xmin) / 1000;
657  // // double distance_p1_p4_incm = 1316.22;
658  // // double mass_deuteron_inMeV = 1875.6;
659  // // double mass_proton_inMeV = 938.272;
660  // // double mass_kaon_inMeV = 493.677;
661  // // double mass_pion_inMeV = 139.570;
662  // // double mass_muon_inMeV = 105.658;
663  // // double mass_electron_inMeV = 0.511;
664  // // for (unsigned int i = 0; i < 1000; i++) {
665  // // *(momentum_for_particle_plots + i) = xmin + (double)i * momentum_step;
666  // // *(tof_for_deuteron_plots + i) =
667  // // calculate_tof_inns(distance_p1_p4_incm, mass_deuteron_inMeV,
668  // // *(momentum_for_particle_plots + i));
669  // // *(tof_for_proton_plots + i) =
670  // // calculate_tof_inns(distance_p1_p4_incm, mass_proton_inMeV,
671  // // *(momentum_for_particle_plots + i));
672  // // *(tof_for_kaon_plots + i) =
673  // // calculate_tof_inns(distance_p1_p4_incm, mass_kaon_inMeV,
674  // // *(momentum_for_particle_plots + i));
675  // // *(tof_for_pion_plots + i) =
676  // // calculate_tof_inns(distance_p1_p4_incm, mass_pion_inMeV,
677  // // *(momentum_for_particle_plots + i));
678  // // *(tof_for_muon_plots + i) =
679  // // calculate_tof_inns(distance_p1_p4_incm, mass_muon_inMeV,
680  // // *(momentum_for_particle_plots + i));
681  // // *(tof_for_electron_plots + i) =
682  // // calculate_tof_inns(distance_p1_p4_incm, mass_electron_inMeV,
683  // // *(momentum_for_particle_plots + i));
684  // // }
685 
686  // // TCanvas* canvas = new TCanvas("canvas", "canvas", 2400, 2400);
687  // pidPlot_noE =
688  // new TGraph(nNotElectronPoints, momentum_array_noE, tof_array_noE);
689  // pidPlot_noE->SetTitle("");
690  // pidPlot_noE->SetMarkerColor(kBlack);
691  // pidPlot_noE->SetMarkerStyle(20);
692  // pidPlot_noE->SetMarkerSize(1);
693  // pidPlot_noE->Draw("AP");
694 
695  // // TGraph* deuteronPlot = new TGraph(1000, momentum_for_particle_plots,
696  // // tof_for_deuteron_plots); TGraph* protonPlot = new TGraph(1000,
697  // // momentum_for_particle_plots, tof_for_proton_plots); TGraph* kaonPlot =
698  // // new TGraph(1000, momentum_for_particle_plots, tof_for_kaon_plots); TGraph*
699  // // pionPlot = new TGraph(1000, momentum_for_particle_plots,
700  // // tof_for_pion_plots); TGraph* muonPlot = new TGraph(1000,
701  // // momentum_for_particle_plots, tof_for_muon_plots); TGraph* electronPlot =
702  // // new TGraph(1000, momentum_for_particle_plots, tof_for_electron_plots);
703 
704  // if (nElectronPoints > 0) {
705  // pidPlot_E = new TGraph(nElectronPoints, momentum_array_E, tof_array_E);
706  // pidPlot_E->SetTitle("");
707  // pidPlot_E->SetMarkerColor(kRed);
708  // pidPlot_E->SetMarkerStyle(20);
709  // pidPlot_E->SetMarkerSize(1);
710  // pidPlot_E->Draw("P SAME");
711  // }
712 
713  // TFile* graphFile = new TFile("graphFile_10pB.root", "RECREATE");
714  // graphFile->cd();
715  // pidPlot_E->SetName("pidPlot_E");
716  // pidPlot_E->Write();
717  // pidPlot_noE->SetName("pidPlot_noE");
718  // pidPlot_noE->Write();
719  // graphFile->Write();
720  // graphFile->Close();
721 
722  // // deuteronPlot->SetLineColor(kMagenta);
723  // // protonPlot->SetLineColor(kRed);
724  // // kaonPlot->SetLineColor(kGreen - 2);
725  // // pionPlot->SetLineColor(kBlue + 1);
726  // // muonPlot->SetLineColor(kOrange - 2);
727  // // electronPlot->SetLineColor(kViolet - 1);
728 
729  // // deuteronPlot->Draw("L SAME");
730  // // protonPlot->Draw("L SAME");
731  // // kaonPlot->Draw("L SAME");
732  // // pionPlot->Draw("L SAME");
733  // // muonPlot->Draw("L SAME");
734  // // electronPlot->Draw("L SAME");
735 
736  // // TLegend* legend = new TLegend(0.5, 0.5, 0.8, 0.87);
737  // // legend->SetBorderSize(0);
738  // // legend->AddEntry(pidPlot_noE, "data", "p");
739  // // if (nElectronPoints > 0) {
740  // // legend->AddEntry(pidPlot_E, "data e", "p");
741  // // }
742  // // legend->AddEntry(deuteronPlot, "D", "l");
743  // // legend->AddEntry(protonPlot, "p", "l");
744  // // legend->AddEntry(kaonPlot, "K", "l");
745  // // legend->AddEntry(pionPlot, "#pi", "l");
746  // // legend->AddEntry(muonPlot, "#mu", "l");
747  // // legend->AddEntry(electronPlot, "e", "l");
748  // // legend->Draw();
749 
750  // // pidPlot_noE->GetXaxis()->SetTitle("Momentum (MeV/c)");
751  // // pidPlot_noE->GetXaxis()->CenterTitle();
752  // // pidPlot_noE->GetXaxis()->SetLimits(xmin, xmax);
753  // // pidPlot_noE->GetYaxis()->SetTitle("TOF (ns)");
754  // // pidPlot_noE->GetYaxis()->CenterTitle();
755  // // pidPlot_noE->GetYaxis()->SetRangeUser(ymin, ymax);
756 
757  // // auto prelim = TLatex(.9, .95, "NO#nuA Preliminary");
758  // // prelim.SetTextColor(4);
759  // // prelim.SetNDC();
760  // // prelim.SetTextSize(1.5 / 30.);
761  // // prelim.SetTextAlign(32);
762  // // prelim.Draw();
763 
764  // // canvas->SaveAs("pid.pdf");
765 }
766 
767 void novatb::TriggerRateAna::decimal2binary(unsigned int decimalNumber,
768  unsigned int* binaryNumber) {
769  for (int index = 31; index >= 0; index--) {
770  unsigned int reference = TMath::Power(2, (int)index);
771  if (decimalNumber >= reference) {
772  *(binaryNumber + index) = 1;
773  decimalNumber = decimalNumber - reference;
774  } else {
775  *(binaryNumber + index) = 0;
776  }
777  }
778 }
779 
780 // double novatb::TriggerRateAna::calculate_tof_inns(double distance_incm,
781 // double mass_inMeV, double momentum_inMeV) {
782 // double gamma = TMath::Sqrt((TMath::Power(momentum_inMeV, 2) +
783 // TMath::Power(mass_inMeV, 2)) / TMath::Power(mass_inMeV, 2)); double beta =
784 // TMath::Sqrt(1 - TMath::Power(1 / gamma, 2)); double velocity_mps =
785 // 299792458 * beta; double tof_s = (distance_incm / 100) / velocity_mps;
786 // return tof_s * 1E9;
787 // }
788 
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)
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
typedef void(XMLCALL *XML_ElementDeclHandler)(void *userData
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.
enum BeamMode string