home *** CD-ROM | disk | FTP | other *** search
/ Monster Media 1994 #1 / monster.zip / monster / MAGAZINE / DDJ9309.ZIP / FPALIB.ZIP / FALIB.C next >
C/C++ Source or Header  |  1993-05-12  |  32KB  |  993 lines
  1. /* Extended IEEE Compatible Floating Point Arithmetic Library
  2. **
  3. ** Version 1.0
  4. ** Copyright (C) 1990, by Fred Motteler
  5. ** All Rights Reserved
  6. **
  7. ** This is a collection of routines to convert arbitrary precision float and
  8. ** integer values to and from ASCII strings.
  9. **
  10. ** Each routine that returns a floating point result returns a condition code
  11. ** value:
  12. **
  13. **    0 0 0 0 N Z I n
  14. **    Where N = negative (FFNEG)
  15. **          Z = zero (FFZERO)
  16. **          I = infinite (FFINF)
  17. **          n = not a number (FFNAN)
  18. **
  19. ** Each routine that returns a integer result returns a condition code
  20. ** value:
  21. **
  22. **    0 0 0 0 Z V S C
  23. **    Where Z = zero (ZERO)
  24. **          V = overflow (OVERFLOW)
  25. **          S = sign (SIGN)
  26. **          C = carry (CARRY)
  27. **
  28. ** Note that these routines recognize zero and infinite input values.  Not a
  29. ** number and denormalized input values are not recognized.  The routines
  30. ** blindly assume that the arguments are valid floating point numbers.
  31. **
  32. ** These routines will return valid zero, infinite, and not a number values.
  33. */
  34. #include <stdio.h>
  35. #ifndef MWC
  36. #include <stdlib.h>
  37. #endif
  38. #include "imlib.h"
  39. #include "ffmlib.h"
  40. #include "fmlib.h"
  41.  
  42. #define SINGLEXP 8
  43. #define SINGLEFRAC 23
  44. #define SINGLETOT 4
  45.  
  46. static unsigned char log10of2AB[4] = {0x3e, 0x9a, 0x20, 0x9b};
  47.  
  48. #ifdef TEST
  49. #define DOUBLEXP 11
  50. #define DOUBLEFRAC 52
  51. #define DOUBLETOT 8
  52. #define EXTENDEXP 15
  53. #define EXTENDFRAC 63
  54. #define EXTENDTOT 10
  55. unsigned char intval1AB[4] = {0x0, 0x0, 0x12, 0x34};
  56. unsigned char intval2AB[4] = {0x12, 0x34, 0x56, 0x78};
  57. unsigned char fltval1AB[4] = {0x3e, 0x9a, 0x20, 0x9b};
  58. unsigned char dblval1AB[8] = {0x52, 0x34, 0x56, 0x78, 0x9a, 0xbc, 0xde, 0xf0};
  59.  
  60. void
  61. main()
  62. {
  63.     unsigned char *strBP, *fltvalBP;
  64.     int i;
  65.  
  66.     strBP = (unsigned char *) malloc(80);
  67.     fltvalBP = (unsigned char *) malloc(16);
  68.  
  69.     printf("LONG to ASCII\n");
  70.     intostr(intval1AB, 4, strBP, 16, 10);
  71.     printf("intval1AB to ASCII: %s  Should be 4660.\n",strBP);
  72.  
  73.     printf("LONG to ASCII\n");
  74.     intostr(intval2AB, 4, strBP, 16, 10);
  75.     printf("intval2AB to ASCII: %s  Should be 305419896.\n",strBP);
  76.  
  77.     printf("SINGLE to ASCII\n");
  78.     fltostr(fltval1AB, SINGLEFRAC, SINGLEXP, strBP, 16);
  79.     printf("fltval1AB to ASCII: %s\n",strBP);
  80.  
  81.     printf("DOUBLE to ASCII\n");
  82.     fltostr(dblval1AB, DOUBLEFRAC, DOUBLEXP, strBP, 32);
  83.     printf("dblval1AB to ASCII: %s\n",strBP);
  84.  
  85.     printf("ASCII to LONG\n");
  86.     strtoint("1024",5, intval1AB, 4, 10);
  87.     printf("Should be 0 0 4 0:");
  88.     for (i = 0; i < 4; i++)
  89.     printf("%2x",intval1AB[i]);
  90.     printf("\n");
  91.     strtoint("305419896",10, intval1AB, 4, 10);
  92.     printf("Should be 12345678:");
  93.     for (i = 0; i < 4; i++)
  94.     printf("%2x",intval1AB[i]);
  95.     printf("\n");
  96.  
  97.     printf("ASCII to FLOAT\n");
  98.     strtoflt("3.459e2", 9, fltvalBP, SINGLEFRAC, SINGLEXP);
  99.     for (i = 0; i < 4; i++)
  100.     printf("%2x",*(fltvalBP + i));
  101.     printf("\n");
  102.     fltostr(fltvalBP, SINGLEFRAC, SINGLEXP, strBP, 16);
  103.     printf("Conversion back to ASCII: %s\n",strBP);
  104.     
  105.     printf("ASCII to DOUBLE\n");
  106.     strtoflt( "3.459e2", 9, fltvalBP, DOUBLEFRAC, DOUBLEXP);
  107.     for (i = 0; i < 8; i++)
  108.     printf("%2x",*(fltvalBP + i));
  109.     printf("\n");
  110.     fltostr(fltvalBP, DOUBLEFRAC, DOUBLEXP, strBP, 32);
  111.     printf("Conversion back ASCII: %s\n",strBP);
  112.  
  113.     exit(0);
  114. }
  115. #endif
  116.  
  117. /* Function:    unsigned char intostr(unsigned char *intvalBP,
  118. **                      int intlenN,
  119. **                      char *strBP,
  120. **                      int slenN,
  121. **                      int radixN)
  122. **
  123. ** This function converts the integer value pointed to intvalBP to a decimal
  124. ** ASCII string representation of the value pointed to by strBP.  intlenN is
  125. ** the length of the integer to convert in bytes.  strBP is the length of the
  126. ** string buffer in bytes.  radixN is the radix of the ASCII string value.
  127. ** Only radixes 2 thru 36 are supported.
  128. */
  129. unsigned char
  130. #ifdef PROTOTYPES
  131. intostr(unsigned char *intvalBP, int intlenN, char *strBP, int slenN,
  132.     int radixN)
  133. #else
  134. intostr(intvalBP, intlenN, strBP, slenN, radixN)
  135. unsigned char *intvalBP;
  136. int intlenN;
  137. char *strBP;
  138. int slenN;
  139. int radixN;
  140. #endif
  141. {
  142.     unsigned char condcodeB, digitB;
  143.     unsigned char *radixBP;
  144.     unsigned char *tempBP;
  145.     unsigned char *digitstackBP;
  146.     int digitindexN;
  147.     int i;
  148.  
  149.     condcodeB = 0;
  150.  
  151.     /* Check if the integer is zero. */
  152.     if (ucheckm(intvalBP, intlenN) == 0)
  153.     {
  154.     *strBP++ = '0';
  155.     *strBP = '\0';
  156.     return(ZERO);
  157.     }
  158.  
  159.     /* Check if integer is negative. */
  160.     if ((*intvalBP & 0x80) != 0)
  161.     {
  162.     /* Change sign of the integer value */
  163.     inegm(intvalBP, intlenN);
  164.     /* Write out a minus (-) sign */
  165.     *strBP++ = '-';
  166.     slenN--;
  167.     condcodeB = SIGN;
  168.     }
  169.  
  170.     /* Convert the radix value to multiple byte integer format. */
  171.     radixBP = FCALLOC(intlenN, 1, "INTOSTR1");
  172.     if (radixN > 36) radixN = 36;
  173.     if (radixN < 2) radixN = 2;
  174.     *(radixBP + intlenN - 1) = (unsigned char) radixN;
  175.  
  176.     /* Allocate space for a digit stack.  This is required because repeated
  177.      * division gives the least significant digit first, however ASCII string
  178.      * representations have the most significant digit first. */
  179.     digitstackBP = FCALLOC(slenN, 1, "INTOSTR2");
  180.     digitindexN = 0;
  181.  
  182.     /* Convert the integer value to a string, one digit at a time.  This 
  183.      * starts with the least significant digits and works towards the most
  184.      * significant. */
  185.     tempBP = (unsigned char *) FMALLOC(intlenN, "INTOSTR3");
  186.     while (ucheckm(intvalBP, intlenN) != 0)
  187.     {
  188.     /* Copy the divisor the the correct place. */
  189.     for (i = 0; i < intlenN; i++)
  190.         *(tempBP + i) = *(radixBP + i);
  191.  
  192.     /* Divide the integer value by the radix. */
  193.     udivm(intvalBP, tempBP, intlenN);
  194.  
  195.     /* The quotient is in tempBP, the remainder is in intvalBP. */
  196.     digitB = *(intvalBP + intlenN - 1);
  197.  
  198.     /* Copy the quotient back to the divisor place. */
  199.     for (i = 0; i < intlenN; i++)
  200.         *(intvalBP + i) = *(tempBP + i);
  201.  
  202.     /* Convert the digit to ASCII */
  203.     if (digitB < 10)
  204.         digitB += 0x30;        /* 0 thru 9 */
  205.     else
  206.         digitB += 0x57;        /* a thru z */
  207.  
  208.     /* Push the digit onto the digit stack */
  209.     *(digitstackBP + digitindexN++) = digitB;
  210.     if( digitindexN >= slenN )
  211.     {
  212.         FFREE(radixBP);
  213.         FFREE(digitstackBP);
  214.         FFREE(tempBP);
  215.         return(OVERFLOW);
  216.     }
  217.     }   
  218.     /* Pop the digits off the stack and copy them to the result string. */ 
  219.     while( digitindexN > 0 )
  220.     {
  221.     *strBP++ = *(digitstackBP + (--digitindexN));
  222.     if( (slenN--) <= 0)
  223.     {
  224.         FFREE(radixBP);
  225.         FFREE(digitstackBP);
  226.         FFREE(tempBP);
  227.         return(OVERFLOW);
  228.     }
  229.     }
  230.     *strBP = '\0';
  231.  
  232.     /* All done, clean up and return. */
  233.     FFREE(radixBP);
  234.     FFREE(digitstackBP);
  235.     FFREE(tempBP);
  236.     return(condcodeB);
  237. }
  238.  
  239. /* Function:    unsigned char fltostr(unsigned char *fltvalBP, int fracbitN,
  240. **                      int expbitN, char *strBP, int slenN)
  241. **
  242. ** This function converts the floating point value pointed to fltvalBP to a
  243. ** decimal ASCII string representation of the value pointed to by strBP.
  244. ** fracbitN is the length of the floating point mantissa in bits.  expbitN
  245. ** is the length of the floating point exponent in bits.  strBP is the
  246. ** length of the string buffer in bytes.
  247. */
  248. unsigned char
  249. #ifdef PROTOTYPES
  250. fltostr(unsigned char *fltvalBP, int fracbitN, int expbitN,
  251.     char *strBP, int slenN)
  252. #else
  253. fltostr(fltvalBP, fracbitN, expbitN, strBP, slenN)
  254. unsigned char *fltvalBP;
  255. int fracbitN;
  256. int expbitN;
  257. char *strBP;
  258. int slenN;
  259. #endif
  260. {
  261.     int expbyteN, fracbyteN;
  262.     unsigned char *expbiasBP;
  263.     unsigned char *exponeBP;
  264.     unsigned char *fltfracBP, *fltexpBP;
  265.     unsigned char *tempfltBP, *tempexpBP;
  266.     unsigned char *ptenfltBP, *onetenthBP, *oneBP;
  267.     unsigned char ptenccB;
  268.     unsigned char fltsignB;
  269.     unsigned char condcodeB;
  270.     int i, totalenN;
  271.     unsigned char minusoneB, zeroB;
  272.     unsigned int mantlenN;
  273.     unsigned char *mantintBP;
  274.     unsigned char *mantlenBP;
  275.  
  276.     minusoneB = 0xff;
  277.     zeroB = 0;
  278.  
  279.     /* Initialize the condition code byte to zero */
  280.     condcodeB = 0;
  281.  
  282.     /* Determine the total byte length of the floating point number */
  283.     totalenN = fftotlen(expbitN, fracbitN);
  284.  
  285.     /* Determine the number of bytes required to hold the mantissa and
  286.      * exponent.  Allocate space for each... */
  287.     expbyteN = ffexplen(expbitN);
  288.     fracbyteN = ffraclen(fracbitN);
  289.  
  290.     fltfracBP = (unsigned char *) FCALLOC(fracbyteN, 1, "FLTOSTR1");
  291.     fltexpBP = (unsigned char *) FMALLOC(expbyteN, "FLTOSTR2");
  292.     expbiasBP = (unsigned char *) FCALLOC(expbyteN, 1, "FLTOSTR3");
  293.     exponeBP = (unsigned char *) FCALLOC(expbyteN, 1, "FLTOSTR4");
  294.  
  295.     /* Isolate the mantissas, exponents, and signs.  Calling ffextall()
  296.      * also sets the implied ms bit of the mantissa to 1. */
  297.     ffextall(fltvalBP, totalenN, fracbitN, fracbyteN, expbitN, expbyteN,
  298.          fltfracBP, fltexpBP, &fltsignB);
  299.  
  300.     /* Write sign bit and decimal point to the output string. */
  301.     if (fltsignB == 0)
  302.     {
  303.     *strBP++ = '+';
  304.     }
  305.     else
  306.     {
  307.     *strBP++ = '-';
  308.     condcodeB = FFNEG;
  309.     ffbitclr(fltvalBP, totalenN, (fracbitN + expbitN));
  310.     }
  311.     *strBP++ = '.';
  312.  
  313.     *(exponeBP + (expbyteN - 1)) = 1;
  314.  
  315.     /* Check the type of floating point number that we have: zero, infinity,
  316.      * or Not-A-Number.  First check for zero value. */
  317.     if (ffchkzero(fltexpBP, expbyteN) == 0)
  318.     {
  319.     /* The exponent value is zero.  Check if the mantissa is also zero.
  320.      * First clear the implied most significant mantissa bit. */
  321.     ffbitclr(fltfracBP, fracbyteN, fracbitN);
  322.     if (ffchkzero(fltfracBP, fracbyteN) == 0)
  323.     {
  324.         /* The mantissa value is also zero, we have a zero value. */
  325.         strcpy((char *) strBP, "0e+0");
  326.         FFREE(fltfracBP);
  327.         FFREE(fltexpBP);
  328.         FFREE(expbiasBP);
  329.         FFREE(exponeBP);
  330.         return((unsigned char) (condcodeB | FFZERO));
  331.     }
  332.     /* Set the implied most significant mantissa bit to 1 */
  333.     ffbitset(fltfracBP, fracbyteN, fracbitN);
  334.     }
  335.  
  336.     /* Having the exponent bias is useful for the next step. */
  337.     ffgenbias(expbiasBP, expbyteN, expbitN);
  338.  
  339.     /* Check if the exponent value is set to the maximum possible
  340.      * value.  This is done by making a copy of the exponent bias,
  341.      * shifting it left once, and then setting the least significant
  342.      * bit to one.  The result is compared with the exponent value. */
  343.     tempexpBP = (unsigned char *) FMALLOC(expbyteN, "FLTOSTR5");
  344.     for (i = 0; i < expbyteN; i++)
  345.     *(tempexpBP + i) = *(expbiasBP + i);
  346.     ushftlm(tempexpBP, expbyteN);
  347.     *(tempexpBP + expbyteN - 1) |= 1;
  348.     if (ucmpm(tempexpBP, fltexpBP, expbyteN) == 0)
  349.     {
  350.     /* The exponent value is set to its maximum value.
  351.      * First clear the implied most significant mantissa bit. */
  352.     ffbitclr(fltfracBP, fracbyteN, fracbitN);
  353.     if (ffchkzero(fltfracBP, fracbyteN) == 0)
  354.     {
  355.         /* The mantissa value is zero, we have an Infinite value. */
  356.         strcpy((char *) (--strBP), "Infinity");
  357.         condcodeB |= FFINF;
  358.     }
  359.     else
  360.     {
  361.         /* The mantissa value is non-zero, we have a Not-A-Number. */
  362.         strcpy((char *) (--strBP), "Not-A-Number");
  363.         condcodeB |= FFNAN;
  364.     }
  365.     FFREE(fltfracBP);
  366.     FFREE(fltexpBP);
  367.     FFREE(expbiasBP);
  368.     FFREE(exponeBP);
  369.     FFREE(tempexpBP);
  370.     return(condcodeB);
  371.     }
  372.  
  373.     /* Ok floating point value. */
  374.     FFREE(tempexpBP);
  375.  
  376.     /* Back to working on the exponent...
  377.      * Subtract off exponent bias, note that if the result is negative, the
  378.      * sign is properly extended.  This is important since it allows
  379.      * exponent overflow and underflow to be detected much more easily. */
  380.     isubm(fltexpBP, expbiasBP, expbyteN);
  381.  
  382.     /* Convert the power of 2 exponent into an approximate power of 10
  383.      * exponent.  This is done by converting the exponent value into a float
  384.      * and then multiplying it by log10(2).  The result is then converter
  385.      * back to a integer. */
  386.     tempfltBP = (unsigned char *) FMALLOC(SINGLETOT, "FLTOSTR6");
  387.     intoflt(fltexpBP, expbyteN, tempfltBP, SINGLEFRAC, SINGLEXP);
  388.     fmultm(tempfltBP, log10of2AB, SINGLEXP, SINGLEFRAC);
  389.     fltoint(tempfltBP, SINGLEFRAC, SINGLEXP, fltexpBP, expbyteN);
  390.  
  391.     /* Add one to the power of 10 exponent */
  392.     iaddm(fltexpBP, exponeBP, expbyteN);
  393.     
  394.     /* Convert the exponent power of ten into a float value */
  395.     ptenfltBP = (unsigned char *) FMALLOC(totalenN, "FLTOSTR7");
  396.     tempexpBP = (unsigned char *) FMALLOC(expbyteN, "FLTOSTR8");
  397.     for (i = 0; i < expbyteN; i++)
  398.     *(tempexpBP + i) = *(fltexpBP + i);
  399.  
  400.     ptenccB = intopten(tempexpBP, expbyteN, ptenfltBP, fracbitN, expbitN);
  401.  
  402.     /* Check if either an overflow or underflow occurred. */
  403.     if (((ptenccB & FFINF) == FFINF) || ((ptenccB & FFZERO) == FFZERO))
  404.     {
  405.     if ((ptenccB & FFINF) == FFINF)
  406.     {
  407.         /* The power of ten is a bit too big... */
  408.         isubm(fltexpBP, exponeBP, expbyteN);
  409.     }
  410.     else
  411.     {
  412.         /* The power of ten is a bit too small... */
  413.         iaddm(fltexpBP, exponeBP, expbyteN);
  414.     }
  415.     /* Re initialize temporary copy of exponent value, and recalculate
  416.      * the modified power of ten. */
  417.     for (i = 0; i < expbyteN; i++)
  418.         *(tempexpBP + i) = *(fltexpBP + i);
  419.     ptenccB = intopten(tempexpBP, expbyteN, ptenfltBP, fracbitN, expbitN);
  420.     }
  421.     FFREE(tempexpBP);
  422.  
  423.     /* Divide the original float by the power of 10. */
  424.     fdivm(fltvalBP, ptenfltBP, expbitN, fracbitN);
  425.  
  426.     /* Check if the result is less than 1/10.  First generate
  427.      * 1/10 in the desired floating point format. */
  428.     onetenthBP = (unsigned char *) FMALLOC(totalenN, "FLTOSTR9");
  429.     intopten(&minusoneB, 1, onetenthBP, fracbitN, expbitN);
  430.     if ((i = fcmpm(fltvalBP, onetenthBP, expbitN, fracbitN)) < 0)
  431.     {
  432.     /* If so, divide the result by 1/10 (same as multiplying it by 10)
  433.     ** and subtract 1 from the power of ten. */
  434.     fdivm(fltvalBP, onetenthBP, expbitN, fracbitN);
  435.     /* Subtract one from the power of 10 exponent */
  436.     isubm(fltexpBP, exponeBP, expbyteN);
  437.     }
  438.     else if (i == 0)
  439.     {
  440.     /* The result is equal to 1/10...  This is a special boundary
  441.      * case that needs to be handled by brute force... */
  442.     strcpy((char *) strBP, "1e+1");
  443.     FFREE(fltfracBP);
  444.     FFREE(fltexpBP);
  445.     FFREE(expbiasBP);
  446.     FFREE(exponeBP);
  447.     FFREE(tempfltBP);
  448.     FFREE(ptenfltBP);
  449.     FFREE(onetenthBP);
  450.     return(condcodeB);
  451.     }
  452.  
  453.     /* Check if the result is greater than 1.  First generate 1 in the
  454.      * desired floating point format. */
  455.     oneBP = (unsigned char *) FMALLOC(totalenN, "FLTOSTR10");
  456.     intopten(&zeroB, 1, oneBP, fracbitN, expbitN);
  457.     if (fcmpm(fltvalBP, oneBP, expbitN, fracbitN) > 0)
  458.     {
  459.     /* If so, multiply the result by 1/10 and add 1 to the power
  460.      * of ten. */
  461.     fmultm(fltvalBP, onetenthBP, expbitN, fracbitN);
  462.     /* Add one to the power of 10 exponent */
  463.     iaddm(fltexpBP, exponeBP, expbyteN);
  464.     }
  465.  
  466.     /* At this point, fltvalBP points to a floating point number between 1/10
  467.      * and 1. */
  468.  
  469.     /* Convert binary length of the mantissa into decimal digit length.
  470.      * This gives the number of significant digits and determines the power
  471.      * of then to multiply the mantissa by to convert it to an integer. */
  472.     mantlenN = (fracbitN + 1);
  473.  
  474.     /* Convert mantissa bit length from an int to a format compatible with
  475.      * two byte integers. */
  476.     mantlenBP = (unsigned char *) FMALLOC(2, "FLTOSTR11");
  477.     *(mantlenBP+1) = (unsigned char) mantlenN;
  478.     *mantlenBP = (unsigned char) (mantlenN >> 8);
  479.  
  480.     /* Now do conversion from binary length to decimal digit length. */
  481.     intoflt(mantlenBP, 2, tempfltBP, SINGLEFRAC, SINGLEXP);
  482.     fmultm(tempfltBP, log10of2AB, SINGLEXP, SINGLEFRAC);
  483.     fltoint(tempfltBP, SINGLEFRAC, SINGLEXP, mantlenBP, 2);
  484.     FFREE(tempfltBP);
  485.  
  486.     /* Convert result back to int */
  487.     mantlenN = (((unsigned int) (*mantlenBP)) << 8) +
  488.     (unsigned int) (*(mantlenBP+1));
  489.     mantlenN++;
  490.     /* And then back to a two byte integer */
  491.     *(mantlenBP+1) = (unsigned char) mantlenN;
  492.     *mantlenBP = (unsigned char) (mantlenN >> 8);
  493.  
  494.     /* Calculate appropriate power of ten to multiply mantissa by and
  495.      * then multiply it.  Convert the result to an int, and then to an
  496.      * ASCII string. */
  497.     tempfltBP = (unsigned char *) FMALLOC(totalenN, "FLTOSTR12");
  498.     intopten(mantlenBP, 2, tempfltBP, fracbitN, expbitN);
  499.     fmultm(fltvalBP, tempfltBP, expbitN, fracbitN);
  500.     mantlenN = mantlenN >> 1;
  501.     mantintBP = (unsigned char *) FMALLOC(mantlenN, "FLTOSTR13");
  502.     fltoint(fltvalBP, fracbitN, expbitN, mantintBP, mantlenN);
  503.     condcodeB = intostr(mantintBP, mantlenN, strBP, (slenN - 2), 10);
  504.  
  505.     /* Free the temporary buffer for the two byte integer mantissa length */
  506.     FFREE(mantlenBP);
  507.  
  508.     /* Add on exponent part */
  509.     strBP += strlen((char *) strBP);
  510.     slenN -= (2 + strlen((char *) strBP));
  511.     *strBP++ = 'e';
  512.     /* Determine the exponent sign. */
  513.     if (((*(fltexpBP + expbyteN - 1)) & 0x80) != 0x80)
  514.     *strBP++ = '+';
  515.     condcodeB = intostr(fltexpBP, expbyteN, strBP, slenN, 10);
  516.     
  517.     FFREE(fltfracBP);
  518.     FFREE(fltexpBP);
  519.     FFREE(expbiasBP);
  520.     FFREE(exponeBP);
  521.     FFREE(tempfltBP);
  522.     FFREE(ptenfltBP);
  523.     FFREE(onetenthBP);
  524.     FFREE(oneBP);
  525.     FFREE(mantintBP);
  526.     return(condcodeB);
  527. }
  528.  
  529. /* Function:    unsigned char *intopten(unsigned char *fltexpBP, int expbyteN,
  530. **                    unsigned char *ptenfltBP, int fracbitN,
  531. **                    int expbitN)
  532. **
  533. ** This function calulates 10 to the integer value pointed to by fltexpBP.
  534. ** fltexpBP points to an integer value that is expbyteN bytes long.
  535. **
  536. ** The float result is written to the memory buffer pointed to by ptenfltBP.
  537. ** fracbitN and expbitN give the number of bits in the floating point's
  538. ** mantissa and exponent.
  539. */
  540. unsigned char
  541. #ifdef PROTOTYPES
  542. intopten(unsigned char *fltexpBP, int expbyteN, unsigned char *ptenfltBP,
  543.      int fracbitN, int expbitN)
  544. #else
  545. intopten(fltexpBP, expbyteN, ptenfltBP, fracbitN, expbitN)
  546. unsigned char *fltexpBP;
  547. int expbyteN;
  548. unsigned char *ptenfltBP;
  549. int fracbitN;
  550. int expbitN;
  551. #endif
  552. {
  553.     unsigned char condcodeB;
  554.     unsigned char *ptenBP;
  555.     unsigned char tenB;
  556.     unsigned char oneB;
  557.     int totalenN;
  558.  
  559.     /* Initialize one and ten values */
  560.     oneB = 1;
  561.     tenB = 10;
  562.  
  563.     /* Initialize the condition code byte to zero */
  564.     condcodeB = 0;
  565.  
  566.     /* Determine the total byte length of the floating point number */
  567.     totalenN = fftotlen(expbitN, fracbitN);
  568.  
  569.     /* Allocate space for the power of ten multiplier/divisor and initialize
  570.      * it to 10. */
  571.     ptenBP = (unsigned char *) FMALLOC(totalenN, "INTOPTEN");
  572.     intoflt(&tenB, 1, ptenBP, fracbitN, expbitN);
  573.     /* Initialize the result to 1. */
  574.     intoflt(&oneB, 1, ptenfltBP, fracbitN, expbitN);
  575.  
  576.     /* Check the sign bit of the exponent */
  577.     if (((*fltexpBP) & 0x80) == 0)
  578.     {
  579.     /* Exponent is positive.  Multiply result value by binary power of 10
  580.      * for each corresponding bit set in the exponent.  Loop until the
  581.      * exponent value is zero. */
  582.     while (ucheckm(fltexpBP, expbyteN) != 0)
  583.     {
  584.         /* Check if the least significant bit is one. */
  585.         if (((*(fltexpBP + expbyteN - 1)) & 1) != 0)
  586.         {
  587.         /* If so, multiply in the power of ten. */
  588.         condcodeB = fmultm(ptenfltBP, ptenBP, expbitN, fracbitN);
  589.         /* Check if an overflow has occurred */
  590.         if ((condcodeB & FFINF) == FFINF)
  591.         {
  592.             FFREE(ptenBP);
  593.             return(condcodeB);
  594.         }
  595.         }
  596.         ushftrm(fltexpBP, expbyteN);
  597.         /* Generate the next binary power of 10 by squaring the previous
  598.          * result. */
  599.         fmultm(ptenBP, ptenBP, expbitN, fracbitN);
  600.     }
  601.     }
  602.     else
  603.     {
  604.     /* Exponent is negative.  Change the sign of the exponent value and
  605.      * then divide the result value by binary power of 10 for each
  606.      * corresponding bit set in the exponent.  Loop until the exponent
  607.      * value is zero. */
  608.     inegm(fltexpBP, expbyteN);
  609.     while (ucheckm(fltexpBP, expbyteN) != 0)
  610.     {
  611.         /* Check if the least significant bit is one. */
  612.         if (((*(fltexpBP + expbyteN - 1)) & 1) != 0)
  613.         {
  614.         /* If so, then divide the previous result by power of ten. */
  615.             fdivm(ptenfltBP, ptenBP, expbitN, fracbitN);
  616.         /* Check if an underflow has occurred */
  617.         if ((condcodeB & FFZERO) == FFZERO)
  618.         {
  619.             FFREE(ptenBP);
  620.             return(condcodeB);
  621.         }
  622.         }
  623.         ushftrm(fltexpBP, expbyteN);
  624.         /* Generate the next binary power of 10 by squaring the previous
  625.          * result. */
  626.         fmultm(ptenBP, ptenBP, expbitN, fracbitN);
  627.     }
  628.     }
  629.  
  630.     /* Clean up and return. */
  631.     FFREE(ptenBP);
  632.     return(condcodeB);
  633. }
  634.  
  635. /* Function:    unsigned char *strtoint(char *strBP, int slenN,
  636. **                        unsigned char *intvalBP, int intlenN,
  637. **                    int radixN)
  638. **
  639. ** This function converts the ASCII number string representation pointed to
  640. ** by strBP to a binary integer value.  intlenN is the length of the integer
  641. ** result in bytes.  strBP is the length of the string buffer in bytes.
  642. ** radixN is the radix of the input number.  Any radix from 2 to 36 is
  643. ** valid.  Valid input digits are 0 thru (radix - 1).  For radixes greater
  644. ** than 10, alphabet characters are used starting with 'a'.  Upper and lower
  645. ** case are treated the same.
  646. */
  647. unsigned char
  648. #ifdef PROTOTYPES
  649. strtoint(char *strBP, int slenN, unsigned char *intvalBP, int intlenN,
  650.      int radixN)
  651. #else
  652. strtoint(strBP, slenN, intvalBP, intlenN, radixN)
  653. char *strBP;
  654. int slenN;
  655. unsigned char *intvalBP;
  656. int intlenN;
  657. int radixN;
  658. #endif
  659. {
  660.     unsigned char *radixBP;
  661.     unsigned char *indigitBP;
  662.     int i;
  663.     unsigned char tempB, signB;
  664.     unsigned char condcodeB;
  665.  
  666.     condcodeB = 0;            /* Assume positive value */
  667.  
  668.     /* Convert the radix value to multiple byte integer format. */
  669.     radixBP = FCALLOC(intlenN, 1, "STRTOINT1");
  670.     if (radixN > 36) radixN = 36;
  671.     if (radixN < 2) radixN = 2;
  672.     *(radixBP + intlenN - 1) = (unsigned char)radixN;
  673.  
  674.     /* Get buffer to read each digit into. */
  675.     indigitBP = FCALLOC(intlenN, 1, "STRTOINT2");
  676.  
  677.     /* Initialize the output value to zero. */
  678.     for (i = 0; i < intlenN; i++)
  679.     *(intvalBP + i) = 0;
  680.  
  681.     /* Check for leading +/- character */
  682.     tempB = *strBP;
  683.     if (tempB == '+')
  684.     {
  685.     strBP++;
  686.     slenN--;
  687.     }
  688.     else if (tempB == '-')
  689.     {
  690.     strBP++;
  691.     slenN--;
  692.     condcodeB = SIGN;
  693.     }
  694.  
  695.     /* Loop thru digits until either the end of the string, the output value
  696.      * overflows, or a non digit character is encountered. */
  697.     while ((slenN-- > 0) && (*strBP != '\0'))
  698.     {
  699.     /* Convert digits from ASCII to binary. */
  700.     tempB = *strBP++;
  701.     if ((tempB >= '0') && (tempB <= '9'))
  702.     {
  703.         tempB -= '0';
  704.     }
  705.     else if ((tempB >= 'A') && (tempB <= 'Z'))
  706.     {
  707.         tempB -= ('A' - (char)10);
  708.     }
  709.     else if ((tempB >= 'a') && (tempB <= 'z'))
  710.     {
  711.         tempB -= ('a' - (char)10);
  712.     }
  713.     else
  714.     {
  715.         FFREE(radixBP);
  716.         FFREE(indigitBP);
  717.         return(condcodeB);
  718.     }
  719.  
  720.     /* If the converted digit is greater than the selected radix, return */
  721.     if (tempB >= radixN)
  722.     {
  723.         FFREE(radixBP);
  724.         FFREE(indigitBP);
  725.         return(condcodeB);
  726.     }
  727.  
  728.     /* Convert the digit value to multiple byte integer format. */
  729.     *(indigitBP + intlenN - 1) = (unsigned char) tempB;
  730.  
  731.         /* Multiply the result value by the radix and add in the new digit. */
  732.     if (((umultm(intvalBP, radixBP, intlenN) & (OVERFLOW | SIGN)) != 0) ||
  733.         ((iaddm(intvalBP, indigitBP, intlenN) & OVERFLOW) != 0))
  734.         {
  735.         /* Overflow has occurred. */
  736.         FFREE(radixBP);
  737.         FFREE(indigitBP);
  738.         if (condcodeB == SIGN)
  739.         {
  740.         /* Negative overflow, return the most negative possible
  741.          * integer. */
  742.         if (*strBP == '\0')
  743.         {
  744.             /* Check if overflow is real, not just due to using
  745.              * unsigned arithmetic.  Unsigned arithmetic causes
  746.              * problems with the minimum negative integer possible. */
  747.             for (i = 1; i < intlenN; i++)
  748.             {
  749.             if (*(intvalBP + i) != 0)
  750.                 break;
  751.             }
  752.             if ((i == intlenN) && (*intvalBP == 0x80))
  753.             return( SIGN );
  754.         }
  755.         *intvalBP = 0x80;
  756.         for (i = 1; i < intlenN; i++)
  757.             *(intvalBP + i) = 0;
  758.         }
  759.         else
  760.         {
  761.         /* Positive overflow, return the most positive possible
  762.          * integer. */
  763.         *intvalBP = 0x7f;
  764.         for (i = 1; i < intlenN; i++)
  765.             *(intvalBP + i) = 0xff;
  766.         }
  767.         return((unsigned char) (condcodeB | OVERFLOW));
  768.     }
  769.     /* Reset the radix value after doing the multiply. */
  770.     *(radixBP + intlenN - 1) = (unsigned char)radixN;
  771.     }
  772.     FFREE(radixBP);
  773.     FFREE(indigitBP);
  774.  
  775.     /* Change sign of result if the input value non-zero and was preceeded
  776.      * by a '-'. */
  777.     if (ucheckm(intvalBP, intlenN) == 0)
  778.         return(ZERO);
  779.     if ((condcodeB & SIGN) != 0)
  780.     inegm(intvalBP, intlenN);
  781.     return(condcodeB);
  782. }
  783.  
  784. /* Function:    unsigned char strtoflt(char *strBP, int slenN,
  785. **                       unsigned char *fltvalBP, int fracbitN,
  786. **                       int expbitN)
  787. **
  788. ** This function converts the decimal ASCII string representation pointed
  789. ** to by strBP to a binary floating point value pointed to fltvalBP
  790. ** fracbitN is the length of the floating point mantissa in bits.  expbitN
  791. ** is the length of the floating point exponent in bits.  strBP is the
  792. ** length of the string buffer in bytes.
  793. */
  794. unsigned char
  795. #ifdef PROTOTYPES
  796. strtoflt(char *strBP, int slenN, unsigned char *fltvalBP, int fracbitN,
  797.      int expbitN)
  798. #else
  799. strtoflt(strBP, slenN, fltvalBP, fracbitN, expbitN)
  800. char *strBP;
  801. int slenN;
  802. unsigned char *fltvalBP;
  803. int fracbitN;
  804. int expbitN;
  805. #endif
  806. {
  807.     int i;
  808.     unsigned char tempB, signB;
  809.     unsigned char dpflagB;
  810.     unsigned char *tenBP, *tempBP;
  811.     unsigned char tenB;
  812.     int totalenN;
  813.     int expbyteN;
  814.     unsigned char *expintBP;
  815.     unsigned char *dpfactorBP;
  816.     unsigned char *ptenfltBP;
  817.     unsigned char condcodeB;
  818.  
  819.     signB = 0;            /* Assume positive value */
  820.     dpflagB = 0xff;        /* Flag set to 0 when decimal point read in */
  821.     tenB = 10;            /* Integer value of ten */
  822.  
  823.     /* Initialize the condition code byte to zero */
  824.     condcodeB = 0;
  825.  
  826.     /* Determine the total byte length of the floating point number */
  827.     totalenN = fftotlen(expbitN, fracbitN);
  828.  
  829.     /* Allocate space for the power of ten multiplier/divisor and initialize
  830.      * it to 10. */
  831.     tenBP = (unsigned char *) FMALLOC(totalenN, "STRTOFLT1");
  832.     intoflt(&tenB, 1, tenBP, fracbitN, expbitN);
  833.  
  834.     /* Initialize the result to 0. */
  835.     for (i = 0; i < totalenN; i++)
  836.     *(fltvalBP + i) = 0;
  837.  
  838.     /* Check for leading +/- character */
  839.     tempB = *strBP;
  840.     if (tempB == '+')
  841.     {
  842.     strBP++;
  843.     slenN--;
  844.     }
  845.     else if (tempB == '-')
  846.     {
  847.     strBP++;
  848.     slenN--;
  849.     signB = 1;
  850.     condcodeB = FFNEG;
  851.     }
  852.  
  853.     /* Loop thru digits until either the end of the string, the output value
  854.      * overflows, or a non digit character is encountered.  This routine
  855.      * handles decimal '.' point characters. */
  856.     tempBP = (unsigned char *) FMALLOC(totalenN, "STRTOFLT2");
  857.     while ((slenN-- > 0) && (*strBP != '\0'))
  858.     {
  859.     /* Convert digits from ASCII to binary. */
  860.     tempB = *strBP++;
  861.     if ((tempB >= '0') && (tempB <= '9'))
  862.     {
  863.         /* Convert the digit from ASCII to binary. */
  864.         tempB -= '0';
  865.  
  866.         /* If a decimal point has been encountered, use the flag to
  867.          * count the number of digits behind the decimal point. */
  868.         if (dpflagB != 0xff)
  869.         dpflagB++;
  870.  
  871.         for (i = 0; i < totalenN; i++)
  872.         *(tempBP + i) = *(tenBP + i);
  873.  
  874.         /* Multiply the previous result by 10, note that if the previous
  875.          * result overflows, then return an overflow.  This should be
  876.          * a rather rare occurence, since even for single precision, the
  877.          * mantissa must more than 38 digits long! */
  878.         if ((fmultm(fltvalBP, tempBP, expbitN, fracbitN) & FFINF) != 0)
  879.         {
  880.         FFREE(tempBP);
  881.         FFREE(tenBP);
  882.         return((unsigned char) (condcodeB | FFINF));
  883.         }
  884.         /* Convert the digit to a float. */
  885.         intoflt(&tempB, 1, tempBP, fracbitN, expbitN);
  886.         /* Add the digit into the result */
  887.         if ((faddm(fltvalBP, tempBP, expbitN, fracbitN ) & FFINF) != 0)
  888.         {
  889.         FFREE(tempBP);
  890.         FFREE(tenBP);
  891.         return((unsigned char) (condcodeB | FFINF));
  892.         }
  893.     }
  894.     else if (tempB == '.')
  895.         dpflagB = 0;        /* Flag we have hit the dec. pt. */
  896.     else if ((tempB == 'e') || (tempB == 'E'))
  897.         break;            /* Done with mantissa */
  898.     else
  899.     {
  900.         FFREE(tempBP);
  901.         FFREE(tenBP);
  902.         return(condcodeB);
  903.     }
  904.     }
  905.     if (dpflagB == 0xff)
  906.     dpflagB = 0;        /* Reset dec. pt. flag for no adjustment. */
  907.  
  908.     /* The ASCII mantissa has been read in.  Note that the above routine
  909.      * continue to multiply the previous mantissa value by 10, even for
  910.      * digits behind the decimal point.  As a result, the mantissa value
  911.      * in fltvalPB probably too big.  If so, it is too big by a factor
  912.      * of 10^(dpflagB).  For now, read in the exponent (if any) first.
  913.      * The "decimal point" factor will be included when the exponent power
  914.      * of ten is multiplied in. (dbflagB will be subtracted from the value
  915.      * of the exponent read in. */
  916.  
  917.     /* Determine the number of bytes required to hold the exponent.
  918.      *  Allocate space for it and zero it it. */
  919.     expbyteN = ffexplen(expbitN);
  920.     expintBP = (unsigned char *) FCALLOC(expbyteN, 1, "STRTOFLT3");
  921.     
  922.     if ((slenN > 0) && (*strBP != '\0'))
  923.     {
  924.     strtoint(strBP, slenN--, expintBP, expbyteN, 10);
  925.     }
  926.  
  927.     /* The ASCII exponent has been read in (if any).  Include the mantissa
  928.      * decimal point factor in the exponent. */
  929.     dpfactorBP = (unsigned char *) FCALLOC(expbyteN, 1, "STRTOFLT4");
  930.     *(dpfactorBP + expbyteN - 1) = dpflagB;
  931.  
  932.     /* It is highly possible that the exponent value with the decimal point
  933.      * factor subtracted off will result in an underflow, yet the entire
  934.      * number is a valid floating point number. */
  935.     if ((isubm(expintBP, dpfactorBP, expbyteN) & SIGN) != 0)
  936.     {
  937.     /* The exponent is negative, and an underflow is possible.  To
  938.      * avoid an unecessary underflow, divide the exponent value by
  939.      * two first, and then multiply the two halves into the mantissa
  940.      * value, one after the other.  If an underflow occurs then, it
  941.      * is legitimate.
  942.      *
  943.      * Check if the exponent is odd, if so, then multiply the mantissa
  944.      * by 10 first. */
  945.     if ((*(expintBP + (expbyteN - 1)) & 1) != 0)
  946.         fmultm(fltvalBP, tenBP, expbitN, fracbitN);
  947.  
  948.     /* Now divide the negative exponent by two. */
  949.     ushftrm(expintBP, expbyteN);
  950.     /* Remember to set the sign bit... */
  951.     *expintBP |= 0x80;
  952.  
  953.     /* Convert the power of 10 exponent (divided by two) into a binary
  954.      * floating point value. */
  955.     ptenfltBP = (unsigned char *) FMALLOC(totalenN, "STRTOFLT5");
  956.     intopten(expintBP, expbyteN, ptenfltBP, fracbitN, expbitN);
  957.  
  958.     /* Multiply the exponent float value times the mantissa float
  959.      * value to get a partial result. */
  960.     fmultm(fltvalBP, ptenfltBP, expbitN, fracbitN);
  961.     }
  962.     else
  963.     {
  964.     /* The exponent is positive, convert the power of 10 exponent into
  965.      * a binary floating point value. */
  966.     ptenfltBP = (unsigned char *) FMALLOC(totalenN, "STRTOFLT6");
  967.     condcodeB =
  968.         intopten(expintBP, expbyteN, ptenfltBP, fracbitN, expbitN);
  969.     }
  970.  
  971.     if ((condcodeB & FFINF) == 0)
  972.     {
  973.     /* Multiply the exponent float value times the mantissa float value
  974.     * to get the final result, only if there has been no overflow. */
  975.     condcodeB = fmultm(fltvalBP, ptenfltBP, expbitN, fracbitN);
  976.     }
  977.  
  978.     /* Set the sign bit if necessary. */
  979.     if (signB != 0)
  980.     {
  981.     ffbitset(fltvalBP, totalenN, (fracbitN + expbitN));
  982.     condcodeB |= FFNEG;
  983.     }
  984.  
  985.     FFREE(ptenfltBP);
  986.     FFREE(dpfactorBP);
  987.     FFREE(expintBP);
  988.     FFREE(tempBP);
  989.     FFREE(tenBP);
  990.     return(condcodeB);
  991. }
  992.  
  993.