strtold.c

Algorithms for Image Processing and Computer Vision Source Code

/*
** Copyright (C) 1993 DJ Delorie, 24 Kirsten Ave, Rochester NH 03867-2954
**
** This file is distributed under the terms listed in the document
** "copying.dj", available from DJ Delorie at the address above.
** A copy of "copying.dj" should accompany this file; if not, a copy
** should be available from where this file was obtained.  This file
** may not be distributed without a verbatim copy of "copying.dj".
**
** This file is distributed WITHOUT ANY WARRANTY; without even the implied
** warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
**
** changed strtod to _strtold (long double version)
** speed up and better accuracy for large exponents 
** D. Buerssner
*/

#include <stdlib.h>
#include <ctype.h>

static long double powten[] =
{
    1e1L, 1e2L, 1e4L, 1e8L, 1e16L, 1e32L, 1e64L, 1e128L, 1e256L,
    1e512L, 1e1024L, 1e2048L, 1e4096L
};

long double _strtold(const char *s, char **sret)
{
  long double r;	/* result */
  int e, ne;	/* exponent */
  int sign;	/* +- 1.0 */
  int esign;
  int flags=0;
  int l2powm1;

  r = 0.0L;
  sign = 1;
  e = ne = 0;
  esign = 1;

  while(*s && isspace(*s))
    s++;

  if (*s == '+')
    s++;
  else if (*s == '-')
  {
    sign = -1;
    s++;
  }

  while ((*s >= '0') && (*s <= '9'))
  {
    flags |= 1;
    r *= 10.0L;
    r += *s - '0';
    s++;
  }

  if (*s == '.')
  {
    s++;
    while ((*s >= '0') && (*s <= '9'))
    {
      flags |= 2;
      r *= 10.0L;
      r += *s - '0';
      s++;
      ne++;
    }
  }
  if (flags == 0)
  {
    if (sret) *sret = (char *)s;
    return 0.0L;
  }

  if ((*s == 'e') || (*s == 'E'))
  {
    s++;
    if (*s == '+')
      s++;
    else if (*s == '-')
    {
      s++;
      esign = -1;
    }
    while ((*s >= '0') && (*s <= '9'))
    {
      e *= 10;
      e += *s - '0';
      s++;
    }
  }
  if (esign < 0)
  {
      esign = -esign;
      e = -e;
  }
  e = e - ne;
  if (e < -4096)
  {
      /* possibly subnormal number, 10^e would overflow */
      r *= 1.0e-2048L;
      e += 2048;
  }
  if (e < 0)
  {
      e = -e;
      esign = -esign;
  }
  if (e >= 8192)
      e = 8191;
  if (e)
  {
    long double d = 1.0L;
    l2powm1 = 0;
    while (e)
    {
	 if (e & 1)
	   d *= powten[l2powm1];
	 e >>= 1;
	 l2powm1++;
    }
    if (esign > 0)
	r *= d;
    else
	r /= d;
  }
  if(sret) *sret = (char *)s;
  return r * sign;
}