Hi, I have been using PostGIS to build and query a DEM, while testing I noticed that I was not getting the expected results from st_value using resample => 'nearest'. Instead of returning the value of the nearest pixel the x-coordinate is always simply rounded down, while the y-coordinate ist always rounded up. The following query showcases this behaviour:
WITH r AS ( SELECT ST_SetValues( ST_AddBand( ST_MakeEmptyRaster( width => 2, height => 2, upperleftx => 0, upperlefty => 2, scalex => 1.0, scaley => -1.0, skewx => 0, skewy => 0, srid => 4326), index => 1, pixeltype => '16BSI', initialvalue => 0, nodataval => -999), 1,1,1, newvalueset =>ARRAY[ARRAY[10.0::float8, 50.0::float8], ARRAY[40.0::float8, 20.0::float8]]) AS rast ) SELECT 'Test 5', round(ST_Value(rast, 1, 'SRID=4326;POINT(0 2)'::geometry, resample => 'nearest')) as top_left, round(ST_Value(rast, 1, 'SRID=4326;POINT(1 2)'::geometry, resample => 'nearest')) as top_right, round(ST_Value(rast, 1, 'SRID=4326;POINT(0 1)'::geometry, resample => 'nearest')) as bottom_left, round(ST_Value(rast, 1, 'SRID=4326;POINT(1 1)'::geometry, resample => 'nearest')) as bottom_right, round(ST_Value(rast, 1, 'SRID=4326;POINT(0.9 1)'::geometry, resample => 'nearest')) as should_be_bottom_right_but_is_bottom_left, round(ST_Value(rast, 1, 'SRID=4326;POINT(1 1.1)'::geometry, resample => 'nearest')) as should_be_bottom_right_but_is_top_right, round(ST_Value(rast, 1, 'SRID=4326;POINT(0.9 1.1)'::geometry, resample => 'nearest')) as should_be_bottom_right_but_is_top_left FROM r; Resulting in: ?column? | top_left | top_right | bottom_left | bottom_right | should_be_bottom_right_but_is_bottom_left | should_be_bottom_right_but_is_top_right | should_be_bottom_right_but_is_top_left ----------+----------+-----------+-------------+--------------+-------------------------------------------+-----------------------------------------+---------------------------------------- Test 5 | 10 | 50 | 40 | 20 | 40 | 50 | 10 When looking through the source code on git I found that the method 'rt_band_get_pixel_resample' that is eventually called when using st_value simply rounds down the given raster coordinates instead of actually using nearest-neighbor: else if (resample == RT_NEAREST) { return rt_band_get_pixel( band, floor(xr), floor(yr), r_value, r_nodata ); } So my question is, is this the intended behaviour? Is there maybe another way to get a value from a raster using actual nearest-neighbor resampling? I'm using Postgres 17.2 and postgis 3.6.0dev if that is relevant. Thanks in advance
smime.p7s
Description: S/MIME cryptographic signature