postgis - ST_MapAlgebra and combination of four rasters

I am trying to create combination of four rasters, but I am lost with ST_MapAlgebra callback function. Official documentation (http://postgis.net/docs/RT_ST_MapAlgebra.html) gives an example, that return only 0 but does not read any values from input.

How should I create, for simplicity, callback function that will sum together four input rasters?

Answer

You may not need ST_MapAlgebra at all. People often reach for ST_MapAlgebra when ST_Union often does the same thing but faster.

If all you want to do is sum the values of pixels to get a new raster where each pixel is the sum of the pixel values of intersected pixels, then ST_Union is what you want. http://postgis.net/docs/manual-2.3/RT_ST_Union.html

Instead of taking the default behavior of LAST, you'd specify SUM for the behavior.

Just to add to my answer, a map algebra callback function returns one pixel value. The reason is that it is only computing the new pixel value of the center pixel of the neighborhood - the pixel value of the new raster for that is located in the center pixel.

So the ST_MapAlgebra gives the callback a neighborhood of pixels of distance X (so if you have distance of 1 and one raster band that means you have a 3x3x1 matrix (because 1 to right and 1 left of the center pixel, and since you have just one band last dimension is one) , whatever your callback function does becomes the new value of the center pixel.

If you have more than one raster, then your mapalgebra call back matrix becomes something like 3x3x2, but it still returns just one value. Because whatever operation you do on an overlay of 2 matrices results in an answer for the center pixel (of a new raster).

So if you just want to sum pixels across rasters for example, you would use a distance = 0. Each matrix would then have exactly 1 element, but you would have as many as you have rasters going in (e.g. 1x1x4 if you have 4 rasters)