// c++ standard library #include #include #include "JSupport/JSupport.hh" #include "JSupport/JSupportToolkit.hh" #include "JSupport/JMeta.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JAutoTreeScanner.hh" #include "JROOT/JROOTClassSelector.hh" #include "JLang/JObjectMultiplexer.hh" #include "JROOT/JManager.hh" #include "Jeep/JParser.hh" #include "km3net-dataformat/online/JDAQPMTIdentifier.hh" #include "km3net-dataformat/online/JDAQ.hh" #include "JDAQ/JDAQTimesliceIO.hh" #include "TH1D.h" #include "TMath.h" using namespace KM3NETDAQ; using namespace JLANG; using namespace JPP; using namespace std; /* * Rebins a histogram with constant bin width in a linear scale, so that the bins will have constant width in log10 scale. * * \param h The histogram. */ inline void BinLogX (TH1* h){ TAxis *axis = h->GetXaxis(); int bins = axis->GetNbins(); Axis_t min = axis->GetXmin(); Axis_t max = axis->GetXmax(); Axis_t width = (max - min) / bins; Axis_t *new_bins = new Axis_t[bins + 1]; for (int i = 0; i <= bins; i++) { new_bins[i] = TMath::Power(10, min + i * width); } axis->Set(bins, new_bins); delete new_bins; } int main(int argc, char **argv) { JMultipleFileScanner inputFile; JLimit_t& nSlices = inputFile.getLimit(); JROOTClassSelector selector; string outputFile; string detectorFile; try { JParser<> zap; zap['f'] = make_field(inputFile , "Path to input file "); zap['o'] = make_field(outputFile , "Path to output file") = "out.root"; zap['C'] = make_field(selector , "timeslice selector, e.g. JDAQTimesliceL1.") = getROOTClassSelection(); zap['s'] = make_field(nSlices , "number of slices" ) = JLimit::max(); zap(argc,argv); } catch(const exception &error) { ERROR(error.what() << endl); } const double factor = 1./(1.0e-6 * getFrameTime()); // [kHz] TH1D* h1 = new TH1D("rate_%", NULL, 100, -2, 2); TH1D* h2 = new TH1D("time_%", NULL, 1, -0.5, 0.5); JAutoTreeScanner map = JType(); JTreeScannerInterface* s = map[selector]; s->configure(inputFile); if (s->hasNext()){ JFrameIndexRange range (s->begin()->getFrameIndex() , s->rbegin()->getFrameIndex()); h2->SetBins(range.second - range.first + 1, range.first - 0.5, range.second + 0.5); } BinLogX(h1); JManager r (h1); JManager t (h2); JObjectMultiplexer scanner(inputFile, selector); counter_type counter = 0; for ( ; scanner.hasNext() && counter != inputFile.getLimit(); ++counter) { JDAQTimeslice* slice = scanner.next(); for(JDAQTimeslice::const_iterator frame = slice->begin() ; frame != slice->end() ; ++frame) { vector pmt_hit_count (NUMBER_OF_PMTS , 0) ; for (JDAQSuperFrame::const_iterator hit = frame->begin() ; hit != frame->end() ; ++hit){ pmt_hit_count[hit->getPMT()]++; } for (int pmt = 0 ; pmt != NUMBER_OF_PMTS ; ++pmt) { r[JDAQPMTIdentifier(frame->getModuleID(),pmt)]->Fill(pmt_hit_count[pmt]*factor); t[JDAQPMTIdentifier(frame->getModuleID(),pmt)]->Fill(slice->getFrameIndex(), pmt_hit_count[pmt]*factor); } } } for (typename JManager < JDAQPMTIdentifier , TH1D >::const_iterator i = r.begin() ; i != r.end() ; ++i){ for(int j=1 ; j < i -> second -> GetNbinsX() ; ++j){ double width = i->second->GetXaxis()->GetBinWidth(j); i -> second->SetBinContent(j ,i->second->GetBinContent(j)/width); i -> second->SetBinError (j ,i->second->GetBinError (j)/width); } } TFile out(outputFile.c_str(), "recreate"); putObject(&out, JMeta(argc, argv)); r.Write(out); t.Write(out); }