#include #include #include #include "km3net-dataformat/offline/Head.hh" #include "km3net-dataformat/offline/MultiHead.hh" #include "km3net-dataformat/offline/Evt.hh" #include "km3net-dataformat/tools/time_converter.hh" #include "JAAnet/JHead.hh" #include "JAAnet/JHeadToolkit.hh" #include "JAAnet/JAAnetToolkit.hh" #include "JDAQ/JDAQEventIO.hh" #include "JDAQ/JDAQTimesliceIO.hh" #include "JDAQ/JDAQSummarysliceIO.hh" #include "JTrigger/JTriggerParameters.hh" #include "JSupport/JTriggeredFileScanner.hh" #include "JSupport/JFileRecorder.hh" #include "JSupport/JMonteCarloFileSupportkit.hh" #include "JSupport/JSupport.hh" #include "JSupport/JMeta.hh" #include "JPhysics/JGeanz.hh" #include "JReconstruction/JEvt.hh" #include "JReconstruction/JEvtToolkit.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \file * * Auxiliary program to store Monte Carlo information from a neutrino or the primary electron in JFIT::JEvt format. * By default the neutrino is considered, a bool variable allows to select the primary electron instead. * * \author mdejong, adomi */ int main(int argc, char **argv) { using namespace std; using namespace JPP; using namespace KM3NETDAQ; typedef JTYPELIST::typelist typelist; JTriggeredFileScanner<> inputFile; JFileRecorder outputFile; JLimit_t& numberOfEvents = inputFile.getLimit(); int debug; bool take_electron; try { JParser<> zap("Auxiliary program to store Monte Carlo true shower in format for subsequent fits."); zap['f'] = make_field(inputFile, "output of JTriggerEfficiency[.sh]"); zap['o'] = make_field(outputFile) = "mcevt.root"; zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['d'] = make_field(debug) = 2; zap['e'] = make_field(take_electron); zap(argc, argv); } catch(const exception& error) { FATAL(error.what() << endl); } const Vec offset = getOffset(getHeader(inputFile)); outputFile.open(); outputFile.put(JMeta(argc, argv)); while (inputFile.hasNext()) { STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl); JTriggeredFileScanner<>::multi_pointer_type ps = inputFile.next(); const JDAQEvent* tev = ps; const Evt* event = ps; JEvt out; const time_converter converter(*event, *tev); vector::const_iterator fermion; if(!take_electron) fermion = find_if(event->mc_trks.begin(), event->mc_trks.end(), is_neutrino); else fermion = find_if(event->mc_trks.begin(), event->mc_trks.end(), is_electron); if (fermion != event->mc_trks.end()) { // if fermion = neutrino => fermion position = neutrino interaction vertex double time = converter.putTime(); JVector3D position = getPosition(*fermion); // if fermion = electron => fermion position and time = EM shower brightest point if(take_electron){ JVector3D electron_dir = getDirection(*fermion); double shower_elongation = geanz.getMaximum(fermion->E); electron_dir *= shower_elongation; position.add(electron_dir); time += (shower_elongation * getInverseSpeedOfLight()); } JTrack3D td(position, getDirection(*fermion), time); JTrack3E ta(td, fermion->E); ta.add(getPosition(offset)); out.push_back(getFit(JMCEVT, ta, 0, 0.0, ta.getE())); } outputFile.put(out); } STATUS(endl); JMultipleFileScanner::typelist> io(inputFile); io >> outputFile; outputFile.close(); }