neo/sound/OggVorbis/vorbissrc/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.29 2002/10/11 07:44: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 "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   /* two different remappings go on here.
 333
 334      First, we collapse the likely sparse codebook down only to
 335      actually represented values/words.  This collapsing needs to be
 336      indexed as map-valueless books are used to encode original entry
 337      positions as integers.
 338
 339      Second, we reorder all vectors, including the entry index above,
 340      by sorted bitreversed codeword to allow treeless decode. */
 341
 342   {
 343     /* perform sort */
 344     ogg_uint32_t *codes=_make_words(s->lengthlist,s->entries,c->used_entries);
 345     ogg_uint32_t **codep=alloca(sizeof(*codep)*n);
 346
 347     if(codes==NULL)goto err_out;
 348
 349     for(i=0;i<n;i++){
 350       codes[i]=bitreverse(codes[i]);
 351       codep[i]=codes+i;
 352     }
 353
 354     qsort(codep,n,sizeof(*codep),sort32a);
 355
 356     sortindex=alloca(n*sizeof(*sortindex));
 357     c->codelist=_ogg_malloc(n*sizeof(*c->codelist));
 358     /* the index is a reverse index */
 359     for(i=0;i<n;i++){
 360       int position=codep[i]-codes;
 361       sortindex[position]=i;
 362     }
 363
 364     for(i=0;i<n;i++)
 365       c->codelist[sortindex[i]]=codes[i];
 366     _ogg_free(codes);
 367   }
 368
 369   c->valuelist=_book_unquantize(s,n,sortindex);
 370   c->dec_index=_ogg_malloc(n*sizeof(*c->dec_index));
 371
 372   for(n=0,i=0;i<s->entries;i++)
 373     if(s->lengthlist[i]>0)
 374       c->dec_index[sortindex[n++]]=i;
 375
 376   c->dec_codelengths=_ogg_malloc(n*sizeof(*c->dec_codelengths));
 377   for(n=0,i=0;i<s->entries;i++)
 378     if(s->lengthlist[i]>0)
 379       c->dec_codelengths[sortindex[n++]]=s->lengthlist[i];
 380
 381   c->dec_firsttablen=_ilog(c->used_entries)-4; /* this is magic */
 382   if(c->dec_firsttablen<5)c->dec_firsttablen=5;
 383   if(c->dec_firsttablen>8)c->dec_firsttablen=8;
 384
 385   tabn=1<<c->dec_firsttablen;
 386   c->dec_firsttable=_ogg_calloc(tabn,sizeof(*c->dec_firsttable));
 387   c->dec_maxlength=0;
 388
 389   for(i=0;i<n;i++){
 390     if(c->dec_maxlength<c->dec_codelengths[i])
 391       c->dec_maxlength=c->dec_codelengths[i];
 392     if(c->dec_codelengths[i]<=c->dec_firsttablen){
 393       ogg_uint32_t orig=bitreverse(c->codelist[i]);
 394       for(j=0;j<(1<<(c->dec_firsttablen-c->dec_codelengths[i]));j++)
 395         c->dec_firsttable[orig|(j<<c->dec_codelengths[i])]=i+1;
 396     }
 397   }
 398
 399   /* now fill in 'unused' entries in the firsttable with hi/lo search
 400      hints for the non-direct-hits */
 401   {
 402     ogg_uint32_t mask=0xfffffffeUL<<(31-c->dec_firsttablen);
 403     long lo=0,hi=0;
 404
 405     for(i=0;i<tabn;i++){
 406       ogg_uint32_t word=i<<(32-c->dec_firsttablen);
 407       if(c->dec_firsttable[bitreverse(word)]==0){
 408         while((lo+1)<n && c->codelist[lo+1]<=word)lo++;
 409         while(    hi<n && word>=(c->codelist[hi]&mask))hi++;
 410
 411         /* we only actually have 15 bits per hint to play with here.
 412            In order to overflow gracefully (nothing breaks, efficiency
 413            just drops), encode as the difference from the extremes. */
 414         {
 415           unsigned long loval=lo;
 416           unsigned long hival=n-hi;
 417
 418           if(loval>0x7fff)loval=0x7fff;
 419           if(hival>0x7fff)hival=0x7fff;
 420           c->dec_firsttable[bitreverse(word)]=
 421             0x80000000UL | (loval<<15) | hival;
 422         }
 423       }
 424     }
 425   }
 426
 427
 428   return(0);
 429  err_out:
 430   vorbis_book_clear(c);
 431   return(-1);
 432 }
 433
 434 static float _dist(int el,float *ref, float *b,int step){
 435   int i;
 436   float acc=0.f;
 437   for(i=0;i<el;i++){
 438     float val=(ref[i]-b[i*step]);
 439     acc+=val*val;
 440   }
 441   return(acc);
 442 }
 443
 444 int _best(codebook *book, float *a, int step){
 445   encode_aux_threshmatch *tt=book->c->thresh_tree;
 446
 447 #if 0
 448   encode_aux_nearestmatch *nt=book->c->nearest_tree;
 449   encode_aux_pigeonhole *pt=book->c->pigeon_tree;
 450 #endif
 451   int dim=book->dim;
 452   int k,o;
 453   /*int savebest=-1;
 454     float saverr;*/
 455
 456   /* do we have a threshhold encode hint? */
 457   if(tt){
 458     int index=0,i;
 459     /* find the quant val of each scalar */
 460     for(k=0,o=step*(dim-1);k<dim;k++,o-=step){
 461
 462       i=tt->threshvals>>1;
 463       if(a[o]<tt->quantthresh[i]){
 464
 465         for(;i>0;i--)
 466           if(a[o]>=tt->quantthresh[i-1])
 467             break;
 468
 469       }else{
 470
 471         for(i++;i<tt->threshvals-1;i++)
 472           if(a[o]<tt->quantthresh[i])break;
 473
 474       }
 475
 476       index=(index*tt->quantvals)+tt->quantmap[i];
 477     }
 478     /* regular lattices are easy :-) */
 479     if(book->c->lengthlist[index]>0) /* is this unused?  If so, we'll
 480                                         use a decision tree after all
 481                                         and fall through*/
 482       return(index);
 483   }
 484
 485 #if 0
 486   /* do we have a pigeonhole encode hint? */
 487   if(pt){
 488     const static_codebook *c=book->c;
 489     int i,besti=-1;
 490     float best=0.f;
 491     int entry=0;
 492
 493     /* dealing with sequentialness is a pain in the ass */
 494     if(c->q_sequencep){
 495       int pv;
 496       long mul=1;
 497       float qlast=0;
 498       for(k=0,o=0;k<dim;k++,o+=step){
 499         pv=(int)((a[o]-qlast-pt->min)/pt->del);
 500         if(pv<0 || pv>=pt->mapentries)break;
 501         entry+=pt->pigeonmap[pv]*mul;
 502         mul*=pt->quantvals;
 503         qlast+=pv*pt->del+pt->min;
 504       }
 505     }else{
 506       for(k=0,o=step*(dim-1);k<dim;k++,o-=step){
 507         int pv=(int)((a[o]-pt->min)/pt->del);
 508         if(pv<0 || pv>=pt->mapentries)break;
 509         entry=entry*pt->quantvals+pt->pigeonmap[pv];
 510       }
 511     }
 512
 513     /* must be within the pigeonholable range; if we quant outside (or
 514        in an entry that we define no list for), brute force it */
 515     if(k==dim && pt->fitlength[entry]){
 516       /* search the abbreviated list */
 517       long *list=pt->fitlist+pt->fitmap[entry];
 518       for(i=0;i<pt->fitlength[entry];i++){
 519         float this=_dist(dim,book->valuelist+list[i]*dim,a,step);
 520         if(besti==-1 || this<best){
 521           best=this;
 522           besti=list[i];
 523         }
 524       }
 525
 526       return(besti);
 527     }
 528   }
 529
 530   if(nt){
 531     /* optimized using the decision tree */
 532     while(1){
 533       float c=0.f;
 534       float *p=book->valuelist+nt->p[ptr];
 535       float *q=book->valuelist+nt->q[ptr];
 536
 537       for(k=0,o=0;k<dim;k++,o+=step)
 538         c+=(p[k]-q[k])*(a[o]-(p[k]+q[k])*.5);
 539
 540       if(c>0.f) /* in A */
 541         ptr= -nt->ptr0[ptr];
 542       else     /* in B */
 543         ptr= -nt->ptr1[ptr];
 544       if(ptr<=0)break;
 545     }
 546     return(-ptr);
 547   }
 548 #endif
 549
 550   /* brute force it! */
 551   {
 552     const static_codebook *c=book->c;
 553     int i,besti=-1;
 554     float best=0.f;
 555     float *e=book->valuelist;
 556     for(i=0;i<book->entries;i++){
 557       if(c->lengthlist[i]>0){
 558         float this=_dist(dim,e,a,step);
 559         if(besti==-1 || this<best){
 560           best=this;
 561           besti=i;
 562         }
 563       }
 564       e+=dim;
 565     }
 566
 567     /*if(savebest!=-1 && savebest!=besti){
 568       fprintf(stderr,"brute force/pigeonhole disagreement:\n"
 569               "original:");
 570       for(i=0;i<dim*step;i+=step)fprintf(stderr,"%g,",a[i]);
 571       fprintf(stderr,"\n"
 572               "pigeonhole (entry %d, err %g):",savebest,saverr);
 573       for(i=0;i<dim;i++)fprintf(stderr,"%g,",
 574                                 (book->valuelist+savebest*dim)[i]);
 575       fprintf(stderr,"\n"
 576               "bruteforce (entry %d, err %g):",besti,best);
 577       for(i=0;i<dim;i++)fprintf(stderr,"%g,",
 578                                 (book->valuelist+besti*dim)[i]);
 579       fprintf(stderr,"\n");
 580       }*/
 581     return(besti);
 582   }
 583 }
 584
 585 long vorbis_book_codeword(codebook *book,int entry){
 586   if(book->c) /* only use with encode; decode optimizations are
 587                  allowed to break this */
 588     return book->codelist[entry];
 589   return -1;
 590 }
 591
 592 long vorbis_book_codelen(codebook *book,int entry){
 593   if(book->c) /* only use with encode; decode optimizations are
 594                  allowed to break this */
 595     return book->c->lengthlist[entry];
 596   return -1;
 597 }
 598
 599 #ifdef _V_SELFTEST
 600
 601 /* Unit tests of the dequantizer; this stuff will be OK
 602    cross-platform, I simply want to be sure that special mapping cases
 603    actually work properly; a bug could go unnoticed for a while */
 604
 605 #include <stdio.h>
 606
 607 /* cases:
 608
 609    no mapping
 610    full, explicit mapping
 611    algorithmic mapping
 612
 613    nonsequential
 614    sequential
 615 */
 616
 617 static long full_quantlist1[]={0,1,2,3,    4,5,6,7, 8,3,6,1};
 618 static long partial_quantlist1[]={0,7,2};
 619
 620 /* no mapping */
 621 static_codebook test1={
 622   4,16,
 623   NULL,
 624   0,
 625   0,0,0,0,
 626   NULL,
 627   NULL,NULL
 628 };
 629 static float *test1_result=NULL;
 630
 631 /* linear, full mapping, nonsequential */
 632 static_codebook test2={
 633   4,3,
 634   NULL,
 635   2,
 636   -533200896,1611661312,4,0,
 637   full_quantlist1,
 638   NULL,NULL
 639 };
 640 static float test2_result[]={-3,-2,-1,0, 1,2,3,4, 5,0,3,-2};
 641
 642 /* linear, full mapping, sequential */
 643 static_codebook test3={
 644   4,3,
 645   NULL,
 646   2,
 647   -533200896,1611661312,4,1,
 648   full_quantlist1,
 649   NULL,NULL
 650 };
 651 static float test3_result[]={-3,-5,-6,-6, 1,3,6,10, 5,5,8,6};
 652
 653 /* linear, algorithmic mapping, nonsequential */
 654 static_codebook test4={
 655   3,27,
 656   NULL,
 657   1,
 658   -533200896,1611661312,4,0,
 659   partial_quantlist1,
 660   NULL,NULL
 661 };
 662 static float test4_result[]={-3,-3,-3, 4,-3,-3, -1,-3,-3,
 663                               -3, 4,-3, 4, 4,-3, -1, 4,-3,
 664                               -3,-1,-3, 4,-1,-3, -1,-1,-3,
 665                               -3,-3, 4, 4,-3, 4, -1,-3, 4,
 666                               -3, 4, 4, 4, 4, 4, -1, 4, 4,
 667                               -3,-1, 4, 4,-1, 4, -1,-1, 4,
 668                               -3,-3,-1, 4,-3,-1, -1,-3,-1,
 669                               -3, 4,-1, 4, 4,-1, -1, 4,-1,
 670                               -3,-1,-1, 4,-1,-1, -1,-1,-1};
 671
 672 /* linear, algorithmic mapping, sequential */
 673 static_codebook test5={
 674   3,27,
 675   NULL,
 676   1,
 677   -533200896,1611661312,4,1,
 678   partial_quantlist1,
 679   NULL,NULL
 680 };
 681 static float test5_result[]={-3,-6,-9, 4, 1,-2, -1,-4,-7,
 682                               -3, 1,-2, 4, 8, 5, -1, 3, 0,
 683                               -3,-4,-7, 4, 3, 0, -1,-2,-5,
 684                               -3,-6,-2, 4, 1, 5, -1,-4, 0,
 685                               -3, 1, 5, 4, 8,12, -1, 3, 7,
 686                               -3,-4, 0, 4, 3, 7, -1,-2, 2,
 687                               -3,-6,-7, 4, 1, 0, -1,-4,-5,
 688                               -3, 1, 0, 4, 8, 7, -1, 3, 2,
 689                               -3,-4,-5, 4, 3, 2, -1,-2,-3};
 690
 691 void run_test(static_codebook *b,float *comp){
 692   float *out=_book_unquantize(b,b->entries,NULL);
 693   int i;
 694
 695   if(comp){
 696     if(!out){
 697       fprintf(stderr,"_book_unquantize incorrectly returned NULL\n");
 698       exit(1);
 699     }
 700
 701     for(i=0;i<b->entries*b->dim;i++)
 702       if(fabs(out[i]-comp[i])>.0001){
 703         fprintf(stderr,"disagreement in unquantized and reference data:\n"
 704                 "position %d, %g != %g\n",i,out[i],comp[i]);
 705         exit(1);
 706       }
 707
 708   }else{
 709     if(out){
 710       fprintf(stderr,"_book_unquantize returned a value array: \n"
 711               " correct result should have been NULL\n");
 712       exit(1);
 713     }
 714   }
 715 }
 716
 717 int main(){
 718   /* run the nine dequant tests, and compare to the hand-rolled results */
 719   fprintf(stderr,"Dequant test 1... ");
 720   run_test(&test1,test1_result);
 721   fprintf(stderr,"OK\nDequant test 2... ");
 722   run_test(&test2,test2_result);
 723   fprintf(stderr,"OK\nDequant test 3... ");
 724   run_test(&test3,test3_result);
 725   fprintf(stderr,"OK\nDequant test 4... ");
 726   run_test(&test4,test4_result);
 727   fprintf(stderr,"OK\nDequant test 5... ");
 728   run_test(&test5,test5_result);
 729   fprintf(stderr,"OK\n\n");
 730
 731   return(0);
 732 }
 733
 734 #endif