#include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TH2D.h" #include "km3net-dataformat/online/JDAQ.hh" #include "km3net-dataformat/online/JDAQHeader.hh" #include "JDAQ/JDAQTimesliceIO.hh" #include "JDetector/JDetector.hh" #include "JDetector/JDetectorToolkit.hh" #include "JDetector/JModuleRouter.hh" #include "JTrigger/JHitR0.hh" #include "JTrigger/JMatchL0.hh" #include "JTrigger/JSuperFrame1D.hh" #include "JTrigger/JSuperFrame2D.hh" #include "JTrigger/JPreprocessor.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JSupport.hh" #include "JSupport/JMeta.hh" #include "JROOT/JManager.hh" #include "JMath/JMathToolkit.hh" #include "JROOT/JROOTClassSelector.hh" #include "JROOT/JRootToolkit.hh" #include "JROOT/JRootFileWriter.hh" #include "JLang/JObjectMultiplexer.hh" #include "JTrigger/JTriggerParameters.hh" #include "JSupport/JTriggerParametersSupportkit.hh" #include "JCalibrate/JCalibrateToT.hh" #include "Jeep/JPrint.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \file * * Monitoring of PMT time-over-threshold distributions. * \author mkarel */ int main(int argc, char **argv) { using namespace std; using namespace JPP; using namespace KM3NETDAQ; JMultipleFileScanner inputFile; JLimit_t& numberOfEvents = inputFile.getLimit(); string outputFile; string detectorFile; JRange T_ns; double ctMax; JRange multiplicity; double deadTime_us; JROOTClassSelector selector; JPreprocessor option; int debug; try { JParser<> zap("Monitoring of PMT time-over-threshold distributions."); zap['f'] = make_field(inputFile, "input file."); zap['o'] = make_field(outputFile, "output file.") = "calibrate_tot.root"; zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['a'] = make_field(detectorFile, "detector file."); zap['T'] = make_field(T_ns, "time window [ns].") = JRange(10.0, 25.0); zap['c'] = make_field(ctMax, "maximal cosine space angle between PMT axes.") = 0.0; zap['M'] = make_field(multiplicity, "multiplicity range of hits on DOM.") = JRange(1, 2); zap['D'] = make_field(deadTime_us, "L0 dead time (us)") = 10; zap['C'] = make_field(selector, "timeslice selector, e.g. JDAQTimesliceL1.") = getROOTClassSelection(); zap['O'] = make_field(option, "hit pre-processing option.") = JPreprocessor::getOptions(); zap['d'] = make_field(debug, "debug flag.") = 1; zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } //----------------------------------------------------------- // check the input parameters //----------------------------------------------------------- if (!T_ns.is_valid()) { FATAL("Invalid time window [ns] " << T_ns << endl); } if (selector == JDAQTimeslice ::Class_Name() || selector == JDAQTimesliceL1::Class_Name()) { JTriggerParameters parameters; try { parameters = getTriggerParameters(inputFile); } catch(const JException& error) { FATAL("No trigger parameters from input." << endl); } if ((selector == JDAQTimeslice ::Class_Name() && parameters.writeL1.prescale > 0) || (selector == JDAQTimesliceL1::Class_Name())) { if (parameters.TMaxLocal_ns < T_ns.getUpperLimit()) { FATAL("Option -T = " << T_ns.getUpperLimit() << " is larger than in the trigger " << parameters.TMaxLocal_ns << endl); } } } if (!multiplicity.is_valid()) { FATAL("Invalid multiplicity " << multiplicity << endl); } if ( multiplicity.getLowerLimit() < 1) { FATAL("Invalid multiplicity " << multiplicity << endl); } //----------------------------------------------------------- // load detector file //----------------------------------------------------------- JDetector detector; try { load(detectorFile, detector); } catch(const JException& error) { FATAL(error); } if (detector.empty()) { FATAL("Empty detector." << endl); } const JModuleRouter router(detector); //----------------------------------------------------------- // initialise histograms //----------------------------------------------------------- const double zmin = -0.5; // [ns] const double zmax = 256.5; // [ns] const int nz = (int) ((zmax-zmin) / 1.0); JManager manager(new TH2D(MAKE_CSTRING("%" << _2SToT), NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS-0.5, nz, zmin, zmax)); typedef JHitR0 hit_type; typedef JSuperFrame2D JSuperFrame2D_t; typedef JSuperFrame1D JSuperFrame1D_t; const JMatchL0 match(T_ns.getUpperLimit()); // time window self-coincidences [ns] const double deadTime_ns = deadTime_us * 1e3; JObjectMultiplexer in(inputFile, selector); TFile out(outputFile.c_str(), "recreate"); putObject(out, JMeta(argc, argv)); JDAQHeader header; for (counter_type counter = 0; in.hasNext() && counter != inputFile.getLimit(); ++counter) { STATUS("event: " << setw(10) << counter << '\r'); DEBUG(endl); const JDAQTimeslice* timeslice = in.next(); if (header.getDetectorID() != timeslice->getDetectorID() || header.getRunNumber () != timeslice->getRunNumber ()) { header = timeslice->getDAQHeader(); putObject(out, header); } for (JDAQTimeslice::const_iterator frame = timeslice->begin(); frame != timeslice->end(); ++frame) { if (router.hasModule(frame->getModuleID())) { TH2D* h2 = manager[frame->getModuleID()]; const JModule& module = router.getModule(frame->getModuleID()); JSuperFrame2D_t& buffer = JSuperFrame2D_t::demultiplex(*frame, module); buffer.preprocess(option, match); JSuperFrame1D_t& data = JSuperFrame1D_t::multiplex(buffer); vector t0(NUMBER_OF_PMTS, numeric_limits::lowest()); for (JSuperFrame1D_t::const_iterator p = data.begin(); p != data.end(); ) { JSuperFrame1D_t::const_iterator q = p; // reject genuine coincidences double ct_max = -1.0; double dt_min = numeric_limits::max(); while (++q != data.end() && q->getT() - p->getT() < T_ns.getUpperLimit()) { const double ct = getDot(module.getPMT(p->getPMT()), module.getPMT(q->getPMT())); const double dt = q->getT() - p->getT(); if (ct > ct_max) { ct_max = ct; } if (dt < dt_min) { dt_min = dt; } } if (multiplicity(distance(p,q)) && ct_max < ctMax && dt_min > T_ns.getLowerLimit()) { for (JSuperFrame1D_t::const_iterator i = p; i != q; ++i) { if (i->getT() > t0[i->getPMT()] + deadTime_ns) { h2->Fill(i->getPMT(), i->getToT()); } t0[i->getPMT()] = i->getT(); } } p = q; } } } } STATUS(endl); out << manager; out.Close(); out.Write(); }