workbench/libs/mathtrans/spatan.c

   1 /*
   2     Copyright © 1995-2003, The AROS Development Team. All rights reserved.
   3     $Id$
   4 */
   5
   6 #include "mathtrans_intern.h"
   7
   8 /*****************************************************************************
   9
  10     NAME */
  11
  12         AROS_LH1(float, SPAtan,
  13
  14 /*  SYNOPSIS */
  15         AROS_LHA(float, fnum1, D0),
  16
  17 /*  LOCATION */
  18         struct Library *, MathTransBase, 5, MathTrans)
  19
  20 /*  FUNCTION
  21         Calculates the angle of a given number representing the tangent
  22         of that angle. The angle will be in radians.
  23
  24     INPUTS
  25
  26     RESULT
  27         Motorola fast floating point number
  28
  29     BUGS
  30
  31     INTERNALS
  32
  33 *****************************************************************************/
  34 {
  35     AROS_LIBFUNC_INIT
  36
  37     LONG fnumabs = fnum1 & (FFPMantisse_Mask | FFPExponent_Mask);
  38     LONG fnumsquared, fnumcubed;
  39
  40     /* check for +- infinity -> output: +-pi/2 */
  41     if ((LONG)FFP_Pinfty == fnumabs ) return (pio2 | (fnum1 & FFPSign_Mask));
  42
  43     /* atan(0) = 0 */
  44     if (0 == fnumabs)
  45     {
  46         SetSR(Zero_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  47         return 0;
  48     }
  49
  50     /* atan(x>= 128) = pi/2 - 1/x */
  51     if ((BYTE) fnumabs >= 0x48)
  52     {
  53         if (fnumabs == fnum1) /* arg has ppositive sign */
  54         {
  55             SetSR(0, Zero_Bit | Negative_Bit | Overflow_Bit);
  56             return SPSub(SPDiv(fnumabs,one),pio2);
  57         }
  58         else
  59         {
  60             LONG res;
  61             SetSR(Negative_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  62
  63             /* workaround a bug in egcs 1.0.3: complains about illegal operands to | */
  64             res = SPSub(SPDiv(fnumabs,one),pio2);
  65             return res | FFPSign_Mask;
  66         }
  67     }
  68
  69     /* atan(x >= 64) = pi/2 - 1/x +1/(3*x^3) */
  70
  71     fnumsquared = SPMul(fnumabs, fnumabs);
  72
  73     if((BYTE) fnumabs >= 0x47)
  74     {
  75         fnumcubed   = SPMul(fnumabs, fnumsquared);
  76
  77         /* pi/2 - 1/x + 1/(3*x^3) = pi/2 + (1-3*x^2)/(3*x^3)*/
  78         if (fnumabs == fnum1) /* arg has positive sign */
  79         {
  80             LONG res;
  81             SetSR(0, Zero_Bit | Negative_Bit | Overflow_Bit);
  82
  83             /* workaround a bug in egcs 1.0.3: complains about illegal operands to | */
  84             res = SPMul(three, fnumsquared);
  85             return SPAdd
  86             (
  87                 pio2, SPDiv
  88                 (
  89                     SPMul(three, fnumcubed),
  90                     SPAdd(res | FFPSign_Mask, one)
  91                 )
  92             );
  93         }
  94         else
  95         {
  96             LONG res;
  97             SetSR(Negative_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  98
  99             /* workaround a bug in egcs 1.0.3: complains about illegal operands to | */
 100             res = SPMul(three, fnumsquared);
 101             res = SPAdd
 102             (
 103                 pio2, SPDiv
 104                 (
 105                     SPMul(three, fnumcubed),
 106                     SPAdd( res | FFPSign_Mask, one )
 107                 )
 108             );
 109             return res | FFPSign_Mask;
 110         }
 111     }
 112
 113     /* atan(x <= 64) */
 114     return SPAsin(SPDiv(SPSqrt(SPAdd(one,fnumsquared)),fnum1));
 115
 116     AROS_LIBFUNC_EXIT
 117 }