[GRASS-SVN] r57839 - grass/trunk/lib/gis

[svn_grass at osgeo.org]

Author: mmetz
Date: 2013-09-25 06:15:32 -0700 (Wed, 25 Sep 2013)
New Revision: 57839

Modified:
   grass/trunk/lib/gis/plot.c
Log:
libgis: fix FPE tolerance for plotting

Modified: grass/trunk/lib/gis/plot.c
===================================================================
--- grass/trunk/lib/gis/plot.c	2013-09-25 10:39:29 UTC (rev 57838)
+++ grass/trunk/lib/gis/plot.c	2013-09-25 13:15:32 UTC (rev 57839)
@@ -449,16 +449,20 @@
     }
 
     /* check if perimeter has odd number of points */
-    if (st->np % 2)
+    if (st->np & 1) {
+	G_debug(1, "Weird internal error: perimeter has odd number of points");
 	return OUT_OF_SYNC;
+    }
 
     /* sort the edge points by col(x) and then by row(y) */
     qsort(st->P, st->np, sizeof(POINT), edge_order);
 
     /* plot */
     for (i = 1; i < st->np; i += 2) {
-	if (st->P[i].y != st->P[i - 1].y)
+	if (st->P[i].y != st->P[i - 1].y) {
+	    G_debug(1, "Weird internal error: edge leaves row");
 	    return OUT_OF_SYNC;
+	}
 	st->row_fill(st->P[i].y, st->P[i - 1].x + shift1, st->P[i].x + shift1);
     }
     if (st->window.proj == PROJECTION_LL) {	/* now do wrap-around, part 2 */
@@ -581,8 +585,10 @@
     }				/* for() */
 
     /* check if perimeter has odd number of points */
-    if (st->np % 2)
+    if (st->np & 1) {
+	G_debug(1, "Weird internal error: perimeter has odd number of points");
 	return OUT_OF_SYNC;
+    }
 
     /* sort the edge points by col(x) and then by row(y) */
     qsort(st->P, st->np, sizeof(POINT), &edge_order);
@@ -590,8 +596,10 @@
     /* plot */
     for (j = 0; j < rings; j++) {
 	for (i = 1; i < st->np; i += 2) {
-	    if (st->P[i].y != st->P[i - 1].y)
+	    if (st->P[i].y != st->P[i - 1].y) {
+		G_debug(1, "Weird internal error: edge leaves row");
 		return OUT_OF_SYNC;
+	    }
 	    st->row_fill(st->P[i].y, st->P[i - 1].x + shift1[j], st->P[i].x + shift1[j]);
 	}
 	if (st->window.proj == PROJECTION_LL) {	/* now do wrap-around, part 2 */
@@ -619,19 +627,27 @@
 
 static int edge(double x0, double y0, double x1, double y1)
 {
-    double m;
-    double dx, dy, x;
+    double m, d;
+    double x;
     int ystart, ystop;
+    int exp;
 
+    /* tolerance to avoid FPE */
+    d = GRASS_EPSILON;
+    if (y0 != y1) {
+	if (fabs(y0) > fabs(y1))
+	    d = fabs(y0);
+	else
+	    d = fabs(y1);
 
-    /* tolerance to avoid FPE */
-    dx = x0 - x1;
-    dy = y0 - y1;
-    if (fabs(dy) < 1e-10)
+	d = frexp(d, &exp);
+	exp -= 52;
+	d = ldexp(d, exp);
+    }
+
+    if (fabs(y0 - y1) < d)
 	return 1;
 
-    m = (x0 - x1) / dy;
-
     if (y0 < y1) {
 	ystart = iceil(y0);
 	ystop = ifloor(y1);
@@ -648,12 +664,14 @@
     if (ystart > ystop)
 	return 1;		/* does not cross center line of row */
 
+    m = (x0 - x1) / (y0 - y1);
     x = m * (ystart - y0) + x0;
     while (ystart <= ystop) {
 	if (!edge_point(x, ystart++))
 	    return 0;
 	x += m;
     }
+
     return 1;
 }