Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Classes | Enumerations | Functions | Variables
JAANET Namespace Reference

Extensions to AAnet data format. More...

Classes

struct  quality_sorter
 Reconstruction type dependent comparison of track quality. More...
 
struct  has_history
 Auxiliary class to test whether given track has specified history. More...
 
struct  JEvtEvaluator
 Auxiliary class to determine value of Evt objects. More...
 
struct  start_run
 Start of run record. More...
 
struct  String
 General purpose string class. More...
 
struct  detector
 Detector file. More...
 
struct  muon_desc_file
 Muon descriptor file. More...
 
struct  target
 Target. More...
 
struct  XSecFile
 Neutrino cross section file. More...
 
struct  cut
 General purpose class of phase space generation. More...
 
struct  cut_primary
 Phase space of primary particle. More...
 
struct  cut_seamuon
 Phase space of atmospheric muon generation. More...
 
struct  cut_in
 Phase space of incoming particle. More...
 
struct  cut_nu
 Phase space of incident neutrino. More...
 
struct  generator
 Description of Monte Carlo event generation applications. More...
 
struct  physics
 Generator for neutrino interaction. More...
 
struct  simul
 Generator for ? More...
 
struct  spectrum
 Neutrino energy spectrum. More...
 
struct  can
 The cylinder used for photon tracking. More...
 
struct  genvol
 Neutrino vertex volume. More...
 
struct  coord_origin
 Coordinate origin. More...
 
struct  genhencut
 Phase space for incident neutrino. More...
 
struct  norma
 Normlisation of CORSIKA events. More...
 
struct  livetime
 Normalisation of MUPAGE events. More...
 
struct  flux
 Neutrino flux. More...
 
struct  seabottom
 The bottom of the sea. More...
 
struct  DAQ
 Normalisation of DAQ data. More...
 
struct  end_event
 End of event record. More...
 
struct  JHead
 Monte Carlo run header. More...
 
class  JAAnetDictionary
 Simple wrapper class around JROOT::JRootDictionary so that other classes could be included in this dictionary if necessary. More...
 
struct  JHit_t
 Auxiliary class to set-up Hit. More...
 
struct  JHits_t
 Auxiliary data structure for list of hits with hit merging capability. More...
 
struct  JParticle
 Auxiliary class to handle particle name, codes and mass. More...
 
struct  JPDB
 Collection of particles. More...
 
struct  JTrk_t
 Auxiliary class to set-up Trk. More...
 

Enumerations

enum  JHitType_t {
  HIT_TYPE_MUON_DIRECT = +5, HIT_TYPE_MUON_SCATTERED = -5, HIT_TYPE_DELTARAYS_DIRECT = +4, HIT_TYPE_DELTARAYS_SCATTERED = -4,
  HIT_TYPE_BREMSSTRAHLUNG_DIRECT = +3, HIT_TYPE_BREMSSTRAHLUNG_SCATTERED = -3, HIT_TYPE_SHOWER_DIRECT = +99, HIT_TYPE_SHOWER_SCATTERED = -99,
  HIT_TYPE_NOISE = -1, HIT_TYPE_UNKNOWN = 0
}
 Enumeration of hit types based on km3 codes. More...
 
enum  JGeant4Type_t {
  GEANT4_TYPE_PHOTON = 1, GEANT4_TYPE_ANTIELECTRON = 2, GEANT4_TYPE_ELECTRON = 3, GEANT4_TYPE_NEUTRINO = 4,
  GEANT4_TYPE_ANTIMUON = 5, GEANT4_TYPE_MUON = 6, GEANT4_TYPE_NEUTRAL_PION = 7, GEANT4_TYPE_PION_PLUS = 8,
  GEANT4_TYPE_PION_MINUS = 9, GEANT4_TYPE_KAON_LONG = 10, GEANT4_TYPE_KAON_PLUS = 11, GEANT4_TYPE_KAON_MINUS = 12,
  GEANT4_TYPE_NEUTRON = 13, GEANT4_TYPE_PROTON = 14, GEANT4_TYPE_ANTIPROTON = 15, GEANT4_TYPE_KAON_SHORT = 16,
  GEANT4_TYPE_ETA = 17, GEANT4_TYPE_LAMBDA = 18, GEANT4_TYPE_SIGMA_PLUS = 19, GEANT4_TYPE_NEUTRAL_SIGMA = 20,
  GEANT4_TYPE_SIGMA_MINUS = 21, GEANT4_TYPE_NEUTRAL_XI = 22, GEANT4_TYPE_XI_MINUS = 23, GEANT4_TYPE_OMEGA_MINUS = 24,
  GEANT4_TYPE_ANTINEUTRON = 25, GEANT4_TYPE_ANTILAMBDA = 26, GEANT4_TYPE_ANTISIGMA_MINUS = 27, GEANT4_TYPE_NEUTRAL_ANTISIGMA = 28,
  GEANT4_TYPE_ANTISIGMA_PLUS = 29, GEANT4_TYPE_NEUTRAL_ANTIXI = 30, GEANT4_TYPE_ANTIXI_PLUS = 31, GEANT4_TYPE_ANTIOMEGA_PLUS = 32,
  GEANT4_TYPE_DEUTERON = 45, GEANT4_TYPE_TRITON = 46, GEANT4_TYPE_ALPHA = 47, GEANT4_TYPE_GEANTINO = 48,
  GEANT4_TYPE_HE3 = 49, GEANT4_TYPE_ANTITAU = 33, GEANT4_TYPE_TAU = 34
}
 Enumeration of hit types based on Geant4 codes, for compatbility with KM3Sim. More...
 
enum  JTrackType_t {
  TRACK_TYPE_ELECTRON = 11, TRACK_TYPE_NUE = 12, TRACK_TYPE_MUON = 13, TRACK_TYPE_NUMU = 14,
  TRACK_TYPE_TAU = 15, TRACK_TYPE_NUTAU = 16, TRACK_TYPE_PHOTON = 22, TRACK_TYPE_NEUTRAL_PION = 111,
  TRACK_TYPE_CHARGED_PION_PLUS = 211, TRACK_TYPE_PION_PLUS = 211, TRACK_TYPE_K_LONG = 130, TRACK_TYPE_K_SHORT = 310,
  TRACK_TYPE_K_PLUS = 321, TRACK_TYPE_PROTON = 2212, TRACK_TYPE_NEUTRON = 2112, TRACK_TYPE_LAMBDA = 3122,
  TRACK_TYPE_SIGMA_PLUS = 3222, TRACK_TYPE_NEUTRAL_SIGMA = 3212, TRACK_TYPE_SIGMA_MINUS = 3112, TRACK_TYPE_NEUTRAL_XI = 3322,
  TRACK_TYPE_XI_MINUS = 3312, TRACK_TYPE_OMEGA_MINUS = 3334, TRACK_TYPE_ANTIELECTRON = -11, TRACK_TYPE_ANTINUE = -12,
  TRACK_TYPE_ANTIMUON = -13, TRACK_TYPE_ANTINUMU = -14, TRACK_TYPE_ANTITAU = -15, TRACK_TYPE_ANTINUTAU = -16,
  TRACK_TYPE_NEUTRAL_ANTIPION = -111, TRACK_TYPE_CHARGED_PION_MINUS = -211, TRACK_TYPE_PION_MINUS = -211, TRACK_TYPE_ANTIK_LONG = -130,
  TRACK_TYPE_ANTIK_SHORT = -310, TRACK_TYPE_K_MINUS = -321, TRACK_TYPE_ANTIPROTON = -2212, TRACK_TYPE_ANTINEUTRON = -2112,
  TRACK_TYPE_ANTILAMBDA = -3122, TRACK_TYPE_ANTISIGMA_PLUS = -3222, TRACK_TYPE_ANTINEUTRAL_SIGMA = -3212, TRACK_TYPE_ANTISIGMA_MINUS = -3112,
  TRACK_TYPE_ANTINEUTRAL_XI = -3322, TRACK_TYPE_ANTIXI_MINUS = -3312, TRACK_TYPE_ANTIOMEGA_MINUS = -3334
}
 Enumeration of track types based on PDG codes. More...
 

Functions

double getTime (const Hit &hit)
 Get true time of hit. More...
 
double getNPE (const Hit &hit)
 Get true charge of hit. More...
 
bool is_noise (const Hit &hit)
 Verify hit origin. More...
 
JTimeRange getTimeRange (const Evt &event)
 Get time range (i.e. More...
 
JTimeRange getTimeRange (const Evt &event, const JTimeRange &T_ns)
 Get time range (i.e. More...
 
JPosition3D getPosition (const Vec &v)
 Get position. More...
 
JPosition3D getPosition (const Trk &track)
 Get position. More...
 
JDirection3D getDirection (const Vec &v)
 Get direction. More...
 
JDirection3D getDirection (const Trk &track)
 Get direction. More...
 
JAxis3D getAxis (const Trk &track)
 Get axis. More...
 
JTrack3E getTrack (const Trk &track)
 Get track. More...
 
JVertex3D getVertex (const Trk &track)
 Get vertex. More...
 
JTransformation3D getTransformation (const Trk &track)
 Get transformation. More...
 
JHitL0 getHit (const Hit &hit)
 Get transformation. More...
 
JHitL0 getHit (const Hit &hit, const JModuleRouter &router)
 Get transformation. More...
 
JHitL0 getHit (const Hit &hit, const JPMTRouter &router)
 Get transformation. More...
 
bool is_photon (const Trk &track)
 Test whether given track is a photon or neutral pion. More...
 
bool is_neutrino (const Trk &track)
 Test whether given track is a neutrino. More...
 
bool is_electron (const Trk &track)
 Test whether given track is a (anti-)electron. More...
 
bool is_muon (const Trk &track)
 Test whether given track is a (anti-)muon. More...
 
bool is_tau (const Trk &track)
 Test whether given track is a (anti-)tau. More...
 
bool is_pion (const Trk &track)
 Test whether given track is a charged pion. More...
 
bool is_proton (const Trk &track)
 Test whether given track is a (anti-)proton. More...
 
bool is_neutron (const Trk &track)
 Test whether given track is a (anti-)neutron. More...
 
bool is_hadron (const Trk &track)
 Test whether given track is a hadron. More...
 
bool has_particleID (const Trk &track, const int type)
 Test whether given track corresponds to given particle type. More...
 
bool has_neutrino (const Evt &evt)
 Test whether given event has an incoming neutrino. More...
 
const Trk & get_neutrino (const Evt &evt)
 Get incoming neutrino. More...
 
int count_electrons (const Evt &evt)
 Count the number of electrons in a given event. More...
 
bool has_electron (const Evt &evt)
 Test whether given event has an electron. More...
 
const Trk & get_electron (const Evt &evt)
 Get first electron from the event tracklist. More...
 
bool from_electron (const Hit &hit)
 Test whether given hit was produced by an electron. More...
 
int count_muons (const Evt &evt)
 Count the number of muons in a given event. More...
 
bool has_muon (const Evt &evt)
 Test whether given event has a muon. More...
 
const Trk & get_muon (const Evt &evt)
 Get first muon from the event tracklist. More...
 
bool from_muon (const Hit &hit)
 Test whether given hit was produced by a muon. More...
 
int count_taus (const Evt &evt)
 Count the number of taus in a given event. More...
 
bool has_tau (const Evt &evt)
 Test whether given event has a tau. More...
 
const Trk & get_tau (const Evt &evt)
 Get first tau from the event tracklist. More...
 
bool from_tau (const Hit &hit)
 Test whether given hit was produced by a tau. More...
 
int count_hadrons (const Evt &evt)
 Count the number of hadrons in a given event (not including neutral pions) More...
 
bool from_hadron (const Hit &hit)
 Test whether given hit was produced by a hadronic shower. More...
 
Trk get_hadronic_cascade (const Evt &evt)
 Auxiliary function to get average direction and total energy of a hadronic cascade. More...
 
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 defined by the track "t" attribute; this could change in the future if the global attribute mc_t is assigned to this purpose. More...
 
bool has_muon_prefit (const Trk &track)
 Test whether given track has muon prefit in history. More...
 
bool has_muon_simplex (const Trk &track)
 Test whether given track has muon simplex fit in history. More...
 
bool has_muon_gandalf (const Trk &track)
 Test whether given track has muon gandalf fit in history. More...
 
bool has_muon_energy (const Trk &track)
 Test whether given track has muon energy fit in history. More...
 
bool has_muon_start (const Trk &track)
 Test whether given track has muon start fit in history. More...
 
bool has_muon_fit (const Trk &track)
 Test whether given track has default muon fit in history. More...
 
bool has_shower_prefit (const Trk &track)
 Test whether given track has shower prefit in history. More...
 
bool has_shower_positionfit (const Trk &track)
 Test whether given track has shower position fit in history. More...
 
bool has_shower_completefit (const Trk &track)
 Test whether given track has shower complete fit in history. More...
 
bool has_shower_fit (const Trk &track)
 Test whether given track has default shower fit in history. More...
 
bool has_aashower_fit (const Trk &track)
 Test whether given track has default shower fit in history. More...
 
template<class JTrackSelector_t >
bool has_reconstructed_track (const Evt &evt, JTrackSelector_t selector)
 Test whether given event has a track according selection. More...
 
template<int reconstruction_type>
bool has_reconstructed_track (const Evt &evt, const JRange< int > range=JRange< int >())
 Test whether given event has a track according selection. More...
 
bool has_reconstructed_muon (const Evt &evt)
 Test whether given event has a track with muon reconstruction. More...
 
bool has_reconstructed_shower (const Evt &evt)
 Test whether given event has a track with shower reconstruction. More...
 
bool has_reconstructed_aashower (const Evt &evt)
 Test whether given event has a track with aashower reconstruction. More...
 
template<class JTrackSelector_t , class JQualitySorter_t >
const Trk & get_best_reconstructed_track (const Evt &evt, JTrackSelector_t selector, JQualitySorter_t comparator)
 Get best reconstructed track. More...
 
template<int reconstruction_type>
const Trk & get_best_reconstructed_track (const Evt &evt, const JRange< int > range=JRange< int >())
 Get best reconstructed track. More...
 
const Trk & get_best_reconstructed_muon (const Evt &evt)
 Get best reconstructed muon. More...
 
const Trk & get_best_reconstructed_shower (const Evt &evt)
 Get best reconstructed shower. More...
 
const Trk & get_best_reconstructed_aashower (const Evt &evt)
 Get best reconstructed aashower. More...
 
void copy (const Head &from, JHead &to)
 Copy header from from to to. More...
 
void copy (const JHead &from, Head &to)
 Copy header from from to to. More...
 
std::string getTag (const std::string &tag)
 Get tag without aanet extension "_<counter>" for identical tags. More...
 
template<class T >
bool is_valid (const T &value)
 Check validity of given value. More...
 
bool operator== (const Head &first, const Head &second)
 Equal operator. More...
 
bool operator< (const Head &first, const Head &second)
 Less than operator. More...
 
template<class T >
void push (JHead &header, T JHead::*p)
 Push data of JHead for subsequent copy to Head. More...
 
bool is_physics (const JHead &header, const std::string &generator)
 Check for generator. More...
 
bool is_sirene (const JHead &header)
 Check for generator. More...
 
bool is_km3 (const JHead &header)
 Check for generator. More...
 
bool is_km3sim (const JHead &header)
 Check for generator. More...
 
bool is_genhen (const JHead &header)
 Check for generator. More...
 
bool is_genie (const JHead &header)
 Check for generator. More...
 
bool is_mupage (const JHead &header)
 Check for generator. More...
 
template<class T >
get (const JHead &head)
 Get object from header. More...
 
template<>
Vec get (const JHead &header)
 Get position offset of detector due to generator. More...
 
template<>
JPosition3D get (const JHead &header)
 Get position offset of detector due to generator. More...
 
template<>
JCylinder3D get (const JHead &head)
 Get cylinder corresponding to can. More...
 
const JEquationParametersgetEquationParameters ()
 Get equation parameters corresponding to Monte Carlo ASCII format, i.e: More...
 
std::istream & operator>> (std::istream &in, JHead &header)
 Read header from input. More...
 
std::ostream & operator<< (std::ostream &out, const JHead &header)
 Write header to output. More...
 
bool operator< (const JHit_t &first, const JHit_t &second)
 Less than operator for hits. More...
 

Variables

static const int AASHOWER_RECONSTRUCTION_TYPE = 101
 AAshower reconstruction type for AAnet. More...
 
static const JEvtEvaluator getEvtValue
 Function object for evaluation of DAQ objects. More...
 
static const char AANET_TAG_SEPARATOR = '_'
 Separator for AAnet tag extension for multiple tags ("_<counter>"). More...
 

Detailed Description

Extensions to AAnet data format.

Author
mlincetto
mdejong

Enumeration Type Documentation

Enumeration of hit types based on km3 codes.

Enumerator
HIT_TYPE_MUON_DIRECT 

Direct light from muon.

HIT_TYPE_MUON_SCATTERED 

Scattered light from muon.

HIT_TYPE_DELTARAYS_DIRECT 

Direct light from delta-rays.

HIT_TYPE_DELTARAYS_SCATTERED 

Scattered light from delta-rays.

HIT_TYPE_BREMSSTRAHLUNG_DIRECT 

Direct light from Bremsstrahlung.

HIT_TYPE_BREMSSTRAHLUNG_SCATTERED 

Scattered light from Bremsstrahlung.

HIT_TYPE_SHOWER_DIRECT 

Direct light from primary shower.

HIT_TYPE_SHOWER_SCATTERED 

Scattered light from primary shower.

HIT_TYPE_NOISE 

Random noise.

HIT_TYPE_UNKNOWN 

Unknown source.

Definition at line 67 of file JAAnetToolkit.hh.

67  {
68  HIT_TYPE_MUON_DIRECT = +5, //!< Direct light from muon
69  HIT_TYPE_MUON_SCATTERED = -5, //!< Scattered light from muon
70  HIT_TYPE_DELTARAYS_DIRECT = +4, //!< Direct light from delta-rays
71  HIT_TYPE_DELTARAYS_SCATTERED = -4, //!< Scattered light from delta-rays
72  HIT_TYPE_BREMSSTRAHLUNG_DIRECT = +3, //!< Direct light from Bremsstrahlung
73  HIT_TYPE_BREMSSTRAHLUNG_SCATTERED = -3, //!< Scattered light from Bremsstrahlung
74  HIT_TYPE_SHOWER_DIRECT = +99, //!< Direct light from primary shower
75  HIT_TYPE_SHOWER_SCATTERED = -99, //!< Scattered light from primary shower
76  HIT_TYPE_NOISE = -1, //!< Random noise
77  HIT_TYPE_UNKNOWN = 0 //!< Unknown source
78  };
Scattered light from primary shower.
Direct light from delta-rays.
Scattered light from muon.
Scattered light from Bremsstrahlung.
Direct light from Bremsstrahlung.
Direct light from primary shower.
Direct light from muon.
Scattered light from delta-rays.

Enumeration of hit types based on Geant4 codes, for compatbility with KM3Sim.

Enumerator
GEANT4_TYPE_PHOTON 
GEANT4_TYPE_ANTIELECTRON 
GEANT4_TYPE_ELECTRON 
GEANT4_TYPE_NEUTRINO 
GEANT4_TYPE_ANTIMUON 
GEANT4_TYPE_MUON 
GEANT4_TYPE_NEUTRAL_PION 
GEANT4_TYPE_PION_PLUS 
GEANT4_TYPE_PION_MINUS 
GEANT4_TYPE_KAON_LONG 
GEANT4_TYPE_KAON_PLUS 
GEANT4_TYPE_KAON_MINUS 
GEANT4_TYPE_NEUTRON 
GEANT4_TYPE_PROTON 
GEANT4_TYPE_ANTIPROTON 
GEANT4_TYPE_KAON_SHORT 
GEANT4_TYPE_ETA 
GEANT4_TYPE_LAMBDA 
GEANT4_TYPE_SIGMA_PLUS 
GEANT4_TYPE_NEUTRAL_SIGMA 
GEANT4_TYPE_SIGMA_MINUS 
GEANT4_TYPE_NEUTRAL_XI 
GEANT4_TYPE_XI_MINUS 
GEANT4_TYPE_OMEGA_MINUS 
GEANT4_TYPE_ANTINEUTRON 
GEANT4_TYPE_ANTILAMBDA 
GEANT4_TYPE_ANTISIGMA_MINUS 
GEANT4_TYPE_NEUTRAL_ANTISIGMA 
GEANT4_TYPE_ANTISIGMA_PLUS 
GEANT4_TYPE_NEUTRAL_ANTIXI 
GEANT4_TYPE_ANTIXI_PLUS 
GEANT4_TYPE_ANTIOMEGA_PLUS 
GEANT4_TYPE_DEUTERON 
GEANT4_TYPE_TRITON 
GEANT4_TYPE_ALPHA 
GEANT4_TYPE_GEANTINO 
GEANT4_TYPE_HE3 
GEANT4_TYPE_ANTITAU 
GEANT4_TYPE_TAU 

Definition at line 18 of file JParticleTypes.hh.

18  { GEANT4_TYPE_PHOTON = 1,
23  GEANT4_TYPE_MUON = 6,
31  GEANT4_TYPE_PROTON = 14,
34  GEANT4_TYPE_ETA = 17,
35  GEANT4_TYPE_LAMBDA = 18,
51  GEANT4_TYPE_TRITON = 46,
52  GEANT4_TYPE_ALPHA = 47,
54  GEANT4_TYPE_HE3 = 49,
55  //
56  // KM3NeT specific codes
57  //
59  GEANT4_TYPE_TAU = 34 };

Enumeration of track types based on PDG codes.

Enumerator
TRACK_TYPE_ELECTRON 
TRACK_TYPE_NUE 
TRACK_TYPE_MUON 
TRACK_TYPE_NUMU 
TRACK_TYPE_TAU 
TRACK_TYPE_NUTAU 
TRACK_TYPE_PHOTON 
TRACK_TYPE_NEUTRAL_PION 
TRACK_TYPE_CHARGED_PION_PLUS 
TRACK_TYPE_PION_PLUS 
TRACK_TYPE_K_LONG 
TRACK_TYPE_K_SHORT 
TRACK_TYPE_K_PLUS 
TRACK_TYPE_PROTON 
TRACK_TYPE_NEUTRON 
TRACK_TYPE_LAMBDA 
TRACK_TYPE_SIGMA_PLUS 
TRACK_TYPE_NEUTRAL_SIGMA 
TRACK_TYPE_SIGMA_MINUS 
TRACK_TYPE_NEUTRAL_XI 
TRACK_TYPE_XI_MINUS 
TRACK_TYPE_OMEGA_MINUS 
TRACK_TYPE_ANTIELECTRON 
TRACK_TYPE_ANTINUE 
TRACK_TYPE_ANTIMUON 
TRACK_TYPE_ANTINUMU 
TRACK_TYPE_ANTITAU 
TRACK_TYPE_ANTINUTAU 
TRACK_TYPE_NEUTRAL_ANTIPION 
TRACK_TYPE_CHARGED_PION_MINUS 
TRACK_TYPE_PION_MINUS 
TRACK_TYPE_ANTIK_LONG 
TRACK_TYPE_ANTIK_SHORT 
TRACK_TYPE_K_MINUS 
TRACK_TYPE_ANTIPROTON 
TRACK_TYPE_ANTINEUTRON 
TRACK_TYPE_ANTILAMBDA 
TRACK_TYPE_ANTISIGMA_PLUS 
TRACK_TYPE_ANTINEUTRAL_SIGMA 
TRACK_TYPE_ANTISIGMA_MINUS 
TRACK_TYPE_ANTINEUTRAL_XI 
TRACK_TYPE_ANTIXI_MINUS 
TRACK_TYPE_ANTIOMEGA_MINUS 

Definition at line 64 of file JParticleTypes.hh.

64  { TRACK_TYPE_ELECTRON = 11,
65  TRACK_TYPE_NUE = 12,
66  TRACK_TYPE_MUON = 13,
67  TRACK_TYPE_NUMU = 14,
68  TRACK_TYPE_TAU = 15,
69  TRACK_TYPE_NUTAU = 16,
70  TRACK_TYPE_PHOTON = 22,
74  TRACK_TYPE_K_LONG = 130,
75  TRACK_TYPE_K_SHORT = 310,
76  TRACK_TYPE_K_PLUS = 321,
77  TRACK_TYPE_PROTON = 2212,
78  TRACK_TYPE_NEUTRON = 2112,
79  TRACK_TYPE_LAMBDA = 3122,
80  TRACK_TYPE_SIGMA_PLUS = 3222,
83  TRACK_TYPE_NEUTRAL_XI = 3322,
84  TRACK_TYPE_XI_MINUS = 3312,
86 
88  TRACK_TYPE_ANTINUE = -12,
89  TRACK_TYPE_ANTIMUON = -13,
90  TRACK_TYPE_ANTINUMU = -14,
91  TRACK_TYPE_ANTITAU = -15,
95  TRACK_TYPE_PION_MINUS = -211,
96  TRACK_TYPE_ANTIK_LONG = -130,
98  TRACK_TYPE_K_MINUS = -321,
99  TRACK_TYPE_ANTIPROTON = -2212,
100  TRACK_TYPE_ANTINEUTRON = -2112,
101  TRACK_TYPE_ANTILAMBDA = -3122,
106  TRACK_TYPE_ANTIXI_MINUS = -3312,
107  TRACK_TYPE_ANTIOMEGA_MINUS = -3334 };

Function Documentation

double JAANET::getTime ( const Hit &  hit)
inline

Get true time of hit.

Parameters
hithit
Returns
true time of hit

Definition at line 87 of file JAAnetToolkit.hh.

88  {
89  // if IGAIN = 0 in km3 program, then pure_xx values are not set.
90 
91  if (hit.pure_a == 0.0)
92  return hit.t;
93  else
94  return hit.pure_t;
95  }
double JAANET::getNPE ( const Hit &  hit)
inline

Get true charge of hit.

Parameters
hithit
Returns
true charge of hit

Definition at line 104 of file JAAnetToolkit.hh.

105  {
106  // if IGAIN = 0 in km3 program, then pure_xx values are not set.
107 
108  if (hit.pure_a == 0.0)
109  return hit.a;
110  else
111  return hit.pure_a;
112  }
bool JAANET::is_noise ( const Hit &  hit)
inline

Verify hit origin.

Parameters
hithit
Returns
true if noise; else false

Definition at line 121 of file JAAnetToolkit.hh.

122  {
123  return hit.type == HIT_TYPE_NOISE;
124  }
JTimeRange JAANET::getTimeRange ( const Evt &  event)
inline

Get time range (i.e.

earlist and latest hit time) of Monte Carlo event. Note that the global event time in not taken into account.

Parameters
eventevent
Returns
time range

Definition at line 134 of file JAAnetToolkit.hh.

135  {
136  JTimeRange timeRange(JTimeRange::DEFAULT_RANGE);
137 
138  for (std::vector<Hit>::const_iterator hit = event.mc_hits.begin(); hit != event.mc_hits.end(); ++hit) {
139  if (!is_noise(*hit)) {
140  timeRange.include(getTime(*hit));
141  }
142  }
143 
144  return timeRange;
145  }
Structure to store the ToT mean and standard deviation of the hits produced by a nanobeacon in a sour...
double getTime(const Hit &hit)
Get true time of hit.
JRange< double > JTimeRange
Type definition for time range.
bool is_noise(const Hit &hit)
Verify hit origin.
JTimeRange JAANET::getTimeRange ( const Evt &  event,
const JTimeRange &  T_ns 
)
inline

Get time range (i.e.

earlist and latest hit time) of Monte Carlo event. The given time window is offset so that the lower limit is equal to the earliest hit time. The time window is then applied to all hits. Note that the global event time in not taken into account.

Parameters
eventevent
T_nstime window [ns]
Returns
time range

Definition at line 158 of file JAAnetToolkit.hh.

159  {
160  using namespace std;
161 
162  JTimeRange timeRange(JTimeRange::DEFAULT_RANGE);
163 
164  double t0 = numeric_limits<double>::max();
165 
166  for (std::vector<Hit>::const_iterator hit = event.mc_hits.begin(); hit != event.mc_hits.end(); ++hit) {
167  if (!is_noise(*hit)) {
168  if (getTime(*hit) < t0) {
169  t0 = getTime(*hit);
170  }
171  }
172  }
173 
174  if (t0 != numeric_limits<double>::max()) {
175 
176  JTimeRange window(T_ns);
177 
178  window.fixLowerLimit(t0);
179 
180  for (std::vector<Hit>::const_iterator hit = event.mc_hits.begin(); hit != event.mc_hits.end(); ++hit) {
181  if (!is_noise(*hit) && window(getTime(*hit))) {
182  timeRange.include(getTime(*hit));
183  }
184  }
185  }
186 
187  return timeRange;
188  }
Structure to store the ToT mean and standard deviation of the hits produced by a nanobeacon in a sour...
double getTime(const Hit &hit)
Get true time of hit.
JRange< double > JTimeRange
Type definition for time range.
bool is_noise(const Hit &hit)
Verify hit origin.
JPosition3D JAANET::getPosition ( const Vec &  v)
inline

Get position.

Parameters
vvector
Returns
position

Definition at line 197 of file JAAnetToolkit.hh.

198  {
199  return JPosition3D(v.x, v.y, v.z);
200  }
JPosition3D JAANET::getPosition ( const Trk &  track)
inline

Get position.

Parameters
tracktrack
Returns
position

Definition at line 209 of file JAAnetToolkit.hh.

210  {
211  return getPosition(track.pos);
212  }
JPosition3D getPosition(const Vec &v)
Get position.
JDirection3D JAANET::getDirection ( const Vec &  v)
inline

Get direction.

Parameters
vvector
Returns
direction

Definition at line 221 of file JAAnetToolkit.hh.

222  {
223  return JDirection3D(v.x, v.y, v.z);
224  }
JDirection3D JAANET::getDirection ( const Trk &  track)
inline

Get direction.

Parameters
tracktrack
Returns
direction

Definition at line 232 of file JAAnetToolkit.hh.

233  {
234  return getDirection(track.dir);
235  }
JDirection3D getDirection(const Vec &v)
Get direction.
JAxis3D JAANET::getAxis ( const Trk &  track)
inline

Get axis.

Parameters
tracktrack
Returns
axis

Definition at line 244 of file JAAnetToolkit.hh.

245  {
246  return JAxis3D(getPosition(track), getDirection(track));
247  }
JDirection3D getDirection(const Vec &v)
Get direction.
JPosition3D getPosition(const Vec &v)
Get position.
JTrack3E JAANET::getTrack ( const Trk &  track)
inline

Get track.

Parameters
tracktrack
Returns
track

Definition at line 256 of file JAAnetToolkit.hh.

257  {
258  return JTrack3E(JTrack3D(getAxis(track), track.t), track.E);
259  }
JAxis3D getAxis(const Trk &track)
Get axis.
JVertex3D JAANET::getVertex ( const Trk &  track)
inline

Get vertex.

Parameters
tracktrack
Returns
vertex

Definition at line 268 of file JAAnetToolkit.hh.

269  {
270  return JVertex3D(getPosition(track), track.t);
271  }
JPosition3D getPosition(const Vec &v)
Get position.
JTransformation3D JAANET::getTransformation ( const Trk &  track)
inline

Get transformation.

Parameters
tracktrack
Returns
transformation

Definition at line 280 of file JAAnetToolkit.hh.

281  {
282  return JTransformation3D(getPosition(track), getDirection(track));
283  }
JDirection3D getDirection(const Vec &v)
Get direction.
JPosition3D getPosition(const Vec &v)
Get position.
JHitL0 JAANET::getHit ( const Hit &  hit)
inline

Get transformation.

Parameters
hithit
Returns
hit

Definition at line 292 of file JAAnetToolkit.hh.

293  {
294  using namespace JPP;
295 
296  return JHitL0(JDAQPMTIdentifier(hit.dom_id, hit.channel_id), JAxis3D(getPosition(hit.pos), getDirection(hit.dir)), JHit(hit.t, hit.tot) );
297  }
JDirection3D getDirection(const Vec &v)
Get direction.
JPosition3D getPosition(const Vec &v)
Get position.
JHitL0 JAANET::getHit ( const Hit &  hit,
const JModuleRouter &  router 
)
inline

Get transformation.

Parameters
hithit
routermodule router
Returns
hit

Definition at line 307 of file JAAnetToolkit.hh.

308  {
309  using namespace JPP;
310 
311  const JPMT& pmt = router.getModule(hit.dom_id).getPMT(hit.channel_id);
312 
313  return JHitL0(JDAQPMTIdentifier(hit.dom_id, hit.channel_id), pmt.getAxis(), JHit(getTime(hit.tdc, pmt.getCalibration()), getToT(hit.tot, pmt.getCalibration())));
314  }
double getTime(const Hit &hit)
Get true time of hit.
double getToT(const T &tot, const JCalibration &cal)
Get calibrated time-over-threshold of hit.
JHitL0 JAANET::getHit ( const Hit &  hit,
const JPMTRouter &  router 
)
inline

Get transformation.

Parameters
hithit
routerpmt router
Returns
hit

Definition at line 324 of file JAAnetToolkit.hh.

325  {
326  using namespace JPP;
327 
328  const JPMT& pmt = router.getPMT(hit.pmt_id);
329 
330  return JHitL0(JDAQPMTIdentifier(hit.dom_id, hit.channel_id), pmt.getAxis(), JHit(getTime(hit.tdc, pmt.getCalibration()), getToT(hit.tot, pmt.getCalibration())));
331  }
double getTime(const Hit &hit)
Get true time of hit.
double getToT(const T &tot, const JCalibration &cal)
Get calibrated time-over-threshold of hit.
bool JAANET::is_photon ( const Trk &  track)
inline

Test whether given track is a photon or neutral pion.

Parameters
tracktrack
Returns
true if photon or neutral pion; else false

Definition at line 340 of file JAAnetToolkit.hh.

340  { return ( track.type == TRACK_TYPE_PHOTON ||
341  abs(track.type) == TRACK_TYPE_NEUTRAL_PION); }
bool JAANET::is_neutrino ( const Trk &  track)
inline

Test whether given track is a neutrino.

Parameters
tracktrack
Returns
true if neutrino; else false

Definition at line 349 of file JAAnetToolkit.hh.

349  { return (abs(track.type) == TRACK_TYPE_NUE ||
350  abs(track.type) == TRACK_TYPE_NUMU ||
351  abs(track.type) == TRACK_TYPE_NUTAU); }
bool JAANET::is_electron ( const Trk &  track)
inline

Test whether given track is a (anti-)electron.

Parameters
tracktrack
Returns
true if (anti-)electron; else false

Definition at line 359 of file JAAnetToolkit.hh.

359 { return (abs(track.type) == TRACK_TYPE_ELECTRON); }
bool JAANET::is_muon ( const Trk &  track)
inline

Test whether given track is a (anti-)muon.

Parameters
tracktrack
Returns
true if (anti-)muon; else false

Definition at line 367 of file JAAnetToolkit.hh.

367 { return (abs(track.type) == TRACK_TYPE_MUON); }
bool JAANET::is_tau ( const Trk &  track)
inline

Test whether given track is a (anti-)tau.

Parameters
tracktrack
Returns
true if (anti-)tau; else false

Definition at line 375 of file JAAnetToolkit.hh.

375 { return (abs(track.type) == TRACK_TYPE_TAU); }
bool JAANET::is_pion ( const Trk &  track)
inline

Test whether given track is a charged pion.

Parameters
tracktrack
Returns
true if charged pion; else false

Definition at line 383 of file JAAnetToolkit.hh.

bool JAANET::is_proton ( const Trk &  track)
inline

Test whether given track is a (anti-)proton.

Parameters
tracktrack
Returns
true if (anti-)proton; else false

Definition at line 391 of file JAAnetToolkit.hh.

391 { return (abs(track.type) == TRACK_TYPE_PROTON); }
bool JAANET::is_neutron ( const Trk &  track)
inline

Test whether given track is a (anti-)neutron.

Parameters
tracktrack
Returns
true if (anti-)neutron; else false

Definition at line 399 of file JAAnetToolkit.hh.

399 { return (abs(track.type) == TRACK_TYPE_NEUTRON); }
bool JAANET::is_hadron ( const Trk &  track)
inline

Test whether given track is a hadron.

Parameters
tracktrack
Returns
true if hadron; else fails

Definition at line 407 of file JAAnetToolkit.hh.

407  { return (!is_neutrino(track) &&
408  !is_electron(track) &&
409  !is_muon (track) &&
410  !is_tau (track)); }
bool is_electron(const Trk &track)
Test whether given track is a (anti-)electron.
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
bool is_neutrino(const Trk &track)
Test whether given track is a neutrino.
bool is_tau(const Trk &track)
Test whether given track is a (anti-)tau.
bool JAANET::has_particleID ( const Trk &  track,
const int  type 
)
inline

Test whether given track corresponds to given particle type.

Parameters
tracktrack
typeparticle type
Returns
true if matches the corresponding particle; else false

Definition at line 419 of file JAAnetToolkit.hh.

419 { return (track.type == type); }
bool JAANET::has_neutrino ( const Evt &  evt)
inline

Test whether given event has an incoming neutrino.

Parameters
evtevent
Returns
true if neutrino; else false

Definition at line 427 of file JAAnetToolkit.hh.

428  {
429  return !evt.mc_trks.empty() && is_neutrino(evt.mc_trks[0]);
430  }
bool is_neutrino(const Trk &track)
Test whether given track is a neutrino.
const Trk& JAANET::get_neutrino ( const Evt &  evt)
inline

Get incoming neutrino.

Parameters
evtevent
Returns
neutrino

Definition at line 438 of file JAAnetToolkit.hh.

439  {
440  if (has_neutrino(evt))
441  return evt.mc_trks[0];
442  else
443  THROW(JIndexOutOfRange, "This event has no neutrino.");
444  }
bool has_neutrino(const Evt &evt)
Test whether given event has an incoming neutrino.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:633
int JAANET::count_electrons ( const Evt &  evt)
inline

Count the number of electrons in a given event.

Parameters
evtevent
Returns
number of electrons

Definition at line 452 of file JAAnetToolkit.hh.

453  {
454  return std::count_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_electron);
455  }
bool is_electron(const Trk &track)
Test whether given track is a (anti-)electron.
bool JAANET::has_electron ( const Evt &  evt)
inline

Test whether given event has an electron.

Parameters
evtevent
Returns
true if event has electron; else false

Definition at line 463 of file JAAnetToolkit.hh.

464  {
465  return count_electrons(evt) != 0;
466  }
int count_electrons(const Evt &evt)
Count the number of electrons in a given event.
const Trk& JAANET::get_electron ( const Evt &  evt)
inline

Get first electron from the event tracklist.

Parameters
evtevent
Returns
electron

Definition at line 474 of file JAAnetToolkit.hh.

475  {
476  if (count_electrons(evt) > 0)
477  return *std::find_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_electron);
478  else
479  THROW(JIndexOutOfRange, "This event has no electron.");
480  }
bool is_electron(const Trk &track)
Test whether given track is a (anti-)electron.
int count_electrons(const Evt &evt)
Count the number of electrons in a given event.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:633
bool JAANET::from_electron ( const Hit &  hit)
inline

Test whether given hit was produced by an electron.

Warning: This method will only work with the output of KM3Sim. For JSirene or KM3, you should check if track.id == hit.origin instead.

Parameters
hithit
Returns
true if hit was produced by an electron; else false

Definition at line 491 of file JAAnetToolkit.hh.

int JAANET::count_muons ( const Evt &  evt)
inline

Count the number of muons in a given event.

Parameters
evtevent
Returns
number of muons

Definition at line 502 of file JAAnetToolkit.hh.

503  {
504  return std::count_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_muon);
505  }
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
bool JAANET::has_muon ( const Evt &  evt)
inline

Test whether given event has a muon.

Parameters
evtevent
Returns
true if event has muon; else false

Definition at line 513 of file JAAnetToolkit.hh.

514  {
515  return count_muons(evt) != 0;
516  }
int count_muons(const Evt &evt)
Count the number of muons in a given event.
const Trk& JAANET::get_muon ( const Evt &  evt)
inline

Get first muon from the event tracklist.

Parameters
evtevent
Returns
muon

Definition at line 524 of file JAAnetToolkit.hh.

525  {
526  if (count_muons(evt) > 0)
527  return *std::find_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_muon);
528  else
529  THROW(JIndexOutOfRange, "This event has no muon.");
530  }
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
int count_muons(const Evt &evt)
Count the number of muons in a given event.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:633
bool JAANET::from_muon ( const Hit &  hit)
inline

Test whether given hit was produced by a muon.

Warning: This method will only work with the output of KM3Sim. For JSirene or KM3, you should check if track.id == hit.origin instead.

Parameters
hithit
Returns
true if hit was produced by a muon; else false

Definition at line 541 of file JAAnetToolkit.hh.

542  {
543  return hit.type == GEANT4_TYPE_MUON || hit.type == GEANT4_TYPE_ANTIMUON;
544  }
int JAANET::count_taus ( const Evt &  evt)
inline

Count the number of taus in a given event.

Parameters
evtevent
Returns
number of taus

Definition at line 552 of file JAAnetToolkit.hh.

553  {
554  return std::count_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_tau);
555  }
bool is_tau(const Trk &track)
Test whether given track is a (anti-)tau.
bool JAANET::has_tau ( const Evt &  evt)
inline

Test whether given event has a tau.

Parameters
evtevent
Returns
true if event has tau; else false

Definition at line 563 of file JAAnetToolkit.hh.

564  {
565  return count_taus(evt) != 0;
566  }
int count_taus(const Evt &evt)
Count the number of taus in a given event.
const Trk& JAANET::get_tau ( const Evt &  evt)
inline

Get first tau from the event tracklist.

Parameters
evtevent
Returns
tau

Definition at line 574 of file JAAnetToolkit.hh.

575  {
576  if (count_taus(evt) > 0)
577  return *std::find_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_tau);
578  else
579  THROW(JIndexOutOfRange, "This event has no tau.");
580  }
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:633
int count_taus(const Evt &evt)
Count the number of taus in a given event.
bool is_tau(const Trk &track)
Test whether given track is a (anti-)tau.
bool JAANET::from_tau ( const Hit &  hit)
inline

Test whether given hit was produced by a tau.

Warning: This method will only work with the output of KM3Sim. For JSirene or KM3, you should check if track.id == hit.origin instead.

Parameters
hithit
Returns
true if hit was produced by a tau; else false

Definition at line 591 of file JAAnetToolkit.hh.

592  {
593  return hit.type == GEANT4_TYPE_TAU || hit.type == GEANT4_TYPE_ANTITAU;
594  }
int JAANET::count_hadrons ( const Evt &  evt)
inline

Count the number of hadrons in a given event (not including neutral pions)

Parameters
evtevent
Returns
number of hadrons

Definition at line 602 of file JAAnetToolkit.hh.

603  {
604  return std::count_if(evt.mc_trks.begin(), evt.mc_trks.end(), is_hadron);
605  }
bool is_hadron(const Trk &track)
Test whether given track is a hadron.
bool JAANET::from_hadron ( const Hit &  hit)
inline

Test whether given hit was produced by a hadronic shower.

Warning: This method will only work with the output of KM3Sim. For JSirene or KM3, you should check if track.id == hit.origin instead.

Parameters
hithit
Returns
true if hit was produced by a hadron; else false

Definition at line 616 of file JAAnetToolkit.hh.

617  {
618  return (!from_electron(hit) && !from_muon(hit) && !from_tau(hit));
619  }
bool from_tau(const Hit &hit)
Test whether given hit was produced by a tau.
bool from_electron(const Hit &hit)
Test whether given hit was produced by an electron.
bool from_muon(const Hit &hit)
Test whether given hit was produced by a muon.
Trk JAANET::get_hadronic_cascade ( const Evt &  evt)
inline

Auxiliary function to get average direction and total energy of a hadronic cascade.

Parameters
evtevent
Returns
cascade (default if there are no hadrons)

Definition at line 627 of file JAAnetToolkit.hh.

628  {
629  Trk cascade;
630 
631  for (vector<Trk>::const_iterator track = evt.mc_trks.begin(); track != evt.mc_trks.end(); ++track) {
632  if (is_hadron(*track)) {
633  cascade.E += track->E;
634  cascade.dir += track->dir * track->E;
635  cascade.t = track->t;
636  cascade.pos = track->pos;
637  }
638  }
639 
640  cascade.dir.normalize();
641 
642  return cascade;
643  }
bool is_hadron(const Trk &track)
Test whether given track is a hadron.
void JAANET::add_time_offset ( Evt &  evt,
const double  tOff 
)
inline

Add time offset to mc event; according to current definition, the absolute time of the event is defined by the track "t" attribute; this could change in the future if the global attribute mc_t is assigned to this purpose.

Parameters
evtevent
tOfftime offset [ns]

Definition at line 653 of file JAAnetToolkit.hh.

654  {
655  for (vector<Trk>::iterator p = evt.mc_trks.begin(); p != evt.mc_trks.end(); p++) {
656  p->t += tOff;
657  }
658  }
bool JAANET::has_muon_prefit ( const Trk &  track)
inline

Test whether given track has muon prefit in history.

Parameters
tracktrack
Returns
true if muon prefit in history; else false

Definition at line 727 of file JAAnetToolkit.hh.

728  {
730  }
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_muon_simplex ( const Trk &  track)
inline

Test whether given track has muon simplex fit in history.

Parameters
tracktrack
Returns
true if muon simplex fit in history; else false

Definition at line 738 of file JAAnetToolkit.hh.

739  {
741  }
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_muon_gandalf ( const Trk &  track)
inline

Test whether given track has muon gandalf fit in history.

Parameters
tracktrack
Returns
true if muon gandalf fit in history; else false

Definition at line 749 of file JAAnetToolkit.hh.

750  {
752  }
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_muon_energy ( const Trk &  track)
inline

Test whether given track has muon energy fit in history.

Parameters
tracktrack
Returns
true if muon energy fit in history; else false

Definition at line 760 of file JAAnetToolkit.hh.

761  {
763  }
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_muon_start ( const Trk &  track)
inline

Test whether given track has muon start fit in history.

Parameters
tracktrack
Returns
true if muon start fit in history; else false

Definition at line 771 of file JAAnetToolkit.hh.

772  {
774  }
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_muon_fit ( const Trk &  track)
inline

Test whether given track has default muon fit in history.

Parameters
tracktrack
Returns
true if muon fit in history; else false

Definition at line 782 of file JAAnetToolkit.hh.

783  {
785  }
End of muon fit applications.
Start of muon fit applications.
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_shower_prefit ( const Trk &  track)
inline

Test whether given track has shower prefit in history.

Parameters
tracktrack
Returns
true if shower prefit in history; else false

Definition at line 793 of file JAAnetToolkit.hh.

794  {
796  }
JShowerPrefit.cc.
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_shower_positionfit ( const Trk &  track)
inline

Test whether given track has shower position fit in history.

Parameters
tracktrack
Returns
true if shower position fit in history; else false

Definition at line 804 of file JAAnetToolkit.hh.

805  {
807  }
JShowerPositionFit.cc.
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_shower_completefit ( const Trk &  track)
inline

Test whether given track has shower complete fit in history.

Parameters
tracktrack
Returns
true if shower complete fit in history; else false

Definition at line 815 of file JAAnetToolkit.hh.

816  {
818  }
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_shower_fit ( const Trk &  track)
inline

Test whether given track has default shower fit in history.

Parameters
tracktrack
Returns
true if shower fit in history; else false

Definition at line 826 of file JAAnetToolkit.hh.

827  {
829  }
Start of shower fit applications.
End of shower fit applications.
static const int JPP_RECONSTRUCTION_TYPE
Jpp reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_aashower_fit ( const Trk &  track)
inline

Test whether given track has default shower fit in history.

Parameters
tracktrack
Returns
true if shower fit in history; else false

Definition at line 837 of file JAAnetToolkit.hh.

838  {
839  return has_history(AASHOWER_RECONSTRUCTION_TYPE, JRange<int>())(track);
840  }
static const int AASHOWER_RECONSTRUCTION_TYPE
AAshower reconstruction type for AAnet.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
template<class JTrackSelector_t >
bool JAANET::has_reconstructed_track ( const Evt &  evt,
JTrackSelector_t  selector 
)
inline

Test whether given event has a track according selection.


The track selector corresponds to the function operator bool selector(const Trk&);.

Parameters
evtevent
selectortrack selector
Returns
true if at least one corresponding track; else false

Definition at line 851 of file JAAnetToolkit.hh.

852  {
853  return std::find_if(evt.trks.begin(), evt.trks.end(), selector) != evt.trks.end();
854  }
template<int reconstruction_type>
bool JAANET::has_reconstructed_track ( const Evt &  evt,
const JRange< int >  range = JRange<int>() 
)
inline

Test whether given event has a track according selection.

Parameters
evtevent
rangerange of application types
Returns
true if at least one corresponding track; else false

Definition at line 864 of file JAAnetToolkit.hh.

865  {
866  return has_reconstructed_track(evt, has_history(reconstruction_type, range));
867  }
bool has_reconstructed_track(const Evt &evt, JTrackSelector_t selector)
Test whether given event has a track according selection.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
bool JAANET::has_reconstructed_muon ( const Evt &  evt)
inline

Test whether given event has a track with muon reconstruction.

Parameters
evtevent
Returns
true if at least one reconstructed muon; else false

Definition at line 875 of file JAAnetToolkit.hh.

876  {
878  }
bool has_reconstructed_track(const Evt &evt, JTrackSelector_t selector)
Test whether given event has a track according selection.
bool has_muon_fit(const Trk &track)
Test whether given track has default muon fit in history.
bool JAANET::has_reconstructed_shower ( const Evt &  evt)
inline

Test whether given event has a track with shower reconstruction.

Parameters
evtevent
Returns
true if at least one reconstructed shower; else false

Definition at line 886 of file JAAnetToolkit.hh.

887  {
889  }
bool has_reconstructed_track(const Evt &evt, JTrackSelector_t selector)
Test whether given event has a track according selection.
bool has_shower_fit(const Trk &track)
Test whether given track has default shower fit in history.
bool JAANET::has_reconstructed_aashower ( const Evt &  evt)
inline

Test whether given event has a track with aashower reconstruction.

Parameters
evtevent
Returns
true if at least one reconstructed shower; else false

Definition at line 897 of file JAAnetToolkit.hh.

898  {
900  }
bool has_reconstructed_track(const Evt &evt, JTrackSelector_t selector)
Test whether given event has a track according selection.
bool has_aashower_fit(const Trk &track)
Test whether given track has default shower fit in history.
template<class JTrackSelector_t , class JQualitySorter_t >
const Trk& JAANET::get_best_reconstructed_track ( const Evt &  evt,
JTrackSelector_t  selector,
JQualitySorter_t  comparator 
)
inline

Get best reconstructed track.


The track selector corresponds to the function operator bool selector(const Trk&); and the track comparator to bool comparator(const Trk&, const Trk&);.

Parameters
evtevent
selectortrack selector
comparatortrack comparator
Returns
track

Definition at line 913 of file JAAnetToolkit.hh.

916  {
917  std::vector<Trk>::const_iterator p = std::find_if(evt.trks.begin(), evt.trks.end(), selector);
918 
919  for (std::vector<Trk>::const_iterator i = p; i != evt.trks.end(); ++i) {
920  if (selector(*i) && comparator(*i, *p)) {
921  p = i;
922  }
923  }
924 
925  if (p != evt.trks.end())
926  return *p;
927  else
928  THROW(JIndexOutOfRange, "This event has no reconstructed track with given selector.");
929  }
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:633
template<int reconstruction_type>
const Trk& JAANET::get_best_reconstructed_track ( const Evt &  evt,
const JRange< int >  range = JRange<int>() 
)
inline

Get best reconstructed track.

Parameters
evtevent
rangerange of application types
Returns
track

Definition at line 939 of file JAAnetToolkit.hh.

940  {
941  return get_best_reconstructed_track(evt, has_history(reconstruction_type, range), quality_sorter<reconstruction_type>());
942  }
Reconstruction type dependent comparison of track quality.
const Trk & get_best_reconstructed_track(const Evt &evt, JTrackSelector_t selector, JQualitySorter_t comparator)
Get best reconstructed track.
bool has_history(const JFit &fit, const int type)
Test whether given fit has specified history.
Definition: JEvtToolkit.hh:262
const Trk& JAANET::get_best_reconstructed_muon ( const Evt &  evt)
inline

Get best reconstructed muon.

Parameters
evtevent
Returns
track

Definition at line 951 of file JAAnetToolkit.hh.

952  {
954  }
Reconstruction type dependent comparison of track quality.
bool has_muon_fit(const Trk &track)
Test whether given track has default muon fit in history.
const Trk & get_best_reconstructed_track(const Evt &evt, JTrackSelector_t selector, JQualitySorter_t comparator)
Get best reconstructed track.
const Trk& JAANET::get_best_reconstructed_shower ( const Evt &  evt)
inline

Get best reconstructed shower.

Parameters
evtevent
Returns
track

Definition at line 962 of file JAAnetToolkit.hh.

963  {
965  }
Reconstruction type dependent comparison of track quality.
const Trk & get_best_reconstructed_track(const Evt &evt, JTrackSelector_t selector, JQualitySorter_t comparator)
Get best reconstructed track.
bool has_shower_fit(const Trk &track)
Test whether given track has default shower fit in history.
const Trk& JAANET::get_best_reconstructed_aashower ( const Evt &  evt)
inline

Get best reconstructed aashower.

Parameters
evtevent
Returns
track

Definition at line 973 of file JAAnetToolkit.hh.

974  {
976  }
Reconstruction type dependent comparison of track quality.
const Trk & get_best_reconstructed_track(const Evt &evt, JTrackSelector_t selector, JQualitySorter_t comparator)
Get best reconstructed track.
bool has_aashower_fit(const Trk &track)
Test whether given track has default shower fit in history.
void JAANET::copy ( const Head &  from,
JHead &  to 
)

Copy header from from to to.

Parameters
fromheader
toheader

Definition at line 40 of file JHead.cc.

41  {
42  using namespace std;
43  using namespace JPP;
44 
46  JRootClassReader cls(to);
47 
48  for (Head::const_iterator i = from.begin(); i != from.end(); ++i) {
49 
50  const JRootClassReader& abc = cls.find(getTag(i->first).c_str());
51 
52  if (abc.is_valid()) {
53 
54  JRedirectString redirect(reader, i->second);
55 
56  reader.getObject(abc);
57  }
58 
59  to.insert(*i); // keep track of parsed tags
60  }
61  }
const char * c_str() const
C-string.
Definition: JTag.hh:195
JNET::JTag getTag(JLANG::JType< KM3NETDAQ::JDAQTimeslice >)
Definition: JDAQTags.hh:80
T & getInstance(const T &object)
Get static instance from temporary object.
Definition: JObject.hh:73
static JNullStream null
Null I/O stream.
Definition: JNullStream.hh:51
const JEquationParameters & getEquationParameters()
Get equation parameters corresponding to Monte Carlo ASCII format, i.e:
void JAANET::copy ( const JHead &  from,
Head &  to 
)

Copy header from from to to.

Parameters
fromheader
toheader

Definition at line 70 of file JHead.cc.

71  {
72  using namespace std;
73  using namespace JPP;
74 
75  to.clear();
76 
77  ostringstream os;
78 
79  JRootWriter writer(os, getEquationParameters(), JAAnetDictionary::getInstance());
80  JRootClassWriter cls(from);
81 
82  TIterator* i = cls.getClass()->GetListOfDataMembers()->MakeIterator();
83 
84  for (const TDataMember* p; (p = (const TDataMember*) i->Next()) != NULL; os.str("")) {
85 
86  if (JRootClass::is_class(*p) && strcmp(p->GetName(), end_event::Class_Name()) != 0) {
87 
88  if (from.find(p->GetName()) != from.end()) { // only copy data of parsed tags
89 
90  writer.putObject(cls.get(*p));
91 
92  to.insert(make_pair(p->GetName(), os.str()));
93  }
94  }
95  }
96 
97  // copy pending data
98 
99  for (JHead::const_iterator i = from.begin(); i != from.end(); ++i) {
100  if (to.find(i->first) == to.end()) {
101  to.insert(*i);
102  }
103  }
104  }
T & getInstance(const T &object)
Get static instance from temporary object.
Definition: JObject.hh:73
const JEquationParameters & getEquationParameters()
Get equation parameters corresponding to Monte Carlo ASCII format, i.e:
std::string JAANET::getTag ( const std::string &  tag)
inline

Get tag without aanet extension "_<counter>" for identical tags.

Parameters
tagtag
Returns
tag

Definition at line 57 of file JHead.hh.

58  {
59  using namespace std;
60 
61  const string::size_type pos = tag.find(AANET_TAG_SEPARATOR);
62 
63  if (pos != string::npos) {
64 
65  for (string::size_type i = pos + 1; i != tag.size(); ++i) {
66  if (!isdigit(tag[i])) {
67  return tag;
68  }
69  }
70 
71  return tag.substr(0, pos);
72  }
73 
74  return tag;
75  }
static const char AANET_TAG_SEPARATOR
Separator for AAnet tag extension for multiple tags (&quot;_&lt;counter&gt;&quot;).
Definition: JHead.hh:28
template<class T >
bool JAANET::is_valid ( const T &  value)
inline

Check validity of given value.

Parameters
valuevalue
Returns
true if given value not-equal to default value; else false

Definition at line 711 of file JHead.hh.

712  {
713  return !value.equals(T());
714  }
bool JAANET::operator== ( const Head &  first,
const Head &  second 
)
inline

Equal operator.

Parameters
firstfirst header
secondsecond header
Returns
true if two headers are equal; else false

Definition at line 885 of file JHead.hh.

887  {
888  return JHead(first).equals(JHead(second));
889  }
Monte Carlo run header.
Definition: JHead.hh:727
bool equals(const JHead &header) const
Test compatibility of headers.
Definition: JHead.hh:784
bool JAANET::operator< ( const Head &  first,
const Head &  second 
)
inline

Less than operator.

Parameters
firstfirst header
secondsecond header
Returns
true if first header is less than second header; else false

Definition at line 899 of file JHead.hh.

901  {
902  return JHead(first).less(JHead(second));
903  }
bool less(const JHead &header) const
Comparison of headers.
Definition: JHead.hh:772
Monte Carlo run header.
Definition: JHead.hh:727
template<class T >
void JAANET::push ( JHead &  header,
T JHead::*  p 
)
inline

Push data of JHead for subsequent copy to Head.

Parameters
headerheader
ppointer to data member to be pushed

Definition at line 102 of file JHeadToolkit.hh.

103  {
104  using namespace JPP;
105 
106  const int offset = (const char*) &(header.*p) - (const char*) &header;
107 
108  const JRootClass cls(getType<JHead>());
109 
110  TIterator* i = cls.getClass()->GetListOfDataMembers()->MakeIterator();
111 
112  for (const TDataMember* p; (p = (const TDataMember*) i->Next()) != NULL; ) {
113 
114  const JRootClass abc = cls.get(*p);
115 
116  if (abc == JRootClass(JType<T>()) && offset == abc.getOffset()) {
117  header[p->GetName()] = "";
118  }
119  }
120  }
int getOffset() const
Get offset of this class with respect to parent class.
Definition: JRootClass.hh:156
Auxiliary class for a type holder.
Definition: JType.hh:19
JRootClass get(const TDataMember &data_member) const
Get ROOT class of given data member.
Definition: JRootClass.hh:266
Auxiliary class to manage access to base classes and data members of ROOT class.
Definition: JRootClass.hh:38
bool JAANET::is_physics ( const JHead &  header,
const std::string &  generator 
)
inline

Check for generator.

Parameters
headerheader
generatorgenerator
Returns
true if this Monte Carlo is produced by given generator; else false

Definition at line 130 of file JHeadToolkit.hh.

131  {
132  for (std::vector<physics>::const_iterator i = header.physics.begin(); i != header.physics.end(); ++i) {
133  if (i->program == generator) {
134  return true;
135  }
136  }
137 
138  return false;
139  }
bool JAANET::is_sirene ( const JHead &  header)
inline

Check for generator.

Parameters
headerheader
Returns
true if this Monte Carlo is produced by JSirene; else false

Definition at line 148 of file JHeadToolkit.hh.

149  {
150  return is_physics(header, JHead::JSIRENE);
151  }
bool is_physics(const JHead &header, const std::string &generator)
Check for generator.
bool JAANET::is_km3 ( const JHead &  header)
inline

Check for generator.

Parameters
headerheader
Returns
true if this Monte Carlo is produced by km3; else false

Definition at line 160 of file JHeadToolkit.hh.

161  {
162  return is_physics(header, JHead::KM3);
163  }
bool is_physics(const JHead &header, const std::string &generator)
Check for generator.
bool JAANET::is_km3sim ( const JHead &  header)
inline

Check for generator.

Parameters
headerheader
Returns
true if this Monte Carlo is produced by KM3Sim; else false

Definition at line 172 of file JHeadToolkit.hh.

173  {
174  return !is_km3(header) && !is_sirene(header);
175  }
bool is_sirene(const JHead &header)
Check for generator.
bool is_km3(const JHead &header)
Check for generator.
bool JAANET::is_genhen ( const JHead &  header)
inline

Check for generator.

Parameters
headerheader
Returns
true if this Monte Carlo is produced by genhen; else false

Definition at line 184 of file JHeadToolkit.hh.

185  {
186  return is_physics(header, JHead::GENHEN);
187  }
bool is_physics(const JHead &header, const std::string &generator)
Check for generator.
bool JAANET::is_genie ( const JHead &  header)
inline

Check for generator.

Parameters
headerheader
Returns
true if this Monte Carlo is produced by genie; else false

Definition at line 196 of file JHeadToolkit.hh.

197  {
198  return is_physics(header, JHead::GENIE) || is_physics(header, JHead::GSEAGEN);
199  }
bool is_physics(const JHead &header, const std::string &generator)
Check for generator.
bool JAANET::is_mupage ( const JHead &  header)
inline

Check for generator.

Parameters
headerheader
Returns
true if this Monte Carlo is produced by mupage; else false

Definition at line 208 of file JHeadToolkit.hh.

209  {
210  return is_physics(header, JHead::MUPAGE);
211  }
bool is_physics(const JHead &header, const std::string &generator)
Check for generator.
template<class T >
T JAANET::get ( const JHead &  head)
inline

Get object from header.

Parameters
headheader
Returns
object

Get object from header.

Parameters
headerheader
Returns
position

Get object from header.

Parameters
headheader
Returns
cylinder

Definition at line 231 of file JHeadToolkit.hh.

232  {
233  if (is_sirene(header) || is_km3(header)) {
234 
235  return header.coord_origin;
236 
237  } else {
238 
239  if (is_valid(header.can))
240  return Vec(0.0, 0.0, -header.can.zmin);
241  else if (is_valid(header.coord_origin))
242  return header.coord_origin;
243  else
244  return Vec(0.0, 0.0, 0.0);
245  }
246  }
bool is_sirene(const JHead &header)
Check for generator.
bool is_valid(const T &value)
Check validity of given value.
Definition: JHead.hh:711
bool is_km3(const JHead &header)
Check for generator.
template<>
Vec JAANET::get ( const JHead &  header)
inline

Get position offset of detector due to generator.

Get object from header.

Parameters
headerheader
Returns
position

Definition at line 231 of file JHeadToolkit.hh.

232  {
233  if (is_sirene(header) || is_km3(header)) {
234 
235  return header.coord_origin;
236 
237  } else {
238 
239  if (is_valid(header.can))
240  return Vec(0.0, 0.0, -header.can.zmin);
241  else if (is_valid(header.coord_origin))
242  return header.coord_origin;
243  else
244  return Vec(0.0, 0.0, 0.0);
245  }
246  }
bool is_sirene(const JHead &header)
Check for generator.
bool is_valid(const T &value)
Check validity of given value.
Definition: JHead.hh:711
bool is_km3(const JHead &header)
Check for generator.
template<>
JPosition3D JAANET::get ( const JHead &  header)
inline

Get position offset of detector due to generator.

Get object from header.

Parameters
headerheader
Returns
position

Definition at line 256 of file JHeadToolkit.hh.

257  {
258  const Vec pos = get<Vec>(header);
259 
260  return JPosition3D(pos.x, pos.y, pos.z);
261  }
template<>
JCylinder3D JAANET::get ( const JHead &  head)
inline

Get cylinder corresponding to can.

Get object from header.

Parameters
headheader
Returns
cylinder

Definition at line 271 of file JHeadToolkit.hh.

272  {
273  using namespace JGEOMETRY2D;
274 
275  return JCylinder3D(JCircle2D(JVector2D(), head.can.r),
276  head.can.zmin,
277  head.can.zmax);
278  }
Data structure for vector in two dimensions.
Definition: JVector2D.hh:30
Data structure for circle in two dimensions.
Definition: JCircle2D.hh:29
const JEquationParameters& JAANET::getEquationParameters ( )
inline

Get equation parameters corresponding to Monte Carlo ASCII format, i.e:

   <key>: <value [<value>]*
   <key>: <value [<value>]*
Returns
equation parameters

Definition at line 325 of file JHeadToolkit.hh.

326  {
327  static const JEquationParameters parameters(":", "\n", "", "");
328 
329  return parameters;
330  }
Simple data structure to support I/O of equations (see class JLANG::JEquation).
std::istream& JAANET::operator>> ( std::istream &  in,
JHead &  header 
)
inline

Read header from input.

Parameters
ininput stream
headerheader
Returns
input stream

Definition at line 340 of file JHeadToolkit.hh.

341  {
342  using namespace JLANG;
343  using namespace JROOT;
344  using namespace std;
345 
347 
348  JRootClassReader cls(header);
349 
350  for (JEquation equation; reader >> equation && equation.getKey() != end_event::Class_Name(); ) {
351 
352  JRedirectString redirect(reader, equation.getValue());
353 
354  const JRootClassReader abc = cls.find(equation.getKey().c_str());
355 
356  if (abc.is_valid()) {
357  reader.getObject(abc);
358  }
359  }
360 
361  return in;
362  }
const std::string & getKey() const
Get key.
Definition: JEquation.hh:150
General purpose equation class.
Definition: JEquation.hh:47
Implementation for ASCII input of objects with ROOT dictionary.
T & getInstance(const T &object)
Get static instance from temporary object.
Definition: JObject.hh:73
ROOT class for reading object.
Definition: JRootClass.hh:478
JRootAddressableClass find(const char *name) const
Find addressable base class or data member with given name within current class.
Definition: JRootClass.hh:433
This class can be used to temporarily redirect an input stream to an input string.
const JEquationParameters & getEquationParameters()
Get equation parameters corresponding to Monte Carlo ASCII format, i.e:
std::ostream& JAANET::operator<< ( std::ostream &  out,
const JHead &  header 
)
inline

Write header to output.

Parameters
outoutput stream
headerheader
Returns
output stream

Definition at line 372 of file JHeadToolkit.hh.

373  {
374  using namespace JLANG;
375  using namespace JROOT;
376 
378 
379  JRootClassWriter cls(header);
380 
381  TIterator* i = cls.getClass()->GetListOfDataMembers()->MakeIterator();
382 
383  for (const TDataMember* p; (p = (const TDataMember*) i->Next()) != NULL; ) {
384  if (!JRootClass::is_static(*p)) {
385  writer.put(p->GetName(), cls.get(*p), true);
386  }
387  }
388 
389  return out;
390  }
Implementation for ASCII output of objects with ROOT dictionary.
ROOT class for writing object.
Definition: JRootClass.hh:568
T & getInstance(const T &object)
Get static instance from temporary object.
Definition: JObject.hh:73
const JEquationParameters & getEquationParameters()
Get equation parameters corresponding to Monte Carlo ASCII format, i.e:
bool JAANET::operator< ( const JHit_t first,
const JHit_t second 
)
inline

Less than operator for hits.

First hit is defined as:

  1. smallest PMT identifier
  2. earliest time if same PMT identifier
Parameters
firstfirst hit
secondsecond hit
Returns
true if first hit earlier than second hit; else false

Definition at line 92 of file JHit_t.hh.

93  {
94  if (first.pmt_id == second.pmt_id)
95  return first.t < second.t;
96  else
97  return first.pmt_id < second.pmt_id;
98  }
double t
Definition: JHit_t.hh:76

Variable Documentation

const int JAANET::AASHOWER_RECONSTRUCTION_TYPE = 101
static

AAshower reconstruction type for AAnet.

Definition at line 62 of file JAAnetToolkit.hh.

const JEvtEvaluator JAANET::getEvtValue
static

Function object for evaluation of DAQ objects.

Definition at line 48 of file JEvtEvaluator.hh.

const char JAANET::AANET_TAG_SEPARATOR = '_'
static

Separator for AAnet tag extension for multiple tags ("_<counter>").

Definition at line 28 of file JHead.hh.