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