Search code examples
pythonlatitude-longitude

Calculate distance between a point and a line segment in latitude and longitude


I have a line segments defined with a start and an end point:

A: 
x1 = 10.7196405787775
y1 = 59.9050401935882

B: 
x2 = 10.7109989561813
y2 = 59.9018650448204

where x defines longitude and y defines latitude.

I also have a point:

P:
x0 = 10.6542116666667
y0 = 59.429105

How do I compute the shortest distance between the line segment and the point? I know how to do this in Cartesian coordinates, but not in long/lat coordinates.


Solution

  • Using the helpful Python geocoding library geopy, and the formula for the midpoint of a great circle from Chris Veness's geodesy formulae, we can find the distance between a great circle arc and a given point:

    from math import sin, cos, atan2, sqrt, degrees, radians, pi
    from geopy.distance import great_circle as distance
    from geopy.point import Point
    
    
    def midpoint(a, b):
        a_lat, a_lon = radians(a.latitude), radians(a.longitude)
        b_lat, b_lon = radians(b.latitude), radians(b.longitude)
        delta_lon = b_lon - a_lon
        B_x = cos(b_lat) * cos(delta_lon)
        B_y = cos(b_lat) * sin(delta_lon)
        mid_lat = atan2(
            sin(a_lat) + sin(b_lat),
            sqrt(((cos(a_lat) + B_x)**2 + B_y**2))
        )
        mid_lon = a_lon + atan2(B_y, cos(a_lat) + B_x)
        # Normalise
        mid_lon = (mid_lon + 3*pi) % (2*pi) - pi
        return Point(latitude=degrees(mid_lat), longitude=degrees(mid_lon))
    

    Which in this example gives:

    # Example:
    a = Point(latitude=59.9050401935882, longitude=10.7196405787775)
    b = Point(latitude=59.9018650448204, longitude=10.7109989561813)
    p = Point(latitude=59.429105, longitude=10.6542116666667)
    
    d = distance(midpoint(a, b), p)
    print d.km
    # 52.8714586903