| @@ -2,27 +2,27 @@ | | | @@ -2,27 +2,27 @@ |
2 | /* | | 2 | /* |
3 | * ==================================================== | | 3 | * ==================================================== |
4 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | | 4 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. |
5 | * | | 5 | * |
6 | * Developed at SunPro, a Sun Microsystems, Inc. business. | | 6 | * Developed at SunPro, a Sun Microsystems, Inc. business. |
7 | * Permission to use, copy, modify, and distribute this | | 7 | * Permission to use, copy, modify, and distribute this |
8 | * software is freely granted, provided that this notice | | 8 | * software is freely granted, provided that this notice |
9 | * is preserved. | | 9 | * is preserved. |
10 | * ==================================================== | | 10 | * ==================================================== |
11 | */ | | 11 | */ |
12 | | | 12 | |
13 | #include <sys/cdefs.h> | | 13 | #include <sys/cdefs.h> |
14 | #if defined(LIBM_SCCS) && !defined(lint) | | 14 | #if defined(LIBM_SCCS) && !defined(lint) |
15 | __RCSID("$NetBSD: e_atan2.c,v 1.12 2002/05/26 22:01:48 wiz Exp $"); | | 15 | __RCSID("$NetBSD: e_atan2.c,v 1.13 2018/03/10 09:44:47 eadler Exp $"); |
16 | #endif | | 16 | #endif |
17 | | | 17 | |
18 | /* __ieee754_atan2(y,x) | | 18 | /* __ieee754_atan2(y,x) |
19 | * Method : | | 19 | * Method : |
20 | * 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x). | | 20 | * 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x). |
21 | * 2. Reduce x to positive by (if x and y are unexceptional): | | 21 | * 2. Reduce x to positive by (if x and y are unexceptional): |
22 | * ARG (x+iy) = arctan(y/x) ... if x > 0, | | 22 | * ARG (x+iy) = arctan(y/x) ... if x > 0, |
23 | * ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0, | | 23 | * ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0, |
24 | * | | 24 | * |
25 | * Special cases: | | 25 | * Special cases: |
26 | * | | 26 | * |
27 | * ATAN2((anything), NaN ) is NaN; | | 27 | * ATAN2((anything), NaN ) is NaN; |
28 | * ATAN2(NAN , (anything) ) is NaN; | | 28 | * ATAN2(NAN , (anything) ) is NaN; |
| @@ -57,27 +57,27 @@ double | | | @@ -57,27 +57,27 @@ double |
57 | __ieee754_atan2(double y, double x) | | 57 | __ieee754_atan2(double y, double x) |
58 | { | | 58 | { |
59 | double z; | | 59 | double z; |
60 | int32_t k,m,hx,hy,ix,iy; | | 60 | int32_t k,m,hx,hy,ix,iy; |
61 | u_int32_t lx,ly; | | 61 | u_int32_t lx,ly; |
62 | | | 62 | |
63 | EXTRACT_WORDS(hx,lx,x); | | 63 | EXTRACT_WORDS(hx,lx,x); |
64 | ix = hx&0x7fffffff; | | 64 | ix = hx&0x7fffffff; |
65 | EXTRACT_WORDS(hy,ly,y); | | 65 | EXTRACT_WORDS(hy,ly,y); |
66 | iy = hy&0x7fffffff; | | 66 | iy = hy&0x7fffffff; |
67 | if(((ix|((lx|-lx)>>31))>0x7ff00000)|| | | 67 | if(((ix|((lx|-lx)>>31))>0x7ff00000)|| |
68 | ((iy|((ly|-ly)>>31))>0x7ff00000)) /* x or y is NaN */ | | 68 | ((iy|((ly|-ly)>>31))>0x7ff00000)) /* x or y is NaN */ |
69 | return x+y; | | 69 | return x+y; |
70 | if(((hx-0x3ff00000)|lx)==0) return atan(y); /* x=1.0 */ | | 70 | if(hx==0x3ff00000&&lx==0) return atan(y); /* x=1.0 */ |
71 | m = ((hy>>31)&1)|((hx>>30)&2); /* 2*sign(x)+sign(y) */ | | 71 | m = ((hy>>31)&1)|((hx>>30)&2); /* 2*sign(x)+sign(y) */ |
72 | | | 72 | |
73 | /* when y = 0 */ | | 73 | /* when y = 0 */ |
74 | if((iy|ly)==0) { | | 74 | if((iy|ly)==0) { |
75 | switch(m) { | | 75 | switch(m) { |
76 | case 0: | | 76 | case 0: |
77 | case 1: return y; /* atan(+-0,+anything)=+-0 */ | | 77 | case 1: return y; /* atan(+-0,+anything)=+-0 */ |
78 | case 2: return pi+tiny;/* atan(+0,-anything) = pi */ | | 78 | case 2: return pi+tiny;/* atan(+0,-anything) = pi */ |
79 | case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */ | | 79 | case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */ |
80 | } | | 80 | } |
81 | } | | 81 | } |
82 | /* when x = 0 */ | | 82 | /* when x = 0 */ |
83 | if((ix|lx)==0) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny; | | 83 | if((ix|lx)==0) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny; |