Public Member Functions | Private Attributes | List of all members
AnalysisExample::RawEVDDP Class Reference
Inheritance diagram for AnalysisExample::RawEVDDP:
art::EDAnalyzer art::detail::Analyzer art::detail::LegacyModule art::Observer art::ModuleBase

Public Member Functions

 RawEVDDP (fhicl::ParameterSet const &pset)
 
virtual ~RawEVDDP ()
 
void beginJob ()
 
void beginRun (const art::Run &run)
 
void reconfigure (fhicl::ParameterSet const &pset)
 
void analyze (const art::Event &evt)
 
- Public Member Functions inherited from art::EDAnalyzer
 EDAnalyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDAnalyzer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Analyzer
virtual ~Analyzer () noexcept
 
 Analyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 Analyzer (Table< Config > const &config)
 
void doBeginJob ()
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Observer
 ~Observer () noexcept
 
 Observer (Observer const &)=delete
 
 Observer (Observer &&)=delete
 
Observeroperator= (Observer const &)=delete
 
Observeroperator= (Observer &&)=delete
 
void registerProducts (ProductDescriptions &, ModuleDescription const &)
 
void fillDescriptions (ModuleDescription const &)
 
std::string const & processName () const
 
bool wantAllEvents () const
 
bool wantEvent (Event const &e)
 
fhicl::ParameterSetID selectorConfig () const
 
Handle< TriggerResultsgetTriggerResults (Event const &e) const
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Attributes

std::string fRawDigitLabel
 
unsigned int fNticks
 
unsigned int fChPerView
 
art::ServiceHandle< geo::GeometryfGeom
 
detinfo::DetectorProperties const * fDetProp = nullptr
 
std::vector< TH2I * > fTimeChanU
 
std::vector< TH2I * > fTimeChanV
 
std::vector< TH2I * > fTimeChanThumbU
 
std::vector< TH2I * > fTimeChanThumbV
 
std::vector< TH1I * > fADCMaxDistU
 
std::vector< TH1I * > fADCMaxDistV
 

Additional Inherited Members

- Public Types inherited from art::EDAnalyzer
using WorkerType = WorkerT< EDAnalyzer >
 
using ModuleType = EDAnalyzer
 
- Protected Member Functions inherited from art::Observer
 Observer (fhicl::ParameterSet const &config)
 
 Observer (std::vector< std::string > const &paths, fhicl::ParameterSet const &config)
 
detail::ProcessAndEventSelectorsprocessAndEventSelectors ()
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollectorconsumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 47 of file RawEVDDP_module.cc.

Constructor & Destructor Documentation

AnalysisExample::RawEVDDP::RawEVDDP ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 103 of file RawEVDDP_module.cc.

104  : EDAnalyzer(parameterSet)
105  {
106  this->reconfigure(parameterSet);
107  }
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:27
void reconfigure(fhicl::ParameterSet const &pset)
AnalysisExample::RawEVDDP::~RawEVDDP ( )
virtual

Definition at line 122 of file RawEVDDP_module.cc.

122  {
123  }

Member Function Documentation

void AnalysisExample::RawEVDDP::analyze ( const art::Event evt)

Definition at line 227 of file RawEVDDP_module.cc.

228  {
229 
230  unsigned int tpcid;//, cryoid;
231  //std::stringstream thumbnameZ0, thumbnameZ1;
232 
233  // get the objects holding all of the raw data information
235  event.getByLabel(fRawDigitLabel, Raw);
236 
237  // put it in a more easily usable form
238  std::vector< art::Ptr<raw::RawDigit> > Digits;
239  art::fill_ptr_vector(Digits, Raw);
240 
241  //loop through all RawDigits (over entire channels)
242 
243  for(size_t d = 0; d < Digits.size(); d++)
244  {
246  digit=Digits.at(d);
247 
248  // get the channel number for this digit
249  uint32_t chan = digit->Channel();
250  tpcid = fGeom->ChannelToWire(chan)[0].TPC;
251  //cryoid = fGeom->ChannelToWire(chan)[0].Cryostat;
252 
253  // ok this loop does not really work for more than 1 cryostat
254 
255  std::vector<short> uncompressed(digit->Samples());
256  raw::Uncompress(digit->ADCs(), uncompressed, digit->Compression());
257 
258  short maxadc = 0;
259  if( fGeom->View(chan) == geo::kU )
260  {
261  for(unsigned int l=0;l<uncompressed.size();l++)
262  {
263  short adcval = uncompressed.at(l);
264  if(uncompressed.at(l) > maxadc)
265  maxadc = adcval;
266 
267  fTimeChanU[tpcid]->Fill(chan - 2*tpcid*fChPerView, l, adcval);
268  if(adcval>0)
269  fTimeChanThumbU[tpcid]->Fill(chan - 2*tpcid*fChPerView,l, adcval);
270  }
271  if(maxadc>0) fADCMaxDistU[tpcid]->Fill(maxadc);
272  }
273  else if( fGeom->View(chan) == geo::kV )
274  {
275  for(unsigned int l=0;l<uncompressed.size();l++)
276  {
277  short adcval = uncompressed.at(l);
278  if(uncompressed.at(l) > maxadc)
279  maxadc = adcval;
280 
281  fTimeChanV[tpcid]->Fill(chan - 2*tpcid*fChPerView-fChPerView, l, adcval);
282  if(adcval>0)
283  fTimeChanThumbV[tpcid]->Fill(chan - 2*tpcid*fChPerView-fChPerView,l, adcval);
284  }
285  if(maxadc > 0) fADCMaxDistV[tpcid]->Fill(maxadc);
286  }
287 
288  } // end RawDigit loop
289 
290  return;
291  }
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
Planes which measure V.
Definition: geo_types.h:126
ChannelID_t Channel() const
DAQ channel this raw data was read from.
Definition: RawDigit.h:212
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
std::vector< TH2I * > fTimeChanU
static QStrList * l
Definition: config.cpp:1044
std::vector< TH2I * > fTimeChanV
Planes which measure U.
Definition: geo_types.h:125
std::vector< TH2I * > fTimeChanThumbU
unsigned int uint32_t
Definition: stdint.h:126
std::vector< TH1I * > fADCMaxDistU
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
raw::Compress_t Compression() const
Compression algorithm used to store the ADC counts.
Definition: RawDigit.h:216
art::ServiceHandle< geo::Geometry > fGeom
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
std::vector< TH2I * > fTimeChanThumbV
std::vector< TH1I * > fADCMaxDistV
void AnalysisExample::RawEVDDP::beginJob ( )
virtual

Reimplemented from art::EDAnalyzer.

Definition at line 127 of file RawEVDDP_module.cc.

127  {
128  // Access ART's TFileService, which will handle creating and writing
129  // histograms and n-tuples for us.
131 
132  //Histogram names and titles
133  std::stringstream name, title;
134 
135  TH2I* TempHisto2;
136  TH1I* TempHisto1;
137 
138  //ofstream outfile;
139  //outfile.open("msglog.txt");
140 
141 
143  mf::LogInfo("RawEVDP")<< "Number of channels per view is "<<fChPerView;
144 
145  if(fGeom->Nviews() != 2 )
146  throw cet::exception("RawEVDDP") << "For DUNE DP expected to have only 2 views ";
147 
148 
149  //
150  unsigned int minT = 0;
151  unsigned int maxT = fNticks;
152  unsigned int binT = (maxT-minT);
153 
154  float minADC = 0.0;
155  float maxADC = 4096;
156  int nBins = (int)((maxADC - minADC)/2.0);
157 
158  for(unsigned int i=0;i<fGeom->NTPC();i++)
159  {
160 
161  name.str("");
162  name << "fTimeChanU";
163  name << i;
164  title.str("");
165  title << "Time vs Channel(Plane U, CRM";
166  title << i<<")";
167  TempHisto2 = tfs->make<TH2I>(name.str().c_str(),title.str().c_str(), fChPerView, 0, fChPerView, binT, minT, maxT);
168  fTimeChanU.push_back(TempHisto2);
169 
170  name.str("");
171  name << "fTimeChanThumbU";
172  name << i;
173  title.str("");
174  title << "Time vs Channel(Plane U, CRM";
175  title << i<<")";
176  TempHisto2 = tfs->make<TH2I>(name.str().c_str(),title.str().c_str(), 32, 0, fChPerView, 32, minT, maxT);
177  fTimeChanThumbU.push_back(TempHisto2);
178 
179  name.str("");
180  name << "fADCMaxU";
181  name << i;
182  title.str("");
183  title << "Max ADC per channel (Plane U, CRM";
184  title << i<<")";
185  TempHisto1 = tfs->make<TH1I>(name.str().c_str(),title.str().c_str(), nBins, minADC, maxADC);
186  fADCMaxDistU.push_back(TempHisto1);
187 
188  //
189  name.str("");
190  name << "fTimeChanV";
191  name << i;
192  title.str("");
193  title << "Time vs Channel(Plane V, CRM";
194  title << i<<")";
195  TempHisto2 = tfs->make<TH2I>(name.str().c_str(),title.str().c_str(), fChPerView, 0, fChPerView, binT, minT, maxT);
196  fTimeChanV.push_back(TempHisto2);
197 
198  name.str("");
199  name << "fTimeChanThumbV";
200  name << i;
201  title.str("");
202  title << "Time vs Channel(Plane V, CRM";
203  title << i<<")";
204  TempHisto2 = tfs->make<TH2I>(name.str().c_str(),title.str().c_str(), 32, 0, fChPerView, 32, minT, maxT);
205  fTimeChanThumbV.push_back(TempHisto2);
206 
207  name.str("");
208  name << "fADCMaxV";
209  name << i;
210  title.str("");
211  title << "Max ADC per channel (Plane V, CRM";
212  title << i<<")";
213  TempHisto1 = tfs->make<TH1I>(name.str().c_str(),title.str().c_str(), nBins, minADC, maxADC);
214  fADCMaxDistV.push_back(TempHisto1);
215  }
216 
217  }
static QCString name
Definition: declinfo.cpp:673
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
std::vector< TH2I * > fTimeChanU
unsigned int Ncryostats() const
Returns the number of cryostats in the detector.
unsigned int Nchannels() const
Returns the number of TPC readout channels in the detector.
std::vector< TH2I * > fTimeChanV
std::vector< TH2I * > fTimeChanThumbU
std::vector< TH1I * > fADCMaxDistU
unsigned int NTPC(unsigned int cstat=0) const
Returns the total number of TPCs in the specified cryostat.
art::ServiceHandle< geo::Geometry > fGeom
unsigned int Nviews() const
Returns the number of views (different wire orientations)
std::vector< TH2I * > fTimeChanThumbV
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
std::vector< TH1I * > fADCMaxDistV
void AnalysisExample::RawEVDDP::beginRun ( const art::Run run)

Definition at line 220 of file RawEVDDP_module.cc.

220  {
221 
222  }
void AnalysisExample::RawEVDDP::reconfigure ( fhicl::ParameterSet const &  pset)

Definition at line 112 of file RawEVDDP_module.cc.

112  {
113  fDetProp = lar::providerFrom<detinfo::DetectorPropertiesService>();
114  fRawDigitLabel = p.get< std::string >("RawDigitLabel");
116  return;
117  }
detinfo::DetectorProperties const * fDetProp
std::string string
Definition: nybbler.cc:12
virtual unsigned int NumberTimeSamples() const =0
p
Definition: test.py:223

Member Data Documentation

std::vector<TH1I*> AnalysisExample::RawEVDDP::fADCMaxDistU
private

Definition at line 95 of file RawEVDDP_module.cc.

std::vector<TH1I*> AnalysisExample::RawEVDDP::fADCMaxDistV
private

Definition at line 96 of file RawEVDDP_module.cc.

unsigned int AnalysisExample::RawEVDDP::fChPerView
private

Definition at line 68 of file RawEVDDP_module.cc.

detinfo::DetectorProperties const* AnalysisExample::RawEVDDP::fDetProp = nullptr
private

Definition at line 87 of file RawEVDDP_module.cc.

art::ServiceHandle<geo::Geometry> AnalysisExample::RawEVDDP::fGeom
private

Definition at line 86 of file RawEVDDP_module.cc.

unsigned int AnalysisExample::RawEVDDP::fNticks
private

Definition at line 67 of file RawEVDDP_module.cc.

std::string AnalysisExample::RawEVDDP::fRawDigitLabel
private

Definition at line 63 of file RawEVDDP_module.cc.

std::vector<TH2I*> AnalysisExample::RawEVDDP::fTimeChanThumbU
private

Definition at line 92 of file RawEVDDP_module.cc.

std::vector<TH2I*> AnalysisExample::RawEVDDP::fTimeChanThumbV
private

Definition at line 93 of file RawEVDDP_module.cc.

std::vector<TH2I*> AnalysisExample::RawEVDDP::fTimeChanU
private

Definition at line 89 of file RawEVDDP_module.cc.

std::vector<TH2I*> AnalysisExample::RawEVDDP::fTimeChanV
private

Definition at line 90 of file RawEVDDP_module.cc.


The documentation for this class was generated from the following file: