El uso de Rasterio o GDAL apilar varias bandas, sin el uso de subproceso comandos

Question

¿Alguien tiene una elocuente forma de apilar múltiples .tif archivos en múltiples banda de la pila utilizando Rasterio y/o GDAL?

Estoy buscando una manera de evitar el uso de un subproceso de comando como gdal_merge.py y no lo tiene como parte de mi secuencia de comandos de python.

Sé que ambos Rasterio y GDAL leer el .archivos tif como matrices, pero ¿cómo puedo pila de dichas matrices y escribe el resultado en forma de bandas independientes?

Answer 1

El uso de rasterio usted podría hacer

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):
        with rasterio.open(layer) as src1:
            dst.write_band(id + 1, src1.read(1))

Ella, por supuesto, asume que sus capas de entrada ya comparten la misma medida, la resolución y el tipo de datos

Answer 2

Si el uso de GDAL 2.1+ es tan simple como gdal.BuildVRT entonces 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)