xx.c

speech signal process tools

/*
  stats - compute statistics on SPS files				     
  Joseph T. Buck							    
 									   
  stats considers an SPS file to be a series of vectors. It computes mean, 
  standard deviation, maximum, and minimum values on an element-by-element 
  basis. The number of vectors used may be controlled with the -o switch.  
  The -n switch suppresses the header.					   
  									   
  From SDS version 1.3 of 10/28/85					   

  Joe Buck started converting it to SPS; Alan Parker finsished it.

  This material contains proprietary software of Entropic Processing, Inc.   
  Any reproduction, distribution, or publication without the the prior	   
  written permission of Entropic Processing, Inc. is strictly prohibited.
  Any public distribution of copies of this work authorized in writing by
  Entropic Processing, Inc. must bear the notice			
 								
               "Copyright 1986 Entropic Processing, Inc."

*/
#ifdef SCCS
	static char *sccs_id = "%W% %G%";
#endif
#include <stdio.h>
#include <sps/sps.h>
#define MAXL 100
#define MAX_FLOAT 1.0e38
#define MIN_FLOAT -MAX_FLOAT
#define SYNTAX USAGE ("stats [-o rec_range] [-e el_range] [-n] file")
main (argc, argv)
int   argc;
char  **argv;
{
    extern  optind;
    extern char *optarg;
    int     c, s_rec, e_rec, n_rec, l_rec, i_rec, i, hflag = 1;
    long    tag;
    char   *prange = NULL, *erange = NULL, *filename = NULL;
    double  data [MAXL], mins[MAXL], maxs[MAXL], avgs[MAXL], devs[MAXL];
    int     minl [MAXL], maxl[MAXL];
    int     s_ele, e_ele;
    double  scale, bias, sqrt ();
    FILE   *istrm;
    struct header *h;

    while ((c = getopt (argc, argv, "no:e:")) != EOF) {
	switch (c) {
	    case 'n':
		hflag = 0;
		break;
	    case 'o': 
		prange = optarg;
		break;
	    case 'e':
	        erange = optarg;
		break;
	    default:
		SYNTAX;
	}
    }
    if (argc == optind + 1)
	filename = argv[optind];
    read_params(NULL,SC_CHECK_FILE,filename);
    if (filename == NULL) {
	if(symtype("filename") == ST_UNDEF) {
	  fprintf (stderr, "stats: no input file.\n");
	  exit(1);
    	}
	filename = getsym_s("filename");
    }

    TRYOPEN ("stats", filename, "r", istrm);
    if (!(h = read_header (istrm)))
	NOTSPS ("stats", filename);
    l_rec = get_rec_len(h);
    if (l_rec > MAXL) {
	fprintf (stderr, "stats: record length (%d) too long\n", l_rec);
	exit (1);
    }
    s_rec = 1;
    e_rec = h -> common.ndrec;
    if(symtype("start") != ST_UNDEF) 
	s_rec = getsym_i("start");
    if(symtype("nan") != ST_UNDEF) 
	e_rec = s_rec + getsym_i("nan") - 1;
    range_switch (prange, &s_rec, &e_rec, 1);
    (void) putsym_s("filename",filename);
    (void) putsym_s("prog","stats");
    (void) putsym_i("start",s_rec);
    (void) putsym_i("nan",e_rec-s_rec+1);
    if (hflag) printf("stats: filename: %s, start: %d, end: %d\n",
	filename, s_rec, e_rec);

    if (e_rec > h -> common.ndrec) {
	fprintf (stderr, "stats: only %d records in file\n", h->common.ndrec);
	exit (1);
    }

    s_ele = 1;
    e_ele = l_rec;
    range_switch (erange, &s_ele, &e_ele, 1);
    if (e_ele > l_rec) {
	fprintf (stderr, "stats: only %d elements per record\n", l_rec);
	exit (1);
    }
/* Initialize statistics */
    for (i = s_ele - 1; i < e_ele; i++) {
	avgs[i] = devs[i] = 0.0;
	mins[i] = MAX_FLOAT;
	maxs[i] = MIN_FLOAT;
    }
    skiprec (istrm, s_rec - 1, size_rec (h));
    n_rec = e_rec - s_rec + 1;
    for (i_rec = s_rec; i_rec <= e_rec; i_rec++) {
	if(get_gen_recd(data,&tag,h,istrm) == EOF) {
		fprintf(stderr,"stats: only read %d records\n",i_rec-s_rec);
		break;
	}
	for (i = s_ele - 1; i < e_ele; i++) {
	    if (data[i] < mins[i]) {
		mins[i] = data[i];
		minl[i] = i_rec;
	    }
	    if (data[i] > maxs[i]) {
		maxs[i] = data[i];
		maxl[i] = i_rec;
	    }
	    avgs[i] += data[i];
	    devs[i] += data[i] * data[i];
	}
    }
    scale = 1.0 / (double) n_rec;
/* bias is the factor that converts the biased estimate of the variance to
   the unbiased estimate: n_rec/(n_rec-1).
*/
    bias = 0.0;
    if (n_rec > 1)
	bias = (double) n_rec / (double) (n_rec - 1);
    for (i = s_ele - 1; i < e_ele; i++) {
	avgs[i] *= scale;
/*
 * for some cases this expression yields a negative number very close
 * zero that should be treated as zero.
 */
	if ((bias * (devs[i] * scale - avgs[i] * avgs[i])) <= 0)
	  devs[i] = 0;
	 else
	  devs[i] = sqrt (bias * (devs[i] * scale - avgs[i] * avgs[i]));
    }
    if (hflag)
	printf ("Element   Minimum    at     Maximum    at      Mean      Std. Dev\n");
    for (i = s_ele - 1; i < e_ele; i++) {
	printf ("%5d%12g%6d%12g%6d%12g%12g\n", i + 1, mins[i], minl[i],
		maxs[i], maxl[i], avgs[i], devs[i]);
    }
    return 0;
}