src/libm/e_atan2.c
author Ryan C. Gordon <icculus@icculus.org>
Thu, 21 Apr 2016 03:16:44 -0400
changeset 11729 d1ce8396c356
parent 11685 1788f503254d
permissions -rw-r--r--
Initial shot at a renderer target for Apple's Metal API.

This isn't complete, but is enough to run testsprite2. It's currently
Mac-only; with a little work to figure out how to properly glue in a Metal
layer to a UIView, this will likely work on iOS, too.

This is only wired up to the configure script right now, and disabled by
default. CMake and Xcode still need their bits filled in as appropriate.
slouken@4870
     1
/*
slouken@4870
     2
 * ====================================================
slouken@4870
     3
 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
slouken@4870
     4
 *
slouken@4870
     5
 * Developed at SunPro, a Sun Microsystems, Inc. business.
slouken@4870
     6
 * Permission to use, copy, modify, and distribute this
slouken@4870
     7
 * software is freely granted, provided that this notice
slouken@4870
     8
 * is preserved.
slouken@4870
     9
 * ====================================================
slouken@4870
    10
 */
slouken@4870
    11
slouken@4870
    12
/* __ieee754_atan2(y,x)
slouken@4870
    13
 * Method :
slouken@4870
    14
 *	1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
slouken@4870
    15
 *	2. Reduce x to positive by (if x and y are unexceptional):
slouken@4870
    16
 *		ARG (x+iy) = arctan(y/x)   	   ... if x > 0,
slouken@4870
    17
 *		ARG (x+iy) = pi - arctan[y/(-x)]   ... if x < 0,
slouken@4870
    18
 *
slouken@4870
    19
 * Special cases:
slouken@4870
    20
 *
slouken@4870
    21
 *	ATAN2((anything), NaN ) is NaN;
slouken@4870
    22
 *	ATAN2(NAN , (anything) ) is NaN;
slouken@4870
    23
 *	ATAN2(+-0, +(anything but NaN)) is +-0  ;
slouken@4870
    24
 *	ATAN2(+-0, -(anything but NaN)) is +-pi ;
slouken@4870
    25
 *	ATAN2(+-(anything but 0 and NaN), 0) is +-pi/2;
slouken@4870
    26
 *	ATAN2(+-(anything but INF and NaN), +INF) is +-0 ;
slouken@4870
    27
 *	ATAN2(+-(anything but INF and NaN), -INF) is +-pi;
slouken@4870
    28
 *	ATAN2(+-INF,+INF ) is +-pi/4 ;
slouken@4870
    29
 *	ATAN2(+-INF,-INF ) is +-3pi/4;
slouken@4870
    30
 *	ATAN2(+-INF, (anything but,0,NaN, and INF)) is +-pi/2;
slouken@4870
    31
 *
slouken@4870
    32
 * Constants:
slouken@4870
    33
 * The hexadecimal values are the intended ones for the following
slouken@4870
    34
 * constants. The decimal values may be used, provided that the
slouken@4870
    35
 * compiler will convert from decimal to binary accurately enough
slouken@4870
    36
 * to produce the hexadecimal values shown.
slouken@4870
    37
 */
slouken@4870
    38
slouken@6044
    39
#include "math_libm.h"
slouken@4870
    40
#include "math_private.h"
slouken@4870
    41
slouken@4870
    42
static const double
slouken@4870
    43
tiny  = 1.0e-300,
slouken@4870
    44
zero  = 0.0,
slouken@4870
    45
pi_o_4  = 7.8539816339744827900E-01, /* 0x3FE921FB, 0x54442D18 */
slouken@4870
    46
pi_o_2  = 1.5707963267948965580E+00, /* 0x3FF921FB, 0x54442D18 */
slouken@4870
    47
pi      = 3.1415926535897931160E+00, /* 0x400921FB, 0x54442D18 */
slouken@4870
    48
pi_lo   = 1.2246467991473531772E-16; /* 0x3CA1A626, 0x33145C07 */
slouken@4870
    49
slouken@4870
    50
double attribute_hidden __ieee754_atan2(double y, double x)
slouken@4870
    51
{
slouken@4870
    52
	double z;
slouken@4870
    53
	int32_t k,m,hx,hy,ix,iy;
slouken@4870
    54
	u_int32_t lx,ly;
slouken@4870
    55
slouken@4870
    56
	EXTRACT_WORDS(hx,lx,x);
slouken@4870
    57
	ix = hx&0x7fffffff;
slouken@4870
    58
	EXTRACT_WORDS(hy,ly,y);
slouken@4870
    59
	iy = hy&0x7fffffff;
slouken@11685
    60
	if(((ix|((lx|-(int32_t)lx)>>31))>0x7ff00000)||
slouken@11685
    61
	   ((iy|((ly|-(int32_t)ly)>>31))>0x7ff00000))	/* x or y is NaN */
slouken@4870
    62
	   return x+y;
slouken@4870
    63
	if(((hx-0x3ff00000)|lx)==0) return atan(y);   /* x=1.0 */
slouken@4870
    64
	m = ((hy>>31)&1)|((hx>>30)&2);	/* 2*sign(x)+sign(y) */
slouken@4870
    65
slouken@4870
    66
    /* when y = 0 */
slouken@4870
    67
	if((iy|ly)==0) {
slouken@4870
    68
	    switch(m) {
slouken@4870
    69
		case 0:
slouken@4870
    70
		case 1: return y; 	/* atan(+-0,+anything)=+-0 */
slouken@4870
    71
		case 2: return  pi+tiny;/* atan(+0,-anything) = pi */
slouken@4870
    72
		case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */
slouken@4870
    73
	    }
slouken@4870
    74
	}
slouken@4870
    75
    /* when x = 0 */
slouken@4870
    76
	if((ix|lx)==0) return (hy<0)?  -pi_o_2-tiny: pi_o_2+tiny;
slouken@4870
    77
slouken@4870
    78
    /* when x is INF */
slouken@4870
    79
	if(ix==0x7ff00000) {
slouken@4870
    80
	    if(iy==0x7ff00000) {
slouken@4870
    81
		switch(m) {
slouken@4870
    82
		    case 0: return  pi_o_4+tiny;/* atan(+INF,+INF) */
slouken@4870
    83
		    case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */
slouken@11683
    84
		    case 2: return  3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/
slouken@11683
    85
		    case 3: return -3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/
slouken@4870
    86
		}
slouken@4870
    87
	    } else {
slouken@4870
    88
		switch(m) {
slouken@4870
    89
		    case 0: return  zero  ;	/* atan(+...,+INF) */
slouken@4870
    90
		    case 1: return -zero  ;	/* atan(-...,+INF) */
slouken@4870
    91
		    case 2: return  pi+tiny  ;	/* atan(+...,-INF) */
slouken@4870
    92
		    case 3: return -pi-tiny  ;	/* atan(-...,-INF) */
slouken@4870
    93
		}
slouken@4870
    94
	    }
slouken@4870
    95
	}
slouken@4870
    96
    /* when y is INF */
slouken@4870
    97
	if(iy==0x7ff00000) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
slouken@4870
    98
slouken@4870
    99
    /* compute y/x */
slouken@4870
   100
	k = (iy-ix)>>20;
slouken@4870
   101
	if(k > 60) z=pi_o_2+0.5*pi_lo; 	/* |y/x| >  2**60 */
slouken@4870
   102
	else if(hx<0&&k<-60) z=0.0; 	/* |y|/x < -2**60 */
slouken@4870
   103
	else z=atan(fabs(y/x));		/* safe to do y/x */
slouken@4870
   104
	switch (m) {
slouken@4870
   105
	    case 0: return       z  ;	/* atan(+,+) */
slouken@4870
   106
	    case 1: {
slouken@4870
   107
	    	      u_int32_t zh;
slouken@4870
   108
		      GET_HIGH_WORD(zh,z);
slouken@4870
   109
		      SET_HIGH_WORD(z,zh ^ 0x80000000);
slouken@4870
   110
		    }
slouken@4870
   111
		    return       z  ;	/* atan(-,+) */
slouken@4870
   112
	    case 2: return  pi-(z-pi_lo);/* atan(+,-) */
slouken@4870
   113
	    default: /* case 3 */
slouken@4870
   114
	    	    return  (z-pi_lo)-pi;/* atan(-,-) */
slouken@4870
   115
	}
slouken@4870
   116
}
slouken@11683
   117
slouken@11683
   118
/*
slouken@11683
   119
 * wrapper atan2(y,x)
slouken@11683
   120
 */
slouken@11683
   121
#ifndef _IEEE_LIBM
slouken@11683
   122
double atan2(double y, double x)
slouken@11683
   123
{
slouken@11683
   124
	double z = __ieee754_atan2(y, x);
slouken@11683
   125
	if (_LIB_VERSION == _IEEE_ || isnan(x) || isnan(y))
slouken@11683
   126
		return z;
slouken@11683
   127
	if (x == 0.0 && y == 0.0)
slouken@11683
   128
		return __kernel_standard(y,x,3); /* atan2(+-0,+-0) */
slouken@11683
   129
	return z;
slouken@11683
   130
}
slouken@11683
   131
#else
slouken@11683
   132
strong_alias(__ieee754_atan2, atan2)
slouken@11683
   133
#endif
slouken@11683
   134
libm_hidden_def(atan2)