Jpp  17.2.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JVisibleEnergyToolkit.hh
Go to the documentation of this file.
1 #ifndef __JSIRENE__JVISIBLEENERGYTOOLKIT__
2 #define __JSIRENE__JVISIBLEENERGYTOOLKIT__
3 
7 
10 #include "JAAnet/JPDB.hh"
11 
12 #include "JGeometry2D/JVector2D.hh"
13 #include "JGeometry2D/JCircle2D.hh"
14 
16 
17 #include "JPhysics/JConstants.hh"
18 #include "JPhysics/JGeane.hh"
19 
20 #include "JSirene/pythia.hh"
21 
22 
23 /**
24  * \file
25  *
26  * Auxiliary methods for evaluating visible energies.
27  * \author bjung
28  */
29 namespace JSIRENE {}
30 namespace JPP { using namespace JSIRENE; }
31 
32 namespace JSIRENE {
33 
34  using JAANET::JPDB;
35  using JAANET::is_muon;
36  using JAANET::getAxis;
37  using JAANET::getPosition;
39 
42 
44 
45  using JPHYSICS::geanc;
46  using JPHYSICS::gWater;
47  using JPHYSICS::MASS_MUON;
50 
51 
52  /**
53  * Get the visible energy of a track.\n
54  * This method accounts for muon radiative energy losses.
55  *
56  * \param track track
57  * \param can detector can
58  * \return visible energy [GeV]
59  */
60  inline Vec getVisibleEnergy(const Trk& track,
61  const JCylinder3D& can) {
62  using namespace std;
63  using namespace JPP;
64 
65  double Evis = 0.0;
66 
67  if (track.is_finalstate()) {
68  if (is_muon(track)) {
69 
70  // Determine muon pathlength inside detector [m]
71 
72  const JCylinder3D::intersection_type& intersection = can.getIntersection(getAxis(track));
73 
74  const double Lmuon = gWater.getX(track.E, MASS_MUON / getSinThetaC());
75  const double Leff = (intersection.first < 0.0 ?
76  min(Lmuon, intersection.second) :
77  min(Lmuon, intersection.second) - intersection.first);
78 
79  // Determine visible energy deposition [GeV]
80 
81  const double dEb = gWater.getEb(track.E, Leff);
82  const double dEc = Leff / geanc();
83 
84  Evis = dEb + dEc;
85 
86  } else if (!is_neutrino(track) && JPDB::getInstance().hasPDG(track.type)) {
87 
88  Evis = pythia(track.type, getKineticEnergy(track));
89  }
90  }
91 
92  return Evis * track.dir / track.dir.len();
93  }
94 
95 
96  /**
97  * Get the visible energy of a track, assuming an infinite detector volume.\n
98  * This method accounts for muon radiative energy losses.
99  *
100  * \param track track
101  */
102  inline Vec getVisibleEnergy(const Trk& track)
103  {
104  static const double Rearth = 6386 * 1e3; // Radius of the Earth [m]
105 
106  static const JCylinder3D can(JCircle2D(JVector2D(), Rearth), Rearth, Rearth);
107 
108  return getVisibleEnergy(track, can);
109  }
110 
111 
112  /**
113  * Get the visible energy of an event.\n
114  * This method accounts for muon radiative energy losses.
115  *
116  * \param evt event
117  * \param can detector can
118  * \return visible energy [GeV]
119  */
120  inline Vec getVisibleEnergy(const Evt& evt,
121  const JCylinder3D& can)
122  {
123  Vec Evis(0.0, 0.0, 0.0);
124 
125  for (std::vector<Trk>::const_iterator track = evt.mc_trks.begin(); track != evt.mc_trks.end(); ++track) {
126  Evis += getVisibleEnergy(*track, can);
127  }
128 
129  return Evis;
130  }
131 
132 
133  /**
134  * Get the visible energy of an event, assuming an infinite detector volume.\n
135  * This method accounts for muon radiative energy losses.
136  *
137  * \param evt event
138  * \return visible energy [GeV]
139  */
140  inline Vec getVisibleEnergy(const Evt& evt)
141  {
142  Vec Evis(0.0, 0.0, 0.0);
143 
144  for (std::vector<Trk>::const_iterator track = evt.mc_trks.cbegin(); track != evt.mc_trks.cend(); ++track) {
145  Evis += getVisibleEnergy(*track);
146  }
147 
148  return Evis;
149  }
150 }
151 
152 #endif
Data structure for vector in two dimensions.
Definition: JVector2D.hh:32
Data structure for circle in two dimensions.
Definition: JCircle2D.hh:33
double geanc()
Equivalent muon track length per unit shower energy.
Definition: JGeane.hh:28
Energy loss of muon.
This file contains converted Fortran code from km3.
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
Vec dir
track direction
Definition: Trk.hh:18
static const double MASS_MUON
muon mass [GeV]
static const JGeaneWater gWater
Function object for energy loss of muon in sea water.
Definition: JGeane.hh:381
virtual double getX(const double E0, const double E1) const override
Get distance traveled by muon.
Definition: JGeane.hh:349
static const JPythia pythia
Function object for relative light yield as a function of GEANT particle code.
Definition: JPythia.hh:96
intersection_type getIntersection(const JAxis3D &axis) const
Get intersection points of axis with cylinder.
Definition: JCylinder3D.hh:278
double E
Energy [GeV] (either MC truth or reconstructed)
Definition: Trk.hh:20
Collection of particles.
Definition: JPDB.hh:109
double getEb(const double E, const double dx) const
Get energy loss due to pair production and bremsstrahlung.
Definition: JGeane.hh:334
double len() const
Get length.
Definition: Vec.hh:145
static const JPDB & getInstance()
Get particle data book.
Definition: JPDB.hh:125
bool is_neutrino(const Trk &track)
Test whether given track is a neutrino.
The Vec class is a straightforward 3-d vector, which also works in pyroot.
Definition: Vec.hh:12
Cylinder object.
Definition: JCylinder3D.hh:39
double getKineticEnergy(const double E, const double m)
Get kinetic energy of particle with given mass.
JAxis3D getAxis(const Trk &track)
Get axis.
JPosition3D getPosition(const Vec &pos)
Get position.
Definition of hit and track types and auxiliary methods for handling Monte Carlo data.
Physics constants.
int type
MC: particle type in PDG encoding.
Definition: Trk.hh:24
bool is_finalstate() const
Test whether given particle is a final state inside the detector.
Definition: Trk.hh:88
Definition of particle types.
Function object for energy dependent energy loss of the muon.
Definition: JGeane.hh:206
Vec getVisibleEnergy(const Trk &track, const JCylinder3D &can)
Get the visible energy of a track.
The Trk class represents a Monte Carlo (MC) particle as well as a reconstructed track/shower.
Definition: Trk.hh:14
std::vector< Trk > mc_trks
MC: list of MC truth tracks.
Definition: Evt.hh:49
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:20
double getSinThetaC()
Get average sine of Cherenkov angle of water corresponding to group velocity.