Is it possible to fit a coordinate to a street in

Question

OSMnx provides solution to calculate the shortest path between two nodes, but I would like to the same with points on streets (I have GPS coordinates recorded from vehicles). I know there is also a method to get the closest node, but I have two question for this problem of mine.

i) When closest node computed is the street where the point is also taken into consideration? (I assume not) ii) If I wanted to implement something like this, I like to know how a street (edge) is represented as a curve (Bézier curve maybe?). Is it possible to get the curve (or the equation of the curve) of an edge?

I asked this question here, because the guidelines for contributing of OSMnx asked it.

Answer 1

Streets and node in OSMnx are shapely.geometry.LineString, and shapely.geometry.Point objects, so there is no curve, only sequence of coordinates. The technical term for what you described is Map Matching. There are different ways of map matching, the simplest one being geometric map matching in which you find nearest geometry (node or edge) to the GPS point. point to point map matching can be easily achieved using built-in osmnx function ox.get_nearest_node(). If you have a luxury of dense GPS tracks, this approach could work reasonably good. For point to line map matching you have to use shapely functions. The problem with this approach is that it is very slow. you can speed up the algorithm using spatial index, but still, it will not be fast enough for most purposes. Note that geometric map matching are least accurate among all approaches. I wrote a function a few weeks ago that does simple point to line map matching using edge GeoDataFrame and node GeoDataFrame that you can get from OSMnx. I abandoned this idea and now I am working on a new algorithm (hopefully much faster), which I will publish on GitHub upon completion. Meanwhile, this may be of some help for you or someone else, so I post it here. This is an early version of abandoned code, not tested enough and not optimized. give it a try and let me know if it works for you.

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