lib/sharedbook.c

   1 /********************************************************************
   2  *                                                                  *
   3  * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE.   *
   4  * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
   5  * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
   6  * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
   7  *                                                                  *
   8  * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2002             *
   9  * by the XIPHOPHORUS Company http://www.xiph.org/                  *
  10  *                                                                  *
  11  ********************************************************************
  12
  13  function: basic shared codebook operations
  14  last mod: $Id: sharedbook.c,v 1.25 2002/01/21 20:51:28 xiphmont Exp $
  15
  16  ********************************************************************/
  17
  18 #include <stdlib.h>
  19 #include <math.h>
  20 #include <string.h>
  21 #include <ogg/ogg.h>
  22 #include "os.h"
  23 #include "vorbis/codec.h"
  24 #include "codebook.h"
  25 #include "scales.h"
  26
  27 /**** pack/unpack helpers ******************************************/
  28 int _ilog(unsigned int v){
  29   int ret=0;
  30   while(v){
  31     ret++;
  32     v>>=1;
  33   }
  34   return(ret);
  35 }
  36
  37 /* 32 bit float (not IEEE; nonnormalized mantissa +
  38    biased exponent) : neeeeeee eeemmmmm mmmmmmmm mmmmmmmm
  39    Why not IEEE?  It's just not that important here. */
  40
  41 #define VQ_FEXP 10
  42 #define VQ_FMAN 21
  43 #define VQ_FEXP_BIAS 768 /* bias toward values smaller than 1. */
  44
  45 /* doesn't currently guard under/overflow */
  46 long _float32_pack(float val){
  47   int sign=0;
  48   long exp;
  49   long mant;
  50   if(val<0){
  51     sign=0x80000000;
  52     val= -val;
  53   }
  54   exp= floor(log(val)/log(2.f));
  55   mant=rint(ldexp(val,(VQ_FMAN-1)-exp));
  56   exp=(exp+VQ_FEXP_BIAS)<<VQ_FMAN;
  57
  58   return(sign|exp|mant);
  59 }
  60
  61 float _float32_unpack(long val){
  62   double mant=val&0x1fffff;
  63   int    sign=val&0x80000000;
  64   long   exp =(val&0x7fe00000L)>>VQ_FMAN;
  65   if(sign)mant= -mant;
  66   return(ldexp(mant,exp-(VQ_FMAN-1)-VQ_FEXP_BIAS));
  67 }
  68
  69 /* given a list of word lengths, generate a list of codewords.  Works
  70    for length ordered or unordered, always assigns the lowest valued
  71    codewords first.  Extended to handle unused entries (length 0) */
  72 ogg_uint32_t *_make_words(long *l,long n,long sparsecount){
  73   long i,j,count=0;
  74   ogg_uint32_t marker[33];
  75   ogg_uint32_t *r=_ogg_malloc((sparsecount?sparsecount:n)*sizeof(*r));
  76   memset(marker,0,sizeof(marker));
  77
  78   for(i=0;i<n;i++){
  79     long length=l[i];
  80     if(length>0){
  81       ogg_uint32_t entry=marker[length];
  82
  83       /* when we claim a node for an entry, we also claim the nodes
  84          below it (pruning off the imagined tree that may have dangled
  85          from it) as well as blocking the use of any nodes directly
  86          above for leaves */
  87
  88       /* update ourself */
  89       if(length<32 && (entry>>length)){
  90         /* error condition; the lengths must specify an overpopulated tree */
  91         _ogg_free(r);
  92         return(NULL);
  93       }
  94       r[count++]=entry;
  95
  96       /* Look to see if the next shorter marker points to the node
  97          above. if so, update it and repeat.  */
  98       {
  99         for(j=length;j>0;j--){
 100
 101           if(marker[j]&1){
 102             /* have to jump branches */
 103             if(j==1)
 104               marker[1]++;
 105             else
 106               marker[j]=marker[j-1]<<1;
 107             break; /* invariant says next upper marker would already
 108                       have been moved if it was on the same path */
 109           }
 110           marker[j]++;
 111         }
 112       }
 113
 114       /* prune the tree; the implicit invariant says all the longer
 115          markers were dangling from our just-taken node.  Dangle them
 116          from our *new* node. */
 117       for(j=length+1;j<33;j++)
 118         if((marker[j]>>1) == entry){
 119           entry=marker[j];
 120           marker[j]=marker[j-1]<<1;
 121         }else
 122           break;
 123     }else
 124       if(sparsecount==0)count++;
 125   }
 126
 127   /* bitreverse the words because our bitwise packer/unpacker is LSb
 128      endian */
 129   for(i=0,count=0;i<n;i++){
 130     ogg_uint32_t temp=0;
 131     for(j=0;j<l[i];j++){
 132       temp<<=1;
 133       temp|=(r[count]>>j)&1;
 134     }
 135
 136     if(sparsecount){
 137       if(l[i])
 138         r[count++]=temp;
 139     }else
 140       r[count++]=temp;
 141   }
 142
 143   return(r);
 144 }
 145
 146 /* there might be a straightforward one-line way to do the below
 147    that's portable and totally safe against roundoff, but I haven't
 148    thought of it.  Therefore, we opt on the side of caution */
 149 long _book_maptype1_quantvals(const static_codebook *b){
 150   long vals=floor(pow((float)b->entries,1.f/b->dim));
 151
 152   /* the above *should* be reliable, but we'll not assume that FP is
 153      ever reliable when bitstream sync is at stake; verify via integer
 154      means that vals really is the greatest value of dim for which
 155      vals^b->bim <= b->entries */
 156   /* treat the above as an initial guess */
 157   while(1){
 158     long acc=1;
 159     long acc1=1;
 160     int i;
 161     for(i=0;i<b->dim;i++){
 162       acc*=vals;
 163       acc1*=vals+1;
 164     }
 165     if(acc<=b->entries && acc1>b->entries){
 166       return(vals);
 167     }else{
 168       if(acc>b->entries){
 169         vals--;
 170       }else{
 171         vals++;
 172       }
 173     }
 174   }
 175 }
 176
 177 /* unpack the quantized list of values for encode/decode ***********/
 178 /* we need to deal with two map types: in map type 1, the values are
 179    generated algorithmically (each column of the vector counts through
 180    the values in the quant vector). in map type 2, all the values came
 181    in in an explicit list.  Both value lists must be unpacked */
 182 float *_book_unquantize(const static_codebook *b,int n,int *sparsemap){
 183   long j,k,count=0;
 184   if(b->maptype==1 || b->maptype==2){
 185     int quantvals;
 186     float mindel=_float32_unpack(b->q_min);
 187     float delta=_float32_unpack(b->q_delta);
 188     float *r=_ogg_calloc(n*b->dim,sizeof(*r));
 189
 190     /* maptype 1 and 2 both use a quantized value vector, but
 191        different sizes */
 192     switch(b->maptype){
 193     case 1:
 194       /* most of the time, entries%dimensions == 0, but we need to be
 195          well defined.  We define that the possible vales at each
 196          scalar is values == entries/dim.  If entries%dim != 0, we'll
 197          have 'too few' values (values*dim<entries), which means that
 198          we'll have 'left over' entries; left over entries use zeroed
 199          values (and are wasted).  So don't generate codebooks like
 200          that */
 201       quantvals=_book_maptype1_quantvals(b);
 202       for(j=0;j<b->entries;j++){
 203         if((sparsemap && b->lengthlist[j]) || !sparsemap){
 204           float last=0.f;
 205           int indexdiv=1;
 206           for(k=0;k<b->dim;k++){
 207             int index= (j/indexdiv)%quantvals;
 208             float val=b->quantlist[index];
 209             val=fabs(val)*delta+mindel+last;
 210             if(b->q_sequencep)last=val;
 211             if(sparsemap)
 212               r[sparsemap[count]*b->dim+k]=val;
 213             else
 214               r[count*b->dim+k]=val;
 215             indexdiv*=quantvals;
 216           }
 217           count++;
 218         }
 219
 220       }
 221       break;
 222     case 2:
 223       for(j=0;j<b->entries;j++){
 224         if((sparsemap && b->lengthlist[j]) || !sparsemap){
 225           float last=0.f;
 226
 227           for(k=0;k<b->dim;k++){
 228             float val=b->quantlist[j*b->dim+k];
 229             val=fabs(val)*delta+mindel+last;
 230             if(b->q_sequencep)last=val;
 231             if(sparsemap)
 232               r[sparsemap[count]*b->dim+k]=val;
 233             else
 234               r[count*b->dim+k]=val;
 235           }
 236           count++;
 237         }
 238       }
 239       break;
 240     }
 241
 242     return(r);
 243   }
 244   return(NULL);
 245 }
 246
 247 void vorbis_staticbook_clear(static_codebook *b){
 248   if(b->allocedp){
 249     if(b->quantlist)_ogg_free(b->quantlist);
 250     if(b->lengthlist)_ogg_free(b->lengthlist);
 251     if(b->nearest_tree){
 252       _ogg_free(b->nearest_tree->ptr0);
 253       _ogg_free(b->nearest_tree->ptr1);
 254       _ogg_free(b->nearest_tree->p);
 255       _ogg_free(b->nearest_tree->q);
 256       memset(b->nearest_tree,0,sizeof(*b->nearest_tree));
 257       _ogg_free(b->nearest_tree);
 258     }
 259     if(b->thresh_tree){
 260       _ogg_free(b->thresh_tree->quantthresh);
 261       _ogg_free(b->thresh_tree->quantmap);
 262       memset(b->thresh_tree,0,sizeof(*b->thresh_tree));
 263       _ogg_free(b->thresh_tree);
 264     }
 265
 266     memset(b,0,sizeof(*b));
 267   }
 268 }
 269
 270 void vorbis_staticbook_destroy(static_codebook *b){
 271   if(b->allocedp){
 272     vorbis_staticbook_clear(b);
 273     _ogg_free(b);
 274   }
 275 }
 276
 277 void vorbis_book_clear(codebook *b){
 278   /* static book is not cleared; we're likely called on the lookup and
 279      the static codebook belongs to the info struct */
 280   if(b->valuelist)_ogg_free(b->valuelist);
 281   if(b->codelist)_ogg_free(b->codelist);
 282
 283   if(b->dec_index)_ogg_free(b->dec_index);
 284   if(b->dec_codelengths)_ogg_free(b->dec_codelengths);
 285   if(b->dec_firsttable)_ogg_free(b->dec_firsttable);
 286
 287   memset(b,0,sizeof(*b));
 288 }
 289
 290 int vorbis_book_init_encode(codebook *c,const static_codebook *s){
 291
 292   memset(c,0,sizeof(*c));
 293   c->c=s;
 294   c->entries=s->entries;
 295   c->used_entries=s->entries;
 296   c->dim=s->dim;
 297   c->codelist=_make_words(s->lengthlist,s->entries,0);
 298   c->valuelist=_book_unquantize(s,s->entries,NULL);
 299
 300   return(0);
 301 }
 302
 303 static ogg_uint32_t bitreverse(ogg_uint32_t x){
 304   x=    ((x>>16)&0x0000ffffUL) | ((x<<16)&0xffff0000UL);
 305   x=    ((x>> 8)&0x00ff00ffUL) | ((x<< 8)&0xff00ff00UL);
 306   x=    ((x>> 4)&0x0f0f0f0fUL) | ((x<< 4)&0xf0f0f0f0UL);
 307   x=    ((x>> 2)&0x33333333UL) | ((x<< 2)&0xccccccccUL);
 308   return((x>> 1)&0x55555555UL) | ((x<< 1)&0xaaaaaaaaUL);
 309 }
 310
 311 static int sort32a(const void *a,const void *b){
 312   return ( (**(ogg_uint32_t **)a>**(ogg_uint32_t **)b)<<1)-1;
 313 }
 314
 315 /* decode codebook arrangement is more heavily optimized than encode */
 316 int vorbis_book_init_decode(codebook *c,const static_codebook *s){
 317   int i,j,n=0,tabn;
 318   int *sortindex;
 319   memset(c,0,sizeof(*c));
 320
 321   /* count actually used entries */
 322   for(i=0;i<s->entries;i++)
 323     if(s->lengthlist[i]>0)
 324       n++;
 325
 326   c->entries=s->entries;
 327   c->used_entries=n;
 328   c->dim=s->dim;
 329
 330   /* two different remappings go on here.
 331
 332      First, we collapse the likely sparse codebook down only to
 333      actually represented values/words.  This collapsing needs to be
 334      indexed as map-valueless books are used to encode original entry
 335      positions as integers.
 336
 337      Second, we reorder all vectors, including the entry index above,
 338      by sorted bitreversed codeword to allow treeless decode. */
 339
 340   {
 341     /* perform sort */
 342     ogg_uint32_t *codes=_make_words(s->lengthlist,s->entries,c->used_entries);
 343     ogg_uint32_t **codep=alloca(sizeof(*codep)*n);
 344
 345     if(codes==NULL)goto err_out;
 346
 347     for(i=0;i<n;i++){
 348       codes[i]=bitreverse(codes[i]);
 349       codep[i]=codes+i;
 350     }
 351
 352     qsort(codep,n,sizeof(*codep),sort32a);
 353
 354     sortindex=alloca(n*sizeof(*sortindex));
 355     c->codelist=_ogg_malloc(n*sizeof(*c->codelist));
 356     /* the index is a reverse index */
 357     for(i=0;i<n;i++){
 358       int position=codep[i]-codes;
 359       sortindex[position]=i;
 360     }
 361
 362     for(i=0;i<n;i++)
 363       c->codelist[sortindex[i]]=codes[i];
 364     _ogg_free(codes);
 365   }
 366
 367   c->valuelist=_book_unquantize(s,n,sortindex);
 368   c->dec_index=_ogg_malloc(n*sizeof(*c->dec_index));
 369
 370   for(n=0,i=0;i<s->entries;i++)
 371     if(s->lengthlist[i]>0)
 372       c->dec_index[sortindex[n++]]=i;
 373
 374   c->dec_codelengths=_ogg_malloc(n*sizeof(*c->dec_codelengths));
 375   for(n=0,i=0;i<s->entries;i++)
 376     if(s->lengthlist[i]>0)
 377       c->dec_codelengths[sortindex[n++]]=s->lengthlist[i];
 378
 379   c->dec_firsttablen=_ilog(c->used_entries)-4; /* this is magic */
 380   if(c->dec_firsttablen<5)c->dec_firsttablen=5;
 381   if(c->dec_firsttablen>8)c->dec_firsttablen=8;
 382
 383   tabn=1<<c->dec_firsttablen;
 384   c->dec_firsttable=_ogg_calloc(tabn,sizeof(*c->dec_firsttable));
 385   c->dec_maxlength=0;
 386
 387   for(i=0;i<n;i++){
 388     if(c->dec_maxlength<c->dec_codelengths[i])
 389       c->dec_maxlength=c->dec_codelengths[i];
 390     if(c->dec_codelengths[i]<=c->dec_firsttablen){
 391       ogg_uint32_t orig=bitreverse(c->codelist[i]);
 392       for(j=0;j<(1<<(c->dec_firsttablen-c->dec_codelengths[i]));j++)
 393         c->dec_firsttable[orig|(j<<c->dec_codelengths[i])]=i+1;
 394     }
 395   }
 396
 397   /* now fill in 'unused' entries in the firsttable with hi/lo search
 398      hints for the non-direct-hits */
 399   {
 400     ogg_uint32_t mask=0xfffffffeUL<<(31-c->dec_firsttablen);
 401
 402     for(i=0;i<tabn;i++){
 403       if(c->dec_firsttable[i]==0){
 404         ogg_uint32_t testword=bitreverse(i);
 405         long lo=0,hi=0;
 406         while((lo+1)<n && c->codelist[lo+1]<=testword)lo++;
 407         while(    hi<n && testword>=(c->codelist[hi]&mask))hi++;
 408
 409         /* we only actually have 15 bits per hint to play with here.
 410            In order to overflow gracefully (nothing breaks, efficiency
 411            just drops), encode as the difference from the extremes. */
 412         {
 413           unsigned long loval=lo;
 414           unsigned long hival=n-hi;
 415
 416           if(loval>0x7fff)loval=0x7fff;
 417           if(hival>0x7fff)hival=0x7fff;
 418           c->dec_firsttable[i]=0x80000000UL | (loval<<15) | hival;
 419         }
 420       }
 421     }
 422   }
 423
 424
 425   return(0);
 426  err_out:
 427   vorbis_book_clear(c);
 428   return(-1);
 429 }
 430
 431 static float _dist(int el,float *ref, float *b,int step){
 432   int i;
 433   float acc=0.f;
 434   for(i=0;i<el;i++){
 435     float val=(ref[i]-b[i*step]);
 436     acc+=val*val;
 437   }
 438   return(acc);
 439 }
 440
 441 int _best(codebook *book, float *a, int step){
 442   encode_aux_nearestmatch *nt=book->c->nearest_tree;
 443   encode_aux_threshmatch *tt=book->c->thresh_tree;
 444   encode_aux_pigeonhole *pt=book->c->pigeon_tree;
 445   int dim=book->dim;
 446   int ptr=0,k,o;
 447   /*int savebest=-1;
 448     float saverr;*/
 449
 450   /* do we have a threshhold encode hint? */
 451   if(tt){
 452     int index=0;
 453     /* find the quant val of each scalar */
 454     for(k=0,o=step*(dim-1);k<dim;k++,o-=step){
 455       int i;
 456       /* linear search the quant list for now; it's small and although
 457          with > ~8 entries, it would be faster to bisect, this would be
 458          a misplaced optimization for now */
 459       for(i=0;i<tt->threshvals-1;i++)
 460         if(a[o]<tt->quantthresh[i])break;
 461
 462       index=(index*tt->quantvals)+tt->quantmap[i];
 463     }
 464     /* regular lattices are easy :-) */
 465     if(book->c->lengthlist[index]>0) /* is this unused?  If so, we'll
 466                                         use a decision tree after all
 467                                         and fall through*/
 468       return(index);
 469   }
 470
 471   /* do we have a pigeonhole encode hint? */
 472   if(pt){
 473     const static_codebook *c=book->c;
 474     int i,besti=-1;
 475     float best=0.f;
 476     int entry=0;
 477
 478     /* dealing with sequentialness is a pain in the ass */
 479     if(c->q_sequencep){
 480       int pv;
 481       long mul=1;
 482       float qlast=0;
 483       for(k=0,o=0;k<dim;k++,o+=step){
 484         pv=(int)((a[o]-qlast-pt->min)/pt->del);
 485         if(pv<0 || pv>=pt->mapentries)break;
 486         entry+=pt->pigeonmap[pv]*mul;
 487         mul*=pt->quantvals;
 488         qlast+=pv*pt->del+pt->min;
 489       }
 490     }else{
 491       for(k=0,o=step*(dim-1);k<dim;k++,o-=step){
 492         int pv=(int)((a[o]-pt->min)/pt->del);
 493         if(pv<0 || pv>=pt->mapentries)break;
 494         entry=entry*pt->quantvals+pt->pigeonmap[pv];
 495       }
 496     }
 497
 498     /* must be within the pigeonholable range; if we quant outside (or
 499        in an entry that we define no list for), brute force it */
 500     if(k==dim && pt->fitlength[entry]){
 501       /* search the abbreviated list */
 502       long *list=pt->fitlist+pt->fitmap[entry];
 503       for(i=0;i<pt->fitlength[entry];i++){
 504         float this=_dist(dim,book->valuelist+list[i]*dim,a,step);
 505         if(besti==-1 || this<best){
 506           best=this;
 507           besti=list[i];
 508         }
 509       }
 510
 511       return(besti);
 512     }
 513   }
 514
 515   if(nt){
 516     /* optimized using the decision tree */
 517     while(1){
 518       float c=0.f;
 519       float *p=book->valuelist+nt->p[ptr];
 520       float *q=book->valuelist+nt->q[ptr];
 521
 522       for(k=0,o=0;k<dim;k++,o+=step)
 523         c+=(p[k]-q[k])*(a[o]-(p[k]+q[k])*.5);
 524
 525       if(c>0.f) /* in A */
 526         ptr= -nt->ptr0[ptr];
 527       else     /* in B */
 528         ptr= -nt->ptr1[ptr];
 529       if(ptr<=0)break;
 530     }
 531     return(-ptr);
 532   }
 533
 534   /* brute force it! */
 535   {
 536     const static_codebook *c=book->c;
 537     int i,besti=-1;
 538     float best=0.f;
 539     float *e=book->valuelist;
 540     for(i=0;i<book->entries;i++){
 541       if(c->lengthlist[i]>0){
 542         float this=_dist(dim,e,a,step);
 543         if(besti==-1 || this<best){
 544           best=this;
 545           besti=i;
 546         }
 547       }
 548       e+=dim;
 549     }
 550
 551     /*if(savebest!=-1 && savebest!=besti){
 552       fprintf(stderr,"brute force/pigeonhole disagreement:\n"
 553               "original:");
 554       for(i=0;i<dim*step;i+=step)fprintf(stderr,"%g,",a[i]);
 555       fprintf(stderr,"\n"
 556               "pigeonhole (entry %d, err %g):",savebest,saverr);
 557       for(i=0;i<dim;i++)fprintf(stderr,"%g,",
 558                                 (book->valuelist+savebest*dim)[i]);
 559       fprintf(stderr,"\n"
 560               "bruteforce (entry %d, err %g):",besti,best);
 561       for(i=0;i<dim;i++)fprintf(stderr,"%g,",
 562                                 (book->valuelist+besti*dim)[i]);
 563       fprintf(stderr,"\n");
 564       }*/
 565     return(besti);
 566   }
 567 }
 568
 569 /* returns the entry number and *modifies a* to the remainder value ********/
 570 int vorbis_book_besterror(codebook *book,float *a,int step,int addmul){
 571   int dim=book->dim,i,o;
 572   int best=_best(book,a,step);
 573   switch(addmul){
 574   case 0:
 575     for(i=0,o=0;i<dim;i++,o+=step)
 576       a[o]-=(book->valuelist+best*dim)[i];
 577     break;
 578   case 1:
 579     for(i=0,o=0;i<dim;i++,o+=step){
 580       float val=(book->valuelist+best*dim)[i];
 581       if(val==0){
 582         a[o]=0;
 583       }else{
 584         a[o]/=val;
 585       }
 586     }
 587     break;
 588   }
 589   return(best);
 590 }
 591
 592 long vorbis_book_codeword(codebook *book,int entry){
 593   if(book->c) /* only use with encode; decode optimizations are
 594                  allowed to break this */
 595     return book->codelist[entry];
 596   return -1;
 597 }
 598
 599 long vorbis_book_codelen(codebook *book,int entry){
 600   if(book->c) /* only use with encode; decode optimizations are
 601                  allowed to break this */
 602     return book->c->lengthlist[entry];
 603   return -1;
 604 }
 605
 606 #ifdef _V_SELFTEST
 607
 608 /* Unit tests of the dequantizer; this stuff will be OK
 609    cross-platform, I simply want to be sure that special mapping cases
 610    actually work properly; a bug could go unnoticed for a while */
 611
 612 #include <stdio.h>
 613
 614 /* cases:
 615
 616    no mapping
 617    full, explicit mapping
 618    algorithmic mapping
 619
 620    nonsequential
 621    sequential
 622 */
 623
 624 static long full_quantlist1[]={0,1,2,3,    4,5,6,7, 8,3,6,1};
 625 static long partial_quantlist1[]={0,7,2};
 626
 627 /* no mapping */
 628 static_codebook test1={
 629   4,16,
 630   NULL,
 631   0,
 632   0,0,0,0,
 633   NULL,
 634   NULL,NULL
 635 };
 636 static float *test1_result=NULL;
 637
 638 /* linear, full mapping, nonsequential */
 639 static_codebook test2={
 640   4,3,
 641   NULL,
 642   2,
 643   -533200896,1611661312,4,0,
 644   full_quantlist1,
 645   NULL,NULL
 646 };
 647 static float test2_result[]={-3,-2,-1,0, 1,2,3,4, 5,0,3,-2};
 648
 649 /* linear, full mapping, sequential */
 650 static_codebook test3={
 651   4,3,
 652   NULL,
 653   2,
 654   -533200896,1611661312,4,1,
 655   full_quantlist1,
 656   NULL,NULL
 657 };
 658 static float test3_result[]={-3,-5,-6,-6, 1,3,6,10, 5,5,8,6};
 659
 660 /* linear, algorithmic mapping, nonsequential */
 661 static_codebook test4={
 662   3,27,
 663   NULL,
 664   1,
 665   -533200896,1611661312,4,0,
 666   partial_quantlist1,
 667   NULL,NULL
 668 };
 669 static float test4_result[]={-3,-3,-3, 4,-3,-3, -1,-3,-3,
 670                               -3, 4,-3, 4, 4,-3, -1, 4,-3,
 671                               -3,-1,-3, 4,-1,-3, -1,-1,-3,
 672                               -3,-3, 4, 4,-3, 4, -1,-3, 4,
 673                               -3, 4, 4, 4, 4, 4, -1, 4, 4,
 674                               -3,-1, 4, 4,-1, 4, -1,-1, 4,
 675                               -3,-3,-1, 4,-3,-1, -1,-3,-1,
 676                               -3, 4,-1, 4, 4,-1, -1, 4,-1,
 677                               -3,-1,-1, 4,-1,-1, -1,-1,-1};
 678
 679 /* linear, algorithmic mapping, sequential */
 680 static_codebook test5={
 681   3,27,
 682   NULL,
 683   1,
 684   -533200896,1611661312,4,1,
 685   partial_quantlist1,
 686   NULL,NULL
 687 };
 688 static float test5_result[]={-3,-6,-9, 4, 1,-2, -1,-4,-7,
 689                               -3, 1,-2, 4, 8, 5, -1, 3, 0,
 690                               -3,-4,-7, 4, 3, 0, -1,-2,-5,
 691                               -3,-6,-2, 4, 1, 5, -1,-4, 0,
 692                               -3, 1, 5, 4, 8,12, -1, 3, 7,
 693                               -3,-4, 0, 4, 3, 7, -1,-2, 2,
 694                               -3,-6,-7, 4, 1, 0, -1,-4,-5,
 695                               -3, 1, 0, 4, 8, 7, -1, 3, 2,
 696                               -3,-4,-5, 4, 3, 2, -1,-2,-3};
 697
 698 void run_test(static_codebook *b,float *comp){
 699   float *out=_book_unquantize(b,b->entries,NULL);
 700   int i;
 701
 702   if(comp){
 703     if(!out){
 704       fprintf(stderr,"_book_unquantize incorrectly returned NULL\n");
 705       exit(1);
 706     }
 707
 708     for(i=0;i<b->entries*b->dim;i++)
 709       if(fabs(out[i]-comp[i])>.0001){
 710         fprintf(stderr,"disagreement in unquantized and reference data:\n"
 711                 "position %d, %g != %g\n",i,out[i],comp[i]);
 712         exit(1);
 713       }
 714
 715   }else{
 716     if(out){
 717       fprintf(stderr,"_book_unquantize returned a value array: \n"
 718               " correct result should have been NULL\n");
 719       exit(1);
 720     }
 721   }
 722 }
 723
 724 int main(){
 725   /* run the nine dequant tests, and compare to the hand-rolled results */
 726   fprintf(stderr,"Dequant test 1... ");
 727   run_test(&test1,test1_result);
 728   fprintf(stderr,"OK\nDequant test 2... ");
 729   run_test(&test2,test2_result);
 730   fprintf(stderr,"OK\nDequant test 3... ");
 731   run_test(&test3,test3_result);
 732   fprintf(stderr,"OK\nDequant test 4... ");
 733   run_test(&test4,test4_result);
 734   fprintf(stderr,"OK\nDequant test 5... ");
 735   run_test(&test5,test5_result);
 736   fprintf(stderr,"OK\n\n");
 737
 738   return(0);
 739 }
 740
 741 #endif