artdaq_demo  3.13.00
DemoViewer_module.cc
1 // P.Murat: a 0-th order example of a DQM client using ROOT web-based graphics
3 // cloned from WFViewer_module.cc w/o much thought
4 // the only purpose is to demonstrate the use of the web-based GUI
5 // - creates two canvases with the following URLs:
6 // http://127.0.0.1:8877/win1/
7 // http://127.0.0.1:8877/win2/
8 // key points:
9 // - create a TApplication running in batch more and using a certain URL
10 // - call gSystem->ProcessEvents() once per event
11 //
12 // for illustration only, save histograms once in the end of the run
13 // (an online application has to do it periodically during the run)
15 #include "TRACE/tracemf.h"
16 #define TRACE_NAME "DemoViewer"
17 
18 #include "art/Framework/Core/EDAnalyzer.h"
19 #include "art/Framework/Core/ModuleMacros.h"
20 #include "art/Framework/Principal/Event.h"
21 #include "art/Framework/Principal/Handle.h"
22 #include "art/Framework/Principal/Run.h"
23 #include "canvas/Utilities/InputTag.h"
24 #include "cetlib_except/exception.h"
25 
26 #include "artdaq-core/Data/ContainerFragment.hh"
27 #include "artdaq-core/Data/Fragment.hh"
28 
29 #include "artdaq-core-demo/Overlays/FragmentType.hh"
30 #include "artdaq-core-demo/Overlays/ToyFragment.hh"
31 
32 #include <TApplication.h>
33 #include <TSystem.h>
34 #include <TAxis.h>
35 #include <TCanvas.h>
36 #include <TFile.h>
37 #include <TGraph.h>
38 #include <TH1D.h>
39 #include <TRootCanvas.h>
40 #include <TStyle.h>
41 
42 #include <algorithm>
43 #include <functional>
44 #include <initializer_list>
45 #include <iostream>
46 #include <limits>
47 #include <memory>
48 #include <numeric>
49 #include <sstream>
50 #include <unordered_map>
51 #include <vector>
52 
53 #include "tracemf.h"
54 #define TRACE_NAME "DemoViewer"
55 
56 namespace demo {
60 class DemoViewer : public art::EDAnalyzer {
61 public:
77  explicit DemoViewer(fhicl::ParameterSet const& p);
78 
79  ~DemoViewer() override;
80 
81  void analyze(art::Event const& e) override;
82 
83  void beginJob() override;
84  void beginRun(art::Run const& e) override;
85  void endRun (art::Run const& e) override;
86 
87 private:
88  DemoViewer(DemoViewer const&) = delete;
89  DemoViewer(DemoViewer&&) = delete;
90  DemoViewer& operator=(DemoViewer const&) = delete;
91  DemoViewer& operator=(DemoViewer&&) = delete;
92 
93  TCanvas* _hCanvas;
94  TCanvas* _gCanvas;
95  std::vector<Double_t> x_;
96  int prescale_;
97  art::RunNumber_t current_run_;
98 
99  size_t max_num_x_plots_;
100  size_t max_num_y_plots_;
101  std::size_t num_x_plots_;
102  std::size_t num_y_plots_;
103 
104  std::string raw_data_label_;
105 
106  std::unordered_map<artdaq::Fragment::fragment_id_t, TGraph*> graphs_;
107  std::unordered_map<artdaq::Fragment::fragment_id_t, TH1D*> histograms_;
108 
109  std::map<artdaq::Fragment::fragment_id_t, std::size_t> id_to_index_;
110 
111  std::string outputFileName_;
112  TFile* fFile_;
113  bool writeOutput_;
114  bool newCanvas_;
115  bool dynamicMode_;
116 
117  TApplication* _app;
118 
119 
120  void getXYDims_ ();
121  void bookCanvas_();
122 };
123 
124  //-----------------------------------------------------------------------------
125 DemoViewer::DemoViewer(fhicl::ParameterSet const& ps)
126  : art::EDAnalyzer (ps)
127  , prescale_ (ps.get<int> ("prescale"))
128  , current_run_ (0)
129  , max_num_x_plots_ (ps.get<std::size_t>("num_x_plots", std::numeric_limits<std::size_t>::max()))
130  , max_num_y_plots_ (ps.get<std::size_t>("num_y_plots", std::numeric_limits<std::size_t>::max()))
131  , num_x_plots_ (0)
132  , num_y_plots_ (0)
133  , raw_data_label_ (ps.get<std::string>("raw_data_label", "daq"))
134  , graphs_ ()
135  , histograms_ ()
136  , outputFileName_ (ps.get<std::string>("fileName", "artdaqdemo_onmon.root"))
137  , writeOutput_ (ps.get<bool> ("write_to_file", false))
138  , newCanvas_ (true)
139  , dynamicMode_ (ps.get<bool> ("dynamic_mode", true))
140 {
141  gStyle->SetOptStat("irm");
142  gStyle->SetMarkerStyle(22);
143  gStyle->SetMarkerColor(4);
144 
145  if (ps.has_key("fragment_ids")) {
146  auto fragment_ids = ps.get<std::vector<artdaq::Fragment::fragment_id_t>>("fragment_ids");
147  for (auto& id : fragment_ids) {
148  auto index = id_to_index_.size();
149  id_to_index_[id] = index;
150  }
151  }
152 }
153  //-----------------------------------------------------------------------------
154 void DemoViewer::beginJob() {
155 
156  int tmp_argc(2);
157  char** tmp_argv;
158 
159  tmp_argv = new char*[2];
160  tmp_argv[0] = new char[100];
161  tmp_argv[1] = new char[100];
162 
163  strcpy(tmp_argv[0],"-b");
164  strcpy(tmp_argv[1],"--web=server:8877");
165 
166  _app = new TApplication("DemoViewer", &tmp_argc, tmp_argv);
167 
168  // app->Run()
169  // app_->Run(true);
170  delete [] tmp_argv;
171 }
172 
173 //-----------------------------------------------------------------------------
174 void DemoViewer::getXYDims_() {
175  // Enforce positive maxes
176  if (max_num_x_plots_ == 0) max_num_x_plots_ = std::numeric_limits<size_t>::max();
177  if (max_num_y_plots_ == 0) max_num_y_plots_ = std::numeric_limits<size_t>::max();
178 
179  num_x_plots_ = num_y_plots_ = static_cast<std::size_t>(ceil(sqrt(id_to_index_.size())));
180 
181  // Do trivial check first to avoid multipling max * max -> undefined
182  if (id_to_index_.size() > max_num_x_plots_ && id_to_index_.size() > max_num_x_plots_ * max_num_y_plots_) {
183  num_x_plots_ = max_num_x_plots_;
184  num_y_plots_ = max_num_y_plots_;
185  auto max = num_x_plots_ * num_y_plots_;
186  auto it = id_to_index_.begin();
187  while (it != id_to_index_.end()) {
188  if (it->second >= max) it = id_to_index_.erase(it);
189  else ++it;
190  }
191  }
192 
193  // Some predefined "nice looking" plotscapes...
194 
195  if (max_num_x_plots_ >= 4 && max_num_y_plots_ >= 2)
196  {
197  switch (id_to_index_.size())
198  {
199  case 1:
200  num_x_plots_ = num_y_plots_ = 1;
201  break;
202  case 2:
203  num_x_plots_ = 2;
204  num_y_plots_ = 1;
205  break;
206  case 3:
207  case 4:
208  num_x_plots_ = 2;
209  num_y_plots_ = 2;
210  break;
211  case 5:
212  case 6:
213  num_x_plots_ = 3;
214  num_y_plots_ = 2;
215  break;
216  case 7:
217  case 8:
218  num_x_plots_ = 4;
219  num_y_plots_ = 2;
220  break;
221  default:
222  break;
223  }
224  }
225  else
226  {
227  // Make sure we fit within specifications
228  while (num_x_plots_ > max_num_x_plots_)
229  {
230  num_x_plots_--;
231  num_y_plots_ = static_cast<size_t>(ceil(id_to_index_.size() / num_x_plots_));
232  }
233  }
234  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_;
235 }
236 
237 //-----------------------------------------------------------------------------
238 //
239 //-----------------------------------------------------------------------------
240 void DemoViewer::bookCanvas_() {
241  newCanvas_ = false;
242  getXYDims_();
243 
244  _hCanvas = new TCanvas("wf0");
245  _hCanvas->Divide(num_x_plots_, num_y_plots_);
246  _hCanvas->SetTitle("ADC Value Distribution");
247  _hCanvas->Update();
248 
249  _gCanvas = new TCanvas("wf1");
250  _gCanvas->Divide(num_x_plots_, num_y_plots_);
251  _gCanvas->SetTitle("ADC Values, Event Snapshot");
252  _gCanvas->Update();
253 }
254 
255 //-----------------------------------------------------------------------------
256 //
257 //-----------------------------------------------------------------------------
258 DemoViewer::~DemoViewer() {
259  // We're going to let ROOT's own garbage collection deal with histograms and Canvases...
260  for (auto& histogram : histograms_) {
261  histogram.second = nullptr;
262  }
263  histograms_.clear();
264  for (auto& graph : graphs_) {
265  graph.second = nullptr;
266  }
267  graphs_.clear();
268 
269  _hCanvas = nullptr;
270  _gCanvas = nullptr;
271  fFile_ = nullptr;
272 }
273 
274 
275 //-----------------------------------------------------------------------------
276 //
277 //-----------------------------------------------------------------------------
278 void DemoViewer::analyze(art::Event const& e) {
279  static std::size_t evt_cntr = -1;
280  evt_cntr++;
281 
282  // John F., 1/22/14 -- there's probably a more elegant way of
283  // collecting fragments of various types using ART interface code;
284  // will investigate. Right now, we're actually re-creating the
285  // fragments locally
286 
287  artdaq::Fragments fragments;
288  artdaq::FragmentPtrs containerFragments;
289 
290  std::vector<art::Handle<artdaq::Fragments>> fragmentHandles;
291  fragmentHandles = e.getMany<std::vector<artdaq::Fragment>>();
292 
293  for (const auto& handle : fragmentHandles) {
294  if (!handle.isValid() || handle->empty()) {
295  continue;
296  }
297 
298  if (handle->front().type() == artdaq::Fragment::ContainerFragmentType) {
299  for (const auto& cont : *handle) {
300  artdaq::ContainerFragment contf(cont);
301  auto ftype = contf.fragment_type();
302  if (ftype != FragmentType::TOY1 && ftype != FragmentType::TOY2) break;
303 
304  for (size_t ii = 0; ii < contf.block_count(); ++ii) {
305  containerFragments.push_back(contf[ii]);
306  fragments.push_back(*containerFragments.back());
307  if (newCanvas_ && !id_to_index_.count(fragments.back().fragmentID())) {
308  auto index = id_to_index_.size();
309  id_to_index_[fragments.back().fragmentID()] = index;
310  }
311  }
312  }
313  }
314  else {
315  for (auto frag : *handle) {
316  fragments.emplace_back(frag);
317  if (newCanvas_ && !id_to_index_.count(fragments.back().fragmentID())) {
318  auto index = id_to_index_.size();
319  id_to_index_[fragments.back().fragmentID()] = index;
320  }
321  }
322  }
323  }
324 
325  if (newCanvas_) {
326  bookCanvas_();
327  }
328 
329  // John F., 1/5/14
330 
331  // Here, we loop over the fragments passed to the analyze
332  // function. A warning is flashed if either (A) the fragments aren't
333  // all from the same event, or (B) there's an unexpected number of
334  // fragments given the number of boardreaders and the number of
335  // fragments per board
336 
337  // For every Nth event, where N is the "prescale" setting, plot the
338  // distribution of ADC counts from each board_id / fragment_id
339  // combo.
340  // also , plot, for the Nth event, a graph of the ADC values
341  // across all channels in each board_id / fragment_id combo
342 
343  artdaq::Fragment::sequence_id_t expected_sequence_id = std::numeric_limits<artdaq::Fragment::sequence_id_t>::max();
344 
345  for (const auto& frag : fragments) {
346  std::unique_ptr<ToyFragment> toyPtr;
347 
348  if (expected_sequence_id == std::numeric_limits<artdaq::Fragment::sequence_id_t>::max()) {
349  expected_sequence_id = frag.sequenceID();
350  }
351 
352  if (expected_sequence_id != frag.sequenceID()) {
353  TLOG(TLVL_WARNING) << "Warning in DemoViewer: expected fragment with sequence ID " << expected_sequence_id
354  << ", received one with sequence ID " << frag.sequenceID();
355  }
356 
357  auto fragtype = static_cast<FragmentType>(frag.type());
358  std::size_t max_adc_count = std::numeric_limits<std::size_t>::max();
359  std::size_t total_adc_values = std::numeric_limits<std::size_t>::max();
360 
361  switch (fragtype) {
362  case FragmentType::TOY1:
363  toyPtr = std::make_unique<ToyFragment>(frag);
364  total_adc_values = toyPtr->total_adc_values();
365  max_adc_count = static_cast<size_t>(pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1);
366  break;
367  case FragmentType::TOY2:
368  toyPtr = std::make_unique<ToyFragment>(frag);
369  total_adc_values = toyPtr->total_adc_values();
370  max_adc_count = static_cast<size_t>(pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1);
371  break;
372  default:
373  throw cet::exception("Error in DemoViewer: unknown fragment type supplied");
374  }
375 
376  artdaq::Fragment::fragment_id_t fid = frag.fragmentID();
377  if (id_to_index_.count(fid) == 0u) {
378  TLOG(TLVL_WARNING) << "Warning in DemoViewer: unexpected Fragment with fragment id " << std::to_string(fid)
379  << " encountered!";
380  continue;
381  }
382 
383  std::size_t ind = id_to_index_[fid];
384 
385  // If a histogram doesn't exist for this board_id / fragment_id combo, create it
386 
387  if (histograms_.count(fid) == 0 || histograms_[fid] == nullptr) {
388  histograms_[fid] = new TH1D(Form("Fragment_%d_hist",fid),"",max_adc_count+1,-0.5,max_adc_count+0.5);
389  histograms_[fid]->SetTitle(Form("Frag %d, Type %s",fid, fragmentTypeToString(fragtype).c_str()));
390  histograms_[fid]->GetXaxis()->SetTitle("ADC value");
391 
392  _hCanvas->cd(ind + 1);
393  histograms_[fid]->Draw();
394  }
395 
396  switch (fragtype) {
397  case FragmentType::TOY1:
398  case FragmentType::TOY2:
399  for (auto val = toyPtr->dataBeginADCs(); val != toyPtr->dataEndADCs(); ++val) {
400  histograms_[fid]->Fill(*val);
401  }
402  break;
403 
404  default:
405  TLOG(TLVL_ERROR) << "Error in DemoViewer: unknown fragment type supplied";
406  throw cet::exception("Error in DemoViewer: unknown fragment type supplied");
407  }
408 
409  if (((evt_cntr % prescale_ - 1) != 0u) && prescale_ > 1) continue;
410 
411  // Create the graph's x-axis
412 
413  if (x_.size() != total_adc_values) {
414  x_.resize(total_adc_values);
415  std::iota(x_.begin(), x_.end(), 0);
416  }
417 
418  // If the graph doesn't exist, create it. Not sure whether to
419  // make it an error if the total_adc_values is new
420 
421  if (graphs_.count(fid) == 0 || graphs_[fid] == nullptr ||
422  static_cast<std::size_t>(graphs_[fid]->GetN()) != total_adc_values) {
423 
424  graphs_[fid] = new TGraph(total_adc_values);
425  graphs_[fid]->SetName(Form("Fragment_%d_graph", fid));
426  graphs_[fid]->SetLineColor(4);
427  std::copy(x_.begin(), x_.end(), graphs_[fid]->GetX());
428 
429  _gCanvas->cd(ind + 1);
430  graphs_[fid]->Draw("ALP");
431  }
432 
433  // Get the data from the fragment
434  switch (fragtype) {
435  case FragmentType::TOY1:
436  case FragmentType::TOY2: {
437  std::copy(toyPtr->dataBeginADCs(), toyPtr->dataBeginADCs() + total_adc_values, graphs_[fid]->GetY());
438  }
439  break;
440 
441  default:
442  TLOG(TLVL_ERROR) << "Error in DemoViewer: unknown fragment type supplied";
443  throw cet::exception("Error in DemoViewer: unknown fragment type supplied"); // NOLINT(cert-err60-cpp)
444  }
445  }
446 
447  gSystem->ProcessEvents();
448 }
449 
450 //-----------------------------------------------------------------------------
451 void DemoViewer::beginRun(art::Run const& e) {
452  if (e.run() == current_run_) return;
453  current_run_ = e.run();
454 
455  _hCanvas = nullptr;
456  _gCanvas = nullptr;
457  for (auto& x : graphs_ ) x.second = nullptr;
458  for (auto& x : histograms_) x.second = nullptr;
459 
460  newCanvas_ = true;
461  if (!dynamicMode_) bookCanvas_();
462 }
463 
464 //-----------------------------------------------------------------------------
465 void DemoViewer::endRun(art::Run const& e) {
466  if (e.run() == current_run_) return;
467  current_run_ = e.run();
468 
469  if (writeOutput_) {
470  fFile_ = new TFile(outputFileName_.c_str(), "RECREATE");
471 
472  _hCanvas->Write("wf0", TObject::kOverwrite);
473  if (_gCanvas != nullptr) {
474  _gCanvas->Write("wf1", TObject::kOverwrite);
475  }
476  fFile_->Write();
477  }
478 }
479 
480 DEFINE_ART_MODULE(DemoViewer) // NOLINT(performance-unnecessary-value-param)
481 } // namespace demo
An example art analysis module which plots events both as histograms and event snapshots (plot of ADC...
DemoViewer(fhicl::ParameterSet const &p)
DemoViewer Constructor.