В поисках питонического способа расчета длины линии WKT

13

Я был совершенно недоволен вычислением длины строк в WGS84 в милях . Это заставило меня задуматься, есть ли более удобный, Pythonic способ вычислить длину строки WKT в соответствии с данным SRID.

Я имею в виду что-то вроде:

srid="WGS84"
line="LINESTRING(3.0 4.0, 3.1 4.1)"
print length(line, srid)

Я ищу точный ответ, а не sin\cosприближения.

Есть идеи?

Адам Матан
источник
Tomkralidis, это веб-сайт ГИС. Ваш ответ не учитывает, что это расстояние между геопространственными координатами (посмотрите SRID). Сам по себе стройный не может вычислить геопространственные расстояния, так как не имеет представления о проекции карты.

Ответы:

18

Модуль геопсии предоставляет формулу Винсенти , которая обеспечивает точные эллипсоидальные расстояния. wktСоедините это с загрузкой в ​​Shapely, и вы получите достаточно простой код:

from geopy import distance
from shapely.wkt import loads

line_wkt="LINESTRING(3.0 4.0, 3.1 4.1)"

# a number of other elipsoids are supported
distance.VincentyDistance.ELLIPSOID = 'WGS-84'
d = distance.distance

line = loads(line_wkt)

# convert the coordinates to xy array elements, compute the distance
dist = d(line.xy[0], line.xy[1])

print dist.meters
SCW
источник
1
+1, было бы +10, если бы я мог. Спасла моя команда часы программирования.
Адам Матан
Отличается ли этот подход от ответа @tomkralidis, если входные координаты уже есть в WGS-84?
LarsVegas
1
@LarsVegas да, Shapely обрабатывает только плоские координаты - поэтому он будет точно измерять расстояния в проецируемом пространстве, но не географически (например, WGS-1984).
Scw
4

Вы также можете использовать свойство длины Shapely , то есть:

from shapely.wkt import loads

l=loads('LINESTRING(3.0 4.0, 3.1 4.1)')
print l.length
tomkralidis
источник
Обратите внимание, что длина для этого конкретного примера будет бессмысленной, так как это географическая система координат (WGS84).
Майк Т
2

Поздно на вечеринку, но с полезным вкладом. Основываясь на ответе scw с использованием geopy, я написал небольшую функцию, которая выполняет расчет для стройного объекта LineString с произвольным количеством координат. Он использует pairsитератор из Stackoverflow.

Основная особенность: строки документов намного длиннее фрагментов.

def line_length(line):
    """Calculate length of a line in meters, given in geographic coordinates.
    Args:
        line: a shapely LineString object with WGS 84 coordinates
    Returns:
        Length of line in meters
    """
    # Swap shapely (lonlat) to geopy (latlon) points
    latlon = lambda lonlat: (lonlat[1], lonlat[0])
    total_length = sum(distance(latlon(a), latlon(b)).meters
                       for (a, b) in pairs(line.coords))
    return round(total_length, 0)


def pairs(lst):
    """Iterate over a list in overlapping pairs without wrap-around.

    Args:
        lst: an iterable/list

    Returns:
        Yields a pair of consecutive elements (lst[k], lst[k+1]) of lst. Last 
        call yields the last two elements.

    Example:
        lst = [4, 7, 11, 2]
        pairs(lst) yields (4, 7), (7, 11), (11, 2)

    Source:
        /programming/1257413/1257446#1257446
    """
    i = iter(lst)
    prev = i.next()
    for item in i:
        yield prev, item
        prev = item
ojdo
источник
1
Это неверно: geopy.distance.distance принимает координаты в (y, x), но стройная линейная строка представляет собой «упорядоченную последовательность из 2 или более (x, y [, z])», поэтому необходимо использовать вспомогательную функцию geopy lonlat () ,
Мартин Берч
@MartinBurch: ой, ты прав. Противная вещь даже не та [, z], но обмен аргументами (y, x)на (x, y)это необходим. Спасибо, что заметили это. Можете ли вы заметить, выглядит ли это редактирование менее подверженным ошибкам?
19