From 14dfdbfd6bf4fd509177d4f08fe0287c1c8b7a7e Mon Sep 17 00:00:00 2001 From: Bborie Park Date: Fri, 16 Nov 2012 01:31:55 +0000 Subject: [PATCH] Added variants of ST_Slope, ST_Aspect and ST_Hillshade to provide support for tiles in a coverage. Ticket is #2078 git-svn-id: http://svn.osgeo.org/postgis/trunk@10688 b70326c6-7e19-0410-871a-916f4a2858ee --- NEWS | 6 +- doc/reference_raster.xml | 94 ++++++++++++++------- raster/rt_pg/rtpostgis.sql.in.c | 119 ++++++++++++++++++++++----- raster/rt_pg/rtpostgis_drop.sql.in.c | 2 + 4 files changed, 168 insertions(+), 53 deletions(-) diff --git a/NEWS b/NEWS index ebff2595e..01329e831 100644 --- a/NEWS +++ b/NEWS @@ -73,11 +73,13 @@ PostGIS 2.1.0 of missing TopoGeometry objects from their natural layer. - ST_HillShade(), ST_Aspect() and ST_Slope() have one new optional parameter to interpolate NODATA pixels before running the operation. - - point variant of ST_SetValue(raster) is now a wrapper around geomval + - Point variant of ST_SetValue(raster) is now a wrapper around geomval variant of ST_SetValues(rast). - - proper support for raster band's isnodata flag in core API and loader. + - Proper support for raster band's isnodata flag in core API and loader. - #1655, Additional default values for parameters of ST_Slope - Additional default values for parameters of ST_Aspect and ST_HillShade + - #2078, New variants of ST_Slope, ST_Aspect and ST_HillShade to provide + solution to handling tiles in a coverage * Fixes * diff --git a/doc/reference_raster.xml b/doc/reference_raster.xml index 7c7b5c27e..c1da8b8f9 100644 --- a/doc/reference_raster.xml +++ b/doc/reference_raster.xml @@ -7696,15 +7696,28 @@ FROM dummy_rast; - raster ST_HillShade - raster rast - integer band=1 - text pixeltype=32BF - double precision azimuth=315 - double precision altitude=45 - double precision max_bright=255 - double precision scale=1.0 - boolean interpolate_nodata=FALSE + raster ST_HillShade + raster rast + integer band=1 + text pixeltype=32BF + double precision azimuth=315 + double precision altitude=45 + double precision max_bright=255 + double precision scale=1.0 + boolean interpolate_nodata=FALSE + + + + raster ST_HillShade + raster rast + integer band + raster customextent + text pixeltype=32BF + double precision azimuth=315 + double precision altitude=45 + double precision max_bright=255 + double precision scale=1.0 + boolean interpolate_nodata=FALSE @@ -7747,7 +7760,7 @@ FROM dummy_rast; - Examples + Examples: Variant 1 WITH foo AS ( SELECT ST_SetValues( @@ -7793,12 +7806,22 @@ FROM foo - raster ST_Aspect - raster rast - integer band=1 - text pixeltype=32BF - text measurement=DEGREES - boolean interpolate_nodata=FALSE + raster ST_Aspect + raster rast + integer band=1 + text pixeltype=32BF + text units=DEGREES + boolean interpolate_nodata=FALSE + + + + raster ST_Aspect + raster rast + integer band + raster customextent + text pixeltype=32BF + text units=DEGREES + boolean interpolate_nodata=FALSE @@ -7809,15 +7832,15 @@ FROM foo Returns the aspect (in degrees by default) of an elevation raster band. Utilizes map algebra and applies the aspect equation to neighboring pixels. - measurement indicates the units of the aspect. Possible values are: RADIANS, DEGREES (default). + units indicates the units of the aspect. Possible values are: RADIANS, DEGREES (default). - When measurement = RADIANS, values are between 0 and 2π radians measured clockwise from North. + When units = RADIANS, values are between 0 and 2π radians measured clockwise from North. - When measurement = DEGREES, values are between 0 and 360 degrees measured clockwise from North. + When units = DEGREES, values are between 0 and 360 degrees measured clockwise from North. @@ -7836,7 +7859,7 @@ FROM foo - Examples + Examples: Variant 1 WITH foo AS ( SELECT ST_SetValues( @@ -7883,13 +7906,24 @@ FROM foo - raster ST_Slope - raster rast - integer nband=1 - text pixeltype=32BF - text measurement=DEGREES - double precision scale=1.0 - boolean interpolate_nodata=FALSE + raster ST_Slope + raster rast + integer nband=1 + text pixeltype=32BF + text units=DEGREES + double precision scale=1.0 + boolean interpolate_nodata=FALSE + + + + raster ST_Slope + raster rast + integer nband + raster customextent + text pixeltype=32BF + text units=DEGREES + double precision scale=1.0 + boolean interpolate_nodata=FALSE @@ -7900,7 +7934,7 @@ FROM foo Returns the slope (in degrees by default) of an elevation raster band. Utilizes map algebra and applies the slope equation to neighboring pixels. - measurement indicates the units of the slope. Possible values are: RADIANS, DEGREES (default), PERCENT. + units indicates the units of the slope. Possible values are: RADIANS, DEGREES (default), PERCENT. @@ -7918,13 +7952,13 @@ FROM foo Availability: 2.0.0 - Enhanced: 2.1.0 Uses ST_MapAlgebra() and added optional measurement, scale, interpolate_nodata function parameters + Enhanced: 2.1.0 Uses ST_MapAlgebra() and added optional units, scale, interpolate_nodata function parameters Changed: 2.1.0 In prior versions, return values were in radians. Now, return values default to degrees - Examples + Examples: Variant 1 WITH foo AS ( SELECT ST_SetValues( diff --git a/raster/rt_pg/rtpostgis.sql.in.c b/raster/rt_pg/rtpostgis.sql.in.c index 89156db19..e611e8670 100644 --- a/raster/rt_pg/rtpostgis.sql.in.c +++ b/raster/rt_pg/rtpostgis.sql.in.c @@ -3327,7 +3327,7 @@ CREATE OR REPLACE FUNCTION _st_slope4ma(value double precision[][][], pos intege _pixheight double precision; _width double precision; _height double precision; - _measure text; + _units text; _scale double precision; dz_dx double precision; @@ -3370,7 +3370,7 @@ CREATE OR REPLACE FUNCTION _st_slope4ma(value double precision[][][], pos intege _pixheight := userargs[2]::double precision; _width := userargs[3]::double precision; _height := userargs[4]::double precision; - _measure := userargs[5]; + _units := userargs[5]; _scale := userargs[6]::double precision; /* ArcGIS returns values for edge pixels @@ -3401,7 +3401,7 @@ CREATE OR REPLACE FUNCTION _st_slope4ma(value double precision[][][], pos intege slope := sqrt(dz_dx * dz_dx + dz_dy * dz_dy) / (8 * _scale); -- output depends on user preference - CASE substring(upper(trim(leading from _measure)) for 3) + CASE substring(upper(trim(leading from _units)) for 3) -- percentages WHEN 'PER' THEN slope := 100.0 * slope; @@ -3418,8 +3418,9 @@ CREATE OR REPLACE FUNCTION _st_slope4ma(value double precision[][][], pos intege $$ LANGUAGE 'plpgsql' IMMUTABLE; CREATE OR REPLACE FUNCTION st_slope( - rast raster, nband integer DEFAULT 1, - pixeltype text DEFAULT '32BF', measurement text DEFAULT 'DEGREES', + rast raster, nband integer, + customextent raster, + pixeltype text DEFAULT '32BF', units text DEFAULT 'DEGREES', scale double precision DEFAULT 1.0, interpolate_nodata boolean DEFAULT FALSE ) RETURNS raster @@ -3432,9 +3433,24 @@ CREATE OR REPLACE FUNCTION st_slope( _pixheight double precision; _width integer; _height integer; + _customextent raster; + _extenttype text; BEGIN + _customextent := customextent; + IF _customextent IS NULL THEN + _extenttype := 'FIRST'; + ELSE + _extenttype := 'CUSTOM'; + END IF; + IF interpolate_nodata IS TRUE THEN - _rast := ST_MapAlgebra(ARRAY[ROW($1, $2)]::rastbandarg[], 'st_invdistweight4ma(double precision[][][], integer[][], text[])'::regprocedure, $3, 'FIRST', NULL, 1, 1); + _rast := ST_MapAlgebra( + ARRAY[ROW(rast, nband)]::rastbandarg[], + 'st_invdistweight4ma(double precision[][][], integer[][], text[])'::regprocedure, + pixeltype, + 'FIRST', NULL, + 1, 1 + ); _nband := 1; _pixtype := NULL; ELSE @@ -3452,15 +3468,24 @@ CREATE OR REPLACE FUNCTION st_slope( ARRAY[ROW(_rast, _nband)]::rastbandarg[], '_st_slope4ma(double precision[][][], integer[][], text[])'::regprocedure, _pixtype, - 'FIRST', NULL, + _extenttype, _customextent, 1, 1, _pixwidth::text, _pixheight::text, _width::text, _height::text, - $4::text, $5::text + units::text, scale::text ); END; $$ LANGUAGE 'plpgsql' IMMUTABLE; +CREATE OR REPLACE FUNCTION st_slope( + rast raster, nband integer DEFAULT 1, + pixeltype text DEFAULT '32BF', units text DEFAULT 'DEGREES', + scale double precision DEFAULT 1.0, interpolate_nodata boolean DEFAULT FALSE +) + RETURNS raster + AS $$ SELECT st_slope($1, $2, NULL::raster, $3, $4, $5, $6) $$ + LANGUAGE 'sql' IMMUTABLE; + ----------------------------------------------------------------------- -- ST_Aspect -- http://webhelp.esri.com/arcgisdesktop/9.3/index.cfm?TopicName=How%20Hillshade%20works @@ -3476,7 +3501,7 @@ CREATE OR REPLACE FUNCTION _st_aspect4ma(value double precision[][][], pos integ _width double precision; _height double precision; - _measure text; + _units text; dz_dx double precision; dz_dy double precision; @@ -3515,7 +3540,7 @@ CREATE OR REPLACE FUNCTION _st_aspect4ma(value double precision[][][], pos integ _width := userargs[1]::double precision; _height := userargs[2]::double precision; - _measure := userargs[3]; + _units := userargs[3]; /* ArcGIS returns values for edge pixels -- check that pixel is not edge pixel @@ -3563,7 +3588,7 @@ CREATE OR REPLACE FUNCTION _st_aspect4ma(value double precision[][][], pos integ END IF; -- output depends on user preference - CASE substring(upper(trim(leading from _measure)) for 3) + CASE substring(upper(trim(leading from _units)) for 3) -- radians WHEN 'rad' THEN RETURN aspect; @@ -3576,8 +3601,9 @@ CREATE OR REPLACE FUNCTION _st_aspect4ma(value double precision[][][], pos integ $$ LANGUAGE 'plpgsql' IMMUTABLE; CREATE OR REPLACE FUNCTION st_aspect( - rast raster, nband integer DEFAULT 1, - pixeltype text DEFAULT '32BF', measurement text DEFAULT 'DEGREES', + rast raster, nband integer, + customextent raster, + pixeltype text DEFAULT '32BF', units text DEFAULT 'DEGREES', interpolate_nodata boolean DEFAULT FALSE ) RETURNS raster @@ -3588,9 +3614,24 @@ CREATE OR REPLACE FUNCTION st_aspect( _pixtype text; _width integer; _height integer; + _customextent raster; + _extenttype text; BEGIN + _customextent := customextent; + IF _customextent IS NULL THEN + _extenttype := 'FIRST'; + ELSE + _extenttype := 'CUSTOM'; + END IF; + IF interpolate_nodata IS TRUE THEN - _rast := ST_MapAlgebra(ARRAY[ROW($1, $2)]::rastbandarg[], 'st_invdistweight4ma(double precision[][][], integer[][], text[])'::regprocedure, $3, 'FIRST', NULL, 1, 1); + _rast := ST_MapAlgebra( + ARRAY[ROW(rast, nband)]::rastbandarg[], + 'st_invdistweight4ma(double precision[][][], integer[][], text[])'::regprocedure, + pixeltype, + 'FIRST', NULL, + 1, 1 + ); _nband := 1; _pixtype := NULL; ELSE @@ -3606,14 +3647,23 @@ CREATE OR REPLACE FUNCTION st_aspect( ARRAY[ROW(_rast, _nband)]::rastbandarg[], '_st_aspect4ma(double precision[][][], integer[][], text[])'::regprocedure, _pixtype, - 'FIRST', NULL, + _extenttype, _customextent, 1, 1, _width::text, _height::text, - $4::text + units::text ); END; $$ LANGUAGE 'plpgsql' IMMUTABLE; +CREATE OR REPLACE FUNCTION st_aspect( + rast raster, nband integer DEFAULT 1, + pixeltype text DEFAULT '32BF', units text DEFAULT 'DEGREES', + interpolate_nodata boolean DEFAULT FALSE +) + RETURNS raster + AS $$ SELECT st_aspect($1, $2, NULL::raster, $3, $4, $5) $$ + LANGUAGE 'sql' IMMUTABLE; + ----------------------------------------------------------------------- -- ST_HillShade -- http://webhelp.esri.com/arcgisdesktop/9.3/index.cfm?TopicName=How%20Hillshade%20works @@ -3750,7 +3800,8 @@ CREATE OR REPLACE FUNCTION _st_hillshade4ma(value double precision[][][], pos in $$ LANGUAGE 'plpgsql' IMMUTABLE; CREATE OR REPLACE FUNCTION st_hillshade( - rast raster, nband integer DEFAULT 1, + rast raster, nband integer, + customextent raster, pixeltype text DEFAULT '32BF', azimuth double precision DEFAULT 315.0, altitude double precision DEFAULT 45.0, max_bright double precision DEFAULT 255.0, scale double precision DEFAULT 1.0, @@ -3766,9 +3817,24 @@ CREATE OR REPLACE FUNCTION st_hillshade( _pixheight double precision; _width integer; _height integer; + _customextent raster; + _extenttype text; BEGIN + _customextent := customextent; + IF _customextent IS NULL THEN + _extenttype := 'FIRST'; + ELSE + _extenttype := 'CUSTOM'; + END IF; + IF interpolate_nodata IS TRUE THEN - _rast := ST_MapAlgebra(ARRAY[ROW($1, $2)]::rastbandarg[], 'st_invdistweight4ma(double precision[][][], integer[][], text[])'::regprocedure, $3, 'FIRST', NULL, 1, 1); + _rast := ST_MapAlgebra( + ARRAY[ROW(rast, nband)]::rastbandarg[], + 'st_invdistweight4ma(double precision[][][], integer[][], text[])'::regprocedure, + pixeltype, + 'FIRST', NULL, + 1, 1 + ); _nband := 1; _pixtype := NULL; ELSE @@ -3786,16 +3852,27 @@ CREATE OR REPLACE FUNCTION st_hillshade( ARRAY[ROW(_rast, _nband)]::rastbandarg[], '_st_hillshade4ma(double precision[][][], integer[][], text[])'::regprocedure, _pixtype, - 'FIRST', NULL, + _extenttype, _customextent, 1, 1, _pixwidth::text, _pixheight::text, _width::text, _height::text, - $4::text, $5::text, - $6::text, $7::text + $5::text, $6::text, + $7::text, $8::text ); END; $$ LANGUAGE 'plpgsql' IMMUTABLE; +CREATE OR REPLACE FUNCTION st_hillshade( + rast raster, nband integer DEFAULT 1, + pixeltype text DEFAULT '32BF', + azimuth double precision DEFAULT 315.0, altitude double precision DEFAULT 45.0, + max_bright double precision DEFAULT 255.0, scale double precision DEFAULT 1.0, + interpolate_nodata boolean DEFAULT FALSE +) + RETURNS RASTER + AS $$ SELECT st_hillshade($1, $2, NULL::raster, $3, $4, $5, $6, $7, $8) $$ + LANGUAGE 'sql' IMMUTABLE; + ----------------------------------------------------------------------- -- Get information about the raster ----------------------------------------------------------------------- diff --git a/raster/rt_pg/rtpostgis_drop.sql.in.c b/raster/rt_pg/rtpostgis_drop.sql.in.c index 0a0373cd1..ef8d56d7d 100644 --- a/raster/rt_pg/rtpostgis_drop.sql.in.c +++ b/raster/rt_pg/rtpostgis_drop.sql.in.c @@ -433,8 +433,10 @@ DROP FUNCTION IF EXISTS _st_contains(geometry, raster, integer); DROP FUNCTION IF EXISTS st_addband(raster, raster[], integer); -- function signatures changed +DROP FUNCTION IF EXISTS st_slope(raster, integer, text, text, double precision, boolean); DROP FUNCTION IF EXISTS st_slope(raster, integer, text, boolean); DROP FUNCTION IF EXISTS st_slope(raster, integer, text); +DROP FUNCTION IF EXISTS st_aspect(raster, integer, text, text, boolean); DROP FUNCTION IF EXISTS st_aspect(raster, integer, text, boolean); DROP FUNCTION IF EXISTS st_aspect(raster, integer, text); DROP FUNCTION IF EXISTS st_hillshade(raster, integer, text, double precision, double precision, double precision, double precision, boolean); -- 2.50.1