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"
14 #include "artdaq-core/Data/Fragment.hh"
15 #include "artdaq-core/Data/ContainerFragment.hh"
26 #if ART_HEX_VERSION < 0x30200
27 # include "art/Framework/Services/Optional/TFileService.h"
29 # include "art_root_io/TFileService.h"
36 class MissingDataCheck;
68 void analyze(art::Event
const& e)
override;
76 std::string raw_data_label_;
77 int expected_n_fragments_;
86 unsigned int timeHigh_;
87 unsigned int timeLow_;
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_;
97 unsigned int frag_id_;
99 unsigned int frag_type_;
100 int contained_frags_;
101 unsigned int frag_data_size_;
109 , raw_data_label_(pset.get<std::string>(
"raw_data_label",
"daq"))
110 , expected_n_fragments_(pset.get<int>(
"expected_n_fragments",-1))
111 , verbosity_(pset.get<int>(
"verbosity",0))
112 , make_tree_(pset.get<bool>(
"make_tree",true))
115 art::ServiceHandle<art::TFileService> tfs;
116 evtree_ = tfs->make<TTree>(
"evtree",
"MissingDataCheck Event Tree");
117 evtree_->Branch(
"run",&run_,
"run/i");
118 evtree_->Branch(
"subrun",&subrun_,
"subrun/i");
119 evtree_->Branch(
"event",&event_,
"event/i");
120 evtree_->Branch(
"timeHigh",&timeHigh_,
"timeHigh/i");
121 evtree_->Branch(
"timeLow",&timeLow_,
"timeLow/i");
123 evtree_->Branch(
"n_frag_exp",&expected_n_fragments_,
"nfrag_exp/I");
124 evtree_->Branch(
"n_frag",&total_n_frags_,
"n_frag/i");
125 evtree_->Branch(
"data_size",&total_data_size_,
"data_size/i");
126 evtree_->Branch(
"n_cont_frag",&total_n_CFs_,
"n_cont_frag/i");
127 evtree_->Branch(
"n_miss_data",&total_n_CFs_missing_,
"n_miss_data/i");
128 evtree_->Branch(
"n_frag_in_cont",&total_n_frags_in_CFs_,
"n_frag_in_cont/i");
129 evtree_->Branch(
"n_frag_broken",&total_n_frags_broken_,
"n_frag_broken/i");
131 fragtree_ = tfs->make<TTree>(
"fragtree",
"MissingDataCheck Fragment Tree");
132 fragtree_->Branch(
"run",&run_,
"run/i");
133 fragtree_->Branch(
"subrun",&subrun_,
"subrun/i");
134 fragtree_->Branch(
"event",&event_,
"event/i");
135 fragtree_->Branch(
"timeHigh",&timeHigh_,
"timeHigh/i");
136 fragtree_->Branch(
"timeLow",&timeLow_,
"timeLow/i");
138 fragtree_->Branch(
"frag_id",&frag_id_,
"frag_id/i");
139 fragtree_->Branch(
"seq_id",&seq_id_,
"seq_id/i");
140 fragtree_->Branch(
"frag_type",&frag_type_,
"frag_type/i");
141 fragtree_->Branch(
"frag_data_size",&frag_data_size_,
"frag_data_size/i");
142 fragtree_->Branch(
"contained_frags",&contained_frags_,
"contained_frags/I");
143 fragtree_->Branch(
"missing_data",&missing_data_,
"missing_data/I");
151 subrun_ = e.subRun();
153 timeHigh_ = e.time().timeHigh();
154 timeLow_ = e.time().timeLow();
158 std::cout <<
"Processing:"
159 <<
" Run " << e.run()
160 <<
", Subrun " << e.subRun()
161 <<
", Event " << e.event()
162 <<
" (Time=" << e.time().timeHigh() <<
" " << e.time().timeLow() <<
")"
166 std::vector< art::Handle< std::vector<artdaq::Fragment> > > fragmentHandles;
167 e.getManyByType(fragmentHandles);
171 std::cout <<
"\tFound " << fragmentHandles.size() <<
" fragment collections." << std::endl;
173 for(
auto const& h : fragmentHandles)
174 std::cout <<
"\t\tCollection " << h.provenance()->productInstanceName()
175 <<
":\t" << h->size() <<
" fragments." << std::endl;
182 for(
auto const& h : fragmentHandles){
183 total_n_frags_ += h->size();
184 for(
auto const& f : *h)
185 total_data_size_ += f.dataSizeBytes();
189 if(expected_n_fragments_==-1)
190 expected_n_fragments_=total_n_frags_;
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;
200 total_n_CFs_missing_=0;
201 total_n_frags_in_CFs_=0;
202 total_n_frags_broken_=0;
203 for(
auto const& h : fragmentHandles){
204 std::string instance_name = h.provenance()->productInstanceName();
206 if(instance_name.compare(0,9,
"Container")==0){
207 total_n_CFs_ += h->size();
209 for(
auto const& f : *h){
210 artdaq::ContainerFragment cf(f);
211 if(cf.missing_data()){
212 ++total_n_CFs_missing_;
214 total_n_frags_in_CFs_ += cf.block_count();
215 frag_id_ = f.fragmentID();
216 seq_id_ = f.sequenceID();
217 frag_type_ = f.type();
218 contained_frags_ = cf.block_count();
219 missing_data_ = cf.missing_data()?1:0;
220 frag_data_size_ = f.dataSizeBytes();
227 for(
auto const& f : *h){
229 if(instance_name.compare(0,4,
"Data")==0 ||
230 instance_name.compare(0,5,
"Error")==0 ||
231 instance_name.compare(0,6,
"Broken")==0)
232 ++total_n_frags_broken_;
234 frag_id_ = f.fragmentID();
235 seq_id_ = f.sequenceID();
236 frag_type_ = f.type();
237 contained_frags_ = -1;
239 frag_data_size_ = f.dataSizeBytes();
248 std::cout <<
"\tTotal container fragments = " << total_n_CFs_
249 <<
" of which " << total_n_CFs_missing_ <<
" have missing data."
251 std::cout <<
"\tTotal fragments in containers = " << total_n_frags_in_CFs_ << std::endl;
263 std::cout <<
"---------------------------------------------------" << std::endl;
264 std::cout <<
"----------- MISSING DATA CHECK SUMMARY ------------" << std::endl;
265 std::cout <<
"---------------------------------------------------" << std::endl;
267 std::cout <<
"Total events processed: " << evtree_->GetEntries()
268 <<
", expected fragments: " << expected_n_fragments_ << std::endl;
270 unsigned int run; evtree_->SetBranchAddress(
"run",&run);
271 unsigned int subrun; evtree_->SetBranchAddress(
"subrun",&subrun);
272 unsigned int event; evtree_->SetBranchAddress(
"event",&event);
273 int n_frag_exp; evtree_->SetBranchAddress(
"n_frag_exp",&n_frag_exp);
274 unsigned int n_frag; evtree_->SetBranchAddress(
"n_frag",&n_frag);
275 unsigned int n_frag_broken; evtree_->SetBranchAddress(
"n_frag_broken",&n_frag_broken);
276 unsigned int n_miss_data; evtree_->SetBranchAddress(
"n_miss_data",&n_miss_data);
278 std::cout <<
"Events missing fragments:\t\t"
279 << evtree_->GetEntries(
"n_frag<n_frag_exp")
280 <<
" / " << evtree_->GetEntries() << std::endl;
282 if(evtree_->GetEntries(
"n_frag<n_frag_exp")>0){
283 for(
int i=0; i<evtree_->GetEntries(); ++i)
285 evtree_->GetEvent(i);
286 if((
int)n_frag<n_frag_exp){
287 std::cout <<
"\tEvent (" << run <<
"," << subrun <<
"," <<
event <<
")"
288 <<
" is missing " << n_frag_exp-n_frag <<
" fragments."
294 std::cout <<
"Events with broken fragments:\t\t"
295 << evtree_->GetEntries(
"n_frag_broken>0")
296 <<
" / " << evtree_->GetEntries() << std::endl;
297 if(evtree_->GetEntries(
"n_frag_broken>0")>0){
298 for(
int i=0; i<evtree_->GetEntries(); ++i)
300 evtree_->GetEvent(i);
301 if(n_frag_broken > 0){
302 std::cout <<
"\tEvent (" << run <<
"," << subrun <<
"," <<
event <<
")"
303 <<
" has " << n_frag_broken <<
" fragments."
308 std::cout <<
"Events missing data in fragments:\t"
309 << evtree_->GetEntries(
"n_miss_data>0")
310 <<
" / " << evtree_->GetEntries() << std::endl;
312 if(evtree_->GetEntries(
"n_miss_data>0")>0){
313 for(
int i=0; i<evtree_->GetEntries(); ++i)
315 evtree_->GetEvent(i);
317 std::cout <<
"\tEvent (" << run <<
"," << subrun <<
"," <<
event <<
")"
318 <<
" has " << n_miss_data <<
" fragments missing data."
325 std::cout <<
"---------------------------------------------------" << std::endl;
326 std::cout <<
"---------------------------------------------------" << std::endl;
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.