I have several DEM rasters in different projections and at different resolutions. The rasters are not rectangular, but they have appropriate NoData values. They overlap enough that I should be able to make a seamless dataset across the entire area, but is there a way I can use the GDAL tools to warp/mosaic them to get that seamless dataset?
I'm currently using the following command line to warp and mosaic the datasets:
gdalwarp -t_srs epsg:3857 -tr <calculated-x> <calculated-y> -r cubicspline output.tif input1.tif input2.tif ...
The problem is that the edges of the raster datasets don't always match each other, so I get effective cliffs at the edges. Is there a way to average them together for a few pixels near the edges of data for the dataset on top?
I'd really prefer to be able to do this in batch via command line or Python bindings.
Here's an (extreme) example of what I'm trying to avoid:
I'm not hugely fond of the ridges at the edge there.
I'm using the USGS NED 1/3 data as a seamless starting point, so I'm guaranteed to have something at every pixel (and the higher-resolution datasets that don't cover the entire area will always have something to be blended with).