Look at the documentation for driver.Create. Your current code is trying to write a 3D array to a single band. That is what raise ValueError("expected array of dim 2") is telling you. You need to create a tiff with bands equal to the number of "stacked images". Then write a single image to each band. You cannot write all the bands at once.
Here is the same function you wrote, updated to reflect this change. There are also modifications because NewFileName isn't the parameter name....
def CreateGeoTiff(Name, Array, driver, NDV,
xsize, ysize, GeoT, Projection, DataType):
Array[np.isnan(Array)] = NDV
DataSet = driver.Create(Name, xsize, ysize, Array.shape[0], DataType)
DataSet.SetGeoTransform(GeoT)
DataSet.SetProjection( Projection.ExportToWkt() )
for i, image in enumerate(Array, 1):
DataSet.GetRasterBand(i).WriteArray( image )
DataSet.GetRasterBand(i).SetNoDataValue(NDV)
DataSet.FlushCache()
return Name
You can make this less verbose because your array will tell you the xsize and ysize.
def CreateGeoTiff(Name, Array, driver, NDV,
GeoT, Projection, DataType):
Array[np.isnan(Array)] = NDV
DataSet = driver.Create(Name, Array.shape[2], Array.shape[1], Array.shape[0], DataType)
DataSet.SetGeoTransform(GeoT)
DataSet.SetProjection( Projection.ExportToWkt() )
for i, image in enumerate(Array, 1):
DataSet.GetRasterBand(i).WriteArray( image )
DataSet.GetRasterBand(i).SetNoDataValue(NDV)
DataSet.FlushCache()
return Name
You can read specific bands in a single call using rasterio by passing a list/tuple of band numbers (Following the GDAL convention, bands are indexed from 1):
import rasterio
rasterio.__version__
'1.0a8'
dataset = rasterio.open('multiband.tif')
dataset.count
4
dataset.read((1,2)) #read 1st two bands into an array.
array([[[ 85, 98, 75, ..., 53, 55, 55],
[ 84, 94, 76, ..., 54, 55, 54],
[ 68, 60, 55, ..., 53, 54, 53],
...,
[ 67, 67, 66, ..., 63, 63, 62],
[255, 255, 255, ..., 255, 255, 255],
[255, 255, 255, ..., 255, 255, 255]],
[[ 78, 88, 65, ..., 41, 43, 45],
[ 77, 84, 66, ..., 42, 43, 44],
[ 61, 51, 46, ..., 41, 42, 43],
...,
[ 77, 77, 77, ..., 71, 70, 69],
[255, 255, 255, ..., 255, 207, 255],
[255, 255, 255, ..., 191, 0, 135]]], dtype=uint8)
Best Answer
You can write to selected bands using the
indexeskeyword.If you want your array to contain only zeros you can also create a new numpy array. This way you don't need to read any data at all.