workbench/libs/mathieeesingtrans/ieeespcos.c

   1 /*
   2     Copyright © 1995-2003, The AROS Development Team. All rights reserved.
   3     $Id$
   4 */
   5
   6 #include "mathieeesingtrans_intern.h"
   7
   8 /*
   9     FUNCTION
  10       Calculate the cosine of a given IEEE single precision number in radians
  11
  12     RESULT
  13       IEEE single precision floating point number
  14
  15       flags:
  16         zero     : result is zero
  17         negative : result is negative
  18         overflow : 0
  19
  20     NOTES
  21
  22     EXAMPLE
  23
  24     BUGS
  25
  26     SEE ALSO
  27
  28     INTERNALS
  29         Algorithm for Calculation of cos(y):
  30         <code>
  31          z    = floor ( |y| / pi );
  32          y_1  = |y| - z * pi;        => 0 <= y_1 < pi
  33
  34          if (y_1 > pi/2 ) then y_1 = pi - y_1;
  35
  36          => 0 <= y_1 < pi/2
  37
  38          Res = 1 - y^2/2! + y^4/4! - y^6/6! + y^8/8! - y^10/10! =
  39              = 1 -(y^2(-1/2!+y^2(1/4!+y^2(-1/6!+y^2(1/8!-1/10!y^2)))));
  40
  41          if (z was an odd number)
  42            Res = -Res;
  43
  44          if (y_1 was greater than pi/2 in the test above)
  45            Res = -Res;
  46       </code>
  47
  48
  49     HISTORY
  50 */
  51
  52 AROS_LH1(float, IEEESPCos,
  53     AROS_LHA(float, y, D0),
  54     struct Library *, MathIeeeSingTransBase, 7, MathIeeeSingTrans
  55 )
  56 {
  57     AROS_LIBFUNC_INIT
  58
  59     LONG z,Res,ysquared,yabs,tmp;
  60     yabs = y & (IEEESPMantisse_Mask + IEEESPExponent_Mask);
  61
  62     if (IEEESP_Pinfty == yabs)
  63     {
  64         SetSR(Overflow_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  65         return IEEESP_NAN;
  66     }
  67
  68     z = IEEESPFloor(IEEESPDiv(yabs, pi));
  69     tmp  = IEEESPMul(z,pi);
  70     tmp |= IEEESPSign_Mask; /* tmp = -tmp; */
  71     yabs = IEEESPAdd(yabs, tmp);
  72     if (yabs > pio2)
  73     {
  74         yabs |= IEEESPSign_Mask;
  75         yabs  =IEEESPAdd(pi, yabs);
  76         tmp = TRUE;
  77     }
  78     else
  79     {
  80         tmp = FALSE;
  81     }
  82
  83     ysquared = IEEESPMul(yabs,yabs);
  84     Res = IEEESPAdd(cosf1,
  85           IEEESPMul(ysquared,
  86           IEEESPAdd(cosf2,
  87           IEEESPMul(ysquared,
  88           IEEESPAdd(cosf3,
  89           IEEESPMul(ysquared,
  90           IEEESPAdd(cosf4,
  91           IEEESPMul(ysquared,
  92           IEEESPAdd(cosf5,
  93           IEEESPMul(ysquared, cosf6))))))))));
  94
  95     if (0 == Res)
  96     {
  97         SetSR(Zero_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  98         return 0;
  99     }
 100
 101     if (TRUE == intern_IEEESPisodd(z)) Res ^= IEEESPSign_Mask;
 102     if (TRUE == tmp)                   Res ^= IEEESPSign_Mask;
 103
 104     if (Res < 0) SetSR(Negative_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
 105
 106     return Res;
 107
 108     AROS_LIBFUNC_EXIT
 109 }