00001
00002
00003
00004
00005
00007
00008 #include "art/Framework/Core/EDAnalyzer.h"
00009 #include "art/Framework/Core/ModuleMacros.h"
00010 #include "art/Framework/Principal/Event.h"
00011 #include "art/Framework/Principal/Handle.h"
00012 #include "canvas/Utilities/Exception.h"
00013
00014 #include "artdaq-core/Data/Fragment.hh"
00015 #include "artdaq-core/Data/ContainerFragment.hh"
00016
00017 #include <algorithm>
00018 #include <cassert>
00019 #include <cmath>
00020 #include <fstream>
00021 #include <iomanip>
00022 #include <vector>
00023 #include <set>
00024 #include <iostream>
00025
00026 #include "art/Framework/Services/Optional/TFileService.h"
00027 #include "TTree.h"
00028
00029 namespace artdaq
00030 {
00031 class MissingDataCheck;
00032 }
00033
00037 class artdaq::MissingDataCheck : public art::EDAnalyzer
00038 {
00039 public:
00051 explicit MissingDataCheck(fhicl::ParameterSet const& pset);
00052
00056 virtual ~MissingDataCheck() = default;
00057
00063 void analyze(art::Event const& e) override;
00064
00065 void endJob() override;
00066
00067 private:
00068 std::string raw_data_label_;
00069 int expected_n_fragments_;
00070 int verbosity_;
00071 bool make_tree_;
00072
00073 TTree* evtree_;
00074 TTree* fragtree_;
00075 unsigned int run_;
00076 unsigned int subrun_;
00077 unsigned int event_;
00078 unsigned int timeHigh_;
00079 unsigned int timeLow_;
00080
00081 unsigned int total_n_frags_;
00082 unsigned int total_data_size_;
00083 unsigned int total_n_CFs_;
00084 unsigned int total_n_CFs_missing_;
00085 unsigned int total_n_frags_in_CFs_;
00086 unsigned int total_n_frags_broken_;
00087
00088
00089 unsigned int frag_id_;
00090 unsigned int seq_id_;
00091 unsigned int frag_type_;
00092 int contained_frags_;
00093 unsigned int frag_data_size_;
00094 int missing_data_;
00095
00096 };
00097
00098
00099 artdaq::MissingDataCheck::MissingDataCheck(fhicl::ParameterSet const& pset)
00100 : EDAnalyzer(pset)
00101 , raw_data_label_(pset.get<std::string>("raw_data_label", "daq"))
00102 , expected_n_fragments_(pset.get<int>("expected_n_fragments",-1))
00103 , verbosity_(pset.get<int>("verbosity",0))
00104 , make_tree_(pset.get<bool>("make_tree",true))
00105 {
00106 if(make_tree_){
00107 art::ServiceHandle<art::TFileService> tfs;
00108 evtree_ = tfs->make<TTree>("evtree","MissingDataCheck Event Tree");
00109 evtree_->Branch("run",&run_,"run/i");
00110 evtree_->Branch("subrun",&subrun_,"subrun/i");
00111 evtree_->Branch("event",&event_,"event/i");
00112 evtree_->Branch("timeHigh",&timeHigh_,"timeHigh/i");
00113 evtree_->Branch("timeLow",&timeLow_,"timeLow/i");
00114
00115 evtree_->Branch("n_frag_exp",&expected_n_fragments_,"nfrag_exp/I");
00116 evtree_->Branch("n_frag",&total_n_frags_,"n_frag/i");
00117 evtree_->Branch("data_size",&total_data_size_,"data_size/i");
00118 evtree_->Branch("n_cont_frag",&total_n_CFs_,"n_cont_frag/i");
00119 evtree_->Branch("n_miss_data",&total_n_CFs_missing_,"n_miss_data/i");
00120 evtree_->Branch("n_frag_in_cont",&total_n_frags_in_CFs_,"n_frag_in_cont/i");
00121 evtree_->Branch("n_frag_broken",&total_n_frags_broken_,"n_frag_broken/i");
00122
00123 fragtree_ = tfs->make<TTree>("fragtree","MissingDataCheck Fragment Tree");
00124 fragtree_->Branch("run",&run_,"run/i");
00125 fragtree_->Branch("subrun",&subrun_,"subrun/i");
00126 fragtree_->Branch("event",&event_,"event/i");
00127 fragtree_->Branch("timeHigh",&timeHigh_,"timeHigh/i");
00128 fragtree_->Branch("timeLow",&timeLow_,"timeLow/i");
00129
00130 fragtree_->Branch("frag_id",&frag_id_,"frag_id/i");
00131 fragtree_->Branch("seq_id",&seq_id_,"seq_id/i");
00132 fragtree_->Branch("frag_type",&frag_type_,"frag_type/i");
00133 fragtree_->Branch("frag_data_size",&frag_data_size_,"frag_data_size/i");
00134 fragtree_->Branch("contained_frags",&contained_frags_,"contained_frags/I");
00135 fragtree_->Branch("missing_data",&missing_data_,"missing_data/I");
00136 }
00137 }
00138
00139
00140 void artdaq::MissingDataCheck::analyze(art::Event const& e)
00141 {
00142 run_ = e.run();
00143 subrun_ = e.subRun();
00144 event_ = e.event();
00145 timeHigh_ = e.time().timeHigh();
00146 timeLow_ = e.time().timeLow();
00147
00148
00149 if(verbosity_>2)
00150 std::cout << "Processing:"
00151 << " Run " << e.run()
00152 << ", Subrun " << e.subRun()
00153 << ", Event " << e.event()
00154 << " (Time=" << e.time().timeHigh() << " " << e.time().timeLow() << ")"
00155 << std::endl;
00156
00157
00158 std::vector< art::Handle< std::vector<artdaq::Fragment> > > fragmentHandles;
00159 e.getManyByType(fragmentHandles);
00160
00161
00162 if(verbosity_>2){
00163 std::cout << "\tFound " << fragmentHandles.size() << " fragment collections." << std::endl;
00164 if(verbosity_>2){
00165 for(auto const& h : fragmentHandles)
00166 std::cout << "\t\tCollection " << h.provenance()->productInstanceName()
00167 << ":\t" << h->size() << " fragments." << std::endl;
00168 }
00169 }
00170
00171
00172 total_n_frags_=0;
00173 total_data_size_=0;
00174 for(auto const& h : fragmentHandles){
00175 total_n_frags_ += h->size();
00176 for(auto const& f : *h)
00177 total_data_size_ += f.dataSizeBytes();
00178 }
00179
00180
00181 if(expected_n_fragments_==-1)
00182 expected_n_fragments_=total_n_frags_;
00183
00184 if(verbosity_>2){
00185 std::cout << "\tTotal fragments = " << total_n_frags_
00186 << " / " << expected_n_fragments_ << std::endl;
00187 std::cout << "\tTotal data size in fragments = " << total_data_size_ << std::endl;
00188
00189 }
00190
00191 total_n_CFs_=0;
00192 total_n_CFs_missing_=0;
00193 total_n_frags_in_CFs_=0;
00194 total_n_frags_broken_=0;
00195 for(auto const& h : fragmentHandles){
00196 std::string instance_name = h.provenance()->productInstanceName();
00197
00198 if(instance_name.compare(0,9,"Container")==0){
00199 total_n_CFs_ += h->size();
00200
00201 for(auto const& f : *h){
00202 artdaq::ContainerFragment cf(f);
00203 if(cf.missing_data()){
00204 ++total_n_CFs_missing_;
00205 }
00206 total_n_frags_in_CFs_ += cf.block_count();
00207 frag_id_ = f.fragmentID();
00208 seq_id_ = f.sequenceID();
00209 frag_type_ = f.type();
00210 contained_frags_ = cf.block_count();
00211 missing_data_ = cf.missing_data()?1:0;
00212 frag_data_size_ = f.dataSizeBytes();
00213 fragtree_->Fill();
00214 }
00215
00216 }
00217
00218 else{
00219 for(auto const& f : *h){
00220
00221 if(instance_name.compare(0,4,"Data")==0 ||
00222 instance_name.compare(0,5,"Error")==0 ||
00223 instance_name.compare(0,6,"Broken")==0)
00224 ++total_n_frags_broken_;
00225
00226 frag_id_ = f.fragmentID();
00227 seq_id_ = f.sequenceID();
00228 frag_type_ = f.type();
00229 contained_frags_ = -1;
00230 missing_data_ = -1;
00231 frag_data_size_ = f.dataSizeBytes();
00232 fragtree_->Fill();
00233 }
00234
00235 }
00236
00237 }
00238
00239 if(verbosity_>2){
00240 std::cout << "\tTotal container fragments = " << total_n_CFs_
00241 << " of which " << total_n_CFs_missing_ << " have missing data."
00242 << std::endl;
00243 std::cout << "\tTotal fragments in containers = " << total_n_frags_in_CFs_ << std::endl;
00244 }
00245
00246 if(make_tree_){
00247 evtree_->Fill();
00248 }
00249
00250 }
00251
00252 void artdaq::MissingDataCheck::endJob()
00253 {
00254 if(verbosity_>0){
00255 std::cout << "---------------------------------------------------" << std::endl;
00256 std::cout << "----------- MISSING DATA CHECK SUMMARY ------------" << std::endl;
00257 std::cout << "---------------------------------------------------" << std::endl;
00258
00259 std::cout << "Total events processed: " << evtree_->GetEntries()
00260 << ", expected fragments: " << expected_n_fragments_ << std::endl;
00261
00262 unsigned int run; evtree_->SetBranchAddress("run",&run);
00263 unsigned int subrun; evtree_->SetBranchAddress("subrun",&subrun);
00264 unsigned int event; evtree_->SetBranchAddress("event",&event);
00265 int n_frag_exp; evtree_->SetBranchAddress("n_frag_exp",&n_frag_exp);
00266 unsigned int n_frag; evtree_->SetBranchAddress("n_frag",&n_frag);
00267 unsigned int n_frag_broken; evtree_->SetBranchAddress("n_frag_broken",&n_frag_broken);
00268 unsigned int n_miss_data; evtree_->SetBranchAddress("n_miss_data",&n_miss_data);
00269
00270 std::cout << "Events missing fragments:\t\t"
00271 << evtree_->GetEntries("n_frag<n_frag_exp")
00272 << " / " << evtree_->GetEntries() << std::endl;
00273 if(verbosity_>1){
00274 if(evtree_->GetEntries("n_frag<n_frag_exp")>0){
00275 for(int i=0; i<evtree_->GetEntries(); ++i)
00276 {
00277 evtree_->GetEvent(i);
00278 if((int)n_frag<n_frag_exp){
00279 std::cout << "\tEvent (" << run << "," << subrun << "," << event << ")"
00280 << " is missing " << n_frag_exp-n_frag << " fragments."
00281 << std::endl;
00282 }
00283 }
00284 }
00285 }
00286 std::cout << "Events with broken fragments:\t\t"
00287 << evtree_->GetEntries("n_frag_broken>0")
00288 << " / " << evtree_->GetEntries() << std::endl;
00289 if(evtree_->GetEntries("n_frag_broken>0")>0){
00290 for(int i=0; i<evtree_->GetEntries(); ++i)
00291 {
00292 evtree_->GetEvent(i);
00293 if(n_frag_broken > 0){
00294 std::cout << "\tEvent (" << run << "," << subrun << "," << event << ")"
00295 << " has " << n_frag_broken << " fragments."
00296 << std::endl;
00297 }
00298 }
00299 }
00300 std::cout << "Events missing data in fragments:\t"
00301 << evtree_->GetEntries("n_miss_data>0")
00302 << " / " << evtree_->GetEntries() << std::endl;
00303 if(verbosity_>1){
00304 if(evtree_->GetEntries("n_miss_data>0")>0){
00305 for(int i=0; i<evtree_->GetEntries(); ++i)
00306 {
00307 evtree_->GetEvent(i);
00308 if(n_miss_data>0){
00309 std::cout << "\tEvent (" << run << "," << subrun << "," << event << ")"
00310 << " has " << n_miss_data << " fragments missing data."
00311 << std::endl;
00312 }
00313 }
00314 }
00315 }
00316
00317 std::cout << "---------------------------------------------------" << std::endl;
00318 std::cout << "---------------------------------------------------" << std::endl;
00319 }
00320 }
00321
00322 DEFINE_ART_MODULE(artdaq::MissingDataCheck)