Jpp
JEventTimesliceWriter.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <fstream>
4 
5 #include "TRandom3.h"
6 
7 #include "evt/Head.hh"
8 #include "evt/Evt.hh"
9 #include "evt/Hit.hh"
10 
11 #include "JAAnet/JAAnetToolkit.hh"
12 #include "JAAnet/JEvtEvaluator.hh"
13 
15 #include "JDetector/JDetector.hh"
21 
22 #include "JDAQ/JDAQSummaryslice.hh"
23 #include "JDAQ/JDAQTimeslice.hh"
26 
28 #include "JSupport/JMeta.hh"
31 #include "JSupport/JSupport.hh"
32 #include "JSupport/JTreeScanner.hh"
34 
35 #include "Jeep/JParser.hh"
36 #include "Jeep/JMessage.hh"
37 
38 
39 /**
40  * \file
41  * Auxiliary program to convert multiple Monte Carlo events to DAQ timeslices.
42  * Events are timed according to a given rate (if > 0), otherwise they are timed according to Evt::mc_t.
43  * \author mdejong, mlincett
44  */
45 int main(int argc, char **argv)
46 {
47  using namespace std;
48  using namespace JPP;
49  using namespace KM3NETDAQ;
50 
51 
52  JMultipleFileScanner<JAAnetTypes_t> inputFile;
53  JFileRecorder <JTYPELIST<JAAnetTypes_t, JDAQTimesliceL0, JMeta, JRootTypes_t>::typelist> outputFile;
54  JLimit_t& numberOfEvents = inputFile.getLimit();
55  string detectorFile;
56  int run;
57  JPMTParametersMap pmtParameters;
58  JK40Rates rates_Hz;
59  double eventRate_Hz;
60  UInt_t seed;
61  int debug;
62 
63  try {
64 
65  JParser<> zap("Auxiliary program to convert multiple Monte Carlo events to time slices.");
66 
67  zap['f'] = make_field(inputFile);
68  zap['o'] = make_field(outputFile) = "timeslice.root";
69  zap['n'] = make_field(numberOfEvents) = JLimit::max();
70  zap['a'] = make_field(detectorFile);
71  zap['R'] = make_field(run) = 1;
72  zap['P'] = make_field(pmtParameters) = JPARSER::initialised();
73  zap['B'] = make_field(rates_Hz) = JPARSER::initialised();
74  zap['r'] = make_field(eventRate_Hz);
75  zap['S'] = make_field(seed) = 0;
76  zap['d'] = make_field(debug) = 0;
77 
78  zap(argc, argv);
79  }
80  catch(const exception &error) {
81  FATAL(error.what() << endl);
82  }
83 
84  gRandom->SetSeed(seed);
85 
86  cout.tie(&cerr);
87 
89 
90  if (pmtParameters.getQE() != 1.0) {
91 
92  WARNING("Correct background rates with global efficiency " << pmtParameters.getQE() << endl);
93 
94  rates_Hz.correct(pmtParameters.getQE());
95  }
96 
97  // parameters
98 
99  DEBUG("PMT paramaters: " << endl << pmtParameters << endl);
100  DEBUG("K40 rates: " << endl << rates_Hz << endl);
101 
102  if (eventRate_Hz < 0.0) {
103  FATAL("Invalid event rate " << eventRate_Hz << "; consider using JRandomTimesliceWriter." << endl);
104  }
105 
106  // detector
107 
108  JDetector detector;
109 
110  try {
111  load(detectorFile, detector);
112  }
113  catch(const JException& error) {
114  FATAL(error);
115  }
116 
117  // header
118 
119  Head header;
120 
121  try {
122  header = getHeader(inputFile);
123  }
124  catch (const JLANG::JNullPointerException error) {
125  FATAL("MC header is invalid");
126  }
127 
128  double liveTimeMC = JHead(header).livetime.numberOfSeconds;
129 
130  if (liveTimeMC < 0) {
131  FATAL("MC live time is negative; input file may be corrupted.");
132  }
133 
134  // background simulator
135 
136  JPMTParametersMap::Throw(false);
137 
138  JDetectorSimulator simbad(detector);
139 
140  simbad.reset(new JPMTDefaultSimulator(pmtParameters, detector));
141  simbad.reset(new JK40DefaultSimulator(rates_Hz));
142  simbad.reset(new JCLBDefaultSimulator());
143 
144  // output file
145 
146  outputFile.open();
147 
148  if (!outputFile.is_open()) FATAL("Error opening file " << outputFile << endl);
149 
150  outputFile.put(*gRandom);
151  outputFile.put(JMeta(argc, argv));
152  outputFile.put(header);
153 
154  // input stream;
155 
156  bool absTime = false;
157 
158  JTreeScannerInterface<Evt>* scan;
159 
160  if (eventRate_Hz == 0.0 && liveTimeMC == 0.0) {
161 
162  NOTICE("Event will be timed according to absolute MC time." << endl);
163 
164  absTime = true;
165 
166  // sorted access;
167  scan = new JTreeScanner<Evt, JEvtEvaluator>(inputFile);
169 
170  } else {
171 
172  // unsorted access;
173  scan = new JTreeScanner<Evt>(inputFile);
175 
176  if (eventRate_Hz == 0.0 && liveTimeMC > 0.0) {
177  int nEntries = scan->getEntries();
178  eventRate_Hz = nEntries / liveTimeMC;
179  DEBUG(nEntries << " events to be written." << endl);
180  NOTICE("MC live time is " << liveTimeMC << " s" << endl);
181  NOTICE("Event rate set to " << eventRate_Hz << " Hz from MC header." << endl);
182  }
183  }
184 
185  // begin timeslice generation;
186 
187  // state variables;
188  bool pendingEvt = false;
189  bool pendingTimeslice = false;
190 
191  Evt* event = NULL; // pending event;
192 
193  JTimeRange timeRange;
194  double tDAQ = 0; // DAQ event time;
195  double tOff = 0; // hit production time offset;
196 
197  int frame_index = 0;
198  JDAQChronometer chronometer;
199  JRandomTimeslice timeslice;
200 
201  counter_type evtCount = 0;
202 
203  while (pendingEvt || scan->hasNext()) {
204 
205  if (!pendingTimeslice) {
206  // update background timeslice;
207  frame_index++;
208  chronometer = JDAQChronometer(detector.getID(), run, frame_index, JDAQUTCExtended(getTimeOfFrame(frame_index)));
209  timeslice = JRandomTimeslice(chronometer, simbad);
210  DEBUG("evt count: " << setw(10) << evtCount << endl);
211  STATUS("frame index: " << setw(10) << frame_index << " | evt count: " << setw(10) << evtCount << "\r"); DEBUG(endl);
212  pendingTimeslice = true;
213  }
214 
215  if (!pendingEvt) {
216  // get new event;
217  event = scan->next();
218  timeRange = getTimeRange(*event);
219  tOff = timeRange.is_valid() ? timeRange.getLowerLimit() : 0;
220  // update DAQ time according to absolute or random;
221  if (absTime)
222  tDAQ = getEvtValue(*event) + tOff;
223  else
224  tDAQ += gRandom->Exp(1.0e9 / eventRate_Hz);
225 
226  DEBUG("event time [s] " << setprecision(5) << tDAQ * 1.0e-9 << endl);
227  pendingEvt = true;
228  }
229 
230  // DEBUG("tDAQ : " << tDAQ << " - frame index: " << frame_index << " - DAQ frame index " << getFrameIndex(tDAQ) << endl);
231  if (getFrameIndex(tDAQ) == frame_index) {
232 
233  // process event;
234  if (timeRange.is_valid()) {
235  DEBUG(*event << endl);
236  // real start time of the hits will be always tDAQ;
237  event->mc_t = tDAQ - tOff;
238  timeslice.add(JEventTimeslice(chronometer, simbad, *event));
239  evtCount++;
240  }
241  // event processed;
242  outputFile.put(*event);
243  pendingEvt = false;
244  } else {
245  // event time is past end of timeslice;
246  DEBUG(timeslice << endl);
247  outputFile.put(timeslice);
248  pendingTimeslice = false;
249  }
250 
251  }
252 
253  if (pendingTimeslice) {
254  DEBUG(timeslice << endl);
255  outputFile.put(timeslice);
256  pendingTimeslice = false;
257  }
258 
259  STATUS(endl);
260 
261  outputFile.put(*gRandom);
262  outputFile.close();
263 
264  NOTICE(evtCount << " events written over " << frame_index << " timeslices. " << endl);
265 }
JPMTParametersMap.hh
JMeta.hh
JK40DefaultSimulator.hh
KM3NETDAQ::JEventTimeslice
Timeslice with Monte Carlo event.
Definition: JEventTimeslice.hh:33
JFileRecorder.hh
JTreeScannerInterface.hh
JAANET::getTimeRange
JTimeRange getTimeRange(const Evt &event)
Get time range (i.e.
Definition: JAAnetToolkit.hh:134
JMessage.hh
JAANET::JHead::livetime
JAANET::livetime livetime
Definition: JHead.hh:1089
KM3NETDAQ::JDAQUTCExtended
Data structure for UTC time.
Definition: JDAQUTCExtended.hh:27
JRandomTimeslice.hh
JPARSER::initialised
Empty structure for specification of parser element that is initialised (i.e.
Definition: JParser.hh:63
JDETECTOR::load
void load(const JString &file_name, JDetector &detector)
Load detector from input file.
Definition: JDetectorToolkit.hh:456
JLANG::JNullPointerException
Exception for null pointer operation.
Definition: JException.hh:216
JSUPPORT::JLimit_t
JLimit JLimit_t
Type definition of limit.
Definition: JLimit.hh:215
JPARSER::JParser
Utility class to parse command line options.
Definition: JParser.hh:1493
NOTICE
#define NOTICE(A)
Definition: JMessage.hh:64
JDAQTimeslice.hh
JTreeScanner.hh
JDAQSummaryslice.hh
JAANET::JHead
Monte Carlo run header.
Definition: JHead.hh:839
JTOOLS::JTimeRange
JRange< double > JTimeRange
Type definition for time range.
Definition: JTools/JTimeRange.hh:19
JSupport.hh
JPP
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Definition: JAAnetToolkit.hh:37
JSUPPORT::getHeader
Head getHeader(const JMultipleFileScanner_t &file_list)
Get Monte Carlo header.
Definition: JMonteCarloFileSupportkit.hh:425
JDetectorSimulator.hh
WARNING
#define WARNING(A)
Definition: JMessage.hh:65
JCLBDefaultSimulator.hh
debug
int debug
debug level
Definition: JSirene.cc:59
JSUPPORT::JLimit::getLimit
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
JAAnetToolkit.hh
JAANET::livetime::numberOfSeconds
double numberOfSeconds
Live time [s].
Definition: JHead.hh:644
KM3NETDAQ::JRandomTimeslice
Timeslice with random data.
Definition: JRandomTimeslice.hh:29
JMultipleFileScanner.hh
STATUS
#define STATUS(A)
Definition: JMessage.hh:63
JParser.hh
JDetectorToolkit.hh
KM3NETDAQ::JDAQTimeslice::add
JDAQTimeslice & add(const JDAQTimeslice &timeslice)
Add another timeslice.
Definition: JDAQTimeslice.hh:134
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1954
JSUPPORT::JMeta
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:71
DEBUG
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
std
Definition: jaanetDictionary.h:36
JPMTDefaultSimulator.hh
KM3NETDAQ
KM3NeT DAQ data structures and auxiliaries.
Definition: DataQueue.cc:39
KM3NETDAQ::setDAQLongprint
void setDAQLongprint(const bool option)
Set DAQ print option.
Definition: JDAQPrint.hh:28
JEventTimeslice.hh
JMonteCarloFileSupportkit.hh
JAANET::getEvtValue
static const JEvtEvaluator getEvtValue
Function object for evaluation of DAQ objects.
Definition: JEvtEvaluator.hh:48
JDetector.hh
FATAL
#define FATAL(A)
Definition: JMessage.hh:67
outputFile
string outputFile
Definition: JDAQTimesliceSelector.cc:37
JROOT::counter_type
Long64_t counter_type
Type definition for counter.
Definition: JCounter.hh:24
JEEP::debug_t
debug
Definition: JMessage.hh:29
main
int main(int argc, char **argv)
Definition: JEventTimesliceWriter.cc:45
KM3NETDAQ::getTimeOfFrame
double getTimeOfFrame(const int frame_index)
Get start time of frame in ns since start of run for a given frame index.
Definition: JDAQClock.hh:185
KM3NETDAQ::JDAQChronometer
DAQ chronometer.
Definition: JDAQChronometer.hh:26
JEvtEvaluator.hh
KM3NETDAQ::getFrameIndex
int getFrameIndex(const double t_ns)
Get frame index for a given time in ns.
Definition: JDAQClock.hh:251