1#ifndef __JAANET__JAANETTOOLKIT__
2#define __JAANET__JAANETTOOLKIT__
137 for (std::vector<Hit>::const_iterator hit = event.
mc_hits.begin(); hit != event.
mc_hits.end(); ++hit) {
288 inline double getW(
const Trk& track,
const size_t index,
const double value)
290 if (index < track.
fitinf.size())
291 return track.
fitinf[index];
409 ((
int)(track.
type / 1000)) % 10 == 0); }
418 ((
int)(track.
type / 1000)) % 10 != 0); }
430 const string typestr = to_string(abs(track.
type));
432 return (
is_baryon(track) || (typestr.size() == 10 && typestr.substr(0,2) ==
"10"));
497 vector<Trk>::const_iterator i = find_if(evt.
mc_trks.cbegin(), evt.
mc_trks.cend(),
499 return i != evt.
mc_trks.cend();
519 vector<Trk>::const_iterator i = find_if(evt.
mc_trks.cbegin(), evt.
mc_trks.cend(),
522 if (i != evt.
mc_trks.cend()) {
return *i; }
544 vector<Trk>::const_iterator i = find_if(evt.
mc_trks.cbegin(), evt.
mc_trks.cend(),
547 return i != evt.
mc_trks.cend();
566 vector<Trk>::const_iterator i = find_if(evt.
mc_trks.cbegin(), evt.
mc_trks.cend(),
572 THROW(
JIndexOutOfRange,
"JAANET::get_target_nucleus(): Cannot retrieve target nucleus for event " << evt.
id <<
".");
588 for (vector<Trk>::const_iterator i = evt.
mc_trks.cbegin(); i != evt.
mc_trks.cend(); ++i) {
609 for (vector<Trk>::const_iterator i = evt.
mc_trks.cbegin(); i != evt.
mc_trks.cend(); ++i) {
634 for (vector<Trk>::const_iterator i = event.
mc_trks.cbegin(); i != event.
mc_trks.cend(); ++i) {
643 THROW(
JValueOutOfRange,
"get_leading_lepton(): Given event does not correspond to a neutrino interaction.");
676 }
else if (!event.
mc_trks.empty()) {
680 if (i != event.
mc_trks.end()) {
686 }
else if (!event.
trks.empty()) {
849 for (vector<Trk>::const_iterator i = event.
mc_trks.cbegin(); i != event.
mc_trks.cend(); ++i) {
850 if (
is_muon(*i) && i->mother_id == tau.
id) {
873 for (vector<Trk>::const_iterator i = event.
mc_trks.cbegin(); i != event.
mc_trks.cend(); ++i) {
911 for (vector<Trk>::const_iterator i = event.
mc_trks.cbegin(); i != event.
mc_trks.cend(); ++i) {
985 const double energy = trk.
E;
1001 using namespace std;
1005 double E0 = neutrino.
E;
1007 for (vector<Trk>::const_iterator track = evt.
mc_trks.cbegin(); track != evt.
mc_trks.end(); ++track) {
1047 using namespace std;
1048 using namespace JPP;
1054 for (vector<Trk>::const_iterator track = evt.
mc_trks.cbegin(); track != evt.
mc_trks.cend(); ++track) {
1061 evt.
mc_trks[track->mother_id] : neutrino );
1063 const double distance = (track->pos - mother.
pos).len();
1088 return neutrino.
E * neutrino.
dir;
1101 using namespace std;
1102 using namespace JPP;
1108 for (vector<Trk>::const_iterator track = evt.
mc_trks.cbegin(); track != evt.
mc_trks.cend(); ++track) {
1112 double kineticEnergy = 0.0;
1117 evt.
mc_trks[track->mother_id] : neutrino );
1119 const double distance = (track->pos - mother.
pos).len();
1120 const double energy = gWater(track->E, -
distance);
1129 P1 += kineticEnergy * track->dir;
1144 using namespace std;
1145 using namespace JPP;
1149 for (vector<Trk>::const_iterator track = event.
mc_trks.cbegin(); track != event.
mc_trks.cend(); ++track) {
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition of particle types.
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
Data structure for direction in three dimensions.
Data structure for position in three dimensions.
double getY() const
Get y position.
double getZ() const
Get z position.
double getX() const
Get x position.
double getDY() const
Get y direction.
double getDX() const
Get x direction.
double getDZ() const
Get z direction.
Exception for accessing an index in a collection that is outside of its range.
Exception for accessing a value in a collection that is outside of its range.
Extensions to Evt data format.
const Trk & get_electron(const Evt &evt)
Get first electron from the event tracklist.
JAxis3D getAxis(const Trk &track)
Get axis.
const Trk & get_tau(const Evt &evt)
Get first tau from the event tracklist.
bool is_neutron(const Trk &track)
Test whether given track is a (anti-)neutron.
bool is_charged_lepton(const Trk &track)
Test whether given track is a charged lepton.
bool is_photon(const Trk &track)
Test whether given track is a photon or neutral pion.
bool has_primary(const Evt &evt)
Auxiliary function to check if an event contains a primary track.
JDirection3D getDirection(const Vec &dir)
Get direction.
int count_taudecay_prongs(const Evt &event)
Count the number of tau-lepton decay prongs in a given event.
int count_muons(const Evt &evt)
Count the number of muons in a given event.
bool from_tau(const Hit &hit)
Test whether given hit was produced by a tau.
bool has_target_nucleus(const Evt &evt)
Test whether given event has well-defined target nucleus.
JTrack3E getTrack(const Trk &track)
Get track.
bool is_electron(const Trk &track)
Test whether given track is a (anti-)electron.
bool has_tau(const Evt &evt)
Test whether given event has a tau.
bool is_nucleus(const Trk &track)
Test whether given track is a nucleus (PDG code corresponds to a baryon or has 10 digits,...
bool from_muon(const Hit &hit)
Test whether given hit was produced by a muon.
void add_time_offset(Evt &evt, const double tOff)
Add time offset to mc event; according to current definition, the absolute time of the event is defin...
Vec getP0(const Evt &evt)
Get momentum vector of the initial state of a neutrino interaction.
double getInvariantMass(const Evt &event)
Get final state invariant mass.
double getKineticEnergy(const Trk &trk)
Get track kinetic energy.
bool from_hadron(const Hit &hit)
Test whether given hit was produced by a hadronic shower.
const Trk & get_target_nucleus(const Evt &evt)
Get target nucleus.
JVertex3D getVertex(const Trk &track)
Get vertex.
int count_taus(const Evt &evt)
Count the number of taus in a given event.
JTimeRange getTimeRange(const Evt &event)
Get time range (i.e. time between earliest and latest hit) of Monte Carlo event.
const Trk & get_leading_lepton(const Evt &event)
Auxiliary function to retrieve the leading lepton of a neutrino interaction.
bool is_muonbundle(const Trk &track)
Test whether given track is a (anti-)muon.
JTransformation3D getTransformation(const Trk &track)
Get transformation.
bool has_muon(const Evt &evt)
Test whether given event has a muon.
bool is_pion(const Trk &track)
Test whether given track is a charged pion.
const Trk & get_primary(const Evt &evt)
Auxiliary function to retrieve the primary track of an event.
static const int AASHOWER_RECONSTRUCTION_TYPE
AAnet shower fit reconstruction type.
bool has_neutrino(const Evt &evt)
Test whether given event has an incoming neutrino.
bool is_finalstate(const Trk &track)
Test whether given track corresponds to a final state particle.
bool has_electron(const Evt &evt)
Test whether given event has an electron.
bool is_neutrino(const Trk &track)
Test whether given track is a neutrino.
double getE1(const Evt &evt)
Get final state energy of a neutrino interaction.
double getTime(const Hit &hit)
Get true time of hit.
double getE0(const Evt &evt)
Get initial state energy of a neutrino interaction.
const Trk & get_muon(const Evt &evt)
Get first muon from the event tracklist.
bool is_initialstate(const Trk &track)
Test whether given track corresponds to an initial state particle.
double getW(const Trk &track, const size_t index, const double value)
Get track information.
int count_electrons(const Evt &evt)
Count the number of electrons in a given event.
JPosition3D getPosition(const Vec &pos)
Get position.
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
JHitType_t
Enumeration of hit types based on km3 codes.
@ HIT_TYPE_UNKNOWN
Unknown source.
@ HIT_TYPE_NOISE
Random noise.
@ HIT_TYPE_DELTARAYS_DIRECT
Direct light from delta-rays.
@ HIT_TYPE_MUON_DIRECT
Direct light from muon.
@ HIT_TYPE_DELTARAYS_SCATTERED
Scattered light from delta-rays.
@ HIT_TYPE_MUON_SCATTERED
Scattered light from muon.
@ HIT_TYPE_SHOWER_DIRECT
Direct light from primary shower.
@ HIT_TYPE_BREMSSTRAHLUNG_SCATTERED
Scattered light from Bremsstrahlung.
@ HIT_TYPE_SHOWER_SCATTERED
Scattered light from primary shower.
@ HIT_TYPE_BREMSSTRAHLUNG_DIRECT
Direct light from Bremsstrahlung.
bool is_noise(const Hit &hit)
Verify hit origin.
bool is_lepton(const Trk &track)
Test whether given track is a lepton.
bool from_electron(const Hit &hit)
Test whether given hit was produced by an electron.
bool is_proton(const Trk &track)
Test whether given track is a (anti-)proton.
bool is_baryon(const Trk &track)
Test whether given track is a baryon (is hadron and third digit of PDG code is not zero)
@ GEANT4_TYPE_ANTIELECTRON
double getNPE(const Hit &hit)
Get true charge of hit.
Vec getP1(const Evt &evt)
Get momentum vector of the final state of a neutrino interaction.
bool is_meson(const Trk &track)
Test whether given track is a meson (is hadron and third digit of PDG code is zero)
bool has_particleID(const Trk &track, const int type)
Test whether given track corresponds to given particle type.
bool has_leptonic_taudecay(const Evt &event)
Test whether given event contains a leptonic tau-decay.
bool has_hadronic_taudecay(const Evt &event)
Test whether given event contains a hadronic tau-decay.
const Trk & get_neutrino(const Evt &evt)
Get incoming neutrino.
bool has_muonic_taudecay(const Evt &event)
Test whether given event contains a tau-lepton decay into a muon.
bool is_hadron(const Trk &track)
Test whether given track is a hadron.
@ TRACK_TYPE_CHARGED_PION_PLUS
@ TRACK_TYPE_NEUTRAL_PION
@ TRACK_TYPE_NEUTRAL_SIGMA
bool is_tau(const Trk &track)
Test whether given track is a (anti-)tau.
int count_hadrons(const Evt &evt)
Count the number of hadrons in a given event (not including neutral pions)
JPredicate< JResult_t T::*, JComparison::eq > make_predicate(JResult_t T::*member, const JResult_t value)
Helper method to create predicate for data member.
static const double MASS_NEUTRON
neutron mass [GeV]
double getKineticEnergy(const double E, const double m)
Get kinetic energy of particle with given energy and mass.
static const double MASS_PROTON
proton mass [GeV]
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
std::vector< Hit > mc_hits
MC: list of MC truth hits.
std::vector< Trk > mc_trks
MC: list of MC truth tracks.
int id
offline event identifier
static int ROOT_IO_VERSION
Streamer version as obtained from ROOT file.
std::vector< Trk > trks
list of reconstructed tracks (can be several because of prefits,showers, etc).
double a
hit amplitude (in p.e.)
int type
particle type or parametrisation used for hit (mc only)
double t
hit time (from tdc+calibration or MC truth)
static const JPDB & getInstance()
Get particle data book.
const JParticle & getPDG(const int pdg) const
Get particle corresponding to given PDG code.
double mass
mass of particle [GeV]
int charge
charge of particle [|e|/3]
The Trk class represents a Monte Carlo (MC) particle as well as a reconstructed track/shower.
int status
MC status code, see km3net-dataformat/definitions/trkmembers.csv for values.
int type
MC: particle type in PDG encoding.
std::vector< double > fitinf
place to store additional fit info, see km3net-dataformat/definitions/fitparameters....
double E
Energy [GeV] (either MC truth or reconstructed)
double t
track time [ns] (when the particle is at pos )
int mother_id
MC id of the parent particle.
Vec pos
postion [m] of the track at time t
The Vec class is a straightforward 3-d vector, which also works in pyroot.
static const int PDG_MUONBUNDLE
muon bundle reached the can level (mupage)
static const int TRK_ST_FINALSTATE
for MC: the particle must be processed by detector simulation ('track_in' tag in evt files)....
static const int TRK_MOTHER_UNDEFINED
KM3NeT Data Definitions v3.6.0 https://git.km3net.de/common/km3net-dataformat.
static const int TRK_ST_MUONBUNDLE
initial state muon bundle (mupage)
static const int TRK_ST_PRIMARYCOSMIC
initial state cosmic ray ('track_primary' tag in evt files from corant).
static const int TRK_ST_PRIMARYNEUTRINO
initial state neutrino ('neutrino' tag in evt files from gseagen and genhen).
static const int TRK_ST_ININUCLEI
Initial state nuclei (gseagen)
static const int TRK_MOTHER_NONE
mother id of a particle if it has no parent