Haversine vs Vincenty: Measuring Distance Between Coordinates
How the haversine formula works on a sphere, why the WGS-84 ellipsoid changes the answer by up to 0.5%, how Vincenty's iteration converges or fails, bearings, midpoints and DMS.
· 10 min read
Two GPS coordinates define a distance, but which distance depends on the model of the Earth you use. Treat the Earth as a sphere and the haversine formula gives the great-circle distance in one line of code. Treat it as the flattened ellipsoid that GPS actually uses and Vincenty's formulae give the geodesic distance to a fraction of a millimetre. This post explains both, shows how far apart they can be, and covers bearings, midpoints and coordinate formats. The examples run in the Coordinate Distance Calculator.
The haversine formula on a sphere
On a sphere of radius R, the shortest path between two points is an arc of a great circle, and its length is R times the central angle between them. With latitudes φ and longitudes λ in radians:
a = sin²(Δφ / 2) + cos φ1 x cos φ2 x sin²(Δλ / 2)
d = 2R x asin(√a)
The name comes from the haversine function, hav(θ) = sin²(θ / 2), which navigators tabulated in the 19th century. The spherical law of cosines gives the same answer in theory, but for short distances it takes the arccosine of a number very close to 1, where floating-point rounding destroys precision. The haversine form stays well-conditioned for small angles, which is why it became the standard.
The radius is a choice. 6,371 km is the usual value; the International Union of Geodesy and Geophysics defines the mean radius of the WGS-84 ellipsoid as (2a + b) / 3 = 6,371,008.8 m, which is what the calculator uses. With that radius, one degree of longitude along the equator is 111.195 km.
Why the sphere is not quite right
The Earth bulges at the equator. WGS-84, the reference ellipsoid for GPS, has an equatorial radius a = 6,378,137 m and a flattening f = 1/298.257223563, which makes the polar radius b = 6,356,752.3 m, about 21 km shorter.
On the ellipsoid, the curvature along a meridian is tighter near the equator and flatter near the poles, so a degree of latitude is about 110.6 km at the equator and 111.7 km at the poles. A sphere uses one value everywhere. The result is that haversine distances can be off by up to about 0.5%, too long on some routes and too short on others. London to New York comes out 15 km short on the sphere (5,570.23 km against 5,585.23 km, -0.27%), while New Delhi to Mumbai, a north-south route nearer the equator, comes out 3.6 km long (1,148.10 against 1,144.53 km). For routing a delivery van that does not matter. For surveying, aviation fuel planning over long sectors, or checking a GIS system, it does.
Vincenty's inverse formula
Thaddeus Vincenty published "Direct and Inverse Solutions of Geodesics on the Ellipsoid with Application of Nested Equations" in Survey Review in 1975. The inverse problem, finding the distance and azimuths between two known points, is solved by mapping both points to an auxiliary sphere using reduced latitudes, tan U = (1 - f) tan φ, and iterating:
- Start with λ equal to the longitude difference L.
- Compute the angular distance σ on the auxiliary sphere and the azimuth α of the geodesic at the equator.
- Update λ = L + (1 - C) f sin α (σ + C sin σ (...)), a correction for the flattening.
- Repeat until λ changes by less than about 10⁻¹² radians, which is well under a millimetre.
The distance then comes from a series in u² = cos²α (a² - b²) / b²: s = b A (σ - Δσ). The iteration usually converges in three to six steps. Geoscience Australia publishes a standard test line from Flinders Peak (37°57′03.72030″ S, 144°25′29.52440″ E) to Buninyong (37°39′10.15610″ S, 143°55′35.38390″ E): 54,972.271 m, forward azimuth 306°52′05.37″ and reverse azimuth 127°10′25.07″. The calculator reproduces those values to the millimetre and the hundredth of a second, which is a good check that an implementation is correct.
Vincenty's method has one weakness. For nearly antipodal points, where the geodesic could pass almost anywhere around the globe, λ can oscillate and never settle. The calculator detects that, says so and falls back to haversine. Charles Karney's 2013 algorithms, published in the Journal of Geodesy and implemented in GeographicLib, converge for every pair of points and are now the reference in survey and GIS software.
Bearings and midpoints
The initial bearing is the direction you set off in, measured clockwise from true north:
θ = atan2(sin Δλ x cos φ2, cos φ1 x sin φ2 - sin φ1 x cos φ2 x cos Δλ)
Along a great circle the bearing changes continuously, except on meridians and the equator. Leaving London for New York, you head 288° (west-north-west) and arrive heading 231° (south-west), which is why great-circle routes look curved on a flat map. The final bearing equals the reverse of the initial bearing from B to A, plus 180°. On the ellipsoid, Vincenty returns both azimuths directly.
The great-circle midpoint is not the average of the coordinates. London to New York has a midpoint near 52.4° N, 41.3° W, north of both cities, because the great circle bows toward the pole. Averaging the latitudes would give 46° N, hundreds of kilometres off.
Bearings here are true bearings. A compass reads magnetic bearings, which differ by the local magnetic declination, currently a few degrees to well over ten degrees depending on location and year.
Coordinate formats and datums
Coordinates come as decimal degrees (40.4461, -79.9822) or degrees, minutes and seconds (40° 26′ 46″ N, 79° 58′ 56″ W). Converting is simple: decimal = degrees + minutes / 60 + seconds / 3600, negative for south and west. The calculator accepts either, with symbols, spaces or colons. Precision is worth thinking about: one degree of latitude is about 111 km, so the fifth decimal place is about 1.1 m and the sixth about 0.11 m. Typing more digits than your GPS can deliver adds nothing.
Coordinates also carry a datum. GPS reports WGS-84, and modern national datums such as NAD83 and GDA2020 agree with it to within a metre or two. Older datums such as NAD27, ED50 and many local grids can differ by tens to hundreds of metres, so convert them before computing distances.
When to use which
Use haversine for most software features: nearest-store search, approximate trip lengths, clustering and anything where 0.5% is acceptable. Use Vincenty or Karney for survey checks, aviation and marine planning, long baselines, and any time you are validating another system. Neither includes height: the distances are on the ellipsoid surface, and terrain or altitude make the real path longer.
For travel times from a distance, use the Speed, Distance and Time Calculator, and the Unit Converter handles any other length units. On construction sites, spot levels taken with GPS feed the grid method in Cut and Fill Earthwork Volumes, and the Cut and Fill Calculator turns them into cut and fill.