[Gdal-dev] Python Problem

[Christopher Condit]

Hi Mateusz-
Thanks for your reply.  The code actually works fine as long as I don't put the bands into an array.  As soon as I add more than one band to the bands structure, I get the exception when attempting to ReadAsArray.  All of my rasters are exactly the same size (720x330).  
Here's a better piece of code illustrating the problem I'm having:

bands=[]
for timeSlice in range(60):
    file = "C:\\develop\\brewin\\vel-%i.tif" % timeSlice
    dataset = gdal.Open( file, GA_ReadOnly )
    bands.append(dataset.GetRasterBand(1))
for iter in range(len(bands)):    
    #The following line prints correctly (720x330)
    print "Dimension: (%i, %i)" % (bands[iter].XSize, bands[iter].YSize)
    #The following line throws 
    #TypeError: Access window out of range in RasterIO().  Requested
    #(0,0) of size 720x1 on raster of 65x16. 
    line = bands[iter].ReadAsArray( 0, 0, bands[iter].XSize, 1 )   

Thanks for any help,

Chris


-----Original Message-----
From: gdal-dev-bounces at lists.maptools.org [mailto:gdal-dev-bounces at lists.maptools.org] On Behalf Of Mateusz Loskot
Sent: Friday, May 19, 2006 8:10 AM
To: gdal-dev at lists.maptools.org
Subject: Re: [Gdal-dev] Python Problem

Christopher Condit wrote:
> Hi All-
> I know little to nothing about Python, but am using it to perform some
> simple math on a collection of raster images.  In this case, calculating
> mean and variance cellwise on a collection of 40 raster files...
> 
> Here's what my code looks like as I collect all of the input bands into
> an array:
>     magnitudeBands = []
>     thetaBands = []
>     for depth in range(40):
>         inFile = "C:\\develop\\data\\l%i\\vel-%i.tif" % (depth, month)
>         inDataset = gdal.Open( inFile, GA_ReadOnly )                
>         magnitudeBands.append(inDataset.GetRasterBand(1))
>         thetaBands.append(inDataset.GetRasterBand(2))        
> 
> ... At this point when I examine magnitudeBands, everything looks
> correct.
> 
>     for y in range(magnitudeBands[0].YSize):        
>         for band in range(len(magnitudeBands)):            
>             line = magnitudeBands[band].ReadAsArray( 0, y,
> magnitudeBands[band].XSize, 1, magnitudeBands[band].XSize, 1 )

You're misusing ReadAsArray method.
Here is small snippet that should will help you to understand how to use
ReadAsArray to read horizontal line:

#########################################################
import gdal
from gdalconst import *
import struct

inFile =  "/home/mloskot/dev/gdal/data/tester.tif"
inDataset = gdal.Open( inFile, GA_ReadOnly )
band = inDataset.GetRasterBand(1)

for y in range(band.YSize):
    # Read horizontal scanline
    line = band.ReadAsArray(0, y, band.XSize, 1)

    # Unpack block of bytes to python tuple
    tpl = struct.unpack('B' * band.XSize, line)

    print "Line: %d, length: %d" % (y, len(tpl))

    # Print scanline values
    #print str(tpl)
#########################################################

It should be easy to add nested loop to iterated through bands of raster.

> But at this point, I get an exception:
> TypeError: Access window out of range in RasterIO().  Requested
> (0,0) of size 720x1 on raster of 9479520x0.

This error says that you passed wrong arguments in 3rd and 4th and
latter positions.

Check RasterIO documentation to learn what is xoff, yoff, xsize, ysize.

http://www.gdal.org/classGDALRasterBand.html#1eef7a21cf4809425c3edced99aa088f

ReadAsArray uses RasterIO function internally.


Cheers

-- 
Mateusz Łoskot
http://mateusz.loskot.net
_______________________________________________
Gdal-dev mailing list
Gdal-dev at lists.maptools.org
http://lists.maptools.org/mailman/listinfo/gdal-dev