lib/codebook.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: basic codebook pack/unpack/code/decode operations
  15  last mod: $Id: codebook.c,v 1.11 2000/02/23 09:24:26 xiphmont Exp $
  16
  17  ********************************************************************/
  18
  19 #include <stdlib.h>
  20 #include <math.h>
  21 #include "vorbis/codec.h"
  22 #include "vorbis/codebook.h"
  23 #include "bitwise.h"
  24 #include "bookinternal.h"
  25
  26 /**** pack/unpack helpers ******************************************/
  27 static int ilog(unsigned int v){
  28   int ret=0;
  29   while(v){
  30     ret++;
  31     v>>=1;
  32   }
  33   return(ret);
  34 }
  35
  36 /* code that packs the 24 bit float can be found in vq/bookutil.c */
  37
  38 static double _float24_unpack(long val){
  39   double mant=val&0x3ffff;
  40   double sign=val&0x800000;
  41   double exp =(val&0x7c0000)>>18;
  42   if(sign)mant= -mant;
  43   return(ldexp(mant,exp-17-VQ_FEXP_BIAS));
  44 }
  45
  46 /* given a list of word lengths, generate a list of codewords.  Works
  47    for length ordered or unordered, always assigns the lowest valued
  48    codewords first */
  49 long *_make_words(long *l,long n){
  50   long i,j;
  51   long marker[33];
  52   long *r=malloc(n*sizeof(long));
  53   memset(marker,0,sizeof(marker));
  54
  55   for(i=0;i<n;i++){
  56     long length=l[i];
  57     long entry=marker[length];
  58
  59     /* when we claim a node for an entry, we also claim the nodes
  60        below it (pruning off the imagined tree that may have dangled
  61        from it) as well as blocking the use of any nodes directly
  62        above for leaves */
  63
  64     /* update ourself */
  65     if(length<32 && (entry>>length)){
  66       /* error condition; the lengths must specify an overpopulated tree */
  67       free(r);
  68       return(NULL);
  69     }
  70     r[i]=entry;
  71
  72     /* Look to see if the next shorter marker points to the node
  73        above. if so, update it and repeat.  */
  74     {
  75       for(j=length;j>0;j--){
  76
  77         if(marker[j]&1){
  78           /* have to jump branches */
  79           if(j==1)
  80             marker[1]++;
  81           else
  82             marker[j]=marker[j-1]<<1;
  83           break; /* invariant says next upper marker would already
  84                     have been moved if it was on the same path */
  85         }
  86         marker[j]++;
  87       }
  88     }
  89
  90     /* prune the tree; the implicit invariant says all the longer
  91        markers were dangling from our just-taken node.  Dangle them
  92        from our *new* node. */
  93     for(j=length+1;j<33;j++)
  94       if((marker[j]>>1) == entry){
  95         entry=marker[j];
  96         marker[j]=marker[j-1]<<1;
  97       }else
  98         break;
  99   }
 100
 101   /* bitreverse the words because our bitwise packer/unpacker is LSb
 102      endian */
 103   for(i=0;i<n;i++){
 104     long temp=0;
 105     for(j=0;j<l[i];j++){
 106       temp<<=1;
 107       temp|=(r[i]>>j)&1;
 108     }
 109     r[i]=temp;
 110   }
 111
 112   return(r);
 113 }
 114
 115 /* build the decode helper tree from the codewords */
 116 decode_aux *_make_decode_tree(codebook *c){
 117   const static_codebook *s=c->c;
 118   long top=0,i,j;
 119   decode_aux *t=malloc(sizeof(decode_aux));
 120   long *ptr0=t->ptr0=calloc(c->entries*2,sizeof(long));
 121   long *ptr1=t->ptr1=calloc(c->entries*2,sizeof(long));
 122   long *codelist=_make_words(s->lengthlist,s->entries);
 123
 124   if(codelist==NULL)return(NULL);
 125   t->aux=c->entries*2;
 126
 127   for(i=0;i<c->entries;i++){
 128     long ptr=0;
 129     for(j=0;j<s->lengthlist[i]-1;j++){
 130       int bit=(codelist[i]>>j)&1;
 131       if(!bit){
 132         if(!ptr0[ptr])
 133           ptr0[ptr]= ++top;
 134         ptr=ptr0[ptr];
 135       }else{
 136         if(!ptr1[ptr])
 137           ptr1[ptr]= ++top;
 138         ptr=ptr1[ptr];
 139       }
 140     }
 141     if(!((codelist[i]>>j)&1))
 142       ptr0[ptr]=-i;
 143     else
 144       ptr1[ptr]=-i;
 145   }
 146   free(codelist);
 147   return(t);
 148 }
 149
 150 /* unpack the quantized list of values for encode/decode ***********/
 151 static double *_book_unquantize(const static_codebook *b){
 152   long j,k;
 153   if(b->quantlist){
 154     double mindel=_float24_unpack(b->q_min);
 155     double delta=_float24_unpack(b->q_delta);
 156     double *r=malloc(sizeof(double)*b->entries*b->dim);
 157
 158     for(j=0;j<b->entries;j++){
 159       double last=0.;
 160       for(k=0;k<b->dim;k++){
 161         double val=b->quantlist[j*b->dim+k]*delta+last+mindel;
 162         r[j*b->dim+k]=val;
 163         if(b->q_sequencep)last=val;
 164       }
 165     }
 166     return(r);
 167   }else
 168     return(NULL);
 169 }
 170
 171 void vorbis_staticbook_clear(static_codebook *b){
 172   if(b->quantlist)free(b->quantlist);
 173   if(b->lengthlist)free(b->lengthlist);
 174   if(b->encode_tree){
 175     free(b->encode_tree->ptr0);
 176     free(b->encode_tree->ptr1);
 177     free(b->encode_tree->p);
 178     free(b->encode_tree->q);
 179     memset(b->encode_tree,0,sizeof(encode_aux));
 180     free(b->encode_tree);
 181   }
 182   memset(b,0,sizeof(static_codebook));
 183 }
 184
 185 void vorbis_book_clear(codebook *b){
 186   /* static book is not cleared; we're likely called on the lookup and
 187      the static codebook belongs to the info struct */
 188   if(b->decode_tree){
 189     free(b->decode_tree->ptr0);
 190     free(b->decode_tree->ptr1);
 191     memset(b->decode_tree,0,sizeof(decode_aux));
 192     free(b->decode_tree);
 193   }
 194   if(b->valuelist)free(b->valuelist);
 195   if(b->codelist)free(b->codelist);
 196   memset(b,0,sizeof(codebook));
 197 }
 198
 199 int vorbis_book_init_encode(codebook *c,const static_codebook *s){
 200   memset(c,0,sizeof(codebook));
 201   c->c=s;
 202   c->entries=s->entries;
 203   c->dim=s->dim;
 204   c->codelist=_make_words(s->lengthlist,s->entries);
 205   c->valuelist=_book_unquantize(s);
 206   return(0);
 207 }
 208
 209 int vorbis_book_init_decode(codebook *c,const static_codebook *s){
 210   memset(c,0,sizeof(codebook));
 211   c->c=s;
 212   c->entries=s->entries;
 213   c->dim=s->dim;
 214   c->valuelist=_book_unquantize(s);
 215   c->decode_tree=_make_decode_tree(c);
 216   if(c->decode_tree==NULL)goto err_out;
 217   return(0);
 218  err_out:
 219   vorbis_book_clear(c);
 220   return(-1);
 221 }
 222
 223 /* packs the given codebook into the bitstream **************************/
 224
 225 int vorbis_staticbook_pack(const static_codebook *c,oggpack_buffer *opb){
 226   long i,j;
 227   int ordered=0;
 228
 229   /* first the basic parameters */
 230   _oggpack_write(opb,0x564342,24);
 231   _oggpack_write(opb,c->dim,16);
 232   _oggpack_write(opb,c->entries,24);
 233
 234   /* pack the codewords.  There are two packings; length ordered and
 235      length random.  Decide between the two now. */
 236
 237   for(i=1;i<c->entries;i++)
 238     if(c->lengthlist[i]<c->lengthlist[i-1])break;
 239   if(i==c->entries)ordered=1;
 240
 241   if(ordered){
 242     /* length ordered.  We only need to say how many codewords of
 243        each length.  The actual codewords are generated
 244        deterministically */
 245
 246     long count=0;
 247     _oggpack_write(opb,1,1);  /* ordered */
 248     _oggpack_write(opb,c->lengthlist[0]-1,5); /* 1 to 32 */
 249
 250     for(i=1;i<c->entries;i++){
 251       long this=c->lengthlist[i];
 252       long last=c->lengthlist[i-1];
 253       if(this>last){
 254         for(j=last;j<this;j++){
 255           _oggpack_write(opb,i-count,ilog(c->entries-count));
 256           count=i;
 257         }
 258       }
 259     }
 260     _oggpack_write(opb,i-count,ilog(c->entries-count));
 261
 262   }else{
 263     /* length random.  Again, we don't code the codeword itself, just
 264        the length.  This time, though, we have to encode each length */
 265     _oggpack_write(opb,0,1);   /* unordered */
 266     for(i=0;i<c->entries;i++)
 267       _oggpack_write(opb,c->lengthlist[i]-1,5);
 268   }
 269
 270   /* is the entry number the desired return value, or do we have a
 271      mapping? */
 272   if(c->quantlist){
 273     /* we have a mapping.  bundle it out. */
 274     _oggpack_write(opb,1,1);
 275
 276     /* values that define the dequantization */
 277     _oggpack_write(opb,c->q_min,24);
 278     _oggpack_write(opb,c->q_delta,24);
 279     _oggpack_write(opb,c->q_quant-1,4);
 280     _oggpack_write(opb,c->q_sequencep,1);
 281
 282     /* quantized values */
 283     for(i=0;i<c->entries*c->dim;i++)
 284       _oggpack_write(opb,c->quantlist[i],c->q_quant);
 285
 286   }else{
 287     /* no mapping. */
 288     _oggpack_write(opb,0,1);
 289   }
 290
 291   return(0);
 292 }
 293
 294 /* unpacks a codebook from the packet buffer into the codebook struct,
 295    readies the codebook auxiliary structures for decode *************/
 296 int vorbis_staticbook_unpack(oggpack_buffer *opb,static_codebook *s){
 297   long i,j;
 298   memset(s,0,sizeof(static_codebook));
 299
 300   /* make sure alignment is correct */
 301   if(_oggpack_read(opb,24)!=0x564342)goto _eofout;
 302
 303   /* first the basic parameters */
 304   s->dim=_oggpack_read(opb,16);
 305   s->entries=_oggpack_read(opb,24);
 306   if(s->entries==-1)goto _eofout;
 307
 308   /* codeword ordering.... length ordered or unordered? */
 309   switch(_oggpack_read(opb,1)){
 310   case 0:
 311     /* unordered */
 312     s->lengthlist=malloc(sizeof(long)*s->entries);
 313     for(i=0;i<s->entries;i++){
 314       long num=_oggpack_read(opb,5);
 315       if(num==-1)goto _eofout;
 316       s->lengthlist[i]=num+1;
 317     }
 318
 319     break;
 320   case 1:
 321     /* ordered */
 322     {
 323       long length=_oggpack_read(opb,5)+1;
 324       s->lengthlist=malloc(sizeof(long)*s->entries);
 325
 326       for(i=0;i<s->entries;){
 327         long num=_oggpack_read(opb,ilog(s->entries-i));
 328         if(num==-1)goto _eofout;
 329         for(j=0;j<num;j++,i++)
 330           s->lengthlist[i]=length;
 331         length++;
 332       }
 333     }
 334     break;
 335   default:
 336     /* EOF */
 337     return(-1);
 338   }
 339
 340   /* Do we have a mapping to unpack? */
 341   if(_oggpack_read(opb,1)){
 342
 343     /* values that define the dequantization */
 344     s->q_min=_oggpack_read(opb,24);
 345     s->q_delta=_oggpack_read(opb,24);
 346     s->q_quant=_oggpack_read(opb,4)+1;
 347     s->q_sequencep=_oggpack_read(opb,1);
 348
 349     /* quantized values */
 350     s->quantlist=malloc(sizeof(double)*s->entries*s->dim);
 351     for(i=0;i<s->entries*s->dim;i++)
 352       s->quantlist[i]=_oggpack_read(opb,s->q_quant);
 353     if(s->quantlist[i-1]==-1)goto _eofout;
 354   }
 355
 356   /* all set */
 357   return(0);
 358
 359  _errout:
 360  _eofout:
 361   vorbis_staticbook_clear(s);
 362   return(-1);
 363 }
 364
 365 /* returns the number of bits ***************************************/
 366 int vorbis_book_encode(codebook *book, int a, oggpack_buffer *b){
 367   _oggpack_write(b,book->codelist[a],book->c->lengthlist[a]);
 368   return(book->c->lengthlist[a]);
 369 }
 370
 371 static int  _best(codebook *book, double *a){
 372   encode_aux *t=book->c->encode_tree;
 373   int dim=book->dim;
 374   int ptr=0,k;
 375   /* optimized, using the decision tree */
 376   while(1){
 377     double c=0.;
 378     double *p=book->valuelist+t->p[ptr];
 379     double *q=book->valuelist+t->q[ptr];
 380
 381     for(k=0;k<dim;k++)
 382       c+=(p[k]-q[k])*(a[k]-(p[k]+q[k])*.5);
 383
 384     if(c>0.) /* in A */
 385       ptr= -t->ptr0[ptr];
 386     else     /* in B */
 387       ptr= -t->ptr1[ptr];
 388     if(ptr<=0)break;
 389   }
 390   return(-ptr);
 391 }
 392
 393 /* returns the number of bits and *modifies a* to the quantization value *****/
 394 int vorbis_book_encodev(codebook *book, double *a, oggpack_buffer *b){
 395   int dim=book->dim;
 396   int best=_best(book,a);
 397   memcpy(a,book->valuelist+best*dim,dim*sizeof(double));
 398   return(vorbis_book_encode(book,best,b));}
 399
 400 /* returns the number of bits and *modifies a* to the quantization error *****/
 401 int vorbis_book_encodevE(codebook *book, double *a, oggpack_buffer *b){
 402   int dim=book->dim,k;
 403   int best=_best(book,a);
 404   for(k=0;k<dim;k++)
 405     a[k]-=(book->valuelist+best*dim)[k];
 406   return(vorbis_book_encode(book,best,b));
 407 }
 408
 409 /* returns the total squared quantization error for best match and sets each
 410    element of a to local error ***************/
 411 double vorbis_book_vE(codebook *book, double *a){
 412   int dim=book->dim,k;
 413   int best=_best(book,a);
 414   double acc=0.;
 415   for(k=0;k<dim;k++){
 416     double val=(book->valuelist+best*dim)[k];
 417     a[k]-=val;
 418     acc+=a[k]*a[k];
 419   }
 420   return(acc);
 421 }
 422
 423 /* returns the entry number or -1 on eof *************************************/
 424 long vorbis_book_decode(codebook *book, oggpack_buffer *b){
 425   long ptr=0;
 426   decode_aux *t=book->decode_tree;
 427   do{
 428     switch(_oggpack_read1(b)){
 429     case 0:
 430       ptr=t->ptr0[ptr];
 431       break;
 432     case 1:
 433       ptr=t->ptr1[ptr];
 434       break;
 435     case -1:
 436       return(-1);
 437     }
 438   }while(ptr>0);
 439   return(-ptr);
 440 }
 441
 442 /* returns the entry number or -1 on eof *************************************/
 443 long vorbis_book_decodev(codebook *book, double *a, oggpack_buffer *b){
 444   long entry=vorbis_book_decode(book,b);
 445   int i;
 446   if(entry==-1)return(-1);
 447   for(i=0;i<book->dim;i++)a[i]+=(book->valuelist+entry*book->dim)[i];
 448   return(entry);
 449 }
 450
 451 #ifdef _V_SELFTEST
 452
 453 /* Simple enough; pack a few candidate codebooks, unpack them.  Code a
 454    number of vectors through (keeping track of the quantized values),
 455    and decode using the unpacked book.  quantized version of in should
 456    exactly equal out */
 457
 458 #include <stdio.h>
 459 #include "vorbis/book/lsp20_0.vqh"
 460 #include "vorbis/book/lsp32_0.vqh"
 461 #define TESTSIZE 40
 462 #define TESTDIM 4
 463
 464 double test1[40]={
 465   0.105939,
 466   0.215373,
 467   0.429117,
 468   0.587974,
 469
 470   0.181173,
 471   0.296583,
 472   0.515707,
 473   0.715261,
 474
 475   0.162327,
 476   0.263834,
 477   0.342876,
 478   0.406025,
 479
 480   0.103571,
 481   0.223561,
 482   0.368513,
 483   0.540313,
 484
 485   0.136672,
 486   0.395882,
 487   0.587183,
 488   0.652476,
 489
 490   0.114338,
 491   0.417300,
 492   0.525486,
 493   0.698679,
 494
 495   0.147492,
 496   0.324481,
 497   0.643089,
 498   0.757582,
 499
 500   0.139556,
 501   0.215795,
 502   0.324559,
 503   0.399387,
 504
 505   0.120236,
 506   0.267420,
 507   0.446940,
 508   0.608760,
 509
 510   0.115587,
 511   0.287234,
 512   0.571081,
 513   0.708603,
 514 };
 515
 516 double test2[40]={
 517   0.088654,
 518   0.165742,
 519   0.279013,
 520   0.395894,
 521
 522   0.110812,
 523   0.218422,
 524   0.283423,
 525   0.371719,
 526
 527   0.136985,
 528   0.186066,
 529   0.309814,
 530   0.381521,
 531
 532   0.123925,
 533   0.211707,
 534   0.314771,
 535   0.433026,
 536
 537   0.088619,
 538   0.192276,
 539   0.277568,
 540   0.343509,
 541
 542   0.068400,
 543   0.132901,
 544   0.223999,
 545   0.302538,
 546
 547   0.202159,
 548   0.306131,
 549   0.360362,
 550   0.416066,
 551
 552   0.072591,
 553   0.178019,
 554   0.304315,
 555   0.376516,
 556
 557   0.094336,
 558   0.188401,
 559   0.325119,
 560   0.390264,
 561
 562   0.091636,
 563   0.223099,
 564   0.282899,
 565   0.375124,
 566 };
 567
 568 static_codebook *testlist[]={&_vq_book_lsp20_0,&_vq_book_lsp32_0,NULL};
 569 double   *testvec[]={test1,test2};
 570
 571 int main(){
 572   oggpack_buffer write;
 573   oggpack_buffer read;
 574   long ptr=0,i;
 575   _oggpack_writeinit(&write);
 576
 577   fprintf(stderr,"Testing codebook abstraction...:\n");
 578
 579   while(testlist[ptr]){
 580     codebook c;
 581     static_codebook s;
 582     double *qv=alloca(sizeof(double)*TESTSIZE);
 583     double *iv=alloca(sizeof(double)*TESTSIZE);
 584     memcpy(qv,testvec[ptr],sizeof(double)*TESTSIZE);
 585     memset(iv,0,sizeof(double)*TESTSIZE);
 586
 587     fprintf(stderr,"\tpacking/coding %ld... ",ptr);
 588
 589     /* pack the codebook, write the testvector */
 590     _oggpack_reset(&write);
 591     vorbis_book_init_encode(&c,testlist[ptr]); /* get it into memory
 592                                                   we can write */
 593     vorbis_staticbook_pack(testlist[ptr],&write);
 594     fprintf(stderr,"Codebook size %ld bytes... ",_oggpack_bytes(&write));
 595     for(i=0;i<TESTSIZE;i+=TESTDIM)
 596       vorbis_book_encodev(&c,qv+i,&write);
 597     vorbis_book_clear(&c);
 598
 599     fprintf(stderr,"OK.\n");
 600     fprintf(stderr,"\tunpacking/decoding %ld... ",ptr);
 601
 602     /* transfer the write data to a read buffer and unpack/read */
 603     _oggpack_readinit(&read,_oggpack_buffer(&write),_oggpack_bytes(&write));
 604     if(vorbis_staticbook_unpack(&read,&s)){
 605       fprintf(stderr,"Error unpacking codebook.\n");
 606       exit(1);
 607     }
 608     if(vorbis_book_init_decode(&c,&s)){
 609       fprintf(stderr,"Error initializing codebook.\n");
 610       exit(1);
 611     }
 612
 613     for(i=0;i<TESTSIZE;i+=TESTDIM)
 614       if(vorbis_book_decodev(&c,iv+i,&read)==-1){
 615         fprintf(stderr,"Error reading codebook test data (EOP).\n");
 616         exit(1);
 617       }
 618     for(i=0;i<TESTSIZE;i++)
 619       if(fabs(qv[i]-iv[i])>.000001){
 620         fprintf(stderr,"input (%g) != output (%g) at position (%ld)\n",
 621                 iv[i],testvec[ptr][i]-qv[i],i);
 622         exit(1);
 623       }
 624
 625     fprintf(stderr,"OK\n");
 626     ptr++;
 627   }
 628   exit(0);
 629 }
 630
 631 #endif