pggeodb.nancy.inra.fr/db_robot - db_robot on pggeodb.nancy.inra.fr

Function: st_invdistweight4ma

Schema

public

Owner

albenard

Descriptions

There is no description for function st_invdistweight4ma

Options

Definition

CREATE OR REPLACE FUNCTION public.st_invdistweight4ma (
 value double precision [],
 pos integer [],
 variadic userargs text [] = NULL::text[]
)
RETURNS double precision AS
$span$
DECLARE
_value double precision[][][];
ndims int;

k double precision DEFAULT 1.;
_k double precision DEFAULT 1.;
z double precision[];
d double precision[];
_d double precision;
z0 double precision;

_z integer;
x integer;
y integer;

cx integer;
cy integer;
cv double precision;
cw double precision DEFAULT NULL;

w integer;
h integer;
max_dx double precision;
max_dy double precision;
BEGIN
-- RAISE NOTICE 'value = %', value;
-- RAISE NOTICE 'userargs = %', userargs;

ndims := array_ndims(value);
-- add a third dimension if 2-dimension
IF ndims = 2 THEN
_value := public._ST_convertarray4ma(value);
ELSEIF ndims != 3 THEN
RAISE EXCEPTION 'First parameter of function must be a 3-dimension array';
ELSE
_value := value;
END IF;

-- only use the first raster passed to this function
IF array_length(_value, 1) > 1 THEN
RAISE NOTICE 'Only using the values from the first raster';
END IF;
_z := array_lower(_value, 1);

-- width and height (0-based)
h := array_upper(_value, 2) - array_lower(_value, 2);
w := array_upper(_value, 3) - array_lower(_value, 3);

-- max distance from center pixel
max_dx := w / 2;
max_dy := h / 2;
-- RAISE NOTICE 'max_dx, max_dy = %, %', max_dx, max_dy;

-- correct width and height (1-based)
w := w + 1;
h := h + 1;
-- RAISE NOTICE 'w, h = %, %', w, h;

-- width and height should be odd numbers
IF w % 2. != 1 THEN
RAISE EXCEPTION 'Width of neighborhood array does not permit for a center pixel';
END IF;
IF h % 2. != 1 THEN
RAISE EXCEPTION 'Height of neighborhood array does not permit for a center pixel';
END IF;

-- center pixel's coordinates
cy := max_dy + array_lower(_value, 2);
cx := max_dx + array_lower(_value, 3);
-- RAISE NOTICE 'cx, cy = %, %', cx, cy;

-- if userargs provided, only use the first two args
IF userargs IS NOT NULL AND array_ndims(userargs) = 1 THEN
-- first arg is power factor
k := userargs[array_lower(userargs, 1)]::double precision;
IF k IS NULL THEN
k := _k;
ELSEIF k < 0. THEN
RAISE NOTICE 'Power factor (< 0) must be between 0 and 1.  Defaulting to 0';
k := 0.;
ELSEIF k > 1. THEN
RAISE NOTICE 'Power factor (> 1) must be between 0 and 1.  Defaulting to 1';
k := 1.;
END IF;

-- second arg is what to do if center pixel has a value
-- this will be a weight to apply for the center pixel
IF array_length(userargs, 1) > 1 THEN
cw := abs(userargs[array_lower(userargs, 1) + 1]::double precision);
IF cw IS NOT NULL THEN
IF cw < 0. THEN
RAISE NOTICE 'Weight (< 0) of center pixel value must be between 0 and 1.  Defaulting to 0';
cw := 0.;
ELSEIF cw > 1 THEN
RAISE NOTICE 'Weight (> 1) of center pixel value must be between 0 and 1.  Defaulting to 1';
cw := 1.;
END IF;
END IF;
END IF;
END IF;
-- RAISE NOTICE 'k = %', k;
k = abs(k) * -1;

-- center pixel value
cv := _value[_z][cy][cx];

-- check to see if center pixel has value
-- RAISE NOTICE 'cw = %', cw;
IF cw IS NULL AND cv IS NOT NULL THEN
RETURN cv;
END IF;

FOR y IN array_lower(_value, 2)..array_upper(_value, 2) LOOP
FOR x IN array_lower(_value, 3)..array_upper(_value, 3) LOOP
-- RAISE NOTICE 'value[%][%][%] = %', _z, y, x, _value[_z][y][x];

-- skip NODATA values and center pixel
IF _value[_z][y][x] IS NULL OR (x = cx AND y = cy) THEN
CONTINUE;
END IF;

z := z || _value[_z][y][x];

-- use pythagorean theorem
_d := sqrt(power(cx - x, 2) + power(cy - y, 2));
-- RAISE NOTICE 'distance = %', _d;

d := d || _d;
END LOOP;
END LOOP;
-- RAISE NOTICE 'z = %', z;
-- RAISE NOTICE 'd = %', d;

-- neighborhood is NODATA
IF z IS NULL OR array_length(z, 1) < 1 THEN
-- center pixel has value
IF cv IS NOT NULL THEN
RETURN cv;
ELSE
RETURN NULL;
END IF;
END IF;

z0 := 0;
_d := 0;
FOR x IN array_lower(z, 1)..array_upper(z, 1) LOOP
d[x] := power(d[x], k);
z[x] := z[x] * d[x];
_d := _d + d[x];
z0 := z0 + z[x];
END LOOP;
z0 := z0 / _d;
-- RAISE NOTICE 'z0 = %', z0;

-- apply weight for center pixel if center pixel has value
IF cv IS NOT NULL THEN
z0 := (cw * cv) + ((1 - cw) * z0);
-- RAISE NOTICE '*z0 = %', z0;
END IF;

RETURN z0;
END;
$span$
LANGUAGE 'plpgsql'
IMMUTABLE
CALLED ON NULL INPUT
SECURITY INVOKER
COST 100;

This file was generated with SQL Manager for PostgreSQL (www.pgsqlmanager.com) at 26/02/2014 11:51