artdaq_demo  v2_10_00
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator
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  std::vector<std::string> fragment_type_labels{"TOY1","TOY2","Container"};
170 
171  for (auto label: fragment_type_labels)
172  {
173  art::Handle<artdaq::Fragments> fragments_with_label;
174  e.getByLabel("daq", label, fragments_with_label);
175 
176  if (!fragments_with_label.isValid()) continue;
177  // for (int i_l = 0; i_l < static_cast<int>(fragments_with_label->size()); ++i_l) {
178  // fragments.emplace_back( (*fragments_with_label)[i_l] );
179  // }
180 
181  if (label == "Container")
182  {
183  for (auto cont : *fragments_with_label)
184  {
185  artdaq::ContainerFragment contf(cont);
186  for (size_t ii = 0; ii < contf.block_count(); ++ii)
187  {
188  size_t fragSize = contf.fragSize(ii);
189  artdaq::Fragment thisfrag;
190  thisfrag.resizeBytes(fragSize);
191 
192  //mf::LogDebug("WFViewer") << "Copying " << fragSize << " bytes from " << contf.at(ii) << " to " << thisfrag.headerAddress();
193  memcpy(thisfrag.headerAddress(), contf.at(ii), fragSize);
194 
195  //mf::LogDebug("WFViewer") << "Putting new fragment into output vector";
196  fragments.push_back(thisfrag);
197  }
198  }
199  }
200  else
201  {
202  for (auto frag : *fragments_with_label)
203  {
204  fragments.emplace_back(frag);
205  }
206  }
207  }
208 
209  // John F., 1/5/14
210 
211  // Here, we loop over the fragments passed to the analyze
212  // function. A warning is flashed if either (A) the fragments aren't
213  // all from the same event, or (B) there's an unexpected number of
214  // fragments given the number of boardreaders and the number of
215  // fragments per board
216 
217  // For every Nth event, where N is the "prescale" setting, plot the
218  // distribution of ADC counts from each board_id / fragment_id
219  // combo. Also, if "digital_sum_only" is set to false in the FHiCL
220  // string, then plot, for the Nth event, a graph of the ADC values
221  // across all channels in each board_id / fragment_id combo
222 
223  artdaq::Fragment::sequence_id_t expected_sequence_id = std::numeric_limits<artdaq::Fragment::sequence_id_t>::max();
224 
225  // for (std::size_t i = 0; i < fragments.size(); ++i) {
226  for (const auto& frag : fragments)
227  {
228  // Pointers to the types of fragment overlays WFViewer can handle;
229  // only one will be used per fragment, of course
230 
231  std::unique_ptr<ToyFragment> toyPtr;
232 
233  // const auto& frag( fragments[i] ); // Basically a shorthand
234 
235  // if (i == 0)
236  if (expected_sequence_id == std::numeric_limits<artdaq::Fragment::sequence_id_t>::max())
237  {
238  expected_sequence_id = frag.sequenceID();
239  }
240 
241  if (expected_sequence_id != frag.sequenceID())
242  {
243  cerr << "Warning in WFViewer: expected fragment with sequence ID " << expected_sequence_id << ", received one with sequence ID " << frag.sequenceID() << endl;
244  }
245 
246  FragmentType fragtype = static_cast<FragmentType>(frag.type());
247  std::size_t max_adc_count = std::numeric_limits<std::size_t>::max();
248  std::size_t total_adc_values = std::numeric_limits<std::size_t>::max();
249 
250  switch (fragtype)
251  {
252  case FragmentType::TOY1:
253  toyPtr.reset(new ToyFragment(frag));
254  total_adc_values = toyPtr->total_adc_values();
255  max_adc_count = pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1;
256  break;
257  case FragmentType::TOY2:
258  toyPtr.reset(new ToyFragment(frag));
259  total_adc_values = toyPtr->total_adc_values();
260  max_adc_count = pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1;
261  break;
262  default:
263  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
264  }
265 
266  artdaq::Fragment::fragment_id_t fragment_id = frag.fragmentID();
267  std::size_t ind = id_to_index_[fragment_id];
268 
269 
270  // If a histogram doesn't exist for this board_id / fragment_id combo, create it
271 
272  if (!histograms_[ind])
273  {
274  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));
275 
276  histograms_[ind]->SetTitle(Form("Frag %d, Type %s", fragment_id,
277  fragmentTypeToString(fragtype).c_str()));
278  histograms_[ind]->GetXaxis()->SetTitle("ADC value");
279  }
280 
281  // For every event, fill the histogram (prescale is ignored here)
282 
283  // Is there some way to templatize an ART module? If not, we're
284  // stuck with this switch code...
285 
286  switch (fragtype)
287  {
288  case FragmentType::TOY1:
289  case FragmentType::TOY2:
290  for (auto val = toyPtr->dataBeginADCs(); val != toyPtr->dataEndADCs(); ++val)
291  histograms_[ind]->Fill(*val);
292  break;
293 
294  default:
295  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
296  }
297 
298  if (evt_cntr % prescale_ - 1 && prescale_ > 1)
299  {
300  continue;
301  }
302 
303  // If we pass the prescale, then if we're not going with
304  // digital_sum_only, plot the ADC counts for this particular event/board/fragment_id
305 
306  if (!digital_sum_only_)
307  {
308  // Create the graph's x-axis
309 
310  if (x_.size() != total_adc_values)
311  {
312  x_.resize(total_adc_values);
313 
314  std::iota(x_.begin(), x_.end(), 0);
315  }
316 
317  // If the graph doesn't exist, create it. Not sure whether to
318  // make it an error if the total_adc_values is new
319 
320  if (!graphs_[ind] || static_cast<std::size_t>(graphs_[ind]->GetN()) != total_adc_values)
321  {
322  graphs_[ind] = std::unique_ptr<TGraph>(new TGraph(total_adc_values));
323  graphs_[ind]->SetName(Form("Fragment_%d_graph", fragment_id));
324  graphs_[ind]->SetLineColor(4);
325  std::copy(x_.begin(), x_.end(), graphs_[ind]->GetX());
326  }
327 
328  // Get the data from the fragment
329 
330  // Is there some way to templatize an ART module? If not, we're stuck with this awkward switch code...
331 
332  switch (fragtype)
333  {
334  case FragmentType::TOY1:
335  case FragmentType::TOY2:
336  {
337  std::copy(toyPtr->dataBeginADCs(), toyPtr->dataBeginADCs() + total_adc_values, graphs_[ind]->GetY());
338  }
339  break;
340 
341  default:
342  throw cet::exception("Error in WFViewer: unknown fragment type supplied");
343  }
344 
345 
346  // And now prepare the graphics without actually drawing anything yet
347 
348  canvas_[1]->cd(ind + 1);
349  TVirtualPad* pad = static_cast<TVirtualPad*>(canvas_[1]->GetPad(ind + 1));
350 
351  Double_t lo_x, hi_x, lo_y, hi_y, dummy;
352 
353  graphs_[ind]->GetPoint(0, lo_x, dummy);
354  graphs_[ind]->GetPoint(graphs_[ind]->GetN() - 1, hi_x, dummy);
355 
356  lo_x -= 0.5;
357  hi_x += 0.5;
358 
359  lo_y = -0.5;
360  hi_y = max_adc_count + 0.5;
361 
362 
363  TH1F* padframe = static_cast<TH1F*>(pad->DrawFrame(lo_x, lo_y, hi_x, hi_y));
364  padframe->SetTitle(Form("Frag %d, Type %s, SeqID %d", static_cast<int>(fragment_id),
365  fragmentTypeToString(fragtype).c_str(),
366  static_cast<int>(expected_sequence_id)));
367  padframe->GetXaxis()->SetTitle("ADC #");
368  pad->SetGrid();
369  padframe->Draw("SAME");
370  }
371 
372  // Draw the histogram
373 
374  canvas_[0]->cd(ind + 1);
375  histograms_[ind]->Draw();
376 
377  canvas_[0]->Modified();
378  canvas_[0]->Update();
379 
380  // And, if desired, the Nth event's ADC counts
381 
382  if (!digital_sum_only_)
383  {
384  canvas_[1]->cd(ind + 1);
385 
386  graphs_[ind]->Draw("PSAME");
387 
388  canvas_[1]->Modified();
389  canvas_[1]->Update();
390  }
391 
392  if (writeOutput_)
393  {
394  canvas_[0]->Write("wf0", TObject::kOverwrite);
395  canvas_[1]->Write("wf1", TObject::kOverwrite);
396  fFile_->Write();
397  }
398  }
399 }
400 
401 void demo::WFViewer::beginRun(art::Run const& e)
402 {
403  if (e.run() == current_run_) return;
404  current_run_ = e.run();
405 
406  if (writeOutput_)
407  {
408  fFile_ = new TFile(outputFileName_.c_str(), "RECREATE");
409  fFile_->cd();
410  }
411 
412  for (int i = 0; i < 2; i++) canvas_[i] = 0;
413  for (auto& x: graphs_) x = 0;
414  for (auto& x: histograms_) x = 0;
415 
416  for (int i = 0; (i < 2 && !digital_sum_only_) || i < 1; i++)
417  {
418  canvas_[i] = std::unique_ptr<TCanvas>(new TCanvas(Form("wf%d", i)));
419  canvas_[i]->Divide(num_x_plots_, num_y_plots_);
420  canvas_[i]->Update();
421  ((TRootCanvas*)canvas_[i]->GetCanvasImp())->DontCallClose();
422  }
423 
424  canvas_[0]->SetTitle("ADC Value Distribution");
425 
426  if (! digital_sum_only_)
427  {
428  canvas_[1]->SetTitle("ADC Values, Event Snapshot");
429  }
430 
431  if (writeOutput_)
432  {
433  canvas_[0]->Write();
434  canvas_[1]->Write();
435  }
436 }
437 
438 
439 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.