Opened 9 years ago
Closed 9 years ago
#6072 closed defect (fixed)
Raster is initialized to 0 instead of no data value in memory file in Python
Reported by: | aly | Owned by: | warmerdam |
---|---|---|---|
Priority: | highest | Milestone: | |
Component: | default | Version: | 1.10.0 |
Severity: | critical | Keywords: | memory driver, ubuntu 12.04 LTS, precise, python 2.7, gdal |
Cc: |
Description
I am trying to reproject and resample (1000m ->> 30m) a raster image (shape 238X406) using gdal.ReprojectImage() function in Python. The input file has -9999 cells as no data value. The result is an array (shape 10834x15055). The data type is float32.
When I write the result to a geotiff file, everything is expected.No data value is set to -9999 and output array has -9999 cells. However, when I write the result to a gdal memory file to save some time (25 seconds shorter processing time, down from 105 seconds to 80 seconds), this time all -9999s (no data value) become 0 (zero) in the output array of memory file. Both results are exactly same, but GeoTIFF file has -9999, memory file has 0. Even though memory file has no data value of -9999, but the output array is initialized to 0.0 instead of -9999.
I am using the same code to produce the results and the only difference is that I call memory driver when I want to write result to memory file (driver=gdal.GetDriverByName("MEM")), and GeoTIFF driver when I want to write result to GeoTIFF file (driver=gdal.GetDriverByName("GTiff"))
My gdal version is '1100000' from gdal.VersionInfo(). My operating system us ubuntu 12.04 LTS "precise".
# Create a file
outFileRead=driver.Create(outFilePath,X,Y,1,dataType,options)
print inFileRead.GetRasterBand(1).GetNoDataValue()
print inFileRead.GetRasterBand(1).ReadAsArray()
print inFileRead.GetRasterBand(1).ReadAsArray().shape
# Reproject gdal.ReprojectImage(inFileRead,outFileRead,inProjection,outProjection,reSamplingType)
print outFileRead.GetRasterBand(1).GetNoDataValue()
print outFileRead.GetRasterBand(1).ReadAsArray()
print outFileRead.GetRasterBand(1).ReadAsArray().shape
the result of prints
-9999.0
[[-9999. -9999. -9999. ..., -9999. -9999. -9999.]
[-9999. -9999. -9999. ..., -9999. -9999. -9999.]
[-9999. -9999. -9999. ..., -9999. -9999. -9999.]
...,
[-9999. -9999. -9999. ..., -9999. -9999. -9999.]
[-9999. -9999. -9999. ..., -9999. -9999. -9999.]
[-9999. -9999. -9999. ..., -9999. -9999. -9999.]]
(406, 238)
-9999.0
[[ 0. 0. 0. ..., 0. 0. 0.]
[ 0. 0. 0. ..., 0. 0. 0.]
[ 0. 0. 0. ..., 0. 0. 0.]
...,
[ 0. 0. 0. ..., 0. 0. 0.]
[ 0. 0. 0. ..., 0. 0. 0.]
[ 0. 0. 0. ..., 0. 0. 0.]]
(15055, 10834)
This shows even if no data value is set to -9999 in memory file, the array is initialized to 0.0 instead of no data value. Array in GeoTIFF is initialized to no data value, -9999, as expected. I think this is a bug in the memory file.
Change History (3)
comment:2 by , 9 years ago
I posted that discussion. The current gdal release for Ubuntu 12.04 LTS "precise" is 1.10.0-1~precise1.
see it here https://launchpad.net/~ubuntugis/+archive/ubuntu/ppa?field.series_filter=precise
comment:3 by , 9 years ago
Resolution: | → fixed |
---|---|
Status: | new → closed |
I assume you manually set the nodata value to -9999 on the output dataset ? The GTiff driver might indeed initialize unset blocks to the nodata value, but this is more an exception than the general case. Starting with GDAL 2.0, gdal.ReprojectImage() will take into account the target nodata value if it has been set before. But you will generally have to initialize the output band with band.Fill(-9999) priorly to gdal.ReprojectImage().
Closing as fixed
See discussion in http://gis.stackexchange.com/questions/158503/9999-no-data-value-becomes-0-when-writing-array-to-gdal-memory-file
I am not sure if the analysis was right, but if it was, all you need to do is to update to recent GDAL version which is 2.0. You can ask for a backport to 1.11 series and it may or may not happen. GDAL 1.10 is out of life.
Thank you for a high quality report, even in this case it would have been better to write first a mail to gdal-dev and ask what developers think about it.