return (float)fmin(fmax((value), -FLT_MAX), FLT_MAX);
}
+/* quicksort */
+#define swap(x, y) { double t; t = x; x = y; y = t; }
+#define order(x, y) if (x > y) swap(x, y)
+
+static double pivot(double *left, double *right) {
+ double l, m, r, *p;
+
+ l = *left;
+ m = *(left + (right - left) / 2);
+ r = *right;
+
+ /* order */
+ order(l, m);
+ order(l, r);
+ order(m, r);
+
+ /* pivot is higher of two values */
+ if (l < m) return m;
+ if (m < r) return r;
+
+ /* find pivot that isn't left */
+ for (p = left + 1; p <= right; ++p) {
+ if (*p != *left)
+ return (*p < *left) ? *left : *p;
+ }
+
+ /* all values are same */
+ return -1;
+}
+
+static double *partition(double *left, double *right, double pivot) {
+ while (left <= right) {
+ while (*left < pivot) ++left;
+ while (*right >= pivot) --right;
+
+ if (left < right) {
+ swap(*left, *right);
+ ++left;
+ --right;
+ }
+ }
+
+ return left;
+}
+
+static void quicksort(double *left, double *right) {
+ double p = pivot(left, right);
+ double *pos;
+
+ if (p != -1) {
+ pos = partition(left, right, p);
+ quicksort(left, pos - 1);
+ quicksort(pos, right);
+ }
+}
/*- rt_context -------------------------------------------------------*/
return bins;
}
+struct rt_quantile_t {
+ double quantile;
+ double value;
+};
+
+/**
+ * Compute the default set of or requested quantiles for a set of data
+ * the quantile formula used is same as Excel and R default method
+ *
+ * @param stats: a populated stats struct for processing
+ * @param quantiles: the quantiles to be computed
+ * @param quantiles_count: the number of quantiles to be computed
+ * @param rtn_count: set to the number of quantiles being returned
+ *
+ * @return the default set of or requested quantiles for a band
+ */
+rt_quantile
+rt_band_get_quantiles(rt_bandstats stats,
+ double *quantiles, int quantiles_count, int *rtn_count) {
+ rt_quantile rtn;
+ int init_quantiles = 0;
+ int i = 0;
+ double h;
+ int hl;
+
+#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_quantile: rt_bandstats object has no value");
+ return NULL;
+ }
+
+ /* quantiles not provided */
+ if (NULL == quantiles) {
+ /* quantile count not specified, default to quartiles */
+ if (quantiles_count < 2)
+ quantiles_count = 5;
+
+ quantiles = rtalloc(sizeof(double) * quantiles_count);
+ init_quantiles = 1;
+ if (NULL == quantiles) {
+ rterror("rt_util_get_quantile: Unable to allocate memory for quantile input");
+ return NULL;
+ }
+
+ quantiles_count--;
+ for (i = 0; i <= quantiles_count; i++)
+ quantiles[i] = ((double) i) / quantiles_count;
+ quantiles_count++;
+ }
+
+ /* check quantiles */
+ for (i = 0; i < quantiles_count; i++) {
+ if (quantiles[i] < 0. || quantiles[i] > 1.) {
+ rterror("rt_util_get_quantile: Quantile value not between 0 and 1");
+ if (init_quantiles) rtdealloc(quantiles);
+ return NULL;
+ }
+ }
+ quicksort(quantiles, quantiles + quantiles_count - 1);
+
+ /* initialize rt_quantile */
+ rtn = rtalloc(sizeof(struct rt_quantile_t) * quantiles_count);
+ if (NULL == rtn) {
+ rterror("rt_util_get_quantile: Unable to allocate memory for quantile output");
+ if (init_quantiles) rtdealloc(quantiles);
+ return NULL;
+ }
+
+ /* sort values */
+ if (!stats->sorted) {
+ quicksort(stats->values, stats->values + stats->count - 1);
+ stats->sorted = 1;
+ }
+
+ /*
+ make quantiles
+
+ formula is that used in R (method 7) and Excel from
+ http://en.wikipedia.org/wiki/Quantile
+ */
+ for (i = 0; i < quantiles_count; i++) {
+ rtn[i].quantile = quantiles[i];
+
+ h = ((stats->count - 1.) * quantiles[i]) + 1.;
+ hl = floor(h);
+
+ /* h greater than hl, do full equation */
+ if (h > hl)
+ rtn[i].value = stats->values[hl - 1] + ((h - hl) * (stats->values[hl] - stats->values[hl - 1]));
+ /* shortcut as second part of equation is zero */
+ else
+ rtn[i].value = stats->values[hl - 1];
+ }
+
+#if POSTGIS_DEBUG_LEVEL > 0
+ stop = clock();
+ elapsed = ((double) (stop - start)) / CLOCKS_PER_SEC;
+ RASTER_DEBUGF(3, "elapsed time = %0.4f", elapsed);
+#endif
+
+ *rtn_count = quantiles_count;
+ RASTER_DEBUG(3, "done");
+ return rtn;
+}
+
/*- rt_raster --------------------------------------------------------*/
struct rt_raster_serialized_t {
int bin_count, double *bin_widths, int bin_widths_count,
int right, int *rtn_count);
+/**
+ * Compute the default set of or requested quantiles for a set of data
+ * the quantile formula used is same as Excel and R default method
+ *
+ * @param stats: a populated stats struct for processing
+ * @param quantiles: the quantiles to be computed
+ * @param quantiles_count: the number of quantiles to be computed
+ * @param rtn_count: the number of quantiles being returned
+ *
+ * @return the default set of or requested quantiles for a band
+ */
+typedef struct rt_quantile_t* rt_quantile;
+rt_quantile rt_band_get_quantiles(rt_bandstats stats,
+ double *quantiles, int quantiles_count, int *rtn_count);
+
/*- rt_raster --------------------------------------------------------*/
/* get histogram */
Datum RASTER_histogram(PG_FUNCTION_ARGS);
+/* get quantiles */
+Datum RASTER_quantile(PG_FUNCTION_ARGS);
+
/* Replace function taken from
* http://ubuntuforums.org/showthread.php?s=aa6f015109fd7e4c7e30d2fd8b717497&t=141670&page=3
}
}
+/* get quantiles */
+struct rt_quantile_t {
+ double quantile;
+ double value;
+};
+
+/**
+ * Returns quantiles for a band
+ */
+PG_FUNCTION_INFO_V1(RASTER_quantile);
+Datum RASTER_quantile(PG_FUNCTION_ARGS)
+{
+ FuncCallContext *funcctx;
+ TupleDesc tupdesc;
+ AttInMetadata *attinmeta;
+
+ int count;
+ rt_quantile quant;
+ rt_quantile quant2;
+ 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;
+ double *quantiles = NULL;
+ int quantiles_count = 0;
+ double quantile;
+ 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_quantile: 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;
+
+ /* quantiles */
+ if (!PG_ARGISNULL(4)) {
+ array = PG_GETARG_ARRAYTYPE_P(4);
+ etype = ARR_ELEMTYPE(array);
+ get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
+
+ switch (etype) {
+ case FLOAT4OID:
+ case FLOAT8OID:
+ break;
+ default:
+ elog(ERROR, "RASTER_quantile: Invalid data type for quanitiles");
+ 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);
+
+ quantiles = palloc(sizeof(double) * n);
+ for (i = 0, j = 0; i < n; i++) {
+ if (nulls[i]) continue;
+
+ switch (etype) {
+ case FLOAT4OID:
+ quantile = (double) DatumGetFloat4(e[i]);
+ break;
+ case FLOAT8OID:
+ quantile = (double) DatumGetFloat8(e[i]);
+ break;
+ }
+
+ if (quantile < 0 || quantile > 1) {
+ elog(NOTICE, "Invalid value for quantile (must be between 0 and 1). Returning NULL");
+ pfree(quantiles);
+ rt_raster_destroy(raster);
+ SRF_RETURN_DONE(funcctx);
+ }
+
+ quantiles[j] = quantile;
+ POSTGIS_RT_DEBUGF(5, "quantiles[%d] = %f", j, quantiles[j]);
+ j++;
+ }
+ quantiles_count = j;
+
+ if (j < 1) {
+ pfree(quantiles);
+ quantiles = NULL;
+ }
+ }
+
+ /* 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 quantiles */
+ quant = rt_band_get_quantiles(stats, quantiles, quantiles_count, &count);
+ if (quantiles_count) pfree(quantiles);
+ pfree(stats);
+ if (NULL == quant || !count) {
+ elog(NOTICE, "Could not retrieve quantiles of raster band of index %d", bandindex);
+ SRF_RETURN_DONE(funcctx);
+ }
+
+ POSTGIS_RT_DEBUGF(3, "%d quantiles returned", count);
+
+ /* Store needed information */
+ funcctx->user_fctx = quant;
+
+ /* 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;
+ quant2 = 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(2 * sizeof(char *));
+
+ values[0] = (char *) palloc(
+ (snprintf(NULL, 0, "%f", quant2[call_cntr].quantile) + 1) * sizeof(char)
+ );
+ values[1] = (char *) palloc(
+ (snprintf(NULL, 0, "%f", quant2[call_cntr].value) + 1) * sizeof(char)
+ );
+
+ snprintf(
+ values[0],
+ (snprintf(NULL, 0, "%f", quant2[call_cntr].quantile) + 1) * sizeof(char),
+ "%f",
+ quant2[call_cntr].quantile
+ );
+ snprintf(
+ values[1],
+ (snprintf(NULL, 0, "%f", quant2[call_cntr].value) + 1) * sizeof(char),
+ "%f",
+ quant2[call_cntr].value
+ );
+
+ /* 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[1]);
+ pfree(values[0]);
+ pfree(values);
+
+ SRF_RETURN_NEXT(funcctx, result);
+ }
+ /* do when there is no more left */
+ else {
+ pfree(quant2);
+ SRF_RETURN_DONE(funcctx);
+ }
+}
+
/* ---------------------------------------------------------------- */
/* Memory allocation / error reporting hooks */
/* ---------------------------------------------------------------- */
AS $$ SELECT min, max, count, proportion FROM _st_histogram($1, $2, FALSE, $3, 0, NULL, FALSE) $$
LANGUAGE 'sql' IMMUTABLE STRICT;
+-----------------------------------------------------------------------
+-- ST_Quantile and ST_ApproxQuantile
+-----------------------------------------------------------------------
+CREATE TYPE quantile AS (
+ quantile double precision,
+ value double precision
+);
+
+-- Cannot be strict as "quantile" can be NULL
+CREATE OR REPLACE FUNCTION _st_quantile(rast raster, nband int, hasnodata boolean, sample_percent double precision, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS 'MODULE_PATHNAME','RASTER_quantile'
+ LANGUAGE 'C' IMMUTABLE;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, nband int, hasnodata boolean, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, $3, 1, $4) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, nband int, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, 1, $3) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, nband int, hasnodata boolean)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, $3, 1, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, nband int)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, 1, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, 1, $2) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, 1, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, nband int, hasnodata boolean, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, $3, 1, ARRAY[$4]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, nband int, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, 1, ARRAY[$3]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, hasnodata boolean, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, $2, 1, ARRAY[$3]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_quantile(rast raster, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, 1, ARRAY[$2]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, nband int, hasnodata boolean, sample_percent double precision, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, $3, $4, $5) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, nband int, sample_percent double precision, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, $3, $4) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, nband int, sample_percent double precision)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, $3, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, sample_percent double precision, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, $2, $3) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, sample_percent double precision)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, $2, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, quantiles double precision[])
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, 0.1, $2) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, nband int)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, 0.1, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster)
+ RETURNS SETOF quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, 0.1, NULL) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, nband int, hasnodata boolean, sample_percent double precision, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, $3, $4, ARRAY[$5]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, nband int, sample_percent double precision, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, $2, FALSE, $3, ARRAY[$4]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, sample_percent double precision, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, FALSE, $2, ARRAY[$3]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
+CREATE OR REPLACE FUNCTION st_approxquantile(rast raster, hasnodata boolean, quantile double precision)
+ RETURNS quantile
+ AS $$ SELECT quantile, value FROM _st_quantile($1, 1, $2, 1, ARRAY[$3]::double precision[]) $$
+ LANGUAGE 'sql' IMMUTABLE STRICT;
+
-----------------------------------------------------------------------
-- MapAlgebra
-----------------------------------------------------------------------
int inc_min;
int inc_max;
};
+struct rt_quantile_t {
+ double quantile;
+ double value;
+};
static void testBandStats() {
rt_bandstats stats = NULL;
rt_histogram histogram = NULL;
double bin_width[] = {10000};
+ double quantiles[] = {0.1, 0.3, 0.5, 0.7, 0.9};
+ rt_quantile quantile = NULL;
int count = 0;
rt_raster raster;
CHECK_EQUALS(stats->min, 1);
CHECK_EQUALS(stats->max, 19998);
+ quantile = (rt_quantile) rt_band_get_quantiles(stats, NULL, 0, &count);
+ CHECK(quantile);
+ rtdealloc(quantile);
+
histogram = (rt_histogram) rt_band_get_histogram(stats, 0, NULL, 0, 0, &count);
CHECK(histogram);
rtdealloc(histogram);
stats = (rt_bandstats) rt_band_get_summary_stats(band, 0, 0.1, 1);
CHECK(stats);
+ quantile = (rt_quantile) rt_band_get_quantiles(stats, NULL, 0, &count);
+ CHECK(quantile);
+ rtdealloc(quantile);
+
+ quantile = (rt_quantile) rt_band_get_quantiles(stats, quantiles, 5, &count);
+ CHECK(quantile);
+ CHECK((count == 5));
+ rtdealloc(quantile);
+
histogram = (rt_histogram) rt_band_get_histogram(stats, 0, NULL, 0, 0, &count);
CHECK(histogram);
rtdealloc(histogram);
CHECK_EQUALS(stats->min, 0);
CHECK_EQUALS(stats->max, 19998);
+ quantile = (rt_quantile) rt_band_get_quantiles(stats, NULL, 0, &count);
+ CHECK(quantile);
+ rtdealloc(quantile);
+
rtdealloc(stats->values);
rtdealloc(stats);
stats = (rt_bandstats) rt_band_get_summary_stats(band, 1, 0.1, 1);
CHECK(stats);
+ quantile = (rt_quantile) rt_band_get_quantiles(stats, NULL, 0, &count);
+ CHECK(quantile);
+ rtdealloc(quantile);
+
rtdealloc(stats->values);
rtdealloc(stats);
rt_stddev.sql \
rt_minmax.sql \
rt_histogram.sql \
+ rt_quantile.sql \
$(NULL)
TEST_PIXEL = \
--- /dev/null
+SELECT * FROM ST_Quantile(
+ 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, ARRAY[0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1]::double precision[]
+);
+SELECT * FROM ST_Quantile(
+ 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, ARRAY[0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1]::double precision[]
+);
+SELECT * FROM ST_Quantile(
+ 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_Quantile(
+ 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_Quantile(
+ 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
+ ),
+ ARRAY[0.05, 0.95]::double precision[]
+);
+SELECT * FROM ST_Quantile(
+ 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_Quantile(
+ 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.05
+);
+SELECT * FROM ST_Quantile(
+ 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.95
+);
+SELECT * FROM ST_Quantile(
+ 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
+ ),
+ TRUE, 0.7
+);
+SELECT * FROM ST_Quantile(
+ 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
+ ),
+ 0.45
+);
--- /dev/null
+0|-10
+0.1|0
+0.2|0
+0.3|0
+0.4|0
+0.5|0
+0.6|0
+0.7|0
+0.8|0
+0.9|0
+1|3.14159
+0|-10
+0.1|-8.685841
+0.2|-7.371682
+0.3|-6.057523
+0.4|-4.743364
+0.5|-3.429205
+0.6|-2.115046
+0.7|-0.800887
+0.8|0.513272
+0.9|1.827431
+1|3.14159
+0|-10
+0.25|0
+0.5|0
+0.75|0
+1|3.14159
+0|-10
+0.25|-6.714602
+0.5|-3.429205
+0.75|-0.143807
+1|3.14159
+0.05|-9.34292
+0.95|2.484511
+0|-10
+0.25|-6.714602
+0.5|-3.429205
+0.75|-0.143807
+1|3.14159
+0.05|0
+0.95|2.484511
+0.7|0
+0.45|-4.086285
projects / postgis / commitdiff