Functions to calculate APPROXIMATE range, bearing, and XY offset between two lat/long positions

latLongCalcs.py

# Functions to calculate APPROXIMATE range, bearing, and XY offset between two lat/long positions
# Uses the Haversine Formula 
# See http://www.movable-type.co.uk/scripts/latlong.html for more details

# Math library needed for cos, sine, atan2, sqrt, degrees, and radians
import math 

# This is the mean radius of the earth in metres
EARTH_RADIUS = 6371000.0

# Calculate the distance in metres between lat1,long1 and lat2,long2
# Lat/long must be in decimal degrees format
def range_m(lat1, long1, lat2, long2):

    lat1_r = math.radians(lat1)
    lat2_r = math.radians(lat2)
    
    delta_lat = math.radians(lat2-lat1)
    delta_long = math.radians(long2-long1)    
    
    a = math.sin(delta_lat/2.0) ** 2 + math.cos(lat1_r) * math.cos(lat2_r) * math.sin(delta_long/2.0) ** 2
    b = 2.0 * math.atan2(math.sqrt(a), math.sqrt(1.0-a))
    range_m = EARTH_RADIUS * b
    
    return range_m

# Calculate the bearing (from north) between lat1,long1 and lat2,long2
# Lat/long must be in decimal degrees format
# Bearing output in radians between -pi and +pi
def bearing_rad(lat1, long1, lat2, long2):
    
    lat1_r = math.radians(lat1)
    lat2_r = math.radians(lat2)

    delta_long = math.radians(long2-long1)

    y = math.sin(delta_long) * math.cos(lat2_r)
    x = math.cos(lat1_r)*math.sin(lat2_r) - math.sin(lat1_r)*math.cos(lat2_r)*math.cos(delta_long)
    return math.atan2(y,x)    

# Calculate the bearing (from north) between lat1,long1 and lat2,long2
# Lat/long must be in decimal degrees format
# Bearing output in degrees (0 to 360, with 0 being north)
def bearing_deg(lat1, long1, lat2, long2):
    
    bearing_r = bearing_rad(lat1, long1, lat2, long2)    
    
    return (math.degrees(bearing_r) + 360.0) % 360.0


# Calculate the XY offset in metres between latRef,longRef and lat,lng
# Lat/long must be in decimal degrees format
# For x +ive is North, -ive is South
# For y +ive is East, -ive is West
def offset(latRef, longRef, lat, lng):

    distance = range_m(latRef, longRef, lat, lng)
    bearing = bearing_rad(latRef, longRef, lat, lng)    

    x = math.cos(bearing) * distance
    y = math.sin(bearing) * distance
    
    return x,y