RawEVD35tTree_module.cc
Go to the documentation of this file.
1 //Use this module to create a Tree for raw data display for 35t detector
2 //Mar. 20, 2014, Seongtae Park
3 //Mar. 2015: Added functionality for photon detectors, M Wallbank (wallbank@fnal.gov)
4 
5 #ifndef RawEVD35tTree_module
6 #define RawEVD35tTree_module
7 
8 // LArSoft includes
16 
17 // Data type includes
18 #include "lardataobj/RawData/raw.h"
25 
26 // Framework includes
31 #include "art_root_io/TFileService.h"
33 #include "canvas/Persistency/Common/FindManyP.h"
35 #include "fhiclcpp/ParameterSet.h"
36 
37 // ROOT includes.
38 #include "TH1.h"
39 #include "TH2.h"
40 #include "TTree.h"
41 #include "TLorentzVector.h"
42 #include "TVector3.h"
43 #include "TCanvas.h"
44 #include "TPad.h"
45 #include "TFile.h"
46 
47 // C++ Includes
48 #include <map>
49 #include <vector>
50 #include <algorithm>
51 #include <fstream>
52 #include <iostream>
53 #include <string>
54 #include <sstream>
55 #include <cmath>
56 
57 #ifdef __MAKECINT__
58 #pragma link C++ class vector<vector<int> >+;
59 #endif
60 
61 namespace AnalysisExample{
62 
64  public:
65 
66  explicit RawEVD35tTree(fhicl::ParameterSet const& pset);
67  virtual ~RawEVD35tTree();
68 
69  void beginJob();
70 
71  void beginRun(const art::Run& run);
72 
73  void reconfigure(fhicl::ParameterSet const& pset);
74 
75  void analyze(const art::Event& evt);
76 
77  void analyzeRCE(const std::vector< art::Ptr<raw::RawDigit> > RawDigits);
78 
79  void analyzeSSP(const std::vector< art::Ptr<raw::OpDetPulse> > RawPulses);
80 
81  private:
82 
83  // Parameters in .fcl file
89 
90  // Branch variables for tree
91  unsigned int fEvent;
92  unsigned int fRun;
93  unsigned int fSubRun;
94 
95  // TPC
96  unsigned int fNUCh;
97  unsigned int fNVCh;
98  unsigned int fNZ0Ch;
99  unsigned int fNZ1Ch;
100  // find channel boundaries for each view
101  unsigned int fUChanMin;
102  unsigned int fUChanMax;
103  unsigned int fVChanMin;
104  unsigned int fVChanMax;
105  unsigned int fZ0ChanMin;
106  unsigned int fZ0ChanMax;
107  unsigned int fZ1ChanMin;
108  unsigned int fZ1ChanMax;
109  unsigned int fNticks;
110  unsigned int fNofAPA;
111  unsigned int fChansPerAPA;
112  std::vector<unsigned int> fChan;
113  std::vector<unsigned int> fAPA;
114  std::vector<unsigned int> fPlane;
115  std::vector<std::vector<int> > fADC;
116  std::vector<std::vector<int> > fTDC;
117 
118  // SSP
119  std::vector<int> fOpChannel;
120  std::vector<int> fPMTFrame;
121  std::vector<std::vector<int> > fWaveform;
122 
124 
125  // Boolean to hold which data is present
126  bool fIsRCE = true;
127  bool fIsSSP = true;
128 
129  // TTree to hold information
130  TTree *tRD;
131 
132  }; // class RawEVD35tTree
133 
134  //-----------------------------------------------------------------------
135 
136  RawEVD35tTree::RawEVD35tTree(fhicl::ParameterSet const& parameterSet): EDAnalyzer(parameterSet) {
137  this->reconfigure(parameterSet);
138  }
139 
140  //-----------------------------------------------------------------------
141 
142  // read in the parameters from the .fcl file
144  fTPCInput = p.get< std::string >("TPCInputModule");
145  fTPCInstance = p.get< std::string >("TPCInstanceName");
146  fSSPInput = p.get< std::string >("SSPInputModule");
147  fSSPInstance = p.get< std::string >("SSPInstanceName");
148  auto const *fDetProp = lar::providerFrom<detinfo::DetectorPropertiesService>();
149  fNticks = fDetProp->NumberTimeSamples();
150  return;
151  }
152 
153 
154  //-----------------------------------------------------------------------
155 
157 
158  //-----------------------------------------------------------------------
159 
161 
163 
164  // Accquiring geometry data
167 
168  // Defining a Tree
169  tRD = tfs->make<TTree>("RawData","Raw Data Display");
170  tRD->Branch("Run",&fRun,"Run/i");
171  tRD->Branch("SubRun",&fSubRun,"SubRun/i");
172  tRD->Branch("Event",&fEvent,"Event/i");
173 
174  // TPC info
175  tRD->Branch("ADC",&fADC);
176  tRD->Branch("TDC",&fTDC);
177  tRD->Branch("Nticks",&fNticks,"Nticks/i");
178  tRD->Branch("NofUChan",&fNUCh,"NofUChan/i");
179  tRD->Branch("NofVChan",&fNVCh,"NofVChan/i");
180  tRD->Branch("NofZ0Chan",&fNZ0Ch,"NofZ0Chan/i");
181  tRD->Branch("NofZ1Chan",&fNZ1Ch,"NofZ1Chan/i");
182  tRD->Branch("Chan",&fChan);
183  tRD->Branch("APA",&fAPA);
184  tRD->Branch("Plane",&fPlane);
185 
186  // SSP info
187  tRD->Branch("Waveform",&fWaveform);
188  tRD->Branch("OpChan",&fOpChannel);
189  tRD->Branch("PMTFrame",&fPMTFrame);
190 
191  // loop through channels in the first APA to find the channel boundaries for each view
192  // will adjust for desired APA after
193  fUChanMin = 0;
196  fZ1ChanMax = fChansPerAPA - 1;
197  for ( unsigned int c = fUChanMin + 1; c < fZ1ChanMax; c++ ){
198  if ( fGeom->View(c) == geo::kV && fGeom->View(c-1) == geo::kU ){
199  fVChanMin = c;
200  fUChanMax = c - 1;
201  }
202  if ( fGeom->View(c) == geo::kZ && fGeom->View(c-1) == geo::kV ){
203  fZ0ChanMin = c;
204  fVChanMax = c-1;
205  }
206  if ( fGeom->View(c) == geo::kZ && fGeom->ChannelToWire(c)[0].TPC == fGeom->ChannelToWire(c-1)[0].TPC + 1 ){
207  fZ1ChanMin = c;
208  fZ0ChanMax = c-1;
209  }
210  }
211 
215  fNZ1Ch=fZ1ChanMax-fZ1ChanMin+1;
216 
217  // ofstream outfile;
218  // outfile.open("msglog.txt");
219  //outfile<<fNUCh<<" "<<fNVCh<<" "<<fNZ0Ch<<" "<<fNZ1Ch<<std::endl;
220  }
221 
222  //-----------------------------------------------------------------------
223 
225  }
226 
227 
228  //-----------------------------------------------------------------------
229 
231 
232  fEvent = event.id().event();
233  fRun = event.run();
234  fSubRun = event.subRun();
235  // unsigned int tpcid, cryoid;
236 
237  // Get the objects holding raw information: RawDigit for TPC data, OpDetPulse for SSP data
239  event.getByLabel(fTPCInput, fTPCInstance, RawTPC);
241  event.getByLabel(fSSPInput, fSSPInstance, RawSSP);
242 
243  // Check to see which data is present
244  try { RawTPC->size(); }
245  catch(std::exception const&) { fIsRCE = false; }
246  try { RawSSP->size(); }
247  catch(std::exception const&) { fIsSSP = false; }
248 
249  // Fill pointer vectors - more useful form for the raw data
250  std::vector< art::Ptr<raw::RawDigit> > RawDigits;
251  if (fIsRCE) art::fill_ptr_vector(RawDigits, RawTPC);
252  std::vector< art::Ptr<raw::OpDetPulse> > RawPulses;
253  if (fIsSSP) art::fill_ptr_vector(RawPulses, RawSSP);
254 
255  fADC.clear();
256  fTDC.clear();
257  fChan.clear();
258  fAPA.clear();
259  fPlane.clear();
260  fWaveform.clear();
261 
262  if (fIsRCE) analyzeRCE(RawDigits);
263  if (fIsSSP) analyzeSSP(RawPulses);
264 
265  tRD->Fill();
266  return;
267 
268  }
269 
271 
272  // Loop over RawDigits (TPC channels)
273  for(size_t d = 0; d < RawDigits.size(); d++) {
274 
276  digit = RawDigits.at(d);
277 
278  // Get the channel number for this digit
279  uint32_t chan = digit->Channel();
280  int nSamples = digit->Samples();
281  unsigned int apa = std::floor( chan/fChansPerAPA );
282  // tpcid=fGeom->ChannelToWire(chan)[0].TPC;
283  // cryoid=fGeom->ChannelToWire(chan)[0].Cryostat;
284  int Plane=0;
285  if( fGeom->View(chan) == geo::kU) Plane=0;
286  if( fGeom->View(chan) == geo::kV) Plane=1;
287  if ( fGeom->View(chan) == geo::kZ && fGeom->ChannelToWire(chan)[0].TPC % 2 == 0) Plane=2;
288  if ( fGeom->View(chan) == geo::kZ && fGeom->ChannelToWire(chan)[0].TPC % 2 == 1) Plane=3;
289 
290  std::vector<short> uncompressed(digit->Samples());
291  raw::Uncompress(digit->ADCs(), uncompressed, digit->Compression());
292 
293  std::vector<int> tmpADC;
294  std::vector<int> tmpTDC;
295  for (int i=0; i<nSamples; i++) {
296  short adc = uncompressed[i];
297  if(adc!=0) {
298  tmpADC.push_back(int(adc));
299  tmpTDC.push_back(int(i));
300  }
301  }
302 
303  fADC.push_back(tmpADC);
304  fTDC.push_back(tmpTDC);
305  fChan.push_back(chan);
306  fAPA.push_back(apa);
307  fPlane.push_back(Plane);
308 
309  } // End RawDigit function
310 
311  }
312 
314 
315  // Loop over OpDetPulses (SSP channels)
316  for(size_t p = 0; p < RawPulses.size(); p++) {
317 
318  art::Ptr<raw::OpDetPulse> pulsePtr;
319  pulsePtr = RawPulses.at(p);
320  raw::OpDetPulse pulse = *pulsePtr;
321 
322  unsigned short opchan = pulse.OpChannel();
323  unsigned short nsamples = pulse.Samples();
324  std::vector<short> waveform = pulse.Waveform();
325  unsigned int PMTFrame = pulse.PMTFrame();
326 
327  std::vector<int> tmpWaveform;
328  for (int i=0; i<nsamples; i++) {
329  short adc = waveform[i];
330  tmpWaveform.push_back(int(adc));
331  }
332 
333  fWaveform.push_back(tmpWaveform);
334  fOpChannel.push_back(opchan);
335  fPMTFrame.push_back(PMTFrame);
336 
337  } // End OpDetPulse function
338 
339  }
340 
342 
343 } // namespace
344 
345 #endif // RawEVD35tTree_module
Store parameters for running LArG4.
const ADCvector_t & ADCs() const
Reference to the compressed ADC count vector.
Definition: RawDigit.h:210
ULong64_t Samples() const
Number of samples in the uncompressed ADC data.
Definition: RawDigit.h:213
std::vector< std::vector< int > > fWaveform
RawEVD35tTree(fhicl::ParameterSet const &pset)
std::string string
Definition: nybbler.cc:12
Planes which measure V.
Definition: geo_types.h:126
std::vector< unsigned int > fPlane
Declaration of signal hit object.
unsigned short OpChannel() const
Definition: OpDetPulse.h:61
ChannelID_t Channel() const
DAQ channel this raw data was read from.
Definition: RawDigit.h:212
struct vector vector
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
Definition of basic raw digits.
Planes which measure Z direction.
Definition: geo_types.h:128
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:27
unsigned int Ncryostats() const
Returns the number of cryostats in the detector.
Particle class.
const std::vector< short > & Waveform() const
Definition: OpDetPulse.h:59
Definition: Run.h:21
unsigned int Nchannels() const
Returns the number of TPC readout channels in the detector.
void reconfigure(fhicl::ParameterSet const &pset)
unsigned short Samples() const
Definition: OpDetPulse.h:62
Planes which measure U.
Definition: geo_types.h:125
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:69
unsigned int uint32_t
Definition: stdint.h:126
Collect all the RawData header files together.
T get(std::string const &key) const
Definition: ParameterSet.h:231
void analyzeRCE(const std::vector< art::Ptr< raw::RawDigit > > RawDigits)
Declaration of cluster object.
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
Definition of data types for geometry description.
unsigned int NTPC(unsigned int cstat=0) const
Returns the total number of TPCs in the specified cryostat.
raw::Compress_t Compression() const
Compression algorithm used to store the ADC counts.
Definition: RawDigit.h:216
p
Definition: test.py:223
std::vector< std::vector< int > > fADC
void analyzeSSP(const std::vector< art::Ptr< raw::OpDetPulse > > RawPulses)
object containing MC truth information necessary for making RawDigits and doing back tracking ...
void analyze(const art::Event &evt)
std::vector< unsigned int > fChan
std::vector< std::vector< int > > fTDC
std::vector< unsigned int > fAPA
TCEvent evt
Definition: DataStructs.cxx:7
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:291
void Uncompress(const std::vector< short > &adc, std::vector< short > &uncompressed, raw::Compress_t compress)
Uncompresses a raw data buffer.
Definition: raw.cxx:755
recob::tracking::Plane Plane
Definition: TrackState.h:17
art::ServiceHandle< geo::Geometry > fGeom
unsigned int PMTFrame() const
Definition: OpDetPulse.h:63
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
Event finding and building.
unsigned int run
void beginRun(const art::Run &run)