The Perl Toolchain Summit needs more sponsors. If your company depends on Perl, please support this very important event.

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), from the OSGB grid to latitude and longitude.

Used on their own, these modules will allow you convert accurately between a grid reference and lat/lon coordinates based on the OSGB Airy 1830 geoid model (the traditional model used for maps in the UK for the last 180 years) and last amended by the OS in 1936. This model is sometimes referred to as OSGB36.

OSGB36 fits the British Isles very well, but is rather different from the WGS84 model that has rapidly become the de facto universal standard model thanks to the popularity of GPS devices and maps on the Internet. So, if you are trying to translate from a OSGB grid reference to lat/lon coordinates that can be used in in Google Earth, Wikipedia, or some other Internet based tool, you will need to do two transformations. First translate your grid ref into OSGB lat/lon, then nudge the result into WGS84. Routines are provided to do both of these operations, but they are only approximate. The inaccuracy of the approximation varies according to where you are in the country but may be as much as several metres in some areas.

To get really accurate results you need to combine this module with its companion Geo::Coordinates::OSTN02 which implements the transformation (known as OSTN02) that now defines the relationship between GPS survey data based on WGS84 and the British National Grid. Using this module you should be able to get results that are accurate to within a few centimetres, but it is a little bit slower and requires a bit more memory to run.

Version: 2.01

SYNOPSIS

  use Geo::Coordinates::OSGB qw(ll_to_grid grid_to_ll);

  # basic conversion routines
  ($easting,$northing) = ll_to_grid($lat,$lon);
  ($lat,$long) = grid_to_ll($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

  # latitude and longitude are returned according to the OSGB36 (as printed on OS maps)
  # if you want them to work in Google Earth or some other tool that uses WGS84 then
  # you need to adjust the results
  ($lat, $long, $elevation) = shift_ll_into_WGS84($lat, $long, $elevation);

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". There are some key concepts explained in that section that you need to know in order to use these modules successfully, so you are recommended to at least skim through it now.

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:

    ll_to_grid
    grid_to_ll

    parse_grid
    parse_trad_grid
    parse_GPS_grid
    parse_landranger_grid
    format_grid_trad
    format_grid_GPS
    format_grid_landranger

    parse_ll
    format_ll_dms
    format_ll_ISO

    shift_ll_into_WGS84
    shift_ll_from_WGS84

None of these is exported by default.

This code is fine tuned to the British national grid system. You could use it elsewhere but you will need to adapt it. Some starting points for doing this are explained in detail in the "Examples" section below.

ll_to_grid(lat,lon)

When called in a void context, or with no arguments ll_to_grid does nothing. When called in a list context, ll_to_grid 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) = ll_to_grid('+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 of the British Grid (which 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 ll_to_grid in a scalar context, it will automatically call format_grid_trad. For example:

    my $gridref = ll_to_grid('+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 get the grid reference formatted differently, just wrap it in the appropriate format routine, like this:

    $gridref = format_grid_GPS(ll_to_grid('+5120','-00025'));
    $gridref = format_grid_landranger(ll_to_grid('+5120','-00025'));

It is not needed for any normal work, but ll_to_grid() also takes an optional third argument that sets the ellipsoid model to use. This normally defaults to "OSGB36", the name of the normal model for working with British maps. If you are working with the highly accurate OSTN02 conversions supplied in the companion module in this distribution, then you will need to produce pseudo-grid references as input to those routines. For these purposes you should call ll_to_grid() like this:

    my $pseudo_gridref = ll_to_grid('+5120','-00025', 'WGS84');

and then transform this to a real grid reference using ETRS89_to_OSGB36() from the companion module.

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 5m even with perfect reception. Most GPS devices can display the accuracy of the current fix you are getting, but you should be aware that all normal consumer-level GPS devices can only ever produce an approximation of an OS grid reference, no matter what level of accuracy they may display. The reasons for this are discussed below in the section on Theory.

format_grid_landranger(e,n)

This routine 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. It is easy to work out the coordinates of the other corners, because all OS Landranger maps cover a 40km square (if you don't count insets or the occasional sheet that includes extra details outside the formal margin).

parse_grid('string')

Attempts to match a grid reference some form or other in the input string and will then call the appropriate grid parsing routine from those defined above. In particular it will parse strings in the form '176-345210' meaning grid ref 345 210 on sheet 176, as well as 'TQ345210' and 'TQ 34500 21000' etc.

grid_to_ll(e,n) or grid_to_ll(grid_ref)

When called in list context grid_to_ll 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 must be an (easting, northing) pair representing the absolute grid reference in metres from the point of origin. You can get these from a grid reference string by calling parse_grid() first.

An optional third argument defines the geoid model to use just as it does for ll_to_grid(). This is only necessary is you are working with the pseudo-grid references produced by the OSTN02 routines. See Theory for more discussion.

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.ordnancesurvey.co.uk/oswebsite/gps/information/index.html

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 British National Grid

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.

Geoid models

This section explains the fundamental problems of mapping a spherical earth onto a flat piece of paper (or computer screen). A basic understanding of this material will help you use these routines more effectively. It will also provide you with a good store of ammunition if you ever get into an argument with someone from the Flat Earth Society.

It is a direct consequence of Newton's law of universal gravitation (and in particular the bit that states that the gravitational attraction between two objects varies inversely as the square of the distance between them) that all planets are roughly spherical. (If they were any other shape gravity would tend to pull them into a sphere). Most useful surfaces for displaying maps (such as pieces of paper or screens) on the other hand are flat. There is therefore a fundamental problem in making any maps of the earth that its curved surface being mapped must be distorted at least slightly in order to get it to fit onto the flat map.

This module sets out to solve the corresponding problem of converting latitude and longitude coordinates (designed for a spherical surface) to and from a rectangular grid (for a flat surface). This projection is in itself is a fairly lengthy bit of maths, but what makes it complicated is that the earth is not quite a sphere. Because our planet is also spinning about its axis, it tends to bulge out slightly in the middle, so it is more of an oblate spheroid than a sphere. This makes the maths even longer, but the real problem is that the earth is not a regular oblate spheroid either, but an irregular lump that closely resembles an oblate spheroid and which is constantly (if slowly) being rearranged by plate tectonics. So the best we can do is to pick an imaginary regular oblate spheroid that provides a good fit for the region of the earth that we are interested in mapping. The British Ordnance Survey did this back in 1830 and have used it ever since (with revisions in 1936) as the base on which the National Grid for Great Britain is constructed. You can also call an oblate spheroid an ellipsoid if you like. The ellipsoid model that the OS defined is called OSGB36 for short, and it's parameters are built into these modules.

The general idea is that you can establish your latitude and longitude by careful observation of the sun, the moon, the planets, or your GPS handset, and that you then do some clever maths to work out the corresponding grid reference, using a suitable idealised ellipsoid model of the earth (which is generally known as a "geoid"). These modules let you do the clever maths, and the model they use is the OSGB36 one. Unfortunately, while this model suits Britain very well, it is less useful in the rest of the world, and there are many other models in use in other countries. In the mid-1980s a new standard geoid model was defined to use with the fledgling global positioning system (GPS). This model is known as WGS84, and is designed to be a compromise model that works equally well for all parts of the globe (or equally poorly depending on your point of view). WGS84 has grown in importance as GPS systems have become consumer items and useful global mapping tools (such as Google Earth) have become freely available through the Internet. Most latitude and longitude coordinates quoted on the Internet (for example in Wikipedia) are WGS84 coordinates. All of this means that there is no such thing as an accurate set of coordinates for every unique spot on earth. There are only approximations based on one or other of the accepted geoid models, however for most practical purposes good approximations are all you need. In Europe the official definition of WGS84 is sometime referred to as ETRS89. For all practical purposes, you can regard ETRS89 as identical to WGS84.

So, if you are working exclusively with British OS maps and you merely want to convert from the grid to the latitude and longitude coordinates prined (as faint blue crosses) on those maps, then all you need from these modules are the plain grid_to_ll() and ll_to_grid() routines. On the other hand if you want to produce latitude and longitude coordinates suitable for Google Earth or Wikipedia from a British grid reference, then you need an extra step. Convert your grid reference using grid_to_ll() and then shift it from the OSGB36 model to the WGS84 model using shift_ll_into_WGS84(). To go the other way round, shift your WGS84 lat/lon coordinated into OSGB36, using shift_ll_from_WGS84(), before you convert them using ll_to_grid().

If you have a requirement for really accurate work (say to within a millimetre or two) then the above routines may not satisfy you, as they are only really accurate to within a metre or two. Instead you will need to use the OS's newly published transformation matrix called OSTN02. This monumental work re-defines the British grid interms of offsets from WGS84 to allow really accurate grid references to be determined from really accurate GPS readings (the sort you get from professional fixed base stations, not from your car's sat nav or your hand held device). The problem with it is that it defines the grid in terms of a deviation in three dimensions from a pseudo-grid based on WGS84 and it does this separately for every square metre of the country, so the data set is huge and takes several seconds to load even on a fast machine. Nevertheless a Perl version of OSTN02 is included as a seperate module in this distribution just in case you really need it (but you don't need it for any "normal" work). Because of the way it is defined, it works differently from the approximate routines described above.

Starting with a really accurate lat/lon reading in WGS84 terms, you need to transform it into a pseudo-grid reference using ll_to_grid() with an optional argument to tell it to use the WGS84 geoid parameters instead of the default OSGB36 parameters. Geo::Coordinates::OSTN02 then provides a routine called ETRS89_to_OSGB36() which will shift this pseudo-grid reference into an accurate OSGB grid reference. To go back the other way, you use OSGB36_to_ETRS89() to make a pseudo-grid reference, and then call grid_to_ll() with the WGS84 parameter to get WGS84 lat/long coordinates.

BUGS

The conversions are only approximate. So after

  ($a1,$b1) = grid_to_ll(ll_to_grid($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 the grid system.

When using ll_to_grid 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 any other ellipsoid model that is suitable for your local area of the earth.

Once you have defined the ellipsoid to use in terms of its major and minor diameters you have to define the projection that defines your local grid and the reference points for the grid. 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.

In Britain the True Point of Origin is a point (near the northern coast of France) at 49° North and 2° West (in the OSGB36 model), and this is given the grid coordinated (400000,-100000) so that the False point of origin, point (0,0) is just to the SW of the Scilly Isles.

AUTHOR

Toby Thurston --- 6 Sep 2007

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.

1 POD Error

The following errors were encountered while parsing the POD:

Around line 1269:

Non-ASCII character seen before =encoding in '49°'. Assuming CP1252