如何使用postgis st_setvalue更新磁贴中的所有单元格？

Question

我在支持postgis的数据库中有两个小栅格:一个DEM和一个streams/ponds图层。两者都存储为10x10瓦片，并且具有相同的srid和单元大小。我想把小溪/池塘“烧”到DEM (即栅格计算器)中。我尝试了ST_MapAlgebra，但不知道如何忽略流/池塘层中的nodata单元。现在，我正在尝试使用ST_SetValue简单地更新DEM。下面的代码似乎只更新了每个平铺的一个单元格(见下图)。如果不包含st_intersects()子句，它将永远运行(即使我只传入一个点来更新单个单元格)。我认为如果栅格没有存储为瓦片，它可能会起作用，但我不希望这样做。

update dem_clip d set rast = st_setvalue(d.rast, geom, st_value(d.rast,geom) - 1000)
from (select ((st_pixelascentroids(rast)).geom) from stream_clip) s
where st_intersects(d.rast,s.geom);

what I want (streams burned into DEM)

what I get

Answer 1

这个函数解决了我的问题。虽然很慢，但我并不着急。我相信它是可以改进的。

CREATE OR REPLACE FUNCTION update_vals(g1 geometry, v1 integer)
 RETURNS void as $$
begin 
update dem_clip_test set rast = 
st_setvalue(rast, g1, 
st_value(rast,g1) + v1)
where st_intersects(rast,g1);
end;
$$ LANGUAGE plpgsql;

select update_vals(geom, -1000) from 
(select ((st_pixelascentroids(rast)).geom) from stream_clip) s;