static char *linebuffer=NULL;
static int lbufsize=0;
-static char *rline(FILE *in,FILE *out,int pass){
+static char *rline(FILE *in,FILE *out){
long sofar=0;
if(feof(in))return NULL;
}
if(linebuffer[0]=='#'){
- if(pass)fprintf(out,"%s\n",linebuffer);
sofar=0;
}else{
return(linebuffer);
vqgen v;
vqbook b;
quant_meta q;
- int *quantlist=NULL;
+ long *quantlist=NULL;
- int entries=-1,dim=-1,dummy;
+ int entries=-1,dim=-1,aux=-1;
FILE *out=NULL;
FILE *in=NULL;
char *line,*name;
exit(1);
}
- ptr=strrchr(filename,'.');
+ ptr=strrchr(filename,'-');
if(ptr){
*ptr='\0';
- name=strdup(ptr);
- sprintf(ptr,".h");
+ name=strdup(filename);
+ sprintf(ptr,".vqh");
}else{
name=strdup(filename);
- strcat(filename,".h");
+ strcat(filename,".vqh");
}
out=fopen(filename,"w");
/* suck in the trained book */
/* read book type, but it doesn't matter */
- line=rline(in,out,1);
+ line=rline(in,out);
- line=rline(in,out,1);
- if(sscanf(line,"%d %d %d",&entries,&dim,&dummy)!=2){
+ line=rline(in,out);
+ if(sscanf(line,"%d %d %d",&entries,&dim,&aux)!=3){
fprintf(stderr,"Syntax error reading book file\n");
exit(1);
}
vqgen_init(&v,dim,0,entries,NULL,NULL);
/* quant */
- line=rline(in,out,1);
+ line=rline(in,out);
if(sscanf(line,"%ld %ld %d %d",&q.min,&q.delta,
&q.quant,&q.sequencep)!=4){
fprintf(stderr,"Syntax error reading book file\n");
/* save quant data; we don't want to requantize later as our method
is currently imperfect wrt repeated application */
i=0;
- quantlist=malloc(sizeof(int)*v.elements*v.entries);
+ quantlist=malloc(sizeof(long)*v.elements*v.entries);
for(j=0;j<entries;j++){
double a;
for(k=0;k<dim;k++){
- line=rline(in,out,0);
+ line=rline(in,out);
sscanf(line,"%lf",&a);
v.entrylist[i]=a;
quantlist[i++]=rint(a);
}
/* ignore bias */
- for(j=0;j<entries;j++)line=rline(in,out,0);
+ for(j=0;j<entries;j++)line=rline(in,out);
free(v.bias);
v.bias=NULL;
/* training points */
{
- double b[80];
+ double *b=alloca(sizeof(double)*(dim+aux));
i=0;
v.entries=0; /* hack to avoid reseeding */
while(1){
- for(k=0;k<dim && k<80;k++){
- line=rline(in,out,0);
+ for(k=0;k<dim+aux;k++){
+ line=rline(in,out);
if(!line)break;
sscanf(line,"%lf",b+k);
}
vqgen_unquantize(&v,&q);
/* build the book */
- vqsp_book(&v,&b);
+ vqsp_book(&v,&b,quantlist);
/* save the book in C header form */
+ fprintf(out,
+ "/********************************************************************\n"
+ " * *\n"
+ " * THIS FILE IS PART OF THE Ogg Vorbis SOFTWARE CODEC SOURCE CODE. *\n"
+ " * USE, DISTRIBUTION AND REPRODUCTION OF THIS SOURCE IS GOVERNED BY *\n"
+ " * THE GNU PUBLIC LICENSE 2, WHICH IS INCLUDED WITH THIS SOURCE. *\n"
+ " * PLEASE READ THESE TERMS DISTRIBUTING. *\n"
+ " * *\n"
+ " * THE OggSQUISH SOURCE CODE IS (C) COPYRIGHT 1994-1999 *\n"
+ " * by 1999 Monty <monty@xiph.org> and The XIPHOPHORUS Company *\n"
+ " * http://www.xiph.org/ *\n"
+ " * *\n"
+ " ********************************************************************\n"
+ "\n"
+ " function: static codebook autogenerated by vq/vqbuild\n"
+ "\n"
+ " ********************************************************************/\n\n");
+
+ fprintf(out,"#ifndef _V_%s_VQH_\n#define _V_%s_VQH_\n",name,name);
+ fprintf(out,"#include \"vqgen.h\"\n\n");
+
+ /* first, the static vectors, then the book structure to tie it together. */
+ /* quantlist */
+ fprintf(out,"static long _vq_quantlist_%s[] = {\n",name);
+ i=0;
+ for(j=0;j<entries;j++){
+ fprintf(out,"\t");
+ for(k=0;k<dim;k++)
+ fprintf(out,"%ld, ",b.quantlist[i++]);
+ fprintf(out,"\n");
+ }
+ fprintf(out,"};\n\n");
+ /* codelist */
+ fprintf(out,"static long _vq_codelist_%s[] = {\n",name);
+ for(j=0;j<entries;j++){
+ fprintf(out,"\t%ld,\n",b.codelist[j]);
+ }
+ fprintf(out,"};\n\n");
+ /* lengthlist */
+ fprintf(out,"static long _vq_lengthlist_%s[] = {\n",name);
+ for(j=0;j<entries;j++){
+ fprintf(out,"\t%ld,\n",b.lengthlist[j]);
+ }
+ fprintf(out,"};\n\n");
+
+ /* ptr0 */
+ fprintf(out,"static long _vq_ptr0_%s[] = {\n",name);
+ for(j=0;j<b.aux;j++){
+ fprintf(out,"\t%ld,\n",b.ptr0[j]);
+ }
+ fprintf(out,"};\n\n");
+
+ /* ptr1 */
+ fprintf(out,"static long _vq_ptr1_%s[] = {\n",name);
+ for(j=0;j<b.aux;j++){
+ fprintf(out,"\t%ld,\n",b.ptr1[j]);
+ }
+ fprintf(out,"};\n\n");
+
+ /* n */
+ fprintf(out,"static double _vq_n_%s[] = {\n",name);
+ i=0;
+ for(j=0;j<b.aux;j++){
+ fprintf(out,"\t");
+ for(k=0;k<dim;k++)
+ fprintf(out,"%g,",b.n[i++]);
+ fprintf(out,"\n");
+ }
+ fprintf(out,"};\n\n");
+
+ /* c */
+ fprintf(out,"static double _vq_c_%s[] = {\n",name);
+ for(j=0;j<b.aux;j++){
+ fprintf(out,"\t%g,\n",b.c[j]);
+ }
+ fprintf(out,"};\n\n");
+ /* tie it all together */
+ fprintf(out,"static vqbook _vq_book_%s = {\n",name);
+ fprintf(out,"\t%ld, %ld, %ld, %ld, %d, %d,\n",
+ b.dim,b.entries,q.min,q.delta,q.quant,q.sequencep);
+ fprintf(out,"\t0,\n"); /* valuelist */
+ fprintf(out,"\t_vq_quantlist_%s,\n",name);
+ fprintf(out,"\t_vq_codelist_%s,\n",name);
+ fprintf(out,"\t_vq_lengthlist_%s,\n",name);
+ fprintf(out,"\t_vq_ptr0_%s,\n",name);
+ fprintf(out,"\t_vq_ptr1_%s,\n",name);
+ fprintf(out,"\t_vq_n_%s,\n",name);
+ fprintf(out,"\t_vq_c_%s,\n",name);
+ fprintf(out,"\t%ld, %ld };\n\n",b.aux,b.aux);
+ fprintf(out,"\n#endif\n");
+ fclose(out);
exit(0);
}
features. */
double global_maxdel=M_PI;
-#define FUDGE ((global_maxdel*2.0)-weight[i])
+#define FUDGE (global_maxdel-weight[i])
double *weight=NULL;
double *vqext_weight(vqgen *v,double *p){
}
}
+ global_maxdel*=1.1;
weight=malloc(sizeof(double)*v->elements);
}
}
{
- double *b=alloca(dim*sizeof(double));
+ double *b=alloca((dim+vqext_aux)*sizeof(double));
i=0;
v.entries=0; /* hack to avoid reseeding */
while(1){
if(v->points==v->entries)_vqgen_seed(v);
}
+int directdsort(const void *a, const void *b){
+ double av=*((double *)a);
+ double bv=*((double *)b);
+ if(av>bv)return(-1);
+ return(1);
+}
+
double vqgen_iterate(vqgen *v){
long i,j,k;
double fdesired=(double)v->points/v->entries;
long desired=fdesired;
+ long desired2=desired*2;
double asserror=0.;
double meterror=0.;
double *new=malloc(sizeof(double)*v->entries*v->elements);
long *nearcount=malloc(v->entries*sizeof(long));
- double *nearbias=malloc(v->entries*desired*sizeof(double));
+ double *nearbias=malloc(v->entries*desired2*sizeof(double));
#ifdef NOISY
char buff[80];
for(j=0;j<v->entries;j++){
double thismetric;
- double *nearbiasptr=nearbias+desired*j;
- long k=nearcount[j]-1;
+ double *nearbiasptr=nearbias+desired2*j;
+ long k=nearcount[j];
/* 'thismetric' is to be the bias value necessary in the current
arrangement for entry j to capture point i */
/* use the primary entry as the threshhold */
thismetric=firstmetric-v->metric_func(v,_now(v,j),ppt);
}
-
- if(k>=0 && thismetric>nearbiasptr[k]){
-
- /* start at the end and search backward for where this entry
- belongs */
-
- for(;k>0;k--) if(nearbiasptr[k-1]>=thismetric)break;
-
- /* insert at k. Shift and inject. */
- memmove(nearbiasptr+k+1,nearbiasptr+k,(desired-k-1)*sizeof(double));
+
+ /* a cute two-stage delayed sorting hack */
+ if(k<desired){
nearbiasptr[k]=thismetric;
+ k++;
+ if(k==desired)
+ qsort(nearbiasptr,desired,sizeof(double),directdsort);
- if(nearcount[j]<desired)nearcount[j]++;
-
- }else{
- if(nearcount[j]<desired){
- /* we checked the thresh earlier. We know this is the
- last entry */
- nearbiasptr[nearcount[j]++]=thismetric;
+ }else if(thismetric>nearbiasptr[desired-1]){
+ nearbiasptr[k]=thismetric;
+ k++;
+ if(k==desired2){
+ qsort(nearbiasptr,desired2,sizeof(double),directdsort);
+ k=desired;
}
}
+ nearcount[j]=k;
}
}
/* inflate/deflate */
- for(i=0;i<v->entries;i++)
- v->bias[i]=nearbias[(i+1)*desired-1];
+ for(i=0;i<v->entries;i++){
+ double *nearbiasptr=nearbias+desired2*i;
+
+ /* due to the delayed sorting, we likely need to finish it off....*/
+ if(nearcount[i]>desired)
+ qsort(nearbiasptr,nearcount[i],sizeof(double),directdsort);
+
+ v->bias[i]=nearbiasptr[desired-1];
+ }
/* assign midpoints */
int quant; /* 0 < quant <= 16 */
int sequencep; /* bitflag */
- double *valuelist; /* list of dim*entries quant/actual entry values */
+ double *valuelist; /* list of dim*entries actual entry values */
+ long *quantlist; /* list of dim*entries quantized entry values */
long *codelist; /* list of bitstream codewords for each entry */
long *lengthlist; /* codeword lengths in bits */
extern double vqgen_iterate(vqgen *v);
extern void vqgen_unquantize(vqgen *v,quant_meta *q);
extern void vqgen_quantize(vqgen *v,quant_meta *q);
+extern long float24_pack(double val);
+extern double float24_unpack(long val);
-extern void vqsp_book(vqgen *v,vqbook *b);
+extern void vqsp_book(vqgen *v,vqbook *b,long *quantlist);
extern int vqenc_entry(vqbook *b,double *val);
#endif
extern double _dist_sq(vqgen *v,double *a, double *b);
/* goes through the split, but just counts it and returns a metric*/
-void vqsp_count(vqgen *v,long *entryindex,long entries,
+void vqsp_count(vqgen *v,long *membership,
+ long *entryindex,long entries,
long *pointindex,long points,
long *entryA,long *entryB,
- double *n, double c, double slack,
+ double *n, double c,
long *entriesA,long *entriesB,long *entriesC){
long i,j,k;
long A=0,B=0,C=0;
memset(entryA,0,sizeof(long)*entries);
memset(entryB,0,sizeof(long)*entries);
+ /* Do the points belonging to this cell occur on sideA, sideB or
+ both? */
+
for(i=0;i<points;i++){
double *ppt=_point(v,pointindex[i]);
- long firstentry=0;
- double firstmetric=_dist_sq(v,_now(v,entryindex[0]),ppt);
+ long firstentry=membership[i];
double position=-c;
-
- for(j=1;j<entries;j++){
- double thismetric=_dist_sq(v,_now(v,entryindex[j]),ppt);
- if(thismetric<firstmetric){
- firstmetric=thismetric;
- firstentry=j;
- }
- }
- /* count point split */
- for(k=0;k<v->elements;k++)
+ if(!entryA[firstentry] || !entryB[firstentry]){
+ for(k=0;k<v->elements;k++)
position+=ppt[k]*n[k];
- if(position>0.){
- entryA[firstentry]++;
- }else{
- entryB[firstentry]++;
+ if(position>0.)
+ entryA[firstentry]=1;
+ else
+ entryB[firstentry]=1;
}
}
- /* look to see if entries are in the slack zone */
/* The entry splitting isn't total, so that storage has to be
allocated for recursion. Reuse the entryA/entryB vectors */
for(j=0;j<entries;j++){
- long total=entryA[j]+entryB[j];
- if((double)entryA[j]/total<slack){
- entryA[j]=0;
- }else if((double)entryB[j]/total<slack){
- entryB[j]=0;
- }
if(entryA[j] && entryB[j])C++;
if(entryA[j])entryA[A++]=entryindex[j];
if(entryB[j])entryB[B++]=entryindex[j];
int lp_split(vqgen *v,vqbook *b,
long *entryindex,long entries,
long *pointindex,long points,
- long depth,double slack){
+ long depth){
/* The encoder, regardless of book, will be using a straight
euclidian distance-to-point metric to determine closest point.
long pointsA=0;
long i,j,k;
+ long *membership=malloc(sizeof(long)*points);
+
+ /* which cells do points belong to? Do this before n^2 best pair chooser. */
+
+ for(i=0;i<points;i++){
+ double *ppt=_point(v,pointindex[i]);
+ long firstentry=0;
+ double firstmetric=_dist_sq(v,_now(v,entryindex[0]),ppt);
+
+ if(points*entries>64*1024)spinnit();
+
+ for(j=1;j<entries;j++){
+ double thismetric=_dist_sq(v,_now(v,entryindex[j]),ppt);
+ if(thismetric<firstmetric){
+ firstmetric=thismetric;
+ firstentry=j;
+ }
+ }
+
+ membership[i]=firstentry;
+ }
+
p=alloca(sizeof(double)*v->elements);
q=alloca(sizeof(double)*v->elements);
n=alloca(sizeof(double)*v->elements);
/* more than one way to do this part. For small sets, we can brute
force it. */
- if(entries<32){
+ {
/* try every pair possibility */
double best=0;
long besti=0;
for(j=i+1;j<entries;j++){
spinnit();
pq_in_out(v,n,&c,_now(v,entryindex[i]),_now(v,entryindex[j]));
- vqsp_count(v,entryindex,entries,
- pointindex,points,
- entryA,entryB,
- n, c, slack,
- &entriesA,&entriesB,&entriesC);
+ vqsp_count(v,membership,
+ entryindex,entries,
+ pointindex,points,
+ entryA,entryB,
+ n, c,
+ &entriesA,&entriesB,&entriesC);
this=(entriesA-entriesC)*(entriesB-entriesC);
if(this>best){
}
}
pq_in_out(v,n,&c,_now(v,entryindex[besti]),_now(v,entryindex[bestj]));
- }else{
- double best=0.;
- long bestj=0;
-
- /* try COG/normal and furthest pairs */
- /* medianpoint */
- for(k=0;k<v->elements;k++){
- spinnit();
-
- /* just sort the index array */
- sortvals=v->entrylist+k;
- els=v->elements;
- qsort(entryindex,entries,sizeof(long),dascsort);
- if(entries&0x1){
- p[k]=v->entrylist[(entryindex[entries/2])*v->elements+k];
- }else{
- p[k]=(v->entrylist[(entryindex[entries/2])*v->elements+k]+
- v->entrylist[(entryindex[entries/2-1])*v->elements+k])/2.;
- }
- }
-
- spinnit();
-
- /* try every normal, but just for distance */
- for(j=0;j<entries;j++){
- double *ppj=_now(v,entryindex[j]);
- double this=_dist_sq(v,p,ppj);
- if(this>best){
- best=this;
- bestj=j;
- }
- }
-
- pq_center_out(v,n,&c,p,_point(v,pointindex[bestj]));
-
-
}
/* find cells enclosing points */
/* count A/B points */
- vqsp_count(v,entryindex,entries,
- pointindex,points,
- entryA,entryB,
- n, c, slack,
- &entriesA,&entriesB,&entriesC);
+ vqsp_count(v,membership,
+ entryindex,entries,
+ pointindex,points,
+ entryA,entryB,
+ n, c,
+ &entriesA,&entriesB,&entriesC);
- /* the point index is split evenly, so we do an Order n
- rearrangement into A first/B last and just pass it on */
+ free(membership);
+
+ /* the point index is split, so we do an Order n rearrangement into
+ A first/B last and just pass it on */
{
long Aptr=0;
long Bptr=points-1;
b->ptr0[thisaux]=entryA[0];
}else{
b->ptr0[thisaux]= -b->aux;
- ret=lp_split(v,b,entryA,entriesA,pointindex,pointsA,depth+1,slack);
+ ret=lp_split(v,b,entryA,entriesA,pointindex,pointsA,depth+1);
}
if(entriesB==1){
ret++;
}else{
b->ptr1[thisaux]= -b->aux;
ret+=lp_split(v,b,entryB,entriesB,pointindex+pointsA,points-pointsA,
- depth+1,slack);
+ depth+1);
}
}
free(entryA);
return(ret);
}
-void vqsp_book(vqgen *v, vqbook *b){
- long *entryindex=malloc(sizeof(double)*v->entries);
- long *pointindex=malloc(sizeof(double)*v->points);
+void vqsp_book(vqgen *v, vqbook *b, long *quantlist){
+ long *entryindex=malloc(sizeof(long)*v->entries);
+ long *pointindex=malloc(sizeof(long)*v->points);
+ long *membership=malloc(sizeof(long)*v->points);
long i,j;
memset(b,0,sizeof(vqbook));
b->c=malloc(sizeof(double)*b->alloc);
b->lengthlist=calloc(b->entries,sizeof(long));
+ /* which cells do points belong to? Only do this once. */
+
+ for(i=0;i<v->points;i++){
+ double *ppt=_point(v,i);
+ long firstentry=0;
+ double firstmetric=_dist_sq(v,_now(v,0),ppt);
+
+ for(j=1;j<v->entries;j++){
+ double thismetric=_dist_sq(v,_now(v,j),ppt);
+ if(thismetric<firstmetric){
+ firstmetric=thismetric;
+ firstentry=j;
+ }
+ }
+
+ membership[i]=firstentry;
+ }
+
+
/* first, generate the encoding decision heirarchy */
fprintf(stderr,"Total leaves: %d\n",
- lp_split(v,b,entryindex,v->entries, pointindex,v->points,0,0));
+ lp_split(v,b,entryindex,v->entries,
+ pointindex,v->points,0));
free(entryindex);
free(pointindex);
/* n and c stay unchanged */
for(i=0;i<b->entries;i++)revindex[index[i]]=i;
for(i=0;i<b->aux;i++){
- if(b->ptr0[i]>=0)b->ptr0[i]=revindex[i];
- if(b->ptr1[i]>=0)b->ptr1[i]=revindex[i];
+ if(b->ptr0[i]>=0)b->ptr0[i]=revindex[b->ptr0[i]];
+ if(b->ptr1[i]>=0)b->ptr1[i]=revindex[b->ptr1[i]];
}
free(revindex);
/* map lengthlist and vallist with index */
b->lengthlist=calloc(b->entries,sizeof(long));
b->valuelist=malloc(sizeof(double)*b->entries*b->dim);
+ b->quantlist=malloc(sizeof(long)*b->entries*b->dim);
for(i=0;i<b->entries;i++){
long e=index[i];
- for(k=0;k<b->dim;k++)
+ for(k=0;k<b->dim;k++){
b->valuelist[i*b->dim+k]=v->entrylist[e*b->dim+k];
+ b->quantlist[i*b->dim+k]=quantlist[e*b->dim+k];
+ }
b->lengthlist[i]=wordlen[e];
}