I am starting with a PostGIS raster table containing anywhere from 6-12 bands. In a nutshell, the goal is to take each of these bands, perform a 3x3 neighborhood average on each pixel in the band, and output a new 12-band raster with the neighborhood-averaged bands and the original geometry intact. I could likely do this with pure GDAL, but the rasters at this stage in the process have already been imported into the database and the original GeoTIFF files deleted, so using GDAL would likely require an infrastructure change.
My current concept uses ST_MapAlgebra to take the bands from the source raster in one-by-one and perform the 3x3 average to return the single averaged band. I then would use ST_Union on a setof these averaged bands to create a single raster. I then need to somehow generate a new table from the output of the ST_Union. I am trying to do this in as few steps as possible, so ideally the process would take 1-2 queries.
My ST_MapAlgebra callback function is this:
CREATE OR REPLACE FUNCTION
callback_fn(pixel float[], pos integer[], variadic userargs text[])
returns float
language plpgsql
immutable
as $$
declare
avg float;
x integer;
y integer;
p integer[];
begin
x = pos[1][1];
y = pos[1][2];
p = pixel[1];
avg := (p[x-1][y-1] + p[x][y-1] + p[x+1][y-1]
+ p[x-1][y] + p[x][y] + p[x+1][y]
+ p[x-1][y+1] + p[x][y+1] + p[x+1][y+1]) / 9.0;
return avg;
end
$$;
My current attempt at an execution query is definitely not syntactically correct, but it gets the point across of what I'm trying to do:
WITH rast AS (
SELECT rast FROM "raster2"
)
CREATE TABLE "raster3" as (select ST_Union(
ARRAY[st_mapalgebra(rast, 1, 'callback_fn(float[], integer[], text[])'::regprocedure),
st_mapalgebra(rast, 2, 'callback_fn(float[], integer[], text[])'::regprocedure)]
));
It is likely I am way over-complicating the process, but I am fairly new to PostGIS. One of my concerns is that an approach like this will not preserve the underlying geometry of the raster, which is a must. How can I change my queries to be functional in the way I need. Aside from that, am I even approaching this in the right way, or is there a much cleaner way?