tremor/sharedbook.c

   1 /********************************************************************
   2  *                                                                  *
   3  * THIS FILE IS PART OF THE OggVorbis 'TREMOR' CODEC SOURCE CODE.   *
   4  *                                                                  *
   5  * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
   6  * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
   7  * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
   8  *                                                                  *
   9  * THE OggVorbis 'TREMOR' SOURCE CODE IS (C) COPYRIGHT 1994-2002    *
  10  * BY THE Xiph.Org FOUNDATION http://www.xiph.org/                  *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: basic shared codebook operations
  15
  16  ********************************************************************/
  17
  18 #include <stdlib.h>
  19 #include <math.h>
  20 #include <string.h>
  21 #include <alloca.h>
  22 #include "ogg.h"
  23 #include "os.h"
  24 #include "misc.h"
  25 #include "ivorbiscodec.h"
  26 #include "codebook.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 static ogg_int32_t _float32_unpack(long val,int *point){
  47   long   mant=val&0x1fffff;
  48   int    sign=val&0x80000000;
  49   long   exp =(val&0x7fe00000L)>>VQ_FMAN;
  50
  51   exp-=(VQ_FMAN-1)+VQ_FEXP_BIAS;
  52
  53   if(mant){
  54     while(!(mant&0x40000000)){
  55       mant<<=1;
  56       exp-=1;
  57     }
  58
  59     if(sign)mant= -mant;
  60   }else{
  61     sign=0;
  62     exp=-9999;
  63   }
  64
  65   *point=exp;
  66   return mant;
  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_uint32_t *)_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   /* get us a starting hint, we'll polish it below */
 151   int bits=_ilog(b->entries);
 152   int vals=b->entries>>((bits-1)*(b->dim-1)/b->dim);
 153
 154   while(1){
 155     long acc=1;
 156     long acc1=1;
 157     int i;
 158     for(i=0;i<b->dim;i++){
 159       acc*=vals;
 160       acc1*=vals+1;
 161     }
 162     if(acc<=b->entries && acc1>b->entries){
 163       return(vals);
 164     }else{
 165       if(acc>b->entries){
 166         vals--;
 167       }else{
 168         vals++;
 169       }
 170     }
 171   }
 172 }
 173
 174 /* different than what _book_unquantize does for mainline:
 175    we repack the book in a fixed point format that shares the same
 176    binary point.  Upon first use, we can shift point if needed */
 177
 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
 183 ogg_int32_t *_book_unquantize(const static_codebook *b,int n,int *sparsemap,
 184                               int *maxpoint){
 185   long j,k,count=0;
 186   if(b->maptype==1 || b->maptype==2){
 187     int quantvals;
 188     int minpoint,delpoint;
 189     ogg_int32_t mindel=_float32_unpack(b->q_min,&minpoint);
 190     ogg_int32_t delta=_float32_unpack(b->q_delta,&delpoint);
 191     ogg_int32_t *r=(ogg_int32_t *)_ogg_calloc(n*b->dim,sizeof(*r));
 192     int *rp=(int *)_ogg_calloc(n*b->dim,sizeof(*rp));
 193
 194     *maxpoint=minpoint;
 195
 196     /* maptype 1 and 2 both use a quantized value vector, but
 197        different sizes */
 198     switch(b->maptype){
 199     case 1:
 200       /* most of the time, entries%dimensions == 0, but we need to be
 201          well defined.  We define that the possible vales at each
 202          scalar is values == entries/dim.  If entries%dim != 0, we'll
 203          have 'too few' values (values*dim<entries), which means that
 204          we'll have 'left over' entries; left over entries use zeroed
 205          values (and are wasted).  So don't generate codebooks like
 206          that */
 207       quantvals=_book_maptype1_quantvals(b);
 208       for(j=0;j<b->entries;j++){
 209         if((sparsemap && b->lengthlist[j]) || !sparsemap){
 210           ogg_int32_t last=0;
 211           int lastpoint=0;
 212           int indexdiv=1;
 213           for(k=0;k<b->dim;k++){
 214             int index= (j/indexdiv)%quantvals;
 215             ogg_int32_t point;
 216             int val=VFLOAT_MULTI(delta,delpoint,
 217                                  abs(b->quantlist[index]),&point);
 218
 219             val=VFLOAT_ADD(mindel,minpoint,val,point,&point);
 220             val=VFLOAT_ADD(last,lastpoint,val,point,&point);
 221
 222             if(b->q_sequencep){
 223               last=val;
 224               lastpoint=point;
 225             }
 226
 227             if(sparsemap){
 228               r[sparsemap[count]*b->dim+k]=val;
 229               rp[sparsemap[count]*b->dim+k]=point;
 230             }else{
 231               r[count*b->dim+k]=val;
 232               rp[count*b->dim+k]=point;
 233             }
 234             if(*maxpoint<point)*maxpoint=point;
 235             indexdiv*=quantvals;
 236           }
 237           count++;
 238         }
 239
 240       }
 241       break;
 242     case 2:
 243       for(j=0;j<b->entries;j++){
 244         if((sparsemap && b->lengthlist[j]) || !sparsemap){
 245           ogg_int32_t last=0;
 246           int         lastpoint=0;
 247
 248           for(k=0;k<b->dim;k++){
 249             ogg_int32_t point;
 250             int val=VFLOAT_MULTI(delta,delpoint,
 251                                  abs(b->quantlist[j*b->dim+k]),&point);
 252
 253             val=VFLOAT_ADD(mindel,minpoint,val,point,&point);
 254             val=VFLOAT_ADD(last,lastpoint,val,point,&point);
 255
 256             if(b->q_sequencep){
 257               last=val;
 258               lastpoint=point;
 259             }
 260
 261             if(sparsemap){
 262               r[sparsemap[count]*b->dim+k]=val;
 263               rp[sparsemap[count]*b->dim+k]=point;
 264             }else{
 265               r[count*b->dim+k]=val;
 266               rp[count*b->dim+k]=point;
 267             }
 268             if(*maxpoint<point)*maxpoint=point;
 269           }
 270           count++;
 271         }
 272       }
 273       break;
 274     }
 275
 276     for(j=0;j<n*b->dim;j++)
 277       if(rp[j]<*maxpoint)
 278         r[j]>>=*maxpoint-rp[j];
 279
 280     _ogg_free(rp);
 281     return(r);
 282   }
 283   return(NULL);
 284 }
 285
 286 void vorbis_staticbook_clear(static_codebook *b){
 287   if(b->quantlist)_ogg_free(b->quantlist);
 288   if(b->lengthlist)_ogg_free(b->lengthlist);
 289   memset(b,0,sizeof(*b));
 290
 291 }
 292
 293 void vorbis_staticbook_destroy(static_codebook *b){
 294   vorbis_staticbook_clear(b);
 295   _ogg_free(b);
 296 }
 297
 298 void vorbis_book_clear(codebook *b){
 299   /* static book is not cleared; we're likely called on the lookup and
 300      the static codebook belongs to the info struct */
 301   if(b->valuelist)_ogg_free(b->valuelist);
 302   if(b->codelist)_ogg_free(b->codelist);
 303
 304   if(b->dec_index)_ogg_free(b->dec_index);
 305   if(b->dec_codelengths)_ogg_free(b->dec_codelengths);
 306   if(b->dec_firsttable)_ogg_free(b->dec_firsttable);
 307
 308   memset(b,0,sizeof(*b));
 309 }
 310
 311 static ogg_uint32_t bitreverse(ogg_uint32_t x){
 312   x=    ((x>>16)&0x0000ffffUL) | ((x<<16)&0xffff0000UL);
 313   x=    ((x>> 8)&0x00ff00ffUL) | ((x<< 8)&0xff00ff00UL);
 314   x=    ((x>> 4)&0x0f0f0f0fUL) | ((x<< 4)&0xf0f0f0f0UL);
 315   x=    ((x>> 2)&0x33333333UL) | ((x<< 2)&0xccccccccUL);
 316   return((x>> 1)&0x55555555UL) | ((x<< 1)&0xaaaaaaaaUL);
 317 }
 318
 319 static int sort32a(const void *a,const void *b){
 320   return (**(ogg_uint32_t **)a>**(ogg_uint32_t **)b)-
 321     (**(ogg_uint32_t **)a<**(ogg_uint32_t **)b);
 322 }
 323
 324 /* decode codebook arrangement is more heavily optimized than encode */
 325 int vorbis_book_init_decode(codebook *c,const static_codebook *s){
 326   int i,j,n=0,tabn;
 327   int *sortindex;
 328   memset(c,0,sizeof(*c));
 329
 330   /* count actually used entries */
 331   for(i=0;i<s->entries;i++)
 332     if(s->lengthlist[i]>0)
 333       n++;
 334
 335   c->entries=s->entries;
 336   c->used_entries=n;
 337   c->dim=s->dim;
 338
 339   c->q_min=s->q_min;
 340   c->q_delta=s->q_delta;
 341
 342   /* two different remappings go on here.
 343
 344      First, we collapse the likely sparse codebook down only to
 345      actually represented values/words.  This collapsing needs to be
 346      indexed as map-valueless books are used to encode original entry
 347      positions as integers.
 348
 349      Second, we reorder all vectors, including the entry index above,
 350      by sorted bitreversed codeword to allow treeless decode. */
 351
 352   {
 353     /* perform sort */
 354     ogg_uint32_t *codes=_make_words(s->lengthlist,s->entries,c->used_entries);
 355     ogg_uint32_t **codep=(ogg_uint32_t **)alloca(sizeof(*codep)*n);
 356
 357     if(codes==NULL)goto err_out;
 358
 359     for(i=0;i<n;i++){
 360       codes[i]=bitreverse(codes[i]);
 361       codep[i]=codes+i;
 362     }
 363
 364     qsort(codep,n,sizeof(*codep),sort32a);
 365
 366     sortindex=(int *)alloca(n*sizeof(*sortindex));
 367     c->codelist=(ogg_uint32_t *)_ogg_malloc(n*sizeof(*c->codelist));
 368     /* the index is a reverse index */
 369     for(i=0;i<n;i++){
 370       int position=codep[i]-codes;
 371       sortindex[position]=i;
 372     }
 373
 374     for(i=0;i<n;i++)
 375       c->codelist[sortindex[i]]=codes[i];
 376     _ogg_free(codes);
 377   }
 378
 379
 380   c->valuelist=_book_unquantize(s,n,sortindex,&c->binarypoint);
 381   c->dec_index=(int *)_ogg_malloc(n*sizeof(*c->dec_index));
 382
 383   for(n=0,i=0;i<s->entries;i++)
 384     if(s->lengthlist[i]>0)
 385       c->dec_index[sortindex[n++]]=i;
 386
 387   c->dec_codelengths=(char *)_ogg_malloc(n*sizeof(*c->dec_codelengths));
 388   for(n=0,i=0;i<s->entries;i++)
 389     if(s->lengthlist[i]>0)
 390       c->dec_codelengths[sortindex[n++]]=s->lengthlist[i];
 391
 392   c->dec_firsttablen=_ilog(c->used_entries)-4; /* this is magic */
 393   if(c->dec_firsttablen<5)c->dec_firsttablen=5;
 394   if(c->dec_firsttablen>8)c->dec_firsttablen=8;
 395
 396   tabn=1<<c->dec_firsttablen;
 397   c->dec_firsttable=(ogg_uint32_t *)_ogg_calloc(tabn,sizeof(*c->dec_firsttable));
 398   c->dec_maxlength=0;
 399
 400   for(i=0;i<n;i++){
 401     if(c->dec_maxlength<c->dec_codelengths[i])
 402       c->dec_maxlength=c->dec_codelengths[i];
 403     if(c->dec_codelengths[i]<=c->dec_firsttablen){
 404       ogg_uint32_t orig=bitreverse(c->codelist[i]);
 405       for(j=0;j<(1<<(c->dec_firsttablen-c->dec_codelengths[i]));j++)
 406         c->dec_firsttable[orig|(j<<c->dec_codelengths[i])]=i+1;
 407     }
 408   }
 409
 410   /* now fill in 'unused' entries in the firsttable with hi/lo search
 411      hints for the non-direct-hits */
 412   {
 413     ogg_uint32_t mask=0xfffffffeUL<<(31-c->dec_firsttablen);
 414     long lo=0,hi=0;
 415
 416     for(i=0;i<tabn;i++){
 417       ogg_uint32_t word=i<<(32-c->dec_firsttablen);
 418       if(c->dec_firsttable[bitreverse(word)]==0){
 419         while((lo+1)<n && c->codelist[lo+1]<=word)lo++;
 420         while(    hi<n && word>=(c->codelist[hi]&mask))hi++;
 421
 422         /* we only actually have 15 bits per hint to play with here.
 423            In order to overflow gracefully (nothing breaks, efficiency
 424            just drops), encode as the difference from the extremes. */
 425         {
 426           unsigned long loval=lo;
 427           unsigned long hival=n-hi;
 428
 429           if(loval>0x7fff)loval=0x7fff;
 430           if(hival>0x7fff)hival=0x7fff;
 431           c->dec_firsttable[bitreverse(word)]=
 432             0x80000000UL | (loval<<15) | hival;
 433         }
 434       }
 435     }
 436   }
 437
 438
 439   return(0);
 440  err_out:
 441   vorbis_book_clear(c);
 442   return(-1);
 443 }
 444