Distance

This script calculates great-circle distances between the two points – that is, the shortest distance over the earth’s surface – using the ‘Haversine’ formula.

It assumes a spherical earth, ignoring ellipsoidal effects – which is accurate enough^* for most purposes… – giving an ‘as-the-crow-flies’ distance between the two points (ignoring any hills!).

Enter the co-ordinates into the text boxes to try it out. It accepts a variety of formats:

• deg-min-sec suffixed with N/S/E/W (e.g. 40°44′55″N, 73 59 11W), or
• signed decimal degrees without compass direction, where negative indicates west/south (e.g. 40.7486, -73.9864):

Lat 1: Long 1:

Lat 2: Long 2:

And you can see it on a map (thanks to the nice guys at Google Maps)

            	var R = 6371; // km var dLat = (lat2-lat1).toRad();
            	var dLon = (lon2-lon1).toRad(); 
            	var a = Math.sin(dLat/2) * Math.sin(dLat/2) + Math.cos(lat1.toRad()) * Math.cos(lat2.toRad()) * Math.sin(dLon/2) * Math.sin(dLon/2); 
            	var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a)); 
            	var d = R * c;

	            	var R = 6371; // km
	            	var d = Math.acos(Math.sin(lat1)*Math.sin(lat2) + Math.cos(lat1)*Math.cos(lat2) *Math.cos(lon2-lon1)) * R;
	            

Bearing

var y = Math.sin(dLon) * Math.cos(lat2);
var x = Math.cos(lat1)*Math.sin(lat2) -
        Math.sin(lat1)*Math.cos(lat2)*Math.cos(dLon);
var brng = Math.atan2(y, x).toBrng();

Since atan2 returns values in the range -π ... +π (that is, -180° ... +180°), to normalise the result to a compass bearing (in the range 0° ... 360°, with -ve values transformed into the range 180° ... 360°), convert to degrees and then use (θ+360) % 360, where % is modulo.

This is the initial bearing which if followed in a straight line along a great-circle arc (orthodrome) will take you from the start point to the end point; in general, the bearing you are following will have varied by the time you get to the end point (if you were to go from say 35°N,45°E (Baghdad) to 35°N,135°E (Osaka), you would start on a bearing of 60° and end up on a bearing of 120°!).

For final bearing, simply take the initial bearing from the end point to the start point and reverse it (using θ = (θ+180) % 360).

Midpoint

var Bx = Math.cos(lat2) * Math.cos(dLon);
var By = Math.cos(lat2) * Math.sin(dLon);
var lat3 = Math.atan2(Math.sin(lat1)+Math.sin(lat2),
                      Math.sqrt((Math.cos(lat1)+Bx)*(Math.cos(lat1)+Bx) + 
                             By*By ) ); 
var lon3 = lon1.toRad() + Math.atan2(By, Math.cos(lat1) + Bx);

Just as the initial bearing may vary from the final bearing, the midpoint may not be located half-way between latitudes/longitudes; the midpoint between 35°N,45°E and 35°N,135°E is around 45°N,90°E.

Destination point given distance and bearing from start point

This page is steadily growing! Given a start point, initial bearing, and distance, this will calculate the destination point and final bearing travelling along a (shortest distance) great circle arc.

var lat2 = Math.asin( Math.sin(lat1)*Math.cos(d/R) + 
                      Math.cos(lat1)*Math.sin(d/R)*Math.cos(brng) );
var lon2 = lon1 + Math.atan2(Math.sin(brng)*Math.sin(d/R)*Math.cos(lat1), 
                             Math.cos(d/R)-Math.sin(lat1)*Math.sin(lat2));

For final bearing, simply take the initial bearing from the end point to the start point and reverse it (using θ = (θ+180) % 360).

Cross-track distance

Here’s a new one: I’ve sometimes been asked about distance of a point from a great-circle path (sometimes called cross track error).

Formula:	dxt = asin(sin(d13/R)*sin(θ13−θ12)) * R
where	d13 is distance from start point to third point θ13 is (initial) bearing from start point to third point θ12 is (initial) bearing from start point to end point R is the earth’s radius
JavaScript:	var dXt = Math.asin(Math.sin(d13/R)*Math.sin(brng13-brng12)) * R;

Here, the great-circle path is identified by a start point and an end point – depending on what initial data you’re working from, you can use the formulae above to obtain the relevant distance and bearings. The sign of d_xt tells you which side of the path the third point is on.

The along-track distance, from the start point to the closest point on the path to the third point, is

var dAt = Math.acos(Math.cos(d13/R)/Math.cos(dXt/R)) * R;

Closest point to the poles

Also new: ‘Clairaut’s formula’ will give you the maximum latitude of a great circle path, given a bearing and latitude on the great circle:

Formula:	latmax = acos(abs(sin(θ)*cos(lat)))
JavaScript:	var latMax = Math.acos(Math.abs(Math.sin(brng)*Math.cos(lat)));

Rhumb lines

A ‘rhumb line’ (or loxodrome) is a path of constant bearing, which crosses all meridians at the same angle.

Sailors used to (and sometimes still) navigate along rhumb lines since it is easier to follow a constant compass bearing than to continually adjust the bearing as is needed to follow a great circle, though they are normally longer than great-circle (orthodrome) routes. Rhumb lines are straight lines on a Mercator Projection map.

var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
var q = (Math.abs(dLat) > 1e-10) ? dLat/dPhi : Math.cos(lat1);
// if dLon over 180° take shorter rhumb across 180° meridian:
if (Math.abs(dLon) > Math.PI) {
  dLon = dLon>0 ? -(2*Math.PI-dLon) : (2*Math.PI+dLon);
}
var d = Math.sqrt(dLat*dLat + q*q*dLon*dLon) * R;
var brng = Math.atan2(dLon, dPhi);

Given a start point and a distance d along constant bearing θ, this will calculate the destination point. If you maintain a constant bearing along a rhumb line, you will gradually spiral in towards one of the poles.

lat2 = lat1 + d*Math.cos(brng);
var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
var q = (Math.abs(lat2-lat1) > 1e-10) ? (lat2-lat1)/dPhi : Math.cos(lat1);
var dLon = d*Math.sin(brng)/q;
// check for some daft bugger going past the pole, normalise latitude if so
if (Math.abs(lat2) > Math.PI/2) lat2 = lat2>0 ? Math.PI-lat2 : -(Math.PI-lat2);
lon2 = (lon1+dLon+Math.PI)%(2*Math.PI) - Math.PI;

If you use Ordnance Survey Grid References, I have implemented a script for converting between Lat/Long & OS Grid References.

Convert between degrees-minutes-seconds & decimal degrees

No, I’ve not included decimal minutes: a decimal system is easy, a sexagesimal system has merits, but mixing the two is a complete sow’s ear. Switch off the option on your GPS!