#include #include #include #include #include #include "km3net-dataformat/offline/Head.hh" #include "km3net-dataformat/offline/MultiHead.hh" #include "km3net-dataformat/offline/Evt.hh" #include "km3net-dataformat/definitions/fitparameters.hh" #include "JDAQ/JDAQEventIO.hh" #include "JDAQ/JDAQTimesliceIO.hh" #include "JDAQ/JDAQSummarysliceIO.hh" #include "JDetector/JDetector.hh" #include "JDetector/JDetectorToolkit.hh" #include "JDetector/JModuleRouter.hh" #include "JTrigger/JHit.hh" #include "JTrigger/JFrame.hh" #include "JTrigger/JTimeslice.hh" #include "JTrigger/JHitL0.hh" #include "JTrigger/JHitL1.hh" #include "JTrigger/JHitR1.hh" #include "JTrigger/JBuildL0.hh" #include "JTrigger/JTriggerParameters.hh" #include "JSupport/JSingleFileScanner.hh" #include "JSupport/JParallelFileScanner.hh" #include "JSupport/JFileRecorder.hh" #include "JSupport/JSupport.hh" #include "JSupport/JMeta.hh" #include "JFit/JLine1Z.hh" #include "JFit/JFitToolkit.hh" #include "JFit/JTimeRange.hh" #include "JReconstruction/JEvt.hh" #include "JReconstruction/JEvtToolkit.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \file * * Program to determine veto of muon trajectory. * \author mdejong */ int main(int argc, char **argv) { using namespace std; using namespace JPP; using namespace KM3NETDAQ; typedef vector JVeto_t; typedef JParallelFileScanner< JTypeList > JParallelFileScanner_t; typedef JParallelFileScanner_t::multi_pointer_type multi_pointer_type; typedef JTYPELIST::typelist typelist; JParallelFileScanner_t inputFile; JFileRecorder outputFile; JLimit_t& numberOfEvents = inputFile.getLimit(); string detectorFile; double roadWidth_m; double R_Hz; size_t numberOfPrefits; JVeto_t veto; bool reprocess; int debug; try { JParser<> zap("Program to determine veto of muon trajectory."); zap['f'] = make_field(inputFile); zap['o'] = make_field(outputFile) = "veto.root"; zap['a'] = make_field(detectorFile); zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['R'] = make_field(roadWidth_m) = numeric_limits::max(); zap['B'] = make_field(R_Hz) = 6.0e3; zap['N'] = make_field(numberOfPrefits) = 1; zap['V'] = make_field(veto); zap['r'] = make_field(reprocess); zap['d'] = make_field(debug) = 1; zap(argc, argv); } catch(const exception& error) { FATAL(error.what() << endl); } JDetector detector; try { load(detectorFile, detector); } catch(const JException& error) { FATAL(error); } const JModuleRouter router(detector); typedef vector JDataL0_t; const JBuildL0 buildL0; double VS = 0.0; // veto [s] for (vector::const_iterator vs = veto.begin(); vs != veto.end(); ++vs) { VS += vs->getLength() * 1.0e-9; } outputFile.open(); outputFile.put(JMeta(argc, argv)); while (inputFile.hasNext()) { STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl); multi_pointer_type ps = inputFile.next(); const JDAQEvent* tev = ps; const JEvt* in = ps; JEvt cp = *in; JEvt out; cp.select(numberOfPrefits, qualitySorter); JDataL0_t dataL0; buildL0(*tev, router, true, back_inserter(dataL0)); for (JEvt::iterator track = cp.begin(); track != cp.end(); ++track) { const JRotation3D R (getDirection(*track)); const JLine1Z tz(getPosition (*track).rotate(R), track->getT()); // hit selection based on start value multiset top; for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) { JHitL0 hit(*i); hit.rotate(R); bool is_veto = false; if (tz.getDistance(hit) <= roadWidth_m) { const double t1 = hit.getT() - tz.getT(hit); for (vector::const_iterator vs = veto.begin(); vs != veto.end() && !is_veto; ++vs) { is_veto = (*vs)(t1); } } if (is_veto) { top.insert(hit.getPMTIdentifier()); } } double npe = 0.0; int count = 0; for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) { JPosition3D pos(module->getPosition()); pos.rotate(R); if (tz.getDistance(pos) <= roadWidth_m) { for (unsigned int i = 0; i != module->size(); ++i) { const JDAQPMTIdentifier id(module->getID(), i); npe += R_Hz * VS; count += top.count(id); } } } out.push_back(JFit(*track).add(JMUONVETO)); out.rbegin()->setW(track->getW()); out.rbegin()->setW(JVETO_NPE, npe); out.rbegin()->setW(JVETO_NUMBER_OF_HITS, count); } // apply default sorter sort(out.begin(), out.end(), qualitySorter); outputFile.put(out); } STATUS(endl); JSingleFileScanner::typelist> io(inputFile); io >> outputFile; outputFile.close(); }