PostGIS 3.0.6dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches

◆ 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.count
int count
Definition genraster.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:67
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 quantile_llist::count, rt_valuecount_t::count, FALSE, FLT_EQ, free(), rt_valuecount_t::percent, POSTGIS_RT_DEBUG, POSTGIS_RT_DEBUGF, rt_band_destroy(), rt_band_get_value_count(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_band(), rt_raster_get_num_bands(), text_to_cstring(), TRUE, rt_valuecount_t::value, and VALUES_LENGTH.

Here is the call graph for this function: