[geos-commits] [SCM] GEOS branch main updated. 2f9e7531804e3a4f09823d8f1b9536cb04d313bf
git at osgeo.org
git at osgeo.org
Tue Jul 18 11:47:38 PDT 2023
This is an automated email from the git hooks/post-receive script. It was
generated because a ref change was pushed to the repository containing
the project "GEOS".
The branch, main has been updated
via 2f9e7531804e3a4f09823d8f1b9536cb04d313bf (commit)
from 356598a9de13d6964773546611e944ff851f8fab (commit)
Those revisions listed above that are new to this repository have
not appeared on any other notification email; so we list those
revisions in full, below.
- Log -----------------------------------------------------------------
commit 2f9e7531804e3a4f09823d8f1b9536cb04d313bf
Author: Martin Davis <mtnclimb at gmail.com>
Date: Tue Jul 18 11:46:19 2023 -0700
Add LargestEmptyCircle handling of polygon obstacles (#939)
diff --git a/capi/geos_c.h.in b/capi/geos_c.h.in
index 781ff7c13..3466cef60 100644
--- a/capi/geos_c.h.in
+++ b/capi/geos_c.h.in
@@ -4129,8 +4129,7 @@ extern GEOSGeometry GEOS_DLL *GEOSMaximumInscribedCircle(
* Constructs the "largest empty circle" (LEC) for a set of obstacle geometries
* and within a polygonal boundary,
* with accuracy to to a specified distance tolerance.
-* The obstacles are point and line geometries.
-* Polygonal obstacles are treated as linear features.
+* The obstacles may be any collection of points, lines and polygons.
* The LEC is the largest circle whose interior does not intersect with any obstacle.
* and which has its **center** inside the given boundary.
* If no boundary is provided, the
@@ -4139,15 +4138,21 @@ extern GEOSGeometry GEOS_DLL *GEOSMaximumInscribedCircle(
* the obstacles (up to the given distance tolerance).
* The LEC is determined by the center point and a point indicating the circle radius
* (which will lie on an obstacle).
+*
+* To compute an LEC which lies **wholly** within a polygonal boundary,
+* include the boundary of the polygon(s) as a linear obstacle.
+*
* The implementation uses a successive-approximation technique over a grid of square cells covering the obstacles and boundary.
* The grid is refined using a branch-and-bound algorithm. Point containment and distance are computed in a performant
* way by using spatial indexes.
+*
* Returns the LEC radius as a two-point linestring, with the start point at the center of the inscribed circle and the end
* on the boundary of the circle.
+*
* \param obstacles The geometries that the LEC must not cross
* \param boundary The area to contain the LEC center (may be null or empty)
* \param tolerance Stop the algorithm when the search area is smaller than this tolerance
-* \return A newly allocated geometry of the LEC radius. NULL on exception.
+* \return A newly allocated geometry of the LEC radius line. NULL on exception.
* Caller is responsible for freeing with GEOSGeom_destroy().
* \see geos::algorithm::construct::LargestEmptyCircle
*
diff --git a/include/geos/algorithm/construct/IndexedDistanceToPoint.h b/include/geos/algorithm/construct/IndexedDistanceToPoint.h
new file mode 100644
index 000000000..68d845245
--- /dev/null
+++ b/include/geos/algorithm/construct/IndexedDistanceToPoint.h
@@ -0,0 +1,88 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2023 Martin Davis <mtnclimb at gmail.com>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************
+ *
+ * Last port: algorithm/construct/IndexedDistanceToPoint.java
+ * https://github.com/locationtech/jts/commit/d92f783163d9440fcc10c729143787bf7b9fe8f9
+ *
+ **********************************************************************/
+
+#pragma once
+
+#include <geos/geom/CoordinateSequence.h>
+#include <geos/geom/Geometry.h>
+#include <geos/geom/Point.h>
+#include <geos/algorithm/construct/IndexedPointInPolygonsLocator.h>
+#include <geos/operation/distance/IndexedFacetDistance.h>
+
+using geos::geom::Geometry;
+using geos::geom::Point;
+using geos::operation::distance::IndexedFacetDistance;
+
+namespace geos {
+namespace algorithm { // geos::algorithm
+namespace construct { // geos::algorithm::construct
+
+/**
+ * \brief Computes the distance between a point and a geometry
+ * (which may be a collection containing any type of geometry).
+ *
+ * Also computes the pair of points containing the input
+ * point and the nearest point on the geometry.
+ *
+ * \author Martin Davis
+ */
+class GEOS_DLL IndexedDistanceToPoint {
+
+public:
+ /**
+ * \brief Creates an instance to find the distance from points to a geometry.
+ *
+ * \param geom the geometry to compute distances to
+ */
+ IndexedDistanceToPoint(const Geometry& geom);
+
+ /**
+ * \brief Computes the distance from the base geometry to the given point.
+ *
+ * \param pt the point to compute the distance to
+ *
+ * \return the computed distance
+ */
+ double distance(const Point& pt);
+
+ /**
+ * \brief Computes the nearest point on the geometry to the point.
+ *
+ * The first location lies on the geometry,
+ * and the second location is the provided point.
+ *
+ * \param pt the point to compute the nearest point for
+ *
+ * \return the points that are nearest
+ */
+ std::unique_ptr<geom::CoordinateSequence> nearestPoints(const Point& pt);
+
+private:
+ void init();
+
+ bool isInArea(const Point& pt);
+
+ //-- members
+ const Geometry& targetGeometry;
+ std::unique_ptr<operation::distance::IndexedFacetDistance> facetDistance;
+ std::unique_ptr<IndexedPointInPolygonsLocator> ptLocator;
+
+};
+
+}}}
\ No newline at end of file
diff --git a/include/geos/algorithm/construct/IndexedPointInPolygonsLocator.h b/include/geos/algorithm/construct/IndexedPointInPolygonsLocator.h
new file mode 100644
index 000000000..dd82d55a4
--- /dev/null
+++ b/include/geos/algorithm/construct/IndexedPointInPolygonsLocator.h
@@ -0,0 +1,79 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2023 Martin Davis <mtnclimb at gmail.com>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************
+ *
+ * Last port: algorithm/construct/IndexedDistanceToPoint.java
+ * https://github.com/locationtech/jts/commit/d92f783163d9440fcc10c729143787bf7b9fe8f9
+ *
+ **********************************************************************/
+
+#pragma once
+
+#include <geos/geom/Coordinate.h>
+#include <geos/geom/Geometry.h>
+#include <geos/geom/Location.h>
+#include <geos/index/strtree/TemplateSTRtree.h>
+#include <geos/algorithm/locate/IndexedPointInAreaLocator.h>
+
+using geos::geom::Geometry;
+using geos::geom::CoordinateXY;
+using geos::geom::Location;
+using geos::index::strtree::TemplateSTRtree;
+using geos::algorithm::locate::IndexedPointInAreaLocator;
+
+namespace geos {
+namespace algorithm { // geos::algorithm
+namespace construct { // geos::algorithm::construct
+
+/**
+ * \brief Determines the location of a point in the polygonal elements of a geometry.
+ *
+ * Uses spatial indexing to provide efficient performance.
+ *
+ * \author Martin Davis
+ */
+class GEOS_DLL IndexedPointInPolygonsLocator {
+
+public:
+ /**
+ * \brief Creates an instance to locate a point in polygonal elements.
+ *
+ * \param geom the geometry to locate in
+ */
+ IndexedPointInPolygonsLocator(const Geometry& geom);
+
+ /** \brief
+ * Determines the [Location](@ref geom::Location) of a point in
+ * the polygonal elements of a
+ * [Geometry](@ref geom::Geometry).
+ *
+ * @param p the point to test
+ * @return the location of the point in the geometry
+ */
+ Location locate(const CoordinateXY* /*const*/ p);
+
+private:
+ void init();
+
+ // Declare type as noncopyable
+ IndexedPointInPolygonsLocator(const IndexedPointInPolygonsLocator& other) = delete;
+ IndexedPointInPolygonsLocator& operator=(const IndexedPointInPolygonsLocator& rhs) = delete;
+
+ //-- members
+ const Geometry& geom;
+ bool isInitialized;
+ TemplateSTRtree<IndexedPointInAreaLocator*> index;
+ std::vector<std::unique_ptr<IndexedPointInAreaLocator>> locators;
+};
+
+}}}
\ No newline at end of file
diff --git a/include/geos/algorithm/construct/LargestEmptyCircle.h b/include/geos/algorithm/construct/LargestEmptyCircle.h
index b84ccb3f4..f19931488 100644
--- a/include/geos/algorithm/construct/LargestEmptyCircle.h
+++ b/include/geos/algorithm/construct/LargestEmptyCircle.h
@@ -22,14 +22,13 @@
#include <geos/geom/Coordinate.h>
#include <geos/geom/Point.h>
#include <geos/geom/Envelope.h>
+#include <geos/algorithm/construct/IndexedDistanceToPoint.h>
#include <geos/algorithm/locate/IndexedPointInAreaLocator.h>
#include <geos/operation/distance/IndexedFacetDistance.h>
#include <memory>
#include <queue>
-
-
namespace geos {
namespace geom {
class Coordinate;
@@ -39,13 +38,9 @@ class GeometryFactory;
class LineString;
class Point;
}
-namespace operation {
-namespace distance {
-class IndexedFacetDistance;
-}
-}
}
+using geos::operation::distance::IndexedFacetDistance;
namespace geos {
namespace algorithm { // geos::algorithm
@@ -53,15 +48,23 @@ namespace construct { // geos::algorithm::construct
/**
* Constructs the Largest Empty Circle for a set of obstacle geometries,
-* up to a specified tolerance. The obstacles are point and line geometries.
+* up to a specified tolerance.
+* The obstacles may be any combination of point, linear and polygonal geometries.
*
-* The Largest Empty Circle is the largest circle which has its center
-* in the convex hull of the obstacles (the boundary), and whose
-* interior does not intersect with any obstacle. The circle center
-* is the point in the interior of the boundary which has the
-* farthest distance from the obstacles (up to tolerance).
-* The circle is determined by the center point and a point lying
-* on an obstacle indicating the circle radius.
+* The Largest Empty Circle (LEC) is the largest circle
+* whose interior does not intersect with any obstacle
+* and whose center lies within a polygonal boundary.
+* The circle center is the point in the interior of the boundary
+* which has the farthest distance from the obstacles
+* (up to the accuracy of the distance tolerance).
+* The circle itself is determined by the center point
+* and a point lying on an obstacle determining the circle radius.
+*
+* The polygonal boundary may be supplied explicitly.
+* If it is not specified the convex hull of the obstacles is used as the boundary.
+*
+* To compute an LEC which lies wholly within
+* a polygonal boundary, include the boundary of the polygon(s) as an obstacle.
*
* The implementation uses a successive-approximation technique
* over a grid of square cells covering the obstacles and boundary.
@@ -77,19 +80,36 @@ public:
/**
* Creates a new instance of a Largest Empty Circle construction.
+ * The obstacles may be any collection of points, lines and polygons.
+ * The constructed circle center lies within the convex hull of the obstacles.
*
- * @param p_obstacles a geometry representing the obstacles (points and lines)
+ * @param p_obstacles a geometry representing the obstacles
* @param p_tolerance the distance tolerance for computing the circle center point
*/
LargestEmptyCircle(const geom::Geometry* p_obstacles, double p_tolerance);
+
+ /**
+ * Creates a new instance of a Largest Empty Circle construction,
+ * interior-disjoint to a set of obstacle geometries
+ * and having its center within a polygonal boundary.
+ * The obstacles may be any collection of points, lines and polygons.
+ * If the boundary is null or empty the convex hull
+ * of the obstacles is used as the boundary.
+ *
+ * @param p_obstacles a geometry representing the obstacles
+ * @param p_boundary a polygonal geometry to contain the LEC center
+ * @param p_tolerance the distance tolerance for computing the circle center point
+ */
LargestEmptyCircle(const geom::Geometry* p_obstacles, const geom::Geometry* p_boundary, double p_tolerance);
+
~LargestEmptyCircle() = default;
/**
* Computes the center point of the Largest Empty Circle
* within a set of obstacles, up to a given tolerance distance.
+ * The obstacles may be any collection of points, lines and polygons.
*
- * @param p_obstacles a geometry representing the obstacles (points and lines)
+ * @param p_obstacles a geometry representing the obstacles
* @param p_tolerance the distance tolerance for computing the center point
* @return the center point of the Largest Empty Circle
*/
@@ -98,8 +118,9 @@ public:
/**
* Computes a radius line of the Largest Empty Circle
* within a set of obstacles, up to a given distance tolerance.
+ * The obstacles may be any collection of points, lines and polygons.
*
- * @param p_obstacles a geometry representing the obstacles (points and lines)
+ * @param p_obstacles a geometry representing the obstacles
* @param p_tolerance the distance tolerance for computing the center point
* @return a line from the center of the circle to a point on the edge
*/
@@ -118,10 +139,10 @@ private:
std::unique_ptr<geom::Geometry> boundary;
const geom::GeometryFactory* factory;
geom::Envelope gridEnv;
- operation::distance::IndexedFacetDistance obstacleDistance;
bool done;
- std::unique_ptr<algorithm::locate::IndexedPointInAreaLocator> ptLocator;
- std::unique_ptr<operation::distance::IndexedFacetDistance> boundaryDistance;
+ std::unique_ptr<algorithm::locate::IndexedPointInAreaLocator> boundaryPtLocater;
+ IndexedDistanceToPoint obstacleDistance;
+ std::unique_ptr<IndexedFacetDistance> boundaryDistance;
geom::CoordinateXY centerPt;
geom::CoordinateXY radiusPt;
diff --git a/include/geos/geom/util/PolygonalExtracter.h b/include/geos/geom/util/PolygonalExtracter.h
new file mode 100644
index 000000000..9dfba65d5
--- /dev/null
+++ b/include/geos/geom/util/PolygonalExtracter.h
@@ -0,0 +1,54 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2001-2002 Vivid Solutions Inc.
+ * Copyright (C) 2006 Refractions Research Inc.
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#pragma once
+
+#include <geos/export.h>
+#include <geos/geom/Geometry.h>
+#include <vector>
+
+namespace geos {
+namespace geom { // geos.geom
+namespace util { // geos.geom.util
+
+/**
+ * \brief Extracts the polygonal (Polygon and MultiPolygon) elements from a Geometry.
+ */
+class GEOS_DLL PolygonalExtracter {
+
+public:
+
+ /**
+ * Pushes the polygonal (Polygon and MultiPolygon) elements from a geometry into
+ * the provided vector.
+ *
+ * @param geom the geometry to extract from
+ * @param polys the vector to add the polygonal elements to
+ */
+ static void getPolygonals(const Geometry& geom, std::vector<const Geometry*>& polys);
+
+private:
+
+ static void getPolygonals(const Geometry* geom, std::vector<const Geometry*>& polys);
+
+ // Declare type as noncopyable
+ PolygonalExtracter(const PolygonalExtracter& other) = delete;
+ PolygonalExtracter& operator=(const PolygonalExtracter& rhs) = delete;
+};
+
+} // namespace geos.geom.util
+} // namespace geos.geom
+} // namespace geos
+
diff --git a/src/algorithm/construct/IndexedDistanceToPoint.cpp b/src/algorithm/construct/IndexedDistanceToPoint.cpp
new file mode 100644
index 000000000..92c703ef8
--- /dev/null
+++ b/src/algorithm/construct/IndexedDistanceToPoint.cpp
@@ -0,0 +1,68 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2023 Martin Davis <mtnclimb at gmail.com>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************
+ *
+ * Last port: algorithm/construct/IndexedDistanceToPoint.java
+ * https://github.com/locationtech/jts/commit/d92f783163d9440fcc10c729143787bf7b9fe8f9
+ *
+ **********************************************************************/
+
+#include <geos/algorithm/construct/IndexedDistanceToPoint.h>
+#include <geos/geom/Location.h>
+
+//using namespace geos::geom;
+using geos::geom::Location;
+
+namespace geos {
+namespace algorithm { // geos.algorithm
+namespace construct { // geos.algorithm.construct
+
+IndexedDistanceToPoint::IndexedDistanceToPoint(const Geometry& geom)
+ : targetGeometry(geom)
+{
+}
+
+/* private */
+void IndexedDistanceToPoint::init()
+{
+ if (facetDistance != nullptr)
+ return;
+ ptLocator.reset(new IndexedPointInPolygonsLocator(targetGeometry));
+ facetDistance.reset(new IndexedFacetDistance(&targetGeometry));
+}
+
+/* public */
+double IndexedDistanceToPoint::distance(const Point& pt)
+{
+ init();
+ //-- distance is 0 if point is inside a target polygon
+ if (isInArea(pt)) {
+ return 0;
+ }
+ return facetDistance->distance(&pt);
+}
+
+/* private */
+bool IndexedDistanceToPoint::isInArea(const Point& pt)
+{
+ return Location::EXTERIOR != ptLocator->locate(pt.getCoordinate());
+}
+
+/* public */
+std::unique_ptr<geom::CoordinateSequence>
+IndexedDistanceToPoint::nearestPoints(const Point& pt)
+{
+ return facetDistance->nearestPoints(&pt);
+}
+
+}}}
\ No newline at end of file
diff --git a/src/algorithm/construct/IndexedPointInPolygonsLocator.cpp b/src/algorithm/construct/IndexedPointInPolygonsLocator.cpp
new file mode 100644
index 000000000..bb9e51abe
--- /dev/null
+++ b/src/algorithm/construct/IndexedPointInPolygonsLocator.cpp
@@ -0,0 +1,70 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2023 Martin Davis <mtnclimb at gmail.com>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************
+ *
+ * Last port: algorithm/construct/IndexedDistanceToPoint.java
+ * https://github.com/locationtech/jts/commit/d92f783163d9440fcc10c729143787bf7b9fe8f9
+ *
+ **********************************************************************/
+
+#include <geos/algorithm/construct/IndexedPointInPolygonsLocator.h>
+#include <geos/geom/Envelope.h>
+#include <geos/geom/util/PolygonalExtracter.h>
+#include <vector>
+
+using geos::algorithm::locate::IndexedPointInAreaLocator;
+using geos::geom::Envelope;
+using geos::geom::util::PolygonalExtracter;
+
+namespace geos {
+namespace algorithm { // geos.algorithm
+namespace construct { // geos.algorithm.construct
+
+/* public */
+IndexedPointInPolygonsLocator::IndexedPointInPolygonsLocator(const Geometry& g)
+ : geom(g), isInitialized(false)
+{
+}
+
+/* private */
+void IndexedPointInPolygonsLocator::init()
+{
+ if (isInitialized) {
+ return;
+ }
+ isInitialized = true;
+ std::vector<const Geometry*> polys;
+ PolygonalExtracter::getPolygonals(geom, polys);
+ for (const Geometry* poly : polys) {
+ IndexedPointInAreaLocator* ptLocator = new IndexedPointInAreaLocator(*poly);
+ locators.emplace_back(ptLocator);
+ index.insert(poly->getEnvelopeInternal(), ptLocator);
+ }
+}
+
+/* public */
+Location IndexedPointInPolygonsLocator::locate(const CoordinateXY* /*const*/ pt)
+{
+ init();
+ Envelope queryEnv(*pt);
+ std::vector<IndexedPointInAreaLocator*> result;
+ index.query(queryEnv, result);
+ for (IndexedPointInAreaLocator* ptLocater : result) {
+ Location loc = ptLocater->locate(pt);
+ if (loc != Location::EXTERIOR)
+ return loc;
+ }
+ return Location::EXTERIOR;
+}
+
+}}}
\ No newline at end of file
diff --git a/src/algorithm/construct/LargestEmptyCircle.cpp b/src/algorithm/construct/LargestEmptyCircle.cpp
index 68f50223a..c4fc40f09 100644
--- a/src/algorithm/construct/LargestEmptyCircle.cpp
+++ b/src/algorithm/construct/LargestEmptyCircle.cpp
@@ -49,8 +49,8 @@ LargestEmptyCircle::LargestEmptyCircle(const Geometry* p_obstacles, const Geomet
: tolerance(p_tolerance)
, obstacles(p_obstacles)
, factory(p_obstacles->getFactory())
- , obstacleDistance(p_obstacles)
, done(false)
+ , obstacleDistance(*p_obstacles)
{
if (obstacles->isEmpty()) {
throw util::IllegalArgumentException("Empty obstacles geometry is not supported");
@@ -164,14 +164,14 @@ LargestEmptyCircle::mayContainCircleCenter(const Cell& cell, const Cell& farthes
double
LargestEmptyCircle::distanceToConstraints(const Coordinate& c)
{
- bool isOutside = ptLocator && (Location::EXTERIOR == ptLocator->locate(&c));
+ bool isOutside = boundaryPtLocater && (Location::EXTERIOR == boundaryPtLocater->locate(&c));
std::unique_ptr<Point> pt(factory->createPoint(c));
if (isOutside) {
double boundaryDist = boundaryDistance->distance(pt.get());
return -boundaryDist;
}
- double dist = obstacleDistance.distance(pt.get());
+ double dist = obstacleDistance.distance(*(pt.get()));
return dist;
}
@@ -200,8 +200,8 @@ LargestEmptyCircle::initBoundary()
gridEnv = *(boundary->getEnvelopeInternal());
// if boundary does not enclose an area cannot create a ptLocator
if (boundary->getDimension() >= 2) {
- ptLocator.reset(new algorithm::locate::IndexedPointInAreaLocator(*(boundary.get())));
- boundaryDistance.reset(new operation::distance::IndexedFacetDistance(boundary.get()));
+ boundaryPtLocater = detail::make_unique<algorithm::locate::IndexedPointInAreaLocator>(*(boundary.get()));
+ boundaryDistance = detail::make_unique<operation::distance::IndexedFacetDistance>(boundary.get());
}
}
@@ -214,7 +214,7 @@ LargestEmptyCircle::compute()
initBoundary();
// if ptLocator is not present then result is degenerate (represented as zero-radius circle)
- if (!ptLocator) {
+ if (!boundaryPtLocater) {
const CoordinateXY* pt = obstacles->getCoordinate();
centerPt = *pt;
radiusPt = *pt;
@@ -273,7 +273,7 @@ LargestEmptyCircle::compute()
// compute radius point
std::unique_ptr<Point> centerPoint(factory->createPoint(centerPt));
- const auto& nearestPts = obstacleDistance.nearestPoints(centerPoint.get());
+ const auto& nearestPts = obstacleDistance.nearestPoints(*(centerPoint.get()));
radiusPt = nearestPts->getAt(0);
// flag computation
diff --git a/src/geom/util/PolygonalExtracter.cpp b/src/geom/util/PolygonalExtracter.cpp
new file mode 100644
index 000000000..009b10969
--- /dev/null
+++ b/src/geom/util/PolygonalExtracter.cpp
@@ -0,0 +1,49 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2001-2002 Vivid Solutions Inc.
+ * Copyright (C) 2006 Refractions Research Inc.
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#include <geos/export.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/MultiPolygon.h>
+#include <vector>
+
+#include <geos/geom/util/PolygonalExtracter.h>
+
+namespace geos {
+namespace geom { // geos.geom
+namespace util { // geos.geom.util
+
+void
+PolygonalExtracter::getPolygonals(const Geometry& geom, std::vector<const Geometry*>& polys)
+{
+ getPolygonals(&geom, polys);
+}
+
+void
+PolygonalExtracter::getPolygonals(const Geometry* geom, std::vector<const Geometry*>& polys)
+{
+ if (dynamic_cast<const Polygon*>(geom) != nullptr
+ || dynamic_cast<const MultiPolygon*>(geom) != nullptr ) {
+ polys.push_back(geom);
+ }
+ else if (dynamic_cast<const GeometryCollection*>(geom) != nullptr) {
+ for (std::size_t i = 0; i < geom->getNumGeometries(); i++) {
+ getPolygonals(geom->getGeometryN(i), polys);
+ }
+ }
+}
+
+}
+}
+}
diff --git a/tests/unit/algorithm/construct/LargestEmptyCircleTest.cpp b/tests/unit/algorithm/construct/LargestEmptyCircleTest.cpp
index 959beb4a0..b0110fa0f 100644
--- a/tests/unit/algorithm/construct/LargestEmptyCircleTest.cpp
+++ b/tests/unit/algorithm/construct/LargestEmptyCircleTest.cpp
@@ -1,30 +1,21 @@
//
// Test Suite for geos::algorithm::construct::LargestEmptyCircle
-
#include <tut/tut.hpp>
// geos
-#include <geos/operation/distance/IndexedFacetDistance.h>
#include <geos/algorithm/construct/LargestEmptyCircle.h>
#include <geos/geom/Coordinate.h>
-#include <geos/geom/Dimension.h>
#include <geos/geom/Geometry.h>
#include <geos/geom/GeometryFactory.h>
-#include <geos/geom/LineString.h>
#include <geos/geom/PrecisionModel.h>
#include <geos/io/WKTReader.h>
-#include <geos/io/WKTWriter.h>
#include <geos/constants.h>
// std
#include <sstream>
#include <string>
#include <memory>
-
-
-using namespace geos;
using namespace geos::geom;
-
using geos::algorithm::construct::LargestEmptyCircle;
namespace tut {
@@ -32,30 +23,24 @@ namespace tut {
// Test Group
//
-// dummy data, not used
-struct test_lec_data {
- geos::geom::Geometry* geom_;
- geos::geom::PrecisionModel pm_;
- geos::geom::GeometryFactory::Ptr factory_;
+struct test_data_LargestEmptyCircle {
+ PrecisionModel pm_;
+ GeometryFactory::Ptr factory_;
geos::io::WKTReader reader_;
- geos::io::WKTWriter writer_;
- test_lec_data():
- geom_(nullptr),
- pm_(geos::geom::PrecisionModel::FLOATING),
+ test_data_LargestEmptyCircle():
+ pm_(PrecisionModel::FLOATING),
factory_(GeometryFactory::create(&pm_, 0)),
reader_(factory_.get())
{}
- ~test_lec_data()
+ ~test_data_LargestEmptyCircle()
{
- factory_->destroyGeometry(geom_);
- geom_ = nullptr;
}
void
- ensure_equals_coordinate(const geos::geom::Coordinate &lhs,
- const geos::geom::Coordinate &rhs, double tolerance)
+ ensure_equals_coordinate(const Coordinate &lhs,
+ const Coordinate &rhs, double tolerance)
{
ensure_equals("x coordinate does not match", lhs.x, rhs.x, tolerance);
ensure_equals("y coordinate does not match", lhs.y, rhs.y, tolerance);
@@ -145,117 +130,88 @@ struct test_lec_data {
};
-typedef test_group<test_lec_data> group;
+typedef test_group<test_data_LargestEmptyCircle> group;
typedef group::object object;
-group test_lec_group("geos::algorithm::construct::LargestEmptyCircle");
+group test_group_LargestEmptyCircle("geos::algorithm::construct::LargestEmptyCircle");
-//
// testPointsSquare
-//
+template<>
template<>
-template<>
-void object::test<1>
-()
+void object::test<1>()
{
-
checkCircle("MULTIPOINT ((100 100), (100 200), (200 200), (200 100))",
- 0.01, 150, 150, 70.71 );}
+ 0.01, 150, 150, 70.71 );
+}
-//
// testPointsTriangleOnHull
-//
template<>
template<>
-void object::test<2>
-()
+void object::test<2>()
{
checkCircle("MULTIPOINT ((100 100), (300 100), (150 50))",
0.01, 216.66, 99.99, 83.33 );
}
-//
// testPointsTriangleInterior
-//
template<>
template<>
-void object::test<3>
-()
+void object::test<3>()
{
checkCircle("MULTIPOINT ((100 100), (300 100), (200 250))",
0.01, 200.00, 141.66, 108.33 );
}
-//
// testLinesOpenDiamond
-//
template<>
template<>
-void object::test<4>
-()
+void object::test<4>()
{
checkCircle("MULTILINESTRING ((50 100, 150 50), (250 50, 350 100), (350 150, 250 200), (50 150, 150 200))",
- 0.01, 200, 125, 90.13 );}
+ 0.01, 200, 125, 90.13 );
+}
-//
// testLinesCrossed
-//
template<>
template<>
-void object::test<5>
-()
+void object::test<5>()
{
checkCircle("MULTILINESTRING ((100 100, 300 300), (100 200, 300 0))",
0.01, 299.99, 150.00, 106.05 );
}
-
-//
// testLinesZigzag
-//
template<>
template<>
-void object::test<6>
-()
+void object::test<6>()
{
-
checkCircle("MULTILINESTRING ((100 100, 200 150, 100 200, 250 250, 100 300, 300 350, 100 400), (70 380, 0 350, 50 300, 0 250, 50 200, 0 150, 50 120))",
0.01, 77.52, 249.99, 54.81 );
}
-//
// testPointsLinesTriangle
-//
template<>
template<>
-void object::test<7>
-()
+void object::test<7>()
{
checkCircle("GEOMETRYCOLLECTION (LINESTRING (100 100, 300 100), POINT (250 200))",
0.01, 196.49, 164.31, 64.31 );
}
-
-//
// testPointsLinesTriangle
-//
template<>
template<>
-void object::test<8>
-()
+void object::test<8>()
{
checkCircleZeroRadius("POINT (100 100)",
0.01 );
}
-//
// testLineFlat
-//
template<>
template<>
-void object::test<9>
-()
+void object::test<9>()
{
checkCircleZeroRadius("LINESTRING (0 0, 50 50)",
0.01 );
@@ -264,18 +220,47 @@ void object::test<9>
// testThinExtent
template<>
template<>
-void object::test<10>
-()
+void object::test<10>()
{
checkCircle("MULTIPOINT ((100 100), (300 100), (200 100.1))",
0.01 );
}
+//------------ Polygon Obstacles -----------------
+
+// testPolygonConcave
+template<>
+template<>
+void object::test<11>()
+{
+ checkCircle("POLYGON ((1 9, 9 6, 6 5, 5 3, 8 3, 9 4, 9 1, 1 1, 1 9))",
+ 0.01, 7.495, 4.216, 1.21);
+}
+
+// testPolygonsBoxes
+template<>
+template<>
+void object::test<12>()
+{
+ checkCircle("MULTIPOLYGON (((1 6, 6 6, 6 1, 1 1, 1 6)), ((6 7, 4 7, 4 9, 6 9, 6 7)))",
+ 0.01, 2.50, 7.50, 1.50);
+}
+
+// testPolygonLines
+template<>
+template<>
+void object::test<13>()
+{
+ checkCircle("GEOMETRYCOLLECTION (POLYGON ((1 6, 6 6, 6 1, 1 1, 1 6)), LINESTRING (6 7, 3 9), LINESTRING (1 7, 3 8))",
+ 0.01, 3.74, 7.14, 1.14);
+}
+
+//---------- Obstacles and Boundary --------------
+
// testBoundaryEmpty
template<>
template<>
-void object::test<11>
-()
+void object::test<14>()
{
checkCircle("MULTIPOINT ((2 2), (8 8), (7 5))",
"POLYGON EMPTY",
@@ -285,8 +270,7 @@ void object::test<11>
// testBoundarySquare
template<>
template<>
-void object::test<12>
-()
+void object::test<15>()
{
checkCircle("MULTIPOINT ((2 2), (6 4), (8 8))",
"POLYGON ((1 9, 9 9, 9 1, 1 1, 1 9))",
@@ -296,8 +280,7 @@ void object::test<12>
//testBoundarySquareObstaclesOutside
template<>
template<>
-void object::test<13>
-()
+void object::test<16>()
{
checkCircle("MULTIPOINT ((10 10), (10 0))",
"POLYGON ((1 9, 9 9, 9 1, 1 1, 1 9))",
@@ -307,8 +290,7 @@ void object::test<13>
// testBoundaryMultiSquares
template<>
template<>
-void object::test<14>
-()
+void object::test<17>()
{
checkCircle("MULTIPOINT ((10 10), (10 0), (5 5))",
"MULTIPOLYGON (((1 9, 9 9, 9 1, 1 1, 1 9)), ((15 20, 20 20, 20 15, 15 15, 15 20)))",
@@ -318,12 +300,20 @@ void object::test<14>
// testBoundaryAsObstacle
template<>
template<>
-void object::test<15>
-()
+void object::test<18>()
{
- checkCircle("GEOMETRYCOLLECTION (POLYGON ((1 9, 9 9, 9 1, 1 1, 1 9)), POINT (4 3), POINT (7 6))",
+ checkCircle("GEOMETRYCOLLECTION (LINESTRING (1 9, 9 9, 9 1, 1 1, 1 9), POINT (4 3), POINT (7 6))",
"POLYGON ((1 9, 9 9, 9 1, 1 1, 1 9))",
0.01, 4, 6, 3 );
}
+// testObstacleEmptyElement
+template<>
+template<>
+void object::test<19>()
+{
+ checkCircle("GEOMETRYCOLLECTION (LINESTRING EMPTY, POINT (4 3), POINT (7 6), POINT (4 6))",
+ 0.01, 5.5, 4.5, 2.12 );
+}
+
} // namespace tut
diff --git a/tests/unit/capi/GEOSLargestEmptyCircleTest.cpp b/tests/unit/capi/GEOSLargestEmptyCircleTest.cpp
new file mode 100644
index 000000000..aff84e2ec
--- /dev/null
+++ b/tests/unit/capi/GEOSLargestEmptyCircleTest.cpp
@@ -0,0 +1,69 @@
+// $Id$
+//
+// Test Suite for C-API GEOSLargestEmptyCircle
+
+#include <tut/tut.hpp>
+// geos
+#include <geos_c.h>
+// std
+#include <cstdarg>
+#include <cstdio>
+#include <cstdlib>
+#include <cstring>
+
+#include "capi_test_utils.h"
+
+namespace tut {
+//
+// Test Group
+//
+
+// Common data used in test cases.
+struct test_data_GEOSLargestEmptyCircle : public capitest::utility {
+
+ test_data_GEOSLargestEmptyCircle()
+ {
+ }
+
+ ~test_data_GEOSLargestEmptyCircle()
+ {
+ }
+
+};
+
+typedef test_group<test_data_GEOSLargestEmptyCircle> group;
+typedef group::object object;
+
+group test_capi_largestemptycircle_group("capi::GEOSLargestEmptyCircle");
+
+//
+// Test Cases
+//
+
+// Points of a square
+template<>
+template<>
+void object::test<1>
+()
+{
+ input_ = GEOSGeomFromWKT("MULTIPOINT ((100 100), (100 200), (200 200), (200 100))");
+ result_ = GEOSLargestEmptyCircle(input_, nullptr, 0.001);
+ expected_ = GEOSGeomFromWKT("LINESTRING (150 150, 100 100)");
+ ensure_geometry_equals_exact(result_, expected_, 0.0001);
+}
+
+// Line obstacles with square boundary
+template<>
+template<>
+void object::test<2>
+()
+{
+ input_ = GEOSGeomFromWKT("MULTILINESTRING ((40 90, 90 60), (90 40, 40 10))");
+ geom2_ = GEOSGeomFromWKT("POLYGON ((0 100, 100 100, 100 0, 0 0, 0 100))");
+ result_ = GEOSLargestEmptyCircle(input_, geom2_, 0.001);
+ expected_ = GEOSGeomFromWKT("LINESTRING (0.00038147 49.99961853, 40 10)");
+ ensure_geometry_equals_exact(result_, expected_, 0.0001);
+}
+
+} // namespace tut
+
diff --git a/tests/unit/capi/GEOSMaximumInscribedCircleTest.cpp b/tests/unit/capi/GEOSMaximumInscribedCircleTest.cpp
index 556c9a993..60ef3f10b 100644
--- a/tests/unit/capi/GEOSMaximumInscribedCircleTest.cpp
+++ b/tests/unit/capi/GEOSMaximumInscribedCircleTest.cpp
@@ -1,6 +1,6 @@
// $Id$
//
-// Test Suite for C-API GEOSGetCentroid
+// Test Suite for C-API GEOSMaximumInscribedCircle
#include <tut/tut.hpp>
// geos
@@ -44,7 +44,7 @@ group test_capimaximuminscribedcircle_group("capi::GEOSMaximumInscribedCircle");
// Test Cases
//
-// Single point
+// Square
template<>
template<>
void object::test<1>
@@ -60,28 +60,11 @@ void object::test<1>
ensure_equals(std::string(wkt_), std::string("LINESTRING (150 150, 150 200)"));
}
-// Single point
-template<>
-template<>
-void object::test<2>
-()
-{
- geom1_ = GEOSGeomFromWKT("MULTIPOINT ((100 100), (100 200), (200 200), (200 100))");
- ensure(nullptr != geom1_);
- geom2_ = GEOSLargestEmptyCircle(geom1_, nullptr, 0.001);
- ensure(nullptr != geom2_);
-
- wkt_ = GEOSWKTWriter_write(wktw_, geom2_);
-
- ensure_equals(std::string(wkt_), std::string("LINESTRING (150 150, 100 100)"));
-}
-
-
// Crash with Inf coords
// https://github.com/libgeos/geos/issues/821
template<>
template<>
-void object::test<3>
+void object::test<2>
()
{
std::string wkb("0106000020E61000000100000001030000000100000005000000000000000000F07F000000000000F07F000000000000F07F000000000000F07F000000000000F07F000000000000F07F000000000000F07F000000000000F07F000000000000F07F000000000000F07F");
diff --git a/tests/unit/capi/capi_test_utils.h b/tests/unit/capi/capi_test_utils.h
index 11305e956..d6a9dc668 100644
--- a/tests/unit/capi/capi_test_utils.h
+++ b/tests/unit/capi/capi_test_utils.h
@@ -93,17 +93,30 @@ namespace capitest {
void
ensure_geometry_equals(GEOSGeometry* g1, GEOSGeometry* g2, double tolerance)
{
- GEOSNormalize(g1);
- GEOSNormalize(g2);
- char rslt = GEOSEqualsExact(g1, g2, tolerance);
- if (rslt != 1)
- {
- char* wkt1 = GEOSWKTWriter_write(wktw_, g1);
- char* wkt2 = GEOSGeomToWKT(g2);
- std::fprintf(stdout, "\n%s != %s\n", wkt1, wkt2);
- GEOSFree(wkt1);
- GEOSFree(wkt2);
+ char rslt;
+ if (g1 == nullptr || g2 == nullptr) {
+ rslt = (g1 == nullptr && g2 == nullptr) ? 1 : 0;
+ }
+ else {
+ GEOSNormalize(g1);
+ GEOSNormalize(g2);
+ rslt = GEOSEqualsExact(g1, g2, tolerance);
+ }
+ report_not_equal("ensure_equals_norm", g1, g2, tolerance, rslt);
+ tut::ensure_equals("GEOSEqualsExact(g1, g2, tolerance)", rslt, 1);
+ }
+
+ void
+ ensure_geometry_equals_exact(GEOSGeometry* g1, GEOSGeometry* g2, double tolerance)
+ {
+ char rslt;
+ if (g1 == nullptr || g2 == nullptr) {
+ rslt = (g1 == nullptr && g2 == nullptr) ? 1 : 0;
+ }
+ else {
+ rslt = GEOSEqualsExact(g1, g2, tolerance);
}
+ report_not_equal("ensure_equals_exact", g1, g2, tolerance, rslt);
tut::ensure_equals("GEOSEqualsExact(g1, g2, tolerance)", rslt, 1);
}
@@ -133,7 +146,24 @@ namespace capitest {
tut::ensure_equals("GEOSEqualsExact(g1, g2, 1e-12)", rslt, 1);
}
-
+ void
+ report_not_equal(const char* tag, GEOSGeometry* g1, GEOSGeometry* g2, double tolerance, char rslt)
+ {
+ if (rslt == 1) return;
+ //TODO: handle rslt exception value
+
+ char* wkt1 = nullptr;
+ char* wkt2 = nullptr;
+ if (g1 != nullptr)
+ wkt1 = GEOSWKTWriter_write(wktw_, g1);
+ if (g2 != nullptr)
+ wkt2 = GEOSWKTWriter_write(wktw_, g2);
+ const char* val1 = (g1 == nullptr) ? "null" : wkt1;
+ const char* val2 = (g2 == nullptr) ? "null" : wkt2;
+ std::fprintf(stdout, "\n%s : %s != %s (tol = %f)\n", tag, val1, val2, tolerance);
+ GEOSFree(wkt1);
+ GEOSFree(wkt2);
+ }
};
-----------------------------------------------------------------------
Summary of changes:
capi/geos_c.h.in | 11 +-
.../algorithm/construct/IndexedDistanceToPoint.h | 88 ++++++++++++
.../construct/IndexedPointInPolygonsLocator.h | 79 +++++++++++
.../geos/algorithm/construct/LargestEmptyCircle.h | 63 ++++++---
include/geos/geom/util/PolygonalExtracter.h | 54 ++++++++
src/algorithm/construct/IndexedDistanceToPoint.cpp | 68 +++++++++
.../construct/IndexedPointInPolygonsLocator.cpp | 70 ++++++++++
src/algorithm/construct/LargestEmptyCircle.cpp | 14 +-
.../{PointExtracter.cpp => PolygonalExtracter.cpp} | 42 ++----
.../algorithm/construct/LargestEmptyCircleTest.cpp | 152 ++++++++++-----------
tests/unit/capi/GEOSLargestEmptyCircleTest.cpp | 69 ++++++++++
tests/unit/capi/GEOSMaximumInscribedCircleTest.cpp | 23 +---
tests/unit/capi/capi_test_utils.h | 52 +++++--
13 files changed, 615 insertions(+), 170 deletions(-)
create mode 100644 include/geos/algorithm/construct/IndexedDistanceToPoint.h
create mode 100644 include/geos/algorithm/construct/IndexedPointInPolygonsLocator.h
create mode 100644 include/geos/geom/util/PolygonalExtracter.h
create mode 100644 src/algorithm/construct/IndexedDistanceToPoint.cpp
create mode 100644 src/algorithm/construct/IndexedPointInPolygonsLocator.cpp
copy src/geom/util/{PointExtracter.cpp => PolygonalExtracter.cpp} (50%)
create mode 100644 tests/unit/capi/GEOSLargestEmptyCircleTest.cpp
hooks/post-receive
--
GEOS
More information about the geos-commits
mailing list