Repository: uscuni/neatnet

Issue:
rough edge in `sgeop.nodes.remove_false_nodes` when multilinestrings passed
When multilinestrings are passed in to `sgeop.nodes.remove_false_nodes()` some strange behavior is observed due to no `.explode()` being called like in the counterpart [`momepy/preprocessing.remove_false_nodes()`](https://github.com/pysal/momepy/blob/f1965cccd48759faf6a2a37739bda310d851124d/momepy/preprocessing.py#L189)

* xref #59

Discussion:

Solution patch:
diff --git a/sgeop/nodes.py b/sgeop/nodes.py
index 55ce2f8..1749c93 100644
--- a/sgeop/nodes.py
+++ b/sgeop/nodes.py
@@ -172,25 +172,44 @@ def weld_edges(
     ).tolist()
 
 
-def remove_false_nodes(gdf, aggfunc="first", **kwargs):
-    """Reimplementation of momepy.remove_false_nodes that preserves attributes
+def remove_false_nodes(
+    gdf: gpd.GeoSeries | gpd.GeoDataFrame, aggfunc: str = "first", **kwargs
+):
+    """Reimplementation of ``momepy.remove_false_nodes()`` that preserves attributes.
 
     Parameters
     ----------
-    gdf : _type_
-        _description_
+    gdf : gpd.GeoSeries | gpd.GeoDataFrame
+        Input edgelines process. If any edges are ``MultiLineString`` they
+        will be exploded into constituent ``LineString`` components.
+    aggfunc : str = "first"
+        Aggregate function for processing non-spatial component.
+    **kwargs
+        Keyword arguments for ``aggfunc``.
 
     Returns
     -------
-    _type_
-        _description_
+    gpd.GeoSeries | gpd.GeoDataFrame
+       The original input ``gdf`` if only 1 edgeline, otherwise the processed
+       edgeline without interstitial nodes.
+
+    Notes
+    -----
+    Any 3D geometries are (potentially) downcast in loops.
     """
+
+    def merge_geometries(block: gpd.GeoSeries) -> shapely.LineString:
+        """Helper in processing the spatial component."""
+        return shapely.line_merge(shapely.GeometryCollection(block.values))
+
     if len(gdf) < 2:
         return gdf
 
     if isinstance(gdf, gpd.GeoSeries):
         gdf = gdf.to_frame("geometry")
 
+    gdf = gdf.explode(ignore_index=True)
+
     labels = get_components(gdf.geometry)
 
     # Process non-spatial component
@@ -199,14 +218,11 @@ def remove_false_nodes(gdf, aggfunc="first", **kwargs):
     aggregated_data.columns = aggregated_data.columns.to_flat_index()
 
     # Process spatial component
-    def merge_geometries(block):
-        merged_geom = shapely.line_merge(shapely.GeometryCollection(block.values))
-        return merged_geom
-
     g = gdf.groupby(group_keys=False, by=labels)[gdf.geometry.name].agg(
         merge_geometries
     )
     aggregated_geometry = gpd.GeoDataFrame(g, geometry=gdf.geometry.name, crs=gdf.crs)
+
     # Recombine
     aggregated = aggregated_geometry.join(aggregated_data)
 
@@ -225,7 +241,7 @@ def merge_geometries(block):
         target_nodes = nodes.geometry.iloc[node_ix[loop_ix == ix]]
         if len(target_nodes) == 2:
             node_coords = shapely.get_coordinates(target_nodes)
-            coords = np.array(loop_geom.coords)
+            coords = shapely.get_coordinates(loop_geom)
             new_start = (
                 node_coords[0]
                 if (node_coords[0] != coords[0]).all()
diff --git a/sgeop/simplify.py b/sgeop/simplify.py
index 43a2ba4..82dd347 100644
--- a/sgeop/simplify.py
+++ b/sgeop/simplify.py
@@ -173,65 +173,6 @@ def simplify_singletons(
         return cleaned_roads
 
 
-def simplify_clusters(
-    artifacts,
-    roads,
-    max_segment_length=1,
-    eps=1e-4,
-    simplification_factor=2,
-    min_dangle_length=20,
-    consolidation_tolerance=10,
-):
-    # Get nodes from the network.
-    nodes = momepy.nx_to_gdf(momepy.node_degree(momepy.gdf_to_nx(roads)), lines=False)
-
-    # collect changes
-    to_drop = []
-    to_add = []
-
-    for _, artifact in artifacts.groupby("comp"):
-        # get artifact cluster polygon
-        cluster_geom = artifact.union_all()
-        # get edges relevant for an artifact
-        edges = roads.iloc[
-            roads.sindex.query(cluster_geom, predicate="intersects")
-        ].copy()
-
-        nx_gx_cluster(
-            edges=edges,
-            cluster_geom=cluster_geom,
-            nodes=nodes,
-            to_drop=to_drop,
-            to_add=to_add,
-            eps=eps,
-            max_segment_length=max_segment_length,
-            min_dangle_length=min_dangle_length,
-            consolidation_tolerance=consolidation_tolerance,
-        )
-
-    cleaned_roads = roads.drop(to_drop)
-
-    # create new roads with fixed geometry. Note that to_add and to_drop lists shall be
-    # global and this step should happen only once, not for every artifact
-    new = gpd.GeoDataFrame(geometry=to_add, crs=roads.crs)
-    new["_status"] = "new"
-    new["geometry"] = new.line_merge().simplify(
-        max_segment_length * simplification_factor
-    )
-    new_roads = pd.concat(
-        [
-            cleaned_roads,
-            new,
-        ],
-        ignore_index=True,
-    ).explode()
-    new_roads = remove_false_nodes(
-        new_roads[~new_roads.is_empty], aggfunc={"_status": _status}
-    ).drop_duplicates("geometry")
-
-    return new_roads
-
-
 def simplify_pairs(
     artifacts,
     roads,
@@ -375,6 +316,65 @@ def simplify_pairs(
     return roads_cleaned
 
 
+def simplify_clusters(
+    artifacts,
+    roads,
+    max_segment_length=1,
+    eps=1e-4,
+    simplification_factor=2,
+    min_dangle_length=20,
+    consolidation_tolerance=10,
+):
+    # Get nodes from the network.
+    nodes = momepy.nx_to_gdf(momepy.node_degree(momepy.gdf_to_nx(roads)), lines=False)
+
+    # collect changes
+    to_drop = []
+    to_add = []
+
+    for _, artifact in artifacts.groupby("comp"):
+        # get artifact cluster polygon
+        cluster_geom = artifact.union_all()
+        # get edges relevant for an artifact
+        edges = roads.iloc[
+            roads.sindex.query(cluster_geom, predicate="intersects")
+        ].copy()
+
+        nx_gx_cluster(
+            edges=edges,
+            cluster_geom=cluster_geom,
+            nodes=nodes,
+            to_drop=to_drop,
+            to_add=to_add,
+            eps=eps,
+            max_segment_length=max_segment_length,
+            min_dangle_length=min_dangle_length,
+            consolidation_tolerance=consolidation_tolerance,
+        )
+
+    cleaned_roads = roads.drop(to_drop)
+
+    # create new roads with fixed geometry. Note that to_add and to_drop lists shall be
+    # global and this step should happen only once, not for every artifact
+    new = gpd.GeoDataFrame(geometry=to_add, crs=roads.crs)
+    new["_status"] = "new"
+    new["geometry"] = new.line_merge().simplify(
+        max_segment_length * simplification_factor
+    )
+    new_roads = pd.concat(
+        [
+            cleaned_roads,
+            new,
+        ],
+        ignore_index=True,
+    ).explode()
+    new_roads = remove_false_nodes(
+        new_roads[~new_roads.is_empty], aggfunc={"_status": _status}
+    ).drop_duplicates("geometry")
+
+    return new_roads
+
+
 def get_type(edges, shared_edge):
     if (  # roundabout special case
         edges.coins_group.nunique() == 1 and edges.shape[0] == edges.coins_count.iloc[0]


Test patch:
diff --git a/sgeop/tests/test_nodes.py b/sgeop/tests/test_nodes.py
index b667d3d..2abee04 100644
--- a/sgeop/tests/test_nodes.py
+++ b/sgeop/tests/test_nodes.py
@@ -2,6 +2,7 @@
 import itertools
 
 import geopandas.testing
+import momepy
 import numpy
 import pandas
 import pytest
@@ -304,3 +305,134 @@ def test_get_components(edgelines, ignore, known):
 def test_weld_edges(edgelines, ignore, known):
     observed = sgeop.nodes.weld_edges(edgelines, ignore=ignore)
     numpy.testing.assert_array_equal(observed, known)
+
+
+class TestRemoveFalseNodes:
+    def setup_method(self):
+        p10 = shapely.Point(1, 0)
+        p20 = shapely.Point(2, 0)
+        p30 = shapely.Point(3, 0)
+        p40 = shapely.Point(4, 0)
+        p50 = shapely.Point(5, 0)
+        p21 = shapely.Point(2, 1)
+        p32 = shapely.Point(3, 2)
+        p41 = shapely.Point(4, 1)
+
+        self.line1020 = shapely.LineString((p10, p20))
+        self.line2030 = shapely.LineString((p20, p30))
+        self.line3040 = shapely.LineString((p30, p40))
+        self.line4050 = shapely.LineString((p40, p50))
+        self.line3021 = shapely.LineString((p30, p21))
+        self.line2132 = shapely.LineString((p21, p32))
+        self.line4132 = shapely.LineString((p41, p32))
+        self.line3041 = shapely.LineString((p30, p41))
+        self.attrs = ["cat"] * 3 + ["dog"] * 3 + ["eel"] * 2
+
+        self.series = geopandas.GeoSeries(
+            [
+                self.line1020,
+                self.line2030,
+                self.line3040,
+                self.line4050,
+                self.line3021,
+                self.line2132,
+                self.line4132,
+                self.line3041,
+            ]
+        )
+
+        self.line102030 = shapely.LineString((p10, p20, p30))
+        self.line304050 = shapely.LineString((p30, p40, p50))
+        self.line3041323130 = shapely.LineString((p30, p41, p32, p21, p30))
+
+        self.known_geoms = [
+            self.line102030,
+            self.line3041323130,
+            self.line304050,
+        ]
+
+    def test_single_series(self):
+        one_in_series = self.series[:0].copy()
+        known = one_in_series
+        observed = sgeop.nodes.remove_false_nodes(one_in_series)
+        geopandas.testing.assert_geoseries_equal(observed, known)
+
+    def test_series(self):
+        known = geopandas.GeoDataFrame(geometry=self.known_geoms)
+        observed = sgeop.nodes.remove_false_nodes(self.series)
+        geopandas.testing.assert_geodataframe_equal(observed, known)
+
+    def test_frame(self):
+        known = geopandas.GeoDataFrame(geometry=self.known_geoms)
+        observed = sgeop.nodes.remove_false_nodes(
+            geopandas.GeoDataFrame(geometry=self.series)
+        )
+        geopandas.testing.assert_geodataframe_equal(observed, known)
+
+    def test_frame_attrs_first(self):
+        known = geopandas.GeoDataFrame(
+            {"animal": ["cat", "dog", "cat"]},
+            geometry=self.known_geoms,
+            columns=["geometry", "animal"],
+        )
+        observed = sgeop.nodes.remove_false_nodes(
+            geopandas.GeoDataFrame({"animal": self.attrs}, geometry=self.series)
+        )
+        geopandas.testing.assert_geodataframe_equal(observed, known)
+
+    def test_frame_attrs_last(self):
+        known = geopandas.GeoDataFrame(
+            {"animal": ["cat", "eel", "dog"]},
+            geometry=self.known_geoms,
+            columns=["geometry", "animal"],
+        )
+        observed = sgeop.nodes.remove_false_nodes(
+            geopandas.GeoDataFrame({"animal": self.attrs}, geometry=self.series),
+            aggfunc="last",
+        )
+        geopandas.testing.assert_geodataframe_equal(observed, known)
+
+    def test_momepy_suite(self):
+        false_network = geopandas.read_file(
+            momepy.datasets.get_path("tests"), layer="network"
+        )
+        false_network["vals"] = range(len(false_network))
+        fixed = sgeop.remove_false_nodes(false_network).reset_index(drop=True)
+        assert len(fixed) == 56
+        assert isinstance(fixed, geopandas.GeoDataFrame)
+        assert false_network.crs.equals(fixed.crs)
+        assert sorted(false_network.columns) == sorted(fixed.columns)
+
+        # check loop order
+        expected = numpy.array(
+            [
+                [-727238.49292668, -1052817.28071986],
+                [-727253.1752498, -1052827.47329062],
+                [-727223.93217677, -1052829.47624082],
+                [-727238.49292668, -1052817.28071986],
+            ]
+        )
+        numpy.testing.assert_almost_equal(
+            numpy.array(fixed.loc[55].geometry.coords), expected
+        )
+
+        fixed_series = sgeop.nodes.remove_false_nodes(
+            false_network.geometry
+        ).reset_index(drop=True)
+        assert len(fixed_series) == 56
+        assert isinstance(fixed_series, geopandas.GeoDataFrame)
+        assert false_network.crs.equals(fixed_series.crs)
+
+        multiindex = false_network.explode(index_parts=True)
+        fixed_multiindex = sgeop.nodes.remove_false_nodes(multiindex)
+        assert len(fixed_multiindex) == 56
+        assert isinstance(fixed, geopandas.GeoDataFrame)
+        assert sorted(false_network.columns) == sorted(fixed.columns)
+
+        # no node of a degree 2
+        df_streets = geopandas.read_file(
+            momepy.datasets.get_path("bubenec"), layer="streets"
+        )
+        known = df_streets.drop([4, 7, 17, 22]).reset_index(drop=True)
+        observed = sgeop.nodes.remove_false_nodes(known).reset_index(drop=True)
+        geopandas.testing.assert_geodataframe_equal(observed, known)