[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