Added ST_TileAsGeom() returning only the extent of the planned tiles as polygons

git-svn-id: http://svn.osgeo.org/postgis/trunk@9069 b70326c6-7e19-0410-871a-916f4a2858ee

diff --git a/raster/scripts/plpgsql/st_tile.sql b/raster/scripts/plpgsql/st_tile.sql
index ee96a9ce90f91f94df3dd6eeed40e615021d950c..d7eee98ed63f63b5fb38566e75487dd1c901156d 100644 (file)
--- a/raster/scripts/plpgsql/st_tile.sql
+++ b/raster/scripts/plpgsql/st_tile.sql
@@ -1,4 +1,4 @@
-----------------------------------------------------------------------\r
+----------------------------------------------------------------------\r
 --\r
 -- $Id: st_tile.sql 8255 2011-11-29 16:34:48Z pracine $\r
 --\r
@@ -6,7 +6,11 @@
 --\r
 ----------------------------------------------------------------------\r
 -- ST_Tile\r
--- Split a raster into a set of raster tiles, one tile per row returned.\r
+-- Split a raster into a set of raster tiles, one tile per row returned. \r
+-- Works on multiband rasters. There is no way to specify the upper left \r
+-- corner of the new tiled grid. The grid start at the upperleft corner \r
+-- of the provided raster.\r
+--\r
 -- rast   - Raster to be tiled.\r
 -- width  - Width of the tiles.\r
 -- height - Height of the tiles\r
@@ -17,7 +21,7 @@
 -- nodatavalue   - nodata value to use to pad the outbound tiles when the provided \r
 --                 raster do not have a nodata value defined. If not provided and \r
 --                 the raster do not have a nodata value defined \r
---                 ST_MinPossibleValue(ST_BandPixelType(rast, band)) is used.\r
+--                 ST_MinPossibleValue(ST_BandPixelType(rast, band)) is used for each band.\r
 --\r
 -- Example producing 120 x 120 pixel tiles\r
 --\r
@@ -25,7 +29,7 @@
 -- SELECT ST_Tile(rast, 120, 120, true), generate_series(1, 3600) rid\r
 -- FROM srtm_22_03;\r
 ----------------------------------------------------------------------------------------------------------------------\r
-DROP FUNCTION ST_Tile(rast raster, width integer, height integer, padwithnodata boolean, nodatavalue double precision);\r
+DROP FUNCTION IF EXISTS ST_Tile(rast raster, width integer, height integer, padwithnodata boolean, nodatavalue double precision);\r
 CREATE OR REPLACE FUNCTION ST_Tile(rast raster, width integer, height integer, padwithnodata boolean DEFAULT FALSE, nodatavalue double precision DEFAULT NULL) \r
     RETURNS SETOF raster AS \r
     $$\r
@@ -38,7 +42,7 @@ CREATE OR REPLACE FUNCTION ST_Tile(rast raster, width integer, height integer, p
         newnodata double precision;\r
         newpixtype text;\r
         nbband integer;\r
-        bnadi integer;\r
+        bandi integer;\r
         newrast raster;\r
         initvalue double precision;\r
         grid record;\r
@@ -98,6 +102,75 @@ CREATE OR REPLACE FUNCTION ST_Tile(rast raster, width integer, height integer, p
     END;\r
     $$\r
     LANGUAGE 'plpgsql';\r
+    \r
+----------------------------------------------------------------------\r
+-- ST_TileAsGeom\r
+-- Split a raster into a set of raster tiles, returning only the geometry \r
+-- corresponding to each tile. \r
+-- There is no way to specify the upper left corner of the new tiled grid. \r
+-- The grid start at the upperleft corner of the provided raster.\r
+--\r
+-- rast   - Raster to be tiled.\r
+-- width  - Width of the tiles.\r
+-- height - Height of the tiles\r
+-- padwithnodata - If TRUE, the produced tiles are strictly width x heigth pixels. \r
+--                 Pixels outside the extent of the passed raster are filled with \r
+--                 nodata value. When FALSE out of bound tiles are clipped to the \r
+--                 extent of the raster. Default to FALSE.\r
+--\r
+-- Example producing 120 x 120 pixel tiles\r
+--\r
+-- SELECT ST_TileAsGeom(rast, 130, 130, true)\r
+-- FROM srtm_22_03;\r
+----------------------------------------------------------------------------------------------------------------------\r
+DROP FUNCTION IF EXISTS ST_TileAsGeom(rast raster, width integer, height integer, padwithnodata boolean);\r
+CREATE OR REPLACE FUNCTION ST_TileAsGeom(rast raster, width integer, height integer, padwithnodata boolean DEFAULT FALSE) \r
+    RETURNS SETOF geometry AS \r
+    $$\r
+    DECLARE\r
+        gridrast raster;\r
+        rastwidth integer;\r
+        rastheight integer;\r
+        gridwidth integer;\r
+        gridheight integer;\r
+        nbband integer;\r
+        rastextent geometry;\r
+    BEGIN\r
+        IF rast IS NULL THEN\r
+            RETURN;\r
+        END IF;\r
+    \r
+        nbband := ST_Numbands(rast);\r
+        IF nbband < 1 THEN\r
+                RAISE NOTICE 'Raster do not have band %. Returning null', band;\r
+                RETURN;\r
+        END IF;\r
+\r
+        rastwidth := ST_Width(rast);\r
+        IF width IS NULL THEN\r
+            width := rastwidth;\r
+        END IF;\r
+\r
+        rastheight := ST_Height(rast);\r
+        IF height IS NULL THEN\r
+            height := rastheight;\r
+        END IF;\r
+\r
+        gridwidth := (rastwidth / width) + CASE WHEN rastwidth % width > 0 THEN 1 ELSE 0 END;\r
+        gridheight := (rastheight / height) + CASE WHEN rastheight % height > 0 THEN 1 ELSE 0 END;\r
+\r
+        gridrast := ST_AddBand(ST_MakeEmptyRaster(gridwidth, gridheight, ST_UpperLeftX(rast), ST_UpperLeftY(rast), ST_ScaleX(rast) * width, ST_ScaleY(rast) * height, ST_SkewX(rast), ST_SkewY(rast), ST_SRID(rast)), '8BUI'::text, 1, 0);\r
+        IF padwithnodata THEN\r
+            RETURN QUERY SELECT (ST_PixelAsPolygons(gridrast)).geom geom;\r
+        ELSE\r
+            rastextent := rast::geometry;\r
+            RETURN QUERY SELECT ST_Intersection(rastextent, (ST_PixelAsPolygons(gridrast)).geom) geom;\r
+        END IF;\r
+        RETURN;\r
+    END;\r
+    $$\r
+    LANGUAGE 'plpgsql';\r
+\r
 \r
 -- Redefine ST_TestRaster()\r
 CREATE OR REPLACE FUNCTION ST_TestRaster(ulx float8, uly float8, val float8) \r
@@ -118,10 +191,6 @@ FROM (SELECT ST_PixelAsPolygons(ST_TestRaster(0, 0, 1)) gv) foo;
 -- Tile it to 10x10 tiles\r
 SELECT ST_Tile(ST_TestRaster(0, 0, 1), 10, 10)\r
 \r
--- Old test for the extent of the tile grid\r
-SELECT ST_AsBinary((gv).geom) geom,  (gv).val\r
-FROM (SELECT ST_PixelAsPolygons(ST_Rescale(ST_TestRaster(0, 0, 1), 0.003, 0.004)) gv) foo;\r
-\r
 -- Display the result of the tile function\r
 SELECT ST_AsBinary((gv).geom) geom,  (gv).val\r
 FROM (SELECT ST_PixelAsPolygons(ST_Tile(ST_TestRaster(0, 0, 1), 3, 4)) gv) foo;\r
@@ -146,6 +215,10 @@ SELECT ST_AsBinary((gv).geom) geom,  (gv).val, rid
 FROM (SELECT ST_PixelAsPolygons(rast) gv, rid\r
       FROM (SELECT ST_Tile(ST_TestRaster(0, 0, 1), 10, 10, true) rast, generate_series(1, 35) rid) foo WHERE rid = 2) foo2;\r
 \r
+-- Test ST_TileAsGeom\r
+SELECT ST_TileAsGeom(ST_TestRaster(0, 0, 1), 10, 10);\r
+\r
+\r
 -- Other tests\r
 SELECT ST_Tile(ST_AddBand(ST_MakeEmptyRaster(48, 63, 0, 0, 0.001, -0.001, 0, 0, 4269), '8BSI'::text, 0, -2), 10, 10, true, -10);\r
 \r