artdaq  v3_09_00
MissingDataCheck_module.cc
1 // Class: MissingDataCheck
3 // Module Type: analyzer
4 // File: MissingDataCheck_module.cc
5 // Description: Prints out information about each event.
7 
8 #include "art/Framework/Core/EDAnalyzer.h"
9 #include "art/Framework/Core/ModuleMacros.h"
10 #include "art/Framework/Principal/Event.h"
11 #include "art/Framework/Principal/Handle.h"
12 #include "canvas/Utilities/Exception.h"
13 
14 #include "artdaq-core/Data/ContainerFragment.hh"
15 #include "artdaq-core/Data/Fragment.hh"
16 
17 #include <algorithm>
18 #include <cassert>
19 #include <cmath>
20 #include <fstream>
21 #include <iomanip>
22 #include <iostream>
23 #include <set>
24 #include <vector>
25 
26 #if ART_HEX_VERSION < 0x30200
27 #include "art/Framework/Services/Optional/TFileService.h"
28 #else
29 #include "art_root_io/TFileService.h"
30 #endif
31 
32 #include "TTree.h"
33 
34 namespace artdaq {
35 class MissingDataCheck;
36 }
37 
41 class artdaq::MissingDataCheck : public art::EDAnalyzer
42 {
43 public:
55  explicit MissingDataCheck(fhicl::ParameterSet const& pset);
56 
60  ~MissingDataCheck() override = default;
61 
67  void analyze(art::Event const& e) override;
68 
72  void endJob() override;
73 
74 private:
75  MissingDataCheck(MissingDataCheck const&) = delete;
77  MissingDataCheck& operator=(MissingDataCheck const&) = delete;
78  MissingDataCheck& operator=(MissingDataCheck&&) = delete;
79 
80  std::string raw_data_label_;
81  int expected_n_fragments_;
82  int verbosity_;
83  bool make_tree_;
84 
85  TTree* evtree_;
86  TTree* fragtree_;
87  unsigned int run_;
88  unsigned int subrun_;
89  unsigned int event_;
90  unsigned int timeHigh_;
91  unsigned int timeLow_;
92 
93  unsigned int total_n_frags_;
94  unsigned int total_data_size_;
95  unsigned int total_n_CFs_;
96  unsigned int total_n_CFs_missing_;
97  unsigned int total_n_frags_in_CFs_;
98  unsigned int total_n_frags_broken_;
99 
100  unsigned int frag_id_;
101  unsigned int seq_id_;
102  unsigned int frag_type_;
103  int contained_frags_;
104  unsigned int frag_data_size_;
105  int missing_data_;
106 };
107 
108 artdaq::MissingDataCheck::MissingDataCheck(fhicl::ParameterSet const& pset)
109  : EDAnalyzer(pset)
110  , raw_data_label_(pset.get<std::string>("raw_data_label", "daq"))
111  , expected_n_fragments_(pset.get<int>("expected_n_fragments", -1))
112  , verbosity_(pset.get<int>("verbosity", 0))
113  , make_tree_(pset.get<bool>("make_tree", true))
114 {
115  if (make_tree_)
116  {
117  art::ServiceHandle<art::TFileService> tfs;
118  evtree_ = tfs->make<TTree>("evtree", "MissingDataCheck Event Tree");
119  evtree_->Branch("run", &run_, "run/i");
120  evtree_->Branch("subrun", &subrun_, "subrun/i");
121  evtree_->Branch("event", &event_, "event/i");
122  evtree_->Branch("timeHigh", &timeHigh_, "timeHigh/i");
123  evtree_->Branch("timeLow", &timeLow_, "timeLow/i");
124 
125  evtree_->Branch("n_frag_exp", &expected_n_fragments_, "nfrag_exp/I");
126  evtree_->Branch("n_frag", &total_n_frags_, "n_frag/i");
127  evtree_->Branch("data_size", &total_data_size_, "data_size/i");
128  evtree_->Branch("n_cont_frag", &total_n_CFs_, "n_cont_frag/i");
129  evtree_->Branch("n_miss_data", &total_n_CFs_missing_, "n_miss_data/i");
130  evtree_->Branch("n_frag_in_cont", &total_n_frags_in_CFs_, "n_frag_in_cont/i");
131  evtree_->Branch("n_frag_broken", &total_n_frags_broken_, "n_frag_broken/i");
132 
133  fragtree_ = tfs->make<TTree>("fragtree", "MissingDataCheck Fragment Tree");
134  fragtree_->Branch("run", &run_, "run/i");
135  fragtree_->Branch("subrun", &subrun_, "subrun/i");
136  fragtree_->Branch("event", &event_, "event/i");
137  fragtree_->Branch("timeHigh", &timeHigh_, "timeHigh/i");
138  fragtree_->Branch("timeLow", &timeLow_, "timeLow/i");
139 
140  fragtree_->Branch("frag_id", &frag_id_, "frag_id/i");
141  fragtree_->Branch("seq_id", &seq_id_, "seq_id/i");
142  fragtree_->Branch("frag_type", &frag_type_, "frag_type/i");
143  fragtree_->Branch("frag_data_size", &frag_data_size_, "frag_data_size/i");
144  fragtree_->Branch("contained_frags", &contained_frags_, "contained_frags/I");
145  fragtree_->Branch("missing_data", &missing_data_, "missing_data/I");
146  }
147 }
148 
149 void artdaq::MissingDataCheck::analyze(art::Event const& e)
150 {
151  run_ = e.run();
152  subrun_ = e.subRun();
153  event_ = e.event();
154  timeHigh_ = e.time().timeHigh();
155  timeLow_ = e.time().timeLow();
156 
157  //print basic run info
158  if (verbosity_ > 2)
159  {
160  std::cout << "Processing:"
161  << " Run " << e.run()
162  << ", Subrun " << e.subRun()
163  << ", Event " << e.event()
164  << " (Time=" << e.time().timeHigh() << " " << e.time().timeLow() << ")"
165  << std::endl;
166  }
167 
168  //get all the artdaq fragment collections in the event.
169  std::vector<art::Handle<std::vector<artdaq::Fragment> > > fragmentHandles;
170  e.getManyByType(fragmentHandles);
171 
172  //print basic fragment number info
173  if (verbosity_ > 2)
174  {
175  std::cout << "\tFound " << fragmentHandles.size() << " fragment collections." << std::endl;
176  if (verbosity_ > 2)
177  {
178  for (auto const& h : fragmentHandles)
179  {
180  std::cout << "\t\tCollection " << h.provenance()->productInstanceName()
181  << ":\t" << h->size() << " fragments." << std::endl;
182  }
183  }
184  }
185 
186  //count total fragments
187  total_n_frags_ = 0;
188  total_data_size_ = 0;
189  for (auto const& h : fragmentHandles)
190  {
191  total_n_frags_ += h->size();
192  for (auto const& f : *h)
193  {
194  total_data_size_ += f.dataSizeBytes();
195  }
196  }
197 
198  //first time through, if this is -1, set it to total fragments seen
199  if (expected_n_fragments_ == -1)
200  {
201  expected_n_fragments_ = total_n_frags_;
202  }
203 
204  if (verbosity_ > 2)
205  {
206  std::cout << "\tTotal fragments = " << total_n_frags_
207  << " / " << expected_n_fragments_ << std::endl;
208  std::cout << "\tTotal data size in fragments = " << total_data_size_ << std::endl;
209  }
210  //count number of container fragments
211  total_n_CFs_ = 0;
212  total_n_CFs_missing_ = 0;
213  total_n_frags_in_CFs_ = 0;
214  total_n_frags_broken_ = 0;
215  for (auto const& h : fragmentHandles)
216  {
217  std::string instance_name = h.provenance()->productInstanceName();
218 
219  if (instance_name.compare(0, 9, "Container") == 0)
220  {
221  total_n_CFs_ += h->size();
222 
223  for (auto const& f : *h)
224  {
225  artdaq::ContainerFragment cf(f);
226  if (cf.missing_data())
227  {
228  ++total_n_CFs_missing_;
229  }
230  total_n_frags_in_CFs_ += cf.block_count();
231  frag_id_ = f.fragmentID();
232  seq_id_ = f.sequenceID();
233  frag_type_ = f.type();
234  contained_frags_ = cf.block_count();
235  missing_data_ = cf.missing_data() ? 1 : 0;
236  frag_data_size_ = f.dataSizeBytes();
237  fragtree_->Fill();
238  }
239  }
240 
241  else
242  {
243  for (auto const& f : *h)
244  {
245  if (instance_name.compare(0, 4, "Data") == 0 ||
246  instance_name.compare(0, 5, "Error") == 0 ||
247  instance_name.compare(0, 6, "Broken") == 0)
248  {
249  ++total_n_frags_broken_;
250  }
251 
252  frag_id_ = f.fragmentID();
253  seq_id_ = f.sequenceID();
254  frag_type_ = f.type();
255  contained_frags_ = -1;
256  missing_data_ = -1;
257  frag_data_size_ = f.dataSizeBytes();
258  fragtree_->Fill();
259  }
260  }
261  }
262 
263  if (verbosity_ > 2)
264  {
265  std::cout << "\tTotal container fragments = " << total_n_CFs_
266  << " of which " << total_n_CFs_missing_ << " have missing data."
267  << std::endl;
268  std::cout << "\tTotal fragments in containers = " << total_n_frags_in_CFs_ << std::endl;
269  }
270 
271  if (make_tree_)
272  {
273  evtree_->Fill();
274  }
275 }
276 
278 {
279  if (verbosity_ > 0)
280  {
281  std::cout << "---------------------------------------------------" << std::endl;
282  std::cout << "----------- MISSING DATA CHECK SUMMARY ------------" << std::endl;
283  std::cout << "---------------------------------------------------" << std::endl;
284 
285  std::cout << "Total events processed: " << evtree_->GetEntries()
286  << ", expected fragments: " << expected_n_fragments_ << std::endl;
287 
288  unsigned int run;
289  evtree_->SetBranchAddress("run", &run);
290  unsigned int subrun;
291  evtree_->SetBranchAddress("subrun", &subrun);
292  unsigned int event;
293  evtree_->SetBranchAddress("event", &event);
294  int n_frag_exp;
295  evtree_->SetBranchAddress("n_frag_exp", &n_frag_exp);
296  unsigned int n_frag;
297  evtree_->SetBranchAddress("n_frag", &n_frag);
298  unsigned int n_frag_broken;
299  evtree_->SetBranchAddress("n_frag_broken", &n_frag_broken);
300  unsigned int n_miss_data;
301  evtree_->SetBranchAddress("n_miss_data", &n_miss_data);
302 
303  std::cout << "Events missing fragments:\t\t"
304  << evtree_->GetEntries("n_frag<n_frag_exp")
305  << " / " << evtree_->GetEntries() << std::endl;
306  if (verbosity_ > 1)
307  {
308  if (evtree_->GetEntries("n_frag<n_frag_exp") > 0)
309  {
310  for (int i = 0; i < evtree_->GetEntries(); ++i)
311  {
312  evtree_->GetEvent(i);
313  if (static_cast<int>(n_frag) < n_frag_exp)
314  {
315  std::cout << "\tEvent (" << run << "," << subrun << "," << event << ")"
316  << " is missing " << n_frag_exp - n_frag << " fragments."
317  << std::endl;
318  }
319  }
320  }
321  }
322  std::cout << "Events with broken fragments:\t\t"
323  << evtree_->GetEntries("n_frag_broken>0")
324  << " / " << evtree_->GetEntries() << std::endl;
325  if (evtree_->GetEntries("n_frag_broken>0") > 0)
326  {
327  for (int i = 0; i < evtree_->GetEntries(); ++i)
328  {
329  evtree_->GetEvent(i);
330  if (n_frag_broken > 0)
331  {
332  std::cout << "\tEvent (" << run << "," << subrun << "," << event << ")"
333  << " has " << n_frag_broken << " fragments."
334  << std::endl;
335  }
336  }
337  }
338  std::cout << "Events missing data in fragments:\t"
339  << evtree_->GetEntries("n_miss_data>0")
340  << " / " << evtree_->GetEntries() << std::endl;
341  if (verbosity_ > 1)
342  {
343  if (evtree_->GetEntries("n_miss_data>0") > 0)
344  {
345  for (int i = 0; i < evtree_->GetEntries(); ++i)
346  {
347  evtree_->GetEvent(i);
348  if (n_miss_data > 0)
349  {
350  std::cout << "\tEvent (" << run << "," << subrun << "," << event << ")"
351  << " has " << n_miss_data << " fragments missing data."
352  << std::endl;
353  }
354  }
355  }
356  }
357 
358  std::cout << "---------------------------------------------------" << std::endl;
359  std::cout << "---------------------------------------------------" << std::endl;
360  }
361 }
362 
363 DEFINE_ART_MODULE(artdaq::MissingDataCheck)// NOLINT(performance-unnecessary-value-param)
void endJob() override
This method is called at the end of the job, used to print a summary.
Check art::Event for potentially missing data.
~MissingDataCheck() override=default
Default virtual Destructor.
MissingDataCheck(fhicl::ParameterSet const &pset)
MissingDataCheck Constructor.
void analyze(art::Event const &e) override
This method is called for each art::Event in a file or run.