#include #include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TTree.h" #include "JDAQ/JDAQEvaluator.hh" #include "km3net-dataformat/offline/Evt.hh" #include "JDetector/JDetector.hh" #include "JDetector/JDetectorToolkit.hh" #include "JROOT/JROOTClassSelector.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JTreeScanner.hh" #include "JSupport/JAutoTreeScanner.hh" #include "JSupport/JSupport.hh" #include "Jeep/JPrint.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" #include "JTrigger/JDAQHitToTSelector.hh" #include "JSupernova.hh" using namespace std; using namespace JPP; using namespace KM3NETDAQ; using namespace JSUPERNOVA; using namespace JROOT; /** * \authors selhedri,iagoos * \file * Application to analyse buffered L1 timeslices for the online CCSN analysis */ /** * Initialize the output ROOT file with CCSN candidate information * \param obs Single-DOM observables * \param frame_number Timeslice frame number * \param tr Output ROOT tree */ inline int initialize_root_tree(JCoincidenceSN &obs, int *frame_number, TTree *tr){ tr->Branch("dom_id", &obs.moduleID, "dom_id/I"); tr->Branch("vetoIndex", &obs.vetoIndex, "vetoIndex/I"); tr->Branch("event_time", &obs.time, "event_time/D"); tr->Branch("timeslice_time_s", &obs.timeslice_time, "timeslice_time_s/D"); tr->Branch("frame_number", frame_number, "frame_number/I"); tr->Branch("multiplicity", &(obs.multiplicity), "multiplicity/I"); return 0; } int main(int argc, char **argv) { typedef JRange JRange_t; typedef JRange JToTRange_t; JMultipleFileScanner<> inputFile; JLimit_t& numberOfEvents = inputFile.getLimit(); string outputFile; string detectorFile; double TMax_ns; double TVeto_ns; JToTRange_t totRange_ns; JDAQHitToTSelector totSelector_ns; JROOTClassSelector selector; int debug; int preTriggerThreshold; // Check input parameters try { JParser<> zap("Example application to study supernova detection background"); zap['f'] = make_field(inputFile); zap['o'] = make_field(outputFile) = "kaboom.root"; zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['a'] = make_field(detectorFile); zap['T'] = make_field(TMax_ns) = 10.0; zap['V'] = make_field(TVeto_ns) = 1000.0; zap['C'] = make_field(selector) = getROOTClassSelection(); zap['d'] = make_field(debug) = 1; zap['P'] = make_field(preTriggerThreshold) = 4; zap['t'] = make_field(totSelector_ns) = JDAQHitToTSelector(4, 255); zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } cout.tie(&cerr); JDetector detector; try { load(detectorFile, detector); } catch(const JException& error) { FATAL(error); } TTree *tr = new TTree("singledom", "SN coincidence info"); JCoincidenceSN SNdata(0,0,0); int frame_number; initialize_root_tree(SNdata, &frame_number, tr); // ----------------------------------- // STEP 0: prepare // ----------------------------------- // configure input streams JTreeScanner<>::debug = debug; JAutoTreeScanner zmap = JType(); JTreeScannerInterface* pts = zmap[selector]; pts->configure(inputFile); // detect input stream size int fEnd = pts->rbegin()->getFrameIndex(); int fStart = pts->begin( )->getFrameIndex(); // crop histogram if -n is specified if (fEnd > inputFile.getUpperLimit()) { fEnd = fStart + inputFile.getUpperLimit(); } int fLength = 1 + fEnd - fStart; NOTICE("begin | end | length = " << fStart << " | " << fEnd << " | " << fLength << endl); // detector routers const JModuleRouter moduleRouter(detector); const JDAQHitRouter hitRouter(detector); // ----------------------------------- // STEP 1: building vetoes from events // ----------------------------------- int runNumber = 0; TH1D* h_vtr = new TH1D("VetoTimeRange","VetoTimeRange", 10000, 0, 10000); JTreeScanner evIn(inputFile); map triggeredEvents; for (; evIn.hasNext(); ) { STATUS("event: " << setw(10) << evIn.getCounter() << '\r'); DEBUG(endl); JDAQEvent* event = evIn.next(); if (!runNumber) { runNumber = event->getRunNumber(); } JVeto vp(*event, hitRouter); triggeredEvents[event->getFrameIndex()].push_back(vp); h_vtr->Fill(vp.getLength()); } STATUS(triggeredEvents.size() << " JDAQEvents loaded in veto buffer." << endl); //----------------------------- // STEP 2: timeslice processing // ---------------------------- counter_type counter = 0; TFile *out = new TFile(outputFile.c_str(), "RECREATE"); tr->SetDirectory(out); for ( ; pts->hasNext() && counter != inputFile.getLimit(); ++counter) { STATUS("timeslice: " << setw(10) << counter << '\r'); DEBUG(endl); const JDAQTimeslice* timeslice = pts->next(); frame_number = timeslice->getFrameIndex(); if (counter == 0) { NOTICE("Start frame index = " << frame_number << " @ T " << timeslice->getTimesliceStart() << endl); } JVetoSet veto; if (triggeredEvents.count(frame_number)) { veto = triggeredEvents.at(frame_number); } // STANDARD PROCESSING JDataSN preTrigger(TMax_ns, preTriggerThreshold);// totSelector_ns); preTrigger(timeslice, moduleRouter, totSelector_ns); JTriggerSN trigger(TVeto_ns); trigger(preTrigger); // generate and configure event-based veto // count trigger JRange_t A = JRange_t(preTriggerThreshold, NUMBER_OF_PMTS); // select single-module clusters with max multiplicity in range A JSNFilterM trgA1(A, 1); // select non-vetoed clusters with max multiplicity in range A JSNFilterMV trgAV(A, veto); JDAQUTCExtended tstart = timeslice->getTimesliceStart(); double event_time = tstart.getTimeNanoSecond() * 1e-9; for(auto &p: trigger){ SNdata = p.getPeak(); SNdata.vetoIndex = trgA1(p) ? (trgAV(p) ? 0 : 2) : (trgAV(p) ? 1 : 3); // 0 if no veto, 1 if only DOM correlation veto, 2 if only muon veto, 3 if both SNdata.timeslice_time = event_time; tr->Fill(); } } //------------------------------------- // STEP 3: Writing data // ------------------------------------ tr->Write(); out->Close(); }