Commit 8e6d642e authored by Sandro Santilli's avatar Sandro Santilli

UnaryUnion: Drop assumption about union not moving vertices

The assumption resulted in invalid geometries being output
by CascadedPolygonUnion intermediary results and then fed
as input by further union calls, which in turn would fail
on first robustness issue.

Fixes #837 in 3.6 branch
Includes automated XML test
parent 459ae2f5
Pipeline #9874591 passed with stage
in 25 minutes and 50 seconds
......@@ -2,6 +2,7 @@ Changes in 3.6.2dev
2017-MM-DD
- Bug fixes / improvements
- Fix exception in UnaryUnion of collection of touching polygons (#837)
- Allow building against python 3 (#774)
- Fix build with android-ndk and other compilers (#799)
- Allows compiling with -Wpointer-bool-conversion (#638)
......
......@@ -14,6 +14,7 @@
**********************************************************************
*
* Last port: operation/union/CascadedPolygonUnion.java r487 (JTS-1.12+)
* Includes custom code to deal with https://trac.osgeo.org/geos/ticket/837
*
**********************************************************************/
......@@ -226,6 +227,16 @@ private:
geom::Geometry* extractByEnvelope(geom::Envelope const& env,
geom::Geometry* geom, std::vector<geom::Geometry*>& disjointGeoms);
void extractByEnvelope(geom::Envelope const& env,
geom::Geometry* geom,
std::vector<geom::Geometry*>& intersectingGeoms,
std::vector<geom::Geometry*>& disjointGeoms);
void extractByEnvelope(geom::Envelope const& env,
std::vector<geom::Geometry*>& sourceGeoms,
std::vector<geom::Geometry*>& intersectingGeoms,
std::vector<geom::Geometry*>& disjointGeoms);
/**
* Encapsulates the actual unioning of two polygonal geometries.
*
......
......@@ -20,7 +20,7 @@
namespace geos {
namespace geom {
class Geometry;
}
}
}
namespace geos {
......
......@@ -8,12 +8,13 @@
*
* 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.
* by the Free Software Foundation.
* See the COPYING file for more information.
*
**********************************************************************
*
* Last port: operation/union/CascadedPolygonUnion.java r487 (JTS-1.12+)
* Includes custom code to deal with https://trac.osgeo.org/geos/ticket/837
*
**********************************************************************/
......@@ -30,6 +31,68 @@
#include <cstddef>
#include <memory>
#include <vector>
#include <sstream>
#include <geos/operation/valid/IsValidOp.h>
#include <geos/operation/IsSimpleOp.h>
#include <geos/algorithm/BoundaryNodeRule.h>
#include <geos/util/TopologyException.h>
#include <string>
#include <iomanip>
//#define GEOS_DEBUG_CASCADED_UNION 1
//#define GEOS_DEBUG_CASCADED_UNION_PRINT_INVALID 1
namespace {
inline bool
check_valid(const geos::geom::Geometry& g, const std::string& label, bool doThrow=false, bool validOnly=false)
{
using namespace geos;
if ( dynamic_cast<const geos::geom::Lineal*>(&g) ) {
if ( ! validOnly ) {
operation::IsSimpleOp sop(g, algorithm::BoundaryNodeRule::getBoundaryEndPoint());
if ( ! sop.isSimple() )
{
if ( doThrow ) {
throw geos::util::TopologyException(
label + " is not simple");
}
return false;
}
}
} else {
operation::valid::IsValidOp ivo(&g);
if ( ! ivo.isValid() )
{
using operation::valid::TopologyValidationError;
TopologyValidationError* err = ivo.getValidationError();
#ifdef GEOS_DEBUG_CASCADED_UNION
std::cerr << label << " is INVALID: "
<< err->toString()
<< " (" << std::setprecision(20)
<< err->getCoordinate() << ")"
<< std::endl
#ifdef GEOS_DEBUG_CASCADED_UNION_PRINT_INVALID
<< "<A>" << std::endl
<< g.toString()
<< std::endl
#endif
;
#endif // GEOS_DEBUG_CASCADED_UNION
if ( doThrow ) {
throw geos::util::TopologyException(
label + " is invalid: " + err->toString(),
err->getCoordinate());
}
return false;
}
}
return true;
}
} // anonymous namespace
namespace geos {
namespace operation { // geos.operation
......@@ -51,7 +114,7 @@ geom::Geometry* CascadedPolygonUnion::Union(std::vector<geom::Polygon*>* polys)
geom::Geometry* CascadedPolygonUnion::Union(const geom::MultiPolygon* multipoly)
{
std::vector<geom::Polygon*> polys;
typedef geom::MultiPolygon::const_iterator iterator;
iterator end = multipoly->end();
for (iterator i = multipoly->begin(); i != end; ++i)
......@@ -71,7 +134,7 @@ geom::Geometry* CascadedPolygonUnion::Union()
/**
* A spatial index to organize the collection
* into groups of close geometries.
* This makes unioning more efficient, since vertices are more likely
* This makes unioning more efficient, since vertices are more likely
* to be eliminated on each round.
*/
index::strtree::STRtree index(STRTREE_NODE_CAPACITY);
......@@ -104,7 +167,7 @@ geom::Geometry* CascadedPolygonUnion::binaryUnion(GeometryListHolder* geoms)
return binaryUnion(geoms, 0, geoms->size());
}
geom::Geometry* CascadedPolygonUnion::binaryUnion(GeometryListHolder* geoms,
geom::Geometry* CascadedPolygonUnion::binaryUnion(GeometryListHolder* geoms,
std::size_t start, std::size_t end)
{
if (end - start <= 1) {
......@@ -122,7 +185,7 @@ geom::Geometry* CascadedPolygonUnion::binaryUnion(GeometryListHolder* geoms,
}
}
GeometryListHolder*
GeometryListHolder*
CascadedPolygonUnion::reduceToGeometries(index::strtree::ItemsList* geomTree)
{
std::auto_ptr<GeometryListHolder> geoms (new GeometryListHolder());
......@@ -146,7 +209,7 @@ CascadedPolygonUnion::reduceToGeometries(index::strtree::ItemsList* geomTree)
return geoms.release();
}
geom::Geometry*
geom::Geometry*
CascadedPolygonUnion::unionSafe(geom::Geometry* g0, geom::Geometry* g1)
{
if (g0 == NULL && g1 == NULL)
......@@ -160,7 +223,7 @@ CascadedPolygonUnion::unionSafe(geom::Geometry* g0, geom::Geometry* g1)
return unionOptimized(g0, g1);
}
geom::Geometry*
geom::Geometry*
CascadedPolygonUnion::unionOptimized(geom::Geometry* g0, geom::Geometry* g1)
{
geom::Envelope const* g0Env = g0->getEnvelopeInternal();
......@@ -172,44 +235,124 @@ CascadedPolygonUnion::unionOptimized(geom::Geometry* g0, geom::Geometry* g1)
if (g0->getNumGeometries() <= 1 && g1->getNumGeometries() <= 1)
return unionActual(g0, g1);
geom::Envelope commonEnv;
geom::Envelope commonEnv;
g0Env->intersection(*g1Env, commonEnv);
return unionUsingEnvelopeIntersection(g0, g1, commonEnv);
}
geom::Geometry*
CascadedPolygonUnion::unionUsingEnvelopeIntersection(geom::Geometry* g0,
/* private */
geom::Geometry*
CascadedPolygonUnion::unionUsingEnvelopeIntersection(geom::Geometry* g0,
geom::Geometry* g1, geom::Envelope const& common)
{
std::vector<geom::Geometry*> disjointPolys;
#if GEOS_DEBUG_CASCADED_UNION
check_valid(*g0, "unionUsingEnvelopeIntersection g0");
check_valid(*g1, "unionUsingEnvelopeIntersection g1");
#endif
std::auto_ptr<geom::Geometry> g0Int(extractByEnvelope(common, g0, disjointPolys));
std::auto_ptr<geom::Geometry> g1Int(extractByEnvelope(common, g1, disjointPolys));
#if GEOS_DEBUG_CASCADED_UNION
check_valid(*g0Int, "unionUsingEnvelopeIntersection g0Int");
check_valid(*g1Int, "unionUsingEnvelopeIntersection g1Int");
#endif
std::auto_ptr<geom::Geometry> u(unionActual(g0Int.get(), g1Int.get()));
disjointPolys.push_back(u.get());
return geom::util::GeometryCombiner::combine(disjointPolys);
#if GEOS_DEBUG_CASCADED_UNION
check_valid(*u, "unionUsingEnvelopeIntersection unionActual return");
#endif
if ( disjointPolys.empty() ) return u.release();
#if GEOS_DEBUG_CASCADED_UNION
for ( size_t i=0; i<disjointPolys.size(); ++i )
{
std::ostringstream os; os << "dp"<< i;
check_valid(*(disjointPolys[i]), os.str());
}
#endif
// TODO: find, in disjointPolys, those which now have their
// environment intersect the environment of the union "u"
// and collect them in another vector to be unioned
std::vector<geom::Geometry*> polysOn;
std::vector<geom::Geometry*> polysOff;
geom::Envelope const* uEnv = u->getEnvelopeInternal(); // TODO: check for EMPTY ?
extractByEnvelope(*uEnv, disjointPolys, polysOn, polysOff);
#if GEOS_DEBUG_CASCADED_UNION
std::cerr << "unionUsingEnvelopeIntersection: " << polysOn.size() << "/" << disjointPolys.size() << " polys intersect union of final thing" << std::endl;
#endif
std::auto_ptr<geom::Geometry> ret;
if ( polysOn.empty() ) {
disjointPolys.push_back(u.get());
ret.reset( geom::util::GeometryCombiner::combine(disjointPolys));
} else {
// TODO: could be further tweaked to only union with polysOn
// and combine with polysOff, but then it'll need again to
// recurse in the check for disjoint/intersecting
ret.reset( geom::util::GeometryCombiner::combine(disjointPolys) );
ret.reset( unionActual(ret.get(), u.get()) );
}
#if GEOS_DEBUG_CASCADED_UNION
check_valid(*ret, "unionUsingEnvelopeIntersection returned geom");
#endif
return ret.release();
}
geom::Geometry*
CascadedPolygonUnion::extractByEnvelope(geom::Envelope const& env,
geom::Geometry* geom, std::vector<geom::Geometry*>& disjointGeoms)
/* private */
void
CascadedPolygonUnion::extractByEnvelope(geom::Envelope const& env,
std::vector<geom::Geometry*>& sourceGeoms,
std::vector<geom::Geometry*>& intersectingGeoms,
std::vector<geom::Geometry*>& disjointGeoms)
{
std::vector<geom::Geometry*> intersectingGeoms;
for (std::vector<geom::Geometry*>::iterator i=sourceGeoms.begin(),
e=sourceGeoms.end(); i!=e; ++i)
{
geom::Geometry* elem = *i;
if (elem->getEnvelopeInternal()->intersects(env))
intersectingGeoms.push_back(elem);
else
disjointGeoms.push_back(elem);
}
}
for (std::size_t i = 0; i < geom->getNumGeometries(); i++) {
/* private */
void
CascadedPolygonUnion::extractByEnvelope(geom::Envelope const& env,
geom::Geometry* geom,
std::vector<geom::Geometry*>& intersectingGeoms,
std::vector<geom::Geometry*>& disjointGeoms)
{
for (std::size_t i = 0; i < geom->getNumGeometries(); i++) {
geom::Geometry* elem = const_cast<geom::Geometry*>(geom->getGeometryN(i));
if (elem->getEnvelopeInternal()->intersects(env))
intersectingGeoms.push_back(elem);
else
disjointGeoms.push_back(elem);
}
}
/* private */
geom::Geometry*
CascadedPolygonUnion::extractByEnvelope(geom::Envelope const& env,
geom::Geometry* geom, std::vector<geom::Geometry*>& disjointGeoms)
{
std::vector<geom::Geometry*> intersectingGeoms;
extractByEnvelope(env, geom, intersectingGeoms, disjointGeoms);
return geomFactory->buildGeometry(intersectingGeoms);
}
geom::Geometry*
geom::Geometry*
CascadedPolygonUnion::unionActual(geom::Geometry* g0, geom::Geometry* g1)
{
return restrictToPolygons(std::auto_ptr<geom::Geometry>(g0->Union(g1))).release();
......@@ -226,7 +369,7 @@ CascadedPolygonUnion::restrictToPolygons(std::auto_ptr<geom::Geometry> g)
}
Polygon::ConstVect polygons;
util::PolygonExtracter::getPolygons(*g, polygons);
geom::util::PolygonExtracter::getPolygons(*g, polygons);
if (polygons.size() == 1)
return std::auto_ptr<Geometry>(polygons[0]->clone());
......
#
# This file is part of project GEOS (http://trac.osgeo.org/geos/)
# This file is part of project GEOS (http://trac.osgeo.org/geos/)
#
#prefix=@prefix@
#top_srcdir=@top_srcdir@
......@@ -9,7 +9,7 @@ AUTOMAKE_OPTIONS = subdir-objects
TESTS = testrunner
CLEANFILES = testrunner
CLEANFILES = testrunner
EXTRA_DIST = testrunner.sh CMakeLists.txt
......@@ -43,6 +43,7 @@ SAFE_XMLTESTS=$(srcdir)/tests/testLeaksBig.xml \
$(srcdir)/tests/ticket/bug605.xml \
$(srcdir)/tests/ticket/bug615.xml \
$(srcdir)/tests/ticket/bug716.xml \
$(srcdir)/tests/ticket/bug837.xml \
$(srcdir)/tests/general/TestBoundary.xml \
$(srcdir)/tests/general/TestBuffer.xml \
$(srcdir)/tests/general/TestBufferMitredJoin.xml \
......@@ -93,11 +94,11 @@ SAFE_XMLTESTS=$(srcdir)/tests/testLeaksBig.xml \
$(srcdir)/tests/safe/16596.xml \
$(srcdir)/tests/safe/TestBufferJagged.xml
INVALID_OUTPUT_XMLTESTS =
INVALID_OUTPUT_XMLTESTS =
FAILING_XMLTESTS = \
$(srcdir)/tests/failure/TestOverlay.xml \
$(srcdir)/tests/ticket/bug488.xml
$(srcdir)/tests/ticket/bug488.xml
XMLTESTS=$(SAFE_XMLTESTS) $(INVALID_OUTPUT_XMLTESTS) $(FAILING_XMLTESTS)
......
This source diff could not be displayed because it is too large. You can view the blob instead.
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment