return stats;
}
+struct rt_histogram_t {
+ uint32_t count;
+ double proportion;
+
+ double min;
+ double max;
+
+ int inc_min;
+ int inc_max;
+};
+
+/**
+ * Count the distribution of data
+ *
+ * @param stats: a populated stats struct for processing
+ * @param bin_count: the number of bins to group the data by
+ * @param bin_width: the width of each bin as an array
+ * @param bin_width_count: number of values in bin_width
+ * @param right: evaluate bins by (a,b] rather than default [a,b)
+ * @param rtn_count: set to the number of bins being returned
+ *
+ * @return the histogram of the data
+ */
+rt_histogram
+rt_band_get_histogram(rt_bandstats stats,
+ int bin_count, double *bin_width, int bin_width_count,
+ int right, int *rtn_count) {
+ rt_histogram bins = NULL;
+ int init_width = 0;
+ int i;
+ int j;
+ double tmp;
+ double value;
+ int sum = 0;
+
+#if POSTGIS_DEBUG_LEVEL > 0
+ clock_t start, stop;
+ double elapsed = 0;
+#endif
+
+ RASTER_DEBUG(3, "starting");
+#if POSTGIS_DEBUG_LEVEL > 0
+ start = clock();
+#endif
+
+ assert(NULL != stats);
+
+ if (stats->count < 1 || NULL == stats->values) {
+ rterror("rt_util_get_histogram: rt_bandstats object has no value");
+ return NULL;
+ }
+
+ /* bin width must be positive numbers and not zero */
+ if (NULL != bin_width && bin_width_count > 0) {
+ for (i = 0; i < bin_width_count; i++) {
+ if (bin_width[i] <= 0.) {
+ rterror("rt_util_get_histogram: bin_width element is less than or equal to zero");
+ return NULL;
+ }
+ }
+ }
+
+ /* # of bins not provided */
+ if (bin_count <= 0) {
+ /*
+ determine # of bins
+ http://en.wikipedia.org/wiki/Histogram
+
+ all computed bins are assumed to have equal width
+ */
+ /* Square-root choice for stats->count < 30 */
+ if (stats->count < 30)
+ bin_count = ceil(sqrt(stats->count));
+ /* Sturges' formula for stats->count >= 30 */
+ else
+ bin_count = ceil(log2((double) stats->count) + 1.);
+
+ /* bin_width_count provided and bin_width has value */
+ if (bin_width_count > 0 && NULL != bin_width) {
+ /* user has defined something specific */
+ if (bin_width_count > bin_count)
+ bin_count = bin_width_count;
+ else if (bin_width_count > 1) {
+ tmp = 0;
+ for (i = 0; i < bin_width_count; i++) tmp += bin_width[i];
+ bin_count = ceil((stats->max - stats->min) / tmp) * bin_width_count;
+ }
+ else
+ bin_count = ceil((stats->max - stats->min) / bin_width[0]);
+ }
+ /* set bin width count to zero so that one can be calculated */
+ else {
+ bin_width_count = 0;
+ }
+ }
+
+ /* min and max the same */
+ if (stats->min == stats->max)
+ bin_count = 1;
+
+ RASTER_DEBUGF(3, "bin_count = %d", bin_count);
+
+ /* bin count = 1, all values are in one bin */
+ if (bin_count < 2) {
+ bins = rtalloc(sizeof(struct rt_histogram_t));
+ if (NULL == bins) {
+ rterror("rt_util_get_histogram: Unable to allocate memory for histogram");
+ return NULL;
+ }
+
+ bins->count = stats->count;
+ bins->proportion = -1;
+ bins->min = stats->min;
+ bins->max = stats->max;
+ bins->inc_min = bins->inc_max = 1;
+
+ *rtn_count = bin_count;
+ return bins;
+ }
+
+ /* establish bin width */
+ if (bin_width_count == 0) {
+ bin_width_count = 1;
+
+ /* bin_width unallocated */
+ if (NULL == bin_width) {
+ bin_width = rtalloc(sizeof(double));
+ if (NULL == bin_width) {
+ rterror("rt_util_get_histogram: Unable to allocate memory for bin widths");
+ return NULL;
+ }
+ init_width = 1;
+ }
+
+ bin_width[0] = (stats->max - stats->min) / bin_count;
+ }
+
+ /* initialize bins */
+ bins = rtalloc(bin_count * sizeof(struct rt_histogram_t));
+ if (NULL == bins) {
+ rterror("rt_util_get_histogram: Unable to allocate memory for histogram");
+ if (init_width) rtdealloc(bin_width);
+ return NULL;
+ }
+ if (!right)
+ tmp = stats->min;
+ else
+ tmp = stats->max;
+ for (i = 0; i < bin_count;) {
+ for (j = 0; j < bin_width_count; j++) {
+ bins[i].count = 0;
+ bins->proportion = -1;
+
+ if (!right) {
+ bins[i].min = tmp;
+ tmp += bin_width[j];
+ bins[i].max = tmp;
+
+ bins[i].inc_min = 1;
+ bins[i].inc_max = 0;
+ }
+ else {
+ bins[i].max = tmp;
+ tmp -= bin_width[j];
+ bins[i].min = tmp;
+
+ bins[i].inc_min = 0;
+ bins[i].inc_max = 1;
+ }
+
+ i++;
+ }
+ }
+ if (!right) {
+ bins[bin_count - 1].inc_max = 1;
+
+ /* align last bin to the max value */
+ if (bins[bin_count - 1].max < stats->max)
+ bins[bin_count - 1].max = stats->max;
+ }
+ else {
+ bins[bin_count - 1].inc_min = 1;
+
+ /* align first bin to the min value */
+ if (bins[bin_count - 1].min > stats->min)
+ bins[bin_count - 1].min = stats->min;
+ }
+
+ /* process the values */
+ for (i = 0; i < stats->count; i++) {
+ value = stats->values[i];
+
+ /* default, [a, b) */
+ if (!right) {
+ for (j = 0; j < bin_count; j++) {
+ if (
+ (!bins[j].inc_max && value < bins[j].max) ||
+ (bins[j].inc_max && value <= bins[j].max)
+ ) {
+ bins[j].count++;
+ sum++;
+ break;
+ }
+ }
+ }
+ else {
+ for (j = 0; j < bin_count; j++) {
+ if (
+ (!bins[j].inc_min && value > bins[j].min) ||
+ (bins[j].inc_min && value >= bins[j].min)
+ ) {
+ bins[j].count++;
+ sum++;
+ break;
+ }
+ }
+ }
+ }
+
+ /* proportions */
+ for (i = 0; i < bin_count; i++) {
+ bins[i].proportion = ((double) bins[i].count) / sum;
+ if (bin_width_count > 1)
+ bins[i].proportion /= (bins[i].max - bins[i].min);
+ }
+
+#if POSTGIS_DEBUG_LEVEL > 0
+ stop = clock();
+ elapsed = ((double) (stop - start)) / CLOCKS_PER_SEC;
+ RASTER_DEBUGF(3, "elapsed time = %0.4f", elapsed);
+
+ for (j = 0; j < bin_count; j++) {
+ RASTER_DEBUGF(5, "(min, max, inc_min, inc_max, count, sum, proportion) = (%f, %f, %d, %d, %d, %d, %f)",
+ bins[j].min, bins[j].max, bins[j].inc_min, bins[j].inc_max, bins[j].count, sum, bins[j].proportion);
+ }
+#endif
+
+ if (init_width) rtdealloc(bin_width);
+ *rtn_count = bin_count;
+ RASTER_DEBUG(3, "done");
+ return bins;
+}
+
/*- rt_raster --------------------------------------------------------*/
struct rt_raster_serialized_t {
rt_bandstats rt_band_get_summary_stats(rt_band band, int hasnodata,
double sample, int inc_vals);
+/**
+ * Count the distribution of data
+ *
+ * @param stats: a populated stats struct for processing
+ * @param bin_count: the number of bins to group the data by
+ * @param bin_width: the width of each bin as an array
+ * @param bin_width_count: number of values in bin_width
+ * @param right: evaluate bins by (a,b] rather than default [a,b)
+ * @param rtn_count: set to the number of bins being returned
+ *
+ * @return the histogram of the data
+ */
+typedef struct rt_histogram_t* rt_histogram;
+rt_histogram rt_band_get_histogram(rt_bandstats stats,
+ int bin_count, double *bin_widths, int bin_widths_count,
+ int right, int *rtn_count);
+
/*- rt_raster --------------------------------------------------------*/
/* Get summary stats */
Datum RASTER_summaryStats(PG_FUNCTION_ARGS);
+/* get histogram */
+Datum RASTER_histogram(PG_FUNCTION_ARGS);
+
/* Replace function taken from
* http://ubuntuforums.org/showthread.php?s=aa6f015109fd7e4c7e30d2fd8b717497&t=141670&page=3
}
}
+/* get histogram */
+struct rt_histogram_t {
+ uint32_t count;
+ double proportion;
+
+ double min;
+ double max;
+
+ int inc_min;
+ int inc_max;
+};
+
+/**
+ * Returns histogram for a band
+ */
+PG_FUNCTION_INFO_V1(RASTER_histogram);
+Datum RASTER_histogram(PG_FUNCTION_ARGS)
+{
+ FuncCallContext *funcctx;
+ TupleDesc tupdesc;
+ AttInMetadata *attinmeta;
+
+ int count;
+ rt_histogram hist;
+ rt_histogram hist2;
+ int call_cntr;
+ int max_calls;
+
+ /* first call of function */
+ if (SRF_IS_FIRSTCALL()) {
+ MemoryContext oldcontext;
+
+ rt_pgraster *pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ rt_raster raster = NULL;
+ rt_band band = NULL;
+ int32_t bandindex = 0;
+ int num_bands = 0;
+ bool hasnodata = FALSE;
+ double sample = 0;
+ int bin_count = 0;
+ double *bin_width = NULL;
+ int bin_width_count = 0;
+ double width;
+ bool right = FALSE;
+ rt_bandstats stats = NULL;
+
+ int i;
+ int j;
+ int n;
+
+ ArrayType *array;
+ Oid etype;
+ Datum *e;
+ bool *nulls;
+ int16 typlen;
+ bool typbyval;
+ char typalign;
+ int ndims = 1;
+ int *dims;
+ int *lbs;
+
+ /* create a function context for cross-call persistence */
+ funcctx = SRF_FIRSTCALL_INIT();
+
+ /* switch to memory context appropriate for multiple function calls */
+ oldcontext = MemoryContextSwitchTo(funcctx->multi_call_memory_ctx);
+
+ raster = rt_raster_deserialize(pgraster);
+ if (!raster) {
+ elog(ERROR, "RASTER_histogram: Could not deserialize raster");
+ PG_RETURN_NULL();
+ }
+
+ /* band index is 1-based */
+ bandindex = PG_GETARG_INT32(1);
+ num_bands = rt_raster_get_num_bands(raster);
+ if (bandindex < 1 || bandindex > num_bands) {
+ elog(NOTICE, "Invalid band index (must use 1-based). Returning NULL");
+ rt_raster_destroy(raster);
+ SRF_RETURN_DONE(funcctx);
+ }
+ assert(0 <= (bandindex - 1));
+
+ /* hasnodata flag */
+ if (!PG_ARGISNULL(2))
+ hasnodata = PG_GETARG_BOOL(2);
+
+ /* sample % */
+ if (!PG_ARGISNULL(3)) {
+ sample = PG_GETARG_FLOAT8(3);
+ if (sample < 0 || sample > 1) {
+ elog(NOTICE, "Invalid sample percentage (must be between 0 and 1). Returning NULL");
+ rt_raster_destroy(raster);
+ SRF_RETURN_DONE(funcctx);
+ }
+ else if (sample == 0)
+ sample = 1;
+ }
+ else
+ sample = 1;
+
+ /* bin_count */
+ if (!PG_ARGISNULL(4)) {
+ bin_count = PG_GETARG_INT32(4);
+ if (bin_count < 1) bin_count = 0;
+ }
+
+ /* bin_width */
+ if (!PG_ARGISNULL(5)) {
+ array = PG_GETARG_ARRAYTYPE_P(5);
+ etype = ARR_ELEMTYPE(array);
+ get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
+
+ switch (etype) {
+ case FLOAT4OID:
+ case FLOAT8OID:
+ break;
+ default:
+ elog(ERROR, "RASTER_histogram: Invalid data type for width");
+ rt_raster_destroy(raster);
+ PG_RETURN_NULL();
+ break;
+ }
+
+ ndims = ARR_NDIM(array);
+ dims = ARR_DIMS(array);
+ lbs = ARR_LBOUND(array);
+
+ deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
+ &nulls, &n);
+
+ bin_width = palloc(sizeof(double) * n);
+ for (i = 0, j = 0; i < n; i++) {
+ if (nulls[i]) continue;
+
+ switch (etype) {
+ case FLOAT4OID:
+ width = (double) DatumGetFloat4(e[i]);
+ break;
+ case FLOAT8OID:
+ width = (double) DatumGetFloat8(e[i]);
+ break;
+ }
+
+ if (width < 0) {
+ elog(NOTICE, "Invalid value for width (must be greater than 0). Returning NULL");
+ pfree(bin_width);
+ rt_raster_destroy(raster);
+ SRF_RETURN_DONE(funcctx);
+ }
+
+ bin_width[j] = width;
+ POSTGIS_RT_DEBUGF(5, "bin_width[%d] = %f", j, bin_width[j]);
+ j++;
+ }
+ bin_width_count = j;
+
+ if (j < 1) {
+ pfree(bin_width);
+ bin_width = NULL;
+ }
+ }
+
+ /* right */
+ if (!PG_ARGISNULL(6))
+ right = PG_GETARG_BOOL(6);
+
+ /* get band */
+ band = rt_raster_get_band(raster, bandindex - 1);
+ if (!band) {
+ elog(NOTICE, "Could not find raster band of index %d. Returning NULL", bandindex);
+ rt_raster_destroy(raster);
+ SRF_RETURN_DONE(funcctx);
+ }
+
+ /* get stats */
+ stats = rt_band_get_summary_stats(band, (int) hasnodata, sample, 1);
+ rt_band_destroy(band);
+ rt_raster_destroy(raster);
+ if (NULL == stats || NULL == stats->values) {
+ elog(NOTICE, "Could not retrieve summary statistics of raster band of index %d", bandindex);
+ SRF_RETURN_DONE(funcctx);
+ }
+
+ /* get histogram */
+ hist = rt_band_get_histogram(stats, bin_count, bin_width, bin_width_count, right, &count);
+ if (bin_width_count) pfree(bin_width);
+ pfree(stats);
+ if (NULL == hist || !count) {
+ elog(NOTICE, "Could not retrieve histogram of raster band of index %d", bandindex);
+ SRF_RETURN_DONE(funcctx);
+ }
+
+ POSTGIS_RT_DEBUGF(3, "%d bins returned", count);
+
+ /* Store needed information */
+ funcctx->user_fctx = hist;
+
+ /* total number of tuples to be returned */
+ funcctx->max_calls = count;
+
+ /* Build a tuple descriptor for our result type */
+ if (get_call_result_type(fcinfo, NULL, &tupdesc) != TYPEFUNC_COMPOSITE) {
+ ereport(ERROR, (
+ errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
+ errmsg(
+ "function returning record called in context "
+ "that cannot accept type record"
+ )
+ ));
+ }
+
+ /*
+ * generate attribute metadata needed later to produce tuples from raw
+ * C strings
+ */
+ attinmeta = TupleDescGetAttInMetadata(tupdesc);
+ funcctx->attinmeta = attinmeta;
+ MemoryContextSwitchTo(oldcontext);
+ }
+
+ /* stuff done on every call of the function */
+ funcctx = SRF_PERCALL_SETUP();
+
+ call_cntr = funcctx->call_cntr;
+ max_calls = funcctx->max_calls;
+ attinmeta = funcctx->attinmeta;
+ hist2 = funcctx->user_fctx;
+
+ /* do when there is more left to send */
+ if (call_cntr < max_calls) {
+ char **values;
+ HeapTuple tuple;
+ Datum result;
+
+ POSTGIS_RT_DEBUGF(3, "Result %d", call_cntr);
+
+ /*
+ * Prepare a values array for building the returned tuple.
+ * This should be an array of C strings which will
+ * be processed later by the type input functions.
+ */
+ values = (char **) palloc(4 * sizeof(char *));
+
+ values[0] = (char *) palloc(
+ (snprintf(NULL, 0, "%f", hist2[call_cntr].min) + 1) * sizeof(char)
+ );
+ values[1] = (char *) palloc(
+ (snprintf(NULL, 0, "%f", hist2[call_cntr].max) + 1) * sizeof(char)
+ );
+ values[2] = (char *) palloc(
+ (snprintf(NULL, 0, "%d", hist2[call_cntr].count) + 1) * sizeof(char)
+ );
+ values[3] = (char *) palloc(
+ (snprintf(NULL, 0, "%f", hist2[call_cntr].proportion) + 1) * sizeof(char)
+ );
+
+ snprintf(
+ values[0],
+ (snprintf(NULL, 0, "%f", hist2[call_cntr].min) + 1) * sizeof(char),
+ "%f",
+ hist2[call_cntr].min
+ );
+ snprintf(
+ values[1],
+ (snprintf(NULL, 0, "%f", hist2[call_cntr].max) + 1) * sizeof(char),
+ "%f",
+ hist2[call_cntr].max
+ );
+ snprintf(
+ values[2],
+ (snprintf(NULL, 0, "%d", hist2[call_cntr].count) + 1) * sizeof(char),
+ "%d",
+ hist2[call_cntr].count
+ );
+ snprintf(
+ values[3],
+ (snprintf(NULL, 0, "%f", hist2[call_cntr].proportion) + 1) * sizeof(char),
+ "%f",
+ hist2[call_cntr].proportion
+ );
+
+ /* build a tuple */
+ tuple = BuildTupleFromCStrings(attinmeta, values);
+
+ /* make the tuple into a datum */
+ result = HeapTupleGetDatum(tuple);
+
+ /* clean up (this is not really necessary) */
+ pfree(values[3]);
+ pfree(values[2]);
+ pfree(values[1]);
+ pfree(values[0]);
+ pfree(values);
+
+ SRF_RETURN_NEXT(funcctx, result);
+ }
+ /* do when there is no more left */
+ else {
+ pfree(hist2);
+ SRF_RETURN_DONE(funcctx);
+ }
+}
+
/* ---------------------------------------------------------------- */
/* Memory allocation / error reporting hooks */
/* ---------------------------------------------------------------- */
AS $$ SELECT min, max FROM _st_summarystats($1, $2, 1, FALSE, 0.1) $$
LANGUAGE 'SQL' IMMUTABLE STRICT;
+-----------------------------------------------------------------------
+-- ST_Histogram and ST_ApproxHistogram
+-----------------------------------------------------------------------
+CREATE TYPE histogram AS (
+ min double precision,
+ max double precision,
+ count integer,
+ proportion double precision
+);
+
+-- Cannot be strict as "width" can be NULL
+CREATE OR REPLACE FUNCTION _st_histogram(rast raster, nband int, hasnodata boolean, sample_percent double precision, bins int, width double precision[], right boolean)
+ RETURNS SETOF histogram
+ AS 'MODULE_PATHNAME','RASTER_histogram'
+ LANGUAGE 'C' IMMUTABLE;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, hasnodata boolean, bins int, width double precision[], right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, 1, $4, $5, $6) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, hasnodata boolean, bins int, right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, 1, $4, NULL, $5) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, hasnodata boolean, bins int)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, 1, $4, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, hasnodata boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, 1, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, 1, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, 1, FALSE, 1, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, bins int, width double precision[], right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, 1, $3, $4, $5) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, bins int, right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, 1, $3, NULL, $4) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_histogram(rast raster, nband int, bins int)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, 1, $3, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, hasnodata boolean, sample_percent double precision, bins int, width double precision[], right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, $4, $5, $6, $7) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, hasnodata boolean, sample_percent double precision, bins int, right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, $4, $5, NULL, $6) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, hasnodata boolean, sample_percent double precision, bins int)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, $4, $5, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, hasnodata boolean, sample_percent double precision)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, $3, $4, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, sample_percent double precision)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, $3, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, 0.1, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, sample_percent double precision)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, 1, FALSE, $2, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, 1, FALSE, 0.1, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, sample_percent double precision, bins int, width double precision[], right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, $3, $4, $5, $6) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, sample_percent double precision, bins int, right boolean)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, $3, $4, NULL, $5) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, sample_percent double precision, bins int)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, $3, $4, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxhistogram(rast raster, nband int, sample_percent double precision)
+ RETURNS SETOF histogram
+ AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, $3, 0, NULL, FALSE) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
-----------------------------------------------------------------------
-- MapAlgebra
-----------------------------------------------------------------------
double *values;
int sorted;
};
+struct rt_histogram_t {
+ uint32_t count;
+ double proportion;
+
+ double min;
+ double max;
+
+ int inc_min;
+ int inc_max;
+};
static void testBandStats() {
rt_bandstats stats = NULL;
+ rt_histogram histogram = NULL;
+ double bin_width[] = {10000};
+ int count = 0;
rt_raster raster;
rt_band band;
CHECK_EQUALS(stats->min, 1);
CHECK_EQUALS(stats->max, 19998);
+ histogram = (rt_histogram) rt_band_get_histogram(stats, 0, NULL, 0, 0, &count);
+ CHECK(histogram);
+ rtdealloc(histogram);
+
+ histogram = (rt_histogram) rt_band_get_histogram(stats, 0, NULL, 0, 1, &count);
+ CHECK(histogram);
+ rtdealloc(histogram);
+
+ histogram = (rt_histogram) rt_band_get_histogram(stats, 0, bin_width, 1, 0, &count);
+ CHECK(histogram);
+ rtdealloc(histogram);
+
rtdealloc(stats->values);
rtdealloc(stats);
stats = (rt_bandstats) rt_band_get_summary_stats(band, 0, 0.1, 1);
CHECK(stats);
+ histogram = (rt_histogram) rt_band_get_histogram(stats, 0, NULL, 0, 0, &count);
+ CHECK(histogram);
+ rtdealloc(histogram);
+
rtdealloc(stats->values);
rtdealloc(stats);
rt_mean.sql \
rt_stddev.sql \
rt_minmax.sql \
+ rt_histogram.sql \
$(NULL)
TEST_PIXEL = \
--- /dev/null
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ )
+);
+SELECT ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ )
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, TRUE, 0, ARRAY[]::double precision[], FALSE
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, TRUE, 1, FALSE
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, TRUE, 5
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, TRUE
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, 0, ARRAY[5]::double precision[], FALSE
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, 3, FALSE
+);
+SELECT * FROM ST_Histogram(
+ ST_SetValue(
+ ST_SetValue(
+ ST_SetValue(
+ ST_AddBand(
+ ST_MakeEmptyRaster(10, 10, 10, 10, 2, 2, 0, 0,-1)
+ , 1, '64BF', 0, 0
+ )
+ , 1, 1, 1, -10
+ )
+ , 1, 5, 4, 0
+ )
+ , 1, 5, 5, 3.14159
+ ),
+ 1, 5
+);
--- /dev/null
+-10|-3.429205|1|0.5
+-3.429205|3.14159|1|0.5
+(-10,-3.429205,1,0.5)
+(-3.429205,3.14159,1,0.5)
+-10|-8.357301|1|0.01
+-8.357301|-6.714602|0|0
+-6.714602|-5.071904|0|0
+-5.071904|-3.429205|0|0
+-3.429205|-1.786506|0|0
+-1.786506|-0.143807|0|0
+-0.143807|1.498891|98|0.98
+1.498891|3.14159|1|0.01
+-10|3.14159|100|-1
+-10|-7.371682|1|0.01
+-7.371682|-4.743364|0|0
+-4.743364|-2.115046|0|0
+-2.115046|0.513272|98|0.98
+0.513272|3.14159|1|0.01
+-10|-8.357301|1|0.01
+-8.357301|-6.714602|0|0
+-6.714602|-5.071904|0|0
+-5.071904|-3.429205|0|0
+-3.429205|-1.786506|0|0
+-1.786506|-0.143807|0|0
+-0.143807|1.498891|98|0.98
+1.498891|3.14159|1|0.01
+-10|-3.429205|1|0.5
+-3.429205|3.14159|1|0.5
+-10|-5|1|0.5
+-5|0|0|0
+0|5|1|0.5
+-10|-5.61947|1|0.5
+-5.61947|-1.23894|0|0
+-1.23894|3.14159|1|0.5
+-10|-7.371682|1|0.5
+-7.371682|-4.743364|0|0
+-4.743364|-2.115046|0|0
+-2.115046|0.513272|0|0
+0.513272|3.14159|1|0.5