[postgis-users] Efficient ways to query PostGIS Raster pixel values using Intersects

[Bborie Park]

Andrew,

Have you considered using 2-raster MapAlgebra instead?

I haven't tested this but for your example...

{{{
WITH foo AS (
	SELECT
		ST_AsRaster(
			GEOMETRYFROMTEXT('MULTIPOLYGON(((10 10, 10 0, 0 0, 0 10, 10 10)))',4326),
			rast
		) AS rast
	FROM marso
	LIMIT 1
)
SELECT
	ST_MapAlgebra(m.rast, f.rast, 'rast1')
FROM marso m
JOIN foo f
	ON ST_Intersects(m.rast, f.rast);
}}}

That'll return a set of rasters who intersect the provided geometry. 
You can then use the aggregate function ST_Union (I haven't played with 
it yet) to bring the rasters together into one raster, which can be 
passed off to ST_AsPNG.

I apologize in advance if this takes you off the beaten path...

-bborie

On 12/08/2011 09:25 AM, Andrew Hill wrote:
> Efficient ways to query PostGIS Raster pixel values using Intersects
>
> Hi,
>
> I've been working on some projects where we would like to find the values
> for pixels in a defined area. From what I can find, these sorts of queries
> have a couple of options. First, finding all tiles that intersect your
> polygon, then converting all pixels to geoms and performing a second
> ST_Intersects. Another (and from my attempts the faster of the two), is to
> perform the first intersects and then doing an XY walk of the pixels. I
> have an example here,
>
> https://gist.github.com/1353427
>
> My problem is that we still aren't finding good methods for doing these
> queries fast. It gets more complex. Take for example,
>
> SELECT rid,ST_AsText(ST_Centroid(rast)) FROM marso WHERE
> ST_Intersects(rast, GEOMETRYFROMTEXT('MULTIPOLYGON(((10 10, 10 0, 0 0, 0
> 10, 10 10)))',4326))
>
> (marso is an indexed 3600x1800 4326 global raster tiled at 100x100)
>
> As expected, I get back a row for every tile my polygon intersects. Now, if
> I want to use those results for display, I can swap in an ST_AsPNG. Here is
> where it gets silly, I still get my 4 tiles, all fully converted to PNG
> (not cropped). So, I have to clip them. Right now, and from the ST_Clip I
> have seen, this isn't going to be very efficient. It will require either XY
> walks or complete pixel->geom conversions to resample the tile. Also, if I
> want to use the result for display, the mechanisms for downsampling are
> nascent. Downsampling a tile for example in this scenario would take both
> an ST_Resample and an ST_Clip, and from what I've done, resulting in a
> really inefficient query.
>
> Here is an example of the SQL to do so and then ST_Clip function used,
>
> sql: https://gist.github.com/1447684
> clip: https://gist.github.com/1447686
>
> So. I'm just curious what sorts of optimizations either (a) I'm missing for
> indexing, reshaping, and returning pixel values and pixels-in-area etc., or
> (b) are in the pipeline? What can we expect from the raster libraries in
> the future? (I'll save questions related to the merging of result tiles for
> a later thread). For web development, speed is directly tied to user
> experience, so it is really important to our plans for using the Raster
> libraries in future work and we'd love to help find some solutions.
>
>
> cheers,
>
> Andrew
>
>
>
>
>
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users

-- 
Bborie Park
Programmer
Center for Vectorborne Diseases
UC Davis
530-752-8380
bkpark at ucdavis.edu