[postgis-tickets] r14595 - Avoid any drift of cutter point on lines split

[Sandro Santilli]

Author: strk
Date: 2016-01-13 03:35:15 -0800 (Wed, 13 Jan 2016)
New Revision: 14595

Modified:
   branches/2.2/NEWS
   branches/2.2/liblwgeom/cunit/cu_tester.h
   branches/2.2/liblwgeom/lwgeom_geos_split.c
Log:
Avoid any drift of cutter point on lines split

Should fix splitting operations on at least arm64, ppc64el and s390x.
See #3422 and #3401 (for 2.2 branch)

Also turn ASSERT_DOUBLE_EQUAL back to a tolerance-free check
(better use a different name for tolerance-aware check,
so caller can decide)

Modified: branches/2.2/NEWS
===================================================================
--- branches/2.2/NEWS	2016-01-13 11:09:24 UTC (rev 14594)
+++ branches/2.2/NEWS	2016-01-13 11:35:15 UTC (rev 14595)
@@ -3,7 +3,8 @@
 
  * Bug Fixes *
 
-  -
+  - #3422, Improve ST_Split robustness on standard precision double
+           systems (arm64, ppc64el, s390c, powerpc, ...)
 
 PostGIS 2.2.1
 2016/01/06

Modified: branches/2.2/liblwgeom/cunit/cu_tester.h
===================================================================
--- branches/2.2/liblwgeom/cunit/cu_tester.h	2016-01-13 11:09:24 UTC (rev 14594)
+++ branches/2.2/liblwgeom/cunit/cu_tester.h	2016-01-13 11:35:15 UTC (rev 14595)
@@ -23,11 +23,10 @@
 /* Our internal callback to register Suites with the main tester */
 typedef void (*PG_SuiteSetup)(void);
 
-#define ASSERT_DOUBLE_EQUAL_TOLERANCE 10e-8
 #define ASSERT_DOUBLE_EQUAL(o,e) do { \
-  if ( fabs((double)o-(double)e) > ASSERT_DOUBLE_EQUAL_TOLERANCE ) \
+  if ( o != e ) \
     fprintf(stderr, "[%s:%d]\n Expected: %g\n Obtained: %g\n", __FILE__, __LINE__, (double)(e), (o)); \
-  CU_ASSERT_DOUBLE_EQUAL((double)o,(double)e,ASSERT_DOUBLE_EQUAL_TOLERANCE); \
+  CU_ASSERT_EQUAL(o,(double)e); \
 } while (0);
 
 #define ASSERT_INT_EQUAL(o,e) do { \

Modified: branches/2.2/liblwgeom/lwgeom_geos_split.c
===================================================================
--- branches/2.2/liblwgeom/lwgeom_geos_split.c	2016-01-13 11:09:24 UTC (rev 14594)
+++ branches/2.2/liblwgeom/lwgeom_geos_split.c	2016-01-13 11:35:15 UTC (rev 14595)
@@ -258,6 +258,7 @@
 		{
 			mindist = dist;
 			seg=i;
+			if ( mindist == 0.0 ) break; /* can't be closer than ON line */
 		}
 		p1 = p2;
 	}
@@ -273,11 +274,18 @@
 
 	/*
 	 * We need to project the
-	 * point on the closest segment.
+	 * point on the closest segment,
+	 * to interpolate Z and M if needed
 	 */
 	getPoint4d_p(ipa, seg, &p1);
 	getPoint4d_p(ipa, seg+1, &p2);
 	closest_point_on_segment(&pt, &p1, &p2, &pt_projected);
+	/* But X and Y we want the ones of the input point,
+	 * as on some architectures the interpolation math moves the
+	 * coordinates (see #3422)
+	 */
+	pt_projected.x = pt.x;
+	pt_projected.y = pt.y;
 
 	LWDEBUGF(3, "Projected point:(%g %g), seg:%d, p1:(%g %g), p2:(%g %g)", pt_projected.x, pt_projected.y, seg, p1.x, p1.y, p2.x, p2.y);