Найти длину между двумя точками на разных линиях - PullRequest
Новая
       12

Найти длину между двумя точками на разных линиях

0 голосов
/ 22 января 2020

У меня есть набор железных дорог как LineString объекты и железнодорожные остановки как Point объекты. Я хотел бы найти расстояние вдоль пути между двумя остановками, которые находятся в двух разных сегментах рельса.

В идеале алгоритм должен выглядеть следующим образом:

  1. найти, к чему LineString s точки принадлежат
  2. найти расстояние маршрута

Геометрия данных: 

LINESTRING (39.6578321 46.9014975, 39.6570969 46.9001473, 39.6559263 46.8989562, 39.6547141 46.8979608, 39.6334559 46.8814632, 39.631423 46.8798763)
LINESTRING (39.5993901 46.8638073, 39.5908007 46.8593408, 39.5779403 46.8525964, 39.5768336 46.8520311, 39.5751525 46.8511401, 39.5746447 46.8508623, 39.5738707 46.8504168, 39.5730826 46.849873, 39.5724067 46.8493054, 39.5719514 46.8488622, 39.5714961 46.848352, 39.571218 46.8480119, 39.5707563 46.8471837, 39.5705042 46.8467756, 39.5704733 46.8467123, 39.5701367 46.8460454, 39.5699599 46.8454946, 39.5697689 46.8447907, 39.5696772 46.844289, 39.5695147 46.8433865, 39.5690113 46.8372209, 39.5684811 46.830653, 39.5678923 46.8235858, 39.567313 46.8164827, 39.5664922 46.8059292, 39.5658968 46.7983574, 39.5648453 46.7849665, 39.5644591 46.7808009, 39.5641158 46.7759223, 39.5616257 46.7440695, 39.5614031 46.7412215, 39.5612064 46.7386438, 39.5609742 46.7360506, 39.5604696 46.7299663, 39.5604382 46.7295303, 39.5603741 46.728642, 39.5603751 46.7277074, 39.5604071 46.7263769, 39.5604906 46.7251804, 39.5734779 46.676293, 39.573848 46.6748994, 39.5790219 46.6550449)
LINESTRING (39.6312641 46.8797614, 39.6304511 46.8792256, 39.625757 46.8768124, 39.6199837 46.8739702)
LINESTRING (39.6198448 46.8739025, 39.614424 46.8712966, 39.6034935 46.8657577, 39.5995486 46.8638936)
POINT (39.6547141 46.8979608)
POINT (39.5609742 46.7360506)

Графически данные выглядят как

enter image description here

Что я должен сделать, чтобы найти расстояние между нижней точкой и верхней точкой? Они l ie на разных LineString с, и между ними есть другие LineString с.

Если лгать на том же LineString, я бы просто использовал shapely.split с каждой точкой, а затем вычел бы результат line.length от общей длины.

Но в таком случае я понятия не имею.

Ответы [ 2 ]

1 голос
/ 23 января 2020

В идеале, если все LineString имеют общие конечные точки, вы можете просто использовать linemerge, чтобы объединить их в одну LineString и найти расстояние вдоль трассы, как обычно, но в вашем В этом случае линии не соприкасаются друг с другом:

enter image description here

В этом случае я предлагаю поискать ближайшие пары LineString и присоединиться к ним. вместе, например, используя следующий код: 

from shapely.geometry import LineString
from shapely.wkt import loads


def merge_nontouching(lines):
    """Makes a LineString out of a list of LineString's by joining the closest pairs"""
    line, *candidates = lines
    while True:
        if not candidates:
            return line
        closest_candidate = min(candidates, key=line.distance)
        line = join_two(line, closest_candidate)
        candidates.pop(candidates.index(closest_candidate))


def join_two(line, other):
    """Joins two nontouching LineString's in one"""
    start, end = line.boundary
    if start.distance(other) > end.distance(other):
        return LineString([*line.coords, *other.coords])
    return LineString([*line.coords[::-1], *other.coords])


l1 = loads('LINESTRING (39.6578321 46.9014975, 39.6570969 46.9001473, 39.6559263 46.8989562, 39.6547141 46.8979608, 39.6334559 46.8814632, 39.631423 46.8798763)')
l2 = loads('LINESTRING (39.5993901 46.8638073, 39.5908007 46.8593408, 39.5779403 46.8525964, 39.5768336 46.8520311, 39.5751525 46.8511401, 39.5746447 46.8508623, 39.5738707 46.8504168, 39.5730826 46.849873, 39.5724067 46.8493054, 39.5719514 46.8488622, 39.5714961 46.848352, 39.571218 46.8480119, 39.5707563 46.8471837, 39.5705042 46.8467756, 39.5704733 46.8467123, 39.5701367 46.8460454, 39.5699599 46.8454946, 39.5697689 46.8447907, 39.5696772 46.844289, 39.5695147 46.8433865, 39.5690113 46.8372209, 39.5684811 46.830653, 39.5678923 46.8235858, 39.567313 46.8164827, 39.5664922 46.8059292, 39.5658968 46.7983574, 39.5648453 46.7849665, 39.5644591 46.7808009, 39.5641158 46.7759223, 39.5616257 46.7440695, 39.5614031 46.7412215, 39.5612064 46.7386438, 39.5609742 46.7360506, 39.5604696 46.7299663, 39.5604382 46.7295303, 39.5603741 46.728642, 39.5603751 46.7277074, 39.5604071 46.7263769, 39.5604906 46.7251804, 39.5734779 46.676293, 39.573848 46.6748994, 39.5790219 46.6550449)')
l3 = loads('LINESTRING (39.6312641 46.8797614, 39.6304511 46.8792256, 39.625757 46.8768124, 39.6199837 46.8739702)')
l4 = loads('LINESTRING (39.6198448 46.8739025, 39.614424 46.8712966, 39.6034935 46.8657577, 39.5995486 46.8638936)')
merged = merge_nontouching([l1, l2, l3, l4])
print(merged.wkt)
# LINESTRING (39.6578321 46.9014975, 39.6570969 46.9001473, 39.6559263 46.8989562, 39.6547141 46.8979608, 39.6334559 46.8814632, 39.631423 46.8798763, 39.6312641 46.8797614, 39.6304511 46.8792256, 39.625757 46.8768124, 39.6199837 46.8739702, 39.6198448 46.8739025, 39.614424 46.8712966, 39.6034935 46.8657577, 39.5995486 46.8638936, 39.5993901 46.8638073, 39.5908007 46.8593408, 39.5779403 46.8525964, 39.5768336 46.8520311, 39.5751525 46.8511401, 39.5746447 46.8508623, 39.5738707 46.8504168, 39.5730826 46.849873, 39.5724067 46.8493054, 39.5719514 46.8488622, 39.5714961 46.848352, 39.571218 46.8480119, 39.5707563 46.8471837, 39.5705042 46.8467756, 39.5704733 46.8467123, 39.5701367 46.8460454, 39.5699599 46.8454946, 39.5697689 46.8447907, 39.5696772 46.844289, 39.5695147 46.8433865, 39.5690113 46.8372209, 39.5684811 46.830653, 39.5678923 46.8235858, 39.567313 46.8164827, 39.5664922 46.8059292, 39.5658968 46.7983574, 39.5648453 46.7849665, 39.5644591 46.7808009, 39.5641158 46.7759223, 39.5616257 46.7440695, 39.5614031 46.7412215, 39.5612064 46.7386438, 39.5609742 46.7360506, 39.5604696 46.7299663, 39.5604382 46.7295303, 39.5603741 46.728642, 39.5603751 46.7277074, 39.5604071 46.7263769, 39.5604906 46.7251804, 39.5734779 46.676293, 39.573848 46.6748994, 39.5790219 46.6550449)

Отсюда вы можете просто использовать функцию project, чтобы получить расстояние: 

p1 = loads('POINT (39.6547141 46.8979608)')
p2 = loads('POINT (39.5609742 46.7360506)')
distance_along_track = abs(merged.project(p1) - merged.project(p2))
print(distance_along_track)
# 0.21041206903426987

Просто две заметки:
1) Это расстояние не в метрах или километрах. На каком-то шаге вам нужно будет выполнить преобразование из координат широты / долготы. К сожалению, я не знаю, как это сделать. Вы можете понять это сами, или я посмотрю и отредактирую это позже:)

2) Это может быть не самый эффективный подход, если у вас есть тысячи сегментов, но для данного примера этого будет достаточно.

0 голосов
/ 23 января 2020

Начиная с предоставленных данных. Вычисление происходит так: 

import shapely.wkt as wkt
import pyproj

wkt1 = 'LINESTRING (39.6578321 46.9014975, 39.6570969 46.9001473, 39.6559263 46.8989562, 39.6547141 46.8979608, 39.6334559 46.8814632, 39.631423 46.8798763)'
wkt2 = 'LINESTRING (39.5993901 46.8638073, 39.5908007 46.8593408, 39.5779403 46.8525964, 39.5768336 46.8520311, 39.5751525 46.8511401, 39.5746447 46.8508623, 39.5738707 46.8504168, 39.5730826 46.849873, 39.5724067 46.8493054, 39.5719514 46.8488622, 39.5714961 46.848352, 39.571218 46.8480119, 39.5707563 46.8471837, 39.5705042 46.8467756, 39.5704733 46.8467123, 39.5701367 46.8460454, 39.5699599 46.8454946, 39.5697689 46.8447907, 39.5696772 46.844289, 39.5695147 46.8433865, 39.5690113 46.8372209, 39.5684811 46.830653, 39.5678923 46.8235858, 39.567313 46.8164827, 39.5664922 46.8059292, 39.5658968 46.7983574, 39.5648453 46.7849665, 39.5644591 46.7808009, 39.5641158 46.7759223, 39.5616257 46.7440695, 39.5614031 46.7412215, 39.5612064 46.7386438, 39.5609742 46.7360506, 39.5604696 46.7299663, 39.5604382 46.7295303, 39.5603741 46.728642, 39.5603751 46.7277074, 39.5604071 46.7263769, 39.5604906 46.7251804, 39.5734779 46.676293, 39.573848 46.6748994, 39.5790219 46.6550449)'

# Create line objects from WKT code
line_1 = wkt.loads( wkt1 )
line_2 = wkt.loads( wkt2 )

# Get coordinates from specific vertices of lines
lon0, lat0 = line_1.coords[3]
lon1, lat1 = line_2.coords[32]

# Prep geodetic ref for computation
geod = pyproj.Geod(ellps='WGS84')

# Compute geodetic inverse problem
forward, back, dist = geod.inv(lon0, lat0, lon1, lat1)
print('Forward bearing: {}\nBackward bearing: {}\nDistance (m): {}'.format(forward, back, dist))

Выход: 

Forward bearing: -158.29035951647245
Backward bearing: 21.64128791894413
Distance (m): 19368.65144050848
...