I am working on a computational model of the abundance of wild pollinators across a landscape. The model itself is complete, and I am now struggling with a post-processing step.
I have my GDAL pollinator supply raster that looks something like this (lighter colors mean higher pollinator visitation to a pixel):
And I have an OGR shapefile of points representing sample locations on the landscape:
I'm trying to run some analysis on the pixels under these points, but to do so, I need to be able to extract the value of a pixel under a point.
Is it possible to extract the value of a pixel under a point using only OGR and GDAL through Python? I would prefer to avoid reading the entire raster into memory through ReadAsArray()
, as my output rasters are very, very large (too large to fit into memory).
I noticed this post, which is similar, but requires a command-line call.
Best Answer
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:
Alternatively, since the reason you gave for not using
numpy
was to avoid reading the entire array in usingReadAsArray()
, below is an example that usesnumpy
and does not read the entire raster in. It uses the built-in gdal.ApplyGeoTransform() function in order to deal with axes rotations.