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