PostGIS  2.5.7dev-r@@SVN_REVISION@@

◆ RASTER_valueCountCoverage()

Datum RASTER_valueCountCoverage ( PG_FUNCTION_ARGS  )

Definition at line 2610 of file rtpg_statistics.c.

2610  {
2611  FuncCallContext *funcctx;
2612  TupleDesc tupdesc;
2613 
2614  uint32_t i;
2615  uint64_t covcount = 0;
2616  uint64_t covtotal = 0;
2617  rt_valuecount covvcnts = NULL;
2618  rt_valuecount covvcnts2;
2619  int call_cntr;
2620  int max_calls;
2621 
2622  POSTGIS_RT_DEBUG(3, "RASTER_valueCountCoverage: Starting");
2623 
2624  /* first call of function */
2625  if (SRF_IS_FIRSTCALL()) {
2626  MemoryContext oldcontext;
2627 
2628  text *tablenametext = NULL;
2629  char *tablename = NULL;
2630  text *colnametext = NULL;
2631  char *colname = NULL;
2632  int32_t bandindex = 1;
2633  bool exclude_nodata_value = TRUE;
2634  double *search_values = NULL;
2635  uint32_t search_values_count = 0;
2636  double roundto = 0;
2637 
2638  int len = 0;
2639  char *sql = NULL;
2640  int spi_result;
2641  Portal portal;
2642  HeapTuple tuple;
2643  Datum datum;
2644  bool isNull = FALSE;
2645  rt_pgraster *pgraster = NULL;
2646  rt_raster raster = NULL;
2647  rt_band band = NULL;
2648  int num_bands = 0;
2649  uint32_t count;
2650  uint32_t total;
2651  rt_valuecount vcnts = NULL;
2652  int exists = 0;
2653 
2654  uint32_t j;
2655  int n;
2656 
2657  ArrayType *array;
2658  Oid etype;
2659  Datum *e;
2660  bool *nulls;
2661  int16 typlen;
2662  bool typbyval;
2663  char typalign;
2664 
2665  /* create a function context for cross-call persistence */
2666  funcctx = SRF_FIRSTCALL_INIT();
2667 
2668  /* switch to memory context appropriate for multiple function calls */
2669  oldcontext = MemoryContextSwitchTo(funcctx->multi_call_memory_ctx);
2670 
2671  /* tablename is null, return null */
2672  if (PG_ARGISNULL(0)) {
2673  elog(NOTICE, "Table name must be provided");
2674  MemoryContextSwitchTo(oldcontext);
2675  SRF_RETURN_DONE(funcctx);
2676  }
2677  tablenametext = PG_GETARG_TEXT_P(0);
2678  tablename = text_to_cstring(tablenametext);
2679  if (!strlen(tablename)) {
2680  elog(NOTICE, "Table name must be provided");
2681  MemoryContextSwitchTo(oldcontext);
2682  SRF_RETURN_DONE(funcctx);
2683  }
2684  POSTGIS_RT_DEBUGF(3, "tablename = %s", tablename);
2685 
2686  /* column name is null, return null */
2687  if (PG_ARGISNULL(1)) {
2688  elog(NOTICE, "Column name must be provided");
2689  MemoryContextSwitchTo(oldcontext);
2690  SRF_RETURN_DONE(funcctx);
2691  }
2692  colnametext = PG_GETARG_TEXT_P(1);
2693  colname = text_to_cstring(colnametext);
2694  if (!strlen(colname)) {
2695  elog(NOTICE, "Column name must be provided");
2696  MemoryContextSwitchTo(oldcontext);
2697  SRF_RETURN_DONE(funcctx);
2698  }
2699  POSTGIS_RT_DEBUGF(3, "colname = %s", colname);
2700 
2701  /* band index is 1-based */
2702  if (!PG_ARGISNULL(2))
2703  bandindex = PG_GETARG_INT32(2);
2704 
2705  /* exclude_nodata_value flag */
2706  if (!PG_ARGISNULL(3))
2707  exclude_nodata_value = PG_GETARG_BOOL(3);
2708 
2709  /* search values */
2710  if (!PG_ARGISNULL(4)) {
2711  array = PG_GETARG_ARRAYTYPE_P(4);
2712  etype = ARR_ELEMTYPE(array);
2713  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
2714 
2715  switch (etype) {
2716  case FLOAT4OID:
2717  case FLOAT8OID:
2718  break;
2719  default:
2720  MemoryContextSwitchTo(oldcontext);
2721  elog(ERROR, "RASTER_valueCountCoverage: Invalid data type for values");
2722  SRF_RETURN_DONE(funcctx);
2723  break;
2724  }
2725 
2726  deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
2727  &nulls, &n);
2728 
2729  search_values = palloc(sizeof(double) * n);
2730  for (i = 0, j = 0; i < (uint32_t) n; i++) {
2731  if (nulls[i]) continue;
2732 
2733  switch (etype) {
2734  case FLOAT4OID:
2735  search_values[j] = (double) DatumGetFloat4(e[i]);
2736  break;
2737  case FLOAT8OID:
2738  search_values[j] = (double) DatumGetFloat8(e[i]);
2739  break;
2740  }
2741 
2742  POSTGIS_RT_DEBUGF(5, "search_values[%d] = %f", j, search_values[j]);
2743  j++;
2744  }
2745  search_values_count = j;
2746 
2747  if (j < 1) {
2748  pfree(search_values);
2749  search_values = NULL;
2750  }
2751  }
2752 
2753  /* roundto */
2754  if (!PG_ARGISNULL(5)) {
2755  roundto = PG_GETARG_FLOAT8(5);
2756  if (roundto < 0.) roundto = 0;
2757  }
2758 
2759  /* iterate through rasters of coverage */
2760  /* connect to database */
2761  spi_result = SPI_connect();
2762  if (spi_result != SPI_OK_CONNECT) {
2763 
2764  if (search_values_count) pfree(search_values);
2765 
2766  MemoryContextSwitchTo(oldcontext);
2767  elog(ERROR, "RASTER_valueCountCoverage: Cannot connect to database using SPI");
2768  SRF_RETURN_DONE(funcctx);
2769  }
2770 
2771  /* create sql */
2772  len = sizeof(char) * (strlen("SELECT \"\" FROM \"\" WHERE \"\" IS NOT NULL") + (strlen(colname) * 2) + strlen(tablename) + 1);
2773  sql = (char *) palloc(len);
2774  if (NULL == sql) {
2775 
2776  if (SPI_tuptable) SPI_freetuptable(SPI_tuptable);
2777  SPI_finish();
2778 
2779  if (search_values_count) pfree(search_values);
2780 
2781  MemoryContextSwitchTo(oldcontext);
2782  elog(ERROR, "RASTER_valueCountCoverage: Cannot allocate memory for sql");
2783  SRF_RETURN_DONE(funcctx);
2784  }
2785 
2786  /* get cursor */
2787  snprintf(sql, len, "SELECT \"%s\" FROM \"%s\" WHERE \"%s\" IS NOT NULL", colname, tablename, colname);
2788  POSTGIS_RT_DEBUGF(3, "RASTER_valueCountCoverage: %s", sql);
2789  portal = SPI_cursor_open_with_args(
2790  "coverage",
2791  sql,
2792  0, NULL,
2793  NULL, NULL,
2794  TRUE, 0
2795  );
2796  pfree(sql);
2797 
2798  /* process resultset */
2799  SPI_cursor_fetch(portal, TRUE, 1);
2800  while (SPI_processed == 1 && SPI_tuptable != NULL) {
2801  tupdesc = SPI_tuptable->tupdesc;
2802  tuple = SPI_tuptable->vals[0];
2803 
2804  datum = SPI_getbinval(tuple, tupdesc, 1, &isNull);
2805  if (SPI_result == SPI_ERROR_NOATTRIBUTE) {
2806 
2807  SPI_freetuptable(SPI_tuptable);
2808  SPI_cursor_close(portal);
2809  SPI_finish();
2810 
2811  if (NULL != covvcnts) pfree(covvcnts);
2812  if (search_values_count) pfree(search_values);
2813 
2814  MemoryContextSwitchTo(oldcontext);
2815  elog(ERROR, "RASTER_valueCountCoverage: Cannot get raster of coverage");
2816  SRF_RETURN_DONE(funcctx);
2817  }
2818  else if (isNull) {
2819  SPI_cursor_fetch(portal, TRUE, 1);
2820  continue;
2821  }
2822 
2823  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(datum);
2824 
2825  raster = rt_raster_deserialize(pgraster, FALSE);
2826  if (!raster) {
2827 
2828  SPI_freetuptable(SPI_tuptable);
2829  SPI_cursor_close(portal);
2830  SPI_finish();
2831 
2832  if (NULL != covvcnts) pfree(covvcnts);
2833  if (search_values_count) pfree(search_values);
2834 
2835  MemoryContextSwitchTo(oldcontext);
2836  elog(ERROR, "RASTER_valueCountCoverage: Cannot deserialize raster");
2837  SRF_RETURN_DONE(funcctx);
2838  }
2839 
2840  /* inspect number of bands*/
2841  num_bands = rt_raster_get_num_bands(raster);
2842  if (bandindex < 1 || bandindex > num_bands) {
2843  elog(NOTICE, "Invalid band index (must use 1-based). Returning NULL");
2844 
2845  rt_raster_destroy(raster);
2846 
2847  SPI_freetuptable(SPI_tuptable);
2848  SPI_cursor_close(portal);
2849  SPI_finish();
2850 
2851  if (NULL != covvcnts) pfree(covvcnts);
2852  if (search_values_count) pfree(search_values);
2853 
2854  MemoryContextSwitchTo(oldcontext);
2855  SRF_RETURN_DONE(funcctx);
2856  }
2857 
2858  /* get band */
2859  band = rt_raster_get_band(raster, bandindex - 1);
2860  if (!band) {
2861  elog(NOTICE, "Cannot find band at index %d. Returning NULL", bandindex);
2862 
2863  rt_raster_destroy(raster);
2864 
2865  SPI_freetuptable(SPI_tuptable);
2866  SPI_cursor_close(portal);
2867  SPI_finish();
2868 
2869  if (NULL != covvcnts) pfree(covvcnts);
2870  if (search_values_count) pfree(search_values);
2871 
2872  MemoryContextSwitchTo(oldcontext);
2873  SRF_RETURN_DONE(funcctx);
2874  }
2875 
2876  /* get counts of values */
2877  vcnts = rt_band_get_value_count(band, (int) exclude_nodata_value, search_values, search_values_count, roundto, &total, &count);
2878  rt_band_destroy(band);
2879  rt_raster_destroy(raster);
2880  if (NULL == vcnts || !count) {
2881  elog(NOTICE, "Cannot count the values for band at index %d", bandindex);
2882 
2883  SPI_freetuptable(SPI_tuptable);
2884  SPI_cursor_close(portal);
2885  SPI_finish();
2886 
2887  if (NULL != covvcnts) free(covvcnts);
2888  if (search_values_count) pfree(search_values);
2889 
2890  MemoryContextSwitchTo(oldcontext);
2891  SRF_RETURN_DONE(funcctx);
2892  }
2893 
2894  POSTGIS_RT_DEBUGF(3, "%d value counts returned", count);
2895 
2896  if (NULL == covvcnts) {
2897  covvcnts = (rt_valuecount) SPI_palloc(sizeof(struct rt_valuecount_t) * count);
2898  if (NULL == covvcnts) {
2899 
2900  SPI_freetuptable(SPI_tuptable);
2901  SPI_cursor_close(portal);
2902  SPI_finish();
2903 
2904  if (search_values_count) pfree(search_values);
2905 
2906  MemoryContextSwitchTo(oldcontext);
2907  elog(ERROR, "RASTER_valueCountCoverage: Cannot allocate memory for value counts of coverage");
2908  SRF_RETURN_DONE(funcctx);
2909  }
2910 
2911  for (i = 0; i < count; i++) {
2912  covvcnts[i].value = vcnts[i].value;
2913  covvcnts[i].count = vcnts[i].count;
2914  covvcnts[i].percent = -1;
2915  }
2916 
2917  covcount = count;
2918  }
2919  else {
2920  for (i = 0; i < count; i++) {
2921  exists = 0;
2922 
2923  for (j = 0; j < covcount; j++) {
2924  if (FLT_EQ(vcnts[i].value, covvcnts[j].value)) {
2925  exists = 1;
2926  break;
2927  }
2928  }
2929 
2930  if (exists) {
2931  covvcnts[j].count += vcnts[i].count;
2932  }
2933  else {
2934  covcount++;
2935  covvcnts = SPI_repalloc(covvcnts, sizeof(struct rt_valuecount_t) * covcount);
2936  if (!covvcnts) {
2937  SPI_freetuptable(SPI_tuptable);
2938  SPI_cursor_close(portal);
2939  SPI_finish();
2940 
2941  if (search_values_count) pfree(search_values);
2942 
2943  MemoryContextSwitchTo(oldcontext);
2944  elog(ERROR, "RASTER_valueCountCoverage: Cannot change allocated memory for value counts of coverage");
2945  SRF_RETURN_DONE(funcctx);
2946  }
2947 
2948  covvcnts[covcount - 1].value = vcnts[i].value;
2949  covvcnts[covcount - 1].count = vcnts[i].count;
2950  covvcnts[covcount - 1].percent = -1;
2951  }
2952  }
2953  }
2954 
2955  covtotal += total;
2956 
2957  pfree(vcnts);
2958 
2959  /* next record */
2960  SPI_cursor_fetch(portal, TRUE, 1);
2961  }
2962 
2963  if (SPI_tuptable) SPI_freetuptable(SPI_tuptable);
2964  SPI_cursor_close(portal);
2965  SPI_finish();
2966 
2967  if (search_values_count) pfree(search_values);
2968 
2969  /* compute percentages */
2970  for (i = 0; i < covcount; i++) {
2971  covvcnts[i].percent = (double) covvcnts[i].count / covtotal;
2972  }
2973 
2974  /* Store needed information */
2975  funcctx->user_fctx = covvcnts;
2976 
2977  /* total number of tuples to be returned */
2978  funcctx->max_calls = covcount;
2979 
2980  /* Build a tuple descriptor for our result type */
2981  if (get_call_result_type(fcinfo, NULL, &tupdesc) != TYPEFUNC_COMPOSITE) {
2982  ereport(ERROR, (
2983  errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
2984  errmsg(
2985  "function returning record called in context "
2986  "that cannot accept type record"
2987  )
2988  ));
2989  }
2990 
2991  BlessTupleDesc(tupdesc);
2992  funcctx->tuple_desc = tupdesc;
2993 
2994  MemoryContextSwitchTo(oldcontext);
2995  }
2996 
2997  /* stuff done on every call of the function */
2998  funcctx = SRF_PERCALL_SETUP();
2999 
3000  call_cntr = funcctx->call_cntr;
3001  max_calls = funcctx->max_calls;
3002  tupdesc = funcctx->tuple_desc;
3003  covvcnts2 = funcctx->user_fctx;
3004 
3005  /* do when there is more left to send */
3006  if (call_cntr < max_calls) {
3007  Datum values[VALUES_LENGTH];
3008  bool nulls[VALUES_LENGTH];
3009  HeapTuple tuple;
3010  Datum result;
3011 
3012  POSTGIS_RT_DEBUGF(3, "Result %d", call_cntr);
3013 
3014  memset(nulls, FALSE, sizeof(bool) * VALUES_LENGTH);
3015 
3016  values[0] = Float8GetDatum(covvcnts2[call_cntr].value);
3017  values[1] = UInt32GetDatum(covvcnts2[call_cntr].count);
3018  values[2] = Float8GetDatum(covvcnts2[call_cntr].percent);
3019 
3020  /* build a tuple */
3021  tuple = heap_form_tuple(tupdesc, values, nulls);
3022 
3023  /* make the tuple into a datum */
3024  result = HeapTupleGetDatum(tuple);
3025 
3026  SRF_RETURN_NEXT(funcctx, result);
3027  }
3028  /* do when there is no more left */
3029  else {
3030  pfree(covvcnts2);
3031  SRF_RETURN_DONE(funcctx);
3032  }
3033 }
TRUE
#define TRUE
Definition: dbfopen.c:169
FALSE
#define FALSE
Definition: dbfopen.c:168
rt_raster_destroy
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:82
FLT_EQ
#define FLT_EQ(x, y)
Definition: librtcore.h:2234
rt_band_destroy
void rt_band_destroy(rt_band band)
Destroy a raster band.
Definition: rt_band.c:340
rt_raster_get_num_bands
uint16_t rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:372
rt_valuecount
struct rt_valuecount_t * rt_valuecount
Definition: librtcore.h:153
rt_band_get_value_count
rt_valuecount rt_band_get_value_count(rt_band band, int exclude_nodata_value, double *search_values, uint32_t search_values_count, double roundto, uint32_t *rtn_total, uint32_t *rtn_count)
Count the number of times provided value(s) occur in the band.
Definition: rt_statistics.c:1608
rt_raster_deserialize
rt_raster rt_raster_deserialize(void *serialized, int header_only)
Return a raster from a serialized form.
Definition: rt_serialize.c:725
rt_raster_get_band
rt_band rt_raster_get_band(rt_raster raster, int bandNum)
Return Nth band, or NULL if unavailable.
Definition: rt_raster.c:381
free
void free(void *)
genraster.value
int value
Definition: genraster.py:61
genraster.count
int count
Definition: genraster.py:56
ovdump.band
band
Definition: ovdump.py:57
rtrowdump.raster
raster
Be careful!! Zeros function's input parameter can be a (height x width) array, not (width x height): ...
Definition: rtrowdump.py:121
text_to_cstring
char * text_to_cstring(const text *textptr)
VALUES_LENGTH
#define VALUES_LENGTH
Definition: rtpg_statistics.c:2396
POSTGIS_RT_DEBUG
#define POSTGIS_RT_DEBUG(level, msg)
Definition: rtpostgis.h:61
POSTGIS_RT_DEBUGF
#define POSTGIS_RT_DEBUGF(level, msg,...)
Definition: rtpostgis.h:65
rt_band_t
Definition: librtcore.h:2313
rt_raster_serialized_t
Struct definitions.
Definition: librtcore.h:2250
rt_raster_t
Definition: librtcore.h:2284
rt_valuecount_t::count
uint32_t count
Definition: librtcore.h:2425
rt_valuecount_t::percent
double percent
Definition: librtcore.h:2426
rt_valuecount_t::value
double value
Definition: librtcore.h:2424
rt_valuecount_t
Definition: librtcore.h:2423
uint32_t
unsigned int uint32_t
Definition: uthash.h:78

References ovdump::band, rt_valuecount_t::count, genraster::count, FALSE, FLT_EQ, free(), rt_valuecount_t::percent, POSTGIS_RT_DEBUG, POSTGIS_RT_DEBUGF, rtrowdump::raster, rt_band_destroy(), rt_band_get_value_count(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_band(), rt_raster_get_num_bands(), rtgdalraster::sql, text_to_cstring(), TRUE, rt_valuecount_t::value, genraster::value, and VALUES_LENGTH.

Here is the call graph for this function: