diff --git a/doc/source/docs/reference/geoseries.rst b/doc/source/docs/reference/geoseries.rst index 0173a8c..1965d11 100644 --- a/doc/source/docs/reference/geoseries.rst +++ b/doc/source/docs/reference/geoseries.rst @@ -69,6 +69,7 @@ Set-theoretic Methods .. autosummary:: :toctree: api/ + GeoSeries.clip_by_rect GeoSeries.difference GeoSeries.intersection GeoSeries.symmetric_difference diff --git a/geopandas/_vectorized.py b/geopandas/_vectorized.py index 457a81f..159e578 100644 --- a/geopandas/_vectorized.py +++ b/geopandas/_vectorized.py @@ -11,6 +11,7 @@ import pandas as pd import shapely.geometry import shapely.geos +import shapely.ops import shapely.wkb import shapely.wkt @@ -728,6 +729,20 @@ def almost_equals(self, other, decimal): # +def clip_by_rect(data, xmin, ymin, xmax, ymax): + if compat.USE_PYGEOS: + return pygeos.clip_by_rect(data, xmin, ymin, xmax, ymax) + else: + clipped_geometries = np.empty(len(data), dtype=object) + clipped_geometries[:] = [ + shapely.ops.clip_by_rect(s, xmin, ymin, xmax, ymax) + if s is not None + else None + for s in data + ] + return clipped_geometries + + def difference(data, other): if compat.USE_PYGEOS: return _binary_method("difference", data, other) diff --git a/geopandas/array.py b/geopandas/array.py index 5f3f24b..9b47d6a 100644 --- a/geopandas/array.py +++ b/geopandas/array.py @@ -564,6 +564,11 @@ class GeometryArray(ExtensionArray): # Binary operations that return new geometries # + def clip_by_rect(self, xmin, ymin, xmax, ymax): + return GeometryArray( + vectorized.clip_by_rect(self.data, xmin, ymin, xmax, ymax), crs=self.crs + ) + def difference(self, other): return GeometryArray( self._binary_method("difference", self, other), crs=self.crs diff --git a/geopandas/base.py b/geopandas/base.py index 2edf647..dab414b 100644 --- a/geopandas/base.py +++ b/geopandas/base.py @@ -2536,6 +2536,75 @@ GeometryCollection """ return _binary_geo("intersection", self, other, align) + def clip_by_rect(self, xmin, ymin, xmax, ymax): + """Returns a ``GeoSeries`` of the portions of geometry within the given + rectangle. + + Note that the results are not exactly equal to + :meth:`~GeoSeries.intersection()`. E.g. in edge cases, + :meth:`~GeoSeries.clip_by_rect()` will not return a point just touching the + rectangle. Check the examples section below for some of these exceptions. + + The geometry is clipped in a fast but possibly dirty way. The output is not + guaranteed to be valid. No exceptions will be raised for topological errors. + + Note: empty geometries or geometries that do not overlap with the specified + bounds will result in ``GEOMETRYCOLLECTION EMPTY``. + + Parameters + ---------- + xmin: float + Minimum x value of the rectangle + ymin: float + Minimum y value of the rectangle + xmax: float + Maximum x value of the rectangle + ymax: float + Maximum y value of the rectangle + + Returns + ------- + GeoSeries + + Examples + -------- + >>> from shapely.geometry import Polygon, LineString, Point + >>> s = geopandas.GeoSeries( + ... [ + ... Polygon([(0, 0), (2, 2), (0, 2)]), + ... Polygon([(0, 0), (2, 2), (0, 2)]), + ... LineString([(0, 0), (2, 2)]), + ... LineString([(2, 0), (0, 2)]), + ... Point(0, 1), + ... ], + ... crs=3857, + ... ) + >>> bounds = (0, 0, 1, 1) + >>> s + 0 POLYGON ((0.000 0.000, 2.000 2.000, 0.000 2.00... + 1 POLYGON ((0.000 0.000, 2.000 2.000, 0.000 2.00... + 2 LINESTRING (0.000 0.000, 2.000 2.000) + 3 LINESTRING (2.000 0.000, 0.000 2.000) + 4 POINT (0.000 1.000) + dtype: geometry + >>> s.clip_by_rect(*bounds) + 0 POLYGON ((0.000 0.000, 0.000 1.000, 1.000 1.00... + 1 POLYGON ((0.000 0.000, 0.000 1.000, 1.000 1.00... + 2 LINESTRING (0.000 0.000, 1.000 1.000) + 3 GEOMETRYCOLLECTION EMPTY + 4 GEOMETRYCOLLECTION EMPTY + dtype: geometry + + See also + -------- + GeoSeries.intersection + """ + from .geoseries import GeoSeries + + geometry_array = GeometryArray(self.geometry.values) + clipped_geometry = geometry_array.clip_by_rect(xmin, ymin, xmax, ymax) + return GeoSeries(clipped_geometry.data, index=self.index, crs=self.crs) + # # Other operations # diff --git a/geopandas/tests/test_geom_methods.py b/geopandas/tests/test_geom_methods.py index 05475a9..05e3d5c 100644 --- a/geopandas/tests/test_geom_methods.py +++ b/geopandas/tests/test_geom_methods.py @@ -85,8 +85,12 @@ class TestGeomMethods: self.g8 = GeoSeries([self.t1, self.t5]) self.empty = GeoSeries([]) self.all_none = GeoSeries([None, None]) + self.all_geometry_collection_empty = GeoSeries( + [GeometryCollection([]), GeometryCollection([])] + ) self.empty_poly = Polygon() self.g9 = GeoSeries(self.g0, index=range(1, 8)) + self.g10 = GeoSeries([self.t1, self.t4]) # Crossed lines self.l3 = LineString([(0, 0), (1, 1)]) @@ -254,6 +258,15 @@ class TestGeomMethods: assert len(self.g0.intersection(self.g9, align=True) == 8) assert len(self.g0.intersection(self.g9, align=False) == 7) + def test_clip_by_rect(self): + self._test_binary_topological( + "clip_by_rect", self.g1, self.g10, *self.sq.bounds + ) + # self.g1 and self.t3.bounds do not intersect + self._test_binary_topological( + "clip_by_rect", self.all_geometry_collection_empty, self.g1, *self.t3.bounds + ) + def test_union_series(self): self._test_binary_topological("union", self.sq, self.g1, self.g2)