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

Опитах се да внедря тази формула: http://andrew.hedges.name/experiments/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. Защо?


person gwaramadze    schedule 16.10.2013    source източник


Отговори (7)


Редактиране: Само като бележка, ако просто се нуждаете от бърз и лесен начин за намиране на разстоянието между две точки, силно препоръчвам да използвате подхода, описан в Отговорът на Kurt по-долу, вместо повторно внедряване на 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 км.

person Michael0x2a    schedule 16.10.2013
comment
това е вярно във всеки език за програмиране, а също и в диференциалното смятане. използването на степени е изключение и се използва само в човешка реч. - person bluesmoon; 17.10.2013
comment
Думата на мъдрите, тази формула изисква всички степени да са положителни. radians(abs(52.123)) трябва да свърши работа... - person Richard Dunn; 04.07.2017
comment
Сигурен ли си, че всички градуси (ъгли?) са положителни? Мисля, че това е грешно. Помислете дали lat1, lon1 = 10, 10 (градуси) и lat2, lon2 = -10, -10 (градуси). Чрез добавяне на abs() около градусите, разстоянието ще бъде нула, което е неправилно. Може би сте искали да вземете абсолютната стойност на dlon и/или dlat, но ако погледнете стойностите на dlon, dlat при изчисляването на a, синусът е четна функция, а косинусът на квадрат е четна функция, така че не вижте някаква полза от вземането на абсолютна стойност на dlat или dlon. - person Dave LeCompte; 24.05.2020
comment
Просто се чудя дали разстоянието по-горе е разстоянието на дъгата или разстоянието в равнината между две места? - person Punita Ojha; 07.07.2021

Актуализация: 04/2018: Обърнете внимание, че разстоянието Vincenty е отхвърлено от версията на GeoPy 1.13 - вместо това трябва да използвате geopy.distance.distance()!


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

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 или една от няколко други единици за разстояние).

person Kurt Peek    schedule 04.04.2017
comment
Благодаря. Можете ли да актуализирате отговора си с въпросните координати, които предоставих, вместо Нюпорт и Кливланд. Това ще даде по-добро разбиране на бъдещите читатели. - person gwaramadze; 05.04.2017
comment
Произволните местоположения на Нюпорт и Кливланд идват от примерната документация за geopy в списъка на PyPI: pypi.python.org/ pypi/geopy - person Jason Parham; 29.05.2017
comment
Трябваше да променя отговора на Kurt Peek на това: Изисква се капитализация: print geopy.distance.VincentyDistance(coords_1, coords_2).km 279.352901604 - person Jim; 14.06.2017
comment
Вероятно трябва да използвате geopy.distance.distance(…) в кода, който е псевдоним на най-добрата в момента (= най-точна) формула за разстояние. (Винсънти в момента.) - person mbirth; 13.01.2018
comment
Вероятно трябва да използвате geopy.distance.geodesic, тъй като разстоянието ще повиши ValueError, ако c1=(1, 179) и c2=(0,0). - person ramwin; 04.07.2018
comment
Това дава ли разстояние или изместване? Имам предвид взема ли предвид някаква вертикална ос? - person deadcode; 10.12.2018
comment
Използване на geopy.distance.vincenty в geopy-1.18.1 резултати: Vincenty е отхвърлен и ще бъде премахнат в geopy 2.0. Вместо това използвайте geopy.distance.geodesic (или стандартното geopy.distance.distance), което е по-точно и винаги се сближава. - person juanmah; 08.03.2019
comment
geopy.distance.great_circle ще работи два пъти по-бързо - person nda; 22.02.2020
comment
DeprecationWarning: Vincenty е отхвърлен и ще бъде премахнат в geopy 2.0. Вместо това използвайте geopy.distance.geodesic (или стандартното geopy.distance.distance), което е по-точно и винаги се сближава. #FYI - person えるまる; 23.04.2020
comment
Актуализирайте geopy 2.0.0: from geopy import distance и след това използвайте dist = distance.distance(coords_1, coords_2).km, което по подразбиране получава геодезичното разстояние. - person Colonel_Old; 23.09.2020

За хора (като мен), които идват тук чрез търсачката и просто търсят решение, което работи веднага, препоръчвам да инсталират mpu. Инсталирайте го чрез pip install mpu --user и го използвайте по този начин, за да получите разстоянието по хаверсинус:

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

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])
person Martin Thoma    schedule 04.07.2016
comment
Има ли начин да се промени зададената височина на една от точките? - person ; 06.12.2019
comment
Можете просто да добавите разликата във височината към разстоянието. Аз обаче не бих го направил. - person Martin Thoma; 05.05.2020
comment
Лион, Париж, 392.2172595594006 км, уау, последната цифра дори не е с размера на атом водород. Много точно! - person mins; 26.10.2020

Стигнах до много по-просто и стабилно решение, което използва 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

person Ramy M. Mousa    schedule 10.08.2019

Има няколко начина за изчисляване на разстоянието въз основа на координатите, т.е. географска ширина и дължина

Инсталирайте и импортирайте

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]

Използване на хаверсинус

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)

Използване на хаверсин със 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)

Използване на geopy

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
person Patel Romil    schedule 26.07.2020

Можете да използвате функцията H3,point_dist() на Uber, за да изчислите сферичното разстояние между две точки (lat, lng) . Можем да зададем единица за връщане („km“, „m“ или „rads“). Мерната единица по подразбиране е Km.

Пример:

import H3

coords_1 = (52.2296756, 21.0122287)
coords_2 = (52.406374, 16.9251681)
distance = h3.point_dist(coords_1,coords_2) #278.4584889328128

Надявам се това да е полезно!

person Ransaka Ravihara    schedule 23.02.2021

person    schedule
comment
Здравейте мислите ли, че има начин да се направи изчислението при получаване на данни директно от шаблона? - person Louis; 03.11.2020