genie::utils::prem Namespace Reference

Preliminary Earth Model. More...

## Functions

double Density (double r)

## Detailed Description

Preliminary Earth Model.

August 07, 2009

## Function Documentation

 double genie::utils::prem::Density ( double r )

Definition at line 26 of file PREM.cxx.

References genie::units::g_cm3, Munits::km, genie::constants::kREarth, and submit_syst::x.

Referenced by genie::NievesQELCCPXSec::CNCTCLimUcalc(), and main().

27 {
28 // Return the Earth density according to the PREM model
29 // Inputs: r, Distance from the centre of the Earth (in std GENIE units)
30 // Outputs: rho, Earth density (in std GENIE units)
31 //
32
33  r = TMath::Max(0., r/units::km); // convert to km
34
35  double rE = constants::kREarth/units::km;
36  double rho = 0.;
37  double x = r / rE;
38
39  if (r <= 1221.5 )
40  {
41  rho = 13.0885 - 8.8381*x*x;
42  }
43  else if (r > 1221.5 && r <= 3480.0 )
44  {
45  rho = 12.5815 - 1.2638*x - 3.6426*x*x - 5.5281*x*x*x;
46  }
47  else if (r > 3480.0 && r <= 5701.0 )
48  {
49  rho = 7.9565 - 6.4761*x + 5.5283*x*x - 3.0807*x*x*x;
50  }
51  else if (r > 5701.0 && r <= 5771.0 )
52  {
53  rho = 5.3197 - 1.4836*x;
54  }
55  else if (r > 5771.0 && r <= 5971.0 )
56  {
57  rho = 11.2494 - 8.0298*x;
58  }
59  else if (r > 5971.0 && r <= 6151.0 )
60  {
61  rho = 7.1089 - 3.8045*x;
62  }
63  else if (r > 6151.0 && r <= 6346.6 )
64  {
65  rho = 2.691 + 0.6924*x;
66  }
67  else if (r > 6346.6 && r <= 6356.0 )
68  {
69  rho = 2.90;
70  }
71  else if (r > 6356.0 && r <= 6368.0 )
72  {
73  rho = 2.60;
74  }
75  else if (r > 6368.0 && r <= rE)
76  {
77  rho = 1.02;
78  }
79
80  rho = rho * units::g_cm3;
81
82  return rho;
83 }
static const double g_cm3
Definition: Units.h:154
static constexpr Double_t km
Definition: Munits.h:148
static const double kREarth
Definition: Constants.h:111
TRandom3 r(0)