Tuesday, 3 October 2017

PostGIS: Faster way to ST_Clip() rasters with a polygon?


I want to clip an image with a polygon from postGIS rasters.
What I do:
- Import raster into postGIS(raster2pgsql, tile-size: 500x500, each tile contains 12 bands landsat 8)

- Import a shapefile into postGIS with SRID = 4326.
I have 21309 tiles in my raster table, and 60 pylygons in my shape table. Now, I want to clip a polygon:


SELECT oid,
lowrite(lo_open(oid, 131072), png) AS num_bytes
FROM (
VALUES (lo_create(0),
ST_AsPNG(
(SELECT ST_Clip(ST_Transform(
(
SELECT ST_Union(foo.rast)

FROM
(SELECT rast
FROM ldcm15, vnm_adm2
WHERE gid=18
AND ST_Intersects(rast_geom_4326, geom)) AS foo),4326), vnm_adm2.geom,true)
FROM vnm_adm2
WHERE gid = 18)))) AS v(oid,png);

rast_geom_4326 is where I stored rast's polygon. Because, I don't want to use ST_Transform(rast, 4326) in the query, it's pain to wait for re-projection.
However, it still takes a long time to get the result.



I'm using PostGIS 2.1 on PostgresSQL 9.3.



Answer



It's generally faster to clip first and then union the clips.


Try changing to


(SELECT ST_Union(ST_Clip(ldcm15.rast, ST_Transform(geom, ST_SRID(ldcm15.rast) ) ) )
FROM ldcm15, vnm_adm2
WHERE gid=18 AND ST_Intersects(rast_geom_4326, geom) )

No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...