This will be a very quick post today. We wanted to use some historical USGS quads for analysis. We had a ready source, used ArcGIS to georeference them to a grid, but then wanted something approaching a seamless product.
First step, create an image to put the mosaic into:
gdal_merge -o usgs_merge.tif -createonly input1.tif input2.tif ...
Then mosaic into that. Since we have a shapefile with the shape of the bounds of each of the images in question, we can use a cutline to perform the mosaic, e.g.
gdalwarp -cutline quadrangle_index_24k.shp -cwhere "quadname_ = West_View" "West_View.tif" usgs_merge.tif
In all honesty, this was scripted in Excel for expediency, but could easily be done in Windows Command with for-in-do, or BASH at a UNIX machine.