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