Verwenden von Rasterio oder GDAL zum Stapeln mehrerer Bänder ohne Verwendung von Unterprozessbefehlen

Hat jemand eine eloquente Möglichkeit, mehrere .tif-Dateien mit Rasterio und / oder GDAL in einen Mehrbandstapel zu stapeln?

Ich suche nach einer Möglichkeit, die Verwendung eines Unterprozessbefehls wie gdal_merge.py zu vermeiden und ihn lieber als Teil meines Python-Skripts zu verwenden.

Ich weiß, dass sowohl Rasterio als auch GDAL die .tif-Dateien als Arrays lesen, aber wie staple ich diese Arrays und schreibe das Ergebnis als separate Bänder aus?

Verwenden rasterioSie könnten tun

import rasterio

file_list = ['file1.tif', 'file2.tif', 'file3.tif']

# Read metadata of first file
with rasterio.open(file_list[0]) as src0:
    meta = src0.meta

# Update meta to reflect the number of layers
meta.update(count = len(file_list))

# Read each layer and write it to stack
with rasterio.open('stack.tif', 'w', **meta) as dst:
    for id, layer in enumerate(file_list, start=1):
        with rasterio.open(layer) as src1:
            dst.write_band(id, src1.read(1))

Es wird natürlich davon ausgegangen, dass Ihre Eingabeebenen bereits den gleichen Umfang, die gleiche Auflösung und den gleichen Datentyp aufweisen

Wenn Sie GDAL 2.1+ verwenden, ist es so einfach wie gdal.BuildVRTdamals gdal.Translate:

from osgeo import gdal
outvrt = '/vsimem/stacked.vrt' #/vsimem is special in-memory virtual "directory"
outtif = '/tmp/stacked.tif'
tifs = ['a.tif', 'b.tif', 'c.tif', 'd.tif'] 
#or for all tifs in a dir
#import glob
#tifs = glob.glob('dir/*.tif')

outds = gdal.BuildVRT(outvrt, tifs, separate=True)
outds = gdal.Translate(outtif, outds)