Haversine Formula in Python (Bearing and Distance Between Two Gps Points)

Haversine Formula in Python (Bearing and Distance between two GPS points)

Here's a Python version:

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

def haversine(lon1, lat1, lon2, lat2):
    """
    Calculate the great circle distance in kilometers between two points 
    on the earth (specified in decimal degrees)
    """
    # convert decimal degrees to radians 
    lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])

    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
    c = 2 * asin(sqrt(a)) 
    r = 6371 # Radius of earth in kilometers. Use 3956 for miles. Determines return value units.
    return c * r

Getting distance between two points based on latitude/longitude

Edit: Just as a note, if you just need a quick and easy way of finding the distance between two points, I strongly recommend using the approach described in Kurt's answer below instead of re-implementing Haversine -- see his post for rationale.

This answer focuses just on answering the specific bug OP ran into.

---

It's because in Python, all the trig functions use radians, not degrees.

You can either convert the numbers manually to radians, or use the radians function from the math module:

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")

The distance is now returning the correct value of 278.545589351 km.

Calculate distance between two latitude-longitude points? (Haversine formula)

This link might be helpful to you, as it details the use of the Haversine formula to calculate the distance.

Excerpt:

This script [in Javascript] calculates great-circle distances between the two points –
that is, the shortest distance over the earth’s surface – using the
‘Haversine’ formula.

function getDistanceFromLatLonInKm(lat1,lon1,lat2,lon2) {
  var R = 6371; // Radius of the earth in km
  var dLat = deg2rad(lat2-lat1);  // deg2rad below
  var dLon = deg2rad(lon2-lon1); 
  var a = 
    Math.sin(dLat/2) * Math.sin(dLat/2) +
    Math.cos(deg2rad(lat1)) * Math.cos(deg2rad(lat2)) * 
    Math.sin(dLon/2) * Math.sin(dLon/2)
    ; 
  var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a)); 
  var d = R * c; // Distance in km
  return d;
}

function deg2rad(deg) {
  return deg * (Math.PI/180)
}

calculate distance between two geocodes in a dataframe

If you refer to this Python's implementation of haversine distance:

df["distance"] = df[["geo1", "geo2"]].apply(lambda x: haversine(*x.geo1, *x.geo2), axis="columns")

>>> df
  Name                      geo1                      geo2    distance
0  ABC  (52.2296756, 21.0122287)  (51.3490756, 23.0922287)  248.451222
1  XYZ  (52.3490756, 23.0922287)  (51.2296756, 21.0122287)  258.456800

Get New GPS with Known start GPS, bearing and distance in feet

Perhaps you could do something like this:

import geopy
import geopy.distance

lat = 30.456025341663068
lon = -86.41408883615411
distance_ft = 86
bearing = 0

start_point = geopy.Point(lat, lon)
end_point = geopy.distance.geodesic(feet=distance_ft).destination(start_point, bearing)

print(end_point.latitude, end_point.longitude)

This should output something like:

30.456261790886277 -86.41408883615411

You can then also use the geodesic method to calculate the distance between the points:

print(geopy.distance.geodesic(start_point, end_point).feet)

And get something like:

86.0000000020017

How can I quickly estimate the distance between two (latitude, longitude) points?

The answers to Haversine Formula in Python (Bearing and Distance between two GPS points) provide Python implementations that answer your question.

Using the implementation below I performed 100,000 iterations in less than 1 second on an older laptop. I think for your purposes this should be sufficient. However, you should profile anything before you optimize for performance.

from math import radians, cos, sin, asin, sqrt
def haversine(lon1, lat1, lon2, lat2):
    """
    Calculate the great circle distance between two points 
    on the earth (specified in decimal degrees)
    """
    # convert decimal degrees to radians 
    lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])
    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
    c = 2 * asin(sqrt(a)) 
    # Radius of earth in kilometers is 6371
    km = 6371* c
    return km

To underestimate haversine(lat1, long1, lat2, long2) * 0.90 or whatever factor you want. I don't see how introducing error to your underestimation is useful.

Efficiently compute pairwise haversine distances between two datasets - NumPy / Python

Pairwise haversine distances

Here's a vectorized way with broadcasting based on this post -

def convert_to_arrays(df1, df2):
    d1 = np.array(df1['coordinates'].tolist())
    d2 = np.array(df2['coordinates'].tolist())
    return d1,d2

def broadcasting_based_lng_lat(data1, data2):
    # data1, data2 are the data arrays with 2 cols and they hold
    # lat., lng. values in those cols respectively
    data1 = np.deg2rad(data1)                     
    data2 = np.deg2rad(data2)                     

    lat1 = data1[:,0]                     
    lng1 = data1[:,1]         

    lat2 = data2[:,0]                     
    lng2 = data2[:,1]         

    diff_lat = lat1[:,None] - lat2
    diff_lng = lng1[:,None] - lng2
    d = np.sin(diff_lat/2)**2 + np.cos(lat1[:,None])*np.cos(lat2) * np.sin(diff_lng/2)**2
    return 2 * 6371 * np.arcsin(np.sqrt(d))

Hence, to solve your case to get all pairwise haversine distances, it would be -

broadcasting_based_lng_lat(*convert_to_arrays(df1,df2))

---

Elementwise haversine distances

For element-wise haversine distance computations between two data, such that each data holds latitude and longitude in two columns each or lists of two elements each, we would skip some of the extensions to 2D and end up with something like this -

def broadcasting_based_lng_lat_elementwise(data1, data2):
    # data1, data2 are the data arrays with 2 cols and they hold
    # lat., lng. values in those cols respectively
    data1 = np.deg2rad(data1)                     
    data2 = np.deg2rad(data2)                     

    lat1 = data1[:,0]                     
    lng1 = data1[:,1]         

    lat2 = data2[:,0]                     
    lng2 = data2[:,1]         

    diff_lat = lat1 - lat2
    diff_lng = lng1 - lng2
    d = np.sin(diff_lat/2)**2 + np.cos(lat1)*np.cos(lat2) * np.sin(diff_lng/2)**2
    return 2 * 6371 * np.arcsin(np.sqrt(d))

Sample run with a dataframe holding the two data in two columns -

In [42]: np.random.seed(0)
    ...: a = np.random.randint(10,100,(5,2)).tolist()
    ...: b = np.random.randint(10,100,(5,2)).tolist()
    ...: df = pd.DataFrame({'A':a,'B':b})

In [43]: df
Out[43]: 
          A         B
0  [54, 57]  [80, 98]
1  [74, 77]  [98, 22]
2  [77, 19]  [68, 75]
3  [93, 31]  [49, 97]
4  [46, 97]  [56, 98]

In [44]: from haversine import haversine

In [45]: [haversine(i,j) for (i,j) in zip(df.A,df.B)]
Out[45]: 
[3235.9659882513424,
 2399.6124657290075,
 2012.0851666001824,
 4702.8069773315865,
 1114.1193334220534]

In [46]: broadcasting_based_lng_lat_elementwise(np.vstack(df.A), np.vstack(df.B))
Out[46]: 
array([3235.96151855, 2399.60915125, 2012.08238739, 4702.80048155,
       1114.11779454])

Those slight differences are largely because haversine library assumes 6371.0088 as the earth radius, while we are taking it as 6371 here.

GPS - Is the haversine formula accurate for distance between two nearby gps points?

The "accurate" distance is subjective here.

Let me explain...

Do you mean the distance over the road? Do you go by bike, or by car? Or by plane?
And do you mean the mathematically shortest distance as a straight line? The shortest possible distance over the surface of the earth?

And do you realize that the smaller the accuracy becomes, the longer the distance will be?
Read more about that here: https://en.wikipedia.org/wiki/Coastline_paradox
It's an interesting read, and if you get it, it'll make you look at your own question in a different way.

Now, if you forget about the real 3d representation, and just assume that the earth is a bit of an elipse-like sphere, things become a lot easier.

In that case, using the haversine formula is probably best for distances of about a mile, especially if you up multiple shorter distances to get a longer one.

If you're talking about GPS data that's measured every second, and you only care about the distance between two points, it's good enough to assume that the earth is flat. But the same goes here: it's fine as long as you don't do it over long distances, or if you add up multiple short distances.

---

Related Topics