home *** CD-ROM | disk | FTP | other *** search
/ Garbo / Garbo.cdr / pc / graphics / gl.pak / TRIG.C < prev    next >
Encoding:
C/C++ Source or Header  |  1989-03-31  |  6.8 KB  |  244 lines
  1. #ifndef lint
  2. static char sccsid[] = "@(#) trig.c 5.1 89/02/20";
  3. #endif
  4.  
  5. /*
  6.  *    Copyright (c) David T. Lewis 1988
  7.  *    All rights reserved.
  8.  *
  9.  *    Permission is granted to use this for any personal noncommercial use.
  10.  *    You may not distribute source or executable code for profit, nor
  11.  *    may you distribute it with a commercial product without the written
  12.  *    consent of the author.  Please send modifications to the author for
  13.  *    inclusion in updates to the program.  Thanks.
  14.  */
  15.  
  16. /* Integer approximations to simple trig functions.  Results are scaled    */
  17. /* by a factor of TRIG_SCALE.                        */
  18.  
  19. /* The macro TRIG_SCALE is defined in graphics.h so that other routines    */
  20. /* can use it.  F_TRIG_SCALE is the same value, but in floating point    */
  21. /* format.  The defines should look like this:                */
  22. /*                                    */
  23. /*    #define TRIG_SCALE 32768                    */
  24. /*    #define F_TRIG_SCALE 32768.0                    */
  25. /*                                    */
  26. /* The macro PTS_PER_QUADRANT defines the size of the sin_table and    */
  27. /* cos_table arrays.  It is defined in graphics.h and should look    */
  28. /* like this:                                */
  29. /*                                    */
  30. /*    #define PTS_PER_QUADRANT 64                    */
  31. /*                                    */
  32.  
  33. #include "config.h"
  34. #if MIX_C
  35. #else
  36. #include <math.h>
  37. #endif /* MIX_C */
  38. #include "bitmaps.h"
  39. #include "graphics.h"
  40.  
  41. #define L_PI 102944L        /* Value of pi times TRIG_SCALE.    */
  42. #define L_PI_2 51472L        /* Value of pi/2 times TRIG_SCALE.    */
  43. #define L_PI_4 25736L        /* Value of pi/4 times TRIG_SCALE.    */
  44.  
  45. #define CIRCLE_SIZE 205887L    /* Value of pi * 2 * TRIG_SCALE.    */
  46.  
  47. #define TABLESLICE 804L        /* (L_PI_2 / PTS_PER_QUADRANT)        */
  48.  
  49. #define XINDEX 0        /* For specifying array indices.    */
  50. #define YINDEX 1
  51.  
  52. /* The following two arrays are tables of values for cos and sin    */
  53. /* for one quadrant of a circle, divided into PTS_PER_QUADRANT slices.    */
  54.  
  55. long cos_table[PTS_PER_QUADRANT] = {
  56.     32758L, 32729L, 32679L, 32610L, 32522L, 32413L, 32286L, 32138L, 
  57.     31972L, 31786L, 31581L, 31357L, 31114L, 30853L, 30572L, 30274L, 
  58.     29957L, 29622L, 29269L, 28899L, 28511L, 28106L, 27684L, 27246L, 
  59.     26791L, 26320L, 25833L, 25330L, 24812L, 24279L, 23732L, 23170L, 
  60.     22595L, 22006L, 21403L, 20788L, 20160L, 19520L, 18868L, 18205L, 
  61.     17531L, 16846L, 16151L, 15447L, 14733L, 14010L, 13279L, 12540L, 
  62.     11793L, 11039L, 10279L, 9512L, 8740L, 7962L, 7180L, 6393L, 
  63.     5602L, 4808L, 4011L, 3212L, 2411L, 1608L, 804L, 0L
  64. };
  65.  
  66. long sin_table[PTS_PER_QUADRANT] = {
  67.     804L, 1608L, 2411L, 3212L, 4011L, 4808L, 5602L, 6393L, 
  68.     7180L, 7962L, 8740L, 9512L, 10279L, 11039L, 11793L, 12540L, 
  69.     13279L, 14010L, 14733L, 15447L, 16151L, 16846L, 17531L, 18205L, 
  70.     18868L, 19520L, 20160L, 20788L, 21403L, 22006L, 22595L, 23170L, 
  71.     23732L, 24279L, 24812L, 25330L, 25833L, 26320L, 26791L, 27246L, 
  72.     27684L, 28106L, 28511L, 28899L, 29269L, 29622L, 29957L, 30274L, 
  73.     30572L, 30853L, 31114L, 31357L, 31581L, 31786L, 31972L, 32138L, 
  74.     32286L, 32413L, 32522L, 32610L, 32679L, 32729L, 32758L, 32768L
  75. };
  76.  
  77. /* Estimate sin by interpolation.  Return long integer value.        */
  78.  
  79. long longsin(theta)
  80. float theta;
  81. {
  82.     long longtheta;
  83.     int quadrant;
  84.     long value1, value2;
  85.     long value;
  86.     long intpart;
  87.     long fractionpart;
  88.  
  89.     /* We have a floating point value for theta.  Put theta     */
  90.     /* into a long integer for the rest of the processing.        */
  91.     /* Keep everything scaled by a factor of TRIG_SCALE.        */
  92.  
  93.     longtheta = (long)(theta * F_TRIG_SCALE);
  94.  
  95.     /* Make angle positive.                        */
  96.  
  97.     while (longtheta < 0) longtheta += CIRCLE_SIZE;
  98.  
  99.     /* Figure quadrant while making it less than PI/2.        */
  100.  
  101.     quadrant = 1;
  102.  
  103.     while (longtheta > L_PI_2)  {
  104.         longtheta -= L_PI_2;
  105.         quadrant++;
  106.         if (quadrant > 4) quadrant = 1;
  107.     }
  108.  
  109.     /* If Quadrant 2 or 4 we need to use the reverse of the        */
  110.     /* interpolation table.                        */
  111.  
  112.     if (quadrant == 2 || quadrant == 4) longtheta = L_PI_2 - longtheta;
  113.  
  114.     /* Express angle in terms of number of table increments.    */
  115.  
  116.     intpart = longtheta / TABLESLICE;
  117.  
  118.     fractionpart = longtheta % TABLESLICE;
  119.  
  120.     /* This looks like a kludge, but it isn't.  When TABLESLICE is    */
  121.     /* set to TABLESLICE, we have the smallest error in the result.    */
  122.     /* This just rounds back down when we overshoot the table.    */
  123.  
  124.     if (intpart >= PTS_PER_QUADRANT)  {
  125.         intpart = PTS_PER_QUADRANT - 1;
  126.         fractionpart = TABLESLICE;
  127.     }
  128.  
  129.     /* Interpolate.                            */
  130.  
  131.     if (intpart <= 0)  {
  132.         value1 = 0;
  133.         value2 = sin_table[0];
  134.     }
  135.     else  {
  136.         value2 = sin_table[(int)(intpart--)];
  137.         value1 = sin_table[(int)intpart];
  138.     }
  139.  
  140.     value = value1 + (fractionpart * (value2 - value1)) / TABLESLICE;
  141.  
  142.     /* If Quadrant 3 or 4 the result should be negative.        */
  143.  
  144.     if (quadrant == 3 || quadrant == 4) return(-value);
  145.     else return(value);
  146. }
  147.  
  148.  
  149. /* Estimate cos by interpolation.  Return long integer value.        */
  150.  
  151. long longcos(theta)
  152. float theta;
  153. {
  154.     long longtheta;
  155.     int quadrant;
  156.     long value1, value2;
  157.     long value;
  158.     long intpart;
  159.     long fractionpart;
  160.  
  161.     /* We have a floating point value for theta.  Put theta     */
  162.     /* into a long integer for the rest of the processing.        */
  163.     /* Keep everything scaled by a factor of TRIG_SCALE.        */
  164.  
  165.     longtheta = (long)(theta * F_TRIG_SCALE);
  166.  
  167.     /* Make angle positive.                        */
  168.  
  169.     while (longtheta < 0) longtheta += CIRCLE_SIZE;
  170.  
  171.     /* Figure quadrant while making it less than PI/2.        */
  172.  
  173.     quadrant = 1;
  174.  
  175.     while (longtheta > L_PI_2)  {
  176.         longtheta -= L_PI_2;
  177.         quadrant++;
  178.         if (quadrant > 4) quadrant = 1;
  179.     }
  180.  
  181.     /* If Quadrant 2 or 4 we need to use the reverse of the        */
  182.     /* interpolation table.                        */
  183.  
  184.     if (quadrant == 2 || quadrant == 4) longtheta = L_PI_2 - longtheta;
  185.  
  186.     /* Express angle in terms of number of table increments.    */
  187.  
  188.     intpart = longtheta / TABLESLICE;
  189.  
  190.     fractionpart = longtheta % TABLESLICE;
  191.  
  192.     /* This looks like a kludge, but it isn't.  When TABLESLICE is    */
  193.     /* set to TABLESLICE, we have the smallest error in the result.    */
  194.     /* This just rounds back down when we overshoot the table.    */
  195.  
  196.     if (intpart >= PTS_PER_QUADRANT)  {
  197.         intpart = PTS_PER_QUADRANT - 1;
  198.         fractionpart = TABLESLICE;
  199.     }
  200.  
  201.     /* Interpolate.                            */
  202.  
  203.     if (intpart <= 0)  {
  204.         value1 = 32768L;
  205.         value2 = cos_table[0];
  206.     }
  207.     else  {
  208.         value2 = cos_table[(int)(intpart--)];
  209.         value1 = cos_table[(int)intpart];
  210.     }
  211.  
  212.     value = value1 + (fractionpart * (value2 - value1)) / TABLESLICE;
  213.  
  214.     /* If Quadrant 2 or 3 the result should be negative.        */
  215.  
  216.     if (quadrant == 2 || quadrant == 3) return(-value);
  217.     else return(value);
  218. }
  219.  
  220. /* Test program:
  221. main()  {
  222.     float theta;
  223.  
  224. float sinval;
  225. float lsinval;
  226. float cosval;
  227. float lcosval;
  228.  
  229.     for (theta = -8.3; theta < 100.0; theta += 0.1)  {
  230.  
  231.         sinval = sin((double)theta);
  232.         lsinval = longsin(theta) / F_TRIG_SCALE;
  233.         cosval = cos((double)theta);
  234.         lcosval = longcos(theta) / F_TRIG_SCALE;
  235.  
  236.         printf("SIN: angle %g    %lg    %g    err %g\n",theta,sinval,
  237.         lsinval,sinval-lsinval);
  238.         printf("COS: angle %g    %lg    %g    err %g\n",theta,cosval,
  239.         lcosval,cosval-lcosval);
  240.     }
  241. }
  242. */
  243.  
  244.