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