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.
131 constexpr Double_t CUncertainty()
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  return 6.673e-11;
141 }
142 
143 ////////////////////////////////////////////////////////////////////////////////
144 /// \f$ cm^{3} g^{-1} s^{-2} \f$
145 constexpr Double_t Gcgs()
146 {
147  return G() / 1000.0;
148 }
149 
150 ////////////////////////////////////////////////////////////////////////////////
151 /// Gravitational constant uncertainty.
152 constexpr Double_t GUncertainty()
153 {
154  return 0.010e-11;
155 }
156 
157 ////////////////////////////////////////////////////////////////////////////////
158 /// \f$ \frac{G}{\hbar C} \f$ in \f$ (GeV/c^{2})^{-2} \f$
159 constexpr Double_t GhbarC()
160 {
161  return 6.707e-39;
162 }
163 
164 ////////////////////////////////////////////////////////////////////////////////
165 /// \f$ \frac{G}{\hbar C} \f$ uncertainty.
166 constexpr Double_t GhbarCUncertainty()
167 {
168  return 0.010e-39;
169 }
170 
171 ////////////////////////////////////////////////////////////////////////////////
172 /// Standard acceleration of gravity in \f$ m s^{-2} \f$
173 constexpr Double_t Gn()
174 {
175  return 9.80665;
176 }
177 
178 ////////////////////////////////////////////////////////////////////////////////
179 /// Standard acceleration of gravity uncertainty.
180 constexpr Double_t GnUncertainty()
181 {
182  return 0.0;
183 }
184 
185 ////////////////////////////////////////////////////////////////////////////////
186 /// Planck's constant in \f$ J s \f$
187 /// \f[ h \f]
188 constexpr Double_t H()
189 {
190  return 6.62606876e-34;
191 }
192 
193 ////////////////////////////////////////////////////////////////////////////////
194 /// \f$ erg s \f$
195 constexpr Double_t Hcgs()
196 {
197  return 1.0e7 * H();
198 }
199 
200 ////////////////////////////////////////////////////////////////////////////////
201 /// Planck's constant uncertainty.
202 constexpr Double_t HUncertainty()
203 {
204  return 0.00000052e-34;
205 }
206 
207 ////////////////////////////////////////////////////////////////////////////////
208 /// \f$ \hbar \f$ in \f$ J s \f$
209 /// \f[ \hbar = \frac{h}{2\pi} \f]
210 constexpr Double_t Hbar()
211 {
212  return 1.054571596e-34;
213 }
214 
215 ////////////////////////////////////////////////////////////////////////////////
216 /// \f$ erg s \f$
217 constexpr Double_t Hbarcgs()
218 {
219  return 1.0e7 * Hbar();
220 }
221 
222 ////////////////////////////////////////////////////////////////////////////////
223 /// \f$ \hbar \f$ uncertainty.
224 constexpr Double_t HbarUncertainty()
225 {
226  return 0.000000082e-34;
227 }
228 
229 ////////////////////////////////////////////////////////////////////////////////
230 /// \f$ hc \f$ in \f$ J m \f$
231 constexpr Double_t HC()
232 {
233  return H() * C();
234 }
235 
236 ////////////////////////////////////////////////////////////////////////////////
237 /// \f$ erg cm \f$
238 constexpr Double_t HCcgs()
239 {
240  return Hcgs() * Ccgs();
241 }
242 
243 ////////////////////////////////////////////////////////////////////////////////
244 /// Boltzmann's constant in \f$ J K^{-1} \f$
245 /// \f[ k \f]
246 constexpr Double_t K()
247 {
248  return 1.3806503e-23;
249 }
250 
251 ////////////////////////////////////////////////////////////////////////////////
252 /// \f$ erg K^{-1} \f$
253 constexpr Double_t Kcgs()
254 {
255  return 1.0e7 * K();
256 }
257 
258 ////////////////////////////////////////////////////////////////////////////////
259 /// Boltzmann's constant uncertainty.
260 constexpr Double_t KUncertainty()
261 {
262  return 0.0000024e-23;
263 }
264 
265 ////////////////////////////////////////////////////////////////////////////////
266 /// Stefan-Boltzmann constant in \f$ W m^{-2} K^{-4}\f$
267 /// \f[ \sigma \f]
268 constexpr Double_t Sigma()
269 {
270  return 5.6704e-8;
271 }
272 
273 ////////////////////////////////////////////////////////////////////////////////
274 /// Stefan-Boltzmann constant uncertainty.
275 constexpr Double_t SigmaUncertainty()
276 {
277  return 0.000040e-8;
278 }
279 
280 ////////////////////////////////////////////////////////////////////////////////
281 /// Avogadro constant (Avogadro's Number) in \f$ mol^{-1} \f$
282 constexpr Double_t Na()
283 {
284  return 6.02214199e+23;
285 }
286 
287 ////////////////////////////////////////////////////////////////////////////////
288 /// Avogadro constant (Avogadro's Number) uncertainty.
289 constexpr Double_t NaUncertainty()
290 {
291  return 0.00000047e+23;
292 }
293 
294 ////////////////////////////////////////////////////////////////////////////////
295 /// [Universal gas constant](http://scienceworld.wolfram.com/physics/UniversalGasConstant.html)
296 /// (\f$ Na K \f$) in \f$ J K^{-1} mol^{-1} \f$
297 //
298 constexpr Double_t R()
299 {
300  return K() * Na();
301 }
302 
303 ////////////////////////////////////////////////////////////////////////////////
304 /// Universal gas constant uncertainty.
305 constexpr Double_t RUncertainty()
306 {
307  return R() * ((KUncertainty() / K()) + (NaUncertainty() / Na()));
308 }
309 
310 ////////////////////////////////////////////////////////////////////////////////
311 /// [Molecular weight of dry air 1976 US Standard Atmosphere](http://atmos.nmsu.edu/jsdap/encyclopediawork.html)
312 /// in \f$ kg kmol^{-1} \f$ or \f$ gm mol^{-1} \f$
313 constexpr Double_t MWair()
314 {
315  return 28.9644;
316 }
317 
318 ////////////////////////////////////////////////////////////////////////////////
319 /// [Dry Air Gas Constant (R / MWair)](http://atmos.nmsu.edu/education_and_outreach/encyclopedia/gas_constant.htm)
320 /// in \f$ J kg^{-1} K^{-1} \f$
321 constexpr Double_t Rgair()
322 {
323  return (1000.0 * R()) / MWair();
324 }
325 
326 ////////////////////////////////////////////////////////////////////////////////
327 /// Euler-Mascheroni Constant.
328 constexpr Double_t EulerGamma()
329 {
330  return 0.577215664901532860606512090082402431042;
331 }
332 
333 ////////////////////////////////////////////////////////////////////////////////
334 /// Elementary charge in \f$ C \f$ .
335 constexpr Double_t Qe()
336 {
337  return 1.602176462e-19;
338 }
339 
340 ////////////////////////////////////////////////////////////////////////////////
341 /// Elementary charge uncertainty.
342 constexpr Double_t QeUncertainty()
343 {
344  return 0.000000063e-19;
345 }
346 
347 ////////////////////////////////////////////////////////////////////////////////
348 // Mathematical Functions
349 
350 ////////////////////////////////////////////////////////////////////////////////
351 // Trigonometrical Functions
352 
353 inline Double_t Sin(Double_t);
354 inline Double_t Cos(Double_t);
355 inline Double_t Tan(Double_t);
356 inline Double_t SinH(Double_t);
357 inline Double_t CosH(Double_t);
358 inline Double_t TanH(Double_t);
359 inline Double_t ASin(Double_t);
360 inline Double_t ACos(Double_t);
361 inline Double_t ATan(Double_t);
362 inline Double_t ATan2(Double_t y, Double_t x);
367 
368 ////////////////////////////////////////////////////////////////////////////////
369 // Elementary Functions
370 
371 inline Double_t Ceil(Double_t x);
372 inline Int_t CeilNint(Double_t x);
373 inline Double_t Floor(Double_t x);
374 inline Int_t FloorNint(Double_t x);
375 template <typename T>
376 inline Int_t Nint(T x);
377 
378 inline Double_t Sq(Double_t x);
379 inline Double_t Sqrt(Double_t x);
380 inline Double_t Exp(Double_t x);
381 inline Double_t Ldexp(Double_t x, Int_t exp);
386 inline Double_t Power(Double_t x, Double_t y);
387 inline Double_t Power(Double_t x, Int_t y);
388 inline Double_t Log(Double_t x);
390 inline Double_t Log10(Double_t x);
391 inline Int_t Finite(Double_t x);
392 inline Int_t Finite(Float_t x);
393 inline Bool_t IsNaN(Double_t x);
394 inline Bool_t IsNaN(Float_t x);
395 
396 inline Double_t QuietNaN();
397 inline Double_t SignalingNaN();
398 inline Double_t Infinity();
399 
400 template <typename T>
401 struct Limits {
402  inline static T Min();
403  inline static T Max();
404  inline static T Epsilon();
405  };
406 
407  // Some integer math
408  Long_t Hypot(Long_t x, Long_t y); // sqrt(px*px + py*py)
409 
410  // Comparing floating points
412  //return kTRUE if absolute difference between af and bf is less than epsilon
413  return TMath::Abs(af-bf) < epsilon ||
414  TMath::Abs(af - bf) < Limits<Double_t>::Min(); // handle 0 < 0 case
415 
416  }
417  inline Bool_t AreEqualRel(Double_t af, Double_t bf, Double_t relPrec) {
418  //return kTRUE if relative difference between af and bf is less than relPrec
419  return TMath::Abs(af - bf) <= 0.5 * relPrec * (TMath::Abs(af) + TMath::Abs(bf)) ||
420  TMath::Abs(af - bf) < Limits<Double_t>::Min(); // handle denormals
421  }
422 
423  /////////////////////////////////////////////////////////////////////////////
424  // Array Algorithms
425 
426  // Min, Max of an array
427  template <typename T> T MinElement(Long64_t n, const T *a);
428  template <typename T> T MaxElement(Long64_t n, const T *a);
429 
430  // Locate Min, Max element number in an array
431  template <typename T> Long64_t LocMin(Long64_t n, const T *a);
432  template <typename Iterator> Iterator LocMin(Iterator first, Iterator last);
433  template <typename T> Long64_t LocMax(Long64_t n, const T *a);
434  template <typename Iterator> Iterator LocMax(Iterator first, Iterator last);
435 
436  // Hashing
437  ULong_t Hash(const void *txt, Int_t ntxt);
438  ULong_t Hash(const char *str);
439 
440  void BubbleHigh(Int_t Narr, Double_t *arr1, Int_t *arr2);
441  void BubbleLow (Int_t Narr, Double_t *arr1, Int_t *arr2);
442 
443  Bool_t Permute(Int_t n, Int_t *a); // Find permutations
444 
445  /////////////////////////////////////////////////////////////////////////////
446  // Geometrical Functions
447 
448  //Sample quantiles
449  void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob,
450  Bool_t isSorted=kTRUE, Int_t *index = 0, Int_t type=7);
451 
452  // IsInside
453  template <typename T> Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y);
454 
455  // Calculate the Cross Product of two vectors
456  template <typename T> T *Cross(const T v1[3],const T v2[3], T out[3]);
457 
458  Float_t Normalize(Float_t v[3]); // Normalize a vector
459  Double_t Normalize(Double_t v[3]); // Normalize a vector
460 
461  //Calculate the Normalized Cross Product of two vectors
462  template <typename T> inline T NormCross(const T v1[3],const T v2[3],T out[3]);
463 
464  // Calculate a normal vector of a plane
465  template <typename T> T *Normal2Plane(const T v1[3],const T v2[3],const T v3[3], T normal[3]);
466 
467  /////////////////////////////////////////////////////////////////////////////
468  // Polynomial Functions
469 
470  Bool_t RootsCubic(const Double_t coef[4],Double_t &a, Double_t &b, Double_t &c);
471 
472  /////////////////////////////////////////////////////////////////////////////
473  // Statistic Functions
474 
475  Double_t Binomial(Int_t n,Int_t k); // Calculate the binomial coefficient n over k
484  Double_t KolmogorovTest(Int_t na, const Double_t *a, Int_t nb, const Double_t *b, Option_t *option);
493  Double_t Prob(Double_t chi2,Int_t ndf);
497  Double_t Vavilov(Double_t x, Double_t kappa, Double_t beta2);
498  Double_t VavilovI(Double_t x, Double_t kappa, Double_t beta2);
500 
501  /////////////////////////////////////////////////////////////////////////////
502  // Statistics over arrays
503 
504  //Mean, Geometric Mean, Median, RMS(sigma)
505 
506  template <typename T> Double_t Mean(Long64_t n, const T *a, const Double_t *w=0);
507  template <typename Iterator> Double_t Mean(Iterator first, Iterator last);
508  template <typename Iterator, typename WeightIterator> Double_t Mean(Iterator first, Iterator last, WeightIterator wfirst);
509 
510  template <typename T> Double_t GeomMean(Long64_t n, const T *a);
511  template <typename Iterator> Double_t GeomMean(Iterator first, Iterator last);
512 
513  template <typename T> Double_t RMS(Long64_t n, const T *a, const Double_t *w=0);
514  template <typename Iterator> Double_t RMS(Iterator first, Iterator last);
515  template <typename Iterator, typename WeightIterator> Double_t RMS(Iterator first, Iterator last, WeightIterator wfirst);
516 
517  template <typename T> Double_t StdDev(Long64_t n, const T *a, const Double_t * w = 0) { return RMS<T>(n,a,w); }
518  template <typename Iterator> Double_t StdDev(Iterator first, Iterator last) { return RMS<Iterator>(first,last); }
519  template <typename Iterator, typename WeightIterator> Double_t StdDev(Iterator first, Iterator last, WeightIterator wfirst) { return RMS<Iterator,WeightIterator>(first,last,wfirst); }
520 
521  template <typename T> Double_t Median(Long64_t n, const T *a, const Double_t *w=0, Long64_t *work=0);
522 
523  //k-th order statistic
524  template <class Element, typename Size> Element KOrdStat(Size n, const Element *a, Size k, Size *work = 0);
525 
526  /////////////////////////////////////////////////////////////////////////////
527  // Special Functions
528 
534 
535  // Bessel functions
536  Double_t BesselI(Int_t n,Double_t x); /// integer order modified Bessel function I_n(x)
537  Double_t BesselK(Int_t n,Double_t x); /// integer order modified Bessel function K_n(x)
538  Double_t BesselI0(Double_t x); /// modified Bessel function I_0(x)
539  Double_t BesselK0(Double_t x); /// modified Bessel function K_0(x)
540  Double_t BesselI1(Double_t x); /// modified Bessel function I_1(x)
541  Double_t BesselK1(Double_t x); /// modified Bessel function K_1(x)
542  Double_t BesselJ0(Double_t x); /// Bessel function J0(x) for any real x
543  Double_t BesselJ1(Double_t x); /// Bessel function J1(x) for any real x
544  Double_t BesselY0(Double_t x); /// Bessel function Y0(x) for positive x
545  Double_t BesselY1(Double_t x); /// Bessel function Y1(x) for positive x
546  Double_t StruveH0(Double_t x); /// Struve functions of order 0
547  Double_t StruveH1(Double_t x); /// Struve functions of order 1
548  Double_t StruveL0(Double_t x); /// Modified Struve functions of order 0
549  Double_t StruveL1(Double_t x); /// Modified Struve functions of order 1
550 
561 }
562 
563 ////////////////////////////////////////////////////////////////////////////////
564 // Trig and other functions
565 
566 #include <float.h>
567 
568 #if defined(R__WIN32) && !defined(__CINT__)
569 # ifndef finite
570 # define finite _finite
571 # endif
572 #endif
573 #if defined(R__AIX) || defined(R__SOLARIS_CC50) || \
574  defined(R__HPUX11) || defined(R__GLIBC) || \
575  (defined(R__MACOSX) )
576 // math functions are defined inline so we have to include them here
577 # include <math.h>
578 # ifdef R__SOLARIS_CC50
579  extern "C" { int finite(double); }
580 # endif
581 // # if defined(R__GLIBC) && defined(__STRICT_ANSI__)
582 // # ifndef finite
583 // # define finite __finite
584 // # endif
585 // # ifndef isnan
586 // # define isnan __isnan
587 // # endif
588 // # endif
589 #else
590 // don't want to include complete <math.h>
591 extern "C" {
592  extern double sin(double);
593  extern double cos(double);
594  extern double tan(double);
595  extern double sinh(double);
596  extern double cosh(double);
597  extern double tanh(double);
598  extern double asin(double);
599  extern double acos(double);
600  extern double atan(double);
601  extern double atan2(double, double);
602  extern double sqrt(double);
603  extern double exp(double);
604  extern double pow(double, double);
605  extern double log(double);
606  extern double log10(double);
607 #ifndef R__WIN32
608 # if !defined(finite)
609  extern int finite(double);
610 # endif
611 # if !defined(isnan)
612  extern int isnan(double);
613 # endif
614  extern double ldexp(double, int);
615  extern double ceil(double);
616  extern double floor(double);
617 #else
618  _CRTIMP double ldexp(double, int);
619  _CRTIMP double ceil(double);
620  _CRTIMP double floor(double);
621 #endif
622 }
623 #endif
624 
625 ////////////////////////////////////////////////////////////////////////////////
627  { return sin(x); }
628 
629 ////////////////////////////////////////////////////////////////////////////////
631  { return cos(x); }
632 
633 ////////////////////////////////////////////////////////////////////////////////
635  { return tan(x); }
636 
637 ////////////////////////////////////////////////////////////////////////////////
639  { return sinh(x); }
640 
641 ////////////////////////////////////////////////////////////////////////////////
643  { return cosh(x); }
644 
645 ////////////////////////////////////////////////////////////////////////////////
647  { return tanh(x); }
648 
649 ////////////////////////////////////////////////////////////////////////////////
651  { if (x < -1.) return -TMath::Pi()/2;
652  if (x > 1.) return TMath::Pi()/2;
653  return asin(x);
654  }
655 
656 ////////////////////////////////////////////////////////////////////////////////
658  { if (x < -1.) return TMath::Pi();
659  if (x > 1.) return 0;
660  return acos(x);
661  }
662 
663 ////////////////////////////////////////////////////////////////////////////////
665  { return atan(x); }
666 
667 ////////////////////////////////////////////////////////////////////////////////
669  { if (x != 0) return atan2(y, x);
670  if (y == 0) return 0;
671  if (y > 0) return Pi()/2;
672  else return -Pi()/2;
673  }
674 
675 ////////////////////////////////////////////////////////////////////////////////
677  { return x*x; }
678 
679 ////////////////////////////////////////////////////////////////////////////////
681  { return sqrt(x); }
682 
683 ////////////////////////////////////////////////////////////////////////////////
685  { return ceil(x); }
686 
687 ////////////////////////////////////////////////////////////////////////////////
689  { return TMath::Nint(ceil(x)); }
690 
691 ////////////////////////////////////////////////////////////////////////////////
693  { return floor(x); }
694 
695 ////////////////////////////////////////////////////////////////////////////////
697  { return TMath::Nint(floor(x)); }
698 
699 ////////////////////////////////////////////////////////////////////////////////
700 /// Round to nearest integer. Rounds half integers to the nearest even integer.
701 template<typename T>
703 {
704  int i;
705  if (x >= 0) {
706  i = int(x + 0.5);
707  if ( i & 1 && x + 0.5 == T(i) ) i--;
708  } else {
709  i = int(x - 0.5);
710  if ( i & 1 && x - 0.5 == T(i) ) i++;
711  }
712  return i;
713 }
714 
715 ////////////////////////////////////////////////////////////////////////////////
717  { return exp(x); }
718 
719 ////////////////////////////////////////////////////////////////////////////////
721  { return ldexp(x, exp); }
722 
723 ////////////////////////////////////////////////////////////////////////////////
725  { return std::pow(x,y); }
726 
727 ////////////////////////////////////////////////////////////////////////////////
729  { return std::pow(x,(LongDouble_t)y); }
730 
731 ////////////////////////////////////////////////////////////////////////////////
733  { return std::pow(x,y); }
734 
735 ////////////////////////////////////////////////////////////////////////////////
737  { return pow(x, y); }
738 
739 ////////////////////////////////////////////////////////////////////////////////
741 #ifdef R__ANSISTREAM
742  return std::pow(x, y);
743 #else
744  return pow(x, (Double_t) y);
745 #endif
746 }
747 
748 ////////////////////////////////////////////////////////////////////////////////
750  { return log(x); }
751 
752 ////////////////////////////////////////////////////////////////////////////////
754  { return log10(x); }
755 
756 ////////////////////////////////////////////////////////////////////////////////
757 /// Check if it is finite with a mask in order to be consistent in presence of
758 /// fast math.
759 /// Inspired from the CMSSW FWCore/Utilities package
761 #if defined(R__FAST_MATH)
762 
763 {
764  const unsigned long long mask = 0x7FF0000000000000LL;
765  union { unsigned long long l; double d;} v;
766  v.d =x;
767  return (v.l&mask)!=mask;
768 }
769 #else
770 # if defined(R__HPUX11)
771  { return isfinite(x); }
772 # elif defined(R__MACOSX)
773 # ifdef isfinite
774  // from math.h
775  { return isfinite(x); }
776 # else
777  // from cmath
778  { return std::isfinite(x); }
779 # endif
780 # else
781  { return finite(x); }
782 # endif
783 #endif
784 
785 ////////////////////////////////////////////////////////////////////////////////
786 /// Check if it is finite with a mask in order to be consistent in presence of
787 /// fast math.
788 /// Inspired from the CMSSW FWCore/Utilities package
790 #if defined(R__FAST_MATH)
791 
792 {
793  const unsigned int mask = 0x7f800000;
794  union { unsigned int l; float d;} v;
795  v.d =x;
796  return (v.l&mask)!=mask;
797 }
798 #else
799 { return std::isfinite(x); }
800 #endif
801 
802 // This namespace provides all the routines necessary for checking if a number
803 // is a NaN also in presence of optimisations affecting the behaviour of the
804 // floating point calculations.
805 // Inspired from the CMSSW FWCore/Utilities package
806 
807 #if defined (R__FAST_MATH)
808 namespace ROOT {
809 namespace Internal {
810 namespace Math {
811 // abridged from GNU libc 2.6.1 - in detail from
812 // math/math_private.h
813 // sysdeps/ieee754/ldbl-96/math_ldbl.h
814 
815 // part of ths file:
816  /*
817  * ====================================================
818  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
819  *
820  * Developed at SunPro, a Sun Microsystems, Inc. business.
821  * Permission to use, copy, modify, and distribute this
822  * software is freely granted, provided that this notice
823  * is preserved.
824  * ====================================================
825  */
826 
827  // A union which permits us to convert between a double and two 32 bit ints.
828  typedef union {
829  Double_t value;
830  struct {
831  UInt_t lsw;
832  UInt_t msw;
833  } parts;
834  } ieee_double_shape_type;
835 
836 #define EXTRACT_WORDS(ix0,ix1,d) \
837  do { \
838  ieee_double_shape_type ew_u; \
839  ew_u.value = (d); \
840  (ix0) = ew_u.parts.msw; \
841  (ix1) = ew_u.parts.lsw; \
842  } while (0)
843 
844  inline Bool_t IsNaN(Double_t x)
845  {
846  UInt_t hx, lx;
847 
848  EXTRACT_WORDS(hx, lx, x);
849 
850  lx |= hx & 0xfffff;
851  hx &= 0x7ff00000;
852  return (hx == 0x7ff00000) && (lx != 0);
853  }
854 
855  typedef union {
856  Float_t value;
857  UInt_t word;
858  } ieee_float_shape_type;
859 
860 #define GET_FLOAT_WORD(i,d) \
861  do { \
862  ieee_float_shape_type gf_u; \
863  gf_u.value = (d); \
864  (i) = gf_u.word; \
865  } while (0)
866 
867  inline Bool_t IsNaN(Float_t x)
868  {
869  UInt_t wx;
870  GET_FLOAT_WORD (wx, x);
871  wx &= 0x7fffffff;
872  return (Bool_t)(wx > 0x7f800000);
873  }
874 } } } // end NS ROOT::Internal::Math
875 #endif // End R__FAST_MATH
876 
877 #if defined(R__FAST_MATH)
880 #else
881  inline Bool_t TMath::IsNaN(Double_t x) { return std::isnan(x); }
882  inline Bool_t TMath::IsNaN(Float_t x) { return std::isnan(x); }
883 #endif
884 
885 ////////////////////////////////////////////////////////////////////////////////
886 // Wrapper to numeric_limits
887 
888 ////////////////////////////////////////////////////////////////////////////////
889 /// Returns a quiet NaN as [defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Quiet_NaN)
891 
892  return std::numeric_limits<Double_t>::quiet_NaN();
893 }
894 
895 ////////////////////////////////////////////////////////////////////////////////
896 /// Returns a signaling NaN as defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Signaling_NaN)
898  return std::numeric_limits<Double_t>::signaling_NaN();
899 }
900 
901 ////////////////////////////////////////////////////////////////////////////////
902 /// Returns an infinity as defined by the IEEE standard
904  return std::numeric_limits<Double_t>::infinity();
905 }
906 
907 ////////////////////////////////////////////////////////////////////////////////
908 /// Returns maximum representation for type T
909 template<typename T>
911  return (std::numeric_limits<T>::min)(); //N.B. use this signature to avoid class with macro min() on Windows
912 }
913 
914 ////////////////////////////////////////////////////////////////////////////////
915 /// Returns minimum double representation
916 template<typename T>
918  return (std::numeric_limits<T>::max)(); //N.B. use this signature to avoid class with macro max() on Windows
919 }
920 
921 ////////////////////////////////////////////////////////////////////////////////
922 /// Returns minimum double representation
923 template<typename T>
926 }
927 
928 ////////////////////////////////////////////////////////////////////////////////
929 // Advanced.
930 
931 ////////////////////////////////////////////////////////////////////////////////
932 /// Calculate the Normalized Cross Product of two vectors
933 template <typename T> inline T TMath::NormCross(const T v1[3],const T v2[3],T out[3])
934 {
935  return Normalize(Cross(v1,v2,out));
936 }
937 
938 ////////////////////////////////////////////////////////////////////////////////
939 /// Return minimum of array a of length n.
940 template <typename T>
942  return *std::min_element(a,a+n);
943 }
944 
945 ////////////////////////////////////////////////////////////////////////////////
946 /// Return maximum of array a of length n.
947 template <typename T>
949  return *std::max_element(a,a+n);
950 }
951 
952 ////////////////////////////////////////////////////////////////////////////////
953 /// Return index of array with the minimum element.
954 /// If more than one element is minimum returns first found.
955 ///
956 /// Implement here since this one is found to be faster (mainly on 64 bit machines)
957 /// than stl generic implementation.
958 /// When performing the comparison, the STL implementation needs to de-reference both the array iterator
959 /// and the iterator pointing to the resulting minimum location
960 template <typename T>
962  if (n <= 0 || !a) return -1;
963  T xmin = a[0];
964  Long64_t loc = 0;
965  for (Long64_t i = 1; i < n; i++) {
966  if (xmin > a[i]) {
967  xmin = a[i];
968  loc = i;
969  }
970  }
971  return loc;
972 }
973 
974 ////////////////////////////////////////////////////////////////////////////////
975 /// Return index of array with the minimum element.
976 /// If more than one element is minimum returns first found.
977 template <typename Iterator>
978 Iterator TMath::LocMin(Iterator first, Iterator last) {
979 
980  return std::min_element(first, last);
981 }
982 
983 ////////////////////////////////////////////////////////////////////////////////
984 /// Return index of array with the maximum element.
985 /// If more than one element is maximum returns first found.
986 ///
987 /// Implement here since it is faster (see comment in LocMin function)
988 template <typename T>
990  if (n <= 0 || !a) return -1;
991  T xmax = a[0];
992  Long64_t loc = 0;
993  for (Long64_t i = 1; i < n; i++) {
994  if (xmax < a[i]) {
995  xmax = a[i];
996  loc = i;
997  }
998  }
999  return loc;
1000 }
1001 
1002 ////////////////////////////////////////////////////////////////////////////////
1003 /// Return index of array with the maximum element.
1004 /// If more than one element is maximum returns first found.
1005 template <typename Iterator>
1006 Iterator TMath::LocMax(Iterator first, Iterator last)
1007 {
1008 
1009  return std::max_element(first, last);
1010 }
1011 
1012 ////////////////////////////////////////////////////////////////////////////////
1013 /// Return the weighted mean of an array defined by the iterators.
1014 template <typename Iterator>
1015 Double_t TMath::Mean(Iterator first, Iterator last)
1016 {
1017  Double_t sum = 0;
1018  Double_t sumw = 0;
1019  while ( first != last )
1020  {
1021  sum += *first;
1022  sumw += 1;
1023  first++;
1024  }
1025 
1026  return sum/sumw;
1027 }
1028 
1029 ////////////////////////////////////////////////////////////////////////////////
1030 /// Return the weighted mean of an array defined by the first and
1031 /// last iterators. The w iterator should point to the first element
1032 /// of a vector of weights of the same size as the main array.
1033 template <typename Iterator, typename WeightIterator>
1034 Double_t TMath::Mean(Iterator first, Iterator last, WeightIterator w)
1035 {
1036 
1037  Double_t sum = 0;
1038  Double_t sumw = 0;
1039  int i = 0;
1040  while ( first != last ) {
1041  if ( *w < 0) {
1042  ::Error("TMath::Mean","w[%d] = %.4e < 0 ?!",i,*w);
1043  return 0;
1044  }
1045  sum += (*w) * (*first);
1046  sumw += (*w) ;
1047  ++w;
1048  ++first;
1049  ++i;
1050  }
1051  if (sumw <= 0) {
1052  ::Error("TMath::Mean","sum of weights == 0 ?!");
1053  return 0;
1054  }
1055 
1056  return sum/sumw;
1057 }
1058 
1059 ////////////////////////////////////////////////////////////////////////////////
1060 /// Return the weighted mean of an array a with length n.
1061 template <typename T>
1063 {
1064  if (w) {
1065  return TMath::Mean(a, a+n, w);
1066  } else {
1067  return TMath::Mean(a, a+n);
1068  }
1069 }
1070 
1071 ////////////////////////////////////////////////////////////////////////////////
1072 /// Return the geometric mean of an array defined by the iterators.
1073 /// \f[ GeomMean = (\prod_{i=0}^{n-1} |a[i]|)^{1/n} \f]
1074 template <typename Iterator>
1075 Double_t TMath::GeomMean(Iterator first, Iterator last)
1076 {
1077  Double_t logsum = 0.;
1078  Long64_t n = 0;
1079  while ( first != last ) {
1080  if (*first == 0) return 0.;
1081  Double_t absa = (Double_t) TMath::Abs(*first);
1082  logsum += TMath::Log(absa);
1083  ++first;
1084  ++n;
1085  }
1086 
1087  return TMath::Exp(logsum/n);
1088 }
1089 
1090 ////////////////////////////////////////////////////////////////////////////////
1091 /// Return the geometric mean of an array a of size n.
1092 /// \f[ GeomMean = (\prod_{i=0}^{n-1} |a[i]|)^{1/n} \f]
1093 template <typename T>
1095 {
1096  return TMath::GeomMean(a, a+n);
1097 }
1098 
1099 ////////////////////////////////////////////////////////////////////////////////
1100 /// Return the Standard Deviation of an array defined by the iterators.
1101 /// Note that this function returns the sigma(standard deviation) and
1102 /// not the root mean square of the array.
1103 ///
1104 /// Use the two pass algorithm, which is slower (! a factor of 2) but much more
1105 /// precise. Since we have a vector the 2 pass algorithm is still faster than the
1106 /// Welford algorithm. (See also ROOT-5545)
1107 template <typename Iterator>
1108 Double_t TMath::RMS(Iterator first, Iterator last)
1109 {
1110 
1111  Double_t n = 0;
1112 
1113  Double_t tot = 0;
1114  Double_t mean = TMath::Mean(first,last);
1115  while ( first != last ) {
1116  Double_t x = Double_t(*first);
1117  tot += (x - mean)*(x - mean);
1118  ++first;
1119  ++n;
1120  }
1121  Double_t rms = (n > 1) ? TMath::Sqrt(tot/(n-1)) : 0.0;
1122  return rms;
1123 }
1124 
1125 ////////////////////////////////////////////////////////////////////////////////
1126 /// Return the weighted Standard Deviation of an array defined by the iterators.
1127 /// Note that this function returns the sigma(standard deviation) and
1128 /// not the root mean square of the array.
1129 ///
1130 /// As in the unweighted case use the two pass algorithm
1131 template <typename Iterator, typename WeightIterator>
1132 Double_t TMath::RMS(Iterator first, Iterator last, WeightIterator w)
1133 {
1134  Double_t tot = 0;
1135  Double_t sumw = 0;
1136  Double_t sumw2 = 0;
1137  Double_t mean = TMath::Mean(first,last,w);
1138  while ( first != last ) {
1139  Double_t x = Double_t(*first);
1140  sumw += *w;
1141  sumw2 += (*w) * (*w);
1142  tot += (*w) * (x - mean)*(x - mean);
1143  ++first;
1144  ++w;
1145  }
1146  // use the correction neff/(neff -1) for the unbiased formula
1147  Double_t rms = TMath::Sqrt(tot * sumw/ (sumw*sumw - sumw2) );
1148  return rms;
1149 }
1150 
1151 ////////////////////////////////////////////////////////////////////////////////
1152 /// Return the Standard Deviation of an array a with length n.
1153 /// Note that this function returns the sigma(standard deviation) and
1154 /// not the root mean square of the array.
1155 template <typename T>
1156 Double_t TMath::RMS(Long64_t n, const T *a, const Double_t * w)
1157 {
1158  return (w) ? TMath::RMS(a, a+n, w) : TMath::RMS(a, a+n);
1159 }
1160 
1161 ////////////////////////////////////////////////////////////////////////////////
1162 /// Calculate the Cross Product of two vectors:
1163 /// out = [v1 x v2]
1164 template <typename T> T *TMath::Cross(const T v1[3],const T v2[3], T out[3])
1165 {
1166  out[0] = v1[1] * v2[2] - v1[2] * v2[1];
1167  out[1] = v1[2] * v2[0] - v1[0] * v2[2];
1168  out[2] = v1[0] * v2[1] - v1[1] * v2[0];
1169 
1170  return out;
1171 }
1172 
1173 ////////////////////////////////////////////////////////////////////////////////
1174 /// Calculate a normal vector of a plane.
1175 ///
1176 /// \param[in] p1, p2,p3 3 3D points belonged the plane to define it.
1177 /// \param[out] normal Pointer to 3D normal vector (normalized)
1178 template <typename T> T * TMath::Normal2Plane(const T p1[3],const T p2[3],const T p3[3], T normal[3])
1179 {
1180  T v1[3], v2[3];
1181 
1182  v1[0] = p2[0] - p1[0];
1183  v1[1] = p2[1] - p1[1];
1184  v1[2] = p2[2] - p1[2];
1185 
1186  v2[0] = p3[0] - p1[0];
1187  v2[1] = p3[1] - p1[1];
1188  v2[2] = p3[2] - p1[2];
1189 
1190  NormCross(v1,v2,normal);
1191  return normal;
1192 }
1193 
1194 ////////////////////////////////////////////////////////////////////////////////
1195 /// Function which returns kTRUE if point xp,yp lies inside the
1196 /// polygon defined by the np points in arrays x and y, kFALSE otherwise.
1197 /// Note that the polygon may be open or closed.
1198 template <typename T> Bool_t TMath::IsInside(T xp, T yp, Int_t np, T *x, T *y)
1199 {
1200  Int_t i, j = np-1 ;
1201  Bool_t oddNodes = kFALSE;
1202 
1203  for (i=0; i<np; i++) {
1204  if ((y[i]<yp && y[j]>=yp) || (y[j]<yp && y[i]>=yp)) {
1205  if (x[i]+(yp-y[i])/(y[j]-y[i])*(x[j]-x[i])<xp) {
1206  oddNodes = !oddNodes;
1207  }
1208  }
1209  j=i;
1210  }
1211 
1212  return oddNodes;
1213 }
1214 
1215 ////////////////////////////////////////////////////////////////////////////////
1216 /// Return the median of the array a where each entry i has weight w[i] .
1217 /// Both arrays have a length of at least n . The median is a number obtained
1218 /// from the sorted array a through
1219 ///
1220 /// median = (a[jl]+a[jh])/2. where (using also the sorted index on the array w)
1221 ///
1222 /// sum_i=0,jl w[i] <= sumTot/2
1223 /// sum_i=0,jh w[i] >= sumTot/2
1224 /// sumTot = sum_i=0,n w[i]
1225 ///
1226 /// If w=0, the algorithm defaults to the median definition where it is
1227 /// a number that divides the sorted sequence into 2 halves.
1228 /// When n is odd or n > 1000, the median is kth element k = (n + 1) / 2.
1229 /// when n is even and n < 1000the median is a mean of the elements k = n/2 and k = n/2 + 1.
1230 ///
1231 /// If the weights are supplied (w not 0) all weights must be >= 0
1232 ///
1233 /// If work is supplied, it is used to store the sorting index and assumed to be
1234 /// >= n . If work=0, local storage is used, either on the stack if n < kWorkMax
1235 /// or on the heap for n >= kWorkMax .
1236 template <typename T> Double_t TMath::Median(Long64_t n, const T *a, const Double_t *w, Long64_t *work)
1237 {
1238 
1239  const Int_t kWorkMax = 100;
1240 
1241  if (n <= 0 || !a) return 0;
1242  Bool_t isAllocated = kFALSE;
1243  Double_t median;
1244  Long64_t *ind;
1245  Long64_t workLocal[kWorkMax];
1246 
1247  if (work) {
1248  ind = work;
1249  } else {
1250  ind = workLocal;
1251  if (n > kWorkMax) {
1252  isAllocated = kTRUE;
1253  ind = new Long64_t[n];
1254  }
1255  }
1256 
1257  if (w) {
1258  Double_t sumTot2 = 0;
1259  for (Int_t j = 0; j < n; j++) {
1260  if (w[j] < 0) {
1261  ::Error("TMath::Median","w[%d] = %.4e < 0 ?!",j,w[j]);
1262  if (isAllocated) delete [] ind;
1263  return 0;
1264  }
1265  sumTot2 += w[j];
1266  }
1267 
1268  sumTot2 /= 2.;
1269 
1270  Sort(n, a, ind, kFALSE);
1271 
1272  Double_t sum = 0.;
1273  Int_t jl;
1274  for (jl = 0; jl < n; jl++) {
1275  sum += w[ind[jl]];
1276  if (sum >= sumTot2) break;
1277  }
1278 
1279  Int_t jh;
1280  sum = 2.*sumTot2;
1281  for (jh = n-1; jh >= 0; jh--) {
1282  sum -= w[ind[jh]];
1283  if (sum <= sumTot2) break;
1284  }
1285 
1286  median = 0.5*(a[ind[jl]]+a[ind[jh]]);
1287 
1288  } else {
1289 
1290  if (n%2 == 1)
1291  median = KOrdStat(n, a,n/2, ind);
1292  else {
1293  median = 0.5*(KOrdStat(n, a, n/2 -1, ind)+KOrdStat(n, a, n/2, ind));
1294  }
1295  }
1296 
1297  if (isAllocated)
1298  delete [] ind;
1299  return median;
1300 }
1301 
1302 ////////////////////////////////////////////////////////////////////////////////
1303 /// Returns k_th order statistic of the array a of size n
1304 /// (k_th smallest element out of n elements).
1305 ///
1306 /// C-convention is used for array indexing, so if you want
1307 /// the second smallest element, call KOrdStat(n, a, 1).
1308 ///
1309 /// If work is supplied, it is used to store the sorting index and
1310 /// assumed to be >= n. If work=0, local storage is used, either on
1311 /// the stack if n < kWorkMax or on the heap for n >= kWorkMax.
1312 /// Note that the work index array will not contain the sorted indices but
1313 /// all indices of the smaller element in arbitrary order in work[0,...,k-1] and
1314 /// all indices of the larger element in arbitrary order in work[k+1,..,n-1]
1315 /// work[k] will contain instead the index of the returned element.
1316 ///
1317 /// Taken from "Numerical Recipes in C++" without the index array
1318 /// implemented by Anna Khreshuk.
1319 ///
1320 /// See also the declarations at the top of this file
1321 template <class Element, typename Size>
1322 Element TMath::KOrdStat(Size n, const Element *a, Size k, Size *work)
1323 {
1324 
1325  const Int_t kWorkMax = 100;
1326 
1327  typedef Size Index;
1328 
1329  Bool_t isAllocated = kFALSE;
1330  Size i, ir, j, l, mid;
1331  Index arr;
1332  Index *ind;
1333  Index workLocal[kWorkMax];
1334  Index temp;
1335 
1336  if (work) {
1337  ind = work;
1338  } else {
1339  ind = workLocal;
1340  if (n > kWorkMax) {
1341  isAllocated = kTRUE;
1342  ind = new Index[n];
1343  }
1344  }
1345 
1346  for (Size ii=0; ii<n; ii++) {
1347  ind[ii]=ii;
1348  }
1349  Size rk = k;
1350  l=0;
1351  ir = n-1;
1352  for(;;) {
1353  if (ir<=l+1) { //active partition contains 1 or 2 elements
1354  if (ir == l+1 && a[ind[ir]]<a[ind[l]])
1355  {temp = ind[l]; ind[l]=ind[ir]; ind[ir]=temp;}
1356  Element tmp = a[ind[rk]];
1357  if (isAllocated)
1358  delete [] ind;
1359  return tmp;
1360  } else {
1361  mid = (l+ir) >> 1; //choose median of left, center and right
1362  {temp = ind[mid]; ind[mid]=ind[l+1]; ind[l+1]=temp;}//elements as partitioning element arr.
1363  if (a[ind[l]]>a[ind[ir]]) //also rearrange so that a[l]<=a[l+1]
1364  {temp = ind[l]; ind[l]=ind[ir]; ind[ir]=temp;}
1365 
1366  if (a[ind[l+1]]>a[ind[ir]])
1367  {temp=ind[l+1]; ind[l+1]=ind[ir]; ind[ir]=temp;}
1368 
1369  if (a[ind[l]]>a[ind[l+1]])
1370  {temp = ind[l]; ind[l]=ind[l+1]; ind[l+1]=temp;}
1371 
1372  i=l+1; //initialize pointers for partitioning
1373  j=ir;
1374  arr = ind[l+1];
1375  for (;;){
1376  do i++; while (a[ind[i]]<a[arr]);
1377  do j--; while (a[ind[j]]>a[arr]);
1378  if (j<i) break; //pointers crossed, partitioning complete
1379  {temp=ind[i]; ind[i]=ind[j]; ind[j]=temp;}
1380  }
1381  ind[l+1]=ind[j];
1382  ind[j]=arr;
1383  if (j>=rk) ir = j-1; //keep active the partition that
1384  if (j<=rk) l=i; //contains the k_th element
1385  }
1386  }
1387 }
1388 
1389 #endif
c
#define c(i)
Definition: RSha256.hxx:119
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:423
TMath::MaxElement
T MaxElement(Long64_t n, const T *a)
Return maximum of array a of length n.
Definition: TMath.h:948
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:1062
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:668
TMath::Kcgs
constexpr Double_t Kcgs()
Definition: TMath.h:259
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:1094
TMath::RadToDeg
constexpr Double_t RadToDeg()
Conversion from radian to degree:
Definition: TMath.h:79
TMath::NaUncertainty
constexpr Double_t NaUncertainty()
Avogadro constant (Avogadro's Number) uncertainty.
Definition: TMath.h:295
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:890
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:630
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:523
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:116
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:749
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:165
TMath::GUncertainty
constexpr Double_t GUncertainty()
Gravitational constant uncertainty.
Definition: TMath.h:158
TMath::Sqrt
Double_t Sqrt(Double_t x)
Definition: TMath.h:680
exp
double exp(double)
TMath::CeilNint
Int_t CeilNint(Double_t x)
Definition: TMath.h:688
TMath::DegToRad
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition: TMath.h:87
TMath::Tan
Double_t Tan(Double_t)
Definition: TMath.h:634
Float_t
float Float_t
Definition: RtypesCore.h:57
log
double log(double)
TMath::Exp
Double_t Exp(Double_t x)
Definition: TMath.h:716
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:168
tanh
double tanh(double)
TMath::PiOver2
constexpr Double_t PiOver2()
Definition: TMath.h:57
Int_t
int Int_t
Definition: RtypesCore.h:45
TMath::Qe
constexpr Double_t Qe()
Elementary charge in .
Definition: TMath.h:341
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:71
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:266
TMath::RUncertainty
constexpr Double_t RUncertainty()
Universal gas constant uncertainty.
Definition: TMath.h:311
TMath::Gn
constexpr Double_t Gn()
Standard acceleration of gravity in .
Definition: TMath.h:179
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:123
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:1198
TMath::LocMax
Long64_t LocMax(Long64_t n, const T *a)
Return index of array with the maximum element.
Definition: TMath.h:989
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:118
TMath::HUncertainty
constexpr Double_t HUncertainty()
Planck's constant uncertainty.
Definition: TMath.h:208
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:1164
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:696
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:274
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:961
TMath::Pi
constexpr Double_t Pi()
Definition: TMath.h:43
TMath::Log2
Double_t Log2(Double_t x)
Definition: TMath.cxx:101
Option_t
const typedef char Option_t
Definition: RtypesCore.h:66
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:702
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:753
TMath::MWair
constexpr Double_t MWair()
Molecular weight of dry air 1976 US Standard Atmosphere in or
Definition: TMath.h:319
TMath::Limits::Min
static T Min()
Returns maximum representation for type T.
Definition: TMath.h:910
xmin
float xmin
Definition: THbookFile.cxx:95
TMath::IsNaN
Bool_t IsNaN(Double_t x)
Definition: TMath.h:881
epsilon
REAL epsilon
Definition: triangle.c:617
TMath::Floor
Double_t Floor(Double_t x)
Definition: TMath.h:692
a
auto * a
Definition: textangle.C:12
TMath::CUncertainty
constexpr Double_t CUncertainty()
Speed of light uncertainty.
Definition: TMath.h:137
TMath::MinElement
T MinElement(Long64_t n, const T *a)
Return minimum of array a of length n.
Definition: TMath.h:941
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TMath::Limits
Definition: TMath.h:407
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:1178
Long_t
long Long_t
Definition: RtypesCore.h:54
TMath::Power
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:724
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:151
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:64
TMath::Sin
Double_t Sin(Double_t)
Definition: TMath.h:626
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:1156
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:903
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:144
TMath::Rgair
constexpr Double_t Rgair()
Dry Air Gas Constant (R / MWair) in
Definition: TMath.h:327
TMath::SinH
Double_t SinH(Double_t)
Definition: TMath.h:638
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:897
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:201
TMath::Limits::Max
static T Max()
Returns minimum double representation.
Definition: TMath.h:917
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:417
TMath::IsNaN
Bool_t IsNaN(Float_t x)
Definition: TMath.h:882
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:50
TMath::R
constexpr Double_t R()
Universal gas constant ( ) in
Definition: TMath.h:304
TMath::CosH
Double_t CosH(Double_t)
Definition: TMath.h:642
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:172
TMath::Ccgs
constexpr Double_t Ccgs()
Definition: TMath.h:130
xmlio::Size
const char * Size
Definition: TXMLSetup.cxx:62
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:186
TMath::ACos
Double_t ACos(Double_t)
Definition: TMath.h:657
TMath::ATan
Double_t ATan(Double_t)
Definition: TMath.h:664
TMath::QeUncertainty
constexpr Double_t QeUncertainty()
Elementary charge uncertainty.
Definition: TMath.h:348
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:720
TMath::H
constexpr Double_t H()
Planck's constant in .
Definition: TMath.h:194
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:109
TMath::Limits::Epsilon
static T Epsilon()
Returns minimum double representation.
Definition: TMath.h:924
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:676
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:933
TMath::TanH
Double_t TanH(Double_t)
Definition: TMath.h:646
TMath::EulerGamma
constexpr Double_t EulerGamma()
Euler-Mascheroni Constant.
Definition: TMath.h:334
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:94
ROOT::Math::Chebyshev::T
double T(double x)
Definition: ChebyshevPol.h:52
d
#define d(i)
Definition: RSha256.hxx:120
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:1322
cosh
double cosh(double)
TMathBase.h
TMath::HC
constexpr Double_t HC()
in
Definition: TMath.h:237
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:252
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:244
pow
double pow(double, double)
asin
double asin(double)
TMath::SigmaUncertainty
constexpr Double_t SigmaUncertainty()
Stefan-Boltzmann constant uncertainty.
Definition: TMath.h:281
TMath::Na
constexpr Double_t Na()
Avogadro constant (Avogadro's Number) in .
Definition: TMath.h:288
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:102
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:760
Math
TMath::Ceil
Double_t Ceil(Double_t x)
Definition: TMath.h:684
TMath::ASin
Double_t ASin(Double_t)
Definition: TMath.h:650
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:188
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:1236
TError.h