GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Classes | Functions
genie::utils::math Namespace Reference

Simple mathematical utilities not found in ROOT's TMath. More...

Classes

class  LongLorentzVector
 

Functions

TMatrixD CholeskyDecomposition (const TMatrixD &cov)
 
TVectorD CholeskyGenerateCorrelatedParams (const TMatrixD &Lch, TVectorD &mean)
 
TVectorD CholeskyGenerateCorrelatedParams (const TMatrixD &Lch, TVectorD &mean, TVectorD &g_uncorrelated)
 
TVectorD CholeskyGenerateCorrelatedParamVariations (const TMatrixD &Lch)
 
TVectorD CholeskyCalculateCorrelatedParamVariations (const TMatrixD &Lch, TVectorD &g_uncorrelated)
 
double KahanSummation (double x[], unsigned int n)
 
double KahanSummation (const vector< double > &x)
 
bool AreEqual (double x1, double x2)
 
bool AreEqual (float x1, float x2)
 
bool IsWithinLimits (double x, Range1D_t range)
 
bool IsWithinLimits (float x, Range1F_t range)
 
bool IsWithinLimits (int i, Range1I_t range)
 
double NonNegative (double x)
 
double NonNegative (float x)
 

Detailed Description

Simple mathematical utilities not found in ROOT's TMath.

Author
Costas Andreopoulos <c.andreopoulos cern.ch> University of Liverpool
Created:
May 06, 2004
License:
Copyright (c) 2003-2024, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org

Function Documentation

bool genie::utils::math::AreEqual ( double  x1,
double  x2 
)

Definition at line 236 of file MathUtils.cxx.

References LOG, and pINFO.

Referenced by genie::PathLengthList::AreAllZero(), genie::Spline::ClosestKnotValueIsZero(), genie::AxialFormFactor::Compare(), genie::ELFormFactors::Compare(), genie::QELFormFactors::Compare(), and genie::DISStructureFunc::Compare().

237 {
238  double err = 0.001*DBL_EPSILON;
239  double dx = TMath::Abs(x1-x2);
240  if(dx<err) {
241  LOG("Math", pINFO) << x1 << " := " << x2;
242  return true;
243  }
244  return false;;
245 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
#define pINFO
Definition: Messenger.h:62
bool genie::utils::math::AreEqual ( float  x1,
float  x2 
)

Definition at line 247 of file MathUtils.cxx.

References LOG, and pINFO.

248 {
249  float err = FLT_EPSILON;
250  float dx = TMath::Abs(x1-x2);
251  if(dx<err) {
252  LOG("Math", pINFO) << x1 << " := " << x2;
253  return true;
254  }
255  return false;;
256 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
#define pINFO
Definition: Messenger.h:62
TVectorD genie::utils::math::CholeskyCalculateCorrelatedParamVariations ( const TMatrixD &  Lch,
TVectorD &  g_uncorrelated 
)

Definition at line 186 of file MathUtils.cxx.

References genie::units::g.

188 {
189  int ncols = cholesky_triangular.GetNcols();
190  int nrows = cholesky_triangular.GetNrows();
191  int npars = g_uncorrelated.GetNrows();
192 
193  assert(ncols==nrows);
194  assert(npars==nrows);
195 
196  // create a vector of unit Gaussian variables
197  // and multiply by Lt to introduce the appropriate correlations
198  TVectorD g(g_uncorrelated);
199  g *= cholesky_triangular;
200 
201  return g;
202 }
static constexpr double g
Definition: Units.h:144
TMatrixD genie::utils::math::CholeskyDecomposition ( const TMatrixD &  cov)

Definition at line 20 of file MathUtils.cxx.

References epsilon, LOG, pERROR, and pINFO.

21 {
22 // Perform a Cholesky decomposition of the input covariance matrix and
23 // return the lower triangular matrix
24 //
25  const double epsilon = 1E-12;
26 
27  int ncols = cov_matrix.GetNcols();
28  int nrows = cov_matrix.GetNrows();
29 
30  assert(ncols==nrows);
31 
32  int n = nrows;
33 
34  TMatrixD L(n, n);
35 
36  for (int i = 0; i < n; ++i) {
37 
38  // calculate the diagonal term first
39  L(i,i) = cov_matrix(i,i);
40  for (int k = 0; k < i; ++k) {
41  double tmp = L(k,i);
42  L(i,i) -= tmp*tmp;
43  }//k
44 
45  if(L(i,i) <= 0) {
46  if(fabs(L(i,i)) < epsilon){
47  L(i,i)=epsilon;
48  LOG("Cholesky", pINFO)
49  << "Changed element (" << i << ", " << i << ") to " << L(i,i);
50  }
51  else{
52  LOG("Cholesky", pERROR)
53  << "Decomposed covariance matrix not positive-definite";
54  LOG("Cholesky", pERROR)
55  << "L(" << i << "," << i << ") = " << L(i,i);
56  exit(1);
57  }
58  }
59  L(i,i) = TMath::Sqrt(L(i,i));
60  // then the off-diagonal terms
61  for (int j = i+1; j < n; ++j) {
62  L(i,j) = cov_matrix(i,j);
63  for (int k = 0; k < i; ++k) {
64  L(i,j) -= L(k,i)*L(k,j);
65  }
66  L(i,j) /= L(i,i);
67  }//j
68  }//i
69 
70  // create the transpose of L
71  TMatrixD LT(TMatrixD::kTransposed,L);
72 
73  return LT;
74 }
#define pERROR
Definition: Messenger.h:59
const double epsilon
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
#define pINFO
Definition: Messenger.h:62
TVectorD genie::utils::math::CholeskyGenerateCorrelatedParams ( const TMatrixD &  Lch,
TVectorD &  mean 
)

Definition at line 76 of file MathUtils.cxx.

References genie::units::g, LOG, and pERROR.

78 {
79 // Generate a vector of correlated params
80 
81  int ncols = cholesky_triangular.GetNcols();
82  int nrows = cholesky_triangular.GetNrows();
83  int npars = mean_params.GetNrows();
84 
85  if(ncols != nrows) {
86  LOG("Cholesky", pERROR)
87  << "Mismatch between number of columns (" << ncols
88  << ") & rows (" << nrows << ")";
89  exit(1);
90  }
91  if(npars != nrows) {
92  LOG("Cholesky", pERROR)
93  << "Mismatch between number of parameters (" << npars
94  << ") & array size (" << nrows << ")";
95  exit(1);
96  }
97 
98  int n = nrows;
99 
100  // create a vector of unit Gaussian variables
101  // and multiply by Lt to introduce the appropriate correlations
102  TVectorD g(n);
103  for (int k = 0; k < n; ++k) {
104  g(k) = RandomGen::Instance()->RndNum().Gaus();
105  }
106  g *= cholesky_triangular;
107 
108  // add the mean value offsets and store the results
109  TVectorD correlated_params(n);
110  for (int i = 0; i < n; ++i) {
111  double v = mean_params[i];
112  v += g(i);
113  correlated_params[i] = v;
114  }
115 
116  return correlated_params;
117 }
static constexpr double g
Definition: Units.h:144
#define pERROR
Definition: Messenger.h:59
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TVectorD genie::utils::math::CholeskyGenerateCorrelatedParams ( const TMatrixD &  Lch,
TVectorD &  mean,
TVectorD &  g_uncorrelated 
)

Definition at line 119 of file MathUtils.cxx.

References LOG, and pERROR.

121 {
122 // Generate a vector of correlated params
123 
124  int ncols = cholesky_triangular.GetNcols();
125  int nrows = cholesky_triangular.GetNrows();
126  int npars = mean_params.GetNrows();
127  int nunco = g_uncorrelated.GetNrows();
128 
129  if(ncols != nrows) {
130  LOG("Cholesky", pERROR)
131  << "Mismatch between number of columns (" << ncols
132  << ") & rows (" << nrows << ")";
133  exit(1);
134  }
135  if(npars != nrows) {
136  LOG("Cholesky", pERROR)
137  << "Mismatch between number of parameters (" << npars
138  << ") & array size (" << nrows << ")";
139  exit(1);
140  }
141  if(nunco != nrows) {
142  LOG("Cholesky", pERROR)
143  << "Mismatch between size of uncorrelated parameter vector (" << nunco
144  << ") & array size (" << nrows << ")";
145  exit(1);
146  }
147 
148  int n = nrows;
149 
150  // create a vector of unit Gaussian variables
151  // and multiply by Lt to introduce the appropriate correlations
152  g_uncorrelated *= cholesky_triangular;
153 
154  // add the mean value offsets and store the results
155  TVectorD correlated_params(n);
156  for (int i = 0; i < n; ++i) {
157  double v = mean_params[i];
158  v += g_uncorrelated(i);
159  correlated_params[i] = v;
160  }
161 
162  return correlated_params;
163 }
#define pERROR
Definition: Messenger.h:59
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TVectorD genie::utils::math::CholeskyGenerateCorrelatedParamVariations ( const TMatrixD &  Lch)

Definition at line 165 of file MathUtils.cxx.

References genie::units::g.

167 {
168  int ncols = cholesky_triangular.GetNcols();
169  int nrows = cholesky_triangular.GetNrows();
170 
171  assert(ncols==nrows);
172 
173  int n = nrows;
174 
175  // create a vector of unit Gaussian variables
176  // and multiply by Lt to introduce the appropriate correlations
177  TVectorD g(n);
178  for (int k = 0; k < n; ++k) {
179  g(k) = RandomGen::Instance()->RndNum().Gaus();
180  }
181  g *= cholesky_triangular;
182 
183  return g;
184 }
static constexpr double g
Definition: Units.h:144
bool genie::utils::math::IsWithinLimits ( double  x,
Range1D_t  range 
)
bool genie::utils::math::IsWithinLimits ( float  x,
Range1F_t  range 
)

Definition at line 263 of file MathUtils.cxx.

References genie::Range1F_t::max, and genie::Range1F_t::min.

264 {
265  return ( x >= range.min && x <= range.max );
266 }
bool genie::utils::math::IsWithinLimits ( int  i,
Range1I_t  range 
)

Definition at line 268 of file MathUtils.cxx.

References genie::Range1I_t::max, and genie::Range1I_t::min.

269 {
270  return ( i >= range.min && i <= range.max );
271 }
double genie::utils::math::KahanSummation ( double  x[],
unsigned int  n 
)

Definition at line 204 of file MathUtils.cxx.

205 {
206 // the Kahan summation algorithm - minimizes the error when adding a sequence
207 // of finite precision floating point numbers (compensated summation)
208 
209  double sum = x[0];
210  double c = 0.0;
211  for(unsigned int i=1; i<n; i++) {
212  double y = x[i]-c;
213  double t = sum+y;
214  c = (t-sum) - y;
215  sum = t;
216  }
217  return sum;
218 }
double genie::utils::math::KahanSummation ( const vector< double > &  x)

Definition at line 220 of file MathUtils.cxx.

221 {
222 // the Kahan summation algorithm - minimizes the error when adding a sequence
223 // of finite precision floating point numbers (compensated summation)
224 
225  double sum = x[0];
226  double c = 0.0;
227  for(unsigned int i=1; i<x.size(); i++) {
228  double y = x[i]-c;
229  double t = sum+y;
230  c = (t-sum) - y;
231  sum = t;
232  }
233  return sum;
234 }
double genie::utils::math::NonNegative ( double  x)

Definition at line 273 of file MathUtils.cxx.

Referenced by genie::PrimaryLeptonGenerator::AddToEventRecord(), genie::OutgoingDarkGenerator::AddToEventRecord(), genie::QELEventGeneratorSuSA::GenerateNucleon(), and genie::NucBindEnergyAggregator::ProcessEventRecord().

274 {
275 // this is used to handle very small numbers in sqrts
276 
277  return TMath::Max(0., x);
278 }
double genie::utils::math::NonNegative ( float  x)

Definition at line 280 of file MathUtils.cxx.

281 {
282 // this is used to handle very small numbers in sqrts
283 
284  return TMath::Max( (float)0., x);
285 }