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: floor backend 1 implementation
14 last mod: $Id: floor1.c 19457 2015-03-03 00:15:29Z giles $
15
16 ********************************************************************/
17
18 #include <stdlib.h>
19 #include <string.h>
20 #include <math.h>
21 #include <ogg/ogg.h>
22 #include "vorbis/codec.h"
23 #include "codec_internal.h"
24 #include "registry.h"
25 #include "codebook.h"
26 #include "misc.h"
27 #include "scales.h"
28
29 #include <stdio.h>
30
31 #define floor1_rangedB 140 /* floor 1 fixed at -140dB to 0dB range */
32
33 typedef struct lsfit_acc{
34 int x0;
35 int x1;
36
37 int xa;
38 int ya;
39 int x2a;
40 int y2a;
41 int xya;
42 int an;
43
44 int xb;
45 int yb;
46 int x2b;
47 int y2b;
48 int xyb;
49 int bn;
50 } lsfit_acc;
51
52 /***********************************************/
53
floor1_free_info(vorbis_info_floor * i)54 static void floor1_free_info(vorbis_info_floor *i){
55 vorbis_info_floor1 *info=(vorbis_info_floor1 *)i;
56 if(info){
57 memset(info,0,sizeof(*info));
58 _ogg_free(info);
59 }
60 }
61
floor1_free_look(vorbis_look_floor * i)62 static void floor1_free_look(vorbis_look_floor *i){
63 vorbis_look_floor1 *look=(vorbis_look_floor1 *)i;
64 if(look){
65 /*fprintf(stderr,"floor 1 bit usage %f:%f (%f total)\n",
66 (float)look->phrasebits/look->frames,
67 (float)look->postbits/look->frames,
68 (float)(look->postbits+look->phrasebits)/look->frames);*/
69
70 memset(look,0,sizeof(*look));
71 _ogg_free(look);
72 }
73 }
74
floor1_pack(vorbis_info_floor * i,oggpack_buffer * opb)75 static void floor1_pack (vorbis_info_floor *i,oggpack_buffer *opb){
76 vorbis_info_floor1 *info=(vorbis_info_floor1 *)i;
77 int j,k;
78 int count=0;
79 int rangebits;
80 int maxposit=info->postlist[1];
81 int maxclass=-1;
82
83 /* save out partitions */
84 oggpack_write(opb,info->partitions,5); /* only 0 to 31 legal */
85 for(j=0;j<info->partitions;j++){
86 oggpack_write(opb,info->partitionclass[j],4); /* only 0 to 15 legal */
87 if(maxclass<info->partitionclass[j])maxclass=info->partitionclass[j];
88 }
89
90 /* save out partition classes */
91 for(j=0;j<maxclass+1;j++){
92 oggpack_write(opb,info->class_dim[j]-1,3); /* 1 to 8 */
93 oggpack_write(opb,info->class_subs[j],2); /* 0 to 3 */
94 if(info->class_subs[j])oggpack_write(opb,info->class_book[j],8);
95 for(k=0;k<(1<<info->class_subs[j]);k++)
96 oggpack_write(opb,info->class_subbook[j][k]+1,8);
97 }
98
99 /* save out the post list */
100 oggpack_write(opb,info->mult-1,2); /* only 1,2,3,4 legal now */
101 /* maxposit cannot legally be less than 1; this is encode-side, we
102 can assume our setup is OK */
103 oggpack_write(opb,ov_ilog(maxposit-1),4);
104 rangebits=ov_ilog(maxposit-1);
105
106 for(j=0,k=0;j<info->partitions;j++){
107 count+=info->class_dim[info->partitionclass[j]];
108 for(;k<count;k++)
109 oggpack_write(opb,info->postlist[k+2],rangebits);
110 }
111 }
112
icomp(const void * a,const void * b)113 static int icomp(const void *a,const void *b){
114 return(**(int **)a-**(int **)b);
115 }
116
floor1_unpack(vorbis_info * vi,oggpack_buffer * opb)117 static vorbis_info_floor *floor1_unpack (vorbis_info *vi,oggpack_buffer *opb){
118 codec_setup_info *ci=vi->codec_setup;
119 int j,k,count=0,maxclass=-1,rangebits;
120
121 vorbis_info_floor1 *info=_ogg_calloc(1,sizeof(*info));
122 /* read partitions */
123 info->partitions=oggpack_read(opb,5); /* only 0 to 31 legal */
124 for(j=0;j<info->partitions;j++){
125 info->partitionclass[j]=oggpack_read(opb,4); /* only 0 to 15 legal */
126 if(info->partitionclass[j]<0)goto err_out;
127 if(maxclass<info->partitionclass[j])maxclass=info->partitionclass[j];
128 }
129
130 /* read partition classes */
131 for(j=0;j<maxclass+1;j++){
132 info->class_dim[j]=oggpack_read(opb,3)+1; /* 1 to 8 */
133 info->class_subs[j]=oggpack_read(opb,2); /* 0,1,2,3 bits */
134 if(info->class_subs[j]<0)
135 goto err_out;
136 if(info->class_subs[j])info->class_book[j]=oggpack_read(opb,8);
137 if(info->class_book[j]<0 || info->class_book[j]>=ci->books)
138 goto err_out;
139 for(k=0;k<(1<<info->class_subs[j]);k++){
140 info->class_subbook[j][k]=oggpack_read(opb,8)-1;
141 if(info->class_subbook[j][k]<-1 || info->class_subbook[j][k]>=ci->books)
142 goto err_out;
143 }
144 }
145
146 /* read the post list */
147 info->mult=oggpack_read(opb,2)+1; /* only 1,2,3,4 legal now */
148 rangebits=oggpack_read(opb,4);
149 if(rangebits<0)goto err_out;
150
151 for(j=0,k=0;j<info->partitions;j++){
152 count+=info->class_dim[info->partitionclass[j]];
153 if(count>VIF_POSIT) goto err_out;
154 for(;k<count;k++){
155 int t=info->postlist[k+2]=oggpack_read(opb,rangebits);
156 if(t<0 || t>=(1<<rangebits))
157 goto err_out;
158 }
159 }
160 info->postlist[0]=0;
161 info->postlist[1]=1<<rangebits;
162
163 /* don't allow repeated values in post list as they'd result in
164 zero-length segments */
165 {
166 int *sortpointer[VIF_POSIT+2];
167 for(j=0;j<count+2;j++)sortpointer[j]=info->postlist+j;
168 qsort(sortpointer,count+2,sizeof(*sortpointer),icomp);
169
170 for(j=1;j<count+2;j++)
171 if(*sortpointer[j-1]==*sortpointer[j])goto err_out;
172 }
173
174 return(info);
175
176 err_out:
177 floor1_free_info(info);
178 return(NULL);
179 }
180
floor1_look(vorbis_dsp_state * vd,vorbis_info_floor * in)181 static vorbis_look_floor *floor1_look(vorbis_dsp_state *vd,
182 vorbis_info_floor *in){
183
184 int *sortpointer[VIF_POSIT+2];
185 vorbis_info_floor1 *info=(vorbis_info_floor1 *)in;
186 vorbis_look_floor1 *look=_ogg_calloc(1,sizeof(*look));
187 int i,j,n=0;
188
189 (void)vd;
190
191 look->vi=info;
192 look->n=info->postlist[1];
193
194 /* we drop each position value in-between already decoded values,
195 and use linear interpolation to predict each new value past the
196 edges. The positions are read in the order of the position
197 list... we precompute the bounding positions in the lookup. Of
198 course, the neighbors can change (if a position is declined), but
199 this is an initial mapping */
200
201 for(i=0;i<info->partitions;i++)n+=info->class_dim[info->partitionclass[i]];
202 n+=2;
203 look->posts=n;
204
205 /* also store a sorted position index */
206 for(i=0;i<n;i++)sortpointer[i]=info->postlist+i;
207 qsort(sortpointer,n,sizeof(*sortpointer),icomp);
208
209 /* points from sort order back to range number */
210 for(i=0;i<n;i++)look->forward_index[i]=sortpointer[i]-info->postlist;
211 /* points from range order to sorted position */
212 for(i=0;i<n;i++)look->reverse_index[look->forward_index[i]]=i;
213 /* we actually need the post values too */
214 for(i=0;i<n;i++)look->sorted_index[i]=info->postlist[look->forward_index[i]];
215
216 /* quantize values to multiplier spec */
217 switch(info->mult){
218 case 1: /* 1024 -> 256 */
219 look->quant_q=256;
220 break;
221 case 2: /* 1024 -> 128 */
222 look->quant_q=128;
223 break;
224 case 3: /* 1024 -> 86 */
225 look->quant_q=86;
226 break;
227 case 4: /* 1024 -> 64 */
228 look->quant_q=64;
229 break;
230 }
231
232 /* discover our neighbors for decode where we don't use fit flags
233 (that would push the neighbors outward) */
234 for(i=0;i<n-2;i++){
235 int lo=0;
236 int hi=1;
237 int lx=0;
238 int hx=look->n;
239 int currentx=info->postlist[i+2];
240 for(j=0;j<i+2;j++){
241 int x=info->postlist[j];
242 if(x>lx && x<currentx){
243 lo=j;
244 lx=x;
245 }
246 if(x<hx && x>currentx){
247 hi=j;
248 hx=x;
249 }
250 }
251 look->loneighbor[i]=lo;
252 look->hineighbor[i]=hi;
253 }
254
255 return(look);
256 }
257
render_point(int x0,int x1,int y0,int y1,int x)258 static int render_point(int x0,int x1,int y0,int y1,int x){
259 y0&=0x7fff; /* mask off flag */
260 y1&=0x7fff;
261
262 {
263 int dy=y1-y0;
264 int adx=x1-x0;
265 int ady=abs(dy);
266 int err=ady*(x-x0);
267
268 int off=err/adx;
269 if(dy<0)return(y0-off);
270 return(y0+off);
271 }
272 }
273
vorbis_dBquant(const float * x)274 static int vorbis_dBquant(const float *x){
275 int i= *x*7.3142857f+1023.5f;
276 if(i>1023)return(1023);
277 if(i<0)return(0);
278 return i;
279 }
280
281 static const float FLOOR1_fromdB_LOOKUP[256]={
282 1.0649863e-07F, 1.1341951e-07F, 1.2079015e-07F, 1.2863978e-07F,
283 1.3699951e-07F, 1.4590251e-07F, 1.5538408e-07F, 1.6548181e-07F,
284 1.7623575e-07F, 1.8768855e-07F, 1.9988561e-07F, 2.128753e-07F,
285 2.2670913e-07F, 2.4144197e-07F, 2.5713223e-07F, 2.7384213e-07F,
286 2.9163793e-07F, 3.1059021e-07F, 3.3077411e-07F, 3.5226968e-07F,
287 3.7516214e-07F, 3.9954229e-07F, 4.2550680e-07F, 4.5315863e-07F,
288 4.8260743e-07F, 5.1396998e-07F, 5.4737065e-07F, 5.8294187e-07F,
289 6.2082472e-07F, 6.6116941e-07F, 7.0413592e-07F, 7.4989464e-07F,
290 7.9862701e-07F, 8.5052630e-07F, 9.0579828e-07F, 9.6466216e-07F,
291 1.0273513e-06F, 1.0941144e-06F, 1.1652161e-06F, 1.2409384e-06F,
292 1.3215816e-06F, 1.4074654e-06F, 1.4989305e-06F, 1.5963394e-06F,
293 1.7000785e-06F, 1.8105592e-06F, 1.9282195e-06F, 2.0535261e-06F,
294 2.1869758e-06F, 2.3290978e-06F, 2.4804557e-06F, 2.6416497e-06F,
295 2.8133190e-06F, 2.9961443e-06F, 3.1908506e-06F, 3.3982101e-06F,
296 3.6190449e-06F, 3.8542308e-06F, 4.1047004e-06F, 4.3714470e-06F,
297 4.6555282e-06F, 4.9580707e-06F, 5.2802740e-06F, 5.6234160e-06F,
298 5.9888572e-06F, 6.3780469e-06F, 6.7925283e-06F, 7.2339451e-06F,
299 7.7040476e-06F, 8.2047000e-06F, 8.7378876e-06F, 9.3057248e-06F,
300 9.9104632e-06F, 1.0554501e-05F, 1.1240392e-05F, 1.1970856e-05F,
301 1.2748789e-05F, 1.3577278e-05F, 1.4459606e-05F, 1.5399272e-05F,
302 1.6400004e-05F, 1.7465768e-05F, 1.8600792e-05F, 1.9809576e-05F,
303 2.1096914e-05F, 2.2467911e-05F, 2.3928002e-05F, 2.5482978e-05F,
304 2.7139006e-05F, 2.8902651e-05F, 3.0780908e-05F, 3.2781225e-05F,
305 3.4911534e-05F, 3.7180282e-05F, 3.9596466e-05F, 4.2169667e-05F,
306 4.4910090e-05F, 4.7828601e-05F, 5.0936773e-05F, 5.4246931e-05F,
307 5.7772202e-05F, 6.1526565e-05F, 6.5524908e-05F, 6.9783085e-05F,
308 7.4317983e-05F, 7.9147585e-05F, 8.4291040e-05F, 8.9768747e-05F,
309 9.5602426e-05F, 0.00010181521F, 0.00010843174F, 0.00011547824F,
310 0.00012298267F, 0.00013097477F, 0.00013948625F, 0.00014855085F,
311 0.00015820453F, 0.00016848555F, 0.00017943469F, 0.00019109536F,
312 0.00020351382F, 0.00021673929F, 0.00023082423F, 0.00024582449F,
313 0.00026179955F, 0.00027881276F, 0.00029693158F, 0.00031622787F,
314 0.00033677814F, 0.00035866388F, 0.00038197188F, 0.00040679456F,
315 0.00043323036F, 0.00046138411F, 0.00049136745F, 0.00052329927F,
316 0.00055730621F, 0.00059352311F, 0.00063209358F, 0.00067317058F,
317 0.00071691700F, 0.00076350630F, 0.00081312324F, 0.00086596457F,
318 0.00092223983F, 0.00098217216F, 0.0010459992F, 0.0011139742F,
319 0.0011863665F, 0.0012634633F, 0.0013455702F, 0.0014330129F,
320 0.0015261382F, 0.0016253153F, 0.0017309374F, 0.0018434235F,
321 0.0019632195F, 0.0020908006F, 0.0022266726F, 0.0023713743F,
322 0.0025254795F, 0.0026895994F, 0.0028643847F, 0.0030505286F,
323 0.0032487691F, 0.0034598925F, 0.0036847358F, 0.0039241906F,
324 0.0041792066F, 0.0044507950F, 0.0047400328F, 0.0050480668F,
325 0.0053761186F, 0.0057254891F, 0.0060975636F, 0.0064938176F,
326 0.0069158225F, 0.0073652516F, 0.0078438871F, 0.0083536271F,
327 0.0088964928F, 0.009474637F, 0.010090352F, 0.010746080F,
328 0.011444421F, 0.012188144F, 0.012980198F, 0.013823725F,
329 0.014722068F, 0.015678791F, 0.016697687F, 0.017782797F,
330 0.018938423F, 0.020169149F, 0.021479854F, 0.022875735F,
331 0.024362330F, 0.025945531F, 0.027631618F, 0.029427276F,
332 0.031339626F, 0.033376252F, 0.035545228F, 0.037855157F,
333 0.040315199F, 0.042935108F, 0.045725273F, 0.048696758F,
334 0.051861348F, 0.055231591F, 0.058820850F, 0.062643361F,
335 0.066714279F, 0.071049749F, 0.075666962F, 0.080584227F,
336 0.085821044F, 0.091398179F, 0.097337747F, 0.10366330F,
337 0.11039993F, 0.11757434F, 0.12521498F, 0.13335215F,
338 0.14201813F, 0.15124727F, 0.16107617F, 0.17154380F,
339 0.18269168F, 0.19456402F, 0.20720788F, 0.22067342F,
340 0.23501402F, 0.25028656F, 0.26655159F, 0.28387361F,
341 0.30232132F, 0.32196786F, 0.34289114F, 0.36517414F,
342 0.38890521F, 0.41417847F, 0.44109412F, 0.46975890F,
343 0.50028648F, 0.53279791F, 0.56742212F, 0.60429640F,
344 0.64356699F, 0.68538959F, 0.72993007F, 0.77736504F,
345 0.82788260F, 0.88168307F, 0.9389798F, 1.F,
346 };
347
render_line(int n,int x0,int x1,int y0,int y1,float * d)348 static void render_line(int n, int x0,int x1,int y0,int y1,float *d){
349 int dy=y1-y0;
350 int adx=x1-x0;
351 int ady=abs(dy);
352 int base=dy/adx;
353 int sy=(dy<0?base-1:base+1);
354 int x=x0;
355 int y=y0;
356 int err=0;
357
358 ady-=abs(base*adx);
359
360 if(n>x1)n=x1;
361
362 if(x<n)
363 d[x]*=FLOOR1_fromdB_LOOKUP[y];
364
365 while(++x<n){
366 err=err+ady;
367 if(err>=adx){
368 err-=adx;
369 y+=sy;
370 }else{
371 y+=base;
372 }
373 d[x]*=FLOOR1_fromdB_LOOKUP[y];
374 }
375 }
376
render_line0(int n,int x0,int x1,int y0,int y1,int * d)377 static void render_line0(int n, int x0,int x1,int y0,int y1,int *d){
378 int dy=y1-y0;
379 int adx=x1-x0;
380 int ady=abs(dy);
381 int base=dy/adx;
382 int sy=(dy<0?base-1:base+1);
383 int x=x0;
384 int y=y0;
385 int err=0;
386
387 ady-=abs(base*adx);
388
389 if(n>x1)n=x1;
390
391 if(x<n)
392 d[x]=y;
393
394 while(++x<n){
395 err=err+ady;
396 if(err>=adx){
397 err-=adx;
398 y+=sy;
399 }else{
400 y+=base;
401 }
402 d[x]=y;
403 }
404 }
405
406 /* the floor has already been filtered to only include relevant sections */
accumulate_fit(const float * flr,const float * mdct,int x0,int x1,lsfit_acc * a,int n,vorbis_info_floor1 * info)407 static int accumulate_fit(const float *flr,const float *mdct,
408 int x0, int x1,lsfit_acc *a,
409 int n,vorbis_info_floor1 *info){
410 long i;
411
412 int xa=0,ya=0,x2a=0,y2a=0,xya=0,na=0, xb=0,yb=0,x2b=0,y2b=0,xyb=0,nb=0;
413
414 memset(a,0,sizeof(*a));
415 a->x0=x0;
416 a->x1=x1;
417 if(x1>=n)x1=n-1;
418
419 for(i=x0;i<=x1;i++){
420 int quantized=vorbis_dBquant(flr+i);
421 if(quantized){
422 if(mdct[i]+info->twofitatten>=flr[i]){
423 xa += i;
424 ya += quantized;
425 x2a += i*i;
426 y2a += quantized*quantized;
427 xya += i*quantized;
428 na++;
429 }else{
430 xb += i;
431 yb += quantized;
432 x2b += i*i;
433 y2b += quantized*quantized;
434 xyb += i*quantized;
435 nb++;
436 }
437 }
438 }
439
440 a->xa=xa;
441 a->ya=ya;
442 a->x2a=x2a;
443 a->y2a=y2a;
444 a->xya=xya;
445 a->an=na;
446
447 a->xb=xb;
448 a->yb=yb;
449 a->x2b=x2b;
450 a->y2b=y2b;
451 a->xyb=xyb;
452 a->bn=nb;
453
454 return(na);
455 }
456
fit_line(lsfit_acc * a,int fits,int * y0,int * y1,vorbis_info_floor1 * info)457 static int fit_line(lsfit_acc *a,int fits,int *y0,int *y1,
458 vorbis_info_floor1 *info){
459 double xb=0,yb=0,x2b=0,y2b=0,xyb=0,bn=0;
460 int i;
461 int x0=a[0].x0;
462 int x1=a[fits-1].x1;
463
464 for(i=0;i<fits;i++){
465 double weight = (a[i].bn+a[i].an)*info->twofitweight/(a[i].an+1)+1.;
466
467 xb+=a[i].xb + a[i].xa * weight;
468 yb+=a[i].yb + a[i].ya * weight;
469 x2b+=a[i].x2b + a[i].x2a * weight;
470 y2b+=a[i].y2b + a[i].y2a * weight;
471 xyb+=a[i].xyb + a[i].xya * weight;
472 bn+=a[i].bn + a[i].an * weight;
473 }
474
475 if(*y0>=0){
476 xb+= x0;
477 yb+= *y0;
478 x2b+= x0 * x0;
479 y2b+= *y0 * *y0;
480 xyb+= *y0 * x0;
481 bn++;
482 }
483
484 if(*y1>=0){
485 xb+= x1;
486 yb+= *y1;
487 x2b+= x1 * x1;
488 y2b+= *y1 * *y1;
489 xyb+= *y1 * x1;
490 bn++;
491 }
492
493 {
494 double denom=(bn*x2b-xb*xb);
495
496 if(denom>0.){
497 double a=(yb*x2b-xyb*xb)/denom;
498 double b=(bn*xyb-xb*yb)/denom;
499 *y0=rint(a+b*x0);
500 *y1=rint(a+b*x1);
501
502 /* limit to our range! */
503 if(*y0>1023)*y0=1023;
504 if(*y1>1023)*y1=1023;
505 if(*y0<0)*y0=0;
506 if(*y1<0)*y1=0;
507
508 return 0;
509 }else{
510 *y0=0;
511 *y1=0;
512 return 1;
513 }
514 }
515 }
516
inspect_error(int x0,int x1,int y0,int y1,const float * mask,const float * mdct,vorbis_info_floor1 * info)517 static int inspect_error(int x0,int x1,int y0,int y1,const float *mask,
518 const float *mdct,
519 vorbis_info_floor1 *info){
520 int dy=y1-y0;
521 int adx=x1-x0;
522 int ady=abs(dy);
523 int base=dy/adx;
524 int sy=(dy<0?base-1:base+1);
525 int x=x0;
526 int y=y0;
527 int err=0;
528 int val=vorbis_dBquant(mask+x);
529 int mse=0;
530 int n=0;
531
532 ady-=abs(base*adx);
533
534 mse=(y-val);
535 mse*=mse;
536 n++;
537 if(mdct[x]+info->twofitatten>=mask[x]){
538 if(y+info->maxover<val)return(1);
539 if(y-info->maxunder>val)return(1);
540 }
541
542 while(++x<x1){
543 err=err+ady;
544 if(err>=adx){
545 err-=adx;
546 y+=sy;
547 }else{
548 y+=base;
549 }
550
551 val=vorbis_dBquant(mask+x);
552 mse+=((y-val)*(y-val));
553 n++;
554 if(mdct[x]+info->twofitatten>=mask[x]){
555 if(val){
556 if(y+info->maxover<val)return(1);
557 if(y-info->maxunder>val)return(1);
558 }
559 }
560 }
561
562 if(info->maxover*info->maxover/n>info->maxerr)return(0);
563 if(info->maxunder*info->maxunder/n>info->maxerr)return(0);
564 if(mse/n>info->maxerr)return(1);
565 return(0);
566 }
567
post_Y(int * A,int * B,int pos)568 static int post_Y(int *A,int *B,int pos){
569 if(A[pos]<0)
570 return B[pos];
571 if(B[pos]<0)
572 return A[pos];
573
574 return (A[pos]+B[pos])>>1;
575 }
576
floor1_fit(vorbis_block * vb,vorbis_look_floor1 * look,const float * logmdct,const float * logmask)577 int *floor1_fit(vorbis_block *vb,vorbis_look_floor1 *look,
578 const float *logmdct, /* in */
579 const float *logmask){
580 long i,j;
581 vorbis_info_floor1 *info=look->vi;
582 long n=look->n;
583 long posts=look->posts;
584 long nonzero=0;
585 lsfit_acc fits[VIF_POSIT+1];
586 int fit_valueA[VIF_POSIT+2]; /* index by range list position */
587 int fit_valueB[VIF_POSIT+2]; /* index by range list position */
588
589 int loneighbor[VIF_POSIT+2]; /* sorted index of range list position (+2) */
590 int hineighbor[VIF_POSIT+2];
591 int *output=NULL;
592 int memo[VIF_POSIT+2];
593
594 for(i=0;i<posts;i++)fit_valueA[i]=-200; /* mark all unused */
595 for(i=0;i<posts;i++)fit_valueB[i]=-200; /* mark all unused */
596 for(i=0;i<posts;i++)loneighbor[i]=0; /* 0 for the implicit 0 post */
597 for(i=0;i<posts;i++)hineighbor[i]=1; /* 1 for the implicit post at n */
598 for(i=0;i<posts;i++)memo[i]=-1; /* no neighbor yet */
599
600 /* quantize the relevant floor points and collect them into line fit
601 structures (one per minimal division) at the same time */
602 if(posts==0){
603 nonzero+=accumulate_fit(logmask,logmdct,0,n,fits,n,info);
604 }else{
605 for(i=0;i<posts-1;i++)
606 nonzero+=accumulate_fit(logmask,logmdct,look->sorted_index[i],
607 look->sorted_index[i+1],fits+i,
608 n,info);
609 }
610
611 if(nonzero){
612 /* start by fitting the implicit base case.... */
613 int y0=-200;
614 int y1=-200;
615 fit_line(fits,posts-1,&y0,&y1,info);
616
617 fit_valueA[0]=y0;
618 fit_valueB[0]=y0;
619 fit_valueB[1]=y1;
620 fit_valueA[1]=y1;
621
622 /* Non degenerate case */
623 /* start progressive splitting. This is a greedy, non-optimal
624 algorithm, but simple and close enough to the best
625 answer. */
626 for(i=2;i<posts;i++){
627 int sortpos=look->reverse_index[i];
628 int ln=loneighbor[sortpos];
629 int hn=hineighbor[sortpos];
630
631 /* eliminate repeat searches of a particular range with a memo */
632 if(memo[ln]!=hn){
633 /* haven't performed this error search yet */
634 int lsortpos=look->reverse_index[ln];
635 int hsortpos=look->reverse_index[hn];
636 memo[ln]=hn;
637
638 {
639 /* A note: we want to bound/minimize *local*, not global, error */
640 int lx=info->postlist[ln];
641 int hx=info->postlist[hn];
642 int ly=post_Y(fit_valueA,fit_valueB,ln);
643 int hy=post_Y(fit_valueA,fit_valueB,hn);
644
645 if(ly==-1 || hy==-1){
646 exit(1);
647 }
648
649 if(inspect_error(lx,hx,ly,hy,logmask,logmdct,info)){
650 /* outside error bounds/begin search area. Split it. */
651 int ly0=-200;
652 int ly1=-200;
653 int hy0=-200;
654 int hy1=-200;
655 int ret0=fit_line(fits+lsortpos,sortpos-lsortpos,&ly0,&ly1,info);
656 int ret1=fit_line(fits+sortpos,hsortpos-sortpos,&hy0,&hy1,info);
657
658 if(ret0){
659 ly0=ly;
660 ly1=hy0;
661 }
662 if(ret1){
663 hy0=ly1;
664 hy1=hy;
665 }
666
667 if(ret0 && ret1){
668 fit_valueA[i]=-200;
669 fit_valueB[i]=-200;
670 }else{
671 /* store new edge values */
672 fit_valueB[ln]=ly0;
673 if(ln==0)fit_valueA[ln]=ly0;
674 fit_valueA[i]=ly1;
675 fit_valueB[i]=hy0;
676 fit_valueA[hn]=hy1;
677 if(hn==1)fit_valueB[hn]=hy1;
678
679 if(ly1>=0 || hy0>=0){
680 /* store new neighbor values */
681 for(j=sortpos-1;j>=0;j--)
682 if(hineighbor[j]==hn)
683 hineighbor[j]=i;
684 else
685 break;
686 for(j=sortpos+1;j<posts;j++)
687 if(loneighbor[j]==ln)
688 loneighbor[j]=i;
689 else
690 break;
691 }
692 }
693 }else{
694 fit_valueA[i]=-200;
695 fit_valueB[i]=-200;
696 }
697 }
698 }
699 }
700
701 output=_vorbis_block_alloc(vb,sizeof(*output)*posts);
702
703 output[0]=post_Y(fit_valueA,fit_valueB,0);
704 output[1]=post_Y(fit_valueA,fit_valueB,1);
705
706 /* fill in posts marked as not using a fit; we will zero
707 back out to 'unused' when encoding them so long as curve
708 interpolation doesn't force them into use */
709 for(i=2;i<posts;i++){
710 int ln=look->loneighbor[i-2];
711 int hn=look->hineighbor[i-2];
712 int x0=info->postlist[ln];
713 int x1=info->postlist[hn];
714 int y0=output[ln];
715 int y1=output[hn];
716
717 int predicted=render_point(x0,x1,y0,y1,info->postlist[i]);
718 int vx=post_Y(fit_valueA,fit_valueB,i);
719
720 if(vx>=0 && predicted!=vx){
721 output[i]=vx;
722 }else{
723 output[i]= predicted|0x8000;
724 }
725 }
726 }
727
728 return(output);
729
730 }
731
floor1_interpolate_fit(vorbis_block * vb,vorbis_look_floor1 * look,int * A,int * B,int del)732 int *floor1_interpolate_fit(vorbis_block *vb,vorbis_look_floor1 *look,
733 int *A,int *B,
734 int del){
735
736 long i;
737 long posts=look->posts;
738 int *output=NULL;
739
740 if(A && B){
741 output=_vorbis_block_alloc(vb,sizeof(*output)*posts);
742
743 /* overly simpleminded--- look again post 1.2 */
744 for(i=0;i<posts;i++){
745 output[i]=((65536-del)*(A[i]&0x7fff)+del*(B[i]&0x7fff)+32768)>>16;
746 if(A[i]&0x8000 && B[i]&0x8000)output[i]|=0x8000;
747 }
748 }
749
750 return(output);
751 }
752
753
floor1_encode(oggpack_buffer * opb,vorbis_block * vb,vorbis_look_floor1 * look,int * post,int * ilogmask)754 int floor1_encode(oggpack_buffer *opb,vorbis_block *vb,
755 vorbis_look_floor1 *look,
756 int *post,int *ilogmask){
757
758 long i,j;
759 vorbis_info_floor1 *info=look->vi;
760 long posts=look->posts;
761 codec_setup_info *ci=vb->vd->vi->codec_setup;
762 int out[VIF_POSIT+2];
763 static_codebook **sbooks=ci->book_param;
764 codebook *books=ci->fullbooks;
765
766 /* quantize values to multiplier spec */
767 if(post){
768 for(i=0;i<posts;i++){
769 int val=post[i]&0x7fff;
770 switch(info->mult){
771 case 1: /* 1024 -> 256 */
772 val>>=2;
773 break;
774 case 2: /* 1024 -> 128 */
775 val>>=3;
776 break;
777 case 3: /* 1024 -> 86 */
778 val/=12;
779 break;
780 case 4: /* 1024 -> 64 */
781 val>>=4;
782 break;
783 }
784 post[i]=val | (post[i]&0x8000);
785 }
786
787 out[0]=post[0];
788 out[1]=post[1];
789
790 /* find prediction values for each post and subtract them */
791 for(i=2;i<posts;i++){
792 int ln=look->loneighbor[i-2];
793 int hn=look->hineighbor[i-2];
794 int x0=info->postlist[ln];
795 int x1=info->postlist[hn];
796 int y0=post[ln];
797 int y1=post[hn];
798
799 int predicted=render_point(x0,x1,y0,y1,info->postlist[i]);
800
801 if((post[i]&0x8000) || (predicted==post[i])){
802 post[i]=predicted|0x8000; /* in case there was roundoff jitter
803 in interpolation */
804 out[i]=0;
805 }else{
806 int headroom=(look->quant_q-predicted<predicted?
807 look->quant_q-predicted:predicted);
808
809 int val=post[i]-predicted;
810
811 /* at this point the 'deviation' value is in the range +/- max
812 range, but the real, unique range can always be mapped to
813 only [0-maxrange). So we want to wrap the deviation into
814 this limited range, but do it in the way that least screws
815 an essentially gaussian probability distribution. */
816
817 if(val<0)
818 if(val<-headroom)
819 val=headroom-val-1;
820 else
821 val=-1-(val<<1);
822 else
823 if(val>=headroom)
824 val= val+headroom;
825 else
826 val<<=1;
827
828 out[i]=val;
829 post[ln]&=0x7fff;
830 post[hn]&=0x7fff;
831 }
832 }
833
834 /* we have everything we need. pack it out */
835 /* mark nontrivial floor */
836 oggpack_write(opb,1,1);
837
838 /* beginning/end post */
839 look->frames++;
840 look->postbits+=ov_ilog(look->quant_q-1)*2;
841 oggpack_write(opb,out[0],ov_ilog(look->quant_q-1));
842 oggpack_write(opb,out[1],ov_ilog(look->quant_q-1));
843
844
845 /* partition by partition */
846 for(i=0,j=2;i<info->partitions;i++){
847 int class=info->partitionclass[i];
848 int cdim=info->class_dim[class];
849 int csubbits=info->class_subs[class];
850 int csub=1<<csubbits;
851 int bookas[8]={0,0,0,0,0,0,0,0};
852 int cval=0;
853 int cshift=0;
854 int k,l;
855
856 /* generate the partition's first stage cascade value */
857 if(csubbits){
858 int maxval[8]={0,0,0,0,0,0,0,0}; /* gcc's static analysis
859 issues a warning without
860 initialization */
861 for(k=0;k<csub;k++){
862 int booknum=info->class_subbook[class][k];
863 if(booknum<0){
864 maxval[k]=1;
865 }else{
866 maxval[k]=sbooks[info->class_subbook[class][k]]->entries;
867 }
868 }
869 for(k=0;k<cdim;k++){
870 for(l=0;l<csub;l++){
871 int val=out[j+k];
872 if(val<maxval[l]){
873 bookas[k]=l;
874 break;
875 }
876 }
877 cval|= bookas[k]<<cshift;
878 cshift+=csubbits;
879 }
880 /* write it */
881 look->phrasebits+=
882 vorbis_book_encode(books+info->class_book[class],cval,opb);
883
884 #ifdef TRAIN_FLOOR1
885 {
886 FILE *of;
887 char buffer[80];
888 sprintf(buffer,"line_%dx%ld_class%d.vqd",
889 vb->pcmend/2,posts-2,class);
890 of=fopen(buffer,"a");
891 fprintf(of,"%d\n",cval);
892 fclose(of);
893 }
894 #endif
895 }
896
897 /* write post values */
898 for(k=0;k<cdim;k++){
899 int book=info->class_subbook[class][bookas[k]];
900 if(book>=0){
901 /* hack to allow training with 'bad' books */
902 if(out[j+k]<(books+book)->entries)
903 look->postbits+=vorbis_book_encode(books+book,
904 out[j+k],opb);
905 /*else
906 fprintf(stderr,"+!");*/
907
908 #ifdef TRAIN_FLOOR1
909 {
910 FILE *of;
911 char buffer[80];
912 sprintf(buffer,"line_%dx%ld_%dsub%d.vqd",
913 vb->pcmend/2,posts-2,class,bookas[k]);
914 of=fopen(buffer,"a");
915 fprintf(of,"%d\n",out[j+k]);
916 fclose(of);
917 }
918 #endif
919 }
920 }
921 j+=cdim;
922 }
923
924 {
925 /* generate quantized floor equivalent to what we'd unpack in decode */
926 /* render the lines */
927 int hx=0;
928 int lx=0;
929 int ly=post[0]*info->mult;
930 int n=ci->blocksizes[vb->W]/2;
931
932 for(j=1;j<look->posts;j++){
933 int current=look->forward_index[j];
934 int hy=post[current]&0x7fff;
935 if(hy==post[current]){
936
937 hy*=info->mult;
938 hx=info->postlist[current];
939
940 render_line0(n,lx,hx,ly,hy,ilogmask);
941
942 lx=hx;
943 ly=hy;
944 }
945 }
946 for(j=hx;j<vb->pcmend/2;j++)ilogmask[j]=ly; /* be certain */
947 return(1);
948 }
949 }else{
950 oggpack_write(opb,0,1);
951 memset(ilogmask,0,vb->pcmend/2*sizeof(*ilogmask));
952 return(0);
953 }
954 }
955
floor1_inverse1(vorbis_block * vb,vorbis_look_floor * in)956 static void *floor1_inverse1(vorbis_block *vb,vorbis_look_floor *in){
957 vorbis_look_floor1 *look=(vorbis_look_floor1 *)in;
958 vorbis_info_floor1 *info=look->vi;
959 codec_setup_info *ci=vb->vd->vi->codec_setup;
960
961 int i,j,k;
962 codebook *books=ci->fullbooks;
963
964 /* unpack wrapped/predicted values from stream */
965 if(oggpack_read(&vb->opb,1)==1){
966 int *fit_value=_vorbis_block_alloc(vb,(look->posts)*sizeof(*fit_value));
967
968 fit_value[0]=oggpack_read(&vb->opb,ov_ilog(look->quant_q-1));
969 fit_value[1]=oggpack_read(&vb->opb,ov_ilog(look->quant_q-1));
970
971 /* partition by partition */
972 for(i=0,j=2;i<info->partitions;i++){
973 int class=info->partitionclass[i];
974 int cdim=info->class_dim[class];
975 int csubbits=info->class_subs[class];
976 int csub=1<<csubbits;
977 int cval=0;
978
979 /* decode the partition's first stage cascade value */
980 if(csubbits){
981 cval=vorbis_book_decode(books+info->class_book[class],&vb->opb);
982
983 if(cval==-1)goto eop;
984 }
985
986 for(k=0;k<cdim;k++){
987 int book=info->class_subbook[class][cval&(csub-1)];
988 cval>>=csubbits;
989 if(book>=0){
990 if((fit_value[j+k]=vorbis_book_decode(books+book,&vb->opb))==-1)
991 goto eop;
992 }else{
993 fit_value[j+k]=0;
994 }
995 }
996 j+=cdim;
997 }
998
999 /* unwrap positive values and reconsitute via linear interpolation */
1000 for(i=2;i<look->posts;i++){
1001 int predicted=render_point(info->postlist[look->loneighbor[i-2]],
1002 info->postlist[look->hineighbor[i-2]],
1003 fit_value[look->loneighbor[i-2]],
1004 fit_value[look->hineighbor[i-2]],
1005 info->postlist[i]);
1006 int hiroom=look->quant_q-predicted;
1007 int loroom=predicted;
1008 int room=(hiroom<loroom?hiroom:loroom)<<1;
1009 int val=fit_value[i];
1010
1011 if(val){
1012 if(val>=room){
1013 if(hiroom>loroom){
1014 val = val-loroom;
1015 }else{
1016 val = -1-(val-hiroom);
1017 }
1018 }else{
1019 if(val&1){
1020 val= -((val+1)>>1);
1021 }else{
1022 val>>=1;
1023 }
1024 }
1025
1026 fit_value[i]=(val+predicted)&0x7fff;
1027 fit_value[look->loneighbor[i-2]]&=0x7fff;
1028 fit_value[look->hineighbor[i-2]]&=0x7fff;
1029
1030 }else{
1031 fit_value[i]=predicted|0x8000;
1032 }
1033
1034 }
1035
1036 return(fit_value);
1037 }
1038 eop:
1039 return(NULL);
1040 }
1041
floor1_inverse2(vorbis_block * vb,vorbis_look_floor * in,void * memo,float * out)1042 static int floor1_inverse2(vorbis_block *vb,vorbis_look_floor *in,void *memo,
1043 float *out){
1044 vorbis_look_floor1 *look=(vorbis_look_floor1 *)in;
1045 vorbis_info_floor1 *info=look->vi;
1046
1047 codec_setup_info *ci=vb->vd->vi->codec_setup;
1048 int n=ci->blocksizes[vb->W]/2;
1049 int j;
1050
1051 if(memo){
1052 /* render the lines */
1053 int *fit_value=(int *)memo;
1054 int hx=0;
1055 int lx=0;
1056 int ly=fit_value[0]*info->mult;
1057 /* guard lookup against out-of-range values */
1058 ly=(ly<0?0:ly>255?255:ly);
1059
1060 for(j=1;j<look->posts;j++){
1061 int current=look->forward_index[j];
1062 int hy=fit_value[current]&0x7fff;
1063 if(hy==fit_value[current]){
1064
1065 hx=info->postlist[current];
1066 hy*=info->mult;
1067 /* guard lookup against out-of-range values */
1068 hy=(hy<0?0:hy>255?255:hy);
1069
1070 render_line(n,lx,hx,ly,hy,out);
1071
1072 lx=hx;
1073 ly=hy;
1074 }
1075 }
1076 for(j=hx;j<n;j++)out[j]*=FLOOR1_fromdB_LOOKUP[ly]; /* be certain */
1077 return(1);
1078 }
1079 memset(out,0,sizeof(*out)*n);
1080 return(0);
1081 }
1082
1083 /* export hooks */
1084 const vorbis_func_floor floor1_exportbundle={
1085 &floor1_pack,&floor1_unpack,&floor1_look,&floor1_free_info,
1086 &floor1_free_look,&floor1_inverse1,&floor1_inverse2
1087 };
1088