Incremental commit of untested code changes to codebook training

(added density limiting via a quantization mesh)

Monty

svn path=/branches/monty_branch_20000524/vorbis/; revision=436

3 files changed, 184 insertions, 74 deletions

diff --git a/vq/train.c b/vq/train.c
index 735991cb..f68b5fef 100644
--- a/vq/train.c
+++ b/vq/train.c
@@ -12,7 +12,7 @@
  ********************************************************************
 
  function: utility main for training codebooks
- last mod: $Id: train.c,v 1.17.2.1 2000/06/01 12:03:04 xiphmont Exp $
+ last mod: $Id: train.c,v 1.17.2.2 2000/06/09 03:59:24 xiphmont Exp $
 
  ********************************************************************/
 
@@ -53,8 +53,10 @@ static void usage(void){
 	  "         -s[ubvector]  <start[,num]>\n"
 	  "         -e[rror]      <desired_error>\n"
 	  "         -i[terations] <maxiterations>\n"
-	  "         -d[istance]   desired minimum cell radius from midpoint\n"
+	  "         -d[istance]   quantization mesh spacing for density limitation\n"
 	  "         -b <dummy>    eliminate cell size biasing; use normal LBG\n\n"
+	  "         -c <dummy>    Use centroid (not median) midpoints\n"
+
 	  "examples:\n"
 	  "   train a new codebook to 1%% tolerance on datafile 'foo':\n"
 	  "      xxxvqtrain book -p 256,6,8 -e .01 foo\n"
@@ -80,6 +82,7 @@ int main(int argc,char *argv[]){
   double desired=.05,mindist=0.;
   int iter=1000;
   int biasp=1;
+  int centroid=0;
 
   FILE *out=NULL;
   char *line;
@@ -132,7 +135,7 @@ int main(int argc,char *argv[]){
       }
       
       vqgen_init(&v,dim,vqext_aux,entries,mindist,
-		 vqext_metric,vqext_weight);
+		 vqext_metric,vqext_weight,centroid);
       init=1;
       
       /* quant setup */
@@ -154,7 +157,7 @@ int main(int argc,char *argv[]){
       }      
       vqgen_unquantize(&v,&q);
 
-      /* bias, points */
+      /* bias */
       i=0;
       for(j=0;j<entries;j++){
 	line=rline(in,out,0);
@@ -162,10 +165,10 @@ int main(int argc,char *argv[]){
 	v.bias[i++]=a;
       }
       
+      v.seed=1;
       {
 	double *b=alloca((dim+vqext_aux)*sizeof(double));
 	i=0;
-	v.entries=0; /* hack to avoid reseeding */
 	while(1){
 	  for(k=0;k<dim+vqext_aux;k++){
 	    line=rline(in,out,0);
@@ -175,7 +178,6 @@ int main(int argc,char *argv[]){
 	  if(feof(in))break;
 	  vqgen_addpoint(&v,b,b+dim);
 	}
-	v.entries=entries;
       }
       
       fclose(in);
@@ -219,6 +221,9 @@ int main(int argc,char *argv[]){
       case 'b':
 	biasp=0;
 	break;
+      case 'c':
+	centroid=1;
+	break;
       default:
 	fprintf(stderr,"Unknown option %s\n",argv[0]);
 	exit(1);
@@ -236,7 +241,7 @@ int main(int argc,char *argv[]){
 	  exit(1);
 	}
 	vqgen_init(&v,dim,vqext_aux,entries,mindist,
-		   vqext_metric,vqext_weight);
+		   vqext_metric,vqext_weight,centroid);
 	init=1;
       }
 
@@ -331,13 +336,21 @@ int main(int argc,char *argv[]){
   for(j=0;j<entries;j++)
     fprintf(out,"%f\n",v.bias[i++]);
 
+  /* we may have done the density limiting mesh trick; refetch the
+     training points from the temp file */
+
+  rewind(v.asciipoints);
   fprintf(out,"# points---\n");
-  i=0;
-  for(j=0;j<v.points;j++)
-    for(k=0;k<dim+vqext_aux;k++)
-      fprintf(out,"%f\n",v.pointlist[i++]);
+  {
+    /* sloppy, no error handling */
+    long bytes;
+    char buff[4096];
+    while((bytes=fread(buff,1,4096,v.asciipoints)))
+      while(bytes)bytes-=fwrite(buff,1,bytes,out);
+  }
 
   fclose(out);
+  fclose(v.asciipoints);
 
   vqgen_unquantize(&v,&q);
   vqgen_cellmetric(&v);
diff --git a/vq/vqgen.c b/vq/vqgen.c
index 797081c0..1aece36e 100644
--- a/vq/vqgen.c
+++ b/vq/vqgen.c
@@ -12,7 +12,7 @@
  ********************************************************************
 
  function: train a VQ codebook 
- last mod: $Id: vqgen.c,v 1.31.2.2 2000/06/09 00:28:33 xiphmont Exp $
+ last mod: $Id: vqgen.c,v 1.31.2.3 2000/06/09 03:59:24 xiphmont Exp $
 
  ********************************************************************/
 
@@ -25,9 +25,6 @@
    considering the undertaking is almost withering.  For now, we brute
    force it all */
 
-#define  CENTROID 1
-#undef   MEDIAN  
-
 #include <stdlib.h>
 #include <stdio.h>
 #include <math.h>
@@ -82,6 +79,7 @@ void _vqgen_seed(vqgen *v){
   long i;
   for(i=0;i<v->entries;i++)
     memcpy(_now(v,i),_point(v,i),sizeof(double)*v->elements);
+  v->seed=1;
 }
 
 int directdsort(const void *a, const void *b){
@@ -244,9 +242,10 @@ void vqgen_unquantize(vqgen *v,quant_meta *q){
 
 void vqgen_init(vqgen *v,int elements,int aux,int entries,double mindist,
 		double  (*metric)(vqgen *,double *, double *),
-		double *(*weight)(vqgen *,double *)){
+		double *(*weight)(vqgen *,double *),int centroid){
   memset(v,0,sizeof(vqgen));
 
+  v->centroid=centroid;
   v->elements=elements;
   v->aux=aux;
   v->mindist=mindist;
@@ -267,25 +266,51 @@ void vqgen_init(vqgen *v,int elements,int aux,int entries,double mindist,
   else
     v->weight_func=_weight_null;
 
+  v->asciipoints=tmpfile();
+
 }
 
 void vqgen_addpoint(vqgen *v, double *p,double *a){
+  int k;
+  for(k=0;k<v->elements;k++)
+    fprintf(v->asciipoints,"%f\n",p[k]);
+  for(k=0;k<v->aux;k++)
+    fprintf(v->asciipoints,"%f\n",a[k]);
+
   if(v->points>=v->allocated){
     v->allocated*=2;
     v->pointlist=realloc(v->pointlist,v->allocated*(v->elements+v->aux)*
 			 sizeof(double));
   }
-  
+
   memcpy(_point(v,v->points),p,sizeof(double)*v->elements);
   if(v->aux)memcpy(_point(v,v->points)+v->elements,a,sizeof(double)*v->aux);
+ 
+  /* deal with the density mesh if it's selected */
+  if(v->mindist>0.){
+    
+    /* quantize to the mesh */
+    for(k=0;k<v->elements+v->aux;k++)
+      _point(v,v->points)[k]=
+	rint(_point(v,v->points)[k]/v->mindist)*v->mindist;
+
+    /* look for prexisting member of this mesh cell */
+    /* just skip out in the event we discard this one */
+    for(k=0;k<v->points;k++)
+      if(!memcmp(_point(v,v->points),
+		 _point(v,k),
+		 (v->elements+v->aux)*sizeof(double))){
+	spinnit("loading... ",v->points);
+	return();
+      }
+  }
   v->points++;
-  if(v->points==v->entries)_vqgen_seed(v);
+  if(v->points==v->entries && !v->seed)_vqgen_seed(v);
   if(!(v->points&0xff))spinnit("loading... ",v->points);
 }
 
 double vqgen_iterate(vqgen *v,int biasp){
   long   i,j,k;
-  long   biasable;
 
   double fdesired=(double)v->points/v->entries;
   long  desired=fdesired;
@@ -320,7 +345,6 @@ double vqgen_iterate(vqgen *v,int biasp){
   /*memset(v->bias,0,sizeof(double)*v->entries);*/
   memset(nearcount,0,sizeof(long)*v->entries);
   memset(v->assigned,0,sizeof(long)*v->entries);
-  biasable=0;
   if(biasp){
     for(i=0;i<v->points;i++){
       double *ppt=v->weight_func(v,_point(v,i));
@@ -377,31 +401,26 @@ double vqgen_iterate(vqgen *v,int biasp){
 	   cells in a codebook to be roughly some minimum size (as with
 	   the low resolution residue books) */
 	
-	if(localmetric>=v->mindist){
+	/* a cute two-stage delayed sorting hack */
+	if(k<desired){
+	  nearbiasptr[k]=thismetric;
+	  k++;
+	  if(k==desired){
+	    spinnit("biasing... ",v->points+v->points+v->entries-i);
+	    qsort(nearbiasptr,desired,sizeof(double),directdsort);
+	  }
 	  
-	  /* a cute two-stage delayed sorting hack */
-	  if(k<desired){
-	    nearbiasptr[k]=thismetric;
-	    k++;
-	    if(k==desired){
-	      spinnit("biasing... ",v->points+v->points+v->entries-i);
-	      qsort(nearbiasptr,desired,sizeof(double),directdsort);
-	    }
-	    
-	  }else if(thismetric>nearbiasptr[desired-1]){
-	    nearbiasptr[k]=thismetric;
-	    k++;
-	    if(k==desired2){
-	      spinnit("biasing... ",v->points+v->points+v->entries-i);
-	      qsort(nearbiasptr,desired2,sizeof(double),directdsort);
-	      k=desired;
-	    }
+	}else if(thismetric>nearbiasptr[desired-1]){
+	  nearbiasptr[k]=thismetric;
+	  k++;
+	  if(k==desired2){
+	    spinnit("biasing... ",v->points+v->points+v->entries-i);
+	    qsort(nearbiasptr,desired2,sizeof(double),directdsort);
+	    k=desired;
 	  }
-	  nearcount[j]=k;
-	}else
-	  biasflag=0;
+	}
+	nearcount[j]=k;
       }
-      biasable+=biasflag;
     }
     
     /* inflate/deflate */
@@ -415,16 +434,8 @@ double vqgen_iterate(vqgen *v,int biasp){
       if(nearcount[i]>desired)
 	qsort(nearbiasptr,nearcount[i],sizeof(double),directdsort);
       
-      /* desired is the *maximum* bias queue size.  If we're using
-	 minimum distance, we're not interested in the max size... we're
-	 interested in the biasable number of points */
-      {
-	long localdesired=(double)biasable/v->entries;
-	if(localdesired)
-	  v->bias[i]=nearbiasptr[localdesired-1];
-	else
-	  v->bias[i]=nearbiasptr[0];
-      }
+      v->bias[i]=nearbiasptr[desired-1];
+
     }
   }else{ 
     memset(v->bias,0,v->entries*sizeof(double));
@@ -456,33 +467,35 @@ double vqgen_iterate(vqgen *v,int biasp){
 
     firstmetric-=v->bias[firstentry];
     meterror+=firstmetric;
-#ifdef MEDIAN
-    /* set up midpoints for next iter */
-    if(v->assigned[j]++){
-      for(k=0;k<v->elements;k++)
-	vN(new,j)[k]+=ppt[k];
-      if(firstmetric>v->max[firstentry])v->max[firstentry]=firstmetric;
-    }else{
-      for(k=0;k<v->elements;k++)
-	vN(new,j)[k]=ppt[k];
-      v->max[firstentry]=firstmetric;
-    }
-#else /* centroid */
-    if(v->assigned[j]++){
-      for(k=0;k<v->elements;k++){
-	if(vN(new,j)[k]>ppt[k])vN(new,j)[k]=ppt[k];
-	if(vN(new2,j)[k]<ppt[k])vN(new2,j)[k]=ppt[k];
+
+    if(v->centroid==0){
+      /* set up midpoints for next iter */
+      if(v->assigned[j]++){
+	for(k=0;k<v->elements;k++)
+	  vN(new,j)[k]+=ppt[k];
+	if(firstmetric>v->max[firstentry])v->max[firstentry]=firstmetric;
+      }else{
+	for(k=0;k<v->elements;k++)
+	  vN(new,j)[k]=ppt[k];
+	v->max[firstentry]=firstmetric;
       }
-      if(firstmetric>v->max[firstentry])v->max[firstentry]=firstmetric;
     }else{
-      for(k=0;k<v->elements;k++){
-	vN(new,j)[k]=ppt[k];
-	vN(new2,j)[k]=ppt[k];
+      /* centroid */
+      if(v->assigned[j]++){
+	for(k=0;k<v->elements;k++){
+	  if(vN(new,j)[k]>ppt[k])vN(new,j)[k]=ppt[k];
+	  if(vN(new2,j)[k]<ppt[k])vN(new2,j)[k]=ppt[k];
+	}
+	if(firstmetric>v->max[firstentry])v->max[firstentry]=firstmetric;
+      }else{
+	for(k=0;k<v->elements;k++){
+	  vN(new,j)[k]=ppt[k];
+	  vN(new2,j)[k]=ppt[k];
+	}
+	v->max[firstentry]=firstmetric;
       }
-      v->max[firstentry]=firstmetric;
     }
 
-
 #endif
   }
 
diff --git a/vq/vqgen.h b/vq/vqgen.h
new file mode 100644
index 00000000..c8afecb1
--- /dev/null
+++ b/vq/vqgen.h
@@ -0,0 +1,84 @@
+/********************************************************************
+ *                                                                  *
+ * THIS FILE IS PART OF THE Ogg Vorbis SOFTWARE CODEC SOURCE CODE.  *
+ * USE, DISTRIBUTION AND REPRODUCTION OF THIS SOURCE IS GOVERNED BY *
+ * THE GNU PUBLIC LICENSE 2, WHICH IS INCLUDED WITH THIS SOURCE.    *
+ * PLEASE READ THESE TERMS DISTRIBUTING.                            *
+ *                                                                  *
+ * THE OggSQUISH SOURCE CODE IS (C) COPYRIGHT 1994-2000             *
+ * by Monty <monty@xiph.org> and The XIPHOPHORUS Company            *
+ * http://www.xiph.org/                                             *
+ *                                                                  *
+ ********************************************************************
+
+ function: build a VQ codebook 
+ last mod: $Id: vqgen.h,v 1.13.2.1 2000/06/09 03:59:25 xiphmont Exp $
+
+ ********************************************************************/
+
+#ifndef _VQGEN_H_
+#define _VQGEN_H_
+
+typedef struct vqgen{
+  int seed;
+  int it;
+  int elements;
+
+  int aux;
+  double mindist;
+  int centroid;
+
+  /* point cache */
+  double *pointlist; 
+  long   points;
+  long   allocated;
+
+  /* entries */
+  double *entrylist;
+  long   *assigned;
+  double *bias;
+  long   entries;
+  double *max;
+  
+  double  (*metric_func) (struct vqgen *v,double *entry,double *point);
+  double *(*weight_func) (struct vqgen *v,double *point);
+
+  FILE *asciipoints;
+} vqgen;
+
+typedef struct {
+  long   min;       /* packed 24 bit float */       
+  long   delta;     /* packed 24 bit float */       
+  int    quant;     /* 0 < quant <= 16 */
+  int    sequencep; /* bitflag */
+} quant_meta;
+
+static inline double *_point(vqgen *v,long ptr){
+  return v->pointlist+((v->elements+v->aux)*ptr);
+}
+
+static inline double *_aux(vqgen *v,long ptr){
+  return _point(v,ptr)+v->aux;
+}
+
+static inline double *_now(vqgen *v,long ptr){
+  return v->entrylist+(v->elements*ptr);
+}
+
+extern void vqgen_init(vqgen *v,
+		       int elements,int aux,int entries,double mindist,
+		       double  (*metric)(vqgen *,double *, double *),
+		       double *(*weight)(vqgen *,double *),int centroid);
+extern void vqgen_addpoint(vqgen *v, double *p,double *aux);
+
+extern double vqgen_iterate(vqgen *v,int biasp);
+extern void vqgen_unquantize(vqgen *v,quant_meta *q);
+extern void vqgen_quantize(vqgen *v,quant_meta *q);
+extern void vqgen_cellmetric(vqgen *v);
+
+#endif
+
+
+
+
+