Repository: geopandas/geopandas (python)

Issue:
BUG: Unexpected behavior in Overlay with keep_geom_type since Version 1.0.0
- [x] I have checked that this issue has not already been reported.

- [x] I have confirmed this bug exists on the latest version of geopandas.

- [x] (optional) I have confirmed this bug exists on the main branch of geopandas.

---

#### Code Sample

```python
import geopandas
from shapely import Polygon

invalid_polygon = Polygon([(0, 0), (1, 1), (1, 2), (1, 1), (0, 0)])
square = Polygon(((0,0),(1,0),(1,1),(0,1),(0,0)))
df1 = geopandas.GeoDataFrame({'A':[1,2]}, geometry=[invalid_polygon, square])
df2 = geopandas.GeoDataFrame({'B':[5,6]}, geometry=[square, square])

print(df1.geom_type[0])
print(df2.geom_type[0])
df1_df2 = df1.overlay(df2, keep_geom_type=True)
print(df1_df2.geom_type[0])

```

#### Problem description
Overlay switches the geometry type of the resulting GeoDataFrame entirely. The argument 'keep_geom_type=True' should keep the geometry type of the Geodataframe as Polygon, but in the example above, it changes it completely to Linestring. This Problem occurs since version 1.0.0, when overlay switched from buffer(0) to make_valid (make_valid does change the geometry type of the invalid Polygon, but buffer(0) does not). This results the first geometry of df1 to be of type Linestring, which seems to cause the behavior.

#### Expected Output
df1_df2 should have same geometry type as df1.

Pull request description:
BUG: Fix regression where `overlay` with `keep_geom_type` returns wrong results if the input contains invalid geometries
@m-richards 

fixes #3389 

Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 53933780..66d9621b 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -21,6 +21,8 @@ Bug fixes:
   (Geo)Series, but instead should be interpreted as column names, like Enums (#3384).
 - Fix regression where constructing a GeoSeries from a pd.Series with GeometryDtype values 
   failed when `crs` was provided (#3383).
+- Fix regression where `overlay` with `keep_geom_type` returns wrong results if the
+  input contains invalid geometries (#3395).
 - Fix the dtype of the GeometryArray backing data being incorrect for zero length 
   GeoDataFrames causing errors in `overlay` (3424).
 
diff --git a/geopandas/array.py b/geopandas/array.py
index 2bb3920d..6bb88ef7 100644
--- a/geopandas/array.py
+++ b/geopandas/array.py
@@ -40,6 +40,9 @@ _names = {
     "GEOMETRYCOLLECTION": "GeometryCollection",
 }
 
+POLYGON_GEOM_TYPES = {"Polygon", "MultiPolygon"}
+LINE_GEOM_TYPES = {"LineString", "MultiLineString", "LinearRing"}
+POINT_GEOM_TYPES = {"Point", "MultiPoint"}
 
 type_mapping = {p.value: _names[p.name] for p in shapely.GeometryType}
 geometry_type_ids = list(type_mapping.keys())
diff --git a/geopandas/tools/_random.py b/geopandas/tools/_random.py
index 007d02fa..210c44e3 100644
--- a/geopandas/tools/_random.py
+++ b/geopandas/tools/_random.py
@@ -4,7 +4,12 @@ import numpy
 
 from shapely.geometry import MultiPoint
 
-from geopandas.array import from_shapely, points_from_xy
+from geopandas.array import (
+    LINE_GEOM_TYPES,
+    POLYGON_GEOM_TYPES,
+    from_shapely,
+    points_from_xy,
+)
 from geopandas.geoseries import GeoSeries
 
 
@@ -45,10 +50,10 @@ def uniform(geom, size, rng=None):
     if geom is None or geom.is_empty:
         return MultiPoint()
 
-    if geom.geom_type in ("Polygon", "MultiPolygon"):
+    if geom.geom_type in POLYGON_GEOM_TYPES:
         return _uniform_polygon(geom, size=size, generator=generator)
 
-    if geom.geom_type in ("LineString", "MultiLineString"):
+    if geom.geom_type in LINE_GEOM_TYPES:
         return _uniform_line(geom, size=size, generator=generator)
 
     warn(
diff --git a/geopandas/tools/clip.py b/geopandas/tools/clip.py
index 641417a5..d1386f30 100644
--- a/geopandas/tools/clip.py
+++ b/geopandas/tools/clip.py
@@ -14,7 +14,13 @@ import pandas.api.types
 from shapely.geometry import MultiPolygon, Polygon, box
 
 from geopandas import GeoDataFrame, GeoSeries
-from geopandas.array import _check_crs, _crs_mismatch_warn
+from geopandas.array import (
+    LINE_GEOM_TYPES,
+    POINT_GEOM_TYPES,
+    POLYGON_GEOM_TYPES,
+    _check_crs,
+    _crs_mismatch_warn,
+)
 
 
 def _mask_is_list_like_rectangle(mask):
@@ -213,25 +219,21 @@ def clip(gdf, mask, keep_geom_type=False, sort=False):
                 stacklevel=2,
             )
         else:
-            polys = ["Polygon", "MultiPolygon"]
-            lines = ["LineString", "MultiLineString", "LinearRing"]
-            points = ["Point", "MultiPoint"]
-
             # Check that the gdf for multiple geom types (points, lines and/or polys)
             orig_types_total = sum(
                 [
-                    gdf.geom_type.isin(polys).any(),
-                    gdf.geom_type.isin(lines).any(),
-                    gdf.geom_type.isin(points).any(),
+                    gdf.geom_type.isin(POLYGON_GEOM_TYPES).any(),
+                    gdf.geom_type.isin(LINE_GEOM_TYPES).any(),
+                    gdf.geom_type.isin(POINT_GEOM_TYPES).any(),
                 ]
             )
 
             # Check how many geometry types are in the clipped GeoDataFrame
             clip_types_total = sum(
                 [
-                    clipped.geom_type.isin(polys).any(),
-                    clipped.geom_type.isin(lines).any(),
-                    clipped.geom_type.isin(points).any(),
+                    clipped.geom_type.isin(POLYGON_GEOM_TYPES).any(),
+                    clipped.geom_type.isin(LINE_GEOM_TYPES).any(),
+                    clipped.geom_type.isin(POINT_GEOM_TYPES).any(),
                 ]
             )
 
@@ -247,9 +249,9 @@ def clip(gdf, mask, keep_geom_type=False, sort=False):
                 orig_type = gdf.geom_type.iloc[0]
                 if new_collection:
                     clipped = clipped.explode(index_parts=False)
-                if orig_type in polys:
-                    clipped = clipped.loc[clipped.geom_type.isin(polys)]
-                elif orig_type in lines:
-                    clipped = clipped.loc[clipped.geom_type.isin(lines)]
+                if orig_type in POLYGON_GEOM_TYPES:
+                    clipped = clipped.loc[clipped.geom_type.isin(POLYGON_GEOM_TYPES)]
+                elif orig_type in LINE_GEOM_TYPES:
+                    clipped = clipped.loc[clipped.geom_type.isin(LINE_GEOM_TYPES)]
 
     return clipped
diff --git a/geopandas/tools/overlay.py b/geopandas/tools/overlay.py
index 053ee587..6332d425 100644
--- a/geopandas/tools/overlay.py
+++ b/geopandas/tools/overlay.py
@@ -6,7 +6,13 @@ import pandas as pd
 
 from geopandas import GeoDataFrame, GeoSeries
 from geopandas._compat import PANDAS_GE_30
-from geopandas.array import _check_crs, _crs_mismatch_warn
+from geopandas.array import (
+    LINE_GEOM_TYPES,
+    POINT_GEOM_TYPES,
+    POLYGON_GEOM_TYPES,
+    _check_crs,
+    _crs_mismatch_warn,
+)
 
 
 def _ensure_geometry_column(df):
@@ -39,7 +45,7 @@ def _overlay_intersection(df1, df2):
         right = df2.geometry.take(idx2)
         right.reset_index(drop=True, inplace=True)
         intersections = left.intersection(right)
-        poly_ix = intersections.geom_type.isin(["Polygon", "MultiPolygon"])
+        poly_ix = intersections.geom_type.isin(POLYGON_GEOM_TYPES)
         intersections.loc[poly_ix] = intersections[poly_ix].make_valid()
 
         # only keep actual intersecting geometries
@@ -96,7 +102,7 @@ def _overlay_difference(df1, df2):
         )
         new_g.append(new)
     differences = GeoSeries(new_g, index=df1.index, crs=df1.crs)
-    poly_ix = differences.geom_type.isin(["Polygon", "MultiPolygon"])
+    poly_ix = differences.geom_type.isin(POLYGON_GEOM_TYPES)
     differences.loc[poly_ix] = differences[poly_ix].make_valid()
     geom_diff = differences[~differences.is_empty].copy()
     dfdiff = df1[~differences.is_empty].copy()
@@ -261,13 +267,10 @@ def overlay(df1, df2, how="intersection", keep_geom_type=None, make_valid=True):
     else:
         keep_geom_type_warning = False
 
-    polys = ["Polygon", "MultiPolygon"]
-    lines = ["LineString", "MultiLineString", "LinearRing"]
-    points = ["Point", "MultiPoint"]
     for i, df in enumerate([df1, df2]):
-        poly_check = df.geom_type.isin(polys).any()
-        lines_check = df.geom_type.isin(lines).any()
-        points_check = df.geom_type.isin(points).any()
+        poly_check = df.geom_type.isin(POLYGON_GEOM_TYPES).any()
+        lines_check = df.geom_type.isin(LINE_GEOM_TYPES).any()
+        points_check = df.geom_type.isin(POINT_GEOM_TYPES).any()
         if sum([poly_check, lines_check, points_check]) > 1:
             raise NotImplementedError(f"df{i + 1} contains mixed geometry types.")
 
@@ -292,11 +295,17 @@ def overlay(df1, df2, how="intersection", keep_geom_type=None, make_valid=True):
     # Computations
     def _make_valid(df):
         df = df.copy()
-        if df.geom_type.isin(polys).all():
+        if df.geom_type.isin(POLYGON_GEOM_TYPES).all():
             mask = ~df.geometry.is_valid
             col = df._geometry_column_name
             if make_valid:
                 df.loc[mask, col] = df.loc[mask, col].make_valid()
+                # Extract only the input geometry type, as make_valid may change it
+                if mask.any():
+                    df = _collection_extract(
+                        df, geom_type="Polygon", keep_geom_type_warning=False
+                    )
+
             elif mask.any():
                 raise ValueError(
                     "You have passed make_valid=False along with "
@@ -306,6 +315,10 @@ def overlay(df1, df2, how="intersection", keep_geom_type=None, make_valid=True):
                 )
         return df
 
+    # Determine the geometry type before make_valid, as make_valid may change it
+    if keep_geom_type:
+        geom_type = df1.geom_type.iloc[0]
+
     df1 = _make_valid(df1)
     df2 = _make_valid(df2)
 
@@ -327,67 +340,69 @@ def overlay(df1, df2, how="intersection", keep_geom_type=None, make_valid=True):
             result.drop(["__idx1", "__idx2"], axis=1, inplace=True)
 
     if keep_geom_type:
-        geom_type = df1.geom_type.iloc[0]
+        result = _collection_extract(result, geom_type, keep_geom_type_warning)
 
-        # First we filter the geometry types inside GeometryCollections objects
-        # (e.g. GeometryCollection([polygon, point]) -> polygon)
-        # we do this separately on only the relevant rows, as this is an expensive
-        # operation (an expensive no-op for geometry types other than collections)
-        is_collection = result.geom_type == "GeometryCollection"
-        if is_collection.any():
-            geom_col = result._geometry_column_name
-            collections = result[[geom_col]][is_collection]
-
-            exploded = collections.reset_index(drop=True).explode(index_parts=True)
-            exploded = exploded.reset_index(level=0)
-
-            orig_num_geoms_exploded = exploded.shape[0]
-            if geom_type in polys:
-                exploded.loc[~exploded.geom_type.isin(polys), geom_col] = None
-            elif geom_type in lines:
-                exploded.loc[~exploded.geom_type.isin(lines), geom_col] = None
-            elif geom_type in points:
-                exploded.loc[~exploded.geom_type.isin(points), geom_col] = None
-            else:
-                raise TypeError(f"`keep_geom_type` does not support {geom_type}.")
-            num_dropped_collection = (
-                orig_num_geoms_exploded - exploded.geometry.isna().sum()
-            )
+    result.reset_index(drop=True, inplace=True)
+    return result
 
-            # level_0 created with above reset_index operation
-            # and represents the original geometry collections
-            # TODO avoiding dissolve to call union_all in this case could further
-            # improve performance (we only need to collect geometries in their
-            # respective Multi version)
-            dissolved = exploded.dissolve(by="level_0")
-            result.loc[is_collection, geom_col] = dissolved[geom_col].values
-        else:
-            num_dropped_collection = 0
-
-        # Now we filter all geometries (in theory we don't need to do this
-        # again for the rows handled above for GeometryCollections, but filtering
-        # them out is probably more expensive as simply including them when this
-        # is typically about only a few rows)
-        orig_num_geoms = result.shape[0]
-        if geom_type in polys:
-            result = result.loc[result.geom_type.isin(polys)]
-        elif geom_type in lines:
-            result = result.loc[result.geom_type.isin(lines)]
-        elif geom_type in points:
-            result = result.loc[result.geom_type.isin(points)]
-        else:
-            raise TypeError(f"`keep_geom_type` does not support {geom_type}.")
-        num_dropped = orig_num_geoms - result.shape[0]
-
-        if (num_dropped > 0 or num_dropped_collection > 0) and keep_geom_type_warning:
-            warnings.warn(
-                "`keep_geom_type=True` in overlay resulted in "
-                f"{num_dropped + num_dropped_collection} dropped geometries of "
-                "different geometry types than df1 has. Set `keep_geom_type=False` to "
-                "retain all geometries",
-                UserWarning,
-                stacklevel=2,
-            )
 
-    result.reset_index(drop=True, inplace=True)
+def _collection_extract(df, geom_type, keep_geom_type_warning):
+    # Check input
+    if geom_type in POLYGON_GEOM_TYPES:
+        geom_types = POLYGON_GEOM_TYPES
+    elif geom_type in LINE_GEOM_TYPES:
+        geom_types = LINE_GEOM_TYPES
+    elif geom_type in POINT_GEOM_TYPES:
+        geom_types = POINT_GEOM_TYPES
+    else:
+        raise TypeError(f"`geom_type` does not support {geom_type}.")
+
+    result = df.copy()
+
+    # First we filter the geometry types inside GeometryCollections objects
+    # (e.g. GeometryCollection([polygon, point]) -> polygon)
+    # we do this separately on only the relevant rows, as this is an expensive
+    # operation (an expensive no-op for geometry types other than collections)
+    is_collection = result.geom_type == "GeometryCollection"
+    if is_collection.any():
+        geom_col = result._geometry_column_name
+        collections = result.loc[is_collection, [geom_col]]
+
+        exploded = collections.reset_index(drop=True).explode(index_parts=True)
+        exploded = exploded.reset_index(level=0)
+
+        orig_num_geoms_exploded = exploded.shape[0]
+        exploded.loc[~exploded.geom_type.isin(geom_types), geom_col] = None
+        num_dropped_collection = (
+            orig_num_geoms_exploded - exploded.geometry.isna().sum()
+        )
+
+        # level_0 created with above reset_index operation
+        # and represents the original geometry collections
+        # TODO avoiding dissolve to call union_all in this case could further
+        # improve performance (we only need to collect geometries in their
+        # respective Multi version)
+        dissolved = exploded.dissolve(by="level_0")
+        result.loc[is_collection, geom_col] = dissolved[geom_col].values
+    else:
+        num_dropped_collection = 0
+
+    # Now we filter all geometries (in theory we don't need to do this
+    # again for the rows handled above for GeometryCollections, but filtering
+    # them out is probably more expensive as simply including them when this
+    # is typically about only a few rows)
+    orig_num_geoms = result.shape[0]
+    result = result.loc[result.geom_type.isin(geom_types)]
+    num_dropped = orig_num_geoms - result.shape[0]
+
+    if (num_dropped > 0 or num_dropped_collection > 0) and keep_geom_type_warning:
+        warnings.warn(
+            "`keep_geom_type=True` in overlay resulted in "
+            f"{num_dropped + num_dropped_collection} dropped geometries of different "
+            "geometry types than df1 has. Set `keep_geom_type=False` to retain all "
+            "geometries",
+            UserWarning,
+            stacklevel=2,
+        )
+
     return result


Test patch:
diff --git a/geopandas/tests/test_overlay.py b/geopandas/tests/test_overlay.py
index 482e4384..ed2f4c1d 100644
--- a/geopandas/tests/test_overlay.py
+++ b/geopandas/tests/test_overlay.py
@@ -110,6 +110,22 @@ def test_overlay(dfs_index, how):
         assert_geodataframe_equal(result, expected, check_column_type=False)
 
 
+@pytest.mark.parametrize("keep_geom_type", [True, False])
+def test_overlay_invalid_input(keep_geom_type):
+    invalid_polygon = Polygon([(0, 0), (1, 1), (1, 2), (1, 1), (0, 0)])
+    square = Polygon(((0, 0), (1, 0), (1, 1), (0, 1), (0, 0)))
+    df1 = geopandas.GeoDataFrame({"A": [1, 2]}, geometry=[invalid_polygon, square])
+    df2 = geopandas.GeoDataFrame({"B": [5, 6]}, geometry=[square, square])
+
+    df1_df2 = df1.overlay(df2, keep_geom_type=keep_geom_type)
+    # As overlay only supports inputs with uniform geometry types, `make_valid` of
+    # `invalid_polygon` should always result in POLYGON EMPTY, so the number of results
+    # will should always be 2.
+    assert len(df1_df2) == 2
+    if keep_geom_type:
+        assert df1_df2.geom_type.isin(["Polygon"]).all()
+
+
 @pytest.mark.filterwarnings("ignore:GeoSeries crs mismatch:UserWarning")
 def test_overlay_nybb(how, nybb_filename):
     polydf = read_file(nybb_filename)
diff --git a/geopandas/tools/tests/test_clip.py b/geopandas/tools/tests/test_clip.py
index dbdbfe6e..bf465e30 100644
--- a/geopandas/tools/tests/test_clip.py
+++ b/geopandas/tools/tests/test_clip.py
@@ -17,6 +17,7 @@ from shapely.geometry import (
 import geopandas
 from geopandas import GeoDataFrame, GeoSeries, clip
 from geopandas._compat import HAS_PYPROJ
+from geopandas.array import POLYGON_GEOM_TYPES
 from geopandas.tools.clip import _mask_is_list_like_rectangle
 
 import pytest
@@ -304,7 +305,7 @@ class TestClipWithSingleRectangleGdf:
         )
         assert np.array_equal(clipped.total_bounds, expected_bounds)
         # Assert returned data is a not geometry collection
-        assert (clipped.geom_type.isin(["Polygon", "MultiPolygon"])).all()
+        assert (clipped.geom_type.isin(POLYGON_GEOM_TYPES)).all()
 
     def test_clip_multiline(self, multi_line, mask):
         """Test that clipping a multiline feature with a poly returns expected