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