Writing a raster file from multi-dimension array using Python

Question

I've been attempting to create a raster file that has multiple bands (10 or more) from a numpy array with multiple dimension(i.e. 10 or more). The script I'm currently using exports the TIFF file I want to produce but the file turns out to be black (blank) image.

The shape of the data has a dimension of (1120,872,12) (rows, cols, bands) and I wanted to create a raster file that has 12 bands.

input_array.shape
 (1120, 872, 12)

---

def CreateTiff(output, array, driver, noData, GeoT, Proj, DataType):
 rows = array.shape[0]
 cols = array.shape[1]
 band = array.shape[2]
 array[np.isnan(array)] = noData
 driver = gdal.GetDriverByName('GTiff')
 DataSet = driver.Create(output, rows, cols, band, gdal.GDT_Float32)
 DataSet.SetGeoTransform(GeoT)
 DataSet.SetProjection(Proj)
 for i, image in range(array.shape[2], 1):
 DataSet.GetRasterBand(i).WriteArray(image)
 DataSet.GetRasterBand(i).SetNoDataValue(NoData)
 DataSet.FlushCache()
 return Name

with this function,

# the raster output layer
output_file = r'rasterfile.tiff'
# Writes raster
CreateTiff(output_file, input_array, driver, noData, GeoTransform, Projection, DataType)

As I mentioned above, this function creates the raster file I was hoping to get but when I inspect it in QGIS, the file is a black(blank) image and the values are zeros.

Any thoughts?

Answer 1

Your range loop is empty:

list(range(array.shape[2], 1)
[]

Try:

for i in range(band):
 DataSet.GetRasterBand(i+1).WriteArray(array[:, :, i])
 DataSet.GetRasterBand(i+1).SetNoDataValue(NoData)

• Destaneon: hey @Luke, thanks for point out the issue with the range loop I used. when I used your suggestion, this error returns: AttributeError: 'NoneType' object has no attribute 'WriteArray'
• ReutKeller: @Destaneon i'm facing similar problem of saving ndarray as raster with many bands, I have tried your function but i'm not sure what is the parameter "DataType" , could you maybe elaborate what should be this parameter? :) thanks
• user2856: @Reut the DataType in Destaneon's is probably supposed to be one of the gdal.GDT_* data type constants, but the function doesn't actually use this argument, gdal.GDT_Float32 is hardcoded instead.