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