NAME

Geo::Coordinates::OSGB --- Convert Coordinates from Lat/Long to UK Grid

A UK-specific implementation of co-ordinate conversion, following formulae from the Ordnance Survey of Great Britain (hence the name).

Version: 1.07

SYNOPSIS

use Geo::Coordinates::OSGB qw(ll2grid grid2ll);

# basic conversion routines
($easting,$northing) = ll2grid($lat,$lon);
($lat,$long) = grid2ll($easting,$northing);

# format full easting and northing into traditional formats
$trad_gr       = format_grid_trad($easting,$northing);  # TQ 234 098
$GPS_gr        = format_grid_GPS($easting,$northing);   # TQ 23451 09893
$landranger_gr = format_grid_landranger($easting, $northing) # 234098 on Sheet 176

# you can call these in list context to get the individual parts
# add "=~ s/\s//g" to the result to remove the spaces

# and there are corresponding parse routines to convert from these formats to full e,n
($e,$n) = parse_trad_grid('TQ 234 098'); # spaces optional, can give list as well
($e,$n) = parse_GPS_grid('TQ 23451 09893'); # spaces optional, can give list as well
($e,$n) = parse_landranger_grid($sheet, $gre, $grn); # gre/grn must be 3 or 5 digits long
($e,$n) = parse_landranger_grid($sheet); # this gives you the SW corner of the sheet

# some convenience routines that bundle these up for you
map2ll();
map2grid();

# set parameters
set_ellipsoid(6377563.396,6356256.91);                 # not needed for UK use
set_projection(49, -2, 400000, -100000, 0.9996012717); # not needed for UK use

DESCRIPTION

This module provides a collection of routines to convert between longitude & latitude and map grid references, using the formulae given in the British Ordnance Survey's excellent information leaflet, referenced below in "Theory".

The module is implemented purely in Perl, and should run on any Perl platform. In this description `OS' means `the Ordnance Survey of Great Britain': the UK government agency that produces the standard maps of England, Wales, and Scotland. Any mention of `sheets' or `maps' refers to one or more of the 204 sheets in the 1:50,000 scale `Landranger' series of OS maps.

FUNCTIONS

The following functions can be exported from the Geo::Coordinates::OSGB module:

ll2grid
grid2ll

format_grid_trad
format_grid_GPS
format_grid_landranger

parse_trad_grid
parse_GPS_grid
parse_landranger_grid

map2ll
map2grid

set_ellipsoid
set_projection

None is exported by default.

This code is fine tuned to the British national grid system. You can use it elsewhere but you will need to adapt it. This is explained in some detail in the "Examples" section below.

The default values for ellipsoid and projection are suitable for mapping between GPS longitude and latitude data and the UK National Grid.

ll2grid(lat,lon)

When called in a void context, or with no arguments ll2grid does nothing. When called in a list context, ll2grid returns two numbers that represent the easting and the northing corresponding to the latitude and longitude supplied.

The parameters can be supplied as real numbers representing degrees or in ISO `degrees, minutes, and seconds' form. That is 'sdd[mm[ss]]' for lat and 'sddd[mm[ss]]' for long. The magnitude of the arguments is used to infer which form is being used. Note the leading s is the sign + (North,East) or - (South,West). If you use the ISO format be sure to 'quote' the arguments, otherwise Perl will think that they are numbers, and strip leading 0s and + signs which may give you unexpected results. For example:

my ($e,$n) = ll2grid('+5120','-00025');

If you have trouble remembering the order of the arguments, note that latitude comes before longitude in the alphabet too.

The easting and northing will be returned as a whole number of metres from the point of origin defined by the projection you have set. In the case of the Britain this is a point a little way to the south-west of the Scilly Isles. If you want the grid presented in a more traditional format you should pass the results to one of the grid formatting routines, which are described below.

If you call ll2grid in a scalar context, it will automatically call format_grid_trad. For example:

my $gridref = ll2grid('+5120','-00025');

In this case the string returned represents the `full national grid reference' with two letters and two sets of three numbers, like this `TQ 102 606'. If you want to remove the spaces, just apply s/\s//g to it.

To force it to call one of the other grid formatting routines, try one of these:

$gridref = ll2grid('+5120','-00025','Trad');
$gridref = ll2grid('+5120','-00025','GPS');
$gridref = ll2grid('+5120','-00025','Landranger');
format_grid_trad(e,n)

Formats an (easting, northing) pair into traditional `full national grid reference' with two letters and two sets of three numbers, like this `TQ 102 606'. If you want to remove the spaces, just apply s/\s//g to it. If you want the individual components call it in a list context.

format_grid_GPS(e,n)

Users who have bought a GPS receiver may initially have been puzzled by the unfamiliar format used to present coordinates in the British national grid format. On my Garmin Legend C it shows this sort of thing in the display.

 TQ 23918
bng 00972

and in the track logs the references look like this

TQ 23918 00972

These are just the same as the references described on the OS sheets, except that the units are metres rather than hectometres, so you get five digits in each of the easting and northings instead of three. format_grid_GPS returns a string representing this format, or a list of the square, the truncated easting, and the truncated northing if you call it in a list context.

Note that, at least until WAAS is working in Europe, the results from your GPS are unlikely to be more accurate than plus or minus 10m even with perfect reception.

format_grid_landranger(e,n)

This does the same as format_grid_trad, but it appends the number of the relevant OS Landranger 1:50,000 scale map to the traditional grid reference. Note that there may be several or no sheets returned. This is because many (most) of the Landranger sheets overlap, and many other valid grid references are not on any of the sheets (because they are in the sea or a remote island. This module does not cope with the detached insets on some sheets (yet).

In a list context you will get back a list like this: (square, easting, northing, sheet) or (square, easting, northing, sheet1, sheet2) etc. There are a few places where three sheets overlap, and one corner of Herefordshire which appears on four maps (sheets 137, 138, 148, and 149). If the GR is not on any sheet, then the list of sheets will be empty.

In a scalar context you will get back the same information in a helpful string form like this "NN 241 738 on OS Sheet 44". Note that the easting and northing will have been truncated to the normal truncated `hectometre' three digit form.

parse_trad_grid(grid_ref)

Turns a traditional grid reference into a full easting and northing pair in metres from the point of origin. The grid_ref can be a string like `TQ203604' or `SW 452 004', or a list like this ('TV', '435904') or a list like this ('NN', '345', '208').

parse_GPS_grid(grid_ref)

Does the same as parse_trad_grid but is looking for five digit numbers like `SW 45202 00421', or a list like this ('NN', '34592', '20804').

parse_landranger_grid($sheet, $e, $n)

This converts an OS Landranger sheet number and a local grid reference into a full easting and northing pair in metres from the point of origin.

The OS Landranger sheet number should be between 1 and 204 inclusive (but I may extend this when I support insets). You can supply (e,n) as 3-digit hectometre numbers or 5-digit metre numbers. In either case if you supply any leading zeros you should 'quote' the numbers to stop Perl thinking that they are octal constants.

This module will croak at you if you give it an undefined sheet number, or if the grid reference that you supply does not exist on the sheet.

In order to get just the coordinates of the SW corner of the sheet, just call it with the sheet number.

grid2ll(e,n) or grid2ll(grid_ref)

When called in list context grid2ll returns a pair of numbers representing longitude and latitude coordinates, as real numbers. Following convention, positive numbers are North and East, negative numbers are South and West. The fractional parts of the results represent fractions of degrees.

When called in scalar context it returns a string in ISO longitude and latitude form, such as '+5025-00403' with the result rounded to the nearest minute (the formulae are not much more accurate than this). In a void context it does nothing.

The arguments can be either an (easting, northing) pair of integers representing the absolute grid reference in metres from the point of origin, or a single string that represents a full grid reference in traditional hectometre form, such as 'NH868943' or 'NH 868 943'.

To force it to read GPS format grid references you can try grid2ll(grid_ref,'GPS'). This should then recognize strings like 'NH 87612 27623'.

map2grid()

Shorthand for format_grid_trad(parse_landranger_grid()).

map2ll()

Shorthand for grid2ll(parse_landranger_grid()).

set_ellipsoid(a,b)

Defines the ellipsoid used to interpret the longitude and latitude values. The arguments a and b are the lengths (in metres) of the semi-major axes of the ellipsoid used to represent the earth's surface. Values used in the UK are given in Annex A of the paper referenced below in "Theory".

You should call set_ellipsoid() before doing anything else, unless you are converting data for the UK National Grid. It will default to the values for the `Airy 1830' ellipsoid that is used with the UK National Grid.

set_projection(lat,long,E,N,F)

Defines the projection used to interpret the grid references. The projection is a `Transverse Mercator projection'. The first two arguments define the longitude and latitude of the true origin of the grid to be used, and the second two are the grid coordinates of this position. Note the order that they are given. Latitude then longitude, followed by easting, then northing. This may seem illogical but it does conform to normal practice.

The fifth argument is the scale factor on the central meridian of the grid area. See the paper referenced below in "Theory" for a table of values suitable for the UK and a full explanation of the theory.

You don't need to call this if you just want to use the normal OS grid.

THEORY

The algorithms and theory for these conversion routines are all from A Guide to Coordinate Systems in Great Britain published by the Ordnance Survey, April 1999 and available at http://www.gps.gov.uk/info.asp

You may also like to read some of the other introductory material there. Should you be hoping to adapt this code to your own custom Mercator projection, you will find the paper called Surveying with the National GPS Network, especially useful.

The true point of origin of the British Grid is the point 49N 2W (near the Channel Islands). If you look at the appropriate OS maps you will notice that the 2W meridian is parallel to all the vertical grid lines. To avoid negative numbers in grid references the (0,0) point on the grid is offset 400 km west and 100 km north of this point. This is called the `false point of origin' and all grid references are measured in metres from this point. The easting is always given before the northing.

For everyday use, the OS suggest that grid references are given in units of 100m (hectometres), and that the country should divided into a series of 100km squares. Within each of these large squares, we need only be concerned with the last three digits of the full national grid reference. If we combine the easting and northing we get the familiar traditional six figure grid reference. Each of these grid references is repeated in each of the large 100km squares but for local use, this does not usually matter. Where it does matter, the OS suggest that the six figure reference is prefixed with the identifier of the large grid square to give a `full national grid reference'. This system is described in the notes of in the corner of every Landranger 1:50,000 scale map.

Modern GPS receivers can all display coordinates in the OS grid system. You just need to set the display units to be `British National Grid' or whatever similar name is used on your unit. Most units display the coordinates as two groups of five digits and a grid square identifier. The units are metres within the grid square (although beware that the GPS fix is unlikely to be accurate down to the last metre).

Each of the large squares is identified in pair of letters: TQ, SU, ND, etc. The grid of the big squares actually used is something like this:

            HP
            HU
         HY
NA NB NC ND
NF NG NH NJ NK
NL NM NN NO NP
   NR NS NT NU
   NW NX NY NZ
      SC SD SE TA
      SH SJ SK TF TG
   SM SN SO SP TL TM
   SR SS ST SU TQ TR
SV SW SX SY SZ TV

with SW covering most of Cornwall, TQ London, and HU the Shetlands. Clearly it could extend much further in each direction. Note that it has the neat feature that N and S are directly above each other, so that most Sx squares are in the south and most Nx squares are in the north.

BUGS

The conversions are only approximate. So after

($a1,$b1) = grid2ll(ll2grid($a,$b));

neither $a==$a1 nor $b==$b1. However abs($a-$a1) and abs($b-$b1) should be less than 0.00001 which will give you accuracy to within a few centimetres. Note that the error increases the further away you are from the reference point of your grid system.

When using ll2grid in scalar mode to get a "TQ999999" type of grid reference OSGB tends to round to the *nearest* 100m grid intersection rather than to the one to the left and below, this may cause your grid references to be to off by one in the last digit, but probably gives more consistent results overall.

The conversion of lat/long or grid to map sheets does not take account of inset areas on the sheets. So if you use ll2map() with the coordinates of the Scilly Isles, it will tell you that they are not on any Landranger sheet, whereas in fact the Scilly Isles are on an inset in the SW corner of Sheet 203. There is nothing in the design that prevents me adding the insets, they just need to be added as extra sheets with names like "Sheet 2003 Inset 1" with their own reference points and special sheet sizes. Collecting the data is another matter.

Not enough testing has been done.

EXAMPLES

This module is intended for use in the UK with the Ordnance Survey's National Grid, however the conversion routines are written in an entirely generic way that you can adapt to provide you with a Transverse Mercator Projection for any location on the the globe (or any other approximately spherical planet for that matter). What you need to do is to define the size of an ellipsoid that closely approximates your planet, and then define a suitable projection on top of it.

For the non-mathematicians among us, an ellipsoid is a sort of squashed ball that approximately represents the shape of the Earth. Using an approximation greatly simplifies working out latitudes and longitudes. Over the years cartographers have used a variety of different approximations that suit their local conditions and the tools they have at their disposal. You need to choose an appropriate approximation for your work. For example coordinates given by a GPS device will be based on one approximation, while coordinates you read from an 18th century chart might use another. You can find the numbers that define frequently used historical and modern approximations in a table such as the one given in the OS paper referenced above.

Choosing the exact ellipsoid model it a matter of some theoretical discussion, but for most practical purposes on earth the default values used by the UK ordnance survey will probably work for your data. However you may get better results setting the ellipsoid to the standard WGS84 shape, which is what the GPS satellite network uses and is appropriate for working with GPS-derived longitude and latitude coordinates.

To use WGS84, add the "set_ellipsoid" function to the list that you import from Geo::Coordinates::OSGB, and call it like this in your script:

set_ellipsoid(6378137.0,6356752.3141)

Next you have to define the projection that defines your local grid and the reference points for the grid. This can be done with the "set_projection" function, which you will also need to import. To get the parameters for "set_projection" you need to start by choosing an point in the middle of your area and finding the longitude and latitude. This is known as the True Origin of the Projection.

You then need to know (or invent) the grid reference for this point. This is entirely arbitrary unless you are trying to conform to someone else's grid. The simplest thing is to set the reference to 0,0. However you will then get negative coordinates for points south and west of your chosen centre point. It is more convenient to choose grid coordinates that make all of your results positive. To do this the eastings and northings should be a point just to the south of the southerly limit of your survey and just to the west of the westerly limit. This point will have grid coordinates (0, 0), and all your other grid coordinates will be positive. This point is known as the False Origin of the Projection.

So assuming your reference point is 150 degrees west and 47 degrees north, and you want this point to be grid reference 100000 east, 200000 north (in metres) then you should set

set_projection(47.0,-150.0,100000,200000,1)

The fifth parameter is the scale factor on the central meridian of the projection: for small areas just set this equal to 1. Making this parameter slightly less than 1 can reduce the scale distortion on the far east and west sides of the grid.

AUTHOR

Toby Thurston --- 7 Sep 2005

web: http://www.wildfire.dircon.co.uk

SEE ALSO

The UK Ordnance Survey's theory paper referenced above in "Theory".

See Geo::Coordinates::Convert for a general approach (not based on the above paper).

See Geo::Coordinates::Lambert for a French approach.