2

I have an osmnx graph G and some co-ordinate points P. Now I want to find all the nodes which are within 1km Euclidean distance of at least one co-ordinate point of P. My code is as follows:

    import osmnx as ox
    G = ox.graph.graph_from_place('New York City, NewYork, United States', network_type="all", retain_all=True)
    P = [(40.718797266, -73.753347584), (40.713511106, -73.759968316), ...]

There are like 44k items in P Is there any efficient way that I can achieve my results?

Sihan Tawsik
  • 423
  • 3
  • 13

1 Answers1

1

This will find all the nodes within 1km of each of 44,000 points. After creating the graph, the node search completes in approx 16 seconds on my laptop and finds 57,134,000 nodes. The resulting GeoDataFrame is (non-uniquely) indexed by the position of the point in the original array of coordinates, and its values are the OSM IDs of all the graph nodes within 1km of that point.

import geopandas as gpd
import numpy as np
import osmnx as ox
ox.settings.log_console = True

# create graph and extract node geometries
G = ox.graph.graph_from_place('New York City, New York, United States', network_type="all", retain_all=True)
gdf_nodes = ox.graph_to_gdfs(G, edges=False)[['geometry']]

# create array of 44,000 points
coords = np.array([(40.718797266, -73.753347584), (40.713511106, -73.759968316)] * 22000)

# get 1km buffers around each point
points = gpd.GeoDataFrame(crs='epsg:4326', geometry=gpd.points_from_xy(coords[:, 1], coords[:, 0]))
buffers = ox.project_gdf(ox.project_gdf(points).buffer(1000), to_latlong=True)
gdf_buffers = gpd.GeoDataFrame(geometry=buffers)

# find all the nodes within the buffer of each point
result = gpd.sjoin(gdf_buffers, gdf_nodes, how='left', predicate='intersects')['index_right']
gboeing
  • 5,691
  • 2
  • 15
  • 41