PostGIS  3.0.6dev-r@@SVN_REVISION@@

◆ RASTER_valueCountCoverage()

Datum RASTER_valueCountCoverage ( PG_FUNCTION_ARGS  )

Definition at line 2614 of file rtpg_statistics.c.

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

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: