PostGIS  2.5.0dev-r@@SVN_REVISION@@

◆ rt_band_get_histogram()

rt_histogram rt_band_get_histogram ( rt_bandstats  stats,
uint32_t  bin_count,
double *  bin_width,
uint32_t  bin_width_count,
int  right,
double  min,
double  max,
uint32_t rtn_count 
)

Count the distribution of data.

Parameters
stats: a populated stats struct for processing
bin_count: the number of bins to group the data by
bin_width: the width of each bin as an array
bin_width_count: number of values in bin_width
right: evaluate bins by (a,b] rather than default [a,b)
min: user-defined minimum value of the histogram a value less than the minimum value is not counted in any bins if min = max, min and max are not used
max: user-defined maximum value of the histogram a value greater than the max value is not counted in any bins if min = max, min and max are not used
rtn_count: set to the number of bins being returned
Returns
the histogram of the data or NULL

Definition at line 412 of file rt_statistics.c.

References genraster::count, rt_bandstats_t::count, rt_histogram_t::count, FLT_EQ, rt_histogram_t::inc_max, rt_histogram_t::inc_min, rt_bandstats_t::max, rt_histogram_t::max, rt_bandstats_t::min, rt_histogram_t::min, rt_histogram_t::percent, RASTER_DEBUG, RASTER_DEBUGF, rtalloc(), rtdealloc(), rterror(), genraster::value, and rt_bandstats_t::values.

Referenced by RASTER_histogram(), RASTER_histogramCoverage(), and test_band_stats().

417  {
418  rt_histogram bins = NULL;
419  int init_width = 0;
420  uint32_t i;
421  uint32_t j;
422  double tmp;
423  double value;
424  int sum = 0;
425  double qmin;
426  double qmax;
427 
428 #if POSTGIS_DEBUG_LEVEL > 0
429  clock_t start, stop;
430  double elapsed = 0;
431 #endif
432 
433  RASTER_DEBUG(3, "starting");
434 #if POSTGIS_DEBUG_LEVEL > 0
435  start = clock();
436 #endif
437 
438  assert(NULL != stats);
439  assert(NULL != rtn_count);
440 
441  if (stats->count < 1 || NULL == stats->values) {
442  rterror("rt_util_get_histogram: rt_bandstats object has no value");
443  return NULL;
444  }
445 
446  /* bin width must be positive numbers and not zero */
447  if (NULL != bin_width && bin_width_count > 0) {
448  for (i = 0; i < bin_width_count; i++) {
449  if (bin_width[i] < 0 || FLT_EQ(bin_width[i], 0.0)) {
450  rterror("rt_util_get_histogram: bin_width element is less than or equal to zero");
451  return NULL;
452  }
453  }
454  }
455 
456  /* ignore min and max parameters */
457  if (FLT_EQ(max, min)) {
458  qmin = stats->min;
459  qmax = stats->max;
460  }
461  else {
462  qmin = min;
463  qmax = max;
464  if (qmin > qmax) {
465  qmin = max;
466  qmax = min;
467  }
468  }
469 
470  /* # of bins not provided */
471  if (bin_count <= 0) {
472  /*
473  determine # of bins
474  http://en.wikipedia.org/wiki/Histogram
475 
476  all computed bins are assumed to have equal width
477  */
478  /* Square-root choice for stats->count < 30 */
479  if (stats->count < 30)
480  bin_count = ceil(sqrt(stats->count));
481  /* Sturges' formula for stats->count >= 30 */
482  else
483  bin_count = ceil(log2((double) stats->count) + 1.);
484 
485  /* bin_width_count provided and bin_width has value */
486  if (bin_width_count > 0 && NULL != bin_width) {
487  /* user has defined something specific */
488  if (bin_width_count > bin_count)
489  bin_count = bin_width_count;
490  else if (bin_width_count > 1) {
491  tmp = 0;
492  for (i = 0; i < bin_width_count; i++) tmp += bin_width[i];
493  bin_count = ceil((qmax - qmin) / tmp) * bin_width_count;
494  }
495  else
496  bin_count = ceil((qmax - qmin) / bin_width[0]);
497  }
498  /* set bin width count to zero so that one can be calculated */
499  else {
500  bin_width_count = 0;
501  }
502  }
503 
504  /* min and max the same */
505  if (FLT_EQ(qmax, qmin)) bin_count = 1;
506 
507  RASTER_DEBUGF(3, "bin_count = %d", bin_count);
508 
509  /* bin count = 1, all values are in one bin */
510  if (bin_count < 2) {
511  bins = rtalloc(sizeof(struct rt_histogram_t));
512  if (NULL == bins) {
513  rterror("rt_util_get_histogram: Could not allocate memory for histogram");
514  return NULL;
515  }
516 
517  bins->count = stats->count;
518  bins->percent = -1;
519  bins->min = qmin;
520  bins->max = qmax;
521  bins->inc_min = bins->inc_max = 1;
522 
523  *rtn_count = bin_count;
524  return bins;
525  }
526 
527  /* establish bin width */
528  if (bin_width_count == 0) {
529  bin_width_count = 1;
530 
531  /* bin_width unallocated */
532  if (NULL == bin_width) {
533  bin_width = rtalloc(sizeof(double));
534  if (NULL == bin_width) {
535  rterror("rt_util_get_histogram: Could not allocate memory for bin widths");
536  return NULL;
537  }
538  init_width = 1;
539  }
540 
541  bin_width[0] = (qmax - qmin) / bin_count;
542  }
543 
544  /* initialize bins */
545  bins = rtalloc(bin_count * sizeof(struct rt_histogram_t));
546  if (NULL == bins) {
547  rterror("rt_util_get_histogram: Could not allocate memory for histogram");
548  if (init_width) rtdealloc(bin_width);
549  return NULL;
550  }
551  if (!right)
552  tmp = qmin;
553  else
554  tmp = qmax;
555  for (i = 0; i < bin_count;) {
556  for (j = 0; j < bin_width_count; j++) {
557  bins[i].count = 0;
558  bins->percent = -1;
559 
560  if (!right) {
561  bins[i].min = tmp;
562  tmp += bin_width[j];
563  bins[i].max = tmp;
564 
565  bins[i].inc_min = 1;
566  bins[i].inc_max = 0;
567  }
568  else {
569  bins[i].max = tmp;
570  tmp -= bin_width[j];
571  bins[i].min = tmp;
572 
573  bins[i].inc_min = 0;
574  bins[i].inc_max = 1;
575  }
576 
577  i++;
578  }
579  }
580  if (!right) {
581  bins[bin_count - 1].inc_max = 1;
582 
583  /* align last bin to the max value */
584  if (bins[bin_count - 1].max < qmax)
585  bins[bin_count - 1].max = qmax;
586  }
587  else {
588  bins[bin_count - 1].inc_min = 1;
589 
590  /* align first bin to the min value */
591  if (bins[bin_count - 1].min > qmin)
592  bins[bin_count - 1].min = qmin;
593  }
594 
595  /* process the values */
596  for (i = 0; i < stats->count; i++) {
597  value = stats->values[i];
598 
599  /* default, [a, b) */
600  if (!right) {
601  for (j = 0; j < bin_count; j++) {
602  if (
603  (!bins[j].inc_max && value < bins[j].max) || (
604  bins[j].inc_max && (
605  (value < bins[j].max) ||
606  FLT_EQ(value, bins[j].max)
607  )
608  )
609  ) {
610  bins[j].count++;
611  sum++;
612  break;
613  }
614  }
615  }
616  /* (a, b] */
617  else {
618  for (j = 0; j < bin_count; j++) {
619  if (
620  (!bins[j].inc_min && value > bins[j].min) || (
621  bins[j].inc_min && (
622  (value > bins[j].min) ||
623  FLT_EQ(value, bins[j].min)
624  )
625  )
626  ) {
627  bins[j].count++;
628  sum++;
629  break;
630  }
631  }
632  }
633  }
634 
635  for (i = 0; i < bin_count; i++) {
636  bins[i].percent = ((double) bins[i].count) / sum;
637  }
638 
639 #if POSTGIS_DEBUG_LEVEL > 0
640  stop = clock();
641  elapsed = ((double) (stop - start)) / CLOCKS_PER_SEC;
642  RASTER_DEBUGF(3, "elapsed time = %0.4f", elapsed);
643 
644  for (j = 0; j < bin_count; j++) {
645  RASTER_DEBUGF(5, "(min, max, inc_min, inc_max, count, sum, percent) = (%f, %f, %d, %d, %d, %d, %f)",
646  bins[j].min, bins[j].max, bins[j].inc_min, bins[j].inc_max, bins[j].count, sum, bins[j].percent);
647  }
648 #endif
649 
650  if (init_width) rtdealloc(bin_width);
651  *rtn_count = bin_count;
652  RASTER_DEBUG(3, "done");
653  return bins;
654 }
uint32_t count
Definition: librtcore.h:2360
#define FLT_EQ(x, y)
Definition: librtcore.h:2234
void rterror(const char *fmt,...)
Wrappers used for reporting errors and info.
Definition: rt_context.c:199
void * rtalloc(size_t size)
Wrappers used for managing memory.
Definition: rt_context.c:171
unsigned int uint32_t
Definition: uthash.h:78
double percent
Definition: librtcore.h:2375
int count
Definition: genraster.py:56
double * values
Definition: librtcore.h:2368
#define RASTER_DEBUGF(level, msg,...)
Definition: librtcore.h:299
void rtdealloc(void *mem)
Definition: rt_context.c:186
#define RASTER_DEBUG(level, msg)
Definition: librtcore.h:295
int value
Definition: genraster.py:61
uint32_t count
Definition: librtcore.h:2374
Here is the call graph for this function:
Here is the caller graph for this function: