是否可以在OSMnx中将坐标拟合到街道？

Question

OSMnx提供了解决方案来计算两个节点之间的最短路径，但是我想对街道上的点也是如此（我从车辆记录了GPS坐标）。我知道还有一种获取最接近节点的方法，但是我的这个问题有两个问题。

i）在计算最近节点时，是否也考虑了该点的街道？（我假设不是）ii）如果要实现类似的功能，我想知道一条街道（边缘）是如何表示为一条曲线的（也许是贝塞尔曲线）。是否可以获取边缘的曲线（或曲线方程）？

我在这里问这个问题，是因为OSMnx的编写准则对此提出了要求。

Answer 1

OSMnx中的街道和节点是shapely.geometry.LineString, 和shapely.geometry.Point对象，因此没有曲线，只有坐标序列。您所描述的技术术语是Map Matching。地图匹配有不同的方法，最简单的一种是几何地图匹配，您可以在其中找到距离 GPS 点最近的几何图形（节点或边缘）。point to point使用内置的osmnx函数可以轻松实现地图匹配ox.get_nearest_node()。如果您拥有密集的 GPS 轨迹，那么这种方法的效果相当好。对于point to line地图匹配，您必须使用形状函数。这种方法的问题是速度非常慢。您可以使用空间索引来加速算法，但对于大多数用途来说，它仍然不够快。请注意，几何地图匹配是所有方法中最不准确的。几周前我编写了一个函数，该函数使用可以从 OSMnx 获取的边缘 GeoDataFrame 和节点 GeoDataFrame 进行简单的点到线地图匹配。我放弃了这个想法，现在正在研究一种新算法（希望更快），完成后我将在 GitHub 上发布。同时，这可能对您或其他人有一些帮助，所以我将其发布在这里。这是废弃代码的早期版本，没有经过足够的测试，也没有优化。尝试一下，让我知道它是否适合您。

def GeoMM(traj, gdfn, gdfe):
"""
performs map matching on a given sequence of points

Parameters
----------

Returns
-------
list of tuples each containing timestamp, projected point to the line, the edge to which GPS point has been projected, the geometry of the edge))

"""

traj = pd.DataFrame(traj, columns=['timestamp', 'xy'])
traj['geom'] = traj.apply(lambda row: Point(row.xy), axis=1)
traj = gpd.GeoDataFrame(traj, geometry=traj['geom'], crs=EPSG3740)
traj.drop('geom', axis=1, inplace=True)

n_sindex = gdfn.sindex

res = []
for gps in traj.itertuples():
    tm = gps[1]
    p = gps[3]
    circle = p.buffer(150)
    possible_matches_index = list(n_sindex.intersection(circle.bounds))
    possible_matches = gdfn.iloc[possible_matches_index]
    precise_matches = possible_matches[possible_matches.intersects(circle)]
    candidate_nodes = list(precise_matches.index)

    candidate_edges = []
    for nid in candidate_nodes:
        candidate_edges.append(G.in_edges(nid))
        candidate_edges.append(G.out_edges(nid))

    candidate_edges = [item for sublist in candidate_edges for item in sublist]
    dist = []
    for edge in candidate_edges:
        # get the geometry
        ls = gdfe[(gdfe.u == edge[0]) & (gdfe.v == edge[1])].geometry
        dist.append([ls.distance(p), edge, ls])

    dist.sort()
    true_edge = dist[0][1]
    true_edge_geom = dist[0][2].item()
    pp = true_edge_geom.interpolate(true_edge_geom.project(p)) # projected point
    res.append((tm, pp, true_edge, true_edge_geom))


    return res