如何高效找出OSMNX图中距坐标点集P1km内的所有节点？

直接暴力遍历所有节点与P中每个点计算距离的方法（O(N*M)时间复杂度）效率极低，针对44000个点的规模，推荐用KDTree空间索引实现批量范围查询，能将时间复杂度降到O(M log N)，大幅提升速度。

核心思路 #

1. 将图和坐标点转换为UTM投影（米为单位），确保欧氏距离等价于实际地面距离，避免经纬度直接计算的误差。
2. 基于图中所有节点的UTM坐标构建KDTree索引，支持快速的邻域查询。
3. 批量查询P中每个点1km范围内的节点，最后去重得到结果。

代码实现 #

import osmnx as ox
import numpy as np
from scipy.spatial import KDTree

# 假设已加载图G和坐标集合P
# 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), ...]

# 1. 将图投影到UTM坐标系（米为单位）
G_proj = ox.project_graph(G)

# 2. 提取节点的UTM坐标与对应的节点ID
nodes_gdf = ox.graph_to_gdfs(G_proj, edges=False)
node_coords = nodes_gdf[['x', 'y']].values
node_ids = nodes_gdf.index.values

# 3. 构建节点坐标的KDTree索引
kdtree = KDTree(node_coords)

# 4. 将P中的经纬度坐标转换为同UTM投影
crs = G_proj.graph['crs']
# 先把P转为GeoDataFrame再投影
p_gdf = ox.geocode_to_gdf(P, crs='EPSG:4326')
p_proj_gdf = ox.projection.project_gdf(p_gdf, to_crs=crs)
p_coords = p_proj_gdf[['x', 'y']].values

# 5. 批量查询每个P点1000米范围内的节点索引
neighbor_indices = kdtree.query_ball_point(p_coords, r=1000)

# 6. 去重并获取最终节点ID
unique_indices = set()
for indices in neighbor_indices:
    unique_indices.update(indices)
target_node_ids = node_ids[list(unique_indices)]

# 可选：获取目标节点的原始属性信息
target_nodes = {node_id: G.nodes[node_id] for node_id in target_node_ids}

关键说明 #

• UTM投影的必要性：经纬度是球面坐标，直接计算欧氏距离无法反映实际地面距离，转换为UTM后，坐标单位为米，欧氏距离就是真实的地面距离，保证1km范围查询的准确性。
• KDTree的优势：相对于暴力遍历，KDTree通过空间划分实现高效的邻域查询，当节点数量较大时（如纽约市的路网节点数可达数十万），效率提升极为明显。
• 去重处理：不同的P点可能会命中同一个节点，通过集合去重可避免重复结果。

内容的提问来源于stack exchange，提问作者Sihan Tawsik