1

qgis screenshot of output GeoTiffI would like to write an array to a geotiff with Python script using the GDAL module. When I try to do something simple like copying an existing GeoTiff via reading it to an array, then, writing to a new array, QGIS reads the integer value 1 as "nan". I have verified that the actual values that QGIS presents as "nan" are the integer value 1 in the GeoTiff. Is there some additional information that I need to write to the GeoTiff file so that QGIS will read my integer values of 1 as 1? My code is copied below:

import os
import gdal
os.chdir(r'C:\Users\pipi\Documents\Rogaine\Tarlo\gpx')

in_ds = gdal.Open('0011aa0011aa.tif')
in_band = in_ds.GetRasterBand(1)
noData = in_band.GetNoDataValue()
if noData:
    pass
else:
    noData = 0

cols = in_band.XSize
rows = in_band.YSize

# Creates an output file with the same dimensions and data type.
out_ds = gdal.GetDriverByName('GTiff').Create(
'master.tif', cols, rows, 1, gdal.GDT_Byte)
out_ds.SetProjection(in_ds.GetProjection())
out_ds.SetGeoTransform(in_ds.GetGeoTransform())
out_band = out_ds.GetRasterBand(1)

data = in_band.ReadAsArray()
out_band.WriteArray(data)
out_band.SetNoDataValue(noData)

del out_ds
4
  • It sounds like the 'noData' variable being used to set the no data value on the output tif is 1. Have you tried either not setting the nodata value for 'out_band' at all or changing it to a different value (e.g. out_band.SetNoDataValue(0) )? ...Changed that value to 0 because the datatype is a byte
    – Josh K
    Jun 27, 2016 at 3:41
  • Currently I am setting the NoDataValue to 0 as the "input dataset" does not have a NoDataValue. If I don't set the NoDataValue (hence noData = None), then, qgis does not show any contrast, with a legend containing values of 0 and 0 instead of 'nan' and 1. Jun 27, 2016 at 3:50
  • Try forcing the computation of stats on the dataset, I've had display issues in the past.
    – user10353
    Jun 27, 2016 at 4:10
  • When I use the statement out_band.ComputeStatistics(False) to the script it does not change the outcome. Is this what you mean by "try forcing the computation of stats"? Jun 27, 2016 at 4:17

2 Answers 2

2

There is no fault in the data or python script. QGIS's default for opening a raster dataset looks for a minimum and maximum value to provide band rendering. However, this dataset contains a single value. Within QGIS's layer properties style menu by setting both the Min and Max values of the dataset to the single dataset value and selecting "No enhancement" one can obtain a plot with a suitable legend.

1

Please try to set array nan's to your no data value before writing to band.

import numpy
...
data = in_band.ReadAsArray()

data[data==numpy.nan] = noData

out_band.WriteArray(data)
...

If the first dataset do not have a no data value set, set it first.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.