CNDIST - Calculate the linear distance in nautical miles between two points
Description: Plane, spherical and geodetic distances are calculated.
Inputs:
| Lat1 | The latitude of point 1 |
| Lon1 | The longitude of point 1 |
| Lat2 | The latitude of point 2 |
| Lon2 | The longitude of point 2 |
Outputs:
| Plane distance | The plane distance in nautical miles (Pythagoras) |
| Spherical distance | The spherical distance in nautical miles |
| Geodetic distance | The geodetic distance in nautical miles (Andoyer) |
| Azimuth | The azimuth of point 2 from point 1 |
Sample execution:
Lat1? 25,30 Lon1? -76,25 Lat2? 26,41 Lon2? -77,15 Plane Distance: 84.1 NM Spherical Distance: 84.0 NM Geodetic distance: 83.9 NM Azimuth 327° 52.2'
(For explanation of notation, conventions etc, see python-programs).
Copyright (C) 2024 Ian Staniforth
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see
https://www.gnu.org/licenses/.
# CNDIST - Calculate the distance in nautical miles between 2 points
from math import *
from CN_LIB import *
lat1 = stod("Lat1? ")
lon1 = stod("Lon1? ")
lat2 = stod("Lat2? ")
lon2 = stod("Lon2? ")
# Plane distance
dlat = abs(lat2-lat1)
dlon = abs(lon2-lon1)*Cos(lat1)
d = sqrt(dlat*dlat+dlon*dlon)
dnm = d*60
print("")
print("Plane Distance: "+"%.1f NM"%dnm)
# Spherical distance
d = aCos(Sin(lat1)*Sin(lat2)+Cos(lat1)*Cos(lat2)*Cos(lon1-lon2))
dnm = d*60
print("Spherical Distance: "+"%.1f NM"%dnm)
# Geodetic distance (Andoyer's method - Meeus p.85)
f = 1/298.25722
a = 6378.14
p1 = lat1
p2 = lat2
l1 = lon1
l2 = lon2
F = (p1+p2)/2
G = (p1-p2)/2
l = (l1-l2)/2
S = Sin(G)*Sin(G)*Cos(l)*Cos(l) + Cos(F)*Cos(F)*Sin(l)*Sin(l)
C = Cos(G)*Cos(G)*Cos(l)*Cos(l) + Sin(F)*Sin(F)*Sin(l)*Sin(l)
w = atan(sqrt(S/C))
R = sqrt(S*C)/w
D = 2*w*a
H1 = (3*R-1)/(2*C)
H2 = (3*R+1)/(2*S)
s = D*(1+f*H1*Sin(F)*Sin(F)*Cos(G)*Cos(G)-f*H2*Cos(F)*Cos(F)*Sin(G)*Sin(G))
snm = s/1.852
print("Geodetic distance: "+"%.1f NM"%snm)
# Azimuth
z = aCos((Sin(lat2)-Sin(lat1)*Cos(d))/(Cos(lat1)*Sin(d)))
zn = 360-z if lon2 < lon1 else z
print("Azimuth "+dtos(zn))