cannam@86: /********************************************************************
cannam@86:  *                                                                  *
cannam@86:  * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE.   *
cannam@86:  * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
cannam@86:  * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
cannam@86:  * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
cannam@86:  *                                                                  *
cannam@86:  * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001             *
cannam@86:  * by the Xiph.Org Foundation http://www.xiph.org/                  *
cannam@86:  *                                                                  *
cannam@86:  ********************************************************************
cannam@86: 
cannam@86:  function: function calls to collect codebook metrics
cannam@86:  last mod: $Id: metrics.c 16037 2009-05-26 21:10:58Z xiphmont $
cannam@86: 
cannam@86:  ********************************************************************/
cannam@86: 
cannam@86: 
cannam@86: #include <stdlib.h>
cannam@86: #include <unistd.h>
cannam@86: #include <math.h>
cannam@86: #include "bookutil.h"
cannam@86: 
cannam@86: /* collect the following metrics:
cannam@86: 
cannam@86:    mean and mean squared amplitude
cannam@86:    mean and mean squared error 
cannam@86:    mean and mean squared error (per sample) by entry
cannam@86:    worst case fit by entry
cannam@86:    entry cell size
cannam@86:    hits by entry
cannam@86:    total bits
cannam@86:    total samples
cannam@86:    (average bits per sample)*/
cannam@86:    
cannam@86: 
cannam@86: /* set up metrics */
cannam@86: 
cannam@86: float meanamplitude_acc=0.f;
cannam@86: float meanamplitudesq_acc=0.f;
cannam@86: float meanerror_acc=0.f;
cannam@86: float meanerrorsq_acc=0.f;
cannam@86: 
cannam@86: float **histogram=NULL;
cannam@86: float **histogram_error=NULL;
cannam@86: float **histogram_errorsq=NULL;
cannam@86: float **histogram_hi=NULL;
cannam@86: float **histogram_lo=NULL;
cannam@86: float bits=0.f;
cannam@86: float count=0.f;
cannam@86: 
cannam@86: static float *_now(codebook *c, int i){
cannam@86:   return c->valuelist+i*c->c->dim;
cannam@86: }
cannam@86: 
cannam@86: int books=0;
cannam@86: 
cannam@86: void process_preprocess(codebook **bs,char *basename){
cannam@86:   int i;
cannam@86:   while(bs[books])books++;
cannam@86:   
cannam@86:   if(books){
cannam@86:     histogram=_ogg_calloc(books,sizeof(float *));
cannam@86:     histogram_error=_ogg_calloc(books,sizeof(float *));
cannam@86:     histogram_errorsq=_ogg_calloc(books,sizeof(float *));
cannam@86:     histogram_hi=_ogg_calloc(books,sizeof(float *));
cannam@86:     histogram_lo=_ogg_calloc(books,sizeof(float *));
cannam@86:   }else{
cannam@86:     fprintf(stderr,"Specify at least one codebook\n");
cannam@86:     exit(1);
cannam@86:   }
cannam@86: 
cannam@86:   for(i=0;i<books;i++){
cannam@86:     codebook *b=bs[i];
cannam@86:     histogram[i]=_ogg_calloc(b->entries,sizeof(float));
cannam@86:     histogram_error[i]=_ogg_calloc(b->entries*b->dim,sizeof(float));
cannam@86:     histogram_errorsq[i]=_ogg_calloc(b->entries*b->dim,sizeof(float));
cannam@86:     histogram_hi[i]=_ogg_calloc(b->entries*b->dim,sizeof(float));
cannam@86:     histogram_lo[i]=_ogg_calloc(b->entries*b->dim,sizeof(float));
cannam@86:   }
cannam@86: }
cannam@86: 
cannam@86: static float _dist(int el,float *a, float *b){
cannam@86:   int i;
cannam@86:   float acc=0.f;
cannam@86:   for(i=0;i<el;i++){
cannam@86:     float val=(a[i]-b[i]);
cannam@86:     acc+=val*val;
cannam@86:   }
cannam@86:   return acc;
cannam@86: }
cannam@86: 
cannam@86: void cell_spacing(codebook *c){
cannam@86:   int j,k;
cannam@86:   float min=-1.f,max=-1.f,mean=0.f,meansq=0.f;
cannam@86:   long total=0;
cannam@86: 
cannam@86:   /* minimum, maximum, mean, ms cell spacing */
cannam@86:   for(j=0;j<c->c->entries;j++){
cannam@86:     if(c->c->lengthlist[j]>0){
cannam@86:       float localmin=-1.;
cannam@86:       for(k=0;k<c->c->entries;k++){
cannam@86:         if(c->c->lengthlist[k]>0){
cannam@86:           float this=_dist(c->c->dim,_now(c,j),_now(c,k));
cannam@86:           if(j!=k &&
cannam@86:              (localmin==-1 || this<localmin))
cannam@86:             localmin=this;
cannam@86:         }
cannam@86:       }
cannam@86:       
cannam@86:       if(min==-1 || localmin<min)min=localmin;
cannam@86:       if(max==-1 || localmin>max)max=localmin;
cannam@86:       mean+=sqrt(localmin);
cannam@86:       meansq+=localmin;
cannam@86:       total++;
cannam@86:     }
cannam@86:   }
cannam@86:   
cannam@86:   fprintf(stderr,"\tminimum cell spacing (closest side): %g\n",sqrt(min));
cannam@86:   fprintf(stderr,"\tmaximum cell spacing (closest side): %g\n",sqrt(max));
cannam@86:   fprintf(stderr,"\tmean closest side spacing: %g\n",mean/total);
cannam@86:   fprintf(stderr,"\tmean sq closest side spacing: %g\n",sqrt(meansq/total));
cannam@86: }
cannam@86: 
cannam@86: void process_postprocess(codebook **bs,char *basename){
cannam@86:   int i,k,book;
cannam@86:   char *buffer=alloca(strlen(basename)+80);
cannam@86: 
cannam@86:   fprintf(stderr,"Done.  Processed %ld data points:\n\n",
cannam@86:           (long)count);
cannam@86: 
cannam@86:   fprintf(stderr,"Global statistics:******************\n\n");
cannam@86: 
cannam@86:   fprintf(stderr,"\ttotal samples: %ld\n",(long)count);
cannam@86:   fprintf(stderr,"\ttotal bits required to code: %ld\n",(long)bits);
cannam@86:   fprintf(stderr,"\taverage bits per sample: %g\n\n",bits/count);
cannam@86: 
cannam@86:   fprintf(stderr,"\tmean sample amplitude: %g\n",
cannam@86:           meanamplitude_acc/count);
cannam@86:   fprintf(stderr,"\tmean squared sample amplitude: %g\n\n",
cannam@86:           sqrt(meanamplitudesq_acc/count));
cannam@86: 
cannam@86:   fprintf(stderr,"\tmean code error: %g\n",
cannam@86:           meanerror_acc/count);
cannam@86:   fprintf(stderr,"\tmean squared code error: %g\n\n",
cannam@86:           sqrt(meanerrorsq_acc/count));
cannam@86: 
cannam@86:   for(book=0;book<books;book++){
cannam@86:     FILE *out;
cannam@86:     codebook *b=bs[book];
cannam@86:     int n=b->c->entries;
cannam@86:     int dim=b->c->dim;
cannam@86: 
cannam@86:     fprintf(stderr,"Book %d statistics:------------------\n",book);
cannam@86: 
cannam@86:     cell_spacing(b);
cannam@86: 
cannam@86:     sprintf(buffer,"%s-%d-mse.m",basename,book);
cannam@86:     out=fopen(buffer,"w");
cannam@86:     if(!out){
cannam@86:       fprintf(stderr,"Could not open file %s for writing\n",buffer);
cannam@86:       exit(1);
cannam@86:     }
cannam@86:     
cannam@86:     for(i=0;i<n;i++){
cannam@86:       for(k=0;k<dim;k++){
cannam@86:         fprintf(out,"%d, %g, %g\n",
cannam@86:                 i*dim+k,(b->valuelist+i*dim)[k],
cannam@86:                 sqrt((histogram_errorsq[book]+i*dim)[k]/histogram[book][i]));
cannam@86:       }
cannam@86:     }
cannam@86:     fclose(out);
cannam@86:       
cannam@86:     sprintf(buffer,"%s-%d-me.m",basename,book);
cannam@86:     out=fopen(buffer,"w");
cannam@86:     if(!out){
cannam@86:       fprintf(stderr,"Could not open file %s for writing\n",buffer);
cannam@86:       exit(1);
cannam@86:     }
cannam@86:     
cannam@86:     for(i=0;i<n;i++){
cannam@86:       for(k=0;k<dim;k++){
cannam@86:         fprintf(out,"%d, %g, %g\n",
cannam@86:                 i*dim+k,(b->valuelist+i*dim)[k],
cannam@86:                 (histogram_error[book]+i*dim)[k]/histogram[book][i]);
cannam@86:       }
cannam@86:     }
cannam@86:     fclose(out);
cannam@86: 
cannam@86:     sprintf(buffer,"%s-%d-worst.m",basename,book);
cannam@86:     out=fopen(buffer,"w");
cannam@86:     if(!out){
cannam@86:       fprintf(stderr,"Could not open file %s for writing\n",buffer);
cannam@86:       exit(1);
cannam@86:     }
cannam@86:     
cannam@86:     for(i=0;i<n;i++){
cannam@86:       for(k=0;k<dim;k++){
cannam@86:         fprintf(out,"%d, %g, %g, %g\n",
cannam@86:                 i*dim+k,(b->valuelist+i*dim)[k],
cannam@86:                 (b->valuelist+i*dim)[k]+(histogram_lo[book]+i*dim)[k],
cannam@86:                 (b->valuelist+i*dim)[k]+(histogram_hi[book]+i*dim)[k]);
cannam@86:       }
cannam@86:     }
cannam@86:     fclose(out);
cannam@86:   }
cannam@86: }
cannam@86: 
cannam@86: float process_one(codebook *b,int book,float *a,int dim,int step,int addmul,
cannam@86:                    float base){
cannam@86:   int j,entry;
cannam@86:   float amplitude=0.f;
cannam@86: 
cannam@86:   if(book==0){
cannam@86:     float last=base;
cannam@86:     for(j=0;j<dim;j++){
cannam@86:       amplitude=a[j*step]-(b->c->q_sequencep?last:0);
cannam@86:       meanamplitude_acc+=fabs(amplitude);
cannam@86:       meanamplitudesq_acc+=amplitude*amplitude;
cannam@86:       count++;
cannam@86:       last=a[j*step];
cannam@86:     }
cannam@86:   }
cannam@86: 
cannam@86:   if(b->c->q_sequencep){
cannam@86:     float temp;
cannam@86:     for(j=0;j<dim;j++){
cannam@86:       temp=a[j*step];
cannam@86:       a[j*step]-=base;
cannam@86:     }
cannam@86:     base=temp;
cannam@86:   }
cannam@86: 
cannam@86:   entry=vorbis_book_besterror(b,a,step,addmul);
cannam@86: 
cannam@86:   if(entry==-1){
cannam@86:     fprintf(stderr,"Internal error: _best returned -1.\n");
cannam@86:     exit(1);
cannam@86:   }
cannam@86:   
cannam@86:   histogram[book][entry]++;  
cannam@86:   bits+=vorbis_book_codelen(b,entry);
cannam@86:           
cannam@86:   for(j=0;j<dim;j++){
cannam@86:     float error=a[j*step];
cannam@86: 
cannam@86:     if(book==books-1){
cannam@86:       meanerror_acc+=fabs(error);
cannam@86:       meanerrorsq_acc+=error*error;
cannam@86:     }
cannam@86:     histogram_errorsq[book][entry*dim+j]+=error*error;
cannam@86:     histogram_error[book][entry*dim+j]+=fabs(error);
cannam@86:     if(histogram[book][entry]==0 || histogram_hi[book][entry*dim+j]<error)
cannam@86:       histogram_hi[book][entry*dim+j]=error;
cannam@86:     if(histogram[book][entry]==0 || histogram_lo[book][entry*dim+j]>error)
cannam@86:       histogram_lo[book][entry*dim+j]=error;
cannam@86:   }
cannam@86:   return base;
cannam@86: }
cannam@86: 
cannam@86: 
cannam@86: void process_vector(codebook **bs,int *addmul,int inter,float *a,int n){
cannam@86:   int bi;
cannam@86:   int i;
cannam@86: 
cannam@86:   for(bi=0;bi<books;bi++){
cannam@86:     codebook *b=bs[bi];
cannam@86:     int dim=b->dim;
cannam@86:     float base=0.f;
cannam@86: 
cannam@86:     if(inter){
cannam@86:       for(i=0;i<n/dim;i++)
cannam@86:         base=process_one(b,bi,a+i,dim,n/dim,addmul[bi],base);
cannam@86:     }else{
cannam@86:       for(i=0;i<=n-dim;i+=dim)
cannam@86:         base=process_one(b,bi,a+i,dim,1,addmul[bi],base);
cannam@86:     }
cannam@86:   }
cannam@86:   
cannam@86:   if((long)(count)%100)spinnit("working.... samples: ",count);
cannam@86: }
cannam@86: 
cannam@86: void process_usage(void){
cannam@86:   fprintf(stderr,
cannam@86:           "usage: vqmetrics [-i] +|*<codebook>.vqh [ +|*<codebook.vqh> ]... \n"
cannam@86:           "                 datafile.vqd [datafile.vqd]...\n\n"
cannam@86:           "       data can be taken on stdin.  -i indicates interleaved coding.\n"
cannam@86:           "       Output goes to output files:\n"
cannam@86:           "       basename-me.m:       gnuplot: mean error by entry value\n"
cannam@86:           "       basename-mse.m:      gnuplot: mean square error by entry value\n"
cannam@86:           "       basename-worst.m:    gnuplot: worst error by entry value\n"
cannam@86:           "       basename-distance.m: gnuplot file showing distance probability\n"
cannam@86:           "\n");
cannam@86: 
cannam@86: }