35#if !defined(R__SOLARIS) && !defined(R__ACC) && !defined(R__FBSD)
67 if(
x==0.0)
return 0.0;
69 return log(
x+ax*sqrt(1.+1./(ax*ax)));
80 if(
x==0.0)
return 0.0;
82 return log(
x+ax*sqrt(1.-1./(ax*ax)));
93 return log((1+
x)/(1-
x))/2;
103 return log(
x)/log(2.0);
118 const Double_t c[20] = {0.42996693560813697, 0.40975987533077106,
119 -0.01858843665014592, 0.00145751084062268,-0.00014304184442340,
120 0.00001588415541880,-0.00000190784959387, 0.00000024195180854,
121 -0.00000003193341274, 0.00000000434545063,-0.00000000060578480,
122 0.00000000008612098,-0.00000000001244332, 0.00000000000182256,
123 -0.00000000000027007, 0.00000000000004042,-0.00000000000000610,
124 0.00000000000000093,-0.00000000000000014, 0.00000000000000002};
127 t=
h=
y=s=
a=alfa=b1=b2=b0=0.;
131 }
else if (
x == -1) {
140 a = -pi3+hf*(b1*b1-b2*b2);
146 }
else if (t <= -0.5) {
170 for (
Int_t i=19;i>=0;i--){
171 b0 =
c[i] + alfa*b1-b2;
175 h = -(s*(b0-
h*b2)+
a);
186 return ::ROOT::Math::erf(
x);
196 return ::ROOT::Math::erfc(
x);
207 Double_t kConst = 0.8862269254527579;
212 Double_t erfi, derfi, y0,y1,dy0,dy1;
217 for (
Int_t iter=0; iter<kMaxit; iter++) {
220 if (
TMath::Abs(dy1) < kEps) {
if (
x < 0)
return -erfi;
else return erfi;}
225 if(
TMath::Abs(derfi/erfi) < kEps) {
if (
x < 0)
return -erfi;
else return erfi;}
249 if (
n <= 0)
return 1.;
268 const Double_t w2 = 1.41421356237309505;
270 const Double_t p10 = 2.4266795523053175e+2, q10 = 2.1505887586986120e+2,
271 p11 = 2.1979261618294152e+1, q11 = 9.1164905404514901e+1,
272 p12 = 6.9963834886191355e+0, q12 = 1.5082797630407787e+1,
273 p13 =-3.5609843701815385e-2, q13 = 1;
275 const Double_t p20 = 3.00459261020161601e+2, q20 = 3.00459260956983293e+2,
276 p21 = 4.51918953711872942e+2, q21 = 7.90950925327898027e+2,
277 p22 = 3.39320816734343687e+2, q22 = 9.31354094850609621e+2,
278 p23 = 1.52989285046940404e+2, q23 = 6.38980264465631167e+2,
279 p24 = 4.31622272220567353e+1, q24 = 2.77585444743987643e+2,
280 p25 = 7.21175825088309366e+0, q25 = 7.70001529352294730e+1,
281 p26 = 5.64195517478973971e-1, q26 = 1.27827273196294235e+1,
282 p27 =-1.36864857382716707e-7, q27 = 1;
284 const Double_t p30 =-2.99610707703542174e-3, q30 = 1.06209230528467918e-2,
285 p31 =-4.94730910623250734e-2, q31 = 1.91308926107829841e-1,
286 p32 =-2.26956593539686930e-1, q32 = 1.05167510706793207e+0,
287 p33 =-2.78661308609647788e-1, q33 = 1.98733201817135256e+0,
288 p34 =-2.23192459734184686e-2, q34 = 1;
339 if (
x > 0)
return 0.5 +0.5*
h;
350 return ::ROOT::Math::tgamma(z);
366 return ::ROOT::Math::inc_gamma(
a,
x);
381 if (
a <= 0 ||
x <= 0)
return 0;
389 for (
Int_t i=1; i<=itmax; i++) {
393 if (
Abs(
d) < fpmin)
d = fpmin;
395 if (
Abs(
c) < fpmin)
c = fpmin;
399 if (
Abs(del-1) < eps)
break;
417 if (
a <= 0 ||
x <= 0)
return 0;
439 Double_t bw = gamma/((
x-mean)*(
x-mean) + gamma*gamma/4);
450 if (
sigma == 0)
return 1.e30;
453 if (arg < -39.0 || arg > 39.0)
return 0.0;
455 if (!norm)
return res;
456 return res/(2.50662827463100024*
sigma);
471 if (
sigma <= 0)
return 0;
473 if (!norm)
return den;
488 return ::ROOT::Math::lgamma(z);
520 if( av0 >= av1 && av0 >= av2 ) {
526 else if (av1 >= av0 && av1 >= av2) {
541 Double_t foofrac = foo/amax, barfrac = bar/amax;
542 Double_t d = amax * Sqrt(1.+foofrac*foofrac+barfrac*barfrac);
595 return Poisson(ix,par);
616 if (ndf <= 0)
return 0;
619 if (chi2 < 0)
return 0;
623 return ::ROOT::Math::chisquared_cdf_c(chi2,ndf);
669 }
else if (u < 0.755) {
672 }
else if (u < 6.8116) {
678 for (
Int_t j=0; j<maxj;j++) {
681 p = 2*(
r[0] -
r[1] +
r[2] -
r[3]);
791 if (!
a || !
b || na <= 2 || nb <= 2) {
792 Error(
"KolmogorovTest",
"Sets must have more than 2 points");
808 for (
Int_t i=0;i<na+nb;i++) {
812 if (ia >= na) {ok =
kTRUE;
break;}
813 }
else if (
a[ia] >
b[ib]) {
816 if (ib >= nb) {ok =
kTRUE;
break;}
820 while(ia < na &&
a[ia] ==
x) {
824 while(ib < nb &&
b[ib] ==
x) {
828 if (ia >= na) {ok =
kTRUE;
break;}
829 if (ib >= nb) {ok =
kTRUE;
break;}
843 printf(
" Kolmogorov Probability = %g, Max Dist = %g\n",prob,rdmax);
877 if ((
sigma < 0 || lg < 0) || (
sigma==0 && lg==0)) {
882 return lg * 0.159154943 / (xx*xx + lg*lg /4);
890 x = xx /
sigma / 1.41421356;
891 y = lg / 2 /
sigma / 1.41421356;
910 Double_t c[6] = { 1.0117281, -0.75197147, 0.012557727, 0.010022008, -0.00024206814, 0.00000050084806};
911 Double_t s[6] = { 1.393237, 0.23115241, -0.15535147, 0.0062183662, 0.000091908299, -0.00000062752596};
912 Double_t t[6] = { 0.31424038, 0.94778839, 1.5976826, 2.2795071, 3.0206370, 3.8897249};
919 Double_t xlim0, xlim1, xlim2, xlim3, xlim4;
920 Double_t a0=0, d0=0, d2=0, e0=0, e2=0, e4=0, h0=0, h2=0, h4=0, h6=0;
921 Double_t p0=0, p2=0, p4=0, p6=0, p8=0, z0=0, z2=0, z4=0, z6=0, z8=0;
922 Double_t xp[6], xm[6], yp[6], ym[6];
923 Double_t mq[6], pq[6], mf[6], pf[6];
938 xlim3 = 3.097 *
y - 0.45;
940 xlim4 = 18.1 *
y + 1.65;
949 k = yrrtpi / (xq + yq);
950 }
else if ( abx > xlim1 ) {
957 d = rrtpi / (d0 + xq*(d2 + xq));
958 k =
d *
y * (a0 + xq);
959 }
else if ( abx > xlim2 ) {
962 h0 = 0.5625 + yq * (4.5 + yq * (10.5 + yq * (6.0 + yq)));
964 h2 = -4.5 + yq * (9.0 + yq * ( 6.0 + yq * 4.0));
965 h4 = 10.5 - yq * (6.0 - yq * 6.0);
966 h6 = -6.0 + yq * 4.0;
967 e0 = 1.875 + yq * (8.25 + yq * (5.5 + yq));
968 e2 = 5.25 + yq * (1.0 + yq * 3.0);
971 d = rrtpi / (h0 + xq * (h2 + xq * (h4 + xq * (h6 + xq))));
972 k =
d *
y * (e0 + xq * (e2 + xq * (e4 + xq)));
973 }
else if ( abx < xlim3 ) {
976 z0 = 272.1014 +
y * (1280.829 +
y *
986 z2 = 211.678 +
y * (902.3066 +
y *
994 z4 = 78.86585 +
y * (308.1852 +
y *
998 (80.39278 +
y * 10.0)
1000 z6 = 22.03523 +
y * (55.02933 +
y *
1002 (53.59518 +
y * 10.0)
1004 z8 = 1.496460 +
y * (13.39880 +
y * 5.0);
1005 p0 = 153.5168 +
y * (549.3954 +
y *
1012 (4.264678 +
y * 0.3183291)
1014 p2 = -34.16955 +
y * (-1.322256+
y *
1019 (12.79458 +
y * 1.2733163)
1021 p4 = 2.584042 +
y * (10.46332 +
y *
1024 (12.79568 +
y * 1.9099744)
1026 p6 = -0.07272979 +
y * (0.9377051 +
y *
1027 (4.266322 +
y * 1.273316));
1028 p8 = 0.0005480304 +
y * 0.3183291;
1030 d = 1.7724538 / (z0 + xq * (z2 + xq * (z4 + xq * (z6 + xq * (z8 + xq)))));
1031 k =
d * (p0 + xq * (p2 + xq * (p4 + xq * (p6 + xq * p8))));
1034 ypy0q = ypy0 * ypy0;
1036 for (j = 0; j <= 5; j++) {
1039 mf[j] = 1.0 / (mq[j] + ypy0q);
1041 ym[j] = mf[j] * ypy0;
1044 pf[j] = 1.0 / (pq[j] + ypy0q);
1046 yp[j] = pf[j] * ypy0;
1048 if ( abx <= xlim4 ) {
1049 for (j = 0; j <= 5; j++) {
1050 k = k +
c[j]*(ym[j]+yp[j]) - s[j]*(xm[j]-xp[j]) ;
1054 for ( j = 0; j <= 5; j++) {
1056 (mq[j] * mf[j] - y0 * ym[j])
1057 + s[j] * yf * xm[j]) / (mq[j]+y0q)
1058 + (
c[j] * (pq[j] * pf[j] - y0 * yp[j])
1059 - s[j] * yf * xp[j]) / (pq[j]+y0q);
1061 k =
y * k + exp( -xq );
1064 return k / 2.506628 /
sigma;
1087 Double_t r,s,t,p,
q,
d,ps3,ps33,qs2,u,
v,tmp,lnu,lnv,su,sv,y1,y2,y3;
1091 if (coef[3] == 0)
return complex;
1092 r = coef[2]/coef[3];
1093 s = coef[1]/coef[3];
1094 t = coef[0]/coef[3];
1097 q = (2*
r*
r*
r)/27.0 - (
r*s)/3 + t;
1186 if (type<1 || type>9){
1187 printf(
"illegal value of type\n");
1193 if (index) ind = index;
1196 isAllocated =
kTRUE;
1205 for (
Int_t i=0; i<nprob; i++){
1215 nppm =
n*prob[i]-0.5;
1239 if (
type == 4) {
a = 0;
b = 1; }
1240 else if (
type == 5) {
a = 0.5;
b = 0.5; }
1241 else if (
type == 6) {
a = 0.;
b = 0.; }
1242 else if (
type == 7) {
a = 1.;
b = 1.; }
1243 else if (
type == 8) {
a = 1./3.;
b =
a; }
1244 else if (
type == 9) {
a = 3./8.;
b =
a; }
1248 nppm =
a + prob[i] * (
n + 1 -
a -
b);
1251 if (gamma < eps) gamma = 0;
1256 int first = (j > 0 && j <=
n) ? j-1 : ( j <=0 ) ? 0 :
n-1;
1257 int second = (j > 0 && j <
n) ? j : ( j <=0 ) ? 0 :
n-1;
1267 quantiles[i] = (1-gamma)*xj + gamma*xjj;
1291 if (Narr <= 0)
return;
1292 double *localArr1 =
new double[Narr];
1293 int *localArr2 =
new int[Narr];
1297 for(iEl = 0; iEl < Narr; iEl++) {
1298 localArr1[iEl] = arr1[iEl];
1299 localArr2[iEl] = iEl;
1302 for (iEl = 0; iEl < Narr; iEl++) {
1303 for (iEl2 = Narr-1; iEl2 > iEl; --iEl2) {
1304 if (localArr1[iEl2-1] < localArr1[iEl2]) {
1305 double tmp = localArr1[iEl2-1];
1306 localArr1[iEl2-1] = localArr1[iEl2];
1307 localArr1[iEl2] = tmp;
1309 int tmp2 = localArr2[iEl2-1];
1310 localArr2[iEl2-1] = localArr2[iEl2];
1311 localArr2[iEl2] = tmp2;
1316 for (iEl = 0; iEl < Narr; iEl++) {
1317 arr2[iEl] = localArr2[iEl];
1319 delete [] localArr2;
1320 delete [] localArr1;
1330 if (Narr <= 0)
return;
1331 double *localArr1 =
new double[Narr];
1332 int *localArr2 =
new int[Narr];
1336 for (iEl = 0; iEl < Narr; iEl++) {
1337 localArr1[iEl] = arr1[iEl];
1338 localArr2[iEl] = iEl;
1341 for (iEl = 0; iEl < Narr; iEl++) {
1342 for (iEl2 = Narr-1; iEl2 > iEl; --iEl2) {
1343 if (localArr1[iEl2-1] > localArr1[iEl2]) {
1344 double tmp = localArr1[iEl2-1];
1345 localArr1[iEl2-1] = localArr1[iEl2];
1346 localArr1[iEl2] = tmp;
1348 int tmp2 = localArr2[iEl2-1];
1349 localArr2[iEl2-1] = localArr2[iEl2];
1350 localArr2[iEl2] = tmp2;
1355 for (iEl = 0; iEl < Narr; iEl++) {
1356 arr2[iEl] = localArr2[iEl];
1358 delete [] localArr2;
1359 delete [] localArr1;
1404 const Double_t p1=1.0, p2=3.5156229, p3=3.0899424,
1405 p4=1.2067492, p5=0.2659732, p6=3.60768e-2, p7=4.5813e-3;
1407 const Double_t q1= 0.39894228, q2= 1.328592e-2, q3= 2.25319e-3,
1408 q4=-1.57565e-3, q5= 9.16281e-3, q6=-2.057706e-2,
1409 q7= 2.635537e-2, q8=-1.647633e-2, q9= 3.92377e-3;
1419 result = p1+
y*(p2+
y*(p3+
y*(p4+
y*(p5+
y*(p6+
y*p7)))));
1438 const Double_t p1=-0.57721566, p2=0.42278420, p3=0.23069756,
1439 p4= 3.488590e-2, p5=2.62698e-3, p6=1.0750e-4, p7=7.4e-6;
1441 const Double_t q1= 1.25331414, q2=-7.832358e-2, q3= 2.189568e-2,
1442 q4=-1.062446e-2, q5= 5.87872e-3, q6=-2.51540e-3, q7=5.3208e-4;
1445 Error(
"TMath::BesselK0",
"*K0* Invalid argument x = %g\n",
x);
1456 result = (exp(-
x)/sqrt(
x))*(q1+
y*(q2+
y*(q3+
y*(q4+
y*(q5+
y*(q6+
y*q7))))));
1472 const Double_t p1=0.5, p2=0.87890594, p3=0.51498869,
1473 p4=0.15084934, p5=2.658733e-2, p6=3.01532e-3, p7=3.2411e-4;
1475 const Double_t q1= 0.39894228, q2=-3.988024e-2, q3=-3.62018e-3,
1476 q4= 1.63801e-3, q5=-1.031555e-2, q6= 2.282967e-2,
1477 q7=-2.895312e-2, q8= 1.787654e-2, q9=-4.20059e-3;
1487 result =
x*(p1+
y*(p2+
y*(p3+
y*(p4+
y*(p5+
y*(p6+
y*p7))))));
1490 result = (exp(ax)/sqrt(ax))*(q1+
y*(q2+
y*(q3+
y*(q4+
y*(q5+
y*(q6+
y*(q7+
y*(q8+
y*q9))))))));
1491 if (
x < 0) result = -result;
1507 const Double_t p1= 1., p2= 0.15443144, p3=-0.67278579,
1508 p4=-0.18156897, p5=-1.919402e-2, p6=-1.10404e-3, p7=-4.686e-5;
1510 const Double_t q1= 1.25331414, q2= 0.23498619, q3=-3.655620e-2,
1511 q4= 1.504268e-2, q5=-7.80353e-3, q6= 3.25614e-3, q7=-6.8245e-4;
1514 Error(
"TMath::BesselK1",
"*K1* Invalid argument x = %g\n",
x);
1525 result = (exp(-
x)/sqrt(
x))*(q1+
y*(q2+
y*(q3+
y*(q4+
y*(q5+
y*(q6+
y*q7))))));
1538 if (
x <= 0 ||
n < 0) {
1539 Error(
"TMath::BesselK",
"*K* Invalid argument(s) (n,x) = (%d, %g)\n",
n,
x);
1551 for (
Int_t j=1; j<
n; j++) {
1568 const Double_t kBigPositive = 1.e10;
1569 const Double_t kBigNegative = 1.e-10;
1572 Error(
"TMath::BesselI",
"*I* Invalid argument (n,x) = (%d, %g)\n",
n,
x);
1579 if (
x == 0)
return 0;
1587 for (
Int_t j=
m; j>=1; j--) {
1593 result *= kBigNegative;
1595 bip *= kBigNegative;
1597 if (j==
n) result=bip;
1601 if ((
x < 0) && (
n%2 == 1)) result = -result;
1613 const Double_t p1 = 57568490574.0, p2 = -13362590354.0, p3 = 651619640.7;
1614 const Double_t p4 = -11214424.18, p5 = 77392.33017, p6 = -184.9052456;
1615 const Double_t p7 = 57568490411.0, p8 = 1029532985.0, p9 = 9494680.718;
1616 const Double_t p10 = 59272.64853, p11 = 267.8532712;
1619 const Double_t q2 = -0.1098628627e-2, q3 = 0.2734510407e-4;
1620 const Double_t q4 = -0.2073370639e-5, q5 = 0.2093887211e-6;
1621 const Double_t q6 = -0.1562499995e-1, q7 = 0.1430488765e-3;
1622 const Double_t q8 = -0.6911147651e-5, q9 = 0.7621095161e-6;
1623 const Double_t q10 = 0.934935152e-7, q11 = 0.636619772;
1625 if ((ax=fabs(
x)) < 8) {
1627 result1 = p1 +
y*(p2 +
y*(p3 +
y*(p4 +
y*(p5 +
y*p6))));
1628 result2 = p7 +
y*(p8 +
y*(p9 +
y*(p10 +
y*(p11 +
y))));
1629 result = result1/result2;
1634 result1 = 1 +
y*(q2 +
y*(q3 +
y*(q4 +
y*q5)));
1635 result2 = q6 +
y*(q7 +
y*(q8 +
y*(q9 -
y*q10)));
1636 result = sqrt(q11/ax)*(cos(xx)*result1-z*sin(xx)*result2);
1648 const Double_t p1 = 72362614232.0, p2 = -7895059235.0, p3 = 242396853.1;
1649 const Double_t p4 = -2972611.439, p5 = 15704.48260, p6 = -30.16036606;
1650 const Double_t p7 = 144725228442.0, p8 = 2300535178.0, p9 = 18583304.74;
1651 const Double_t p10 = 99447.43394, p11 = 376.9991397;
1654 const Double_t q2 = 0.183105e-2, q3 = -0.3516396496e-4;
1655 const Double_t q4 = 0.2457520174e-5, q5 = -0.240337019e-6;
1656 const Double_t q6 = 0.04687499995, q7 = -0.2002690873e-3;
1657 const Double_t q8 = 0.8449199096e-5, q9 = -0.88228987e-6;
1658 const Double_t q10 = 0.105787412e-6, q11 = 0.636619772;
1660 if ((ax=fabs(
x)) < 8) {
1662 result1 =
x*(p1 +
y*(p2 +
y*(p3 +
y*(p4 +
y*(p5 +
y*p6)))));
1663 result2 = p7 +
y*(p8 +
y*(p9 +
y*(p10 +
y*(p11 +
y))));
1664 result = result1/result2;
1669 result1 = 1 +
y*(q2 +
y*(q3 +
y*(q4 +
y*q5)));
1670 result2 = q6 +
y*(q7 +
y*(q8 +
y*(q9 +
y*q10)));
1671 result = sqrt(q11/ax)*(cos(xx)*result1-z*sin(xx)*result2);
1672 if (
x < 0) result = -result;
1683 const Double_t p1 = -2957821389., p2 = 7062834065.0, p3 = -512359803.6;
1684 const Double_t p4 = 10879881.29, p5 = -86327.92757, p6 = 228.4622733;
1685 const Double_t p7 = 40076544269., p8 = 745249964.8, p9 = 7189466.438;
1686 const Double_t p10 = 47447.26470, p11 = 226.1030244, p12 = 0.636619772;
1689 const Double_t q2 = -0.1098628627e-2, q3 = 0.2734510407e-4;
1690 const Double_t q4 = -0.2073370639e-5, q5 = 0.2093887211e-6;
1691 const Double_t q6 = -0.1562499995e-1, q7 = 0.1430488765e-3;
1692 const Double_t q8 = -0.6911147651e-5, q9 = 0.7621095161e-6;
1693 const Double_t q10 = -0.934945152e-7, q11 = 0.636619772;
1697 result1 = p1 +
y*(p2 +
y*(p3 +
y*(p4 +
y*(p5 +
y*p6))));
1698 result2 = p7 +
y*(p8 +
y*(p9 +
y*(p10 +
y*(p11 +
y))));
1704 result1 = 1 +
y*(q2 +
y*(q3 +
y*(q4 +
y*q5)));
1705 result2 = q6 +
y*(q7 +
y*(q8 +
y*(q9 +
y*q10)));
1706 result = sqrt(q11/
x)*(sin(xx)*result1+z*cos(xx)*result2);
1717 const Double_t p1 = -0.4900604943e13, p2 = 0.1275274390e13;
1718 const Double_t p3 = -0.5153438139e11, p4 = 0.7349264551e9;
1719 const Double_t p5 = -0.4237922726e7, p6 = 0.8511937935e4;
1720 const Double_t p7 = 0.2499580570e14, p8 = 0.4244419664e12;
1721 const Double_t p9 = 0.3733650367e10, p10 = 0.2245904002e8;
1722 const Double_t p11 = 0.1020426050e6, p12 = 0.3549632885e3;
1725 const Double_t q2 = 0.183105e-2, q3 = -0.3516396496e-4;
1726 const Double_t q4 = 0.2457520174e-5, q5 = -0.240337019e-6;
1727 const Double_t q6 = 0.04687499995, q7 = -0.2002690873e-3;
1728 const Double_t q8 = 0.8449199096e-5, q9 = -0.88228987e-6;
1729 const Double_t q10 = 0.105787412e-6, q11 = 0.636619772;
1733 result1 =
x*(p1 +
y*(p2 +
y*(p3 +
y*(p4 +
y*(p5 +
y*p6)))));
1734 result2 = p7 +
y*(p8 +
y*(p9 +
y*(p10 +
y*(p11 +
y*(p12+
y)))));
1740 result1 = 1 +
y*(q2 +
y*(q3 +
y*(q4 +
y*q5)));
1741 result2 = q6 +
y*(q7 +
y*(q8 +
y*(q9 +
y*q10)));
1742 result = sqrt(q11/
x)*(sin(xx)*result1+z*cos(xx)*result2);
1754 const Int_t n1 = 15;
1755 const Int_t n2 = 25;
1756 const Double_t c1[16] = { 1.00215845609911981, -1.63969292681309147,
1757 1.50236939618292819, -.72485115302121872,
1758 .18955327371093136, -.03067052022988,
1759 .00337561447375194, -2.6965014312602e-4,
1760 1.637461692612e-5, -7.8244408508e-7,
1761 3.021593188e-8, -9.6326645e-10,
1762 2.579337e-11, -5.8854e-13,
1763 1.158e-14, -2
e-16 };
1764 const Double_t c2[26] = { .99283727576423943, -.00696891281138625,
1765 1.8205103787037e-4, -1.063258252844e-5,
1766 9.8198294287e-7, -1.2250645445e-7,
1767 1.894083312e-8, -3.44358226e-9,
1768 7.1119102e-10, -1.6288744e-10,
1769 4.065681e-11, -1.091505e-11,
1770 3.12005e-12, -9.4202e-13,
1771 2.9848e-13, -9.872e-14,
1772 3.394e-14, -1.208e-14,
1773 4.44e-15, -1.68e-15,
1792 for (i = n1; i >= 0; --i) {
1793 b0 =
c1[i] + alfa*b1 - b2;
1805 for (i = n2; i >= 0; --i) {
1806 b0 =
c2[i] + alfa*b1 - b2;
1823 const Int_t n3 = 16;
1824 const Int_t n4 = 22;
1825 const Double_t c3[17] = { .5578891446481605, -.11188325726569816,
1826 -.16337958125200939, .32256932072405902,
1827 -.14581632367244242, .03292677399374035,
1828 -.00460372142093573, 4.434706163314e-4,
1829 -3.142099529341e-5, 1.7123719938e-6,
1830 -7.416987005e-8, 2.61837671e-9,
1831 -7.685839e-11, 1.9067e-12,
1832 -4.052e-14, 7.5e-16,
1834 const Double_t c4[23] = { 1.00757647293865641, .00750316051248257,
1835 -7.043933264519e-5, 2.66205393382e-6,
1836 -1.8841157753e-7, 1.949014958e-8,
1837 -2.6126199e-9, 4.236269e-10,
1838 -7.955156e-11, 1.679973e-11,
1839 -3.9072e-12, 9.8543e-13,
1840 -2.6636e-13, 7.645e-14,
1841 -2.313e-14, 7.33e-15,
1844 -4
e-17, 2
e-17,-1
e-17 };
1855 }
else if (
v <= 0.3) {
1860 for (i = 1; i <= i1; ++i) {
1861 h = -
h*
y / ((2*i+ 1)*(2*i + 3));
1871 for (i = n3; i >= 0; --i) {
1872 b0 =
c3[i] + alfa*b1 - b2;
1883 for (i = n4; i >= 0; --i) {
1884 b0 = c4[i] + alfa*b1 - b2;
1910 for (i=1; i<=60;i++) {
1911 r*=(
x/(2*i+1))*(
x/(2*i+1));
1919 for (i=1; i<=km; i++) {
1920 r*=(2*i-1)*(2*i-1)/
x/
x;
1927 for (i=1; i<=16; i++) {
1928 r=0.125*
r*(2.0*i-1.0)*(2.0*i-1.0)/(i*
x);
1934 sl0=-2.0/(pi*
x)*s+bi0;
1952 for (i=1; i<=60;i++) {
1953 r*=
x*
x/(4.0*i*i-1.0);
1962 for (i=1; i<=km; i++) {
1963 r*=(2*i+3)*(2*i+1)/
x/
x;
1967 sl1=2.0/pi*(-1.0+1.0/(
x*
x)+3.0*s/(
x*
x*
x*
x));
1971 for (i=1; i<=16; i++) {
1972 r=-0.125*
r*(4.0-(2.0*i-1.0)*(2.0*i-1.0))/(i*
x);
1986 return ::ROOT::Math::beta(p,
q);
2006 d = 1.0 - qab*
x/qap;
2010 for (
m=1;
m<=itmax;
m++) {
2012 aa =
m*(
b-
m)*
x/((qam+ m2)*(
a+m2));
2019 aa = -(
a+
m)*(qab +
m)*
x/((
a+m2)*(qap+m2));
2030 Info(
"TMath::BetaCf",
"a or b too big, or itmax too small, a=%g, b=%g, x=%g, h=%g, itmax=%d",
2046 if ((
x<0) || (
x>1) || (p<=0) || (
q<=0)){
2047 Error(
"TMath::BetaDist",
"parameter value outside allowed range");
2064 if ((
x<0) || (
x>1) || (p<=0) || (
q<=0)){
2065 Error(
"TMath::BetaDistI",
"parameter value outside allowed range");
2077 return ::ROOT::Math::inc_beta(
x,
a,
b);
2086 if (k==0 ||
n==k)
return 1;
2110 if(k <= 0)
return 1.0;
2111 if(k >
n)
return 0.0;
2159 Double_t c[]={0, 0.01, 0.222222, 0.32, 0.4, 1.24, 2.2,
2160 4.67, 6.66, 6.73, 13.32, 60.0, 70.0,
2161 84.0, 105.0, 120.0, 127.0, 140.0, 175.0,
2162 210.0, 252.0, 264.0, 294.0, 346.0, 420.0,
2163 462.0, 606.0, 672.0, 707.0, 735.0, 889.0,
2164 932.0, 966.0, 1141.0, 1182.0, 1278.0, 1740.0,
2172 if (ndf <= 0)
return 0;
2185 if (ch >
c[6]*ndf + 6)
2192 p1 = 1 + ch * (
c[7]+ch);
2193 p2 = ch * (
c[9] + ch * (
c[8] + ch));
2194 t = -0.5 + (
c[7] + 2 * ch) / p1 - (
c[9] + ch * (
c[10] + 3 * ch)) / p2;
2195 ch = ch - (1 -
TMath::Exp(
a +
g + 0.5 * ch + cp * aa) *p2 / p1) / t;
2200 if (ch <
e)
return ch;
2203 for (
Int_t i=0; i<maxit; i++){
2210 a = 0.5 * t -
b * cp;
2211 s1 = (
c[19] +
a * (
c[17] +
a * (
c[14] +
a * (
c[13] +
a * (
c[12] +
c[11] *
a))))) /
c[24];
2212 s2 = (
c[24] +
a * (
c[29] +
a * (
c[32] +
a * (
c[33] +
c[35] *
a)))) /
c[37];
2213 s3 = (
c[19] +
a * (
c[25] +
a * (
c[28] +
c[31] *
a))) /
c[37];
2214 s4 = (
c[20] +
a * (
c[27] +
c[34] *
a) + cp * (
c[22] +
a * (
c[30] +
c[36] *
a))) /
c[38];
2215 s5 = (
c[13] +
c[21] *
a + cp * (
c[18] +
c[26] *
a)) /
c[37];
2216 s6 = (
c[15] + cp * (
c[23] +
c[16] * cp)) /
c[38];
2217 ch = ch + t * (1 + 0.5 * t *
s1 -
b * cp * (
s1 -
b * (s2 -
b * (s3 -
b * (s4 -
b * (s5 -
b * s6))))));
2242 return ::ROOT::Math::fdistribution_pdf(
F,
N,M);
2310 if ((
x<mu) || (gamma<=0) || (beta <=0)) {
2311 Error(
"TMath::GammaDist",
"illegal parameter values x = %f , gamma = %f beta = %f",
x,gamma,beta);
2314 return ::ROOT::Math::gamma_pdf(
x, gamma, beta, mu);
2399 if ((
x<theta) || (
sigma<=0) || (
m<=0)) {
2400 Error(
"TMath::Lognormal",
"illegal parameter values");
2418 if ((p<=0)||(p>=1)) {
2419 Error(
"TMath::NormQuantile",
"probability outside (0, 1)");
2423 Double_t a0 = 3.3871328727963666080e0;
2424 Double_t a1 = 1.3314166789178437745e+2;
2425 Double_t a2 = 1.9715909503065514427e+3;
2426 Double_t a3 = 1.3731693765509461125e+4;
2427 Double_t a4 = 4.5921953931549871457e+4;
2428 Double_t a5 = 6.7265770927008700853e+4;
2429 Double_t a6 = 3.3430575583588128105e+4;
2430 Double_t a7 = 2.5090809287301226727e+3;
2431 Double_t b1 = 4.2313330701600911252e+1;
2432 Double_t b2 = 6.8718700749205790830e+2;
2433 Double_t b3 = 5.3941960214247511077e+3;
2434 Double_t b4 = 2.1213794301586595867e+4;
2435 Double_t b5 = 3.9307895800092710610e+4;
2436 Double_t b6 = 2.8729085735721942674e+4;
2437 Double_t b7 = 5.2264952788528545610e+3;
2438 Double_t c0 = 1.42343711074968357734e0;
2442 Double_t c4 = 1.27045825245236838258e0;
2443 Double_t c5 = 2.41780725177450611770e-1;
2444 Double_t c6 = 2.27238449892691845833e-2;
2445 Double_t c7 = 7.74545014278341407640e-4;
2446 Double_t d1 = 2.05319162663775882187e0;
2447 Double_t d2 = 1.67638483018380384940e0;
2448 Double_t d3 = 6.89767334985100004550e-1;
2449 Double_t d4 = 1.48103976427480074590e-1;
2450 Double_t d5 = 1.51986665636164571966e-2;
2451 Double_t d6 = 5.47593808499534494600e-4;
2452 Double_t d7 = 1.05075007164441684324e-9;
2453 Double_t e0 = 6.65790464350110377720e0;
2454 Double_t e1 = 5.46378491116411436990e0;
2455 Double_t e2 = 1.78482653991729133580e0;
2456 Double_t e3 = 2.96560571828504891230e-1;
2457 Double_t e4 = 2.65321895265761230930e-2;
2458 Double_t e5 = 1.24266094738807843860e-3;
2459 Double_t e6 = 2.71155556874348757815e-5;
2460 Double_t e7 = 2.01033439929228813265e-7;
2462 Double_t f2 = 1.36929880922735805310e-1;
2463 Double_t f3 = 1.48753612908506148525e-2;
2464 Double_t f4 = 7.86869131145613259100e-4;
2465 Double_t f5 = 1.84631831751005468180e-5;
2466 Double_t f6 = 1.42151175831644588870e-7;
2467 Double_t f7 = 2.04426310338993978564e-15;
2478 quantile =
q* (((((((a7 *
r + a6) *
r + a5) *
r + a4) *
r + a3)
2479 *
r + a2) *
r + a1) *
r + a0) /
2480 (((((((b7 *
r + b6) *
r + b5) *
r + b4) *
r + b3)
2481 *
r + b2) *
r + b1) *
r + 1.);
2492 quantile=(((((((c7 *
r + c6) *
r + c5) *
r + c4) *
r +
c3)
2493 *
r +
c2) *
r +
c1) *
r + c0) /
2494 (((((((d7 *
r + d6) *
r + d5) *
r + d4) *
r + d3)
2495 *
r + d2) *
r + d1) *
r + 1);
2498 quantile=(((((((e7 *
r + e6) *
r + e5) *
r + e4) *
r + e3)
2499 *
r + e2) *
r + e1) *
r + e0) /
2500 (((((((f7 *
r + f6) *
r + f5) *
r + f4) *
r + f3)
2501 *
r + f2) *
r +
f1) *
r + 1);
2503 if (
q<0) quantile=-quantile;
2523 for(i=
n-2; i>-1; i--) {
2530 if(i1==-1)
return kFALSE;
2534 for(i=
n-1;i>i1;i--) {
2545 for(i=0;i<(
n-i1-1)/2;i++) {
2639 if (ndf<1 || p>=1 || p<=0) {
2640 Error(
"TMath::StudentQuantile",
"illegal parameter values");
2643 if ((lower_tail && p>0.5)||(!lower_tail && p<0.5)){
2645 q=2*(lower_tail ? (1-p) : p);
2648 q=2*(lower_tail? p : (1-p));
2668 if (ndf<5)
c+=0.3*(ndf-4.5)*(
x+0.6);
2669 c+=(((0.05*
d*
x-5.)*
x-7.)*
x-2.)*
x +
b;
2670 y=(((((0.4*
y+6.3)*
y+36.)*
y+94.5)/
c -
y-3.)/
b+1)*
x;
2675 y=((1./(((ndf+6.)/(ndf*
y)-0.089*
d-0.822)*(ndf+2.)*3)+0.5/(ndf+4.))*
y-1.)*
2676 (ndf+1.)/(ndf+2.)+1/
y;
2681 if(neg) quantile=-quantile;
2777 if (
x < ac[0])
v = 0;
2778 else if (
x >=ac[8])
v = 1;
2781 k =
Int_t(xx*ac[10]);
2782 v =
TMath::Min(wcm[k] + (xx - k*ac[9])*(wcm[k+1]-wcm[k])*ac[10], 1.);
2798 return ::ROOT::Math::landau_cdf(
x);
2808 Double_t BKMNX1 = 0.02, BKMNY1 = 0.05, BKMNX2 = 0.12, BKMNY2 = 0.05,
2809 BKMNX3 = 0.22, BKMNY3 = 0.05, BKMXX1 = 0.1 , BKMXY1 = 1,
2810 BKMXX2 = 0.2 , BKMXY2 = 1 , BKMXX3 = 0.3 , BKMXY3 = 1;
2812 Double_t FBKX1 = 2/(BKMXX1-BKMNX1), FBKX2 = 2/(BKMXX2-BKMNX2),
2813 FBKX3 = 2/(BKMXX3-BKMNX3), FBKY1 = 2/(BKMXY1-BKMNY1),
2814 FBKY2 = 2/(BKMXY2-BKMNY2), FBKY3 = 2/(BKMXY3-BKMNY3);
2816 Double_t FNINV[] = {0, 1, 0.5, 0.33333333, 0.25, 0.2};
2818 Double_t EDGEC[]= {0, 0, 0.16666667e+0, 0.41666667e-1, 0.83333333e-2,
2819 0.13888889e-1, 0.69444444e-2, 0.77160493e-3};
2821 Double_t U1[] = {0, 0.25850868e+0, 0.32477982e-1, -0.59020496e-2,
2822 0. , 0.24880692e-1, 0.47404356e-2,
2823 -0.74445130e-3, 0.73225731e-2, 0. ,
2824 0.11668284e-2, 0. , -0.15727318e-2,-0.11210142e-2};
2826 Double_t U2[] = {0, 0.43142611e+0, 0.40797543e-1, -0.91490215e-2,
2827 0. , 0.42127077e-1, 0.73167928e-2,
2828 -0.14026047e-2, 0.16195241e-1, 0.24714789e-2,
2829 0.20751278e-2, 0. , -0.25141668e-2,-0.14064022e-2};
2831 Double_t U3[] = {0, 0.25225955e+0, 0.64820468e-1, -0.23615759e-1,
2832 0. , 0.23834176e-1, 0.21624675e-2,
2833 -0.26865597e-2, -0.54891384e-2, 0.39800522e-2,
2834 0.48447456e-2, -0.89439554e-2, -0.62756944e-2,-0.24655436e-2};
2836 Double_t U4[] = {0, 0.12593231e+1, -0.20374501e+0, 0.95055662e-1,
2837 -0.20771531e-1, -0.46865180e-1, -0.77222986e-2,
2838 0.32241039e-2, 0.89882920e-2, -0.67167236e-2,
2839 -0.13049241e-1, 0.18786468e-1, 0.14484097e-1};
2841 Double_t U5[] = {0, -0.24864376e-1, -0.10368495e-2, 0.14330117e-2,
2842 0.20052730e-3, 0.18751903e-2, 0.12668869e-2,
2843 0.48736023e-3, 0.34850854e-2, 0. ,
2844 -0.36597173e-3, 0.19372124e-2, 0.70761825e-3, 0.46898375e-3};
2846 Double_t U6[] = {0, 0.35855696e-1, -0.27542114e-1, 0.12631023e-1,
2847 -0.30188807e-2, -0.84479939e-3, 0. ,
2848 0.45675843e-3, -0.69836141e-2, 0.39876546e-2,
2849 -0.36055679e-2, 0. , 0.15298434e-2, 0.19247256e-2};
2851 Double_t U7[] = {0, 0.10234691e+2, -0.35619655e+1, 0.69387764e+0,
2852 -0.14047599e+0, -0.19952390e+1, -0.45679694e+0,
2853 0. , 0.50505298e+0};
2854 Double_t U8[] = {0, 0.21487518e+2, -0.11825253e+2, 0.43133087e+1,
2855 -0.14500543e+1, -0.34343169e+1, -0.11063164e+1,
2856 -0.21000819e+0, 0.17891643e+1, -0.89601916e+0,
2857 0.39120793e+0, 0.73410606e+0, 0. ,-0.32454506e+0};
2859 Double_t V1[] = {0, 0.27827257e+0, -0.14227603e-2, 0.24848327e-2,
2860 0. , 0.45091424e-1, 0.80559636e-2,
2861 -0.38974523e-2, 0. , -0.30634124e-2,
2862 0.75633702e-3, 0.54730726e-2, 0.19792507e-2};
2864 Double_t V2[] = {0, 0.41421789e+0, -0.30061649e-1, 0.52249697e-2,
2865 0. , 0.12693873e+0, 0.22999801e-1,
2866 -0.86792801e-2, 0.31875584e-1, -0.61757928e-2,
2867 0. , 0.19716857e-1, 0.32596742e-2};
2869 Double_t V3[] = {0, 0.20191056e+0, -0.46831422e-1, 0.96777473e-2,
2870 -0.17995317e-2, 0.53921588e-1, 0.35068740e-2,
2871 -0.12621494e-1, -0.54996531e-2, -0.90029985e-2,
2872 0.34958743e-2, 0.18513506e-1, 0.68332334e-2,-0.12940502e-2};
2874 Double_t V4[] = {0, 0.13206081e+1, 0.10036618e+0, -0.22015201e-1,
2875 0.61667091e-2, -0.14986093e+0, -0.12720568e-1,
2876 0.24972042e-1, -0.97751962e-2, 0.26087455e-1,
2877 -0.11399062e-1, -0.48282515e-1, -0.98552378e-2};
2879 Double_t V5[] = {0, 0.16435243e-1, 0.36051400e-1, 0.23036520e-2,
2880 -0.61666343e-3, -0.10775802e-1, 0.51476061e-2,
2881 0.56856517e-2, -0.13438433e-1, 0. ,
2882 0. , -0.25421507e-2, 0.20169108e-2,-0.15144931e-2};
2884 Double_t V6[] = {0, 0.33432405e-1, 0.60583916e-2, -0.23381379e-2,
2885 0.83846081e-3, -0.13346861e-1, -0.17402116e-2,
2886 0.21052496e-2, 0.15528195e-2, 0.21900670e-2,
2887 -0.13202847e-2, -0.45124157e-2, -0.15629454e-2, 0.22499176e-3};
2889 Double_t V7[] = {0, 0.54529572e+1, -0.90906096e+0, 0.86122438e-1,
2890 0. , -0.12218009e+1, -0.32324120e+0,
2891 -0.27373591e-1, 0.12173464e+0, 0. ,
2892 0. , 0.40917471e-1};
2894 Double_t V8[] = {0, 0.93841352e+1, -0.16276904e+1, 0.16571423e+0,
2895 0. , -0.18160479e+1, -0.50919193e+0,
2896 -0.51384654e-1, 0.21413992e+0, 0. ,
2897 0. , 0.66596366e-1};
2899 Double_t W1[] = {0, 0.29712951e+0, 0.97572934e-2, 0. ,
2900 -0.15291686e-2, 0.35707399e-1, 0.96221631e-2,
2901 -0.18402821e-2, -0.49821585e-2, 0.18831112e-2,
2902 0.43541673e-2, 0.20301312e-2, -0.18723311e-2,-0.73403108e-3};
2904 Double_t W2[] = {0, 0.40882635e+0, 0.14474912e-1, 0.25023704e-2,
2905 -0.37707379e-2, 0.18719727e+0, 0.56954987e-1,
2906 0. , 0.23020158e-1, 0.50574313e-2,
2907 0.94550140e-2, 0.19300232e-1};
2909 Double_t W3[] = {0, 0.16861629e+0, 0. , 0.36317285e-2,
2910 -0.43657818e-2, 0.30144338e-1, 0.13891826e-1,
2911 -0.58030495e-2, -0.38717547e-2, 0.85359607e-2,
2912 0.14507659e-1, 0.82387775e-2, -0.10116105e-1,-0.55135670e-2};
2914 Double_t W4[] = {0, 0.13493891e+1, -0.26863185e-2, -0.35216040e-2,
2915 0.24434909e-1, -0.83447911e-1, -0.48061360e-1,
2916 0.76473951e-2, 0.24494430e-1, -0.16209200e-1,
2917 -0.37768479e-1, -0.47890063e-1, 0.17778596e-1, 0.13179324e-1};
2919 Double_t W5[] = {0, 0.10264945e+0, 0.32738857e-1, 0. ,
2920 0.43608779e-2, -0.43097757e-1, -0.22647176e-2,
2921 0.94531290e-2, -0.12442571e-1, -0.32283517e-2,
2922 -0.75640352e-2, -0.88293329e-2, 0.52537299e-2, 0.13340546e-2};
2924 Double_t W6[] = {0, 0.29568177e-1, -0.16300060e-2, -0.21119745e-3,
2925 0.23599053e-2, -0.48515387e-2, -0.40797531e-2,
2926 0.40403265e-3, 0.18200105e-2, -0.14346306e-2,
2927 -0.39165276e-2, -0.37432073e-2, 0.19950380e-2, 0.12222675e-2};
2929 Double_t W8[] = {0, 0.66184645e+1, -0.73866379e+0, 0.44693973e-1,
2930 0. , -0.14540925e+1, -0.39529833e+0,
2931 -0.44293243e-1, 0.88741049e-1};
2934 if (rkappa <0.01 || rkappa >12) {
2935 Error(
"Vavilov distribution",
"illegal value of kappa");
2943 Double_t x,
y, xx, yy,
x2,
x3, y2, y3,
xy, p2, p3, q2, q3, pq;
2944 if (rkappa >= 0.29) {
2949 AC[0] = (-0.032227*beta2-0.074275)*rkappa + (0.24533*beta2+0.070152)*wk + (-0.55610*beta2-3.1579);
2950 AC[8] = (-0.013483*beta2-0.048801)*rkappa + (-1.6921*beta2+8.3656)*wk + (-0.73275*beta2-3.5226);
2953 for (j=1; j<=4; j++) {
2954 DRK[j+1] = DRK[1]*DRK[j];
2955 DSIGM[j+1] = DSIGM[1]*DSIGM[j];
2956 ALFA[j+1] = (FNINV[j]-beta2*FNINV[j+1])*DRK[j];
2960 HC[2]=ALFA[3]*DSIGM[3];
2961 HC[3]=(3*ALFA[2]*ALFA[2] + ALFA[4])*DSIGM[4]-3;
2962 HC[4]=(10*ALFA[2]*ALFA[3]+ALFA[5])*DSIGM[5]-10*
HC[2];
2966 for (j=2; j<=7; j++)
2968 HC[8]=0.39894228*
HC[1];
2970 else if (rkappa >=0.22) {
2973 x = 1+(rkappa-BKMXX3)*FBKX3;
2987 AC[1] = W1[1] + W1[2]*
x + W1[4]*
x3 + W1[5]*
y + W1[6]*y2 + W1[7]*y3 +
2988 W1[8]*
xy + W1[9]*p2 + W1[10]*p3 + W1[11]*q2 + W1[12]*q3 + W1[13]*pq;
2989 AC[2] = W2[1] + W2[2]*
x + W2[3]*
x2 + W2[4]*
x3 + W2[5]*
y + W2[6]*y2 +
2990 W2[8]*
xy + W2[9]*p2 + W2[10]*p3 + W2[11]*q2;
2991 AC[3] = W3[1] + W3[3]*
x2 + W3[4]*
x3 + W3[5]*
y + W3[6]*y2 + W3[7]*y3 +
2992 W3[8]*
xy + W3[9]*p2 + W3[10]*p3 + W3[11]*q2 + W3[12]*q3 + W3[13]*pq;
2993 AC[4] = W4[1] + W4[2]*
x + W4[3]*
x2 + W4[4]*
x3 + W4[5]*
y + W4[6]*y2 + W4[7]*y3 +
2994 W4[8]*
xy + W4[9]*p2 + W4[10]*p3 + W4[11]*q2 + W4[12]*q3 + W4[13]*pq;
2995 AC[5] = W5[1] + W5[2]*
x + W5[4]*
x3 + W5[5]*
y + W5[6]*y2 + W5[7]*y3 +
2996 W5[8]*
xy + W5[9]*p2 + W5[10]*p3 + W5[11]*q2 + W5[12]*q3 + W5[13]*pq;
2997 AC[6] = W6[1] + W6[2]*
x + W6[3]*
x2 + W6[4]*
x3 + W6[5]*
y + W6[6]*y2 + W6[7]*y3 +
2998 W6[8]*
xy + W6[9]*p2 + W6[10]*p3 + W6[11]*q2 + W6[12]*q3 + W6[13]*pq;
2999 AC[8] = W8[1] + W8[2]*
x + W8[3]*
x2 + W8[5]*
y + W8[6]*y2 + W8[7]*y3 + W8[8]*
xy;
3001 }
else if (rkappa >= 0.12) {
3004 x = 1 + (rkappa-BKMXX2)*FBKX2;
3018 AC[1] = V1[1] + V1[2]*
x + V1[3]*
x2 + V1[5]*
y + V1[6]*y2 + V1[7]*y3 +
3019 V1[9]*p2 + V1[10]*p3 + V1[11]*q2 + V1[12]*q3;
3020 AC[2] = V2[1] + V2[2]*
x + V2[3]*
x2 + V2[5]*
y + V2[6]*y2 + V2[7]*y3 +
3021 V2[8]*
xy + V2[9]*p2 + V2[11]*q2 + V2[12]*q3;
3022 AC[3] = V3[1] + V3[2]*
x + V3[3]*
x2 + V3[4]*
x3 + V3[5]*
y + V3[6]*y2 + V3[7]*y3 +
3023 V3[8]*
xy + V3[9]*p2 + V3[10]*p3 + V3[11]*q2 + V3[12]*q3 + V3[13]*pq;
3024 AC[4] = V4[1] + V4[2]*
x + V4[3]*
x2 + V4[4]*
x3 + V4[5]*
y + V4[6]*y2 + V4[7]*y3 +
3025 V4[8]*
xy + V4[9]*p2 + V4[10]*p3 + V4[11]*q2 + V4[12]*q3;
3026 AC[5] = V5[1] + V5[2]*
x + V5[3]*
x2 + V5[4]*
x3 + V5[5]*
y + V5[6]*y2 + V5[7]*y3 +
3027 V5[8]*
xy + V5[11]*q2 + V5[12]*q3 + V5[13]*pq;
3028 AC[6] = V6[1] + V6[2]*
x + V6[3]*
x2 + V6[4]*
x3 + V6[5]*
y + V6[6]*y2 + V6[7]*y3 +
3029 V6[8]*
xy + V6[9]*p2 + V6[10]*p3 + V6[11]*q2 + V6[12]*q3 + V6[13]*pq;
3030 AC[7] = V7[1] + V7[2]*
x + V7[3]*
x2 + V7[5]*
y + V7[6]*y2 + V7[7]*y3 +
3031 V7[8]*
xy + V7[11]*q2;
3032 AC[8] = V8[1] + V8[2]*
x + V8[3]*
x2 + V8[5]*
y + V8[6]*y2 + V8[7]*y3 +
3033 V8[8]*
xy + V8[11]*q2;
3037 if (rkappa >=0.02) itype = 3;
3039 x = 1+(rkappa-BKMXX1)*FBKX1;
3054 AC[1] = U1[1] + U1[2]*
x + U1[3]*
x2 + U1[5]*
y + U1[6]*y2 + U1[7]*y3 +
3055 U1[8]*
xy + U1[10]*p3 + U1[12]*q3 + U1[13]*pq;
3056 AC[2] = U2[1] + U2[2]*
x + U2[3]*
x2 + U2[5]*
y + U2[6]*y2 + U2[7]*y3 +
3057 U2[8]*
xy + U2[9]*p2 + U2[10]*p3 + U2[12]*q3 + U2[13]*pq;
3058 AC[3] = U3[1] + U3[2]*
x + U3[3]*
x2 + U3[5]*
y + U3[6]*y2 + U3[7]*y3 +
3059 U3[8]*
xy + U3[9]*p2 + U3[10]*p3 + U3[11]*q2 + U3[12]*q3 + U3[13]*pq;
3060 AC[4] = U4[1] + U4[2]*
x + U4[3]*
x2 + U4[4]*
x3 + U4[5]*
y + U4[6]*y2 + U4[7]*y3 +
3061 U4[8]*
xy + U4[9]*p2 + U4[10]*p3 + U4[11]*q2 + U4[12]*q3;
3062 AC[5] = U5[1] + U5[2]*
x + U5[3]*
x2 + U5[4]*
x3 + U5[5]*
y + U5[6]*y2 + U5[7]*y3 +
3063 U5[8]*
xy + U5[10]*p3 + U5[11]*q2 + U5[12]*q3 + U5[13]*pq;
3064 AC[6] = U6[1] + U6[2]*
x + U6[3]*
x2 + U6[4]*
x3 + U6[5]*
y + U6[7]*y3 +
3065 U6[8]*
xy + U6[9]*p2 + U6[10]*p3 + U6[12]*q3 + U6[13]*pq;
3066 AC[7] = U7[1] + U7[2]*
x + U7[3]*
x2 + U7[4]*
x3 + U7[5]*
y + U7[6]*y2 + U7[8]*
xy;
3068 AC[8] = U8[1] + U8[2]*
x + U8[3]*
x2 + U8[4]*
x3 + U8[5]*
y + U8[6]*y2 + U8[7]*y3 +
3069 U8[8]*
xy + U8[9]*p2 + U8[10]*p3 + U8[11]*q2 + U8[13]*pq;
3073 AC[9] = (AC[8] - AC[0])/npt;
3076 x = (AC[7]-AC[8])/(AC[7]*AC[8]);
3079 AC[11] = p2*(AC[1]*
TMath::Exp(-AC[2]*(AC[7]+AC[5]*p2)-
3080 AC[3]*
TMath::Exp(-AC[4]*(AC[7]+AC[6]*p2)))-0.045*
y/AC[7])/(1+
x*
y*AC[7]);
3081 AC[12] = (0.045+
x*AC[11])*
y;
3083 if (itype == 4) AC[13] = 0.995/
LandauI(AC[8]);
3085 if (mode==0)
return;
3092 for (k=1; k<=npt; k++) {
3095 WCM[k] = WCM[k-1] + fl + fu;
3099 for (k=1; k<=npt; k++)
3109 if (rlam < AC[0] || rlam > AC[8])
3116 x = (rlam +
HC[0])*
HC[1];
3119 for (k=2; k<=8; k++) {
3121 h[k+1] =
x*
h[k]-fn*
h[k-1];
3124 for (k=2; k<=6; k++)
3128 else if (itype == 2) {
3132 else if (itype == 3) {
3138 v = (AC[11]*
x + AC[12])*
x;
3141 else if (itype == 4) {
3149#ifdef R__HAS_VECCORE
static const double x2[5]
static const double x3[11]
#define NamespaceImp(name)
void Info(const char *location, const char *msgfmt,...)
Use this function for informational messages.
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
UInt_t Hash(const TString &s)
void ToUpper()
Change string to upper case.
UInt_t Hash(ECaseCompare cmp=kExact) const
Return hash value.
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
double landau_pdf(double x, double xi=1, double x0=0)
Probability density function of the Landau distribution:
double fdistribution_cdf(double x, double n, double m, double x0=0)
Cumulative distribution function of the F-distribution (lower tail).
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...
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.
Double_t DiLog(Double_t x)
Modified Struve functions of order 1.
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...
Double_t GamSer(Double_t a, Double_t x)
Computation of the incomplete gamma function P(a,x) via its series representation.
Double_t VavilovDenEval(Double_t rlam, Double_t *AC, Double_t *HC, Int_t itype)
Internal function, called by Vavilov and VavilovSet.
Double_t BesselI(Int_t n, Double_t x)
Compute the Integer Order Modified Bessel function I_n(x) for n=0,1,2,... and any real x.
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).
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.
Double_t Factorial(Int_t i)
Compute factorial(n).
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 ...
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
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...
Short_t Max(Short_t a, Short_t b)
Double_t Vavilov(Double_t x, Double_t kappa, Double_t beta2)
Returns the value of the Vavilov density function.
Double_t Binomial(Int_t n, Int_t k)
Calculate the binomial coefficient n over k.
Float_t Normalize(Float_t v[3])
Normalize a vector v in place.
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...
Double_t Log2(Double_t x)
Double_t BesselK1(Double_t x)
modified Bessel function I_1(x)
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...
Double_t BesselI1(Double_t x)
modified Bessel function K_0(x)
Double_t Erf(Double_t x)
Computation of the error function erf(x).
Bool_t Permute(Int_t n, Int_t *a)
Simple recursive algorithm to find the permutations of n natural numbers, not necessarily all distinc...
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754
Double_t PoissonI(Double_t x, Double_t par)
Compute the Discrete Poisson distribution function for (x,par).
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...
Double_t StruveL1(Double_t x)
Modified Struve functions of order 0.
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 ...
ULong_t Hash(const void *txt, Int_t ntxt)
Calculates hash index from any char string.
Double_t BreitWigner(Double_t x, Double_t mean=0, Double_t gamma=1)
Calculate a Breit Wigner function with mean and gamma.
Double_t Landau(Double_t x, Double_t mpv=0, Double_t sigma=1, Bool_t norm=kFALSE)
The LANDAU function.
Double_t Voigt(Double_t x, Double_t sigma, Double_t lg, Int_t r=4)
Computation of Voigt function (normalised).
Double_t Student(Double_t T, Double_t ndf)
Computes density function for Student's t- distribution (the probability function (integral of densit...
constexpr Double_t PiOver2()
Double_t BetaDistI(Double_t x, Double_t p, Double_t q)
Computes the distribution function of the Beta distribution.
Int_t FloorNint(Double_t x)
Double_t BesselK0(Double_t x)
modified Bessel function I_0(x)
Double_t BesselY0(Double_t x)
Bessel function J1(x) for any real x.
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...
Double_t ErfInverse(Double_t x)
returns the inverse error function x must be <-1<x<1
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...
Double_t Erfc(Double_t x)
Compute the complementary error function erfc(x).
Double_t VavilovI(Double_t x, Double_t kappa, Double_t beta2)
Returns the value of the Vavilov distribution function.
Double_t Beta(Double_t p, Double_t q)
Calculates Beta-function Gamma(p)*Gamma(q)/Gamma(p+q).
Double_t Poisson(Double_t x, Double_t par)
Compute the Poisson distribution function for (x,par).
Double_t Sqrt(Double_t x)
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Double_t BesselJ0(Double_t x)
modified Bessel function K_1(x)
Double_t Gamma(Double_t z)
Computation of gamma(z) for all z.
Short_t Min(Short_t a, Short_t b)
Double_t StruveL0(Double_t x)
Struve functions of order 1.
Double_t NormQuantile(Double_t p)
Computes quantiles for standard normal distribution N(0, 1) at probability p.
Double_t GamCf(Double_t a, Double_t x)
Computation of the incomplete gamma function P(a,x) via its continued fraction representation.
Double_t Hypot(Double_t x, Double_t y)
Double_t StruveH0(Double_t x)
Bessel function Y1(x) for positive x.
Double_t LnGamma(Double_t z)
Computation of ln[gamma(z)] for all z.
Double_t KolmogorovProb(Double_t z)
Calculates the Kolmogorov distribution function,.
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.
Double_t ChisquareQuantile(Double_t p, Double_t ndf)
Evaluate the quantiles of the chi-squared probability distribution function.
Double_t FDist(Double_t F, Double_t N, Double_t M)
Computes the density function of F-distribution (probability function, integral of density,...
Double_t SignalingNaN()
Returns a signaling NaN as defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Signaling_NaN)
void BubbleLow(Int_t Narr, Double_t *arr1, Int_t *arr2)
Opposite ordering of the array arr2[] to that of BubbleHigh.
Double_t BesselK(Int_t n, Double_t x)
integer order modified Bessel function I_n(x)
Double_t BesselJ1(Double_t x)
Bessel function J0(x) for any real x.
Double_t BetaIncomplete(Double_t x, Double_t a, Double_t b)
Calculates the incomplete Beta-function.
Double_t StruveH1(Double_t x)
Struve functions of order 0.
Double_t Freq(Double_t x)
Computation of the normal frequency function freq(x).
Double_t LandauI(Double_t x)
Returns the value of the Landau distribution function at point x.
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.
Double_t BesselI0(Double_t x)
integer order modified Bessel function K_n(x)
void VavilovSet(Double_t rkappa, Double_t beta2, Bool_t mode, Double_t *WCM, Double_t *AC, Double_t *HC, Int_t &itype, Int_t &npt)
Internal function, called by Vavilov and VavilovI.
Double_t Log10(Double_t x)
Double_t StudentI(Double_t T, Double_t ndf)
Calculates the cumulative distribution function of Student's t-distribution second parameter stands f...
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,...
Double_t BesselY1(Double_t x)
Bessel function Y0(x) for positive x.
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.
constexpr Double_t HC()
in
Double_t ErfcInverse(Double_t x)
returns the inverse of the complementary error function x must be 0<x<2 implement using the quantile ...
static T Epsilon()
Returns minimum double representation.
static uint64_t sum(uint64_t i)