From b88337d4978812158e35dde58aec334a5fec5302 Mon Sep 17 00:00:00 2001 From: Bborie Park Date: Fri, 28 Sep 2012 23:09:31 +0000 Subject: [PATCH] Addition of C-based ST_Union(raster) aggregate function (ticket #1364). Renamed low level function rt_raster_has_no_band() to raster_has_band(). Updated docs and NEWS for ST_Union. git-svn-id: http://svn.osgeo.org/postgis/trunk@10344 b70326c6-7e19-0410-871a-916f4a2858ee --- NEWS | 1 + doc/reference_raster.xml | 52 +-- raster/rt_core/rt_api.c | 48 +-- raster/rt_core/rt_api.h | 8 +- raster/rt_pg/rt_pg.c | 551 +++++++++++++++++++++++++- raster/rt_pg/rtpostgis.sql.in.c | 156 ++------ raster/test/core/testapi.c | 26 +- raster/test/regress/rt_union.sql | 90 ++++- raster/test/regress/rt_union_expected | 77 +++- 9 files changed, 788 insertions(+), 221 deletions(-) diff --git a/NEWS b/NEWS index 55b42c4d6..d784b4261 100644 --- a/NEWS +++ b/NEWS @@ -33,6 +33,7 @@ PostGIS 2.1.0 * Enhancements * - #823, tiger geocoder: Make loader_generate_script download portion less greedy - #1363, ST_AddBand(raster, ...) array version rewritten in C + - #1364, ST_Union(raster, ...) aggregate function rewritten in C - #1661, Add aggregate variant of ST_SameAlignment - #1719, Add support for Point and GeometryCollection ST_MakeValid inputs - #1796, Big performance boost for distance calculations in geography diff --git a/doc/reference_raster.xml b/doc/reference_raster.xml index e2f2a4012..72e4e3d7c 100644 --- a/doc/reference_raster.xml +++ b/doc/reference_raster.xml @@ -8614,37 +8614,36 @@ UPDATE wind ST_Union - Returns the union of a set of raster tiles into a single raster composed of 1 band. If no band is specified for unioning, band num 1 is assumed. The resulting raster's extent is the extent of the whole set. In the case of intersection, the resulting value is defined by p_expression which is one of the following: LAST - the default when none is specified, MEAN, SUM, FIRST, MAX, MIN. + Returns the union of a set of raster tiles into a single raster composed of 1 band. - raster ST_Union - setof raster rast + raster ST_Union + setof raster rast - raster ST_Union - raster set rast - integer band_num + raster ST_Union + raster set rast + integer nband - raster ST_Union - raster set rast - text p_expression + raster ST_Union + raster set rast + text uniontype - raster ST_Union - raster set rast - integer band_num - text p_expression - + raster ST_Union + raster set rast + integer nband + text uniontype @@ -8652,19 +8651,16 @@ UPDATE wind Description - Returns the union of a set of raster tiles into a single raster composed of 1 band. If no band is specified for unioning, band num 1 is assumed. The resulting raster's extent is the extent of the whole set. In the case of intersection, the resulting value is defined by p_expression which is one of the following: LAST - the default when none is specified, MEAN, SUM, FIRST, MAX, MIN + Returns the union of a set of raster tiles into a single raster composed of 1 band. If nband is not specified, band 1 is assumed. The resulting raster's extent is the extent of the whole set. In the case of intersection, the resulting value is defined by uniontype which is one of the following: LAST (default), FIRST, MIN, MAX, COUNT, SUM, MEAN. - There are several other variants of this function not installed by default in PostGIS 2.0.0 -- these can be found in the raster/scripts/plpgsql/st_union.sql file of postgis source code. - - The ST_Union function in 2.0.0 is currently implemented predominantly in plpgsql. Because of the memory copying needed to copy between the C and plpgsql layer, this function is much much slower than it needs to be. - Future 2.0 releases will have this function implemented in C, so you should witness significant improvements in speed when that happens. As a general rule of thumb you want to minimize the size of the rasters, that ST_Union works with. - One approach is to clip first and then union the clipped versions. Refer to select parcels example in . That example if unioning is done before clipping takes about 4 times longer. With the higher res imagery the timing the ratio between is even higher. Availability: 2.0.0 + Enhanced: 2.1.0 Improved Speed (fully C-Based). Examples: Reconstitute a single band chunked raster tile - -- this creates a single band from first band of raster tiles + +-- this creates a single band from first band of raster tiles -- that form the original file system tile SELECT filename, ST_Union(rast) As file_rast FROM sometable WHERE filename IN('dem01', 'dem02') GROUP BY filename; @@ -8672,20 +8668,24 @@ FROM sometable WHERE filename IN('dem01', 'dem02') GROUP BY filename; Examples: Return a multi-band raster that is the union of tiles intersecting geometry - -- this creates a multi band raster collecting all the tiles that intersect a line + +-- this creates a multi band raster collecting all the tiles that intersect a line SELECT ST_AddBand(NULL,ARRAY[ST_Union(rast,1), ST_Union(rast,2), ST_Union(rast,3) ]) FROM aerials.boston WHERE ST_Intersects(rast, ST_GeomFromText('LINESTRING(230486 887771, 230500 88772)',26986) ); - - + See Also - , , , + + , + , + + - + Raster Processing Builtin Functions diff --git a/raster/rt_core/rt_api.c b/raster/rt_core/rt_api.c index 2f2aecfc8..6ef0aacbc 100644 --- a/raster/rt_core/rt_api.c +++ b/raster/rt_core/rt_api.c @@ -982,25 +982,25 @@ rt_pixtype_index_from_name(const char* pixname) { if (strcmp(pixname, "1BB") == 0) return PT_1BB; - if (strcmp(pixname, "2BUI") == 0) + else if (strcmp(pixname, "2BUI") == 0) return PT_2BUI; - if (strcmp(pixname, "4BUI") == 0) + else if (strcmp(pixname, "4BUI") == 0) return PT_4BUI; - if (strcmp(pixname, "8BSI") == 0) + else if (strcmp(pixname, "8BSI") == 0) return PT_8BSI; - if (strcmp(pixname, "8BUI") == 0) + else if (strcmp(pixname, "8BUI") == 0) return PT_8BUI; - if (strcmp(pixname, "16BSI") == 0) + else if (strcmp(pixname, "16BSI") == 0) return PT_16BSI; - if (strcmp(pixname, "16BUI") == 0) + else if (strcmp(pixname, "16BUI") == 0) return PT_16BUI; - if (strcmp(pixname, "32BSI") == 0) + else if (strcmp(pixname, "32BSI") == 0) return PT_32BSI; - if (strcmp(pixname, "32BUI") == 0) + else if (strcmp(pixname, "32BUI") == 0) return PT_32BUI; - if (strcmp(pixname, "32BF") == 0) + else if (strcmp(pixname, "32BF") == 0) return PT_32BF; - if (strcmp(pixname, "64BF") == 0) + else if (strcmp(pixname, "64BF") == 0) return PT_64BF; return PT_END; @@ -8088,16 +8088,16 @@ rt_raster_is_empty(rt_raster raster) { } /** - * Return TRUE if the raster do not have a band of this number. + * Return TRUE if the raster has a band of this number. * * @param raster: the raster to get info from * @param nband: the band number. 0-based * - * @return TRUE if the raster do not have a band of this number, FALSE otherwise + * @return TRUE if the raster has a band of this number, FALSE otherwise */ int -rt_raster_has_no_band(rt_raster raster, int nband) { - return (NULL == raster || nband >= raster->numBands || nband < 0); +rt_raster_has_band(rt_raster raster, int nband) { + return !(NULL == raster || nband >= raster->numBands || nband < 0); } /** @@ -12230,12 +12230,6 @@ rt_raster_from_two_rasters( *noerr = 0; - /* rasters must have same srid */ - if (rt_raster_get_srid(rast1) != rt_raster_get_srid(rast2)) { - rterror("rt_raster_from_two_rasters: The two rasters provided do not have the same SRID"); - return NULL; - } - /* rasters must be aligned */ if (!rt_raster_same_alignment(rast1, rast2, &aligned)) { rterror("rt_raster_from_two_rasters: Unable to test for alignment on the two rasters"); @@ -12984,8 +12978,8 @@ _rti_param_populate( _param->band[i] = NULL; _param->offset[i] = NULL; - /* set isnull and isempty */ - if (itrset[i].rast == NULL) { + /* set isempty */ + if (itrset[i].raster == NULL) { _param->isempty[i] = 1; (*allnull)++; @@ -12993,7 +12987,7 @@ _rti_param_populate( continue; } - else if (rt_raster_is_empty(itrset[i].rast)) { + else if (rt_raster_is_empty(itrset[i].raster)) { _param->isempty[i] = 1; (*allempty)++; @@ -13002,13 +12996,13 @@ _rti_param_populate( } /* check band number */ - if (rt_raster_has_no_band(itrset[i].rast, itrset[i].nband)) { + if (!rt_raster_has_band(itrset[i].raster, itrset[i].nband)) { rterror("_rti_param_populate: Band %d not found for raster %d", itrset[i].nband, i); return 0; } - _param->raster[i] = itrset[i].rast; - _param->band[i] = rt_raster_get_band(itrset[i].rast, itrset[i].nband); + _param->raster[i] = itrset[i].raster; + _param->band[i] = rt_raster_get_band(itrset[i].raster, itrset[i].nband); if (_param->band[i] == NULL) { rterror("_rti_param_populate: Unable to get band %d for raster %d", itrset[i].nband, i); return 0; @@ -13684,12 +13678,14 @@ rt_raster_iterator( } /* no pixels in neighborhood */ + /* if (!status) { RASTER_DEBUG(3, "no pixels in neighborhood, using empty"); _param->arg->values[i] = _param->empty.values; _param->arg->nodata[i] = _param->empty.nodata; continue; } + */ /* convert set of rt_pixel to 2D array */ status = rt_pixel_set_to_array( diff --git a/raster/rt_core/rt_api.h b/raster/rt_core/rt_api.h index f134d117a..7d354b663 100644 --- a/raster/rt_core/rt_api.h +++ b/raster/rt_core/rt_api.h @@ -1300,14 +1300,14 @@ rt_raster rt_raster_deserialize(void* serialized, int header_only); int rt_raster_is_empty(rt_raster raster); /** - * Return TRUE if the raster do not have a band of this number. + * Return TRUE if the raster has a band of this number. * * @param raster: the raster to get info from * @param nband: the band number. 0-based * - * @return TRUE if the raster do not have a band of this number, FALSE otherwise + * @return TRUE if the raster has a band of this number, FALSE otherwise */ -int rt_raster_has_no_band(rt_raster raster, int nband); +int rt_raster_has_band(rt_raster raster, int nband); /** * Copy one band from one raster to another. Bands are duplicated from @@ -2138,7 +2138,7 @@ struct rt_reclassexpr_t { /* raster iterator */ struct rt_iterator_t { - rt_raster rast; + rt_raster raster; uint16_t nband; /* 0-based */ }; diff --git a/raster/rt_pg/rt_pg.c b/raster/rt_pg/rt_pg.c index 738c3af79..89781714a 100644 --- a/raster/rt_pg/rt_pg.c +++ b/raster/rt_pg/rt_pg.c @@ -352,6 +352,10 @@ Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS); /* one-raster neighborhood MapAlgebra */ Datum RASTER_mapAlgebraFctNgb(PG_FUNCTION_ARGS); +/* raster union aggregate */ +Datum RASTER_union_transfn(PG_FUNCTION_ARGS); +Datum RASTER_union_finalfn(PG_FUNCTION_ARGS); + /* string replacement function taken from * http://ubuntuforums.org/showthread.php?s=aa6f015109fd7e4c7e30d2fd8b717497&t=141670&page=3 */ @@ -4574,7 +4578,7 @@ Datum RASTER_hasNoBand(PG_FUNCTION_ARGS) /* Get band number */ bandindex = PG_GETARG_INT32(1); - hasnoband = rt_raster_has_no_band(raster, bandindex - 1); + hasnoband = !rt_raster_has_band(raster, bandindex - 1); rt_raster_destroy(raster); PG_FREE_IF_COPY(pgraster, 0); @@ -4718,7 +4722,7 @@ Datum RASTER_mapAlgebraExpr(PG_FUNCTION_ARGS) * Check if the raster has the required band. Otherwise, return a raster * without band **/ - if (rt_raster_has_no_band(raster, nband - 1)) { + if (!rt_raster_has_band(raster, nband - 1)) { elog(NOTICE, "Raster does not have the required band. Returning a raster " "without a band"); rt_raster_destroy(raster); @@ -5379,7 +5383,7 @@ Datum RASTER_mapAlgebraFct(PG_FUNCTION_ARGS) * Check if the raster has the required band. Otherwise, return a raster * without band **/ - if (rt_raster_has_no_band(raster, nband - 1)) { + if (!rt_raster_has_band(raster, nband - 1)) { elog(NOTICE, "Raster does not have the required band. Returning a raster " "without a band"); rt_raster_destroy(raster); @@ -12011,6 +12015,7 @@ Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS) } /* SRID must match */ + /* if (rt_raster_get_srid(_rast[0]) != rt_raster_get_srid(_rast[1])) { elog(NOTICE, "The two rasters provided have different SRIDs. Returning NULL"); for (k = 0; k < set_count; k++) { @@ -12021,6 +12026,7 @@ Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS) } PG_RETURN_NULL(); } + */ /* same alignment */ if (!rt_raster_same_alignment(_rast[0], _rast[1], &aligned)) { @@ -12122,7 +12128,7 @@ Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS) } /* specified band not found */ - if (rt_raster_has_no_band(_rast[i], bandindex[i] - 1)) { + if (!rt_raster_has_band(_rast[i], bandindex[i] - 1)) { elog(NOTICE, "The %s raster does not have the band at index %d. Returning no band raster of correct extent", (i != 1 ? "FIRST" : "SECOND"), bandindex[i] ); @@ -12187,12 +12193,23 @@ Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS) PG_RETURN_POINTER(pgrtn); } break; + case ET_LAST: + case ET_CUSTOM: + elog(ERROR, "RASTER_mapAlgebra2: ET_LAST and ET_CUSTOM are not implemented"); + for (k = 0; k < set_count; k++) { + if (_rast[k] != NULL) + rt_raster_destroy(_rast[k]); + if (pgrastpos[k] != -1) + PG_FREE_IF_COPY(pgrast[k], pgrastpos[k]); + } + PG_RETURN_NULL(); + break; } /* both rasters do not have specified bands */ if ( - (!_isempty[0] && rt_raster_has_no_band(_rast[0], bandindex[0] - 1)) && - (!_isempty[1] && rt_raster_has_no_band(_rast[1], bandindex[1] - 1)) + (!_isempty[0] && !rt_raster_has_band(_rast[0], bandindex[0] - 1)) && + (!_isempty[1] && !rt_raster_has_band(_rast[1], bandindex[1] - 1)) ) { elog(NOTICE, "The two rasters provided do not have the respectively specified band indices. Returning no band raster of correct extent"); @@ -12213,7 +12230,7 @@ Datum RASTER_mapAlgebra2(PG_FUNCTION_ARGS) /* get bands */ for (i = 0; i < set_count; i++) { - if (_isempty[i] || rt_raster_has_no_band(_rast[i], bandindex[i] - 1)) { + if (_isempty[i] || !rt_raster_has_band(_rast[i], bandindex[i] - 1)) { _hasnodata[i] = 1; _nodataval[i] = 0; @@ -13058,7 +13075,7 @@ Datum RASTER_mapAlgebraFctNgb(PG_FUNCTION_ARGS) * Check if the raster has the required band. Otherwise, return a raster * without band **/ - if (rt_raster_has_no_band(raster, nband - 1)) { + if (!rt_raster_has_band(raster, nband - 1)) { elog(NOTICE, "Raster does not have the required band. Returning a raster " "without a band"); rt_raster_destroy(raster); @@ -13522,6 +13539,524 @@ Datum RASTER_mapAlgebraFctNgb(PG_FUNCTION_ARGS) PG_RETURN_POINTER(pgrtn); } +/* ---------------------------------------------------------------- */ +/* ST_Union aggregate functions */ +/* ---------------------------------------------------------------- */ + +typedef enum { + UT_LAST = 0, + UT_FIRST, + UT_MIN, + UT_MAX, + UT_COUNT, + UT_SUM, + UT_MEAN +} rtpg_union_type; + +typedef struct { + int numraster; + rt_raster *raster; + + rtpg_union_type uniontype; +} rtpg_union_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); + + 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; +} + +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; + return 1; + } + /* second NODATA */ + else if (!arg->nodata[0][0][0] && arg->nodata[1][0][0]) { + *value = arg->values[0][0][0]; + return 1; + } + /* first NODATA */ + else if (arg->nodata[0][0][0] && !arg->nodata[1][0][0]) { + *value = arg->values[1][0][0]; + 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; + } + + 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; + + 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; + } + + return 1; +} + +/* UNION aggregate transition function */ +PG_FUNCTION_INFO_V1(RASTER_union_transfn); +Datum RASTER_union_transfn(PG_FUNCTION_ARGS) +{ + MemoryContext aggcontext; + MemoryContext oldcontext; + rtpg_union_arg *iwr; + + rt_pgraster *pgraster = NULL; + rt_raster raster = NULL; + rt_raster _raster = NULL; + rt_band _band = NULL; + int nband = 1; + int noerr = 1; + int isempty[2] = {0}; + int hasband[2] = {0}; + + int i = 0; + + rt_iterator itrset; + char *utypename = NULL; + rtpg_union_type utype = UT_LAST; + rt_pixtype pixtype = PT_END; + int hasnodata = 1; + double nodataval = 0; + + /* cannot be called directly as this is exclusive aggregate function */ + if (!AggCheckCallContext(fcinfo, &aggcontext)) { + elog(ERROR, "RASTER_union_transfn: Cannot be called in a non-aggregate context"); + PG_RETURN_NULL(); + } + + /* switch to aggcontext */ + oldcontext = MemoryContextSwitchTo(aggcontext); + + if (PG_ARGISNULL(0)) { + /* allocate container in aggcontext */ + iwr = (rtpg_union_arg *) palloc(sizeof(rtpg_union_arg)); + + iwr->numraster = 1; + iwr->raster = NULL; + iwr->uniontype = UT_LAST; + } + else + iwr = (rtpg_union_arg *) PG_GETARG_POINTER(0); + + /* raster arg is NOT NULL */ + if (!PG_ARGISNULL(1)) { + /* deserialize raster */ + pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); + + /* Get raster object */ + raster = rt_raster_deserialize(pgraster, FALSE); + if (raster == NULL) { + elog(ERROR, "RASTER_union_transfn: Could not deserialize raster"); + + pfree(iwr->raster); + pfree(iwr); + PG_FREE_IF_COPY(pgraster, 1); + + MemoryContextSwitchTo(oldcontext); + PG_RETURN_NULL(); + } + } + + /* if more than 2 arguments, determine the type of argument 3 */ + /* band number or UNION type */ + if (PG_NARGS() > 2 && !PG_ARGISNULL(2)) { + Oid calltype = get_fn_expr_argtype(fcinfo->flinfo, 3); + + switch (calltype) { + /* UNION type */ + case TEXTOID: + utypename = text_to_cstring(PG_GETARG_TEXT_P(3)); + utype = rtpg_uniontype_index_from_name(rtpg_strtoupper(utypename)); + iwr->uniontype = utype; + break; + case INT2OID: + case INT4OID: + nband = PG_GETARG_INT32(2); + if (nband < 1) { + elog(ERROR, "RASTER_union_transfn: Band number must be greater than zero (1-based)"); + + pfree(iwr->raster); + pfree(iwr); + if (raster != NULL) { + rt_raster_destroy(raster); + PG_FREE_IF_COPY(pgraster, 1); + } + + MemoryContextSwitchTo(oldcontext); + PG_RETURN_NULL(); + } + break; + } + } + + /* UNION type */ + if (PG_NARGS() > 3 && !PG_ARGISNULL(3)) { + utypename = text_to_cstring(PG_GETARG_TEXT_P(3)); + utype = rtpg_uniontype_index_from_name(rtpg_strtoupper(utypename)); + iwr->uniontype = utype; + } + if (utype == UT_MEAN) + iwr->numraster = 2; + + /* allocate iwr->raster */ + if (iwr->raster == NULL) { + iwr->raster = (rt_raster *) rtalloc(sizeof(rt_raster) * iwr->numraster); + if (iwr->raster == NULL) { + elog(ERROR, "RASTER_union_transfn: Unable to allocate memory for pointer to raster"); + + pfree(iwr); + if (raster != NULL) { + rt_raster_destroy(raster); + PG_FREE_IF_COPY(pgraster, 1); + } + + MemoryContextSwitchTo(oldcontext); + PG_RETURN_NULL(); + } + + memset(iwr->raster, 0, sizeof(rt_raster) * iwr->numraster); + } + + for (i = 0; i < iwr->numraster; i++) { + /* raster flags */ + isempty[0] = rt_raster_is_empty(iwr->raster[i]); + isempty[1] = rt_raster_is_empty(raster); + if (!isempty[0]) + hasband[0] = rt_raster_has_band(iwr->raster[i], 0); + if (!isempty[1]) + hasband[1] = rt_raster_has_band(raster, nband - 1); + + /* determine pixtype, hasnodata and nodataval */ + _band = NULL; + if (!isempty[0] && hasband[0]) + _band = rt_raster_get_band(iwr->raster[i], 0); + else if (!isempty[1] && hasband[1]) + _band = rt_raster_get_band(raster, nband - 1); + else { + pixtype = PT_64BF; + hasnodata = 1; + nodataval = rt_pixtype_get_min_value(pixtype); + } + if (_band != NULL) { + pixtype = rt_band_get_pixtype(_band); + hasnodata = rt_band_get_hasnodata_flag(_band); + if (hasnodata) + nodataval = rt_band_get_nodata(_band); + } + + /* UT_MEAN requires two passes, first for UT_COUNT and second for UT_SUM */ + if (iwr->uniontype == UT_MEAN) { + /* first pass, UT_COUNT */ + if (i < 1) + utype = UT_COUNT; + else + utype = UT_SUM; + } + + /* force band settings for UT_COUNT */ + if (utype == UT_COUNT) { + pixtype = PT_32BUI; + hasnodata = 0; + nodataval = 0; + } + + /* build itrset */ + itrset = palloc(sizeof(struct rt_iterator_t) * 2); + if (itrset == NULL) { + elog(ERROR, "RASTER_union_transfn: Unable to allocate memory for iterator arguments"); + + pfree(iwr->raster); + pfree(iwr); + if (raster != NULL) { + rt_raster_destroy(raster); + PG_FREE_IF_COPY(pgraster, 1); + } + + MemoryContextSwitchTo(oldcontext); + PG_RETURN_NULL(); + } + + itrset[0].raster = iwr->raster[i]; + itrset[0].nband = 0; + itrset[1].raster = raster; + itrset[1].nband = nband - 1; + + /* pass everything to iterator */ + _raster = rt_raster_iterator( + itrset, 2, + ET_UNION, NULL, + pixtype, + hasnodata, nodataval, + 0, 0, + &utype, + rtpg_union_callback, + &noerr + ); + + /* proactive cleanup */ + pfree(itrset); + + if (!noerr) { + elog(ERROR, "RASTER_union_transfn: Unable to run raster interator function"); + + pfree(iwr->raster); + pfree(iwr); + if (raster != NULL) { + rt_raster_destroy(raster); + PG_FREE_IF_COPY(pgraster, 1); + } + + MemoryContextSwitchTo(oldcontext); + PG_RETURN_NULL(); + } + + /* replace working raster */ + if (iwr->raster[i] != NULL) + rt_raster_destroy(iwr->raster[i]); + iwr->raster[i] = _raster; + } + + if (raster != NULL) { + rt_raster_destroy(raster); + PG_FREE_IF_COPY(pgraster, 1); + } + + /* switch back to local context */ + MemoryContextSwitchTo(oldcontext); + + PG_RETURN_POINTER(iwr); +} + +/* UNION aggregate final function */ +PG_FUNCTION_INFO_V1(RASTER_union_finalfn); +Datum RASTER_union_finalfn(PG_FUNCTION_ARGS) +{ + rtpg_union_arg *iwr; + rt_pgraster *pgraster = NULL; + + /* cannot be called directly as this is exclusive aggregate function */ + if (!AggCheckCallContext(fcinfo, NULL)) { + elog(ERROR, "RASTER_union_finalfn: Cannot be called in a non-aggregate context"); + PG_RETURN_NULL(); + } + + /* NULL, return null */ + if (PG_ARGISNULL(0)) + PG_RETURN_NULL(); + + iwr = (rtpg_union_arg *) PG_GETARG_POINTER(0); + + /* all UNION types except UT_MEAN have one raster */ + if (iwr->numraster < 2) { + /* return NULL if raster is NULL */ + if (iwr->raster[0] == NULL) { + pfree(iwr->raster); + pfree(iwr); + PG_RETURN_NULL(); + } + + pgraster = rt_raster_serialize(iwr->raster[0]); + rt_raster_destroy(iwr->raster[0]); + pfree(iwr->raster); + pfree(iwr); + } + /* special case for UT_MEAN */ + else { + int i; + + /* return NULL if either raster is NULL */ + if (iwr->raster[0] == NULL || iwr->raster[1] == NULL) { + for (i = 0; i < iwr->numraster; i++) { + if (iwr->raster[i] != NULL) + rt_raster_destroy(iwr->raster[i]); + } + pfree(iwr->raster); + pfree(iwr); + PG_RETURN_NULL(); + } + + /* if second raster is empty, return empty raster */ + if (rt_raster_is_empty(iwr->raster[1])) { + pgraster = rt_raster_serialize(iwr->raster[1]); + + for (i = 0; i < iwr->numraster; i++) + rt_raster_destroy(iwr->raster[i]); + pfree(iwr->raster); + pfree(iwr); + } + else { + rt_iterator itrset; + rt_raster _raster = NULL; + rt_band _band = NULL; + rt_pixtype pixtype; + int hasnodata = 1; + double nodataval = 0; + int noerr = 1; + + _band = rt_raster_get_band(iwr->raster[1], 0); + pixtype = rt_band_get_pixtype(_band); + hasnodata = rt_band_get_hasnodata_flag(_band); + if (hasnodata) + nodataval = rt_band_get_nodata(_band); + + /* build itrset */ + itrset = palloc(sizeof(struct rt_iterator_t) * 2); + if (itrset == NULL) { + elog(ERROR, "RASTER_union_finalfn: Unable to allocate memory for iterator arguments"); + + for (i = 0; i < iwr->numraster; i++) + rt_raster_destroy(iwr->raster[i]); + pfree(iwr->raster); + pfree(iwr); + PG_RETURN_NULL(); + } + + itrset[0].raster = iwr->raster[0]; + itrset[0].nband = 0; + itrset[1].raster = iwr->raster[1]; + itrset[1].nband = 0; + + /* pass everything to iterator */ + _raster = rt_raster_iterator( + itrset, 2, + ET_UNION, NULL, + pixtype, + hasnodata, nodataval, + 0, 0, + NULL, + rtpg_union_mean_callback, + &noerr + ); + + /* proactive cleanup */ + pfree(itrset); + for (i = 0; i < iwr->numraster; i++) + rt_raster_destroy(iwr->raster[i]); + pfree(iwr->raster); + pfree(iwr); + + if (!noerr) { + elog(ERROR, "RASTER_union_finalfn: Unable to run raster interator function"); + 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); +} + /* ---------------------------------------------------------------- */ /* Memory allocation / error reporting hooks */ /* ---------------------------------------------------------------- */ diff --git a/raster/rt_pg/rtpostgis.sql.in.c b/raster/rt_pg/rtpostgis.sql.in.c index 065ffeab7..c95bbef71 100644 --- a/raster/rt_pg/rtpostgis.sql.in.c +++ b/raster/rt_pg/rtpostgis.sql.in.c @@ -4398,131 +4398,55 @@ CREATE OR REPLACE FUNCTION st_intersection( LANGUAGE 'sql' STABLE; ----------------------------------------------------------------------- --- st_union aggregate +-- ***NEW*** ST_Union aggregate ----------------------------------------------------------------------- --- Main state function -CREATE OR REPLACE FUNCTION _ST_MapAlgebra4UnionState(rast1 raster, rast2 raster, p_expression text, p_nodata1expr text, p_nodata2expr text, p_nodatanodataval double precision,t_expression text,t_nodata1expr text, t_nodata2expr text,t_nodatanodataval double precision) - RETURNS raster AS - $$ - DECLARE - t_raster raster; - p_raster raster; - BEGIN - -- With the new ST_MapAlgebraExpr we must split the main expression in three expressions: expression, nodata1expr, nodata2expr and a nodatanodataval - -- ST_MapAlgebraExpr(rast1 raster, band1 integer, rast2 raster, band2 integer, expression text, pixeltype text, extentexpr text, nodata1expr text, nodata2expr text, nodatanodatadaval double precision) - -- We must make sure that when NULL is passed as the first raster to ST_MapAlgebraExpr, ST_MapAlgebraExpr resolve the nodata1expr - -- Note: rast2 is always a single band raster since it is the accumulated raster thus far - -- There we always set that to band 1 regardless of what band num is requested - IF upper(p_expression) = 'LAST' THEN - --RAISE NOTICE 'last asked for '; - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, '[rast2.val]'::text, NULL::text, 'UNION'::text, '[rast2.val]'::text, '[rast1.val]'::text, NULL::double precision); - ELSIF upper(p_expression) = 'FIRST' THEN - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, '[rast1.val]'::text, NULL::text, 'UNION'::text, '[rast2.val]'::text, '[rast1.val]'::text, NULL::double precision); - ELSIF upper(p_expression) = 'MIN' THEN - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, 'LEAST([rast1.val], [rast2.val])'::text, NULL::text, 'UNION'::text, '[rast2.val]'::text, '[rast1.val]'::text, NULL::double precision); - ELSIF upper(p_expression) = 'MAX' THEN - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, 'GREATEST([rast1.val], [rast2.val])'::text, NULL::text, 'UNION'::text, '[rast2.val]'::text, '[rast1.val]'::text, NULL::double precision); - ELSIF upper(p_expression) = 'COUNT' THEN - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, '[rast1.val] + 1'::text, NULL::text, 'UNION'::text, '1'::text, '[rast1.val]'::text, 0::double precision); - ELSIF upper(p_expression) = 'SUM' THEN - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, '[rast1.val] + [rast2.val]'::text, NULL::text, 'UNION'::text, '[rast2.val]'::text, '[rast1.val]'::text, NULL::double precision); - ELSIF upper(p_expression) = 'RANGE' THEN - -- have no idea what this is - t_raster = ST_MapAlgebraExpr(rast1, 2, rast2, 1, 'LEAST([rast1.val], [rast2.val])'::text, NULL::text, 'UNION'::text, '[rast2.val]'::text, '[rast1.val]'::text, NULL::double precision); - p_raster := _ST_MapAlgebra4UnionState(rast1, rast2, 'MAX'::text, NULL::text, NULL::text, NULL::double precision, NULL::text, NULL::text, NULL::text, NULL::double precision); - RETURN ST_AddBand(p_raster, t_raster, 1, 2); - ELSIF upper(p_expression) = 'MEAN' THEN - -- looks like t_raster is used to keep track of accumulated count - -- and p_raster is there to keep track of accumulated sum and final state function - -- would then do a final map to divide them. This one is currently broken because - -- have not reworked it so it can do without a final function - t_raster = ST_MapAlgebraExpr(rast1, 2, rast2, 1, '[rast1.val] + 1'::text, NULL::text, 'UNION'::text, '1'::text, '[rast1.val]'::text, 0::double precision); - p_raster := _ST_MapAlgebra4UnionState(rast1, rast2, 'SUM'::text, NULL::text, NULL::text, NULL::double precision, NULL::text, NULL::text, NULL::text, NULL::double precision); - RETURN ST_AddBand(p_raster, t_raster, 1, 2); - ELSE - IF t_expression NOTNULL AND t_expression != '' THEN - t_raster = ST_MapAlgebraExpr(rast1, 2, rast2, 1, t_expression, NULL::text, 'UNION'::text, t_nodata1expr, t_nodata2expr, t_nodatanodataval::double precision); - p_raster = ST_MapAlgebraExpr(rast1, 1, rast2, 1, p_expression, NULL::text, 'UNION'::text, p_nodata1expr, p_nodata2expr, p_nodatanodataval::double precision); - RETURN ST_AddBand(p_raster, t_raster, 1, 2); - END IF; - RETURN ST_MapAlgebraExpr(rast1, 1, rast2, 1, p_expression, NULL, 'UNION'::text, NULL::text, NULL::text, NULL::double precision); - END IF; - END; - $$ - LANGUAGE 'plpgsql'; +CREATE OR REPLACE FUNCTION _st_union_transfn(internal, raster, integer, text) + RETURNS internal + AS 'MODULE_PATHNAME', 'RASTER_union_transfn' + LANGUAGE 'c' IMMUTABLE; --- State function when there is a primary expression, band number and no alternative nodata expressions -CREATE OR REPLACE FUNCTION _ST_MapAlgebra4UnionState(rast1 raster,rast2 raster,bandnum integer, p_expression text) - RETURNS raster - AS $$ - SELECT _ST_MapAlgebra4UnionState($1, ST_Band($2,$3), $4, NULL, NULL, NULL, NULL, NULL, NULL, NULL) - $$ LANGUAGE 'sql'; +CREATE OR REPLACE FUNCTION _st_union_finalfn(internal) + RETURNS raster + AS 'MODULE_PATHNAME', 'RASTER_union_finalfn' + LANGUAGE 'c' IMMUTABLE; --- State function when there is no expressions but allows specifying band -CREATE OR REPLACE FUNCTION _ST_MapAlgebra4UnionState(rast1 raster,rast2 raster, bandnum integer) - RETURNS raster - AS $$ - SELECT _ST_MapAlgebra4UnionState($1,ST_Band($2,$3), 'LAST', NULL, NULL, NULL, NULL, NULL, NULL, NULL) - $$ LANGUAGE 'sql'; - --- State function when there is no expressions and assumes band 1 -CREATE OR REPLACE FUNCTION _ST_MapAlgebra4UnionState(rast1 raster,rast2 raster) - RETURNS raster - AS $$ - SELECT _ST_MapAlgebra4UnionState($1,$2, 'LAST', NULL, NULL, NULL, NULL, NULL, NULL, NULL) - $$ LANGUAGE 'sql'; - --- State function when there isan expressions and assumes band 1 -CREATE OR REPLACE FUNCTION _ST_MapAlgebra4UnionState(rast1 raster,rast2 raster, p_expression text) - RETURNS raster - AS $$ - SELECT _ST_MapAlgebra4UnionState($1,$2, $3, NULL, NULL, NULL, NULL, NULL, NULL, NULL) - $$ LANGUAGE 'sql'; - --- Final function with only the primary expression -CREATE OR REPLACE FUNCTION _ST_MapAlgebra4UnionFinal1(rast raster) - RETURNS raster AS - $$ - DECLARE - BEGIN - -- NOTE: I have to sacrifice RANGE. Sorry RANGE. Any 2 banded raster is going to be treated - -- as a MEAN - IF ST_NumBands(rast) = 2 THEN - RETURN ST_MapAlgebraExpr(rast, 1, rast, 2, 'CASE WHEN [rast2.val] > 0 THEN [rast1.val] / [rast2.val]::float8 ELSE NULL END'::text, NULL::text, 'UNION'::text, NULL::text, NULL::text, NULL::double precision); - ELSE - RETURN rast; - END IF; - END; - $$ - LANGUAGE 'plpgsql'; - --- Variant with primary expression defaulting to 'LAST' and working on first band -CREATE AGGREGATE ST_Union(raster) ( - SFUNC = _ST_MapAlgebra4UnionState, - STYPE = raster, - FINALFUNC = _ST_MapAlgebra4UnionFinal1 +CREATE AGGREGATE st_union(raster, integer, text) ( + SFUNC = _st_union_transfn, + STYPE = internal, + FINALFUNC = _st_union_finalfn ); --- Variant with primary expression defaulting to 'LAST' and working on specified band -CREATE AGGREGATE ST_Union(raster, integer) ( - SFUNC = _ST_MapAlgebra4UnionState, - STYPE = raster, - FINALFUNC = _ST_MapAlgebra4UnionFinal1 +CREATE OR REPLACE FUNCTION _st_union_transfn(internal, raster, integer) + RETURNS internal + AS 'MODULE_PATHNAME', 'RASTER_union_transfn' + LANGUAGE 'c' IMMUTABLE; + +CREATE AGGREGATE st_union(raster, integer) ( + SFUNC = _st_union_transfn, + STYPE = internal, + FINALFUNC = _st_union_finalfn ); --- Variant with simple primary expressions but without alternative nodata, temporary and final expressions --- and working on first band --- supports LAST, MIN,MAX,MEAN,FIRST,SUM -CREATE AGGREGATE ST_Union(raster, text) ( - SFUNC = _ST_MapAlgebra4UnionState, - STYPE = raster, - FINALFUNC = _ST_MapAlgebra4UnionFinal1 +CREATE OR REPLACE FUNCTION _st_union_transfn(internal, raster) + RETURNS internal + AS 'MODULE_PATHNAME', 'RASTER_union_transfn' + LANGUAGE 'c' IMMUTABLE; + +CREATE AGGREGATE st_union(raster) ( + SFUNC = _st_union_transfn, + STYPE = internal, + FINALFUNC = _st_union_finalfn ); -CREATE AGGREGATE ST_Union(raster, integer, text) ( - SFUNC = _ST_MapAlgebra4UnionState, - STYPE = raster, - FINALFUNC = _ST_MapAlgebra4UnionFinal1 +CREATE OR REPLACE FUNCTION _st_union_transfn(internal, raster, text) + RETURNS internal + AS 'MODULE_PATHNAME', 'RASTER_union_transfn' + LANGUAGE 'c' IMMUTABLE; + +CREATE AGGREGATE st_union(raster, text) ( + SFUNC = _st_union_transfn, + STYPE = internal, + FINALFUNC = _st_union_finalfn ); ------------------------------------------------------------------- diff --git a/raster/test/core/testapi.c b/raster/test/core/testapi.c index 65ca4162a..9263270f2 100644 --- a/raster/test/core/testapi.c +++ b/raster/test/core/testapi.c @@ -119,7 +119,7 @@ static void testGDALPolygonize() { char *wkt = NULL; rt = fillRasterToPolygonize(1, -1.0); - CHECK(!rt_raster_has_no_band(rt, 0)); + CHECK(rt_raster_has_band(rt, 0)); nPols = 0; gv = rt_raster_gdal_polygonize(rt, 0, TRUE, &nPols); @@ -173,8 +173,8 @@ static void testGDALPolygonize() { rt = fillRasterToPolygonize(1, 2.0); #endif - /* We can check rt_raster_has_no_band here too */ - CHECK(!rt_raster_has_no_band(rt, 0)); + /* We can check rt_raster_has_band here too */ + CHECK(rt_raster_has_band(rt, 0)); nPols = 0; gv = rt_raster_gdal_polygonize(rt, 0, TRUE, &nPols); @@ -230,8 +230,8 @@ static void testGDALPolygonize() { rt = fillRasterToPolygonize(1, 3.0); #endif - /* We can check rt_raster_has_no_band here too */ - CHECK(!rt_raster_has_no_band(rt, 0)); + /* We can check rt_raster_has_band here too */ + CHECK(rt_raster_has_band(rt, 0)); nPols = 0; gv = rt_raster_gdal_polygonize(rt, 0, TRUE, &nPols); @@ -275,8 +275,8 @@ static void testGDALPolygonize() { /* Fourth test: NODATA value = 0 */ rt = fillRasterToPolygonize(1, 0.0); - /* We can check rt_raster_has_no_band here too */ - CHECK(!rt_raster_has_no_band(rt, 0)); + /* We can check rt_raster_has_band here too */ + CHECK(rt_raster_has_band(rt, 0)); nPols = 0; gv = rt_raster_gdal_polygonize(rt, 0, TRUE, &nPols); @@ -315,8 +315,8 @@ static void testGDALPolygonize() { /* Last test: There is no NODATA value (all values are valid) */ rt = fillRasterToPolygonize(0, 0.0); - /* We can check rt_raster_has_no_band here too */ - CHECK(!rt_raster_has_no_band(rt, 0)); + /* We can check rt_raster_has_band here too */ + CHECK(rt_raster_has_band(rt, 0)); nPols = 0; gv = rt_raster_gdal_polygonize(rt, 0, TRUE, &nPols); @@ -1202,7 +1202,7 @@ static void testRasterFromBand() { CHECK(rast); CHECK(!rt_raster_is_empty(rast)); - CHECK(!rt_raster_has_no_band(rast, 1)); + CHECK(rt_raster_has_band(rast, 1)); deepRelease(rast); deepRelease(raster); @@ -7827,9 +7827,9 @@ static void testRasterIterator() { /* allocate itrset */ itrset = rtalloc(sizeof(struct rt_iterator_t) * num); - itrset[0].rast = rast1; + itrset[0].raster = rast1; itrset[0].nband = 0; - itrset[1].rast = rast2; + itrset[1].raster = rast2; itrset[1].nband = 0; /* 1 raster, 0 distance, FIRST or SECOND or LAST or UNION or INTERSECTION */ @@ -8053,7 +8053,7 @@ main() rt_raster_destroy(emptyraster); /* Once we add a band to this raster, we'll check the opposite */ - CHECK(rt_raster_has_no_band(raster, 1)); + CHECK(!rt_raster_has_band(raster, 1)); } diff --git a/raster/test/regress/rt_union.sql b/raster/test/regress/rt_union.sql index 3a63c6a8f..b55fd6d8c 100644 --- a/raster/test/regress/rt_union.sql +++ b/raster/test/regress/rt_union.sql @@ -1,16 +1,74 @@ -SELECT '#1 ' as run, ST_AsText((gval).geom), (gval).val::text - FROM (SELECT ST_DumpAsPolygons(ST_Union(rast) ) As gval - FROM ( - SELECT i As rid, ST_AddBand( - ST_MakeEmptyRaster(10, 10, 10*i, 10*i, 2, 2, 0, 0, ST_SRID(ST_Point(0,0) )), - '8BUI'::text, 10*i, 0) As rast - FROM generate_series(0,10) As i ) As foo ) As foofoo ORDER BY (gval).val; - -SELECT '#2 ' As run, ST_AsText((gval).geom), (gval).val::text - FROM (SELECT ST_DumpAsPolygons(ST_Union(rast,'MEAN') ) As gval - FROM ( - SELECT i As rid, ST_AddBand( - ST_MakeEmptyRaster(10, 10, 10*i, 10*i, 2, 2, 0, 0, ST_SRID(ST_Point(0,0) )), - '8BUI'::text, 10*i, 0) As rast - FROM generate_series(0,10) As i ) As foo ) As foofoo - ORDER BY (gval).val LIMIT 2; +DROP TABLE IF EXISTS raster_union_in; +CREATE TABLE raster_union_in ( + rid integer, + rast raster +); +DROP TABLE IF EXISTS raster_union_out; +CREATE TABLE raster_union_out ( + uniontype text, + rast raster +); + +INSERT INTO raster_union_in + SELECT 0, NULL::raster AS rast UNION ALL + SELECT 1, ST_AddBand(ST_MakeEmptyRaster(2, 2, 0, 0, 1, -1, 0, 0, 0), 1, '8BUI', 1, 0) AS rast UNION ALL + SELECT 2, ST_AddBand(ST_MakeEmptyRaster(2, 2, 1, -1, 1, -1, 0, 0, 0), 1, '8BUI', 2, 0) AS rast +; + +INSERT INTO raster_union_out + SELECT + 'LAST', + ST_Union(rast, 1) AS rast + FROM raster_union_in; + +INSERT INTO raster_union_out + SELECT + 'FIRST', + ST_Union(rast, 1, 'FIRST') AS rast + FROM raster_union_in; + +INSERT INTO raster_union_out + SELECT + 'MIN', + ST_Union(rast, 1, 'MIN') AS rast + FROM raster_union_in; + +INSERT INTO raster_union_out + SELECT + 'MAX', + ST_Union(rast, 1, 'MAX') AS rast + FROM raster_union_in; + +INSERT INTO raster_union_out + SELECT + 'COUNT', + ST_Union(rast, 1, 'COUNT') AS rast + FROM raster_union_in; + +INSERT INTO raster_union_out + SELECT + 'SUM', + ST_Union(rast, 1, 'SUM') AS rast + FROM raster_union_in; + +INSERT INTO raster_union_out + SELECT + 'MEAN', + ST_Union(rast, 'mean') AS rast + FROM raster_union_in; + +SELECT + uniontype, + x, + y, + val +FROM ( + SELECT + uniontype, + (ST_PixelAsPoints(rast)).* + FROM raster_union_out +) foo +ORDER BY uniontype, y, x; + +DROP TABLE IF EXISTS raster_union_in; +DROP TABLE IF EXISTS raster_union_out; diff --git a/raster/test/regress/rt_union_expected b/raster/test/regress/rt_union_expected index fd18af417..c3caa1856 100644 --- a/raster/test/regress/rt_union_expected +++ b/raster/test/regress/rt_union_expected @@ -1,12 +1,65 @@ -#1 |POLYGON((10 10,10 30,20 30,20 20,28 20,30 20,30 10,10 10))|10 -#1 |POLYGON((20 20,20 40,30 40,30 30,38 30,40 30,40 20,20 20))|20 -#1 |POLYGON((30 30,30 50,40 50,40 40,48 40,50 40,50 30,30 30))|30 -#1 |POLYGON((40 40,40 60,50 60,50 50,58 50,60 50,60 40,40 40))|40 -#1 |POLYGON((50 50,50 70,60 70,60 60,68 60,70 60,70 50,50 50))|50 -#1 |POLYGON((60 60,60 80,70 80,70 70,78 70,80 70,80 60,60 60))|60 -#1 |POLYGON((70 70,70 90,80 90,80 80,88 80,90 80,90 70,70 70))|70 -#1 |POLYGON((80 80,80 100,90 100,90 90,98 90,100 90,100 80,80 80))|80 -#1 |POLYGON((90 90,90 110,100 110,100 100,108 100,110 100,110 90,90 90))|90 -#1 |POLYGON((100 100,100 120,120 120,120 100,100 100))|100 -#2 |POLYGON((10 10,10 30,20 30,20 20,28 20,30 20,30 10,10 10))|10 -#2 |POLYGON((20 20,20 30,30 30,30 20,20 20))|15 +NOTICE: table "raster_union_in" does not exist, skipping +NOTICE: table "raster_union_out" does not exist, skipping +COUNT|1|1|1 +COUNT|2|1|1 +COUNT|3|1|0 +COUNT|1|2|1 +COUNT|2|2|2 +COUNT|3|2|1 +COUNT|1|3|0 +COUNT|2|3|1 +COUNT|3|3|1 +FIRST|1|1|1 +FIRST|2|1|1 +FIRST|3|1| +FIRST|1|2|1 +FIRST|2|2|1 +FIRST|3|2|2 +FIRST|1|3| +FIRST|2|3|2 +FIRST|3|3|2 +LAST|1|1|1 +LAST|2|1|1 +LAST|3|1| +LAST|1|2|1 +LAST|2|2|2 +LAST|3|2|2 +LAST|1|3| +LAST|2|3|2 +LAST|3|3|2 +MAX|1|1|1 +MAX|2|1|1 +MAX|3|1| +MAX|1|2|1 +MAX|2|2|2 +MAX|3|2|2 +MAX|1|3| +MAX|2|3|2 +MAX|3|3|2 +MEAN|1|1|1 +MEAN|2|1|1 +MEAN|3|1| +MEAN|1|2|1 +MEAN|2|2|2 +MEAN|3|2|2 +MEAN|1|3| +MEAN|2|3|2 +MEAN|3|3|2 +MIN|1|1|1 +MIN|2|1|1 +MIN|3|1| +MIN|1|2|1 +MIN|2|2|1 +MIN|3|2|2 +MIN|1|3| +MIN|2|3|2 +MIN|3|3|2 +SUM|1|1|1 +SUM|2|1|1 +SUM|3|1| +SUM|1|2|1 +SUM|2|2|3 +SUM|3|2|2 +SUM|1|3| +SUM|2|3|2 +SUM|3|3|2 -- 2.40.0