artdaq_demo  v3_04_00
WFViewer_module.cc
1 #include "art/Framework/Core/EDAnalyzer.h"
2 #include "art/Framework/Principal/Handle.h"
3 #include "art/Framework/Principal/Event.h"
4 #include "art/Framework/Principal/Run.h"
5 #include "art/Framework/Core/ModuleMacros.h"
6 #include "canvas/Utilities/InputTag.h"
7 
8 #include "artdaq-core/Data/Fragment.hh"
9 #include "artdaq-core/Data/ContainerFragment.hh"
10 
11 #include "artdaq-core-demo/Overlays/FragmentType.hh"
12 #include "artdaq-core-demo/Overlays/ToyFragment.hh"
13 
14 #include "cetlib_except/exception.h"
15 
16 #include <TFile.h>
17 #include <TRootCanvas.h>
18 #include <TCanvas.h>
19 #include <TGraph.h>
20 #include <TAxis.h>
21 #include <TH1D.h>
22 #include <TStyle.h>
23 
24 #include <numeric>
25 #include <vector>
26 #include <functional>
27 #include <algorithm>
28 #include <iostream>
29 #include <sstream>
30 #include <initializer_list>
31 #include <limits>
32 
33 using std::cout;
34 using std::cerr;
35 using std::endl;
36 
37 namespace demo
38 {
42  class WFViewer : public art::EDAnalyzer
43  {
44  public:
61  explicit WFViewer(fhicl::ParameterSet const& p);
62 
66  virtual ~WFViewer() = default;
67 
72  void analyze(art::Event const& e) override;
73 
78  void beginRun(art::Run const&) override;
79 
80  private:
81 
82  std::unique_ptr<TCanvas> canvas_[2];
83  std::vector<Double_t> x_;
84  int prescale_;
85  bool digital_sum_only_;
86  art::RunNumber_t current_run_;
87 
88  std::size_t num_x_plots_;
89  std::size_t num_y_plots_;
90 
91  std::string raw_data_label_;
92  std::vector<artdaq::Fragment::fragment_id_t> fragment_ids_;
93 
94  std::vector<std::unique_ptr<TGraph>> graphs_;
95  std::vector<std::unique_ptr<TH1D>> histograms_;
96 
97  std::map<artdaq::Fragment::fragment_id_t, std::size_t> id_to_index_;
98  std::string outputFileName_;
99  TFile* fFile_;
100  bool writeOutput_;
101  };
102 }
103 
104 demo::WFViewer::WFViewer(fhicl::ParameterSet const& ps):
105  art::EDAnalyzer(ps)
106  , prescale_(ps.get<int>("prescale"))
107  , digital_sum_only_(ps.get<bool>("digital_sum_only", false))
108  , current_run_(0)
109  , num_x_plots_(ps.get<std::size_t>("num_x_plots", std::numeric_limits<std::size_t>::max()))
110  , num_y_plots_(ps.get<std::size_t>("num_y_plots", std::numeric_limits<std::size_t>::max()))
111  , raw_data_label_(ps.get<std::string>("raw_data_label", "daq"))
112  , fragment_ids_(ps.get<std::vector<artdaq::Fragment::fragment_id_t>>("fragment_ids"))
113  , graphs_(fragment_ids_.size())
114  , histograms_(fragment_ids_.size())
115  , outputFileName_(ps.get<std::string>("fileName", "artdaqdemo_onmon.root"))
116  , writeOutput_(ps.get<bool>("write_to_file", false))
117 {
118  if (num_x_plots_ == std::numeric_limits<std::size_t>::max() ||
119  num_y_plots_ == std::numeric_limits<std::size_t>::max())
120  {
121  switch (fragment_ids_.size())
122  {
123  case 1: num_x_plots_ = num_y_plots_ = 1;
124  break;
125  case 2: num_x_plots_ = 2;
126  num_y_plots_ = 1;
127  break;
128  case 3:
129  case 4: num_x_plots_ = 2;
130  num_y_plots_ = 2;
131  break;
132  case 5:
133  case 6: num_x_plots_ = 3;
134  num_y_plots_ = 2;
135  break;
136  case 7:
137  case 8: num_x_plots_ = 4;
138  num_y_plots_ = 2;
139  break;
140  default:
141  num_x_plots_ = num_y_plots_ = static_cast<std::size_t>(ceil(sqrt(fragment_ids_.size())));
142  }
143  }
144 
145  // id_to_index_ will translate between a fragment's ID and where in
146  // the vector of graphs and histograms it's located
147 
148  for (std::size_t i_f = 0; i_f < fragment_ids_.size(); ++i_f)
149  {
150  id_to_index_[fragment_ids_[i_f]] = i_f;
151  }
152 
153  gStyle->SetOptStat("irm");
154  gStyle->SetMarkerStyle(22);
155  gStyle->SetMarkerColor(4);
156 }
157 
158 void demo::WFViewer::analyze(art::Event const& e)
159 {
160  static std::size_t evt_cntr = -1;
161  evt_cntr++;
162 
163  // John F., 1/22/14 -- there's probably a more elegant way of
164  // collecting fragments of various types using ART interface code;
165  // will investigate. Right now, we're actually re-creating the
166  // fragments locally
167 
168  artdaq::Fragments fragments;
169  artdaq::FragmentPtrs containerFragments;
170  std::vector<std::string> fragment_type_labels{"TOY1","TOY2","ContainerTOY1", "ContainerTOY2"};
171 
172  for (auto label: fragment_type_labels)
173  {
174  art::Handle<artdaq::Fragments> fragments_with_label;
175  e.getByLabel("daq", label, fragments_with_label);
176 
177  if (!fragments_with_label.isValid()) continue;
178  // for (int i_l = 0; i_l < static_cast<int>(fragments_with_label->size()); ++i_l) {
179  // fragments.emplace_back( (*fragments_with_label)[i_l] );
180  // }
181 
182  if (label == "Container" || label == "ContainerTOY1" || label == "ContainerTOY2")
183  {
184  for (auto cont : *fragments_with_label)
185  {
186  artdaq::ContainerFragment contf(cont);
187  for (size_t ii = 0; ii < contf.block_count(); ++ii)
188  {
189  containerFragments.push_back(contf[ii]);
190  fragments.push_back(*containerFragments.back());
191  }
192  }
193  }
194  else
195  {
196  for (auto frag : *fragments_with_label)
197  {
198  fragments.emplace_back(frag);
199  }
200  }
201  }
202 
203  // John F., 1/5/14
204 
205  // Here, we loop over the fragments passed to the analyze
206  // function. A warning is flashed if either (A) the fragments aren't
207  // all from the same event, or (B) there's an unexpected number of
208  // fragments given the number of boardreaders and the number of
209  // fragments per board
210 
211  // For every Nth event, where N is the "prescale" setting, plot the
212  // distribution of ADC counts from each board_id / fragment_id
213  // combo. Also, if "digital_sum_only" is set to false in the FHiCL
214  // string, then plot, for the Nth event, a graph of the ADC values
215  // across all channels in each board_id / fragment_id combo
216 
217  artdaq::Fragment::sequence_id_t expected_sequence_id = std::numeric_limits<artdaq::Fragment::sequence_id_t>::max();
218 
219  // for (std::size_t i = 0; i < fragments.size(); ++i) {
220  for (const auto& frag : fragments)
221  {
222  // Pointers to the types of fragment overlays WFViewer can handle;
223  // only one will be used per fragment, of course
224 
225  std::unique_ptr<ToyFragment> toyPtr;
226 
227  // const auto& frag( fragments[i] ); // Basically a shorthand
228 
229  // if (i == 0)
230  if (expected_sequence_id == std::numeric_limits<artdaq::Fragment::sequence_id_t>::max())
231  {
232  expected_sequence_id = frag.sequenceID();
233  }
234 
235  if (expected_sequence_id != frag.sequenceID())
236  {
237  cerr << "Warning in WFViewer: expected fragment with sequence ID " << expected_sequence_id << ", received one with sequence ID " << frag.sequenceID() << endl;
238  }
239 
240  FragmentType fragtype = static_cast<FragmentType>(frag.type());
241  std::size_t max_adc_count = std::numeric_limits<std::size_t>::max();
242  std::size_t total_adc_values = std::numeric_limits<std::size_t>::max();
243 
244  switch (fragtype)
245  {
246  case FragmentType::TOY1:
247  toyPtr.reset(new ToyFragment(frag));
248  total_adc_values = toyPtr->total_adc_values();
249  max_adc_count = pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1;
250  break;
251  case FragmentType::TOY2:
252  toyPtr.reset(new ToyFragment(frag));
253  total_adc_values = toyPtr->total_adc_values();
254  max_adc_count = pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1;
255  break;
256  default:
257  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
258  }
259 
260  artdaq::Fragment::fragment_id_t fragment_id = frag.fragmentID();
261  if (!id_to_index_.count(fragment_id))
262  {
263  cerr << "Warning in WFViewer: unexpected Fragment with fragment_id " << std::to_string(fragment_id) << " encountered!";
264  continue;
265  }
266  std::size_t ind = id_to_index_[fragment_id];
267 
268 
269  // If a histogram doesn't exist for this board_id / fragment_id combo, create it
270 
271  if (!histograms_[ind])
272  {
273  histograms_[ind] = std::unique_ptr<TH1D>(new TH1D(Form("Fragment_%d_hist", fragment_id), "", max_adc_count + 1, -0.5, max_adc_count + 0.5));
274 
275  histograms_[ind]->SetTitle(Form("Frag %d, Type %s", fragment_id,
276  fragmentTypeToString(fragtype).c_str()));
277  histograms_[ind]->GetXaxis()->SetTitle("ADC value");
278  }
279 
280  // For every event, fill the histogram (prescale is ignored here)
281 
282  // Is there some way to templatize an ART module? If not, we're
283  // stuck with this switch code...
284 
285  switch (fragtype)
286  {
287  case FragmentType::TOY1:
288  case FragmentType::TOY2:
289  for (auto val = toyPtr->dataBeginADCs(); val != toyPtr->dataEndADCs(); ++val)
290  histograms_[ind]->Fill(*val);
291  break;
292 
293  default:
294  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
295  }
296 
297  if (evt_cntr % prescale_ - 1 && prescale_ > 1)
298  {
299  continue;
300  }
301 
302  // If we pass the prescale, then if we're not going with
303  // digital_sum_only, plot the ADC counts for this particular event/board/fragment_id
304 
305  if (!digital_sum_only_)
306  {
307  // Create the graph's x-axis
308 
309  if (x_.size() != total_adc_values)
310  {
311  x_.resize(total_adc_values);
312 
313  std::iota(x_.begin(), x_.end(), 0);
314  }
315 
316  // If the graph doesn't exist, create it. Not sure whether to
317  // make it an error if the total_adc_values is new
318 
319  if (!graphs_[ind] || static_cast<std::size_t>(graphs_[ind]->GetN()) != total_adc_values)
320  {
321  graphs_[ind] = std::unique_ptr<TGraph>(new TGraph(total_adc_values));
322  graphs_[ind]->SetName(Form("Fragment_%d_graph", fragment_id));
323  graphs_[ind]->SetLineColor(4);
324  std::copy(x_.begin(), x_.end(), graphs_[ind]->GetX());
325  }
326 
327  // Get the data from the fragment
328 
329  // Is there some way to templatize an ART module? If not, we're stuck with this awkward switch code...
330 
331  switch (fragtype)
332  {
333  case FragmentType::TOY1:
334  case FragmentType::TOY2:
335  {
336  std::copy(toyPtr->dataBeginADCs(), toyPtr->dataBeginADCs() + total_adc_values, graphs_[ind]->GetY());
337  }
338  break;
339 
340  default:
341  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
342  }
343 
344 
345  // And now prepare the graphics without actually drawing anything yet
346 
347  canvas_[1]->cd(ind + 1);
348  TVirtualPad* pad = static_cast<TVirtualPad*>(canvas_[1]->GetPad(ind + 1));
349 
350  Double_t lo_x, hi_x, lo_y, hi_y, dummy;
351 
352  graphs_[ind]->GetPoint(0, lo_x, dummy);
353  graphs_[ind]->GetPoint(graphs_[ind]->GetN() - 1, hi_x, dummy);
354 
355  lo_x -= 0.5;
356  hi_x += 0.5;
357 
358  lo_y = -0.5;
359  hi_y = max_adc_count + 0.5;
360 
361 
362  TH1F* padframe = static_cast<TH1F*>(pad->DrawFrame(lo_x, lo_y, hi_x, hi_y));
363  padframe->SetTitle(Form("Frag %d, Type %s, SeqID %d", static_cast<int>(fragment_id),
364  fragmentTypeToString(fragtype).c_str(),
365  static_cast<int>(expected_sequence_id)));
366  padframe->GetXaxis()->SetTitle("ADC #");
367  pad->SetGrid();
368  padframe->Draw("SAME");
369  }
370 
371  // Draw the histogram
372 
373  canvas_[0]->cd(ind + 1);
374  histograms_[ind]->Draw();
375 
376  canvas_[0]->Modified();
377  canvas_[0]->Update();
378 
379  // And, if desired, the Nth event's ADC counts
380 
381  if (!digital_sum_only_)
382  {
383  canvas_[1]->cd(ind + 1);
384 
385  graphs_[ind]->Draw("PSAME");
386 
387  canvas_[1]->Modified();
388  canvas_[1]->Update();
389  }
390 
391  if (writeOutput_)
392  {
393  canvas_[0]->Write("wf0", TObject::kOverwrite);
394  canvas_[1]->Write("wf1", TObject::kOverwrite);
395  fFile_->Write();
396  }
397  }
398 }
399 
400 void demo::WFViewer::beginRun(art::Run const& e)
401 {
402  if (e.run() == current_run_) return;
403  current_run_ = e.run();
404 
405  if (writeOutput_)
406  {
407  fFile_ = new TFile(outputFileName_.c_str(), "RECREATE");
408  fFile_->cd();
409  }
410 
411  for (int i = 0; i < 2; i++) canvas_[i] = 0;
412  for (auto& x: graphs_) x = 0;
413  for (auto& x: histograms_) x = 0;
414 
415  for (int i = 0; (i < 2 && !digital_sum_only_) || i < 1; i++)
416  {
417  canvas_[i] = std::unique_ptr<TCanvas>(new TCanvas(Form("wf%d", i)));
418  canvas_[i]->Divide(num_x_plots_, num_y_plots_);
419  canvas_[i]->Update();
420  ((TRootCanvas*)canvas_[i]->GetCanvasImp())->DontCallClose();
421  }
422 
423  canvas_[0]->SetTitle("ADC Value Distribution");
424 
425  if (! digital_sum_only_)
426  {
427  canvas_[1]->SetTitle("ADC Values, Event Snapshot");
428  }
429 
430  if (writeOutput_)
431  {
432  canvas_[0]->Write();
433  canvas_[1]->Write();
434  }
435 }
436 
437 
438 DEFINE_ART_MODULE(demo::WFViewer)
virtual ~WFViewer()=default
WFViewer default Destructor.
void analyze(art::Event const &e) override
Analyze an event. Called by art for each event in run (based on command line options) ...
void beginRun(art::Run const &) override
Art calls this function at the beginning of the run. Used for set-up of ROOT histogram objects and to...
An example art analysis module which plots events both as histograms and event snapshots (plot of ADC...
WFViewer(fhicl::ParameterSet const &p)
WFViewer Constructor.