From b29c926c9249de154a307f00f1bc6da50d9abf35 Mon Sep 17 00:00:00 2001
From: Lorenzo Moneta <Lorenzo.Moneta@cern.ch>
Date: Tue, 6 May 2008 07:58:46 +0000
Subject: [PATCH] - merge changes with development branch (up to rev 23650).
Main changes are:
TMath:
from David: The following methods now have additional interface working with Iterators as well:
MinElement
MaxElement
LocMin
LocMax
Mean
GeomMean
RMS
BinarySearch
Sort -> (changed to SortItr)
- use now std in binary search algorithm. Old algorithm failed some tests and it is found to be a
little bit slower
- do not use std instead in LocMin and LocMax since generic std algorithm is slower on 64 bit arch.
- TMath::PoissonI
fix a problem for large values, by using directly the TMath::Poisson implementation
- TMath::Landau and TMath::Landau
remove duplications by using directly ROOT::Math::landau_pdf and ROOT::Math::landau_cdf
- TMath::Prob
use ROOT::Math::chisquared_cdf_c (gives also less rounding errors for small values)
- TMath::GammaDist
use ROOT::Math::gamma_pdf , implementation now returns correct value, 0, for large values of
shape parameter instead of nan
- TMath::LogNormal:
use ROOT::Math::lognormal_pdf
- compile TMath dictionary with -O2
- fix a warning in -Wshadow
- TMath.h is now reordered by functionality and alphabetical order
test:
various improvement from David, see http://root.cern.ch/viewcvs?rev=22880&root=root&view=rev
- add dictionary for the enumeration used in the Integration class
- from David: some clean up of Minimizer class (remove typedef definitions)
git-svn-id: http://root.cern.ch/svn/root/trunk@23655 27541ba8-7e3a-0410-8455-c3a389f83636
---
math/mathcore/Module.mk | 6 +-
math/mathcore/inc/LinkDef.h | 7 -
math/mathcore/inc/Math/LinkDef.h | 6 +
math/mathcore/inc/Math/Minimizer.h | 14 +-
math/mathcore/inc/TMath.h | 514 +++++++++++-------
math/mathcore/src/ProbFuncMathCore.cxx | 3 +-
math/mathcore/src/TMath.cxx | 170 +-----
math/mathcore/test/Makefile | 33 +-
math/mathcore/test/binarySearchTime.cxx | 1 -
math/mathcore/test/stdsort.cxx | 85 ++-
math/mathcore/test/stressTMath.cxx | 182 +++++--
math/mathcore/test/testIntegration.cxx | 40 +-
.../mathcore/test/testIntegrationMultiDim.cxx | 2 +-
math/mathcore/test/testRootFinder.cxx | 19 +-
math/mathcore/test/testSortOrder.cxx | 71 +--
math/mathcore/test/testSpecFuncBeta.cxx | 92 +++-
math/mathcore/test/testSpecFuncBetaI.cxx | 109 ++--
math/mathcore/test/testSpecFuncErf.cxx | 166 ++++--
math/mathcore/test/testSpecFuncGamma.cxx | 142 +++--
math/mathcore/test/testTMath.cxx | 260 ++++-----
20 files changed, 1110 insertions(+), 812 deletions(-)
diff --git a/math/mathcore/Module.mk b/math/mathcore/Module.mk
index 253d3d6c700..e210e5af19f 100644
--- a/math/mathcore/Module.mk
+++ b/math/mathcore/Module.mk
@@ -79,6 +79,7 @@ ALLMAPS += $(MATHCOREMAP)
# include all dependency files
INCLUDEFILES += $(MATHCOREDEP)
+
##### local rules #####
.PHONY: all-$(MODNAME) clean-$(MODNAME) distclean-$(MODNAME) \
test-$(MODNAME)
@@ -142,4 +143,7 @@ test-$(MODNAME): all-$(MODNAME)
##### extra rules ######
$(MATHCOREO): CXXFLAGS += -DUSE_ROOT_ERROR
-$(MATHCOREDO): CXXFLAGS += -DUSE_ROOT_ERROR
+$(MATHCOREDO): CXXFLAGS += -DUSE_ROOT_ERROR
+# add optimization to G__Math compilation
+$(MATHCOREDO1) : NOOPT = $(OPT)
+
diff --git a/math/mathcore/inc/LinkDef.h b/math/mathcore/inc/LinkDef.h
index 19fecb2e85a..1c808810107 100644
--- a/math/mathcore/inc/LinkDef.h
+++ b/math/mathcore/inc/LinkDef.h
@@ -114,13 +114,6 @@
#pragma link C++ function TMath::IsInside(Float_t, Float_t, Int_t, Float_t*, Float_t*);
#pragma link C++ function TMath::IsInside(Int_t, Int_t, Int_t, Int_t*, Int_t*);
-#pragma link C++ function TMath::Sort(Int_t, const Short_t*, Int_t*, Bool_t);
-#pragma link C++ function TMath::Sort(Int_t, const Int_t*, Int_t*, Bool_t);
-#pragma link C++ function TMath::Sort(Int_t, const Float_t*, Int_t*, Bool_t);
-#pragma link C++ function TMath::Sort(Int_t, const Double_t*, Int_t*, Bool_t);
-#pragma link C++ function TMath::Sort(Int_t, const Long_t*, Int_t*, Bool_t);
-#pragma link C++ function TMath::Sort(Int_t, const Long64_t*, Int_t*, Bool_t);
-
#pragma link C++ function TMath::Sort(Long64_t, const Short_t*, Long64_t*, Bool_t);
#pragma link C++ function TMath::Sort(Long64_t, const Int_t*, Long64_t*, Bool_t);
#pragma link C++ function TMath::Sort(Long64_t, const Float_t*, Long64_t*, Bool_t);
diff --git a/math/mathcore/inc/Math/LinkDef.h b/math/mathcore/inc/Math/LinkDef.h
index 0800fd8b669..3265e506470 100644
--- a/math/mathcore/inc/Math/LinkDef.h
+++ b/math/mathcore/inc/Math/LinkDef.h
@@ -60,6 +60,12 @@
#pragma link C++ class ROOT::Math::AdaptiveIntegratorMultiDim+;
#pragma link C++ typedef ROOT::Math::Integrator;
+#pragma link C++ namespace ROOT::Math::IntegrationOneDim;
+#pragma link C++ enum ROOT::Math::IntegrationOneDim::Type;
+#pragma link C++ namespace ROOT::Math::IntegrationMultiDim;
+#pragma link C++ enum ROOT::Math::IntegrationMultiDim::Type;
+
+
#pragma link C++ class ROOT::Math::BasicFitMethodFunction<ROOT::Math::IBaseFunctionMultiDim>+;
#pragma link C++ class ROOT::Math::BasicFitMethodFunction<ROOT::Math::IGradientFunctionMultiDim>+;
diff --git a/math/mathcore/inc/Math/Minimizer.h b/math/mathcore/inc/Math/Minimizer.h
index e46bbfea51c..dd59611a733 100644
--- a/math/mathcore/inc/Math/Minimizer.h
+++ b/math/mathcore/inc/Math/Minimizer.h
@@ -13,8 +13,8 @@
#ifndef ROOT_Math_Minimizer
#define ROOT_Math_Minimizer
-#ifndef ROOT_Math_IFunctionfwd
-#include "Math/IFunctionfwd.h"
+#ifndef ROOT_Math_IFunction
+#include "Math/IFunction.h"
#endif
@@ -73,9 +73,6 @@ class Minimizer {
public:
- typedef ROOT::Math::IMultiGenFunction IObjFunction;
- typedef ROOT::Math::IMultiGradFunction IGradObjFunction;
-
/**
Default constructor
*/
@@ -123,10 +120,13 @@ public:
virtual void Clear() {}
/// set the function to minimize
- virtual void SetFunction(const IObjFunction & func) = 0;
+ virtual void SetFunction(const ROOT::Math::IMultiGenFunction & func) = 0;
/// set a function to minimize using gradient
- virtual void SetFunction(const IGradObjFunction & func) = 0;
+ virtual void SetFunction(const ROOT::Math::IMultiGradFunction & func)
+ {
+ SetFunction(static_cast<const ROOT::Math::IMultiGenFunction &> (func));
+ }
/// add variables . Return number of variables succesfully added
diff --git a/math/mathcore/inc/TMath.h b/math/mathcore/inc/TMath.h
index 64d1b71536d..bc4575a15fe 100644
--- a/math/mathcore/inc/TMath.h
+++ b/math/mathcore/inc/TMath.h
@@ -35,7 +35,10 @@
namespace TMath {
- // Fundamental constants
+ /* ************************* */
+ /* * Fundamental constants * */
+ /* ************************* */
+
inline Double_t Pi() { return 3.14159265358979323846; }
inline Double_t TwoPi() { return 2.0 * Pi(); }
inline Double_t PiOver2() { return Pi() / 2.0; }
@@ -120,7 +123,13 @@ namespace TMath {
inline Double_t Qe() { return 1.602176462e-19; } // C
inline Double_t QeUncertainty() { return 0.000000063e-19; }
- // Trigo
+ /* ************************** */
+ /* * Mathematical Functions * */
+ /* ************************** */
+
+ /* ***************************** */
+ /* * Trigonometrical Functions * */
+ /* ***************************** */
inline Double_t Sin(Double_t);
inline Double_t Cos(Double_t);
inline Double_t Tan(Double_t);
@@ -136,7 +145,10 @@ namespace TMath {
Double_t ATanH(Double_t);
Double_t Hypot(Double_t x, Double_t y);
- // Misc
+
+ /* ************************ */
+ /* * Elementary Functions * */
+ /* ************************ */
inline Double_t Sqrt(Double_t x);
inline Double_t Ceil(Double_t x);
inline Int_t CeilNint(Double_t x);
@@ -157,113 +169,156 @@ namespace TMath {
// Some integer math
Long_t Hypot(Long_t x, Long_t y); // sqrt(px*px + py*py)
+ /* ******************** */
+ /* * Array Algorithms * */
+ /* ******************** */
+
// Min, Max of an array
template <typename T> T MinElement(Long64_t n, const T *a);
template <typename T> T MaxElement(Long64_t n, const T *a);
// Locate Min, Max element number in an array
template <typename T> Long64_t LocMin(Long64_t n, const T *a);
+ template <typename Iterator> Iterator LocMin(Iterator first, Iterator last);
template <typename T> Long64_t LocMax(Long64_t n, const T *a);
-
- //Mean, Geometric Mean, Median, RMS
-
- template <typename T> Double_t Mean(Long64_t n, const T *a, const Double_t *w=0);
- template <typename T> Double_t GeomMean(Long64_t n, const T *a);
- template <typename T> Double_t RMS(Long64_t n, const T *a);
-
- template <class Element, class Index, class Size> Double_t MedianImp(Size n, const Element *a, const Double_t *w=0, Index *work=0);
- template <typename T> Double_t Median(Long64_t n, const T *a, const Double_t *w=0, Long64_t *work=0);
-
- //k-th order statistic
- template <class Element, typename Size> Element KOrdStat(Size n, const Element *a, Size k, Size *work = 0);
-
- //Sample quantiles
- void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob,
- Bool_t isSorted=kTRUE, Int_t *index = 0, Int_t type=7);
+ template <typename Iterator> Iterator LocMax(Iterator first, Iterator last);
// Binary search
template <typename T> Long64_t BinarySearch(Long64_t n, const T *array, T value);
template <typename T> Long64_t BinarySearch(Long64_t n, const T **array, T value);
+ template <typename Iterator, typename Element> Iterator BinarySearch(Iterator first, Iterator last, Element value);
// Hashing
ULong_t Hash(const void *txt, Int_t ntxt);
ULong_t Hash(const char *str);
- // IsInside
- template <typename T> Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y);
-
// Sorting
- template <class Element, class Index, class Size> void SortImp(Size n, const Element*, Index* index, Bool_t down=kTRUE);
- template <typename Element, typename Index, typename Size> void Sort(Size n, const Element* a, Index* index, Bool_t down=kTRUE);
+ template <typename Element, typename Index>
+ void Sort(Long64_t n, const Element* a, Index* index, Bool_t down=kTRUE);
+ template <typename Iterator, typename IndexIterator>
+ void SortItr(Iterator first, Iterator last, IndexIterator index, Bool_t down=kTRUE);
+
void BubbleHigh(Int_t Narr, Double_t *arr1, Int_t *arr2);
void BubbleLow (Int_t Narr, Double_t *arr1, Int_t *arr2);
- // Advanced
- template <typename T> T *Cross(const T v1[3],const T v2[3], T out[3]); // Calculate the Cross Product of two vectors
- Float_t Normalize(Float_t v[3]); // Normalize a vector
- Double_t Normalize(Double_t v[3]); // Normalize a vector
- template <typename T> inline T NormCross(const T v1[3],const T v2[3],T out[3]); //Calculate the Normalized Cross Product of two vectors
- template <typename T> T *Normal2Plane(const T v1[3],const T v2[3],const T v3[3], T normal[3]); // Calculate a normal vector of a plane
+ Bool_t Permute(Int_t n, Int_t *a); // Find permutations
+
+ /* ************************* */
+ /* * Geometrical Functions * */
+ /* ************************* */
+
+ //Sample quantiles
+ void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob,
+ Bool_t isSorted=kTRUE, Int_t *index = 0, Int_t type=7);
+
+ // IsInside
+ template <typename T> Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y);
+
+ // Calculate the Cross Product of two vectors
+ template <typename T> T *Cross(const T v1[3],const T v2[3], T out[3]);
+
+ Float_t Normalize(Float_t v[3]); // Normalize a vector
+ Double_t Normalize(Double_t v[3]); // Normalize a vector
+
+ //Calculate the Normalized Cross Product of two vectors
+ template <typename T> inline T NormCross(const T v1[3],const T v2[3],T out[3]);
+
+ // Calculate a normal vector of a plane
+ template <typename T> T *Normal2Plane(const T v1[3],const T v2[3],const T v3[3], T normal[3]);
+
+ /* ************************ */
+ /* * Polinomial Functions * */
+ /* ************************ */
+
Bool_t RootsCubic(const Double_t coef[4],Double_t &a, Double_t &b, Double_t &c);
- Double_t BreitWigner(Double_t x, Double_t mean=0, Double_t gamma=1);
- Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE);
- Double_t Landau(Double_t x, Double_t mpv=0, Double_t sigma=1, Bool_t norm=kFALSE);
- Double_t Voigt(Double_t x, Double_t sigma, Double_t lg, Int_t R = 4);
+
+ /* *********************** */
+ /* * Statistic Functions * */
+ /* *********************** */
+
+ Double_t Binomial(Int_t n,Int_t k); // Calculate the binomial coefficient n over k
+ Double_t BinomialI(Double_t p, Int_t n, Int_t k);
+ Double_t BreitWigner(Double_t x, Double_t mean=0, Double_t gamma=1);
+ Double_t CauchyDist(Double_t x, Double_t t=0, Double_t s=1);
+ Double_t ChisquareQuantile(Double_t p, Double_t ndf);
+ Double_t FDist(Double_t F, Double_t N, Double_t M);
+ Double_t FDistI(Double_t F, Double_t N, Double_t M);
+ Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE);
+ Double_t KolmogorovProb(Double_t z);
+ Double_t KolmogorovTest(Int_t na, const Double_t *a, Int_t nb, const Double_t *b, Option_t *option);
+ Double_t Landau(Double_t x, Double_t mpv=0, Double_t sigma=1, Bool_t norm=kFALSE);
+ Double_t LandauI(Double_t x);
+ Double_t LaplaceDist(Double_t x, Double_t alpha=0, Double_t beta=1);
+ Double_t LaplaceDistI(Double_t x, Double_t alpha=0, Double_t beta=1);
+ Double_t LogNormal(Double_t x, Double_t sigma, Double_t theta=0, Double_t m=1);
+ Double_t NormQuantile(Double_t p);
+ Double_t Poisson(Double_t x, Double_t par);
+ Double_t PoissonI(Double_t x, Double_t par);
+ Double_t Prob(Double_t chi2,Int_t ndf);
+ Double_t Student(Double_t T, Double_t ndf);
+ Double_t StudentI(Double_t T, Double_t ndf);
+ Double_t StudentQuantile(Double_t p, Double_t ndf, Bool_t lower_tail=kTRUE);
+ Double_t Vavilov(Double_t x, Double_t kappa, Double_t beta2);
+ Double_t VavilovI(Double_t x, Double_t kappa, Double_t beta2);
+ Double_t Voigt(Double_t x, Double_t sigma, Double_t lg, Int_t R = 4);
+
+ /* ************************** */
+ /* * Statistics over arrays * */
+ /* ************************** */
+
+ //Mean, Geometric Mean, Median, RMS
+
+ template <typename T> Double_t Mean(Long64_t n, const T *a, const Double_t *w=0);
+ template <typename Iterator> Double_t Mean(Iterator first, Iterator last);
+ template <typename Iterator, typename WeightIterator> Double_t Mean(Iterator first, Iterator last, WeightIterator w);
+
+ template <typename T> Double_t GeomMean(Long64_t n, const T *a);
+ template <typename Iterator> Double_t GeomMean(Iterator first, Iterator last);
+
+ template <typename T> Double_t RMS(Long64_t n, const T *a);
+ template <typename Iterator> Double_t RMS(Iterator first, Iterator last);
+
+ template <typename T> Double_t Median(Long64_t n, const T *a, const Double_t *w=0, Long64_t *work=0);
+
+ //k-th order statistic
+ template <class Element, typename Size> Element KOrdStat(Size n, const Element *a, Size k, Size *work = 0);
+
+ /* ******************* */
+ /* * Special Functions */
+ /* ******************* */
+
+ Double_t Beta(Double_t p, Double_t q);
+ Double_t BetaCf(Double_t x, Double_t a, Double_t b);
+ Double_t BetaDist(Double_t x, Double_t p, Double_t q);
+ Double_t BetaDistI(Double_t x, Double_t p, Double_t q);
+ Double_t BetaIncomplete(Double_t x, Double_t a, Double_t b);
// Bessel functions
- Double_t BesselI(Int_t n,Double_t x); // integer order modified Bessel function I_n(x)
- Double_t BesselK(Int_t n,Double_t x); // integer order modified Bessel function K_n(x)
- Double_t BesselI0(Double_t x); // modified Bessel function I_0(x)
- Double_t BesselK0(Double_t x); // modified Bessel function K_0(x)
- Double_t BesselI1(Double_t x); // modified Bessel function I_1(x)
- Double_t BesselK1(Double_t x); // modified Bessel function K_1(x)
- Double_t BesselJ0(Double_t x); // Bessel function J0(x) for any real x
- Double_t BesselJ1(Double_t x); // Bessel function J1(x) for any real x
- Double_t BesselY0(Double_t x); // Bessel function Y0(x) for positive x
- Double_t BesselY1(Double_t x); // Bessel function Y1(x) for positive x
- Double_t StruveH0(Double_t x); // Struve functions of order 0
- Double_t StruveH1(Double_t x); // Struve functions of order 1
- Double_t StruveL0(Double_t x); // Modified Struve functions of order 0
- Double_t StruveL1(Double_t x); // Modified Struve functions of order 1
-
- // Statistics
- Double_t Beta(Double_t p, Double_t q);
- Double_t BetaCf(Double_t x, Double_t a, Double_t b);
- Double_t BetaDist(Double_t x, Double_t p, Double_t q);
- Double_t BetaDistI(Double_t x, Double_t p, Double_t q);
- Double_t BetaIncomplete(Double_t x, Double_t a, Double_t b);
- Double_t Binomial(Int_t n,Int_t k); // Calculate the binomial coefficient n over k
- Double_t BinomialI(Double_t p, Int_t n, Int_t k);
- Double_t CauchyDist(Double_t x, Double_t t=0, Double_t s=1);
- Double_t ChisquareQuantile(Double_t p, Double_t ndf);
- Double_t DiLog(Double_t x);
- Double_t Erf(Double_t x);
- Double_t ErfInverse(Double_t x);
- Double_t Erfc(Double_t x);
+ Double_t BesselI(Int_t n,Double_t x); // integer order modified Bessel function I_n(x)
+ Double_t BesselK(Int_t n,Double_t x); // integer order modified Bessel function K_n(x)
+ Double_t BesselI0(Double_t x); // modified Bessel function I_0(x)
+ Double_t BesselK0(Double_t x); // modified Bessel function K_0(x)
+ Double_t BesselI1(Double_t x); // modified Bessel function I_1(x)
+ Double_t BesselK1(Double_t x); // modified Bessel function K_1(x)
+ Double_t BesselJ0(Double_t x); // Bessel function J0(x) for any real x
+ Double_t BesselJ1(Double_t x); // Bessel function J1(x) for any real x
+ Double_t BesselY0(Double_t x); // Bessel function Y0(x) for positive x
+ Double_t BesselY1(Double_t x); // Bessel function Y1(x) for positive x
+ Double_t StruveH0(Double_t x); // Struve functions of order 0
+ Double_t StruveH1(Double_t x); // Struve functions of order 1
+ Double_t StruveL0(Double_t x); // Modified Struve functions of order 0
+ Double_t StruveL1(Double_t x); // Modified Struve functions of order 1
+
+ Double_t DiLog(Double_t x);
+ Double_t Erf(Double_t x);
+ Double_t ErfInverse(Double_t x);
+ Double_t Erfc(Double_t x);
inline Double_t ErfcInverse(Double_t x) {return TMath::ErfInverse(1-x);}
- Double_t FDist(Double_t F, Double_t N, Double_t M);
- Double_t FDistI(Double_t F, Double_t N, Double_t M);
- Double_t Freq(Double_t x);
- Double_t Gamma(Double_t z);
- Double_t Gamma(Double_t a,Double_t x);
- Double_t GammaDist(Double_t x, Double_t gamma, Double_t mu=0, Double_t beta=1);
- Double_t KolmogorovProb(Double_t z);
- Double_t KolmogorovTest(Int_t na, const Double_t *a, Int_t nb, const Double_t *b, Option_t *option);
- Double_t LandauI(Double_t x);
- Double_t LaplaceDist(Double_t x, Double_t alpha=0, Double_t beta=1);
- Double_t LaplaceDistI(Double_t x, Double_t alpha=0, Double_t beta=1);
- Double_t LnGamma(Double_t z);
- Double_t LogNormal(Double_t x, Double_t sigma, Double_t theta=0, Double_t m=1);
- Double_t NormQuantile(Double_t p);
- Bool_t Permute(Int_t n, Int_t *a); // Find permutations
- Double_t Poisson(Double_t x, Double_t par);
- Double_t PoissonI(Double_t x, Double_t par);
- Double_t Prob(Double_t chi2,Int_t ndf);
- Double_t Student(Double_t T, Double_t ndf);
- Double_t StudentI(Double_t T, Double_t ndf);
- Double_t StudentQuantile(Double_t p, Double_t ndf, Bool_t lower_tail=kTRUE);
- Double_t Vavilov(Double_t x, Double_t kappa, Double_t beta2);
- Double_t VavilovI(Double_t x, Double_t kappa, Double_t beta2);
+ Double_t Freq(Double_t x);
+ Double_t Gamma(Double_t z);
+ Double_t Gamma(Double_t a,Double_t x);
+ Double_t GammaDist(Double_t x, Double_t gamma, Double_t mu=0, Double_t beta=1);
+ Double_t LnGamma(Double_t z);
}
@@ -436,7 +491,28 @@ Long64_t TMath::LocMin(Long64_t n, const T *a) {
// Return index of array with the minimum element.
// If more than one element is minimum returns first found.
- return (Long64_t) (std::min_element(a, a+n) - a);
+ // Implement here since this one is found to be faster (mainly on 64 bit machines)
+ // than stl generic implementation.
+ // When performing the comparison, the STL implementation needs to de-reference both the array iterator
+ // and the iterator pointing to the resulting minimum location
+
+ if (n <= 0 || !a) return -1;
+ T xmin = a[0];
+ Long64_t loc = 0;
+ for (Long64_t i = 1; i < n; i++) {
+ if (xmin > a[i]) {
+ xmin = a[i];
+ loc = i;
+ }
+ }
+ return loc;
+}
+
+template <typename Iterator>
+Iterator TMath::LocMin(Iterator first, Iterator last) {
+ // Return index of array with the minimum element.
+ // If more than one element is minimum returns first found.
+ return std::min_element(first, last);
}
template <typename T>
@@ -444,91 +520,146 @@ Long64_t TMath::LocMax(Long64_t n, const T *a) {
// Return index of array with the maximum element.
// If more than one element is maximum returns first found.
- return (Long64_t) (std::max_element(a, a+n) - a);
+ // Implement here since it is faster (see comment in LocMin function)
+
+ if (n <= 0 || !a) return -1;
+ T xmax = a[0];
+ Long64_t loc = 0;
+ for (Long64_t i = 1; i < n; i++) {
+ if (xmax < a[i]) {
+ xmax = a[i];
+ loc = i;
+ }
+ }
+ return loc;
+}
+
+template <typename Iterator>
+Iterator TMath::LocMax(Iterator first, Iterator last)
+{
+ // Return index of array with the maximum element.
+ // If more than one element is maximum returns first found.
+
+ return std::max_element(first, last);
}
template<typename T>
struct CompareDesc {
- CompareDesc(const T * d) : fData(d) {}
+ CompareDesc(T d) : fData(d) {}
bool operator()(int i1, int i2) {
- return fData[i1] > fData[i2];
+ return *(fData + i1) > *(fData + i2);
}
- const T * fData;
+ T fData;
};
template<typename T>
struct CompareAsc {
- CompareAsc(const T * d) : fData(d) {}
+ CompareAsc(T d) : fData(d) {}
bool operator()(int i1, int i2) {
- return fData[i1] < fData[i2];
+ return *(fData + i1) < *(fData + i2);
}
- const T * fData;
+ T fData;
};
-template <typename T> Double_t TMath::Mean(Long64_t n, const T *a, const Double_t *w)
+template <typename Iterator>
+Double_t TMath::Mean(Iterator first, Iterator last)
{
// Return the weighted mean of an array a with length n.
- if (n <= 0 || !a) return 0;
+ Double_t sum = 0;
+ Double_t sumw = 0;
+ while ( first != last )
+ {
+ sum += *first;
+ sumw += 1;
+ first++;
+ }
+
+ return sum/sumw;
+}
+
+template <typename Iterator, typename WeightIterator>
+Double_t TMath::Mean(Iterator first, Iterator last, WeightIterator w)
+{
+ // Return the weighted mean of an array a with length n.
Double_t sum = 0;
Double_t sumw = 0;
- if (w) {
- for (Long64_t i = 0; i < n; i++) {
- if (w[i] < 0) {
- ::Error("TMath::Mean","w[%d] = %.4e < 0 ?!",i,w[i]);
- return 0;
- }
- sum += w[i]*a[i];
- sumw += w[i];
- }
- if (sumw <= 0) {
- ::Error("TMath::Mean","sum of weights == 0 ?!");
+ int i = 0;
+ while ( first != last ) {
+ if ( *w < 0) {
+ ::Error("TMath::Mean","w[%d] = %.4e < 0 ?!",i,*w);
return 0;
}
- } else {
- sumw = n;
- for (Long64_t i = 0; i < n; i++)
- sum += a[i];
+ sum += (*w) * (*first);
+ sumw += (*w) ;
+ ++w;
+ ++first;
+ ++i;
+ }
+ if (sumw <= 0) {
+ ::Error("TMath::Mean","sum of weights == 0 ?!");
+ return 0;
}
- return sum/sumw;
+ return sum/sumw;
}
+template <typename T>
+Double_t TMath::Mean(Long64_t n, const T *a, const Double_t *w)
+{
+ // Return the weighted mean of an array a with length n.
+
+ if (w) {
+ return TMath::Mean(a, a+n, w);
+ } else {
+ return TMath::Mean(a, a+n);
+ }
+}
-template <typename T> Double_t TMath::GeomMean(Long64_t n, const T *a)
+template <typename Iterator>
+Double_t TMath::GeomMean(Iterator first, Iterator last)
{
// Return the geometric mean of an array a with length n.
// geometric_mean = (Prod_i=0,n-1 |a[i]|)^1/n
- if (n <= 0 || !a) return 0;
-
Double_t logsum = 0.;
- for (Long64_t i = 0; i < n; i++) {
- if (a[i] == 0) return 0.;
- Double_t absa = (Double_t) TMath::Abs(a[i]);
+ Long64_t n = 0;
+ while ( first != last ) {
+ if (*first == 0) return 0.;
+ Double_t absa = (Double_t) TMath::Abs(*first);
logsum += TMath::Log(absa);
+ ++first;
+ ++n;
}
return TMath::Exp(logsum/n);
}
-template <typename T> Double_t TMath::RMS(Long64_t n, const T *a)
+template <typename T>
+Double_t TMath::GeomMean(Long64_t n, const T *a)
{
- // Return the RMS of an array a with length n.
+ return TMath::GeomMean(a, a+n);
+}
- if (n <= 0 || !a) return 0;
+template <typename Iterator>
+Double_t TMath::RMS(Iterator first, Iterator last)
+{
+ // Return the RMS of an array a with length n.
+ Double_t n = 0;
Double_t tot = 0, tot2 =0, adouble;
- for (Long64_t i=0;i<n;i++) {
- adouble=Double_t(a[i]);
+ while ( first != last ) {
+ adouble=Double_t(*first);
tot += adouble; tot2 += adouble*adouble;
+ ++first;
+ ++n;
}
Double_t n1 = 1./n;
Double_t mean = tot*n1;
@@ -536,28 +667,28 @@ template <typename T> Double_t TMath::RMS(Long64_t n, const T *a)
return rms;
}
-template <typename T> Double_t TMath::Median(Long64_t n, const T *a, const Double_t *w, Long64_t *work)
+template <typename T>
+Double_t TMath::RMS(Long64_t n, const T *a)
{
- // Return the median of the array a where each entry i has weight w[i] .
- // Both arrays have a length of at least n . The median is a number obtained
- // from the sorted array a through
- //
- // median = (a[jl]+a[jh])/2. where (using also the sorted index on the array w)
- //
- // sum_i=0,jl w[i] <= sumTot/2
- // sum_i=0,jh w[i] >= sumTot/2
- // sumTot = sum_i=0,n w[i]
- //
- // If w=0, the algorithm defaults to the median definition where it is
- // a number that divides the sorted sequence into 2 halves.
- // When n is odd or n > 1000, the median is kth element k = (n + 1) / 2.
- // when n is even and n < 1000the median is a mean of the elements k = n/2 and k = n/2 + 1.
- //
- // If work is supplied, it is used to store the sorting index and assumed to be
- // >= n . If work=0, local storage is used, either on the stack if n < kWorkMax
- // or on the heap for n >= kWorkMax .
+ return TMath::RMS(a, a+n);
+}
- return MedianImp(n, a, w, work);
+template <typename Iterator, typename Element>
+Iterator TMath::BinarySearch(Iterator first, Iterator last, Element value)
+{
+ // Binary search in an array of n values to locate value.
+ //
+ // The values in the iterators range are supposed to be sorted
+ // prior to this call. If match is found, function returns
+ // position of element. If no match found, function gives nearest
+ // element smaller than value.
+
+ Iterator pind;
+ pind = std::lower_bound(first, last, value);
+ if ( (pind != last) && (*pind == value) )
+ return pind;
+ else
+ return ( pind - 1);
}
@@ -569,28 +700,12 @@ template <typename T> Long64_t TMath::BinarySearch(Long64_t n, const T *array,
// If match is found, function returns position of element.
// If no match found, function gives nearest element smaller than value.
-
-#ifdef USE_NEW_STD_IMPL
const T* pind;
pind = std::lower_bound(array, array + n, value);
if ( (pind != array + n) && (*pind == value) )
return (pind - array);
else
return ( pind - array - 1);
-#else
-
- Long64_t nabove, nbelow, middle;
- nabove = n+1;
- nbelow = 0;
- while(nabove-nbelow > 1) {
- middle = (nabove+nbelow)/2;
- if (value == array[middle-1]) return middle-1;
- if (value < array[middle-1]) nabove = middle;
- else nbelow = middle;
- }
- return nbelow-1;
-
-#endif
}
template <typename T> Long64_t TMath::BinarySearch(Long64_t n, const T **array, T value)
@@ -601,39 +716,53 @@ template <typename T> Long64_t TMath::BinarySearch(Long64_t n, const T **array,
// If match is found, function returns position of element.
// If no match found, function gives nearest element smaller than value.
-#ifdef USE_NEW_STD_IMPL
const T* pind;
pind = std::lower_bound(*array, *array + n, value);
if ( (pind != *array + n) && (*pind == value) )
return (pind - *array);
else
return ( pind - *array - 1);
-#else
- Long64_t nabove, nbelow, middle;
- nabove = n+1;
- nbelow = 0;
- while(nabove-nbelow > 1) {
- middle = (nabove+nbelow)/2;
- if (value == *array[middle-1]) return middle-1;
- if (value < *array[middle-1]) nabove = middle;
- else nbelow = middle;
- }
- return nbelow-1;
-#endif
+}
+
+template <typename Iterator, typename IndexIterator>
+void TMath::SortItr(Iterator first, Iterator last, IndexIterator index, Bool_t down)
+{
+ // Sort the n1 elements of the Short_t array a.
+ // In output the array index contains the indices of the sorted array.
+ // If down is false sort in increasing order (default is decreasing order).
+
+ // NOTE that the array index must be created with a length >= n1
+ // before calling this function.
+ int i = 0;
+
+ IndexIterator cindex = index;
+ for ( Iterator cfirst = first; cfirst != last; ++cfirst )
+ {
+ *cindex = i++;
+ ++cindex;
+ }
+
+ if ( down )
+ std::sort(index, cindex, CompareDesc<Iterator>(first) );
+ else
+ std::sort(index, cindex, CompareAsc<Iterator>(first) );
}
-template <typename Element, typename Index, typename Size> void TMath::Sort(Size n, const Element* a, Index* index, Bool_t down)
+template <typename Element, typename Index> void TMath::Sort(Long64_t n, const Element* a, Index* index, Bool_t down)
{
// Sort the n1 elements of the Short_t array a.
// In output the array index contains the indices of the sorted array.
// If down is false sort in increasing order (default is decreasing order).
- // This is a translation of the CERNLIB routine sortzv (M101)
- // based on the quicksort algorithm.
+
// NOTE that the array index must be created with a length >= n1
// before calling this function.
- SortImp(n,a,index,down);
+ for(Long64_t i = 0; i < n; i++) { index[i] = i; }
+ if ( down )
+ std::sort(index, index + n, CompareDesc<const Element*>(a) );
+ else
+ std::sort(index, index + n, CompareAsc<const Element*>(a) );
}
template <typename T> T *TMath::Cross(const T v1[3],const T v2[3], T out[3])
@@ -693,9 +822,7 @@ template <typename T> Bool_t TMath::IsInside(T xp, T yp, Int_t np, T *x, T *y)
return kFALSE;
}
-
-template <class Element, class Index, class Size>
-Double_t TMath::MedianImp(Size n, const Element *a,const Double_t *w, Index *work)
+template <typename T> Double_t TMath::Median(Long64_t n, const T *a, const Double_t *w, Long64_t *work)
{
// Return the median of the array a where each entry i has weight w[i] .
// Both arrays have a length of at least n . The median is a number obtained
@@ -721,8 +848,8 @@ Double_t TMath::MedianImp(Size n, const Element *a,const Double_t *w, Index *wor
if (n <= 0 || !a) return 0;
Bool_t isAllocated = kFALSE;
Double_t median;
- Index *ind;
- Index workLocal[kWorkMax];
+ Long64_t *ind;
+ Long64_t workLocal[kWorkMax];
if (work) {
ind = work;
@@ -730,7 +857,7 @@ Double_t TMath::MedianImp(Size n, const Element *a,const Double_t *w, Index *wor
ind = workLocal;
if (n > kWorkMax) {
isAllocated = kTRUE;
- ind = new Index[n];
+ ind = new Long64_t[n];
}
}
@@ -746,7 +873,7 @@ Double_t TMath::MedianImp(Size n, const Element *a,const Double_t *w, Index *wor
sumTot2 /= 2.;
- SortImp(n, a, ind, kFALSE);
+ Sort(n, a, ind, kFALSE);
Double_t sum = 0.;
Int_t jl;
@@ -778,27 +905,8 @@ Double_t TMath::MedianImp(Size n, const Element *a,const Double_t *w, Index *wor
return median;
}
-template <class Element, class Index, class Size>
-void TMath::SortImp(Size n1, const Element *a,
- Index *index, Bool_t down)
-{
- // Templated version of the Sort.
- //
- // Sort the n1 elements of the array a.of Element
- // In output the array index contains the indices of the sorted array.
- // If down is false sort in increasing order (default is decreasing order).
- //
- // NOTE that the array index must be created with a length >= n1
- // before calling this function.
- //
- // See also the declarations at the top of this file.
- for(Size i = 0; i < n1; i++) { index[i] = i; }
- if ( down )
- std::sort(index, index + n1, CompareDesc<Element>(a) );
- else
- std::sort(index, index + n1, CompareAsc<Element>(a) );
-}
+
template <class Element, typename Size>
Element TMath::KOrdStat(Size n, const Element *a, Size k, Size *work)
diff --git a/math/mathcore/src/ProbFuncMathCore.cxx b/math/mathcore/src/ProbFuncMathCore.cxx
index c126fbc6c77..b6bfbfeba78 100644
--- a/math/mathcore/src/ProbFuncMathCore.cxx
+++ b/math/mathcore/src/ProbFuncMathCore.cxx
@@ -328,8 +328,7 @@ namespace Math {
lan = 1-(a2[1]+(a2[2]+a2[3]*u)*u)*u;
}
return lan;
-
- return 0;
+
}
diff --git a/math/mathcore/src/TMath.cxx b/math/mathcore/src/TMath.cxx
index dc6bd785cb1..5abc8d27e2c 100644
--- a/math/mathcore/src/TMath.cxx
+++ b/math/mathcore/src/TMath.cxx
@@ -26,12 +26,13 @@
#include "TString.h"
#include <Math/SpecFuncMathCore.h>
+#include <Math/PdfFuncMathCore.h>
+#include <Math/ProbFuncMathCore.h>
//const Double_t
// TMath::Pi = 3.14159265358979323846,
// TMath::E = 2.7182818284590452354;
-const Int_t kWorkMax = 100;
// Without this macro the THtml doc for TMath can not be generated
#if !defined(R__ALPHA) && !defined(R__SOLARIS) && !defined(R__ACC) && !defined(R__FBSD)
@@ -501,64 +502,8 @@ Double_t TMath::Landau(Double_t x, Double_t mpv, Double_t sigma, Bool_t norm)
// This function has been adapted from the CERNLIB routine G110 denlan.
// If norm=kTRUE (default is kFALSE) the result is divided by sigma
- static Double_t p1[5] = {0.4259894875,-0.1249762550, 0.03984243700, -0.006298287635, 0.001511162253};
- static Double_t q1[5] = {1.0 ,-0.3388260629, 0.09594393323, -0.01608042283, 0.003778942063};
-
- static Double_t p2[5] = {0.1788541609, 0.1173957403, 0.01488850518, -0.001394989411, 0.0001283617211};
- static Double_t q2[5] = {1.0 , 0.7428795082, 0.3153932961, 0.06694219548, 0.008790609714};
-
- static Double_t p3[5] = {0.1788544503, 0.09359161662,0.006325387654, 0.00006611667319,-0.000002031049101};
- static Double_t q3[5] = {1.0 , 0.6097809921, 0.2560616665, 0.04746722384, 0.006957301675};
-
- static Double_t p4[5] = {0.9874054407, 118.6723273, 849.2794360, -743.7792444, 427.0262186};
- static Double_t q4[5] = {1.0 , 106.8615961, 337.6496214, 2016.712389, 1597.063511};
-
- static Double_t p5[5] = {1.003675074, 167.5702434, 4789.711289, 21217.86767, -22324.94910};
- static Double_t q5[5] = {1.0 , 156.9424537, 3745.310488, 9834.698876, 66924.28357};
-
- static Double_t p6[5] = {1.000827619, 664.9143136, 62972.92665, 475554.6998, -5743609.109};
- static Double_t q6[5] = {1.0 , 651.4101098, 56974.73333, 165917.4725, -2815759.939};
-
- static Double_t a1[3] = {0.04166666667,-0.01996527778, 0.02709538966};
-
- static Double_t a2[2] = {-1.845568670,-4.284640743};
-
- if (sigma <= 0) return 0;
- Double_t v = (x-mpv)/sigma;
- Double_t u, ue, us, den;
- if (v < -5.5) {
- u = TMath::Exp(v+1.0);
- if (u < 1e-10) return 0.0;
- ue = TMath::Exp(-1/u);
- us = TMath::Sqrt(u);
- den = 0.3989422803*(ue/us)*(1+(a1[0]+(a1[1]+a1[2]*u)*u)*u);
- } else if(v < -1) {
- u = TMath::Exp(-v-1);
- den = TMath::Exp(-u)*TMath::Sqrt(u)*
- (p1[0]+(p1[1]+(p1[2]+(p1[3]+p1[4]*v)*v)*v)*v)/
- (q1[0]+(q1[1]+(q1[2]+(q1[3]+q1[4]*v)*v)*v)*v);
- } else if(v < 1) {
- den = (p2[0]+(p2[1]+(p2[2]+(p2[3]+p2[4]*v)*v)*v)*v)/
- (q2[0]+(q2[1]+(q2[2]+(q2[3]+q2[4]*v)*v)*v)*v);
- } else if(v < 5) {
- den = (p3[0]+(p3[1]+(p3[2]+(p3[3]+p3[4]*v)*v)*v)*v)/
- (q3[0]+(q3[1]+(q3[2]+(q3[3]+q3[4]*v)*v)*v)*v);
- } else if(v < 12) {
- u = 1/v;
- den = u*u*(p4[0]+(p4[1]+(p4[2]+(p4[3]+p4[4]*u)*u)*u)*u)/
- (q4[0]+(q4[1]+(q4[2]+(q4[3]+q4[4]*u)*u)*u)*u);
- } else if(v < 50) {
- u = 1/v;
- den = u*u*(p5[0]+(p5[1]+(p5[2]+(p5[3]+p5[4]*u)*u)*u)*u)/
- (q5[0]+(q5[1]+(q5[2]+(q5[3]+q5[4]*u)*u)*u)*u);
- } else if(v < 300) {
- u = 1/v;
- den = u*u*(p6[0]+(p6[1]+(p6[2]+(p6[3]+p6[4]*u)*u)*u)*u)/
- (q6[0]+(q6[1]+(q6[2]+(q6[3]+q6[4]*u)*u)*u)*u);
- } else {
- u = 1/(v-v*TMath::Log(v)/(v+1));
- den = u*u*(1+(a2[0]+a2[1]*u)*u);
- }
+ if (sigma <= 0) return 0;
+ Double_t den = ::ROOT::Math::landau_pdf(x, sigma, mpv);
if (!norm) return den;
return den/sigma;
}
@@ -671,22 +616,16 @@ Double_t TMath::Poisson(Double_t x, Double_t par)
Double_t TMath::PoissonI(Double_t x, Double_t par)
{
// compute the Poisson distribution function for (x,par)
- // This is a non-smooth function
+ // This is a non-smooth function.
+ // This function is equivalent to ROOT::Math::poisson_pdf
//Begin_Html
/*
<img src="gif/PoissonI.gif">
*/
//End_Html
-
- const Double_t kMaxInt = 2e6;
- if(x<0) return 0;
- if(x<1) return TMath::Exp(-par);
- Double_t gam;
Int_t ix = Int_t(x);
- if(x < kMaxInt) gam = TMath::Power(par,ix)/TMath::Gamma(ix+1);
- else gam = TMath::Power(par,x)/TMath::Gamma(x+1);
- return gam/TMath::Exp(par);
+ return Poisson(ix,par);
}
//______________________________________________________________________________
@@ -714,20 +653,7 @@ Double_t TMath::Prob(Double_t chi2,Int_t ndf)
else return 1;
}
- if (ndf==1) {
- Double_t v = 1.-Erf(Sqrt(chi2)/Sqrt(2.));
- return v;
- }
-
- // Gaussian approximation for large ndf
- Double_t q = Sqrt(2*chi2)-Sqrt(Double_t(2*ndf-1));
- if (ndf > 30 && q > 5) {
- Double_t v = 0.5*(1-Erf(q/Sqrt(2.)));
- return v;
- }
-
- // Evaluate the incomplete gamma function
- return (1-Gamma(0.5*ndf,0.5*chi2));
+ return ::ROOT::Math::chisquared_cdf_c(chi2,ndf);
}
//______________________________________________________________________________
@@ -1993,7 +1919,7 @@ Double_t TMath::StruveL0(Double_t x)
if (x<=20.) {
a0=2.0*x/pi;
- for (int i=1; i<=60;i++) {
+ for (i=1; i<=60;i++) {
r*=(x/(2*i+1))*(x/(2*i+1));
s+=r;
if(TMath::Abs(r/s)<1.e-12) break;
@@ -2355,9 +2281,10 @@ Double_t TMath::GammaDist(Double_t x, Double_t gamma, Double_t mu, Double_t beta
// gamma - shape parameter
// mu - location parameter
// beta - scale parameter
- // The formula was taken from "Engineering Statistics Handbook" on site
+ //
+ // The definition can be found in "Engineering Statistics Handbook" on site
// http://www.itl.nist.gov/div898/handbook/eda/section3/eda366b.htm
- // Implementation by Anna Kreshuk.
+ // use now implementation in ROOT::Math::gamma_pdf
//Begin_Html
/*
<img src="gif/gammadist.gif">
@@ -2368,10 +2295,7 @@ Double_t TMath::GammaDist(Double_t x, Double_t gamma, Double_t mu, Double_t beta
Error("TMath::GammaDist", "illegal parameter values");
return 0;
}
- Double_t temp = (x-mu)/beta;
- Double_t temp2 = beta * TMath::Gamma(gamma);
- Double_t result = (TMath::Power(temp, gamma-1) * TMath::Exp(-temp))/temp2;
- return result;
+ return ::ROOT::Math::gamma_pdf(x, gamma, beta, mu);
}
//______________________________________________________________________________
@@ -2418,7 +2342,7 @@ Double_t TMath::LogNormal(Double_t x, Double_t sigma, Double_t theta, Double_t m
// m is the scale parameter
// The formula was taken from "Engineering Statistics Handbook" on site
// http://www.itl.nist.gov/div898/handbook/eda/section3/eda3669.htm
- // Implementation by Anna Kreshuk.
+ // Implementation using ROOT::Math::lognormal_pdf
//Begin_Html
/*
<img src="gif/lognormal.gif">
@@ -2429,11 +2353,11 @@ Double_t TMath::LogNormal(Double_t x, Double_t sigma, Double_t theta, Double_t m
Error("TMath::Lognormal", "illegal parameter values");
return 0;
}
- Double_t templog2 = TMath::Log((x-theta)/m)*TMath::Log((x-theta)/m);
- Double_t temp1 = TMath::Exp(-templog2/(2*sigma*sigma));
- Double_t temp2 = (x-theta)*sigma*TMath::Sqrt(2*TMath::Pi());
+ // lognormal_pdf uses same definition of http://en.wikipedia.org/wiki/Log-normal_distribution
+ // where mu = log(m)
+
+ return ::ROOT::Math::lognormal_pdf(x, TMath::Log(m), sigma, theta);
- return temp1/temp2;
}
//______________________________________________________________________________
@@ -2781,62 +2705,8 @@ Double_t TMath::LandauI(Double_t x)
//The algorithm was taken from the Cernlib function dislan(G110)
//Reference: K.S.Kolbig and B.Schorr, "A program package for the Landau
//distribution", Computer Phys.Comm., 31(1984), 97-111
-
- Double_t p1[] = {0.2514091491e+0,-0.6250580444e-1, 0.1458381230e-1, -0.2108817737e-2, 0.7411247290e-3};
- Double_t q1[] = {1.0 ,-0.5571175625e-2, 0.6225310236e-1, -0.3137378427e-2, 0.1931496439e-2};
-
- Double_t p2[] = {0.2868328584e+0, 0.3564363231e+0, 0.1523518695e+0, 0.2251304883e-1};
- Double_t q2[] = {1.0 , 0.6191136137e+0, 0.1720721448e+0, 0.2278594771e-1};
-
- Double_t p3[] = {0.2868329066e+0, 0.3003828436e+0, 0.9950951941e-1, 0.8733827185e-2};
- Double_t q3[] = {1.0 , 0.4237190502e+0, 0.1095631512e+0, 0.8693851567e-2};
-
- Double_t p4[] = {0.1000351630e+1, 0.4503592498e+1, 0.1085883880e+2, 0.7536052269e+1};
- Double_t q4[] = {1.0 , 0.5539969678e+1, 0.1933581111e+2, 0.2721321508e+2};
-
- Double_t p5[] = {0.1000006517e+1, 0.4909414111e+2, 0.8505544753e+2, 0.1532153455e+3};
- Double_t q5[] = {1.0 , 0.5009928881e+2, 0.1399819104e+3, 0.4200002909e+3};
-
- Double_t p6[] = {0.1000000983e+1, 0.1329868456e+3, 0.9162149244e+3, -0.9605054274e+3};
- Double_t q6[] = {1.0 , 0.1339887843e+3, 0.1055990413e+4, 0.5532224619e+3};
-
- Double_t a1[] = {0, -0.4583333333e+0, 0.6675347222e+0,-0.1641741416e+1};
-
- Double_t a2[] = {0, 1.0 ,-0.4227843351e+0,-0.2043403138e+1};
-
- Double_t u, v;
- Double_t lan;
- v = x;
- if (v < -5.5) {
- u = TMath::Exp(v+1);
- lan = 0.3989422803*TMath::Exp(-1./u)*TMath::Sqrt(u)*(1+(a1[1]+(a1[2]+a1[3]*u)*u)*u);
- }
- else if (v < -1 ) {
- u = TMath::Exp(-v-1);
- lan = (TMath::Exp(-u)/TMath::Sqrt(u))*(p1[0]+(p1[1]+(p1[2]+(p1[3]+p1[4]*v)*v)*v)*v)/
- (q1[0]+(q1[1]+(q1[2]+(q1[3]+q1[4]*v)*v)*v)*v);
- }
- else if (v < 1)
- lan = (p2[0]+(p2[1]+(p2[2]+p2[3]*v)*v)*v)/(q2[0]+(q2[1]+(q2[2]+q2[3]*v)*v)*v);
- else if (v < 4)
- lan = (p3[0]+(p3[1]+(p3[2]+p3[3]*v)*v)*v)/(q3[0]+(q3[1]+(q3[2]+q3[3]*v)*v)*v);
- else if (v < 12) {
- u = 1./v;
- lan = (p4[0]+(p4[1]+(p4[2]+p4[3]*u)*u)*u)/(q4[0]+(q4[1]+(q4[2]+q4[3]*u)*u)*u);
- }
- else if (v < 50) {
- u = 1./v;
- lan = (p5[0]+(p5[1]+(p5[2]+p5[3]*u)*u)*u)/(q5[0]+(q5[1]+(q5[2]+q5[3]*u)*u)*u);
- }
- else if (v < 300) {
- u = 1./v;
- lan = (p6[0]+(p6[1]+(p6[2]+p6[3]*u)*u)*u)/(q6[0]+(q6[1]+(q6[2]+q6[3]*u)*u)*u);
- }
- else {
- u = 1./(v-v*TMath::Log(v)/(v+1));
- lan = 1-(a2[1]+(a2[2]+a2[3]*u)*u)*u;
- }
- return lan;
+
+ return ::ROOT::Math::landau_cdf(x);
}
diff --git a/math/mathcore/test/Makefile b/math/mathcore/test/Makefile
index 62f5193dcd6..d330135d982 100644
--- a/math/mathcore/test/Makefile
+++ b/math/mathcore/test/Makefile
@@ -12,10 +12,9 @@ include $(ROOTSYS)/test/Makefile.arch
#------------------------------------------------------------------------------
ifeq ($(PLATFORM),win32)
-EXTRALIBS = -include:_G__cpp_setupG__MathCore "$(ROOTSYS)/lib/libMathCore.lib" \
-"$(ROOTSYS)/lib/libMathMore.lib"
+EXTRALIBS = -include:_G__cpp_setupG__MathCore "$(ROOTSYS)/lib/libMathMore.lib"
else
-EXTRALIBS = -lMathCore -lMathMore
+EXTRALIBS = -lMathMore
CXXFLAGS += -g
endif
@@ -117,58 +116,56 @@ PROGRAMS =$(SPECFUNBETA) $(SPECFUNBETAI) $(SPECFUNGAMMA) $(SPECFUNERF) $(T
all: $(PROGRAMS)
-
-
$(STRESSTF1): $(STRESSTF1OBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(STRESSTMATH): $(STRESSTMATHOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(TESTTMATH): $(TESTTMATHOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(BSEARCHTIME): $(BSEARCHTIMEOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(TESTBSEARCH): $(TESTBSEARCHOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(TESTSORT): $(TESTSORTOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(TESTSORTORDER): $(TESTSORTORDEROBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(SPECFUNERF): $(SPECFUNERFOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(SPECFUNGAMMA): $(SPECFUNGAMMAOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(SPECFUNBETA): $(SPECFUNBETAOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(SPECFUNBETAI): $(SPECFUNBETAIOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(INTEGRATION): $(INTEGRATIONOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(INTEGRATIONMULTI): $(INTEGRATIONMULTIOBJ)
- $(LD) $(LDFLAGS) $^ $(LIBS) $(EXTRALIBS) $(OutPutOpt)$@
+ $(LD) $(LDFLAGS) $^ $(LIBS) $(OutPutOpt)$@
@echo "$@ done"
$(ROOTFINDER): $(ROOTFINDEROBJ)
diff --git a/math/mathcore/test/binarySearchTime.cxx b/math/mathcore/test/binarySearchTime.cxx
index a06cd91158a..8d53b5a831f 100644
--- a/math/mathcore/test/binarySearchTime.cxx
+++ b/math/mathcore/test/binarySearchTime.cxx
@@ -131,7 +131,6 @@ int main(int argc, char **argv)
if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
{
- cout << "No graphics mode!" << endl;
showGraphics = false;
}
diff --git a/math/mathcore/test/stdsort.cxx b/math/mathcore/test/stdsort.cxx
index 48a6e1ebb33..75a2509d009 100644
--- a/math/mathcore/test/stdsort.cxx
+++ b/math/mathcore/test/stdsort.cxx
@@ -23,6 +23,7 @@ const int maxsize = 500;
const int increment = 10;
const int arraysize = (maxsize-minsize)/10 + 1;
+bool showGraphics = true;
template<typename T>
struct Compare {
@@ -88,41 +89,67 @@ void stdsort()
for ( int i = minsize; i <= maxsize; i += increment)
cout << tM[(i-minsize)/10] << ' ' << tS[(i-minsize)/10] << endl;
- TCanvas* c1 = new TCanvas("c1", "Comparision of Sorting Time", 600, 400);
- TH2F* hpx = new TH2F("hpx", "Comparision of Sorting Time", arraysize, minsize, maxsize, arraysize, 0,tM[arraysize-1]);
- hpx->SetStats(kFALSE);
- hpx->Draw();
-
- TGraph* gM = new TGraph(arraysize, &index[0], &tM[0]);
- gM->SetLineColor(2);
- gM->SetLineWidth(3);
- gM->SetTitle("TMath::Sort()");
- gM->Draw("SAME");
-
- TGraph* gS = new TGraph(arraysize, &index[0], &tS[0]);
- gS->SetLineColor(3);
- gS->SetLineWidth(3);
- gS->SetTitle("std::sort()");
- gS->Draw("SAME");
-
- TLegend* legend = new TLegend(0.15,0.72,0.4,0.86);
- legend->AddEntry(gM, "TMath::Sort()");
- legend->AddEntry(gS, "std::sort()");
- legend->Draw();
-
- hpx->GetXaxis()->SetTitle("Array Size");
- hpx->GetYaxis()->SetTitle("Time");
-
-
- c1->Show();
+ if ( showGraphics )
+ {
+ TCanvas* c1 = new TCanvas("c1", "Comparision of Sorting Time", 600, 400);
+ TH2F* hpx = new TH2F("hpx", "Comparision of Sorting Time", arraysize, minsize, maxsize, arraysize, 0,tM[arraysize-1]);
+ hpx->SetStats(kFALSE);
+ hpx->Draw();
+
+ TGraph* gM = new TGraph(arraysize, &index[0], &tM[0]);
+ gM->SetLineColor(2);
+ gM->SetLineWidth(3);
+ gM->SetTitle("TMath::Sort()");
+ gM->Draw("SAME");
+
+ TGraph* gS = new TGraph(arraysize, &index[0], &tS[0]);
+ gS->SetLineColor(3);
+ gS->SetLineWidth(3);
+ gS->SetTitle("std::sort()");
+ gS->Draw("SAME");
+
+ TLegend* legend = new TLegend(0.15,0.72,0.4,0.86);
+ legend->AddEntry(gM, "TMath::Sort()");
+ legend->AddEntry(gS, "std::sort()");
+ legend->Draw();
+
+ hpx->GetXaxis()->SetTitle("Array Size");
+ hpx->GetYaxis()->SetTitle("Time");
+
+
+ c1->Show();
+ }
}
int main(int argc, char **argv)
{
- TApplication theApp("App",&argc,argv);
+ if ( argc > 1 && argc != 2 )
+ {
+ cerr << "Usage: " << argv[0] << " [-ng]\n";
+ cerr << " where:\n";
+ cerr << " -ng : no graphics mode";
+ cerr << endl;
+ exit(1);
+ }
+
+ if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
+ {
+ showGraphics = false;
+ }
+
+ TApplication* theApp = 0;
+ if ( showGraphics )
+ theApp = new TApplication("App",&argc,argv);
+
stdsort();
- theApp.Run();
+
+ if ( showGraphics )
+ {
+ theApp->Run();
+ delete theApp;
+ theApp = 0;
+ }
return 0;
}
diff --git a/math/mathcore/test/stressTMath.cxx b/math/mathcore/test/stressTMath.cxx
index 9db07d8b0e0..cc905b10a7d 100644
--- a/math/mathcore/test/stressTMath.cxx
+++ b/math/mathcore/test/stressTMath.cxx
@@ -5,112 +5,159 @@
#include <TStopwatch.h>
#include <TMath.h>
-using namespace std;
-using namespace TMath;
+using std::cout;
+using std::endl;
-const unsigned int SIZE = 100000;
const unsigned int NUMTEST = 500;
-//#define DEBUG
+// #define DEBUG
template <typename T> T randD() {
- static TRandom2 r( time( 0 ) );
+ // use default seed to get same sequence
+ static TRandom2 r;
return (T) r.Uniform(-500,500);
}
-template <typename T> double stressVector(const char* type)
+double Time(TStopwatch & w) {
+ //return w.CpuTime();
+ return w.RealTime();
+}
+
+
+template <typename T> double stressVector(unsigned int size, const char* type)
{
cout << "Generating random vector of '"
- << type << "'..." << endl;
+ << type << "' and size " << size << " ..." << endl << endl;
double totalTime = 0;
+ double totalUnitTime = 0;
- T *vector = new T[SIZE];
- generate(vector, &vector[SIZE], randD<T>);
+ T *vector = new T[size];
+ std::generate(vector, &vector[size], randD<T>);
#ifdef DEBUG
- for ( unsigned int i = 0; i < SIZE; ++i )
+ for ( unsigned int i = 0; i < size; ++i )
cout << vector[i] << " " << endl;
#endif
TStopwatch w;
+ std::cout.precision(6);
+ unsigned int ntest = 3 * NUMTEST;
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- MinElement(SIZE, vector);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::MinElement(size, vector);
w.Stop();
- cout << "MinMaxElement() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "MinMaxElement() \tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+ ntest = 3 * NUMTEST;
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- LocMin(SIZE, vector);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::LocMin(size, vector);
w.Stop();
- cout << "LocMin/Max() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "LocMin/Max() \t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+ ntest = 10 * NUMTEST;
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- Mean(SIZE, vector);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::Mean(size, vector);
w.Stop();
- cout << "Mean() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "Mean() \t\t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+ ntest = (unsigned int) ( NUMTEST/2.5 );
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- Median(SIZE, vector);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::Median(size, vector);
w.Stop();
- cout << "Median() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "Median() \t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+ ntest = (unsigned int) ( 10 * NUMTEST );
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- RMS(SIZE, vector);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::RMS(size, vector);
w.Stop();
- cout << "RMS() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "RMS() \t\t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+ ntest = (unsigned int) ( NUMTEST/2.5 );
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- GeomMean(SIZE, vector);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::GeomMean(size, vector);
w.Stop();
- cout << "GeomMean() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "GeomMean() \t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
- Int_t index[SIZE];
+ Int_t * index =new Int_t[size];
+ ntest = NUMTEST/10;
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- Sort(SIZE, vector, index, kFALSE);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::Sort(size, vector, index, kFALSE);
w.Stop();
- cout << "Sort() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "Sort() \t\t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+ std::sort(vector, vector + size);
+#ifdef DEBUG
+ for ( unsigned int i = 0; i < size; ++i )
+ cout << vector[i] << " " << endl;
+#endif
+ ntest = 20000*NUMTEST;
w.Start( kTRUE );
- for ( unsigned int i = 0; i < NUMTEST; ++i )
- BinarySearch(SIZE, vector, vector[1]);
+ for ( unsigned int i = 0; i < ntest; ++i )
+ TMath::BinarySearch(size, vector, vector[ i % size ]);
w.Stop();
- cout << "BinarySearch() Time: " << w.CpuTime()/NUMTEST << endl;
- totalTime += w.CpuTime()/NUMTEST;
+ cout << "BinarySearch() \t\tTotal Time: " << Time(w) << " (s)\t\t"
+ << " Time/call: " << Time(w)/(ntest)*1.E6 << " (microsec)" << endl;
+ totalUnitTime += Time(w)/ntest;
+ totalTime += Time(w);
+
+ cout << "\nTotal Time : " << totalTime << " (s)\n"
+ << "Total Time/call : " << totalUnitTime*1.E3 << " (ms)\n" << endl;
- cout << "Total Time: " << totalTime << "\n" << endl;
+ delete [] vector;
+ delete [] index;
- return totalTime;
+ return totalUnitTime;
}
-void stressTMath()
+void stressTMath(unsigned int size, const char * type)
{
double totalTime = 0;
cout << "Stress Test Start..." << endl;
- totalTime += stressVector<Short_t>("Short_t");
- totalTime += stressVector<Int_t>("Int_t");
- totalTime += stressVector<Float_t>("Float_t");
- totalTime += stressVector<Double_t>("Double_t");
- totalTime += stressVector<Long_t>("Long_t");
- totalTime += stressVector<Long64_t>("Long64_t");
+ if ( strcmp(type, "Short_t") == 0 )
+ totalTime += stressVector<Short_t>(size, type);
+ else if ( strcmp(type, "Int_t") == 0 )
+ totalTime += stressVector<Int_t>(size, type);
+ else if ( strcmp(type, "Float_t") == 0 )
+ totalTime += stressVector<Float_t>(size, type);
+ else if ( strcmp(type, "Long_t") == 0 )
+ totalTime += stressVector<Long_t>(size, type);
+ else if ( strcmp(type, "Long64_t") == 0 )
+ totalTime += stressVector<Long64_t>(size, type);
+ else
+ totalTime += stressVector<Double_t>(size, "Double_t");
- cout << "Total Test Time: " << totalTime << "\n" << endl;
+ //cout << "Total Test Time: " << totalTime << "\n" << endl;
cout << "End of Stress Test..." << endl;
@@ -118,9 +165,34 @@ void stressTMath()
}
-int main()
+int main(int argc, char* argv[])
{
- stressTMath();
+ // Default size and data type
+ unsigned int size = 100000;
+ const char * type = "Double_t";
+
+ if ( argc > 1 ) {
+ if (strcmp(argv[1], "-h") == 0) {
+ cout << "Usage: " << argv[0]
+ << " [TYPE OF ARRAY] [SIZE OF ARRAY]\n\n"
+ << "where [TYPE OF ARRAY] is one of the following:\n"
+ << "\t\"Short_t\"\n"
+ << "\t\"Int_t\"\n"
+ << "\t\"Float_t\"\n"
+ << "\t\"Long_t\"\n"
+ << "\t\"Long64_t\"\n"
+ << "\t \"Double_t\"\n"
+ << endl;
+ return 1;
+ }
+ type = argv[1];
+ }
+
+
+ if ( argc > 2 )
+ size = (unsigned int) atoi(argv[2]);
+
+ stressTMath(size, type);
return 0;
}
diff --git a/math/mathcore/test/testIntegration.cxx b/math/mathcore/test/testIntegration.cxx
index c6e256d2723..e8625abc499 100644
--- a/math/mathcore/test/testIntegration.cxx
+++ b/math/mathcore/test/testIntegration.cxx
@@ -1,9 +1,13 @@
#include "Math/Integrator.h"
#include "Math/IntegratorMultiDim.h"
#include "Math/AllIntegrationTypes.h"
-#include "Math/WrappedFunction.h"
+#include "Math/Functor.h"
#include "Math/GaussIntegrator.h"
+#include <cmath>
+
+const double ERRORLIMIT = 1E-3;
+
double f(double x) {
return x;
}
@@ -13,34 +17,45 @@ double f2(const double * x) {
}
-void testIntegration1D() {
+int testIntegration1D() {
- ROOT::Math::WrappedFunction<> wf(f);
+ const double RESULT = 0.5;
+ int status = 0;
+
+ ROOT::Math::Functor1D wf(&f);
ROOT::Math::Integrator ig(ROOT::Math::IntegrationOneDim::ADAPTIVESINGULAR);
ig.SetFunction(wf);
double val = ig.Integral(0,1);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
ROOT::Math::Integrator ig2(ROOT::Math::IntegrationOneDim::NONADAPTIVE);
ig2.SetFunction(wf);
val = ig2.Integral(0,1);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
ROOT::Math::Integrator ig3(wf, ROOT::Math::IntegrationOneDim::ADAPTIVE);
val = ig3.Integral(0,1);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
//ROOT::Math::GaussIntegrator ig4;
ROOT::Math::Integrator ig4(ROOT::Math::IntegrationOneDim::GAUSS);
ig4.SetFunction(wf);
val = ig4.Integral(0,1);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
+ return status;
}
-void testIntegrationMultiDim() {
+int testIntegrationMultiDim() {
+
+ const double RESULT = 1.0;
+ int status = 0;
- ROOT::Math::WrappedMultiFunction<> wf(f2,2);
+ ROOT::Math::Functor wf(&f2,2);
double a[2] = {0,0};
double b[2] = {1,1};
@@ -48,26 +63,33 @@ void testIntegrationMultiDim() {
ig.SetFunction(wf);
double val = ig.Integral(a,b);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
ROOT::Math::IntegratorMultiDim ig2(ROOT::Math::IntegrationMultiDim::VEGAS);
ig2.SetFunction(wf);
val = ig2.Integral(a,b);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
ROOT::Math::IntegratorMultiDim ig3(wf,ROOT::Math::IntegrationMultiDim::PLAIN);
val = ig3.Integral(a,b);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
ROOT::Math::IntegratorMultiDim ig4(wf,ROOT::Math::IntegrationMultiDim::MISER);
val = ig4.Integral(a,b);
std::cout << "integral result is " << val << std::endl;
+ status += std::fabs(val-RESULT) > ERRORLIMIT;
-
-
+ return status;
}
int main() {
- testIntegration1D();
- testIntegrationMultiDim();
+ int status = 0;
+
+ status += testIntegration1D();
+ status += testIntegrationMultiDim();
+
+ return status;
}
diff --git a/math/mathcore/test/testIntegrationMultiDim.cxx b/math/mathcore/test/testIntegrationMultiDim.cxx
index 514cb0752ac..c3d1494c600 100644
--- a/math/mathcore/test/testIntegrationMultiDim.cxx
+++ b/math/mathcore/test/testIntegrationMultiDim.cxx
@@ -196,7 +196,6 @@ int main(int argc, char **argv)
if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
{
- std::cout << "No graphics mode!" << std::endl;
showGraphics = false;
}
@@ -214,4 +213,5 @@ int main(int argc, char **argv)
}
return 0;
+
}
diff --git a/math/mathcore/test/testRootFinder.cxx b/math/mathcore/test/testRootFinder.cxx
index 25ec4801937..a29b525cecd 100644
--- a/math/mathcore/test/testRootFinder.cxx
+++ b/math/mathcore/test/testRootFinder.cxx
@@ -7,6 +7,7 @@
#include <iostream>
+const double ERRORLIMIT = 1E-8;
const int iterTest = 10000;
int myfuncCalls = 0;
@@ -15,20 +16,23 @@ double myfunc ( double x ) {
return x*x - 5;
}
-void printStats(TStopwatch& timer, double root) {
+int printStats(TStopwatch& timer, double root) {
//std::cout << "Return code: " << status << std::endl;
+ double difference = root - sqrt(5.0);
std::cout << "Result: " << root << std::endl;
- std::cout << "Exact result: " << sqrt(5.0) << " difference: " << root - sqrt(5.0) << std::endl;
+ std::cout << "Exact result: " << sqrt(5.0) << " difference: " << difference << std::endl;
std::cout << "Time: " << timer.RealTime()/(double) iterTest << std::endl;
std::cout << "Number of calls to function: " << myfuncCalls/iterTest << std::endl;
+ return difference > ERRORLIMIT;
}
-void testRootFinder() {
+int testRootFinder() {
TStopwatch timer;
double root;
+ int status = 0;
ROOT::Math::Polynomial polyf(2);
std::vector<double> p(3);
@@ -51,7 +55,7 @@ void testRootFinder() {
}
timer.Stop();
std::cout << "RootFinder Stats:" << std::endl;
- printStats(timer, root);
+ status += printStats(timer, root);
TF1* f1 = new TF1("f1", "x*x - 5", 0, 5);
@@ -63,13 +67,14 @@ void testRootFinder() {
}
timer.Stop();
std::cout << "\nTF1 Stats:" << std::endl;
- printStats(timer, root);
+ status += printStats(timer, root);
+ return status;
}
int main() {
- testRootFinder();
- return 0;
+ int status = testRootFinder();
+ return status;
}
diff --git a/math/mathcore/test/testSortOrder.cxx b/math/mathcore/test/testSortOrder.cxx
index 729780e88e6..69a4a93849b 100644
--- a/math/mathcore/test/testSortOrder.cxx
+++ b/math/mathcore/test/testSortOrder.cxx
@@ -20,33 +20,36 @@ const int arraysize = (maxsize-minsize)/10 + 1;
template<typename T>
struct CompareDesc {
- CompareDesc(const T * d) : fData(d) {}
+ CompareDesc(T d) : fData(d) {}
bool operator()(int i1, int i2) {
- return fData[i1] > fData[i2];
+ return *(fData + i1) > *(fData + i2);
}
- const T * fData;
+ T fData;
};
template<typename T>
struct CompareAsc {
- CompareAsc(const T * d) : fData(d) {}
+ CompareAsc(T d) : fData(d) {}
bool operator()(int i1, int i2) {
- return fData[i1] < fData[i2];
+ return *(fData + i1) < *(fData + i2);
}
- const T * fData;
+ T fData;
};
#endif
-template <typename T> void testSort(const int n)
+template <typename T> bool testSort(const int n)
{
vector<T> k(n);
- vector<T> index(n);
+ vector<T> indexM(n);
+ vector<T> indexS(n);
+
+ bool equals = true;
TRandom2 r( time( 0 ) );
@@ -57,55 +60,61 @@ template <typename T> void testSort(const int n)
}
cout << endl;
-
-
- for(Int_t i = 0; i < n; i++) { index[i] = i; }
- TMath::Sort(n,&k[0],&index[0],kTRUE);
+ for(Int_t i = 0; i < n; i++) { indexM[i] = i; }
+ TMath::Sort(n,&k[0],&indexM[0],kTRUE);
cout << "TMath[kTRUE]\n\tindex = ";
for ( Int_t i = 0; i < n; ++i )
- cout << k[index[i]] << ' ';
+ cout << k[indexM[i]] << ' ';
cout << endl;
-
-
- for(Int_t i = 0; i < n; i++) { index[i] = i; }
- std::sort(&index[0],&index[n], CompareDesc<T>(&k[0]) );
+ for(Int_t i = 0; i < n; i++) { indexS[i] = i; }
+ std::sort(&indexS[0],&indexS[n], CompareDesc<const T*>(&k[0]) );
cout << "std::sort[CompareDesc]\n\tindex = ";
for ( Int_t i = 0; i < n; ++i )
- cout << k[index[i]] << ' ';
+ cout << k[indexS[i]] << ' ';
cout << endl;
+ equals &= std::equal(indexM.begin(), indexM.end(), indexS.begin());
+ cout << "Equals? " << (char*) (equals?"OK":"FAILED") << endl;
+
- for(Int_t i = 0; i < n; i++) { index[i] = i; }
- TMath::Sort(n,&k[0],&index[0],kFALSE);
+ for(Int_t i = 0; i < n; i++) { indexM[i] = i; }
+ TMath::Sort(n,&k[0],&indexM[0],kFALSE);
cout << "TMath[kFALSE]\n\tindex = ";
for ( Int_t i = 0; i < n; ++i )
- cout << k[index[i]] << ' ';
+ cout << k[indexM[i]] << ' ';
cout << endl;
-
-
- for(Int_t i = 0; i < n; i++) { index[i] = i; }
- std::sort(&index[0],&index[n], CompareAsc<T>(&k[0]) );
+ for(Int_t i = 0; i < n; i++) { indexS[i] = i; }
+ std::sort(&indexS[0],&indexS[n], CompareAsc<const T*>(&k[0]) );
cout << "std::sort[CompareAsc]\n\tindex = ";
for ( Int_t i = 0; i < n; ++i )
- cout << k[index[i]] << ' ';
+ cout << k[indexS[i]] << ' ';
cout << endl;
+
+
+ equals &= std::equal(indexM.begin(), indexM.end(), indexS.begin());
+ cout << "Equals? " << (char*) (equals?"OK":"FAILED") << endl;
+
+ return equals;
}
-void stdsort()
+bool stdsort()
{
- testSort<Int_t>(20);
+ return testSort<Int_t>(20);
}
-int main(int argc, char **argv)
+int main(int /* argc */ , char ** /* argv */ )
{
- stdsort();
+ bool equals = stdsort();
- return 0;
+ if ( !equals )
+ return 1;
+ else
+ return 0;
}
diff --git a/math/mathcore/test/testSpecFuncBeta.cxx b/math/mathcore/test/testSpecFuncBeta.cxx
index 54e27c0f2e6..2d73db2ecf2 100644
--- a/math/mathcore/test/testSpecFuncBeta.cxx
+++ b/math/mathcore/test/testSpecFuncBeta.cxx
@@ -2,6 +2,8 @@
#include <fstream>
#include <vector>
+#include <cmath>
+
#include <TMath.h>
#include <Math/SpecFunc.h>
@@ -12,11 +14,13 @@
#include <TGraph.h>
#include <TLegend.h>
+const double ERRORLIMIT = 1E-8;
const double MIN = -2.5;
const double MAX = +2.5;
const double INCREMENT = 0.01;
-const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT);
-inline int arrayindex(double i) { return ARRAYSIZE - (int) ( (MAX - i) / INCREMENT ) -1 ; };
+const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT) + 1;
+
+bool showGraphics = true;
using namespace std;
@@ -31,57 +35,99 @@ TGraph* drawPoints(Double_t x[], Double_t y[], int color, int style = 1)
return g;
}
-void testSpecFuncBeta()
+int testSpecFuncBeta()
{
vector<Double_t> x( ARRAYSIZE );
vector<Double_t> yb( ARRAYSIZE );
vector<Double_t> ymb( ARRAYSIZE );
- TCanvas* c1 = new TCanvas("c1", "BetaFunction", 600, 400);
- TH2F* hpx = new TH2F("hpx", "BetaFunction(p,b)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, 0, 5);
- hpx->SetStats(kFALSE);
- hpx->Draw();
+ int status = 0;
int color = 2;
-
TGraph *gb, *gmb;
+ TCanvas* c1 = new TCanvas("c1", "BetaFunction", 600, 400);
+ TH2F* hpx;
+ {
+ hpx = new TH2F("hpx", "BetaFunction(p,b)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, 0, 5);
+ hpx->SetStats(kFALSE);
+ hpx->Draw();
+ }
+
for ( double b = 0.9; b < 2; b+=0.4)
{
cout << "** b = " << b << " **" << endl;
+ unsigned int index = 0;
for ( double i = MIN; i < MAX; i += INCREMENT )
{
// cout << "i:"; cout.width(5); cout << i
-// << " index: "; cout.width(5); cout << arrayindex(i)
+// << " index: "; cout.width(5); cout << index
// << " TMath::Beta(p,b): "; cout.width(10); cout << TMath::Beta(i,b)
// << " ROOT::Math::beta(p,b): "; cout.width(10); cout << ROOT::Math::beta(i,b)
// << endl;
- x[arrayindex(i)] = i;
- yb[arrayindex(i)] = TMath::Beta(i,b);
- ymb[arrayindex(i)] = ROOT::Math::beta(i,b);
+ x[index] = i;
+ yb[index] = TMath::Beta(i,b);
+ ymb[index] = ROOT::Math::beta(i,b);
+ if ( std::fabs( yb[index] - ymb[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " b " << b
+ << " yb[index] " << yb[index]
+ << " ymb[index] " << ymb[index]
+ << " " << std::fabs( yb[index] - ymb[index] )
+ << endl;
+ status += 1;
+ }
+ index += 1;
}
gb = drawPoints(&x[0], &yb[0], color++);
gmb = drawPoints(&x[0], &ymb[0], color++, 7);
}
- TLegend* legend = new TLegend(0.61,0.72,0.86,0.86);
- legend->AddEntry(gb, "TMath::Beta()");
- legend->AddEntry(gmb, "ROOT::Math::beta()");
- legend->Draw();
-
- c1->Show();
+ if ( showGraphics )
+ {
+ TLegend* legend = new TLegend(0.61,0.72,0.86,0.86);
+ legend->AddEntry(gb, "TMath::Beta()");
+ legend->AddEntry(gmb, "ROOT::Math::beta()");
+ legend->Draw();
+
+ c1->Show();
+ }
cout << "Test Done!" << endl;
- return;
+ return status;
}
int main(int argc, char **argv)
{
- TApplication theApp("App",&argc,argv);
- testSpecFuncBeta();
- theApp.Run();
+ if ( argc > 1 && argc != 2 )
+ {
+ cerr << "Usage: " << argv[0] << " [-ng]\n";
+ cerr << " where:\n";
+ cerr << " -ng : no graphics mode";
+ cerr << endl;
+ exit(1);
+ }
+
+ if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
+ {
+ showGraphics = false;
+ }
+
+ TApplication* theApp = 0;
+ if ( showGraphics )
+ theApp = new TApplication("App",&argc,argv);
+
+ int status = testSpecFuncBeta();
+
+ if ( showGraphics )
+ {
+ theApp->Run();
+ delete theApp;
+ theApp = 0;
+ }
- return 0;
+ return status;
}
diff --git a/math/mathcore/test/testSpecFuncBetaI.cxx b/math/mathcore/test/testSpecFuncBetaI.cxx
index 42ea4462b28..0e50fef5c0f 100644
--- a/math/mathcore/test/testSpecFuncBetaI.cxx
+++ b/math/mathcore/test/testSpecFuncBetaI.cxx
@@ -2,6 +2,8 @@
#include <fstream>
#include <vector>
+#include <cmath>
+
#include <TMath.h>
#include <Math/SpecFunc.h>
@@ -12,11 +14,13 @@
#include <TGraph.h>
#include <TLegend.h>
+const double ERRORLIMIT = 1E-8;
const double MIN = 0;
const double MAX = 1;
const double INCREMENT = 0.01;
-const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT);
-inline int arrayindex(double i) { return ARRAYSIZE - (int) ( (MAX - i) / INCREMENT ) -1 ; };
+const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT) + 1;
+
+bool showGraphics = true;
using namespace std;
@@ -31,57 +35,96 @@ TGraph* drawPoints(Double_t x[], Double_t y[], int color, int style = 1)
return g;
}
-void testSpecFuncBetaI()
+int testSpecFuncBetaI()
{
vector<Double_t> x( ARRAYSIZE );
vector<Double_t> yb( ARRAYSIZE );
vector<Double_t> ymb( ARRAYSIZE );
-// ofstream outputFile ("values.txt");
-
- TCanvas* c1 = new TCanvas("c1", "Two Graphs", 600, 400);
- TH2F* hpx = new TH2F("hpx", "Two Graphs(hpx)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, 0, 5);
- hpx->SetStats(kFALSE);
- hpx->Draw();
+ int status = 0;
- int color = 2;
-
- TGraph *gb, *gmb;
double b = 0.2, a= 0.9;
cout << "** b = " << b << " **" << endl;
+ unsigned int index = 0;
for ( double i = MIN; i < MAX; i += INCREMENT )
{
- cout << "i:"; cout.width(5); cout << i
- << " index: "; cout.width(5); cout << arrayindex(i)
- << " TMath::BetaIncomplete(x,a,b): "; cout.width(10); cout << TMath::BetaIncomplete(i,a,b)
- << " ROOT::Math::inc_beta(a,a,b): "; cout.width(10); cout << ROOT::Math::inc_beta(i,a,b)
- << endl;
+// cout << "i:"; cout.width(5); cout << i
+// << " index: "; cout.width(5); cout << index
+// << " TMath::BetaIncomplete(x,a,b): "; cout.width(10); cout << TMath::BetaIncomplete(i,a,b)
+// << " ROOT::Math::inc_beta(a,a,b): "; cout.width(10); cout << ROOT::Math::inc_beta(i,a,b)
+// << endl;
- x[arrayindex(i)] = i;
- yb[arrayindex(i)] = TMath::BetaIncomplete(i,a,b);
- ymb[arrayindex(i)] = ROOT::Math::inc_beta(i,a,b);
+ x[index] = i;
+
+ yb[index] = TMath::BetaIncomplete(i,a,b);
+ ymb[index] = ROOT::Math::inc_beta(i,a,b);
+ if ( std::fabs( yb[index] - ymb[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yb[index] " << yb[index]
+ << " ymb[index] " << ymb[index]
+ << " " << std::fabs( yb[index] - ymb[index] )
+ << endl;
+ status += 1;
+ }
+ index += 1;
}
-
- gb = drawPoints(&x[0], &yb[0], color++);
- gmb = drawPoints(&x[0], &ymb[0], color++, 7);
- TLegend* legend = new TLegend(0.61,0.72,0.86,0.86);
- legend->AddEntry(gb, "TMath::BetaIncomplete()");
- legend->AddEntry(gmb, "ROOT::Math::inc_beta()");
- legend->Draw();
+ if ( showGraphics )
+ {
- c1->Show();
+ TCanvas* c1 = new TCanvas("c1", "Two Graphs", 600, 400);
+ TH2F* hpx = new TH2F("hpx", "Two Graphs(hpx)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, 0, 5);
+ hpx->SetStats(kFALSE);
+ hpx->Draw();
+
+ int color = 2;
+
+ TGraph *gb, *gmb;
+ gb = drawPoints(&x[0], &yb[0], color++);
+ gmb = drawPoints(&x[0], &ymb[0], color++, 7);
+
+ TLegend* legend = new TLegend(0.61,0.72,0.86,0.86);
+ legend->AddEntry(gb, "TMath::BetaIncomplete()");
+ legend->AddEntry(gmb, "ROOT::Math::inc_beta()");
+ legend->Draw();
+
+ c1->Show();
+ }
cout << "Test Done!" << endl;
- return;
+ return status;
}
int main(int argc, char **argv)
{
- TApplication theApp("App",&argc,argv);
- testSpecFuncBetaI();
- theApp.Run();
+ if ( argc > 1 && argc != 2 )
+ {
+ cerr << "Usage: " << argv[0] << " [-ng]\n";
+ cerr << " where:\n";
+ cerr << " -ng : no graphics mode";
+ cerr << endl;
+ exit(1);
+ }
+
+ if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
+ {
+ showGraphics = false;
+ }
+
+ TApplication* theApp = 0;
+ if ( showGraphics )
+ theApp = new TApplication("App",&argc,argv);
+
+ int status = testSpecFuncBetaI();
+
+ if ( showGraphics )
+ {
+ theApp->Run();
+ delete theApp;
+ theApp = 0;
+ }
- return 0;
+ return status;
}
diff --git a/math/mathcore/test/testSpecFuncErf.cxx b/math/mathcore/test/testSpecFuncErf.cxx
index 604ab27a1dc..e2f852d43de 100644
--- a/math/mathcore/test/testSpecFuncErf.cxx
+++ b/math/mathcore/test/testSpecFuncErf.cxx
@@ -2,6 +2,8 @@
#include <fstream>
#include <vector>
+#include <cmath>
+
#include <TMath.h>
#include <Math/SpecFunc.h>
@@ -14,11 +16,13 @@
#include <TGraph.h>
#include <TLegend.h>
+const double ERRORLIMIT = 1E-8;
const double MIN = -2.5;
const double MAX = +2.5;
const double INCREMENT = 0.01;
-const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT);
-inline int arrayindex(double i) { return ARRAYSIZE - (int) ( (MAX - i) / INCREMENT ) -1 ; };
+const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT) + 1;
+
+bool showGraphics = true;
using namespace std;
@@ -33,7 +37,7 @@ TGraph* drawPoints(Double_t x[], Double_t y[], int color, int style = 1)
return g;
}
-void testSpecFuncErf()
+int testSpecFuncErf()
{
vector<Double_t> x( ARRAYSIZE );
vector<Double_t> yerf( ARRAYSIZE );
@@ -44,64 +48,134 @@ void testSpecFuncErf()
vector<Double_t> yierfc( ARRAYSIZE );
// vector<Double_t> yndtri( ARRAYSIZE );
- ofstream outputFile ("values.txt");
+ int status = 0;
+
+// ofstream outputFile ("values.txt");
+ unsigned int index = 0;
for ( double i = MIN; i < MAX; i += INCREMENT )
{
- outputFile << "i:"; outputFile.width(5); outputFile << i
- << " index: "; outputFile.width(5); outputFile << arrayindex(i)
- << " TMath::Erf(x): "; outputFile.width(10); outputFile << TMath::Erf(i)
- << " ROOT::Math::erf(x): "; outputFile.width(10); outputFile << ROOT::Math::erf(i)
- << " TMath::Erfc(x): "; outputFile.width(10); outputFile << TMath::Erfc(i)
- << " ROOT::Math::erfc(x): "; outputFile.width(10); outputFile << ROOT::Math::erfc(i)
- << " TMath::ErfInverse(x): "; outputFile.width(10); outputFile << TMath::ErfInverse(i)
- << " TMath::ErfcInverse(x): "; outputFile.width(10); outputFile << TMath::ErfcInverse(i)
-// << " ROOT::Math::Cephes::ndtri(x): "; outputFile.width(10); outputFile << ROOT::Math::Cephes::ndtri(i)
- << endl;
-
- x[arrayindex(i)] = i;
- yerf[arrayindex(i)] = TMath::Erf(i);
- ymerf[arrayindex(i)] = ROOT::Math::erf(i);
- yerfc[arrayindex(i)] = TMath::Erfc(i);
- ymerfc[arrayindex(i)] = ROOT::Math::erfc(i);
- yierf[arrayindex(i)] = TMath::ErfInverse(i);
- yierfc[arrayindex(i)] = TMath::ErfcInverse(i);
-// yndtri[arrayindex(i)] = ROOT::Math::Cephes::ndtri(i);
+// outputFile << "i:"; outputFile.width(5); outputFile << i
+// << " index: "; outputFile.width(5); outputFile << index
+// << " TMath::Erf(x): "; outputFile.width(10); outputFile << TMath::Erf(i)
+// << " ROOT::Math::erf(x): "; outputFile.width(10); outputFile << ROOT::Math::erf(i)
+// << " TMath::Erfc(x): "; outputFile.width(10); outputFile << TMath::Erfc(i)
+// << " ROOT::Math::erfc(x): "; outputFile.width(10); outputFile << ROOT::Math::erfc(i)
+// << " TMath::ErfInverse(x): "; outputFile.width(10); outputFile << TMath::ErfInverse(i)
+// << " TMath::ErfcInverse(x): "; outputFile.width(10); outputFile << TMath::ErfcInverse(i)
+// // << " ROOT::Math::Cephes::ndtri(x): "; outputFile.width(10); outputFile << ROOT::Math::Cephes::ndtri(i)
+// << endl;
+
+ x[index] = i;
+
+ yerf[index] = TMath::Erf(i);
+ ymerf[index] = ROOT::Math::erf(i);
+ if ( std::fabs( yerf[index] - ymerf[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yerf[index] " << yerf[index]
+ << " ymerf[index] " << ymerf[index]
+ << " " << std::fabs( yerf[index] - ymerf[index] )
+ << endl;
+ status += 1;
+ }
+
+ yerfc[index] = TMath::Erfc(i);
+ ymerfc[index] = ROOT::Math::erfc(i);
+ if ( std::fabs( yerfc[index] - ymerfc[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yerfc[index] " << yerfc[index]
+ << " ymerfc[index] " << ymerfc[index]
+ << " " << std::fabs( yerfc[index] - ymerfc[index] )
+ << endl;
+ status += 1;
+ }
+
+ yierf[index] = TMath::ErfInverse(yerf[index]);
+ if ( std::fabs( yierf[index] - i ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yierf[index] " << yierf[index]
+ << " " << std::fabs( yierf[index] - i )
+ << endl;
+ status += 1;
+ }
+
+ yierfc[index] = TMath::ErfcInverse(yerfc[index]);
+ if ( std::fabs( yierfc[index] - i ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yierfc[index] " << yierfc[index]
+ << " " << std::fabs( yierfc[index] - i )
+ << endl;
+ status += 1;
+ }
+
+// yndtri[index] = ROOT::Math::Cephes::ndtri(i);
+
+ index += 1;
}
- TCanvas* c1 = new TCanvas("c1", "Two Graphs", 600, 400);
- TH2F* hpx = new TH2F("hpx", "Two Graphs(hpx)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, -1,2);
- hpx->SetStats(kFALSE);
- hpx->Draw();
+ if ( showGraphics )
+ {
- TGraph* gerf = drawPoints(&x[0], &yerf[0], 14);
- TGraph* gmerf = drawPoints(&x[0], &ymerf[0], 5, 7);
- TGraph* gerfc = drawPoints(&x[0], &yerfc[0], 2);
- TGraph* gmerfc = drawPoints(&x[0], &ymerfc[0], 3, 7);
+ TCanvas* c1 = new TCanvas("c1", "Two Graphs", 600, 400);
+ TH2F* hpx = new TH2F("hpx", "Two Graphs(hpx)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, -1,2);
+ hpx->SetStats(kFALSE);
+ hpx->Draw();
+
+ TGraph* gerf = drawPoints(&x[0], &yerf[0], 14);
+ TGraph* gmerf = drawPoints(&x[0], &ymerf[0], 5, 7);
+ TGraph* gerfc = drawPoints(&x[0], &yerfc[0], 2);
+ TGraph* gmerfc = drawPoints(&x[0], &ymerfc[0], 3, 7);
// drawPoints(&x[0], &yierf[0], 21);
// drawPoints(&x[0], &yierfc[0], 28);
// drawPoints(&x[0], &yndtri[0], 9);
- TLegend* legend = new TLegend(0.61,0.62,0.86,0.86);
- legend->AddEntry(gerf, "TMath::Erf()");
- legend->AddEntry(gmerf, "ROOT:Math::erf()");
- legend->AddEntry(gerfc, "TMath::Erfc()");
- legend->AddEntry(gmerfc, "ROOT::Math::erfInverse()");
- legend->Draw();
-
- c1->Show();
-
+ TLegend* legend = new TLegend(0.61,0.62,0.86,0.86);
+ legend->AddEntry(gerf, "TMath::Erf()");
+ legend->AddEntry(gmerf, "ROOT:Math::erf()");
+ legend->AddEntry(gerfc, "TMath::Erfc()");
+ legend->AddEntry(gmerfc, "ROOT::Math::erfInverse()");
+ legend->Draw();
+
+ c1->Show();
+ }
+
cout << "Test Done!" << endl;
- return;
+ return status;
}
-
int main(int argc, char **argv)
{
- TApplication theApp("App",&argc,argv);
- testSpecFuncErf();
- theApp.Run();
+ if ( argc > 1 && argc != 2 )
+ {
+ cerr << "Usage: " << argv[0] << " [-ng]\n";
+ cerr << " where:\n";
+ cerr << " -ng : no graphics mode";
+ cerr << endl;
+ exit(1);
+ }
+
+ if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
+ {
+ showGraphics = false;
+ }
+
+ TApplication* theApp = 0;
+ if ( showGraphics )
+ theApp = new TApplication("App",&argc,argv);
+
+ int status = testSpecFuncErf();
+
+ if ( showGraphics )
+ {
+ theApp->Run();
+ delete theApp;
+ theApp = 0;
+ }
- return 0;
+ return status;
}
diff --git a/math/mathcore/test/testSpecFuncGamma.cxx b/math/mathcore/test/testSpecFuncGamma.cxx
index 18de74b5273..c5329dafcd3 100644
--- a/math/mathcore/test/testSpecFuncGamma.cxx
+++ b/math/mathcore/test/testSpecFuncGamma.cxx
@@ -2,6 +2,8 @@
#include <fstream>
#include <vector>
+#include <cmath>
+
#include <TMath.h>
#include <Math/SpecFuncMathCore.h>
@@ -12,11 +14,13 @@
#include <TGraph.h>
#include <TLegend.h>
+const double ERRORLIMIT = 1E-8;
const double MIN = -2.5;
const double MAX = +2.5;
const double INCREMENT = 0.01;
-const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT);
-inline int arrayindex(double i) { return ARRAYSIZE - (int) ( (MAX - i) / INCREMENT ) -1 ; };
+const int ARRAYSIZE = (int) (( MAX - MIN ) / INCREMENT) + 1;
+
+bool showGraphics = true;
using namespace std;
@@ -31,7 +35,7 @@ TGraph* drawPoints(Double_t x[], Double_t y[], int color, int style = 1)
return g;
}
-void testSpecFuncGamma()
+int testSpecFuncGamma()
{
vector<Double_t> x( ARRAYSIZE );
vector<Double_t> yg( ARRAYSIZE );
@@ -43,12 +47,15 @@ void testSpecFuncGamma()
Double_t a = 0.56;
+ int status = 0;
+
//ofstream cout ("values.txt");
+ unsigned int index = 0;
for ( double i = MIN; i < MAX; i += INCREMENT )
{
// cout << "i:"; cout.width(5); cout << i
-// << " index: "; cout.width(5); cout << arrayindex(i)
+// << " index: "; cout.width(5); cout << index
// << " TMath::Gamma(x): "; cout.width(10); cout << TMath::Gamma(i)
// << " ROOT::Math::tgamma(x): "; cout.width(10); cout << ROOT::Math::tgamma(i)
// << " TMath::Gamma(a, x): "; cout.width(10); cout << TMath::Gamma(a, i)
@@ -57,49 +64,108 @@ void testSpecFuncGamma()
// << " ROOT::Math::lgamma(x): "; cout.width(10); cout << ROOT::Math::lgamma(i)
// << endl;
- x[arrayindex(i)] = i;
- yg[arrayindex(i)] = TMath::Gamma(i);
- ymtg[arrayindex(i)] = ROOT::Math::tgamma(i);
- yga[arrayindex(i)] = TMath::Gamma(a, i);
- ymga[arrayindex(i)] = ROOT::Math::inc_gamma(a, i);
- ylng[arrayindex(i)] = TMath::LnGamma(i);
- ymlng[arrayindex(i)] = ROOT::Math::lgamma(i);
+ x[index] = i;
+ yg[index] = TMath::Gamma(i);
+ ymtg[index] = ROOT::Math::tgamma(i);
+ // take the infinity values out of the error checking!
+ if ( std::fabs(yg[index]) < 1E+12 && std::fabs( yg[index] - ymtg[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yg[index] " << yg[index]
+ << " ymtg[index] " << ymtg[index]
+ << " " << std::fabs( yg[index] - ymtg[index] )
+ << endl;
+ status += 1;
+ }
+
+ yga[index] = TMath::Gamma(a, i);
+ ymga[index] = ROOT::Math::inc_gamma(a, i);
+ if ( std::fabs( yga[index] - ymga[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " yga[index] " << yga[index]
+ << " ymga[index] " << ymga[index]
+ << " " << std::fabs( yga[index] - ymga[index] )
+ << endl;
+ status += 1;
+ }
+
+ ylng[index] = TMath::LnGamma(i);
+ ymlng[index] = ROOT::Math::lgamma(i);
+ if ( std::fabs( ylng[index] - ymlng[index] ) > ERRORLIMIT )
+ {
+ cout << "i " << i
+ << " ylng[index] " << ylng[index]
+ << " ymlng[index] " << ymlng[index]
+ << " " << std::fabs( ylng[index] - ymlng[index] )
+ << endl;
+ status += 1;
+ }
+
+ index += 1;
}
- TCanvas* c1 = new TCanvas("c1", "Two Graphs", 600, 400);
- TH2F* hpx = new TH2F("hpx", "Two Graphs(hpx)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, -1,5);
- hpx->SetStats(kFALSE);
- hpx->Draw();
-
- TGraph* gg = drawPoints(&x[0], &yg[0], 1);
- TGraph* gmtg = drawPoints(&x[0], &ymtg[0], 2, 7);
- TGraph* gga = drawPoints(&x[0], &yga[0], 3);
- TGraph* gmga = drawPoints(&x[0], &ymga[0], 4, 7);
- TGraph* glng = drawPoints(&x[0], &ylng[0], 5);
- TGraph* gmlng = drawPoints(&x[0], &ymlng[0], 6, 7);
-
- TLegend* legend = new TLegend(0.61,0.52,0.86,0.86);
- legend->AddEntry(gg, "TMath::Gamma()");
- legend->AddEntry(gmtg, "ROOT::Math::tgamma()");
- legend->AddEntry(gga, "TMath::GammaI()");
- legend->AddEntry(gmga, "ROOT::Math::inc_gamma()");
- legend->AddEntry(glng, "TMath::LnGamma()");
- legend->AddEntry(gmlng, "ROOT::Math::lgamma()");
- legend->Draw();
-
- c1->Show();
+ if ( showGraphics )
+ {
+
+ TCanvas* c1 = new TCanvas("c1", "Two Graphs", 600, 400);
+ TH2F* hpx = new TH2F("hpx", "Two Graphs(hpx)", ARRAYSIZE, MIN, MAX, ARRAYSIZE, -1,5);
+ hpx->SetStats(kFALSE);
+ hpx->Draw();
+
+ TGraph* gg = drawPoints(&x[0], &yg[0], 1);
+ TGraph* gmtg = drawPoints(&x[0], &ymtg[0], 2, 7);
+ TGraph* gga = drawPoints(&x[0], &yga[0], 3);
+ TGraph* gmga = drawPoints(&x[0], &ymga[0], 4, 7);
+ TGraph* glng = drawPoints(&x[0], &ylng[0], 5);
+ TGraph* gmlng = drawPoints(&x[0], &ymlng[0], 6, 7);
+
+ TLegend* legend = new TLegend(0.61,0.52,0.86,0.86);
+ legend->AddEntry(gg, "TMath::Gamma()");
+ legend->AddEntry(gmtg, "ROOT::Math::tgamma()");
+ legend->AddEntry(gga, "TMath::GammaI()");
+ legend->AddEntry(gmga, "ROOT::Math::inc_gamma()");
+ legend->AddEntry(glng, "TMath::LnGamma()");
+ legend->AddEntry(gmlng, "ROOT::Math::lgamma()");
+ legend->Draw();
+
+ c1->Show();
+ }
cout << "Test Done!" << endl;
- return;
+ return status;
}
int main(int argc, char **argv)
{
- TApplication theApp("App",&argc,argv);
- testSpecFuncGamma();
- theApp.Run();
+ if ( argc > 1 && argc != 2 )
+ {
+ cerr << "Usage: " << argv[0] << " [-ng]\n";
+ cerr << " where:\n";
+ cerr << " -ng : no graphics mode";
+ cerr << endl;
+ exit(1);
+ }
+
+ if ( argc == 2 && strcmp( argv[1], "-ng") == 0 )
+ {
+ showGraphics = false;
+ }
+
+ TApplication* theApp = 0;
+ if ( showGraphics )
+ theApp = new TApplication("App",&argc,argv);
+
+ int status = testSpecFuncGamma();
+
+ if ( showGraphics )
+ {
+ theApp->Run();
+ delete theApp;
+ theApp = 0;
+ }
- return 0;
+ return status;
}
diff --git a/math/mathcore/test/testTMath.cxx b/math/mathcore/test/testTMath.cxx
index ae22d86fe38..c1601479f48 100644
--- a/math/mathcore/test/testTMath.cxx
+++ b/math/mathcore/test/testTMath.cxx
@@ -1,34 +1,42 @@
#include <iostream>
+#include <vector>
+#include <string>
+#include <cstring>
#include <TMath.h>
using namespace std;
using namespace TMath;
+bool showVector = true;
+
+template <typename T>
void testNormCross()
{
-//Float_t NormCross(const Float_t v1[3],const Float_t v2[3],Float_t out[3])
+ //Float_t NormCross(const Float_t v1[3],const Float_t v2[3],Float_t out[3])
- Float_t fv1[] = {1,2,3};
- Float_t fv2[] = {4,5,6};
- Float_t fout[3];
+ T fv1[] = {1,2,3};
+ T fv2[] = {4,5,6};
+ T fout[3];
- TMath::NormCross(fv1,fv2,fout);
- cout << "NormCross(const Float_t v1[3],const Float_t v2[3],Float_t out[3]): out = ["
- << fout[0] << ", " << fout[1] << ", " << fout[2] << "]"
- << endl;
+ std::string type;
+ if ( strcmp( typeid(T).name(), "f" ) == 0 )
+ type = "Float_t";
+ else if ( strcmp( typeid(T).name(), "d" ) == 0 )
+ type = "Double_t";
+ else
+ type = typeid(T).name();
- Double_t dv1[] = {1,2,3};
- Double_t dv2[] = {4,5,6};
- Double_t dout[3];
+ TMath::NormCross(fv1,fv2,fout);
- TMath::NormCross(dv1,dv2,dout);
- cout << "NormCross(const Double_t v1[3],const Double_t v2[3],Double_t out[3]): out = ["
- << dout[0] << ", " << dout[1] << ", " << dout[2] << "]\n"
- << endl;
+ cout << "NormCross(const " << type << " v1[3],const "
+ << type << " v2[3]," << type << " out[3]): out = ["
+ << fout[0] << ", " << fout[1] << ", " << fout[2] << "]"
+ << endl;
}
+template <typename T>
void testArrayFunctions()
{
const Long64_t n = 10;
@@ -36,13 +44,16 @@ void testArrayFunctions()
Int_t index[n];
Long64_t is;
-// Short_t MinElement(Long64_t n, const Short_t *a)
- Short_t sa[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
+ T sa[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
- cout << "Vector a[] = {" << sa[0];
- for ( Int_t i = 1; i < n; ++i )
- cout << ", " << sa[i];
- cout << "}\n" << endl;
+ if ( showVector )
+ {
+ cout << "Vector a[] = {" << sa[0];
+ for ( Int_t i = 1; i < n; ++i )
+ cout << ", " << sa[i];
+ cout << "}\n" << endl;
+ showVector = false;
+ }
cout << "Min: a[" << LocMin(n, sa) << "] = " << MinElement(n, sa)
<< " Max: a[" << LocMax(n, sa) << "] = " << MaxElement(n, sa)
@@ -60,165 +71,112 @@ void testArrayFunctions()
cout << "}" << endl;
sort(sa, sa+n);
- is = BinarySearch(n, sa, (Short_t) 57);
- cout << "BinarySearch(n, a, 57) = " << is << "\n" << endl;
-
-// Int_t MinElement(Long64_t n, const Int_t *a)
- Int_t ia[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
- cout << "Min: a[" << LocMin(n, ia) << "] = " << MinElement(n, ia)
- << " Max: a[" << LocMax(n, ia) << "] = " << MaxElement(n, ia)
- << " Mean: " << Mean(n, ia)
- << " GeomMean: " << GeomMean(n, ia)
- << " RMS: " << RMS(n, ia)
- << " Median: " << Median(n, ia)
- << " KOrdStat(3): " << KOrdStat(n, sa, k)
- << endl;
-
- Sort(n, ia, index, kFALSE);
- cout << "Sorted a[] = {" << ia[index[0]];
- for ( Int_t i = 1; i < n; ++i )
- cout << ", " << ia[index[i]];
- cout << "}" << endl;
-
- sort(sa, sa+n);
- is = BinarySearch(n, sa, (Short_t) 57);
- cout << "BinarySearch(n, a, 57) = " << is << "\n" << endl;
-
-// Float_t MinElement(Long64_t n, const Float_t *a)
- Float_t fa[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
- cout << "Min: a[" << LocMin(n, fa) << "] = " << MinElement(n, fa)
- << " Max: a[" << LocMax(n, fa) << "] = " << MaxElement(n, fa)
- << " Mean: " << Mean(n, fa)
- << " GeomMean: " << GeomMean(n, fa)
- << " RMS: " << RMS(n, fa)
- << " Median: " << Median(n, fa)
- << " KOrdStat(3): " << KOrdStat(n, sa, k)
- << endl;
-
- Sort(n, fa, index, kFALSE);
- cout << "Sorted a[] = {" << fa[index[0]];
- for ( Int_t i = 1; i < n; ++i )
- cout << ", " << fa[index[i]];
- cout << "}" << endl;
-
- sort(sa, sa+n);
- is = BinarySearch(n, sa, (Short_t) 57);
- cout << "BinarySearch(n, a, 57) = " << is << "\n" << endl;
-
-// Double_t MinElement(Long64_t n, const Double_t *a)
- Double_t da[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
- cout << "Min: a[" << LocMin(n, da) << "] = " << MinElement(n, da)
- << " Max: a[" << LocMax(n, da) << "] = " << MaxElement(n, da)
- << " Mean: " << Mean(n, da)
- << " GeomMean: " << GeomMean(n, da)
- << " RMS: " << RMS(n, da)
- << " Median: " << Median(n, da)
- << " KOrdStat(3): " << KOrdStat(n, sa, k)
- << endl;
-
- Sort(n, da, index, kFALSE);
- cout << "Sorted a[] = {" << da[index[0]];
- for ( Int_t i = 1; i < n; ++i )
- cout << ", " << da[index[i]];
- cout << "}" << endl;
-
- sort(sa, sa+n);
- is = BinarySearch(n, sa, (Short_t) 57);
+ is = BinarySearch(n, sa, (T) 57);
cout << "BinarySearch(n, a, 57) = " << is << "\n" << endl;
+}
-// Long_t MinElement(Long64_t n, const Long_t *a)
- Long_t la[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
- cout << "Min: a[" << LocMin(n, la) << "] = " << MinElement(n, la)
- << " Max: a[" << LocMax(n, la) << "] = " << MaxElement(n, la)
- << " Mean: " << Mean(n, la)
- << " GeomMean: " << GeomMean(n, la)
- << " RMS: " << RMS(n, la)
- << " Median: " << Median(n, la)
- << " KOrdStat(3): " << KOrdStat(n, sa, k)
- << endl;
-
- Sort(n, la, index, kFALSE);
- cout << "Sorted a[] = {" << la[index[0]];
- for ( Int_t i = 1; i < n; ++i )
- cout << ", " << la[index[i]];
- cout << "}" << endl;
-
- sort(sa, sa+n);
- is = BinarySearch(n, sa, (Short_t) 57);
- cout << "BinarySearch(n, a, 57) = " << is << "\n" << endl;
+template <typename T>
+void testIteratorFunctions()
+{
+ const Long64_t n = 10;
+ vector<Int_t> index(n);
+ Long64_t is;
-// Long64_t MinElement(Long64_t n, const Long64_t *a)
- Long64_t l64a[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
- cout << "Min: a[" << LocMin(n, l64a) << "] = " << MinElement(n, l64a)
- << " Max: a[" << LocMax(n, l64a) << "] = " << MaxElement(n, l64a)
- << " Mean: " << Mean(n, l64a)
- << " GeomMean: " << GeomMean(n, l64a)
- << " RMS: " << RMS(n, l64a)
- << " Median: " << Median(n, l64a)
- << " KOrdStat(3): " << KOrdStat(n, sa, k)
+ T tsa[n] = { 2, 55 ,23, 57, -9, 24, 6, 82, -4, 10};
+ vector<T> sa(n);
+ for ( int i = 0; i < n; ++i ) sa[i] = tsa[i];
+
+ if ( showVector )
+ {
+ cout << "\nVector a[] = {" << sa[0];
+ for ( Int_t i = 1; i < n; ++i )
+ cout << ", " << sa[i];
+ cout << "}\n" << endl;
+ showVector = false;
+ }
+
+ cout << "Min: " << *LocMin(sa.begin(), sa.end())
+ << " Max: " << *LocMax(sa.begin(), sa.end())
+ << " Mean: " << Mean(sa.begin(), sa.end())
+ << " GeomMean: " << GeomMean(sa.begin(), sa.end())
+ << " RMS: " << RMS(sa.begin(), sa.end())
<< endl;
- Sort(n, l64a, index, kFALSE);
- cout << "Sorted a[] = {" << l64a[index[0]];
+ TMath::SortItr(sa.begin(), sa.end(), index.begin(), kFALSE);
+ cout << "Sorted a[] = {" << sa[ index[0] ];
for ( Int_t i = 1; i < n; ++i )
- cout << ", " << l64a[index[i]];
+ cout << ", " << sa[ index[i] ];
cout << "}" << endl;
- sort(sa, sa+n);
- is = BinarySearch(n, sa, (Short_t) 57);
+ sort(&sa[0], &sa[n]);
+ is = BinarySearch(n, &sa[0], (T) 57);
cout << "BinarySearch(n, a, 57) = " << is << "\n" << endl;
-
}
-void testPlane()
+template <typename T>
+void testPoints(T x, T y)
{
const Int_t n = 4;
- Double_t dx[4] = {0.0, 0.0, 2.0, 2.0};
- Double_t dy[4] = {0.0, 2.0, 2.0, 0.0};
- cout << "Point(1.3,0.5) IsInside?: " << IsInside(1.3, 0.5, n, dx, dy) << endl;
+ T dx[4] = {0, 0, 2, 2};
+ T dy[4] = {0, 2, 2, 0};
+ cout << "Point(" << x << "," << y << ") IsInside?: "
+ << IsInside( x, y, n, dx, dy) << endl;
- Float_t fx[4] = {0.0, 0.0, 2.0, 2.0};
- Float_t fy[4] = {0.0, 2.0, 2.0, 0.0};
- cout << "Point(-0.2f,1.7f) IsInside?: " << IsInside(-0.2f, 1.7f, n, fx, fy) << endl;
-
- Int_t ix[4] = {0, 0, 2, 2};
- Int_t iy[4] = {0, 2, 2, 0};
- cout << "Point(1,1) IsInside?: " << IsInside(1, 1, n, ix, iy) << endl;
+}
- Double_t dp1[3] = {0,0,0};
- Double_t dp2[3] = {1,0,0};
- Double_t dp3[3] = {0,1,0};
- Double_t dn[3];
+template <typename T>
+void testPlane()
+{
+ T dp1[3] = {0,0,0};
+ T dp2[3] = {1,0,0};
+ T dp3[3] = {0,1,0};
+ T dn[3];
Normal2Plane(dp1, dp2, dp3, dn);
cout << "Normal: ("
<< dn[0] << ", "
<< dn[1] << ", "
<< dn[2] << ")"
<< endl;
-
- Double_t fp1[3] = {0,0,0};
- Double_t fp2[3] = {1,0,0};
- Double_t fp3[3] = {0,1,0};
- Double_t fn[3];
- Normal2Plane(fp1, fp2, fp3, fn);
- cout << "Normal: ("
- << fn[0] << ", "
- << fn[1] << ", "
- << fn[2] << ")"
- << endl;
-
}
-
void testTMath()
{
cout << "Starting tests on TMath..." << endl;
- testNormCross();
- testArrayFunctions();
- testPlane();
+ cout << "\nNormCross tests: " << endl;
+
+ testNormCross<Float_t>();
+ testNormCross<Double_t>();
+
+ cout << "\nArray functions tests: " << endl;
+
+ testArrayFunctions<Short_t>();
+ testArrayFunctions<Int_t>();
+ testArrayFunctions<Float_t>();
+ testArrayFunctions<Double_t>();
+ testArrayFunctions<Long_t>();
+ testArrayFunctions<Long64_t>();
+
+ cout << "\nIterator functions tests: " << endl;
+
+ testIteratorFunctions<Short_t>();
+ testIteratorFunctions<Int_t>();
+ testIteratorFunctions<Float_t>();
+ testIteratorFunctions<Double_t>();
+ testIteratorFunctions<Long_t>();
+ testIteratorFunctions<Long64_t>();
+
+ cout << "\nPoint functions tests: " << endl;
+
+ testPoints<Double_t>(1.3, 0.5);
+ testPoints<Float_t>(-0.2, 1.7);
+ testPoints<Int_t>(1, 1);
+
+ cout << "\nPLane functions tests: " << endl;
+ testPlane<Double_t>();
+ testPlane<Float_t>();
}
int main()
--
GitLab