home *** CD-ROM | disk | FTP | other *** search
/ Club Amiga de Montreal - CAM / CAM_CD_1.iso / files / 211.lha / Spiff / float.c < prev    next >
Encoding:
C/C++ Source or Header  |  1996-02-14  |  12.8 KB  |  812 lines
  1. /*                        Copyright (c) 1988 Bellcore
  2. **                            All Rights Reserved
  3. **       Permission is granted to copy or use this program, EXCEPT that it
  4. **       may not be sold for profit, the copyright notice must be reproduced
  5. **       on copies, and credit should be given to Bellcore where it is due.
  6. **       BELLCORE MAKES NO WARRANTY AND ACCEPTS NO LIABILITY FOR THIS PROGRAM.
  7. */
  8.  
  9.  
  10. #ifndef lint
  11. static char rcsid[]= "$Header: float.c,v 1.1 88/09/15 11:33:53 daniel Rel $";
  12. #endif
  13.  
  14. #include <ctype.h>
  15. #include "misc.h"
  16. #include "floatrep.h"
  17. #include "float.h"
  18. #include "strings.h"
  19.  
  20. #define _F_GETEND(x)    (x + (strlen(x)-1)) 
  21.  
  22. /*
  23. int floatcnt = 0;
  24. */
  25. /*
  26. **    routines to convert strings to our internal floating point form
  27. **        isfloat just looks at the string
  28. **        to see if a conversion is reasonable
  29. **            it does look-ahead on when it sees an 'e' and such.
  30. **        atocf actually does the conversion.
  31. **    these two routines could probably be combined
  32. */
  33.  
  34. /*
  35. **    test to see if the string can reasonably
  36. **        be interpreted as floating point number
  37. **    returns 0 if string can't be interpreted as a float
  38. **    otherwise returns the number of digits that will be used in F_atof
  39. */
  40. F_isfloat(str,need_decimal,allow_sign)
  41. char *str;
  42. int need_decimal;    /* if non-zero, require that a decimal point be present
  43.                 otherwise, accept strings like "123" */
  44. int allow_sign;        /* if non-zero, allow + or - to set the sign */
  45. {
  46.     int man_length = 0;    /* length of the fractional part (mantissa) */
  47.     int exp_length = 0;    /* length of the exponential part */
  48.     int got_a_digit = 0;    /* flag to set if we ever see a digit */
  49.  
  50.     /*
  51.     **    look for an optional leading sign marker
  52.     */
  53.     if (allow_sign && ('+' == *str  || '-' == *str))
  54.     {
  55.         str++; man_length++;
  56.     }
  57.     /*
  58.     **    count up the digits on the left hand side
  59.     **         of the decimal point
  60.     */
  61.     while(isdigit(*str))
  62.     {
  63.         got_a_digit = 1;
  64.         str++; man_length++;
  65.     }
  66.  
  67.     /*
  68.     **    check for a decimal point
  69.     */
  70.     if ('.' == *str)
  71.     {
  72.         str++; man_length++;
  73.     }
  74.     else
  75.     {
  76.         if (need_decimal)
  77.         {
  78.             return(0);
  79.         }
  80.     }
  81.  
  82.     /*
  83.     **    collect the digits on the right hand
  84.     **        side of the decimal point
  85.     */
  86.     while(isdigit(*str))
  87.     {
  88.         got_a_digit = 1;
  89.         str++; man_length++;
  90.     }
  91.  
  92.     if (!got_a_digit)
  93.         return(0);
  94.  
  95.     /*
  96.     **    now look ahead for an exponent
  97.     */
  98.     if ('e' == *str ||
  99.         'E' == *str ||
  100.         'd' == *str ||
  101.         'D' == *str)
  102.     {
  103.         str++; exp_length++;
  104.         if ('+' == *str  || '-' == *str)
  105.         {
  106.             str++; exp_length++;
  107.         }
  108.  
  109.         if (!isdigit(*str))
  110.         {
  111.             /*
  112.             **    look ahead went too far,
  113.             **    so return just the length of the mantissa
  114.             */
  115.             return(man_length);
  116.         }
  117.  
  118.         while (isdigit(*str))
  119.         {
  120.             str++; exp_length++;
  121.         }
  122.     }
  123.     return(man_length+exp_length);    /* return the total length */
  124. }
  125.  
  126. /*
  127. **    routine to convert a string to our internal
  128. **    floating point representation
  129. **
  130. **        similar to atof()
  131. */
  132. F_float
  133. F_atof(str,allflag)
  134. char *str;
  135. int allflag;    /* require that exactly all the characters are used */
  136. {
  137.     char *beg = str; /* place holder for beginning of the string */ 
  138.     char man[R_MANMAX];    /* temporary location to build the mantissa */
  139.     int length = 0;    /* length of the mantissa so far */
  140.     int got_a_digit = 0;    /* flag to set if we get a non-zero digit */ 
  141.     int i;
  142.     int resexp;
  143.  
  144.     F_float res;    /* where we build the result */
  145.  
  146. /*
  147. floatcnt++;
  148. */
  149.     res = R_makefloat();
  150.  
  151.     R_setsign(res,R_POSITIVE);
  152.  
  153.     resexp = 0;
  154.     man[0] = '\0';
  155.  
  156.     /*
  157.     **    check for leading sign
  158.     */
  159.     if ('+' == *str)
  160.     {
  161.         /*
  162.         **    sign should already be positive, see above in this
  163.         **        routine, so just skip the plus sign
  164.         */
  165.         str++;
  166.     }
  167.     else
  168.     {
  169.         if ('-' == *str)
  170.         {
  171.             R_setsign(res,R_NEGATIVE);
  172.             str++;
  173.         }
  174.     }
  175.  
  176.     /*
  177.     **    skip any leading zeros
  178.     */
  179.     while('0' == *str)
  180.     {
  181.         str++;
  182.     }
  183.  
  184.     /*
  185.     **    now snarf up the digits on the left hand side
  186.     **        of the decimal point
  187.     */
  188.     while(isdigit(*str))
  189.     {
  190.         got_a_digit = 1;
  191.         man[length++] = *str++;
  192.         man[length] = '\0';
  193.         resexp++;
  194.     }
  195.  
  196.     /*
  197.     **    skip the decimal point if there is one
  198.     */
  199.     if ('.' == *str)
  200.         str++;
  201.  
  202.     /*
  203.     **    trim off any leading zeros (on the right hand side)
  204.     **    if there were no digits in front of the decimal point.
  205.     */
  206.  
  207.     if (!got_a_digit)
  208.     {
  209.         while('0' == *str)
  210.         {
  211.             str++;
  212.             resexp--;
  213.         }
  214.     }
  215.  
  216.     /*
  217.     **    now snarf up the digits on the right hand side
  218.     */
  219.     while(isdigit(*str))
  220.     {
  221.         man[length++] = *str++;
  222.         man[length] = '\0';
  223.     }
  224.  
  225.     if ('e' == *str ||
  226.             'E' == *str ||
  227.             'd' == *str ||
  228.             'D' == *str )
  229.     {
  230.         str++;
  231.         resexp += atoi(str);
  232.     }
  233.  
  234.     if (allflag)
  235.     {
  236.         if ('+' == *str ||
  237.             '-' == *str)
  238.         {
  239.             str++;
  240.         }
  241.         while (isdigit(*str))
  242.         {
  243.             str++;
  244.         }
  245.         if ('\0' != *str)
  246.         {
  247.             (void) sprintf(Z_err_buf,
  248.                     "didn't use up all of %s in atocf",
  249.                     beg);
  250.             Z_fatal(Z_err_buf);
  251.         }
  252.     }
  253.  
  254.     /*
  255.     **    check for special case of all zeros in the mantissa
  256.     */
  257.     for (i=0;i<length;i++)
  258.     {
  259.         if (man[i] != '0')
  260.         {
  261.             /*
  262.             **    the mantissa is non-zero, so return it unchanged
  263.             */
  264.             S_trimzeros(man);
  265.             /*
  266.             **    save a copy of the mantissa
  267.             */
  268.             R_setfrac(res,man);
  269.             R_setexp(res,resexp);
  270.             return(res);
  271.         }
  272.     }
  273.  
  274.     /*
  275.     **    the answer is 0, so . . .
  276.     */
  277.     R_setzero(res);
  278.     return(res);
  279. }
  280.  
  281.  
  282. /*
  283. **    add s2 to s1
  284. */
  285. static
  286. void
  287. _F_stradd(s1,s2)
  288. char *s1,*s2;
  289. {
  290.     char *end1 = s1 + (strlen(s1)-1);
  291.     char *end2 = s2 + (strlen(s2)-1);
  292.  
  293.     static char result[R_MANMAX];
  294.     char *resptr = result+(R_MANMAX-1); /*point to the end of the array */
  295.     int carry = 0;
  296.     int tmp,val1,val2;
  297.  
  298.     *resptr-- = '\0';
  299.  
  300.     while ((end1 >= s1) ||  ( end2 >= s2))
  301.     {
  302.         if (end1 >= s1)
  303.         {
  304.             val1 = *end1 - '0';
  305.             --end1;
  306.         }
  307.         else
  308.         {
  309.             val1 = 0;
  310.         }
  311.  
  312.         if (end2 >= s2)
  313.         {
  314.             val2 = *end2 - '0';
  315.             --end2;
  316.         }
  317.         else
  318.         {
  319.             val2 = 0;
  320.         }
  321.  
  322.         tmp = val1 + val2 + carry;
  323.         if (tmp > 9)
  324.         {
  325.             carry = 1;
  326.             tmp -= 10;
  327.         }
  328.         else
  329.         {
  330.             carry = 0;
  331.         }
  332.  
  333.         *resptr-- = tmp+'0';
  334.     }
  335.     if (carry)
  336.     {
  337.         *resptr =  '1';
  338.     }
  339.     else
  340.     {
  341.         resptr++;
  342.     }
  343.     (void) strcpy(s1,resptr);
  344.     return;
  345. }
  346.  
  347. /*
  348. **    add zero(s) onto the end of a string
  349. */
  350. static void
  351. addzeros(ptr,count)
  352. char *ptr;
  353. int count;
  354. {
  355.     for(;count> 0;count--)
  356.     {
  357.         (void) strcat(ptr,"0");
  358.     }
  359.     return;
  360. }
  361.  
  362. /*
  363. **    subtract two mantissa strings
  364. */
  365. F_float
  366. F_floatsub(p1,p2)
  367. F_float  p1,p2;
  368. {
  369.     static F_float result;
  370.     static needinit = 1;
  371.     static char man1[R_MANMAX],man2[R_MANMAX],diff[R_MANMAX];
  372.     int exp1,exp2;
  373.     char *diffptr,*big,*small;
  374.     int man_cmp_val,i,borrow;
  375.  
  376.     if (needinit)
  377.     {
  378.         result = R_makefloat();
  379.         needinit = 0;
  380.     }
  381.  
  382.     man1[0] = '\0';
  383.     man2[0] = '\0';
  384.  
  385.     exp1 = R_getexp(p1);
  386.     exp2 = R_getexp(p2);
  387.  
  388.     /*
  389.     **    line up the mantissas
  390.     */
  391.     while (exp1 < exp2)
  392.     {
  393.         (void) strcat(man1,"0");
  394.         exp1++;
  395.     }
  396.  
  397.     while(exp1 > exp2)
  398.     {
  399.         (void) strcat(man2,"0");
  400.         exp2++;
  401.     }
  402.  
  403.     if (exp1 != exp2)    /* boiler plate assertion */
  404.     {
  405.         Z_fatal("mantissas didn't get lined up properly in floatsub");
  406.     }
  407.  
  408.     (void) strcat(man1,R_getfrac(p1));
  409.     (void) strcat(man2,R_getfrac(p2));
  410.     
  411.     /*
  412.     **    now that the mantissa are aligned,
  413.     **    if the strings are the same, return 0
  414.     */
  415.     if((man_cmp_val = strcmp(man1,man2)) == 0)
  416.     {
  417.         R_setzero(result);
  418.         return(result);
  419.     }
  420.  
  421.     /*
  422.     **    pad the shorter string with 0's
  423.     **        when this loop finishes, both mantissas should
  424.     **        have the same length
  425.     */
  426.     if (strlen(man1)> strlen(man2))
  427.     {
  428.         addzeros(man2,strlen(man1)-strlen(man2));
  429.     }
  430.     else
  431.     {
  432.         if (strlen(man1)<strlen(man2))
  433.         {
  434.             addzeros(man1,strlen(man2)-strlen(man1));
  435.         }
  436.     }
  437.  
  438.     if (strlen(man1) != strlen(man2))    /* pure boilerplate */
  439.     {
  440.         Z_fatal("lengths not equal in F_floatsub");
  441.     }
  442.  
  443.     if (man_cmp_val < 0)
  444.     {
  445.         big = man2;
  446.         small = man1;
  447.     }
  448.     else
  449.     {
  450.         big = man1;
  451.         small = man2;
  452.     }
  453.  
  454.     /*
  455.     **    find the difference between the mantissas
  456.     */
  457.     for(i=(strlen(big)-1),borrow=0,diff[strlen(big)] = '\0';i>=0;i--)
  458.     {
  459.         char from;
  460.         if (borrow)
  461.         {
  462.             if (big[i] == '0')
  463.             {
  464.                 from = '9';
  465.             }
  466.             else
  467.             {
  468.                 from = big[i]-1;
  469.                 borrow = 0;
  470.             }
  471.         }
  472.         else
  473.         {
  474.             if(big[i]<small[i])
  475.             {
  476.                 from = '9'+1;
  477.                 borrow = 1;
  478.             }
  479.             else
  480.             {
  481.                 from = big[i];
  482.             }
  483.         }
  484.         diff[i] = (from-small[i]) + '0';
  485.     }
  486.  
  487.     /*
  488.     ** trim the leading zeros on the difference
  489.     */
  490.     diffptr = diff;
  491.     while('0' == *diffptr)
  492.     {
  493.         diffptr++;
  494.         exp1--;
  495.     }
  496.  
  497.     R_setexp(result,exp1); /* exponents are equal at the point */
  498.     R_setfrac(result,diffptr);
  499.     R_setsign(result,R_POSITIVE);
  500.     return(result);
  501. }
  502.  
  503. F_floatcmp(f1,f2)
  504. F_float f1,f2;
  505. {
  506.     static char man1[R_MANMAX],man2[R_MANMAX];
  507.  
  508.     /*
  509.     **        special case for zero
  510.     */
  511.     if (R_zerofloat(f1))
  512.     {
  513.         if (R_zerofloat(f2))
  514.         {
  515.             return(0);
  516.         }
  517.         else
  518.         {
  519.             return(-1);
  520.         }
  521.     }
  522.     else
  523.     {
  524.         if (R_zerofloat(f2))
  525.         {
  526.             return(1);
  527.         }
  528.     }
  529.  
  530.     /*
  531.     **    to reach this point, both numbers must be non zeros
  532.     */
  533.     if (R_getexp(f1) < R_getexp(f2))
  534.     {
  535.         return(-1);
  536.     }
  537.  
  538.     if (R_getexp(f1) > R_getexp(f2))
  539.     {
  540.         return(1);
  541.     }
  542.  
  543.     (void) strcpy(man1,R_getfrac(f1));
  544.     S_trimzeros(man1);
  545.  
  546.     (void) strcpy(man2,R_getfrac(f2));
  547.     S_trimzeros(man2);
  548.     return(strcmp(man1,man2));
  549. }
  550.  
  551. F_float
  552. F_floatmul(f1,f2)
  553. F_float f1,f2;
  554. {
  555.     static char prod[R_MANMAX];
  556.     char *end;
  557.     int count1 = 0;
  558.     int count2 = 0;
  559.     int tmp,len;
  560.     char *end1;
  561.     char *end2;
  562.     static char man1[R_MANMAX],man2[R_MANMAX];
  563.     char *bigman,*smallman;
  564.     static F_float result;
  565.     static int needinit = 1;
  566.  
  567.     if (needinit)
  568.     {
  569.         result = R_makefloat();
  570.         needinit = 0;
  571.     }
  572.     /*
  573.     **    special case for a zero result
  574.     */
  575.     if (R_zerofloat(f1) || R_zerofloat(f2))
  576.     {
  577.         R_setzero(result);
  578.         return(result);
  579.     }
  580.  
  581.     (void) strcpy(man1,R_getfrac(f1));
  582.     (void) strcpy(man2,R_getfrac(f2));
  583.  
  584.     end1 = _F_GETEND(man1);
  585.     end2 = _F_GETEND(man2);
  586.  
  587.     /*
  588.     **    decide which number will cause multiplication loop to go
  589.     **    around the least
  590.     */
  591.     while(end1 >= man1)
  592.     {
  593.         count1 += *end1 - '0';
  594.         end1--;
  595.     }
  596.  
  597.     while(end2 >= man2)
  598.     {
  599.         count2 += *end2 - '0';
  600.         end2--;
  601.     }
  602.  
  603.  
  604.     if (count1 > count2)
  605.     {
  606.         bigman = man1;
  607.         smallman = man2;
  608.     }
  609.     else
  610.     {
  611.         bigman = man2;
  612.         smallman = man1;
  613.     }
  614.     S_trimzeros(bigman);
  615.     S_trimzeros(smallman);
  616.     len = strlen(bigman) +  strlen(smallman);
  617.  
  618.     end = _F_GETEND(smallman);
  619.     (void) strcpy(prod,"0");
  620.  
  621.     /*
  622.     **    multiplication by repeated addition
  623.     */
  624.     while(end >= smallman)
  625.     {
  626.         for(tmp = 0;tmp<*end-'0';tmp++)
  627.         {
  628.             _F_stradd(prod,bigman);
  629.         }
  630.         addzeros(bigman,1);
  631.         end--;
  632.     }
  633.  
  634.     R_setfrac(result,prod);
  635.     R_setexp(result,(((R_getexp(f1) + R_getexp(f2)) - len)+ strlen(prod)));
  636.  
  637.     if (R_getsign(f1) == R_getsign(f2))
  638.     {
  639.         R_setsign(result,R_POSITIVE);
  640.     }
  641.     else
  642.     {
  643.         R_setsign(result,R_NEGATIVE);
  644.     }
  645.     return(result);
  646. }
  647.  
  648. _F_xor(x,y)
  649. {
  650.     return(((x) && !(y)) || (!(x) && (y)));
  651. }
  652. #define    _F_SAMESIGN(x,y)    _F_xor((x<0),(y<0))
  653. #define _F_ABSADD(x,y)        (Z_ABS(x) + Z_ABS(y))
  654.  
  655. _F_ABSDIFF(x,y)
  656. {
  657.     if (Z_ABS(x) < Z_ABS(y))
  658.     {
  659.         return(Z_ABS(y) - Z_ABS(x));
  660.     }
  661.     else
  662.     {
  663.         return(Z_ABS(x) - Z_ABS(y));
  664.     }
  665. }
  666. /*
  667. **    add two floats without regard to sign
  668. */
  669. F_float
  670. F_floatmagadd(p1,p2)
  671. F_float p1,p2;
  672. {
  673.     static F_float result;
  674.     static int needinit = 1;
  675.  
  676.     static char  man1[R_MANMAX],man2[R_MANMAX];
  677.  
  678.     int digits;    /* count of the number of digits needed to represent the
  679.                 result */
  680.     int resexp;    /* exponent of the result */
  681.     int len;    /* length of the elements before adding */
  682.     char *diffptr;
  683.  
  684.     if (needinit)
  685.     {
  686.         result = R_makefloat();
  687.         needinit = 0;
  688.     }
  689.     (void) strcpy(man1,"");
  690.     (void) strcpy(man2,"");
  691.  
  692.     /*
  693.     **    find the difference in the exponents number of digits
  694.     */
  695.     if( _F_SAMESIGN(R_getexp(p1),R_getexp(p2)))
  696.     {
  697.         digits =  _F_ABSDIFF(R_getexp(p1),R_getexp(p2));
  698.     }
  699.     else
  700.     {
  701.         digits = _F_ABSADD(R_getexp(p1),R_getexp(p2));
  702.     }
  703.  
  704.     /*
  705.     **    make sure that there is room to store the result
  706.     */
  707.     if (digits>0)
  708.     { 
  709.         if (R_getexp(p1) < R_getexp(p2))
  710.         {
  711.             /*
  712.             **    leave room for terminator
  713.             */
  714.             if (digits+strlen(R_getfrac(p1)) > (R_MANMAX-1))
  715.             {
  716.                 (void) sprintf(Z_err_buf,
  717.                    "numbers differ by too much in magnitude");
  718.                 Z_fatal(Z_err_buf);
  719.             }
  720.         }
  721.         else
  722.         {
  723.             /*
  724.             **    leave room for terminator
  725.             */
  726.             if (digits+strlen(R_getfrac(p2)) > (R_MANMAX-1))
  727.             {
  728.                 (void) sprintf(Z_err_buf,
  729.                    "numbers differ by too much in magnitude");
  730.                 Z_fatal(Z_err_buf);
  731.             }
  732.         }
  733.     }
  734.     else
  735.     {
  736.         /*
  737.         **    leave room for terminator and possible carry
  738.         */
  739.         if (Z_MAX(strlen(R_getfrac(p1)),
  740.             strlen(R_getfrac(p2))) > (R_MANMAX-2))
  741.         {                        
  742.             (void) sprintf(Z_err_buf,
  743.                "numbers differ by too much in magnitude");
  744.             Z_fatal(Z_err_buf);
  745.         }
  746.     }
  747.  
  748.     /*
  749.     **    pad zeroes on the front of the smaller number
  750.     */
  751.     if (R_getexp(p1) < R_getexp(p2))
  752.     {
  753.  
  754.         addzeros(man1,digits);
  755.         resexp = R_getexp(p2);
  756.     }
  757.     else
  758.     {
  759.         addzeros(man2,digits);
  760.         resexp = R_getexp(p1);
  761.     }
  762.     (void) strcat(man1,R_getfrac(p1));
  763.     (void) strcat(man2,R_getfrac(p2));
  764.  
  765.     len = Z_MAX(strlen(man1),strlen(man2));
  766.  
  767.     /*
  768.     **    add the two values
  769.     */
  770.     _F_stradd(man1,man2);
  771.  
  772.     /*
  773.     **    adjust the exponent to account for a
  774.     **        possible carry
  775.     */
  776.     resexp += strlen(man1) - len;
  777.  
  778.  
  779.     /*
  780.     ** trim the leading zeros on the sum
  781.     */
  782.     diffptr = man1;
  783.     while('0' == *diffptr)
  784.     {
  785.         diffptr++;
  786.         resexp--;
  787.     }
  788.  
  789.     R_setfrac(result,diffptr);
  790.     R_setexp(result,resexp);
  791.     R_setsign(result,R_POSITIVE);
  792.  
  793.     return(result);
  794. }
  795.  
  796. /*
  797. **    useful debugging routine. we don't call it in the release,
  798. **    so it is commented out, but we'll leave it for future use
  799. */
  800.  
  801. /*
  802. F_printfloat(fl)
  803. F_float fl;
  804. {
  805.     (void) printf("fraction = :%s: exp = %d sign = %c\n",
  806.             R_getfrac(fl),
  807.             R_getexp(fl),
  808.             ((R_getsign(fl) == R_POSITIVE) ? '+': '-'));
  809.  
  810. }
  811. */
  812.