[geos-commits] r2358 - in trunk/source: algorithm
headers/geos/algorithm
svn_geos at osgeo.org
svn_geos at osgeo.org
Tue Apr 14 09:23:08 EDT 2009
Author: strk
Date: 2009-04-14 09:23:07 -0400 (Tue, 14 Apr 2009)
New Revision: 2358
Modified:
trunk/source/algorithm/CGAlgorithms.cpp
trunk/source/headers/geos/algorithm/CGAlgorithms.h
Log:
Sync CGAlgorithms with JTS-1.9
Modified: trunk/source/algorithm/CGAlgorithms.cpp
===================================================================
--- trunk/source/algorithm/CGAlgorithms.cpp 2009-04-14 12:55:39 UTC (rev 2357)
+++ trunk/source/algorithm/CGAlgorithms.cpp 2009-04-14 13:23:07 UTC (rev 2358)
@@ -14,15 +14,18 @@
*
**********************************************************************
*
- * Last port: algorithm/CGAlgorithms.java rev. 1.34 (JTS-1.7.1)
+ * Last port: algorithm/CGAlgorithms.java rev. 1.46 (JTS-1.9)
*
**********************************************************************/
#include <geos/algorithm/CGAlgorithms.h>
#include <geos/algorithm/RobustDeterminant.h>
#include <geos/algorithm/LineIntersector.h>
+#include <geos/algorithm/RayCrossingCounter.h>
#include <geos/geom/CoordinateSequence.h>
#include <geos/geom/Coordinate.h>
+#include <geos/geom/Location.h>
+#include <geos/util/IllegalArgumentException.h>
#include <algorithm>
//#include <cstdio>
@@ -53,48 +56,7 @@
bool
CGAlgorithms::isPointInRing(const Coordinate& p, const CoordinateSequence* ring)
{
- double xInt; // x intersection of segment with ray
- int crossings = 0; // number of segment/ray crossings
- double x1; // translated coordinates
- double y1;
- double x2;
- double y2;
-
- /*
- * For each segment l = (i-1, i), see if it crosses ray from
- * test point in positive x direction.
- */
- size_t nPts=ring->getSize();
- for(size_t i=1; i<nPts; i++)
- {
- const Coordinate &p1=ring->getAt(i);
- const Coordinate &p2=ring->getAt(i-1);
- x1 = p1.x - p.x;
- y1 = p1.y - p.y;
- x2 = p2.x - p.x;
- y2 = p2.y - p.y;
-
- if (((y1 > 0) && (y2 <= 0)) ||
- ((y2 > 0) && (y1 <= 0)))
- {
- /*
- * segment straddles x axis, so compute intersection.
- */
- xInt = RobustDeterminant::signOfDet2x2(x1, y1, x2, y2)
- / (y2 - y1);
-
- /*
- * crosses ray if strictly positive intersection.
- */
- if (0.0 < xInt) crossings++;
- }
- }
-
- /*
- * p is inside if number of crossings is odd.
- */
- if ((crossings % 2) == 1) return true;
- return false;
+ return locatePointInRing(p, *ring) != Location::EXTERIOR;
}
/*public static*/
@@ -102,47 +64,23 @@
CGAlgorithms::isPointInRing(const Coordinate& p,
const Coordinate::ConstVect& ring)
{
- double xInt; // x intersection of segment with ray
- int crossings = 0; // number of segment/ray crossings
- double x1; // translated coordinates
- double y1;
- double x2;
- double y2;
+ return locatePointInRing(p, ring) != Location::EXTERIOR;
+}
- /*
- * For each segment l = (i-1, i), see if it crosses ray from
- * test point in positive x direction.
- */
- for(size_t i=1, nPts=ring.size(); i<nPts; ++i)
- {
- const Coordinate *p1=ring[i];
- const Coordinate *p2=ring[i-1];
- x1 = p1->x - p.x;
- y1 = p1->y - p.y;
- x2 = p2->x - p.x;
- y2 = p2->y - p.y;
+/*public static*/
+int
+CGAlgorithms::locatePointInRing(const Coordinate& p,
+ const CoordinateSequence& ring)
+{
+ return RayCrossingCounter::locatePointInRing(p, ring);
+}
- if (((y1 > 0) && (y2 <= 0)) ||
- ((y2 > 0) && (y1 <= 0)))
- {
- /*
- * segment straddles x axis, so compute intersection.
- */
- xInt = RobustDeterminant::signOfDet2x2(x1, y1, x2, y2)
- / (y2 - y1);
-
- /*
- * crosses ray if strictly positive intersection.
- */
- if (0.0 < xInt) crossings++;
- }
- }
-
- /*
- * p is inside if number of crossings is odd.
- */
- if ((crossings % 2) == 1) return true;
- return false;
+/*public static*/
+int
+CGAlgorithms::locatePointInRing(const Coordinate& p,
+ const std::vector<const geom::Coordinate*>& ring)
+{
+ return RayCrossingCounter::locatePointInRing(p, ring);
}
/*public static*/
@@ -169,12 +107,18 @@
CGAlgorithms::isCCW(const CoordinateSequence* ring)
{
// # of points without closing endpoint
- const std::size_t nPts=ring->getSize()-1;
+ const std::size_t nPts=ring->getSize()-1;
+ // sanity check
+ if (nPts < 3)
+ {
+ throw util::IllegalArgumentException("Ring has fewer than 3 points, so orientation cannot be determined");
+ }
+
// find highest point
const Coordinate *hiPt=&ring->getAt(0);
int hiIndex=0;
- for (std::size_t i=1; i<=nPts; ++i)
+ for (std::size_t i=1; i<=nPts; ++i)
{
const Coordinate *p=&ring->getAt(i);
if (p->y > hiPt->y) {
@@ -379,17 +323,32 @@
double
CGAlgorithms::length(const CoordinateSequence* pts)
{
+ // optimized for processing CoordinateSequences
+
size_t npts=pts->getSize();
+ if (npts <= 1) return 0.0;
- if (npts<1) return 0.0;
+ double len = 0.0;
- double sum=0.0;
+ const Coordinate& p = pts->getAt(0);
+ double x0 = p.x;
+ double y0 = p.y;
- for(size_t i=1; i<npts; ++i)
+ for(size_t i = 1; i < npts; ++i)
{
- sum+=pts->getAt(i).distance(pts->getAt(i - 1));
+ const Coordinate& p = pts->getAt(i);
+ double x1 = p.x;
+ double y1 = p.y;
+ double dx = x1 - x0;
+ double dy = y1 - y0;
+
+ len += sqrt(dx * dx + dy * dy);
+
+ x0 = x1;
+ y0 = y1;
}
- return sum;
+
+ return len;
}
Modified: trunk/source/headers/geos/algorithm/CGAlgorithms.h
===================================================================
--- trunk/source/headers/geos/algorithm/CGAlgorithms.h 2009-04-14 12:55:39 UTC (rev 2357)
+++ trunk/source/headers/geos/algorithm/CGAlgorithms.h 2009-04-14 13:23:07 UTC (rev 2358)
@@ -14,7 +14,7 @@
*
**********************************************************************
*
- * Last port: algorithm/CGAlgorithms.java rev. 1.34 (JTS-1.7.1)
+ * Last port: algorithm/CGAlgorithms.java rev. 1.46 (JTS-1.9)
*
**********************************************************************/
@@ -62,18 +62,20 @@
CGAlgorithms(){};
/** \brief
- * Test whether a point lies inside a ring.
+ * Tests whether a point lies inside a ring.
*
* The ring may be oriented in either direction.
- * If the point lies on the ring boundary the result
- * of this method is unspecified.
+ * A point lying exactly on the ring boundary is considered
+ * to be inside the ring.
*
- * This algorithm does not attempt to first check the
+ * This algorithm does not first check the
* point against the envelope of the ring.
*
* @param p point to check for ring inclusion
- * @param ring assumed to have first point identical to last point
+ * @param ring is assumed to have first point identical to last point
* @return <code>true</code> if p is inside ring
+ *
+ * @see locatePointInRing
*/
static bool isPointInRing(const geom::Coordinate& p,
const geom::CoordinateSequence* ring);
@@ -83,8 +85,29 @@
const std::vector<const geom::Coordinate*>& ring);
/** \brief
- * Test whether a point lies on a linestring.
+ * Determines whether a point lies in the interior,
+ * on the boundary, or in the exterior of a ring.
*
+ * The ring may be oriented in either direction.
+ *
+ * This method does <i>not</i> first check the point against
+ * the envelope of the ring.
+ *
+ * @param p point to check for ring inclusion
+ * @param ring an array of coordinates representing the ring
+ * (which must have first point identical to last point)
+ * @return the {@link Location} of p relative to the ring
+ */
+ static int locatePointInRing(const geom::Coordinate& p,
+ const geom::CoordinateSequence& ring);
+
+ /// Same as above, but taking a vector of const Coordinates
+ static int locatePointInRing(const geom::Coordinate& p,
+ const std::vector<const geom::Coordinate*>& ring);
+
+ /** \brief
+ * Test whether a point lies on the given line segment
+ *
* @return true true if
* the point is a vertex of the line or lies in the interior of a line
* segment in the linestring
More information about the geos-commits
mailing list