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 | |
32 | typedef 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 | |
53 | static 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 | |
61 | static 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 | |
74 | static 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 | |
112 | static int icomp(const void *a,const void *b){ |
113 | return(**(int **)a-**(int **)b); |
114 | } |
115 | |
116 | static 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 | |
180 | static 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 | |
257 | static 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 | |
273 | static 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 | |
280 | static 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 | |
347 | static 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 | |
376 | static 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 */ |
406 | static 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 | |
456 | static 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 | |
516 | static 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 | |
567 | static 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 | |
576 | int *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 | |
731 | int *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 | |
753 | int 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 | |
955 | static 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 | |
1041 | static 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 */ |
1083 | const 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 | |