#include #include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TH1D.h" #include "TF1.h" #include "JSupernova.hh" #include "km3net-dataformat/online/JDAQClock.hh" #include "JDAQ/JDAQEvaluator.hh" #include "JDetector/JDetector.hh" #include "JDetector/JDetectorToolkit.hh" #include "JROOT/JManager.hh" #include "JROOT/JROOTClassSelector.hh" #include "JTools/JQuantile.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JTreeScanner.hh" #include "JSupport/JAutoTreeScanner.hh" #include "JSupport/JSupport.hh" #include "JSupport/JMeta.hh" #include "JLang/JObjectMultiplexer.hh" #include "JMath/JMathToolkit.hh" #include "Jeep/JPrint.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" #include "JTrigger/JDAQHitToTSelector.hh" using namespace std; using namespace JPP; using namespace KM3NETDAQ; using namespace JSUPERNOVA; /** * \author mlincett * \file * Example application to test supernova trigger on run files. */ int main(int argc, char **argv) { typedef JRange JRange_t; JMultipleFileScanner<> inputFile; JLimit_t& numberOfEvents = inputFile.getLimit(); string outputFile; string detectorFile; double TMax_ns; double TVeto_ns; JROOTClassSelector selector; JDAQHitToTSelector totSelector_ns; int debug; int preTriggerThreshold; JRange_t M; // Check input parameters try { JParser<> zap("Example program test supernova triggers."); zap['f'] = make_field(inputFile); zap['o'] = make_field(outputFile) = "kexing.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['M'] = make_field(M) = JRange(6,10); zap['t'] = make_field(preTriggerThreshold) = 4; zap['t'] = make_field(totSelector_ns) = JDAQHitToTSelector(4, 255); zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } JDetector detector; try { load(detectorFile, detector); } catch(const JException& error) { FATAL(error); } // Configure input streams JTreeScanner<>::debug = debug; JAutoTreeScanner zmap = JType(); JTreeScannerInterface* pts = zmap[selector]; pts->configure(inputFile); // Configure routers const JModuleRouter moduleRouter(detector); const JDAQHitRouter hitRouter(detector); // ----------------------------------- // STEP 1: building vetoes from events // ----------------------------------- 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(); 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 // ------------------------------- // 0 -> no filter // 1 -> count once per track // 2 -> suppress track // 3 -> suppress based on trigger veto // Output data structures typedef JManager JManager_t; JManager_t SNT(new TH1D("SNT_F%", NULL, 100, 0.0, 100)); JManager_t MUL(new TH1D("MUL_F%", NULL, 1 + 31, -0.5, 31 + 1 - 0.5)); const int nStages = 5; vector > trgHistory(nStages); // Loop counter_type counter = 0; for ( ; pts->hasNext() && counter != inputFile.getLimit(); ++counter) { STATUS("timeslice: " << setw(10) << counter << '\r'); DEBUG(endl); const JDAQTimeslice* timeslice = pts->next(); int fIndex = timeslice->getFrameIndex(); JDataSN preTrigger(TMax_ns, preTriggerThreshold); preTrigger(timeslice, moduleRouter, totSelector_ns); JTriggerSN trigger(TVeto_ns); trigger(preTrigger); // generate and configure event-based veto JVetoSet veto; if (triggeredEvents.count(fIndex)) { veto = triggeredEvents.at(fIndex); } // trigger.setVeto(veto); // count trigger JSNFilterM trgF0(M, 0); JSNFilterM trgF1(M, 1); JSNFilterMV trgFV(M, veto); JRange A = JRange(2,31); JSNFilterMV trgAV(A, veto); int rawCount = count_if(preTrigger.begin(), preTrigger.end(), trgF0); int trgCountF0 = count_if(trigger.begin(), trigger.end(), trgF0); int trgCountF1 = count_if(trigger.begin(), trigger.end(), trgF1); int trgCountFV = count_if(trigger.begin(), trigger.end(), trgFV); int domCountF1 = trigger.getModules(trgF1).size(); trgHistory[0].push_back(rawCount); trgHistory[1].push_back(trgCountF0); trgHistory[2].push_back(trgCountF1); trgHistory[3].push_back(trgCountFV); trgHistory[4].push_back(domCountF1); trigger.fill(MUL[1], trgF0); trigger.fill(MUL[2], trgF1); trigger.fill(MUL[3], trgFV); trigger.fill(MUL[4], trgAV); } //------------------------------------- // STEP 3: generating trigger summaries // ------------------------------------ for (int i = 0; i < nStages; i++) { for (unsigned j = 0; j < trgHistory[i].size(); j++) { SNT[i]->Fill(trgHistory[i][j]); } } SNT[0]->SetTitle("M[6,10] count before clustering"); SNT[1]->SetTitle("M[6,10] count after clustering"); SNT[2]->SetTitle("M[6,10] count after track self-veto"); SNT[3]->SetTitle("M[6,10] count after track trigger-veto"); SNT[4]->SetTitle("M[6,10] count after track self-veto, unique modules"); if (outputFile != "") { TFile out(outputFile.c_str(), "RECREATE"); h_vtr->Write(); SNT.Write(out); MUL.Write(out); putObject(out, JMeta(argc,argv)); out.Close(); } }