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-2015 *
9 * by the Xiph.Org Foundation http://www.xiph.org/ *
10 * *
11 ********************************************************************
12
13 function: basic shared codebook operations
14 last mod: $Id: sharedbook.c 19457 2015-03-03 00:15:29Z giles $
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
ov_ilog(ogg_uint32_t v)30 int ov_ilog(ogg_uint32_t v){
31 int ret;
32 for(ret=0;v;ret++)v>>=1;
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 /* doesn't currently guard under/overflow */
_float32_pack(float val)45 long _float32_pack(float val){
46 int sign=0;
47 long exp;
48 long mant;
49 if(val<0){
50 sign=0x80000000;
51 val= -val;
52 }
53 exp= floor(log(val)/log(2.f)+.001); //+epsilon
54 mant=rint(ldexp(val,(VQ_FMAN-1)-exp));
55 exp=(exp+VQ_FEXP_BIAS)<<VQ_FMAN;
56
57 return(sign|exp|mant);
58 }
59
_float32_unpack(long val)60 float _float32_unpack(long val){
61 double mant=val&0x1fffff;
62 int sign=val&0x80000000;
63 long exp =(val&0x7fe00000L)>>VQ_FMAN;
64 if(sign)mant= -mant;
65 return(ldexp(mant,exp-(VQ_FMAN-1)-VQ_FEXP_BIAS));
66 }
67
68 /* given a list of word lengths, generate a list of codewords. Works
69 for length ordered or unordered, always assigns the lowest valued
70 codewords first. Extended to handle unused entries (length 0) */
_make_words(char * l,long n,long sparsecount)71 ogg_uint32_t *_make_words(char *l,long n,long sparsecount){
72 long i,j,count=0;
73 ogg_uint32_t marker[33];
74 ogg_uint32_t *r=_ogg_malloc((sparsecount?sparsecount:n)*sizeof(*r));
75 memset(marker,0,sizeof(marker));
76
77 for(i=0;i<n;i++){
78 long length=l[i];
79 if(length>0){
80 ogg_uint32_t entry=marker[length];
81
82 /* when we claim a node for an entry, we also claim the nodes
83 below it (pruning off the imagined tree that may have dangled
84 from it) as well as blocking the use of any nodes directly
85 above for leaves */
86
87 /* update ourself */
88 if(length<32 && (entry>>length)){
89 /* error condition; the lengths must specify an overpopulated tree */
90 _ogg_free(r);
91 return(NULL);
92 }
93 r[count++]=entry;
94
95 /* Look to see if the next shorter marker points to the node
96 above. if so, update it and repeat. */
97 {
98 for(j=length;j>0;j--){
99
100 if(marker[j]&1){
101 /* have to jump branches */
102 if(j==1)
103 marker[1]++;
104 else
105 marker[j]=marker[j-1]<<1;
106 break; /* invariant says next upper marker would already
107 have been moved if it was on the same path */
108 }
109 marker[j]++;
110 }
111 }
112
113 /* prune the tree; the implicit invariant says all the longer
114 markers were dangling from our just-taken node. Dangle them
115 from our *new* node. */
116 for(j=length+1;j<33;j++)
117 if((marker[j]>>1) == entry){
118 entry=marker[j];
119 marker[j]=marker[j-1]<<1;
120 }else
121 break;
122 }else
123 if(sparsecount==0)count++;
124 }
125
126 /* any underpopulated tree must be rejected. */
127 /* Single-entry codebooks are a retconned extension to the spec.
128 They have a single codeword '0' of length 1 that results in an
129 underpopulated tree. Shield that case from the underformed tree check. */
130 if(!(count==1 && marker[2]==2)){
131 for(i=1;i<33;i++)
132 if(marker[i] & (0xffffffffUL>>(32-i))){
133 _ogg_free(r);
134 return(NULL);
135 }
136 }
137
138 /* bitreverse the words because our bitwise packer/unpacker is LSb
139 endian */
140 for(i=0,count=0;i<n;i++){
141 ogg_uint32_t temp=0;
142 for(j=0;j<l[i];j++){
143 temp<<=1;
144 temp|=(r[count]>>j)&1;
145 }
146
147 if(sparsecount){
148 if(l[i])
149 r[count++]=temp;
150 }else
151 r[count++]=temp;
152 }
153
154 return(r);
155 }
156
157 /* there might be a straightforward one-line way to do the below
158 that's portable and totally safe against roundoff, but I haven't
159 thought of it. Therefore, we opt on the side of caution */
_book_maptype1_quantvals(const static_codebook * b)160 long _book_maptype1_quantvals(const static_codebook *b){
161 long vals=floor(pow((float)b->entries,1.f/b->dim));
162
163 /* the above *should* be reliable, but we'll not assume that FP is
164 ever reliable when bitstream sync is at stake; verify via integer
165 means that vals really is the greatest value of dim for which
166 vals^b->bim <= b->entries */
167 /* treat the above as an initial guess */
168 while(1){
169 long acc=1;
170 long acc1=1;
171 int i;
172 for(i=0;i<b->dim;i++){
173 acc*=vals;
174 acc1*=vals+1;
175 }
176 if(acc<=b->entries && acc1>b->entries){
177 return(vals);
178 }else{
179 if(acc>b->entries){
180 vals--;
181 }else{
182 vals++;
183 }
184 }
185 }
186 }
187
188 /* unpack the quantized list of values for encode/decode ***********/
189 /* we need to deal with two map types: in map type 1, the values are
190 generated algorithmically (each column of the vector counts through
191 the values in the quant vector). in map type 2, all the values came
192 in in an explicit list. Both value lists must be unpacked */
_book_unquantize(const static_codebook * b,int n,int * sparsemap)193 float *_book_unquantize(const static_codebook *b,int n,int *sparsemap){
194 long j,k,count=0;
195 if(b->maptype==1 || b->maptype==2){
196 int quantvals;
197 float mindel=_float32_unpack(b->q_min);
198 float delta=_float32_unpack(b->q_delta);
199 float *r=_ogg_calloc(n*b->dim,sizeof(*r));
200
201 /* maptype 1 and 2 both use a quantized value vector, but
202 different sizes */
203 switch(b->maptype){
204 case 1:
205 /* most of the time, entries%dimensions == 0, but we need to be
206 well defined. We define that the possible vales at each
207 scalar is values == entries/dim. If entries%dim != 0, we'll
208 have 'too few' values (values*dim<entries), which means that
209 we'll have 'left over' entries; left over entries use zeroed
210 values (and are wasted). So don't generate codebooks like
211 that */
212 quantvals=_book_maptype1_quantvals(b);
213 for(j=0;j<b->entries;j++){
214 if((sparsemap && b->lengthlist[j]) || !sparsemap){
215 float last=0.f;
216 int indexdiv=1;
217 for(k=0;k<b->dim;k++){
218 int index= (j/indexdiv)%quantvals;
219 float val=b->quantlist[index];
220 val=fabs(val)*delta+mindel+last;
221 if(b->q_sequencep)last=val;
222 if(sparsemap)
223 r[sparsemap[count]*b->dim+k]=val;
224 else
225 r[count*b->dim+k]=val;
226 indexdiv*=quantvals;
227 }
228 count++;
229 }
230
231 }
232 break;
233 case 2:
234 for(j=0;j<b->entries;j++){
235 if((sparsemap && b->lengthlist[j]) || !sparsemap){
236 float last=0.f;
237
238 for(k=0;k<b->dim;k++){
239 float val=b->quantlist[j*b->dim+k];
240 val=fabs(val)*delta+mindel+last;
241 if(b->q_sequencep)last=val;
242 if(sparsemap)
243 r[sparsemap[count]*b->dim+k]=val;
244 else
245 r[count*b->dim+k]=val;
246 }
247 count++;
248 }
249 }
250 break;
251 }
252
253 return(r);
254 }
255 return(NULL);
256 }
257
vorbis_staticbook_destroy(static_codebook * b)258 void vorbis_staticbook_destroy(static_codebook *b){
259 if(b->allocedp){
260 if(b->quantlist)_ogg_free(b->quantlist);
261 if(b->lengthlist)_ogg_free(b->lengthlist);
262 memset(b,0,sizeof(*b));
263 _ogg_free(b);
264 } /* otherwise, it is in static memory */
265 }
266
vorbis_book_clear(codebook * b)267 void vorbis_book_clear(codebook *b){
268 /* static book is not cleared; we're likely called on the lookup and
269 the static codebook belongs to the info struct */
270 if(b->valuelist)_ogg_free(b->valuelist);
271 if(b->codelist)_ogg_free(b->codelist);
272
273 if(b->dec_index)_ogg_free(b->dec_index);
274 if(b->dec_codelengths)_ogg_free(b->dec_codelengths);
275 if(b->dec_firsttable)_ogg_free(b->dec_firsttable);
276
277 memset(b,0,sizeof(*b));
278 }
279
vorbis_book_init_encode(codebook * c,const static_codebook * s)280 int vorbis_book_init_encode(codebook *c,const static_codebook *s){
281
282 memset(c,0,sizeof(*c));
283 c->c=s;
284 c->entries=s->entries;
285 c->used_entries=s->entries;
286 c->dim=s->dim;
287 c->codelist=_make_words(s->lengthlist,s->entries,0);
288 //c->valuelist=_book_unquantize(s,s->entries,NULL);
289 c->quantvals=_book_maptype1_quantvals(s);
290 c->minval=(int)rint(_float32_unpack(s->q_min));
291 c->delta=(int)rint(_float32_unpack(s->q_delta));
292
293 return(0);
294 }
295
bitreverse(ogg_uint32_t x)296 static ogg_uint32_t bitreverse(ogg_uint32_t x){
297 x= ((x>>16)&0x0000ffffUL) | ((x<<16)&0xffff0000UL);
298 x= ((x>> 8)&0x00ff00ffUL) | ((x<< 8)&0xff00ff00UL);
299 x= ((x>> 4)&0x0f0f0f0fUL) | ((x<< 4)&0xf0f0f0f0UL);
300 x= ((x>> 2)&0x33333333UL) | ((x<< 2)&0xccccccccUL);
301 return((x>> 1)&0x55555555UL) | ((x<< 1)&0xaaaaaaaaUL);
302 }
303
sort32a(const void * a,const void * b)304 static int sort32a(const void *a,const void *b){
305 return ( **(ogg_uint32_t **)a>**(ogg_uint32_t **)b)-
306 ( **(ogg_uint32_t **)a<**(ogg_uint32_t **)b);
307 }
308
309 /* decode codebook arrangement is more heavily optimized than encode */
vorbis_book_init_decode(codebook * c,const static_codebook * s)310 int vorbis_book_init_decode(codebook *c,const static_codebook *s){
311 int i,j,n=0,tabn;
312 int *sortindex;
313
314 memset(c,0,sizeof(*c));
315
316 /* count actually used entries and find max length */
317 for(i=0;i<s->entries;i++)
318 if(s->lengthlist[i]>0)
319 n++;
320
321 c->entries=s->entries;
322 c->used_entries=n;
323 c->dim=s->dim;
324
325 if(n>0){
326 /* two different remappings go on here.
327
328 First, we collapse the likely sparse codebook down only to
329 actually represented values/words. This collapsing needs to be
330 indexed as map-valueless books are used to encode original entry
331 positions as integers.
332
333 Second, we reorder all vectors, including the entry index above,
334 by sorted bitreversed codeword to allow treeless decode. */
335
336 /* perform sort */
337 ogg_uint32_t *codes=_make_words(s->lengthlist,s->entries,c->used_entries);
338 ogg_uint32_t **codep=alloca(sizeof(*codep)*n);
339
340 if(codes==NULL)goto err_out;
341
342 for(i=0;i<n;i++){
343 codes[i]=bitreverse(codes[i]);
344 codep[i]=codes+i;
345 }
346
347 qsort(codep,n,sizeof(*codep),sort32a);
348
349 sortindex=alloca(n*sizeof(*sortindex));
350 c->codelist=_ogg_malloc(n*sizeof(*c->codelist));
351 /* the index is a reverse index */
352 for(i=0;i<n;i++){
353 int position=codep[i]-codes;
354 sortindex[position]=i;
355 }
356
357 for(i=0;i<n;i++)
358 c->codelist[sortindex[i]]=codes[i];
359 _ogg_free(codes);
360
361 c->valuelist=_book_unquantize(s,n,sortindex);
362 c->dec_index=_ogg_malloc(n*sizeof(*c->dec_index));
363
364 for(n=0,i=0;i<s->entries;i++)
365 if(s->lengthlist[i]>0)
366 c->dec_index[sortindex[n++]]=i;
367
368 c->dec_codelengths=_ogg_malloc(n*sizeof(*c->dec_codelengths));
369 c->dec_maxlength=0;
370 for(n=0,i=0;i<s->entries;i++)
371 if(s->lengthlist[i]>0){
372 c->dec_codelengths[sortindex[n++]]=s->lengthlist[i];
373 if(s->lengthlist[i]>c->dec_maxlength)
374 c->dec_maxlength=s->lengthlist[i];
375 }
376
377 if(n==1 && c->dec_maxlength==1){
378 /* special case the 'single entry codebook' with a single bit
379 fastpath table (that always returns entry 0 )in order to use
380 unmodified decode paths. */
381 c->dec_firsttablen=1;
382 c->dec_firsttable=_ogg_calloc(2,sizeof(*c->dec_firsttable));
383 c->dec_firsttable[0]=c->dec_firsttable[1]=1;
384
385 }else{
386 c->dec_firsttablen=ov_ilog(c->used_entries)-4; /* this is magic */
387 if(c->dec_firsttablen<5)c->dec_firsttablen=5;
388 if(c->dec_firsttablen>8)c->dec_firsttablen=8;
389
390 tabn=1<<c->dec_firsttablen;
391 c->dec_firsttable=_ogg_calloc(tabn,sizeof(*c->dec_firsttable));
392
393 for(i=0;i<n;i++){
394 if(c->dec_codelengths[i]<=c->dec_firsttablen){
395 ogg_uint32_t orig=bitreverse(c->codelist[i]);
396 for(j=0;j<(1<<(c->dec_firsttablen-c->dec_codelengths[i]));j++)
397 c->dec_firsttable[orig|(j<<c->dec_codelengths[i])]=i+1;
398 }
399 }
400
401 /* now fill in 'unused' entries in the firsttable with hi/lo search
402 hints for the non-direct-hits */
403 {
404 ogg_uint32_t mask=0xfffffffeUL<<(31-c->dec_firsttablen);
405 long lo=0,hi=0;
406
407 for(i=0;i<tabn;i++){
408 ogg_uint32_t word=i<<(32-c->dec_firsttablen);
409 if(c->dec_firsttable[bitreverse(word)]==0){
410 while((lo+1)<n && c->codelist[lo+1]<=word)lo++;
411 while( hi<n && word>=(c->codelist[hi]&mask))hi++;
412
413 /* we only actually have 15 bits per hint to play with here.
414 In order to overflow gracefully (nothing breaks, efficiency
415 just drops), encode as the difference from the extremes. */
416 {
417 unsigned long loval=lo;
418 unsigned long hival=n-hi;
419
420 if(loval>0x7fff)loval=0x7fff;
421 if(hival>0x7fff)hival=0x7fff;
422 c->dec_firsttable[bitreverse(word)]=
423 0x80000000UL | (loval<<15) | hival;
424 }
425 }
426 }
427 }
428 }
429 }
430
431 return(0);
432 err_out:
433 vorbis_book_clear(c);
434 return(-1);
435 }
436
vorbis_book_codeword(codebook * book,int entry)437 long vorbis_book_codeword(codebook *book,int entry){
438 if(book->c) /* only use with encode; decode optimizations are
439 allowed to break this */
440 return book->codelist[entry];
441 return -1;
442 }
443
vorbis_book_codelen(codebook * book,int entry)444 long vorbis_book_codelen(codebook *book,int entry){
445 if(book->c) /* only use with encode; decode optimizations are
446 allowed to break this */
447 return book->c->lengthlist[entry];
448 return -1;
449 }
450
451 #ifdef _V_SELFTEST
452
453 /* Unit tests of the dequantizer; this stuff will be OK
454 cross-platform, I simply want to be sure that special mapping cases
455 actually work properly; a bug could go unnoticed for a while */
456
457 #include <stdio.h>
458
459 /* cases:
460
461 no mapping
462 full, explicit mapping
463 algorithmic mapping
464
465 nonsequential
466 sequential
467 */
468
469 static long full_quantlist1[]={0,1,2,3, 4,5,6,7, 8,3,6,1};
470 static long partial_quantlist1[]={0,7,2};
471
472 /* no mapping */
473 static_codebook test1={
474 4,16,
475 NULL,
476 0,
477 0,0,0,0,
478 NULL,
479 0
480 };
481 static float *test1_result=NULL;
482
483 /* linear, full mapping, nonsequential */
484 static_codebook test2={
485 4,3,
486 NULL,
487 2,
488 -533200896,1611661312,4,0,
489 full_quantlist1,
490 0
491 };
492 static float test2_result[]={-3,-2,-1,0, 1,2,3,4, 5,0,3,-2};
493
494 /* linear, full mapping, sequential */
495 static_codebook test3={
496 4,3,
497 NULL,
498 2,
499 -533200896,1611661312,4,1,
500 full_quantlist1,
501 0
502 };
503 static float test3_result[]={-3,-5,-6,-6, 1,3,6,10, 5,5,8,6};
504
505 /* linear, algorithmic mapping, nonsequential */
506 static_codebook test4={
507 3,27,
508 NULL,
509 1,
510 -533200896,1611661312,4,0,
511 partial_quantlist1,
512 0
513 };
514 static float test4_result[]={-3,-3,-3, 4,-3,-3, -1,-3,-3,
515 -3, 4,-3, 4, 4,-3, -1, 4,-3,
516 -3,-1,-3, 4,-1,-3, -1,-1,-3,
517 -3,-3, 4, 4,-3, 4, -1,-3, 4,
518 -3, 4, 4, 4, 4, 4, -1, 4, 4,
519 -3,-1, 4, 4,-1, 4, -1,-1, 4,
520 -3,-3,-1, 4,-3,-1, -1,-3,-1,
521 -3, 4,-1, 4, 4,-1, -1, 4,-1,
522 -3,-1,-1, 4,-1,-1, -1,-1,-1};
523
524 /* linear, algorithmic mapping, sequential */
525 static_codebook test5={
526 3,27,
527 NULL,
528 1,
529 -533200896,1611661312,4,1,
530 partial_quantlist1,
531 0
532 };
533 static float test5_result[]={-3,-6,-9, 4, 1,-2, -1,-4,-7,
534 -3, 1,-2, 4, 8, 5, -1, 3, 0,
535 -3,-4,-7, 4, 3, 0, -1,-2,-5,
536 -3,-6,-2, 4, 1, 5, -1,-4, 0,
537 -3, 1, 5, 4, 8,12, -1, 3, 7,
538 -3,-4, 0, 4, 3, 7, -1,-2, 2,
539 -3,-6,-7, 4, 1, 0, -1,-4,-5,
540 -3, 1, 0, 4, 8, 7, -1, 3, 2,
541 -3,-4,-5, 4, 3, 2, -1,-2,-3};
542
run_test(static_codebook * b,float * comp)543 void run_test(static_codebook *b,float *comp){
544 float *out=_book_unquantize(b,b->entries,NULL);
545 int i;
546
547 if(comp){
548 if(!out){
549 fprintf(stderr,"_book_unquantize incorrectly returned NULL\n");
550 exit(1);
551 }
552
553 for(i=0;i<b->entries*b->dim;i++)
554 if(fabs(out[i]-comp[i])>.0001){
555 fprintf(stderr,"disagreement in unquantized and reference data:\n"
556 "position %d, %g != %g\n",i,out[i],comp[i]);
557 exit(1);
558 }
559
560 }else{
561 if(out){
562 fprintf(stderr,"_book_unquantize returned a value array: \n"
563 " correct result should have been NULL\n");
564 exit(1);
565 }
566 }
567 }
568
main()569 int main(){
570 /* run the nine dequant tests, and compare to the hand-rolled results */
571 fprintf(stderr,"Dequant test 1... ");
572 run_test(&test1,test1_result);
573 fprintf(stderr,"OK\nDequant test 2... ");
574 run_test(&test2,test2_result);
575 fprintf(stderr,"OK\nDequant test 3... ");
576 run_test(&test3,test3_result);
577 fprintf(stderr,"OK\nDequant test 4... ");
578 run_test(&test4,test4_result);
579 fprintf(stderr,"OK\nDequant test 5... ");
580 run_test(&test5,test5_result);
581 fprintf(stderr,"OK\n\n");
582
583 return(0);
584 }
585
586 #endif
587