[geos-commits] r3788 - in branches/3.3: . include/geos/algorithm src/algorithm tests/unit/algorithm

svn_geos at osgeo.org svn_geos at osgeo.org
Fri Mar 8 10:01:29 PST 2013


Author: strk
Date: 2013-03-08 10:01:29 -0800 (Fri, 08 Mar 2013)
New Revision: 3788

Modified:
   branches/3.3/NEWS
   branches/3.3/include/geos/algorithm/LineIntersector.h
   branches/3.3/src/algorithm/LineIntersector.cpp
   branches/3.3/tests/unit/algorithm/RobustLineIntersectionTest.cpp
Log:
Fix RobustLineIntersector handling of invalid intersection points (#622)

Adds new testcases. Many of them fail, probably due to the lack of
double double use, but one of them only fails with the old
RobustLineIntersector heuristic handling invalid intersection points

Modified: branches/3.3/NEWS
===================================================================
--- branches/3.3/NEWS	2013-03-08 17:51:21 UTC (rev 3787)
+++ branches/3.3/NEWS	2013-03-08 18:01:29 UTC (rev 3788)
@@ -4,6 +4,7 @@
 - Bug fixes / improvements
     - Fix area boundary return from GEOSPointOnSurface (#623)
     - Speedup GEOSWKBReader_read (#621)
+    - Fix RobustLineIntersector handling of invalid intersection points (#622)
 
 Changes in 3.3.8
 2013-02-28

Modified: branches/3.3/include/geos/algorithm/LineIntersector.h
===================================================================
--- branches/3.3/include/geos/algorithm/LineIntersector.h	2013-03-08 17:51:21 UTC (rev 3787)
+++ branches/3.3/include/geos/algorithm/LineIntersector.h	2013-03-08 18:01:29 UTC (rev 3788)
@@ -14,7 +14,7 @@
  *
  **********************************************************************
  *
- * Last port: algorithm/RobustLineIntersector.java rev. 1.38 (JTS-1.10)
+ * Last port: algorithm/RobustLineIntersector.java r785 (JTS-1.13+)
  *
  **********************************************************************/
 

Modified: branches/3.3/src/algorithm/LineIntersector.cpp
===================================================================
--- branches/3.3/src/algorithm/LineIntersector.cpp	2013-03-08 17:51:21 UTC (rev 3787)
+++ branches/3.3/src/algorithm/LineIntersector.cpp	2013-03-08 18:01:29 UTC (rev 3788)
@@ -14,7 +14,7 @@
  *
  **********************************************************************
  *
- * Last port: algorithm/RobustLineIntersector.java rev. 1.38 (JTS-1.10)
+ * Last port: algorithm/RobustLineIntersector.java r785 (JTS-1.13+)
  *
  **********************************************************************/
 
@@ -22,7 +22,7 @@
 #include <geos/algorithm/CGAlgorithms.h>
 #include <geos/algorithm/HCoordinate.h>
 #include <geos/algorithm/NotRepresentableException.h>
-#include <geos/algorithm/CentralEndpointIntersector.h>
+//#include <geos/algorithm/CentralEndpointIntersector.h>
 #include <geos/geom/Coordinate.h>
 #include <geos/geom/PrecisionModel.h>
 #include <geos/geom/Envelope.h>
@@ -54,6 +54,54 @@
 namespace geos {
 namespace algorithm { // geos.algorithm
 
+namespace { // anonymous
+
+/**
+ * Finds the endpoint of the segments P and Q which
+ * is closest to the other segment.
+ * This is a reasonable surrogate for the true
+ * intersection points in ill-conditioned cases
+ * (e.g. where two segments are nearly coincident,
+ * or where the endpoint of one segment lies almost on the other segment).
+ * <p>
+ * This replaces the older CentralEndpoint heuristic,
+ * which chose the wrong endpoint in some cases
+ * where the segments had very distinct slopes
+ * and one endpoint lay almost on the other segment.
+ *
+ * @param p1 an endpoint of segment P
+ * @param p2 an endpoint of segment P
+ * @param q1 an endpoint of segment Q
+ * @param q2 an endpoint of segment Q
+ * @return the nearest endpoint to the other segment
+ */
+Coordinate nearestEndpoint(const Coordinate& p1, const Coordinate& p2,
+    const Coordinate& q1, const Coordinate& q2)
+{
+  Coordinate nearestPt = p1;
+  double minDist = CGAlgorithms::distancePointLine(p1, q1, q2);
+
+  double dist = CGAlgorithms::distancePointLine(p2, q1, q2);
+  if (dist < minDist) {
+    minDist = dist;
+    nearestPt = p2;
+  }
+  dist = CGAlgorithms::distancePointLine(q1, p1, p2);
+  if (dist < minDist) {
+    minDist = dist;
+    nearestPt = q1;
+  }
+  dist = CGAlgorithms::distancePointLine(q2, p1, p2);
+  if (dist < minDist) {
+    minDist = dist;
+    nearestPt = q2;
+  }
+  return nearestPt;
+}
+
+
+} // anonymous namespace
+
 /*public static*/
 double
 LineIntersector::computeEdgeDistance(const Coordinate& p,const Coordinate& p0,const Coordinate& p1)
@@ -716,7 +764,8 @@
 
 	if (! isInSegmentEnvelopes(intPt))
 	{
-		intPt = CentralEndpointIntersector::getIntersection(p1, p2, q1, q2);
+		//intPt = CentralEndpointIntersector::getIntersection(p1, p2, q1, q2);
+		intPt = nearestEndpoint(p1, p2, q1, q2);
 #if GEOS_DEBUG
 		cerr << "Intersection outside segment envelopes, snapped to "
 		     << intPt.toString() << endl;
@@ -850,34 +899,11 @@
 
 	} catch (const NotRepresentableException& /* e */) {
 		// compute an approximate result
-		intPt = CentralEndpointIntersector::getIntersection(p1, p2, q1, q2);
+		//intPt = CentralEndpointIntersector::getIntersection(p1, p2, q1, q2);
+		intPt = nearestEndpoint(p1, p2, q1, q2);
     	}
 }
 
 } // namespace geos.algorithm
 } // namespace geos
 
-/**********************************************************************
- * $Log$
- * Revision 1.42  2006/04/12 11:57:11  strk
- * Fixed debugging line
- *
- * Revision 1.41  2006/04/07 09:12:57  strk
- * kept isInSegmentEnvelopes() check even when not debugging
- *
- * Revision 1.40  2006/04/06 21:31:40  strk
- * Const correctness for debugging function
- *
- * Revision 1.39  2006/04/04 12:39:08  strk
- * Changed NotRepresentableCoordinate exception handler to throw
- * a TopologyException. This allows further handling.
- *
- * Revision 1.38  2006/03/21 11:12:23  strk
- * Cleanups: headers inclusion and Log section
- *
- * Revision 1.37  2006/03/16 10:38:14  strk
- * Bug #63 - Remove unreferenced local variable warning
- *
- * Revision 1.36  2006/03/09 16:46:45  strk
- * geos::geom namespace definition, first pass at headers split
- **********************************************************************/

Modified: branches/3.3/tests/unit/algorithm/RobustLineIntersectionTest.cpp
===================================================================
--- branches/3.3/tests/unit/algorithm/RobustLineIntersectionTest.cpp	2013-03-08 17:51:21 UTC (rev 3787)
+++ branches/3.3/tests/unit/algorithm/RobustLineIntersectionTest.cpp	2013-03-08 18:01:29 UTC (rev 3788)
@@ -1,6 +1,6 @@
 // $Id$
 // 
-// Ported from JTS junit/algorithm/RobustLineIntersectionTest.java rev. 1.3
+// Ported from JTS junit/algorithm/RobustLineIntersectionTest.java r785
 
 #include <tut.hpp>
 // geos
@@ -44,7 +44,7 @@
 	                   double distanceTolerance)
         {
                 bool isEqual = equals(p, q, distanceTolerance);
-		ensure("testIntPoints", isEqual);
+		ensure("testIntPoints: expected: " + p.toString() + " obtained " + q.toString(), isEqual);
         }
 
 	/**
@@ -176,6 +176,32 @@
 		                    intPt, distanceTolerance);
         }
 
+	void computeIntersectionNone(const std::string& wkt1,
+	                         const std::string& wkt2)
+                // throws ParseException
+        {
+		GeomPtr g1(reader.read(wkt1));
+		GeomPtr g2(reader.read(wkt2));
+
+		LineString* l1ptr = dynamic_cast<LineString*>(g1.get());
+		LineString* l2ptr = dynamic_cast<LineString*>(g2.get());
+
+		ensure(0 != l1ptr);
+		ensure(0 != l2ptr);
+
+		LineString& l1 = *l1ptr;
+		LineString& l2 = *l2ptr;
+
+	        std::vector<Coordinate> pt;
+		pt.push_back(l1.getCoordinateN(0));
+		pt.push_back(l1.getCoordinateN(1));
+		pt.push_back(l2.getCoordinateN(0));
+		pt.push_back(l2.getCoordinateN(1));
+
+	  std::vector<Coordinate> intPt;
+		computeIntersection(pt, 0, intPt, 0);
+        }
+
 	test_robustlineintersection_data()
 		:
 		pm(),
@@ -287,7 +313,99 @@
 	}
 #endif // fails
 
+  /**
+   * Test involving two non-almost-parallel lines.
+   * Does not seem to cause problems with basic line intersection algorithm.
+   *
+   */
+	//     (testLeduc_1)
+	template<>
+	template<>
+	void object::test<6>()
+	{         
+    computeIntersection(
+        "LINESTRING (305690.0434123494 254176.46578338774, 305601.9999843455 254243.19999846347)",
+        "LINESTRING (305689.6153764265 254177.33102743194, 305692.4999844298 254171.4999983967)",
+        1,
+        "POINT (305690.0434123494 254176.46578338774)",
+        0);
+	}
 
+#if 0 // fails: finds an intersection (we don't have DD)
+  /**
+   * Test from Tomas Fa - JTS list 6/13/2012
+   *
+   * Fails using original JTS DeVillers determine orientation test.
+   * Succeeds using DD and Shewchuk orientation
+   *
+   */
+  // testTomasFa_2
+	template<>
+	template<>
+	void object::test<7>()
+	{         
+    computeIntersectionNone(
+        "LINESTRING (-5.9 163.1, 76.1 250.7)",
+        "LINESTRING (14.6 185.0, 96.6 272.6)");
+	}
+#endif // fails
 
+#if 0 // fails: finds an intersection (we don't have DD)
+  /**
+   * Test from Tomas Fa - JTS list 6/13/2012
+   *
+   * Fails using original JTS DeVillers determine orientation test.
+   * Succeeds using DD and Shewchuk orientation
+   *
+   */
+  // testTomasFa_1
+	template<>
+	template<>
+	void object::test<8>()
+	{         
+    computeIntersectionNone(
+        "LINESTRING (-42.0 163.2, 21.2 265.2)",
+        "LINESTRING (-26.2 188.7, 37.0 290.7)");
+	}
+#endif // fails
+
+  /**
+   * Following cases were failures when using the CentralEndpointIntersector heuristic.
+   * This is because one segment lies at a significant angle to the other,
+   * with only one endpoint is close to the other segment.
+   * The CE heuristic chose the wrong endpoint to return.
+   * The fix is to use a new heuristic which out of the 4 endpoints
+   * chooses the one which is closest to the other segment.
+   * This works in all known failure cases.
+   *
+   */
+  // public void testCentralEndpointHeuristicFailure()
+  template<>
+	template<>
+	void object::test<9>()
+  {
+    computeIntersection(
+        "LINESTRING (163.81867067 -211.31840378, 165.9174252 -214.1665075)",
+        "LINESTRING (2.84139601 -57.95412726, 469.59990601 -502.63851732)",
+        1,
+        "POINT (163.81867067 -211.31840378)",
+        0);
+  }
+
+  // public void testCentralEndpointHeuristicFailure2()
+  template<>
+	template<>
+	void object::test<10>()
+  {
+    computeIntersection(
+        "LINESTRING (-58.00593335955 -1.43739086465, -513.86101637525 -457.29247388035)",
+        "LINESTRING (-215.22279674875 -158.65425425385, -218.1208801283 -160.68343590235)",
+        1,
+        "POINT ( -215.22279674875 -158.65425425385 )",
+        0);
+  }
+
+
+
 } // namespace tut
 



More information about the geos-commits mailing list