From cb310dff59d4e8b069706e637e7318fac80bcc41 Mon Sep 17 00:00:00 2001 From: Simone Parmeggiani <70708372+SimoParmeg@users.noreply.github.com> Date: Sat, 22 Oct 2022 15:49:15 +0200 Subject: [PATCH] ENH: Added make_valid method from shapely (#2539) --- CHANGELOG.md | 1 + geopandas/_vectorized.py | 21 ++++++++++++ geopandas/array.py | 3 ++ geopandas/base.py | 39 ++++++++++++++++++++++ geopandas/tests/test_geom_methods.py | 48 +++++++++++++++++++++++++++- 5 files changed, 111 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 611c067..3e6b63f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -18,6 +18,7 @@ for more details. New features and improvements: - Added ``normalize()`` method from shapely to GeoSeries/GeoDataframe (#2537) +- Added ``make_valid()`` method from shapely to GeoSeries/GeoDataframe (#2539) - Add ``where`` filter to ``read_file`` (#2552) Deprecations and compatibility notes: diff --git a/geopandas/_vectorized.py b/geopandas/_vectorized.py index e5736f2..45d383d 100644 --- a/geopandas/_vectorized.py +++ b/geopandas/_vectorized.py @@ -15,6 +15,7 @@ import shapely.geos import shapely.ops import shapely.wkb import shapely.wkt +import shapely.validation from shapely.geometry.base import BaseGeometry @@ -975,6 +976,26 @@ def normalize(data): return out +def make_valid(data): + if compat.USE_SHAPELY_20: + return shapely.make_valid(data) + elif compat.USE_PYGEOS: + return pygeos.make_valid(data) + elif not compat.SHAPELY_GE_18: + raise NotImplementedError( + f"shapely >= 1.8 or PyGEOS is required, " + f"version {shapely.__version__} is installed" + ) + else: + out = np.empty(len(data), dtype=object) + with compat.ignore_shapely2_warnings(): + out[:] = [ + shapely.validation.make_valid(geom) if geom is not None else None + for geom in data + ] + return out + + def project(data, other, normalized=False): if compat.USE_SHAPELY_20: return shapely.line_locate_point(data, other, normalized=normalized) diff --git a/geopandas/array.py b/geopandas/array.py index ee56019..779e729 100644 --- a/geopandas/array.py +++ b/geopandas/array.py @@ -514,6 +514,9 @@ class GeometryArray(ExtensionArray): def normalize(self): return GeometryArray(vectorized.normalize(self.data), crs=self.crs) + def make_valid(self): + return GeometryArray(vectorized.make_valid(self.data), crs=self.crs) + # # Binary predicates # diff --git a/geopandas/base.py b/geopandas/base.py index 437a611..8302b32 100644 --- a/geopandas/base.py +++ b/geopandas/base.py @@ -725,6 +725,45 @@ GeometryCollection """ return _delegate_geo_method("normalize", self) + def make_valid(self): + """ + Repairs invalid geometries. + + Returns a ``GeoSeries`` with valid geometries. + If the input geometry is already valid, then it will be preserved. + In many cases, in order to create a valid geometry, the input + geometry must be split into multiple parts or multiple geometries. + If the geometry must be split into multiple parts of the same type + to be made valid, then a multi-part geometry will be returned + (e.g. a MultiPolygon). + If the geometry must be split into multiple parts of different types + to be made valid, then a GeometryCollection will be returned. + + Examples + -------- + >>> from shapely.geometry import MultiPolygon, Polygon, LineString, Point + >>> s = geopandas.GeoSeries( + ... [ + ... Polygon([(0, 0), (0, 2), (1, 1), (2, 2), (2, 0), (1, 1), (0, 0)]), + ... Polygon([(0, 2), (0, 1), (2, 0), (0, 0), (0, 2)]), + ... LineString([(0, 0), (1, 1), (1, 0)]), + ... ], + ... crs='EPSG:3857', + ... ) + >>> s + 0 POLYGON ((0.000 0.000, 0.000 2.000, 1.000 1.00... + 1 POLYGON ((0.000 2.000, 0.000 1.000, 2.000 0.00... + 2 LINESTRING (0.000 0.000, 1.000 1.000, 1.000 0.... + dtype: geometry + + >>> s.make_valid() + 0 MULTIPOLYGON (((1.000 1.000, 0.000 0.000, 0.00... + 1 GEOMETRYCOLLECTION (POLYGON ((2.000 0.000, 0.0... + 2 LINESTRING (0.000 0.000, 1.000 1.000, 1.000 0.... + dtype: geometry + """ + return _delegate_geo_method("make_valid", self) + # # Reduction operations that return a Shapely geometry # diff --git a/geopandas/tests/test_geom_methods.py b/geopandas/tests/test_geom_methods.py index d29932d..606102e 100644 --- a/geopandas/tests/test_geom_methods.py +++ b/geopandas/tests/test_geom_methods.py @@ -5,7 +5,16 @@ import numpy as np from numpy.testing import assert_array_equal from pandas import DataFrame, Index, MultiIndex, Series -from shapely.geometry import LinearRing, LineString, MultiPoint, Point, Polygon +import shapely + +from shapely.geometry import ( + LinearRing, + LineString, + MultiPoint, + Point, + Polygon, + MultiPolygon, +) from shapely.geometry.collection import GeometryCollection from shapely.ops import unary_union from shapely import wkt @@ -682,6 +691,43 @@ class TestGeomMethods: expected = GeoSeries([polygon2, linestring, point]) assert_geoseries_equal(series.normalize(), expected) + @pytest.mark.skipif( + not compat.SHAPELY_GE_18, + reason="make_valid keyword introduced in shapely 1.8.0", + ) + def test_make_valid(self): + polygon1 = Polygon([(0, 0), (0, 2), (1, 1), (2, 2), (2, 0), (1, 1), (0, 0)]) + polygon2 = Polygon([(0, 2), (0, 1), (2, 0), (0, 0), (0, 2)]) + linestring = LineString([(0, 0), (1, 1), (1, 0)]) + series = GeoSeries([polygon1, polygon2, linestring]) + out_polygon1 = MultiPolygon( + [ + Polygon([(1, 1), (0, 0), (0, 2), (1, 1)]), + Polygon([(2, 0), (1, 1), (2, 2), (2, 0)]), + ] + ) + out_polygon2 = GeometryCollection( + [Polygon([(2, 0), (0, 0), (0, 1), (2, 0)]), LineString([(0, 2), (0, 1)])] + ) + expected = GeoSeries([out_polygon1, out_polygon2, linestring]) + assert not series.is_valid.all() + result = series.make_valid() + assert_geoseries_equal(result, expected) + assert result.is_valid.all() + + @pytest.mark.skipif( + compat.SHAPELY_GE_18, + reason="make_valid keyword introduced in shapely 1.8.0", + ) + def test_make_valid_shapely_pre18(self): + s = GeoSeries([Point(1, 1)]) + with pytest.raises( + NotImplementedError, + match=f"shapely >= 1.8 or PyGEOS is required, " + f"version {shapely.__version__} is installed", + ): + s.make_valid() + def test_convex_hull(self): # the convex hull of a square should be the same as the square squares = GeoSeries([self.sq for i in range(3)])