GRASS GIS 8 Programmer's Manual 8.3.2(2024)-exported
Loading...
Searching...
No Matches
as66.c
Go to the documentation of this file.
1/*-Algorithm AS 66
2 * The Normal Integral, by I. D. Hill, 1973.
3 * Applied Statistics 22(3):424-427.
4 *
5 * Translation to C by James Darrell McCauley, mccauley@ecn.purdue.edu.
6 *
7 * Calculates the upper or lower tail area of the standardized normal
8 * curve corresponding to any given argument.
9 *
10 * x - the argument value
11 * upper: 1 -> the area from x to \infty
12 * 0 -> the area from -\infty to x
13 *
14 * Notes:
15 * The constant LTONE should be set to the value at which the
16 * lower tail area becomes 1.0 to the accuracy of the machine.
17 * LTONE=(n+9)/3 gives the required value accurately enough, for a
18 * machine that produces n decimal digits in its real numbers.
19 *
20 * The constant UTZERO should be set to the value at which the upper
21 * tail area becomes 0.0 to the accuracy of the machine. This may be
22 * taken as the value such that exp(-0.5 * UTZERO * UTZERO) /
23 * (UTZERO * sqrt(2*M_PI)) is just greater than the smallest allowable
24 * real numbers.
25 */
26
27#include <math.h>
28
29#define LTONE 7.0
30#define UTZERO 18.66
31
32double Cdhc_alnorm(double x, int upper)
33{
34 double ret, z, y;
35 int up;
36
37 up = upper;
38 z = x;
39
40 if (x < 0.0) {
41 up = up == 0 ? 1 : 0;
42 z = -x;
43 }
44
45 if (!(z <= LTONE || (up == 1 && z <= UTZERO)))
46 ret = 0.0;
47 else {
48 y = 0.5 * z * z;
49 if (z <= 1.28)
50 ret = 0.5 - z * (0.398942280444 -
51 0.399903438504 * y /
52 (y + 5.75885480458 -
53 29.8213557808 /
54 (y + 2.62433121679 +
55 48.6959930692 / (y + 5.92885724438))));
56 else
57 ret = 0.398942280385 * exp(-y) /
58 (z - 3.8052e-8 +
59 1.00000615302 /
60 (z + 3.98064794e-4 +
61 1.98615381364 /
62 (z - 0.151679116635 +
63 5.29330324926 /
64 (z + 4.8385912808 -
65 15.1508972451 /
66 (z + 0.742380924027 +
67 30.789933034 / (z + 3.99019417011))))));
68 }
69
70 if (up == 0)
71 ret = 1.0 - ret;
72
73 return ret;
74}
double Cdhc_alnorm(double x, int upper)
Definition as66.c:32
#define UTZERO
Definition as66.c:30
#define LTONE
Definition as66.c:29
#define x