PostGIS  2.1.10dev-r@@SVN_REVISION@@
Datum RASTER_nMapAlgebra ( PG_FUNCTION_ARGS  )

Definition at line 16694 of file rt_pg.c.

References ovdump::band, rtpg_nmapalgebra_arg_t::callback, rtpg_nmapalgebra_arg_t::cextent, rtpg_nmapalgebra_arg_t::distance, ES_NONE, ET_CUSTOM, ET_LAST, ET_SECOND, rtpg_nmapalgebra_arg_t::extenttype, FALSE, rtpg_nmapalgebra_arg_t::hasband, rtpg_nmapalgebra_arg_t::hasnodata, rt_iterator_t::nband, rtpg_nmapalgebra_arg_t::nband, rt_iterator_t::nbnodata, rtpg_nmapalgebra_arg_t::nodataval, rtpg_nmapalgebra_arg_t::numraster, rtpg_nmapalgebra_arg_t::pgcextent, rtpg_nmapalgebra_arg_t::pixtype, POSTGIS_RT_DEBUG, POSTGIS_RT_DEBUGF, PT_END, rtrowdump::raster, rt_iterator_t::raster, rtpg_nmapalgebra_arg_t::raster, rt_band_get_hasnodata_flag(), rt_band_get_min_value(), rt_band_get_nodata(), rt_band_get_pixtype(), rt_pixtype_index_from_name(), rt_pixtype_name(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_band(), rt_raster_is_empty(), rt_raster_iterator(), rt_raster_new(), rt_raster_serialize(), rt_util_extent_type(), rtpg_nmapalgebra_arg_destroy(), rtpg_nmapalgebra_arg_init(), rtpg_nmapalgebra_callback(), rtpg_nmapalgebra_rastbandarg_process(), rtpg_strtoupper(), rtpg_trim(), rt_raster_serialized_t::size, TRUE, rtpg_nmapalgebra_callback_arg::ufc_info, rtpg_nmapalgebra_callback_arg::ufc_noid, rtpg_nmapalgebra_callback_arg::ufc_rettype, and rtpg_nmapalgebra_callback_arg::ufl_info.

16695 {
16696  rtpg_nmapalgebra_arg arg = NULL;
16697  rt_iterator itrset;
16698  int i = 0;
16699  int noerr = 0;
16700  int allnull = 0;
16701  int allempty = 0;
16702  int noband = 0;
16703 
16704  rt_raster raster = NULL;
16705  rt_band band = NULL;
16706  rt_pgraster *pgraster = NULL;
16707 
16708  POSTGIS_RT_DEBUG(3, "Starting...");
16709 
16710  if (PG_ARGISNULL(0))
16711  PG_RETURN_NULL();
16712 
16713  /* init argument struct */
16714  arg = rtpg_nmapalgebra_arg_init();
16715  if (arg == NULL) {
16716  elog(ERROR, "RASTER_nMapAlgebra: Could not initialize argument structure");
16717  PG_RETURN_NULL();
16718  }
16719 
16720  /* let helper function process rastbandarg (0) */
16721  if (!rtpg_nmapalgebra_rastbandarg_process(arg, PG_GETARG_ARRAYTYPE_P(0), &allnull, &allempty, &noband)) {
16723  elog(ERROR, "RASTER_nMapAlgebra: Could not process rastbandarg");
16724  PG_RETURN_NULL();
16725  }
16726 
16727  POSTGIS_RT_DEBUGF(4, "allnull, allempty, noband = %d, %d, %d", allnull, allempty, noband);
16728 
16729  /* all rasters are NULL, return NULL */
16730  if (allnull == arg->numraster) {
16731  elog(NOTICE, "All input rasters are NULL. Returning NULL");
16733  PG_RETURN_NULL();
16734  }
16735 
16736  /* pixel type (2) */
16737  if (!PG_ARGISNULL(2)) {
16738  char *pixtypename = text_to_cstring(PG_GETARG_TEXT_P(2));
16739 
16740  /* Get the pixel type index */
16741  arg->pixtype = rt_pixtype_index_from_name(pixtypename);
16742  if (arg->pixtype == PT_END) {
16744  elog(ERROR, "RASTER_nMapAlgebra: Invalid pixel type: %s", pixtypename);
16745  PG_RETURN_NULL();
16746  }
16747  }
16748 
16749  /* distancex (3) */
16750  if (!PG_ARGISNULL(3))
16751  arg->distance[0] = PG_GETARG_INT32(3);
16752  /* distancey (4) */
16753  if (!PG_ARGISNULL(4))
16754  arg->distance[1] = PG_GETARG_INT32(4);
16755 
16756  if (arg->distance[0] < 0 || arg->distance[1] < 0) {
16758  elog(ERROR, "RASTER_nMapAlgebra: Distance for X and Y axis must be greater than or equal to zero");
16759  PG_RETURN_NULL();
16760  }
16761 
16762  /* extent type (5) */
16763  if (!PG_ARGISNULL(5)) {
16764  char *extenttypename = rtpg_strtoupper(rtpg_trim(text_to_cstring(PG_GETARG_TEXT_P(5))));
16765  arg->extenttype = rt_util_extent_type(extenttypename);
16766  }
16767  POSTGIS_RT_DEBUGF(4, "extenttype: %d", arg->extenttype);
16768 
16769  /* custom extent (6) */
16770  if (arg->extenttype == ET_CUSTOM) {
16771  if (PG_ARGISNULL(6)) {
16772  elog(NOTICE, "Custom extent is NULL. Returning NULL");
16774  PG_RETURN_NULL();
16775  }
16776 
16777  arg->pgcextent = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(6));
16778 
16779  /* only need the raster header */
16781  if (arg->cextent == NULL) {
16783  elog(ERROR, "RASTER_nMapAlgebra: Could not deserialize custom extent");
16784  PG_RETURN_NULL();
16785  }
16786  else if (rt_raster_is_empty(arg->cextent)) {
16787  elog(NOTICE, "Custom extent is an empty raster. Returning empty raster");
16789 
16790  raster = rt_raster_new(0, 0);
16791  if (raster == NULL) {
16792  elog(ERROR, "RASTER_nMapAlgebra: Could not create empty raster");
16793  PG_RETURN_NULL();
16794  }
16795 
16796  pgraster = rt_raster_serialize(raster);
16797  rt_raster_destroy(raster);
16798  if (!pgraster) PG_RETURN_NULL();
16799 
16800  SET_VARSIZE(pgraster, pgraster->size);
16801  PG_RETURN_POINTER(pgraster);
16802  }
16803  }
16804 
16805  noerr = 1;
16806  /* all rasters are empty, return empty raster */
16807  if (allempty == arg->numraster) {
16808  elog(NOTICE, "All input rasters are empty. Returning empty raster");
16809  noerr = 0;
16810  }
16811  /* all rasters don't have indicated band, return empty raster */
16812  else if (noband == arg->numraster) {
16813  elog(NOTICE, "All input rasters do not have bands at indicated indexes. Returning empty raster");
16814  noerr = 0;
16815  }
16816  if (!noerr) {
16818 
16819  raster = rt_raster_new(0, 0);
16820  if (raster == NULL) {
16821  elog(ERROR, "RASTER_nMapAlgebra: Could not create empty raster");
16822  PG_RETURN_NULL();
16823  }
16824 
16825  pgraster = rt_raster_serialize(raster);
16826  rt_raster_destroy(raster);
16827  if (!pgraster) PG_RETURN_NULL();
16828 
16829  SET_VARSIZE(pgraster, pgraster->size);
16830  PG_RETURN_POINTER(pgraster);
16831  }
16832 
16833  /* do regprocedure last (1) */
16834  if (!PG_ARGISNULL(1) || get_fn_expr_argtype(fcinfo->flinfo, 1) == REGPROCEDUREOID) {
16835  POSTGIS_RT_DEBUG(4, "processing callbackfunc");
16836  arg->callback.ufc_noid = PG_GETARG_OID(1);
16837 
16838  /* get function info */
16839  fmgr_info(arg->callback.ufc_noid, &(arg->callback.ufl_info));
16840 
16841  /* function cannot return set */
16842  noerr = 0;
16843  if (arg->callback.ufl_info.fn_retset) {
16844  noerr = 1;
16845  }
16846  /* function should have correct # of args */
16847  else if (arg->callback.ufl_info.fn_nargs != 3) {
16848  noerr = 2;
16849  }
16850 
16851  /* check that callback function return type is supported */
16852  if (
16853  get_func_result_type(
16854  arg->callback.ufc_noid,
16855  &(arg->callback.ufc_rettype),
16856  NULL
16857  ) != TYPEFUNC_SCALAR
16858  ) {
16859  noerr = 3;
16860  }
16861 
16862  if (!(
16863  arg->callback.ufc_rettype == FLOAT8OID ||
16864  arg->callback.ufc_rettype == FLOAT4OID ||
16865  arg->callback.ufc_rettype == INT4OID ||
16866  arg->callback.ufc_rettype == INT2OID
16867  )) {
16868  noerr = 4;
16869  }
16870 
16871  /*
16872  TODO: consider adding checks of the userfunction parameters
16873  should be able to use get_fn_expr_argtype() of fmgr.c
16874  */
16875 
16876  if (noerr != 0) {
16878  switch (noerr) {
16879  case 4:
16880  elog(ERROR, "RASTER_nMapAlgebra: Function provided must return a double precision, float, int or smallint");
16881  break;
16882  case 3:
16883  elog(ERROR, "RASTER_nMapAlgebra: Function provided must return scalar (double precision, float, int, smallint)");
16884  break;
16885  case 2:
16886  elog(ERROR, "RASTER_nMapAlgebra: Function provided must have three input parameters");
16887  break;
16888  case 1:
16889  elog(ERROR, "RASTER_nMapAlgebra: Function provided must return double precision, not resultset");
16890  break;
16891  }
16892  PG_RETURN_NULL();
16893  }
16894 
16895  if (func_volatile(arg->callback.ufc_noid) == 'v')
16896  elog(NOTICE, "Function provided is VOLATILE. Unless required and for best performance, function should be IMMUTABLE or STABLE");
16897 
16898  /* prep function call data */
16899 #if POSTGIS_PGSQL_VERSION > 90
16900  InitFunctionCallInfoData(arg->callback.ufc_info, &(arg->callback.ufl_info), arg->callback.ufl_info.fn_nargs, InvalidOid, NULL, NULL);
16901 #else
16902  InitFunctionCallInfoData(arg->callback.ufc_info, &(arg->callback.ufl_info), arg->callback.ufl_info.fn_nargs, NULL, NULL);
16903 #endif
16904  memset(arg->callback.ufc_info.argnull, FALSE, sizeof(bool) * arg->callback.ufl_info.fn_nargs);
16905 
16906  /* userargs (7) */
16907  if (!PG_ARGISNULL(7))
16908  arg->callback.ufc_info.arg[2] = PG_GETARG_DATUM(7);
16909  else {
16910  if (arg->callback.ufl_info.fn_strict) {
16911  /* build and assign an empty TEXT array */
16912  /* TODO: manually free the empty array? */
16913  arg->callback.ufc_info.arg[2] = PointerGetDatum(
16914  construct_empty_array(TEXTOID)
16915  );
16916  arg->callback.ufc_info.argnull[2] = FALSE;
16917  }
16918  else {
16919  arg->callback.ufc_info.arg[2] = (Datum) NULL;
16920  arg->callback.ufc_info.argnull[2] = TRUE;
16921  }
16922  }
16923  }
16924  else {
16926  elog(ERROR, "RASTER_nMapAlgebra: callbackfunc must be provided");
16927  PG_RETURN_NULL();
16928  }
16929 
16930  /* determine nodataval and possibly pixtype */
16931  /* band to check */
16932  switch (arg->extenttype) {
16933  case ET_LAST:
16934  i = arg->numraster - 1;
16935  break;
16936  case ET_SECOND:
16937  if (arg->numraster > 1) {
16938  i = 1;
16939  break;
16940  }
16941  default:
16942  i = 0;
16943  break;
16944  }
16945  /* find first viable band */
16946  if (!arg->hasband[i]) {
16947  for (i = 0; i < arg->numraster; i++) {
16948  if (arg->hasband[i])
16949  break;
16950  }
16951  if (i >= arg->numraster)
16952  i = arg->numraster - 1;
16953  }
16954  band = rt_raster_get_band(arg->raster[i], arg->nband[i]);
16955 
16956  /* set pixel type if PT_END */
16957  if (arg->pixtype == PT_END)
16958  arg->pixtype = rt_band_get_pixtype(band);
16959 
16960  /* set hasnodata and nodataval */
16961  arg->hasnodata = 1;
16962  if (rt_band_get_hasnodata_flag(band))
16963  rt_band_get_nodata(band, &(arg->nodataval));
16964  else
16965  arg->nodataval = rt_band_get_min_value(band);
16966 
16967  POSTGIS_RT_DEBUGF(4, "pixtype, hasnodata, nodataval: %s, %d, %f", rt_pixtype_name(arg->pixtype), arg->hasnodata, arg->nodataval);
16968 
16969  /* init itrset */
16970  itrset = palloc(sizeof(struct rt_iterator_t) * arg->numraster);
16971  if (itrset == NULL) {
16973  elog(ERROR, "RASTER_nMapAlgebra: Could not allocate memory for iterator arguments");
16974  PG_RETURN_NULL();
16975  }
16976 
16977  /* set itrset */
16978  for (i = 0; i < arg->numraster; i++) {
16979  itrset[i].raster = arg->raster[i];
16980  itrset[i].nband = arg->nband[i];
16981  itrset[i].nbnodata = 1;
16982  }
16983 
16984  /* pass everything to iterator */
16985  noerr = rt_raster_iterator(
16986  itrset, arg->numraster,
16987  arg->extenttype, arg->cextent,
16988  arg->pixtype,
16989  arg->hasnodata, arg->nodataval,
16990  arg->distance[0], arg->distance[1],
16991  &(arg->callback),
16993  &raster
16994  );
16995 
16996  /* cleanup */
16997  pfree(itrset);
16999 
17000  if (noerr != ES_NONE) {
17001  elog(ERROR, "RASTER_nMapAlgebra: Could not run raster iterator function");
17002  PG_RETURN_NULL();
17003  }
17004  else if (raster == NULL)
17005  PG_RETURN_NULL();
17006 
17007  pgraster = rt_raster_serialize(raster);
17008  rt_raster_destroy(raster);
17009 
17010  POSTGIS_RT_DEBUG(3, "Finished");
17011 
17012  if (!pgraster)
17013  PG_RETURN_NULL();
17014 
17015  SET_VARSIZE(pgraster, pgraster->size);
17016  PG_RETURN_POINTER(pgraster);
17017 }
static int rtpg_nmapalgebra_callback(rt_iterator_arg arg, void *userarg, double *value, int *nodata)
Definition: rt_pg.c:16539
static char * rtpg_strtoupper(char *str)
Definition: rt_pg.c:730
void * rt_raster_serialize(rt_raster raster)
Return this raster in serialized form.
Definition: rt_api.c:8158
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_api.c:5387
static char * rtpg_trim(const char *input)
Definition: rt_pg.c:856
const char * rt_pixtype_name(rt_pixtype pixtype)
Definition: rt_api.c:1168
Definition: rt_api.h:184
rt_extenttype rt_util_extent_type(const char *name)
Definition: rt_api.c:302
rt_pixtype rt_band_get_pixtype(rt_band band)
Return pixeltype of this band.
Definition: rt_api.c:1900
rt_raster * raster
Definition: rt_pg.c:16283
tuple band
Definition: ovdump.py:57
rt_pgraster * pgcextent
Definition: rt_pg.c:16296
tuple raster
Be careful!! Zeros function's input parameter can be a (height x width) array, not (width x height): ...
Definition: rtrowdump.py:123
static int rtpg_nmapalgebra_rastbandarg_process(rtpg_nmapalgebra_arg arg, ArrayType *array, int *allnull, int *allempty, int *noband)
Definition: rt_pg.c:16360
rtpg_nmapalgebra_callback_arg callback
Definition: rt_pg.c:16299
#define POSTGIS_RT_DEBUG(level, msg)
Definition: rt_pg.h:58
rt_pixtype pixtype
Definition: rt_pg.c:16289
double rt_band_get_min_value(rt_band band)
Returns the minimal possible value for the band according to the pixel type.
Definition: rt_api.c:3073
rt_errorstate rt_band_get_nodata(rt_band band, double *nodata)
Get NODATA value.
Definition: rt_api.c:3058
static rtpg_nmapalgebra_arg rtpg_nmapalgebra_arg_init()
Definition: rt_pg.c:16302
rt_raster raster
Definition: rt_api.h:2360
uint16_t nband
Definition: rt_api.h:2361
rt_band rt_raster_get_band(rt_raster raster, int n)
Return Nth band, or NULL if unavailable.
Definition: rt_api.c:5686
#define FALSE
Definition: dbfopen.c:169
uint8_t nbnodata
Definition: rt_api.h:2362
Struct definitions.
Definition: rt_api.h:2175
#define POSTGIS_RT_DEBUGF(level, msg,...)
Definition: rt_pg.h:62
int rt_raster_is_empty(rt_raster raster)
Return TRUE if the raster is empty.
Definition: rt_api.c:8550
rt_errorstate rt_raster_iterator(rt_iterator itrset, uint16_t itrcount, rt_extenttype extenttype, rt_raster customextent, rt_pixtype pixtype, uint8_t hasnodata, double nodataval, uint16_t distancex, uint16_t distancey, void *userarg, int(*callback)(rt_iterator_arg arg, void *userarg, double *value, int *nodata), rt_raster *rtnraster)
n-raster iterator.
Definition: rt_api.c:13927
static void rtpg_nmapalgebra_arg_destroy(rtpg_nmapalgebra_arg arg)
Definition: rt_pg.c:16336
rt_raster rt_raster_deserialize(void *serialized, int header_only)
Return a raster from a serialized form.
Definition: rt_api.c:8350
rt_pixtype rt_pixtype_index_from_name(const char *pixname)
Definition: rt_api.c:1138
rt_raster rt_raster_new(uint32_t width, uint32_t height)
Construct a raster with given dimensions.
Definition: rt_api.c:5353
rt_extenttype extenttype
Definition: rt_pg.c:16295
#define TRUE
Definition: dbfopen.c:170
int rt_band_get_hasnodata_flag(rt_band band)
Get hasnodata flag value.
Definition: rt_api.c:2002
FunctionCallInfoData ufc_info
Definition: rt_pg.c:16276

Here is the call graph for this function: