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