Merge pull request #313 from statisticsnorway/nwa-speed

undo nwa changes because of floating point precision...

Author: mortewle

pyproject.toml

@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "ssb-sgis"
-version = "1.3.8"
+version = "1.3.9"
 description = "GIS functions used at Statistics Norway."
 authors = ["Morten Letnes <morten.letnes@ssb.no>"]
 license = "MIT"

src/sgis/geopandas_tools/general.py

@@ -915,8 +915,8 @@ def _prepare_make_edge_cols_simple(
     lines: GeoDataFrame,
 ) -> tuple[GeoDataFrame, GeoDataFrame]:
     """Faster version of _prepare_make_edge_cols."""
+    lines = lines[(lines.geometry.notna()) & (~lines.geometry.is_empty)]
     endpoints = lines.geometry.boundary.explode(ignore_index=False)
-
     if len(lines) and len(endpoints) / len(lines) != 2:
         raise ValueError(
             "The lines should have only two endpoints each. "

src/sgis/networkanalysis/_points.py

@@ -84,17 +84,17 @@ def _get_edges_and_weights(
 
         search_factor_multiplier = 1 + rules.search_factor / 100
         distances = distances.loc[
-            lambda df: (df.distance <= rules.search_tolerance)
+            lambda df: (df["distance"] <= rules.search_tolerance)
             & (
-                df.distance
-                <= df.distance_min * search_factor_multiplier + rules.search_factor
+                df["distance"]
+                <= df["distance_min"] * search_factor_multiplier + rules.search_factor
             )
         ]
 
         edges = self._make_edges(distances, from_col=from_col, to_col=to_col)
 
         weighs = self._convert_distance_to_weight(
-            distances=list(distances.distance), rules=rules
+            distances=list(distances["distance"]), rules=rules
         )
 
         return edges, weighs

src/sgis/networkanalysis/_service_area.py

@@ -97,7 +97,6 @@ def _service_area(
             partly_within["remaining_distance"] = break_ - partly_within[weight]
 
             split_lines = _split_lines(partly_within, directed)
-
             split_lines = make_edge_wkt_cols(split_lines)
 
             # select the cutted lines shorter than the cut distance that intersects

src/sgis/networkanalysis/network.py

@@ -139,16 +139,16 @@ def _nodes_are_up_to_date(self) -> bool:
         """
         if not hasattr(self, "_nodes"):
             return False
-        new_or_missing = (~self.gdf.source.isin(self._nodes.node_id)) | (
-            ~self.gdf.target.isin(self._nodes.node_id)
+        new_or_missing = (~self.gdf["source"].isin(self._nodes["node_id"])) | (
+            ~self.gdf["target"].isin(self._nodes["node_id"])
         )
 
         if any(new_or_missing):
             return False
 
         removed = ~(
-            (self._nodes.node_id.isin(self.gdf.source))
-            | (self._nodes.node_id.isin(self.gdf.target))
+            (self._nodes["node_id"].isin(self.gdf["source"]))
+            | (self._nodes["node_id"].isin(self.gdf["target"]))
         )
 
         if any(removed):
@@ -208,9 +208,12 @@ def _warn_if_undirected(self) -> None:
         warnings.warn(mess, stacklevel=2)
 
     @property
-    def percent_bidirectional(self) -> float:
+    def percent_bidirectional(self) -> float | None:
         """The percentage of lines that appear in both directions."""
-        return self._percent_bidirectional
+        try:
+            return self._percent_bidirectional
+        except AttributeError:
+            return None
 
     def copy(self) -> "Network":
         """Returns a shallow copy of the class instance."""
@@ -224,7 +227,7 @@ def __repr__(self) -> str:
         """The print representation."""
         cl = self.__class__.__name__
         km = int(sum(self.gdf.length) / 1000)
-        return f"{cl}({km} km, percent_bidirectional={self._percent_bidirectional})"
+        return f"{cl}({km} km, percent_bidirectional={self.percent_bidirectional})"
 
     def __len__(self) -> int:
         """Number og rows in the GeoDataFrame."""

src/sgis/networkanalysis/networkanalysis.py

@@ -20,7 +20,6 @@
 from shapely import force_2d
 
 from ..geopandas_tools.general import _push_geom_col
-from ..geopandas_tools.sfilter import sfilter_inverse
 from ._get_route import _get_k_routes
 from ._get_route import _get_route
 from ._get_route import _get_route_frequencies
@@ -31,8 +30,6 @@
 from .cutting_lines import split_lines_by_nearest_point
 from .network import Network
 from .networkanalysisrules import NetworkAnalysisRules
-from .nodes import _map_node_ids_from_wkt
-from .nodes import make_node_ids
 
 
 class NetworkAnalysis:
@@ -413,8 +410,8 @@ def od_cost_matrix(
         if lines:
             results = _push_geom_col(results)
 
-        if self.rules.split_lines:
-            self._unsplit_network()
+        # if self.rules.split_lines:
+        #     self._unsplit_network()
 
         if self._log:
             minutes_elapsed = round((perf_counter() - time_) / 60, 1)
@@ -650,8 +647,8 @@ def multiindex_mapper(x: tuple[int, int]) -> tuple[int, int]:
 
         results.geometry = force_2d(results.geometry)
 
-        if self.rules.split_lines:
-            self._unsplit_network()
+        # if self.rules.split_lines:
+        #     self._unsplit_network()
 
         if self._log:
             minutes_elapsed = round((perf_counter() - time_) / 60, 1)
@@ -760,8 +757,8 @@ def get_route(
         results["origin"] = results["origin"].map(self.origins.idx_dict)
         results["destination"] = results["destination"].map(self.destinations.idx_dict)
 
-        if self.rules.split_lines:
-            self._unsplit_network()
+        # if self.rules.split_lines:
+        #     self._unsplit_network()
 
         if self._log:
             minutes_elapsed = round((perf_counter() - time_) / 60, 1)
@@ -928,8 +925,8 @@ def get_k_routes(
         if isinstance(results, GeoDataFrame):
             results = _push_geom_col(results)
 
-        if self.rules.split_lines:
-            self._unsplit_network()
+        # if self.rules.split_lines:
+        #     self._unsplit_network()
 
         if self._log:
             minutes_elapsed = round((perf_counter() - time_) / 60, 1)
@@ -1054,8 +1051,8 @@ def service_area(
 
             results = _push_geom_col(results)
 
-        if self.rules.split_lines:
-            self._unsplit_network()
+        # if self.rules.split_lines:
+        #     self._unsplit_network()
 
         if self._log:
             minutes_elapsed = round((perf_counter() - time_) / 60, 1)
@@ -1183,8 +1180,8 @@ def precice_service_area(
 
             results = _push_geom_col(results)
 
-        if self.rules.split_lines:
-            self._unsplit_network()
+        # if self.rules.split_lines:
+        #     self._unsplit_network()
 
         if self._log:
             minutes_elapsed = round((perf_counter() - time_) / 60, 1)
@@ -1371,6 +1368,8 @@ def _prepare_network_analysis(
                 "rowwise=True"
             )
 
+        self._unsplit_network()
+
         self.network.gdf = self.rules._validate_weight(self.network.gdf)
 
         self.origins = Origins(origins)
@@ -1383,6 +1382,7 @@ def _prepare_network_analysis(
             self.destinations._make_temp_idx(
                 start=max(self.origins.gdf["temp_idx"].astype(int)) + 1
             )
+
         else:
             self.destinations = None
 
@@ -1413,12 +1413,13 @@ def _get_edges_and_weights(
         """
         if self.rules.split_lines:
             self._split_lines()
+            self.network._make_node_ids()
             self.origins._make_temp_idx(
-                start=max(self.network.nodes.node_id.astype(int)) + 1
+                start=max(self.network.nodes["node_id"].astype(int)) + 1
             )
             if self.destinations is not None:
                 self.destinations._make_temp_idx(
-                    start=max(self.origins.gdf.temp_idx.astype(int)) + 1
+                    start=max(self.origins.gdf["temp_idx"].astype(int)) + 1
                 )
 
         edges: list[tuple[str, str]] = self.network.get_edges()
@@ -1427,7 +1428,6 @@ def _get_edges_and_weights(
 
         self.network.gdf["src_tgt_wt"] = self.network._create_edge_ids(edges, weights)
 
-        # add edges between origins+destinations to the network nodes
         edges_start, weights_start = self.origins._get_edges_and_weights(
             nodes=self.network.nodes,
             rules=self.rules,
@@ -1484,7 +1484,7 @@ def _split_lines(self) -> None:
 
         points = points.drop_duplicates(points.geometry.name)
 
-        self.network.gdf["_meters2"] = self.network.gdf.length
+        self.network.gdf["meters_"] = self.network.gdf.length
 
         # create an id from before the split, used to revert the split later
         self.network.gdf["temp_idx__"] = range(len(self.network.gdf))
@@ -1497,39 +1497,16 @@ def _split_lines(self) -> None:
         )
 
         # save the unsplit lines for later
-        splitted = lines.loc[lines["splitted"] == 1]
+        splitted = lines.loc[lines["splitted"] == 1, "temp_idx__"]
         self.network._not_splitted = self.network.gdf.loc[
-            lambda x: x["temp_idx__"].isin(splitted["temp_idx__"])
+            self.network.gdf["temp_idx__"].isin(splitted)
         ]
 
-        new_lines, new_nodes = make_node_ids(splitted)
-        new_nodes = sfilter_inverse(new_nodes, self.network.nodes.buffer(1e-5))
-        new_nodes["node_id"] = (
-            new_nodes["node_id"].astype(int) + len(self.network.nodes) + 1
-        ).astype(str)
-        self.network._new_node_ids = list(new_nodes["node_id"])
-
         # adjust weight to new length
-        new_lines[self.rules.weight] = new_lines[self.rules.weight] * (
-            new_lines.length / new_lines["_meters2"]
-        )
-        self.network._nodes = pd.concat(
-            [self.network._nodes, new_nodes],
-            ignore_index=True,
-        )
-
-        lines = pd.concat(
-            [
-                self.network.gdf.loc[
-                    lambda x: ~x["temp_idx__"].isin(splitted["temp_idx__"])
-                ],
-                new_lines,
-            ],
-            ignore_index=True,
+        lines[self.rules.weight] = lines[self.rules.weight] * (
+            lines.length / lines["meters_"]
         )
 
-        lines = _map_node_ids_from_wkt(lines, self.network._nodes)
-
         self.network.gdf = lines
 
     def _unsplit_network(self):
@@ -1540,10 +1517,7 @@ def _unsplit_network(self):
         self.network.gdf = pd.concat(
             [lines, self.network._not_splitted], ignore_index=True
         ).drop("temp_idx__", axis=1)
-        self.network._nodes = self.network._nodes[
-            lambda x: ~x["node_id"].isin(self.network._new_node_ids)
-        ]
-        del self.network._not_splitted, self.network._new_node_ids
+        del self.network._not_splitted
 
     @staticmethod
     def _make_graph(
@@ -1587,7 +1561,7 @@ def _graph_is_up_to_date(self) -> bool:
         for points in ["origins", "destinations"]:
             if self[points] is None:
                 continue
-            if not hasattr(self, points) or self[points] is None:
+            if points not in self.wkts:
                 return False
             if self._points_have_changed(self[points].gdf, what=points):
                 return False
@@ -1603,7 +1577,7 @@ def _points_have_changed(self, points: GeoDataFrame, what: str) -> bool:
         if not np.array_equal(self.wkts[what], points.geometry.to_wkt().values):
             return True
 
-        if not all(x in self.graph.vs["name"] for x in list(points.temp_idx.values)):
+        if not all(x in self.graph.vs["name"] for x in list(points["temp_idx"].values)):
             return True
 
         return False
@@ -1617,6 +1591,8 @@ def _update_wkts(self) -> None:
         """
         self.wkts = {}
 
+        self.wkts["network"] = self.network.gdf.geometry.to_wkt().values
+
         if not hasattr(self, "origins"):
             return
 

src/sgis/networkanalysis/nodes.py

@@ -95,4 +95,6 @@ def _map_node_ids_from_wkt(lines, nodes, wkt: bool = True) -> GeoDataFrame:
     }
     lines["source"] = lines[geomcol1].map(id_dict)
     lines["target"] = lines[geomcol2].map(id_dict)
+    assert lines["source"].notna().all(), lines.loc[lines["source"].isna(), geomcol1]
+    assert lines["target"].notna().all(), lines.loc[lines["target"].isna(), geomcol2]
     return lines