GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GAstroFlux.h
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*!
3 
4 \class genie::flux::GAstroFlux
5 
6 \brief A base class for the concrete astrophysical neutrino flux drivers.
7 
8  <<<< NOTE: CODE UNDER DEVELOPMENT >>>>
9 
10  COORDINATE SYSTEMS / NEUTRINO GENERATION & PROPAGATION :
11 
12  Neutrinos are generated on a sphere with radius R_{earth}.
13  Especially,
14  - For diffuse fluxes:
15  Neutrinos can be generated anywhere on that surface.
16  - For point sources:
17  Neutrinos are generated at fixed right ascension and declination.
18  Then time is randomized, to account for Earth's rotation, and the
19  Equatorial Coordinates are converted to GEF. So, neutrinos from each
20  point source are generated on circles parallel to the Earth's Equator.
21 
22  Initially, neutrino coordinates are generated in the Geocentric Earth-
23  Fixed (GEF) Coordinate System (later to be converted to the appropriate
24  detector coordinate system - See further below).
25 
26  * Definition:
27  Geocentric Earth-Fixed (GEF) Coordinate System
28  +z: Points to North Pole.
29  xy: Equatorial plane.
30  +x: Points to the Prime Meridian.
31  +y: As needed to make a right-handed coordinate system.
32 
33  Neutrinos are then propagated towards the detector.
34  The Earth opaqueness to ultra high energy neutrinos is taken into
35  account. The Earth density profile is modelled using the PREM
36  (Preliminary Earth Model, The Encyclopedia of Solid Earth Geophysics,
37  David E. James, ed., Van Nostrand Reinhold, New York, 1989, p.331).
38 
39  The detector position is determined in the Spherical/Geographic System
40  by its geographic latitude (angle relative to Equator), its geographic
41  longitude (angle relative to Prime Meridian) and its depth from the
42  surface.
43 
44  The generated flux neutrinos, propagated through the Earth towards the
45  detector) are then positioned on the surface of a sphere with radius Rd
46  which should fully enclose the neutrino detector. The centre of that
47  sphere is taken to be the origin of the detector coordinate system.
48  The transverse coords are appropriately randomized so that neutrinos
49  from any given direction bath the entire sphere enclosing the detector.
50 
51  The final flux neutrino coordinates are given in the detector coordinate
52  system. The default detector coordinate system is the Topocentric Horizontal
53  (THZ) Coordinate System. Alternative user-defined topocentric systems can
54  be defined by specifying the appropriate rotation from THZ.
55 
56  * Definition:
57  Topocentric Horizontal (THZ) Coordinate System
58  (default detector coordinate system)
59  +z: Points towards the local zenith.
60  +x: On same plane as local meridian, pointing south.
61  +y: As needed to make a right-handed coordinate system.
62  origin: detector centre
63 
64  WEIGHTING SCHEMES:
65 
66  For a detector with geometrical cross section ~ 1km^2, the solid
67  angle acceptance changes by ~10 orders of magnitude across the
68  surface of the Earth.
69 
70  The driver supports both weighted and un-weighted flux generation
71  schemes. However, because of the enormous changes in solid angle
72  acceptance and energy, only the weighted scheme is practical.
73 
74  PHYSICS:
75 
76  The relative neutrino population needs to be set by the user using
77  the SetRelNuPopulations() method.
78  If run without arguments, the following relative populations are set:
79  nue:numu:nutau:nuebar:numubar:nutaubar = 1:2:0:1:2:0
80 
81  The energy spectrum is follows a power law. The user needs to
82  specify the power-law index by calling SetEnergyPowLawIdx().
83 
84 \author Costas Andreopoulos <c.andreopoulos \at cern.ch>
85  University of Liverpool
86 
87 \created March 27, 2010
88 
89 \cpright Copyright (c) 2003-2024, The GENIE Collaboration
90  For the full text of the license visit http://copyright.genie-mc.org
91 */
92 //____________________________________________________________________________
93 
94 #ifndef _GASTRO_FLUX_H_
95 #define _GASTRO_FLUX_H_
96 
97 #include <string>
98 #include <map>
99 
100 #include <TLorentzVector.h>
101 #include <TVector3.h>
102 #include <TRotation.h>
103 
106 
107 class TH1D;
108 class TH2D;
109 
110 using std::string;
111 using std::map;
112 
113 namespace genie {
114 namespace flux {
115 
116 const double kAstroDefMaxEv = 1E+20 * units::GeV; ///<
117 const double kAstroDefMinEv = 1E-3 * units::GeV; ///<
118 const int kAstroNlog10EvBins = 1000; ///<
119 const int kAstroNCosThetaBins = 500; ///<
120 const int kAstroNPhiBins = 500; ///<
121 
122 //
123 //
124 //
125 
126 class GAstroFlux: public GFluxI {
127 
128 public :
129  virtual ~GAstroFlux();
130 
131  //
132  // methods implementing the GENIE GFluxI interface
133  //
134  virtual const PDGCodeList & FluxParticles (void) { return *fPdgCList; }
135  virtual double MaxEnergy (void);
136  virtual bool GenerateNext (void);
137  virtual int PdgCode (void) { return fgPdgC; }
138  virtual double Weight (void) { return fgWeight; }
139  virtual const TLorentzVector & Momentum (void) { return fgP4; }
140  virtual const TLorentzVector & Position (void) { return fgX4; }
141  virtual bool End (void) { return false; }
142  virtual long int Index (void) { return -1; }
143  virtual void Clear (Option_t * opt);
144  virtual void GenerateWeighted (bool gen_weighted);
145 
146  //
147  // configuration methods specific to all astrophysical neutrino flux drivers
148  //
149  void ForceMinEnergy (double emin);
150  void ForceMaxEnergy (double emax);
151  void SetDetectorPosition (double latitude, double longitude, double depth, double size);
152  void SetRelNuPopulations (double nnue=1, double nnumu=2, double nnutau=0, double nnuebar=1, double nnumubar=2, double nnutaubar=0);
153  void SetEnergyPowLawIdx (double n);
154  void SetUserCoordSystem (TRotation & rotation); ///< rotation Topocentric Horizontal -> User-defined Topocentric Coord System
155 
156 protected:
157 
158  class NuGenerator;
159  class NuPropagator;
160 
161  // abstract class, ctor hidden
162  GAstroFlux();
163 
164  // protected methods
165  void Initialize (void);
166  void CleanUp (void);
167  void ResetSelection (void);
168 
169  //
170  // protected data members
171  //
172 
173  PDGCodeList * fPdgCList; ///< declared list of neutrino pdg-codes that can be thrown by current instance
174  int fgPdgC; ///< (current) generated nu pdg-code
175  TLorentzVector fgP4; ///< (current) generated nu 4-momentum
176  TLorentzVector fgX4; ///< (current) generated nu 4-position
177  double fgWeight; ///< (current) generated nu weight
178  // configuration properties set by the user
179  double fMaxEvCut; ///< (config) user-defined maximum energy cut
180  double fMinEvCut; ///< (config) user-defined minimum energy cut
181  bool fGenWeighted; ///< (config) generate a weighted or unweighted flux?
182  double fDetGeoLatitude; ///< (config) detector: geographic latitude
183  double fDetGeoLongitude; ///< (config) detector: geographic longitude
184  double fDetGeoDepth; ///< (config) detector: depth from surface
185  double fDetSize; ///< (config) detector: size (detector should be enclosed in sphere of this radius)
186  map<int,double> fRelNuPopulations; ///< (config) relative neutrino populations
187  TRotation fRotGEF2THz; ///< (config) coord. system rotation: GEF translated to detector centre -> THZ
188  TRotation fRotTHz2User; ///< (config) coord. system rotation: THZ -> Topocentric user-defined
189  // internal flags and utility objects
190  TVector3 fDetCenter; ///<
191  TH1D * fEnergySpectrum; ///<
195 
196  //
197  // utility classes
198  //
199  class NuGenerator {
200  public:
203  bool SelectNuPdg (bool weighted, const map<int,double> & nupdgpdf, int & nupdg, double & wght);
204  bool SelectEnergy(bool weighted, TH1D & log10epdf, double log10emin, double log10emax, double & log10e, double & wght);
205  bool SelectOrigin(bool weighted, TH2D & opdf, double & phi, double & costheta, double & wght);
206  };
207  class NuPropagator {
208  public:
209  NuPropagator(double stepsz) : fStepSize(stepsz/units::km) { }
211  bool Go(double phi_start, double costheta_start, const TVector3 & detector_centre, double detector_sz, int nu_pdg, double Ev);
212  int NuPdgAtDetVolBoundary (void) { return fNuPdg; }
213  TVector3 & X3AtDetVolBoundary (void) { return fX3; }
214  TVector3 & P3AtDetVolBoundary (void) { return fP3; }
215  private:
216  double fStepSize;
217  int fNuPdg;
218  TVector3 fX3;
219  TVector3 fP3;
220  };
221 
222 };
223 
224 
225 //
226 // Concrete astrophysical flux drivers
227 //
228 
229 //............................................................................
230 // GENIE diffuse astrophysical neutrino flux driver
231 //
233 public :
236 
237  //
238  //
239 };
240 
241 //............................................................................
242 // GENIE concrete flux driver for astrophysical point neutrino sources
243 //
245 public :
248 
249  //
250  bool GenerateNext (void);
251 
252  //
253  void AddPointSource(string name, double ra, double dec, double rel_intensity);
254 
255 private:
256 
257  bool SelectSource(void);
258 
259  map<int, string> fPntSrcName; ///< point source name
260  map<int, double> fPntSrcRA; ///< right ascension
261  map<int, double> fPntSrcDec; ///< declination
262  map<int, double> fPntSrcRelI; ///< relative intensity
263  double fPntSrcTotI; ///< sum of all relative intensities
264 
265  unsigned int fSelSourceId;
266 };
267 //............................................................................
268 
269 } // flux namespace
270 } // genie namespace
271 
272 #endif // _GASTRO_FLUX_H_
virtual const TLorentzVector & Momentum(void)
returns the flux neutrino 4-momentum
Definition: GAstroFlux.h:139
void ForceMinEnergy(double emin)
Definition: GAstroFlux.cxx:138
TLorentzVector fgP4
(current) generated nu 4-momentum
Definition: GAstroFlux.h:175
virtual void Clear(Option_t *opt)
reset state variables based on opt
Definition: GAstroFlux.cxx:150
double fDetSize
(config) detector: size (detector should be enclosed in sphere of this radius)
Definition: GAstroFlux.h:185
double fMinEvCut
(config) user-defined minimum energy cut
Definition: GAstroFlux.h:180
bool fGenWeighted
(config) generate a weighted or unweighted flux?
Definition: GAstroFlux.h:181
virtual double MaxEnergy(void)
declare the max flux neutrino energy that can be generated (for init. purposes)
Definition: GAstroFlux.cxx:42
void SetEnergyPowLawIdx(double n)
Definition: GAstroFlux.cxx:258
TRotation fRotTHz2User
(config) coord. system rotation: THZ -&gt; Topocentric user-defined
Definition: GAstroFlux.h:188
double fMaxEvCut
(config) user-defined maximum energy cut
Definition: GAstroFlux.h:179
void SetUserCoordSystem(TRotation &rotation)
rotation Topocentric Horizontal -&gt; User-defined Topocentric Coord System
Definition: GAstroFlux.cxx:281
const double kAstroDefMaxEv
Definition: GAstroFlux.h:116
void ForceMaxEnergy(double emax)
Definition: GAstroFlux.cxx:144
map< int, double > fRelNuPopulations
(config) relative neutrino populations
Definition: GAstroFlux.h:186
static constexpr double km
Definition: Units.h:64
map< int, double > fPntSrcRelI
relative intensity
Definition: GAstroFlux.h:262
map< int, string > fPntSrcName
point source name
Definition: GAstroFlux.h:259
const int kAstroNPhiBins
Definition: GAstroFlux.h:120
virtual bool GenerateNext(void)
generate the next flux neutrino (return false in err)
Definition: GAstroFlux.cxx:47
A list of PDG codes.
Definition: PDGCodeList.h:32
bool SelectNuPdg(bool weighted, const map< int, double > &nupdgpdf, int &nupdg, double &wght)
Definition: GAstroFlux.cxx:363
virtual const PDGCodeList & FluxParticles(void)
declare list of flux neutrinos that can be generated (for init. purposes)
Definition: GAstroFlux.h:134
double fPntSrcTotI
sum of all relative intensities
Definition: GAstroFlux.h:263
double fDetGeoDepth
(config) detector: depth from surface
Definition: GAstroFlux.h:184
TRotation fRotGEF2THz
(config) coord. system rotation: GEF translated to detector centre -&gt; THZ
Definition: GAstroFlux.h:187
NuGenerator * fNuGen
Definition: GAstroFlux.h:193
virtual double Weight(void)
returns the flux neutrino weight (if any)
Definition: GAstroFlux.h:138
static constexpr double GeV
Definition: Units.h:28
double fDetGeoLongitude
(config) detector: geographic longitude
Definition: GAstroFlux.h:183
bool SelectEnergy(bool weighted, TH1D &log10epdf, double log10emin, double log10emax, double &log10e, double &wght)
Definition: GAstroFlux.cxx:411
virtual long int Index(void)
returns corresponding index for current flux neutrino (e.g. for a flux ntuple returns the current ent...
Definition: GAstroFlux.h:142
double fDetGeoLatitude
(config) detector: geographic latitude
Definition: GAstroFlux.h:182
void SetDetectorPosition(double latitude, double longitude, double depth, double size)
Definition: GAstroFlux.cxx:162
void AddPointSource(string name, double ra, double dec, double rel_intensity)
Definition: GAstroFlux.cxx:574
PDGCodeList * fPdgCList
declared list of neutrino pdg-codes that can be thrown by current instance
Definition: GAstroFlux.h:173
int fgPdgC
(current) generated nu pdg-code
Definition: GAstroFlux.h:174
virtual const TLorentzVector & Position(void)
returns the flux neutrino 4-position (note: expect SI rather than physical units) ...
Definition: GAstroFlux.h:140
bool SelectOrigin(bool weighted, TH2D &opdf, double &phi, double &costheta, double &wght)
Definition: GAstroFlux.cxx:446
bool Go(double phi_start, double costheta_start, const TVector3 &detector_centre, double detector_sz, int nu_pdg, double Ev)
Definition: GAstroFlux.cxx:473
virtual void GenerateWeighted(bool gen_weighted)
set whether to generate weighted or unweighted neutrinos
Definition: GAstroFlux.cxx:157
map< int, double > fPntSrcDec
declination
Definition: GAstroFlux.h:261
map< int, double > fPntSrcRA
right ascension
Definition: GAstroFlux.h:260
const int kAstroNCosThetaBins
Definition: GAstroFlux.h:119
const int kAstroNlog10EvBins
Definition: GAstroFlux.h:118
TLorentzVector fgX4
(current) generated nu 4-position
Definition: GAstroFlux.h:176
A base class for the concrete astrophysical neutrino flux drivers.
Definition: GAstroFlux.h:126
NuPropagator * fNuPropg
Definition: GAstroFlux.h:194
bool GenerateNext(void)
generate the next flux neutrino (return false in err)
Definition: GAstroFlux.cxx:569
const double kAstroDefMinEv
Definition: GAstroFlux.h:117
virtual bool End(void)
true if no more flux nu&#39;s can be thrown (eg reaching end of beam sim ntuples)
Definition: GAstroFlux.h:141
virtual int PdgCode(void)
returns the flux neutrino pdg code
Definition: GAstroFlux.h:137
void SetRelNuPopulations(double nnue=1, double nnumu=2, double nnutau=0, double nnuebar=1, double nnumubar=2, double nnutaubar=0)
Definition: GAstroFlux.cxx:208
double fgWeight
(current) generated nu weight
Definition: GAstroFlux.h:177
GENIE Interface for user-defined flux classes.
Definition: GFluxI.h:29