PostGIS  3.1.6dev-r@@SVN_REVISION@@

◆ RASTER_neighborhood()

Datum RASTER_neighborhood ( PG_FUNCTION_ARGS  )

Definition at line 2056 of file rtpg_pixel.c.

2057 {
2058  rt_pgraster *pgraster = NULL;
2059  rt_raster raster = NULL;
2060  rt_band band = NULL;
2061  int bandindex = 1;
2062  int num_bands = 0;
2063  int x = 0;
2064  int y = 0;
2065  int _x = 0;
2066  int _y = 0;
2067  int distance[2] = {0};
2068  bool exclude_nodata_value = TRUE;
2069  double pixval;
2070  int isnodata = 0;
2071 
2072  rt_pixel npixels = NULL;
2073  int count;
2074  double **value2D = NULL;
2075  int **nodata2D = NULL;
2076 
2077  int i = 0;
2078  int j = 0;
2079  int k = 0;
2080  Datum *value1D = NULL;
2081  bool *nodata1D = NULL;
2082  int dim[2] = {0};
2083  int lbound[2] = {1, 1};
2084  ArrayType *mdArray = NULL;
2085 
2086  int16 typlen;
2087  bool typbyval;
2088  char typalign;
2089 
2090  /* pgraster is null, return nothing */
2091  if (PG_ARGISNULL(0))
2092  PG_RETURN_NULL();
2093  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
2094 
2095  raster = rt_raster_deserialize(pgraster, FALSE);
2096  if (!raster) {
2097  PG_FREE_IF_COPY(pgraster, 0);
2098  elog(ERROR, "RASTER_neighborhood: Could not deserialize raster");
2099  PG_RETURN_NULL();
2100  }
2101 
2102  /* band index is 1-based */
2103  if (!PG_ARGISNULL(1))
2104  bandindex = PG_GETARG_INT32(1);
2105  num_bands = rt_raster_get_num_bands(raster);
2106  if (bandindex < 1 || bandindex > num_bands) {
2107  elog(NOTICE, "Invalid band index (must use 1-based). Returning NULL");
2109  PG_FREE_IF_COPY(pgraster, 0);
2110  PG_RETURN_NULL();
2111  }
2112 
2113  /* pixel column, 1-based */
2114  x = PG_GETARG_INT32(2);
2115  _x = x - 1;
2116 
2117  /* pixel row, 1-based */
2118  y = PG_GETARG_INT32(3);
2119  _y = y - 1;
2120 
2121  /* distance X axis */
2122  distance[0] = PG_GETARG_INT32(4);
2123  if (distance[0] < 0) {
2124  elog(NOTICE, "Invalid value for distancex (must be >= zero). Returning NULL");
2126  PG_FREE_IF_COPY(pgraster, 0);
2127  PG_RETURN_NULL();
2128  }
2129  distance[0] = (uint16_t) distance[0];
2130 
2131  /* distance Y axis */
2132  distance[1] = PG_GETARG_INT32(5);
2133  if (distance[1] < 0) {
2134  elog(NOTICE, "Invalid value for distancey (must be >= zero). Returning NULL");
2136  PG_FREE_IF_COPY(pgraster, 0);
2137  PG_RETURN_NULL();
2138  }
2139  distance[1] = (uint16_t) distance[1];
2140 
2141  /* exclude_nodata_value flag */
2142  if (!PG_ARGISNULL(6))
2143  exclude_nodata_value = PG_GETARG_BOOL(6);
2144 
2145  /* get band */
2146  band = rt_raster_get_band(raster, bandindex - 1);
2147  if (!band) {
2148  elog(NOTICE, "Could not find band at index %d. Returning NULL", bandindex);
2150  PG_FREE_IF_COPY(pgraster, 0);
2151  PG_RETURN_NULL();
2152  }
2153 
2154  /* get neighborhood */
2155  count = 0;
2156  npixels = NULL;
2157  if (distance[0] > 0 || distance[1] > 0) {
2159  band,
2160  _x, _y,
2161  distance[0], distance[1],
2162  exclude_nodata_value,
2163  &npixels
2164  );
2165  /* error */
2166  if (count < 0) {
2167  elog(NOTICE, "Could not get the pixel's neighborhood for band at index %d", bandindex);
2168 
2171  PG_FREE_IF_COPY(pgraster, 0);
2172 
2173  PG_RETURN_NULL();
2174  }
2175  }
2176 
2177  /* get pixel's value */
2178  if (
2179  (_x >= 0 && _x < rt_band_get_width(band)) &&
2180  (_y >= 0 && _y < rt_band_get_height(band))
2181  ) {
2182  if (rt_band_get_pixel(
2183  band,
2184  _x, _y,
2185  &pixval,
2186  &isnodata
2187  ) != ES_NONE) {
2188  elog(NOTICE, "Could not get the pixel of band at index %d. Returning NULL", bandindex);
2191  PG_FREE_IF_COPY(pgraster, 0);
2192  PG_RETURN_NULL();
2193  }
2194  }
2195  /* outside band extent, set to NODATA */
2196  else {
2197  /* has NODATA, use NODATA */
2200  /* no NODATA, use min possible value */
2201  else
2203  isnodata = 1;
2204  }
2205  POSTGIS_RT_DEBUGF(4, "pixval: %f", pixval);
2206 
2207 
2208  /* add pixel to neighborhood */
2209  count++;
2210  if (count > 1)
2211  npixels = (rt_pixel) repalloc(npixels, sizeof(struct rt_pixel_t) * count);
2212  else
2213  npixels = (rt_pixel) palloc(sizeof(struct rt_pixel_t));
2214  if (npixels == NULL) {
2215 
2218  PG_FREE_IF_COPY(pgraster, 0);
2219 
2220  elog(ERROR, "RASTER_neighborhood: Could not reallocate memory for neighborhood");
2221  PG_RETURN_NULL();
2222  }
2223  npixels[count - 1].x = _x;
2224  npixels[count - 1].y = _y;
2225  npixels[count - 1].nodata = 1;
2226  npixels[count - 1].value = pixval;
2227 
2228  /* set NODATA */
2229  if (!exclude_nodata_value || !isnodata) {
2230  npixels[count - 1].nodata = 0;
2231  }
2232 
2233  /* free unnecessary stuff */
2236  PG_FREE_IF_COPY(pgraster, 0);
2237 
2238  /* convert set of rt_pixel to 2D array */
2239  /* dim is passed with element 0 being Y-axis and element 1 being X-axis */
2241  npixels, count, NULL,
2242  _x, _y,
2243  distance[0], distance[1],
2244  &value2D,
2245  &nodata2D,
2246  &(dim[1]), &(dim[0])
2247  );
2248  pfree(npixels);
2249  if (count != ES_NONE) {
2250  elog(NOTICE, "Could not create 2D array of neighborhood");
2251  PG_RETURN_NULL();
2252  }
2253 
2254  /* 1D arrays for values and nodata from 2D arrays */
2255  value1D = palloc(sizeof(Datum) * dim[0] * dim[1]);
2256  nodata1D = palloc(sizeof(bool) * dim[0] * dim[1]);
2257 
2258  if (value1D == NULL || nodata1D == NULL) {
2259 
2260  for (i = 0; i < dim[0]; i++) {
2261  pfree(value2D[i]);
2262  pfree(nodata2D[i]);
2263  }
2264  pfree(value2D);
2265  pfree(nodata2D);
2266 
2267  elog(ERROR, "RASTER_neighborhood: Could not allocate memory for return 2D array");
2268  PG_RETURN_NULL();
2269  }
2270 
2271  /* copy values from 2D arrays to 1D arrays */
2272  k = 0;
2273  /* Y-axis */
2274  for (i = 0; i < dim[0]; i++) {
2275  /* X-axis */
2276  for (j = 0; j < dim[1]; j++) {
2277  nodata1D[k] = (bool) nodata2D[i][j];
2278  if (!nodata1D[k])
2279  value1D[k] = Float8GetDatum(value2D[i][j]);
2280  else
2281  value1D[k] = PointerGetDatum(NULL);
2282 
2283  k++;
2284  }
2285  }
2286 
2287  /* no more need for 2D arrays */
2288  for (i = 0; i < dim[0]; i++) {
2289  pfree(value2D[i]);
2290  pfree(nodata2D[i]);
2291  }
2292  pfree(value2D);
2293  pfree(nodata2D);
2294 
2295  /* info about the type of item in the multi-dimensional array (float8). */
2296  get_typlenbyvalalign(FLOAT8OID, &typlen, &typbyval, &typalign);
2297 
2298  mdArray = construct_md_array(
2299  value1D, nodata1D,
2300  2, dim, lbound,
2301  FLOAT8OID,
2302  typlen, typbyval, typalign
2303  );
2304 
2305  pfree(value1D);
2306  pfree(nodata1D);
2307 
2308  PG_RETURN_ARRAYTYPE_P(mdArray);
2309 }
#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:1221
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:82
struct rt_pixel_t * rt_pixel
Definition: librtcore.h:147
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:1745
@ ES_NONE
Definition: librtcore.h:180
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:372
rt_errorstate rt_band_get_nodata(rt_band band, double *nodata)
Get NODATA value.
Definition: rt_band.c:1730
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:286
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:1374
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:381
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:65
double value
Definition: librtcore.h:2339
uint8_t nodata
Definition: librtcore.h:2338
Struct definitions.
Definition: librtcore.h:2251

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: