Use gdal to create raster with average values from point data

Question

I have dense points file that I'd like to rasterize using gdal in python. I'm relatively new to code and gdal and find the documentation confusing to navigate. I want to be able to rasterize these points to a 20m resolution raster that takes the average of points within each raster cell. I know ArcGIS Pro has this option through the 'pointstoraster' tool, where you can select the cell assignment type.

What's the open source equivalent of doing this?

• user30184: Try gdal_grid gdal.org/en/stable/programs/gdal_grid.html with "average". But I am not sure if you want to interpolate at all.
• monicaJ: I'm just looking to take the average value of a series of points that fall within the same pixel of a pre-defined resolution (in my case 20m) for a raster that covers the same extent as the points data. Would gdal.grid work that way?

Answer 1

Not sure if there's a simpler process, but you can use gdal_rasterize to calculate sum and count rasters from your points and then gdal_calc to calculate the average.

Given a point dataset with a "value" attribute field I want to average:

# Calculate sum with -a attribute and -add arguments, 
# note Float32 so we don't get caught by integer division or 
# have to cast to float in gdal_calc
gdal_rasterize -a value -add -tr 20 20 -ot Float32 points.shp sum.tif
# Calculate count with -burn 1 and -add arguments
gdal_rasterize -burn 1 -add -tr 20 20 -ot Int16 points.shp count.tif
# Calculate average
gdal_calc -A sum.tif -B count.tif --outfile=average.tif --calc="A/B"

Or in Python using osgeo.gdal.Rasterize and osgeo_utils.gdal_calc.Calc:

from osgeo import gdal
from osgeo_utils import gdal_calc
points = "points.shp"
sum = gdal.Rasterize(
 "sum.tif", points, 
 outputType=gdal.GDT_Float32, xRes=20, yRes=20, 
 attribute="value", add=True)
count = gdal.Rasterize(
 "count.tif", points, 
 outputType=gdal.GDT_Int16, xRes=20, yRes=20, 
 burnValues=[1], add=True)
average = gdal_calc.Calc(calc="A/B", A="sum.tif", B="count.tif", outfile="average.tif")