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