OLD | NEW |
1 /* origin: FreeBSD /usr/src/lib/msun/src/e_atan2f.c */ | 1 /* origin: FreeBSD /usr/src/lib/msun/src/e_atan2f.c */ |
2 /* | 2 /* |
3 * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com. | 3 * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com. |
4 */ | 4 */ |
5 /* | 5 /* |
6 * ==================================================== | 6 * ==================================================== |
7 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | 7 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. |
8 * | 8 * |
9 * Developed at SunPro, a Sun Microsystems, Inc. business. | 9 * Developed at SunPro, a Sun Microsystems, Inc. business. |
10 * Permission to use, copy, modify, and distribute this | 10 * Permission to use, copy, modify, and distribute this |
11 * software is freely granted, provided that this notice | 11 * software is freely granted, provided that this notice |
12 * is preserved. | 12 * is preserved. |
13 * ==================================================== | 13 * ==================================================== |
14 */ | 14 */ |
15 | 15 |
16 #include "libm.h" | 16 #include "libm.h" |
17 | 17 |
18 static const float | 18 static const float pi = 3.1415927410e+00, /* 0x40490fdb */ |
19 pi = 3.1415927410e+00, /* 0x40490fdb */ | 19 pi_lo = -8.7422776573e-08; /* 0xb3bbbd2e */ |
20 pi_lo = -8.7422776573e-08; /* 0xb3bbbd2e */ | |
21 | 20 |
22 float atan2f(float y, float x) | 21 float atan2f(float y, float x) { |
23 { | 22 float z; |
24 » float z; | 23 uint32_t m, ix, iy; |
25 » uint32_t m,ix,iy; | |
26 | 24 |
27 » if (isnan(x) || isnan(y)) | 25 if (isnan(x) || isnan(y)) |
28 » » return x+y; | 26 return x + y; |
29 » GET_FLOAT_WORD(ix, x); | 27 GET_FLOAT_WORD(ix, x); |
30 » GET_FLOAT_WORD(iy, y); | 28 GET_FLOAT_WORD(iy, y); |
31 » if (ix == 0x3f800000) /* x=1.0 */ | 29 if (ix == 0x3f800000) /* x=1.0 */ |
32 » » return atanf(y); | 30 return atanf(y); |
33 » m = ((iy>>31)&1) | ((ix>>30)&2); /* 2*sign(x)+sign(y) */ | 31 m = ((iy >> 31) & 1) | ((ix >> 30) & 2); /* 2*sign(x)+sign(y) */ |
34 » ix &= 0x7fffffff; | 32 ix &= 0x7fffffff; |
35 » iy &= 0x7fffffff; | 33 iy &= 0x7fffffff; |
36 | 34 |
37 » /* when y = 0 */ | 35 /* when y = 0 */ |
38 » if (iy == 0) { | 36 if (iy == 0) { |
39 » » switch (m) { | 37 switch (m) { |
40 » » case 0: | 38 case 0: |
41 » » case 1: return y; /* atan(+-0,+anything)=+-0 */ | 39 case 1: |
42 » » case 2: return pi; /* atan(+0,-anything) = pi */ | 40 return y; /* atan(+-0,+anything)=+-0 */ |
43 » » case 3: return -pi; /* atan(-0,-anything) =-pi */ | 41 case 2: |
44 » » } | 42 return pi; /* atan(+0,-anything) = pi */ |
45 » } | 43 case 3: |
46 » /* when x = 0 */ | 44 return -pi; /* atan(-0,-anything) =-pi */ |
47 » if (ix == 0) | 45 } |
48 » » return m&1 ? -pi/2 : pi/2; | 46 } |
49 » /* when x is INF */ | 47 /* when x = 0 */ |
50 » if (ix == 0x7f800000) { | 48 if (ix == 0) |
51 » » if (iy == 0x7f800000) { | 49 return m & 1 ? -pi / 2 : pi / 2; |
52 » » » switch (m) { | 50 /* when x is INF */ |
53 » » » case 0: return pi/4; /* atan(+INF,+INF) */ | 51 if (ix == 0x7f800000) { |
54 » » » case 1: return -pi/4; /* atan(-INF,+INF) */ | 52 if (iy == 0x7f800000) { |
55 » » » case 2: return 3*pi/4; /*atan(+INF,-INF)*/ | 53 switch (m) { |
56 » » » case 3: return -3*pi/4; /*atan(-INF,-INF)*/ | 54 case 0: |
57 » » » } | 55 return pi / 4; /* atan(+INF,+INF) */ |
58 » » } else { | 56 case 1: |
59 » » » switch (m) { | 57 return -pi / 4; /* atan(-INF,+INF) */ |
60 » » » case 0: return 0.0f; /* atan(+...,+INF) */ | 58 case 2: |
61 » » » case 1: return -0.0f; /* atan(-...,+INF) */ | 59 return 3 * pi / 4; /*atan(+INF,-INF)*/ |
62 » » » case 2: return pi; /* atan(+...,-INF) */ | 60 case 3: |
63 » » » case 3: return -pi; /* atan(-...,-INF) */ | 61 return -3 * pi / 4; /*atan(-INF,-INF)*/ |
64 » » » } | 62 } |
65 » » } | 63 } else { |
66 » } | 64 switch (m) { |
67 » /* |y/x| > 0x1p26 */ | 65 case 0: |
68 » if (ix+(26<<23) < iy || iy == 0x7f800000) | 66 return 0.0f; /* atan(+...,+INF) */ |
69 » » return m&1 ? -pi/2 : pi/2; | 67 case 1: |
| 68 return -0.0f; /* atan(-...,+INF) */ |
| 69 case 2: |
| 70 return pi; /* atan(+...,-INF) */ |
| 71 case 3: |
| 72 return -pi; /* atan(-...,-INF) */ |
| 73 } |
| 74 } |
| 75 } |
| 76 /* |y/x| > 0x1p26 */ |
| 77 if (ix + (26 << 23) < iy || iy == 0x7f800000) |
| 78 return m & 1 ? -pi / 2 : pi / 2; |
70 | 79 |
71 » /* z = atan(|y/x|) with correct underflow */ | 80 /* z = atan(|y/x|) with correct underflow */ |
72 » if ((m&2) && iy+(26<<23) < ix) /*|y/x| < 0x1p-26, x < 0 */ | 81 if ((m & 2) && iy + (26 << 23) < ix) /*|y/x| < 0x1p-26, x < 0 */ |
73 » » z = 0.0; | 82 z = 0.0; |
74 » else | 83 else |
75 » » z = atanf(fabsf(y/x)); | 84 z = atanf(fabsf(y / x)); |
76 » switch (m) { | 85 switch (m) { |
77 » case 0: return z; /* atan(+,+) */ | 86 case 0: |
78 » case 1: return -z; /* atan(-,+) */ | 87 return z; /* atan(+,+) */ |
79 » case 2: return pi - (z-pi_lo); /* atan(+,-) */ | 88 case 1: |
80 » default: /* case 3 */ | 89 return -z; /* atan(-,+) */ |
81 » » return (z-pi_lo) - pi; /* atan(-,-) */ | 90 case 2: |
82 » } | 91 return pi - (z - pi_lo); /* atan(+,-) */ |
| 92 default: /* case 3 */ |
| 93 return (z - pi_lo) - pi; /* atan(-,-) */ |
| 94 } |
83 } | 95 } |
OLD | NEW |