Index: raster/test/core/testapi.c =================================================================== --- raster/test/core/testapi.c (revision 7056) +++ raster/test/core/testapi.c (working copy) @@ -972,6 +972,21 @@ CHECK_EQUALS(flag, 1); } +static void testRasterFromBand(rt_context ctx, rt_raster raster) { + uint32_t bandNums[] = {1,3}; + int lenBandNums = 2; + rt_raster rast; + + rast = rt_raster_from_band(ctx, raster, bandNums, lenBandNums); + assert(rast); + + CHECK(rast); + CHECK(!rt_raster_is_empty(ctx, rast)); + CHECK(!rt_raster_has_no_band(ctx, rast, 1)); + + rt_raster_destroy(ctx, rast); +} + int main() { @@ -1282,6 +1297,10 @@ printf("Testing band hasnodata flag\n"); testBandHasNoData(ctx, band_64BF); + printf("Testing rt_raster_from_band\n"); + testRasterFromBand(ctx, raster); + printf("Successfully tested rt_raster_from_band\n"); + deepRelease(ctx, raster); rt_context_destroy(ctx); Index: raster/test/regress/Makefile.in =================================================================== --- raster/test/regress/Makefile.in (revision 7056) +++ raster/test/regress/Makefile.in (working copy) @@ -41,6 +41,7 @@ create_rt_box2d_test.sql \ rt_box2d.sql \ rt_addband.sql \ + rt_band.sql \ $(NULL) TEST_PROPS = \ Index: raster/test/regress/rt_band_expected =================================================================== --- raster/test/regress/rt_band_expected (revision 0) +++ raster/test/regress/rt_band_expected (revision 0) @@ -0,0 +1,20 @@ +123.4567 +1234.567 +1234.567 +1234.5678 +987.654321 +9876.54321 +1234.5678 +987.654321 +9876.54321 +1234.5678 +1234.5678 +9876.54321 +987.654321 +1234.5678 +4 +3 +2 +1 +1 +1 Index: raster/test/regress/rt_band.sql =================================================================== --- raster/test/regress/rt_band.sql (revision 0) +++ raster/test/regress/rt_band.sql (revision 0) @@ -0,0 +1,256 @@ +SELECT ST_Value(ST_Band(ST_AddBand(ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0, -1), 1, '64BF', 123.4567, NULL), ARRAY[1]), 3, 3); +SELECT ST_Value(ST_Band(ST_AddBand(ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0, -1), 1, '64BF', 1234.567, NULL), 1), 3, 3); +SELECT ST_Value(ST_Band(ST_AddBand(ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0, -1), 1, '64BF', 1234.567, NULL)), 3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1] + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[2] + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[3] + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + 1 + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + 2 + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + 3 + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ) + ), +3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1,3] + ), +1, 3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1,3] + ), +2, 3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[2,3] + ), +1, 3, 3); +SELECT ST_Value( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1,1] + ), +2, 3, 3); +SELECT ST_NumBands( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1,1,3,3] + ) +); +SELECT ST_NumBands( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1,1,3] + ) +); +SELECT ST_NumBands( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[1,2] + ) +); +SELECT ST_NumBands( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + ARRAY[3] + ) +); +SELECT ST_NumBands( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ), + 2 + ) +); +SELECT ST_NumBands( + ST_Band( + ST_AddBand( + ST_AddBand( + ST_AddBand( + ST_MakeEmptyRaster(200, 200, 10, 10, 2, 2, 0, 0,-1) + , 1, '64BF', 1234.5678, NULL + ) + , '64BF', 987.654321, NULL + ) + , '64BF', 9876.54321, NULL + ) + ) +); Index: raster/rt_pg/rtpostgis.sql.in.c =================================================================== --- raster/rt_pg/rtpostgis.sql.in.c (revision 7056) +++ raster/rt_pg/rtpostgis.sql.in.c (working copy) @@ -229,6 +229,34 @@ LANGUAGE 'SQL' IMMUTABLE STRICT; ----------------------------------------------------------------------- +-- Constructor ST_Band +----------------------------------------------------------------------- +CREATE OR REPLACE FUNCTION st_band(rast raster, nbands int[]) + RETURNS RASTER + AS 'MODULE_PATHNAME', 'RASTER_band' + LANGUAGE 'C' IMMUTABLE STRICT; + +CREATE OR REPLACE FUNCTION st_band(rast raster, nband int) + RETURNS RASTER + AS $$ SELECT st_band($1, ARRAY[$2]) $$ + LANGUAGE 'SQL' IMMUTABLE STRICT; + +CREATE OR REPLACE FUNCTION st_band(rast raster, nbands text) + RETURNS RASTER + AS $$ SELECT st_band($1, regexp_split_to_array(regexp_replace($2, '[[:space:]]', '', 'g'), ',')::int[]) $$ + LANGUAGE 'SQL' IMMUTABLE STRICT; + +CREATE OR REPLACE FUNCTION st_band(rast raster, nbands text, delimiter char) + RETURNS RASTER + AS $$ SELECT st_band($1, regexp_split_to_array(regexp_replace($2, '[[:space:]]', '', 'g'), $3)::int[]) $$ + LANGUAGE 'SQL' IMMUTABLE STRICT; + +CREATE OR REPLACE FUNCTION st_band(rast raster) + RETURNS RASTER + AS $$ SELECT st_band($1, ARRAY[1]) $$ + LANGUAGE 'SQL' IMMUTABLE STRICT; + +----------------------------------------------------------------------- -- MapAlgebra ----------------------------------------------------------------------- -- This function can not be STRICT, because nodatavalueexpr can be NULL (could be just '' though) Index: raster/rt_pg/rt_pg.c =================================================================== --- raster/rt_pg/rt_pg.c (revision 7056) +++ raster/rt_pg/rt_pg.c (working copy) @@ -49,6 +49,10 @@ #include "rt_api.h" #include "../raster_config.h" +#include /* for get_typlenbyvalalign */ +#include /* for ArrayType */ +#include /* for INT2OID, INT4OID and TEXTOID */ + #define POSTGIS_RASTER_WARN_ON_TRUNCATION /* @@ -158,7 +162,10 @@ /* Raster analysis */ Datum RASTER_mapAlgebra(PG_FUNCTION_ARGS); +/* create new raster from existing raster's bands */ +Datum RASTER_band(PG_FUNCTION_ARGS); + /* Replace function taken from http://ubuntuforums.org/showthread.php?s=aa6f015109fd7e4c7e30d2fd8b717497&t=141670&page=3 */ /* --------------------------------------------------------------------------- Name : replace - Search & replace a substring by another one. @@ -2677,6 +2684,120 @@ PG_RETURN_POINTER(pgraster); } +/** + * Return new raster from selected bands of existing raster through ST_Band. + * second argument is an array of band numbers (1 based) + */ +PG_FUNCTION_INFO_V1(RASTER_band); +Datum RASTER_band(PG_FUNCTION_ARGS) +{ + rt_pgraster *pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); + rt_pgraster *pgrast; + rt_raster raster; + rt_raster rast; + rt_context ctx = get_rt_context(fcinfo); + + ArrayType *array; + Oid etype; + Datum *e; + bool *nulls; + int16 typlen; + bool typbyval; + char typalign; + int ndims = 1; + int *dims; + int *lbs; + + uint32_t *bandNums; + uint32 idx = 0; + int n; + int i = 0; + int j = 0; + + raster = rt_raster_deserialize(ctx, pgraster); + if (!raster) { + elog(ERROR, "RASTER_band: Could not deserialize raster"); + rt_context_destroy(ctx); + PG_RETURN_NULL(); + } + + /* process bandNums */ + if (PG_ARGISNULL(1)) { + elog(ERROR, "RASTER_band: Band number must be provided"); + rt_raster_destroy(ctx, raster); + rt_context_destroy(ctx); + PG_RETURN_NULL(); + } + + array = PG_GETARG_ARRAYTYPE_P(1); + etype = ARR_ELEMTYPE(array); + get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign); + + switch (etype) { + case INT2OID: + case INT4OID: + break; + default: + elog(ERROR, "RASTER_band: Invalid data type for band numbers"); + rt_raster_destroy(ctx, raster); + rt_context_destroy(ctx); + PG_RETURN_NULL(); + break; + } + + ndims = ARR_NDIM(array); + dims = ARR_DIMS(array); + lbs = ARR_LBOUND(array); + + deconstruct_array(array, etype, typlen, typbyval, typalign, &e, + &nulls, &n); + + bandNums = (uint32_t *) palloc(sizeof(uint32_t) * n); + for (i = 0, j = 0; i < n; i++) { + if (nulls[i]) + continue; + + switch (etype) { + case INT2OID: + idx = (uint32_t) DatumGetInt16(e[i]); + break; + case INT4OID: + idx = (uint32_t) DatumGetInt32(e[i]); + break; + } + + POSTGIS_RT_DEBUGF(3, "band idx (before): %d", idx); + if (idx > pgraster->numBands || idx < 1) { + pfree(bandNums); + elog(ERROR, "RASTER_band: Invalid band number provided"); + rt_raster_destroy(ctx, raster); + rt_context_destroy(ctx); + PG_RETURN_NULL(); + } + + bandNums[j] = idx - 1; + POSTGIS_RT_DEBUGF(3, "bandNums[%d] = %d", j, bandNums[j]); + j++; + } + + rast = rt_raster_from_band(ctx, raster, bandNums, j); + pfree(bandNums); + rt_raster_destroy(ctx, raster); + if (!rast) { + elog(ERROR, "RASTER_band: Could not create new raster"); + rt_context_destroy(ctx); + PG_RETURN_NULL(); + } + + pgrast = rt_raster_serialize(ctx, rast); + rt_context_destroy(ctx); + if (!pgrast) + PG_RETURN_NULL(); + + SET_VARSIZE(pgrast, pgrast->size); + PG_RETURN_POINTER(pgrast); +} + /* ---------------------------------------------------------------- */ /* Memory allocation / error reporting hooks */ /* ---------------------------------------------------------------- */ Index: raster/rt_core/rt_api.c =================================================================== --- raster/rt_core/rt_api.c (revision 7056) +++ raster/rt_core/rt_api.c (working copy) @@ -3618,3 +3618,56 @@ return rt_raster_add_band(ctx, torast, newband, toindex); } +/** + * Construct a new rt_raster from an existing rt_raster and an array + * of band numbers + * + * @param ctx : context, for thread safety + * @param raster : the source raster + * @param bandNums : array of band numbers to extract from source raster + * and add to the new raster (0 based) + * @param count : number of elements in bandNums + * @return a new rt_raster or 0 on error + * + */ +rt_raster +rt_raster_from_band(rt_context ctx, rt_raster raster, + uint32_t *bandNums, int count) { + rt_raster rast = NULL; + int i = 0; + int idx; + int32_t flag; + + assert(NULL != ctx); + assert(NULL != raster); + assert(NULL != bandNums); + + RASTER_DEBUGF(3, "rt_raster_from_band: source raster has %d bands", + rt_raster_get_num_bands(ctx, raster)); + + /* create new raster */ + rast = rt_raster_new(ctx, raster->width, raster->height); + if (!rast) { + ctx->err("rt_raster_from_band: Out of memory allocating new raster\n"); + return 0; + } + + /* copy bands */ + for (i = 0; i < count; i++) { + idx = bandNums[i]; + flag = rt_raster_copy_band(ctx, rast, raster, idx, i); + + if (flag < 0) { + ctx->err("rt_raster_from_band: Unable to copy band\n"); + rt_raster_destroy(ctx, rast); + return 0; + } + + RASTER_DEBUGF(3, "rt_raster_from_band: band created at index %d", + flag); + } + + RASTER_DEBUGF(3, "rt_raster_from_band: new raster has %d bands", + rt_raster_get_num_bands(ctx, rast)); + return rast; +} Index: raster/rt_core/rt_api.h =================================================================== --- raster/rt_core/rt_api.h (revision 7056) +++ raster/rt_core/rt_api.h (working copy) @@ -737,6 +737,21 @@ int32_t rt_raster_copy_band(rt_context ctx, rt_raster torast, rt_raster fromrast, int fromindex, int toindex); +/** + * Construct a new rt_raster from an existing rt_raster and an array + * of band numbers + * + * @param ctx : context, for thread safety + * @param raster : the source raster + * @param bandNums : array of band numbers to extract from source raster + * and add to the new raster (0 based) + * @param count : number of elements in bandNums + * @return a new rt_raster or 0 on error + * + */ +rt_raster rt_raster_from_band(rt_context ctx, rt_raster raster, + uint32_t *bandNums, int count); + /*- utilities -------------------------------------------------------*/ /* Set of functions to clamp double to int of different size