#include <GAstroFlux.h>

Public Member Functions
	NuPropagator (double stepsz)

	~NuPropagator ()

bool	Go (double phi_start, double costheta_start, const TVector3 &detector_centre, double detector_sz, int nu_pdg, double Ev)

int	NuPdgAtDetVolBoundary (void)

TVector3 &	X3AtDetVolBoundary (void)

TVector3 &	P3AtDetVolBoundary (void)

Private Attributes
double	fStepSize

int	fNuPdg

TVector3	fX3

TVector3	fP3

Detailed Description

Definition at line 207 of file GAstroFlux.h.

Constructor & Destructor Documentation

genie::flux::GAstroFlux::NuPropagator::NuPropagator ( double stepsz )

inline

Definition at line 209 of file GAstroFlux.h.

209 : fStepSize(stepsz/units::km) { }

genie::units::km

static constexpr double km

Definition: Units.h:64

genie::flux::GAstroFlux::NuPropagator::fStepSize

double fStepSize

Definition: GAstroFlux.h:216

genie::flux::GAstroFlux::NuPropagator::~NuPropagator ( )

inline

Definition at line 210 of file GAstroFlux.h.

210 { }

Member Function Documentation

bool GAstroFlux::NuPropagator::Go	(	double	phi_start,
		double	costheta_start,
		const TVector3 &	detector_centre,
		double	detector_sz,
		int	nu_pdg,
		double	Ev
	)

Definition at line 473 of file GAstroFlux.cxx.

 {
   // initialize neutrino code
   fNuPdg = nu_pdg;
 
   //
   // initialize neutrino position vector
   //
   double sintheta  = TMath::Sqrt(1-costheta*costheta);
   double cosphi    = TMath::Cos(phi);
   double sinphi    = TMath::Sin(phi);
   double REarth    = constants::kREarth/units::km;
   double zs = REarth * costheta;
   double ys = REarth * sintheta * cosphi;
   double xs = REarth * sintheta * sinphi;
 
   TVector3 start_position(xs,ys,zs);
   fX3 = start_position - detector_centre;
 
   //
   // initialize neutrino momentum 4-vector
   //
   TVector3 direction_unit_vec = -1. * fX3.Unit();
   fP3 = Ev * direction_unit_vec;
 
   //
   // step through the Earth
   //
 
   LOG("Flux", pWARN) << "|dist|    = " << fX3.Mag();
   LOG("Flux", pWARN) << "|detsize| = " << detector_sz;
 
   while(1) {
     double currdist = fX3.Mag();
     if(currdist <= detector_sz - 0.1) break;
 
     double stepsz = (currdist-detector_sz > fStepSize) ?
                          fStepSize : currdist-detector_sz;
     if(stepsz <= 0.) break;
 
 //    LOG("Flux", pWARN) << "Stepping by... |dr| = " << stepsz;
 
     //
     // check earth density at current position, calculate interaction
     // probability over next step size, decide whether it interacts and
     // what happens if it does...
     //
     // ... todo ...
 
     fX3 += (stepsz * direction_unit_vec);
   }
 
   return true;
 }