- GEOMETRYCOLLECTION support for ST_MakeValid (Sandro Santilli / Vizzuality)
- ST_PixelOfValue (Bborie Park / UC Davis)
- Casts to/from PostgreSQL geotypes (point/path/polygon).
+ - #2030, n-raster (and n-band) ST_MapAlgebra (Bborie Park / UC Davis)
* Enhancements *
- #823, tiger geocoder: Make loader_generate_script download portion less greedy
</refsection>
</refentry>
+ <refentry id="rastbandarg">
+ <refnamediv>
+ <refname>rastbandarg</refname>
+ <refpurpose>A composite type for use when needing to express a raster and a band index of that raster.</refpurpose>
+ </refnamediv>
+
+ <refsection>
+ <title>Description</title>
+ <para>
+ A composite type for use when needing to express a raster and a band index of that raster.
+
+ <variablelist>
+
+ <varlistentry>
+ <term>
+ <parameter>rast </parameter>
+ <type>raster</type>
+ </term>
+ <listitem>
+ <para>
+ The raster in question/
+ </para>
+ </listitem>
+ </varlistentry>
+
+ <varlistentry>
+ <term>
+ <parameter>nband </parameter>
+ <type>integer</type>
+ </term>
+ <listitem>
+ <para>
+ 1-based value indicating the band of raster
+ </para>
+ </listitem>
+ </varlistentry>
+
+ </variablelist>
+
+ </para>
+ </refsection>
+
+ <refsection>
+ <title>See Also</title>
+ <para>
+ <xref linkend="RT_ST_MapAlgebra" />
+ </para>
+ </refsection>
+ </refentry>
+
<refentry id="raster">
<refnamediv>
<refname>raster</refname>
</refsection>
</refentry>
+ <refentry id="RT_ST_MapAlgebra">
+ <refnamediv>
+ <refname>ST_MapAlgebra</refname>
+ <refpurpose>
+ Returns a one-band raster given one or more input rasters, band indexes and one user-specified callback function.
+ </refpurpose>
+ </refnamediv>
+
+ <refsynopsisdiv>
+ <funcsynopsis>
+
+ <funcprototype>
+ <funcdef>raster <function>ST_MapAlgebra</function></funcdef>
+ <paramdef><type>rastbandarg[] </type> <parameter>rastbandargset</parameter></paramdef>
+ <paramdef><type>regprocedure </type> <parameter>callbackfunc</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>pixeltype=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>extenttype=INTERSECTION</parameter></paramdef>
+ <paramdef choice="opt"><type>raster </type> <parameter>customextent=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancex=0</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancey=0</parameter></paramdef>
+ <paramdef choice="opt"><type>text[] </type> <parameter>VARIADIC userargs=NULL</parameter></paramdef>
+ </funcprototype>
+
+ <funcprototype>
+ <funcdef>raster <function>ST_MapAlgebra</function></funcdef>
+ <paramdef><type>raster </type> <parameter>rast</parameter></paramdef>
+ <paramdef><type>integer[] </type> <parameter>nband</parameter></paramdef>
+ <paramdef><type>regprocedure </type> <parameter>callbackfunc</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>pixeltype=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>extenttype=FIRST</parameter></paramdef>
+ <paramdef choice="opt"><type>raster </type> <parameter>customextent=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancex=0</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancey=0</parameter></paramdef>
+ <paramdef choice="opt"><type>text[] </type> <parameter>VARIADIC userargs=NULL</parameter></paramdef>
+ </funcprototype>
+
+ <funcprototype>
+ <funcdef>raster <function>ST_MapAlgebra</function></funcdef>
+ <paramdef><type>raster </type> <parameter>rast</parameter></paramdef>
+ <paramdef><type>integer </type> <parameter>nband</parameter></paramdef>
+ <paramdef><type>regprocedure </type> <parameter>callbackfunc</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>pixeltype=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>extenttype=FIRST</parameter></paramdef>
+ <paramdef choice="opt"><type>raster </type> <parameter>customextent=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancex=0</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancey=0</parameter></paramdef>
+ <paramdef choice="opt"><type>text[] </type> <parameter>VARIADIC userargs=NULL</parameter></paramdef>
+ </funcprototype>
+
+ <funcprototype>
+ <funcdef>raster <function>ST_MapAlgebra</function></funcdef>
+ <paramdef><type>raster </type> <parameter>rast1</parameter></paramdef>
+ <paramdef><type>integer </type> <parameter>nband1</parameter></paramdef>
+ <paramdef><type>raster </type> <parameter>rast2</parameter></paramdef>
+ <paramdef><type>integer </type> <parameter>nband2</parameter></paramdef>
+ <paramdef><type>regprocedure </type> <parameter>callbackfunc</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>pixeltype=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>text </type> <parameter>extenttype=INTERSECTION</parameter></paramdef>
+ <paramdef choice="opt"><type>raster </type> <parameter>customextent=NULL</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancex=0</parameter></paramdef>
+ <paramdef choice="opt"><type>integer </type> <parameter>distancey=0</parameter></paramdef>
+ <paramdef choice="opt"><type>text[] </type> <parameter>VARIADIC userargs=NULL</parameter></paramdef>
+ </funcprototype>
+
+ </funcsynopsis>
+ </refsynopsisdiv>
+
+ <refsection>
+ <title>Description</title>
+
+ <para>
+ Returns a one-band raster given one or more input rasters, band indexes and one user-specified callback function.
+ </para>
+
+ <para>
+ The <varname>tworastuserfunc</varname> parameter must be the name and signature of an SQL or PL/pgSQL function, cast to a regprocedure. An example PL/pgSQL function example is:
+ <programlisting>
+CREATE OR REPLACE FUNCTION sample_callbackfunc(value double precision[][][], position integer[][], VARIADIC userargs text[])
+ RETURNS double precision
+ AS $$
+ BEGIN
+ RETURN 0;
+ END;
+ $$ LANGUAGE 'plpgsql' IMMUTABLE;
+ </programlisting>
+
+ The <varname>callbackfunc</varname> must have three arguments: a 3-dimension double precision array, a 2-dimension double precision array and a variadic 1-dimension text array. The first argument <varname>value</varname> is the set of values (as double precision) from all input rasters. The three dimensions (where indexes are 1-based) are: raster #, row y, column x. The second argument <varname>position</varname> is the set of pixel positions from the output raster and input rasters. The outer dimension (where indexes are 0-based) is the raster #. The position at outer dimension index 0 is the output raster's pixel position. For each outer dimension, there are two elements in the inner dimension for X and Y. The third argument <varname>userargs</varname> is for passing through any user-specified arguments.
+ </para>
+
+ <para>
+ Passing a <type>regprocedure</type> argument to a SQL function requires the full function signature to be passed, then cast to a <type>regprocedure</type> type. To pass the above example PL/pgSQL function as an argument, the SQL for the argument is:
+
+ <programlisting>
+'sample_callbackfunc(double precision[], integer[], text[])'::regprocedure
+ </programlisting>
+
+ Note that the argument contains the name of the function, the types of the function arguments, quotes around the name and argument types, and a cast to a <type>regprocedure</type>.
+ </para>
+
+ <para>
+ The third argument to the <varname>callbackfunc</varname> is a <type>variadic text</type> array. All trailing text arguments are passed through to the specified <varname>callbackfunc</varname>, and are contained in the <varname>userargs</varname> argument.
+ </para>
+
+ <note>
+ <para>
+ For more information about the VARIADIC keyword, please refer to the PostgreSQL documentation and the "SQL Functions with Variable Numbers of Arguments" section of <ulink url="http://www.postgresql.org/docs/current/static/xfunc-sql.html">Query Language (SQL) Functions</ulink>.
+ </para>
+ </note>
+
+ <note>
+ <para>
+ The <type>text[]</type> argument to the <varname>callbackfunc</varname> is required, regardless of whether you choose to pass any arguments to the callback function for processing or not.
+ </para>
+ </note>
+
+ <para>
+ If <varname>pixeltype</varname> is passed in, the one band of the new raster will be of that pixeltype. If pixeltype is passed NULL or left out, the new raster band will have the same pixeltype as the specified band of the first raster (for extent types: INTERSECTION, UNION, FIRST, CUSTOM) or the specified band of the appropriate raster (for extent types: SECOND, LAST). If in doubt, always specify <varname>pixeltype</varname>.
+ </para>
+
+ <para>
+ Possible values for <varname>extenttype</varname> are: INTERSECTION (default), UNION, FIRST (default for one raster variants), SECOND, LAST, CUSTOM. If <varname>extentype</varname> is CUSTOM, a raster must be provided for <varname>customextent</varname>. See example 4 of Variant 1.
+ </para>
+
+ <para>
+ Variant 1 accepts an array of <varname>rastbandarg</varname> allowing the use of a map algebra operation on many rasters and/or many bands. See example Variant 1.
+ </para>
+
+ <para>
+ Variants 2 and 3 operate upon one or more bands of one raster. See example Variant 2 and 3.
+ </para>
+
+ <para>
+ Variants 4 operate upon two rasters with one band per raster. See example Variant 4.
+ </para>
+
+<!-- enable once these functions are deprecated
+ <warning>
+ <para>
+ The functions <xref linkend="RT_ST_MapAlgebraExpr" />, <xref linkend="RT_ST_MapAlgebraExpr2" />, <xref linkend="RT_ST_MapAlgebraFct" />, <xref linkend="RT_ST_MapAlgebraFct2" /> and <xref linkend="RT_ST_MapAlgebraFctNgb" /> are deprecated as of PostGIS 2.1.0 and use of these functions are strongly discouraged. These functions will be removed from PostGIS as of 2.2.0.
+ </para>
+ </warning>
+-->
+
+ <para>Availability: 2.1.0</para>
+
+ </refsection>
+
+ <refsection>
+ <title>Examples: Variant 1</title>
+
+ <programlisting>
+-- one raster, one band
+WITH foo AS (
+ SELECT 1 AS rid, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0) AS rast
+)
+SELECT
+ ST_MapAlgebra(
+ ARRAY[ROW(rast, 1)]::rastbandarg[],
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+FROM foo
+ </programlisting>
+
+ <programlisting>
+-- one raster, several bands
+WITH foo AS (
+ SELECT 1 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0), 2, '8BUI', 10, 0), 3, '32BUI', 100, 0) AS rast
+)
+SELECT
+ ST_MapAlgebra(
+ ARRAY[ROW(rast, 3), ROW(rast, 1), ROW(rast, 3), ROW(rast, 2)]::rastbandarg[],
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+FROM foo
+ </programlisting>
+
+ <programlisting>
+-- several rasters, several bands
+WITH foo AS (
+ SELECT 1 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0), 2, '8BUI', 10, 0), 3, '32BUI', 100, 0) AS rast UNION ALL
+ SELECT 2 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 1, 1, -1, 0, 0, 0), 1, '16BUI', 2, 0), 2, '8BUI', 20, 0), 3, '32BUI', 300, 0) AS rast
+)
+SELECT
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 3), ROW(t2.rast, 1), ROW(t2.rast, 3), ROW(t1.rast, 2)]::rastbandarg[],
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+FROM foo t1
+CROSS JOIN foo t2
+WHERE t1.rid = 1
+ AND t2.rid = 2
+ </programlisting>
+
+ <programlisting>
+-- complete example of tiles of a coverage with neighborhood
+WITH foo AS (
+ SELECT 0 AS rid, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0) AS rast UNION ALL
+ SELECT 1, ST_AddBand(ST_MakeEmptyRaster(2, 2, 2, 0, 1, -1, 0, 0, 0), 1, '16BUI', 2, 0) AS rast UNION ALL
+ SELECT 2, ST_AddBand(ST_MakeEmptyRaster(2, 2, 4, 0, 1, -1, 0, 0, 0), 1, '16BUI', 3, 0) AS rast UNION ALL
+
+ SELECT 3, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, -2, 1, -1, 0, 0, 0), 1, '16BUI', 10, 0) AS rast UNION ALL
+ SELECT 4, ST_AddBand(ST_MakeEmptyRaster(2, 2, 2, -2, 1, -1, 0, 0, 0), 1, '16BUI', 20, 0) AS rast UNION ALL
+ SELECT 5, ST_AddBand(ST_MakeEmptyRaster(2, 2, 4, -2, 1, -1, 0, 0, 0), 1, '16BUI', 30, 0) AS rast UNION ALL
+
+ SELECT 6, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, -4, 1, -1, 0, 0, 0), 1, '16BUI', 100, 0) AS rast UNION ALL
+ SELECT 7, ST_AddBand(ST_MakeEmptyRaster(2, 2, 2, -4, 1, -1, 0, 0, 0), 1, '16BUI', 200, 0) AS rast UNION ALL
+ SELECT 8, ST_AddBand(ST_MakeEmptyRaster(2, 2, 4, -4, 1, -1, 0, 0, 0), 1, '16BUI', 300, 0) AS rast
+)
+SELECT
+ t1.rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(ST_Union(t2.rast), 1)]::rastbandarg[],
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure,
+ '32BUI',
+ 'CUSTOM', t1.rast,
+ 1, 1
+ ) AS rast
+FROM raster_nmapalgebra_in t1
+CROSS JOIN raster_nmapalgebra_in t2
+WHERE t1.rid = 4
+ AND t2.rid BETWEEN 0 AND 8
+ AND ST_Intersects(t1.rast, t2.rast)
+GROUP BY t1.rid, t1.rast
+ </programlisting>
+
+ </refsection>
+
+ <refsection>
+ <title>Examples: Variants 2 and 3</title>
+
+ <programlisting>
+-- one raster, several bands
+WITH foo AS (
+ SELECT 1 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0), 2, '8BUI', 10, 0), 3, '32BUI', 100, 0) AS rast
+)
+SELECT
+ ST_MapAlgebra(
+ rast, ARRAY[3, 1, 3, 2]::integer[],
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+FROM foo
+ </programlisting>
+
+ <programlisting>
+-- one raster, one band
+WITH foo AS (
+ SELECT 1 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0), 2, '8BUI', 10, 0), 3, '32BUI', 100, 0) AS rast
+)
+SELECT
+ ST_MapAlgebra(
+ rast, 2,
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+FROM foo
+ </programlisting>
+
+ </refsection>
+
+ <refsection>
+ <title>Examples: Variant 4</title>
+
+ <programlisting>
+-- two rasters, two bands
+WITH foo AS (
+ SELECT 1 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0), 2, '8BUI', 10, 0), 3, '32BUI', 100, 0) AS rast UNION ALL
+ SELECT 2 AS rid, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 1, 1, -1, 0, 0, 0), 1, '16BUI', 2, 0), 2, '8BUI', 20, 0), 3, '32BUI', 300, 0) AS rast
+)
+SELECT
+ ST_MapAlgebra(
+ t1.rast, 2,
+ t2.rast, 1,
+ 'sample_callbackfunc(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+FROM foo t1
+CROSS JOIN foo t2
+WHERE t1.rid = 1
+ AND t2.rid = 2
+ </programlisting>
+
+ </refsection>
+
+ <refsection>
+ <title>See Also</title>
+
+ <para>
+ <xref linkend="rastbandarg" />,
+ <xref linkend="RT_ST_Union" />
+ </para>
+ </refsection>
+
+
+ </refentry>
+
<refentry id="RT_ST_MapAlgebraExpr">
<refnamediv>
<refname>ST_MapAlgebraExpr</refname>
return NULL;
}
+ else if (rt_raster_is_empty(rast)) {
+ rterror("rt_raster_iterator: Computed raster for %s extent is empty",
+ extenttype == ET_UNION ? "union" : "intersection"
+ );
+
+ _rti_param_destroy(_param);
+
+ return NULL;
+ }
rtnrast = rast;
rast = NULL;
}
Datum RASTER_copyBand(PG_FUNCTION_ARGS);
Datum RASTER_addBandRasterArray(PG_FUNCTION_ARGS);
-/* Raster analysis */
-Datum RASTER_mapAlgebraExpr(PG_FUNCTION_ARGS);
-Datum RASTER_mapAlgebraFct(PG_FUNCTION_ARGS);
-
/* create new raster from existing raster's bands */
Datum RASTER_band(PG_FUNCTION_ARGS);
/* determine if two rasters are aligned */
Datum RASTER_sameAlignment(PG_FUNCTION_ARGS);
-/* two-raster MapAlgebra */
-Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS);
+/* one-raster MapAlgebra */
+Datum RASTER_mapAlgebraExpr(PG_FUNCTION_ARGS);
+Datum RASTER_mapAlgebraFct(PG_FUNCTION_ARGS);
/* one-raster neighborhood MapAlgebra */
Datum RASTER_mapAlgebraFctNgb(PG_FUNCTION_ARGS);
+/* two-raster MapAlgebra */
+Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS);
+
+/* n-raster MapAlgebra */
+Datum RASTER_nMapAlgebra(PG_FUNCTION_ARGS);
+
/* raster union aggregate */
Datum RASTER_union_transfn(PG_FUNCTION_ARGS);
Datum RASTER_union_finalfn(PG_FUNCTION_ARGS);
}
/* ---------------------------------------------------------------- */
-/* ST_Union aggregate functions */
+/* n-raster MapAlgebra */
/* ---------------------------------------------------------------- */
-typedef enum {
- UT_LAST = 0,
- UT_FIRST,
- UT_MIN,
- UT_MAX,
- UT_COUNT,
- UT_SUM,
- UT_MEAN
-} rtpg_union_type;
+typedef struct {
+ Oid ufc_noid;
+ FmgrInfo ufl_info;
+ FunctionCallInfoData ufc_info;
+ int ufc_nullcount;
+} rtpg_nmapalgebra_callback_arg;
-/* internal function translating text of UNION type to enum */
-static rtpg_union_type rtpg_uniontype_index_from_name(const char *cutype) {
- assert(cutype && strlen(cutype) > 0);
+typedef struct rtpg_nmapalgebra_arg_t *rtpg_nmapalgebra_arg;
+struct rtpg_nmapalgebra_arg_t {
+ int numraster;
+ rt_pgraster **pgraster;
+ rt_raster *raster;
+ uint8_t *isempty; /* flag indicating if raster is empty */
+ uint8_t *ownsdata; /* is the raster self owned or just a pointer to another raster */
+ int *nband; /* source raster's band index, 0-based */
+ uint8_t *hasband; /* does source raster have band at index nband? */
- if (strcmp(cutype, "LAST") == 0)
- return UT_LAST;
- else if (strcmp(cutype, "FIRST") == 0)
- return UT_FIRST;
- else if (strcmp(cutype, "MIN") == 0)
- return UT_MIN;
- else if (strcmp(cutype, "MAX") == 0)
- return UT_MAX;
- else if (strcmp(cutype, "COUNT") == 0)
- return UT_COUNT;
- else if (strcmp(cutype, "SUM") == 0)
- return UT_SUM;
- else if (strcmp(cutype, "MEAN") == 0)
- return UT_MEAN;
+ rt_pixtype pixtype; /* output raster band's pixel type */
+ int hasnodata; /* NODATA flag */
+ double nodataval; /* NODATA value */
- return UT_LAST;
-}
+ int distance[2]; /* distance in X and Y axis */
-typedef struct rtpg_union_band_arg_t *rtpg_union_band_arg;
-struct rtpg_union_band_arg_t {
- int nband; /* source raster's band index, 0-based */
- rtpg_union_type uniontype;
+ rt_extenttype extenttype; /* oupt raster's extent type */
+ rt_pgraster *pgcextent; /* custom extent of type rt_pgraster */
+ rt_raster cextent; /* custom extent of type rt_raster */
- int numraster;
- rt_raster *raster;
+ rtpg_nmapalgebra_callback_arg callback;
};
-typedef struct rtpg_union_arg_t *rtpg_union_arg;
-struct rtpg_union_arg_t {
- int numband; /* number of bandargs */
- rtpg_union_band_arg bandarg;
-};
+static rtpg_nmapalgebra_arg rtpg_nmapalgebra_arg_init() {
+ rtpg_nmapalgebra_arg arg = NULL;
-static void rtpg_union_arg_destroy(rtpg_union_arg arg) {
+ arg = palloc(sizeof(struct rtpg_nmapalgebra_arg_t));
+ if (arg == NULL) {
+ elog(ERROR, "rtpg_nmapalgebra_arg_init: Unable to allocate memory for arguments");
+ return 0;
+ }
+
+ arg->numraster = 0;
+ arg->pgraster = NULL;
+ arg->raster = NULL;
+ arg->isempty = NULL;
+ arg->ownsdata = NULL;
+ arg->nband = NULL;
+ arg->hasband = NULL;
+
+ arg->pixtype = PT_END;
+ arg->hasnodata = 1;
+ arg->nodataval = 0;
+
+ arg->distance[0] = 0;
+ arg->distance[1] = 0;
+
+ arg->extenttype = ET_INTERSECTION;
+ arg->pgcextent = NULL;
+ arg->cextent = NULL;
+
+ arg->callback.ufc_noid = InvalidOid;
+ arg->callback.ufc_nullcount = 0;
+
+ return arg;
+}
+
+static void rtpg_nmapalgebra_arg_destroy(rtpg_nmapalgebra_arg arg) {
int i = 0;
- int j = 0;
- if (arg->bandarg != NULL) {
- for (i = 0; i < arg->numband; i++) {
- if (!arg->bandarg[i].numraster)
+ if (arg->raster != NULL) {
+ for (i = 0; i < arg->numraster; i++) {
+ if (arg->raster[i] == NULL || !arg->ownsdata[i])
continue;
- for (j = 0; j < arg->bandarg[i].numraster; j++) {
- if (arg->bandarg[i].raster[j] == NULL)
- continue;
- rt_raster_destroy(arg->bandarg[i].raster[j]);
- }
-
- pfree(arg->bandarg[i].raster);
+ rt_raster_destroy(arg->raster[i]);
}
- pfree(arg->bandarg);
+ pfree(arg->raster);
+ pfree(arg->pgraster);
+ pfree(arg->isempty);
+ pfree(arg->ownsdata);
+ pfree(arg->nband);
}
+ if (arg->cextent != NULL)
+ rt_raster_destroy(arg->cextent);
+
pfree(arg);
}
-static int rtpg_union_callback(
- rt_iterator_arg arg, void *userarg,
- double *value, int *nodata
-) {
- rtpg_union_type *utype = (rtpg_union_type *) userarg;
+static int rtpg_nmapalgebra_rastbandarg_process(rtpg_nmapalgebra_arg arg, ArrayType *array, int *allnull, int *allempty, int *noband) {
+ Oid etype;
+ Datum *e;
+ bool *nulls;
+ int16 typlen;
+ bool typbyval;
+ char typalign;
+ int n = 0;
- if (arg == NULL)
+ HeapTupleHeader tup;
+ bool isnull;
+ Datum tupv;
+
+ int i;
+ int j;
+ int nband;
+
+ if (arg == NULL || array == NULL) {
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: NULL values not permitted for parameters");
+ return 0;
+ }
+
+ etype = ARR_ELEMTYPE(array);
+ get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
+
+ deconstruct_array(
+ array,
+ etype,
+ typlen, typbyval, typalign,
+ &e, &nulls, &n
+ );
+
+ if (!n) {
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: Invalid argument for rastbandarg");
return 0;
+ }
+ /* prep arg */
+ arg->numraster = n;
+ arg->pgraster = palloc(sizeof(rt_pgraster *) * arg->numraster);
+ arg->raster = palloc(sizeof(rt_raster) * arg->numraster);
+ arg->isempty = palloc(sizeof(uint8_t) * arg->numraster);
+ arg->ownsdata = palloc(sizeof(uint8_t) * arg->numraster);
+ arg->nband = palloc(sizeof(int) * arg->numraster);
+ arg->hasband = palloc(sizeof(uint8_t) * arg->numraster);
if (
- arg->rasters != 2 ||
- arg->rows != 1 ||
- arg->columns != 1
+ arg->pgraster == NULL ||
+ arg->raster == NULL ||
+ arg->isempty == NULL ||
+ arg->ownsdata == NULL ||
+ arg->nband == NULL ||
+ arg->hasband == NULL
) {
- elog(ERROR, "rtpg_union_callback: Invalid arguments passed to callback");
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: Unable to allocate memory for processing rastbandarg");
return 0;
}
- *value = 0;
- *nodata = 0;
+ *allnull = 0;
+ *allempty = 0;
+ *noband = 0;
- /* handle NODATA situations except for COUNT, which is a special case */
- if (*utype != UT_COUNT) {
- /* both NODATA */
- if (arg->nodata[0][0][0] && arg->nodata[1][0][0]) {
- *nodata = 1;
- POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
- return 1;
+ /* process each element */
+ for (i = 0; i < n; i++) {
+ if (nulls[i]) {
+ arg->numraster--;
+ continue;
}
- /* second NODATA */
- else if (!arg->nodata[0][0][0] && arg->nodata[1][0][0]) {
- *value = arg->values[0][0][0];
- POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
- return 1;
+
+ POSTGIS_RT_DEBUGF(4, "Processing rastbandarg at index %d", i);
+
+ arg->raster[i] = NULL;
+ arg->isempty[i] = 0;
+ arg->ownsdata[i] = 1;
+ arg->nband[i] = 0;
+ arg->hasband[i] = 0;
+
+ /* each element is a tuple */
+ tup = (HeapTupleHeader) DatumGetPointer(e[i]);
+ if (NULL == tup) {
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: Invalid argument for rastbandarg at index %d", i);
+ return 0;
}
- /* first NODATA */
- else if (arg->nodata[0][0][0] && !arg->nodata[1][0][0]) {
- *value = arg->values[1][0][0];
- POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
- return 1;
+
+ /* first element, raster */
+ POSTGIS_RT_DEBUG(4, "Processing first element (raster)");
+ tupv = GetAttributeByName(tup, "rast", &isnull);
+ if (isnull) {
+ elog(NOTICE, "First argument (nband) of rastbandarg at index %d is NULL. Assuming NULL raster", i);
+ arg->isempty[i] = 1;
+ arg->ownsdata[i] = 0;
+
+ (*allnull)++;
+ (*allempty)++;
+ (*noband)++;
+
+ continue;
}
- }
+ else {
+ arg->pgraster[i] = (rt_pgraster *) PG_DETOAST_DATUM(tupv);
+
+ /* see if this is a copy of an existing pgraster */
+ for (j = 0; j < i; j++) {
+ if (arg->pgraster[i] == arg->pgraster[j]) {
+ POSTGIS_RT_DEBUG(4, "raster matching existing same raster found");
+ arg->raster[i] = arg->raster[j];
+ arg->ownsdata[i] = 0;
+ break;
+ }
+ }
- switch (*utype) {
- case UT_FIRST:
- *value = arg->values[0][0][0];
- break;
- case UT_MIN:
- if (arg->values[0][0][0] < arg->values[1][0][0])
- *value = arg->values[0][0][0];
- else
- *value = arg->values[1][0][0];
- break;
- case UT_MAX:
- if (arg->values[0][0][0] > arg->values[1][0][0])
- *value = arg->values[0][0][0];
- else
- *value = arg->values[1][0][0];
- break;
- case UT_COUNT:
- /* both NODATA */
- if (arg->nodata[0][0][0] && arg->nodata[1][0][0])
- *value = 0;
- /* second NODATA */
- else if (!arg->nodata[0][0][0] && arg->nodata[1][0][0])
- *value = arg->values[0][0][0];
- /* first NODATA */
- else if (arg->nodata[0][0][0] && !arg->nodata[1][0][0])
- *value = 1;
- /* has value, increment */
- else
- *value = arg->values[0][0][0] + 1;
- break;
- case UT_SUM:
- *value = arg->values[0][0][0] + arg->values[1][0][0];
- break;
- case UT_MEAN:
- break;
- case UT_LAST:
- default:
- *value = arg->values[1][0][0];
- break;
+ if (arg->ownsdata[i]) {
+ POSTGIS_RT_DEBUG(4, "deserializing raster");
+ arg->raster[i] = rt_raster_deserialize(arg->pgraster[i], FALSE);
+ if (arg->raster[i] == NULL) {
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: Could not deserialize raster at index %d", i);
+ return 0;
+ }
+ }
+ }
+
+ /* is raster empty? */
+ arg->isempty[i] = rt_raster_is_empty(arg->raster[i]);
+ if (arg->isempty[i]) {
+ (*allempty)++;
+ (*noband)++;
+
+ continue;
+ }
+
+ /* second element, nband */
+ POSTGIS_RT_DEBUG(4, "Processing second element (nband)");
+ tupv = GetAttributeByName(tup, "nband", &isnull);
+ if (isnull) {
+ nband = 1;
+ elog(NOTICE, "First argument (nband) of rastbandarg at index %d is NULL. Assuming nband = %d", i, nband);
+ }
+ else
+ nband = DatumGetInt32(tupv);
+
+ if (nband < 1) {
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: Band number provided for rastbandarg at index %d must be greater than zero (1-based)", i);
+ return 0;
+ }
+
+ arg->nband[i] = nband - 1;
+ arg->hasband[i] = rt_raster_has_band(arg->raster[i], arg->nband[i]);
+ if (!arg->hasband[i]) {
+ (*noband)++;
+ POSTGIS_RT_DEBUGF(4, "Band at index %d not found in raster", nband);
+ }
}
- POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
+ if (arg->numraster < n) {
+ arg->pgraster = repalloc(arg->pgraster, sizeof(rt_pgraster *) * arg->numraster);
+ arg->raster = repalloc(arg->raster, sizeof(rt_raster) * arg->numraster);
+ arg->isempty = repalloc(arg->isempty, sizeof(uint8_t) * arg->numraster);
+ arg->ownsdata = repalloc(arg->ownsdata, sizeof(uint8_t) * arg->numraster);
+ arg->nband = repalloc(arg->nband, sizeof(int) * arg->numraster);
+ arg->hasband = repalloc(arg->hasband, sizeof(uint8_t) * arg->numraster);
+ if (
+ arg->pgraster == NULL ||
+ arg->raster == NULL ||
+ arg->isempty == NULL ||
+ arg->ownsdata == NULL ||
+ arg->nband == NULL ||
+ arg->hasband == NULL
+ ) {
+ elog(ERROR, "rtpg_nmapalgebra_rastbandarg_process: Unable to reallocate memory for processed rastbandarg");
+ return 0;
+ }
+ }
+ POSTGIS_RT_DEBUGF(4, "arg->numraster = %d", arg->numraster);
return 1;
}
-static int rtpg_union_mean_callback(
+/*
+ Callback for RASTER_nMapAlgebra
+*/
+static int rtpg_nmapalgebra_callback(
rt_iterator_arg arg, void *userarg,
double *value, int *nodata
) {
- if (arg == NULL)
- return 0;
+ rtpg_nmapalgebra_callback_arg *callback = (rtpg_nmapalgebra_callback_arg *) userarg;
- if (
- arg->rasters != 2 ||
- arg->rows != 1 ||
- arg->columns != 1
- ) {
- elog(ERROR, "rtpg_union_mean_callback: Invalid arguments passed to callback");
+ int16 typlen;
+ bool typbyval;
+ char typalign;
+
+ ArrayType *mdValues = NULL;
+ Datum *_values = NULL;
+ bool *_nodata = NULL;
+
+ ArrayType *mdPos = NULL;
+ Datum *_pos = NULL;
+ bool *_null = NULL;
+
+ int i = 0;
+ int x = 0;
+ int y = 0;
+ int z = 0;
+ int dim[3] = {0};
+ int lbound[3] = {1, 1, 1};
+ Datum datum = (Datum) NULL;
+
+ if (arg == NULL)
return 0;
- }
*value = 0;
- *nodata = 1;
+ *nodata = 0;
- POSTGIS_RT_DEBUGF(4, "rast0: %f %d", arg->values[0][0][0], arg->nodata[0][0][0]);
- POSTGIS_RT_DEBUGF(4, "rast1: %f %d", arg->values[1][0][0], arg->nodata[1][0][0]);
+ dim[0] = arg->rasters;
+ dim[1] = arg->rows;
+ dim[2] = arg->columns;
- if (
- !arg->nodata[0][0][0] &&
- FLT_NEQ(arg->values[0][0][0], 0) &&
- !arg->nodata[1][0][0]
- ) {
- *value = arg->values[1][0][0] / arg->values[0][0][0];
- *nodata = 0;
+ _values = palloc(sizeof(Datum) * arg->rasters * arg->rows * arg->columns);
+ _nodata = palloc(sizeof(bool) * arg->rasters * arg->rows * arg->columns);
+ if (_values == NULL || _nodata == NULL) {
+ elog(ERROR, "rtpg_nmapalgebra_callback: Unable to allocate memory for values array");
+ return 0;
}
- POSTGIS_RT_DEBUGF(4, "value, nodata = (%f, %d)", *value, *nodata);
+ /* build mdValues */
+ i = 0;
+ /* raster */
+ for (z = 0; z < arg->rasters; z++) {
+ /* Y axis */
+ for (y = 0; y < arg->rows; y++) {
+ /* X axis */
+ for (x = 0; x < arg->columns; x++) {
+ POSTGIS_RT_DEBUGF(4, "(z, y ,x) = (%d, %d, %d)", z, y, x);
+ POSTGIS_RT_DEBUGF(4, "(value, nodata) = (%f, %d)", arg->values[z][y][x], arg->nodata[z][y][x]);
+
+ _nodata[i] = (bool) arg->nodata[z][y][x];
+ if (!_nodata[i])
+ _values[i] = Float8GetDatum(arg->values[z][y][x]);
+ else
+ _values[i] = (Datum) NULL;
- return 1;
-}
+ i++;
+ }
+ }
+ }
+
+ /* info about the type of item in the multi-dimensional array (float8). */
+ get_typlenbyvalalign(FLOAT8OID, &typlen, &typbyval, &typalign);
+
+ /* construct mdValues */
+ mdValues = construct_md_array(
+ _values, _nodata,
+ 3, dim, lbound,
+ FLOAT8OID,
+ typlen, typbyval, typalign
+ );
+ pfree(_nodata);
+ pfree(_values);
+
+ _pos = palloc(sizeof(Datum) * (arg->rasters + 1) * 2);
+ _null = palloc(sizeof(bool) * (arg->rasters + 1) * 2);
+ if (_pos == NULL || _null == NULL) {
+ elog(ERROR, "rtpg_nmapalgebra_callback: Unable to allocate memory for position array");
+ pfree(mdValues);
+ return 0;
+ }
+ memset(_null, 0, sizeof(bool) * (arg->rasters + 1) * 2);
+
+ /* build mdPos */
+ i = 0;
+ _pos[i] = arg->dst_pixel[0] + 1;
+ i++;
+ _pos[i] = arg->dst_pixel[1] + 1;
+ i++;
+
+ for (z = 0; z < arg->rasters; z++) {
+ _pos[i] = arg->src_pixel[z][0] + 1;
+ i++;
+
+ _pos[i] = arg->src_pixel[z][1] + 1;
+ i++;
+ }
+
+ /* info about the type of item in the multi-dimensional array (int4). */
+ get_typlenbyvalalign(INT4OID, &typlen, &typbyval, &typalign);
+
+ /* reuse dim and lbound, just tweak to what we need */
+ dim[1] = 2;
+ lbound[0] = 0;
+
+ /* construct mdValues */
+ mdPos = construct_md_array(
+ _pos, _null,
+ 2, dim, lbound,
+ INT4OID,
+ typlen, typbyval, typalign
+ );
+ pfree(_pos);
+ pfree(_null);
+
+ callback->ufc_info.arg[0] = PointerGetDatum(mdValues);
+ callback->ufc_info.arg[1] = PointerGetDatum(mdPos);
+
+ /* function is strict and null parameter is passed */
+ /* http://archives.postgresql.org/pgsql-general/2011-11/msg00424.php */
+ if (callback->ufl_info.fn_strict && callback->ufc_nullcount) {
+ *nodata = 1;
+
+ pfree(mdValues);
+ pfree(mdPos);
+
+ return 1;
+ }
+
+ /* call user callback function */
+ datum = FunctionCallInvoke(&(callback->ufc_info));
+ pfree(mdValues);
+ pfree(mdPos);
+
+ /* result is not null*/
+ if (!callback->ufc_info.isnull)
+ *value = DatumGetFloat8(datum);
+ else
+ *nodata = 1;
+
+ return 1;
+}
+
+/*
+ ST_MapAlgebra for n rasters
+*/
+PG_FUNCTION_INFO_V1(RASTER_nMapAlgebra);
+Datum RASTER_nMapAlgebra(PG_FUNCTION_ARGS)
+{
+ rtpg_nmapalgebra_arg arg = NULL;
+ rt_iterator itrset;
+ int i = 0;
+ int noerr = 0;
+ int allnull = 0;
+ int allempty = 0;
+ int noband = 0;
+
+ rt_raster raster = NULL;
+ rt_band band = NULL;
+ rt_pgraster *pgraster = NULL;
+
+ POSTGIS_RT_DEBUG(3, "Starting...");
+
+ if (PG_ARGISNULL(0))
+ PG_RETURN_NULL();
+
+ /* init argument struct */
+ arg = rtpg_nmapalgebra_arg_init();
+ if (arg == NULL) {
+ elog(ERROR, "RASTER_nMapAlgebra: Unable to initialize argument structure");
+ PG_RETURN_NULL();
+ }
+
+ /* let helper function process rastbandarg (0) */
+ if (!rtpg_nmapalgebra_rastbandarg_process(arg, PG_GETARG_ARRAYTYPE_P(0), &allnull, &allempty, &noband)) {
+ elog(ERROR, "RASTER_nMapAlgebra: Unable to process rastbandarg");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ POSTGIS_RT_DEBUGF(4, "allnull, allempty, noband = %d, %d, %d", allnull, allempty, noband);
+
+ /* all rasters are NULL, return NULL */
+ if (allnull == arg->numraster) {
+ elog(NOTICE, "All input rasters are NULL. Returning NULL");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ /* pixel type (2) */
+ if (!PG_ARGISNULL(2)) {
+ char *pixtypename = text_to_cstring(PG_GETARG_TEXT_P(2));
+
+ /* Get the pixel type index */
+ arg->pixtype = rt_pixtype_index_from_name(pixtypename);
+ if (arg->pixtype == PT_END) {
+ elog(ERROR, "RASTER_nMapAlgebra: Invalid pixel type: %s", pixtypename);
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+ }
+
+ /* distancex (3) */
+ if (!PG_ARGISNULL(3))
+ arg->distance[0] = PG_GETARG_INT32(3);
+ /* distancey (4) */
+ if (!PG_ARGISNULL(4))
+ arg->distance[1] = PG_GETARG_INT32(4);
+
+ if (arg->distance[0] < 0 || arg->distance[1] < 0) {
+ elog(ERROR, "RASTER_nMapAlgebra: Distance for X and Y axis must be greater than or equal to zero");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ /* extent type (5) */
+ if (!PG_ARGISNULL(5)) {
+ char *extenttypename = rtpg_strtoupper(rtpg_trim(text_to_cstring(PG_GETARG_TEXT_P(5))));
+ arg->extenttype = rt_util_extent_type(extenttypename);
+ }
+ POSTGIS_RT_DEBUGF(4, "extenttype: %d", arg->extenttype);
+
+ /* custom extent (6) */
+ if (arg->extenttype == ET_CUSTOM) {
+ if (PG_ARGISNULL(6)) {
+ elog(NOTICE, "Custom extent is NULL. Returning NULL");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ arg->pgcextent = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(6));
+
+ /* only need the raster header */
+ arg->cextent = rt_raster_deserialize(arg->pgcextent, TRUE);
+ if (arg->cextent == NULL) {
+ elog(ERROR, "RASTER_nMapAlgebra: Could not deserialize custom extent");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+ else if (rt_raster_is_empty(arg->cextent)) {
+ elog(ERROR, "RASTER_nMapAlgebra: Custom extent is an empty raster. Returning empty raster");
+ rtpg_nmapalgebra_arg_destroy(arg);
+
+ raster = rt_raster_new(0, 0);
+ if (raster == NULL) {
+ elog(ERROR, "RASTER_nMapAlgebra: Unable to create empty raster");
+ PG_RETURN_NULL();
+ }
+
+ pgraster = rt_raster_serialize(raster);
+ rt_raster_destroy(raster);
+ if (!pgraster) PG_RETURN_NULL();
+
+ SET_VARSIZE(pgraster, pgraster->size);
+ PG_RETURN_POINTER(pgraster);
+ }
+ }
+
+ noerr = 1;
+ /* all rasters are empty, return empty raster */
+ if (allempty == arg->numraster) {
+ elog(NOTICE, "All input rasters are empty. Returning empty raster");
+ noerr = 0;
+ }
+ /* all rasters don't have indicated band, return empty raster */
+ else if (noband == arg->numraster) {
+ elog(NOTICE, "All input rasters do not have bands at indicated indexes. Returning empty raster");
+ noerr = 0;
+ }
+ if (!noerr) {
+ rtpg_nmapalgebra_arg_destroy(arg);
+
+ raster = rt_raster_new(0, 0);
+ if (raster == NULL) {
+ elog(ERROR, "RASTER_nMapAlgebra: Unable to create empty raster");
+ PG_RETURN_NULL();
+ }
+
+ pgraster = rt_raster_serialize(raster);
+ rt_raster_destroy(raster);
+ if (!pgraster) PG_RETURN_NULL();
+
+ SET_VARSIZE(pgraster, pgraster->size);
+ PG_RETURN_POINTER(pgraster);
+ }
+
+ /* do regprocedure last (1) */
+ if (!PG_ARGISNULL(1) || get_fn_expr_argtype(fcinfo->flinfo, 1) == REGPROCEDUREOID) {
+ POSTGIS_RT_DEBUG(4, "processing callbackfunc");
+ arg->callback.ufc_noid = PG_GETARG_OID(1);
+
+ /* get function info */
+ fmgr_info(arg->callback.ufc_noid, &(arg->callback.ufl_info));
+
+ /* function cannot return set */
+ noerr = 1;
+ if (arg->callback.ufl_info.fn_retset) {
+ elog(ERROR, "RASTER_nMapAlgebra: Function provided must return double precision, not resultset");
+ noerr = 0;
+ }
+ /* function should have correct # of args */
+ else if (arg->callback.ufl_info.fn_nargs != 3) {
+ elog(ERROR, "RASTER_nMapAlgebra: Function provided must have three input parameters");
+ noerr = 0;
+ }
+
+ /*
+ TODO: consider adding checks of the userfunction parameters
+ should be able to use get_fn_expr_argtype() of fmgr.c
+ */
+
+ if (!noerr) {
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ if (func_volatile(arg->callback.ufc_noid) == 'v')
+ elog(NOTICE, "Function provided is VOLATILE. Unless required and for best performance, function should be IMMUTABLE or STABLE");
+
+ /* prep function call data */
+ InitFunctionCallInfoData(arg->callback.ufc_info, &(arg->callback.ufl_info), arg->callback.ufl_info.fn_nargs, InvalidOid, NULL, NULL);
+ memset(arg->callback.ufc_info.argnull, FALSE, arg->callback.ufl_info.fn_nargs);
+
+ /* userargs (7) */
+ if (!PG_ARGISNULL(7))
+ arg->callback.ufc_info.arg[2] = PG_GETARG_DATUM(7);
+ else {
+ arg->callback.ufc_info.arg[2] = (Datum) NULL;
+ arg->callback.ufc_info.argnull[2] = TRUE;
+ arg->callback.ufc_nullcount++;
+ }
+ }
+ else {
+ elog(ERROR, "RASTER_nMapAlgebra: callbackfunc must be provided");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ /* determine nodataval and possibly pixtype */
+ /* band to check */
+ switch (arg->extenttype) {
+ case ET_LAST:
+ i = arg->numraster - 1;
+ break;
+ case ET_SECOND:
+ if (arg->numraster > 1) {
+ i = 1;
+ break;
+ }
+ default:
+ i = 0;
+ break;
+ }
+ /* find first viable band */
+ if (!arg->hasband[i]) {
+ for (i = 0; i < arg->numraster; i++) {
+ if (arg->hasband[i])
+ break;
+ }
+ }
+ band = rt_raster_get_band(arg->raster[i], arg->nband[i]);
+
+ /* set pixel type if PT_END */
+ if (arg->pixtype == PT_END)
+ arg->pixtype = rt_band_get_pixtype(band);
+
+ /* set hasnodata and nodataval */
+ arg->hasnodata = 1;
+ if (rt_band_get_hasnodata_flag(band))
+ arg->nodataval = rt_band_get_nodata(band);
+ else
+ arg->nodataval = rt_band_get_min_value(band);
+
+ POSTGIS_RT_DEBUGF(4, "pixtype, hasnodata, nodataval: %s, %d, %f", rt_pixtype_name(arg->pixtype), arg->hasnodata, arg->nodataval);
+
+ /* init itrset */
+ itrset = palloc(sizeof(struct rt_iterator_t) * arg->numraster);
+ if (itrset == NULL) {
+ elog(ERROR, "RASTER_nMapAlgebra: Unable to allocate memory for iterator arguments");
+ rtpg_nmapalgebra_arg_destroy(arg);
+ PG_RETURN_NULL();
+ }
+
+ /* set itrset */
+ for (i = 0; i < arg->numraster; i++) {
+ itrset[i].raster = arg->raster[i];
+ itrset[i].nband = arg->nband[i];
+ itrset[i].nbnodata = 1;
+ }
+
+ /* pass everything to iterator */
+ raster = rt_raster_iterator(
+ itrset, arg->numraster,
+ arg->extenttype, arg->cextent,
+ arg->pixtype,
+ arg->hasnodata, arg->nodataval,
+ arg->distance[0], arg->distance[1],
+ &(arg->callback),
+ rtpg_nmapalgebra_callback,
+ &noerr
+ );
+
+ /* cleanup */
+ pfree(itrset);
+ rtpg_nmapalgebra_arg_destroy(arg);
+
+ if (!noerr) {
+ elog(ERROR, "RASTER_nMapAlgebra: Unable to run raster interator function");
+ PG_RETURN_NULL();
+ }
+
+ pgraster = rt_raster_serialize(raster);
+ rt_raster_destroy(raster);
+
+ POSTGIS_RT_DEBUG(3, "Finished");
+
+ if (!pgraster)
+ PG_RETURN_NULL();
+
+ SET_VARSIZE(pgraster, pgraster->size);
+ PG_RETURN_POINTER(pgraster);
+}
+
+/* ---------------------------------------------------------------- */
+/* ST_Union aggregate functions */
+/* ---------------------------------------------------------------- */
+
+typedef enum {
+ UT_LAST = 0,
+ UT_FIRST,
+ UT_MIN,
+ UT_MAX,
+ UT_COUNT,
+ UT_SUM,
+ UT_MEAN
+} rtpg_union_type;
+
+/* internal function translating text of UNION type to enum */
+static rtpg_union_type rtpg_uniontype_index_from_name(const char *cutype) {
+ assert(cutype && strlen(cutype) > 0);
+
+ if (strcmp(cutype, "LAST") == 0)
+ return UT_LAST;
+ else if (strcmp(cutype, "FIRST") == 0)
+ return UT_FIRST;
+ else if (strcmp(cutype, "MIN") == 0)
+ return UT_MIN;
+ else if (strcmp(cutype, "MAX") == 0)
+ return UT_MAX;
+ else if (strcmp(cutype, "COUNT") == 0)
+ return UT_COUNT;
+ else if (strcmp(cutype, "SUM") == 0)
+ return UT_SUM;
+ else if (strcmp(cutype, "MEAN") == 0)
+ return UT_MEAN;
+
+ return UT_LAST;
+}
+
+typedef struct rtpg_union_band_arg_t *rtpg_union_band_arg;
+struct rtpg_union_band_arg_t {
+ int nband; /* source raster's band index, 0-based */
+ rtpg_union_type uniontype;
+
+ int numraster;
+ rt_raster *raster;
+};
+
+typedef struct rtpg_union_arg_t *rtpg_union_arg;
+struct rtpg_union_arg_t {
+ int numband; /* number of bandargs */
+ rtpg_union_band_arg bandarg;
+};
+
+static void rtpg_union_arg_destroy(rtpg_union_arg arg) {
+ int i = 0;
+ int j = 0;
+
+ if (arg->bandarg != NULL) {
+ for (i = 0; i < arg->numband; i++) {
+ if (!arg->bandarg[i].numraster)
+ continue;
+
+ for (j = 0; j < arg->bandarg[i].numraster; j++) {
+ if (arg->bandarg[i].raster[j] == NULL)
+ continue;
+ rt_raster_destroy(arg->bandarg[i].raster[j]);
+ }
+
+ pfree(arg->bandarg[i].raster);
+ }
+
+ pfree(arg->bandarg);
+ }
+
+ pfree(arg);
+}
+
+static int rtpg_union_callback(
+ rt_iterator_arg arg, void *userarg,
+ double *value, int *nodata
+) {
+ rtpg_union_type *utype = (rtpg_union_type *) userarg;
+
+ if (arg == NULL)
+ return 0;
+
+ if (
+ arg->rasters != 2 ||
+ arg->rows != 1 ||
+ arg->columns != 1
+ ) {
+ elog(ERROR, "rtpg_union_callback: Invalid arguments passed to callback");
+ return 0;
+ }
+
+ *value = 0;
+ *nodata = 0;
+
+ /* handle NODATA situations except for COUNT, which is a special case */
+ if (*utype != UT_COUNT) {
+ /* both NODATA */
+ if (arg->nodata[0][0][0] && arg->nodata[1][0][0]) {
+ *nodata = 1;
+ POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
+ return 1;
+ }
+ /* second NODATA */
+ else if (!arg->nodata[0][0][0] && arg->nodata[1][0][0]) {
+ *value = arg->values[0][0][0];
+ POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
+ return 1;
+ }
+ /* first NODATA */
+ else if (arg->nodata[0][0][0] && !arg->nodata[1][0][0]) {
+ *value = arg->values[1][0][0];
+ POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
+ return 1;
+ }
+ }
+
+ switch (*utype) {
+ case UT_FIRST:
+ *value = arg->values[0][0][0];
+ break;
+ case UT_MIN:
+ if (arg->values[0][0][0] < arg->values[1][0][0])
+ *value = arg->values[0][0][0];
+ else
+ *value = arg->values[1][0][0];
+ break;
+ case UT_MAX:
+ if (arg->values[0][0][0] > arg->values[1][0][0])
+ *value = arg->values[0][0][0];
+ else
+ *value = arg->values[1][0][0];
+ break;
+ case UT_COUNT:
+ /* both NODATA */
+ if (arg->nodata[0][0][0] && arg->nodata[1][0][0])
+ *value = 0;
+ /* second NODATA */
+ else if (!arg->nodata[0][0][0] && arg->nodata[1][0][0])
+ *value = arg->values[0][0][0];
+ /* first NODATA */
+ else if (arg->nodata[0][0][0] && !arg->nodata[1][0][0])
+ *value = 1;
+ /* has value, increment */
+ else
+ *value = arg->values[0][0][0] + 1;
+ break;
+ case UT_SUM:
+ *value = arg->values[0][0][0] + arg->values[1][0][0];
+ break;
+ case UT_MEAN:
+ break;
+ case UT_LAST:
+ default:
+ *value = arg->values[1][0][0];
+ break;
+ }
+
+ POSTGIS_RT_DEBUGF(4, "value, nodata = %f, %d", *value, *nodata);
+
+
+ return 1;
+}
+
+static int rtpg_union_mean_callback(
+ rt_iterator_arg arg, void *userarg,
+ double *value, int *nodata
+) {
+ if (arg == NULL)
+ return 0;
+
+ if (
+ arg->rasters != 2 ||
+ arg->rows != 1 ||
+ arg->columns != 1
+ ) {
+ elog(ERROR, "rtpg_union_mean_callback: Invalid arguments passed to callback");
+ return 0;
+ }
+
+ *value = 0;
+ *nodata = 1;
+
+ POSTGIS_RT_DEBUGF(4, "rast0: %f %d", arg->values[0][0][0], arg->nodata[0][0][0]);
+ POSTGIS_RT_DEBUGF(4, "rast1: %f %d", arg->values[1][0][0], arg->nodata[1][0][0]);
+
+ if (
+ !arg->nodata[0][0][0] &&
+ FLT_NEQ(arg->values[0][0][0], 0) &&
+ !arg->nodata[1][0][0]
+ ) {
+ *value = arg->values[1][0][0] / arg->values[0][0][0];
+ *nodata = 0;
+ }
+
+ POSTGIS_RT_DEBUGF(4, "value, nodata = (%f, %d)", *value, *nodata);
+
+ return 1;
+}
+
+/* called for ST_Union(raster, unionarg[]) */
+static int rtpg_union_unionarg_process(rtpg_union_arg arg, ArrayType *array) {
+ Oid etype;
+ Datum *e;
+ bool *nulls;
+ int16 typlen;
+ bool typbyval;
+ char typalign;
+ int n = 0;
+
+ HeapTupleHeader tup;
+ bool isnull;
+ Datum tupv;
+
+ int i;
+ int nband = 1;
+ char *utypename = NULL;
+ rtpg_union_type utype = UT_LAST;
+
+ etype = ARR_ELEMTYPE(array);
+ get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
+
+ deconstruct_array(
+ array,
+ etype,
+ typlen, typbyval, typalign,
+ &e, &nulls, &n
+ );
+
+ if (!n) {
+ elog(ERROR, "rtpg_union_unionarg_process: Invalid argument for unionarg");
+ return 0;
+ }
+
+ /* prep arg */
+ arg->numband = n;
+ arg->bandarg = palloc(sizeof(struct rtpg_union_band_arg_t) * arg->numband);
+ if (arg->bandarg == NULL) {
+ elog(ERROR, "rtpg_union_unionarg_process: Unable to allocate memory for band information");
+ return 0;
+ }
+
+ /* process each element */
+ for (i = 0; i < n; i++) {
+ if (nulls[i]) {
+ arg->numband--;
+ continue;
+ }
+
+ POSTGIS_RT_DEBUGF(4, "Processing unionarg at index %d", i);
+
+ /* each element is a tuple */
+ tup = (HeapTupleHeader) DatumGetPointer(e[i]);
+ if (NULL == tup) {
+ elog(ERROR, "rtpg_union_unionarg_process: Invalid argument for unionarg");
+ return 0;
+ }
+
+ /* first element, bandnum */
+ tupv = GetAttributeByName(tup, "nband", &isnull);
+ if (isnull) {
+ nband = i + 1;
+ elog(NOTICE, "First argument (nband) of unionarg is NULL. Assuming nband = %d", nband);
+ }
+ else
+ nband = DatumGetInt32(tupv);
+
+ if (nband < 1) {
+ elog(ERROR, "rtpg_union_unionarg_process: Band number must be greater than zero (1-based)");
+ return 0;
+ }
+
+ /* second element, uniontype */
+ tupv = GetAttributeByName(tup, "uniontype", &isnull);
+ if (isnull) {
+ elog(NOTICE, "Second argument (uniontype) of unionarg is NULL. Assuming uniontype = LAST");
+ utype = UT_LAST;
+ }
+ else {
+ utypename = text_to_cstring((text *) DatumGetPointer(tupv));
+ utype = rtpg_uniontype_index_from_name(rtpg_strtoupper(utypename));
+ }
+
+ arg->bandarg[i].uniontype = utype;
+ arg->bandarg[i].nband = nband - 1;
+ arg->bandarg[i].raster = NULL;
+
+ if (utype != UT_MEAN)
+ arg->bandarg[i].numraster = 1;
+ else
+ arg->bandarg[i].numraster = 2;
+ }
+
+ if (arg->numband < n) {
+ arg->bandarg = repalloc(arg->bandarg, sizeof(struct rtpg_union_band_arg_t) * arg->numband);
+ if (arg->bandarg == NULL) {
+ elog(ERROR, "rtpg_union_unionarg_process: Unable to reallocate memory for band information");
+ return 0;
+ }
+ }
+
+ return 1;
+}
+
+/* called for ST_Union(raster) */
+static int rtpg_union_noarg(rtpg_union_arg arg, rt_raster raster) {
+ int numbands;
+ int i;
+
+ if (rt_raster_is_empty(raster))
+ return 1;
+
+ numbands = rt_raster_get_num_bands(raster);
+ if (numbands <= arg->numband)
+ return 1;
+
+ /* more bands to process */
+ POSTGIS_RT_DEBUG(4, "input raster has more bands, adding more bandargs");
+ if (arg->numband)
+ arg->bandarg = repalloc(arg->bandarg, sizeof(struct rtpg_union_band_arg_t) * numbands);
+ else
+ arg->bandarg = palloc(sizeof(struct rtpg_union_band_arg_t) * numbands);
+ if (arg->bandarg == NULL) {
+ elog(ERROR, "rtpg_union_noarg: Unable to reallocate memory for band information");
+ return 0;
+ }
+
+ i = arg->numband;
+ arg->numband = numbands;
+ for (; i < arg->numband; i++) {
+ POSTGIS_RT_DEBUGF(4, "Adding bandarg for band at index %d", i);
+ arg->bandarg[i].uniontype = UT_LAST;
+ arg->bandarg[i].nband = i;
+ arg->bandarg[i].numraster = 1;
+
+ arg->bandarg[i].raster = (rt_raster *) palloc(sizeof(rt_raster) * arg->bandarg[i].numraster);
+ if (arg->bandarg[i].raster == NULL) {
+ elog(ERROR, "rtpg_union_noarg: Unable to allocate memory for working rasters");
+ return 0;
+ }
+ memset(arg->bandarg[i].raster, 0, sizeof(rt_raster) * arg->bandarg[i].numraster);
+
+ /* add new working rt_raster but only if working raster already exists */
+ if (!rt_raster_is_empty(arg->bandarg[0].raster[0])) {
+ arg->bandarg[i].raster[0] = rt_raster_clone(arg->bandarg[0].raster[0], 0);
+ if (arg->bandarg[i].raster[0] == NULL) {
+ elog(ERROR, "rtpg_union_noarg: Unable to create working raster");
+ return 0;
+ }
+ }
+ }
+
+ return 1;
+}
/* UNION aggregate transition function */
PG_FUNCTION_INFO_V1(RASTER_union_transfn);
iwr->bandarg[0].raster = NULL;
break;
/* only other type allowed is unionarg */
- default: {
- ArrayType *array;
- Oid etype;
- Datum *e;
- bool *nulls;
- int16 typlen;
- bool typbyval;
- char typalign;
- int n = 0;
-
- HeapTupleHeader tup;
- bool isnull;
- Datum tupv;
-
+ default:
POSTGIS_RT_DEBUG(4, "Processing arg 3 as unionarg");
+ if (!rtpg_union_unionarg_process(iwr, PG_GETARG_ARRAYTYPE_P(2))) {
+ elog(ERROR, "RASTER_union_transfn: Unable to process unionarg");
- array = PG_GETARG_ARRAYTYPE_P(2);
- etype = ARR_ELEMTYPE(array);
- get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
-
- deconstruct_array(
- array,
- etype,
- typlen, typbyval, typalign,
- &e, &nulls, &n
- );
-
- if (!n) {
- elog(ERROR, "RASTER_union_transfn: Invalid argument for unionarg");
rtpg_union_arg_destroy(iwr);
if (raster != NULL) {
rt_raster_destroy(raster);
PG_RETURN_NULL();
}
- /* prep iwr */
- iwr->numband = n;
- iwr->bandarg = palloc(sizeof(struct rtpg_union_band_arg_t) * iwr->numband);
- if (iwr->bandarg == NULL) {
- elog(ERROR, "RASTER_union_transfn: Unable to allocate memory for band information");
-
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
- }
-
- /* process each element */
- for (i = 0; i < n; i++) {
- if (nulls[i]) {
- iwr->numband--;
- continue;
- }
-
- POSTGIS_RT_DEBUGF(4, "Processing unionarg at index %d", i);
-
- /* each element is a tuple */
- tup = (HeapTupleHeader) DatumGetPointer(e[i]);
- if (NULL == tup) {
- elog(ERROR, "RASTER_union_transfn: Invalid argument for unionarg");
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
- }
-
- /* first element, bandnum */
- tupv = GetAttributeByName(tup, "nband", &isnull);
- if (isnull) {
- nband = i + 1;
- elog(NOTICE, "First argument (nband) of unionarg is NULL. Assuming nband = %d", nband);
- }
- else
- nband = DatumGetInt32(tupv);
-
- if (nband < 1) {
- elog(ERROR, "RASTER_union_transfn: Band number must be greater than zero (1-based)");
-
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
- }
-
- /* second element, uniontype */
- tupv = GetAttributeByName(tup, "uniontype", &isnull);
- if (isnull) {
- elog(NOTICE, "Second argument (uniontype) of unionarg is NULL. Assuming uniontype = LAST");
- utype = UT_LAST;
- }
- else {
- utypename = text_to_cstring((text *) DatumGetPointer(tupv));
- utype = rtpg_uniontype_index_from_name(rtpg_strtoupper(utypename));
- }
-
- iwr->bandarg[i].uniontype = utype;
- iwr->bandarg[i].nband = nband - 1;
- iwr->bandarg[i].raster = NULL;
-
- if (utype != UT_MEAN)
- iwr->bandarg[i].numraster = 1;
- else
- iwr->bandarg[i].numraster = 2;
- }
-
- if (iwr->numband < n) {
- iwr->bandarg = repalloc(iwr->bandarg, sizeof(struct rtpg_union_band_arg_t) * iwr->numband);
- if (iwr->bandarg == NULL) {
- elog(ERROR, "RASTER_union_transfn: Unable to reallocate memory for band information");
-
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
- }
- }
-
break;
- }
}
}
/* only do this if raster isn't empty */
else if (nargs < 3) {
POSTGIS_RT_DEBUG(4, "no additional args, checking input raster");
+ if (!rtpg_union_noarg(iwr, raster)) {
+ elog(ERROR, "RASTER_union_transfn: Unable to check and balance number of bands");
- do {
- int numbands;
-
- if (rt_raster_is_empty(raster))
- break;
-
- numbands = rt_raster_get_num_bands(raster);
- if (numbands <= iwr->numband)
- break;
-
- /* more bands to process */
- POSTGIS_RT_DEBUG(4, "input raster has more bands, adding more bandargs");
- if (iwr->numband)
- iwr->bandarg = repalloc(iwr->bandarg, sizeof(struct rtpg_union_band_arg_t) * numbands);
- else
- iwr->bandarg = palloc(sizeof(struct rtpg_union_band_arg_t) * numbands);
- if (iwr->bandarg == NULL) {
- elog(ERROR, "RASTER_union_transfn: Unable to reallocate memory for band information");
-
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
+ rtpg_union_arg_destroy(iwr);
+ if (raster != NULL) {
+ rt_raster_destroy(raster);
+ PG_FREE_IF_COPY(pgraster, 1);
}
- i = iwr->numband;
- iwr->numband = numbands;
- for (; i < iwr->numband; i++) {
- POSTGIS_RT_DEBUGF(4, "Adding bandarg for band at index %d", i);
- iwr->bandarg[i].uniontype = UT_LAST;
- iwr->bandarg[i].nband = i;
- iwr->bandarg[i].numraster = 1;
-
- iwr->bandarg[i].raster = (rt_raster *) palloc(sizeof(rt_raster) * iwr->bandarg[i].numraster);
- if (iwr->bandarg[i].raster == NULL) {
- elog(ERROR, "RASTER_union_transfn: Unable to allocate memory for working rasters");
-
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
- }
- memset(iwr->bandarg[i].raster, 0, sizeof(rt_raster) * iwr->bandarg[i].numraster);
-
- /* add new working rt_raster but only if working raster already exists */
- if (!rt_raster_is_empty(iwr->bandarg[0].raster[0])) {
- iwr->bandarg[i].raster[0] = rt_raster_clone(iwr->bandarg[0].raster[0], 0);
- if (iwr->bandarg[i].raster[0] == NULL) {
- elog(ERROR, "RASTER_union_transfn: Unable to create working raster");
-
- rtpg_union_arg_destroy(iwr);
- if (raster != NULL) {
- rt_raster_destroy(raster);
- PG_FREE_IF_COPY(pgraster, 1);
- }
-
- MemoryContextSwitchTo(oldcontext);
- PG_RETURN_NULL();
- }
- }
-
- }
+ MemoryContextSwitchTo(oldcontext);
+ PG_RETURN_NULL();
}
- while (0);
}
/* init itrset */
}
if (_band != NULL) {
pixtype = rt_band_get_pixtype(_band);
- hasnodata = rt_band_get_hasnodata_flag(_band);
- if (hasnodata)
+ hasnodata = 1;
+ if (rt_band_get_hasnodata_flag(_band))
nodataval = rt_band_get_nodata(_band);
else
nodataval = rt_band_get_min_value(_band);
AS 'MODULE_PATHNAME', 'RASTER_gdal_version'
LANGUAGE 'c' IMMUTABLE;
+-----------------------------------------------------------------------
+-- generic composite type of a raster and its band index
+-----------------------------------------------------------------------
+
+CREATE TYPE rastbandarg AS (
+ rast raster,
+ nband integer
+);
+
-----------------------------------------------------------------------
-- Raster Accessors
-----------------------------------------------------------------------
END;
$$ LANGUAGE 'plpgsql' IMMUTABLE;
+-----------------------------------------------------------------------
+-- n-Raster ST_MapAlgebra
+-----------------------------------------------------------------------
+
+CREATE OR REPLACE FUNCTION _st_mapalgebra(
+ rastbandargset rastbandarg[],
+ callbackfunc regprocedure,
+ pixeltype text DEFAULT NULL,
+ distancex integer DEFAULT 0, distancey integer DEFAULT 0,
+ extenttype text DEFAULT 'INTERSECTION', customextent raster DEFAULT NULL,
+ VARIADIC userargs text[] DEFAULT NULL
+)
+ RETURNS raster
+ AS 'MODULE_PATHNAME', 'RASTER_nMapAlgebra'
+ LANGUAGE 'c' STABLE;
+
+CREATE OR REPLACE FUNCTION st_mapalgebra(
+ rastbandargset rastbandarg[],
+ callbackfunc regprocedure,
+ pixeltype text DEFAULT NULL,
+ extenttype text DEFAULT 'INTERSECTION', customextent raster DEFAULT NULL,
+ distancex integer DEFAULT 0, distancey integer DEFAULT 0,
+ VARIADIC userargs text[] DEFAULT NULL
+)
+ RETURNS raster
+ AS $$ SELECT _ST_MapAlgebra($1, $2, $3, $6, $7, $4, $5, VARIADIC $8) $$
+ LANGUAGE 'sql' STABLE;
+
+CREATE OR REPLACE FUNCTION st_mapalgebra(
+ rast raster, nband int[],
+ callbackfunc regprocedure,
+ pixeltype text DEFAULT NULL,
+ extenttype text DEFAULT 'FIRST', customextent raster DEFAULT NULL,
+ distancex integer DEFAULT 0, distancey integer DEFAULT 0,
+ VARIADIC userargs text[] DEFAULT NULL
+)
+ RETURNS raster
+ AS $$
+ DECLARE
+ x int;
+ argset rastbandarg[];
+ BEGIN
+ IF $2 IS NULL OR array_ndims($2) < 1 OR array_length($2, 1) < 1 THEN
+ RAISE EXCEPTION 'Populated 1D array must be provided for nband';
+ RETURN NULL;
+ END IF;
+
+ FOR x IN array_lower($2, 1)..array_upper($2, 1) LOOP
+ IF $2[x] IS NULL THEN
+ CONTINUE;
+ END IF;
+
+ argset := argset || ROW($1, $2[x])::rastbandarg;
+ END LOOP;
+
+ IF array_length(argset, 1) < 1 THEN
+ RAISE EXCEPTION 'Populated 1D array must be provided for nband';
+ RETURN NULL;
+ END IF;
+
+ RETURN _ST_MapAlgebra(argset, $3, $4, $7, $8, $5, $6, VARIADIC $9);
+ END;
+ $$ LANGUAGE 'plpgsql' STABLE;
+
+CREATE OR REPLACE FUNCTION st_mapalgebra(
+ rast raster, nband int,
+ callbackfunc regprocedure,
+ pixeltype text DEFAULT NULL,
+ extenttype text DEFAULT 'FIRST', customextent raster DEFAULT NULL,
+ distancex integer DEFAULT 0, distancey integer DEFAULT 0,
+ VARIADIC userargs text[] DEFAULT NULL
+)
+ RETURNS raster
+ AS $$ SELECT _ST_MapAlgebra(ARRAY[ROW($1, $2)]::rastbandarg[], $3, $4, $7, $8, $5, $6, VARIADIC $9) $$
+ LANGUAGE 'sql' STABLE;
+
+CREATE OR REPLACE FUNCTION st_mapalgebra(
+ rast1 raster, nband1 int,
+ rast2 raster, nband2 int,
+ callbackfunc regprocedure,
+ pixeltype text DEFAULT NULL,
+ extenttype text DEFAULT 'INTERSECTION', customextent raster DEFAULT NULL,
+ distancex integer DEFAULT 0, distancey integer DEFAULT 0,
+ VARIADIC userargs text[] DEFAULT NULL
+)
+ RETURNS raster
+ AS $$ SELECT _ST_MapAlgebra(ARRAY[ROW($1, $2), ROW($3, $4)]::rastbandarg[], $5, $6, $9, $10, $7, $8, VARIADIC $11) $$
+ LANGUAGE 'sql' STABLE;
+
-----------------------------------------------------------------------
-- Get information about the raster
-----------------------------------------------------------------------
(nband integer,
uniontype text);
END IF;
+
+ -- create rastbandarg type if it does not exist
+ IF NOT EXISTS(SELECT typname
+ FROM pg_type
+ WHERE typname = 'rastbandarg') THEN
+ CREATE TYPE rastbandarg AS (
+ rast raster,
+ nband integer
+ );
+ END IF;
+
END$$;
-- make geometry cast ASSIGNMENT
rtn = rt_raster_iterator(
itrset, 1,
- ET_LAST, NULL,
+ ET_INTERSECTION, NULL,
+ PT_32BUI,
+ 1, 0,
+ 0, 0,
+ userargs,
+ testRasterIterator1_callback,
+ &noerr
+ );
+ CHECK(noerr);
+ CHECK((rt_raster_get_width(rtn) == 5));
+ CHECK((rt_raster_get_height(rtn) == 5));
+ CHECK((rt_raster_get_x_offset(rtn) == 0));
+ CHECK((rt_raster_get_y_offset(rtn) == 0));
+ CHECK((rt_raster_get_x_scale(rtn) == 1));
+ CHECK((rt_raster_get_y_scale(rtn) == -1));
+ CHECK((rt_raster_get_x_skew(rtn) == 0));
+ CHECK((rt_raster_get_y_skew(rtn) == 0));
+ CHECK((rt_raster_get_srid(rtn) == 0));
+
+ if (rtn != NULL) deepRelease(rtn);
+ rtn = NULL;
+
+ /* 1 raster, 0 distance, FIRST or SECOND or LAST or UNION or INTERSECTION */
+ userargs->rasters = 1;
+ userargs->rows = 1;
+ userargs->columns = 1;
+
+ rtn = rt_raster_iterator(
+ itrset, 1,
+ ET_UNION, NULL,
PT_32BUI,
1, 0,
0, 0,
rt_mapalgebrafctngb_userfunc \
rt_intersection \
rt_clip \
+ rt_mapalgebra \
rt_union \
rt_invdistweight4ma
--- /dev/null
+SET client_min_messages TO warning;
+
+DROP TABLE IF EXISTS raster_nmapalgebra_in;
+CREATE TABLE raster_nmapalgebra_in (
+ rid integer,
+ rast raster
+);
+
+INSERT INTO raster_nmapalgebra_in
+ SELECT 0, NULL::raster AS rast UNION ALL
+ SELECT 1, ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0) AS rast UNION ALL
+ SELECT 2, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '8BUI', 1, 0) AS rast UNION ALL
+ SELECT 3, ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 1, -1, 1, -1, 0, 0, 0), 1, '8BUI', 2, 0), 2, '32BF', 20, 0) AS rast UNION ALL
+ SELECT 4, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 1, -1, 1, -1, 0, 0, 0), 1, '8BUI', 2, 0), 2, '32BF', 20, 0), 3, '16BUI', 200, 0) AS rast
+;
+
+CREATE OR REPLACE FUNCTION raster_nmapalgebra_test(
+ value double precision[][][],
+ pos int[][],
+ VARIADIC userargs text[]
+)
+ RETURNS double precision
+ AS $$
+ BEGIN
+ RAISE NOTICE 'value = %', value;
+ RAISE NOTICE 'pos = %', pos;
+ RAISE NOTICE 'userargs = %', userargs;
+
+ IF userargs IS NULL OR array_length(userargs, 1) < 1 THEN
+ RETURN 255;
+ ELSE
+ RETURN userargs[array_lower(userargs, 1)];
+ END IF;
+ END;
+ $$ LANGUAGE 'plpgsql' IMMUTABLE;
+
+SET client_min_messages TO notice;
+
+SELECT
+ rid,
+ ST_Value(
+ ST_MapAlgebra(
+ ARRAY[ROW(rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ),
+ 1, 1, 1
+ ) = 255
+FROM raster_nmapalgebra_in
+WHERE rid IN (0, 1);
+
+SELECT
+ rid,
+ ST_Value(
+ ST_MapAlgebra(
+ ARRAY[ROW(rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ),
+ 1, 1, 1
+ ) = 255
+FROM raster_nmapalgebra_in
+WHERE rid IN (2,3,4);
+
+SELECT
+ rid,
+ round(ST_Value(
+ ST_MapAlgebra(
+ ARRAY[ROW(rast, 2)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'INTERSECTION', NULL,
+ 0, 0,
+ '3.14'
+ ),
+ 1, 1, 1
+ )::numeric, 2) = 3.14
+FROM raster_nmapalgebra_in
+WHERE rid IN (3,4);
+
+WITH foo AS (
+ SELECT
+ rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(rast, 3)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ '8BUI',
+ 'INTERSECTION', NULL,
+ 1, 1
+ ) AS rast
+ FROM raster_nmapalgebra_in
+ WHERE rid IN (3,4)
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1)),
+ ST_Value(rast, 1, 1, 1)
+FROM foo;
+
+INSERT INTO raster_nmapalgebra_in
+ SELECT 10, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0) AS rast UNION ALL
+ SELECT 11, ST_AddBand(ST_MakeEmptyRaster(2, 2, 2, 0, 1, -1, 0, 0, 0), 1, '16BUI', 2, 0) AS rast UNION ALL
+ SELECT 12, ST_AddBand(ST_MakeEmptyRaster(2, 2, 4, 0, 1, -1, 0, 0, 0), 1, '16BUI', 3, 0) AS rast UNION ALL
+
+ SELECT 13, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, -2, 1, -1, 0, 0, 0), 1, '16BUI', 10, 0) AS rast UNION ALL
+ SELECT 14, ST_AddBand(ST_MakeEmptyRaster(2, 2, 2, -2, 1, -1, 0, 0, 0), 1, '16BUI', 20, 0) AS rast UNION ALL
+ SELECT 15, ST_AddBand(ST_MakeEmptyRaster(2, 2, 4, -2, 1, -1, 0, 0, 0), 1, '16BUI', 30, 0) AS rast UNION ALL
+
+ SELECT 16, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, -4, 1, -1, 0, 0, 0), 1, '16BUI', 100, 0) AS rast UNION ALL
+ SELECT 17, ST_AddBand(ST_MakeEmptyRaster(2, 2, 2, -4, 1, -1, 0, 0, 0), 1, '16BUI', 200, 0) AS rast UNION ALL
+ SELECT 18, ST_AddBand(ST_MakeEmptyRaster(2, 2, 4, -4, 1, -1, 0, 0, 0), 1, '16BUI', 300, 0) AS rast
+;
+
+WITH foo AS (
+ SELECT
+ t1.rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(ST_Union(t2.rast), 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ '32BUI',
+ 'CUSTOM', t1.rast,
+ 1, 1
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 10
+ AND t2.rid BETWEEN 10 AND 18
+ AND ST_Intersects(t1.rast, t2.rast)
+ GROUP BY t1.rid, t1.rast
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1)),
+ ST_Value(rast, 1, 1, 1)
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(ST_Union(t2.rast), 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ '32BUI',
+ 'CUSTOM', t1.rast,
+ 1, 1
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 14
+ AND t2.rid BETWEEN 10 AND 18
+ AND ST_Intersects(t1.rast, t2.rast)
+ GROUP BY t1.rid, t1.rast
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1)),
+ ST_Value(rast, 1, 1, 1)
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(ST_Union(t2.rast), 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ '32BUI',
+ 'CUSTOM', t1.rast,
+ 1, 1,
+ '1000'
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 17
+ AND t2.rid BETWEEN 10 AND 18
+ AND ST_Intersects(t1.rast, t2.rast)
+ GROUP BY t1.rid, t1.rast
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1)),
+ ST_Value(rast, 1, 1, 1)
+FROM foo;
+
+INSERT INTO raster_nmapalgebra_in
+ SELECT 20, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0) AS rast UNION ALL
+ SELECT 21, ST_AddBand(ST_MakeEmptyRaster(2, 2, 1, -1, 1, -1, 0, 0, 0), 1, '16BUI', 2, 0) AS rast UNION ALL
+ SELECT 22, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, -2, 1, -1, 0, 0, 0), 1, '16BUI', 3, 0) AS rast
+;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 20
+ AND t2.rid = 21
+)
+SELECT
+ rid1,
+ rid2,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 20
+ AND t2.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 21
+ AND t2.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'UNION', NULL,
+ 0, 0
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 20
+ AND t2.rid = 21
+)
+SELECT
+ rid1,
+ rid2,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'UNION', NULL,
+ 0, 0
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 20
+ AND t2.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ t3.rid AS rid3,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1), ROW(t3.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'UNION', NULL,
+ 0, 0
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ CROSS JOIN raster_nmapalgebra_in t3
+ WHERE t1.rid = 20
+ AND t2.rid = 21
+ AND t3.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ rid3,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ t3.rid AS rid3,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1), ROW(t3.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'FIRST', NULL,
+ 0, 0
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ CROSS JOIN raster_nmapalgebra_in t3
+ WHERE t1.rid = 20
+ AND t2.rid = 21
+ AND t3.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ rid3,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ t3.rid AS rid3,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1), ROW(t3.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'SECOND', NULL,
+ 0, 0
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ CROSS JOIN raster_nmapalgebra_in t3
+ WHERE t1.rid = 20
+ AND t2.rid = 21
+ AND t3.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ rid3,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ t3.rid AS rid3,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t2.rast, 1), ROW(t3.rast, 1)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ NULL,
+ 'LAST', NULL,
+ 0, 0
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ CROSS JOIN raster_nmapalgebra_in t3
+ WHERE t1.rid = 20
+ AND t2.rid = 21
+ AND t3.rid = 22
+)
+SELECT
+ rid1,
+ rid2,
+ rid3,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+INSERT INTO raster_nmapalgebra_in
+ SELECT 30, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '16BUI', 1, 0), 2, '8BUI', 10, 0), 3, '32BUI', 100, 0) AS rast UNION ALL
+ SELECT 31, ST_AddBand(ST_AddBand(ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 1, 1, -1, 0, 0, 0), 1, '16BUI', 2, 0), 2, '8BUI', 20, 0), 3, '32BUI', 300, 0) AS rast
+;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 1), ROW(t1.rast, 2), ROW(t1.rast, 3)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ WHERE t1.rid = 30
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 3), ROW(t1.rast, 1), ROW(t1.rast, 3)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ WHERE t1.rid = 30
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 2), ROW(t1.rast, 2)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ '16BUI'::text
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ WHERE t1.rid = 31
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid1,
+ t2.rid AS rid2,
+ ST_MapAlgebra(
+ ARRAY[ROW(t1.rast, 2), ROW(t2.rast, 1), ROW(t2.rast, 2)]::rastbandarg[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure,
+ '16BUI'
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ CROSS JOIN raster_nmapalgebra_in t2
+ WHERE t1.rid = 30
+ AND t2.rid = 31
+)
+SELECT
+ rid1,
+ rid2,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid,
+ ST_MapAlgebra(
+ t1.rast, ARRAY[3, 1, 3]::int[],
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ WHERE t1.rid = 30
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+WITH foo AS (
+ SELECT
+ t1.rid AS rid,
+ ST_MapAlgebra(
+ t1.rast, 2,
+ 'raster_nmapalgebra_test(double precision[], int[], text[])'::regprocedure
+ ) AS rast
+ FROM raster_nmapalgebra_in t1
+ WHERE t1.rid = 30
+)
+SELECT
+ rid,
+ (ST_Metadata(rast)),
+ (ST_BandMetadata(rast, 1))
+FROM foo;
+
+DROP FUNCTION IF EXISTS raster_nmapalgebra_test(double precision[], int[], text[]);
+DROP TABLE IF EXISTS raster_nmapalgebra_in;
--- /dev/null
+NOTICE: First argument (nband) of rastbandarg at index 0 is NULL. Assuming NULL raster
+NOTICE: All input rasters are NULL. Returning NULL
+NOTICE: All input rasters do not have bands at indicated indexes. Returning empty raster
+NOTICE: Could not find raster band of index 1 when getting pixel value. Returning NULL
+0|
+1|
+NOTICE: value = {{{1}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+2|t
+3|t
+4|t
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = {3.14}
+NOTICE: value = {{{20}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = {3.14}
+3|t
+4|t
+NOTICE: All input rasters do not have bands at indicated indexes. Returning empty raster
+NOTICE: Raster provided has no bands
+NOTICE: Could not find raster band of index 1 when getting pixel value. Returning NULL
+NOTICE: value = {{{NULL,NULL,NULL},{NULL,200,200},{NULL,200,200}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL,NULL,NULL},{200,200,NULL},{200,200,NULL}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL,200,200},{NULL,200,200},{NULL,NULL,NULL}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{200,200,NULL},{200,200,NULL},{NULL,NULL,NULL}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+3|(0,0,0,0,1,1,0,0,0,0)|(,,,)|
+4|(1,-1,2,2,1,-1,0,0,0,1)|(8BUI,0,f,)|255
+NOTICE: value = {{{NULL,NULL,NULL},{NULL,1,1},{NULL,1,1}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL,NULL,NULL},{1,1,2},{1,1,2}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL,1,1},{NULL,1,1},{NULL,10,10}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1,1,2},{1,1,2},{10,10,20}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+10|(0,0,2,2,1,-1,0,0,0,1)|(32BUI,0,f,)|255
+NOTICE: value = {{{1,2,2},{10,20,20},{10,20,20}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{2,2,3},{20,20,30},{20,20,30}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{10,20,20},{10,20,20},{100,200,200}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{20,20,30},{20,20,30},{200,200,300}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+14|(2,-2,2,2,1,-1,0,0,0,1)|(32BUI,0,f,)|255
+NOTICE: value = {{{10,20,20},{100,200,200},{100,200,200}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = {1000}
+NOTICE: value = {{{20,20,30},{200,200,300},{200,200,300}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = {1000}
+NOTICE: value = {{{100,200,200},{100,200,200},{NULL,NULL,NULL}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = {1000}
+NOTICE: value = {{{200,200,300},{200,200,300},{NULL,NULL,NULL}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = {1000}
+17|(2,-4,2,2,1,-1,0,0,0,1)|(32BUI,0,f,)|1000
+NOTICE: value = {{{1}},{{2}}}
+NOTICE: pos = [0:1][1:2]={{1,1},{2,2}}
+NOTICE: userargs = <NULL>
+20|21|(1,-1,1,1,1,-1,0,0,0,1)|(16BUI,0,f,)
+ERROR: rt_raster_iterator: Computed raster for intersection extent is empty
+NOTICE: value = {{{2}},{{3}}}
+NOTICE: pos = [0:1][1:2]={{1,1},{1,2}}
+NOTICE: userargs = <NULL>
+21|22|(1,-2,1,1,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{1,1},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{2,1},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{3,1},{3,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{1,2},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{2}}}
+NOTICE: pos = [0:1][1:2]={{2,2},{2,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}}}
+NOTICE: pos = [0:1][1:2]={{3,2},{3,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{1,3},{1,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}}}
+NOTICE: pos = [0:1][1:2]={{2,3},{2,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}}}
+NOTICE: pos = [0:1][1:2]={{3,3},{3,3}}
+NOTICE: userargs = <NULL>
+20|21|(0,0,3,3,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{1,1},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{2,1},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{1,2},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}}}
+NOTICE: pos = [0:1][1:2]={{2,2},{2,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{3}}}
+NOTICE: pos = [0:1][1:2]={{1,3},{1,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{3}}}
+NOTICE: pos = [0:1][1:2]={{2,3},{2,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{3}}}
+NOTICE: pos = [0:1][1:2]={{1,4},{1,4}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{3}}}
+NOTICE: pos = [0:1][1:2]={{2,4},{2,4}}
+NOTICE: userargs = <NULL>
+20|22|(0,0,2,4,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{1}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,1},{0,0}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,1},{1,0}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{3,1},{3,1},{2,0}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{1,2},{0,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{2,2},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{3,2},{3,2},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{1,3},{1,3},{0,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{2,3},{2,3},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{3,3},{3,3},{2,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{1,4},{1,4},{0,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{2,4},{2,4},{1,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{3,4},{3,4},{2,3}}
+NOTICE: userargs = <NULL>
+20|21|22|(0,0,3,4,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{1}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,1},{0,0}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,1},{1,0}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{NULL}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{1,2},{0,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{2,2},{1,1}}
+NOTICE: userargs = <NULL>
+20|21|22|(0,0,2,2,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{1}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{2,2},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{3,2},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{2,3},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{NULL}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{3,3},{2,2}}
+NOTICE: userargs = <NULL>
+20|21|22|(1,-1,2,2,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{NULL}},{{NULL}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,3},{0,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{2}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,3},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{1,4},{0,3}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{NULL}},{{NULL}},{{3}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{2,4},{1,3}}
+NOTICE: userargs = <NULL>
+20|21|22|(0,-2,2,2,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{1}},{{10}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,1},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{10}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,1},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{10}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{1,2},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{1}},{{10}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{2,2},{2,2}}
+NOTICE: userargs = <NULL>
+30|(0,0,2,2,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,1},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,1},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{1,2},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{2,2},{2,2}}
+NOTICE: userargs = <NULL>
+30|(0,0,2,2,1,-1,0,0,0,1)|(32BUI,0,f,)
+NOTICE: value = {{{20}},{{20}}}
+NOTICE: pos = [0:1][1:2]={{1,1},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{20}},{{20}}}
+NOTICE: pos = [0:1][1:2]={{2,1},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{20}},{{20}}}
+NOTICE: pos = [0:1][1:2]={{1,2},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{20}},{{20}}}
+NOTICE: pos = [0:1][1:2]={{2,2},{2,2}}
+NOTICE: userargs = <NULL>
+31|(0,1,2,2,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{10}},{{2}},{{20}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,1},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{10}},{{2}},{{20}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,1},{2,2}}
+NOTICE: userargs = <NULL>
+30|31|(0,0,2,1,1,-1,0,0,0,1)|(16BUI,0,f,)
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{1,1},{1,1},{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{2,1},{2,1},{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{1,2},{1,2},{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{100}},{{1}},{{100}}}
+NOTICE: pos = [0:2][1:2]={{2,2},{2,2},{2,2}}
+NOTICE: userargs = <NULL>
+30|(0,0,2,2,1,-1,0,0,0,1)|(32BUI,0,f,)
+NOTICE: value = {{{10}}}
+NOTICE: pos = [0:0][1:2]={{1,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{10}}}
+NOTICE: pos = [0:0][1:2]={{2,1}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{10}}}
+NOTICE: pos = [0:0][1:2]={{1,2}}
+NOTICE: userargs = <NULL>
+NOTICE: value = {{{10}}}
+NOTICE: pos = [0:0][1:2]={{2,2}}
+NOTICE: userargs = <NULL>
+30|(0,0,2,2,1,-1,0,0,0,1)|(8BUI,0,f,)
SUM|3|3|2
LAST|1|1|1
LAST|2|1|1
-LAST|3|1|0
+LAST|3|1|
LAST|1|2|1
LAST|2|2|2
LAST|3|2|2
-LAST|1|3|0
+LAST|1|3|
LAST|2|3|2
LAST|3|3|2
COUNT|1|1|1