PostGIS  3.4.0dev-r@@SVN_REVISION@@

◆ RASTER_neighborhood()

Datum RASTER_neighborhood ( PG_FUNCTION_ARGS  )

Definition at line 2255 of file rtpg_pixel.c.

2256 {
2257  rt_pgraster *pgraster = NULL;
2258  rt_raster raster = NULL;
2259  rt_band band = NULL;
2260  int bandindex = 1;
2261  int num_bands = 0;
2262  int x = 0;
2263  int y = 0;
2264  int _x = 0;
2265  int _y = 0;
2266  int distance[2] = {0};
2267  bool exclude_nodata_value = TRUE;
2268  double pixval;
2269  int isnodata = 0;
2270 
2271  rt_pixel npixels = NULL;
2272  int count;
2273  double **value2D = NULL;
2274  int **nodata2D = NULL;
2275 
2276  int i = 0;
2277  int j = 0;
2278  int k = 0;
2279  Datum *value1D = NULL;
2280  bool *nodata1D = NULL;
2281  int dim[2] = {0};
2282  int lbound[2] = {1, 1};
2283  ArrayType *mdArray = NULL;
2284 
2285  int16 typlen;
2286  bool typbyval;
2287  char typalign;
2288 
2289  /* pgraster is null, return nothing */
2290  if (PG_ARGISNULL(0))
2291  PG_RETURN_NULL();
2292  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
2293 
2294  raster = rt_raster_deserialize(pgraster, FALSE);
2295  if (!raster) {
2296  PG_FREE_IF_COPY(pgraster, 0);
2297  elog(ERROR, "RASTER_neighborhood: Could not deserialize raster");
2298  PG_RETURN_NULL();
2299  }
2300 
2301  /* band index is 1-based */
2302  if (!PG_ARGISNULL(1))
2303  bandindex = PG_GETARG_INT32(1);
2304  num_bands = rt_raster_get_num_bands(raster);
2305  if (bandindex < 1 || bandindex > num_bands) {
2306  elog(NOTICE, "Invalid band index (must use 1-based). Returning NULL");
2308  PG_FREE_IF_COPY(pgraster, 0);
2309  PG_RETURN_NULL();
2310  }
2311 
2312  /* pixel column, 1-based */
2313  x = PG_GETARG_INT32(2);
2314  _x = x - 1;
2315 
2316  /* pixel row, 1-based */
2317  y = PG_GETARG_INT32(3);
2318  _y = y - 1;
2319 
2320  /* distance X axis */
2321  distance[0] = PG_GETARG_INT32(4);
2322  if (distance[0] < 0) {
2323  elog(NOTICE, "Invalid value for distancex (must be >= zero). Returning NULL");
2325  PG_FREE_IF_COPY(pgraster, 0);
2326  PG_RETURN_NULL();
2327  }
2328  distance[0] = (uint16_t) distance[0];
2329 
2330  /* distance Y axis */
2331  distance[1] = PG_GETARG_INT32(5);
2332  if (distance[1] < 0) {
2333  elog(NOTICE, "Invalid value for distancey (must be >= zero). Returning NULL");
2335  PG_FREE_IF_COPY(pgraster, 0);
2336  PG_RETURN_NULL();
2337  }
2338  distance[1] = (uint16_t) distance[1];
2339 
2340  /* exclude_nodata_value flag */
2341  if (!PG_ARGISNULL(6))
2342  exclude_nodata_value = PG_GETARG_BOOL(6);
2343 
2344  /* get band */
2345  band = rt_raster_get_band(raster, bandindex - 1);
2346  if (!band) {
2347  elog(NOTICE, "Could not find band at index %d. Returning NULL", bandindex);
2349  PG_FREE_IF_COPY(pgraster, 0);
2350  PG_RETURN_NULL();
2351  }
2352 
2353  /* get neighborhood */
2354  count = 0;
2355  npixels = NULL;
2356  if (distance[0] > 0 || distance[1] > 0) {
2358  band,
2359  _x, _y,
2360  distance[0], distance[1],
2361  exclude_nodata_value,
2362  &npixels
2363  );
2364  /* error */
2365  if (count < 0) {
2366  elog(NOTICE, "Could not get the pixel's neighborhood for band at index %d", bandindex);
2367 
2370  PG_FREE_IF_COPY(pgraster, 0);
2371 
2372  PG_RETURN_NULL();
2373  }
2374  }
2375 
2376  /* get pixel's value */
2377  if (
2378  (_x >= 0 && _x < rt_band_get_width(band)) &&
2379  (_y >= 0 && _y < rt_band_get_height(band))
2380  ) {
2381  if (rt_band_get_pixel(
2382  band,
2383  _x, _y,
2384  &pixval,
2385  &isnodata
2386  ) != ES_NONE) {
2387  elog(NOTICE, "Could not get the pixel of band at index %d. Returning NULL", bandindex);
2390  PG_FREE_IF_COPY(pgraster, 0);
2391  PG_RETURN_NULL();
2392  }
2393  }
2394  /* outside band extent, set to NODATA */
2395  else {
2396  /* has NODATA, use NODATA */
2399  /* no NODATA, use min possible value */
2400  else
2402  isnodata = 1;
2403  }
2404  POSTGIS_RT_DEBUGF(4, "pixval: %f", pixval);
2405 
2406 
2407  /* add pixel to neighborhood */
2408  count++;
2409  if (count > 1)
2410  npixels = (rt_pixel) repalloc(npixels, sizeof(struct rt_pixel_t) * count);
2411  else
2412  npixels = (rt_pixel) palloc(sizeof(struct rt_pixel_t));
2413  if (npixels == NULL) {
2414 
2417  PG_FREE_IF_COPY(pgraster, 0);
2418 
2419  elog(ERROR, "RASTER_neighborhood: Could not reallocate memory for neighborhood");
2420  PG_RETURN_NULL();
2421  }
2422  npixels[count - 1].x = _x;
2423  npixels[count - 1].y = _y;
2424  npixels[count - 1].nodata = 1;
2425  npixels[count - 1].value = pixval;
2426 
2427  /* set NODATA */
2428  if (!exclude_nodata_value || !isnodata) {
2429  npixels[count - 1].nodata = 0;
2430  }
2431 
2432  /* free unnecessary stuff */
2435  PG_FREE_IF_COPY(pgraster, 0);
2436 
2437  /* convert set of rt_pixel to 2D array */
2438  /* dim is passed with element 0 being Y-axis and element 1 being X-axis */
2440  npixels, count, NULL,
2441  _x, _y,
2442  distance[0], distance[1],
2443  &value2D,
2444  &nodata2D,
2445  &(dim[1]), &(dim[0])
2446  );
2447  pfree(npixels);
2448  if (count != ES_NONE) {
2449  elog(NOTICE, "Could not create 2D array of neighborhood");
2450  PG_RETURN_NULL();
2451  }
2452 
2453  /* 1D arrays for values and nodata from 2D arrays */
2454  value1D = palloc(sizeof(Datum) * dim[0] * dim[1]);
2455  nodata1D = palloc(sizeof(bool) * dim[0] * dim[1]);
2456 
2457  if (value1D == NULL || nodata1D == NULL) {
2458 
2459  for (i = 0; i < dim[0]; i++) {
2460  pfree(value2D[i]);
2461  pfree(nodata2D[i]);
2462  }
2463  pfree(value2D);
2464  pfree(nodata2D);
2465 
2466  elog(ERROR, "RASTER_neighborhood: Could not allocate memory for return 2D array");
2467  PG_RETURN_NULL();
2468  }
2469 
2470  /* copy values from 2D arrays to 1D arrays */
2471  k = 0;
2472  /* Y-axis */
2473  for (i = 0; i < dim[0]; i++) {
2474  /* X-axis */
2475  for (j = 0; j < dim[1]; j++) {
2476  nodata1D[k] = (bool) nodata2D[i][j];
2477  if (!nodata1D[k])
2478  value1D[k] = Float8GetDatum(value2D[i][j]);
2479  else
2480  value1D[k] = PointerGetDatum(NULL);
2481 
2482  k++;
2483  }
2484  }
2485 
2486  /* no more need for 2D arrays */
2487  for (i = 0; i < dim[0]; i++) {
2488  pfree(value2D[i]);
2489  pfree(nodata2D[i]);
2490  }
2491  pfree(value2D);
2492  pfree(nodata2D);
2493 
2494  /* info about the type of item in the multi-dimensional array (float8). */
2495  get_typlenbyvalalign(FLOAT8OID, &typlen, &typbyval, &typalign);
2496 
2497  mdArray = construct_md_array(
2498  value1D, nodata1D,
2499  2, dim, lbound,
2500  FLOAT8OID,
2501  typlen, typbyval, typalign
2502  );
2503 
2504  pfree(value1D);
2505  pfree(nodata1D);
2506 
2507  PG_RETURN_ARRAYTYPE_P(mdArray);
2508 }
#define TRUE
Definition: dbfopen.c:73
#define FALSE
Definition: dbfopen.c:72
uint16_t rt_band_get_width(rt_band band)
Return width of this band.
Definition: rt_band.c:640
int rt_band_get_hasnodata_flag(rt_band band)
Get hasnodata flag value.
Definition: rt_band.c:674
rt_errorstate rt_band_get_pixel(rt_band band, int x, int y, double *value, int *nodata)
Get pixel value.
Definition: rt_band.c:1376
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:86
struct rt_pixel_t * rt_pixel
Definition: librtcore.h:149
double rt_band_get_min_value(rt_band band)
Returns the minimal possible value for the band according to the pixel type.
Definition: rt_band.c:1902
@ ES_NONE
Definition: librtcore.h:182
void rt_band_destroy(rt_band band)
Destroy a raster band.
Definition: rt_band.c:340
uint16_t rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:376
rt_errorstate rt_band_get_nodata(rt_band band, double *nodata)
Get NODATA value.
Definition: rt_band.c:1887
rt_errorstate rt_pixel_set_to_array(rt_pixel npixel, uint32_t count, rt_mask mask, int x, int y, uint16_t distancex, uint16_t distancey, double ***value, int ***nodata, int *dimx, int *dimy)
Definition: rt_pixel.c:288
uint32_t rt_band_get_nearest_pixel(rt_band band, int x, int y, uint16_t distancex, uint16_t distancey, int exclude_nodata_value, rt_pixel *npixels)
Get nearest pixel(s) with value (not NODATA) to specified pixel.
Definition: rt_band.c:1529
uint16_t rt_band_get_height(rt_band band)
Return height of this band.
Definition: rt_band.c:649
rt_raster rt_raster_deserialize(void *serialized, int header_only)
Return a raster from a serialized form.
Definition: rt_serialize.c:725
rt_band rt_raster_get_band(rt_raster raster, int bandNum)
Return Nth band, or NULL if unavailable.
Definition: rt_raster.c:385
static double distance(double x1, double y1, double x2, double y2)
Definition: lwtree.c:1032
int count
Definition: genraster.py:57
band
Definition: ovdump.py:58
pixval
Definition: pixval.py:94
Definition: pixval.py:1
raster
Be careful!! Zeros function's input parameter can be a (height x width) array, not (width x height): ...
Definition: rtrowdump.py:121
#define POSTGIS_RT_DEBUGF(level, msg,...)
Definition: rtpostgis.h:69
double value
Definition: librtcore.h:2491
uint8_t nodata
Definition: librtcore.h:2490
Struct definitions.
Definition: librtcore.h:2403

References ovdump::band, genraster::count, distance(), ES_NONE, FALSE, rt_pixel_t::nodata, pixval::pixval, POSTGIS_RT_DEBUGF, rtrowdump::raster, rt_band_destroy(), rt_band_get_hasnodata_flag(), rt_band_get_height(), rt_band_get_min_value(), rt_band_get_nearest_pixel(), rt_band_get_nodata(), rt_band_get_pixel(), rt_band_get_width(), rt_pixel_set_to_array(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_band(), rt_raster_get_num_bands(), TRUE, rt_pixel_t::value, rt_pixel_t::x, pixval::x, rt_pixel_t::y, and pixval::y.

Here is the call graph for this function: