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