artdaq_demo  v3_00_01
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/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  std::size_t ind = id_to_index_[fragment_id];
262 
263 
264  // If a histogram doesn't exist for this board_id / fragment_id combo, create it
265 
266  if (!histograms_[ind])
267  {
268  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));
269 
270  histograms_[ind]->SetTitle(Form("Frag %d, Type %s", fragment_id,
271  fragmentTypeToString(fragtype).c_str()));
272  histograms_[ind]->GetXaxis()->SetTitle("ADC value");
273  }
274 
275  // For every event, fill the histogram (prescale is ignored here)
276 
277  // Is there some way to templatize an ART module? If not, we're
278  // stuck with this switch code...
279 
280  switch (fragtype)
281  {
282  case FragmentType::TOY1:
283  case FragmentType::TOY2:
284  for (auto val = toyPtr->dataBeginADCs(); val != toyPtr->dataEndADCs(); ++val)
285  histograms_[ind]->Fill(*val);
286  break;
287 
288  default:
289  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
290  }
291 
292  if (evt_cntr % prescale_ - 1 && prescale_ > 1)
293  {
294  continue;
295  }
296 
297  // If we pass the prescale, then if we're not going with
298  // digital_sum_only, plot the ADC counts for this particular event/board/fragment_id
299 
300  if (!digital_sum_only_)
301  {
302  // Create the graph's x-axis
303 
304  if (x_.size() != total_adc_values)
305  {
306  x_.resize(total_adc_values);
307 
308  std::iota(x_.begin(), x_.end(), 0);
309  }
310 
311  // If the graph doesn't exist, create it. Not sure whether to
312  // make it an error if the total_adc_values is new
313 
314  if (!graphs_[ind] || static_cast<std::size_t>(graphs_[ind]->GetN()) != total_adc_values)
315  {
316  graphs_[ind] = std::unique_ptr<TGraph>(new TGraph(total_adc_values));
317  graphs_[ind]->SetName(Form("Fragment_%d_graph", fragment_id));
318  graphs_[ind]->SetLineColor(4);
319  std::copy(x_.begin(), x_.end(), graphs_[ind]->GetX());
320  }
321 
322  // Get the data from the fragment
323 
324  // Is there some way to templatize an ART module? If not, we're stuck with this awkward switch code...
325 
326  switch (fragtype)
327  {
328  case FragmentType::TOY1:
329  case FragmentType::TOY2:
330  {
331  std::copy(toyPtr->dataBeginADCs(), toyPtr->dataBeginADCs() + total_adc_values, graphs_[ind]->GetY());
332  }
333  break;
334 
335  default:
336  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
337  }
338 
339 
340  // And now prepare the graphics without actually drawing anything yet
341 
342  canvas_[1]->cd(ind + 1);
343  TVirtualPad* pad = static_cast<TVirtualPad*>(canvas_[1]->GetPad(ind + 1));
344 
345  Double_t lo_x, hi_x, lo_y, hi_y, dummy;
346 
347  graphs_[ind]->GetPoint(0, lo_x, dummy);
348  graphs_[ind]->GetPoint(graphs_[ind]->GetN() - 1, hi_x, dummy);
349 
350  lo_x -= 0.5;
351  hi_x += 0.5;
352 
353  lo_y = -0.5;
354  hi_y = max_adc_count + 0.5;
355 
356 
357  TH1F* padframe = static_cast<TH1F*>(pad->DrawFrame(lo_x, lo_y, hi_x, hi_y));
358  padframe->SetTitle(Form("Frag %d, Type %s, SeqID %d", static_cast<int>(fragment_id),
359  fragmentTypeToString(fragtype).c_str(),
360  static_cast<int>(expected_sequence_id)));
361  padframe->GetXaxis()->SetTitle("ADC #");
362  pad->SetGrid();
363  padframe->Draw("SAME");
364  }
365 
366  // Draw the histogram
367 
368  canvas_[0]->cd(ind + 1);
369  histograms_[ind]->Draw();
370 
371  canvas_[0]->Modified();
372  canvas_[0]->Update();
373 
374  // And, if desired, the Nth event's ADC counts
375 
376  if (!digital_sum_only_)
377  {
378  canvas_[1]->cd(ind + 1);
379 
380  graphs_[ind]->Draw("PSAME");
381 
382  canvas_[1]->Modified();
383  canvas_[1]->Update();
384  }
385 
386  if (writeOutput_)
387  {
388  canvas_[0]->Write("wf0", TObject::kOverwrite);
389  canvas_[1]->Write("wf1", TObject::kOverwrite);
390  fFile_->Write();
391  }
392  }
393 }
394 
395 void demo::WFViewer::beginRun(art::Run const& e)
396 {
397  if (e.run() == current_run_) return;
398  current_run_ = e.run();
399 
400  if (writeOutput_)
401  {
402  fFile_ = new TFile(outputFileName_.c_str(), "RECREATE");
403  fFile_->cd();
404  }
405 
406  for (int i = 0; i < 2; i++) canvas_[i] = 0;
407  for (auto& x: graphs_) x = 0;
408  for (auto& x: histograms_) x = 0;
409 
410  for (int i = 0; (i < 2 && !digital_sum_only_) || i < 1; i++)
411  {
412  canvas_[i] = std::unique_ptr<TCanvas>(new TCanvas(Form("wf%d", i)));
413  canvas_[i]->Divide(num_x_plots_, num_y_plots_);
414  canvas_[i]->Update();
415  ((TRootCanvas*)canvas_[i]->GetCanvasImp())->DontCallClose();
416  }
417 
418  canvas_[0]->SetTitle("ADC Value Distribution");
419 
420  if (! digital_sum_only_)
421  {
422  canvas_[1]->SetTitle("ADC Values, Event Snapshot");
423  }
424 
425  if (writeOutput_)
426  {
427  canvas_[0]->Write();
428  canvas_[1]->Write();
429  }
430 }
431 
432 
433 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.