I recently went to considerable effort to stitch together NSW 1:25k e-Topo (PDF) maps, for use in OziExplorer. I figured it would be worthwhile to document the steps, with a bit of background. The same process may be useful for stitching other georeferenced maps with collars.
First thing is, it's worth remembering that the UTM system slices the world up into a series of 6 degree zones, which are shaped like sections of an orange. These are flattened, and any grid is overlaid on this. For NSW (and many other places), individual maps are a fixed number of degrees wide. In the case of standard NSW 1:25k maps, they are 0.25 degrees wide. But 6 degrees at the equator is much wider than 6 degrees at 80 degrees south. Or, in the case of NSW, 0.25 degrees at Murwillumbah (24.5km) in the north of the state is wider than 0.25 degrees at Eden (22.2km) in the south.
Each 4cm grid square is still 1km wide, so the Murwillumbah map (9541-2N MURWILLUMBAH) is physically wider than the Eden map (8823-1N EDEN)!
In fact, if you measure really carefully, you'll find that the top border of the maps is longer than the bottom border. Only by about 1.5mm or so, so barely noticeable. But it's important when it comes to stitching together maps to remember that the maps are not actually rectangular. They're more like a trapezium, for those who remember their high school maths.
The other thing that it's worth noting is that the maps are printed to be aligned to the lat/lon “grid”, not the UTM grid. The UTM grid is generally at an angle to the lat/lon. For NSW, this ranges from about +1.5 degrees to -1.5 degrees at the edges of the zones, changing slightly depending on how far north/south you are. So the maps needing to be joined together are all at slightly different angles!
The upshot of all of this is that maps are not nice rectangles that can just have the borders stripped off and dropped into a grid. However, there are some tools in the GDAL toolkit that will help us.
First step is to strip off the border of the PDF and convert to a raster format. We can do this for a single image using gdalwarp, by setting a few options, including:
- -t_srs <target spatial reference> - I've used EPSG:28356, which is MGA94 Zone 56 (eastern NSW). All maps need to be in the same co-ordinate system in order to merge them in the next step.
- -r <resampling_method> - cubic seems to be the recommended option given that we aren't trying to resize
- -cutline - applies a mask to the raster (confusingly, this doesn't actually crop the image - you also need “crop_to_cutline”)
- -cwhere - restrict desired cutline features based on attribute query - the cutline file is all of the NSW map boundaries, which can be handily downloaded from the https://maps.ozultimate.com/wiki/downloads page. The query is to restrict it to the individual map boundary
- -crop_to_cutline - crops to the mask applied under “cutline”
- -dstalpha - sets an output alpha band so that the pixels that are outside the cutline but inside the bounding box are transparent. This is needed when merging in the next step
Sample command in full - you'll need to set your own directories, both for input and output file, and for the location of the NSW map boundaries GeoJSON:
gdalwarp -t_srs EPSG:28356 -r cubic -cutline "E:\data\maps\nsw_map_boundaries\20160506_nsw_map_bounds.geojson" -cwhere "name = '9030-4S SPRINGWOOD'" -crop_to_cutline -dstalpha "E:\data\maps\maps\nsw\9030-4S SPRINGWOOD.pdf" "E:\data\maps\maps\nsw\9030-4S SPRINGWOOD.tif"
Unfortunately this approach turns a 10MB PDF into a 70MB GeoTiff!
You can batch automate the conversion of these in Windows by running a Powershell command eg
Get-ChildItem -Path "E:\data\maps\maps\nsw\*.pdf" | ForEach-Object {
$srcFile = $_.FullName
$destFile = $_.FullName -replace '.pdf', '.tif'
$shortName = $_.BaseName
&gdalwarp -t_srs EPSG:28356 -r cubic -cutline "E:\data\maps\nsw_map_boundaries\20160506_nsw_map_bounds.geojson" -cwhere "name = '$shortName'" -crop_to_cutline -dstalpha $srcFile $destFile
}
You now have a directory of GeoTiffs. In Windows, you can run the following at the command line to create a list of all the files in a text document:
dir /b /s "E:\data\maps\maps\nsw\*.tif" > "E:\data\maps\maps\nsw\tiff_list.txt"
Finally, you can run gdal_merge, using the list from the previous step to create one huge GeoTiff. Note that normal GeoTiffs max out at 4GB, so you need to be selective about how many files you include in the merge. If you try to include too many, you will get errors, and a partially completed output. I got away with 56 1:25k maps, but the actual number will depend on various factors. There is a BIGTIFF=YES option that can be added to the -co option, but I don't know how this works in downstream processes/systems.
gdal_merge -o "E:\data\maps\maps\nsw\sydney_basin.tif" -of GTiff -co COMPRESS=LZW --optfile "E:\data\maps\maps\nsw\tiff_list.txt"
For OziExplorer, the TIF can be used directly (use File→Import Map→Single DRG Map - and you will need the GeoTiff support library DLL), but it is faster and more space-efficient if converted to OZFX4 format. This can be done using Img2Ozf, downloadable from the OziExplorer website. Default settings are reasonable. As an example, a 4GB GeoTiff was converted to a 400MB OZFX4 file - a 90% reduction in size.