Reprojecting and Aggregating Rasters with GDAL

When working with raster datasets of different projections and resolutions, it is often desirable to reproject them to the same projection and align them to the same pixel grid. In this post, we will explore the recently introduced options in the open-source GDAL utility gdalwarp that makes this process much simpler and efficient. In particular, we will be exploring the -r sum (Resample with Sum), -r average (Resample with Average) and -tap (Target Aligned Pixels).

We will take the following 3 raster datasets and clip, resample and align them to a common pixel grid.

  • LandScan Global: A high-quality global population grid that is available at 1km resolution in the geographic CRS WGS84 Lat/Lon (EPSG:4326).
  • GHS Population Grid: A 100m resolution global population dataset that is distributed in the World Mollweide Equal Area Projection (ESRI:54009).
  • NLCD Tree Canopy Cover: A 30m resolution gridded dataset with percent canopy estimate of tree cover in the NAD83 CONUS Alberts Projection (EPSG:5070).

As you can see we have datasets that have widely varying pixel sizes and projections. If we wanted to compare them with each other – we must first harmonize them on a unified pixel grid. We will learn how to reproject, resample and align these to the NAD83 California Albers Projection (EPSG:3311) and at 1km resolution.

This post assumes familiarity with GDAL command-line utilities. You can review my course material for the Mastering GDAL Tools course for installation and usage instructions.

Clipping Rasters

We will first clip the rasters to the boundary of California using a shapefile obtained from US Census Bureau. We specify the polygon shapefile california.shp using the -cutline argument.

gdalwarp -cutline california.shp -crop_to_cutline \
	landscan.tif landscan_clipped.tif
gdalwarp -cutline california.shp -crop_to_cutline \
	ghsl.tif ghsl_clipped.tif
gdalwarp -cutline california.shp -crop_to_cutline \
	nlcd_tcc.tif nlcd_tcc_clipped.tif

Resampling and Aggregating Rasters

When reprojecting rasters to a target projection, we must specify the resampling method. The default resampling method used by GDAL (and many other programs) is Nearest Neighbor. This method is fast and works ok where the source and target pixel sizes are similar. If you are downsampling a 30m raster to 1km – nearest neighbor resampling will give you wrong results. Also when resampling population datasets, we must sum all the pixels within the target pixels to get an accurate value of the target pixel.

GDAL version 3.1 introduced a new resampling method sum. This method is useful for aggregating population datasets where the value of the target pixel is computed by summing the source pixel that fall within the target pixel. This method uses a weighted-sum – meaning the value of the pixel will be weighted by the overlap of the pixel with the target pixel. (i.e. if the pixel is only 50% within the target pixel, the sum will consider only half the value). Another useful option when trying to align pixel grids is the -tap (target aligned pixels) which makes sure that there is no sub-pixel shift between rasters with same projection and resolution.

Let’s apply these to resample and reproject the population rasters to a common projection of EPSG:3311 and pixel grid of 1km x 1km.

gdalwarp -t_srs EPSG:3311 -tr 1000 1000 -r sum  -tap \
	landscan_clipped.tif landscan_resampled.tif
gdalwarp -t_srs EPSG:3311 -tr 1000 1000 -r sum -tap \
	ghsl_clipped.tif ghsl_resampled.tif 

Lastly, we have the NLCD Tree Canopy Cover (TCC) raster. This raster contains pixels indicating the percentage tree cover within each 30m x 30m pixel. When we downsample this data to a 1km x 1km grid, we can use the average resampling method. This will give us the fractional tree cover for the target pixel by averaging all the source pixels that are contained within.

gdalwarp -t_srs EPSG:3311 -tr 1000 1000 -r average \
	nlcd_tcc_clipped.tif nlcd_tcc_resampled.tif -tap

The result of the resampling operation is that now we have scaled and aligned the rasters where they have the same pixel grid and can be compared with each other.

Leave a Reply