Merging GDAL tiles and crop via bounding box?

Question

I am using GDAL.

Assume I have four adjacent GDAL tiles in the formation of a square. I have a bounding box with a corner in each tile so that the rectangle formed spans all four tiles. My goal is to generate a DEM from those four adjacent GDAL files that is cropped to the bounding box. Also as an added bonus, assume these GDAL tiles are in folder of GDAL tiles. In this case, the DEM is a Geotile.

In essence I am asking how to perform three steps here:

1. Find all GDAL tiles covered by the bounding box.
2. Merge the tiles into one giant DEM
3. Crop the DEM to the bounding box and output as GEOTIFF.

The methods I know to do this would be really efficient (ex: open each GDAL sequentially and find all the ones that contain a corner point, merge tiles, then crop, surely GDAL has a more efficient means of doing this).

Answer 1

This can be achieved with the help of GDAL's Virtual Raster Format. With this you can essentially skip the step of creating one giant DEM. The VRT will be handled by GDAL like a giant, merged DEM but is just a small XML file containing the file paths for each tile as well as some metadata. This can then be fed to gdalwarp together with a bounding box or a cutline.

This assumes that your tiles do not overlap. My guess is you are probably working with SRTM tiles for which this is the case.

1. Merge the tiles into one giant, virtual DEM

gdalbuildvrt big_DEM.vrt tile_folder/*.tif

2. Crop the DEM to your liking and write the result to a new file. In this case with a cropline vector dataset.

gdalwarp -of GTiff -cutline region_of_interest.kml -crop_to_cutline big_DEM.vrt small_DEM.tif

Answer 2

The accepted answer is probably superior, but I just chose a quick and dirty method to merge files.

First I build a index file using gdaltileindex. I then use gdal tools to returns a list of files that the bounding box intersects with. I then use gdal_merge.py to merge them the intersecting DEMs into one large dem and used the crop parameter of gdalmerge.py to crop to the bounding box.