e_atan2f.c revision 2117
10SN/A/* e_atan2f.c -- float version of e_atan2.c.
211928Sserb * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
30SN/A */
40SN/A
50SN/A/*
60SN/A * ====================================================
72362SN/A * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
80SN/A *
92362SN/A * Developed at SunPro, a Sun Microsystems, Inc. business.
100SN/A * Permission to use, copy, modify, and distribute this
110SN/A * software is freely granted, provided that this notice
120SN/A * is preserved.
130SN/A * ====================================================
140SN/A */
150SN/A
160SN/A#ifndef lint
170SN/Astatic char rcsid[] = "$Id: e_atan2f.c,v 1.2 1994/08/18 23:05:11 jtc Exp $";
180SN/A#endif
190SN/A
200SN/A#include "math.h"
212362SN/A#include "math_private.h"
222362SN/A
232362SN/A#ifdef __STDC__
240SN/Astatic const float
250SN/A#else
260SN/Astatic float
270SN/A#endif
280SN/Atiny  = 1.0e-30,
290SN/Azero  = 0.0,
300SN/Api_o_4  = 7.8539818525e-01, /* 0x3f490fdb */
310SN/Api_o_2  = 1.5707963705e+00, /* 0x3fc90fdb */
320SN/Api      = 3.1415925026e+00, /* 0x40490fda */
330SN/Api_lo   = 1.5099578832e-07; /* 0x34222168 */
340SN/A
350SN/A#ifdef __STDC__
360SN/A	float __ieee754_atan2f(float y, float x)
375428SN/A#else
380SN/A	float __ieee754_atan2f(y,x)
3913629Savstepan	float  y,x;
400SN/A#endif
410SN/A{
420SN/A	float z;
438926SN/A	int32_t k,m,hx,hy,ix,iy;
440SN/A
450SN/A	GET_FLOAT_WORD(hx,x);
460SN/A	ix = hx&0x7fffffff;
470SN/A	GET_FLOAT_WORD(hy,y);
480SN/A	iy = hy&0x7fffffff;
490SN/A	if((ix>0x7f800000)||
500SN/A	   (iy>0x7f800000))	/* x or y is NaN */
510SN/A	   return x+y;
520SN/A	if(hx==0x3f800000) return atanf(y);   /* x=1.0 */
530SN/A	m = ((hy>>31)&1)|((hx>>30)&2);	/* 2*sign(x)+sign(y) */
540SN/A
558234SN/A    /* when y = 0 */
560SN/A	if(iy==0) {
5713629Savstepan	    switch(m) {
580SN/A		case 0:
590SN/A		case 1: return y; 	/* atan(+-0,+anything)=+-0 */
600SN/A		case 2: return  pi+tiny;/* atan(+0,-anything) = pi */
610SN/A		case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */
6213629Savstepan	    }
6313629Savstepan	}
640SN/A    /* when x = 0 */
6513629Savstepan	if(ix==0) return (hy<0)?  -pi_o_2-tiny: pi_o_2+tiny;
660SN/A
6713629Savstepan    /* when x is INF */
689680SN/A	if(ix==0x7f800000) {
690SN/A	    if(iy==0x7f800000) {
700SN/A		switch(m) {
7110071SN/A		    case 0: return  pi_o_4+tiny;/* atan(+INF,+INF) */
720SN/A		    case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */
730SN/A		    case 2: return  (float)3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/
740SN/A		    case 3: return (float)-3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/
7513629Savstepan		}
7613629Savstepan	    } else {
770SN/A		switch(m) {
780SN/A		    case 0: return  zero  ;	/* atan(+...,+INF) */
790SN/A		    case 1: return -zero  ;	/* atan(-...,+INF) */
800SN/A		    case 2: return  pi+tiny  ;	/* atan(+...,-INF) */
810SN/A		    case 3: return -pi-tiny  ;	/* atan(-...,-INF) */
820SN/A		}
830SN/A	    }
840SN/A	}
856447SN/A    /* when y is INF */
860SN/A	if(iy==0x7f800000) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
870SN/A
8813629Savstepan    /* compute y/x */
890SN/A	k = (iy-ix)>>23;
900SN/A	if(k > 60) z=pi_o_2+(float)0.5*pi_lo; 	/* |y/x| >  2**60 */
910SN/A	else if(hx<0&&k<-60) z=0.0; 	/* |y|/x < -2**60 */
920SN/A	else z=atanf(fabsf(y/x));	/* safe to do y/x */
930SN/A	switch (m) {
940SN/A	    case 0: return       z  ;	/* atan(+,+) */
950SN/A	    case 1: {
960SN/A	    	      u_int32_t zh;
970SN/A		      GET_FLOAT_WORD(zh,z);
980SN/A		      SET_FLOAT_WORD(z,zh ^ 0x80000000);
990SN/A		    }
1000SN/A		    return       z  ;	/* atan(-,+) */
1010SN/A	    case 2: return  pi-(z-pi_lo);/* atan(+,-) */
1020SN/A	    default: /* case 3 */
1030SN/A	    	    return  (z-pi_lo)-pi;/* atan(-,-) */
1045906SN/A	}
1055906SN/A}
1065906SN/A