7

When I want the distance between two points [(117.454361,38.8459879),(117.459880 ,38.846255)] (longitude,latitude) on the earth, I take the GeoSeries.distance method, but the method does not give me the right answer.

The real distance is about 479 meters, however the method give the result far from it, why?

Code:

import geopandas as gpd
from shapely.geometry import Point


geom=[Point(xy) for xy in zip([117.454361,117.459880],[38.8459879,38.846255])]
gdf=gpd.GeoDataFrame(geometry=geom,crs={'init':'epsg:4326'})
dis=gdf.distance(gdf.shift())
print(dis)
9

This is a distance in degrees, the coordinates of your data. I can get this using Pythagoras' theorem from your coordinates:

>>> a = [117.454361, 38.8459879]
>>> b = [117.45988, 38.846255]
>>> math.sqrt((a[0]-b[0])**2 + (a[1]-b[1])**2)
0.005525459565494833

to get a distance in metres, transform to a metric coordinate system, such as the UTM zone for your data.

Possible help here:

Calculate distance between a coordinate and a county in GeoPandas

| improve this answer | |
7
import geopandas as gpd
from shapely.geometry import Point
geom=[Point(xy) for xy in zip([117.454361,117.459880],[38.8459879,38.846255])]
gdf=gpd.GeoDataFrame(geometry=geom,crs={'init':'epsg:4326'})
gdf.to_crs(epsg=3310,inplace=True)
l=gdf.distance(gdf.shift())
print(l)

The result is 479.450134meters.

| improve this answer | |
  • 1
    What is crs 3310? Why are you using it? – larsks Jun 3 '19 at 14:33
  • 1
    I believe it's a conversion to UTM, a rectilinear projection coordinate system where distance can be calculated with pythagorean theorem. But I don't know why it's EPSG 3310 in particular. Is that a general solution for anywhere on Earth? – user2561747 Jul 31 '19 at 20:43
3

Pythagoras only works on a flat plane and not an sphere. The distance between two points on the surface of a sphere is found using great-circle distance:

enter image description here

where φ's are latitude and λ's are longitudes. To convert the distance to meter you need to know the radius of the sphere (6371km for Earth) and multiply it by Δσ in radians. Here is a code that does that:

def haversine(coord1, coord2):
    import math
    # Coordinates in decimal degrees (e.g. 2.89078, 12.79797)
    lon1, lat1 = coord1
    lon2, lat2 = coord2
    R = 6371000  # radius of Earth in meters
    phi_1 = math.radians(lat1)
    phi_2 = math.radians(lat2)

    delta_phi = math.radians(lat2 - lat1)
    delta_lambda = math.radians(lon2 - lon1)

    a = math.sin(delta_phi / 2.0) ** 2 + math.cos(phi_1) * math.cos(phi_2) * math.sin(delta_lambda / 2.0) ** 2

    c = 2 * math.atan2(math.sqrt(a), math.sqrt(1 - a))

    meters = R * c  # output distance in meters
    km = meters / 1000.0  # output distance in kilometers

    meters = round(meters)
    km = round(km, 3)
    print(f"Distance: {meters} m")
    print(f"Distance: {km} km")

and run it for the coordinates you get 479m.

haversine(coord1= (117.454361,38.8459879), coord2= (117.459880,38.846255))

The other way to this is to do what @CJ Xu did and convert the coordinates to UTM.

| improve this answer | |
2

Its pretty easy with Geopy

import geopy.distance
dist = geopy.distance.geodesic((38.8459879,117.454361),(38.846255,117.459880))
dist.meters

enter image description here

| improve this answer | |
  • Are the coordinates indise (lat, lon) or (lon, lat). Thanks! – M.O. Apr 22 at 3:59
  • 1
    its (lat, long) – virtuvious Apr 22 at 14:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.