铁路网中矿区至港口最短路径计算问题求助

核心问题分析 #

你遇到的所有距离返回0的问题，根源在于错误使用了矿区/港口的原始坐标作为路径的起点/终点——这些点并不在你用momepy.gdf_to_nx生成的铁路网图节点中。NetworkX在找不到指定的源/目标节点时，部分场景下会默认返回0（或因节点不存在触发未捕获异常，但你的代码环境可能表现为返回0）。此外还有两个关键问题：

• 你已经计算出的nearest_rail_point并未被使用
• 使用WGS84（EPSG:4326）经纬度计算路径长度，结果是度数而非实际距离（米/公里），毫无实际意义

分步修复 #

1. 转换投影坐标系，确保距离计算正确 #

先将所有地理数据转换为适合澳大利亚的投影坐标系，比如EPSG:3112（GDA94 / Australian Albers），这样计算的长度单位是米：

# 定义投影坐标系
projected_crs = 'EPSG:3112'

# 转换所有GeoDataFrame到投影坐标系
ports = ports.to_crs(projected_crs)
mine_sites = mine_sites.to_crs(projected_crs)
rail_network = rail_network.to_crs(projected_crs)

2. 正确关联最近铁路节点（而非任意铁路点） #

nearest_points返回的是铁路线段上的最近点，但这个点不一定是图的节点。我们需要找到铁路网图中距离该点最近的节点：

# 先将铁路网转为图，提取所有节点坐标（元组格式）
G = momepy.gdf_to_nx(rail_network)
graph_nodes = list(G.nodes())

# 定义函数：找到图中距离目标点最近的节点
def get_nearest_graph_node(point, graph_nodes):
    min_dist = float('inf')
    nearest_node = None
    point_coords = (point.x, point.y)
    for node in graph_nodes:
        dist = ((node[0] - point_coords[0])**2 + (node[1] - point_coords[1])**2)**0.5
        if dist < min_dist:
            min_dist = dist
            nearest_node = node
    return nearest_node

# 为港口和矿区匹配最近的图节点
ports['nearest_rail_node'] = ports['nearest_rail_point'].apply(lambda x: get_nearest_graph_node(x, graph_nodes))
mine_sites['nearest_rail_node'] = mine_sites['nearest_rail_point'].apply(lambda x: get_nearest_graph_node(x, graph_nodes))

3. 使用正确的节点计算最短路径 #

现在用匹配到的铁路网节点作为起点和终点计算路径长度：

shortest_distances = {}

for _, mine_row in mine_sites.iterrows():
    mine_name = mine_row['mine_name']
    mine_start_node = mine_row['nearest_rail_node']
    shortest_distances[mine_name] = {}
    
    for _, port_row in ports.iterrows():
        port_name = port_row['port_name']
        port_end_node = port_row['nearest_rail_node']

        try:
            # 这里的weight='length'是momepy生成图时默认存储的线段长度字段
            shortest_distance = nx.shortest_path_length(G, 
                                                        source=mine_start_node, 
                                                        target=port_end_node, 
                                                        weight='length')
            # 转换为公里（可选）
            shortest_distances[mine_name][port_name] = round(shortest_distance / 1000, 2)
        except nx.NetworkXNoPath:
            shortest_distances[mine_name][port_name] = float('inf')
        except nx.NodeNotFound:
            # 处理节点不存在的情况
            shortest_distances[mine_name][port_name] = None

4. 额外优化：直接将矿区/港口点添加到铁路网图中 #

如果需要更精确的路径（从矿区最近点而非最近节点开始），可以将这些点添加到铁路网图中，连接到最近的铁路线段对应的节点：

# 示例：添加一个矿区点到图中
def add_point_to_graph(G, point, rail_gdf):
    # 找到距离点最近的铁路线段
    nearest_rail = rail_network.distance(point).idxmin()
    rail_segment = rail_network.loc[nearest_rail]
    # 获取线段的两个节点
    u = (rail_segment.geometry.coords[0][0], rail_segment.geometry.coords[0][1])
    v = (rail_segment.geometry.coords[-1][0], rail_segment.geometry.coords[-1][1])
    # 计算点到u和v的距离
    dist_u = point.distance(Point(u))
    dist_v = point.distance(Point(v))
    # 添加新节点
    new_node = (point.x, point.y)
    G.add_node(new_node)
    # 添加到两个端点的边
    G.add_edge(new_node, u, length=dist_u)
    G.add_edge(new_node, v, length=dist_v)
    return new_node

# 为每个矿区添加节点
mine_sites['graph_start_node'] = mine_sites['geometry'].apply(lambda x: add_point_to_graph(G, x, rail_network))
# 港口同理
ports['graph_end_node'] = ports['geometry'].apply(lambda x: add_point_to_graph(G, x, rail_network))

完整修正代码（整合版） #

import geopandas as gpd
from shapely.geometry import Point
from shapely.ops import nearest_points
import momepy
import networkx as nx

# Ports
ports = [
    {'port_name': 'Geraldton', 'geometry': Point(114.59786, -28.77688)},
    {'port_name': 'Bunbury', 'geometry': Point(115.673447, -33.318797)},
    {'port_name': 'Albany', 'geometry': Point(117.895025, -35.032831)}, 
    {'port_name': 'Esperance', 'geometry': Point(121.897114, -33.871834)}
]  

# Mine sites
mines = [
    {'mine_name': 'Gold', 'geometry': Point(117.94568, -34.93467)},
    {'mine_name': 'Silver', 'geometry': Point(115.16923, -29.65613)},
    {'mine_name': 'Bronze', 'geometry': Point(115.11039, -29.51287)}, 
    {'mine_name': 'Platinum', 'geometry': Point(115.11130, -29.42621)}
]

# Convert to GeoDataFrame
crs = 'EPSG:4326'
ports = gpd.GeoDataFrame(ports, crs=crs, geometry='geometry')
mine_sites = gpd.GeoDataFrame(mines, crs=crs, geometry='geometry')

# 读取铁路网数据（请替换为你的文件路径）
rail_network = gpd.read_file("railway-lines.shp")

# 转换到投影坐标系
projected_crs = 'EPSG:3112'
ports = ports.to_crs(projected_crs)
mine_sites = mine_sites.to_crs(projected_crs)
rail_network = rail_network.to_crs(projected_crs)

# 计算最近铁路点
ports['nearest_rail_point'] = ports['geometry'].apply(lambda x: nearest_points(x, rail_network.unary_union)[1])
mine_sites['nearest_rail_point'] = mine_sites['geometry'].apply(lambda x: nearest_points(x, rail_network.unary_union)[1])

# 转换铁路网为图
G = momepy.gdf_to_nx(rail_network)
graph_nodes = list(G.nodes())

# 找到最近图节点
def get_nearest_graph_node(point, graph_nodes):
    min_dist = float('inf')
    nearest_node = None
    point_coords = (point.x, point.y)
    for node in graph_nodes:
        dist = ((node[0] - point_coords[0])**2 + (node[1] - point_coords[1])**2)**0.5
        if dist < min_dist:
            min_dist = dist
            nearest_node = node
    return nearest_node

ports['nearest_rail_node'] = ports['nearest_rail_point'].apply(lambda x: get_nearest_graph_node(x, graph_nodes))
mine_sites['nearest_rail_node'] = mine_sites['nearest_rail_point'].apply(lambda x: get_nearest_graph_node(x, graph_nodes))

# 计算最短路径
shortest_distances = {}

for _, mine_row in mine_sites.iterrows():
    mine_name = mine_row['mine_name']
    mine_start = mine_row['nearest_rail_node']
    shortest_distances[mine_name] = {}
    
    for _, port_row in ports.iterrows():
        port_name = port_row['port_name']
        port_end = port_row['nearest_rail_node']

        try:
            dist = nx.shortest_path_length(G, source=mine_start, target=port_end, weight='length')
            shortest_distances[mine_name][port_name] = round(dist / 1000, 2)  # 转为公里
        except nx.NetworkXNoPath:
            shortest_distances[mine_name][port_name] = float('inf')
        except nx.NodeNotFound:
            shortest_distances[mine_name][port_name] = None

# 打印结果
for mine, port_dist in shortest_distances.items():
    print(f"{mine}矿区到各港口的最短铁路距离（公里）：")
    for port, dist in port_dist.items():
        print(f"- {port}: {dist}")

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