Fix 695467: Add and use fz_ftoa function (like dtoa but with floats).

The dtoa function is for doubles (which is what MuJS uses) but for MuPDF
we only need and want float precision in our output formatting.

1 files changed, 540 insertions, 0 deletions

diff --git a/source/fitz/strtod.c b/source/fitz/strtod.c
new file mode 100644
index 00000000..f548b6f3
--- /dev/null
+++ b/source/fitz/strtod.c
@@ -0,0 +1,540 @@
+/* The authors of this software are Rob Pike and Ken Thompson.
+ * Copyright (c) 2002 by Lucent Technologies.
+ * Permission to use, copy, modify, and distribute this software for any
+ * purpose without fee is hereby granted, provided that this entire notice
+ * is included in all copies of any software which is or includes a copy
+ * or modification of this software and in all copies of the supporting
+ * documentation for such software.
+ * THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
+ * WARRANTY. IN PARTICULAR, NEITHER THE AUTHORS NOR LUCENT TECHNOLOGIES MAKE ANY
+ * REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY
+ * OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
+ */
+
+#include "mupdf/fitz.h"
+
+#include <stdio.h>
+#include <math.h>
+#include <float.h>
+#include <string.h>
+#include <stdlib.h>
+#include <errno.h>
+
+#ifndef INFINITY
+#define INFINITY (DBL_MAX+DBL_MAX)
+#endif
+#ifndef NAN
+#define NAN (INFINITY-INFINITY)
+#endif
+
+typedef unsigned long ulong;
+
+static ulong
+umuldiv(ulong a, ulong b, ulong c)
+{
+	double d;
+
+	d = ((double)a * (double)b) / (double)c;
+	if(d >= 4294967295.)
+		d = 4294967295.;
+	return (ulong)d;
+}
+
+/*
+ * This routine will convert to arbitrary precision
+ * floating point entirely in multi-precision fixed.
+ * The answer is the closest floating point number to
+ * the given decimal number. Exactly half way are
+ * rounded ala ieee rules.
+ * Method is to scale input decimal between .500 and .999...
+ * with external power of 2, then binary search for the
+ * closest mantissa to this decimal number.
+ * Nmant is is the required precision. (53 for ieee dp)
+ * Nbits is the max number of bits/word. (must be <= 28)
+ * Prec is calculated - the number of words of fixed mantissa.
+ */
+enum
+{
+	Nbits	= 28,				/* bits safely represented in a ulong */
+	Nmant	= 53,				/* bits of precision required */
+	Prec	= (Nmant+Nbits+1)/Nbits,	/* words of Nbits each to represent mantissa */
+	Sigbit	= 1<<(Prec*Nbits-Nmant),	/* first significant bit of Prec-th word */
+	Ndig	= 1500,
+	One	= (ulong)(1<<Nbits),
+	Half	= (ulong)(One>>1),
+	Maxe	= 310,
+
+	Fsign	= 1<<0,		/* found - */
+	Fesign	= 1<<1,		/* found e- */
+	Fdpoint	= 1<<2,		/* found . */
+
+	S0	= 0,		/* _		_S0	+S1	#S2	.S3 */
+	S1,			/* _+		#S2	.S3 */
+	S2,			/* _+#		#S2	.S4	eS5 */
+	S3,			/* _+.		#S4 */
+	S4,			/* _+#.#	#S4	eS5 */
+	S5,			/* _+#.#e	+S6	#S7 */
+	S6,			/* _+#.#e+	#S7 */
+	S7			/* _+#.#e+#	#S7 */
+};
+
+static	int	xcmp(char*, char*);
+static	int	fpcmp(char*, ulong*);
+static	void	frnorm(ulong*);
+static	void	divascii(char*, int*, int*, int*);
+static	void	mulascii(char*, int*, int*, int*);
+
+typedef	struct	Tab	Tab;
+struct	Tab
+{
+	int	bp;
+	int	siz;
+	char*	cmp;
+};
+
+double
+fz_strtod(const char *as, char **aas)
+{
+	int na, ex, dp, bp, c, i, flag, state;
+	ulong low[Prec], hig[Prec], mid[Prec];
+	double d;
+	char *s, a[Ndig];
+
+	flag = 0;	/* Fsign, Fesign, Fdpoint */
+	na = 0;		/* number of digits of a[] */
+	dp = 0;		/* na of decimal point */
+	ex = 0;		/* exonent */
+
+	state = S0;
+	for(s=(char*)as;; s++) {
+		c = *s;
+		if(c >= '0' && c <= '9') {
+			switch(state) {
+			case S0:
+			case S1:
+			case S2:
+				state = S2;
+				break;
+			case S3:
+			case S4:
+				state = S4;
+				break;
+
+			case S5:
+			case S6:
+			case S7:
+				state = S7;
+				ex = ex*10 + (c-'0');
+				continue;
+			}
+			if(na == 0 && c == '0') {
+				dp--;
+				continue;
+			}
+			if(na < Ndig-50)
+				a[na++] = c;
+			continue;
+		}
+		switch(c) {
+		case '\t':
+		case '\n':
+		case '\v':
+		case '\f':
+		case '\r':
+		case ' ':
+			if(state == S0)
+				continue;
+			break;
+		case '-':
+			if(state == S0)
+				flag |= Fsign;
+			else
+				flag |= Fesign;
+		case '+':
+			if(state == S0)
+				state = S1;
+			else
+			if(state == S5)
+				state = S6;
+			else
+				break;	/* syntax */
+			continue;
+		case '.':
+			flag |= Fdpoint;
+			dp = na;
+			if(state == S0 || state == S1) {
+				state = S3;
+				continue;
+			}
+			if(state == S2) {
+				state = S4;
+				continue;
+			}
+			break;
+		case 'e':
+		case 'E':
+			if(state == S2 || state == S4) {
+				state = S5;
+				continue;
+			}
+			break;
+		}
+		break;
+	}
+
+	/*
+	 * clean up return char-pointer
+	 */
+	switch(state) {
+	case S0:
+		if(xcmp(s, "nan") == 0) {
+			if(aas != NULL)
+				*aas = s+3;
+			goto retnan;
+		}
+	case S1:
+		if(xcmp(s, "infinity") == 0) {
+			if(aas != NULL)
+				*aas = s+8;
+			goto retinf;
+		}
+		if(xcmp(s, "inf") == 0) {
+			if(aas != NULL)
+				*aas = s+3;
+			goto retinf;
+		}
+	case S3:
+		if(aas != NULL)
+			*aas = (char*)as;
+		goto ret0;	/* no digits found */
+	case S6:
+		s--;		/* back over +- */
+	case S5:
+		s--;		/* back over e */
+		break;
+	}
+	if(aas != NULL)
+		*aas = s;
+
+	if(flag & Fdpoint)
+	while(na > 0 && a[na-1] == '0')
+		na--;
+	if(na == 0)
+		goto ret0;	/* zero */
+	a[na] = 0;
+	if(!(flag & Fdpoint))
+		dp = na;
+	if(flag & Fesign)
+		ex = -ex;
+	dp += ex;
+	if(dp < -Maxe){
+		errno = ERANGE;
+		goto ret0;	/* underflow by exp */
+	} else
+	if(dp > +Maxe)
+		goto retinf;	/* overflow by exp */
+
+	/*
+	 * normalize the decimal ascii number
+	 * to range .[5-9][0-9]* e0
+	 */
+	bp = 0;		/* binary exponent */
+	while(dp > 0)
+		divascii(a, &na, &dp, &bp);
+	while(dp < 0 || a[0] < '5')
+		mulascii(a, &na, &dp, &bp);
+
+	/* close approx by naive conversion */
+	mid[0] = 0;
+	mid[1] = 1;
+	for(i=0; (c=a[i]) != '\0'; i++) {
+		mid[0] = mid[0]*10 + (c-'0');
+		mid[1] = mid[1]*10;
+		if(i >= 8)
+			break;
+	}
+	low[0] = umuldiv(mid[0], One, mid[1]);
+	hig[0] = umuldiv(mid[0]+1, One, mid[1]);
+	for(i=1; i<Prec; i++) {
+		low[i] = 0;
+		hig[i] = One-1;
+	}
+
+	/* binary search for closest mantissa */
+	for(;;) {
+		/* mid = (hig + low) / 2 */
+		c = 0;
+		for(i=0; i<Prec; i++) {
+			mid[i] = hig[i] + low[i];
+			if(c)
+				mid[i] += One;
+			c = mid[i] & 1;
+			mid[i] >>= 1;
+		}
+		frnorm(mid);
+
+		/* compare */
+		c = fpcmp(a, mid);
+		if(c > 0) {
+			c = 1;
+			for(i=0; i<Prec; i++)
+				if(low[i] != mid[i]) {
+					c = 0;
+					low[i] = mid[i];
+				}
+			if(c)
+				break;	/* between mid and hig */
+			continue;
+		}
+		if(c < 0) {
+			for(i=0; i<Prec; i++)
+				hig[i] = mid[i];
+			continue;
+		}
+
+		/* only hard part is if even/odd roundings wants to go up */
+		c = mid[Prec-1] & (Sigbit-1);
+		if(c == Sigbit/2 && (mid[Prec-1]&Sigbit) == 0)
+			mid[Prec-1] -= c;
+		break;	/* exactly mid */
+	}
+
+	/* normal rounding applies */
+	c = mid[Prec-1] & (Sigbit-1);
+	mid[Prec-1] -= c;
+	if(c >= Sigbit/2) {
+		mid[Prec-1] += Sigbit;
+		frnorm(mid);
+	}
+	goto out;
+
+ret0:
+	return 0;
+
+retnan:
+	return NAN;
+
+retinf:
+	/*
+	 * Unix strtod requires these. Plan 9 would return Inf(0) or Inf(-1). */
+	errno = ERANGE;
+	if(flag & Fsign)
+		return -HUGE_VAL;
+	return HUGE_VAL;
+
+out:
+	d = 0;
+	for(i=0; i<Prec; i++)
+		d = d*One + mid[i];
+	if(flag & Fsign)
+		d = -d;
+	d = ldexp(d, bp - Prec*Nbits);
+	if(d == 0){	/* underflow */
+		errno = ERANGE;
+	}
+	return d;
+}
+
+static void
+frnorm(ulong *f)
+{
+	int i, c;
+
+	c = 0;
+	for(i=Prec-1; i>0; i--) {
+		f[i] += c;
+		c = f[i] >> Nbits;
+		f[i] &= One-1;
+	}
+	f[0] += c;
+}
+
+static int
+fpcmp(char *a, ulong* f)
+{
+	ulong tf[Prec];
+	int i, d, c;
+
+	for(i=0; i<Prec; i++)
+		tf[i] = f[i];
+
+	for(;;) {
+		/* tf *= 10 */
+		for(i=0; i<Prec; i++)
+			tf[i] = tf[i]*10;
+		frnorm(tf);
+		d = (tf[0] >> Nbits) + '0';
+		tf[0] &= One-1;
+
+		/* compare next digit */
+		c = *a;
+		if(c == 0) {
+			if('0' < d)
+				return -1;
+			if(tf[0] != 0)
+				goto cont;
+			for(i=1; i<Prec; i++)
+				if(tf[i] != 0)
+					goto cont;
+			return 0;
+		}
+		if(c > d)
+			return +1;
+		if(c < d)
+			return -1;
+		a++;
+	cont:;
+	}
+}
+
+static void
+divby(char *a, int *na, int b)
+{
+	int n, c;
+	char *p;
+
+	p = a;
+	n = 0;
+	while(n>>b == 0) {
+		c = *a++;
+		if(c == 0) {
+			while(n) {
+				c = n*10;
+				if(c>>b)
+					break;
+				n = c;
+			}
+			goto xx;
+		}
+		n = n*10 + c-'0';
+		(*na)--;
+	}
+	for(;;) {
+		c = n>>b;
+		n -= c<<b;
+		*p++ = c + '0';
+		c = *a++;
+		if(c == 0)
+			break;
+		n = n*10 + c-'0';
+	}
+	(*na)++;
+xx:
+	while(n) {
+		n = n*10;
+		c = n>>b;
+		n -= c<<b;
+		*p++ = c + '0';
+		(*na)++;
+	}
+	*p = 0;
+}
+
+static	Tab	tab1[] =
+{
+	{ 1, 0, "" },
+	{ 3, 1, "7" },
+	{ 6, 2, "63" },
+	{ 9, 3, "511" },
+	{ 13, 4, "8191" },
+	{ 16, 5, "65535" },
+	{ 19, 6, "524287" },
+	{ 23, 7, "8388607" },
+	{ 26, 8, "67108863" },
+	{ 27, 9, "134217727" },
+};
+
+static void
+divascii(char *a, int *na, int *dp, int *bp)
+{
+	int b, d;
+	Tab *t;
+
+	d = *dp;
+	if(d >= (int)(nelem(tab1)))
+		d = (int)(nelem(tab1))-1;
+	t = tab1 + d;
+	b = t->bp;
+	if(memcmp(a, t->cmp, t->siz) > 0)
+		d--;
+	*dp -= d;
+	*bp += b;
+	divby(a, na, b);
+}
+
+static void
+mulby(char *a, char *p, char *q, int b)
+{
+	int n, c;
+
+	n = 0;
+	*p = 0;
+	for(;;) {
+		q--;
+		if(q < a)
+			break;
+		c = *q - '0';
+		c = (c<<b) + n;
+		n = c/10;
+		c -= n*10;
+		p--;
+		*p = c + '0';
+	}
+	while(n) {
+		c = n;
+		n = c/10;
+		c -= n*10;
+		p--;
+		*p = c + '0';
+	}
+}
+
+static	Tab	tab2[] =
+{
+	{ 1, 1, "" },				/* dp = 0-0 */
+	{ 3, 3, "125" },
+	{ 6, 5, "15625" },
+	{ 9, 7, "1953125" },
+	{ 13, 10, "1220703125" },
+	{ 16, 12, "152587890625" },
+	{ 19, 14, "19073486328125" },
+	{ 23, 17, "11920928955078125" },
+	{ 26, 19, "1490116119384765625" },
+	{ 27, 19, "7450580596923828125" },		/* dp 8-9 */
+};
+
+static void
+mulascii(char *a, int *na, int *dp, int *bp)
+{
+	char *p;
+	int d, b;
+	Tab *t;
+
+	d = -*dp;
+	if(d >= (int)(nelem(tab2)))
+		d = (int)(nelem(tab2))-1;
+	t = tab2 + d;
+	b = t->bp;
+	if(memcmp(a, t->cmp, t->siz) < 0)
+		d--;
+	p = a + *na;
+	*bp -= b;
+	*dp += d;
+	*na += d;
+	mulby(a, p+d, p, b);
+}
+
+static int
+xcmp(char *a, char *b)
+{
+	int c1, c2;
+
+	while((c1 = *b++) != '\0') {
+		c2 = *a++;
+		if(c2 >= 'A' && c2 <= 'Z')
+			c2 = c2 - 'A' + 'a';
+		if(c1 != c2)
+			return 1;
+	}
+	return 0;
+}