There are various ways to serve rasters on Geoserver, but when dealing with hundreds or thousands of tiles, it is usually best to process imagery using GDAL. Although one could use Geoserver’s ImageMosaic plugin to serve tiles directly, I have not had much success with this approach when dealing larger datasets (e.g. >1000 tiles) . In my experience, using the GDAL command line tools to mosaic geotiff tiles into a single, seamless raster usually works best — and the added bonus is that the same mosaic can be loaded into other programs, like SAGA, PostgreSQL/PostGIS, or QGIS.
In addition to creating a seamless mosaic, there are two other requirements to create an image that serves well with Geoserver, namely build internal tiles and overviews. The whole process takes considerable time to process, especially to build overviews (third command), but in the end the image produced loads quickly in QGIS and serves from Geoserver quickly as a WMS. I do need to look into ways to increase the speed of this process by using multiple cores, but that will be for another post. For now, here is a methodology that works for me.
The Task
In Ontario, the most recent government commissioned aerial photographs are delivered in 1 km tiles with a pixel resolution of ~20cm — do the math, that’s a lot of pixels. Most organizations that manage these data must use hundreds or thousands of tiles. Serving them out as a WMS layer is the ultimate goal, so the tiles need to be turned into a seamless mosaic, internally tiled, and overviews must be produced.
Before you start
To use these commands, one needs to have the GDAL command line tools installed. Check the GDAL site for instructions on how to install the software. If you are a Windows user, note that the QGIS Standalone Installer contains a command line tool that includes GDAL utilities (among others).
Step 1: Mosaic tiles as a Virtual Raster Table
The creation of a virtual raster catalog (VRT format) is the first step. The creation of a VRT takes minutes as the gdalbuildvrt tool simply creates a table containing the spatial extents of and path to each inputted raster — in this case all .tif files within the current directory. However, just because you have created a VRT does not mean it will load in QGIS. Without building at least overviews, QGIS will attempt to load ALL tiles at once and eventually crash.
Navigate to directory that contains the imagery, and run the following command:
gdalbuildvrt mosaic_output.vrt *.tif
Step 2: Convert the VRT into a seamless GeoTIFF
The second step is to convert the virtual raster (mosaic_output.vrt) into a seamless GeoTIFF using gdal_translate. To do so, run the following command:
gdal_translate -b 1 -b 2 -b 3 -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co TILED=YES -co BIGTIFF=YES mosaic_output.vrt mosaic_rgb.tif
This command has many options:
- -b 1 -b 2 -b 3: this takes the first three bands of the imputed imagery, where b1 = red, b2 = blue, and b3 = green. If the imagery ONLY has three spectral bands, there is no need to add this option. You can also use this option to reorder spectral bands. If there is an infrared band 4, one could create a near-infrared image with -b 4-b 2 -b 3, where b4 = IR, b2 = blue, and b3 = green.
- -co COMPRESS=JPEG: this option compresses the output using JPEG compression. It has been stated elsewhere that imagery in Geoserver should not be compressed since it can slow down performance. In my experience, I do not mind sacrificing a little performance for an image that is not upwards of 100GB.
- -co PHOTOMETRIC=YCBCR: option to use the YCbCr color model.
- -co TILED=YES: this creates internal tiles.
- BIGTIFF=YES: if your GeoTIFF will be larger than 4GB, this option is required.
- mosaic_output.vrt: input the vrt file.
- mosaic_rgb.tif: the output file name.
Step 3: Build Overviews
The last step is to build overviews using gdaladdo, a tool that builds and rebuilds overviews.
gdaladdo --config COMPRESS_OVERVIEW JPEG --config PHOTOMETRIC_OVERVIEW YCBCR --config INTERLEAVE_OVERVIEW PIXEL -r average mosaic_rgb.tif 2 4 8 16 32 64 128 256 512 1024
As with the previous command, there are a number of options with this command:
- –config COMPRESS_OVERVIEW JPEG: this compresses the overviews with JPEG compression.
- –config PHOTOMETRIC_OVERVIEW YCBCR: as with the previous command, this option sets the overviews to the YCbCr color model.
- –config INTERLEAVE_OVERVIEW PIXEL: this creates internal overviews that will be appended to the original file, as opposed to being added to an external file.
- -r average: this option sets the resampling method. Nearest is the default if this option is not specified otherwise.
- mosaic_rgb.tif: this is the input raster. Note that there is no output file specified because this tool simply creates internal or external overviews.
- 2 4 8 16 32 64 128 256 512 1024: these are the overview levels.
And that’s it! In the end, I converted over 1000 1km tiles into a ~20GB GeoTIFF image that works nicely on Geoserver. This method was specifically used with Ontario’s SWOOP 2010 and 2015 imagery. Some minor tweaks may be needed for your imagery.
