[geos-commits] [SCM] GEOS branch master updated. 5f11ee083576cfef840981181f33f431a779905b

git at osgeo.org git at osgeo.org
Thu Feb 28 12:38:05 PST 2019


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, master has been updated
       via  5f11ee083576cfef840981181f33f431a779905b (commit)
       via  5cf6893d58d92bbdd7e6bf5cc1c9a893c20160ca (commit)
       via  e04ba73ca7f2a033a1867251c5d5fff46491cfc7 (commit)
       via  589a05f24d08fc755e76270734d397a18b3cbaf2 (commit)
      from  a2fc78f764ec5edfc6b37b1ec9480de21f6bc0ef (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 5f11ee083576cfef840981181f33f431a779905b
Merge: a2fc78f 5cf6893
Author: Daniel Baston <dbaston at gmail.com>
Date:   Thu Feb 28 15:37:26 2019 -0500

    Merge remote-tracking branch 'rouault/add_GEOSMakeValid_cpp'


commit 5cf6893d58d92bbdd7e6bf5cc1c9a893c20160ca
Author: Even Rouault <even.rouault at spatialys.com>
Date:   Thu Feb 28 16:21:11 2019 +0100

    NEWS: mention GEOSBuildArea and GEOSMakeValid

diff --git a/NEWS b/NEWS
index 10bff91..481f3c0 100644
--- a/NEWS
+++ b/NEWS
@@ -1,3 +1,8 @@
+3.8.0 changes
+- New things:
+  - CAPI: GEOSBuildArea (#952, Even Rouault)
+  - CAPI: GEOSMakeValid (#952, Even Rouault)
+
 Changes in 3.7.0rc1
 2018-08-19
 Fixes / enhancements since 3.7.0beta2

commit e04ba73ca7f2a033a1867251c5d5fff46491cfc7
Author: Even Rouault <even.rouault at spatialys.com>
Date:   Mon Feb 25 22:23:56 2019 +0100

    Add operation::valid::MakeValid
    
    and add GEOSMakeValid_r() (fixes #952)

diff --git a/capi/geos_c.cpp b/capi/geos_c.cpp
index ce1986f..78e28c8 100644
--- a/capi/geos_c.cpp
+++ b/capi/geos_c.cpp
@@ -731,6 +731,12 @@ extern "C" {
     }
 
     Geometry*
+    GEOSMakeValid(const Geometry* g)
+    {
+        return GEOSMakeValid_r(handle, g);
+    }
+
+    Geometry*
     GEOSLineMerge(const Geometry* g)
     {
         return GEOSLineMerge_r(handle, g);
diff --git a/capi/geos_c.h.in b/capi/geos_c.h.in
index d8eb68e..6b5fc57 100644
--- a/capi/geos_c.h.in
+++ b/capi/geos_c.h.in
@@ -959,6 +959,9 @@ extern char GEOS_DLL GEOSisValidDetail_r(GEOSContextHandle_t handle,
                                          char** reason,
                                          GEOSGeometry** location);
 
+extern GEOSGeometry GEOS_DLL *GEOSMakeValid_r(GEOSContextHandle_t handle,
+                                              const GEOSGeometry* g);
+
 /************************************************************************
  *
  *  Geometry info
@@ -1937,6 +1940,8 @@ extern char GEOS_DLL GEOSisValidDetail(const GEOSGeometry* g,
                                        int flags,
                                        char** reason, GEOSGeometry** location);
 
+extern GEOSGeometry GEOS_DLL *GEOSMakeValid(const GEOSGeometry* g);
+
 /************************************************************************
  *
  *  Geometry info
diff --git a/capi/geos_ts_c.cpp b/capi/geos_ts_c.cpp
index 2306d75..1beaa76 100644
--- a/capi/geos_ts_c.cpp
+++ b/capi/geos_ts_c.cpp
@@ -68,6 +68,7 @@
 #include <geos/operation/sharedpaths/SharedPathsOp.h>
 #include <geos/operation/union/CascadedPolygonUnion.h>
 #include <geos/operation/valid/IsValidOp.h>
+#include <geos/operation/valid/MakeValid.h>
 #include <geos/precision/GeometryPrecisionReducer.h>
 #include <geos/linearref/LengthIndexedLine.h>
 #include <geos/triangulate/DelaunayTriangulationBuilder.h>
@@ -3187,6 +3188,37 @@ extern "C" {
     }
 
     Geometry*
+    GEOSMakeValid_r(GEOSContextHandle_t extHandle, const Geometry* g)
+    {
+        if(nullptr == extHandle) {
+            return nullptr;
+        }
+
+        GEOSContextHandleInternal_t* handle = nullptr;
+        handle = reinterpret_cast<GEOSContextHandleInternal_t*>(extHandle);
+        if(0 == handle->initialized) {
+            return nullptr;
+        }
+
+        Geometry* out = nullptr;
+
+        try {
+            // BuildArea
+            using geos::operation::valid::MakeValid;
+            MakeValid makeValid;
+            out = makeValid.build(g).release();
+        }
+        catch(const std::exception& e) {
+            handle->ERROR_MESSAGE("%s", e.what());
+        }
+        catch(...) {
+            handle->ERROR_MESSAGE("Unknown exception thrown");
+        }
+
+        return out;
+    }
+
+    Geometry*
     GEOSPolygonizer_getCutEdges_r(GEOSContextHandle_t extHandle, const Geometry* const* g, unsigned int ngeoms)
     {
         if(0 == extHandle) {
diff --git a/include/geos/operation/valid/MakeValid.h b/include/geos/operation/valid/MakeValid.h
new file mode 100644
index 0000000..d3d3c7b
--- /dev/null
+++ b/include/geos/operation/valid/MakeValid.h
@@ -0,0 +1,79 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright 2009-2010 Sandro Santilli <strk at kbt.io>
+ * Copyright (C) 2019 Even Rouault <even.rouault at spatialys.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.
+ **********************************************************************
+ *
+ * Ported from rtgeom_geos.c from
+ *   rttopo - topology library
+ *   http://git.osgeo.org/gitea/rttopo/librttopo
+ * with relicensing from GPL to LGPL with Copyright holder permission.
+ *
+ **********************************************************************/
+
+#ifndef GEOS_OP_VALID_MAKEVALID_H
+#define GEOS_OP_VALID_MAKEVALID_H
+
+#include <geos/export.h>
+
+#include <memory>
+
+#ifdef _MSC_VER
+#pragma warning(push)
+#pragma warning(disable: 4251) // warning C4251: needs to have dll-interface to be used by clients of class
+#endif
+
+// Forward declarations
+namespace geos {
+namespace geom {
+class Geometry;
+}
+}
+
+namespace geos {
+namespace operation { // geos::operation
+namespace valid { // geos::operation::valid
+
+/** \brief  The function attempts to create a valid representation of a given
+ * invalid geometry without losing any of the input vertices.
+ *
+ * Already-valid geometries are returned without further intervention.
+ *
+ * Supported inputs are: POINTS, MULTIPOINTS, LINESTRINGS, MULTILINESTRINGS, POLYGONS, MULTIPOLYGONS and GEOMETRYCOLLECTIONS containing any mix of them.
+ *
+ * In case of full or partial dimensional collapses, the output geometry may be a collection of lower-to-equal dimension geometries or a geometry of lower dimension.
+ *
+ * Single polygons may become multi-geometries in case of self-intersections.
+ */
+class GEOS_DLL MakeValid {
+
+public:
+
+    /** \brief
+     * Create a MakeValid object.
+     */
+    MakeValid() = default;
+
+    ~MakeValid() = default;
+
+    /** \brief Return a valid version of the input geometry. */
+    std::unique_ptr<geom::Geometry> build(const geom::Geometry* geom);
+};
+
+} // namespace geos::operation::valid
+} // namespace geos::operation
+} // namespace geos
+
+#ifdef _MSC_VER
+#pragma warning(pop)
+#endif
+
+#endif // GEOS_OP_VALID_MAKEVALID_H
diff --git a/src/operation/valid/MakeValid.cpp b/src/operation/valid/MakeValid.cpp
new file mode 100644
index 0000000..9064ba8
--- /dev/null
+++ b/src/operation/valid/MakeValid.cpp
@@ -0,0 +1,344 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright 2009-2010 Sandro Santilli <strk at kbt.io>
+ * Copyright (C) 2019 Even Rouault <even.rouault at spatialys.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.
+ **********************************************************************
+ *
+ * Ported from rtgeom_geos.c from
+ *   rttopo - topology library
+ *   http://git.osgeo.org/gitea/rttopo/librttopo
+ * with relicensing from GPL to LGPL with Copyright holder permission.
+ *
+ **********************************************************************/
+
+#include <geos/operation/valid/MakeValid.h>
+#include <geos/operation/valid/IsValidOp.h>
+#include <geos/operation/polygonize/BuildArea.h>
+#include <geos/geom/Geometry.h>
+#include <geos/geom/GeometryCollection.h>
+#include <geos/geom/GeometryFactory.h>
+#include <geos/geom/LineString.h>
+#include <geos/geom/Point.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/MultiLineString.h>
+#include <geos/geom/MultiPolygon.h>
+#include <geos/util/UniqueCoordinateArrayFilter.h>
+#include <geos/util/UnsupportedOperationException.h>
+
+// std
+#include <cassert>
+#include <algorithm>
+#include <utility>
+#include <vector>
+
+#ifdef _MSC_VER
+#pragma warning(disable:4355)
+#endif
+
+using namespace geos::geom;
+
+namespace geos {
+namespace operation { // geos.operation
+namespace valid { // geos.operation.valid
+
+/*
+ * Fully node given linework
+ */
+static std::unique_ptr<geom::Geometry>
+nodeLineWithFirstCoordinate(const geom::Geometry* geom)
+{
+  /*
+   * Union with first geometry point, obtaining full noding
+   * and dissolving of duplicated repeated points
+   *
+   * TODO: substitute this with UnaryUnion?
+   */
+
+  if( geom->isEmpty() )
+      return nullptr;
+
+  const auto geomType = geom->getGeometryTypeId();
+  assert( geomType == GEOS_LINESTRING || geomType == GEOS_MULTILINESTRING );
+
+  geom::Geometry* point;
+  if( geomType == GEOS_LINESTRING ) {
+      auto line = dynamic_cast<const geom::LineString*>(geom);
+      assert(line);
+      point = line->getPointN(0);
+  } else {
+      auto mls = dynamic_cast<const geom::MultiLineString*>(geom);
+      assert(mls);
+      auto line = dynamic_cast<const geom::LineString*>(mls->getGeometryN(0));
+      assert(line);
+      point = line->getPointN(0);
+  }
+  auto noded = geom->Union(point);
+  delete point;
+
+  return std::unique_ptr<geom::Geometry>(noded);
+}
+
+
+static std::unique_ptr<geom::Geometry> MakeValidLine(const geom::LineString* line)
+{
+    return nodeLineWithFirstCoordinate(line);
+}
+
+static std::unique_ptr<geom::Geometry> MakeValidMultiLine(const geom::MultiLineString* mls)
+{
+
+    std::vector<std::unique_ptr<geom::Geometry>> points;
+    std::vector<std::unique_ptr<geom::Geometry>> lines;
+
+    for(const auto subgeom: *mls) {
+        auto line = dynamic_cast<const geom::LineString*>(subgeom);
+        assert(line);
+        auto validSubGeom = MakeValidLine(line);
+        if( !validSubGeom || validSubGeom->isEmpty() ) {
+            continue;
+        }
+        auto validLineType = validSubGeom->getGeometryTypeId();
+        if( validLineType == GEOS_POINT ) {
+            points.emplace_back(std::move(validSubGeom));
+        }
+        else if( validLineType == GEOS_LINESTRING ) {
+            lines.emplace_back(std::move(validSubGeom));
+        } else if( validLineType == GEOS_MULTILINESTRING ) {
+            auto mlsValid = dynamic_cast<const geom::MultiLineString*>(validSubGeom.get());
+            for(const auto subgeomMlsValid: *mlsValid) {
+                lines.emplace_back(subgeomMlsValid->clone());
+            }
+        } else {
+            throw util::UnsupportedOperationException();
+        }
+    }
+
+    std::unique_ptr<geom::Geometry> pointsRet;
+    if( !points.empty() ) {
+        if( points.size() > 1 ) {
+            auto pointsRawPtr = new std::vector<geom::Geometry*>(points.size());
+            for( size_t i = 0; i < points.size(); i++ ) {
+                auto pointMoved = std::move(points[i]);
+                (*pointsRawPtr)[i] = pointMoved.release();
+            }
+            pointsRet.reset(
+                mls->getFactory()->createMultiPoint(pointsRawPtr));
+        } else {
+            pointsRet = std::move(points[0]);
+        }
+    }
+
+    std::unique_ptr<geom::Geometry> linesRet;
+    if( !lines.empty() ) {
+        if( lines.size() > 1 ) {
+            auto linesRawPtr = new std::vector<geom::Geometry*>(lines.size());
+            for( size_t i = 0; i < lines.size(); i++ ) {
+                auto lineMoved = std::move(lines[i]);
+                (*linesRawPtr)[i] = lineMoved.release();
+            }
+            linesRet.reset(
+                mls->getFactory()->createMultiLineString(linesRawPtr));
+        } else {
+            linesRet = std::move(lines[0]);
+        }
+    }
+
+    if( pointsRet && linesRet ) {
+        auto geoms = new std::vector<geom::Geometry*>(2);
+        (*geoms)[0] = pointsRet.release();
+        (*geoms)[1] = linesRet.release();
+        return std::unique_ptr<geom::Geometry>(
+                    mls->getFactory()->createGeometryCollection(geoms));
+    } else if( pointsRet ) {
+        return pointsRet;
+    } else if( linesRet ) {
+        return linesRet;
+    }
+
+    return nullptr;
+}
+
+static std::unique_ptr<geom::Geometry> extractUniquePoints(const geom::Geometry* geom)
+{
+
+    // Code taken from GEOSGeom_extractUniquePoints_r()
+
+    /* 1: extract points */
+    std::vector<const geom::Coordinate*> coords;
+    geos::util::UniqueCoordinateArrayFilter filter(coords);
+    geom->apply_ro(&filter);
+
+    /* 2: for each point, create a geometry and put into a vector */
+    auto points = new std::vector<geom::Geometry*>();
+    points->reserve(coords.size());
+    const GeometryFactory* factory = geom->getFactory();
+    for(std::vector<const geom::Coordinate*>::iterator it = coords.begin(),
+            itE = coords.end();
+            it != itE; ++it) {
+        auto point = factory->createPoint(*(*it));
+        points->push_back(point);
+    }
+
+    /* 3: create a multipoint */
+    return std::unique_ptr<geom::Geometry>(factory->createMultiPoint(points));
+}
+
+static std::unique_ptr<geom::Geometry> MakeValidPoly(const geom::Geometry* geom)
+{
+    assert( geom->getGeometryTypeId() == GEOS_POLYGON ||
+            geom->getGeometryTypeId() == GEOS_MULTIPOLYGON );
+
+    std::unique_ptr<geom::Geometry> bound(geom->getBoundary());
+    if( !bound )
+        return nullptr;
+
+    /* Use noded boundaries as initial "cut" edges */
+    auto cut_edges = nodeLineWithFirstCoordinate(bound.get());
+    if( !cut_edges )
+        return nullptr;
+
+    /* NOTE: the noding process may drop lines collapsing to points.
+    *       We want to retrieve any of those */
+    auto pi = extractUniquePoints(bound.get());
+    auto po = extractUniquePoints(cut_edges.get());
+    std::unique_ptr<geom::Geometry> collapse_points(pi->difference(po.get()));
+    assert(collapse_points);
+    pi.reset();
+    po.reset();
+
+    /* And use an empty geometry as initial "area" */
+    const GeometryFactory* factory = geom->getFactory();
+    std::unique_ptr<geom::Geometry> area(factory->createPolygon());
+    assert(area);
+
+    /*
+    * See if an area can be build with the remaining edges
+    * and if it can, symdifference with the original area.
+    * Iterate this until no more polygons can be created
+    * with left-over edges.
+    */
+    while( cut_edges->getNumGeometries() ) {
+
+        // ASSUMPTION: cut_edges should already be fully noded
+        auto new_area = geos::operation::polygonize::BuildArea().build(cut_edges.get());
+        assert(new_area); // never return nullptr, but exception
+        if( new_area->isEmpty() ) {
+            /* no more rings can be build with thes edges */
+            break;
+        }
+
+        // We succeeded in building a ring !
+        // Save the new ring boundaries first (to compute further cut edges later)
+        std::unique_ptr<geom::Geometry> new_area_bound(new_area->getBoundary());
+        assert(new_area_bound);
+
+        // Now symdif new and old area
+        std::unique_ptr<geom::Geometry> symdif(area->symDifference(new_area.get()));
+        assert(symdif);
+
+        area = std::move(symdif);
+
+        /*
+        * Now let's re-set cut_edges with what's left
+        * from the original boundary.
+        * ASSUMPTION: only the previous cut-edges can be
+        *             left, so we don't need to reconsider
+        *             the whole original boundaries
+        *
+        * NOTE: this is an expensive operation.
+        *
+        */
+        std::unique_ptr<geom::Geometry> new_cut_edges(cut_edges->difference(new_area_bound.get()));
+        assert(new_cut_edges);
+
+        cut_edges = std::move(new_cut_edges);
+    }
+
+    auto vgeoms = new std::vector<Geometry*>(3);
+    unsigned int nvgeoms=0;
+
+    if( !area->isEmpty() ) {
+        (*vgeoms)[nvgeoms++] = area.release();
+    }
+    if( !cut_edges->isEmpty() ) {
+        (*vgeoms)[nvgeoms++] = cut_edges.release();
+    }
+    if( !collapse_points->isEmpty() ) {
+        (*vgeoms)[nvgeoms++] = collapse_points.release();
+    }
+
+    if( nvgeoms == 1 ) {
+        /* Return cut edges */
+        auto ret = std::unique_ptr<geom::Geometry>((*vgeoms)[0]);
+        delete vgeoms;
+        return ret;
+    }
+
+    /* Collect areas and lines (if any line) */
+    vgeoms->resize(nvgeoms);
+    return std::unique_ptr<geom::Geometry>(
+                    factory->createGeometryCollection(vgeoms));
+}
+
+static std::unique_ptr<geom::Geometry> MakeValidCollection(const geom::GeometryCollection* coll)
+{
+    auto validGeoms = new std::vector<geom::Geometry*>();
+    try {
+        for( auto geom: *coll )
+        {
+            validGeoms->push_back(MakeValid().build(geom).release());
+        }
+        return std::unique_ptr<geom::Geometry>(
+            GeometryFactory::create()->createGeometryCollection(validGeoms));
+    }
+    catch( ... ) {
+        for( auto geom: *validGeoms ) {
+            delete geom;
+        }
+        delete validGeoms;
+        throw;
+    }
+}
+
+/** Return a valid version of the input geometry. */
+std::unique_ptr<geom::Geometry> MakeValid::build(const geom::Geometry* geom)
+{
+
+    IsValidOp ivo(geom);
+    if( ivo.getValidationError() == nullptr ) {
+        return std::unique_ptr<geom::Geometry>(geom->clone());
+    }
+
+    auto typeId = geom->getGeometryTypeId();
+    if( typeId == GEOS_LINESTRING ) {
+        auto lineString = dynamic_cast<const LineString*>(geom);
+        return MakeValidLine(lineString);
+    }
+    if( typeId == GEOS_MULTILINESTRING ) {
+        auto mls = dynamic_cast<const MultiLineString*>(geom);
+        return MakeValidMultiLine(mls);
+    }
+    if( typeId == GEOS_POLYGON ||
+        typeId == GEOS_MULTIPOLYGON ) {
+        return MakeValidPoly(geom);
+    }
+    if( typeId == GEOS_GEOMETRYCOLLECTION ) {
+        auto coll = dynamic_cast<const GeometryCollection*>(geom);
+        return MakeValidCollection(coll);
+    }
+
+    throw util::UnsupportedOperationException();
+}
+
+} // namespace geos.operation.valid
+} // namespace geos.operation
+} // namespace geos
+
diff --git a/src/operation/valid/Makefile.am b/src/operation/valid/Makefile.am
index 4647f09..4c082a7 100644
--- a/src/operation/valid/Makefile.am
+++ b/src/operation/valid/Makefile.am
@@ -21,6 +21,7 @@ libopvalid_la_SOURCES = \
     SweeplineNestedRingTester.cpp \
     TopologyValidationError.cpp \
     IndexedNestedRingTester.cpp \
-    IndexedNestedRingTester.h
+    IndexedNestedRingTester.h \
+    MakeValid.cpp
 
 libopvalid_la_LIBADD =
diff --git a/tests/unit/Makefile.am b/tests/unit/Makefile.am
index da53fb9..3a9ceb6 100644
--- a/tests/unit/Makefile.am
+++ b/tests/unit/Makefile.am
@@ -177,7 +177,8 @@ geos_unit_SOURCES = \
 	capi/GEOSisValidDetailTest.cpp \
 	capi/GEOSisClosedTest.cpp \
 	capi/GEOSInterpolateTest.cpp \
-	capi/GEOSBuildArea.cpp
+	capi/GEOSBuildArea.cpp \
+	capi/GEOSMakeValid.cpp
 
 noinst_HEADERS = \
 	utility.h
diff --git a/tests/unit/capi/GEOSMakeValid.cpp b/tests/unit/capi/GEOSMakeValid.cpp
new file mode 100644
index 0000000..6bf8d21
--- /dev/null
+++ b/tests/unit/capi/GEOSMakeValid.cpp
@@ -0,0 +1,81 @@
+// Test Suite for C-API MakeValid
+
+#include <tut/tut.hpp>
+// geos
+#include <geos_c.h>
+// std
+#include <cstdarg>
+#include <cstdio>
+#include <cstdlib>
+#include <cmath>
+
+namespace tut {
+//
+// Test Group
+//
+
+// Common data used in test cases.
+struct test_capimakevalid_data {
+    GEOSWKTWriter* wktw_;
+    GEOSGeometry* geom1_ = nullptr;
+    GEOSGeometry* geom2_ = nullptr;
+
+    static void
+    notice(const char* fmt, ...)
+    {
+        std::fprintf(stdout, "NOTICE: ");
+
+        va_list ap;
+        va_start(ap, fmt);
+        std::vfprintf(stdout, fmt, ap);
+        va_end(ap);
+
+        std::fprintf(stdout, "\n");
+    }
+
+    std::string
+    toWKT(GEOSGeometry* g)
+    {
+        char* wkt = GEOSWKTWriter_write(wktw_, g);
+        std::string ret(wkt);
+        GEOSFree(wkt);
+        return ret;
+    }
+
+    test_capimakevalid_data()
+    {
+        initGEOS(notice, notice);
+        wktw_ = GEOSWKTWriter_create();
+        GEOSWKTWriter_setTrim(wktw_, 1);
+        GEOSWKTWriter_setOutputDimension(wktw_, 3);
+    }
+
+    ~test_capimakevalid_data()
+    {
+        GEOSGeom_destroy(geom1_);
+        GEOSGeom_destroy(geom2_);
+        GEOSWKTWriter_destroy(wktw_);
+        finishGEOS();
+    }
+
+};
+
+typedef test_group<test_capimakevalid_data> group;
+typedef group::object object;
+
+group test_capimakevalid_group("capi::GEOSMakeValid");
+
+//
+// Test Cases
+//
+
+template<>
+template<>
+void object::test<1>
+()
+{
+    geom1_ = GEOSGeomFromWKT("POLYGON((0 0,1 1,0 1,1 0,0 0))");
+    geom2_ = GEOSMakeValid(geom1_);
+    ensure_equals(toWKT(geom2_), std::string("MULTIPOLYGON (((0 0, 0.5 0.5, 1 0, 0 0)), ((0.5 0.5, 0 1, 1 1, 0.5 0.5)))"));
+}
+} // namespace tut
diff --git a/tests/xmltester/CMakeLists.txt b/tests/xmltester/CMakeLists.txt
index fd492ce..a3eed66 100644
--- a/tests/xmltester/CMakeLists.txt
+++ b/tests/xmltester/CMakeLists.txt
@@ -63,6 +63,7 @@ if(GEOS_ENABLE_TESTS)
     ${XMLTESTS_DIR}/test.xml
     ${XMLTESTS_DIR}/testLeaksBig.xml
     ${XMLTESTS_DIR}/buildarea.xml
+    ${XMLTESTS_DIR}/makevalid.xml
     ${XMLTESTS_DIR}/general/TestBoundary.xml
     ${XMLTESTS_DIR}/general/TestBuffer.xml
     ${XMLTESTS_DIR}/general/TestBufferMitredJoin.xml
diff --git a/tests/xmltester/Makefile.am b/tests/xmltester/Makefile.am
index 819fccc..6289e3d 100644
--- a/tests/xmltester/Makefile.am
+++ b/tests/xmltester/Makefile.am
@@ -31,6 +31,7 @@ SAFE_XMLTESTS= \
 	$(srcdir)/tests/test.xml \
 	$(srcdir)/tests/testLeaksBig.xml \
 	$(srcdir)/tests/buildarea.xml \
+	$(srcdir)/tests/makevalid.xml \
 	$(srcdir)/tests/general/TestBoundary.xml \
 	$(srcdir)/tests/general/TestBuffer.xml \
 	$(srcdir)/tests/general/TestBufferMitredJoin.xml \
diff --git a/tests/xmltester/XMLTester.cpp b/tests/xmltester/XMLTester.cpp
index e996002..3b65168 100644
--- a/tests/xmltester/XMLTester.cpp
+++ b/tests/xmltester/XMLTester.cpp
@@ -42,6 +42,7 @@
 #include <geos/operation/buffer/BufferParameters.h>
 #include <geos/operation/buffer/BufferOp.h>
 #include <geos/operation/polygonize/BuildArea.h>
+#include <geos/operation/valid/MakeValid.h>
 #include <geos/precision/MinimumClearance.h>
 #include <geos/util.h>
 #include <geos/util/Interrupt.h>
@@ -1774,6 +1775,32 @@ XMLTester::parseTest(const tinyxml2::XMLNode* node)
             }
         }
 
+        else if (opName == "makevalid")
+        {
+            GeomPtr gExpected(parseGeometry(opRes, "expected"));
+            gExpected->normalize();
+
+            auto gGot = geos::operation::valid::MakeValid().build(gA);
+            if( gGot )
+            {
+                GeomPtr gRealRes(gGot.release());
+                gRealRes->normalize();
+
+                if (gExpected->equals(gRealRes.get())) success=1;
+
+                actual_result=printGeom(gRealRes.get());
+                expected_result=printGeom(gExpected.get());
+                if( actual_result == expected_result ) success=1;
+
+                if ( testValidOutput )
+                    success &= int(testValid(gRealRes.get(), "result"));
+            }
+            else
+            {
+                success = false;
+            }
+        }
+
         else {
             std::cerr << *curr_file << ":";
             std::cerr << " case" << caseCount << ":";
diff --git a/tests/xmltester/tests/makevalid.xml b/tests/xmltester/tests/makevalid.xml
new file mode 100644
index 0000000..15b06ab
--- /dev/null
+++ b/tests/xmltester/tests/makevalid.xml
@@ -0,0 +1,157 @@
+<run>
+  <precisionModel scale="10" offsetx="0.0" offsety="0.0"/>
+
+<case>
+  <desc>makevalid/point/already_valid</desc>
+  <a>POINT(0 0)</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        POINT(0 0)
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/point/empty</desc>
+  <a>POINT EMPTY</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        POINT EMPTY
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/linestring/already_valid</desc>
+  <a>LINESTRING(0 0,1 1)</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        LINESTRING(0 0,1 1)
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/linestring/invalid_result_point</desc>
+  <a>LINESTRING(0 0,0 0)</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        POINT (0.0 0.0)
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/linestring/empty</desc>
+  <a>LINESTRING EMPTY</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        LINESTRING EMPTY
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/mulilinestring/empty</desc>
+  <a>MULTILINESTRING EMPTY</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        MULTILINESTRING EMPTY
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/mulilinestring/case1</desc>
+  <a>MULTILINESTRING((0 0,0 0),(1 1,2 2))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        GEOMETRYCOLLECTION (LINESTRING (1.0 1.0, 2.0 2.0), POINT (0.0 0.0))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/mulilinestring/case2</desc>
+  <a>MULTILINESTRING((0 0,0 0),(1 1,2 2),(2 2,3 3))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        GEOMETRYCOLLECTION (MULTILINESTRING ((2.0 2.0, 3.0 3.0), (1.0 1.0, 2.0 2.0)), POINT (0.0 0.0))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/mulilinestring/case2</desc>
+  <a>MULTILINESTRING((0 0,0 0),(1 1,2 2),(2 2,3 3),(4 4,4 4))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        GEOMETRYCOLLECTION (MULTILINESTRING ((2.0 2.0, 3.0 3.0), (1.0 1.0, 2.0 2.0)), MULTIPOINT (4.0 4.0, 0.0 0.0))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/polygon/already_valid</desc>
+  <a>POLYGON((0 0,0 1,1 1,0 0))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        POLYGON((0 0,0 1,1 1,0 0))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/polygon/crossing</desc>
+  <a>POLYGON((0 0,1 1,0 1,1 0,0 0))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        MULTIPOLYGON (((0.0 1.0, 1.0 1.0, 0.5 0.5, 0.0 1.0)), ((0.0 0.0, 0.5 0.5, 1.0 0.0, 0.0 0.0)))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>makevalid/polygon/hole_touching_two_places</desc>
+  <a>POLYGON((0 0,0 1,1 1,1 0,0 0),(0 0.5,0.5 0.1,1 0.5,0 0.5))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        MULTIPOLYGON (((0.0 0.0, 0.0 0.5, 0.5 0.1, 1.0 0.5, 1.0 0.0, 0.0 0.0)), ((0.0 0.5, 0.0 1.0, 1.0 1.0, 1.0 0.5, 0.0 0.5)))
+    </op>
+  </test>
+</case>
+
+
+<!-- not completely sure I would expect this result -->
+<case>
+  <desc>makevalid/multipolygon/second_part_overlapping</desc>
+  <a>MULTIPOLYGON(((0 0,0 1,1 1,1 0,0 0)),((0.8 0.1,2 0.1,2 0.9,0.8 0.9,0.8 0.1)))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        MULTIPOLYGON (((0.0 0.0, 0.0 1.0, 1.0 1.0, 1.0 0.9, 0.8 0.9, 0.8 0.1, 1.0 0.1, 1.0 0.0, 0.0 0.0)), ((1.0 0.1, 1.0 0.9, 2.0 0.9, 2.0 0.1, 1.0 0.1)))
+    </op>
+  </test>
+</case>
+
+
+<case>
+  <desc>makevalid/multipolygon/first_part_crossing_second_part_overlapping</desc>
+  <a>MULTIPOLYGON(((0 0,1 1,0 1,1 0,0 0)),((0.8 0.1,2 0.1,2 0.9,0.8 0.9,0.8 0.1)))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        MULTIPOLYGON (((0.8 0.2, 0.8 0.8, 0.9 0.9, 2.0 0.9, 2.0 0.1, 0.9 0.1, 0.8 0.2)), ((0.0 1.0, 1.0 1.0, 0.9 0.9, 0.8 0.9, 0.8 0.8, 0.5 0.5, 0.0 1.0)), ((0.0 0.0, 0.5 0.5, 0.8 0.2, 0.8 0.1, 0.9 0.1, 1.0 0.0, 0.0 0.0)))
+    </op>
+  </test>
+</case>
+
+<case>
+    <desc>makevalid/geometry_collection</desc>
+  <a>GEOMETRYCOLLECTION(POINT EMPTY,LINESTRING EMPTY,POLYGON((0 0,0 1,1 1,1 0,0 0),(0 0.5,0.5 0.1,1 0.5,0 0.5)))</a>
+  <test>
+    <op name="makevalid" arg1="a">
+        GEOMETRYCOLLECTION (MULTIPOLYGON (((0.0 0.0, 0.0 0.5, 0.5 0.1, 1.0 0.5, 1.0 0.0, 0.0 0.0)), ((0.0 0.5, 0.0 1.0, 1.0 1.0, 1.0 0.5, 0.0 0.5))),LINESTRING EMPTY, POINT EMPTY)
+    </op>
+  </test>
+</case>
+
+</run>

commit 589a05f24d08fc755e76270734d397a18b3cbaf2
Author: Even Rouault <even.rouault at spatialys.com>
Date:   Mon Feb 25 15:41:55 2019 +0100

    Add operation::polygonize::BuildArea
    
    And GEOSBuildArea_r() (refs #952)

diff --git a/capi/geos_c.cpp b/capi/geos_c.cpp
index 9f72e40..ce1986f 100644
--- a/capi/geos_c.cpp
+++ b/capi/geos_c.cpp
@@ -725,6 +725,12 @@ extern "C" {
     }
 
     Geometry*
+    GEOSBuildArea(const Geometry* g)
+    {
+        return GEOSBuildArea_r(handle, g);
+    }
+
+    Geometry*
     GEOSLineMerge(const Geometry* g)
     {
         return GEOSLineMerge_r(handle, g);
diff --git a/capi/geos_c.h.in b/capi/geos_c.h.in
index 23b1a6d..d8eb68e 100644
--- a/capi/geos_c.h.in
+++ b/capi/geos_c.h.in
@@ -631,6 +631,10 @@ extern GEOSGeometry GEOS_DLL *GEOSPolygonize_full_r(GEOSContextHandle_t handle,
                               const GEOSGeometry* input, GEOSGeometry** cuts,
                               GEOSGeometry** dangles, GEOSGeometry** invalidRings);
 
+extern GEOSGeometry GEOS_DLL *GEOSBuildArea_r(
+                                              GEOSContextHandle_t handle,
+                                              const GEOSGeometry* g);
+
 extern GEOSGeometry GEOS_DLL *GEOSLineMerge_r(GEOSContextHandle_t handle,
                                               const GEOSGeometry* g);
 extern GEOSGeometry GEOS_DLL *GEOSReverse_r(GEOSContextHandle_t handle,
@@ -1631,6 +1635,8 @@ extern GEOSGeometry GEOS_DLL *GEOSPolygonizer_getCutEdges(const GEOSGeometry * c
 extern GEOSGeometry GEOS_DLL *GEOSPolygonize_full(const GEOSGeometry* input,
     GEOSGeometry** cuts, GEOSGeometry** dangles, GEOSGeometry** invalid);
 
+extern GEOSGeometry GEOS_DLL *GEOSBuildArea(const GEOSGeometry* g);
+
 extern GEOSGeometry GEOS_DLL *GEOSLineMerge(const GEOSGeometry* g);
 extern GEOSGeometry GEOS_DLL *GEOSReverse(const GEOSGeometry* g);
 extern GEOSGeometry GEOS_DLL *GEOSSimplify(const GEOSGeometry* g, double tolerance);
diff --git a/capi/geos_ts_c.cpp b/capi/geos_ts_c.cpp
index ca0625a..2306d75 100644
--- a/capi/geos_ts_c.cpp
+++ b/capi/geos_ts_c.cpp
@@ -63,6 +63,7 @@
 #include <geos/operation/intersection/Rectangle.h>
 #include <geos/operation/intersection/RectangleIntersection.h>
 #include <geos/operation/polygonize/Polygonizer.h>
+#include <geos/operation/polygonize/BuildArea.h>
 #include <geos/operation/relate/RelateOp.h>
 #include <geos/operation/sharedpaths/SharedPathsOp.h>
 #include <geos/operation/union/CascadedPolygonUnion.h>
@@ -3156,6 +3157,36 @@ extern "C" {
     }
 
     Geometry*
+    GEOSBuildArea_r(GEOSContextHandle_t extHandle, const Geometry* g)
+    {
+        if(nullptr == extHandle) {
+            return nullptr;
+        }
+
+        GEOSContextHandleInternal_t* handle = nullptr;
+        handle = reinterpret_cast<GEOSContextHandleInternal_t*>(extHandle);
+        if(0 == handle->initialized) {
+            return nullptr;
+        }
+
+        Geometry* out = nullptr;
+
+        try {
+            using geos::operation::polygonize::BuildArea;
+            BuildArea builder;
+            out = builder.build(g).release();
+        }
+        catch(const std::exception& e) {
+            handle->ERROR_MESSAGE("%s", e.what());
+        }
+        catch(...) {
+            handle->ERROR_MESSAGE("Unknown exception thrown");
+        }
+
+        return out;
+    }
+
+    Geometry*
     GEOSPolygonizer_getCutEdges_r(GEOSContextHandle_t extHandle, const Geometry* const* g, unsigned int ngeoms)
     {
         if(0 == extHandle) {
diff --git a/include/geos/operation/polygonize/BuildArea.h b/include/geos/operation/polygonize/BuildArea.h
new file mode 100644
index 0000000..06c1e83
--- /dev/null
+++ b/include/geos/operation/polygonize/BuildArea.h
@@ -0,0 +1,78 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright 2011-2014 Sandro Santilli <strk at kbt.io>
+ * Copyright (C) 2019 Even Rouault <even.rouault at spatialys.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.
+ **********************************************************************
+ *
+ * Ported from rtgeom_geos.c from
+ *   rttopo - topology library
+ *   http://git.osgeo.org/gitea/rttopo/librttopo
+ * with relicensing from GPL to LGPL with Copyright holder permission.
+ *
+ **********************************************************************/
+
+#ifndef GEOS_OP_POLYGONIZE_BUILDAREA_H
+#define GEOS_OP_POLYGONIZE_BUILDAREA_H
+
+#include <geos/export.h>
+
+#include <memory>
+
+#ifdef _MSC_VER
+#pragma warning(push)
+#pragma warning(disable: 4251) // warning C4251: needs to have dll-interface to be used by clients of class
+#endif
+
+// Forward declarations
+namespace geos {
+namespace geom {
+class Geometry;
+}
+}
+
+namespace geos {
+namespace operation { // geos::operation
+namespace polygonize { // geos::operation::polygonize
+
+/** \brief Creates an areal geometry formed by the constituent linework of given geometry.
+ *
+ * The return type can be a Polygon or MultiPolygon, depending on input.
+ * If the input lineworks do not form polygons NULL is returned.
+ * The inputs can be LINESTRINGS, MULTILINESTRINGS, POLYGONS, MULTIPOLYGONS, and GeometryCollections.
+ *
+ * This algorithm will assume all inner geometries represent holes
+ *
+ * This is the equivalent of the PostGIS ST_BuildArea() function.
+ */
+class GEOS_DLL BuildArea {
+
+public:
+
+    /** \brief
+     * Create a BuildArea object.
+     */
+    BuildArea() = default;
+
+    ~BuildArea() = default;
+
+    /** \brief Return the area built fromthe constituent linework of the input geometry. */
+    std::unique_ptr<geom::Geometry> build(const geom::Geometry* geom);
+};
+
+} // namespace geos::operation::polygonize
+} // namespace geos::operation
+} // namespace geos
+
+#ifdef _MSC_VER
+#pragma warning(pop)
+#endif
+
+#endif // GEOS_OP_POLYGONIZE_BUILDAREA_H
diff --git a/include/geos/operation/polygonize/Makefile.am b/include/geos/operation/polygonize/Makefile.am
index 861a2ef..f9e83ce 100644
--- a/include/geos/operation/polygonize/Makefile.am
+++ b/include/geos/operation/polygonize/Makefile.am
@@ -12,4 +12,5 @@ geos_HEADERS = \
     PolygonizeDirectedEdge.h \
     PolygonizeEdge.h \
     PolygonizeGraph.h \
-    Polygonizer.h
+    Polygonizer.h \
+    BuildArea.h
diff --git a/src/operation/polygonize/BuildArea.cpp b/src/operation/polygonize/BuildArea.cpp
new file mode 100644
index 0000000..5b61e4b
--- /dev/null
+++ b/src/operation/polygonize/BuildArea.cpp
@@ -0,0 +1,246 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright 2011-2014 Sandro Santilli <strk at kbt.io>
+ * Copyright (C) 2019 Even Rouault <even.rouault at spatialys.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.
+ **********************************************************************
+ *
+ * Ported from rtgeom_geos.c from
+ *   rttopo - topology library
+ *   http://git.osgeo.org/gitea/rttopo/librttopo
+ * with relicensing from GPL to LGPL with Copyright holder permission.
+ *
+ **********************************************************************/
+
+#include <geos/operation/polygonize/Polygonizer.h>
+#include <geos/operation/polygonize/BuildArea.h>
+#include <geos/operation/union/CascadedPolygonUnion.h>
+#include <geos/geom/Geometry.h>
+#include <geos/geom/GeometryCollection.h>
+#include <geos/geom/GeometryFactory.h>
+#include <geos/geom/LineString.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/MultiPolygon.h>
+
+//#define DUMP_GEOM
+#ifdef DUMP_GEOM
+#include <geos/io/WKTWriter.h>
+#endif
+
+// std
+#include <algorithm>
+#include <vector>
+
+#ifdef _MSC_VER
+#pragma warning(disable:4355)
+#endif
+
+using namespace std;
+using namespace geos::geom;
+
+namespace geos {
+namespace operation { // geos.operation
+namespace polygonize { // geos.operation.polygonize
+
+struct Face {
+    const geom::Polygon* poly = nullptr;
+    std::unique_ptr<geom::Geometry> env; // envelope
+    double envarea = 0.0; // envelope area
+    Face* parent = nullptr; /* if this face is an hole of another one, or NULL */
+
+    size_t countParents() const {
+        const Face* f = this;
+        size_t pcount = 0;
+        while ( f->parent ) {
+            ++pcount;
+            f = f->parent;
+        }
+        return pcount;
+    }
+};
+
+static std::unique_ptr<Face> newFace(const geom::Polygon* p) {
+    auto f = std::unique_ptr<Face>(new Face());
+    f->poly = p;
+    f->env.reset(p->getEnvelope());
+    f->envarea = f->env->getArea();
+    return f;
+}
+
+struct CompareByEnvarea {
+
+    bool operator()(const std::unique_ptr<Face> &a,
+                    const std::unique_ptr<Face> &b) const {
+        return a->envarea > b->envarea;
+    }
+};
+
+/* Find holes of each face */
+static void findFaceHoles(std::vector<std::unique_ptr<Face>>& faces) {
+
+    /* We sort by decreasing envelope area so that we know holes are only
+     * after their shells */
+    std::sort(faces.begin(), faces.end(), CompareByEnvarea());
+
+    const size_t nfaces = faces.size();
+    for( size_t i = 0; i < nfaces; ++i ) {
+        auto& f = faces[i];
+        const size_t nholes = f->poly->getNumInteriorRing();
+        for( size_t h = 0; h < nholes; h++ ) {
+            const auto hole = f->poly->getInteriorRingN(h);
+            for( auto j=i+1; j < nfaces; ++j ) {
+                auto& f2 = faces[j];
+                if( f2->parent ) {
+                    continue; /* hole already assigned */
+                }
+                const auto f2er = f2->poly->getExteriorRing();
+                /* TODO: can be optimized as the ring would have the
+                *       same vertices, possibly in different order.
+                *       maybe comparing number of points could already be
+                *       useful.
+                */
+                if( f2er->equals(hole) ) {
+                    f2->parent = f.get();
+                    break;
+                }
+            }
+        }
+    }
+}
+
+static std::unique_ptr<geom::MultiPolygon> collectFacesWithEvenAncestors(
+    std::vector<std::unique_ptr<Face>>& faces) {
+    std::vector<geom::Geometry*>* geoms = new std::vector<geom::Geometry*>();
+    for( auto& face: faces ) {
+        if( face->countParents() % 2 ) {
+            continue; /* we skip odd parents geoms */
+        }
+        geoms->push_back(face->poly->clone());
+    }
+    return std::unique_ptr<geom::MultiPolygon>(
+        GeometryFactory::create()->createMultiPolygon(geoms));
+}
+
+#ifdef DUMP_GEOM
+
+static void dumpGeometry(const geom::Geometry* geom)
+{
+  geos::io::WKTWriter oWriter;
+  std::cerr << oWriter.write(geom) << std::endl;
+}
+#endif
+
+/** Return the area built from the constituent linework of the input geometry. */
+unique_ptr<geom::Geometry> BuildArea::build(const geom::Geometry* geom) {
+    Polygonizer polygonizer;
+    polygonizer.add(geom);
+    std::vector<geom::Polygon*>* polys = polygonizer.getPolygons();
+    if( !polys )
+        return nullptr;
+
+    // No geometries in collection, early out
+    if( polys->empty() ) {
+        auto emptyGeomCollection = unique_ptr<geom::Geometry>(
+            GeometryFactory::create()->createGeometryCollection());
+        emptyGeomCollection->setSRID(geom->getSRID());
+        delete polys;
+        return emptyGeomCollection;
+    }
+
+    // Return first geometry if we only have one in collection
+    if( polys->size() == 1 ) {
+        auto ret = unique_ptr<geom::Geometry>((*polys)[0]);
+        ret->setSRID(geom->getSRID());
+        delete polys;
+        return ret;
+    }
+
+    /*
+    * Polygonizer returns a polygon for each face in the built topology.
+    *
+    * This means that for any face with holes we'll have other faces
+    * representing each hole. We can imagine a parent-child relationship
+    * between these faces.
+    *
+    * In order to maximize the number of visible rings in output we
+    * only use those faces which have an even number of parents.
+    *
+    * Example:
+    *
+    *   +---------------+
+    *   |     L0        |  L0 has no parents
+    *   |  +---------+  |
+    *   |  |   L1    |  |  L1 is an hole of L0
+    *   |  |  +---+  |  |
+    *   |  |  |L2 |  |  |  L2 is an hole of L1 (which is an hole of L0)
+    *   |  |  |   |  |  |
+    *   |  |  +---+  |  |
+    *   |  +---------+  |
+    *   |               |
+    *   +---------------+
+    *
+    * See http://trac.osgeo.org/postgis/ticket/1806
+    *
+    */
+
+    try {
+        /* Prepare face structures for later analysis */
+        std::vector<std::unique_ptr<Face>> faces;
+        for(auto poly: *polys) {
+            faces.emplace_back(newFace(poly));
+        }
+
+        /* Find faces representing other faces holes */
+        findFaceHoles(faces);
+
+        /* Build a MultiPolygon composed only by faces with an
+        * even number of ancestors */
+        auto tmp = collectFacesWithEvenAncestors(faces);
+
+        for(auto poly: *polys) {
+            delete poly;
+        }
+        delete polys;
+        polys = nullptr;
+
+#ifdef DUMP_GEOM
+        std::cerr << "after collectFacesWithEvenAncestors:" << std::endl;
+        dumpGeometry(tmp.get());
+#endif
+
+        /* Run a single overlay operation to dissolve shared edges */
+        auto shp = std::unique_ptr<geom::Geometry>(
+            geos::operation::geounion::CascadedPolygonUnion::CascadedPolygonUnion::Union(tmp.get()));
+        if( shp ) {
+            shp->setSRID(geom->getSRID());
+        }
+
+#ifdef DUMP_GEOM
+        std::cerr << "after CascadedPolygonUnion:" << std::endl;
+        dumpGeometry(shp.get());
+#endif
+
+        return shp;
+    }
+    catch( ... ) {
+        if( polys ) {
+            for(auto poly: *polys) {
+                delete poly;
+            }
+            delete polys;
+        }
+        throw;
+    }
+}
+
+} // namespace geos.operation.polygonize
+} // namespace geos.operation
+} // namespace geos
+
diff --git a/src/operation/polygonize/Makefile.am b/src/operation/polygonize/Makefile.am
index f670b0c..0e8d417 100644
--- a/src/operation/polygonize/Makefile.am
+++ b/src/operation/polygonize/Makefile.am
@@ -12,6 +12,7 @@ liboppolygonize_la_SOURCES = \
     PolygonizeEdge.cpp \
     PolygonizeGraph.cpp \
     Polygonizer.cpp \
+    BuildArea.cpp \
     EdgeRing.cpp 
 
 liboppolygonize_la_LIBADD = 
diff --git a/tests/unit/Makefile.am b/tests/unit/Makefile.am
index 186dd34..da53fb9 100644
--- a/tests/unit/Makefile.am
+++ b/tests/unit/Makefile.am
@@ -176,7 +176,8 @@ geos_unit_SOURCES = \
 	capi/GEOSUnaryUnionTest.cpp \
 	capi/GEOSisValidDetailTest.cpp \
 	capi/GEOSisClosedTest.cpp \
-	capi/GEOSInterpolateTest.cpp
+	capi/GEOSInterpolateTest.cpp \
+	capi/GEOSBuildArea.cpp
 
 noinst_HEADERS = \
 	utility.h
diff --git a/tests/unit/capi/GEOSBuildArea.cpp b/tests/unit/capi/GEOSBuildArea.cpp
new file mode 100644
index 0000000..c1350d9
--- /dev/null
+++ b/tests/unit/capi/GEOSBuildArea.cpp
@@ -0,0 +1,81 @@
+// Test Suite for C-API BuildArea
+
+#include <tut/tut.hpp>
+// geos
+#include <geos_c.h>
+// std
+#include <cstdarg>
+#include <cstdio>
+#include <cstdlib>
+#include <cmath>
+
+namespace tut {
+//
+// Test Group
+//
+
+// Common data used in test cases.
+struct test_capi_buildarea_data {
+    GEOSWKTWriter* wktw_;
+    GEOSGeometry* geom1_ = nullptr;
+    GEOSGeometry* geom2_ = nullptr;
+
+    static void
+    notice(const char* fmt, ...)
+    {
+        std::fprintf(stdout, "NOTICE: ");
+
+        va_list ap;
+        va_start(ap, fmt);
+        std::vfprintf(stdout, fmt, ap);
+        va_end(ap);
+
+        std::fprintf(stdout, "\n");
+    }
+
+    std::string
+    toWKT(GEOSGeometry* g)
+    {
+        char* wkt = GEOSWKTWriter_write(wktw_, g);
+        std::string ret(wkt);
+        GEOSFree(wkt);
+        return ret;
+    }
+
+    test_capi_buildarea_data()
+    {
+        initGEOS(notice, notice);
+        wktw_ = GEOSWKTWriter_create();
+        GEOSWKTWriter_setTrim(wktw_, 1);
+        GEOSWKTWriter_setOutputDimension(wktw_, 3);
+    }
+
+    ~test_capi_buildarea_data()
+    {
+        GEOSGeom_destroy(geom1_);
+        GEOSGeom_destroy(geom2_);
+        GEOSWKTWriter_destroy(wktw_);
+        finishGEOS();
+    }
+
+};
+
+typedef test_group<test_capi_buildarea_data> group;
+typedef group::object object;
+
+group test_capi_buildarea_group("capi::GEOSBuildArea");
+
+//
+// Test Cases
+//
+
+template<>
+template<>
+void object::test<1>
+()
+{
+    geom1_ = GEOSGeomFromWKT("GEOMETRYCOLLECTION(LINESTRING(0 0,0 1,1 1),LINESTRING(1 1,1 0,0 0))");
+    geom2_ = GEOSBuildArea(geom1_);
+    ensure_equals(toWKT(geom2_), std::string("POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0))"));
+}
+} // namespace tut
diff --git a/tests/xmltester/CMakeLists.txt b/tests/xmltester/CMakeLists.txt
index 1d780e7..fd492ce 100644
--- a/tests/xmltester/CMakeLists.txt
+++ b/tests/xmltester/CMakeLists.txt
@@ -62,6 +62,7 @@ if(GEOS_ENABLE_TESTS)
     ${XMLTESTS_DIR}/split.xml
     ${XMLTESTS_DIR}/test.xml
     ${XMLTESTS_DIR}/testLeaksBig.xml
+    ${XMLTESTS_DIR}/buildarea.xml
     ${XMLTESTS_DIR}/general/TestBoundary.xml
     ${XMLTESTS_DIR}/general/TestBuffer.xml
     ${XMLTESTS_DIR}/general/TestBufferMitredJoin.xml
diff --git a/tests/xmltester/Makefile.am b/tests/xmltester/Makefile.am
index 87b2b4f..819fccc 100644
--- a/tests/xmltester/Makefile.am
+++ b/tests/xmltester/Makefile.am
@@ -30,6 +30,7 @@ SAFE_XMLTESTS= \
 	$(srcdir)/tests/split.xml \
 	$(srcdir)/tests/test.xml \
 	$(srcdir)/tests/testLeaksBig.xml \
+	$(srcdir)/tests/buildarea.xml \
 	$(srcdir)/tests/general/TestBoundary.xml \
 	$(srcdir)/tests/general/TestBuffer.xml \
 	$(srcdir)/tests/general/TestBufferMitredJoin.xml \
diff --git a/tests/xmltester/XMLTester.cpp b/tests/xmltester/XMLTester.cpp
index 085d1e6..e996002 100644
--- a/tests/xmltester/XMLTester.cpp
+++ b/tests/xmltester/XMLTester.cpp
@@ -41,6 +41,7 @@
 #include <geos/operation/buffer/BufferBuilder.h>
 #include <geos/operation/buffer/BufferParameters.h>
 #include <geos/operation/buffer/BufferOp.h>
+#include <geos/operation/polygonize/BuildArea.h>
 #include <geos/precision/MinimumClearance.h>
 #include <geos/util.h>
 #include <geos/util/Interrupt.h>
@@ -1747,6 +1748,32 @@ XMLTester::parseTest(const tinyxml2::XMLNode* node)
             success = lineE.get()->equalsExact(lineO.get(), tol) ? 1 : 0;
         }
 
+        else if (opName == "buildarea")
+        {
+            GeomPtr gExpected(parseGeometry(opRes, "expected"));
+            gExpected->normalize();
+
+            auto gGot = BuildArea().build(gA);
+            if( gGot )
+            {
+                GeomPtr gRealRes(gGot.release());
+                gRealRes->normalize();
+
+                if (gExpected->equals(gRealRes.get())) success=1;
+
+                actual_result=printGeom(gRealRes.get());
+                expected_result=printGeom(gExpected.get());
+                if( actual_result == expected_result ) success=1;
+
+                if ( testValidOutput )
+                    success &= int(testValid(gRealRes.get(), "result"));
+            }
+            else
+            {
+                success = false;
+            }
+        }
+
         else {
             std::cerr << *curr_file << ":";
             std::cerr << " case" << caseCount << ":";
diff --git a/tests/xmltester/tests/buildarea.xml b/tests/xmltester/tests/buildarea.xml
new file mode 100644
index 0000000..fd9af75
--- /dev/null
+++ b/tests/xmltester/tests/buildarea.xml
@@ -0,0 +1,98 @@
+<run>
+  <precisionModel scale="100" offsetx="0.0" offsety="0.0"/>
+
+<case>
+  <desc>buildarea/point</desc>
+  <a>POINT(0 0)</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        GEOMETRYCOLLECTION EMPTY
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/linestring</desc>
+  <a>LINESTRING(0 0,1 1)</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        GEOMETRYCOLLECTION EMPTY
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/polygon</desc>
+  <a>POLYGON((0 0,0 1,1 1,0 0))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        POLYGON((0 0,0 1,1 1,0 0))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/polygon_from_linestrings</desc>
+  <a>GEOMETRYCOLLECTION(LINESTRING(0 0,0 1,1 1),LINESTRING(1 1,1 0,0 0))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        POLYGON((0 0,0 1,1 1,1 0,0 0))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/multipolygon_of_two_polygons_from_linestrings</desc>
+  <a>GEOMETRYCOLLECTION(LINESTRING(0 0,0 1,1 1),LINESTRING(1 1,1 0,0 0),LINESTRING(10 0,10 1,11 1),LINESTRING(11 1,11 0,10 0))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        MULTIPOLYGON (((10.0 0.0, 10.0 1.0, 11.0 1.0, 11.0 0.0, 10.0 0.0)), ((0.0 0.0, 0.0 1.0, 1.0 1.0, 1.0 0.0, 0.0 0.0)))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/polygon_with_hole_from_linestrings</desc>
+  <a>GEOMETRYCOLLECTION(LINESTRING(0 0,0 1,1 1),LINESTRING(1 1,1 0,0 0),LINESTRING(0.25 0.25,0.25 0.75,0.75 0.75),LINESTRING(0.75 0.75,0.75 0.25,0.25 0.25))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        POLYGON ((0.00 0.00, 0.00 1.00, 1.00 1.00, 1.00 0.00, 0.00 0.00), (0.25 0.25, 0.75 0.25, 0.75 0.75, 0.25 0.75, 0.25 0.25))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/polygon_with_hole_and_island_from_linestrings</desc>
+  <a>GEOMETRYCOLLECTION(LINESTRING(0 0,0 1,1 1),LINESTRING(1 1,1 0,0 0),LINESTRING(0.25 0.25,0.25 0.75,0.75 0.75),LINESTRING(0.75 0.75,0.75 0.25,0.25 0.25),LINESTRING(0.3 0.3,0.3 0.7,0.7 0.7),LINESTRING(0.7 0.7,0.7 0.3,0.3 0.3))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        MULTIPOLYGON (((0.00 0.00, 0.00 1.00, 1.00 1.00, 1.00 0.00, 0.00 0.00), (0.25 0.25, 0.75 0.25, 0.75 0.75, 0.25 0.75, 0.25 0.25)),((0.30 0.30, 0.30 0.70, 0.70 0.70, 0.70 0.30, 0.30 0.30)))
+    </op>
+  </test>
+</case>
+
+<case>
+  <desc>buildarea/self_touching_multipolygons</desc>
+  <a>GEOMETRYCOLLECTION (LINESTRING (100 100, 100 400, 400 400, 400 100, 100 100), LINESTRING (150 350, 350 350, 350 150, 100 100, 150 350), LINESTRING (200 250, 250 200, 100 100, 200 250), LINESTRING (200 250, 160 290, 150 250, 200 250), LINESTRING (350 350, 200 300, 200 250, 300 250, 350 350))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        MULTIPOLYGON (((200 250, 250 200, 100 100, 200 250)), ((100 100, 100 400, 400 400, 400 100, 100 100), (150 350, 100 100, 350 150, 350 350, 150 350)), ((200 250, 150 250, 160 290, 200 250)), ((350 350, 300 250, 200 250, 200 300, 350 350)))
+    </op>
+  </test>
+</case>
+
+<!-- Test commented out. See https://github.com/libgeos/geos/pull/151#discussion_r260897618 -->
+<!-- The current result is POLYGON ((0.00 0.00, 0.00 300.00, 50.00 300.00, 50.00 450.00, 200.00 450.00, 200.00 500.00, 500.00 500.00, 500.00 200.00, 450.00 200.00, 450.00 50.00, 300.00 50.00, 300.00 0.00, 0.00 0.00)) -->
+<!--
+<case>
+  <desc>buildarea/checkerboard</desc>
+  <a>MULTILINESTRING ((50 300, 100 300), (50 300, 50 450, 200 450), (300 50, 300 0, 0 0, 0 300, 50 300), (100 150, 50 150, 50 300), (150 150, 100 150), (100 150, 100 300), (150 100, 100 100, 100 150), (100 300, 150 300), (100 300, 100 400, 200 400), (300 100, 150 100), (150 100, 150 150), (300 50, 150 50, 150 100), (300 150, 150 150), (150 150, 150 300), (150 300, 200 300), (150 300, 150 350, 200 350), (200 300, 300 300, 300 200), (200 300, 200 350), (300 200, 200 200, 200 300), (200 350, 350 350), (200 350, 200 400), (200 400, 350 400), (200 400, 200 450), (200 450, 350 450, 350 400), (200 450, 200 500, 500 500, 500 200, 450 200), (450 200, 450 50, 300 50), (300 100, 300 50), (400 200, 400 100, 300 100), (300 150, 300 100), (350 200, 350 150, 300 150), (300 200, 300 150), (350 200, 300 200), (400 200, 350 200), (350 350, 350 200), (350 350, 400 350), (350 400, 350 350), (350 400, 400 400, 400 350), (450 200, 400 200), (400 350, 400 200), (400 350, 450 350, 450 200))</a>
+  <test>
+    <op name="buildarea" arg1="a">
+        MULTIPOLYGON (((150 100, 150 50, 300 50, 300 0, 0 0, 0 300, 50 300, 50 150, 100 150, 100 100, 150 100)), ((200 400, 100 400, 100 300, 50 300, 50 450, 200 450, 200 400)), ((100 150, 100 300, 150 300, 150 150, 100 150)), ((150 150, 300 150, 300 100, 150 100, 150 150)), ((200 350, 200 300, 150 300, 150 350, 200 350)), ((300 200, 200 200, 200 300, 300 300, 300 200)), ((350 400, 350 450, 200 450, 200 500, 500 500, 500 200, 450 200, 450 350, 400 350, 400 400, 350 400)), ((200 350, 200 400, 350 400, 350 350, 200 350)), ((450 200, 450 50, 300 50, 300 100, 400 100, 400 200, 450 200)), ((300 150, 300 200, 350 200, 350 150, 300 150)), ((400 350, 400 200, 350 200, 350 350, 400 350)))
+    </op>
+  </test>
+</case>
+-->
+
+</run>

-----------------------------------------------------------------------

Summary of changes:
 NEWS                                          |   5 +
 capi/geos_c.cpp                               |  12 +
 capi/geos_c.h.in                              |  11 +
 capi/geos_ts_c.cpp                            |  63 +++++
 include/geos/operation/polygonize/BuildArea.h |  78 ++++++
 include/geos/operation/polygonize/Makefile.am |   3 +-
 include/geos/operation/valid/MakeValid.h      |  79 ++++++
 src/operation/polygonize/BuildArea.cpp        | 246 ++++++++++++++++++
 src/operation/polygonize/Makefile.am          |   1 +
 src/operation/valid/MakeValid.cpp             | 344 ++++++++++++++++++++++++++
 src/operation/valid/Makefile.am               |   3 +-
 tests/unit/Makefile.am                        |   4 +-
 tests/unit/capi/GEOSBuildArea.cpp             |  81 ++++++
 tests/unit/capi/GEOSMakeValid.cpp             |  81 ++++++
 tests/xmltester/CMakeLists.txt                |   2 +
 tests/xmltester/Makefile.am                   |   2 +
 tests/xmltester/XMLTester.cpp                 |  54 ++++
 tests/xmltester/tests/buildarea.xml           |  98 ++++++++
 tests/xmltester/tests/makevalid.xml           | 157 ++++++++++++
 19 files changed, 1321 insertions(+), 3 deletions(-)
 create mode 100644 include/geos/operation/polygonize/BuildArea.h
 create mode 100644 include/geos/operation/valid/MakeValid.h
 create mode 100644 src/operation/polygonize/BuildArea.cpp
 create mode 100644 src/operation/valid/MakeValid.cpp
 create mode 100644 tests/unit/capi/GEOSBuildArea.cpp
 create mode 100644 tests/unit/capi/GEOSMakeValid.cpp
 create mode 100644 tests/xmltester/tests/buildarea.xml
 create mode 100644 tests/xmltester/tests/makevalid.xml


hooks/post-receive
-- 
GEOS


More information about the geos-commits mailing list