Commit | Line | Data |
---|---|---|
94759b7b TL |
1 | |
2 | /* | |
3 | floating-point arctangent | |
4 | ||
5 | atan returns the value of the arctangent of its | |
6 | argument in the range [-pi/2,pi/2]. | |
7 | ||
8 | atan2 returns the arctangent of arg1/arg2 | |
9 | in the range [-pi,pi]. | |
10 | ||
11 | there are no error returns. | |
12 | ||
13 | coefficients are #5077 from Hart & Cheney. (19.56D) | |
14 | */ | |
15 | ||
16 | ||
17 | double static sq2p1 =2.414213562373095048802e0; | |
18 | static double sq2m1 = .414213562373095048802e0; | |
19 | static double pio2 =1.570796326794896619231e0; | |
20 | static double pio4 = .785398163397448309615e0; | |
21 | static double p4 = .161536412982230228262e2; | |
22 | static double p3 = .26842548195503973794141e3; | |
23 | static double p2 = .11530293515404850115428136e4; | |
24 | static double p1 = .178040631643319697105464587e4; | |
25 | static double p0 = .89678597403663861959987488e3; | |
26 | static double q4 = .5895697050844462222791e2; | |
27 | static double q3 = .536265374031215315104235e3; | |
28 | static double q2 = .16667838148816337184521798e4; | |
29 | static double q1 = .207933497444540981287275926e4; | |
30 | static double q0 = .89678597403663861962481162e3; | |
31 | ||
32 | ||
33 | /* | |
34 | atan makes its argument positive and | |
35 | calls the inner routine satan. | |
36 | */ | |
37 | ||
38 | double | |
39 | atan(arg) | |
40 | double arg; | |
41 | { | |
42 | double satan(); | |
43 | ||
44 | if(arg>0) | |
45 | return(satan(arg)); | |
46 | else | |
47 | return(-satan(-arg)); | |
48 | } | |
49 | ||
50 | ||
51 | /* | |
52 | atan2 discovers what quadrant the angle | |
53 | is in and calls atan. | |
54 | */ | |
55 | ||
56 | double | |
57 | atan2(arg1,arg2) | |
58 | double arg1,arg2; | |
59 | { | |
60 | double satan(); | |
61 | ||
62 | if((arg1+arg2)==arg1) | |
63 | if(arg1 >= 0.) return(pio2); | |
64 | else return(-pio2); | |
65 | else if(arg2 <0.) | |
66 | if(arg1 >= 0.) | |
67 | return(pio2+pio2 - satan(-arg1/arg2)); | |
68 | else | |
69 | return(-pio2-pio2 + satan(arg1/arg2)); | |
70 | else if(arg1>0) | |
71 | return(satan(arg1/arg2)); | |
72 | else | |
73 | return(-satan(-arg1/arg2)); | |
74 | } | |
75 | ||
76 | /* | |
77 | satan reduces its argument (known to be positive) | |
78 | to the range [0,0.414...] and calls xatan. | |
79 | */ | |
80 | ||
81 | static double | |
82 | satan(arg) | |
83 | double arg; | |
84 | { | |
85 | double xatan(); | |
86 | ||
87 | if(arg < sq2m1) | |
88 | return(xatan(arg)); | |
89 | else if(arg > sq2p1) | |
90 | return(pio2 - xatan(1.0/arg)); | |
91 | else | |
92 | return(pio4 + xatan((arg-1.0)/(arg+1.0))); | |
93 | } | |
94 | ||
95 | /* | |
96 | xatan evaluates a series valid in the | |
97 | range [-0.414...,+0.414...]. | |
98 | */ | |
99 | ||
100 | static double | |
101 | xatan(arg) | |
102 | double arg; | |
103 | { | |
104 | double argsq; | |
105 | double value; | |
106 | ||
107 | argsq = arg*arg; | |
108 | value = ((((p4*argsq + p3)*argsq + p2)*argsq + p1)*argsq + p0); | |
109 | value = value/(((((argsq + q4)*argsq + q3)*argsq + q2)*argsq + q1)*argsq + q0); | |
110 | return(value*arg); | |
111 | } |