Logo ROOT  
Reference Guide
TMath.h
Go to the documentation of this file.
1 // @(#)root/mathcore:$Id$
2 // Authors: Rene Brun, Anna Kreshuk, Eddy Offermann, Fons Rademakers 29/07/95
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2004, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #ifndef ROOT_TMath
13 #define ROOT_TMath
14 
15 #include "TMathBase.h"
16 
17 #include "TError.h"
18 #include <algorithm>
19 #include <limits>
20 #include <cmath>
21 
22 ////////////////////////////////////////////////////////////////////////////////
23 ///
24 /// TMath
25 ///
26 /// Encapsulate most frequently used Math functions.
27 /// NB. The basic functions Min, Max, Abs and Sign are defined
28 /// in TMathBase.
29 
30 namespace TMath {
31 
32 ////////////////////////////////////////////////////////////////////////////////
33 // Fundamental constants
34 
35 ////////////////////////////////////////////////////////////////////////////////
36 /// \f[ \pi\f]
37 constexpr Double_t Pi()
38 {
39  return 3.14159265358979323846;
40 }
41 
42 ////////////////////////////////////////////////////////////////////////////////
43 /// \f[ 2\pi\f]
44 constexpr Double_t TwoPi()
45 {
46  return 2.0 * Pi();
47 }
48 
49 ////////////////////////////////////////////////////////////////////////////////
50 /// \f[ \frac{\pi}{2} \f]
51 constexpr Double_t PiOver2()
52 {
53  return Pi() / 2.0;
54 }
55 
56 ////////////////////////////////////////////////////////////////////////////////
57 /// \f[ \frac{\pi}{4} \f]
58 constexpr Double_t PiOver4()
59 {
60  return Pi() / 4.0;
61 }
62 
63 ////////////////////////////////////////////////////////////////////////////////
64 /// \f$ \frac{1.}{\pi}\f$
65 constexpr Double_t InvPi()
66 {
67  return 1.0 / Pi();
68 }
69 
70 ////////////////////////////////////////////////////////////////////////////////
71 /// Conversion from radian to degree:
72 /// \f[ \frac{180}{\pi} \f]
73 constexpr Double_t RadToDeg()
74 {
75  return 180.0 / Pi();
76 }
77 
78 ////////////////////////////////////////////////////////////////////////////////
79 /// Conversion from degree to radian:
80 /// \f[ \frac{\pi}{180} \f]
81 constexpr Double_t DegToRad()
82 {
83  return Pi() / 180.0;
84 }
85 
86 ////////////////////////////////////////////////////////////////////////////////
87 /// \f[ \sqrt{2} \f]
88 constexpr Double_t Sqrt2()
89 {
90  return 1.4142135623730950488016887242097;
91 }
92 
93 ////////////////////////////////////////////////////////////////////////////////
94 /// Base of natural log:
95 /// \f[ e \f]
96 constexpr Double_t E()
97 {
98  return 2.71828182845904523536;
99 }
100 
101 ////////////////////////////////////////////////////////////////////////////////
102 /// Natural log of 10 (to convert log to ln)
103 constexpr Double_t Ln10()
104 {
105  return 2.30258509299404568402;
106 }
107 
108 ////////////////////////////////////////////////////////////////////////////////
109 /// Base-10 log of e (to convert ln to log)
110 constexpr Double_t LogE()
111 {
112  return 0.43429448190325182765;
113 }
114 
115 ////////////////////////////////////////////////////////////////////////////////
116 /// Velocity of light in \f$ m s^{-1} \f$
117 constexpr Double_t C()
118 {
119  return 2.99792458e8;
120 }
121 
122 ////////////////////////////////////////////////////////////////////////////////
123 /// \f$ cm s^{-1} \f$
124 constexpr Double_t Ccgs()
125 {
126  return 100.0 * C();
127 }
128 
129 ////////////////////////////////////////////////////////////////////////////////
130 /// Speed of light uncertainty.
132 {
133  return 0.0;
134 }
135 
136 ////////////////////////////////////////////////////////////////////////////////
137 /// Gravitational constant in: \f$ m^{3} kg^{-1} s^{-2} \f$
138 constexpr Double_t G()
139 {
140  // use 2018 value from NIST (https://physics.nist.gov/cgi-bin/cuu/Value?bg|search_for=G)
141  return 6.67430e-11;
142 }
143 
144 ////////////////////////////////////////////////////////////////////////////////
145 /// \f$ cm^{3} g^{-1} s^{-2} \f$
146 constexpr Double_t Gcgs()
147 {
148  return G() * 1000.0;
149 }
150 
151 ////////////////////////////////////////////////////////////////////////////////
152 /// Gravitational constant uncertainty.
154 {
155  // use 2018 value from NIST
156  return 0.00015e-11;
157 }
158 
159 ////////////////////////////////////////////////////////////////////////////////
160 /// \f$ \frac{G}{\hbar C} \f$ in \f$ (GeV/c^{2})^{-2} \f$
161 constexpr Double_t GhbarC()
162 {
163  // use new value from NIST (2018)
164  return 6.70883e-39;
165 }
166 
167 ////////////////////////////////////////////////////////////////////////////////
168 /// \f$ \frac{G}{\hbar C} \f$ uncertainty.
170 {
171  // use new value from NIST (2018)
172  return 0.00015e-39;
173 }
174 
175 ////////////////////////////////////////////////////////////////////////////////
176 /// Standard acceleration of gravity in \f$ m s^{-2} \f$
177 constexpr Double_t Gn()
178 {
179  return 9.80665;
180 }
181 
182 ////////////////////////////////////////////////////////////////////////////////
183 /// Standard acceleration of gravity uncertainty.
185 {
186  return 0.0;
187 }
188 
189 ////////////////////////////////////////////////////////////////////////////////
190 /// Planck's constant in \f$ J s \f$
191 /// \f[ h \f]
192 constexpr Double_t H()
193 {
194  return 6.62607015e-34;
195 }
196 
197 ////////////////////////////////////////////////////////////////////////////////
198 /// \f$ erg s \f$
199 constexpr Double_t Hcgs()
200 {
201  return 1.0e7 * H();
202 }
203 
204 ////////////////////////////////////////////////////////////////////////////////
205 /// Planck's constant uncertainty.
207 {
208  // Planck constant is exact according to 2019 redefinition
209  // (https://en.wikipedia.org/wiki/2019_redefinition_of_the_SI_base_units)
210  return 0.0;
211 }
212 
213 ////////////////////////////////////////////////////////////////////////////////
214 /// \f$ \hbar \f$ in \f$ J s \f$
215 /// \f[ \hbar = \frac{h}{2\pi} \f]
216 constexpr Double_t Hbar()
217 {
218  return 1.054571817e-34;
219 }
220 
221 ////////////////////////////////////////////////////////////////////////////////
222 /// \f$ erg s \f$
223 constexpr Double_t Hbarcgs()
224 {
225  return 1.0e7 * Hbar();
226 }
227 
228 ////////////////////////////////////////////////////////////////////////////////
229 /// \f$ \hbar \f$ uncertainty.
231 {
232  // hbar is an exact constant
233  return 0.0;
234 }
235 
236 ////////////////////////////////////////////////////////////////////////////////
237 /// \f$ hc \f$ in \f$ J m \f$
238 constexpr Double_t HC()
239 {
240  return H() * C();
241 }
242 
243 ////////////////////////////////////////////////////////////////////////////////
244 /// \f$ erg cm \f$
245 constexpr Double_t HCcgs()
246 {
247  return Hcgs() * Ccgs();
248 }
249 
250 ////////////////////////////////////////////////////////////////////////////////
251 /// Boltzmann's constant in \f$ J K^{-1} \f$
252 /// \f[ k \f]
253 constexpr Double_t K()
254 {
255  return 1.380649e-23;
256 }
257 
258 ////////////////////////////////////////////////////////////////////////////////
259 /// \f$ erg K^{-1} \f$
260 constexpr Double_t Kcgs()
261 {
262  return 1.0e7 * K();
263 }
264 
265 ////////////////////////////////////////////////////////////////////////////////
266 /// Boltzmann's constant uncertainty.
268 {
269  // constant is exact according to 2019 redefinition
270  // (https://en.wikipedia.org/wiki/2019_redefinition_of_the_SI_base_units)
271  return 0.0;
272 }
273 
274 ////////////////////////////////////////////////////////////////////////////////
275 /// Stefan-Boltzmann constant in \f$ W m^{-2} K^{-4}\f$
276 /// \f[ \sigma \f]
277 constexpr Double_t Sigma()
278 {
279  return 5.670373e-8;
280 }
281 
282 ////////////////////////////////////////////////////////////////////////////////
283 /// Stefan-Boltzmann constant uncertainty.
285 {
286  return 0.0;
287 }
288 
289 ////////////////////////////////////////////////////////////////////////////////
290 /// Avogadro constant (Avogadro's Number) in \f$ mol^{-1} \f$
291 constexpr Double_t Na()
292 {
293  return 6.02214076e+23;
294 }
295 
296 ////////////////////////////////////////////////////////////////////////////////
297 /// Avogadro constant (Avogadro's Number) uncertainty.
299 {
300  // constant is exact according to 2019 redefinition
301  // (https://en.wikipedia.org/wiki/2019_redefinition_of_the_SI_base_units)
302  return 0.0;
303 }
304 
305 ////////////////////////////////////////////////////////////////////////////////
306 /// [Universal gas constant](http://scienceworld.wolfram.com/physics/UniversalGasConstant.html)
307 /// (\f$ Na K \f$) in \f$ J K^{-1} mol^{-1} \f$
308 //
309 constexpr Double_t R()
310 {
311  return K() * Na();
312 }
313 
314 ////////////////////////////////////////////////////////////////////////////////
315 /// Universal gas constant uncertainty.
317 {
318  return R() * ((KUncertainty() / K()) + (NaUncertainty() / Na()));
319 }
320 
321 ////////////////////////////////////////////////////////////////////////////////
322 /// [Molecular weight of dry air 1976 US Standard Atmosphere](http://atmos.nmsu.edu/jsdap/encyclopediawork.html)
323 /// in \f$ kg kmol^{-1} \f$ or \f$ gm mol^{-1} \f$
324 constexpr Double_t MWair()
325 {
326  return 28.9644;
327 }
328 
329 ////////////////////////////////////////////////////////////////////////////////
330 /// [Dry Air Gas Constant (R / MWair)](http://atmos.nmsu.edu/education_and_outreach/encyclopedia/gas_constant.htm)
331 /// in \f$ J kg^{-1} K^{-1} \f$
332 constexpr Double_t Rgair()
333 {
334  return (1000.0 * R()) / MWair();
335 }
336 
337 ////////////////////////////////////////////////////////////////////////////////
338 /// Euler-Mascheroni Constant.
339 constexpr Double_t EulerGamma()
340 {
341  return 0.577215664901532860606512090082402431042;
342 }
343 
344 ////////////////////////////////////////////////////////////////////////////////
345 /// Elementary charge in \f$ C \f$ .
346 constexpr Double_t Qe()
347 {
348  return 1.602176634e-19;
349 }
350 
351 ////////////////////////////////////////////////////////////////////////////////
352 /// Elementary charge uncertainty.
354 {
355  // constant is exact according to 2019 redefinition
356  // (https://en.wikipedia.org/wiki/2019_redefinition_of_the_SI_base_units)
357  return 0.0;
358 }
359 
360 ////////////////////////////////////////////////////////////////////////////////
361 // Mathematical Functions
362 
363 ////////////////////////////////////////////////////////////////////////////////
364 // Trigonometrical Functions
365 
366 inline Double_t Sin(Double_t);
367 inline Double_t Cos(Double_t);
368 inline Double_t Tan(Double_t);
369 inline Double_t SinH(Double_t);
370 inline Double_t CosH(Double_t);
371 inline Double_t TanH(Double_t);
372 inline Double_t ASin(Double_t);
373 inline Double_t ACos(Double_t);
374 inline Double_t ATan(Double_t);
375 inline Double_t ATan2(Double_t y, Double_t x);
380 
381 ////////////////////////////////////////////////////////////////////////////////
382 // Elementary Functions
383 
384 inline Double_t Ceil(Double_t x);
385 inline Int_t CeilNint(Double_t x);
386 inline Double_t Floor(Double_t x);
387 inline Int_t FloorNint(Double_t x);
388 template <typename T>
389 inline Int_t Nint(T x);
390 
391 inline Double_t Sq(Double_t x);
392 inline Double_t Sqrt(Double_t x);
393 inline Double_t Exp(Double_t x);
394 inline Double_t Ldexp(Double_t x, Int_t exp);
399 inline Double_t Power(Double_t x, Double_t y);
400 inline Double_t Power(Double_t x, Int_t y);
401 inline Double_t Log(Double_t x);
403 inline Double_t Log10(Double_t x);
404 inline Int_t Finite(Double_t x);
405 inline Int_t Finite(Float_t x);
406 inline Bool_t IsNaN(Double_t x);
407 inline Bool_t IsNaN(Float_t x);
408 
409 inline Double_t QuietNaN();
410 inline Double_t SignalingNaN();
411 inline Double_t Infinity();
412 
413 template <typename T>
414 struct Limits {
415  inline static T Min();
416  inline static T Max();
417  inline static T Epsilon();
418  };
419 
420  // Some integer math
421  Long_t Hypot(Long_t x, Long_t y); // sqrt(px*px + py*py)
422 
423  // Comparing floating points
425  //return kTRUE if absolute difference between af and bf is less than epsilon
426  return TMath::Abs(af-bf) < epsilon ||
427  TMath::Abs(af - bf) < Limits<Double_t>::Min(); // handle 0 < 0 case
428 
429  }
430  inline Bool_t AreEqualRel(Double_t af, Double_t bf, Double_t relPrec) {
431  //return kTRUE if relative difference between af and bf is less than relPrec
432  return TMath::Abs(af - bf) <= 0.5 * relPrec * (TMath::Abs(af) + TMath::Abs(bf)) ||
433  TMath::Abs(af - bf) < Limits<Double_t>::Min(); // handle denormals
434  }
435 
436  /////////////////////////////////////////////////////////////////////////////
437  // Array Algorithms
438 
439  // Min, Max of an array
440  template <typename T> T MinElement(Long64_t n, const T *a);
441  template <typename T> T MaxElement(Long64_t n, const T *a);
442 
443  // Locate Min, Max element number in an array
444  template <typename T> Long64_t LocMin(Long64_t n, const T *a);
445  template <typename Iterator> Iterator LocMin(Iterator first, Iterator last);
446  template <typename T> Long64_t LocMax(Long64_t n, const T *a);
447  template <typename Iterator> Iterator LocMax(Iterator first, Iterator last);
448 
449  // Hashing
450  ULong_t Hash(const void *txt, Int_t ntxt);
451  ULong_t Hash(const char *str);
452 
453  void BubbleHigh(Int_t Narr, Double_t *arr1, Int_t *arr2);
454  void BubbleLow (Int_t Narr, Double_t *arr1, Int_t *arr2);
455 
456  Bool_t Permute(Int_t n, Int_t *a); // Find permutations
457 
458  /////////////////////////////////////////////////////////////////////////////
459  // Geometrical Functions
460 
461  //Sample quantiles
462  void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob,
463  Bool_t isSorted=kTRUE, Int_t *index = 0, Int_t type=7);
464 
465  // IsInside
466  template <typename T> Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y);
467 
468  // Calculate the Cross Product of two vectors
469  template <typename T> T *Cross(const T v1[3],const T v2[3], T out[3]);
470 
471  Float_t Normalize(Float_t v[3]); // Normalize a vector
472  Double_t Normalize(Double_t v[3]); // Normalize a vector
473 
474  //Calculate the Normalized Cross Product of two vectors
475  template <typename T> inline T NormCross(const T v1[3],const T v2[3],T out[3]);
476 
477  // Calculate a normal vector of a plane
478  template <typename T> T *Normal2Plane(const T v1[3],const T v2[3],const T v3[3], T normal[3]);
479 
480  /////////////////////////////////////////////////////////////////////////////
481  // Polynomial Functions
482 
483  Bool_t RootsCubic(const Double_t coef[4],Double_t &a, Double_t &b, Double_t &c);
484 
485  /////////////////////////////////////////////////////////////////////////////
486  // Statistic Functions
487 
488  Double_t Binomial(Int_t n,Int_t k); // Calculate the binomial coefficient n over k
497  Double_t KolmogorovTest(Int_t na, const Double_t *a, Int_t nb, const Double_t *b, Option_t *option);
506  Double_t Prob(Double_t chi2,Int_t ndf);
510  Double_t Vavilov(Double_t x, Double_t kappa, Double_t beta2);
511  Double_t VavilovI(Double_t x, Double_t kappa, Double_t beta2);
513 
514  /////////////////////////////////////////////////////////////////////////////
515  // Statistics over arrays
516 
517  //Mean, Geometric Mean, Median, RMS(sigma)
518 
519  template <typename T> Double_t Mean(Long64_t n, const T *a, const Double_t *w=0);
520  template <typename Iterator> Double_t Mean(Iterator first, Iterator last);
521  template <typename Iterator, typename WeightIterator> Double_t Mean(Iterator first, Iterator last, WeightIterator wfirst);
522 
523  template <typename T> Double_t GeomMean(Long64_t n, const T *a);
524  template <typename Iterator> Double_t GeomMean(Iterator first, Iterator last);
525 
526  template <typename T> Double_t RMS(Long64_t n, const T *a, const Double_t *w=0);
527  template <typename Iterator> Double_t RMS(Iterator first, Iterator last);
528  template <typename Iterator, typename WeightIterator> Double_t RMS(Iterator first, Iterator last, WeightIterator wfirst);
529 
530  template <typename T> Double_t StdDev(Long64_t n, const T *a, const Double_t * w = 0) { return RMS<T>(n,a,w); }
531  template <typename Iterator> Double_t StdDev(Iterator first, Iterator last) { return RMS<Iterator>(first,last); }
532  template <typename Iterator, typename WeightIterator> Double_t StdDev(Iterator first, Iterator last, WeightIterator wfirst) { return RMS<Iterator,WeightIterator>(first,last,wfirst); }
533 
534  template <typename T> Double_t Median(Long64_t n, const T *a, const Double_t *w=0, Long64_t *work=0);
535 
536  //k-th order statistic
537  template <class Element, typename Size> Element KOrdStat(Size n, const Element *a, Size k, Size *work = 0);
538 
539  /////////////////////////////////////////////////////////////////////////////
540  // Special Functions
541 
547 
548  // Bessel functions
549  Double_t BesselI(Int_t n,Double_t x); /// integer order modified Bessel function I_n(x)
550  Double_t BesselK(Int_t n,Double_t x); /// integer order modified Bessel function K_n(x)
551  Double_t BesselI0(Double_t x); /// modified Bessel function I_0(x)
552  Double_t BesselK0(Double_t x); /// modified Bessel function K_0(x)
553  Double_t BesselI1(Double_t x); /// modified Bessel function I_1(x)
554  Double_t BesselK1(Double_t x); /// modified Bessel function K_1(x)
555  Double_t BesselJ0(Double_t x); /// Bessel function J0(x) for any real x
556  Double_t BesselJ1(Double_t x); /// Bessel function J1(x) for any real x
557  Double_t BesselY0(Double_t x); /// Bessel function Y0(x) for positive x
558  Double_t BesselY1(Double_t x); /// Bessel function Y1(x) for positive x
559  Double_t StruveH0(Double_t x); /// Struve functions of order 0
560  Double_t StruveH1(Double_t x); /// Struve functions of order 1
561  Double_t StruveL0(Double_t x); /// Modified Struve functions of order 0
562  Double_t StruveL1(Double_t x); /// Modified Struve functions of order 1
563 
574 }
575 
576 ////////////////////////////////////////////////////////////////////////////////
577 // Trig and other functions
578 
579 #include <float.h>
580 
581 #if defined(R__WIN32) && !defined(__CINT__)
582 # ifndef finite
583 # define finite _finite
584 # endif
585 #endif
586 #if defined(R__AIX) || defined(R__SOLARIS_CC50) || \
587  defined(R__HPUX11) || defined(R__GLIBC) || \
588  (defined(R__MACOSX) )
589 // math functions are defined inline so we have to include them here
590 # include <math.h>
591 # ifdef R__SOLARIS_CC50
592  extern "C" { int finite(double); }
593 # endif
594 // # if defined(R__GLIBC) && defined(__STRICT_ANSI__)
595 // # ifndef finite
596 // # define finite __finite
597 // # endif
598 // # ifndef isnan
599 // # define isnan __isnan
600 // # endif
601 // # endif
602 #else
603 // don't want to include complete <math.h>
604 extern "C" {
605  extern double sin(double);
606  extern double cos(double);
607  extern double tan(double);
608  extern double sinh(double);
609  extern double cosh(double);
610  extern double tanh(double);
611  extern double asin(double);
612  extern double acos(double);
613  extern double atan(double);
614  extern double atan2(double, double);
615  extern double sqrt(double);
616  extern double exp(double);
617  extern double pow(double, double);
618  extern double log(double);
619  extern double log10(double);
620 #ifndef R__WIN32
621 # if !defined(finite)
622  extern int finite(double);
623 # endif
624 # if !defined(isnan)
625  extern int isnan(double);
626 # endif
627  extern double ldexp(double, int);
628  extern double ceil(double);
629  extern double floor(double);
630 #else
631  _CRTIMP double ldexp(double, int);
632  _CRTIMP double ceil(double);
633  _CRTIMP double floor(double);
634 #endif
635 }
636 #endif
637 
638 ////////////////////////////////////////////////////////////////////////////////
640  { return sin(x); }
641 
642 ////////////////////////////////////////////////////////////////////////////////
644  { return cos(x); }
645 
646 ////////////////////////////////////////////////////////////////////////////////
648  { return tan(x); }
649 
650 ////////////////////////////////////////////////////////////////////////////////
652  { return sinh(x); }
653 
654 ////////////////////////////////////////////////////////////////////////////////
656  { return cosh(x); }
657 
658 ////////////////////////////////////////////////////////////////////////////////
660  { return tanh(x); }
661 
662 ////////////////////////////////////////////////////////////////////////////////
664  {
665  return asin(x);
666  }
667 
668 ////////////////////////////////////////////////////////////////////////////////
670  {
671  return acos(x);
672  }
673 
674 ////////////////////////////////////////////////////////////////////////////////
676  { return atan(x); }
677 
678 ////////////////////////////////////////////////////////////////////////////////
680  { if (x != 0) return atan2(y, x);
681  if (y == 0) return 0;
682  if (y > 0) return Pi()/2;
683  else return -Pi()/2;
684  }
685 
686 ////////////////////////////////////////////////////////////////////////////////
688  { return x*x; }
689 
690 ////////////////////////////////////////////////////////////////////////////////
692  { return sqrt(x); }
693 
694 ////////////////////////////////////////////////////////////////////////////////
696  { return ceil(x); }
697 
698 ////////////////////////////////////////////////////////////////////////////////
700  { return TMath::Nint(ceil(x)); }
701 
702 ////////////////////////////////////////////////////////////////////////////////
704  { return floor(x); }
705 
706 ////////////////////////////////////////////////////////////////////////////////
708  { return TMath::Nint(floor(x)); }
709 
710 ////////////////////////////////////////////////////////////////////////////////
711 /// Round to nearest integer. Rounds half integers to the nearest even integer.
712 template<typename T>
714 {
715  int i;
716  if (x >= 0) {
717  i = int(x + 0.5);
718  if ( i & 1 && x + 0.5 == T(i) ) i--;
719  } else {
720  i = int(x - 0.5);
721  if ( i & 1 && x - 0.5 == T(i) ) i++;
722  }
723  return i;
724 }
725 
726 ////////////////////////////////////////////////////////////////////////////////
728  { return exp(x); }
729 
730 ////////////////////////////////////////////////////////////////////////////////
732  { return ldexp(x, exp); }
733 
734 ////////////////////////////////////////////////////////////////////////////////
736  { return std::pow(x,y); }
737 
738 ////////////////////////////////////////////////////////////////////////////////
740  { return std::pow(x,(LongDouble_t)y); }
741 
742 ////////////////////////////////////////////////////////////////////////////////
744  { return std::pow(x,y); }
745 
746 ////////////////////////////////////////////////////////////////////////////////
748  { return pow(x, y); }
749 
750 ////////////////////////////////////////////////////////////////////////////////
752 #ifdef R__ANSISTREAM
753  return std::pow(x, y);
754 #else
755  return pow(x, (Double_t) y);
756 #endif
757 }
758 
759 ////////////////////////////////////////////////////////////////////////////////
761  { return log(x); }
762 
763 ////////////////////////////////////////////////////////////////////////////////
765  { return log10(x); }
766 
767 ////////////////////////////////////////////////////////////////////////////////
768 /// Check if it is finite with a mask in order to be consistent in presence of
769 /// fast math.
770 /// Inspired from the CMSSW FWCore/Utilities package
772 #if defined(R__FAST_MATH)
773 
774 {
775  const unsigned long long mask = 0x7FF0000000000000LL;
776  union { unsigned long long l; double d;} v;
777  v.d =x;
778  return (v.l&mask)!=mask;
779 }
780 #else
781 # if defined(R__HPUX11)
782  { return isfinite(x); }
783 # elif defined(R__MACOSX)
784 # ifdef isfinite
785  // from math.h
786  { return isfinite(x); }
787 # else
788  // from cmath
789  { return std::isfinite(x); }
790 # endif
791 # else
792  { return finite(x); }
793 # endif
794 #endif
795 
796 ////////////////////////////////////////////////////////////////////////////////
797 /// Check if it is finite with a mask in order to be consistent in presence of
798 /// fast math.
799 /// Inspired from the CMSSW FWCore/Utilities package
801 #if defined(R__FAST_MATH)
802 
803 {
804  const unsigned int mask = 0x7f800000;
805  union { unsigned int l; float d;} v;
806  v.d =x;
807  return (v.l&mask)!=mask;
808 }
809 #else
810 { return std::isfinite(x); }
811 #endif
812 
813 // This namespace provides all the routines necessary for checking if a number
814 // is a NaN also in presence of optimisations affecting the behaviour of the
815 // floating point calculations.
816 // Inspired from the CMSSW FWCore/Utilities package
817 
818 #if defined (R__FAST_MATH)
819 namespace ROOT {
820 namespace Internal {
821 namespace Math {
822 // abridged from GNU libc 2.6.1 - in detail from
823 // math/math_private.h
824 // sysdeps/ieee754/ldbl-96/math_ldbl.h
825 
826 // part of ths file:
827  /*
828  * ====================================================
829  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
830  *
831  * Developed at SunPro, a Sun Microsystems, Inc. business.
832  * Permission to use, copy, modify, and distribute this
833  * software is freely granted, provided that this notice
834  * is preserved.
835  * ====================================================
836  */
837 
838  // A union which permits us to convert between a double and two 32 bit ints.
839  typedef union {
840  Double_t value;
841  struct {
842  UInt_t lsw;
843  UInt_t msw;
844  } parts;
845  } ieee_double_shape_type;
846 
847 #define EXTRACT_WORDS(ix0,ix1,d) \
848  do { \
849  ieee_double_shape_type ew_u; \
850  ew_u.value = (d); \
851  (ix0) = ew_u.parts.msw; \
852  (ix1) = ew_u.parts.lsw; \
853  } while (0)
854 
855  inline Bool_t IsNaN(Double_t x)
856  {
857  UInt_t hx, lx;
858 
859  EXTRACT_WORDS(hx, lx, x);
860 
861  lx |= hx & 0xfffff;
862  hx &= 0x7ff00000;
863  return (hx == 0x7ff00000) && (lx != 0);
864  }
865 
866  typedef union {
867  Float_t value;
868  UInt_t word;
869  } ieee_float_shape_type;
870 
871 #define GET_FLOAT_WORD(i,d) \
872  do { \
873  ieee_float_shape_type gf_u; \
874  gf_u.value = (d); \
875  (i) = gf_u.word; \
876  } while (0)
877 
878  inline Bool_t IsNaN(Float_t x)
879  {
880  UInt_t wx;
881  GET_FLOAT_WORD (wx, x);
882  wx &= 0x7fffffff;
883  return (Bool_t)(wx > 0x7f800000);
884  }
885 } } } // end NS ROOT::Internal::Math
886 #endif // End R__FAST_MATH
887 
888 #if defined(R__FAST_MATH)
891 #else
892  inline Bool_t TMath::IsNaN(Double_t x) { return std::isnan(x); }
893  inline Bool_t TMath::IsNaN(Float_t x) { return std::isnan(x); }
894 #endif
895 
896 ////////////////////////////////////////////////////////////////////////////////
897 // Wrapper to numeric_limits
898 
899 ////////////////////////////////////////////////////////////////////////////////
900 /// Returns a quiet NaN as [defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Quiet_NaN)
902 
903  return std::numeric_limits<Double_t>::quiet_NaN();
904 }
905 
906 ////////////////////////////////////////////////////////////////////////////////
907 /// Returns a signaling NaN as defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Signaling_NaN)
909  return std::numeric_limits<Double_t>::signaling_NaN();
910 }
911 
912 ////////////////////////////////////////////////////////////////////////////////
913 /// Returns an infinity as defined by the IEEE standard
915  return std::numeric_limits<Double_t>::infinity();
916 }
917 
918 ////////////////////////////////////////////////////////////////////////////////
919 /// Returns maximum representation for type T
920 template<typename T>
922  return (std::numeric_limits<T>::min)(); //N.B. use this signature to avoid class with macro min() on Windows
923 }
924 
925 ////////////////////////////////////////////////////////////////////////////////
926 /// Returns minimum double representation
927 template<typename T>
929  return (std::numeric_limits<T>::max)(); //N.B. use this signature to avoid class with macro max() on Windows
930 }
931 
932 ////////////////////////////////////////////////////////////////////////////////
933 /// Returns minimum double representation
934 template<typename T>
937 }
938 
939 ////////////////////////////////////////////////////////////////////////////////
940 // Advanced.
941 
942 ////////////////////////////////////////////////////////////////////////////////
943 /// Calculate the Normalized Cross Product of two vectors
944 template <typename T> inline T TMath::NormCross(const T v1[3],const T v2[3],T out[3])
945 {
946  return Normalize(Cross(v1,v2,out));
947 }
948 
949 ////////////////////////////////////////////////////////////////////////////////
950 /// Return minimum of array a of length n.
951 template <typename T>
953  return *std::min_element(a,a+n);
954 }
955 
956 ////////////////////////////////////////////////////////////////////////////////
957 /// Return maximum of array a of length n.
958 template <typename T>
960  return *std::max_element(a,a+n);
961 }
962 
963 ////////////////////////////////////////////////////////////////////////////////
964 /// Return index of array with the minimum element.
965 /// If more than one element is minimum returns first found.
966 ///
967 /// Implement here since this one is found to be faster (mainly on 64 bit machines)
968 /// than stl generic implementation.
969 /// When performing the comparison, the STL implementation needs to de-reference both the array iterator
970 /// and the iterator pointing to the resulting minimum location
971 template <typename T>
973  if (n <= 0 || !a) return -1;
974  T xmin = a[0];
975  Long64_t loc = 0;
976  for (Long64_t i = 1; i < n; i++) {
977  if (xmin > a[i]) {
978  xmin = a[i];
979  loc = i;
980  }
981  }
982  return loc;
983 }
984 
985 ////////////////////////////////////////////////////////////////////////////////
986 /// Return index of array with the minimum element.
987 /// If more than one element is minimum returns first found.
988 template <typename Iterator>
989 Iterator TMath::LocMin(Iterator first, Iterator last) {
990 
991  return std::min_element(first, last);
992 }
993 
994 ////////////////////////////////////////////////////////////////////////////////
995 /// Return index of array with the maximum element.
996 /// If more than one element is maximum returns first found.
997 ///
998 /// Implement here since it is faster (see comment in LocMin function)
999 template <typename T>
1001  if (n <= 0 || !a) return -1;
1002  T xmax = a[0];
1003  Long64_t loc = 0;
1004  for (Long64_t i = 1; i < n; i++) {
1005  if (xmax < a[i]) {
1006  xmax = a[i];
1007  loc = i;
1008  }
1009  }
1010  return loc;
1011 }
1012 
1013 ////////////////////////////////////////////////////////////////////////////////
1014 /// Return index of array with the maximum element.
1015 /// If more than one element is maximum returns first found.
1016 template <typename Iterator>
1017 Iterator TMath::LocMax(Iterator first, Iterator last)
1018 {
1019 
1020  return std::max_element(first, last);
1021 }
1022 
1023 ////////////////////////////////////////////////////////////////////////////////
1024 /// Return the weighted mean of an array defined by the iterators.
1025 template <typename Iterator>
1026 Double_t TMath::Mean(Iterator first, Iterator last)
1027 {
1028  Double_t sum = 0;
1029  Double_t sumw = 0;
1030  while ( first != last )
1031  {
1032  sum += *first;
1033  sumw += 1;
1034  first++;
1035  }
1036 
1037  return sum/sumw;
1038 }
1039 
1040 ////////////////////////////////////////////////////////////////////////////////
1041 /// Return the weighted mean of an array defined by the first and
1042 /// last iterators. The w iterator should point to the first element
1043 /// of a vector of weights of the same size as the main array.
1044 template <typename Iterator, typename WeightIterator>
1045 Double_t TMath::Mean(Iterator first, Iterator last, WeightIterator w)
1046 {
1047 
1048  Double_t sum = 0;
1049  Double_t sumw = 0;
1050  int i = 0;
1051  while ( first != last ) {
1052  if ( *w < 0) {
1053  ::Error("TMath::Mean","w[%d] = %.4e < 0 ?!",i,*w);
1054  return 0;
1055  }
1056  sum += (*w) * (*first);
1057  sumw += (*w) ;
1058  ++w;
1059  ++first;
1060  ++i;
1061  }
1062  if (sumw <= 0) {
1063  ::Error("TMath::Mean","sum of weights == 0 ?!");
1064  return 0;
1065  }
1066 
1067  return sum/sumw;
1068 }
1069 
1070 ////////////////////////////////////////////////////////////////////////////////
1071 /// Return the weighted mean of an array a with length n.
1072 template <typename T>
1074 {
1075  if (w) {
1076  return TMath::Mean(a, a+n, w);
1077  } else {
1078  return TMath::Mean(a, a+n);
1079  }
1080 }
1081 
1082 ////////////////////////////////////////////////////////////////////////////////
1083 /// Return the geometric mean of an array defined by the iterators.
1084 /// \f[ GeomMean = (\prod_{i=0}^{n-1} |a[i]|)^{1/n} \f]
1085 template <typename Iterator>
1086 Double_t TMath::GeomMean(Iterator first, Iterator last)
1087 {
1088  Double_t logsum = 0.;
1089  Long64_t n = 0;
1090  while ( first != last ) {
1091  if (*first == 0) return 0.;
1092  Double_t absa = (Double_t) TMath::Abs(*first);
1093  logsum += TMath::Log(absa);
1094  ++first;
1095  ++n;
1096  }
1097 
1098  return TMath::Exp(logsum/n);
1099 }
1100 
1101 ////////////////////////////////////////////////////////////////////////////////
1102 /// Return the geometric mean of an array a of size n.
1103 /// \f[ GeomMean = (\prod_{i=0}^{n-1} |a[i]|)^{1/n} \f]
1104 template <typename T>
1106 {
1107  return TMath::GeomMean(a, a+n);
1108 }
1109 
1110 ////////////////////////////////////////////////////////////////////////////////
1111 /// Return the Standard Deviation of an array defined by the iterators.
1112 /// Note that this function returns the sigma(standard deviation) and
1113 /// not the root mean square of the array.
1114 ///
1115 /// Use the two pass algorithm, which is slower (! a factor of 2) but much more
1116 /// precise. Since we have a vector the 2 pass algorithm is still faster than the
1117 /// Welford algorithm. (See also ROOT-5545)
1118 template <typename Iterator>
1119 Double_t TMath::RMS(Iterator first, Iterator last)
1120 {
1121 
1122  Double_t n = 0;
1123 
1124  Double_t tot = 0;
1125  Double_t mean = TMath::Mean(first,last);
1126  while ( first != last ) {
1127  Double_t x = Double_t(*first);
1128  tot += (x - mean)*(x - mean);
1129  ++first;
1130  ++n;
1131  }
1132  Double_t rms = (n > 1) ? TMath::Sqrt(tot/(n-1)) : 0.0;
1133  return rms;
1134 }
1135 
1136 ////////////////////////////////////////////////////////////////////////////////
1137 /// Return the weighted Standard Deviation of an array defined by the iterators.
1138 /// Note that this function returns the sigma(standard deviation) and
1139 /// not the root mean square of the array.
1140 ///
1141 /// As in the unweighted case use the two pass algorithm
1142 template <typename Iterator, typename WeightIterator>
1143 Double_t TMath::RMS(Iterator first, Iterator last, WeightIterator w)
1144 {
1145  Double_t tot = 0;
1146  Double_t sumw = 0;
1147  Double_t sumw2 = 0;
1148  Double_t mean = TMath::Mean(first,last,w);
1149  while ( first != last ) {
1150  Double_t x = Double_t(*first);
1151  sumw += *w;
1152  sumw2 += (*w) * (*w);
1153  tot += (*w) * (x - mean)*(x - mean);
1154  ++first;
1155  ++w;
1156  }
1157  // use the correction neff/(neff -1) for the unbiased formula
1158  Double_t rms = TMath::Sqrt(tot * sumw/ (sumw*sumw - sumw2) );
1159  return rms;
1160 }
1161 
1162 ////////////////////////////////////////////////////////////////////////////////
1163 /// Return the Standard Deviation of an array a with length n.
1164 /// Note that this function returns the sigma(standard deviation) and
1165 /// not the root mean square of the array.
1166 template <typename T>
1167 Double_t TMath::RMS(Long64_t n, const T *a, const Double_t * w)
1168 {
1169  return (w) ? TMath::RMS(a, a+n, w) : TMath::RMS(a, a+n);
1170 }
1171 
1172 ////////////////////////////////////////////////////////////////////////////////
1173 /// Calculate the Cross Product of two vectors:
1174 /// out = [v1 x v2]
1175 template <typename T> T *TMath::Cross(const T v1[3],const T v2[3], T out[3])
1176 {
1177  out[0] = v1[1] * v2[2] - v1[2] * v2[1];
1178  out[1] = v1[2] * v2[0] - v1[0] * v2[2];
1179  out[2] = v1[0] * v2[1] - v1[1] * v2[0];
1180 
1181  return out;
1182 }
1183 
1184 ////////////////////////////////////////////////////////////////////////////////
1185 /// Calculate a normal vector of a plane.
1186 ///
1187 /// \param[in] p1, p2,p3 3 3D points belonged the plane to define it.
1188 /// \param[out] normal Pointer to 3D normal vector (normalized)
1189 template <typename T> T * TMath::Normal2Plane(const T p1[3],const T p2[3],const T p3[3], T normal[3])
1190 {
1191  T v1[3], v2[3];
1192 
1193  v1[0] = p2[0] - p1[0];
1194  v1[1] = p2[1] - p1[1];
1195  v1[2] = p2[2] - p1[2];
1196 
1197  v2[0] = p3[0] - p1[0];
1198  v2[1] = p3[1] - p1[1];
1199  v2[2] = p3[2] - p1[2];
1200 
1201  NormCross(v1,v2,normal);
1202  return normal;
1203 }
1204 
1205 ////////////////////////////////////////////////////////////////////////////////
1206 /// Function which returns kTRUE if point xp,yp lies inside the
1207 /// polygon defined by the np points in arrays x and y, kFALSE otherwise.
1208 /// Note that the polygon may be open or closed.
1209 template <typename T> Bool_t TMath::IsInside(T xp, T yp, Int_t np, T *x, T *y)
1210 {
1211  Int_t i, j = np-1 ;
1212  Bool_t oddNodes = kFALSE;
1213 
1214  for (i=0; i<np; i++) {
1215  if ((y[i]<yp && y[j]>=yp) || (y[j]<yp && y[i]>=yp)) {
1216  if (x[i]+(yp-y[i])/(y[j]-y[i])*(x[j]-x[i])<xp) {
1217  oddNodes = !oddNodes;
1218  }
1219  }
1220  j=i;
1221  }
1222 
1223  return oddNodes;
1224 }
1225 
1226 ////////////////////////////////////////////////////////////////////////////////
1227 /// Return the median of the array a where each entry i has weight w[i] .
1228 /// Both arrays have a length of at least n . The median is a number obtained
1229 /// from the sorted array a through
1230 ///
1231 /// median = (a[jl]+a[jh])/2. where (using also the sorted index on the array w)
1232 ///
1233 /// sum_i=0,jl w[i] <= sumTot/2
1234 /// sum_i=0,jh w[i] >= sumTot/2
1235 /// sumTot = sum_i=0,n w[i]
1236 ///
1237 /// If w=0, the algorithm defaults to the median definition where it is
1238 /// a number that divides the sorted sequence into 2 halves.
1239 /// When n is odd or n > 1000, the median is kth element k = (n + 1) / 2.
1240 /// when n is even and n < 1000the median is a mean of the elements k = n/2 and k = n/2 + 1.
1241 ///
1242 /// If the weights are supplied (w not 0) all weights must be >= 0
1243 ///
1244 /// If work is supplied, it is used to store the sorting index and assumed to be
1245 /// >= n . If work=0, local storage is used, either on the stack if n < kWorkMax
1246 /// or on the heap for n >= kWorkMax .
1247 template <typename T> Double_t TMath::Median(Long64_t n, const T *a, const Double_t *w, Long64_t *work)
1248 {
1249 
1250  const Int_t kWorkMax = 100;
1251 
1252  if (n <= 0 || !a) return 0;
1253  Bool_t isAllocated = kFALSE;
1254  Double_t median;
1255  Long64_t *ind;
1256  Long64_t workLocal[kWorkMax];
1257 
1258  if (work) {
1259  ind = work;
1260  } else {
1261  ind = workLocal;
1262  if (n > kWorkMax) {
1263  isAllocated = kTRUE;
1264  ind = new Long64_t[n];
1265  }
1266  }
1267 
1268  if (w) {
1269  Double_t sumTot2 = 0;
1270  for (Int_t j = 0; j < n; j++) {
1271  if (w[j] < 0) {
1272  ::Error("TMath::Median","w[%d] = %.4e < 0 ?!",j,w[j]);
1273  if (isAllocated) delete [] ind;
1274  return 0;
1275  }
1276  sumTot2 += w[j];
1277  }
1278 
1279  sumTot2 /= 2.;
1280 
1281  Sort(n, a, ind, kFALSE);
1282 
1283  Double_t sum = 0.;
1284  Int_t jl;
1285  for (jl = 0; jl < n; jl++) {
1286  sum += w[ind[jl]];
1287  if (sum >= sumTot2) break;
1288  }
1289 
1290  Int_t jh;
1291  sum = 2.*sumTot2;
1292  for (jh = n-1; jh >= 0; jh--) {
1293  sum -= w[ind[jh]];
1294  if (sum <= sumTot2) break;
1295  }
1296 
1297  median = 0.5*(a[ind[jl]]+a[ind[jh]]);
1298 
1299  } else {
1300 
1301  if (n%2 == 1)
1302  median = KOrdStat(n, a,n/2, ind);
1303  else {
1304  median = 0.5*(KOrdStat(n, a, n/2 -1, ind)+KOrdStat(n, a, n/2, ind));
1305  }
1306  }
1307 
1308  if (isAllocated)
1309  delete [] ind;
1310  return median;
1311 }
1312 
1313 ////////////////////////////////////////////////////////////////////////////////
1314 /// Returns k_th order statistic of the array a of size n
1315 /// (k_th smallest element out of n elements).
1316 ///
1317 /// C-convention is used for array indexing, so if you want
1318 /// the second smallest element, call KOrdStat(n, a, 1).
1319 ///
1320 /// If work is supplied, it is used to store the sorting index and
1321 /// assumed to be >= n. If work=0, local storage is used, either on
1322 /// the stack if n < kWorkMax or on the heap for n >= kWorkMax.
1323 /// Note that the work index array will not contain the sorted indices but
1324 /// all indices of the smaller element in arbitrary order in work[0,...,k-1] and
1325 /// all indices of the larger element in arbitrary order in work[k+1,..,n-1]
1326 /// work[k] will contain instead the index of the returned element.
1327 ///
1328 /// Taken from "Numerical Recipes in C++" without the index array
1329 /// implemented by Anna Khreshuk.
1330 ///
1331 /// See also the declarations at the top of this file
1332 template <class Element, typename Size>
1333 Element TMath::KOrdStat(Size n, const Element *a, Size k, Size *work)
1334 {
1335 
1336  const Int_t kWorkMax = 100;
1337 
1338  typedef Size Index;
1339 
1340  Bool_t isAllocated = kFALSE;
1341  Size i, ir, j, l, mid;
1342  Index arr;
1343  Index *ind;
1344  Index workLocal[kWorkMax];
1345  Index temp;
1346 
1347  if (work) {
1348  ind = work;
1349  } else {
1350  ind = workLocal;
1351  if (n > kWorkMax) {
1352  isAllocated = kTRUE;
1353  ind = new Index[n];
1354  }
1355  }
1356 
1357  for (Size ii=0; ii<n; ii++) {
1358  ind[ii]=ii;
1359  }
1360  Size rk = k;
1361  l=0;
1362  ir = n-1;
1363  for(;;) {
1364  if (ir<=l+1) { //active partition contains 1 or 2 elements
1365  if (ir == l+1 && a[ind[ir]]<a[ind[l]])
1366  {temp = ind[l]; ind[l]=ind[ir]; ind[ir]=temp;}
1367  Element tmp = a[ind[rk]];
1368  if (isAllocated)
1369  delete [] ind;
1370  return tmp;
1371  } else {
1372  mid = (l+ir) >> 1; //choose median of left, center and right
1373  {temp = ind[mid]; ind[mid]=ind[l+1]; ind[l+1]=temp;}//elements as partitioning element arr.
1374  if (a[ind[l]]>a[ind[ir]]) //also rearrange so that a[l]<=a[l+1]
1375  {temp = ind[l]; ind[l]=ind[ir]; ind[ir]=temp;}
1376 
1377  if (a[ind[l+1]]>a[ind[ir]])
1378  {temp=ind[l+1]; ind[l+1]=ind[ir]; ind[ir]=temp;}
1379 
1380  if (a[ind[l]]>a[ind[l+1]])
1381  {temp = ind[l]; ind[l]=ind[l+1]; ind[l+1]=temp;}
1382 
1383  i=l+1; //initialize pointers for partitioning
1384  j=ir;
1385  arr = ind[l+1];
1386  for (;;){
1387  do i++; while (a[ind[i]]<a[arr]);
1388  do j--; while (a[ind[j]]>a[arr]);
1389  if (j<i) break; //pointers crossed, partitioning complete
1390  {temp=ind[i]; ind[i]=ind[j]; ind[j]=temp;}
1391  }
1392  ind[l+1]=ind[j];
1393  ind[j]=arr;
1394  if (j>=rk) ir = j-1; //keep active the partition that
1395  if (j<=rk) l=i; //contains the k_th element
1396  }
1397  }
1398 }
1399 
1400 #endif
c
#define c(i)
Definition: RSha256.hxx:101
l
auto * l
Definition: textangle.C:4
TMath::Normalize
Double_t Normalize(Double_t v[3])
Normalize a vector v in place.
Definition: TMath.cxx:512
TMath::FDist
Double_t FDist(Double_t F, Double_t N, Double_t M)
Computes the density function of F-distribution (probability function, integral of density,...
Definition: TMath.cxx:2240
ROOT::Math::Cephes::gamma
double gamma(double x)
Definition: SpecFuncCephes.cxx:339
m
auto * m
Definition: textangle.C:8
TMath::BesselJ1
Double_t BesselJ1(Double_t x)
Bessel function J0(x) for any real x.
Definition: TMath.cxx:1644
TMath::Voigt
Double_t Voigt(Double_t x, Double_t sigma, Double_t lg, Int_t r=4)
Computation of Voigt function (normalised).
Definition: TMath.cxx:875
TMath::AreEqualRel
Bool_t AreEqualRel(Double_t af, Double_t bf, Double_t relPrec)
Definition: TMath.h:430
TMath::MaxElement
T MaxElement(Long64_t n, const T *a)
Return maximum of array a of length n.
Definition: TMath.h:959
n
const Int_t n
Definition: legend1.C:16
TMath::Hbar
constexpr Double_t Hbar()
in
Definition: TMath.h:216
TMath::Mean
Double_t Mean(Long64_t n, const T *a, const Double_t *w=0)
Return the weighted mean of an array a with length n.
Definition: TMath.h:1073
first
Definition: first.py:1
TMath::ATanH
Double_t ATanH(Double_t)
Definition: TMath.cxx:90
TMath::Hypot
Double_t Hypot(Double_t x, Double_t y)
Definition: TMath.cxx:57
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
TMath::ChisquareQuantile
Double_t ChisquareQuantile(Double_t p, Double_t ndf)
Evaluate the quantiles of the chi-squared probability distribution function.
Definition: TMath.cxx:2157
TMath::Hbarcgs
constexpr Double_t Hbarcgs()
Definition: TMath.h:223
TMath::ATan2
Double_t ATan2(Double_t y, Double_t x)
Definition: TMath.h:679
TMath::Kcgs
constexpr Double_t Kcgs()
Definition: TMath.h:260
TMath::PoissonI
Double_t PoissonI(Double_t x, Double_t par)
Compute the Discrete Poisson distribution function for (x,par).
Definition: TMath.cxx:592
TMath::BesselK0
Double_t BesselK0(Double_t x)
modified Bessel function I_0(x)
Definition: TMath.cxx:1435
TMath::BesselY0
Double_t BesselY0(Double_t x)
Bessel function J1(x) for any real x.
Definition: TMath.cxx:1680
TMath::GeomMean
Double_t GeomMean(Long64_t n, const T *a)
Return the geometric mean of an array a of size n.
Definition: TMath.h:1105
TMath::RadToDeg
constexpr Double_t RadToDeg()
Conversion from radian to degree:
Definition: TMath.h:73
TMath::NaUncertainty
constexpr Double_t NaUncertainty()
Avogadro constant (Avogadro's Number) uncertainty.
Definition: TMath.h:298
Option_t
const char Option_t
Definition: RtypesCore.h:66
TMath::CauchyDist
Double_t CauchyDist(Double_t x, Double_t t=0, Double_t s=1)
Computes the density of Cauchy distribution at point x by default, standard Cauchy distribution is us...
Definition: TMath.cxx:2141
finite
int finite(double)
TMath::QuietNaN
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754
Definition: TMath.h:901
TMath::BesselI0
Double_t BesselI0(Double_t x)
integer order modified Bessel function K_n(x)
Definition: TMath.cxx:1401
TMath::StruveL1
Double_t StruveL1(Double_t x)
Modified Struve functions of order 0.
Definition: TMath.cxx:1943
TMath::Cos
Double_t Cos(Double_t)
Definition: TMath.h:643
tan
double tan(double)
floor
double floor(double)
TMath::KolmogorovProb
Double_t KolmogorovProb(Double_t z)
Calculates the Kolmogorov distribution function,.
Definition: TMath.cxx:656
TMath::StdDev
Double_t StdDev(Long64_t n, const T *a, const Double_t *w=0)
Definition: TMath.h:530
F
#define F(x, y, z)
ldexp
double ldexp(double, int)
TMath::BesselY1
Double_t BesselY1(Double_t x)
Bessel function Y0(x) for positive x.
Definition: TMath.cxx:1714
TMath::LogE
constexpr Double_t LogE()
Base-10 log of e (to convert ln to log)
Definition: TMath.h:110
TMath::BesselJ0
Double_t BesselJ0(Double_t x)
modified Bessel function K_1(x)
Definition: TMath.cxx:1609
TMath::LogNormal
Double_t LogNormal(Double_t x, Double_t sigma, Double_t theta=0, Double_t m=1)
Computes the density of LogNormal distribution at point x.
Definition: TMath.cxx:2397
r
ROOT::R::TRInterface & r
Definition: Object.C:4
xmax
float xmax
Definition: THbookFile.cxx:95
TMath::Poisson
Double_t Poisson(Double_t x, Double_t par)
Compute the Poisson distribution function for (x,par).
Definition: TMath.cxx:564
Long64_t
long long Long64_t
Definition: RtypesCore.h:73
TMath::Log
Double_t Log(Double_t x)
Definition: TMath.h:760
TMath::StudentI
Double_t StudentI(Double_t T, Double_t ndf)
Calculates the cumulative distribution function of Student's t-distribution second parameter stands f...
Definition: TMath.cxx:2605
TMath::Gaus
Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE)
Calculate a gaussian function with mean and sigma.
Definition: TMath.cxx:448
TMath::GhbarC
constexpr Double_t GhbarC()
in
Definition: TMath.h:161
TMath::GUncertainty
constexpr Double_t GUncertainty()
Gravitational constant uncertainty.
Definition: TMath.h:153
TMath::Sqrt
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
exp
double exp(double)
TMath::CeilNint
Int_t CeilNint(Double_t x)
Definition: TMath.h:699
TMath::DegToRad
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition: TMath.h:81
TMath::Tan
Double_t Tan(Double_t)
Definition: TMath.h:647
Float_t
float Float_t
Definition: RtypesCore.h:57
log
double log(double)
TMath::Exp
Double_t Exp(Double_t x)
Definition: TMath.h:727
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:162
tanh
double tanh(double)
TMath::PiOver2
constexpr Double_t PiOver2()
Definition: TMath.h:51
Int_t
int Int_t
Definition: RtypesCore.h:45
TMath::Qe
constexpr Double_t Qe()
Elementary charge in .
Definition: TMath.h:346
N
#define N
TMath::VavilovI
Double_t VavilovI(Double_t x, Double_t kappa, Double_t beta2)
Returns the value of the Vavilov distribution function.
Definition: TMath.cxx:2767
TMath::LaplaceDist
Double_t LaplaceDist(Double_t x, Double_t alpha=0, Double_t beta=1)
Computes the probability density function of Laplace distribution at point x, with location parameter...
Definition: TMath.cxx:2325
TMath::InvPi
constexpr Double_t InvPi()
Definition: TMath.h:65
sin
double sin(double)
x
Double_t x[n]
Definition: legend1.C:17
cos
double cos(double)
TMath::StruveH0
Double_t StruveH0(Double_t x)
Bessel function Y1(x) for positive x.
Definition: TMath.cxx:1752
TMath::Erf
Double_t Erf(Double_t x)
Computation of the error function erf(x).
Definition: TMath.cxx:184
TMath::KUncertainty
constexpr Double_t KUncertainty()
Boltzmann's constant uncertainty.
Definition: TMath.h:267
TMath::RUncertainty
constexpr Double_t RUncertainty()
Universal gas constant uncertainty.
Definition: TMath.h:316
TMath::Gn
constexpr Double_t Gn()
Standard acceleration of gravity in .
Definition: TMath.h:177
TMath::KolmogorovTest
Double_t KolmogorovTest(Int_t na, const Double_t *a, Int_t nb, const Double_t *b, Option_t *option)
Statistical test whether two one-dimensional sets of points are compatible with coming from the same ...
Definition: TMath.cxx:782
TMath::Abs
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
TMath::Binomial
Double_t Binomial(Int_t n, Int_t k)
Calculate the binomial coefficient n over k.
Definition: TMath.cxx:2083
TMath::Prob
Double_t Prob(Double_t chi2, Int_t ndf)
Computation of the probability for a certain Chi-squared (chi2) and number of degrees of freedom (ndf...
Definition: TMath.cxx:614
TMath::Sort
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMathBase.h:362
TMath::Landau
Double_t Landau(Double_t x, Double_t mpv=0, Double_t sigma=1, Bool_t norm=kFALSE)
The LANDAU function.
Definition: TMath.cxx:469
ceil
double ceil(double)
TMath::C
constexpr Double_t C()
Velocity of light in .
Definition: TMath.h:117
TMath::IsInside
Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y)
Function which returns kTRUE if point xp,yp lies inside the polygon defined by the np points in array...
Definition: TMath.h:1209
TMath::LocMax
Long64_t LocMax(Long64_t n, const T *a)
Return index of array with the maximum element.
Definition: TMath.h:1000
Bool_t
bool Bool_t
Definition: RtypesCore.h:63
TMath::LnGamma
Double_t LnGamma(Double_t z)
Computation of ln[gamma(z)] for all z.
Definition: TMath.cxx:486
ROOT::Math::beta
double beta(double x, double y)
Calculates the beta function.
Definition: SpecFuncMathCore.cxx:111
v
@ v
Definition: rootcling_impl.cxx:3635
b
#define b(i)
Definition: RSha256.hxx:100
TMath::HUncertainty
constexpr Double_t HUncertainty()
Planck's constant uncertainty.
Definition: TMath.h:206
bool
TMath::Cross
T * Cross(const T v1[3], const T v2[3], T out[3])
Calculate the Cross Product of two vectors: out = [v1 x v2].
Definition: TMath.h:1175
TMath::BetaCf
Double_t BetaCf(Double_t x, Double_t a, Double_t b)
Continued fraction evaluation by modified Lentz's method used in calculation of incomplete Beta funct...
Definition: TMath.cxx:1993
TMath::Beta
Double_t Beta(Double_t p, Double_t q)
Calculates Beta-function Gamma(p)*Gamma(q)/Gamma(p+q).
Definition: TMath.cxx:1984
q
float * q
Definition: THbookFile.cxx:89
TMath::LaplaceDistI
Double_t LaplaceDistI(Double_t x, Double_t alpha=0, Double_t beta=1)
Computes the distribution function of Laplace distribution at point x, with location parameter alpha ...
Definition: TMath.cxx:2341
TMath::StudentQuantile
Double_t StudentQuantile(Double_t p, Double_t ndf, Bool_t lower_tail=kTRUE)
Computes quantiles of the Student's t-distribution 1st argument is the probability,...
Definition: TMath.cxx:2633
TMath::HbarUncertainty
constexpr Double_t HbarUncertainty()
uncertainty.
Definition: TMath.h:230
TMath::FloorNint
Int_t FloorNint(Double_t x)
Definition: TMath.h:707
TMath::Gamma
Double_t Gamma(Double_t z)
Computation of gamma(z) for all z.
Definition: TMath.cxx:348
TMath::RootsCubic
Bool_t RootsCubic(const Double_t coef[4], Double_t &a, Double_t &b, Double_t &c)
Calculates roots of polynomial of 3rd order a*x^3 + b*x^2 + c*x + d, where.
Definition: TMath.cxx:1084
TMath::NormQuantile
Double_t NormQuantile(Double_t p)
Computes quantiles for standard normal distribution N(0, 1) at probability p.
Definition: TMath.cxx:2416
TMath::Sigma
constexpr Double_t Sigma()
Stefan-Boltzmann constant in .
Definition: TMath.h:277
atan2
double atan2(double, double)
TMath::LocMin
Long64_t LocMin(Long64_t n, const T *a)
Return index of array with the minimum element.
Definition: TMath.h:972
TMath::Pi
constexpr Double_t Pi()
Definition: TMath.h:37
TMath::Log2
Double_t Log2(Double_t x)
Definition: TMath.cxx:101
log10
double log10(double)
TMath::Nint
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
Definition: TMath.h:713
TMath::BetaIncomplete
Double_t BetaIncomplete(Double_t x, Double_t a, Double_t b)
Calculates the incomplete Beta-function.
Definition: TMath.cxx:2075
TMath::Permute
Bool_t Permute(Int_t n, Int_t *a)
Simple recursive algorithm to find the permutations of n natural numbers, not necessarily all distinc...
Definition: TMath.cxx:2517
TMath::StruveL0
Double_t StruveL0(Double_t x)
Struve functions of order 1.
Definition: TMath.cxx:1897
TMath::Log10
Double_t Log10(Double_t x)
Definition: TMath.h:764
TMath::MWair
constexpr Double_t MWair()
Molecular weight of dry air 1976 US Standard Atmosphere in or
Definition: TMath.h:324
TMath::Limits::Min
static T Min()
Returns maximum representation for type T.
Definition: TMath.h:921
xmin
float xmin
Definition: THbookFile.cxx:95
TMath::IsNaN
Bool_t IsNaN(Double_t x)
Definition: TMath.h:892
epsilon
REAL epsilon
Definition: triangle.c:617
TMath::Floor
Double_t Floor(Double_t x)
Definition: TMath.h:703
a
auto * a
Definition: textangle.C:12
TMath::CUncertainty
constexpr Double_t CUncertainty()
Speed of light uncertainty.
Definition: TMath.h:131
TMath::MinElement
T MinElement(Long64_t n, const T *a)
Return minimum of array a of length n.
Definition: TMath.h:952
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TMath::Limits
Definition: TMath.h:414
TMath::Normal2Plane
T * Normal2Plane(const T v1[3], const T v2[3], const T v3[3], T normal[3])
Calculate a normal vector of a plane.
Definition: TMath.h:1189
Long_t
long Long_t
Definition: RtypesCore.h:54
TMath::Power
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:735
TMath::BesselK1
Double_t BesselK1(Double_t x)
modified Bessel function I_1(x)
Definition: TMath.cxx:1504
TMath::DiLog
Double_t DiLog(Double_t x)
Modified Struve functions of order 1.
Definition: TMath.cxx:110
TMath::Factorial
Double_t Factorial(Int_t i)
Compute factorial(n).
Definition: TMath.cxx:247
TMath::Gcgs
constexpr Double_t Gcgs()
Definition: TMath.h:146
TMath::GammaDist
Double_t GammaDist(Double_t x, Double_t gamma, Double_t mu=0, Double_t beta=1)
Computes the density function of Gamma distribution at point x.
Definition: TMath.cxx:2308
TMath::BreitWigner
Double_t BreitWigner(Double_t x, Double_t mean=0, Double_t gamma=1)
Calculate a Breit Wigner function with mean and gamma.
Definition: TMath.cxx:437
TMath::PiOver4
constexpr Double_t PiOver4()
Definition: TMath.h:58
TMath::Sin
Double_t Sin(Double_t)
Definition: TMath.h:639
TMath::Normalize
Float_t Normalize(Float_t v[3])
Normalize a vector v in place.
Definition: TMath.cxx:495
UInt_t
unsigned int UInt_t
Definition: RtypesCore.h:46
isnan
int isnan(double)
TMath::BetaDist
Double_t BetaDist(Double_t x, Double_t p, Double_t q)
Computes the probability density function of the Beta distribution (the distribution function is comp...
Definition: TMath.cxx:2044
TMath::RMS
Double_t RMS(Long64_t n, const T *a, const Double_t *w=0)
Return the Standard Deviation of an array a with length n.
Definition: TMath.h:1167
y
Double_t y[n]
Definition: legend1.C:17
TMath::BinomialI
Double_t BinomialI(Double_t p, Int_t n, Int_t k)
Suppose an event occurs with probability p per trial Then the probability P of its occurring k or mor...
Definition: TMath.cxx:2108
TMath::Infinity
Double_t Infinity()
Returns an infinity as defined by the IEEE standard.
Definition: TMath.h:914
sqrt
double sqrt(double)
ULong_t
unsigned long ULong_t
Definition: RtypesCore.h:55
TMath::G
constexpr Double_t G()
Gravitational constant in: .
Definition: TMath.h:138
TMath::Rgair
constexpr Double_t Rgair()
Dry Air Gas Constant (R / MWair) in
Definition: TMath.h:332
TMath::SinH
Double_t SinH(Double_t)
Definition: TMath.h:651
TMath::SignalingNaN
Double_t SignalingNaN()
Returns a signaling NaN as defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Signaling_NaN)
Definition: TMath.h:908
TMath::BesselI1
Double_t BesselI1(Double_t x)
modified Bessel function K_0(x)
Definition: TMath.cxx:1469
acos
double acos(double)
TMath::Hcgs
constexpr Double_t Hcgs()
Definition: TMath.h:199
TMath::Limits::Max
static T Max()
Returns minimum double representation.
Definition: TMath.h:928
TMath::ASinH
Double_t ASinH(Double_t)
Definition: TMath.cxx:64
TMath::AreEqualAbs
Bool_t AreEqualAbs(Double_t af, Double_t bf, Double_t epsilon)
Definition: TMath.h:424
TMath::IsNaN
Bool_t IsNaN(Float_t x)
Definition: TMath.h:893
TMath::Student
Double_t Student(Double_t T, Double_t ndf)
Computes density function for Student's t- distribution (the probability function (integral of densit...
Definition: TMath.cxx:2583
TMath::TwoPi
constexpr Double_t TwoPi()
Definition: TMath.h:44
TMath::R
constexpr Double_t R()
Universal gas constant ( ) in
Definition: TMath.h:309
TMath::CosH
Double_t CosH(Double_t)
Definition: TMath.h:655
TMath::StruveH1
Double_t StruveH1(Double_t x)
Struve functions of order 0.
Definition: TMath.cxx:1821
TMath::Freq
Double_t Freq(Double_t x)
Computation of the normal frequency function freq(x).
Definition: TMath.cxx:265
TMath::GhbarCUncertainty
constexpr Double_t GhbarCUncertainty()
uncertainty.
Definition: TMath.h:169
TMath::Ccgs
constexpr Double_t Ccgs()
Definition: TMath.h:124
xmlio::Size
const char * Size
Definition: TXMLSetup.cxx:56
TMath::Vavilov
Double_t Vavilov(Double_t x, Double_t kappa, Double_t beta2)
Returns the value of the Vavilov density function.
Definition: TMath.cxx:2734
TMath::LandauI
Double_t LandauI(Double_t x)
Returns the value of the Landau distribution function at point x.
Definition: TMath.cxx:2796
sigma
const Double_t sigma
Definition: h1analysisProxy.h:11
sum
static long int sum(long int i)
Definition: Factory.cxx:2272
TMath::Erfc
Double_t Erfc(Double_t x)
Compute the complementary error function erfc(x).
Definition: TMath.cxx:194
TMath::GnUncertainty
constexpr Double_t GnUncertainty()
Standard acceleration of gravity uncertainty.
Definition: TMath.h:184
TMath::ACos
Double_t ACos(Double_t)
Definition: TMath.h:669
TMath::ATan
Double_t ATan(Double_t)
Definition: TMath.h:675
TMath::QeUncertainty
constexpr Double_t QeUncertainty()
Elementary charge uncertainty.
Definition: TMath.h:353
LongDouble_t
long double LongDouble_t
Definition: RtypesCore.h:61
v1
@ v1
Definition: rootcling_impl.cxx:3637
Double_t
double Double_t
Definition: RtypesCore.h:59
TMath::BubbleLow
void BubbleLow(Int_t Narr, Double_t *arr1, Int_t *arr2)
Opposite ordering of the array arr2[] to that of BubbleHigh.
Definition: TMath.cxx:1328
TMath::Ldexp
Double_t Ldexp(Double_t x, Int_t exp)
Definition: TMath.h:731
TMath::H
constexpr Double_t H()
Planck's constant in .
Definition: TMath.h:192
TMath::ACosH
Double_t ACosH(Double_t)
Definition: TMath.cxx:77
TMath::BesselI
Double_t BesselI(Int_t n, Double_t x)
Compute the Integer Order Modified Bessel function I_n(x) for n=0,1,2,...
Definition: TMath.cxx:1565
TMath::Ln10
constexpr Double_t Ln10()
Natural log of 10 (to convert log to ln)
Definition: TMath.h:103
TMath::Limits::Epsilon
static T Epsilon()
Returns minimum double representation.
Definition: TMath.h:935
v3
@ v3
Definition: rootcling_impl.cxx:3639
TMath::BesselK
Double_t BesselK(Int_t n, Double_t x)
integer order modified Bessel function I_n(x)
Definition: TMath.cxx:1536
TMath::Sq
Double_t Sq(Double_t x)
Definition: TMath.h:687
v2
@ v2
Definition: rootcling_impl.cxx:3638
TMath::FDistI
Double_t FDistI(Double_t F, Double_t N, Double_t M)
Calculates the cumulative distribution function of F-distribution, this function occurs in the statis...
Definition: TMath.cxx:2259
TMath::NormCross
T NormCross(const T v1[3], const T v2[3], T out[3])
Calculate the Normalized Cross Product of two vectors.
Definition: TMath.h:944
TMath::TanH
Double_t TanH(Double_t)
Definition: TMath.h:659
TMath::EulerGamma
constexpr Double_t EulerGamma()
Euler-Mascheroni Constant.
Definition: TMath.h:339
TMath::BetaDistI
Double_t BetaDistI(Double_t x, Double_t p, Double_t q)
Computes the distribution function of the Beta distribution.
Definition: TMath.cxx:2062
TMath::Sqrt2
constexpr Double_t Sqrt2()
Definition: TMath.h:88
ROOT::Math::Chebyshev::T
double T(double x)
Definition: ChebyshevPol.h:34
d
#define d(i)
Definition: RSha256.hxx:102
sinh
double sinh(double)
atan
double atan(double)
TMath::KOrdStat
Element KOrdStat(Size n, const Element *a, Size k, Size *work=0)
Returns k_th order statistic of the array a of size n (k_th smallest element out of n elements).
Definition: TMath.h:1333
cosh
double cosh(double)
TMathBase.h
TMath::HC
constexpr Double_t HC()
in
Definition: TMath.h:238
TMath::ErfInverse
Double_t ErfInverse(Double_t x)
returns the inverse error function x must be <-1<x<1
Definition: TMath.cxx:203
TMath::K
constexpr Double_t K()
Boltzmann's constant in .
Definition: TMath.h:253
TMath::ErfcInverse
Double_t ErfcInverse(Double_t x)
returns the inverse of the complementary error function x must be 0<x<2 implement using the quantile ...
Definition: TMath.cxx:237
type
int type
Definition: TGX11.cxx:121
RooFit::Index
RooCmdArg Index(RooCategory &icat)
Definition: RooGlobalFunc.cxx:97
TMath::Hash
ULong_t Hash(const void *txt, Int_t ntxt)
Calculates hash index from any char string.
Definition: TMath.cxx:1383
TMath::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)
Computes sample quantiles, corresponding to the given probabilities.
Definition: TMath.cxx:1183
TMath::HCcgs
constexpr Double_t HCcgs()
Definition: TMath.h:245
pow
double pow(double, double)
asin
double asin(double)
TMath::SigmaUncertainty
constexpr Double_t SigmaUncertainty()
Stefan-Boltzmann constant uncertainty.
Definition: TMath.h:284
TMath::Na
constexpr Double_t Na()
Avogadro constant (Avogadro's Number) in .
Definition: TMath.h:291
TMath::BubbleHigh
void BubbleHigh(Int_t Narr, Double_t *arr1, Int_t *arr2)
Bubble sort variant to obtain the order of an array's elements into an index in order to do more usef...
Definition: TMath.cxx:1289
ROOT
VSD Structures.
Definition: StringConv.hxx:21
TMath::E
constexpr Double_t E()
Base of natural log:
Definition: TMath.h:96
Cross
TGLVector3 Cross(const TGLVector3 &v1, const TGLVector3 &v2)
Definition: TGLUtil.h:323
TMath::Finite
Int_t Finite(Double_t x)
Check if it is finite with a mask in order to be consistent in presence of fast math.
Definition: TMath.h:771
Math
Namespace for new Math classes and functions.
TMath::Ceil
Double_t Ceil(Double_t x)
Definition: TMath.h:695
TMath::ASin
Double_t ASin(Double_t)
Definition: TMath.h:663
TMath
TMath.
Definition: TMathBase.h:35
int
Error
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition: TError.cxx:187
TMath::Median
Double_t Median(Long64_t n, const T *a, const Double_t *w=0, Long64_t *work=0)
Return the median of the array a where each entry i has weight w[i] .
Definition: TMath.h:1247
TError.h