vq/metrics.c

   1 /********************************************************************
   2  *                                                                  *
   3  * THIS FILE IS PART OF THE Ogg Vorbis SOFTWARE CODEC SOURCE CODE.  *
   4  * USE, DISTRIBUTION AND REPRODUCTION OF THIS SOURCE IS GOVERNED BY *
   5  * THE GNU PUBLIC LICENSE 2, WHICH IS INCLUDED WITH THIS SOURCE.    *
   6  * PLEASE READ THESE TERMS DISTRIBUTING.                            *
   7  *                                                                  *
   8  * THE OggSQUISH SOURCE CODE IS (C) COPYRIGHT 1994-2000             *
   9  * by Monty <monty@xiph.org> and The XIPHOPHORUS Company            *
  10  * http://www.xiph.org/                                             *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: function calls to collect codebook metrics
  15  last mod: $Id: metrics.c,v 1.6 2000/02/16 16:18:37 xiphmont Exp $
  16
  17  ********************************************************************/
  18
  19
  20 #include <stdlib.h>
  21 #include <unistd.h>
  22 #include <math.h>
  23 #include "bookutil.h"
  24
  25 /* set up metrics */
  26
  27 double meanamplitude_acc=0.;
  28 double meanamplitudesq_acc=0.;
  29 double meanerror_acc=0.;
  30 double meanerrorsq_acc=0.;
  31 double meandev_acc=0.;
  32
  33 double *histogram=NULL;
  34 double *histogram_error=NULL;
  35 double *histogram_errorsq=NULL;
  36 double *histogram_distance=NULL;
  37 double *histogram_hi=NULL;
  38 double *histogram_lo=NULL;
  39
  40 double count=0.;
  41
  42 int books=0;
  43 int dim=-1;
  44 double *work;
  45
  46 static double *_now(codebook *c, int i){
  47   return c->valuelist+i*c->c->dim;
  48 }
  49
  50 int histerrsort(const void *a, const void *b){
  51   double av=histogram_distance[*((long *)a)];
  52   double bv=histogram_distance[*((long *)b)];
  53   if(av<bv)return(-1);
  54   return(1);
  55 }
  56
  57 void process_preprocess(codebook **bs,char *basename){
  58   while(bs[books]){
  59     codebook *b=bs[books];
  60     if(dim==-1){
  61       dim=b->c->dim;
  62       work=malloc(sizeof(double)*dim);
  63     }else{
  64       if(dim!=b->c->dim){
  65         fprintf(stderr,"Each codebook in a cascade must have the same dimensional order\n");
  66         exit(1);
  67       }
  68     }
  69     books++;
  70   }
  71
  72   if(books){
  73     const static_codebook *b=bs[books-1]->c;
  74     histogram=calloc(b->entries,sizeof(double));
  75     histogram_distance=calloc(b->entries,sizeof(double));
  76     histogram_errorsq=calloc(b->entries*dim,sizeof(double));
  77     histogram_error=calloc(b->entries*dim,sizeof(double));
  78     histogram_hi=calloc(b->entries*dim,sizeof(double));
  79     histogram_lo=calloc(b->entries*dim,sizeof(double));
  80   }else{
  81     fprintf(stderr,"Specify at least one codebook\n");
  82     exit(1);
  83   }
  84 }
  85
  86 static double _dist(int el,double *a, double *b){
  87   int i;
  88   double acc=0.;
  89   for(i=0;i<el;i++){
  90     double val=(a[i]-b[i]);
  91     acc+=val*val;
  92   }
  93   return acc;
  94 }
  95
  96 void cell_spacing(codebook **b){
  97   int i,j,k;
  98   for(i=0;i<books;i++){
  99     double min,max,mean=0.,meansq=0.;
 100     codebook *c=b[i];
 101
 102     fprintf(stderr,"\nCell spacing for book %d:\n",i);
 103
 104     /* minimum, maximum, mean, ms cell spacing */
 105     for(j=0;j<c->c->entries;j++){
 106       double localmin=-1.;
 107       for(k=0;k<c->c->entries;k++){
 108         double this=_dist(c->c->dim,_now(c,j),_now(c,k));
 109         if(j!=k &&
 110            (localmin==-1 || this<localmin))
 111           localmin=this;
 112       }
 113
 114       if(j==0 || localmin<min)min=localmin;
 115       if(j==0 || localmin>max)max=localmin;
 116       mean+=sqrt(localmin);
 117       meansq+=localmin;
 118     }
 119
 120     fprintf(stderr,"\tminimum cell spacing (closest side): %g\n",sqrt(min));
 121     fprintf(stderr,"\tmaximum cell spacing (closest side): %g\n",sqrt(max));
 122     fprintf(stderr,"\tmean closest side spacing: %g\n",mean/c->c->entries);
 123     fprintf(stderr,"\tmean sq closest side spacing: %g\n",sqrt(meansq/c->c->entries));
 124   }
 125 }
 126
 127 void process_postprocess(codebook **b,char *basename){
 128   int i,j,k;
 129   char *buffer=alloca(strlen(basename)+80);
 130   codebook *bb=b[books-1];
 131
 132   fprintf(stderr,"Done.  Processed %ld data points for %ld entries:\n",
 133           (long)count,bb->c->entries);
 134   fprintf(stderr,"\tglobal mean amplitude: %g\n",
 135           meanamplitude_acc/(count*dim));
 136   fprintf(stderr,"\tglobal mean squared amplitude: %g\n",
 137           sqrt(meanamplitudesq_acc/(count*dim)));
 138
 139   fprintf(stderr,"\tglobal mean error: %g\n",
 140           meanerror_acc/(count*dim));
 141   fprintf(stderr,"\tglobal mean squared error: %g\n",
 142           sqrt(meanerrorsq_acc/(count*dim)));
 143   fprintf(stderr,"\tglobal mean deviation: %g\n",
 144           meandev_acc/(count*dim));
 145
 146   cell_spacing(b);
 147
 148   {
 149     FILE *out;
 150
 151     sprintf(buffer,"%s-mse.m",basename);
 152     out=fopen(buffer,"w");
 153     if(!out){
 154       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 155       exit(1);
 156     }
 157
 158     for(i=0;i<bb->c->entries;i++){
 159       for(k=0;k<bb->c->dim;k++){
 160         fprintf(out,"%ld, %g, %g\n",
 161                 i*bb->c->dim+k,(bb->valuelist+i*bb->c->dim)[k],
 162                 sqrt((histogram_errorsq+i*bb->c->dim)[k]/histogram[i]));
 163       }
 164     }
 165     fclose(out);
 166
 167     sprintf(buffer,"%s-me.m",basename);
 168     out=fopen(buffer,"w");
 169     if(!out){
 170       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 171       exit(1);
 172     }
 173
 174     for(i=0;i<bb->c->entries;i++){
 175       for(k=0;k<bb->c->dim;k++){
 176         fprintf(out,"%ld, %g, %g\n",
 177                 i*bb->c->dim+k,(bb->valuelist+i*bb->c->dim)[k],
 178                 (histogram_error+i*bb->c->dim)[k]/histogram[i]);
 179       }
 180     }
 181     fclose(out);
 182
 183     sprintf(buffer,"%s-worst.m",basename);
 184     out=fopen(buffer,"w");
 185     if(!out){
 186       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 187       exit(1);
 188     }
 189
 190     for(i=0;i<bb->c->entries;i++){
 191       for(k=0;k<bb->c->dim;k++){
 192         fprintf(out,"%ld, %g, %g, %g\n",
 193                 i*bb->c->dim+k,(bb->valuelist+i*bb->c->dim)[k],
 194                 (bb->valuelist+i*bb->c->dim)[k]+(histogram_lo+i*bb->c->dim)[k],
 195                 (bb->valuelist+i*bb->c->dim)[k]+(histogram_hi+i*bb->c->dim)[k]);
 196       }
 197     }
 198     fclose(out);
 199   }
 200
 201   {
 202     FILE *out;
 203     long *index=alloca(sizeof(long)*bb->c->entries);
 204     sprintf(buffer,"%s-distance.m",basename);
 205     out=fopen(buffer,"w");
 206     if(!out){
 207       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 208       exit(1);
 209     }
 210     for(j=0;j<bb->c->entries;j++){
 211       if(histogram[j])histogram_distance[j]/=histogram[j];
 212       index[j]=j;
 213     }
 214
 215     qsort(index,bb->c->entries,sizeof(long),histerrsort);
 216
 217     for(j=0;j<bb->c->entries;j++)
 218       for(k=0;k<histogram[index[j]];k++)
 219         fprintf(out,"%g,\n",histogram_distance[index[j]]);
 220     fclose(out);
 221
 222   }
 223 }
 224
 225 void process_vector(codebook **bs,double *a){
 226   int bi;
 227   int i;
 228   double amplitude=0.;
 229   double distance=0.;
 230   double base=0.;
 231   int entry;
 232   double *e;
 233   memcpy(work,a,sizeof(double)*dim);
 234
 235   for(bi=0;bi<books;bi++){
 236     codebook *b=bs[bi];
 237     entry=codebook_entry(b,work);
 238     e=b->valuelist+b->c->dim*entry;
 239     for(i=0;i<b->c->dim;i++)work[i]-=e[i];
 240   }
 241
 242   for(i=0;i<dim;i++){
 243     double error=work[i];
 244     if(bs[0]->c->q_sequencep){
 245       amplitude=a[i]-base;
 246       base=a[i];
 247     }else
 248       amplitude=a[i];
 249
 250     meanamplitude_acc+=fabs(amplitude);
 251     meanamplitudesq_acc+=amplitude*amplitude;
 252     meanerror_acc+=fabs(error);
 253     meanerrorsq_acc+=error*error;
 254
 255     if(amplitude)
 256       meandev_acc+=fabs(error/amplitude);
 257     else
 258       meandev_acc+=fabs(error); /* yeah, yeah */
 259
 260     histogram_errorsq[entry*dim+i]+=error*error;
 261     histogram_error[entry*dim+i]+=fabs(error);
 262     if(histogram[entry]==0 || histogram_hi[entry*dim+i]<error)
 263       histogram_hi[entry*dim+i]=error;
 264     if(histogram[entry]==0 || histogram_lo[entry*dim+i]>error)
 265       histogram_lo[entry*dim+i]=error;
 266     distance+=error*error;
 267   }
 268
 269   histogram[entry]++;
 270   histogram_distance[entry]+=sqrt(distance);
 271
 272   if((long)(count++)%100)spinnit("working.... lines: ",count);
 273 }
 274
 275 void process_usage(void){
 276   fprintf(stderr,
 277           "usage: vqmetrics <codebook>.vqh datafile.vqd [datafile.vqd]...\n\n"
 278           "       data can be taken on stdin.  Output goes to output files:\n"
 279           "       basename-me.m:       gnuplot: mean error by entry value\n"
 280           "       basename-mse.m:      gnuplot: mean square error by entry value\n"
 281           "       basename-worst.m:    gnuplot: worst error by entry value\n"
 282           "       basename-distance.m: gnuplot file showing distance probability\n"
 283           "\n");
 284
 285 }