DetAna_module.cc
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 /// \brief
3 /// \author messier@indiana.edu
4 ///////////////////////////////////////////////////////////////////////////////
5 
6 // C/C++ includes
7 #include <vector>
8 
9 // ART includes
19 #include "fhiclcpp/ParameterSet.h"
21 
22 // NOvA includes
24 #include "DAQDataFormats/NanoSliceVersionConvention.h"
25 #include "Geometry/Geometry.h"
27 #include "NovaDAQConventions/DAQConventions.h"
28 #include "RawData/DAQHeader.h"
29 #include "RawData/RawDigit.h"
30 
31 // ROOT includes
32 #include "TH1.h"
33 #include "TH2.h"
34 #include "TSystem.h"
35 
36 namespace mcchk {
37  /// A module to check the results from the Monte Carlo generator
38  class DetAna : public art::EDAnalyzer {
39 
40  public:
41  explicit DetAna(fhicl::ParameterSet const& pset);
42  virtual ~DetAna();
43 
44  void analyze(art::Event const& evt);
45  void beginRun(art::Run const& run);
46 
47  private:
48  std::string fRawDataLabel; ///< module label for the raw data
49 
50  double fRawDigit_FDOffset; ///< How far to shift FD Number of digit bins
51 
52  TH1F* fNumDigits; ///< number of digits per event
53  TH1F* fADC; ///< ADC value for each digit
54  TH1F* fADCPerSample; ///< ADC value for each sample
55  TH1F* fADCDiff[3]; ///< Difference between ADC(1/2/3) and ADC(0)
56  TH2F* fHitCellsXView; ///< cells hit in the x view
57  TH2F* fHitCellsYView; ///< cells hit in the y view
58  };
59 }
60 
61 ///////////////////////////////////////////////////////////////////////////////
62 namespace mcchk {
63  static int msg1cnt = 0;
64 
65  //......................................................................
67  : EDAnalyzer(pset),
68  fRawDataLabel(pset.get<std::string>("RawDataLabel", "daq")),
69  fRawDigit_FDOffset(pset.get<double>("RawDigit_FDOffset")),
70  fNumDigits(0) // Flag that all histograms are uninitialized
71  {
72  }
73 
74  //......................................................................
76  {
77  }
78 
79  //......................................................................
80  void DetAna::beginRun(art::Run const& )
81  {
82  if(fNumDigits) { return; } // Only make histograms once
83 
84  // Helpful for spacing/readability
85  // Also helpful in case histograms are made in loops over, say, axis labels
86  char HistoName[200];
87  char HistoTitle[200];
88  char TitleHelper[200];
89 
92 
93  const std::set<unsigned int> planesX = geo->GetPlanesByView(geo::kX);
94  const std::set<unsigned int> planesY = geo->GetPlanesByView(geo::kY);
95  int fPlanes = planesX.size() + planesY.size();
96 
97  int fCellsX = geo->Plane(*planesX.begin())->Ncells();
98  int fCellsY = geo->Plane(*planesY.begin())->Ncells();
99 
100  double edgeL = 0.;
101  double edgeR = 10000.;
102  if(geo->DetId() == novadaq::cnv::DetId::kFARDET) {
103  edgeL += fRawDigit_FDOffset;
104  edgeR += fRawDigit_FDOffset;
105  }
106  sprintf(HistoName, "fNumDigits");
107  sprintf(HistoTitle, "Number of RawDigits per Event;Digits;");
108  fNumDigits = tfs->make<TH1F>(HistoName, HistoTitle,
109  1000, edgeL, edgeR);
110 
111  sprintf(HistoName, "fADC");
112  sprintf(HistoTitle, "ADC Value of the RawDigit;ADC;");
113  fADC = tfs->make<TH1F>(HistoName, HistoTitle,
114  1000, 0., 1000.);
115 
116  sprintf(HistoName, "fADCPerSample");
117  sprintf(HistoTitle, "ADC Value for each Sample Read Out in the RawDigit;ADC;");
118  fADCPerSample = tfs->make<TH1F>(HistoName, HistoTitle,
119  1000, 0., 1000.);
120 
121  sprintf(TitleHelper, ") and ADC(0);#Delta ADC;");
122 
123  for(unsigned int i = 0, n = 3; i < n; ++i) {
124  sprintf(HistoName, "fADCDiff%d", i + 1);
125  sprintf(HistoTitle, "Difference Between ADC(%d%s", i + 1, TitleHelper);
126 
127  fADCDiff[i] = tfs->make<TH1F>(HistoName, HistoTitle, 1000, 0., 1000.);
128  }
129 
130  sprintf(TitleHelper, "Cells Hit in the");
131 
132  sprintf(HistoName, "fHitCellsXView");
133  sprintf(HistoTitle, "%s XZ View;Plane;Cell", TitleHelper);
134  fHitCellsXView = tfs->make<TH2F>(HistoName, HistoTitle,
135  fPlanes, 0., (float)fPlanes,
136  fCellsX, 0., (float)fCellsX);
137 
138  sprintf(HistoName, "fHitCellsYView");
139  sprintf(HistoTitle, "%s YZ View;Plane;Cell", TitleHelper);
140  fHitCellsYView = tfs->make<TH2F>(HistoName, HistoTitle,
141  fPlanes, 0., (float)fPlanes,
142  fCellsY, 0., (float)fCellsY);
143  }
144 
145  //......................................................................
147  {
148  // Get rawdata::RawDigits from the event
149  bool isRawDigitListEmpty = false;
151  try {
152  evt.getByLabel(fRawDataLabel, diglist);
153  if(diglist->empty()) { isRawDigitListEmpty = true; }
154  }
155  catch(...) { isRawDigitListEmpty = true; }
156 
157  if(isRawDigitListEmpty) {
158  if(msg1cnt < 5) {
159  mf::LogWarning("DetAna RawDigit") << "Error retrieving raw digits." << std::endl;
160  ++msg1cnt;
161  }
162  fNumDigits->Fill(0.);
163  return;
164  }
165 
166  std::vector<art::Ptr<rawdata::RawDigit> > rdptrs;
167  art::fill_ptr_vector(rdptrs, diglist);
168 
169  fNumDigits->Fill(rdptrs.size()); // Fill number of digits in the event
170 
174 
176  const unsigned int nSamples = conv.getNSamples(rdptrs[0]->Version());
177  const unsigned int nPretrig = conv.getNPretriggeredSamples(rdptrs[0]->Version());
178 
179  // Loop over RawDigits
180  // rd should be an art::Ptr<rawdata::RawDigit>
181  for(const auto& rd : rdptrs) {
182  fADC->Fill(rd->ADC());
183 
184  // Fill the ADC value of each sample
185  for(unsigned int i_adc = 0, n_adc = rd->NADC(); i_adc < n_adc; ++i_adc) {
186  fADCPerSample->Fill(rd->ADC(i_adc));
187  }
188 
189  if(rd->NADC() == nSamples) {
190  if(nSamples >= 4 && nPretrig >= 3) { // Make sure there are enough samples
191  for(unsigned int i = 0, n = 3; i < n; ++i) {
192  // Fill the difference between the ADC(1/2/3) and ADC(0)
193  fADCDiff[i]->Fill(rd->ADC(i + 1) - rd->ADC(0));
194  }
195  }
196  else {
197  if(msg1cnt < 5) {
198  mf::LogWarning("Multipoint") << "Not enough samples to fill ADC difference plots." << std::endl;
199  ++msg1cnt;
200  }
201  }
202  }
203  else {
204  mf::LogWarning("DetAna") << "The number of ADC samples in the RawDigit did not match"
205  << " the number of samples in the DAQDataFormat." << std::endl;
206  }
207 
208  // Ask channel map for plane and cell number
209  view = geo->Plane(cmap->GetPlane(rd.get()))->View();
210  if(view == geo::kX) {
211  fHitCellsXView->Fill(1.*cmap->GetPlane(rd.get())+0.5, 1.*cmap->GetCell(rd.get())+0.5);
212  }
213  if(view == geo::kY) {
214  fHitCellsYView->Fill(1.*cmap->GetPlane(rd.get())+0.5, 1.*cmap->GetCell(rd.get())+0.5);
215  }
216  } // end of loop over RawDigits
217  }
218 } // end of namespace mcchk
219 
220 ///////////////////////////////////////////////////////////////////////////////
221 namespace mcchk {
223 }
std::string fRawDataLabel
module label for the raw data
double fRawDigit_FDOffset
How far to shift FD Number of digit bins.
Simple module to analyze MC cosmics distributions.
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
Vertical planes which measure X.
Definition: PlaneGeo.h:28
const PlaneGeo * Plane(unsigned int i) const
DEFINE_ART_MODULE(TestTMapFile)
DetAna(fhicl::ParameterSet const &pset)
Horizontal planes which measure Y.
Definition: PlaneGeo.h:29
Definition: Run.h:31
TH1F * fADC
ADC value for each digit.
View_t View() const
Which coordinate does this plane measure.
Definition: PlaneGeo.h:53
void analyze(art::Event const &evt)
uint32_t getNSamples(const version_t ver) const
Get the number of samples.
uint32_t getNPretriggeredSamples(const version_t ver) const
Get number of pretriggered samples.
TH2F * fHitCellsYView
cells hit in the y view
novadaq::cnv::DetId DetId() const
Prefer ds::DetectorService::DetId() instead.
Definition: GeometryBase.h:243
int evt
unsigned short GetPlane(const rawdata::RawDigit *dig)
Definition: CMap.cxx:285
EDAnalyzer(Table< Config > const &config)
Definition: EDAnalyzer.h:100
virtual ~DetAna()
Definition: run.py:1
TH1F * fNumDigits
number of digits per event
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
TH1F * fADCDiff[3]
Difference between ADC(1/2/3) and ADC(0)
T * make(ARGS...args) const
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
TH2F * fHitCellsXView
cells hit in the x view
A module to check the results from the Monte Carlo generator.
const std::set< unsigned int > & GetPlanesByView(View_t v=kXorY) const
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
unsigned short GetCell(const rawdata::RawDigit *dig)
Definition: CMap.cxx:327
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:464
void beginRun(art::Run const &run)
static int msg1cnt
Helper for AttenCurve.
Definition: Path.h:10
TH1F * fADCPerSample
ADC value for each sample.
Encapsulate the geometry of one entire detector (near, far, ndos)