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