1 /* origin: FreeBSD /usr/src/lib/msun/src/e_atan2l.c */
3 * ====================================================
4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
6 * Developed at SunSoft, a Sun Microsystems, Inc. business.
7 * Permission to use, copy, modify, and distribute this
8 * software is freely granted, provided that this notice
10 * ====================================================
14 * See comments in atan2.c.
15 * Converted to long double by David Schultz <das@FreeBSD.ORG>.
20 #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
21 long double atan2l(long double y, long double x)
25 #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
26 #include "__invtrigl.h"
28 long double atan2l(long double y, long double x)
30 union IEEEl2bits ux, uy;
33 uint16_t exptx, expsignx, expty, expsigny;
35 if (isnan(x) || isnan(y))
40 expsigny = uy.xbits.expsign;
41 expty = expsigny & 0x7fff;
43 expsignx = ux.xbits.expsign;
44 exptx = expsignx & 0x7fff;
45 m = ((expsigny>>15)&1) | ((expsignx>>14)&2); /* 2*sign(x)+sign(y) */
50 case 1: return y; /* atan(+-0,+anything)=+-0 */
51 case 2: return 2*pio2_hi+0x1p-120f; /* atan(+0,-anything) = pi */
52 case 3: return -2*pio2_hi-0x1p-120f; /* atan(-0,-anything) =-pi */
56 return m&1 ? -pio2_hi-0x1p-120f : pio2_hi+0x1p-120f;
58 if (exptx == 0x7fff) {
59 if (expty == 0x7fff) {
61 case 0: return pio2_hi*0.5+0x1p-120f; /* atan(+INF,+INF) */
62 case 1: return -pio2_hi*0.5-0x1p-120f; /* atan(-INF,+INF) */
63 case 2: return 1.5*pio2_hi+0x1p-120f; /* atan(+INF,-INF) */
64 case 3: return -1.5*pio2_hi-0x1p-120f; /* atan(-INF,-INF) */
68 case 0: return 0.0; /* atan(+...,+INF) */
69 case 1: return -0.0; /* atan(-...,+INF) */
70 case 2: return 2*pio2_hi+0x1p-120f; /* atan(+...,-INF) */
71 case 3: return -2*pio2_hi-0x1p-120f; /* atan(-...,-INF) */
76 if (exptx+120 < expty || expty == 0x7fff)
77 return m&1 ? -pio2_hi-0x1p-120f : pio2_hi+0x1p-120f;
79 if ((m&2) && expty+120 < exptx) /* |y/x| tiny, x<0 */
82 z = atanl(fabsl(y/x));
84 case 0: return z; /* atan(+,+) */
85 case 1: return -z; /* atan(-,+) */
86 case 2: return 2*pio2_hi-(z-2*pio2_lo); /* atan(+,-) */
88 return (z-2*pio2_lo)-2*pio2_hi; /* atan(-,-) */