[gdal-dev] Raster algebra with Gdal and Python

[Alexander Bruy]

Hi all

I want to write a raster calculator and have two problems with Python
and GDAL.

My first problem: when calculating a difference between two raster bands I
get a wrong results. Raster in ERDAS IMAGINE format (HFA) and bands extracted
correctly. When I print each band and subtraction result I get something
similar to this

BAND 1
[[ 118  121  121 ...,  111  108  107]....]]
BAND 5
[[ 156  162  159 ...,  148  145  142]....]]
DIFFERENCE
[[218 215 218 ..., 219 219 221]....]]

But when calculate difference between two separate TIFFs (each band in
separate file) - result is ok

RASTER 1
[[ 118.  121.  121. ...,  111.  108.  107.]....]]
RASTER 5
[[ 156.  162.  159. ...,  148.  145.  142.]....]]
DIFFERENCE
[[-38. -41. -38. ..., -37. -37. -35.]....]]

Here is my code and sample files http://www.4shared.com/file/212381308/a3e3f140/gdal_numpytar.html.
Run examples as:
 - bands subtraction: test_bands.py landsat-6bands/s1986-06-12-163.img
 - files subtraction: test_files.py landsat-separate-files/band1.tif landsat-separate-files/band5.tif 

Can anyone help? What is wrong in my code?

Second problem: processing large rasters. I get error

File "C:/Documents and Settings/alex/.qgis//python/plugins\rastercalc\rastercalcutils.py", line 162, in layerAsArray
    array = gdalData.ReadAsArray()
File "C:\OSGeo4W\apps\gdal-16\pymod\osgeo\gdal.py", line 603, in ReadAsArray
    return gdalnumeric.DatasetReadAsArray( self, xoff, yoff, xsize, ysize )
File "C:\OSGeo4W\apps\gdal-16\pymod\osgeo\gdal_array.py", line 108, in DatasetReadAsArray
    return numpy.concatenate( array_list )
MemoryError

Large rasters is landsats with 6 bands, 30m/pixel, scene size is aprox. 185x185 km, aprox. 9000x9000 pixels. On disk one file occupied ~500-700 Mb. Some files
can be larger, up to 15000x15000 pix. and more.
Is it possible to work with large files? Maybe it is possible to process
small pieces of this files and then merge this pieces into one large file?
Where I can found more information about this?

Thanks

P.S.: sorry for my English, I'm Ukrainian

-- 
Alexander Bruy
alexander.bruy at gmail.com