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/ContainerFragment.hh"
15 #include "artdaq-core/Data/Fragment.hh"
26 #if ART_HEX_VERSION < 0x30200
27 #include "art/Framework/Services/Optional/TFileService.h"
29 #include "art_root_io/TFileService.h"
35 class MissingDataCheck;
67 void analyze(art::Event
const& e)
override;
80 std::string raw_data_label_;
81 int expected_n_fragments_;
90 unsigned int timeHigh_;
91 unsigned int timeLow_;
93 unsigned int total_n_frags_;
94 unsigned int total_data_size_;
95 unsigned int total_n_CFs_;
96 unsigned int total_n_CFs_missing_;
97 unsigned int total_n_frags_in_CFs_;
98 unsigned int total_n_frags_broken_;
100 unsigned int frag_id_;
101 unsigned int seq_id_;
102 unsigned int frag_type_;
103 int contained_frags_;
104 unsigned int frag_data_size_;
110 , raw_data_label_(pset.get<std::string>(
"raw_data_label",
"daq"))
111 , expected_n_fragments_(pset.get<int>(
"expected_n_fragments", -1))
112 , verbosity_(pset.get<int>(
"verbosity", 0))
113 , make_tree_(pset.get<bool>(
"make_tree", true))
117 art::ServiceHandle<art::TFileService> tfs;
118 evtree_ = tfs->make<TTree>(
"evtree",
"MissingDataCheck Event Tree");
119 evtree_->Branch(
"run", &run_,
"run/i");
120 evtree_->Branch(
"subrun", &subrun_,
"subrun/i");
121 evtree_->Branch(
"event", &event_,
"event/i");
122 evtree_->Branch(
"timeHigh", &timeHigh_,
"timeHigh/i");
123 evtree_->Branch(
"timeLow", &timeLow_,
"timeLow/i");
125 evtree_->Branch(
"n_frag_exp", &expected_n_fragments_,
"nfrag_exp/I");
126 evtree_->Branch(
"n_frag", &total_n_frags_,
"n_frag/i");
127 evtree_->Branch(
"data_size", &total_data_size_,
"data_size/i");
128 evtree_->Branch(
"n_cont_frag", &total_n_CFs_,
"n_cont_frag/i");
129 evtree_->Branch(
"n_miss_data", &total_n_CFs_missing_,
"n_miss_data/i");
130 evtree_->Branch(
"n_frag_in_cont", &total_n_frags_in_CFs_,
"n_frag_in_cont/i");
131 evtree_->Branch(
"n_frag_broken", &total_n_frags_broken_,
"n_frag_broken/i");
133 fragtree_ = tfs->make<TTree>(
"fragtree",
"MissingDataCheck Fragment Tree");
134 fragtree_->Branch(
"run", &run_,
"run/i");
135 fragtree_->Branch(
"subrun", &subrun_,
"subrun/i");
136 fragtree_->Branch(
"event", &event_,
"event/i");
137 fragtree_->Branch(
"timeHigh", &timeHigh_,
"timeHigh/i");
138 fragtree_->Branch(
"timeLow", &timeLow_,
"timeLow/i");
140 fragtree_->Branch(
"frag_id", &frag_id_,
"frag_id/i");
141 fragtree_->Branch(
"seq_id", &seq_id_,
"seq_id/i");
142 fragtree_->Branch(
"frag_type", &frag_type_,
"frag_type/i");
143 fragtree_->Branch(
"frag_data_size", &frag_data_size_,
"frag_data_size/i");
144 fragtree_->Branch(
"contained_frags", &contained_frags_,
"contained_frags/I");
145 fragtree_->Branch(
"missing_data", &missing_data_,
"missing_data/I");
152 subrun_ = e.subRun();
154 timeHigh_ = e.time().timeHigh();
155 timeLow_ = e.time().timeLow();
160 std::cout <<
"Processing:"
161 <<
" Run " << e.run()
162 <<
", Subrun " << e.subRun()
163 <<
", Event " << e.event()
164 <<
" (Time=" << e.time().timeHigh() <<
" " << e.time().timeLow() <<
")"
169 std::vector<art::Handle<std::vector<artdaq::Fragment> > > fragmentHandles;
170 e.getManyByType(fragmentHandles);
175 std::cout <<
"\tFound " << fragmentHandles.size() <<
" fragment collections." << std::endl;
178 for (
auto const& h : fragmentHandles)
180 std::cout <<
"\t\tCollection " << h.provenance()->productInstanceName()
181 <<
":\t" << h->size() <<
" fragments." << std::endl;
188 total_data_size_ = 0;
189 for (
auto const& h : fragmentHandles)
191 total_n_frags_ += h->size();
192 for (
auto const& f : *h)
194 total_data_size_ += f.dataSizeBytes();
199 if (expected_n_fragments_ == -1)
201 expected_n_fragments_ = total_n_frags_;
206 std::cout <<
"\tTotal fragments = " << total_n_frags_
207 <<
" / " << expected_n_fragments_ << std::endl;
208 std::cout <<
"\tTotal data size in fragments = " << total_data_size_ << std::endl;
212 total_n_CFs_missing_ = 0;
213 total_n_frags_in_CFs_ = 0;
214 total_n_frags_broken_ = 0;
215 for (
auto const& h : fragmentHandles)
217 std::string instance_name = h.provenance()->productInstanceName();
219 if (instance_name.compare(0, 9,
"Container") == 0)
221 total_n_CFs_ += h->size();
223 for (
auto const& f : *h)
225 artdaq::ContainerFragment cf(f);
226 if (cf.missing_data())
228 ++total_n_CFs_missing_;
230 total_n_frags_in_CFs_ += cf.block_count();
231 frag_id_ = f.fragmentID();
232 seq_id_ = f.sequenceID();
233 frag_type_ = f.type();
234 contained_frags_ = cf.block_count();
235 missing_data_ = cf.missing_data() ? 1 : 0;
236 frag_data_size_ = f.dataSizeBytes();
243 for (
auto const& f : *h)
245 if (instance_name.compare(0, 4,
"Data") == 0 ||
246 instance_name.compare(0, 5,
"Error") == 0 ||
247 instance_name.compare(0, 6,
"Broken") == 0)
249 ++total_n_frags_broken_;
252 frag_id_ = f.fragmentID();
253 seq_id_ = f.sequenceID();
254 frag_type_ = f.type();
255 contained_frags_ = -1;
257 frag_data_size_ = f.dataSizeBytes();
265 std::cout <<
"\tTotal container fragments = " << total_n_CFs_
266 <<
" of which " << total_n_CFs_missing_ <<
" have missing data."
268 std::cout <<
"\tTotal fragments in containers = " << total_n_frags_in_CFs_ << std::endl;
281 std::cout <<
"---------------------------------------------------" << std::endl;
282 std::cout <<
"----------- MISSING DATA CHECK SUMMARY ------------" << std::endl;
283 std::cout <<
"---------------------------------------------------" << std::endl;
285 std::cout <<
"Total events processed: " << evtree_->GetEntries()
286 <<
", expected fragments: " << expected_n_fragments_ << std::endl;
289 evtree_->SetBranchAddress(
"run", &run);
291 evtree_->SetBranchAddress(
"subrun", &subrun);
293 evtree_->SetBranchAddress(
"event", &event);
295 evtree_->SetBranchAddress(
"n_frag_exp", &n_frag_exp);
297 evtree_->SetBranchAddress(
"n_frag", &n_frag);
298 unsigned int n_frag_broken;
299 evtree_->SetBranchAddress(
"n_frag_broken", &n_frag_broken);
300 unsigned int n_miss_data;
301 evtree_->SetBranchAddress(
"n_miss_data", &n_miss_data);
303 std::cout <<
"Events missing fragments:\t\t"
304 << evtree_->GetEntries(
"n_frag<n_frag_exp")
305 <<
" / " << evtree_->GetEntries() << std::endl;
308 if (evtree_->GetEntries(
"n_frag<n_frag_exp") > 0)
310 for (
int i = 0; i < evtree_->GetEntries(); ++i)
312 evtree_->GetEvent(i);
313 if (static_cast<int>(n_frag) < n_frag_exp)
315 std::cout <<
"\tEvent (" << run <<
"," << subrun <<
"," <<
event <<
")"
316 <<
" is missing " << n_frag_exp - n_frag <<
" fragments."
322 std::cout <<
"Events with broken fragments:\t\t"
323 << evtree_->GetEntries(
"n_frag_broken>0")
324 <<
" / " << evtree_->GetEntries() << std::endl;
325 if (evtree_->GetEntries(
"n_frag_broken>0") > 0)
327 for (
int i = 0; i < evtree_->GetEntries(); ++i)
329 evtree_->GetEvent(i);
330 if (n_frag_broken > 0)
332 std::cout <<
"\tEvent (" << run <<
"," << subrun <<
"," <<
event <<
")"
333 <<
" has " << n_frag_broken <<
" fragments."
338 std::cout <<
"Events missing data in fragments:\t"
339 << evtree_->GetEntries(
"n_miss_data>0")
340 <<
" / " << evtree_->GetEntries() << std::endl;
343 if (evtree_->GetEntries(
"n_miss_data>0") > 0)
345 for (
int i = 0; i < evtree_->GetEntries(); ++i)
347 evtree_->GetEvent(i);
350 std::cout <<
"\tEvent (" << run <<
"," << subrun <<
"," <<
event <<
")"
351 <<
" has " << n_miss_data <<
" fragments missing data."
358 std::cout <<
"---------------------------------------------------" << std::endl;
359 std::cout <<
"---------------------------------------------------" << std::endl;
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.