Jpp  19.1.0
the software that should make you happy
examples/JDynamics/JConvertEvt.cc
Go to the documentation of this file.
1 
2 #include <string>
3 #include <iostream>
4 #include <iomanip>
5 #include <algorithm>
6 #include <memory>
7 
8 #include "TROOT.h"
9 #include "TFile.h"
10 
14 
15 #include "JDetector/JDetector.hh"
17 #include "JDetector/JPMTRouter.hh"
21 
22 #include "JDynamics/JDynamics.hh"
23 
24 #include "JCompass/JSupport.hh"
25 
26 #include "JAcoustics/JSupport.hh"
27 
28 #include "JTrigger/JHit.hh"
29 #include "JTrigger/JHitToolkit.hh"
30 
31 #include "JAAnet/JHead.hh"
32 #include "JAAnet/JHeadToolkit.hh"
33 #include "JAAnet/JAAnetToolkit.hh"
34 
35 #include "JDAQ/JDAQEventIO.hh"
37 
39 #include "JSupport/JTreeScanner.hh"
43 #include "JSupport/JSupport.hh"
44 #include "JSupport/JMeta.hh"
45 
46 #include "Jeep/JeepToolkit.hh"
47 #include "Jeep/JParser.hh"
48 #include "Jeep/JMessage.hh"
49 
50 
51 /**
52  * \file
53  * Auxiliary program to convert trigger files to Evt format.
54  *
55  * \author mdejong, jseneca
56  */
57 int main(int argc, char **argv)
58 {
59  using namespace std;
60  using namespace JPP;
61  using namespace KM3NETDAQ;
62 
64  JACOUSTICS::JEvt>::typelist calibration_types;
65  typedef JMultipleFileScanner<calibration_types> JCalibration_t;
66 
67  JMultipleFileScanner<> inputFile;
69  JLimit_t& numberOfEvents = inputFile.getLimit();
70  string detectorFile;
71  JCalibration_t calibrationFile;
72  double Tmax_s;
73  JPMTParametersMap pmtParameters;
74  int debug;
75 
76  try {
77 
78  JParser<> zap("Auxiliary program to convert fit results to Evt format.\
79  \nThe option -L corresponds to the name of a shared library \
80  \nand function so to rearrange the order of fit results.");
81 
82  zap['f'] = make_field(inputFile);
83  zap['o'] = make_field(outputFile) = "aanet.root";
84  zap['n'] = make_field(numberOfEvents) = JLimit::max();
85  zap['a'] = make_field(detectorFile) = "";
86  zap['+'] = make_field(calibrationFile) = JPARSER::initialised();
87  zap['T'] = make_field(Tmax_s) = 100.0;
88  zap['P'] = make_field(pmtParameters, "PMT simulation data (or corresponding file name)") = JPARSER::initialised();
89  zap['d'] = make_field(debug) = 2;
90 
91  zap(argc, argv);
92  }
93  catch(const exception& error) {
94  FATAL(error.what() << endl);
95  }
96 
97 
99  typedef JSingleFileScanner_t::pointer_type pointer_type;
100 
101 
103 
104  if (detectorFile != "") {
105  try {
106  load(detectorFile, detector);
107  }
108  catch(const JException& error) {
109  FATAL(error);
110  }
111  }
112 
113  unique_ptr<JDynamics> dynamics;
114 
115  try {
116 
117  dynamics.reset(new JDynamics(detector, Tmax_s));
118 
119  dynamics->load(calibrationFile);
120  }
121  catch(const exception& error) {
122  if (!calibrationFile.empty()) {
123  FATAL(error.what());
124  }
125  }
126 
127 
128  const JPMTRouter pmt_router(dynamics ? dynamics->getDetector() : detector);
129  const JModuleRouter mod_router(dynamics ? dynamics->getDetector() : detector);
130 
131 
132  outputFile.open();
133 
134  Head header;
135 
136  try {
137  header = getHeader(inputFile);
138  } catch(const exception& error) {}
139  {
140  JAANET::JHead buffer(header);
141 
142  buffer.DAQ.livetime_s = getLivetime(inputFile.begin(), inputFile.end());
143  buffer.push(&JAANET::JHead::DAQ);
144 
145  copy(buffer, header);
146  }
147 
148  outputFile.put(header);
149 
150 
151  outputFile.put(JMeta(argc, argv));
152 
153 
154  counter_type number_of_events = 0;
155 
156  for (JMultipleFileScanner<>::const_iterator i = inputFile.begin(); i != inputFile.end(); ++i) {
157 
158  STATUS("Processing: " << *i << endl);
159 
160  JSingleFileScanner_t in(*i);
161  JTreeScanner<Evt> mc(*i);
162 
163  in.setLimit(inputFile.getLimit());
164 
165  Vec offset(0,0,0);
166 
167  int mc_run_id = 0;
168 
169  try {
170 
171  const JAANET::JHead head(getHeader(*i));
172 
173  offset = getOffset(head);
174  mc_run_id = head.start_run.run_id;
175 
176  } catch(const exception& error) {}
177 
178 
179  while (in.hasNext()) {
180 
181  STATUS("event: " << setw(10) << in.getCounter() << '\r'); DEBUG(endl);
182 
183  pointer_type ps = in.next();
184 
185  JDAQEvent* tev = ps;
186 
187  if (dynamics) {
188  dynamics->update(*tev); // update detector
189  }
190 
191  Evt out; // output event
192 
193  if (mc.getEntries() != 0) {
194 
195  Evt* event = mc.getEntry(tev->getCounter());
196 
197  if (event != NULL) {
198 
199  out = *event; // import Monte Carlo true data
200 
201  if (out.mc_run_id == 0) {
202  out.mc_run_id = mc_run_id;
203  }
204 
205  for (vector<Trk>::iterator i = out.mc_trks.begin(); i != out.mc_trks.end(); ++i) {
206  i->pos += offset; // offset true positions
207  }
208  }
209  }
210 
211  read(out, *tev); // import DAQ data
212 
213  if (!pmt_router->empty()) { // import calibration data
214 
215  for (vector<Hit>::iterator i = out.mc_hits.begin(); i != out.mc_hits.end(); ++i) {
216 
217  if (pmt_router.hasPMT(i->pmt_id)) {
218 
219  const JPMTAddress address = pmt_router.getAddress(i->pmt_id);
220  const JPMTIdentifier id = pmt_router.getIdentifier(address);
221  const JPMT& pmt = pmt_router.getPMT(address);
222 
223  i->dom_id = id.getID();
224  i->channel_id = id.getPMTAddress();
225  i->pos = getPosition (pmt);
226  i->dir = getDirection(pmt);
227 
228  } else {
229 
230  FATAL("Missing PMT" << i->pmt_id << endl);
231  }
232  }
233  }
234 
235  if (!mod_router->empty()) { // import calibration data
236 
237  for (vector<Hit>::iterator i = out.hits.begin(); i != out.hits.end(); ++i) {
238 
239  if (mod_router.hasModule(i->dom_id)) {
240 
241  const JPMTIdentifier id(i->dom_id, i->channel_id);
242 
243  const JPMTParameters& parameters = pmtParameters.getPMTParameters(id);
244  const JPMTAnalogueSignalProcessor cpu(parameters);
245 
246  const JPMT& pmt = mod_router.getPMT(id);
247 
249 
250  const JTRIGGER::JHit hit(getTime(i->tdc, pmt.getCalibration()), i->tot);
251 
252  i->pmt_id = pmt.getID();
253  i->pos = getPosition (pmt);
254  i->dir = getDirection(pmt);
255  i->t = hit.getT();
256  i->tot = hit.getToT();
257  i->a = cpu.getNPE(i->tot);
258 
259  } else {
260 
261  FATAL("Missing module " << i->dom_id << endl);
262  }
263  }
264  }
265 
266  out.id = ++number_of_events;
267 
268  outputFile.put(out);
269  }
270  }
271 
272  STATUS(endl);
273 
275 
276  io >> outputFile;
277 
278  outputFile.close();
279 }
Definition of hit and track types and auxiliary methods for handling Monte Carlo data.
ROOT TTree parameter settings.
ROOT TTree parameter settings.
string outputFile
Data structure for detector geometry and calibration.
Dynamic detector calibration.
Recording of objects on file according a format that follows from the file name extension.
Tools for handling different hit types.
General purpose messaging.
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
#define STATUS(A)
Definition: JMessage.hh:63
#define FATAL(A)
Definition: JMessage.hh:67
int debug
debug level
Definition: JSirene.cc:69
ROOT I/O of application specific meta data.
Direct access to module in detector data structure.
PMT analogue signal processor.
Direct access to PMT in detector data structure.
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:2142
Scanning of objects from a single file according a format that follows from the extension of each fil...
Support methods.
ROOT TTree parameter settings of various packages.
Basic data structure for time and time over threshold information of hit.
Auxiliary methods for handling file names, type names and environment.
Monte Carlo run header.
Definition: JHead.hh:1236
JAANET::start_run start_run
Definition: JHead.hh:1583
JAANET::DAQ DAQ
Definition: JHead.hh:1607
void push(T JHead::*pd)
Push given data member to Head.
Definition: JHead.hh:1374
const JCalibration & getCalibration() const
Get calibration.
Detector data structure.
Definition: JDetector.hh:96
Router for direct addressing of module data in detector data structure.
bool hasModule(const JObjectID &id) const
Has module.
const JPMT & getPMT(const JPMTIdentifier &id) const
Get PMT parameters.
Address of PMT in detector data structure.
Definition: JPMTAddress.hh:35
Auxiliary class for map of PMT parameters.
const JPMTParameters & getPMTParameters(const JPMTIdentifier &id) const
Get PMT parameters.
Data structure for PMT parameters.
bool slewing
time slewing of analogue signal
Router for direct addressing of PMT data in detector data structure.
Definition: JPMTRouter.hh:37
bool hasPMT(const JObjectID &id) const
Has PMT.
Definition: JPMTRouter.hh:116
const JPMT & getPMT(const JPMTAddress &address) const
Get PMT.
Definition: JPMTRouter.hh:92
JPMTIdentifier getIdentifier(const JPMTAddress &address) const
Get identifier of PMT.
Definition: JPMTRouter.hh:128
const JPMTAddress & getAddress(const JObjectID &id) const
Get address of PMT.
Definition: JPMTRouter.hh:80
Data structure for PMT geometry, calibration and status.
Definition: JPMT.hh:49
General exception.
Definition: JException.hh:24
int getID() const
Get identifier.
Definition: JObjectID.hh:50
Utility class to parse command line options.
Definition: JParser.hh:1698
Object writing to file.
General purpose class for object reading from a list of file names.
Object reading from a list of files.
Hit data structure.
static void setSlewing(const bool slewing)
Set slewing option.
double getToT() const
Get calibrated time over threshold of hit.
double getT() const
Get calibrated time of hit.
JTriggerCounter_t getCounter() const
Get trigger counter.
int main(int argc, char **argv)
JDirection3D getDirection(const Vec &dir)
Get direction.
double getTime(const Hit &hit)
Get true time of hit.
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:162
JPosition3D getPosition(const Vec &pos)
Get position.
Vec getOffset(const JHead &header)
Get offset.
std::istream & read(std::istream &in, JTestSummary &summary, const char delimiter=' ')
Read test summary.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Long64_t counter_type
Type definition for counter.
double getLivetime(const std::string &file_name)
Get data taking live time.
Head getHeader(const JMultipleFileScanner_t &file_list)
Get Monte Carlo header.
KM3NeT DAQ data structures and auxiliaries.
Definition: DataQueue.cc:39
Definition: JSTDTypes.hh:14
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:21
std::vector< Hit > hits
list of hits
Definition: Evt.hh:38
int mc_run_id
MC run identifier.
Definition: Evt.hh:27
std::vector< Hit > mc_hits
MC: list of MC truth hits.
Definition: Evt.hh:48
std::vector< Trk > mc_trks
MC: list of MC truth tracks.
Definition: Evt.hh:49
int id
offline event identifier
Definition: Evt.hh:22
The Head class reflects the header of Monte-Carlo event files, which consists of keys (also referred ...
Definition: Head.hh:65
double livetime_s
Live time [s].
Definition: JHead.hh:1053
Detector file.
Definition: JHead.hh:227
int run_id
MC run number.
Definition: JHead.hh:143
Acoustic event fit.
Orientation of module.
virtual double getNPE(const double tot_ns) const override
Get number of photo-electrons.
Dynamic detector calibration.
Definition: JDynamics.hh:84
Auxiliary class for recursive type list generation.
Definition: JTypeList.hh:351
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:68
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:45
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:72
Auxiliary base class for file name.
The Vec class is a straightforward 3-d vector, which also works in pyroot.
Definition: Vec.hh:13