DetRespDrift_module.cc
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 // Class: DetRespDrift
3 // Module Type: filter
4 // File: DetRespDrift_module.cc
5 //
6 // Generated at Wed Jul 8 2015 by Luke Vinton
7 ////////////////////////////////////////////////////////////////////////
8 
9 
10 // framework includes....
18 #include "fhiclcpp/ParameterSet.h"
22 
23 // novasoft stuff
25 #include "RecoBase/Track.h"
26 #include "RecoBase/CellHit.h"
27 #include "Calibrator/Calibrator.h"
29 
30 // C++
31 #include <memory>
32 #include <string>
33 #include <iostream>
34 #include <vector>
35 
36 // ROOT
37 #include "TH2D.h"
38 #include "TFile.h"
39 #include "TVector3.h"
40 #include "TTree.h"
41 
42 namespace caldp {
43  class PCHit;
44 }
45 
46 
47 
48 namespace calib {
49  class DetRespDrift;
50 
51  class DetRespDrift : public art::EDFilter {
52  public:
53  explicit DetRespDrift(fhicl::ParameterSet const & p);
54  virtual ~DetRespDrift();
55 
56  bool filter(art::Event & e) override;
57  void reconfigure(fhicl::ParameterSet const & p);
58  void beginJob() override;
59  void endJob() override;
60  bool beginRun(art::Run & r) override;
61  bool endRun(art::Run & r) override;
62 
63  private:
64  double getPECorr(art::Event & e,caldp::PCHit const &pchit);
65 
66  TTree* fTree;
67  //declare variables for ttree
68  int run;
69  uint32_t EventTimeHigh;
70  double path;
71  double trueE;
72  double PECorr;
73  double pe;
75 
76 
82 
83  };
84 
85 
86  DetRespDrift::DetRespDrift(fhicl::ParameterSet const & p)
87  {
88  this->reconfigure(p);
89 
90  }
91 
92  DetRespDrift::~DetRespDrift()
93  {
94  // Clean up dynamic memory and other resources here.
95  }
96 
97  bool DetRespDrift::filter(art::Event & e)
98  {
99  /// Get the reconstructed tracks
101  e.getByLabel(fLabel, tracks);
102  if( tracks->size() == 0 ) return true;
103 
104  run = (int)(e.run());
105  art::Timestamp ts = e.time();
106  EventTimeHigh = ts.timeHigh();
107 
108  art::FindManyP<caldp::PCHit> PCHitXYGetter(tracks, e, art::InputTag(fPCHitLabel,fQualXYName));
109  art::FindManyP<caldp::PCHit> PCHitZGetter(tracks, e, art::InputTag(fPCHitLabel,fQualZName));
110  art::FindManyP<caldp::PCHit> PCHitAvgGetter(tracks, e, art::InputTag(fPCHitLabel,fQualAvgName));
111 
112  if( !PCHitAvgGetter.isValid() && !PCHitXYGetter.isValid() && !PCHitZGetter.isValid()) return false;
113 
114  // loop over stopping tracks in event
115  for( unsigned int i=0; i<tracks->size(); ++i ){
116 
117  // loop over PCHits associated to stopping track
118  std::vector<art::Ptr<caldp::PCHit> > hitsXY = PCHitXYGetter.at(i);
119  std::vector<art::Ptr<caldp::PCHit> > hitsZ = PCHitZGetter.at(i);
120  std::vector<art::Ptr<caldp::PCHit> > hitsAvg = PCHitAvgGetter.at(i);
121 
122  // XY hits
123  for( unsigned int j=0; j<hitsXY.size(); ++j ){
124 
125  path = hitsXY[j]->Path();
126  trueE = hitsXY[j]->TrueMeV();
127  PECorr = this->getPECorr(e,*hitsXY[j]);
128  pe = hitsXY[j]->PE();
129  pathQual = "XY";
130 
131  fTree->Fill(); //fill ttree
132  }
133 
134  // Z hits
135  for( unsigned int j=0; j<hitsZ.size(); ++j ){
136 
137  path = hitsZ[j]->Path();
138  trueE = hitsZ[j]->TrueMeV();
139  PECorr = this->getPECorr(e,*hitsZ[j]);
140  pe = hitsZ[j]->PE();
141  pathQual = "Z";
142 
143  fTree->Fill(); //fill ttree
144  }
145 
146  // Avg hits
147  for( unsigned int j=0; j<hitsAvg.size(); ++j ){
148 
149  path = hitsAvg[j]->Path();
150  trueE = hitsAvg[j]->TrueMeV();
151  PECorr = this->getPECorr(e,*hitsAvg[j]);
152  pe = hitsAvg[j]->PE();
153  pathQual = "Avg";
154 
155  fTree->Fill(); //fill ttree
156  }
157 
158  }
159 
160  return true;
161 
162  }
163 
164  double DetRespDrift::getPECorr(art::Event & e, caldp::PCHit const &pchit){
165 
166  //Calibrator necessary to get Atten information
168 
169  //Make mock cellhit
170  rb::CellHit cellhit;
171  cellhit.SetPlane(pchit.Plane());
172  cellhit.SetCell(pchit.Cell());
173  cellhit.SetView(pchit.View());
174  cellhit.SetTNS(pchit.TNS(),pchit.GoodTime());
175  if(!e.isRealData())cellhit.SetMC(); // assume only filled for MC...
176  cellhit.SetPE(pchit.PE());
177 
178  double pecorr = calibrator->GetPECorr(cellhit, pchit.W());
179  return pecorr;
180  }
181 
182  void DetRespDrift::reconfigure(fhicl::ParameterSet const & pset)
183  {
184  fPCHitLabel = pset.get< std::string >("PCHitLabel"); //Label of PCHitList module
185  fQualXYName = pset.get< std::string >("QualXYName"); //Instance label, "XY" quality hits
186  fQualZName = pset.get< std::string >("QualZName"); //Instance label, "Z" quality hits
187  fQualAvgName = pset.get< std::string >("QualAvgName"); //Instance label, "Avg" quality hits
188  fLabel = pset.get<std::string>("Label");
189 
190  }
191 
193  {
194  // histograms for output
196 
197  fTree = tfs->make<TTree>("fTree","tree to hold variables");
198  fTree->Branch("run",&run);
199  fTree->Branch("EventTimeHigh",&EventTimeHigh);
200  fTree->Branch("pe",&pe);
201  fTree->Branch("PECorr",&PECorr);
202  fTree->Branch("trueE",&trueE);
203  fTree->Branch("path",&path);
204  fTree->Branch("pathQual",&pathQual);
205 
206 
207  }
208 
209  bool DetRespDrift::beginRun(art::Run &r){
210  return true;
211  }
212 
213  bool DetRespDrift::endRun(art::Run &r){
214 
216 
217  return true;
218  }
219 
220  void DetRespDrift::endJob()
221  {
222 
223  }
224 }// end namespace
226 
float TNS() const
Return uncorrected hit time.
Definition: PCHit.h:52
void SetTNS(float tns, bool good)
Definition: CellHit.h:56
const char * p
Definition: xmltok.h:285
bool GoodTime() const
Return quality of timing fit for cell.
Definition: PCHit.h:54
constexpr std::uint32_t timeHigh() const
Definition: Timestamp.h:35
int Cell() const
Return cell value.
Definition: PCHit.h:26
bool isRealData() const
Definition: Event.h:83
DEFINE_ART_MODULE(TestTMapFile)
"Pre-calibration hit". Common input to calibration procedures
Definition: PCHit.h:16
Definition: Run.h:31
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:42
void SetCell(unsigned short cell)
Definition: CellHit.h:52
CDPStorage service.
void SetPE(float pe)
Definition: CellHit.h:55
float PE() const
Return PE value.
Definition: PCHit.h:38
void beginJob()
Encapsulate the geometry of one entire detector (near, far, ndos)
T get(std::string const &key) const
Definition: ParameterSet.h:231
const double j
Definition: BetheBloch.cxx:29
int Plane() const
Return plane value.
Definition: PCHit.h:24
geo::View_t View() const
Return view.
Definition: PCHit.h:36
Definition: run.py:1
void SetMC(bool isMC=true)
Definition: RawDigit.h:106
Histograms used by attenuation calibration.
const std::string path
Definition: plot_BEN.C:43
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
A rawdata::RawDigit with channel information decoded.
Definition: CellHit.h:27
T * make(ARGS...args) const
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
TRandom3 r(0)
Timestamp time() const
Definition: Event.h:61
double GetPECorr(rb::CellHit const &cellhit, double w)
Float_t e
Definition: plot.C:35
RunNumber_t run() const
Definition: Event.h:77