[gdal-dev] Getting data values from a first GeoTIFF and do some map algebra type operations with a second GeoTIFF

[Travis Kirstine]

Have you looked at the samples at:
http://svn.osgeo.org/gdal/trunk/gdal/swig/python/samples/

The val_repl.py script show how to utilize numpy with gdal


On 23 November 2011 15:02, jdmorgan <jdmorgan at unca.edu> wrote:
> Hello,
>
> I am new to GDAL and raster programming in general and hoping that I can get
> some guidance. I have GeoTIFF which are projected in WGS84 and a single
> band.  I am trying to start simply with the concept of getting out the data
> values from the GeoTIFF from one file so that I compare them to another.
> However, I am a bit stuck at the point of getting out data values from a
> first GeoTIFF and do some map algebra type operations with a second
> GeoTIFF.  Here is what I have so far:
>
> from osgeo import gdal, gdalconst
>
> import os, sys, time
>
> import numpy
>
> ds = gdal.Open('C:\\Data\\GIS\\python\\gdal\\CMIP3-tb1_2021_2050.tif',
> gdal.GA_ReadOnly)
>
> cols = ds.RasterXSize
>
> rows = ds.RasterYSize
>
> bands = ds.RasterCount
>
> print cols, rows, bands
>
> gdal.AllRegister()
>
> gt = ds.GetGeoTransform()
>
> print gt
>
> xOrigin = gt[0]
>
> yOrigin = gt[3]
>
> pixWidth = gt[1]
>
> print "pixWidth " + str(pixWidth)
>
> pixHeight = gt[5]
>
> print "pixHeight " + str(pixHeight)
>
> for r in range(rows):
>
>               data = ds.GetRasterBand(1).ReadAsArray(0, r, cols, 1)
>
>               print str(data)
>
>
>
> Is this the correct way to be going about this?  Any guidance or best
> practices/algorithms would be greatly appreciated.
>
>
>
> Thanks,
>
> Derek
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>