#include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TH1D.h" #include "JDAQ/JDAQEventIO.hh" #include "JDAQ/JDAQSummarysliceIO.hh" #include "JDAQ/JDAQEvaluator.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JTreeScanner.hh" #include "JSupport/JSupport.hh" #include "JROOT/JRootToolkit.hh" #include "JAcoustics/JEvent.hh" #include "JAcoustics/JSupport.hh" #include "JAcoustics/JAcousticsToolkit.hh" #include "JAcoustics/JEvtToolkit.hh" #include "Jeep/JProperties.hh" #include "Jeep/JPrint.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \file * * Example program to process DAQ and acoustic events. * \author mdejong */ int main(int argc, char **argv) { using namespace std; using namespace JPP; using namespace KM3NETDAQ; JMultipleFileScanner inputFile; JMultipleFileScanner<> acousticsFile; JLimit_t& numberOfEvents = inputFile.getLimit(); string outputFile; double Tmax_s = 300.0; // time window [s] size_t Nmin = 3; // minimum number of emitters int debug; try { JProperties properties; properties.insert(gmake_property(Tmax_s)); properties.insert(gmake_property(Nmin)); JParser<> zap("Example program to process DAQ and acoustic events."); zap['f'] = make_field(inputFile); zap['A'] = make_field(acousticsFile); zap['o'] = make_field(outputFile) = "example.root"; zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['@'] = make_field(properties) = JPARSER::initialised(); zap['d'] = make_field(debug) = 2; zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } TH1D h1("h1", NULL, 101, -0.5, +100.5); TH1D hn("hn", NULL, 101, -0.5, +100.5); TH1D ht("ht", NULL, 100, -1.0e3, +1.0e3); map data; // time-of-emission -> range of times-of-emission typedef JTreeScanner JTreeScanner_t; JTreeScanner_t in(acousticsFile); for (JTreeScanner_t::iterator p = in.begin(), q; p != in.end(); p = q) { STATUS("event " << FIXED(20,6) << p->begin()->getToE() << '\r'); DEBUG(endl); for (q = p; ++q != in.end() && q->begin()->getToE() <= p->rbegin()->getToE() + Tmax_s; ) {} h1.Fill((double) getNumberOfEmitters(p,q)); hn.Fill((double) distance(p,q)); if (getNumberOfEmitters(p,q) >= Nmin) { double range = 0.0; for (JTreeScanner_t::iterator i = p; i != q; ++i) { if (i->rbegin()->getToE() - p->begin()->getToE() > range) { range = i->rbegin()->getToE() - p->begin()->getToE(); } } data[p->begin()->getToE()] = range; } } STATUS(endl); while (inputFile.hasNext()) { const JDAQEvent* evt = inputFile.next(); const double t0 = getUNIXTime(*evt); STATUS("DAQ event " << FIXED(20,6) << t0 << '\r'); DEBUG(endl); map::const_iterator p = data.lower_bound(t0); double t1 = (p != data.end() ? p->first : 0.0); if (p != data.begin()) { --p; if (fabs(t1 - t0) > fabs(p->first - t0)) { t1 = p->first; } if (fabs(t1 - t0) > fabs(p->first + p->second - t0)) { t1 = p->first + p->second; } } ht.Fill(t0 - t1); } STATUS(endl); TFile out(outputFile.c_str(), "recreate"); out << h1 << hn << ht; out.Write(); out.Close(); }