Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
JConvertEvt.cc File Reference

Auxiliary program to convert fit results to AAnet format. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <algorithm>
#include "TROOT.h"
#include "TFile.h"
#include "km3net-dataformat/offline/Head.hh"
#include "km3net-dataformat/offline/Evt.hh"
#include "km3net-dataformat/offline/io_online.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JPMTRouter.hh"
#include "JDetector/JModuleRouter.hh"
#include "JDetector/JPMTParametersMap.hh"
#include "JDetector/JPMTAnalogueSignalProcessor.hh"
#include "JTrigger/JHit.hh"
#include "JTrigger/JHitToolkit.hh"
#include "JAAnet/JHead.hh"
#include "JAAnet/JHeadToolkit.hh"
#include "JAAnet/JAAnetToolkit.hh"
#include "JDAQ/JDAQEventIO.hh"
#include "JDAQ/JDAQSummarysliceIO.hh"
#include "JSupport/JParallelFileScanner.hh"
#include "JSupport/JTreeScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMonteCarloFileSupportkit.hh"
#include "JSupport/JSupportToolkit.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JMeta.hh"
#include "JReconstruction/JEvt.hh"
#include "JReconstruction/JEvtToolkit.hh"
#include "Jeep/JeepToolkit.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Detailed Description

Auxiliary program to convert fit results to AAnet format.

Author
mdejong

Definition in file JConvertEvt.cc.

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 53 of file JConvertEvt.cc.

54 {
55  using namespace std;
56  using namespace JPP;
57  using namespace KM3NETDAQ;
58 
59  JMultipleFileScanner<> inputFile;
61  JLimit_t& numberOfEvents = inputFile.getLimit();
62  string detectorFile;
63  JPMTParametersMap pmtParameters;
65  size_t numberOfFits;
66  int debug;
67 
68  try {
69 
70  JParser<> zap("Auxiliary program to convert fit results to AAnet format.\
71  \nThe option -L corresponds to the name of a shared library \
72  \nand function so to rearrange the order of fit results.");
73 
74  zap['f'] = make_field(inputFile);
75  zap['o'] = make_field(outputFile) = "aanet.root";
76  zap['n'] = make_field(numberOfEvents) = JLimit::max();
77  zap['a'] = make_field(detectorFile) = "";
78  zap['P'] = make_field(pmtParameters, "PMT simulation data (or corresponding file name)") = JPARSER::initialised();
80  zap['N'] = make_field(numberOfFits) = 0;
81  zap['d'] = make_field(debug) = 2;
82 
83  zap(argc, argv);
84  }
85  catch(const exception& error) {
86  FATAL(error.what() << endl);
87  }
88 
89 
90  typedef JParallelFileScanner< JTypeList<JDAQEvent, JEvt> > JParallelFileScanner_t;
91  typedef JParallelFileScanner_t::multi_pointer_type multi_pointer_type;
92 
93 
95 
96  if (detectorFile != "") {
97  try {
98  load(detectorFile, detector);
99  }
100  catch(const JException& error) {
101  FATAL(error);
102  }
103  }
104 
105  const JPMTRouter pmt_router(detector);
106  const JModuleRouter mod_router(detector);
107 
108 
109  outputFile.open();
110 
111  Head header;
112 
113  try {
114  header = getHeader(inputFile);
115  } catch(const exception& error) {}
116  {
117  JHead buffer(header);
118 
119  buffer.DAQ.livetime_s = getLivetime(inputFile.begin(), inputFile.end());
120  buffer.push(&JHead::DAQ);
121 
122  copy(buffer, header);
123  }
124 
125  outputFile.put(header);
126 
127 
128  outputFile.put(JMeta(argc, argv));
129 
130 
131  counter_type number_of_events = 0;
132 
133  for (JMultipleFileScanner<>::const_iterator i = inputFile.begin(); i != inputFile.end(); ++i) {
134 
135  STATUS("Processing: " << *i << endl);
136 
137  JParallelFileScanner_t in(*i);
138  JTreeScanner<Evt> mc(*i);
139 
140  in.setLimit(inputFile.getLimit());
141 
142  Vec center(0,0,0);
143 
144  int mc_run_id = 0;
145 
146  try {
147 
148  const JHead head(getHeader(*i));
149 
150  center = get<Vec>(head);
151  mc_run_id = head.start_run.run_id;
152 
153  } catch(const exception& error) {}
154 
155 
156  while (in.hasNext()) {
157 
158  STATUS("event: " << setw(10) << in.getCounter() << '\r'); DEBUG(endl);
159 
160  multi_pointer_type ps = in.next();
161 
162  JDAQEvent* tev = ps;
163  JEvt* evt = ps;
164 
165  JEvt::iterator __end = evt->end();
166 
167  if (numberOfFits > 0) {
168  advance(__end = evt->begin(), min(numberOfFits, evt->size()));
169  }
170 
171  if (qualitySorter.is_valid()) {
172  partial_sort(evt->begin(), __end, evt->end(), qualitySorter);
173  }
174 
175  Evt out; // output event
176 
177  if (mc.getEntries() != 0) {
178 
179  Evt* event = mc.getEntry(tev->getCounter());
180 
181  if (event != NULL) {
182 
183  out = *event; // import Monte Carlo true data
184 
185  if (out.mc_run_id == 0) {
186  out.mc_run_id = mc_run_id;
187  }
188 
189  for (vector<Trk>::iterator i = out.mc_trks.begin(); i != out.mc_trks.end(); ++i) {
190  i->pos += center; // offset true positions
191  }
192  }
193  }
194 
195  read(out, *tev); // import DAQ data
196 
197  if (!pmt_router->empty()) { // import calibration data
198 
199  for (vector<Hit>::iterator i = out.mc_hits.begin(); i != out.mc_hits.end(); ++i) {
200 
201  if (pmt_router.hasPMT(i->pmt_id)) {
202 
203  const JPMTAddress address = pmt_router.getAddress(i->pmt_id);
204  const JPMTIdentifier id = pmt_router.getIdentifier(address);
205  const JPMT& pmt = pmt_router.getPMT(address);
206 
207  i->dom_id = id.getID();
208  i->channel_id = id.getPMTAddress();
209  i->pos = getPosition (pmt);
210  i->dir = getDirection(pmt);
211 
212  } else {
213 
214  FATAL("Missing PMT" << i->pmt_id << endl);
215  }
216  }
217  }
218 
219  if (!mod_router->empty()) { // import calibration data
220 
221  for (vector<Hit>::iterator i = out.hits.begin(); i != out.hits.end(); ++i) {
222 
223  if (mod_router.hasModule(i->dom_id)) {
224 
225  const JPMTIdentifier id(i->dom_id, i->channel_id);
226  const JPMTAnalogueSignalProcessor cpu(pmtParameters.getPMTParameters(id));
227 
228  const JPMT& pmt = mod_router.getModule(i->dom_id).getPMT(i->channel_id);
229 
230  const JHit hit(getTime(i->tdc, pmt.getCalibration()), i->tot);
231 
232  i->pmt_id = pmt.getID();
233  i->pos = getPosition (pmt);
234  i->dir = getDirection(pmt);
235  i->t = hit.getT();
236  i->tot = hit.getToT();
237  i->a = cpu.getNPE(i->tot);
238 
239  } else {
240 
241  FATAL("Missing module " << i->dom_id << endl);
242  }
243  }
244  }
245 
246  for (JEvt::const_iterator fit = evt->begin(); fit != __end; ++fit) { // import reconstruction data
247 
248  Trk trk;
249 
250  trk.id = out.trks.size() + 1;
251  trk.pos = Vec(fit->getX(), fit->getY(), fit->getZ());
252  trk.dir = Vec(fit->getDX(), fit->getDY(), fit->getDZ());
253  trk.t = fit->getT();
254  trk.E = fit->getE();
255  trk.lik = fit->getQ();
257 
258  for (JHistory::const_iterator i = fit->getHistory().begin(); i != fit->getHistory().end(); ++i) {
259  trk.rec_stages.push_back(i->type);
260  }
261 
262  for (int i = 0; i != fit->getN(); ++i) {
263  trk.fitinf.push_back(fit->getW(i));
264  }
265 
266  out.trks.push_back(trk);
267  }
268 
269  out.id = ++number_of_events;
270 
271  outputFile.put(out);
272  }
273  }
274 
275  STATUS(endl);
276 
278 
279  io >> outputFile;
280 
281  outputFile.close();
282 }
Router for direct addressing of PMT data in detector data structure.
Definition: JPMTRouter.hh:33
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:1500
General exception.
Definition: JException.hh:23
double t
track time [ns] (when the particle is at pos )
Definition: Trk.hh:18
#define STATUS(A)
Definition: JMessage.hh:63
Detector data structure.
Definition: JDetector.hh:80
bool read(Vec &v, std::istream &is)
Read a Vec(tor) from a stream.
Definition: io_ascii.hh:141
Router for direct addressing of module data in detector data structure.
Vec dir
track direction
Definition: Trk.hh:17
General purpose sorter of fit results.
General purpose class for parallel reading of objects from a single file or multiple files...
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:66
Long64_t counter_type
Type definition for counter.
double getTime(const Hit &hit)
Get true time of hit.
string outputFile
double E
Energy [GeV] (either MC truth or reconstructed)
Definition: Trk.hh:19
Acoustics hit.
Head getHeader(const JMultipleFileScanner_t &file_list)
Get Monte Carlo header.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
The Vec class is a straightforward 3-d vector, which also works in pyroot.
Definition: Vec.hh:12
Detector file.
Definition: JHead.hh:196
JDirection3D getDirection(const Vec &dir)
Get direction.
Acoustic event fit.
int mc_run_id
MC run identifier.
Definition: Evt.hh:26
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
int getID() const
Get identifier.
Definition: JObjectID.hh:50
Data structure for PMT geometry and calibration.
Definition: JPMT.hh:47
JPosition3D getPosition(const Vec &pos)
Get position.
std::vector< double > fitinf
place to store additional fit info, for jgandalf, see JFitParameters.hh
Definition: Trk.hh:30
Auxiliary class for map of PMT parameters.
int debug
debug level
Definition: JSirene.cc:63
counter_type advance(counter_type &counter, const counter_type value, const counter_type limit=std::numeric_limits< counter_type >::max())
Advance counter.
Monte Carlo run header.
Definition: JHead.hh:1091
The Head class reflects the header of Monte-Carlo event files, which consists of keys (also referred ...
Definition: Head.hh:67
#define FATAL(A)
Definition: JMessage.hh:67
int id
track identifier
Definition: Trk.hh:15
Vec pos
postion of the track at time t [m]
Definition: Trk.hh:16
Address of PMT in detector data structure.
Definition: JPMTAddress.hh:32
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
std::vector< Trk > trks
list of reconstructed tracks (can be several because of prefits,showers, etc).
Definition: Evt.hh:36
double getLivetime(const std::string &file_name)
Get data taking live time.
std::vector< int > rec_stages
list of identifyers of succesfull fitting stages resulting in this track
Definition: Trk.hh:25
double lik
likelihood or lambda value (for aafit, lambda)
Definition: Trk.hh:22
std::vector< Hit > mc_hits
MC: list of MC truth hits.
Definition: Evt.hh:45
static const int JPP_RECONSTRUCTION_TYPE
KM3NeT Data Definitions v2.0.0 https://git.km3net.de/common/km3net-dataformat.
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:139
int id
offline event identifier
Definition: Evt.hh:21
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
std::vector< Hit > hits
list of hits
Definition: Evt.hh:35
do set_variable DETECTOR_TXT $WORKDIR detector
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 source JAcoustics sh $DETECTOR_ID typeset A TRIPODS get_tripods $WORKDIR tripod txt TRIPODS for EMITTER in
Definition: JCanberra.sh:36
The Trk class represents a Monte Carlo (MC) particle as well as a reconstructed track/shower.
Definition: Trk.hh:13
JTriggerCounter_t getCounter() const
Get trigger counter.
std::vector< Trk > mc_trks
MC: list of MC truth tracks.
Definition: Evt.hh:46
bool qualitySorter(const JRECONSTRUCTION::JFit &first, const JRECONSTRUCTION::JFit &second)
Comparison of fit results.
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:19
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
int rec_type
identifyer for the overall fitting algorithm/chain/strategy
Definition: Trk.hh:24