[gdal-dev] Cannot use in memory raster object with gdal warp (GDAL Python API)
umbertofilippo
umbertofilippo at tiscali.it
Fri Feb 21 07:52:04 PST 2020
Sorry if this is a copy of a question I already asked on gis.stackexchange (
here
<https://gis.stackexchange.com/questions/351577/cannot-use-in-memory-raster-object-with-gdal-warp-gdal-python-api>
), but I thought I might get more chance to receive an answer here.
These are the operations I am doing in my code:
1 - reading an existing raster (tif)
2 - creating a in memory raster object
3 - reading the existing raster object in chunk as array
4 - making some calculations
5 - writing the array to my in memory object
What I want to do now is, with the in memory raster object ready to be used,
feed it to gdal warp, which should also output an intermediate raster object
that I'll be use afterwards.
However, when I try to GetRasterBand() out of this latter raster object, I
receive an error because it is actually equal to None:
AttributeError: 'NoneType' object has no attribute 'GetRasterBand'
I am pretty sure it is about saving/updating/flushing the status of the
raster object after the calculations but I could not find out if this is the
actual problem.
How can I use in memory (intermediate) files in my script?
Here is the code code I am using (the error is at the very last line):
from osgeo import gdal, gdalconst
import numpy as np
in_raster = gdal.Open(r'/my/raster.tif', gdalconst.GA_ReadOnly)
in_b = in_r.GetRasterBand(1)
b1 = in_raster.GetRasterBand(1)
xsize = b1.XSize
ysize = b1.YSize
step = 10 # THIS IS THE CHUNK SIZE!
ystep = ysize / step
ystep = int(ystep)
yresidual = ysize - (ystep * step)
driver = gdal.GetDriverByName("MEM")
temp_calc = driver.Create('',
xsize,
ysize,
1,
gdal.GDT_Float32)
temp_b = temp_calc.GetRasterBand(1)
for i in range(step):
if i != step-1:
arr_part = in_b.ReadAsArray(
0, ystep * i, xsize, ystep).astype('float32')
data = np.nan_to_num(arr_part)
temp_b.WriteArray(data, 0, ystep * i)
else:
arr_part = in_b.ReadAsArray(
0, ystep * i, xsize, ystep + yresidual).astype('float32')
# replace nan with zeros
#
https://docs.scipy.org/doc/numpy/reference/generated/numpy.nan_to_num.html
data = np.nan_to_num(arr_part)
# write array chunk to (in memory) file
# https://gdal.org/python/osgeo.gdal.Band-class.html#WriteArray
temp_b.WriteArray(data, 0, ystep * i)
# get raster resolution
in_gt = in_raster.GetGeoTransform()
in_x = in_gt[1]
in_y = -in_gt[5]
kwargs = {'xRes': in_x,
'yRes': in_y,
'targetAlignedPixels': True,
'resampleAlg': 'average',
'format': 'MEM'
}
ds = gdal.Warp('',
temp_calc,
**kwargs)
b1 = ds.GetRasterBand(1)
# ds is None and the above line throws the exception
--
Sent from: http://osgeo-org.1560.x6.nabble.com/GDAL-Dev-f3742093.html
More information about the gdal-dev
mailing list