hist.c

speech signal process tools

/************************************************************************/
/*    file: hist.c							*/
/*    Desc: General purpose histogram manipulating routines.		*/
/*    Date: Nov 27, 1990						*/
/*									*/
/************************************************************************/
/*
DOC:filename: hist.c
DOC:include:  ../include/hist.h
DOC:package:  numeric histogram package
DOC:purpose:  create, manipulate and de-create histograms
DOC:
*/
#include <util/utillib.h>

#define HIST_C_VERSION "V1.0"

#ifndef M_PI
#define M_PI    3.14159265358979323846
#endif
/***********************************************************/
/* DOC-P:init_hist()*/
/* DOC-PS:allocate and init a HIST structure*/
/* Based on the number of bins, and the range of values to */
/* be within the histograms values, build and return a     */
/* histogram.                                              */
void init_hist(HIST ***hist, int num_bins, float from, float to)
{
    HIST **th;
    int i;
    double dist;

    /* what's the span of possible values */
    dist = (double) (to - from);

    alloc_2dimarr(th,num_bins,1,HIST);
    /* initialize the count, and set up the ranges */
    for (i=0; i<num_bins; i++){
        th[i]->count=0;
        th[i]->from=from+(dist*((double)i/(double)num_bins));
        th[i]->to=from+(dist*((double)(i+1)/(double)num_bins));
    }
    *hist=th;
}

/***********************************************************/
/* DOC-P: free_hist()       */
/* DOC-PS:Free the memory associated with a histogram      */
void free_hist(HIST ***hist, int num_bins)
{
    HIST **th;

    th = *hist;
    free_2dimarr(th,num_bins,HIST);
    *hist = (HIST **)0;
}

/************************************************************/
/* DOC-P:   diff_hist()       */
/* DOC-PS:   subtract the counts of hist 2 from hist 1. */
/* subtract the counts of hist 2 from hist 1.  This may     */
/* create negative count values                             */
void diff_hist(HIST **h1, HIST **h2, HIST **hd, int num_bins)
{
    int i;

    for (i=0; i<num_bins; i++)
        hd[i]->count = h1[i]->count - h2[i]->count;
}

/************************************************************/
/* DOC-P:  subtract_hist()       */
/* DOC-PS: subtract h1 from h2, all negative values are set to zero */
void subtract_hist(HIST **h1, HIST **h2, HIST **hs, int num_bins)
{
    int i;

    for (i=0; i<num_bins; i++){
        hs[i]->count = h2[i]->count - h1[i]->count;
        if (hs[i]->count < 0)
            hs[i]->count = 0;
    }
}

/************************************************************/
/* DOC-P:   float percentile_hist()       */
/* DOC-PS: return the mid-point of the bin with the Nth percentile */
/* return the mid-point of the histogram bin containing the */
/* percentile specified.                                    */
float percentile_hist(HIST **hist, int num_bins, float percentile)
{
    int i, pct_area, area=0;

    pct_area = (int)((float)hist_area(hist,num_bins) * percentile);
   
    for (i=0; (i<num_bins) && (area+hist[i]->count < pct_area); i++){
        area+=hist[i]->count;
    }
    return(hist[i]->from+((hist[i]->to-hist[i]->from)/2.0));
}

/************************************************************/
/* DOC-P:   percentage_left_hist()       */
/* DOC-PS:  return the % of HIST to the left of the bin */
/* return the percentage of the histogram to the left of    */
/* specified bin value, including that bin                  */
float percentage_left_hist(HIST **hist, int num_bins, float value)
{
    int i, area=0, left_area=0;

    area = (float)hist_area(hist,num_bins);
   
    for (i=0; (i<num_bins) &&((hist[i]->to+hist[i]->from)/2.0) <= value;i++){
        left_area+=hist[i]->count;
    }

    return((float)left_area/(float)area*100.0);
}

/***************************************************************/
/* DOC-P:   float do_least_squares()                           */
/* DOC-PS:  returns the square-distance between two histograms */
/*          Used as a distance metric in the direct search.    */
/***************************************************************/

float do_least_squares(HIST **noise, HIST **normal, int num_bins)
{
    double sqr_sum=0.0, sqr, extend_db=5.0;
    int i,end=0;

    i=0;
    while ((i<num_bins) && (normal[i]->count <= 0))
        i++;
    end=i;
    if ((i-=num_bins)<0)
        i=0; 

    for (; end<num_bins && (normal[end]->count > 0); end++)
        ;
    if (end>=num_bins) end=num_bins-1;
    end += ((float)(float)num_bins/(normal[num_bins-1]->to-normal[0]->from)) * extend_db;
    if (end>=num_bins) end=num_bins-1;
    for (; i<end; i++){
        sqr      = (float)(noise[i]->count - normal[i]->count) *
                   (float)(noise[i]->count - normal[i]->count) ;
        if (noise[i]->count == 0)
            sqr_sum += sqr*sqr;
        else
            sqr_sum += sqr;
    }
    return(sqr_sum);
}


/*************************************************************/
/* DOC-P:   hist_copy()                                      */
/* DOC-PS:  copy one histogram to another.                   */
void hist_copy(HIST **from, HIST **to, int num_bins, int start, int end)
{
    int i;
 
    for (i=start; (i<num_bins) && (i<=end); i++)
         to[i]->count = from[i]->count;
}

/*************************************************************/
/* DOC-P:   erase_hist()                                     */
/* DOC-PS:  reset all bin counts to zero.                    */
void erase_hist(HIST **hist, int num_bins)
{
    int i;

    for (i=0; i<num_bins; i++)
        hist[i]->count = 0;
}

/**********************************************************************/
/* DOC-P:   dump_hist()                                               */
/* DOC-PS:  print the contents of the histogram to the file pointer   */
void dump_hist(HIST **hist, int num_bins, FILE *fp)
{
    int i;

    fprintf(fp,"Dump of a histogram\n");
    for (i=0; i<num_bins; i++)
        fprintf(fp,"%5.3f : %5.3f -> %d\n",
                hist[i]->from,hist[i]->to,hist[i]->count);
}

/*************************************************************/
/* DOC-P:   dump_esps_hist()                                 */
/* DOC-PS:  dump a HIST in ESPS format for aplot             */
/* Print the histogram to the file "fname" in a format       */
/* readable by the entropic 'aplot' program                  */
void dump_esps_hist(HIST **hist, int num_bins, char *fname)
{
    int i,height;
    FILE *fp;

/*  printf("Dump of a histogram in ESPS format to %s\n",fname); */
    if ((fp=fopen(fname,"w")) == NULL){
        fprintf(stderr,"Warning: unable to open %s for writing\n",fname);
        return;
    }

    fprintf(fp,"%d\n",num_bins*3); /*each bin has three points */
    fprintf(fp,"1\n");             /*I'm not sure why */
    /* the X axis - beg, end, increment*/
    fprintf(fp,"%5.1f %5.1f 10.0\n",hist[0]->from,hist[num_bins-1]->to);
     /* the Y axis - beg, end, increment*/
    height= (int)(max_hist(hist,num_bins)+10)/10*10;
    fprintf(fp,"%d %d %d\n",0,height,height/10);

    for (i=0; i<num_bins; i++){
        fprintf(fp,"%5.1f %d\n", hist[i]->from,hist[i]->count);
        fprintf(fp,"%5.1f %d\n", hist[i]->to,hist[i]->count);
        fprintf(fp,"%5.1f 0\n", hist[i]->to);
    }
    fflush(fp);
    fclose(fp);
}

/**************************************************************/
/* DOC-P:   read_esps_hist()                                  */
/* DOC-PS:  read an ESPS format histogram                     */
/* Read an entropic histogram file, then create a histogram to*/
/* hold it's information and return that histogram            */
void read_esps_hist(HIST ***hist, int *num_bins, char *fname)
{
    int i;
    FILE *fp;
    HIST **t_hist;
    char buff[400];
    float hist_beg, hist_end, hist_incr;

    printf("Reading of a histogram in ESPS format from %s\n",fname);
    if ((fp=fopen(fname,"r")) == NULL){
        fprintf(stderr,"Warning: unable to open %s for reading\n",fname);
        return;
    }

    fscanf(fp,"%d\n",&i);
    *num_bins = i/3;
    printf("number of bins to read %d\n",*num_bins);

    fgets(buff,400,fp);  /* skip a line */
    
    fscanf(fp,"%f %f %f\n",&hist_beg, &hist_end, &hist_incr);
    printf("hist beg %f,  end %f  incr %f\n",hist_beg,hist_end,hist_incr);
    init_hist(&t_hist,*num_bins,hist_beg,hist_end);

    fgets(buff,400,fp); /* skip a line */

    for (i=0; i<*num_bins; i++){
        fscanf(fp,"%*f %d\n",&(t_hist[i]->count));
        fgets(buff,400,fp);  /* skip a line */
        fgets(buff,400,fp);  /* skip a line */
    }

    fclose(fp);
    *hist = t_hist;
}

/*************************************************************/
/* DOC-P:   dump_gnuplot_hist()                              */
/* DOC-PS:  dump a HIST in gnuplot format                    */
/* Write to the file 'fname' the histogram in a form usable  */
/* by gnuplot                                                */
void dump_gnuplot_hist(HIST **hist, int num_bins, char *fname)
{
    int i;
    FILE *fp;

/*  printf("Dump of a histogram in GNUPLOT format to %s\n",fname); */
    if ((fp=fopen(fname,"w")) == NULL){
        fprintf(stderr,"Warning: unable to open %s for writing\n",fname);
        return;
    }

    for (i=0; i<num_bins; i++){
        fprintf(fp,"%5.3f %d\n", (hist[i]->to-hist[i]->from)/2.0+hist[i]->from,
                hist[i]->count);
    }

    fflush(fp);
    fclose(fp);
}

/*************************************************************/
/* DOC-P:   dump_gnuplot_2hist()                             */
/* DOC-PS:  dump a 2 HIST's in gnuplot format                */
/* Write to the file 'fname' the histogram in a form usable  */
/* by gnuplot                                                */
void dump_gnuplot_2hist(HIST **hist1, HIST **hist2, int num_bins, char *fname)
{
    int i;
    FILE *fp;

/*  printf("Dump of a histogram in GNUPLOT format to %s\n",fname); */
    if ((fp=fopen(fname,"w")) == NULL){