Получение расстояния между двумя точками на основе широты / долготы

158

Я попытался реализовать эту формулу: http://andrew.hedges.name/experiment/haversine/ Апплет хорошо подходит для двух пунктов, которые я тестирую:

введите описание изображения здесь

Все же мой код не работает.

from math import sin, cos, sqrt, atan2

R = 6373.0

lat1 = 52.2296756
lon1 = 21.0122287
lat2 = 52.406374
lon2 = 16.9251681

dlon = lon2 - lon1
dlat = lat2 - lat1
a = (sin(dlat/2))**2 + cos(lat1) * cos(lat2) * (sin(dlon/2))**2
c = 2 * atan2(sqrt(a), sqrt(1-a))
distance = R * c

print "Result", distance
print "Should be", 278.546

Расстояние, которое он возвращает, составляет 5447,05546147 . Зачем?

gwaramadze
источник

Ответы:

207

Редактировать: просто в качестве примечания: если вам просто нужен быстрый и простой способ определения расстояния между двумя точками, я настоятельно рекомендую использовать подход, описанный в ответе Курта ниже, вместо повторной реализации Haversine - см. Его пост для обоснования.

Этот ответ сфокусирован только на ответе на конкретную ошибку, с которой столкнулся OP.


Это потому, что в Python все триггерные функции используют радианы , а не градусы.

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

from math import sin, cos, sqrt, atan2, radians

# approximate radius of earth in km
R = 6373.0

lat1 = radians(52.2296756)
lon1 = radians(21.0122287)
lat2 = radians(52.406374)
lon2 = radians(16.9251681)

dlon = lon2 - lon1
dlat = lat2 - lat1

a = sin(dlat / 2)**2 + cos(lat1) * cos(lat2) * sin(dlon / 2)**2
c = 2 * atan2(sqrt(a), sqrt(1 - a))

distance = R * c

print("Result:", distance)
print("Should be:", 278.546, "km")

Расстояние теперь возвращает правильное значение 278.545589351км.

Michael0x2a
источник
13
это верно для любого языка программирования, а также в дифференциальном исчислении. использование степеней является исключением и используется только в человеческой речи.
Bluesmoon
11
Говоря мудрым, эта формула требует, чтобы все степени были положительными. radians(abs(52.123))должен сделать трюк ...
Ричард Данн
1
Вы уверены, что все градусы (углы?) Положительны? Я думаю, что это неправильно. Рассмотрим, если lat1, lon1 = 10, 10 (градусы) и lat2, lon2 = -10, -10 (градусы). Добавляя abs () вокруг градусов, расстояние будет равно нулю, что неверно. Возможно, вы хотели взять абсолютное значение dlon и / или dlat, но если вы посмотрите на значения dlon, значения dlat в расчете a, синус - четная функция, а квадрат косинуса - четная функция, поэтому я не буду увидеть любую выгоду от принятия абсолютного значения dlat или dlon, либо.
Дейв ЛеКомпт
238

Обновление: 04/2018: обратите внимание, что расстояние Vincenty устарело с версии GeoPy 1.13 - вместо этого вы должны использовать geopy.distance.distance ()!


Ответы выше основаны на формуле Хаверсайна , которая предполагает, что земля является сферой, что приводит к ошибкам примерно до 0,5% (согласно help(geopy.distance)). Расстояние Винсенти использует более точные эллипсоидальные модели, такие как WGS-84 , и реализовано в геопсии . Например,

import geopy.distance

coords_1 = (52.2296756, 21.0122287)
coords_2 = (52.406374, 16.9251681)

print geopy.distance.vincenty(coords_1, coords_2).km

напечатает расстояние в 279.352901604километрах, используя эллипсоид по умолчанию WGS-84. (Вы также можете выбрать .milesодин или несколько других единиц измерения расстояния).

Курт Пик
источник
1
Спасибо. Можете ли вы обновить свой ответ с помощью координат, которые я указывал в вопросе вместо Ньюпорта и Кливленда. Это даст лучшее понимание будущим читателям.
Гварамадзе
1
Произвольные местоположения Ньюпорта и Кливленда взяты из примерной документации по геопее в листинге PyPI: pypi.python.org/pypi/geopy
Джейсон Пархам,
Я должен был изменить ответ Курта Пика на это: Требуется капитализация:print geopy.distance.VincentyDistance(coords_1, coords_2).km 279.352901604
Джим
4
Вам, вероятно, следует использовать geopy.distance.distance(…)в коде, который является псевдонимом наилучшей (= наиболее точной) формулы расстояния в настоящее время. (Винсенти на данный момент.)
Рождение
10
Использование geopy.distance.vincenty в выходных данных geopy-1.18.1: Vincenty устарела и будет удалена в geopy 2.0. Вместо этого используйте geopy.distance.geodesic(или значение по умолчанию geopy.distance.distance), что является более точным и всегда сходится.
Джуанма
88

Для людей (таких как я), которые приходят сюда через поисковик и просто ищут решение, которое работает «из коробки», я рекомендую установить mpu. Установите это через pip install mpu --userи используйте это как это, чтобы получить расстояние haversine :

import mpu

# Point one
lat1 = 52.2296756
lon1 = 21.0122287

# Point two
lat2 = 52.406374
lon2 = 16.9251681

# What you were looking for
dist = mpu.haversine_distance((lat1, lon1), (lat2, lon2))
print(dist)  # gives 278.45817507541943.

Альтернативный пакет есть gpxpy.

Если вы не хотите зависимостей, вы можете использовать:

import math


def distance(origin, destination):
    """
    Calculate the Haversine distance.

    Parameters
    ----------
    origin : tuple of float
        (lat, long)
    destination : tuple of float
        (lat, long)

    Returns
    -------
    distance_in_km : float

    Examples
    --------
    >>> origin = (48.1372, 11.5756)  # Munich
    >>> destination = (52.5186, 13.4083)  # Berlin
    >>> round(distance(origin, destination), 1)
    504.2
    """
    lat1, lon1 = origin
    lat2, lon2 = destination
    radius = 6371  # km

    dlat = math.radians(lat2 - lat1)
    dlon = math.radians(lon2 - lon1)
    a = (math.sin(dlat / 2) * math.sin(dlat / 2) +
         math.cos(math.radians(lat1)) * math.cos(math.radians(lat2)) *
         math.sin(dlon / 2) * math.sin(dlon / 2))
    c = 2 * math.atan2(math.sqrt(a), math.sqrt(1 - a))
    d = radius * c

    return d


if __name__ == '__main__':
    import doctest
    doctest.testmod()

Другой альтернативный пакет [haversine][1]

from haversine import haversine, Unit

lyon = (45.7597, 4.8422) # (lat, lon)
paris = (48.8567, 2.3508)

haversine(lyon, paris)
>> 392.2172595594006  # in kilometers

haversine(lyon, paris, unit=Unit.MILES)
>> 243.71201856934454  # in miles

# you can also use the string abbreviation for units:
haversine(lyon, paris, unit='mi')
>> 243.71201856934454  # in miles

haversine(lyon, paris, unit=Unit.NAUTICAL_MILES)
>> 211.78037755311516  # in nautical miles

Они утверждают, что имеют оптимизацию производительности для расстояний между всеми точками в двух векторах

from haversine import haversine_vector, Unit

lyon = (45.7597, 4.8422) # (lat, lon)
paris = (48.8567, 2.3508)
new_york = (40.7033962, -74.2351462)

haversine_vector([lyon, lyon], [paris, new_york], Unit.KILOMETERS)

>> array([ 392.21725956, 6163.43638211])
Мартин Тома
источник
Есть ли способ изменить данный Highet одного из пунктов?
Йовель Cohen
Вы можете просто добавить разницу в высоте к расстоянию. Я бы не стал этого делать.
Мартин Тома
16

Я прибыл в гораздо более простым и надежным решением , которое использует geodesicиз geopyпакета , так как вы будете весьма вероятно , использовать его в вашем проекте в любом случае так без дополнительной установки пакета необходимо.

Вот мое решение:

from geopy.distance import geodesic


origin = (30.172705, 31.526725)  # (latitude, longitude) don't confuse
dist = (30.288281, 31.732326)

print(geodesic(origin, dist).meters)  # 23576.805481751613
print(geodesic(origin, dist).kilometers)  # 23.576805481751613
print(geodesic(origin, dist).miles)  # 14.64994773134371

geopy

Рами М. Муса
источник
5
import numpy as np


def Haversine(lat1,lon1,lat2,lon2, **kwarg):
    """
    This uses the ‘haversine’ formula to calculate the great-circle distance between two points – that is, 
    the shortest distance over the earth’s surface – giving an ‘as-the-crow-flies’ distance between the points 
    (ignoring any hills they fly over, of course!).
    Haversine
    formula:    a = sin²(Δφ/2) + cos φ1 ⋅ cos φ2 ⋅ sin²(Δλ/2)
    c = 2 ⋅ atan2( √a, √(1−a) )
    d = R ⋅ c
    where   φ is latitude, λ is longitude, R is earth’s radius (mean radius = 6,371km);
    note that angles need to be in radians to pass to trig functions!
    """
    R = 6371.0088
    lat1,lon1,lat2,lon2 = map(np.radians, [lat1,lon1,lat2,lon2])

    dlat = lat2 - lat1
    dlon = lon2 - lon1
    a = np.sin(dlat/2)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2) **2
    c = 2 * np.arctan2(a**0.5, (1-a)**0.5)
    d = R * c
    return round(d,4)
Sudhirln92
источник
0

Есть несколько способов рассчитать расстояние на основе координат, то есть широты и долготы

Установить и импортировать

from geopy import distance
from math import sin, cos, sqrt, atan2, radians
from sklearn.neighbors import DistanceMetric
import osrm
import numpy as np

Определить координаты

lat1, lon1, lat2, lon2, R = 20.9467,72.9520, 21.1702, 72.8311, 6373.0
coordinates_from = [lat1, lon1]
coordinates_to = [lat2, lon2]

Использование haversine

dlon = radians(lon2) - radians(lon1)
dlat = radians(lat2) - radians(lat1)
    
a = sin(dlat / 2)**2 + cos(lat1) * cos(lat2) * sin(dlon / 2)**2
c = 2 * atan2(sqrt(a), sqrt(1 - a))
    
distance_haversine_formula = R * c
print('distance using haversine formula: ', distance_haversine_formula)

Использование haversine со sklearn

dist = DistanceMetric.get_metric('haversine')
    
X = [[radians(lat1), radians(lon1)], [radians(lat2), radians(lon2)]]
distance_sklearn = R * dist.pairwise(X)
print('distance using sklearn: ', np.array(distance_sklearn).item(1))

Использование OSRM

osrm_client = osrm.Client(host='http://router.project-osrm.org')
coordinates_osrm = [[lon1, lat1], [lon2, lat2]] # note that order is lon, lat
    
osrm_response = osrm_client.route(coordinates=coordinates_osrm, overview=osrm.overview.full)
dist_osrm = osrm_response.get('routes')[0].get('distance')/1000 # in km
print('distance using OSRM: ', dist_osrm)

Использование геопии

distance_geopy = distance.distance(coordinates_from, coordinates_to).km
print('distance using geopy: ', distance_geopy)
    
distance_geopy_great_circle = distance.great_circle(coordinates_from, coordinates_to).km 
print('distance using geopy great circle: ', distance_geopy_great_circle)

Вывод

distance using haversine formula:  26.07547017310917
distance using sklearn:  27.847882224769783
distance using OSRM:  33.091699999999996
distance using geopy:  27.7528030550408
distance using geopy great circle:  27.839182219511834
Патель Ромил
источник