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