---

Histogram.cpp

Go to the documentation of this file.

/**************************************************************************
 *
 * Source code for Real Time Image Library (libRTImage)
 *
 * Leeds Vision Group give permission for this code to be copied, modified
 * and distributed within the University of Leeds subject to the following
 * conditions:-
 *
 * - The code is not to be used for commercial gain.
 * - The code and use thereof will be attributed to the author where
 *   appropriate (inluding demonstrations which rely on it's use).
 * - All modified, distributions of the source files will retain this header.
 *
 ****************************************************************************/

#include <cstdio>
#include <cmath>
#include <iostream>
#include <fstream>

using namespace std; 

#include "Histogram.h"

Histogram::Histogram(float        minimum_data_value,
                     float        maximum_data_value,
                     unsigned int bands_per_dimension,
                     unsigned int dimensionality)
 : min_val       (minimum_data_value),
   max_val       (maximum_data_value),
   bands_per_dim (bands_per_dimension),
   dim           (dimensionality),
   max_prob      (0),
   max_freq      (0)
{
    band_size = (maximum_data_value - minimum_data_value) / bands_per_dim;

    no_bins = (unsigned int)pow((double)bands_per_dimension, (double)dim);
       
    frequency_data = new unsigned int[no_bins];   
    probs          = new float[no_bins];   

    reset_histogram() ;
}


Histogram::Histogram(char const* histfile)
 : frequency_data(0), probs(0)
{
    load_histogram(histfile);
}


void Histogram::reset_histogram()
{  
    unsigned int cnt;

    no_data_items = 0;
    for (cnt = 0; cnt < no_bins; cnt++) {
        frequency_data[cnt] = 0;
    }
    max_prob = 0;
    max_freq = 0;
}


void Histogram::add_new_data_item(int *data) 
{
    unsigned int   increment;
    unsigned int   d_cnt;
    unsigned int   b_cnt;
    register float band_min;
    register float band_max;
    unsigned int   bin_no = 0; 
    register float d_val;

    increment = 1; 

    for (d_cnt = 0; d_cnt < dim; d_cnt++) {
    // Loop through dimensions
         band_min = min_val;
         band_max = min_val + band_size;
         d_val = (float)(data[d_cnt]);
         
         for (b_cnt = 0; b_cnt < bands_per_dim; b_cnt++) {
         // Loop through bands for this dimension
             if (d_val >= band_min && d_val <= band_max) break;

             band_min += band_size;
             band_max += band_size;
         }
         // Update bin no. w.r.t this dimension
         bin_no += b_cnt * increment;

         // Decrease size of increment
         increment *= bands_per_dim;
    }

    // Update frequency data
    frequency_data[bin_no]++;
    no_data_items++;
    
    if (frequency_data[bin_no] > max_freq) max_freq = frequency_data[bin_no];
} 




void Histogram::add_new_data_item(float *data) 
{
    unsigned int increment;
    unsigned int d_cnt;
    unsigned int b_cnt;
    float        band_min;
    float        band_max;
    unsigned int bin_no = 0; 

    increment = 1; 

    for (d_cnt = 0; d_cnt < dim; d_cnt++) {
    // Loop through dimensions
         band_min = min_val;
         band_max = min_val + band_size;
         for (b_cnt = 0; b_cnt < bands_per_dim; b_cnt++){
         // Loop through bands for this dimension
             if (data[d_cnt] >= band_min && data[d_cnt] <= band_max) break;

             band_min += band_size;
             band_max += band_size;
         }
         // Update bin no. w.r.t this dimension
         bin_no += b_cnt * increment;

         // Decrease size of increment
         increment *= bands_per_dim;
    }

    // Update frequency data
    frequency_data[bin_no]++;
    no_data_items++;
    
    if (frequency_data[bin_no] > max_freq) max_freq = frequency_data[bin_no];
} 

void Histogram::calculate_probabilities()
{
    unsigned int *f_ptr = frequency_data; 
    float        *p_ptr = probs;
    unsigned int cnt ;

    if(no_data_items > 0){
        max_prob = 0 ;
        for (cnt = 0; cnt < no_bins; cnt++, f_ptr++, p_ptr++){
        // Loop through all bins
            *p_ptr = (float)(*f_ptr) / (float)no_data_items;
            if (*p_ptr > max_prob) max_prob = *p_ptr;
        }
    }
    else{
        for (cnt = 0; cnt < no_bins; cnt++){
        // Loop through all bins
            *p_ptr = 0;
        }
        max_prob = 0;
    }

}

float Histogram::calculate_proportion_threshold(float proportion)
{
    float        *p_ptr ;
    unsigned int cnt ;
    float        min_p = 0 ;
    float        new_min_p ;
    float        p_above_thresh = 1.0 ;
    unsigned int no_with_p ;

    while(p_above_thresh > proportion){
        new_min_p = 1.0 ;
        no_with_p = 0 ;
        for (cnt = 0, p_ptr = probs ; cnt < no_bins; cnt++, p_ptr++){
        // Loop through all bins
            if(*p_ptr > min_p && 
               *p_ptr < new_min_p){
            // Lowest this round
                new_min_p = *p_ptr ;
                no_with_p = 1 ;
            }
            else if(*p_ptr > min_p && 
                    *p_ptr == new_min_p){
            // Joint lowest
                no_with_p++ ;
            }
        }
        p_above_thresh -= new_min_p * no_with_p ;
        if(p_above_thresh == 1.0){
        // Error, everything is 0
            break ;
        }
        min_p = new_min_p ;
    }

    return min_p ;
}

void Histogram::calculate_histogram()
{
    unsigned int *f_ptr = frequency_data;
    float        *p_ptr = probs;

    max_prob = 0 ;
    for (unsigned int cnt = 0; cnt < no_bins; cnt++, f_ptr++, p_ptr++){
    // Loop through all bins
        *p_ptr = (float)(*f_ptr) / (float)max_freq;
        if (*p_ptr > max_prob) max_prob = *p_ptr;
    }
}


float Histogram::get_histogram_value(unsigned int *bin)
{
    unsigned int bin_no = 0;
    unsigned int increment;

    increment = 1;

    for (unsigned int d_cnt = 0; d_cnt < dim; d_cnt++) {
    // Loop through dimensions

         bin_no += bin[d_cnt] * increment;
         increment *= bands_per_dim;
    }

    return probs[bin_no];
}


float Histogram::get_mapped_value(unsigned int *data)
{
#if 0
    for (unsigned int d_cnt = 0; d_cnt < dim; d_cnt++)
        data[d_cnt] /= (int) band_size;
    
    return get_histogram_value(data);
#else
    unsigned int   increment;
    unsigned int   d_cnt;
    unsigned int   b_cnt;
    register float band_min;
    register float band_max;
    unsigned int   bin_no = 0; 
    register float d_val;

    increment = 1; 

    for (d_cnt = 0; d_cnt < dim; d_cnt++) {
    // Loop through dimensions
         band_min = min_val;
         band_max = min_val + band_size;
         d_val = (float)(data[d_cnt]);
         
         for (b_cnt = 0; b_cnt < bands_per_dim; b_cnt++) {
         // Loop through bands for this dimension
             if (d_val >= band_min && d_val <= band_max) break;

             band_min += band_size;
             band_max += band_size;
         }
         // Update bin no. w.r.t this dimension
         bin_no += b_cnt * increment;

         // Decrease size of increment
         increment *= bands_per_dim;
    }
    
    return probs[bin_no];
#endif
}


float Histogram::calculate_difference(Histogram *hist)
{
    unsigned int b_cnt;
    float        diff;
    float        ret_val = 0;
    float        *t_ptr;
    float        *h_ptr;

    t_ptr = probs;
    h_ptr = hist->probs;
    for (b_cnt = 0; b_cnt < no_bins; b_cnt++) {
        diff = *t_ptr - *h_ptr;
        if (diff < 0) diff = -diff;
        ret_val += diff;

        t_ptr++;
        h_ptr++;
    }

    return ret_val;
}


void Histogram::save_histogram(char const* histfile)
{
    FILE* fd = fopen(histfile, "w");

    fprintf(fd, "<?xml version=\"1.0\" ?>\n\n");
    
    fprintf(fd, "<histogram min-value=\"%f\" max-value=\"%f\" " 
                "bands=\"%u\" dimensionality=\"%u\">\n",
                min_val, max_val, bands_per_dim, dim);

    for (unsigned int b_cnt = 0; b_cnt < no_bins; b_cnt++) {
        float p = probs[b_cnt];
        if (p) fprintf(fd, "  <prob index=\"%u\">%f</prob>\n", b_cnt, p);
    }
    
    fprintf(fd, "</histogram>\n");
    fclose(fd);
}


void Histogram::load_histogram(char const* histfile)
{
    delete [] frequency_data;
    delete [] probs;
    
    FILE* fd = fopen(histfile, "r");
    if (fd == 0) throw "Could not open hst file.";

    int n;  // return values for fscanf, set to number of items read

    /* xml version string */
    
    float xmlversion = 0;
   
    n = fscanf(fd, "<?xml version=\"%f\" ?>\n\n", &xmlversion);
    
    if (xmlversion == 0) throw "Incompatible hst file format";


    /* histogram attributes */
    
    n = fscanf(fd, "<histogram min-value=\"%f\" max-value=\"%f\" "
                              "bands=\"%u\" dimensionality=\"%u\">\n",
               &min_val, &max_val, &bands_per_dim, &dim);
    
    if (n != 4) throw "error in hst file, unable to import data.";

    band_size = (max_val - min_val) / bands_per_dim;
    no_bins   = (unsigned int) pow((double)bands_per_dim, (double)dim);
   
    frequency_data = new unsigned int[no_bins];   
    probs          = new float[no_bins];   
   
    for (unsigned int cnt = 0; cnt < no_bins; cnt++) {
        frequency_data[cnt] = 0;
        probs[cnt]          = 0.0;
    }

    
    /* probabilities */
    
    unsigned int bin, i = 0;
    float p = 0;    // probability
    
    n = fscanf(fd, "  <prob index=\"%u\">%f</prob>\n", &bin, &p);
    max_prob = p;
    while (n == 2) {
        probs[bin] = p;
        n = fscanf(fd, "  <prob index=\"%u\">%f</prob>\n", &bin, &p);
        if (p > max_prob) max_prob = p;
        i++;
    };
    
    /* end of file */
    
    char last[2];
    n = fscanf(fd, "/histogram%[>]", last);
    if (n != 1)
        throw "Reading hst file.";
    
    fclose(fd);
}

void Histogram::get_best_bin(unsigned int *bin)
{
    unsigned int best_bin = 0 ; 
    unsigned int b_cnt ; 
    unsigned int increment = 1 ; 
    int          d_cnt ; 
    float        highest_p = -1 ;

    for (b_cnt = 0; b_cnt < no_bins; b_cnt++) {
    // Loop through bins
        if(probs[b_cnt] > highest_p){
            highest_p = probs[b_cnt] ;
            best_bin = b_cnt ;
        }
    }

    for(d_cnt = dim-1 ; d_cnt> 0 ; d_cnt--){
        increment *= bands_per_dim; 
    }

    for(d_cnt = dim-1 ; d_cnt>= 0 ; d_cnt--){
        bin[d_cnt] = 0 ;
        while(best_bin >= increment){
            bin[d_cnt] ++ ;
            best_bin -= increment ;
        }
        increment /= bands_per_dim ;
    }
}

---