How to set a spatial filter with Python/GDAL?

Question

Is there a way with GDAL (Python) to set a spatial filter (bounding-box) on a raster before I process it? I want to read it as an array.

BandReadAsArray(band, xoff=0, yoff=0, win_xsize=None, win_ysize=None, buf_xsize=None, buf_ysize=None, buf_obj=None)

Answer 1

You can do it by converting your bounding box to pixel offset and size:

from osgeo import gdal

def map_to_pixel(mx,my,gt):
    #Convert from map to pixel coordinates.
    #Only works for geotransforms with no rotation.
    px = int((mx - gt[0]) / gt[1]) #x pixel
    py = int((my - gt[3]) / gt[5]) #y pixel

    return px,py

def extent_to_offset(xmin,ymin,xmax,ymax,gt):
    pxmin,pymin = map_to_pixel(xmin,ymin,gt)
    pxmax,pymax = map_to_pixel(xmax,ymax,gt)

    return pxmin,pymin,pxmax-pxmin,pymax-pymin

xmin,ymin,xmax,ymax = [123,234,345,456]
ds=gdal.Open('some.tif')
gt = ds.GetGeoTransform()
xoff, yoff, xsize, ysize = extent_to_offset(xmin,ymin,xmax,ymax,gt)
band=ds.GetRasterBand(1)
array=band.ReadAsArray(xoff, yoff, xsize, ysize)