common/strtod.c

   1 /*
   2  * strtod.c --
   3  *
   4  *      Source code for the "strtod" library procedure.
   5  *
   6  * Copyright (c) 1988-1993 The Regents of the University of California.
   7  * Copyright (c) 1994 Sun Microsystems, Inc.
   8  *
   9  * Permission to use, copy, modify, and distribute this
  10  * software and its documentation for any purpose and without
  11  * fee is hereby granted, provided that the above copyright
  12  * notice appear in all copies.  The University of California
  13  * makes no representations about the suitability of this
  14  * software for any purpose.  It is provided "as is" without
  15  * express or implied warranty.
  16  *
  17  * RCS: @(#) $Id$
  18  */
  19
  20 #include <ctype.h>
  21 #include <errno.h>
  22 extern  int     errno;
  23
  24 #ifndef __STDC__
  25 # ifdef __GNUC__
  26 #  define const __const__
  27 # else
  28 #  define const
  29 # endif
  30 #endif
  31
  32 #ifndef TRUE
  33 #define TRUE 1
  34 #define FALSE 0
  35 #endif
  36 #ifndef NULL
  37 #define NULL 0
  38 #endif
  39
  40 static int maxExponent = 511;   /* Largest possible base 10 exponent.  Any
  41                                  * exponent larger than this will already
  42                                  * produce underflow or overflow, so there's
  43                                  * no need to worry about additional digits.
  44                                  */
  45 static double powersOf10[] = {  /* Table giving binary powers of 10.  Entry */
  46     10.,                        /* is 10^2^i.  Used to convert decimal */
  47     100.,                       /* exponents into floating-point numbers. */
  48     1.0e4,
  49     1.0e8,
  50     1.0e16,
  51     1.0e32,
  52     1.0e64,
  53     1.0e128,
  54     1.0e256
  55 };
  56
  57 /*
  58  *----------------------------------------------------------------------
  59  *
  60  * strtod --
  61  *
  62  *      This procedure converts a floating-point number from an ASCII
  63  *      decimal representation to internal double-precision format.
  64  *
  65  * Results:
  66  *      The return value is the double-precision floating-point
  67  *      representation of the characters in string.  If endPtr isn't
  68  *      NULL, then *endPtr is filled in with the address of the
  69  *      next character after the last one that was part of the
  70  *      floating-point number.
  71  *
  72  * Side effects:
  73  *      None.
  74  *
  75  *----------------------------------------------------------------------
  76  */
  77
  78 double
  79 sc_strtod(
  80     const char *string,         /* A decimal ASCII floating-point number,
  81                                  * optionally preceded by white space.
  82                                  * Must have form "-I.FE-X", where I is the
  83                                  * integer part of the mantissa, F is the
  84                                  * fractional part of the mantissa, and X
  85                                  * is the exponent.  Either of the signs
  86                                  * may be "+", "-", or omitted.  Either I
  87                                  * or F may be omitted, or both.  The decimal
  88                                  * point isn't necessary unless F is present.
  89                                  * The "E" may actually be an "e".  E and X
  90                                  * may both be omitted (but not just one).
  91                                  */
  92     char **endPtr               /* If non-NULL, store terminating character's
  93                                  * address here. */
  94         )
  95 {
  96     int sign, expSign = FALSE;
  97     double fraction, dblExp, *d;
  98     register const char *p;
  99     register int c;
 100     int exp = 0;                /* Exponent read from "EX" field. */
 101     int fracExp = 0;            /* Exponent that derives from the fractional
 102                                  * part.  Under normal circumstatnces, it is
 103                                  * the negative of the number of digits in F.
 104                                  * However, if I is very long, the last digits
 105                                  * of I get dropped (otherwise a long I with a
 106                                  * large negative exponent could cause an
 107                                  * unnecessary overflow on I alone).  In this
 108                                  * case, fracExp is incremented one for each
 109                                  * dropped digit. */
 110     int mantSize;               /* Number of digits in mantissa. */
 111     int decPt;                  /* Number of mantissa digits BEFORE decimal
 112                                  * point. */
 113     const char *pExp;           /* Temporarily holds location of exponent
 114                                  * in string. */
 115
 116     /*
 117      * Strip off leading blanks and check for a sign.
 118      */
 119
 120     p = string;
 121     while (isspace(*p)) {
 122         p += 1;
 123     }
 124     if (*p == '-') {
 125         sign = TRUE;
 126         p += 1;
 127     } else {
 128         if (*p == '+') {
 129             p += 1;
 130         }
 131         sign = FALSE;
 132     }
 133
 134     /*
 135      * Count the number of digits in the mantissa (including the decimal
 136      * point), and also locate the decimal point.
 137      */
 138
 139     decPt = -1;
 140     for (mantSize = 0; ; mantSize += 1)
 141     {
 142         c = *p;
 143         if (!isdigit(c)) {
 144             if ((c != '.') || (decPt >= 0)) {
 145                 break;
 146             }
 147             decPt = mantSize;
 148         }
 149         p += 1;
 150     }
 151
 152     /*
 153      * Now suck up the digits in the mantissa.  Use two integers to
 154      * collect 9 digits each (this is faster than using floating-point).
 155      * If the mantissa has more than 18 digits, ignore the extras, since
 156      * they can't affect the value anyway.
 157      */
 158
 159     pExp  = p;
 160     p -= mantSize;
 161     if (decPt < 0) {
 162         decPt = mantSize;
 163     } else {
 164         mantSize -= 1;                  /* One of the digits was the point. */
 165     }
 166     if (mantSize > 18) {
 167         fracExp = decPt - 18;
 168         mantSize = 18;
 169     } else {
 170         fracExp = decPt - mantSize;
 171     }
 172     if (mantSize == 0) {
 173         fraction = 0.0;
 174         p = string;
 175         goto done;
 176     } else {
 177         int frac1, frac2;
 178         frac1 = 0;
 179         for ( ; mantSize > 9; mantSize -= 1)
 180         {
 181             c = *p;
 182             p += 1;
 183             if (c == '.') {
 184                 c = *p;
 185                 p += 1;
 186             }
 187             frac1 = 10*frac1 + (c - '0');
 188         }
 189         frac2 = 0;
 190         for (; mantSize > 0; mantSize -= 1)
 191         {
 192             c = *p;
 193             p += 1;
 194             if (c == '.') {
 195                 c = *p;
 196                 p += 1;
 197             }
 198             frac2 = 10*frac2 + (c - '0');
 199         }
 200         fraction = (1.0e9 * frac1) + frac2;
 201     }
 202
 203     /*
 204      * Skim off the exponent.
 205      */
 206
 207     p = pExp;
 208     if ((*p == 'E') || (*p == 'e')) {
 209         p += 1;
 210         if (*p == '-') {
 211             expSign = TRUE;
 212             p += 1;
 213         } else {
 214             if (*p == '+') {
 215                 p += 1;
 216             }
 217             expSign = FALSE;
 218         }
 219         while (isdigit(*p)) {
 220             exp = exp * 10 + (*p - '0');
 221             p += 1;
 222         }
 223     }
 224     if (expSign) {
 225         exp = fracExp - exp;
 226     } else {
 227         exp = fracExp + exp;
 228     }
 229
 230     /*
 231      * Generate a floating-point number that represents the exponent.
 232      * Do this by processing the exponent one bit at a time to combine
 233      * many powers of 2 of 10. Then combine the exponent with the
 234      * fraction.
 235      */
 236
 237     if (exp < 0) {
 238         expSign = TRUE;
 239         exp = -exp;
 240     } else {
 241         expSign = FALSE;
 242     }
 243     if (exp > maxExponent) {
 244         exp = maxExponent;
 245         errno = ERANGE;
 246     }
 247     dblExp = 1.0;
 248     for (d = powersOf10; exp != 0; exp >>= 1, d += 1) {
 249         if (exp & 01) {
 250             dblExp *= *d;
 251         }
 252     }
 253     if (expSign) {
 254         fraction /= dblExp;
 255     } else {
 256         fraction *= dblExp;
 257     }
 258
 259 done:
 260     if (endPtr != NULL) {
 261         *endPtr = (char *) p;
 262     }
 263
 264     if (sign) {
 265         return -fraction;
 266     }
 267     return fraction;
 268 }