1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
|
/* @(#)z_atangent.c 1.0 98/08/13 */
/******************************************************************
* The following routines are coded directly from the algorithms
* and coefficients given in "Software Manual for the Elementary
* Functions" by William J. Cody, Jr. and William Waite, Prentice
* Hall, 1980.
******************************************************************/
/*
FUNCTION
<<atan>>, <<atanf>>, <<atan2>>, <<atan2f>>, <<atangent>>, <<atangentf>>---arc tangent
INDEX
atan2
INDEX
atan2f
INDEX
atan
INDEX
atanf
ANSI_SYNOPSIS
#include <math.h>
double atan(double <[x]>);
float atan(float <[x]>);
double atan2(double <[y]>,double <[x]>);
float atan2f(float <[y]>,float <[x]>);
TRAD_SYNOPSIS
#include <math.h>
double atan2(<[y]>,<[x]>);
double <[y]>;
double <[x]>;
float atan2f(<[y]>,<[x]>);
float <[y]>;
float <[x]>;
#include <math.h>
double atan(<[x]>);
double <[x]>;
float atanf(<[x]>);
float <[x]>;
DESCRIPTION
<<atan2>> computes the inverse tangent (arc tangent) of y / x.
<<atan2f>> is identical to <<atan2>>, save that it operates on <<floats>>.
<<atan>> computes the inverse tangent (arc tangent) of the input value.
<<atanf>> is identical to <<atan>>, save that it operates on <<floats>>.
RETURNS
@ifnottex
<<atan>> returns a value in radians, in the range of -pi/2 to pi/2.
<<atan2>> returns a value in radians, in the range of -pi/2 to pi/2.
@end ifnottex
@tex
<<atan>> returns a value in radians, in the range of $-\pi/2$ to $\pi/2$.
<<atan2>> returns a value in radians, in the range of $-\pi/2$ to $\pi/2$.
@end tex
PORTABILITY
<<atan>> is ANSI C. <<atanf>> is an extension.
<<atan2>> is ANSI C. <<atan2f>> is an extension.
*/
/******************************************************************
* Arctangent
*
* Input:
* x - floating point value
*
* Output:
* arctangent of x
*
* Description:
* This routine calculates arctangents.
*
*****************************************************************/
#include <float.h>
#include "fdlibm.h"
#include "zmath.h"
#ifndef _DOUBLE_IS_32BITS
static const double ROOT3 = 1.73205080756887729353;
static const double a[] = { 0.0, 0.52359877559829887308, 1.57079632679489661923,
1.04719755119659774615 };
static const double q[] = { 0.41066306682575781263e+2,
0.86157349597130242515e+2,
0.59578436142597344465e+2,
0.15024001160028576121e+2 };
static const double p[] = { -0.13688768894191926929e+2,
-0.20505855195861651981e+2,
-0.84946240351320683534e+1,
-0.83758299368150059274 };
double
_DEFUN (atangent, (double, double, double, int),
double x _AND
double v _AND
double u _AND
int arctan2)
{
double f, g, R, P, Q, A, res;
int N;
int branch = 0;
int expv, expu;
/* Preparation for calculating arctan2. */
if (arctan2)
{
if (u == 0.0)
if (v == 0.0)
{
errno = ERANGE;
return (z_notanum.d);
}
else
{
branch = 1;
res = __PI_OVER_TWO;
}
if (!branch)
{
int e;
/* Get the exponent values of the inputs. */
g = frexp (v, &expv);
g = frexp (u, &expu);
/* See if a divide will overflow. */
e = expv - expu;
if (e > DBL_MAX_EXP)
{
branch = 1;
res = __PI_OVER_TWO;
}
/* Also check for underflow. */
else if (e < DBL_MIN_EXP)
{
branch = 2;
res = 0.0;
}
}
}
if (!branch)
{
if (arctan2)
f = fabs (v / u);
else
f = fabs (x);
if (f > 1.0)
{
f = 1.0 / f;
N = 2;
}
else
N = 0;
if (f > (2.0 - ROOT3))
{
A = ROOT3 - 1.0;
f = (((A * f - 0.5) - 0.5) + f) / (ROOT3 + f);
N++;
}
/* Check for values that are too small. */
if (-z_rooteps < f && f < z_rooteps)
res = f;
/* Calculate the Taylor series. */
else
{
g = f * f;
P = (((p[3] * g + p[2]) * g + p[1]) * g + p[0]) * g;
Q = (((g + q[3]) * g + q[2]) * g + q[1]) * g + q[0];
R = P / Q;
res = f + f * R;
}
if (N > 1)
res = -res;
res += a[N];
}
if (arctan2)
{
if (u < 0.0 || branch == 2)
res = __PI - res;
if (v < 0.0 || branch == 1)
res = -res;
}
else if (x < 0.0)
{
res = -res;
}
return (res);
}
#endif /* _DOUBLE_IS_32BITS */
|
