diff --git a/aequilibrae/project/network/connector_creation.py b/aequilibrae/project/network/connector_creation.py index 1b3b2fd67..c7545261c 100644 --- a/aequilibrae/project/network/connector_creation.py +++ b/aequilibrae/project/network/connector_creation.py @@ -1,116 +1,61 @@ -from math import pi, sqrt from sqlite3 import Connection from typing import Optional -import numpy as np -from scipy.cluster.vq import kmeans2, whiten -from scipy.spatial.distance import cdist -import shapely.wkb -from shapely.geometry import LineString +from shapely.geometry import LineString, Polygon from aequilibrae.utils.db_utils import commit_and_close INFINITE_CAPACITY = 99999 def connector_creation( - geo, zone_id: int, srid: int, mode_id: str, network, link_types="", connectors=1, conn_: Optional[Connection] = None + zone_id: int, + mode_id: str, + network, + link_types="", + connectors=1, + conn_: Optional[Connection] = None, + delimiting_area: Polygon = None, ): if len(mode_id) > 1: raise Exception("We can only add centroid connectors for one mode at a time") - with conn_ or commit_and_close(network.project.connect()) as conn: + with conn_ or commit_and_close(network.project.path_to_file) as conn: logger = network.project.logger - if conn.execute("select count(*) from nodes where node_id=?", [zone_id]).fetchone() is None: + if sum(conn.execute("select count(*) from nodes where node_id=?", [zone_id]).fetchone()) == 0: logger.warning("This centroid does not exist. Please create it first") return - proj_nodes = network.nodes - node = proj_nodes.get(zone_id) sql = "select count(*) from links where a_node=? and instr(modes,?) > 0" if conn.execute(sql, [zone_id, mode_id]).fetchone()[0] > 0: logger.warning("Mode is already connected") return - if len(link_types) > 0: - lt = f"*[{link_types}]*" - else: - lt = "".join([x[0] for x in conn.execute("Select link_type_id from link_types").fetchall()]) - lt = f"*[{lt}]*" + proj_nodes = network.nodes.data - sql = """select node_id, ST_asBinary(geometry), modes, link_types from nodes where ST_Within(geometry, GeomFromWKB(?, ?)) and - (nodes.rowid in (select rowid from SpatialIndex where f_table_name = 'nodes' and - search_frame = GeomFromWKB(?, ?))) - and link_types glob ? and instr(modes, ?)>0""" + centroid = proj_nodes.query("node_id == @zone_id") # type: gpd.GeoDataFrame + centroid = centroid.rename(columns={"node_id": "zone_id"})[["zone_id", "geometry"]] - # We expand the area by its average radius until it is 20 times - # beginning with a strict search within the zone - buffer = 0 - increase = sqrt(geo.area / pi) - dt = [] - while dt == [] and buffer <= increase * 10: - wkb = geo.buffer(buffer).wkb - dt = conn.execute(sql, [wkb, srid, wkb, srid, lt, mode_id]).fetchall() - buffer += increase + nodes = proj_nodes.query("is_centroid != 1 and modes.str.contains(@mode_id)", engine="python") - if buffer > increase: - msg = f"Could not find node inside zone {zone_id}. Search area was expanded until we found a suitable node" - logger.warning(msg) - if dt == []: - logger.warning( - f"FAILED! Could not find suitable nodes to connect within 5 times the diameter of zone {zone_id}." - ) - return - - coords = [] - nodes = [] - for node_id, wkb, modes, link_types in dt: - geo = shapely.wkb.loads(wkb) - coords.append([geo.x, geo.y]) - nodes.append(node_id) - - num_connectors = connectors - if len(nodes) == 0: - raise Exception("We could not find any candidate nodes that satisfied your criteria") - elif len(nodes) < connectors: - logger.warning( - f"We have fewer possible nodes than required connectors for zone {zone_id}. Will connect all of them." - ) - num_connectors = len(nodes) - - if num_connectors == len(coords): - all_nodes = nodes - else: - features = np.array(coords) - whitened = whiten(features) - centroids, allocation = kmeans2(whitened, num_connectors) + if len(link_types) > 0: + nodes = nodes[nodes.link_types.str.contains("|".join(list(link_types)))] - all_nodes = set() - for i in range(num_connectors): - nds = [x for x, y in zip(nodes, list(allocation)) if y == i] - centr = centroids[i] - positions = [x for x, y in zip(whitened, allocation) if y == i] - if positions: - dist = cdist(np.array([centr]), np.array(positions)).flatten() - node_to_connect = nds[dist.argmin()] - all_nodes.add(node_to_connect) + if delimiting_area is not None: + nodes = nodes[nodes.geometry.within(delimiting_area)] - nds = list(all_nodes) - data = [zone_id] + nds - sql = f'select b_node from links where a_node=? and b_node in ({",".join(["?"] * len(nds))})' - data = [x[0] for x in conn.execute(sql, data).fetchall()] + if nodes.empty: + zone_id = centroid["zone_id"].values[0] + logger.warning(f"No nodes found for centroid {zone_id} (mode {mode_id} and link types {link_types})") + return - if data: - qry = ",".join(["?"] * len(data)) - dt = [mode_id, zone_id] + data - conn.execute(f"Update links set modes=modes || ? where a_node=? and b_node in ({qry})", dt) - nds = [x for x in nds if x not in data] - logger.warning(f"Mode {mode_id} added to {len(data)} existing centroid connectors for zone {zone_id}") + joined = nodes[["node_id", "geometry"]].sjoin_nearest(centroid, distance_col="distance_connector") + joined = joined.nsmallest(connectors, "distance_connector") + centr_geo = centroid.geometry.values[0] links = network.links - for node_to_connect in nds: + for _, rec in joined.iterrows(): link = links.new() - node_to = proj_nodes.get(node_to_connect) - link.geometry = LineString([node.geometry, node_to.geometry]) + link.geometry = LineString([centr_geo, rec.geometry]) link.modes = mode_id link.direction = 0 link.link_type = "centroid_connector" @@ -118,5 +63,5 @@ def connector_creation( link.capacity_ab = INFINITE_CAPACITY link.capacity_ba = INFINITE_CAPACITY link.save() - if nds: - logger.warning(f"{len(nds)} new centroid connectors for mode {mode_id} added for centroid {zone_id}") + if not joined.empty: + logger.warning(f"{joined.shape[0]} new centroid connectors for mode {mode_id} added for centroid {zone_id}") diff --git a/aequilibrae/project/network/node.py b/aequilibrae/project/network/node.py index dc181598f..c4a653147 100644 --- a/aequilibrae/project/network/node.py +++ b/aequilibrae/project/network/node.py @@ -113,7 +113,14 @@ def __save_existing_node(self): sql = f"Update Nodes set {txts}" return data, sql - def connect_mode(self, area: Polygon, mode_id: str, link_types="", connectors=1, conn: Optional[Connection] = None): + def connect_mode( + self, + mode_id: str, + link_types="", + connectors=1, + conn: Optional[Connection] = None, + area: Optional[Polygon] = None, + ): """Adds centroid connectors for the desired mode to the network file Centroid connectors are created by connecting the zone centroid to one or more nodes selected from @@ -129,7 +136,6 @@ def connect_mode(self, area: Polygon, mode_id: str, link_types="", connectors=1, If fewer candidates than required connectors are found, all candidates are connected. :Arguments: - **area** (:obj:`Polygon`): Initial area where AequilibraE will look for nodes to connect **mode_id** (:obj:`str`): Mode ID we are trying to connect @@ -137,20 +143,21 @@ def connect_mode(self, area: Polygon, mode_id: str, link_types="", connectors=1, be considered. eg: yCdR. Defaults to ALL link types **connectors** (:obj:`int`, *Optional*): Number of connectors to add. Defaults to 1 + + **area** (:obj:`Polygon`, *Optional*): Area limiting the search for connectors """ if self.is_centroid != 1 or self.__original__["is_centroid"] != 1: self._logger.warning("Connecting a mode only makes sense for centroids and not for regular nodes") return connector_creation( - area, - self.node_id, - self.__srid__, - mode_id, + zone_id=self.node_id, + mode_id=mode_id, link_types=link_types, connectors=connectors, network=self.project.network, conn_=conn, + delimiting_area=area, ) def __setattr__(self, instance, value) -> None: diff --git a/aequilibrae/project/zone.py b/aequilibrae/project/zone.py index fd308362e..612f9ec7d 100644 --- a/aequilibrae/project/zone.py +++ b/aequilibrae/project/zone.py @@ -90,18 +90,16 @@ def add_centroid(self, point: Point, robust=True) -> None: data = [self.zone_id, point.wkb, self.__srid__] conn.execute(sql, data) - def connect_mode(self, mode_id: str, link_types="", connectors=1, conn: Optional[Connection] = None) -> None: + def connect_mode( + self, mode_id: str, link_types="", connectors=1, conn: Optional[Connection] = None, limit_to_zone=True + ) -> None: """Adds centroid connectors for the desired mode to the network file Centroid connectors are created by connecting the zone centroid to one or more nodes selected from all those that satisfy the mode and link_types criteria and are inside the zone. - The selection of the nodes that will be connected is done simply by computing running the - `KMeans2 `_ - clustering algorithm from SciPy and selecting the nodes closest to each cluster centroid. - - When there are no node candidates inside the zone, the search area is progressively expanded until - at least one candidate is found. + The selection of the nodes that will be connected is done simply by searching for the node closest to the + zone centroid, or the N closest nodes to the centroid. If fewer candidates than required connectors are found, all candidates are connected. @@ -112,16 +110,21 @@ def connect_mode(self, mode_id: str, link_types="", connectors=1, conn: Optional eg: yCdR. Defaults to ALL link types **connectors** (:obj:`int`, *Optional*): Number of connectors to add. Defaults to 1 + + **conn** (:obj:`sqlite3.Connection`, *Optional*): Connection to the database. + + **limit_to_zone** (:obj:`bool`): Limits the search for nodes inside the zone. Defaults to ``True``. """ + + area = self.geometry if limit_to_zone else None connector_creation( - self.geometry, zone_id=self.zone_id, - srid=self.__srid__, mode_id=mode_id, link_types=link_types, connectors=connectors, network=self.project.network, conn_=conn, + delimiting_area=area, ) def disconnect_mode(self, mode_id: str) -> None: diff --git a/docs/source/examples/creating_models/plot_create_zoning.py b/docs/source/examples/creating_models/plot_create_zoning.py index e8f03afad..d2843dea8 100644 --- a/docs/source/examples/creating_models/plot_create_zoning.py +++ b/docs/source/examples/creating_models/plot_create_zoning.py @@ -153,8 +153,7 @@ # %% # When connecting a centroid not associated with a zone, we need to tell AequilibraE what is the initial area around # the centroid that needs to be considered when looking for candidate nodes. -# Distance here is in degrees, so 0.01 is equivalent to roughly 1.1km -airport.connect_mode(airport.geometry.buffer(0.01), mode_id="c", link_types="ytrusP", connectors=1) +airport.connect_mode(mode_id="c", link_types="ytrusP", connectors=1) # %% project.close() diff --git a/tests/aequilibrae/project/test_zone.py b/tests/aequilibrae/project/test_zone.py index 61906728c..6896ede23 100644 --- a/tests/aequilibrae/project/test_zone.py +++ b/tests/aequilibrae/project/test_zone.py @@ -154,9 +154,8 @@ def test_connect_mode(self): self.assertIsNot(0, curr.fetchone()[0], "failed to add connectors for mode t") # Cannot connect a centroid that does not exist - with self.assertRaises(ValueError): - zone2 = zones.get(2) - zone2.connect_mode("c") + zone2 = zones.get(2) + zone2.connect_mode("c", conn=self.proj.conn) def test_disconnect_mode(self): self.__change_project()