|
| 1 | +use super::atan; |
| 2 | +use super::fabs; |
| 3 | + |
| 4 | +const PI: f64 = 3.1415926535897931160E+00; /* 0x400921FB, 0x54442D18 */ |
| 5 | +const PI_LO: f64 = 1.2246467991473531772E-16; /* 0x3CA1A626, 0x33145C07 */ |
| 6 | + |
| 7 | +#[inline] |
| 8 | +pub fn atan2(y: f64, x: f64) -> f64 { |
| 9 | + if x.is_nan() || y.is_nan() { |
| 10 | + return x + y; |
| 11 | + } |
| 12 | + let mut ix = (x.to_bits() >> 32) as u32; |
| 13 | + let lx = x.to_bits() as u32; |
| 14 | + let mut iy = (y.to_bits() >> 32) as u32; |
| 15 | + let ly = y.to_bits() as u32; |
| 16 | + if (ix - 0x3ff00000 | lx) == 0 { |
| 17 | + /* x = 1.0 */ |
| 18 | + return atan(y); |
| 19 | + } |
| 20 | + let m = ((iy >> 31) & 1) | ((ix >> 30) & 2); /* 2*sign(x)+sign(y) */ |
| 21 | + ix &= 0x7fffffff; |
| 22 | + iy &= 0x7fffffff; |
| 23 | + |
| 24 | + /* when y = 0 */ |
| 25 | + if (iy | ly) == 0 { |
| 26 | + return match m { |
| 27 | + 0 | 1 => y, /* atan(+-0,+anything)=+-0 */ |
| 28 | + 2 => PI, /* atan(+0,-anything) = PI */ |
| 29 | + _ => -PI, /* atan(-0,-anything) =-PI */ |
| 30 | + }; |
| 31 | + } |
| 32 | + /* when x = 0 */ |
| 33 | + if (ix | lx) == 0 { |
| 34 | + return if m & 1 != 0 { -PI / 2.0 } else { PI / 2.0 }; |
| 35 | + } |
| 36 | + /* when x is INF */ |
| 37 | + if ix == 0x7ff00000 { |
| 38 | + if iy == 0x7ff00000 { |
| 39 | + return match m { |
| 40 | + 0 => PI / 4.0, /* atan(+INF,+INF) */ |
| 41 | + 1 => -PI / 4.0, /* atan(-INF,+INF) */ |
| 42 | + 2 => 3.0 * PI / 4.0, /* atan(+INF,-INF) */ |
| 43 | + _ => -3.0 * PI / 4.0, /* atan(-INF,-INF) */ |
| 44 | + }; |
| 45 | + } else { |
| 46 | + return match m { |
| 47 | + 0 => 0.0, /* atan(+...,+INF) */ |
| 48 | + 1 => -0.0, /* atan(-...,+INF) */ |
| 49 | + 2 => PI, /* atan(+...,-INF) */ |
| 50 | + _ => -PI, /* atan(-...,-INF) */ |
| 51 | + }; |
| 52 | + } |
| 53 | + } |
| 54 | + /* |y/x| > 0x1p64 */ |
| 55 | + if ix + (64 << 20) < iy || iy == 0x7ff00000 { |
| 56 | + return if m & 1 != 0 { -PI / 2.0 } else { PI / 2.0 }; |
| 57 | + } |
| 58 | + |
| 59 | + /* z = atan(|y/x|) without spurious underflow */ |
| 60 | + let z = if (m & 2 != 0) && iy + (64 << 20) < ix { |
| 61 | + /* |y/x| < 0x1p-64, x<0 */ |
| 62 | + 0.0 |
| 63 | + } else { |
| 64 | + atan(fabs(y / x)) |
| 65 | + }; |
| 66 | + match m { |
| 67 | + 0 => z, /* atan(+,+) */ |
| 68 | + 1 => -z, /* atan(-,+) */ |
| 69 | + 2 => PI - (z - PI_LO), /* atan(+,-) */ |
| 70 | + _ => (z - PI_LO) - PI, /* atan(-,-) */ |
| 71 | + } |
| 72 | +} |
0 commit comments