Jpp  18.3.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
JAAnet.cc File Reference

Example program to analyse track fit results from Evt formatted data. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <cmath>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "km3net-dataformat/offline/Head.hh"
#include "km3net-dataformat/offline/Evt.hh"
#include "km3net-dataformat/online/JDAQClock.hh"
#include "JAAnet/JAAnetToolkit.hh"
#include "JTools/JRange.hh"
#include "JMath/JMathToolkit.hh"
#include "JTrigger/JHitL0.hh"
#include "JPhysics/JConstants.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JSupport.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

Example program to analyse track fit results from Evt formatted data.

Author
mdejong

Definition in file JAAnet.cc.

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 58 of file JAAnet.cc.

59 {
60  using namespace std;
61  using namespace JPP;
62  using namespace KM3NETDAQ;
63 
65 
66  JMultipleFileScanner<Evt> inputFile;
67  JLimit_t& numberOfEvents = inputFile.getLimit();
68  string detectorFile;
69  string outputFile;
70  JRange_t E;
71  int debug;
72 
73  try {
74 
75  JParser<> zap("Example program to analyse track fit results from Evt formatted data.");
76 
77  zap['f'] = make_field(inputFile);
78  zap['n'] = make_field(numberOfEvents) = JLimit::max();
79  zap['a'] = make_field(detectorFile) = "";
80  zap['o'] = make_field(outputFile) = "histogram.root";
81  zap['E'] = make_field(E) = JRange_t();
82  zap['d'] = make_field(debug) = 2;
83 
84  zap(argc, argv);
85  }
86  catch(const exception& error) {
87  FATAL(error.what() << endl);
88  }
89 
90 
92 
93  if (detectorFile != "") {
94  try {
95  load(detectorFile, detector);
96  }
97  catch(const JException& error) {
98  FATAL(error);
99  }
100  }
101 
102  const JModuleRouter router(detector);
103 
104  TFile out(outputFile.c_str(), "recreate");
105 
106  TH1D hx("hx", NULL, 100, -3.0, +2.3); // [log(deg)]
107  TH1D hd("hd", NULL, 100, 0.0, 10.0); // [m]
108  TH1D ht("ht", NULL, 100, -100.0, 100.0); // [ns]
109  TH1D he("he", NULL, 100, -5.0, +5.0); // [log10(E)]
110  TH1D h1("h1", NULL, 100, -50.0, +50.0); // [ns]
111 
112 
113  while (inputFile.hasNext()) {
114 
115  STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl);
116 
117  const Evt* evt = inputFile.next();
118 
119  //if (has_muon(*evt) && has_reconstructed_track<JPP_RECONSTRUCTION_TYPE>(*evt, rec_stages_range(JMUONBEGIN, JMUONEND))) {
120  if (has_muon(*evt) && has_reconstructed_jppmuon(*evt)) {
121 
123 
124  JTrack3E ta = getTrack(get_muon(*evt));
125  //JTrack3E tb = getTrack(get_best_reconstructed_track<JPP_RECONSTRUCTION_TYPE>(*evt, rec_stages_range(JMUONBEGIN, JMUONEND)));
127 
128  ta.add(getTimeSinceRTS(evt->mc_t));
129 
132 
133  if (E(ta.getE())) {
134 
135  hx.Fill(log10(getAngle(ta.getDirection(), tb.getDirection())));
136  hd.Fill((ta.getPosition() - tb.getPosition()).getLength());
137  ht.Fill( ta.getT() - tb.getT());
138  he.Fill(log10(tb.getE()/ta.getE()));
139  }
140 
141  for (vector<Hit>::const_iterator i = evt->hits.begin(); i != evt->hits.end(); ++i) {
142 
143  if (router.hasModule(i->dom_id)) {
144 
145  const JHitL0 hit = getHit(*i, router);
146 
147  h1.Fill(hit.getT() - tb.getT(hit.getPosition()));
148  }
149  }
150  }
151  }
152  STATUS(endl);
153 
154  out.Write();
155  out.Close();
156 }
double getT() const
Get calibrated time of hit.
Utility class to parse command line options.
Definition: JParser.hh:1514
double getAngle(const JQuaternion3D &first, const JQuaternion3D &second)
Get space angle between quanternions.
General exception.
Definition: JException.hh:24
const Trk & get_muon(const Evt &evt)
Get first muon from the event tracklist.
then usage $script< input file >[option[primary[working directory]]] nWhere option can be E
Definition: JMuonPostfit.sh:40
JTrack3E getTrack(const Trk &track)
Get track.
#define STATUS(A)
Definition: JMessage.hh:63
Detector data structure.
Definition: JDetector.hh:89
double getIntersection(const JVector3D &pos) const
Get longitudinal position along axis of position of closest approach with given position.
Definition: JAxis3D.hh:146
Router for direct addressing of module data in detector data structure.
std::string comment
use as you like
Definition: Trk.hh:35
static const JGeaneWater gWater
Function object for energy loss of muon in sea water.
Definition: JGeane.hh:381
string outputFile
bool has_reconstructed_jppmuon(const Evt &evt)
Test whether given event has a track with muon reconstruction.
3D track with energy.
Definition: JTrack3E.hh:30
void move(const double step, const double velocity, const JGeane &geane)
Move vertex along this track with given velocity.
Definition: JTrack3E.hh:117
JVersor3D getDirection(const JVector3D &pos) const
Get photon direction of Cherenkov light on PMT.
Definition: JTrack3D.hh:147
JTime & add(const JTime &value)
Addition operator.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Type definition of range.
Definition: JHead.hh:41
Detector file.
Definition: JHead.hh:226
double getE() const
Get energy.
Definition: JTrack3E.hh:93
const Trk & get_best_reconstructed_jppmuon(const Evt &evt)
Get best reconstructed muon.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1989
set_variable E_E log10(E_{fit}/E_{#mu})"
const JPosition3D & getPosition() const
Get position.
Definition: JPosition3D.hh:130
double mc_t
MC: time where the mc-event was put in the timeslice, since start of run (offset+frameidx*timeslice_d...
Definition: Evt.hh:47
double getT(const JVector3D &pos) const
Get arrival time of Cherenkov light at given position.
Definition: JTrack3D.hh:126
double getTimeSinceRTS(const int frame_index)
Get time in ns since last RTS for a given frame index.
Definition: JDAQClock.hh:263
#define FATAL(A)
Definition: JMessage.hh:67
const double getSpeedOfLight()
Get speed of light.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
bool has_muon(const Evt &evt)
Test whether given event has a muon.
General purpose class for object reading from a list of file names.
Data structure for L0 hit.
Definition: JHitL0.hh:27
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:84
std::vector< Hit > hits
list of hits
Definition: Evt.hh:38
do set_variable DETECTOR_TXT $WORKDIR detector
int debug
debug level
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:20
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62