RTK Torch: Local base station corrections vs Point Perfect

Hi Jim,

I’m from the UK, so your “0.27 inches or 8 cm” jumped out at me. I knocked up some python just to check the numbers. It’s 0.27 feet or 8cm…

Here’s the code I came up with. The standard great-circle distance formula works great for large distances. For small distances, that rounds to zero and haversines are needed.

import math
import argparse
from haversine import haversine, Unit

def greatCircleDistance(lat1, lon1, lat2, lon2):
    # https://en.wikipedia.org/wiki/Great-circle_distance
    lat1_radians = math.radians(lat1)
    lon1_radians = math.radians(lon1)
    lat2_radians = math.radians(lat2)
    lon2_radians = math.radians(lon2)
    delta_lon = lon1_radians - lon2_radians
    cos_delta_lon = math.cos(delta_lon)
    sin_lat1 = math.sin(lat1_radians)
    cos_lat1 = math.cos(lat1_radians)
    sin_lat2 = math.sin(lat2_radians)
    cos_lat2 = math.cos(lat2_radians)
    central_angle = math.acos((sin_lat1 * sin_lat2) + (cos_lat1 * cos_lat2 * cos_delta_lon))
    earth_radius = 6371008.8
    distance = earth_radius * central_angle
    return distance

def greatCircleDistanceHaversine(lat1, lon1, lat2, lon2):
    return haversine((lat1, lon1), (lat2, lon2), unit=Unit.METERS)

if __name__ == '__main__':
    
    parser = argparse.ArgumentParser(
        description='Great Circle Distance Calculator')

    parser.add_argument('-lat1', dest='lat1', default=30.52999699417858, type=float,
                        help='Latitude 1 in Degrees')

    parser.add_argument('-lon1', dest='lon1', default=-97.70260421030356, type=float,
                        help='Longitude 1 in Degrees')

    parser.add_argument('-lat2', dest='lat2', default=30.529997691258917, type=float,
                        help='Latitude 2 in Degrees')

    parser.add_argument('-lon2', dest='lon2', default=-97.70260388406844, type=float,
                        help='Longitude 2 in Degrees')

    args = parser.parse_args()

    print("{:.4f}m".format(greatCircleDistance(args.lat1,args.lon1,args.lat2,args.lon2)))
    print("{:.4f}m".format(greatCircleDistanceHaversine(args.lat1,args.lon1,args.lat2,args.lon2)))

The output is:

0.0000m

0.0836m

Now, this will (probably) spark another discussion about datums and geoids. Looking at the haversine source, they’re using _AVG_EARTH_RADIUS_KM = 6371.0088. Strictly, that radius should be adjusted to match your preferred datum at your position… :smiley:

https://github.com/mapado/haversine/blo … versine.py

Best wishes,

Paul