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.3 2000/01/10 10:42:04 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 int histerrsort(const void *a, const void *b){
  47   double av=histogram_distance[*((long *)a)];
  48   double bv=histogram_distance[*((long *)b)];
  49   if(av<bv)return(-1);
  50   return(1);
  51 }
  52
  53 void process_preprocess(codebook **bs,char *basename){
  54   while(bs[books]){
  55     codebook *b=bs[books];
  56     if(dim==-1){
  57       dim=b->dim;
  58       work=malloc(sizeof(double)*dim);
  59     }else{
  60       if(dim!=b->dim){
  61         fprintf(stderr,"Each codebook in a cascade must have the same dimensional order\n");
  62         exit(1);
  63       }
  64     }
  65     books++;
  66   }
  67
  68   if(books){
  69     codebook *b=bs[books-1];
  70     histogram=calloc(b->entries,sizeof(double));
  71     histogram_distance=calloc(b->entries,sizeof(double));
  72     histogram_errorsq=calloc(b->entries*dim,sizeof(double));
  73     histogram_error=calloc(b->entries*dim,sizeof(double));
  74     histogram_hi=calloc(b->entries*dim,sizeof(double));
  75     histogram_lo=calloc(b->entries*dim,sizeof(double));
  76   }else{
  77     fprintf(stderr,"Specify at least one codebook\n");
  78     exit(1);
  79   }
  80 }
  81
  82 void process_postprocess(codebook **b,char *basename){
  83   int i,j,k;
  84   char *buffer=alloca(strlen(basename)+80);
  85   codebook *bb=b[books-1];
  86
  87   fprintf(stderr,"Done.  Processed %ld data points for %ld entries:\n",
  88           (long)count,bb->entries);
  89   fprintf(stderr,"\tglobal mean amplitude: %g\n",
  90           meanamplitude_acc/(count*dim));
  91   fprintf(stderr,"\tglobal mean squared amplitude: %g\n",
  92           sqrt(meanamplitudesq_acc/(count*dim)));
  93
  94   fprintf(stderr,"\tglobal mean error: %g\n",
  95           meanerror_acc/(count*dim));
  96   fprintf(stderr,"\tglobal mean squared error: %g\n",
  97           sqrt(meanerrorsq_acc/(count*dim)));
  98   fprintf(stderr,"\tglobal mean deviation: %g\n",
  99           meandev_acc/(count*dim));
 100   {
 101     FILE *out;
 102
 103     sprintf(buffer,"%s-mse.m",basename);
 104     out=fopen(buffer,"w");
 105     if(!out){
 106       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 107       exit(1);
 108     }
 109
 110     for(i=0;i<bb->entries;i++){
 111       for(k=0;k<bb->dim;k++){
 112         fprintf(out,"%ld, %g, %g\n",
 113                 i*bb->dim+k,(bb->valuelist+i*bb->dim)[k],
 114                 sqrt((histogram_errorsq+i*bb->dim)[k]/histogram[i]));
 115       }
 116     }
 117     fclose(out);
 118
 119     sprintf(buffer,"%s-me.m",basename);
 120     out=fopen(buffer,"w");
 121     if(!out){
 122       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 123       exit(1);
 124     }
 125
 126     for(i=0;i<bb->entries;i++){
 127       for(k=0;k<bb->dim;k++){
 128         fprintf(out,"%ld, %g, %g\n",
 129                 i*bb->dim+k,(bb->valuelist+i*bb->dim)[k],
 130                 (histogram_error+i*bb->dim)[k]/histogram[i]);
 131       }
 132     }
 133     fclose(out);
 134
 135     sprintf(buffer,"%s-worst.m",basename);
 136     out=fopen(buffer,"w");
 137     if(!out){
 138       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 139       exit(1);
 140     }
 141
 142     for(i=0;i<bb->entries;i++){
 143       for(k=0;k<bb->dim;k++){
 144         fprintf(out,"%ld, %g, %g, %g\n",
 145                 i*bb->dim+k,(bb->valuelist+i*bb->dim)[k],
 146                 (bb->valuelist+i*bb->dim)[k]+(histogram_lo+i*bb->dim)[k],
 147                 (bb->valuelist+i*bb->dim)[k]+(histogram_hi+i*bb->dim)[k]);
 148       }
 149     }
 150     fclose(out);
 151   }
 152
 153   {
 154     FILE *out;
 155     long *index=alloca(sizeof(long)*bb->entries);
 156     sprintf(buffer,"%s-distance.m",basename);
 157     out=fopen(buffer,"w");
 158     if(!out){
 159       fprintf(stderr,"Could not open file %s for writing\n",buffer);
 160       exit(1);
 161     }
 162     for(j=0;j<bb->entries;j++){
 163       if(histogram[j])histogram_distance[j]/=histogram[j];
 164       index[j]=j;
 165     }
 166
 167     qsort(index,bb->entries,sizeof(long),histerrsort);
 168
 169     for(j=0;j<bb->entries;j++)
 170       for(k=0;k<histogram[index[j]];k++)
 171         fprintf(out,"%g,\n",histogram_distance[index[j]]);
 172     fclose(out);
 173
 174   }
 175
 176 }
 177
 178 void process_vector(codebook **bs,double *a){
 179   int bi;
 180   int i;
 181   double amplitude=0.;
 182   double distance=0.;
 183   double base=0.;
 184   int entry;
 185   double *e;
 186   memcpy(work,a,sizeof(double)*dim);
 187
 188   for(bi=0;bi<books;bi++){
 189     codebook *b=bs[bi];
 190     entry=codebook_entry(b,work);
 191     e=b->valuelist+b->dim*entry;
 192     for(i=0;i<b->dim;i++)work[i]-=e[i];
 193   }
 194
 195   for(i=0;i<dim;i++){
 196     double error=work[i];
 197     if(bs[0]->q_sequencep){
 198       amplitude=a[i]-base;
 199       base=a[i];
 200     }else
 201       amplitude=a[i];
 202
 203     meanamplitude_acc+=fabs(amplitude);
 204     meanamplitudesq_acc+=amplitude*amplitude;
 205     meanerror_acc+=fabs(error);
 206     meanerrorsq_acc+=error*error;
 207
 208     if(amplitude)
 209       meandev_acc+=fabs(error/amplitude);
 210     else
 211       meandev_acc+=fabs(error); /* yeah, yeah */
 212
 213     histogram_errorsq[entry*dim+i]+=error*error;
 214     histogram_error[entry*dim+i]+=fabs(error);
 215     if(histogram[entry]==0 || histogram_hi[entry*dim+i]<error)
 216       histogram_hi[entry*dim+i]=error;
 217     if(histogram[entry]==0 || histogram_lo[entry*dim+i]>error)
 218       histogram_lo[entry*dim+i]=error;
 219     distance+=error*error;
 220   }
 221
 222   histogram[entry]++;
 223   histogram_distance[entry]+=sqrt(distance);
 224
 225   if((long)(count++)%100)spinnit("working.... lines: ",count);
 226 }
 227
 228 void process_usage(void){
 229   fprintf(stderr,
 230           "usage: vqmetrics <codebook>.vqh datafile.vqd [datafile.vqd]...\n\n"
 231           "       data can be taken on stdin.  Output goes to output files:\n"
 232           "       basename-me.m:       gnuplot: mean error by entry value\n"
 233           "       basename-mse.m:      gnuplot: mean square error by entry value\n"
 234           "       basename-worst.m:    gnuplot: worst error by entry value\n"
 235           "       basename-distance.m: gnuplot file showing distance probability\n"
 236           "\n");
 237
 238 }