You can use the gdal.Dataset or gdal.Band ReadRaster method. See the GDAL and OGR API tutorials and the example below. ReadRaster does not use/require numpy, the return value is raw binary data and needs to be unpacked using the standard python struct module.
An example:
from math import floor
from osgeo import gdal,ogr
import struct
src_filename = '/tmp/test.tif'
shp_filename = '/tmp/test.shp'
src_ds=gdal.Open(src_filename)
gt_forward=src_ds.GetGeoTransform()
gt_reverse=gdal.InvGeoTransform(gt_forward)
rb=src_ds.GetRasterBand(1)
ds=ogr.Open(shp_filename)
lyr=ds.GetLayer()
for feat in lyr:
geom = feat.GetGeometryRef()
mx,my=geom.GetX(), geom.GetY() #coord in map units
#Convert from map to pixel coordinates.
px, py = gdal.ApplyGeoTransform(gt_reverse, mx, my)
px = floor(px) #x pixel
py = floor(py) #y pixel
structval=rb.ReadRaster(px,py,1,1,buf_type=gdal.GDT_UInt16) #Assumes 16 bit int aka 'short'
intval = struct.unpack('h' , structval) #use the 'short' format code (2 bytes) not int (4 bytes)
print intval[0] #intval is a tuple, length=1 as we only asked for 1 pixel value
Alternatively, since the reason you gave for not using numpy
was to avoid reading the entire array in using ReadAsArray()
, below is an example that does use numpy
and does not read the entire raster into memory, it only reads the raster value at the given point. It uses the built-in gdal.ApplyGeoTransform() function in order to deal with axes rotations.
from math import floor
from osgeo import gdal,ogr
src_filename = '/tmp/test.tif'
shp_filename = '/tmp/test.shp'
src_ds=gdal.Open(src_filename)
gt_forward=src_ds.GetGeoTransform()
gt_reverse=gdal.InvGeoTransform(gt_forward)
rb=src_ds.GetRasterBand(1)
ds=ogr.Open(shp_filename)
lyr=ds.GetLayer()
for feat in lyr:
geom = feat.GetGeometryRef()
mx,my=geom.GetX(), geom.GetY() #coord in map units
#Convert from map to pixel coordinates.
px, py = gdal.ApplyGeoTransform(gt_reverse, mx, my)
px = floor(px) #x pixel
py = floor(py) #y pixel
intval=rb.ReadAsArray(px,py,1,1)
print intval[0] #intval is a numpy array, length=1 as we only asked for 1 pixel value