6

I have 2 pairs of GPS coordinates A = (44.27364400,-121.17116400) and B = (44.27357900, -121.17006800) and a distance of 10m. I'm trying to find the GPS coordinates of the point that's 10m from A toward B. I know I can use haversine to find the distance between A and B coutesy of: https://stackoverflow.com/a/4913653/5369777 is there someway I can add or subtract the distance 10m to return the degrees for the new point?

2
  • geo.sphere.haversine_distance pypi.org/project/geo-py
    – Mapperz
    Feb 3, 2020 at 20:38
  • I'm not seeing exactly how this applies. I can get the midpoint with sphere.haversine_distance/2 and adjust the distance, but as my question stated, how to convert back to gps coordinates?
    – Derek_P
    Feb 4, 2020 at 0:15

2 Answers 2

5

Charles Karney wrote the powerfull GeographicLib program to solve all of this problems.

You can install the geographiclib Python package, it comes as a dependency with geopandas too.

from geographiclib.geodesic import Geodesic

A = (44.27364400, -121.17116400) #Point A (lat, lon)
B = (44.27357900, -121.17006800) #Point B (lat, lon)
s = 10 #Distance (m)

#Define the ellipsoid
geod = Geodesic.WGS84

#Solve the Inverse problem
inv = geod.Inverse(A[0],A[1],B[0],B[1])
azi1 = inv['azi1']
print('Initial Azimuth from A to B = ' + str(azi1))

#Solve the Direct problem
dir = geod.Direct(A[0],A[1],azi1,s)
C = (dir['lat2'],dir['lon2'])
print('C = ' + str(C))

Returns:

Initial Azimuth from A to B = 94.71831670212772
C = (44.27363659722271, -121.17103916871612)
0

I found that I needed compass bearing to get the correct direction for my distance. For bearing I found compassbearing.py @ https://gist.github.com/jeromer/2005586

The inputs of which are pointA, pointB). Output is bearing. From there I can use geopy Point and geopy.distance vincenty to figure out the distance and subsequent new coordinate.

import math from geopy import Point from geopy.distance import vincenty

def calculate_initial_compass_bearing(pointA, pointB):
    """
    Calculates the bearing between two points.
    The formulae used is the following:
    θ = atan2(sin(Δlong).cos(lat2),
              cos(lat1).sin(lat2) − 
    sin(lat1).cos(lat2).cos(Δlong))
    :Parameters:
    - `pointA: The tuple representing the 
    latitude/longitude for the
    first point. Latitude and longitude must be in 
    decimal degrees
    - `pointB: The tuple representing the latitude/longitude for the
    second point. Latitude and longitude must be in decimal degrees
    :Returns:
    The bearing in degrees
    :Returns Type:
    float
    """
    if (type(pointA) != tuple) or (type(pointB) != tuple):
        raise TypeError("Only tuples are supported as arguments")

    lat1 = math.radians(pointA[0])
    lat2 = math.radians(pointB[0])

    diffLong = math.radians(pointB[1] - pointA[1])

    x = math.sin(diffLong) * math.cos(lat2)
    y = math.cos(lat1) * math.sin(lat2) - (math.sin(lat1)
        * math.cos(lat2) * math.cos(diffLong))

    initial_bearing = math.atan2(x, y)

    # Now we have the initial bearing but math.atan2 return values
    # from -180° to + 180° which is not what we want for a compass bearing
    # The solution is to normalize the initial bearing as shown below
    initial_bearing = math.degrees(initial_bearing)
    compass_bearing = (initial_bearing + 360) % 360

    return compass_bearing


bearing = calculate_initial_compass_bearing((44.27364400, -121.17116400),(44.27357900, -121.17006800))

distKm = .01
lat1 = 44.27364400
lon1 = -121.17116400

new_coord = vincenty(kilometers=distKm).destination(Point(lat1, lon1), bearing).format_decimal())

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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