Calculate the great-circle distance between geographic (LL) coordinates. Also calculate the initial bearing of the great-circle arc (at its starting point).
calcGCdist(lon1, lat1, lon2, lat2, R=6371.2)
A list obect containing:
a
-- Haversine \(a\) = square of half the chord length between the points,
c
-- Haversine \(c\) = angular distance in radians,
d
-- Haversine \(d\) = great-circle distance (km) between two points,
d2
-- Law of Cosines \(d\) = great-circle distance (km) between two points,
theta
-- Initial bearing \(\theta\) (degrees) for the start point.
numeric
-- Longitude coordinate (degrees) of the start point.
numeric
-- Latitude coordinate(degrees) of the start point.
numeric
-- Longitude coordinate(degrees) of the end point.
numeric
-- Latitude coordinate(degrees) of the end point.
numeric
-- Mean radius (km) of the Earth.
Rowan Haigh, Program Head -- Offshore Rockfish
Pacific Biological Station (PBS), Fisheries & Oceans Canada (DFO), Nanaimo BC
locus opus: Offsite, Vancouver BC
Last modified Rd: 2023-11-03
The great-circle distance is calculated between two points along a
spherical surface using the shortest distance and disregarding
topography.
Method 1: Haversine Formula
$$a = \sin^2((\phi_2 - \phi_1)/2) + \cos(\phi_1) \cos(\phi_2) \sin^2((\lambda_2 - \lambda_1)/2)$$
$$c = 2~\mathrm{atan2}(\sqrt{a}, \sqrt{1-a})$$
$$d = R c$$
where
\(\phi\) = latitude (in radians),
\(\lambda\) = longitude (in radians),
\(R\) = radius (km) of the Earth,
\(a\) = square of half the chord length between the points,
\(c\) = angular distance in radians,
\(d\) = great-circle distance (km) between two points.
Method 2: Spherical Law of Cosines
$$d = \mathrm{acos}(\sin(\phi_1)\sin(\phi_2) + \cos(\phi_1)\cos(\phi_2)\cos(\lambda_2 - \lambda_1)) R$$
The initial bearing (aka forward azimuth) for the start point can be calculated using:
$$\theta = \mathrm{atan2}(\sin(\lambda_2-\lambda_1)\cos(\phi_2), \cos(\phi_1)\sin(\phi_2) - \sin(\phi_1)\cos(\phi_2)\cos(\lambda_2-\lambda_1))$$
Movable Type Scripts -- Calculate distance, bearing and more between Latitude/Longitude points
In package PBSmapping:
addCompass
,
calcArea
,
calcCentroid
,
calcLength
local(envir=.PBSmapEnv,expr={
#-- Distance between southern BC waters and north geomagnetic pole
print(calcGCdist(-126.5,48.6,-72.7,80.4))
})
Run the code above in your browser using DataLab