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