Before going further, I have to mention that, at a fine level of accuracy, there are different ways of measuring latitude & longitude. The Ordnance Survey uses ‘OSGB-36’, based on an elliptical model of the earth’s surface which is a good fit to the UK. GPS systems generally use the world-wide ‘WGS-84’, based on an elliptical model which is a best approximation to the entire earth. At Greenwich, these differ by about 126m (they coincide somewhere in the Atlantic ocean; there’s more on Wikipedia).
A source of possible confusion to be aware of is that a grid reference identifies a square (with the size of the square determined by the number of digits), whereas a latitude/longitude coordinate identifies a point (with precision implied by the number of digits).
Note also that an alternative way of expressing OS Grid References is as all-numeric eastings and northings. As square TG is six squares across, three squares up within the grid, grid reference TG 5140 1317 can also be expressed as 65140,31317.
Beware that eastings and northings can be used either as a grid reference (identifying a grid square), or as a coordinate (identifying a point). In this example, the coordinate representing the centre of grid reference TG 5140 1317 would be 651405,313175. It is often ambiguous whether eastings and northings are being used as grid references or as coordinates.
The Ordnance Survey grid is a Transverse Mercator projection (with origin at 49°N, 2°W) based on the Airy 1830 ellipsoid using the OSGB36 datum. GPS is based on WGS84/GRS80, which as mentioned can vary from OSGB36 by as much as 120m or 6" or arc (OSGB36/Airy is a better fit to the UK geiod than the geocentric WGS84 which covers the entire world). I have written some separate notes on converting between OSGB-36 & WGS-84.
The JavaScript implementation should be quite simple to translate to other languages, if required. Since JavaScript lacks a power operator, I opted to keep the script easier to read by using temporary variables and multiplication rather than the Math.pow method. Since JavaScript is untyped, expressions such as ‘5/4’ give floating-point results – typed languages may need explicit casts to give floating-point results. Also to keep the scripts simple, I have included minimal error checking, and no user options for specifying precision.
As this script is for converting OS grid references, I have ‘hard-wired’ in the Airy 1830 axes and the National Grid projection origin & scale factor; if you want to convert to other transverse mercator projections, you will need to change these constants. UTM (Universal Transverse Mercator projection) uses a scale factor of 0.9996 and origins at 6° intervals of longitude, based on the WGS84 ellipse (though previously UTM projections were based on International 1924 and Clark 1866 ellipses, among others).
Aside from the transformation maths, the other tricky bit of the script is converting grid letter-pairs to/from numeric eastings & northings. To follow what’s going on, it is worth noting that the letter-pairs define a 5x5 grid of 5x5 sub-grids; the eastings & northings work from a ‘false origin’ at grid square SV, which is displaced from grid square AA by 10 squares E, 19 squares N, with the northing axis inverted; and letter ‘I’ is skipped. OS Grid References apply to the UK only.
Code: Select all
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Convert latitude/longitude <=> OS National Grid Reference points (c) Chris Veness 2005-2009 */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/*
* convert geodesic co-ordinates to OS grid reference
*/
function LatLongToOSGrid(p) {
var lat = p.lat.toRad(), lon = p.lon.toRad();
var a = 6377563.396, b = 6356256.910; // Airy 1830 major & minor semi-axes
var F0 = 0.9996012717; // NatGrid scale factor on central meridian
var lat0 = (49).toRad(), lon0 = (-2).toRad(); // NatGrid true origin
var N0 = -100000, E0 = 400000; // northing & easting of true origin, metres
var e2 = 1 - (b*b)/(a*a); // eccentricity squared
var n = (a-b)/(a+b), n2 = n*n, n3 = n*n*n;
var cosLat = Math.cos(lat), sinLat = Math.sin(lat);
var nu = a*F0/Math.sqrt(1-e2*sinLat*sinLat); // transverse radius of curvature
var rho = a*F0*(1-e2)/Math.pow(1-e2*sinLat*sinLat, 1.5); // meridional radius of curvature
var eta2 = nu/rho-1;
var Ma = (1 + n + (5/4)*n2 + (5/4)*n3) * (lat-lat0);
var Mb = (3*n + 3*n*n + (21/8)*n3) * Math.sin(lat-lat0) * Math.cos(lat+lat0);
var Mc = ((15/8)*n2 + (15/8)*n3) * Math.sin(2*(lat-lat0)) * Math.cos(2*(lat+lat0));
var Md = (35/24)*n3 * Math.sin(3*(lat-lat0)) * Math.cos(3*(lat+lat0));
var M = b * F0 * (Ma - Mb + Mc - Md); // meridional arc
var cos3lat = cosLat*cosLat*cosLat;
var cos5lat = cos3lat*cosLat*cosLat;
var tan2lat = Math.tan(lat)*Math.tan(lat);
var tan4lat = tan2lat*tan2lat;
var I = M + N0;
var II = (nu/2)*sinLat*cosLat;
var III = (nu/24)*sinLat*cos3lat*(5-tan2lat+9*eta2);
var IIIA = (nu/720)*sinLat*cos5lat*(61-58*tan2lat+tan4lat);
var IV = nu*cosLat;
var V = (nu/6)*cos3lat*(nu/rho-tan2lat);
var VI = (nu/120) * cos5lat * (5 - 18*tan2lat + tan4lat + 14*eta2 - 58*tan2lat*eta2);
var dLon = lon-lon0;
var dLon2 = dLon*dLon, dLon3 = dLon2*dLon, dLon4 = dLon3*dLon, dLon5 = dLon4*dLon, dLon6 = dLon5*dLon;
var N = I + II*dLon2 + III*dLon4 + IIIA*dLon6;
var E = E0 + IV*dLon + V*dLon3 + VI*dLon5;
return gridrefNumToLet(E, N, 8);
}
/*
* convert OS grid reference to geodesic co-ordinates
*/
function OSGridToLatLong(gridRef) {
var gr = gridrefLetToNum(gridRef);
var E = gr[0], N = gr[1];
var a = 6377563.396, b = 6356256.910; // Airy 1830 major & minor semi-axes
var F0 = 0.9996012717; // NatGrid scale factor on central meridian
var lat0 = 49*Math.PI/180, lon0 = -2*Math.PI/180; // NatGrid true origin
var N0 = -100000, E0 = 400000; // northing & easting of true origin, metres
var e2 = 1 - (b*b)/(a*a); // eccentricity squared
var n = (a-b)/(a+b), n2 = n*n, n3 = n*n*n;
var lat=lat0, M=0;
do {
lat = (N-N0-M)/(a*F0) + lat;
var Ma = (1 + n + (5/4)*n2 + (5/4)*n3) * (lat-lat0);
var Mb = (3*n + 3*n*n + (21/8)*n3) * Math.sin(lat-lat0) * Math.cos(lat+lat0);
var Mc = ((15/8)*n2 + (15/8)*n3) * Math.sin(2*(lat-lat0)) * Math.cos(2*(lat+lat0));
var Md = (35/24)*n3 * Math.sin(3*(lat-lat0)) * Math.cos(3*(lat+lat0));
M = b * F0 * (Ma - Mb + Mc - Md); // meridional arc
} while (N-N0-M >= 0.00001); // ie until < 0.01mm
var cosLat = Math.cos(lat), sinLat = Math.sin(lat);
var nu = a*F0/Math.sqrt(1-e2*sinLat*sinLat); // transverse radius of curvature
var rho = a*F0*(1-e2)/Math.pow(1-e2*sinLat*sinLat, 1.5); // meridional radius of curvature
var eta2 = nu/rho-1;
var tanLat = Math.tan(lat);
var tan2lat = tanLat*tanLat, tan4lat = tan2lat*tan2lat, tan6lat = tan4lat*tan2lat;
var secLat = 1/cosLat;
var nu3 = nu*nu*nu, nu5 = nu3*nu*nu, nu7 = nu5*nu*nu;
var VII = tanLat/(2*rho*nu);
var VIII = tanLat/(24*rho*nu3)*(5+3*tan2lat+eta2-9*tan2lat*eta2);
var IX = tanLat/(720*rho*nu5)*(61+90*tan2lat+45*tan4lat);
var X = secLat/nu;
var XI = secLat/(6*nu3)*(nu/rho+2*tan2lat);
var XII = secLat/(120*nu5)*(5+28*tan2lat+24*tan4lat);
var XIIA = secLat/(5040*nu7)*(61+662*tan2lat+1320*tan4lat+720*tan6lat);
var dE = (E-E0), dE2 = dE*dE, dE3 = dE2*dE, dE4 = dE2*dE2, dE5 = dE3*dE2, dE6 = dE4*dE2, dE7 = dE5*dE2;
lat = lat - VII*dE2 + VIII*dE4 - IX*dE6;
var lon = lon0 + X*dE - XI*dE3 + XII*dE5 - XIIA*dE7;
return new LatLon(lat.toDeg(), lon.toDeg());
}
/*
* convert standard grid reference ('SU387148') to fully numeric ref ([438700,114800])
* returned co-ordinates are in metres, centred on grid square for conversion to lat/long
*
* note that northern-most grid squares will give 7-digit northings
* no error-checking is done on gridref (bad input will give bad results or NaN)
*/
function gridrefLetToNum(gridref) {
// get numeric values of letter references, mapping A->0, B->1, C->2, etc:
var l1 = gridref.toUpperCase().charCodeAt(0) - 'A'.charCodeAt(0);
var l2 = gridref.toUpperCase().charCodeAt(1) - 'A'.charCodeAt(0);
// shuffle down letters after 'I' since 'I' is not used in grid:
if (l1 > 7) l1--;
if (l2 > 7) l2--;
// convert grid letters into 100km-square indexes from false origin (grid square SV):
var e = ((l1-2)%5)*5 + (l2%5);
var n = (19-Math.floor(l1/5)*5) - Math.floor(l2/5);
// skip grid letters to get numeric part of ref, stripping any spaces:
gridref = gridref.slice(2).replace(/ /g,'');
// append numeric part of references to grid index:
e += gridref.slice(0, gridref.length/2);
n += gridref.slice(gridref.length/2);
// normalise to 1m grid, rounding up to centre of grid square:
switch (gridref.length) {
case 6: e += '50'; n += '50'; break;
case 8: e += '5'; n += '5'; break;
// 10-digit refs are already 1m
}
return [e, n];
}
/*
* convert numeric grid reference (in metres) to standard-form grid ref
*/
function gridrefNumToLet(e, n, digits) {
// get the 100km-grid indices
var e100k = Math.floor(e/100000), n100k = Math.floor(n/100000);
if (e100k<0 || e100k>6 || n100k<0 || n100k>12) return '';
// translate those into numeric equivalents of the grid letters
var l1 = (19-n100k) - (19-n100k)%5 + Math.floor((e100k+10)/5);
var l2 = (19-n100k)*5%25 + e100k%5;
// compensate for skipped 'I' and build grid letter-pairs
if (l1 > 7) l1++;
if (l2 > 7) l2++;
var letPair = String.fromCharCode(l1+'A'.charCodeAt(0), l2+'A'.charCodeAt(0));
// strip 100km-grid indices from easting & northing, and reduce precision
e = Math.floor((e%100000)/Math.pow(10,5-digits/2));
n = Math.floor((n%100000)/Math.pow(10,5-digits/2));
var gridRef = letPair + e.padLZ(digits/2) + n.padLZ(digits/2);
return gridRef;
}
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/*
* extend Number object with methods for converting degrees/radians
*/
Number.prototype.toRad = function() { // convert degrees to radians
return this * Math.PI / 180;
}
Number.prototype.toDeg = function() { // convert radians to degrees (signed)
return this * 180 / Math.PI;
}
/*
* pad a number with sufficient leading zeros to make it w chars wide
*/
Number.prototype.padLZ = function(w) {
var n = this.toString();
for (var i=0; i<w-n.length; i++) n = '0' + n;
return n;
}
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */