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