Repository: pysal/momepy

Issue:
BUG: verify handling of MultiIndex
It seems that when your GeoDataFrame is using MultiIindex, it does affect some parts of `momepy` and will likely break when using Graph. We should test and verify that it behaves as we want, both here and in libpysal. In some cases, raising may be the best option but if we can support it, it is better. I don't think that there's scope in supporting MultiIndex in Graph so it may need to be checked only in functions that do not use graphs.

Discussion:
codecov[bot]: ## [Codecov](https://app.codecov.io/gh/pysal/momepy/pull/622?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal) Report
All modified and coverable lines are covered by tests :white_check_mark:
> Project coverage is 97.7%. Comparing base [(`4037c70`)](https://app.codecov.io/gh/pysal/momepy/commit/4037c703ae35c0258939c4a2f326dbf08538bc20?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal) to head [(`efbb4d5`)](https://app.codecov.io/gh/pysal/momepy/commit/efbb4d506108c97e70314c78ae0c48da063f4ca7?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal).
> Report is 58 commits behind head on main.

<details><summary>Additional details and impacted files</summary>


[![Impacted file tree graph](https://app.codecov.io/gh/pysal/momepy/pull/622/graphs/tree.svg?width=650&height=150&src=pr&token=VNn0WR5JWT&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal)](https://app.codecov.io/gh/pysal/momepy/pull/622?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal)

```diff
@@           Coverage Diff           @@
##            main    #622     +/-   ##
=======================================
+ Coverage   97.4%   97.7%   +0.3%     
=======================================
  Files         26      39     +13     
  Lines       4328    6226   +1898     
=======================================
+ Hits        4214    6082   +1868     
- Misses       114     144     +30     
```

| [Files](https://app.codecov.io/gh/pysal/momepy/pull/622?dropdown=coverage&src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal) | Coverage Δ | |
|---|---|---|
| [momepy/functional/tests/test\_multi\_index.py](https://app.codecov.io/gh/pysal/momepy/pull/622?src=pr&el=tree&filepath=momepy%2Ffunctional%2Ftests%2Ftest_multi_index.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal#diff-bW9tZXB5L2Z1bmN0aW9uYWwvdGVzdHMvdGVzdF9tdWx0aV9pbmRleC5weQ==) | `100.0% <100.0%> (ø)` | |

... and [1 file with indirect coverage changes](https://app.codecov.io/gh/pysal/momepy/pull/622/indirect-changes?src=pr&el=tree-more&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=pysal)

</details>

Patch:
diff --git a/momepy/functional/_elements.py b/momepy/functional/_elements.py
index 1b6b2f3..a30def3 100644
--- a/momepy/functional/_elements.py
+++ b/momepy/functional/_elements.py
@@ -10,7 +10,7 @@ from joblib import Parallel, delayed
 from libpysal.cg import voronoi_frames
 from libpysal.graph import Graph
 from packaging.version import Version
-from pandas import Series
+from pandas import MultiIndex, Series
 
 GPD_GE_013 = Version(gpd.__version__) >= Version("0.13.0")
 GPD_GE_10 = Version(gpd.__version__) >= Version("1.0dev")
@@ -103,6 +103,12 @@ def morphological_tessellation(
     4  POLYGON ((1603084.231 6464104.386, 1603083.773...
 
     """
+
+    if isinstance(geometry.index, MultiIndex):
+        raise ValueError(
+            "MultiIndex is not supported in `momepy.morphological_tessellation`."
+        )
+
     if isinstance(clip, GeoSeries | GeoDataFrame):
         clip = clip.union_all() if GPD_GE_10 else clip.unary_union
 
@@ -218,6 +224,11 @@ def enclosed_tessellation(
     126  POLYGON ((1603528.593 6464221.033, 1603527.796...                0
     """
 
+    if isinstance(geometry.index, MultiIndex):
+        raise ValueError(
+            "MultiIndex is not supported in `momepy.enclosed_tessellation`."
+        )
+
     # convert to GeoDataFrame and add position (we will need it later)
     enclosures = enclosures.geometry.to_frame()
     enclosures["position"] = range(len(enclosures))
@@ -360,6 +371,12 @@ def verify_tessellation(tessellation, geometry):
 
     >>> excluded, multipolygons = momepy.verify_tessellation(tessellation, buildings)
     """
+
+    if isinstance(geometry.index, MultiIndex) or isinstance(
+        tessellation.index, MultiIndex
+    ):
+        raise ValueError("MultiIndex is not supported in `momepy.verify_tessellation`.")
+
     # check against input layer
     ids_original = geometry.index
     ids_generated = tessellation.index
@@ -397,7 +414,7 @@ def get_nearest_street(
     buildings: GeoSeries | GeoDataFrame,
     streets: GeoSeries | GeoDataFrame,
     max_distance: float | None = None,
-) -> np.ndarray:
+) -> Series:
     """Identify the nearest street for each building.
 
     Parameters
@@ -429,30 +446,26 @@ def get_nearest_street(
     Get street index.
 
     >>> momepy.get_nearest_street(buildings, streets)
-    array([ 0., 33., 10.,  8.,  8.,  8.,  8.,  8., 33., 11., 11., 28., 28.,
-           28., 28., 28., 16.,  8.,  8.,  8.,  8.,  8.,  8., 11., 28., 28.,
-           28.,  8.,  8.,  8.,  8., 16., 28., 28., 28., 28., 28.,  1., 21.,
-           21., 21., 21., 21., 12., 12., 12., 26., 26., 26., 19., 19., 19.,
-           19., 21., 21., 21., 32., 32., 32., 32., 32., 26., 26.,  5.,  5.,
-            5.,  5.,  2.,  2.,  2.,  2.,  2.,  2., 25., 25., 25., 19., 19.,
-           19., 19.,  5., 25.,  6., 33., 33., 33., 33., 33., 33., 33., 34.,
-           34., 34., 34., 34., 34., 34., 34.,  6.,  6.,  6.,  6.,  6., 34.,
-           33.,  6., 34., 34., 34., 34.,  0.,  0.,  0.,  0.,  0.,  0., 34.,
-           34., 34.,  0.,  0., 14.,  2.,  2., 25., 24.,  2.,  2.,  2.,  2.,
-           24., 24., 24., 24., 24., 28., 12., 28., 34., 34., 32., 21., 16.,
-           19.], dtype=float32)
+    0       0.0
+    1      33.0
+    2      10.0
+    3       8.0
+    4       8.0
+        ...
+    139    34.0
+    140    32.0
+    141    21.0
+    142    16.0
+    143    19.0
+    Length: 144, dtype: float64
     """
     blg_idx, str_idx = streets.sindex.nearest(
         buildings.geometry, return_all=False, max_distance=max_distance
     )
 
-    if streets.index.dtype == "object":
-        ids = np.empty(len(buildings), dtype=object)
-    else:
-        ids = np.empty(len(buildings), dtype=np.float32)
-        ids[:] = np.nan
+    ids = pd.Series(None, index=buildings.index, dtype=streets.index.dtype)
 
-    ids[blg_idx] = streets.index[str_idx]
+    ids.iloc[blg_idx] = streets.index[str_idx]
     return ids
 
 
@@ -522,6 +535,15 @@ def get_nearest_node(
     143    22.0
     Length: 144, dtype: float64
     """
+
+    if (
+        isinstance(buildings.index, MultiIndex)
+        or isinstance(nearest_edge.index, MultiIndex)
+        or isinstance(nodes.index, MultiIndex)
+        or isinstance(edges.index, MultiIndex)
+    ):
+        raise ValueError("MultiIndex is not supported in `momepy.get_nearest_node`.")
+
     # treat possibly missing edge index
     a = np.empty(len(buildings))
     na_mask = np.isnan(nearest_edge)
@@ -606,6 +628,12 @@ def generate_blocks(
     >>> tessellation["block_id"] = tessellation_id
     """
 
+    if (
+        isinstance(buildings.index, MultiIndex)
+        or isinstance(tessellation.index, MultiIndex)
+        or isinstance(edges.index, MultiIndex)
+    ):
+        raise ValueError("MultiIndex is not supported in `momepy.generate_blocks`.")
     id_name: str = "bID"
 
     # slice the tessellations by the street network
diff --git a/momepy/functional/_shape.py b/momepy/functional/_shape.py
index 8e1d734..a4de964 100644
--- a/momepy/functional/_shape.py
+++ b/momepy/functional/_shape.py
@@ -4,7 +4,7 @@ import shapely
 from geopandas import GeoDataFrame, GeoSeries
 from numpy.typing import NDArray
 from packaging.version import Version
-from pandas import DataFrame, Series
+from pandas import DataFrame, MultiIndex, Series
 
 from momepy.functional import _dimension
 
@@ -724,6 +724,10 @@ def centroid_corner_distance(
             "momepy.centroid_corner_distance requires geopandas 0.13 or later. "
         )
 
+    result_index = geometry.index
+    if isinstance(geometry.index, MultiIndex):
+        geometry = geometry.reset_index(drop=True)
+
     def _ccd(points: DataFrame, eps: float) -> Series:
         centroid = points.values[0, 2:]
         pts = points.values[:-1, :2]
@@ -738,7 +742,7 @@ def centroid_corner_distance(
         coords = geometry.exterior.get_coordinates(index_parts=False)
     coords[["cent_x", "cent_y"]] = geometry.centroid.get_coordinates(index_parts=False)
     ccd = coords.groupby(level=0).apply(_ccd, eps=eps)
-    ccd.index = geometry.index
+    ccd.index = result_index
     return ccd
 
 


Test patch:
diff --git a/momepy/functional/tests/test_elements.py b/momepy/functional/tests/test_elements.py
index f6b0b26..1fc2215 100644
--- a/momepy/functional/tests/test_elements.py
+++ b/momepy/functional/tests/test_elements.py
@@ -222,7 +222,7 @@ class TestElements:
 
         streets.index = streets.index.astype(str)
         nearest = mm.get_nearest_street(self.df_buildings, streets, 10)
-        assert (nearest == None).sum() == 137  # noqa: E711
+        assert pd.isna(nearest).sum() == 137  # noqa: E711
 
     def test_get_nearest_node(self):
         nodes, edges = mm.nx_to_gdf(mm.gdf_to_nx(self.df_streets))
@@ -332,6 +332,37 @@ class TestElements:
         else:
             assert len(blocks.sindex.query_bulk(blocks.geometry, "overlaps")[0]) == 0
 
+    def test_multi_index(self):
+        buildings = self.df_buildings.set_index(["uID", "uID"])
+        with pytest.raises(
+            ValueError,
+            match="MultiIndex is not supported in `momepy.morphological_tessellation`.",
+        ):
+            mm.morphological_tessellation(buildings)
+        with pytest.raises(
+            ValueError,
+            match="MultiIndex is not supported in `momepy.enclosed_tessellation`.",
+        ):
+            mm.enclosed_tessellation(buildings, self.enclosures)
+        with pytest.raises(
+            ValueError,
+            match="MultiIndex is not supported in `momepy.verify_tessellation`.",
+        ):
+            mm.verify_tessellation(buildings, self.enclosures)
+
+        with pytest.raises(
+            ValueError,
+            match="MultiIndex is not supported in `momepy.get_nearest_node`.",
+        ):
+            mm.get_nearest_node(
+                buildings, self.enclosures, self.enclosures, self.enclosures
+            )
+
+        with pytest.raises(
+            ValueError, match="MultiIndex is not supported in `momepy.generate_blocks`"
+        ):
+            mm.generate_blocks(buildings, self.enclosures, self.enclosures)
+
     def test_tess_single_building_edge_case(self):
         tessellations = mm.enclosed_tessellation(
             self.df_buildings, self.enclosures.geometry, n_jobs=-1