GENIEGenerator
Loading...
Searching...
No Matches
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 \at cern.ch> University of Liverpool
Created:\n May 06, 2004
License:\n Copyright (c) 2003-2025, The GENIE Collaboration
For the full text of the license visit http://copyright.genie-mc.org

Function Documentation

◆ AreEqual() [1/2]

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

Definition at line 236 of file MathUtils.cxx.

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}
pINFO
#define pINFO
Definition Messenger.h:62
LOG
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE,...
Definition Messenger.h:96

References LOG, and pINFO.

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

◆ AreEqual() [2/2]

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

Definition at line 247 of file MathUtils.cxx.

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}

References LOG, and pINFO.

◆ CholeskyCalculateCorrelatedParamVariations()

TVectorD genie::utils::math::CholeskyCalculateCorrelatedParamVariations ( const TMatrixD & Lch,
TVectorD & g_uncorrelated )

Definition at line 186 of file MathUtils.cxx.

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}

◆ CholeskyDecomposition()

TMatrixD genie::utils::math::CholeskyDecomposition ( const TMatrixD & cov)

Definition at line 20 of file MathUtils.cxx.

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}
pERROR
#define pERROR
Definition Messenger.h:59
epsilon
const double epsilon
Definition gCalcHEDISDiffXsec.cxx:76

References epsilon, LOG, pERROR, and pINFO.

◆ CholeskyGenerateCorrelatedParams() [1/2]

TVectorD genie::utils::math::CholeskyGenerateCorrelatedParams ( const TMatrixD & Lch,
TVectorD & mean )

Definition at line 76 of file MathUtils.cxx.

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}
genie::RandomGen::Instance
static RandomGen * Instance()
Access instance.
Definition RandomGen.cxx:74
genie::RandomGen::RndNum
TRandom3 & RndNum(void) const
rnd number generator used by MC integrators & other numerical methods
Definition RandomGen.h:77

References genie::RandomGen::Instance(), LOG, and pERROR.

◆ CholeskyGenerateCorrelatedParams() [2/2]

TVectorD genie::utils::math::CholeskyGenerateCorrelatedParams ( const TMatrixD & Lch,
TVectorD & mean,
TVectorD & g_uncorrelated )

Definition at line 119 of file MathUtils.cxx.

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}

References LOG, and pERROR.

◆ CholeskyGenerateCorrelatedParamVariations()

TVectorD genie::utils::math::CholeskyGenerateCorrelatedParamVariations ( const TMatrixD & Lch)

Definition at line 165 of file MathUtils.cxx.

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}

References genie::RandomGen::Instance().

◆ IsWithinLimits() [1/3]

bool genie::utils::math::IsWithinLimits ( double x,
Range1D_t range )

Definition at line 258 of file MathUtils.cxx.

259{
260 return ( x >= range.min && x <= range.max );
261}
genie::Range1D_t::max
double max
Definition Range1.h:53
genie::Range1D_t::min
double min
Definition Range1.h:52

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

Referenced by genie::utils::kinematics::ApplyCutsToKineLimits(), genie::DFRKinematicsGenerator::ComputeMaxXSec(), genie::KPhaseSpace::IsAllowed(), and genie::utils::kinematics::PhaseSpaceVolume().

◆ IsWithinLimits() [2/3]

bool genie::utils::math::IsWithinLimits ( float x,
Range1F_t range )

Definition at line 263 of file MathUtils.cxx.

264{
265 return ( x >= range.min && x <= range.max );
266}
genie::Range1F_t::min
float min
Definition Range1.h:38
genie::Range1F_t::max
float max
Definition Range1.h:39

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

◆ IsWithinLimits() [3/3]

bool genie::utils::math::IsWithinLimits ( int i,
Range1I_t range )

Definition at line 268 of file MathUtils.cxx.

269{
270 return ( i >= range.min && i <= range.max );
271}
genie::Range1I_t::min
int min
Definition Range1.h:66
genie::Range1I_t::max
int max
Definition Range1.h:67

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

◆ KahanSummation() [1/2]

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}

◆ KahanSummation() [2/2]

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}

◆ NonNegative() [1/2]

double genie::utils::math::NonNegative ( double x)

Definition at line 273 of file MathUtils.cxx.

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

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

◆ NonNegative() [2/2]

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}
Generated on for GENIEGenerator by doxygen 1.14.0