wiki:UsersWikiwktrasterfunctions

Version 8 (modified by robe, 15 years ago) ( diff )

--

WKT Raster User contributed Functions

  • This section contains additional functions not currently in WKT Raster. Some of these may become obsolete as WKT Raster becomes more mature
-- Example use case: 
/** SELECT p.gis, p.geom, upgis_ptval(r.rast, 1, p.geom) FROM 
sometableofpts AS p INNER JOIN raster_table As r 
  ON ST_Intersects(p.geom, ST_Envelope(rast));
***/
CREATE OR REPLACE FUNCTION upgis_ptval(param_rast raster, param_bnum integer, param_pt geometry)
 RETURNS float AS
  $$
  DECLARE 
    var_result float = 0; var_rid integer;
    var_cols integer; var_rows integer;var_c integer;  var_r integer;
    var_upx float; var_upy float; var_sizex float; var_sizey float; 
  BEGIN
  -- 1 does point intersects chosen tile
        SELECT ST_Width(param_rast), ST_Height(param_rast), ST_UpperLeftX(param_rast), ST_UpperLeftY(param_rast), ST_PixelSizeX(param_rast), ST_PixelSizeY(param_rast)
      INTO var_cols, var_rows, var_upx, var_upy, var_sizex, var_sizey
    WHERE ST_Intersects(ST_Envelope(param_rast), param_pt) 
    LIMIT 1;
    
   -- 2 if no return null
    IF var_cols IS NULL THEN -- no intersection
      var_result := NULL;
    ELSE
   --  3 find row and column
    SELECT 1 + floor( ( (ST_X(param_pt) - var_upx) / (var_sizex*var_cols) )*var_cols )::integer, 1 + floor( ( (ST_Y(param_pt) - var_upy)/(var_sizey*var_rows) )*var_rows )::integer
      INTO var_c, var_r;
      
   -- 4 edge case
      IF var_c = 0 THEN var_c := 1; END IF;
      IF var_r = 0 THEN var_r := 1; END IF; 
      IF var_c = var_cols + 1 THEN var_c := var_cols; END IF;
      IF var_r = var_rows + 1 THEN var_r := var_rows; END IF;
   -- 5 get cell value
      SELECT ST_Value(param_rast,param_bnum, var_c , var_r )
        INTO var_result;
    END IF;
    --RAISE NOTICE 'col: %, row: %, val: %, pt: %', var_c, var_r, var_result, ST_AsText(param_pt);
    RETURN var_result;
  END
 $$
LANGUAGE 'plpgsql' IMMUTABLE;
COMMENT ON FUNCTION upgis_ptval(raster, integer, geometry) IS 'This function takes as input raster, band_num, point geometry and returns the band pixel value of the cell the point falls in.  It assumes the spatial reference system of the geometry and raster are the same.';
Note: See TracWiki for help on using the wiki.