![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to atan2.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [CSRG BSD Unix] / 43BSDReno / lib / libm / mc68881|
Return to atan2.c CVS log | Up to [CSRG BSD Unix] / 43BSDReno / lib / libm / mc68881 |
1.1 root 1: /*-
2: * Copyright (c) 1990 The Regents of the University of California.
3: * All rights reserved.
4: *
5: * This code is derived from software contributed to Berkeley by
6: * the Systems Programming Group of the University of Utah Computer
7: * Science Department.
8: *
9: * Redistribution and use in source and binary forms are permitted
10: * provided that: (1) source distributions retain this entire copyright
11: * notice and comment, and (2) distributions including binaries display
12: * the following acknowledgement: ``This product includes software
13: * developed by the University of California, Berkeley and its contributors''
14: * in the documentation or other materials provided with the distribution
15: * and in all advertising materials mentioning features or use of this
16: * software. Neither the name of the University nor the names of its
17: * contributors may be used to endorse or promote products derived
18: * from this software without specific prior written permission.
19: * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
20: * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
21: * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
22: */
23:
24: #ifndef lint
25: static char sccsid[] = "@(#)atan2.c 5.1 (Berkeley) 5/17/90";
26: #endif /* not lint */
27:
28: /*
29: * ATAN2(Y,X)
30: * RETURN ARG (X+iY)
31: * DOUBLE PRECISION (IEEE DOUBLE 53 BITS)
32: *
33: * Scaled down version to weed out special cases. "Normal" cases are
34: * handled by calling atan2__A(), an assembly coded support routine in
35: * support.s.
36: *
37: * Required system supported functions :
38: * copysign(x,y)
39: * atan2__A(y,x)
40: *
41: * Method :
42: * 1. Deal with special cases
43: * 2. Call atan2__A() to do the others
44: *
45: * Special cases:
46: * Notations: atan2(y,x) == ARG (x+iy) == ARG(x,y).
47: *
48: * ARG( NAN , (anything) ) is NaN;
49: * ARG( (anything), NaN ) is NaN;
50: * ARG(+(anything but NaN), +-0) is +-0 ;
51: * ARG(-(anything but NaN), +-0) is +-PI ;
52: * ARG( 0, +-(anything but 0 and NaN) ) is +-PI/2;
53: * ARG( +INF,+-(anything but INF and NaN) ) is +-0 ;
54: * ARG( -INF,+-(anything but INF and NaN) ) is +-PI;
55: * ARG( +INF,+-INF ) is +-PI/4 ;
56: * ARG( -INF,+-INF ) is +-3PI/4;
57: * ARG( (anything but,0,NaN, and INF),+-INF ) is +-PI/2;
58: *
59: * Accuracy:
60: * atan2(y,x) returns (PI/pi) * the exact ARG (x+iy) nearly rounded,
61: * where
62: *
63: * in decimal:
64: * pi = 3.141592653589793 23846264338327 .....
65: * 53 bits PI = 3.141592653589793 115997963 ..... ,
66: * 56 bits PI = 3.141592653589793 227020265 ..... ,
67: *
68: * in hexadecimal:
69: * pi = 3.243F6A8885A308D313198A2E....
70: * 53 bits PI = 3.243F6A8885A30 = 2 * 1.921FB54442D18 error=.276ulps
71: * 56 bits PI = 3.243F6A8885A308 = 4 * .C90FDAA22168C2 error=.206ulps
72: *
73: * In a test run with 356,000 random argument on [-1,1] * [-1,1] on a
74: * VAX, the maximum observed error was 1.41 ulps (units of the last place)
75: * compared with (PI/pi)*(the exact ARG(x+iy)).
76: *
77: * Note:
78: * We use machine PI (the true pi rounded) in place of the actual
79: * value of pi for all the trig and inverse trig functions. In general,
80: * if trig is one of sin, cos, tan, then computed trig(y) returns the
81: * exact trig(y*pi/PI) nearly rounded; correspondingly, computed arctrig
82: * returns the exact arctrig(y)*PI/pi nearly rounded. These guarantee the
83: * trig functions have period PI, and trig(arctrig(x)) returns x for
84: * all critical values x.
85: *
86: * Constants:
87: * The hexadecimal values are the intended ones for the following constants.
88: * The decimal values may be used, provided that the compiler will convert
89: * from decimal to binary accurately enough to produce the hexadecimal values
90: * shown.
91: */
92:
93: static double
94: PIo4 = 7.8539816339744827900E-1 , /*Hex 2^ -1 * 1.921FB54442D18 */
95: PIo2 = 1.5707963267948965580E0 , /*Hex 2^ 0 * 1.921FB54442D18 */
96: PI = 3.1415926535897931160E0 ; /*Hex 2^ 1 * 1.921FB54442D18 */
97:
98: double atan2(y,x)
99: double y,x;
100: {
101: static double zero=0, one=1;
102: double copysign(),atan2__A(),signy,signx;
103: int finite();
104:
105: /* if x or y is NAN */
106: if(x!=x) return(x); if(y!=y) return(y);
107:
108: /* copy down the sign of y and x */
109: signy = copysign(one,y);
110: signx = copysign(one,x);
111:
112: /* when y = 0 */
113: if(y==zero) return((signx==one)?y:copysign(PI,signy));
114:
115: /* when x = 0 */
116: if(x==zero) return(copysign(PIo2,signy));
117:
118: /* when x is INF */
119: if(!finite(x))
120: if(!finite(y))
121: return(copysign((signx==one)?PIo4:3*PIo4,signy));
122: else
123: return(copysign((signx==one)?zero:PI,signy));
124:
125: /* when y is INF */
126: if(!finite(y)) return(copysign(PIo2,signy));
127:
128: /* else let atan2__A do the work */
129: return(atan2__A(y,x));
130: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.