Subversion Repositories planix.SVN

Rev

Rev 2 | Details | Compare with Previous | Last modification | View Log | RSS feed

Rev Author Line No. Line
2 - 1 
/*
2 
	floating-point arctangent
3 
 
4 
	atan returns the value of the arctangent of its
5 
	argument in the range [-pi/2,pi/2].
6 
 
7 
	atan2 returns the arctangent of arg1/arg2
8 
	in the range [-pi,pi].
9 
 
10 
	there are no error returns.
11 
 
12 
	coefficients are #5077 from Hart & Cheney. (19.56D)
13 
*/
14 
 
15 
#include <math.h>
16 
 
17 
#define sq2p1 2.414213562373095048802e0
18 
#define sq2m1  .414213562373095048802e0
19 
#define pio2 1.570796326794896619231e0
20 
#define pio4  .785398163397448309615e0
21 
#define p4  .161536412982230228262e2
22 
#define p3  .26842548195503973794141e3
23 
#define p2  .11530293515404850115428136e4
24 
#define p1  .178040631643319697105464587e4
25 
#define p0  .89678597403663861959987488e3
26 
#define q4  .5895697050844462222791e2
27 
#define q3  .536265374031215315104235e3
28 
#define q2  .16667838148816337184521798e4
29 
#define q1  .207933497444540981287275926e4
30 
#define q0  .89678597403663861962481162e3
31 
 
32 
 
33 
/*
34 
	xatan evaluates a series valid in the
35 
	range [-0.414...,+0.414...].
36 
 */
37 
 
38 
static
39 
double
40 
xatan(double arg)
41 
{
42 
	double argsq, value;
43 
 
44 
	/* get denormalized add in following if range arg**10 is much smaller
45 
	    than q1, so check for that case
46 
	*/
47 
	if(-.01 < arg && arg < .01)
48 
		value = p0/q0;
49 
	else {
50 
		argsq = arg*arg;
51 
		value = ((((p4*argsq + p3)*argsq + p2)*argsq + p1)*argsq + p0);
52 
		value = value/(((((argsq + q4)*argsq + q3)*argsq + q2)*argsq + q1)*argsq + q0);
53 
	}
54 
	return value*arg;
55 
}
56 
 
57 
/*
58 
	satan reduces its argument (known to be positive)
59 
	to the range [0,0.414...] and calls xatan.
60 
 */
61 
 
62 
static
63 
double
64 
satan(double arg)
65 
{
66 
 
67 
	if(arg < sq2m1)
68 
		return xatan(arg);
69 
	if(arg > sq2p1)
70 
		return pio2 - xatan(1.0/arg);
71 
	return pio4 + xatan((arg-1.0)/(arg+1.0));
72 
}
73 
 
74 
/*
75 
	atan makes its argument positive and
76 
	calls the inner routine satan.
77 
 */
78 
 
79 
double
80 
atan(double arg)
81 
{
82 
 
83 
	if(arg > 0)
84 
		return satan(arg);
85 
	return -satan(-arg);
86 
}