Jpp  17.1.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JTimesliceReprocessor.cc
Go to the documentation of this file.
1 
2 #include <string>
3 #include <iostream>
4 #include <iomanip>
5 
10 
11 #include "JDAQ/JDAQTimesliceIO.hh"
12 #include "JDAQ/JDAQEventIO.hh"
14 
15 #include "JDetector/JDetector.hh"
18 
19 #include "JTrigger/JHit.hh"
20 #include "JTrigger/JHitToolkit.hh"
23 #include "JTrigger/JTimeslice.hh"
24 #include "JTrigger/JBuildL1.hh"
25 #include "JTrigger/JBuildL2.hh"
26 #include "JTrigger/JTimesliceL1.hh"
29 
32 
36 #include "JSupport/JSupport.hh"
37 #include "JSupport/JMeta.hh"
39 
40 #include "Jeep/JParser.hh"
41 #include "Jeep/JMessage.hh"
42 
43 
44 /**
45  * \file
46  * Auxiliary program to re-process KM3NETDAQ::JDAQTimeslice data.
47  * \author mdejong
48  */
49 int main(int argc, char **argv)
50 {
51  using namespace std;
52  using namespace JPP;
53  using namespace KM3NETDAQ;
54 
55  typedef JAllTypes_t typelist;
56 
59  JLimit_t& numberOfEvents = inputFile.getLimit();
60  string detectorFile;
62  bool reuse_parameters;
63  JROOTClassSelector selector;
64  int debug;
65 
66  try {
67 
68  JParser<> zap("Auxiliary program to re-process time slice data.");
69 
70  zap['f'] = make_field(inputFile);
71  zap['o'] = make_field(outputFile) = "timeslice.root";
72  zap['n'] = make_field(numberOfEvents) = JLimit::max();
73  zap['a'] = make_field(detectorFile);
75  zap['U'] = make_field(reuse_parameters);
76  zap['C'] = make_field(selector, "timeslice selector, e.g. JDAQTimesliceL1.") = getROOTClassSelection<JDAQTimesliceTypes_t>();
77  zap['d'] = make_field(debug) = 0;
78 
79  zap(argc, argv);
80  }
81  catch(const exception &error) {
82  FATAL(error.what() << endl);
83  }
84 
85 
86  cout.tie(&cerr);
87 
89 
90  DEBUG("Frame time [ms] " << getFrameTime() * 1e-6 << endl);
91  DEBUG("Reset time [ms] " << getRTS() * 1e-6 << endl);
92  DEBUG("Trigger" << endl << parameters << endl);
93 
95 
96  try {
97  load(detectorFile, detector);
98  }
99  catch(const JException& error) {
100  FATAL(error);
101  }
102 
103  if (reuse_parameters) {
104 
105  try {
106 
107  parameters = getTriggerParameters(inputFile);
108 
109  NOTICE("Set trigger parameters from input." << endl);
110  }
111  catch(const JException& error) {
112  FATAL("No trigger parameters from input." << endl);
113  }
114  }
115  // detector
116 
117  if (parameters.disableHighRateVeto) {
118 
119  NOTICE("Disabling high-rate veto of all PMTs." << endl);
120 
121  detector.setPMTStatus(HIGH_RATE_VETO_DISABLE);
122  }
123 
125 
126  const JModuleRouter moduleRouter(detector);
127 
128  if (parameters.writeSummary()) { WARNING("Discard writeSummary option during reprocesing of data." << endl); }
129 
130  //typedef JHit hit_type;
131  //typedef int hit_type;
132  typedef double hit_type;
133 
134  typedef JSuperFrame1D<hit_type> JSuperFrame1D_t;
135  typedef JSuperFrame2D<hit_type> JSuperFrame2D_t;
136  typedef JTimeslice <hit_type> JTimeslice_t;
137  typedef JBuildL1 <hit_type> JBuildL1_t;
138  typedef JBuildL2 <hit_type> JBuildL2_t;
139 
140  const JBuildL1_t buildL1(parameters);
141  const JBuildL2_t buildL2(parameters.L2);
142  const JBuildL2_t buildSN(parameters.SN);
143 
144  JTimesliceRouter timesliceRouter(parameters.numberOfBins);
145 
146  outputFile.open();
147 
148  if (!outputFile.is_open()) {
149  FATAL("Error opening file " << outputFile << endl);
150  }
151 
152  outputFile.put(JMeta(argc, argv));
153  outputFile.put(parameters);
154 
156 
157  counter_type counter = 0;
158 
159  for ( ; in.hasNext() && counter != inputFile.getLimit(); ++counter) {
160 
161  STATUS("event: " << setw(10) << counter << '\r'); DEBUG(endl);
162 
163  const JDAQTimeslice* timeslice = in.next();
164 
165  DEBUG(*timeslice << endl);
166 
167  timesliceRouter.configure(*timeslice);
168 
169  JTimeslice_t timesliceL0(timeslice->getDAQChronometer());
170  JTimeslice_t timesliceL1(timeslice->getDAQChronometer());
171  JTimeslice_t timesliceL2(timeslice->getDAQChronometer());
172  JTimeslice_t timesliceSN(timeslice->getDAQChronometer());
173 
174  for (JDAQTimeslice::const_iterator super_frame = timeslice->begin(); super_frame != timeslice->end(); ++super_frame) {
175 
176  if (moduleRouter.hasModule(super_frame->getModuleID())) {
177 
178  // calibration
179 
180  const JModule& module = moduleRouter.getModule(super_frame->getModuleID());
181  JSuperFrame2D_t& buffer = JSuperFrame2D_t::demultiplex(*super_frame, module);
182 
183  // Apply high-rate veto
184 
185  buffer.applyHighRateVeto(parameters.highRateVeto_Hz);
186 
187  // L0
188 
189  timesliceL0.push_back(JSuperFrame1D_t(buffer));
190 
191  // L1
192 
193  timesliceL1.push_back(JSuperFrame1D_t(super_frame->getDAQChronometer(),
194  super_frame->getModuleIdentifier(),
195  module.getPosition()));
196 
197  buildL1(*timesliceL0.rbegin(), back_inserter(*timesliceL1.rbegin()));
198 
199  // L2
200 
201  timesliceL2.push_back(JSuperFrame1D_t(super_frame->getDAQChronometer(),
202  super_frame->getModuleIdentifier(),
203  module.getPosition()));
204 
205  buildL2(buffer, *timesliceL1.rbegin(), back_inserter(*timesliceL2.rbegin()));
206 
207  // SN
208 
209  timesliceSN.push_back(JSuperFrame1D_t(super_frame->getDAQChronometer(),
210  super_frame->getModuleIdentifier(),
211  module.getPosition()));
212 
213  buildSN(buffer, *timesliceL1.rbegin(), back_inserter(*timesliceSN.rbegin()));
214  }
215  }
216 
217  if (parameters.writeL0()) {
218  outputFile.put(*timeslice);
219  }
220 
221  if (parameters.writeL1()) {
222  outputFile.put(JTimesliceL1<JDAQTimesliceL1>(timesliceL1, timesliceRouter, moduleRouter, parameters.TMaxLocal_ns));
223  }
224 
225  if (parameters.writeL2()) {
226  outputFile.put(JTimesliceL1<JDAQTimesliceL2>(timesliceL2, timesliceRouter, moduleRouter, parameters.L2.TMaxLocal_ns));
227  }
228 
229  if (parameters.writeSN()) {
230  outputFile.put(JTimesliceL1<JDAQTimesliceSN>(timesliceSN, timesliceRouter, moduleRouter, parameters.SN.TMaxLocal_ns));
231  }
232  }
233  STATUS(endl);
234 
236 
237  io >> outputFile;
238 
239  outputFile.close();
240 }
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:70
Object writing to file.
Utility class to parse command line options.
Definition: JParser.hh:1517
General exception.
Definition: JException.hh:23
#define WARNING(A)
Definition: JMessage.hh:65
debug
Definition: JMessage.hh:29
int main(int argc, char *argv[])
Definition: Main.cc:15
ROOT TTree parameter settings of various packages.
const JModule & getModule(const JObjectID &id) const
Get module parameters.
Data structure for a composite optical module.
Definition: JModule.hh:68
Router for fast addressing of hits in KM3NETDAQ::JDAQTimeslice data structure as a function of the op...
#define STATUS(A)
Definition: JMessage.hh:63
Detector data structure.
Definition: JDetector.hh:89
void configure(const JDAQTimeslice &timeslice)
Configure.
Auxiliary class to select ROOT class based on class name.
Recording of objects on file according a format that follows from the file name extension.
Router for direct addressing of module data in detector data structure.
static const int HIGH_RATE_VETO_DISABLE
Enable (disable) use of high-rate veto test if this status bit is 0 (1);.
Definition: pmt_status.hh:13
*fatal Wrong number of arguments esac JCookie sh typeset Z DETECTOR typeset Z SOURCE_RUN typeset Z TARGET_RUN set_variable PARAMETERS_FILE $WORKDIR parameters
Definition: diff-Tuna.sh:38
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:83
Long64_t counter_type
Type definition for counter.
Auxiliary class for multiplexing object iterators.
Basic data structure for time and time over threshold information of hit.
string outputFile
Data structure for detector geometry and calibration.
Tools for handling different hit types.
1-dimensional frame with time calibrated data from one optical module.
Type list.
Definition: JTypeList.hh:22
void setDAQLongprint(const bool option)
Set DAQ print option.
Definition: JDAQPrint.hh:28
Scanning of objects from a single file according a format that follows from the extension of each fil...
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Template L2 builder.
Definition: JBuildL2.hh:45
double getMaximalDistance(const JDetector &detector)
Get maximal distance between modules in detector.
Detector file.
Definition: JHead.hh:226
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1993
double getFrameTime()
Get frame time duration.
Definition: JDAQClock.hh:162
ROOT I/O of application specific meta data.
#define NOTICE(A)
Definition: JMessage.hh:64
const JDAQChronometer & getDAQChronometer() const
Get DAQ chronometer.
Data time slice.
virtual const pointer_type & next() override
Get next element.
virtual bool hasNext() override
Check availability of next element.
const JPosition3D & getPosition() const
Get position.
Definition: JPosition3D.hh:130
General purpose messaging.
Template L1 hit builder.
Definition: JBuildL1.hh:85
#define FATAL(A)
Definition: JMessage.hh:67
Direct access to module in detector data structure.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Auxiliary class to build JDAQTimeslice for L1 timeslice.
Definition: JTimesliceL1.hh:36
double getRTS()
Get TDC dynamic range.
Definition: JDAQClock.hh:173
Utility class to parse command line options.
bool hasModule(const JObjectID &id) const
Has module.
2-dimensional frame with time calibrated data from one optical module.
Object reading from a list of files.
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
do set_variable DETECTOR_TXT $WORKDIR detector
KM3NeT DAQ constants, bit handling, etc.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:46
JTriggerParameters getTriggerParameters(const JMultipleFileScanner_t &file_list)
Get trigger parameters.
int debug
debug level
Time slice with calibrated data.
Definition: JTimeslice.hh:26
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62