Wednesday 2 September 2009

Processing Ordnance Survey Rasters with GDAL

Ordnance survey produce raster maps of the UK at several resolutions. They are provided as a set of tiles, geographically chunking the whole of the UK according to the British national grid. As with most things Ordnance survey, they seem to go out of their way to make them as hard to use in a standard GIS application as possible.

The set of raster tiles for the UK at 1:10k comes in at about 26GB. A whopping 10,000+ tiles. They do not come with any geo-referencing and are not made available in any format that includes this information (i.e Geotif). Instead, they are provided as tif, bmp or jpeg files with the option to download the geo-reference files (tfw etc) separately. To add insult to injury, it would seem that, in trying to save disk space, instead of a standard 3 band rgb tif file, what you actually get is a single band tif file with a colour palette. Its not so bad really, just another minor annoyance.

As I have mentioned before, we use a SharpMap based WPF control to render our data. This makes rendering OS tiles non-trivial for the following reasons:

  1. There is a bug in GdalRasterLayer which prevents it from rendering tif files that make use of colour palettes
  2. There is no dedicated layer type for rendering a raster tile index/set unless you want to serve them up over a WMS

Now, I was actually able to fix the problem in (1) and because of this (and with a better understanding of GDAL) I was able to render a raster tile set using nothing but a single GdalRasterLayer. However, maybe that's something I'll come back to in another post. For now, I want to talk about our original solution, which was to convert the tiles to ECW's using GDAL.

For 1:250k, 1:50k and 1:25k we found it was possible to create a single ECW for the whole of the UK. For 1:10k the huge number of tiles involved really put a strain on the memory. I don't know if it is achievable or not but one thing is for sure..it would take a very very long time (somewhere in the region of 1-2 weeks of solid processing). In the end, for the 1:10k stuff, we decided to create a single ECW for each top level square of the British national grid, giving us something like 50 ecw's in all. The time taken to process each of these regions varied between a few minutes and a few hours.

So, here are the steps involved to convert a set ordnance survey rasters into a single ECW:

  1. Copy all the tiles at a particular resolution (say 1:250k) into their own directory. e.g C:\250k\
  2. Download the geo-reference files from ordnance survey and copy them into the same directory as (1)
  3. Using the fwtools shell, cd to the directory with the raster tiles and the geo-reference files in.
  4. Build a virtual raster file(an xml file) of the entire set using:
    gdalbuildvrt 250k.vrt *.tif
  5. Use gdal_translate on the virtual raster file to convert the entire set of tiles into a single ECW:
    gdal_translate -of ECW -co LARGE_OK=YES -expand rgb -a_srs EPSG:27700 250k.vrt 250k.ecw
Issues, workarounds and a handy bit of C#

Generally speaking, this worked well. The only real problem I faced was when dealing with the 10k national set. When trying to create a virtual raster of a top level grid square with a command like:
gdalbuildvrt SU.vrt SU*.tif
occasionally the command would fail. This would happen if the tiles in question had differing colour palettes or a differing number of bands.


The solution was to introduce an intermediate stage whereby the source rasters were 'normalised' to 3 band rgb tif files. I achieved this by using the following command on each file:
gdal_translate -co COMPRESS=DEFLATE -co TFW=YES -expand rgb HP60NE.tif Expanded\HP60NE.tif
Then, when it came to making the final ECW the -expand option had to be left out, due to the palette expansion having already happened in the previous stage.

I wrapped all of this up into a neat little C# program (mainly because I couldn't be bothered to work out the java script or batch file commands) that checked the integrity of the files, 'normalised' their palettes if necessary, built the vrt and finally, converted the raster tiles into a single ECW file, logging all the results into a separate file to record timings and errors.

The resulting tiles2ecw utility can be downloaded from here. You can either compile it into a console application or run it directly from the command line using the rather brillient cs-script.

No comments:

Post a Comment