misc/libtremor/tremor/sharedbook.c
changeset 9431 0f5961910e27
parent 9357 a501f5ec7b34
parent 9429 7a97a554ac80
child 9433 f0a8ac191839
equal deleted inserted replaced
9357:a501f5ec7b34 9431:0f5961910e27
     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 "ogg.h"
       
    22 #include "misc.h"
       
    23 #include "ivorbiscodec.h"
       
    24 #include "codebook.h"
       
    25 
       
    26 /**** pack/unpack helpers ******************************************/
       
    27 int _ilog(unsigned int v){
       
    28   int ret=0;
       
    29   while(v){
       
    30     ret++;
       
    31     v>>=1;
       
    32   }
       
    33   return(ret);
       
    34 }
       
    35 
       
    36 /* 32 bit float (not IEEE; nonnormalized mantissa +
       
    37    biased exponent) : neeeeeee eeemmmmm mmmmmmmm mmmmmmmm 
       
    38    Why not IEEE?  It's just not that important here. */
       
    39 
       
    40 #define VQ_FEXP 10
       
    41 #define VQ_FMAN 21
       
    42 #define VQ_FEXP_BIAS 768 /* bias toward values smaller than 1. */
       
    43 
       
    44 static ogg_int32_t _float32_unpack(long val,int *point){
       
    45   long   mant=val&0x1fffff;
       
    46   int    sign=val&0x80000000;
       
    47   long   exp =(val&0x7fe00000L)>>VQ_FMAN;
       
    48 
       
    49   exp-=(VQ_FMAN-1)+VQ_FEXP_BIAS;
       
    50 
       
    51   if(mant){
       
    52     while(!(mant&0x40000000)){
       
    53       mant<<=1;
       
    54       exp-=1;
       
    55     }
       
    56 
       
    57     if(sign)mant= -mant;
       
    58   }else{
       
    59     sign=0;
       
    60     exp=-9999;
       
    61   }
       
    62 
       
    63   *point=exp;
       
    64   return mant;
       
    65 }
       
    66 
       
    67 /* given a list of word lengths, generate a list of codewords.  Works
       
    68    for length ordered or unordered, always assigns the lowest valued
       
    69    codewords first.  Extended to handle unused entries (length 0) */
       
    70 ogg_uint32_t *_make_words(long *l,long n,long sparsecount){
       
    71   long i,j,count=0;
       
    72   ogg_uint32_t marker[33];
       
    73   ogg_uint32_t *r=(ogg_uint32_t *)_ogg_malloc((sparsecount?sparsecount:n)*sizeof(*r));
       
    74   memset(marker,0,sizeof(marker));
       
    75 
       
    76   for(i=0;i<n;i++){
       
    77     long length=l[i];
       
    78     if(length>0){
       
    79       ogg_uint32_t entry=marker[length];
       
    80       
       
    81       /* when we claim a node for an entry, we also claim the nodes
       
    82 	 below it (pruning off the imagined tree that may have dangled
       
    83 	 from it) as well as blocking the use of any nodes directly
       
    84 	 above for leaves */
       
    85       
       
    86       /* update ourself */
       
    87       if(length<32 && (entry>>length)){
       
    88 	/* error condition; the lengths must specify an overpopulated tree */
       
    89 	_ogg_free(r);
       
    90 	return(NULL);
       
    91       }
       
    92       r[count++]=entry;
       
    93     
       
    94       /* Look to see if the next shorter marker points to the node
       
    95 	 above. if so, update it and repeat.  */
       
    96       {
       
    97 	for(j=length;j>0;j--){
       
    98 	  
       
    99 	  if(marker[j]&1){
       
   100 	    /* have to jump branches */
       
   101 	    if(j==1)
       
   102 	      marker[1]++;
       
   103 	    else
       
   104 	      marker[j]=marker[j-1]<<1;
       
   105 	    break; /* invariant says next upper marker would already
       
   106 		      have been moved if it was on the same path */
       
   107 	  }
       
   108 	  marker[j]++;
       
   109 	}
       
   110       }
       
   111       
       
   112       /* prune the tree; the implicit invariant says all the longer
       
   113 	 markers were dangling from our just-taken node.  Dangle them
       
   114 	 from our *new* node. */
       
   115       for(j=length+1;j<33;j++)
       
   116 	if((marker[j]>>1) == entry){
       
   117 	  entry=marker[j];
       
   118 	  marker[j]=marker[j-1]<<1;
       
   119 	}else
       
   120 	  break;
       
   121     }else
       
   122       if(sparsecount==0)count++;
       
   123   }
       
   124     
       
   125   /* bitreverse the words because our bitwise packer/unpacker is LSb
       
   126      endian */
       
   127   for(i=0,count=0;i<n;i++){
       
   128     ogg_uint32_t temp=0;
       
   129     for(j=0;j<l[i];j++){
       
   130       temp<<=1;
       
   131       temp|=(r[count]>>j)&1;
       
   132     }
       
   133 
       
   134     if(sparsecount){
       
   135       if(l[i])
       
   136 	r[count++]=temp;
       
   137     }else
       
   138       r[count++]=temp;
       
   139   }
       
   140 
       
   141   return(r);
       
   142 }
       
   143 
       
   144 /* there might be a straightforward one-line way to do the below
       
   145    that's portable and totally safe against roundoff, but I haven't
       
   146    thought of it.  Therefore, we opt on the side of caution */
       
   147 long _book_maptype1_quantvals(const static_codebook *b){
       
   148   /* get us a starting hint, we'll polish it below */
       
   149   int bits=_ilog(b->entries);
       
   150   int vals=b->entries>>((bits-1)*(b->dim-1)/b->dim);
       
   151 
       
   152   while(1){
       
   153     long acc=1;
       
   154     long acc1=1;
       
   155     int i;
       
   156     for(i=0;i<b->dim;i++){
       
   157       acc*=vals;
       
   158       acc1*=vals+1;
       
   159     }
       
   160     if(acc<=b->entries && acc1>b->entries){
       
   161       return(vals);
       
   162     }else{
       
   163       if(acc>b->entries){
       
   164 	vals--;
       
   165       }else{
       
   166 	vals++;
       
   167       }
       
   168     }
       
   169   }
       
   170 }
       
   171 
       
   172 /* different than what _book_unquantize does for mainline:
       
   173    we repack the book in a fixed point format that shares the same
       
   174    binary point.  Upon first use, we can shift point if needed */
       
   175 
       
   176 /* we need to deal with two map types: in map type 1, the values are
       
   177    generated algorithmically (each column of the vector counts through
       
   178    the values in the quant vector). in map type 2, all the values came
       
   179    in in an explicit list.  Both value lists must be unpacked */
       
   180 
       
   181 ogg_int32_t *_book_unquantize(const static_codebook *b,int n,int *sparsemap,
       
   182 			      int *maxpoint){
       
   183   long j,k,count=0;
       
   184   if(b->maptype==1 || b->maptype==2){
       
   185     int quantvals;
       
   186     int minpoint,delpoint;
       
   187     ogg_int32_t mindel=_float32_unpack(b->q_min,&minpoint);
       
   188     ogg_int32_t delta=_float32_unpack(b->q_delta,&delpoint);
       
   189     ogg_int32_t *r=(ogg_int32_t *)_ogg_calloc(n*b->dim,sizeof(*r));
       
   190     int *rp=(int *)_ogg_calloc(n*b->dim,sizeof(*rp));
       
   191 
       
   192     *maxpoint=minpoint;
       
   193 
       
   194     /* maptype 1 and 2 both use a quantized value vector, but
       
   195        different sizes */
       
   196     switch(b->maptype){
       
   197     case 1:
       
   198       /* most of the time, entries%dimensions == 0, but we need to be
       
   199 	 well defined.  We define that the possible vales at each
       
   200 	 scalar is values == entries/dim.  If entries%dim != 0, we'll
       
   201 	 have 'too few' values (values*dim<entries), which means that
       
   202 	 we'll have 'left over' entries; left over entries use zeroed
       
   203 	 values (and are wasted).  So don't generate codebooks like
       
   204 	 that */
       
   205       quantvals=_book_maptype1_quantvals(b);
       
   206       for(j=0;j<b->entries;j++){
       
   207 	if((sparsemap && b->lengthlist[j]) || !sparsemap){
       
   208 	  ogg_int32_t last=0;
       
   209 	  int lastpoint=0;
       
   210 	  int indexdiv=1;
       
   211 	  for(k=0;k<b->dim;k++){
       
   212 	    int index= (j/indexdiv)%quantvals;
       
   213 	    int point=0;
       
   214 	    int val=VFLOAT_MULTI(delta,delpoint,
       
   215 				 abs(b->quantlist[index]),&point);
       
   216 
       
   217 	    val=VFLOAT_ADD(mindel,minpoint,val,point,&point);
       
   218 	    val=VFLOAT_ADD(last,lastpoint,val,point,&point);
       
   219 	    
       
   220 	    if(b->q_sequencep){
       
   221 	      last=val;	  
       
   222 	      lastpoint=point;
       
   223 	    }
       
   224 	    
       
   225 	    if(sparsemap){
       
   226 	      r[sparsemap[count]*b->dim+k]=val;
       
   227 	      rp[sparsemap[count]*b->dim+k]=point;
       
   228 	    }else{
       
   229 	      r[count*b->dim+k]=val;
       
   230 	      rp[count*b->dim+k]=point;
       
   231 	    }
       
   232 	    if(*maxpoint<point)*maxpoint=point;
       
   233 	    indexdiv*=quantvals;
       
   234 	  }
       
   235 	  count++;
       
   236 	}
       
   237 
       
   238       }
       
   239       break;
       
   240     case 2:
       
   241       for(j=0;j<b->entries;j++){
       
   242 	if((sparsemap && b->lengthlist[j]) || !sparsemap){
       
   243 	  ogg_int32_t last=0;
       
   244 	  int         lastpoint=0;
       
   245 
       
   246 	  for(k=0;k<b->dim;k++){
       
   247 	    int point=0;
       
   248 	    int val=VFLOAT_MULTI(delta,delpoint,
       
   249 				 abs(b->quantlist[j*b->dim+k]),&point);
       
   250 
       
   251 	    val=VFLOAT_ADD(mindel,minpoint,val,point,&point);
       
   252 	    val=VFLOAT_ADD(last,lastpoint,val,point,&point);
       
   253 	    
       
   254 	    if(b->q_sequencep){
       
   255 	      last=val;	  
       
   256 	      lastpoint=point;
       
   257 	    }
       
   258 
       
   259 	    if(sparsemap){
       
   260 	      r[sparsemap[count]*b->dim+k]=val;
       
   261 	      rp[sparsemap[count]*b->dim+k]=point;
       
   262 	    }else{
       
   263 	      r[count*b->dim+k]=val;
       
   264 	      rp[count*b->dim+k]=point;
       
   265 	    }
       
   266 	    if(*maxpoint<point)*maxpoint=point;
       
   267 	  }
       
   268 	  count++;
       
   269 	}
       
   270       }
       
   271       break;
       
   272     }
       
   273 
       
   274     for(j=0;j<n*b->dim;j++)
       
   275       if(rp[j]<*maxpoint)
       
   276 	r[j]>>=*maxpoint-rp[j];
       
   277 	    
       
   278     _ogg_free(rp);
       
   279     return(r);
       
   280   }
       
   281   return(NULL);
       
   282 }
       
   283 
       
   284 void vorbis_staticbook_clear(static_codebook *b){
       
   285   if(b->quantlist)_ogg_free(b->quantlist);
       
   286   if(b->lengthlist)_ogg_free(b->lengthlist);
       
   287   memset(b,0,sizeof(*b));
       
   288 
       
   289 }
       
   290 
       
   291 void vorbis_staticbook_destroy(static_codebook *b){
       
   292   vorbis_staticbook_clear(b);
       
   293   _ogg_free(b);
       
   294 }
       
   295 
       
   296 void vorbis_book_clear(codebook *b){
       
   297   /* static book is not cleared; we're likely called on the lookup and
       
   298      the static codebook belongs to the info struct */
       
   299   if(b->valuelist)_ogg_free(b->valuelist);
       
   300   if(b->codelist)_ogg_free(b->codelist);
       
   301 
       
   302   if(b->dec_index)_ogg_free(b->dec_index);
       
   303   if(b->dec_codelengths)_ogg_free(b->dec_codelengths);
       
   304   if(b->dec_firsttable)_ogg_free(b->dec_firsttable);
       
   305 
       
   306   memset(b,0,sizeof(*b));
       
   307 }
       
   308 
       
   309 static ogg_uint32_t bitreverse(ogg_uint32_t x){
       
   310   x=    ((x>>16)&0x0000ffffUL) | ((x<<16)&0xffff0000UL);
       
   311   x=    ((x>> 8)&0x00ff00ffUL) | ((x<< 8)&0xff00ff00UL);
       
   312   x=    ((x>> 4)&0x0f0f0f0fUL) | ((x<< 4)&0xf0f0f0f0UL);
       
   313   x=    ((x>> 2)&0x33333333UL) | ((x<< 2)&0xccccccccUL);
       
   314   return((x>> 1)&0x55555555UL) | ((x<< 1)&0xaaaaaaaaUL);
       
   315 }
       
   316 
       
   317 static int sort32a(const void *a,const void *b){
       
   318   return (**(ogg_uint32_t **)a>**(ogg_uint32_t **)b)-
       
   319     (**(ogg_uint32_t **)a<**(ogg_uint32_t **)b);
       
   320 }
       
   321 
       
   322 /* decode codebook arrangement is more heavily optimized than encode */
       
   323 int vorbis_book_init_decode(codebook *c,const static_codebook *s){
       
   324   int i,j,n=0,tabn;
       
   325   int *sortindex;
       
   326   memset(c,0,sizeof(*c));
       
   327   
       
   328   /* count actually used entries */
       
   329   for(i=0;i<s->entries;i++)
       
   330     if(s->lengthlist[i]>0)
       
   331       n++;
       
   332 
       
   333   c->entries=s->entries;
       
   334   c->used_entries=n;
       
   335   c->dim=s->dim;
       
   336 
       
   337   if(n>0){
       
   338     /* two different remappings go on here.  
       
   339        
       
   340        First, we collapse the likely sparse codebook down only to
       
   341        actually represented values/words.  This collapsing needs to be
       
   342        indexed as map-valueless books are used to encode original entry
       
   343        positions as integers.
       
   344        
       
   345        Second, we reorder all vectors, including the entry index above,
       
   346        by sorted bitreversed codeword to allow treeless decode. */
       
   347     
       
   348     /* perform sort */
       
   349     ogg_uint32_t *codes=_make_words(s->lengthlist,s->entries,c->used_entries);
       
   350     ogg_uint32_t **codep=(ogg_uint32_t **)alloca(sizeof(*codep)*n);
       
   351     
       
   352     if(codes==NULL)goto err_out;
       
   353 
       
   354     for(i=0;i<n;i++){
       
   355       codes[i]=bitreverse(codes[i]);
       
   356       codep[i]=codes+i;
       
   357     }
       
   358 
       
   359     qsort(codep,n,sizeof(*codep),sort32a);
       
   360 
       
   361     sortindex=(int *)alloca(n*sizeof(*sortindex));
       
   362     c->codelist=(ogg_uint32_t *)_ogg_malloc(n*sizeof(*c->codelist));
       
   363     /* the index is a reverse index */
       
   364     for(i=0;i<n;i++){
       
   365       int position=codep[i]-codes;
       
   366       sortindex[position]=i;
       
   367     }
       
   368 
       
   369     for(i=0;i<n;i++)
       
   370       c->codelist[sortindex[i]]=codes[i];
       
   371     _ogg_free(codes);
       
   372     
       
   373     
       
   374     
       
   375     c->valuelist=_book_unquantize(s,n,sortindex,&c->binarypoint);
       
   376     c->dec_index=(int *)_ogg_malloc(n*sizeof(*c->dec_index));
       
   377     
       
   378     for(n=0,i=0;i<s->entries;i++)
       
   379       if(s->lengthlist[i]>0)
       
   380 	c->dec_index[sortindex[n++]]=i;
       
   381     
       
   382     c->dec_codelengths=(char *)_ogg_malloc(n*sizeof(*c->dec_codelengths));
       
   383     for(n=0,i=0;i<s->entries;i++)
       
   384       if(s->lengthlist[i]>0)
       
   385 	c->dec_codelengths[sortindex[n++]]=s->lengthlist[i];
       
   386     
       
   387     c->dec_firsttablen=_ilog(c->used_entries)-4; /* this is magic */
       
   388     if(c->dec_firsttablen<5)c->dec_firsttablen=5;
       
   389     if(c->dec_firsttablen>8)c->dec_firsttablen=8;
       
   390     
       
   391     tabn=1<<c->dec_firsttablen;
       
   392     c->dec_firsttable=(ogg_uint32_t *)_ogg_calloc(tabn,sizeof(*c->dec_firsttable));
       
   393     c->dec_maxlength=0;
       
   394     
       
   395     for(i=0;i<n;i++){
       
   396       if(c->dec_maxlength<c->dec_codelengths[i])
       
   397 	c->dec_maxlength=c->dec_codelengths[i];
       
   398       if(c->dec_codelengths[i]<=c->dec_firsttablen){
       
   399 	ogg_uint32_t orig=bitreverse(c->codelist[i]);
       
   400 	for(j=0;j<(1<<(c->dec_firsttablen-c->dec_codelengths[i]));j++)
       
   401 	  c->dec_firsttable[orig|(j<<c->dec_codelengths[i])]=i+1;
       
   402       }
       
   403     }
       
   404     
       
   405     /* now fill in 'unused' entries in the firsttable with hi/lo search
       
   406        hints for the non-direct-hits */
       
   407     {
       
   408       ogg_uint32_t mask=0xfffffffeUL<<(31-c->dec_firsttablen);
       
   409       long lo=0,hi=0;
       
   410       
       
   411       for(i=0;i<tabn;i++){
       
   412 	ogg_uint32_t word=i<<(32-c->dec_firsttablen);
       
   413 	if(c->dec_firsttable[bitreverse(word)]==0){
       
   414 	  while((lo+1)<n && c->codelist[lo+1]<=word)lo++;
       
   415 	  while(    hi<n && word>=(c->codelist[hi]&mask))hi++;
       
   416 	  
       
   417 	  /* we only actually have 15 bits per hint to play with here.
       
   418 	     In order to overflow gracefully (nothing breaks, efficiency
       
   419 	     just drops), encode as the difference from the extremes. */
       
   420 	  {
       
   421 	    unsigned long loval=lo;
       
   422 	    unsigned long hival=n-hi;
       
   423 	    
       
   424 	    if(loval>0x7fff)loval=0x7fff;
       
   425 	    if(hival>0x7fff)hival=0x7fff;
       
   426 	    c->dec_firsttable[bitreverse(word)]=
       
   427 	      0x80000000UL | (loval<<15) | hival;
       
   428 	  }
       
   429 	}
       
   430       }
       
   431     }
       
   432   }
       
   433 
       
   434   return(0);
       
   435  err_out:
       
   436   vorbis_book_clear(c);
       
   437   return(-1);
       
   438 }
       
   439