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.
Edit:
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).