MuonTrackHits_module.cc
Go to the documentation of this file.
1 
2 ////////////////////////////////////////////////////////////////////////
3 // Class: MuonTrackHits //
4 // Module Type: filter //
5 // File: MuonTrackHits_module.cc //
6 // //
7 // Created on Thurs Dec 10 2015 by Diana Patricia Mendez //
8 // Based on MuondEdxAna //
9 // //
10 // Added two new TTrees //
11 // *Eliminated the cut to tricells within tracklength < 200cm //
12 // *fTree stores all the hits (hitId: 0=avghit, 1=zhit, 2=tricell) //
13 // *fTree_CalibTrack saves info about all the muon tracks //
14 ////////////////////////////////////////////////////////////////////////
15 
16 // framework includes....
24 #include "fhiclcpp/ParameterSet.h"
25 #include "art_root_io/TFileService.h"
28 #include "Geometry/LiveGeometry.h"
30 
31 // novasoft stuff
33 #include "RecoBase/Track.h"
34 #include "RecoBase/CellHit.h"
35 #include "Calibrator/Calibrator.h"
36 #include "NovaTimingUtilities/TimingUtilities.h"
37 
38 // C++
39 #include <memory>
40 #include <string>
41 #include <iostream>
42 #include <vector>
43 
44 // ROOT
45 #include "TH3D.h"
46 #include "TH2D.h"
47 #include "TProfile.h"
48 #include "TProfile2D.h"
49 #include "TFile.h"
50 #include "TVector3.h"
51 #include "TTree.h"
52 
53 namespace caldp {
54  class PCHit;
55 }
56 
57 namespace calib {
58  class MuonTrackHits;
59 
60  class MuonTrackHits : public art::EDFilter {
61  public:
62  explicit MuonTrackHits(fhicl::ParameterSet const & p);
63  virtual ~MuonTrackHits();
64 
65  bool filter(art::Event & e) override;
66  void reconfigure(fhicl::ParameterSet const & p);
67  void beginJob() override;
68  void endJob() override;
69  bool beginRun(art::Run & r) override;
70  bool endRun(art::Run & r) override;
71 
72  private:
73  double getPECorr(art::Event & e,caldp::PCHit const &pchit);
74  double getPECorrToGeV(art::Event & e, caldp::PCHit const &pchit);
75 
76  //keep
77  TH2D* fdEdx;
78  TH2D* fdEdx_minus;
79  TH2D* fdEdx_plus;
80  TH2D* fdEdx_true;
81 
82 
83  TTree* fTree;
85 
86  int diblock;
87  int view;
88  double FlightLength;
90  double x;
91  double xCust;
92  double plane;
93  double path;
94  double trueE;
95  double PE;
96  double PECorr;
97  double tns;
98  uint32_t EventTimeHigh;
99  double cell;
100  double w;
101  double PECorrToGeV;
102  double RecoGeV;
104  int run;
105  int subRun;
106  int event;
107  int evt_year;
109  int evt_day;
110  int evt_hour;
111  int evt_min;
112  int evt_sec;
113  int APD;
116  int nhits;
117  int navghits;
118  int nzhits;
120  int nXYHits;
121  int hitId;
122 
123  bool badTNS;
124 
125  double xStart;
126  double yStart;
127  double zStart;
128  double xEnd;
129  double yEnd;
130  double zEnd;
131 
132  double dCosX;
133  double dCosY;
134  double dCosZ;
135  //double zenith;
136  //double azimuth;
138  int xPlanes;
139  int yPlanes;
141 
147 
148  double fSysShift;
149 
150  // *time_t* of start of Nova epoch, 01-Jan-2010 00:00:00, UTC
151  // This is the value subtracted from the UNIX time_t, timeval or timespec
152  // seconds field.
153  // Using unsigned long long rather than uint64_t to ensure consistent behavior
154  // on OSX and LINUX
155  //const uint32_t NOVA_EPOCH = 1262304000;
156  //const unsigned long long NOVA_TIME_FACTOR = 64000000;
157 
158  };
159 
160 
161  MuonTrackHits::MuonTrackHits(fhicl::ParameterSet const & p)
162  : EDFilter(p)
163  {
164  this->reconfigure(p);
165  produces< TH2D, art::InRun >("dEdx");
166  produces< TH2D, art::InRun >("dEdxMinus");
167  produces< TH2D, art::InRun >("dEdxPlus");
168  produces< TH2D, art::InRun >("dEdxTrue");
169  }
170 
172  {
173  // Clean up dynamic memory and other resources here.
174  }
175 
177  {
178 
179  run = (int)(e.run());
180  subRun = (int)(e.subRun());
181  event = (int)(e.event());
182 
183 
184  art::Timestamp ts = e.time();
185  EventTimeHigh = ts.timeHigh();
186  const time_t timeSec = ts.timeHigh();
187  struct tm* timeStruct = localtime(&timeSec);
188  evt_year=timeStruct->tm_year+1900;
189  evt_month=timeStruct->tm_mon+1;
190  evt_day=timeStruct->tm_mday;
191  evt_hour=timeStruct->tm_hour+1;
192  evt_min=timeStruct->tm_min;
193  evt_sec=timeStruct->tm_sec;
194 
195  /// Get the reconstructed tracks
197  e.getByLabel(fStopperLabel, tracks);
198  if( tracks->size() == 0 ) return false;
199 
203  if( !PCHitGetter.isValid() ) return false;
204 
205  // loop over stopping tracks in event
206  for( unsigned int i=0; i<tracks->size(); ++i ){
207 
208  //track level counters:
209  int windowTricells = 0;
210  int windowAvgHits = 0;
211  int windowZHits = 0;
212  int windowPlanes = 0;
213 
214  totalPlanes = 0;
215  xPlanes = 0;
216  yPlanes = 0;
217 
218  std::vector<bool> trackWindowPlanes (1000, false); //vector to prevent double counting planes
219 
220  std::vector<bool> trackPlanes (1000, false); //vector to prevent double counting total number of planes in track
221  std::vector<bool> trackXPlanes (1000, false); //vector to prevent double counting total number of x planes in track
222  std::vector<bool> trackYPlanes (1000, false); //vector to prevent double counting total number of y planes in track
223 
224  bool windowXView = false;
225  bool windowYView = false;
226 
227  //make custom track length to check against inbuilt length
228  TVector3 startXYZ = (*tracks)[i].Start();
229  TVector3 stopXYZ = (*tracks)[i].Stop();
230  double trackLength = (startXYZ-stopXYZ).Mag();
231  xStart = startXYZ.X();
232  yStart = startXYZ.Y();
233  zStart = startXYZ.Z();
234  xEnd = stopXYZ.X();
235  yEnd = stopXYZ.Y();
236  zEnd = stopXYZ.Z();
237  // zenith = startXYZ.Theta();
238  // azimuth = startXYZ.Phi();
239 
240  //cosine giving initial direction of track
241  dCosX = (*tracks)[i].Dir().X()/(*tracks)[i].Dir().Mag();
242  dCosY = (*tracks)[i].Dir().Y()/(*tracks)[i].Dir().Mag();
243  dCosZ = (*tracks)[i].Dir().Z()/(*tracks)[i].Dir().Mag();
244 
245 
246  // needed for track window, but good generally to reject protons
247  // if( (*tracks)[i].TotalLength() < 200. ) continue;
248 
251  int stopPlaneID = 10000;
252  int firstMuCPlane = -1;
253  bool isNearDet = false;
254 
255  //remove tracks that stop outside/ near edge of the active detector. New files won't need this as the selection is in the updated pchitslist_module.
257  distanceToLiveDetEdgeZ = fLiveGeom->DistToEdgeZ((*tracks)[i].Stop());
258 
259  //continue if track ends in the muon catcher
260  if(fGeo->DetId() == novadaq::cnv::kNEARDET){
261  isNearDet = true;
262  cid = fGeo->CellId(stopXYZ[0], stopXYZ[1], stopXYZ[2], 1.5, 1.5, 6., 1.0);
263  if(cid) stopPlaneID = fGeo->getPlaneID(cid);
264  firstMuCPlane = fGeo->FirstPlaneInMuonCatcher();
265  if(stopPlaneID >= firstMuCPlane) continue;
266  }
267 
269 
270  // loop over PCHits associated to stopping track
271  std::vector<art::Ptr<caldp::PCHit> > hits = PCHitGetter.at(i);
272  std::vector<art::Ptr<caldp::PCHit> > avghits = PCHitAvgGetter.at(i);
273  std::vector<art::Ptr<caldp::PCHit> > zhits = PCHitZGetter.at(i);
274 
275  navghits = avghits.size();
276  nzhits = zhits.size();
277  ntricells = hits.size();
278  nhits = hits.size() + avghits.size() + zhits.size();
279 
280  int tplane;
281  int trackView;
282 
283  for(unsigned int k=0; k<hits.size(); ++k){
284 
285  tplane = hits[k]->Plane();
286  trackView = (int)hits[k]->View();
287 
288  if(trackPlanes[tplane] == false) {
289  totalPlanes ++;
290  trackPlanes[tplane] = true;
291  };
292  if(trackView == 0 && trackXPlanes[tplane] == false){
293  xPlanes ++;
294  trackXPlanes[tplane] = true;
295  };
296  if(trackView == 1 && trackYPlanes[tplane] == false){
297  yPlanes ++;
298  trackYPlanes[tplane] = true;
299  };
300 
301  };
302 
303  for(unsigned int k=0; k<avghits.size(); ++k){
304  tplane = avghits[k]->Plane();
305  trackView = (int)avghits[k]->View();
306 
307  if(trackPlanes[tplane] == false){
308  totalPlanes ++;
309  trackPlanes[tplane] = true;
310  };
311  if(trackView == 0 && trackXPlanes[tplane] == false){
312  xPlanes ++;
313  trackXPlanes[tplane] = true;
314  };
315  if(trackView == 1 && trackYPlanes[tplane] == false){
316  yPlanes ++;
317  trackYPlanes[tplane] = true;
318  };
319 
320  };
321 
322  for(unsigned int k=0; k<zhits.size(); ++k){
323  tplane = zhits[k]->Plane();
324  trackView = (int)zhits[k]->View();
325 
326  if(trackPlanes[tplane] == false){
327  totalPlanes ++;
328  trackPlanes[tplane] = true;
329  };
330  if(trackView == 0 && trackXPlanes[tplane] == false){
331  xPlanes ++;
332  trackXPlanes[tplane] = true;
333  };
334  if(trackView == 1 && trackYPlanes[tplane] == false){
335  yPlanes ++;
336  trackYPlanes[tplane] = true;
337  };
338 
339  };
340 
341 
342  /*
343  if(xPlanes==0 || yPlanes==0){
344  std::cout << "***** HEY! *****" << std::endl;
345  std::cout << "run: " << run << ", subRun: " << subRun << ", event: " << event << " xPlanes: " << xPlanes << " yPlanes: " << yPlanes << ", TotalTrackLength: " << TotalTrackLength << std::endl;
346  std::cout << "xStart: " << xStart << ", yStart: " << yStart << ", zStart: " << zStart << std::endl;
347  std::cout << "xEnd: " << xEnd << ", yEnd: " << yEnd << ", zEnd: " << zEnd << std::endl;
348  std::cout << "* * * * * * * * * *" << std::endl;
349  }
350  */
351 
352  int cont;
353  for(unsigned int k=0; k<1000;k++){
354  cont = 0;
355  if(trackPlanes[k] == true){
356  for(unsigned int kk=0;kk<6;kk++){
357  if(trackPlanes[k+kk]==true) cont ++;
358  //else cont = cont;
359  };
360  if (cont>=5) break;
361  };
362  };
363 
364  if(cont>=5) continuity = 1;
365  else continuity = 0;
366 
367 
368  fTree_CalibTracks->Fill();
369 
370  // ---------------- Fill vector containing hits per APD ---------------------- //
371  std::vector<int> APDhits (100,0);
372  std::vector< std::vector<int> > DiblockAPDhits (100, std::vector<int> (100,0));
373  badTNS = false; //temporary handle for showery events in old data
374 
375  //avg. hits
376  if(avghits.size() != 0)
377  for( unsigned int p=0; p<avghits.size(); ++p ) {
378  APDhits[avghits[p]->APD()] += 1;
379  DiblockAPDhits.at(avghits[p]->Diblock()).at(avghits[p]->APD()) += 1;
380  if(avghits[p]->TNS() < 0 || avghits[p]->TNS() > 1400) badTNS = true;
381  FlightLength = avghits[p]->FlightLen();
382  TotalTrackLength = (*tracks)[i].TotalLength();
383  x = trackLength - FlightLength;
384  if(x<100 || x>200) continue;
385  windowAvgHits ++; //track level counters
386  }
387  //z hits
388  if(zhits.size() != 0)
389  for( unsigned int q=0; q<zhits.size(); ++q ){
390  APDhits[zhits[q]->APD()] += 1;
391  DiblockAPDhits.at(zhits[q]->Diblock()).at(zhits[q]->APD()) += 1;
392  if(zhits[q]->TNS() < 0 || zhits[q]->TNS() > 1400) badTNS = true;
393  FlightLength = zhits[q]->FlightLen();
394  TotalTrackLength = (*tracks)[i].TotalLength();
395  x = trackLength - FlightLength;
396  if(x<100 || x>200) continue;
397  windowZHits ++; //track level counters
398  }
399  // tricell hits
400  if(hits.size() != 0)
401  for( unsigned int n=0; n<hits.size(); ++n ){
402  APDhits[hits[n]->APD()] += 1;
403  DiblockAPDhits.at(hits[n]->Diblock()).at(hits[n]->APD()) += 1;
404  if(hits[n]->TNS() < 0 || hits[n]->TNS() > 1400) badTNS = true;
405  FlightLength = hits[n]->FlightLen();
406  TotalTrackLength = (*tracks)[i].TotalLength();
407  x = trackLength - FlightLength;
408  if(x<100 || x>200) continue;
409  windowTricells ++; //track level counters
410  }
411  // ---------------- Filled APD vector -------------------------------------------- //
412 
413 
414  for( unsigned int j=0; j<avghits.size(); ++j) {
415  if(isNearDet)
416  if(avghits[j]->Plane() >= firstMuCPlane) continue;
417 
418  hitId = 0;
419  nXYHits = avghits.size();
420  diblock = avghits[j]->Diblock();
421  view = (int)avghits[j]->View();
422  FlightLength = avghits[j]->FlightLen();
423  TotalTrackLength= (*tracks)[i].TotalLength();
424  x = trackLength - FlightLength;
425  plane = avghits[j]->Plane();
426  path = avghits[j]->Path();
427  trueE = avghits[j]->TrueMeV();
428  PE = avghits[j]->PE();
429  PECorr = this->getPECorr(e,*avghits[j]);
430  PECorrToGeV = this->getPECorrToGeV(e,*avghits[j]);
432  tns = avghits[j]->TNS();
433  cell = avghits[j]->Cell();
434  w = avghits[j]->W();
435  APD = avghits[j]->APD();
436  sameAPDhitsEvent= APDhits[APD];
437  sameAPDhitDiblockEvent= DiblockAPDhits[diblock][APD];
438 
439  fTree->Fill(); //fill ttree
440 
441  }
442 
443 
444  for( unsigned int j=0; j<zhits.size(); ++j) {
445  if(isNearDet)
446  if(zhits[j]->Plane() >= firstMuCPlane) continue;
447 
448  hitId = 1;
449  nXYHits = zhits.size();
450  diblock = zhits[j]->Diblock();
451  view = (int)zhits[j]->View();
452  FlightLength = zhits[j]->FlightLen();
453  TotalTrackLength= (*tracks)[i].TotalLength();
454  x = trackLength - FlightLength;
455  plane = zhits[j]->Plane();
456  path = zhits[j]->Path();
457  trueE = zhits[j]->TrueMeV();
458  PE = zhits[j]->PE();
459  PECorr = this->getPECorr(e,*zhits[j]);
460  PECorrToGeV = this->getPECorrToGeV(e,*zhits[j]);
462  tns = zhits[j]->TNS();
463  cell = zhits[j]->Cell();
464  w = zhits[j]->W();
465  APD = zhits[j]->APD();
466  sameAPDhitsEvent= APDhits[APD];
467  sameAPDhitDiblockEvent= DiblockAPDhits[diblock][APD];
468 
469  fTree->Fill(); //fill ttree
470 
471  }
472 
473 
474  for( unsigned int j=0; j<hits.size(); ++j ) { // loop over number of tricells per track
475 
476  if(isNearDet)
477  if(hits[j]->Plane() >= firstMuCPlane) continue;
478 
479  //test poisson lambda:
480  double poissonLambda = hits[j]->PoissonLambda();
481  std::cout<< "poissonLambda: " << poissonLambda <<std::endl;
482 
483  hitId = 2;
484  nXYHits = hits.size();
485  diblock = hits[j]->Diblock();
486  view = (int)hits[j]->View();
487  FlightLength = hits[j]->FlightLen();
488  TotalTrackLength= (*tracks)[i].TotalLength();
489  x = trackLength - FlightLength;
490  plane = hits[j]->Plane();
491  path = hits[j]->Path();
492  trueE = hits[j]->TrueMeV();
493  PE = hits[j]->PE();
494  PECorr = this->getPECorr(e,*hits[j]);
495  PECorrToGeV = this->getPECorrToGeV(e,*hits[j]);
497  tns = hits[j]->TNS();
498  cell = hits[j]->Cell();
499  w = hits[j]->W();
500  APD = hits[j]->APD();
501  sameAPDhitsEvent= APDhits[APD];
502  sameAPDhitDiblockEvent= DiblockAPDhits[diblock][APD];
503 
504  fTree->Fill(); //fill ttree
505 
506 
507  if(x<100 || x>200) continue;
508 
509  //track level counters
510  if(view == 0) windowXView = true;
511  if(view == 1) windowYView = true;
512  //check if plane already counted
513  if(trackWindowPlanes[plane] == false) windowPlanes ++;
514  trackWindowPlanes[plane] = true;
515 
516 
517  } // end loop over tricells
518 
519 
520  int nViewsInWindow = 0;
521  if(windowXView == true) nViewsInWindow++;
522  if(windowYView == true) nViewsInWindow++;
523 
524 
525 
526  }// loop over tracks
527 
528  return true;
529  }
530 
531 
533 
534  //Calibrator necessary to get Atten information
536 
537  //Make mock cellhit
538  rb::CellHit cellhit;
539  cellhit.SetPlane(pchit.Plane());
540  cellhit.SetCell(pchit.Cell());
541  cellhit.SetView(pchit.View());
542  cellhit.SetTNS(pchit.TNS(),pchit.GoodTime());
543  if(!e.isRealData())cellhit.SetMC(); // assume only filled for MC...
544  cellhit.SetPE(pchit.PE());
545 
546  double pecorr = calibrator->GetPECorr(cellhit, pchit.W());
547  return pecorr;
548  }
549 
551 
552  //Calibrator necessary to get Atten information
554  rb::CellHit cellhit;
555  if(!e.isRealData())cellhit.SetMC(); // assume only filled for MC...
556 
557  double pecorrtogev = calibrator->GetPECorrToGeVScale(cellhit);
558  return pecorrtogev;
559  }
560 
562  {
563  fPCHitLabel = pset.get< std::string >("PCHitLabel"); //Label of PCHitList module
564  fQualXYName = pset.get< std::string >("QualXYName"); //Instance label, "XY" quality hits
565  fQualZName = pset.get< std::string >("QualZName"); //Instance label, "Z" quality hits
566  fQualAvgName = pset.get< std::string >("QualAvgName"); //Instance label, "Avg" quality hits
567  fStopperLabel = pset.get<std::string>("StopperLabel");
568  fSysShift = pset.get<double>("SysShift"); // track end point systematic shift
569  }
570 
572  {
573  // histograms for output
575 
576  fTree = tfs->make<TTree>("fTree","tree to hold muon hits within each track variables");
577  fTree->Branch("nXYHits",&nXYHits);
578  fTree->Branch("nhits",&nhits);
579  fTree->Branch("hitId",&hitId);
580  fTree->Branch("xStart",&xStart);
581  fTree->Branch("yStart",&yStart);
582  fTree->Branch("zStart",&zStart);
583  fTree->Branch("dCosX",&dCosX);
584  fTree->Branch("dCosY",&dCosY);
585  fTree->Branch("dCosZ",&dCosZ);
586  // fTree->Branch("zenith",&zenith);
587  // fTree->Branch("azimuth",&azimuth);
588  fTree->Branch("x",&x);
589  fTree->Branch("badTNS",&badTNS);
590  fTree->Branch("run",&run);
591  fTree->Branch("diblock",&diblock);
592  fTree->Branch("view",&view);
593  fTree->Branch("FlightLength",&FlightLength);
594  fTree->Branch("TotalTrackLength",&TotalTrackLength);
595  fTree->Branch("xCust",&xCust);
596  fTree->Branch("distanceToLiveDetEdgeZ",&distanceToLiveDetEdgeZ);
597  fTree->Branch("plane",&plane);
598  fTree->Branch("path",&path);
599  fTree->Branch("trueE",&trueE);
600  fTree->Branch("PE",&PE);
601  fTree->Branch("PECorr",&PECorr);
602  fTree->Branch("tns",&tns);
603  fTree->Branch("cell",&cell);
604  fTree->Branch("w",&w);
605  fTree->Branch("PECorrToGeV",&PECorrToGeV);
606  fTree->Branch("RecoGeV",&RecoGeV);
607  fTree->Branch("EventTimeHigh",&EventTimeHigh);
608  fTree->Branch("evt_year",&evt_year);
609  fTree->Branch("evt_month",&evt_month);
610  fTree->Branch("evt_day",&evt_day);
611  fTree->Branch("evt_hour",&evt_hour);
612  fTree->Branch("evt_min",&evt_min);
613  fTree->Branch("evt_sec",&evt_sec);
614  fTree->Branch("APD",&APD);
615  fTree->Branch("sameAPDhitsEvent",&sameAPDhitsEvent);
616  fTree->Branch("sameAPDhitDiblockEvent",&sameAPDhitDiblockEvent);
617  fTree->Branch("totalPlanes",&totalPlanes);
618  fTree->Branch("xPlanes",&xPlanes);
619  fTree->Branch("yPlanes",&yPlanes);
620  fTree->Branch("continuity",&continuity);
621 
622  fTree_CalibTracks = tfs->make<TTree>("fTree_CalibTracks","tree to hold info about muon tracks");
623  fTree_CalibTracks->Branch("run",&run);
624  fTree_CalibTracks->Branch("subRun",&subRun);
625  fTree_CalibTracks->Branch("event",&event);
626  fTree_CalibTracks->Branch("nhits",&nhits);
627  fTree_CalibTracks->Branch("totalPlanes",&totalPlanes);
628  fTree_CalibTracks->Branch("xPlanes",&xPlanes);
629  fTree_CalibTracks->Branch("yPlanes",&yPlanes);
630  fTree_CalibTracks->Branch("continuity",&continuity);
631  fTree_CalibTracks->Branch("navghits",&navghits);
632  fTree_CalibTracks->Branch("nzhits",&nzhits);
633  fTree_CalibTracks->Branch("ntricells",&ntricells);
634  fTree_CalibTracks->Branch("totalPlanes",&totalPlanes);
635  fTree_CalibTracks->Branch("xStart",&xStart);
636  fTree_CalibTracks->Branch("yStart",&yStart);
637  fTree_CalibTracks->Branch("zStart",&zStart);
638  fTree_CalibTracks->Branch("xEnd",&xEnd);
639  fTree_CalibTracks->Branch("yEnd",&yEnd);
640  fTree_CalibTracks->Branch("zEnd",&zEnd);
641  fTree_CalibTracks->Branch("dCosX",&dCosX);
642  fTree_CalibTracks->Branch("dCosY",&dCosY);
643  fTree_CalibTracks->Branch("dCosZ",&dCosZ);
644  // fTree_CalibTracks->Branch("zenith",&zenith);
645  // fTree_CalibTracks->Branch("azimuth",&azimuth);
646  fTree_CalibTracks->Branch("TotalTrackLength",&TotalTrackLength);
647 
648 
649  fdEdx = tfs->make<TH2D>("dEdx",";Distance to track end (cm);PECorr / cm;",
650  50,0.,500.,200,0.,200.);
651  fdEdx_minus = tfs->make<TH2D>("dEdx_minus",";Distance to track end (cm);PECorr / cm;",
652  50,0.,500.,200,0.,200.);
653  fdEdx_plus = tfs->make<TH2D>("dEdx_plus",";Distance to track end (cm);PECorr / cm;",
654  50,0.,500.,200,0.,200.);
655  fdEdx_true = tfs->make<TH2D>("dEdx_true",";Distance to track end (cm);True MeV / cm;",
656  50,0.,500.,200,0.,5.);
657 
658  }
659 
661  return true;
662  }
663 
665 
667  std::unique_ptr<TH2D> pdEdx(tfs->make<TH2D>(*fdEdx));
668  std::unique_ptr<TH2D> pdEdx_minus(tfs->make<TH2D>(*fdEdx_minus));
669  std::unique_ptr<TH2D> pdEdx_plus(tfs->make<TH2D>(*fdEdx_plus));
670  std::unique_ptr<TH2D> pdEdx_true(tfs->make<TH2D>(*fdEdx_true));
671 
672  r.put(std::move(pdEdx),"dEdx");
673  r.put(std::move(pdEdx_minus),"dEdxMinus");
674  r.put(std::move(pdEdx_plus),"dEdxPlus");
675  r.put(std::move(pdEdx_true),"dEdxTrue");
676  return true;
677  }
678 
680  {
681 
682  }
683 
684 }// end namespace
void reconfigure(fhicl::ParameterSet const &p)
float TNS() const
Return uncorrected hit time.
Definition: PCHit.h:54
EventNumber_t event() const
double GetPECorrToGeVScale(rb::CellHit const &cellhit)
void SetTNS(float tns, bool good)
Definition: CellHit.h:56
bool beginRun(art::Run &r) override
const char * p
Definition: xmltok.h:285
bool GoodTime() const
Return quality of timing fit for cell.
Definition: PCHit.h:56
constexpr std::uint32_t timeHigh() const
Definition: Timestamp.h:34
int Cell() const
Return cell value.
Definition: PCHit.h:26
DEFINE_ART_MODULE(TestTMapFile)
"Pre-calibration hit". Common input to calibration procedures
Definition: PCHit.h:16
Definition: Run.h:21
Module that kips a configurable number of events between each that it allows through. Note that this module really skips (N-1) events, it uses a simple modular division as its critera. This module will cut down the data sample to 1/N of its original size.
void SetView(geo::View_t view)
Definition: CellHit.h:54
void SetPlane(unsigned short plane)
Definition: CellHit.h:53
float W() const
Return W value.
Definition: PCHit.h:44
void SetCell(unsigned short cell)
Definition: CellHit.h:52
bool isRealData() const
CDPStorage service.
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Definition: DataViewImpl.h:446
void hits()
Definition: readHits.C:15
void SetPE(float pe)
Definition: CellHit.h:55
Timestamp time() const
float PE() const
Return PE value.
Definition: PCHit.h:38
double DistToEdgeZ(TVector3 vertex)
void beginJob()
Encapsulate the geometry of one entire detector (near, far, ndos)
std::void_t< T > n
T get(std::string const &key) const
Definition: ParameterSet.h:231
int getPlaneID(const CellUniqueId &id) const
novadaq::cnv::DetId DetId() const
Prefer ds::DetectorService::DetId() instead.
Definition: GeometryBase.h:243
Near Detector in the NuMI cavern.
const CellUniqueId CellId(const double &x, const double &y, const double &z, double dxds=0., double dyds=0., double dzds=1., double step=0.01) const
SubRunNumber_t subRun() const
const double j
Definition: BetheBloch.cxx:29
RunNumber_t run() const
int Plane() const
Return plane value.
Definition: PCHit.h:24
geo::View_t View() const
Return view.
Definition: PCHit.h:36
bool filter(art::Event &e) override
Definition: run.py:1
OStream cout
Definition: OStream.cxx:6
void SetMC(bool isMC=true)
Definition: RawDigit.h:106
Histograms used by attenuation calibration.
unsigned long long int CellUniqueId
Definition: CellUniqueId.h:15
A rawdata::RawDigit with channel information decoded.
Definition: CellHit.h:27
double getPECorr(art::Event &e, caldp::PCHit const &pchit)
bool endRun(art::Run &r) override
float Mag() const
TRandom3 r(0)
double GetPECorr(rb::CellHit const &cellhit, double w)
double getPECorrToGeV(art::Event &e, caldp::PCHit const &pchit)
Float_t e
Definition: plot.C:35
ProductID put(std::unique_ptr< PROD > &&edp, FullSemantic< Level::Run > const semantic)
Definition: DataViewImpl.h:730
const unsigned int FirstPlaneInMuonCatcher() const
Returns the index of the first plane contained in the muon catcher.
enum BeamMode string