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