00001 #include "art/Framework/Core/EDAnalyzer.h"
00002 #include "art/Framework/Principal/Handle.h"
00003 #include "art/Framework/Principal/Event.h"
00004 #include "art/Framework/Principal/Run.h"
00005 #include "art/Framework/Core/ModuleMacros.h"
00006 #include "canvas/Utilities/InputTag.h"
00007
00008 #include "artdaq-core/Data/Fragment.hh"
00009 #include "artdaq-core/Data/ContainerFragment.hh"
00010
00011 #include "artdaq-core-demo/Overlays/FragmentType.hh"
00012 #include "artdaq-core-demo/Overlays/ToyFragment.hh"
00013
00014 #include "cetlib_except/exception.h"
00015
00016 #include <TFile.h>
00017 #include <TRootCanvas.h>
00018 #include <TCanvas.h>
00019 #include <TGraph.h>
00020 #include <TAxis.h>
00021 #include <TH1D.h>
00022 #include <TStyle.h>
00023
00024 #include <numeric>
00025 #include <vector>
00026 #include <functional>
00027 #include <algorithm>
00028 #include <iostream>
00029 #include <sstream>
00030 #include <initializer_list>
00031 #include <limits>
00032
00033 using std::cout;
00034 using std::cerr;
00035 using std::endl;
00036
00037 namespace demo
00038 {
00042 class WFViewer : public art::EDAnalyzer
00043 {
00044 public:
00061 explicit WFViewer(fhicl::ParameterSet const& p);
00062
00066 virtual ~WFViewer() = default;
00067
00072 void analyze(art::Event const& e) override;
00073
00078 void beginRun(art::Run const&) override;
00079
00080 private:
00081
00082 std::unique_ptr<TCanvas> canvas_[2];
00083 std::vector<Double_t> x_;
00084 int prescale_;
00085 bool digital_sum_only_;
00086 art::RunNumber_t current_run_;
00087
00088 std::size_t num_x_plots_;
00089 std::size_t num_y_plots_;
00090
00091 std::string raw_data_label_;
00092 std::vector<artdaq::Fragment::fragment_id_t> fragment_ids_;
00093
00094 std::vector<std::unique_ptr<TGraph>> graphs_;
00095 std::vector<std::unique_ptr<TH1D>> histograms_;
00096
00097 std::map<artdaq::Fragment::fragment_id_t, std::size_t> id_to_index_;
00098 std::string outputFileName_;
00099 TFile* fFile_;
00100 bool writeOutput_;
00101 };
00102 }
00103
00104 demo::WFViewer::WFViewer(fhicl::ParameterSet const& ps):
00105 art::EDAnalyzer(ps)
00106 , prescale_(ps.get<int>("prescale"))
00107 , digital_sum_only_(ps.get<bool>("digital_sum_only", false))
00108 , current_run_(0)
00109 , num_x_plots_(ps.get<std::size_t>("num_x_plots", std::numeric_limits<std::size_t>::max()))
00110 , num_y_plots_(ps.get<std::size_t>("num_y_plots", std::numeric_limits<std::size_t>::max()))
00111 , raw_data_label_(ps.get<std::string>("raw_data_label", "daq"))
00112 , fragment_ids_(ps.get<std::vector<artdaq::Fragment::fragment_id_t>>("fragment_ids"))
00113 , graphs_(fragment_ids_.size())
00114 , histograms_(fragment_ids_.size())
00115 , outputFileName_(ps.get<std::string>("fileName", "artdaqdemo_onmon.root"))
00116 , writeOutput_(ps.get<bool>("write_to_file", false))
00117 {
00118 if (num_x_plots_ == std::numeric_limits<std::size_t>::max() ||
00119 num_y_plots_ == std::numeric_limits<std::size_t>::max())
00120 {
00121 switch (fragment_ids_.size())
00122 {
00123 case 1: num_x_plots_ = num_y_plots_ = 1;
00124 break;
00125 case 2: num_x_plots_ = 2;
00126 num_y_plots_ = 1;
00127 break;
00128 case 3:
00129 case 4: num_x_plots_ = 2;
00130 num_y_plots_ = 2;
00131 break;
00132 case 5:
00133 case 6: num_x_plots_ = 3;
00134 num_y_plots_ = 2;
00135 break;
00136 case 7:
00137 case 8: num_x_plots_ = 4;
00138 num_y_plots_ = 2;
00139 break;
00140 default:
00141 num_x_plots_ = num_y_plots_ = static_cast<std::size_t>(ceil(sqrt(fragment_ids_.size())));
00142 }
00143 }
00144
00145
00146
00147
00148 for (std::size_t i_f = 0; i_f < fragment_ids_.size(); ++i_f)
00149 {
00150 id_to_index_[fragment_ids_[i_f]] = i_f;
00151 }
00152
00153 gStyle->SetOptStat("irm");
00154 gStyle->SetMarkerStyle(22);
00155 gStyle->SetMarkerColor(4);
00156 }
00157
00158 void demo::WFViewer::analyze(art::Event const& e)
00159 {
00160 static std::size_t evt_cntr = -1;
00161 evt_cntr++;
00162
00163
00164
00165
00166
00167
00168 artdaq::Fragments fragments;
00169 artdaq::FragmentPtrs containerFragments;
00170 std::vector<std::string> fragment_type_labels{"TOY1","TOY2","ContainerTOY1", "ContainerTOY2"};
00171
00172 for (auto label: fragment_type_labels)
00173 {
00174 art::Handle<artdaq::Fragments> fragments_with_label;
00175 e.getByLabel("daq", label, fragments_with_label);
00176
00177 if (!fragments_with_label.isValid()) continue;
00178
00179
00180
00181
00182 if (label == "Container" || label == "ContainerTOY1" || label == "ContainerTOY2")
00183 {
00184 for (auto cont : *fragments_with_label)
00185 {
00186 artdaq::ContainerFragment contf(cont);
00187 for (size_t ii = 0; ii < contf.block_count(); ++ii)
00188 {
00189 containerFragments.push_back(contf[ii]);
00190 fragments.push_back(*containerFragments.back());
00191 }
00192 }
00193 }
00194 else
00195 {
00196 for (auto frag : *fragments_with_label)
00197 {
00198 fragments.emplace_back(frag);
00199 }
00200 }
00201 }
00202
00203
00204
00205
00206
00207
00208
00209
00210
00211
00212
00213
00214
00215
00216
00217 artdaq::Fragment::sequence_id_t expected_sequence_id = std::numeric_limits<artdaq::Fragment::sequence_id_t>::max();
00218
00219
00220 for (const auto& frag : fragments)
00221 {
00222
00223
00224
00225 std::unique_ptr<ToyFragment> toyPtr;
00226
00227
00228
00229
00230 if (expected_sequence_id == std::numeric_limits<artdaq::Fragment::sequence_id_t>::max())
00231 {
00232 expected_sequence_id = frag.sequenceID();
00233 }
00234
00235 if (expected_sequence_id != frag.sequenceID())
00236 {
00237 cerr << "Warning in WFViewer: expected fragment with sequence ID " << expected_sequence_id << ", received one with sequence ID " << frag.sequenceID() << endl;
00238 }
00239
00240 FragmentType fragtype = static_cast<FragmentType>(frag.type());
00241 std::size_t max_adc_count = std::numeric_limits<std::size_t>::max();
00242 std::size_t total_adc_values = std::numeric_limits<std::size_t>::max();
00243
00244 switch (fragtype)
00245 {
00246 case FragmentType::TOY1:
00247 toyPtr.reset(new ToyFragment(frag));
00248 total_adc_values = toyPtr->total_adc_values();
00249 max_adc_count = pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1;
00250 break;
00251 case FragmentType::TOY2:
00252 toyPtr.reset(new ToyFragment(frag));
00253 total_adc_values = toyPtr->total_adc_values();
00254 max_adc_count = pow(2, frag.template metadata<ToyFragment::Metadata>()->num_adc_bits) - 1;
00255 break;
00256 default:
00257 throw cet::exception("Error in WFViewer: unknown fragment type supplied");
00258 }
00259
00260 artdaq::Fragment::fragment_id_t fragment_id = frag.fragmentID();
00261 std::size_t ind = id_to_index_[fragment_id];
00262
00263
00264
00265
00266 if (!histograms_[ind])
00267 {
00268 histograms_[ind] = std::unique_ptr<TH1D>(new TH1D(Form("Fragment_%d_hist", fragment_id), "", max_adc_count + 1, -0.5, max_adc_count + 0.5));
00269
00270 histograms_[ind]->SetTitle(Form("Frag %d, Type %s", fragment_id,
00271 fragmentTypeToString(fragtype).c_str()));
00272 histograms_[ind]->GetXaxis()->SetTitle("ADC value");
00273 }
00274
00275
00276
00277
00278
00279
00280 switch (fragtype)
00281 {
00282 case FragmentType::TOY1:
00283 case FragmentType::TOY2:
00284 for (auto val = toyPtr->dataBeginADCs(); val != toyPtr->dataEndADCs(); ++val)
00285 histograms_[ind]->Fill(*val);
00286 break;
00287
00288 default:
00289 throw cet::exception("Error in WFViewer: unknown fragment type supplied");
00290 }
00291
00292 if (evt_cntr % prescale_ - 1 && prescale_ > 1)
00293 {
00294 continue;
00295 }
00296
00297
00298
00299
00300 if (!digital_sum_only_)
00301 {
00302
00303
00304 if (x_.size() != total_adc_values)
00305 {
00306 x_.resize(total_adc_values);
00307
00308 std::iota(x_.begin(), x_.end(), 0);
00309 }
00310
00311
00312
00313
00314 if (!graphs_[ind] || static_cast<std::size_t>(graphs_[ind]->GetN()) != total_adc_values)
00315 {
00316 graphs_[ind] = std::unique_ptr<TGraph>(new TGraph(total_adc_values));
00317 graphs_[ind]->SetName(Form("Fragment_%d_graph", fragment_id));
00318 graphs_[ind]->SetLineColor(4);
00319 std::copy(x_.begin(), x_.end(), graphs_[ind]->GetX());
00320 }
00321
00322
00323
00324
00325
00326 switch (fragtype)
00327 {
00328 case FragmentType::TOY1:
00329 case FragmentType::TOY2:
00330 {
00331 std::copy(toyPtr->dataBeginADCs(), toyPtr->dataBeginADCs() + total_adc_values, graphs_[ind]->GetY());
00332 }
00333 break;
00334
00335 default:
00336 throw cet::exception("Error in WFViewer: unknown fragment type supplied");
00337 }
00338
00339
00340
00341
00342 canvas_[1]->cd(ind + 1);
00343 TVirtualPad* pad = static_cast<TVirtualPad*>(canvas_[1]->GetPad(ind + 1));
00344
00345 Double_t lo_x, hi_x, lo_y, hi_y, dummy;
00346
00347 graphs_[ind]->GetPoint(0, lo_x, dummy);
00348 graphs_[ind]->GetPoint(graphs_[ind]->GetN() - 1, hi_x, dummy);
00349
00350 lo_x -= 0.5;
00351 hi_x += 0.5;
00352
00353 lo_y = -0.5;
00354 hi_y = max_adc_count + 0.5;
00355
00356
00357 TH1F* padframe = static_cast<TH1F*>(pad->DrawFrame(lo_x, lo_y, hi_x, hi_y));
00358 padframe->SetTitle(Form("Frag %d, Type %s, SeqID %d", static_cast<int>(fragment_id),
00359 fragmentTypeToString(fragtype).c_str(),
00360 static_cast<int>(expected_sequence_id)));
00361 padframe->GetXaxis()->SetTitle("ADC #");
00362 pad->SetGrid();
00363 padframe->Draw("SAME");
00364 }
00365
00366
00367
00368 canvas_[0]->cd(ind + 1);
00369 histograms_[ind]->Draw();
00370
00371 canvas_[0]->Modified();
00372 canvas_[0]->Update();
00373
00374
00375
00376 if (!digital_sum_only_)
00377 {
00378 canvas_[1]->cd(ind + 1);
00379
00380 graphs_[ind]->Draw("PSAME");
00381
00382 canvas_[1]->Modified();
00383 canvas_[1]->Update();
00384 }
00385
00386 if (writeOutput_)
00387 {
00388 canvas_[0]->Write("wf0", TObject::kOverwrite);
00389 canvas_[1]->Write("wf1", TObject::kOverwrite);
00390 fFile_->Write();
00391 }
00392 }
00393 }
00394
00395 void demo::WFViewer::beginRun(art::Run const& e)
00396 {
00397 if (e.run() == current_run_) return;
00398 current_run_ = e.run();
00399
00400 if (writeOutput_)
00401 {
00402 fFile_ = new TFile(outputFileName_.c_str(), "RECREATE");
00403 fFile_->cd();
00404 }
00405
00406 for (int i = 0; i < 2; i++) canvas_[i] = 0;
00407 for (auto& x: graphs_) x = 0;
00408 for (auto& x: histograms_) x = 0;
00409
00410 for (int i = 0; (i < 2 && !digital_sum_only_) || i < 1; i++)
00411 {
00412 canvas_[i] = std::unique_ptr<TCanvas>(new TCanvas(Form("wf%d", i)));
00413 canvas_[i]->Divide(num_x_plots_, num_y_plots_);
00414 canvas_[i]->Update();
00415 ((TRootCanvas*)canvas_[i]->GetCanvasImp())->DontCallClose();
00416 }
00417
00418 canvas_[0]->SetTitle("ADC Value Distribution");
00419
00420 if (! digital_sum_only_)
00421 {
00422 canvas_[1]->SetTitle("ADC Values, Event Snapshot");
00423 }
00424
00425 if (writeOutput_)
00426 {
00427 canvas_[0]->Write();
00428 canvas_[1]->Write();
00429 }
00430 }
00431
00432
00433 DEFINE_ART_MODULE(demo::WFViewer)