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

◆ RASTER_histogram()

Datum RASTER_histogram ( PG_FUNCTION_ARGS  )

Definition at line 921 of file rtpg_statistics.c.

922{
923 FuncCallContext *funcctx;
924 TupleDesc tupdesc;
925
926 int i;
927 rt_histogram hist;
928 rt_histogram hist2;
929 int call_cntr;
930 int max_calls;
931
932 /* first call of function */
933 if (SRF_IS_FIRSTCALL()) {
934 MemoryContext oldcontext;
935
936 rt_pgraster *pgraster = NULL;
937 rt_raster raster = NULL;
938 rt_band band = NULL;
939 int32_t bandindex = 1;
940 int num_bands = 0;
941 bool exclude_nodata_value = TRUE;
942 double sample = 0;
943 uint32_t bin_count = 0;
944 double *bin_width = NULL;
945 uint32_t bin_width_count = 0;
946 double width = 0;
947 bool right = FALSE;
948 double min = 0;
949 double max = 0;
950 rt_bandstats stats = NULL;
951 uint32_t count;
952
953 int j;
954 int n;
955
956 ArrayType *array;
957 Oid etype;
958 Datum *e;
959 bool *nulls;
960 int16 typlen;
961 bool typbyval;
962 char typalign;
963
964 POSTGIS_RT_DEBUG(3, "RASTER_histogram: Starting");
965
966 /* create a function context for cross-call persistence */
967 funcctx = SRF_FIRSTCALL_INIT();
968
969 /* switch to memory context appropriate for multiple function calls */
970 oldcontext = MemoryContextSwitchTo(funcctx->multi_call_memory_ctx);
971
972 /* pgraster is null, return nothing */
973 if (PG_ARGISNULL(0)) {
974 MemoryContextSwitchTo(oldcontext);
975 SRF_RETURN_DONE(funcctx);
976 }
977 pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
978
979 raster = rt_raster_deserialize(pgraster, FALSE);
980 if (!raster) {
981 PG_FREE_IF_COPY(pgraster, 0);
982 MemoryContextSwitchTo(oldcontext);
983 elog(ERROR, "RASTER_histogram: Cannot deserialize raster");
984 SRF_RETURN_DONE(funcctx);
985 }
986
987 /* band index is 1-based */
988 if (!PG_ARGISNULL(1))
989 bandindex = PG_GETARG_INT32(1);
990 num_bands = rt_raster_get_num_bands(raster);
991 if (bandindex < 1 || bandindex > num_bands) {
992 elog(NOTICE, "Invalid band index (must use 1-based). Returning NULL");
993 rt_raster_destroy(raster);
994 PG_FREE_IF_COPY(pgraster, 0);
995 MemoryContextSwitchTo(oldcontext);
996 SRF_RETURN_DONE(funcctx);
997 }
998
999 /* exclude_nodata_value flag */
1000 if (!PG_ARGISNULL(2))
1001 exclude_nodata_value = PG_GETARG_BOOL(2);
1002
1003 /* sample % */
1004 if (!PG_ARGISNULL(3)) {
1005 sample = PG_GETARG_FLOAT8(3);
1006 if (sample < 0 || sample > 1) {
1007 elog(NOTICE, "Invalid sample percentage (must be between 0 and 1). Returning NULL");
1008 rt_raster_destroy(raster);
1009 PG_FREE_IF_COPY(pgraster, 0);
1010 MemoryContextSwitchTo(oldcontext);
1011 SRF_RETURN_DONE(funcctx);
1012 }
1013 else if (FLT_EQ(sample, 0.0))
1014 sample = 1;
1015 }
1016 else
1017 sample = 1;
1018
1019 /* bin_count */
1020 if (!PG_ARGISNULL(4)) {
1021 bin_count = PG_GETARG_INT32(4);
1022 if (bin_count < 1) bin_count = 0;
1023 }
1024
1025 /* bin_width */
1026 if (!PG_ARGISNULL(5)) {
1027 array = PG_GETARG_ARRAYTYPE_P(5);
1028 etype = ARR_ELEMTYPE(array);
1029 get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
1030
1031 switch (etype) {
1032 case FLOAT4OID:
1033 case FLOAT8OID:
1034 break;
1035 default:
1036 rt_raster_destroy(raster);
1037 PG_FREE_IF_COPY(pgraster, 0);
1038 MemoryContextSwitchTo(oldcontext);
1039 elog(ERROR, "RASTER_histogram: Invalid data type for width");
1040 SRF_RETURN_DONE(funcctx);
1041 break;
1042 }
1043
1044 deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
1045 &nulls, &n);
1046
1047 bin_width = palloc(sizeof(double) * n);
1048 for (i = 0, j = 0; i < n; i++) {
1049 if (nulls[i]) continue;
1050
1051 switch (etype) {
1052 case FLOAT4OID:
1053 width = (double) DatumGetFloat4(e[i]);
1054 break;
1055 case FLOAT8OID:
1056 width = (double) DatumGetFloat8(e[i]);
1057 break;
1058 }
1059
1060 if (width < 0 || FLT_EQ(width, 0.0)) {
1061 elog(NOTICE, "Invalid value for width (must be greater than 0). Returning NULL");
1062 pfree(bin_width);
1063 rt_raster_destroy(raster);
1064 PG_FREE_IF_COPY(pgraster, 0);
1065 MemoryContextSwitchTo(oldcontext);
1066 SRF_RETURN_DONE(funcctx);
1067 }
1068
1069 bin_width[j] = width;
1070 POSTGIS_RT_DEBUGF(5, "bin_width[%d] = %f", j, bin_width[j]);
1071 j++;
1072 }
1073 bin_width_count = j;
1074
1075 if (j < 1) {
1076 pfree(bin_width);
1077 bin_width = NULL;
1078 }
1079 }
1080
1081 /* right */
1082 if (!PG_ARGISNULL(6))
1083 right = PG_GETARG_BOOL(6);
1084
1085 /* min */
1086 if (!PG_ARGISNULL(7)) min = PG_GETARG_FLOAT8(7);
1087
1088 /* max */
1089 if (!PG_ARGISNULL(8)) max = PG_GETARG_FLOAT8(8);
1090
1091 /* get band */
1092 band = rt_raster_get_band(raster, bandindex - 1);
1093 if (!band) {
1094 elog(NOTICE, "Cannot find band at index %d. Returning NULL", bandindex);
1095 rt_raster_destroy(raster);
1096 PG_FREE_IF_COPY(pgraster, 0);
1097 MemoryContextSwitchTo(oldcontext);
1098 SRF_RETURN_DONE(funcctx);
1099 }
1100
1101 /* get stats */
1102 stats = rt_band_get_summary_stats(band, (int) exclude_nodata_value, sample, 1, NULL, NULL, NULL);
1103 rt_band_destroy(band);
1104 rt_raster_destroy(raster);
1105 PG_FREE_IF_COPY(pgraster, 0);
1106 if (NULL == stats || NULL == stats->values) {
1107 elog(NOTICE, "Cannot compute summary statistics for band at index %d", bandindex);
1108 MemoryContextSwitchTo(oldcontext);
1109 SRF_RETURN_DONE(funcctx);
1110 }
1111 else if (stats->count < 1) {
1112 elog(NOTICE, "Cannot compute histogram for band at index %d as the band has no values", bandindex);
1113 MemoryContextSwitchTo(oldcontext);
1114 SRF_RETURN_DONE(funcctx);
1115 }
1116
1117 /* get histogram */
1118 hist = rt_band_get_histogram(stats, (uint32_t)bin_count, bin_width, bin_width_count, right, min, max, &count);
1119 if (bin_width_count) pfree(bin_width);
1120 pfree(stats);
1121 if (NULL == hist || !count) {
1122 elog(NOTICE, "Cannot compute histogram for band at index %d", bandindex);
1123 MemoryContextSwitchTo(oldcontext);
1124 SRF_RETURN_DONE(funcctx);
1125 }
1126
1127 POSTGIS_RT_DEBUGF(3, "%d bins returned", count);
1128
1129 /* Store needed information */
1130 funcctx->user_fctx = hist;
1131
1132 /* total number of tuples to be returned */
1133 funcctx->max_calls = count;
1134
1135 /* Build a tuple descriptor for our result type */
1136 if (get_call_result_type(fcinfo, NULL, &tupdesc) != TYPEFUNC_COMPOSITE) {
1137 ereport(ERROR, (
1138 errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
1139 errmsg(
1140 "function returning record called in context "
1141 "that cannot accept type record"
1142 )
1143 ));
1144 }
1145
1146 BlessTupleDesc(tupdesc);
1147 funcctx->tuple_desc = tupdesc;
1148
1149 MemoryContextSwitchTo(oldcontext);
1150 }
1151
1152 /* stuff done on every call of the function */
1153 funcctx = SRF_PERCALL_SETUP();
1154
1155 call_cntr = funcctx->call_cntr;
1156 max_calls = funcctx->max_calls;
1157 tupdesc = funcctx->tuple_desc;
1158 hist2 = funcctx->user_fctx;
1159
1160 /* do when there is more left to send */
1161 if (call_cntr < max_calls) {
1162 Datum values[VALUES_LENGTH];
1163 bool nulls[VALUES_LENGTH];
1164 HeapTuple tuple;
1165 Datum result;
1166
1167 POSTGIS_RT_DEBUGF(3, "Result %d", call_cntr);
1168
1169 memset(nulls, FALSE, sizeof(bool) * VALUES_LENGTH);
1170
1171 values[0] = Float8GetDatum(hist2[call_cntr].min);
1172 values[1] = Float8GetDatum(hist2[call_cntr].max);
1173 values[2] = Int64GetDatum(hist2[call_cntr].count);
1174 values[3] = Float8GetDatum(hist2[call_cntr].percent);
1175
1176 /* build a tuple */
1177 tuple = heap_form_tuple(tupdesc, values, nulls);
1178
1179 /* make the tuple into a datum */
1180 result = HeapTupleGetDatum(tuple);
1181
1182 SRF_RETURN_NEXT(funcctx, result);
1183 }
1184 /* do when there is no more left */
1185 else {
1186 pfree(hist2);
1187 SRF_RETURN_DONE(funcctx);
1188 }
1189}
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_band_get_histogram
rt_histogram rt_band_get_histogram(rt_bandstats stats, uint32_t bin_count, double *bin_widths, uint32_t bin_widths_count, int right, double min, double max, uint32_t *rtn_count)
Count the distribution of data.
Definition rt_statistics.c:412
rt_band_get_summary_stats
rt_bandstats rt_band_get_summary_stats(rt_band band, int exclude_nodata_value, double sample, int inc_vals, uint64_t *cK, double *cM, double *cQ)
Compute summary statistics for a band.
Definition rt_statistics.c:110
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
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
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_bandstats_t::count
uint32_t count
Definition librtcore.h:2361
rt_bandstats_t::values
double * values
Definition librtcore.h:2369
rt_bandstats_t
Definition librtcore.h:2359
rt_histogram_t
Definition librtcore.h:2374
rt_raster_serialized_t
Struct definitions.
Definition librtcore.h:2251
rt_raster_t
Definition librtcore.h:2285

References rt_bandstats_t::count, FALSE, FLT_EQ, POSTGIS_RT_DEBUG, POSTGIS_RT_DEBUGF, rt_band_destroy(), rt_band_get_histogram(), rt_band_get_summary_stats(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_band(), rt_raster_get_num_bands(), TRUE, rt_bandstats_t::values, and VALUES_LENGTH.

Here is the call graph for this function: