ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/tcos.c
(Generate patch)

Comparing ray/src/common/tcos.c (file contents):
Revision 3.1 by gwlarson, Wed Dec 16 18:15:07 1998 UTC vs.
Revision 3.9 by greg, Fri Feb 8 16:39:01 2013 UTC

# Line 1 | Line 1
1 – /* Copyright (c) 1998 Silicon Graphics, Inc. */
2 –
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ SGI";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6 –
4   /*
5   * Table-based cosine approximation.
6   *
# Line 11 | Line 8 | static char SCCSid[] = "$SunId$ SGI";
8   * to avoid conversion and guarantee that tsin(x)^2 + tcos(x)^2 == 1.
9   *
10   * No interpolation in this version.
11 + *
12 + * External symbols declared in rtmath.h
13   */
14  
15 + #include "copyright.h"
16 +
17   #include <math.h>
18  
19 + #include "rtmath.h"
20 +
21   #ifndef NCOSENTRY
22 < #define NCOSENTRY       256
22 > #define NCOSENTRY       1024
23   #endif
24  
22 – #ifdef M_PI
23 – #define PI              ((double)M_PI)
24 – #else
25 – #define PI              3.14159265358979323846
26 – #endif
25  
28 –
26   double
27 < tcos(x)                         /* approximate cosine */
31 < register double x;
27 > tcos(double x)                          /* approximate cosine */
28   {
29          static double   costab[NCOSENTRY+1];
30          register int    i;
# Line 40 | Line 36 | register double        x;
36          if (x < 0.)
37                  x = -x;
38          i = (NCOSENTRY*2./PI) * x  +  0.5;
39 <        if (i >= 4*NCOSENTRY)
40 <                i %= 4*NCOSENTRY;
39 >        while (i >= 4*NCOSENTRY)
40 >                i -= 4*NCOSENTRY;
41          switch (i / NCOSENTRY) {
42          case 0:
43                  return(costab[i]);
# Line 53 | Line 49 | register double        x;
49                  return(costab[(4*NCOSENTRY)-i]);
50          }
51          return(0.);             /* should never be reached */
52 + }
53 +
54 +
55 + /* Fast arctangent approximation due to Rajan et al. 2006 */
56 + double
57 + atan2a(double y, double x)
58 + {
59 +        double  ratio, aratio, val;
60 +
61 +        if (x == 0)
62 +                return (y > 0) ? PI/2. : 3./2.*PI;
63 +
64 +        aratio = (ratio = y/x) >= 0 ? ratio : -ratio;
65 +
66 +        if (aratio > 1.01)
67 +                return PI/2. - atan2a(x, y);
68 +
69 +        val = PI/4.*ratio - ratio*(aratio - 1.)*(0.2447 + 0.0663*aratio);
70 +
71 +        return val + PI*(x < 0);
72   }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines