1 /****************************************************************************
2 Quantizer core functions
3 quality setting, error distribution, etc.
4
5 Copyright (C) 2017 Krzysztof Nikiel
6
7 This program is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 This program is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with this program. If not, see <http://www.gnu.org/licenses/>.
19 ****************************************************************************/
20
21 #include <math.h>
22 #include <stdio.h>
23 #include "quantize.h"
24 #include "huff2.h"
25
26 #ifdef HAVE_IMMINTRIN_H
27 # include <immintrin.h>
28 #endif
29
30 #ifdef __SSE2__
31 # ifdef __GNUC__
32 # include <cpuid.h>
33 # endif
34 #endif
35
36 #ifdef _MSC_VER
37 # include <immintrin.h>
38 # include <intrin.h>
39 # define __SSE2__
40 # define bit_SSE2 (1 << 26)
41 #endif
42
43 #ifdef __GNUC__
44 #define GCC_VERSION (__GNUC__ * 10000 \
45 + __GNUC_MINOR__ * 100 \
46 + __GNUC_PATCHLEVEL__)
47 #endif
48
49 #define MAGIC_NUMBER 0.4054
50 #define NOISEFLOOR 0.4
51
52 // band sound masking
bmask(CoderInfo * coderInfo,double * xr0,double * bandqual,int gnum,double quality)53 static void bmask(CoderInfo *coderInfo, double *xr0, double *bandqual,
54 int gnum, double quality)
55 {
56 int sfb, start, end, cnt;
57 int *cb_offset = coderInfo->sfb_offset;
58 int last;
59 double avgenrg;
60 double powm = 0.4;
61 double totenrg = 0.0;
62 int gsize = coderInfo->groups.len[gnum];
63 double *xr;
64 int win;
65 int enrgcnt = 0;
66
67
68 for (sfb = 0; sfb < coderInfo->sfbn; sfb++)
69 {
70 start = coderInfo->sfb_offset[sfb];
71 end = coderInfo->sfb_offset[sfb + 1];
72
73 xr = xr0;
74 for (win = 0; win < gsize; win++)
75 {
76 for (cnt = start; cnt < end; cnt++)
77 {
78 totenrg += xr[cnt] * xr[cnt];
79 enrgcnt++;
80 }
81
82 xr += BLOCK_LEN_SHORT;
83 }
84 }
85
86 if (totenrg < ((NOISEFLOOR * NOISEFLOOR) * (double)enrgcnt))
87 {
88 for (sfb = 0; sfb < coderInfo->sfbn; sfb++)
89 bandqual[sfb] = 0.0;
90
91 return;
92 }
93
94 for (sfb = 0; sfb < coderInfo->sfbn; sfb++)
95 {
96 double avge, maxe;
97 double target;
98
99 start = cb_offset[sfb];
100 end = cb_offset[sfb + 1];
101
102 avge = 0.0;
103 maxe = 0.0;
104 xr = xr0;
105 for (win = 0; win < gsize; win++)
106 {
107 for (cnt = start; cnt < end; cnt++)
108 {
109 double e = xr[cnt]*xr[cnt];
110 avge += e;
111 if (maxe < e)
112 maxe = e;
113 }
114 xr += BLOCK_LEN_SHORT;
115 }
116 maxe *= gsize;
117
118 #define NOISETONE 0.2
119 if (coderInfo->block_type == ONLY_SHORT_WINDOW)
120 {
121 last = BLOCK_LEN_SHORT;
122 avgenrg = totenrg / last;
123 avgenrg *= end - start;
124
125 target = NOISETONE * pow(avge/avgenrg, powm);
126 target += (1.0 - NOISETONE) * 0.45 * pow(maxe/avgenrg, powm);
127
128 target *= 1.5;
129 }
130 else
131 {
132 last = BLOCK_LEN_LONG;
133 avgenrg = totenrg / last;
134 avgenrg *= end - start;
135
136 target = NOISETONE * pow(avge/avgenrg, powm);
137 target += (1.0 - NOISETONE) * 0.45 * pow(maxe/avgenrg, powm);
138 }
139
140 target *= 10.0 / (1.0 + ((double)(start+end)/last));
141
142 bandqual[sfb] = target * quality;
143 }
144 }
145
146 enum {MAXSHORTBAND = 36};
147 // use band quality levels to quantize a group of windows
qlevel(CoderInfo * coderInfo,const double * xr0,const double * bandqual,int gnum,int pnslevel)148 static void qlevel(CoderInfo *coderInfo,
149 const double *xr0,
150 const double *bandqual,
151 int gnum,
152 int pnslevel
153 )
154 {
155 int sb, cnt;
156 #if !defined(__clang__) && defined(__GNUC__) && (GCC_VERSION >= 40600)
157 /* 2^0.25 (1.50515 dB) step from AAC specs */
158 static const double sfstep = 1.0 / log10(sqrt(sqrt(2.0)));
159 #else
160 static const double sfstep = 20 / 1.50515;
161 #endif
162 int gsize = coderInfo->groups.len[gnum];
163 double pnsthr = 0.1 * pnslevel;
164 #ifdef __SSE2__
165 int cpuid[4];
166 int sse2 = 0;
167
168 cpuid[3] = 0;
169 # ifdef __GNUC__
170 __cpuid(1, cpuid[0], cpuid[1], cpuid[2], cpuid[3]);
171 # endif
172 # ifdef _MSC_VER
173 __cpuid(cpuid, 1);
174 # endif
175 if (cpuid[3] & bit_SSE2)
176 sse2 = 1;
177 #endif
178
179 for (sb = 0; sb < coderInfo->sfbn; sb++)
180 {
181 double sfacfix;
182 int sfac;
183 double rmsx;
184 double etot;
185 int xitab[8 * MAXSHORTBAND];
186 int *xi;
187 int start, end;
188 const double *xr;
189 int win;
190
191 if (coderInfo->book[coderInfo->bandcnt] != HCB_NONE)
192 {
193 coderInfo->bandcnt++;
194 continue;
195 }
196
197 start = coderInfo->sfb_offset[sb];
198 end = coderInfo->sfb_offset[sb+1];
199
200 etot = 0.0;
201 xr = xr0;
202 for (win = 0; win < gsize; win++)
203 {
204 for (cnt = start; cnt < end; cnt++)
205 {
206 double e = xr[cnt] * xr[cnt];
207 etot += e;
208 }
209 xr += BLOCK_LEN_SHORT;
210 }
211 etot /= (double)gsize;
212 rmsx = sqrt(etot / (end - start));
213
214 if ((rmsx < NOISEFLOOR) || (!bandqual[sb]))
215 {
216 coderInfo->book[coderInfo->bandcnt++] = HCB_ZERO;
217 continue;
218 }
219
220 #ifndef DRM
221 if (bandqual[sb] < pnsthr)
222 {
223 coderInfo->book[coderInfo->bandcnt] = HCB_PNS;
224 coderInfo->sf[coderInfo->bandcnt] +=
225 lrint(log10(etot) * (0.5 * sfstep));
226 coderInfo->bandcnt++;
227 continue;
228 }
229 #endif
230
231 sfac = lrint(log10(bandqual[sb] / rmsx) * sfstep);
232 if ((SF_OFFSET - sfac) < 10)
233 sfacfix = 0.0;
234 else
235 sfacfix = pow(10, sfac / sfstep);
236
237 xr = xr0 + start;
238 end -= start;
239 xi = xitab;
240 for (win = 0; win < gsize; win++)
241 {
242 #ifdef __SSE2__
243 if (sse2)
244 {
245 for (cnt = 0; cnt < end; cnt += 4)
246 {
247 __m128 x = {xr[cnt], xr[cnt + 1], xr[cnt + 2], xr[cnt + 3]};
248
249 x = _mm_max_ps(x, _mm_sub_ps((__m128){0, 0, 0, 0}, x));
250 x = _mm_mul_ps(x, (__m128){sfacfix, sfacfix, sfacfix, sfacfix});
251 x = _mm_mul_ps(x, _mm_sqrt_ps(x));
252 x = _mm_sqrt_ps(x);
253 x = _mm_add_ps(x, (__m128){MAGIC_NUMBER, MAGIC_NUMBER, MAGIC_NUMBER, MAGIC_NUMBER});
254
255 *(__m128i*)(xi + cnt) = _mm_cvttps_epi32(x);
256 }
257 for (cnt = 0; cnt < end; cnt++)
258 {
259 if (xr[cnt] < 0)
260 xi[cnt] = -xi[cnt];
261 }
262 xi += cnt;
263 xr += BLOCK_LEN_SHORT;
264 continue;
265 }
266 #endif
267
268 for (cnt = 0; cnt < end; cnt++)
269 {
270 double tmp = fabs(xr[cnt]);
271
272 tmp *= sfacfix;
273 tmp = sqrt(tmp * sqrt(tmp));
274
275 xi[cnt] = (int)(tmp + MAGIC_NUMBER);
276 if (xr[cnt] < 0)
277 xi[cnt] = -xi[cnt];
278 }
279 xi += cnt;
280 xr += BLOCK_LEN_SHORT;
281 }
282 huffbook(coderInfo, xitab, gsize * end);
283 coderInfo->sf[coderInfo->bandcnt++] += SF_OFFSET - sfac;
284 }
285 }
286
BlocQuant(CoderInfo * coder,double * xr,AACQuantCfg * aacquantCfg)287 int BlocQuant(CoderInfo *coder, double *xr, AACQuantCfg *aacquantCfg)
288 {
289 double bandlvl[MAX_SCFAC_BANDS];
290 int cnt;
291 double *gxr;
292
293 coder->global_gain = 0;
294
295 coder->bandcnt = 0;
296 coder->datacnt = 0;
297 #ifdef DRM
298 coder->iLenReordSpData = 0; /* init length of reordered spectral data */
299 coder->iLenLongestCW = 0; /* init length of longest codeword */
300 coder->cur_cw = 0; /* init codeword counter */
301 #endif
302
303 {
304 int lastis;
305 int lastsf;
306
307 gxr = xr;
308 for (cnt = 0; cnt < coder->groups.n; cnt++)
309 {
310 bmask(coder, gxr, bandlvl, cnt,
311 (double)aacquantCfg->quality/DEFQUAL);
312 qlevel(coder, gxr, bandlvl, cnt, aacquantCfg->pnslevel);
313 gxr += coder->groups.len[cnt] * BLOCK_LEN_SHORT;
314 }
315
316 coder->global_gain = 0;
317 for (cnt = 0; cnt < coder->bandcnt; cnt++)
318 {
319 int book = coder->book[cnt];
320 if (!book)
321 continue;
322 if ((book != HCB_INTENSITY) && (book != HCB_INTENSITY2))
323 {
324 coder->global_gain = coder->sf[cnt];
325 break;
326 }
327 }
328
329 lastsf = coder->global_gain;
330 lastis = 0;
331 // fixme: move SF range check to quantizer
332 for (cnt = 0; cnt < coder->bandcnt; cnt++)
333 {
334 int book = coder->book[cnt];
335 if ((book == HCB_INTENSITY) || (book == HCB_INTENSITY2))
336 {
337 int diff = coder->sf[cnt] - lastis;
338
339 if (diff < -60)
340 diff = -60;
341 if (diff > 60)
342 diff = 60;
343
344 lastis += diff;
345 coder->sf[cnt] = lastis;
346 }
347 else if (book == HCB_ESC)
348 {
349 int diff = coder->sf[cnt] - lastsf;
350
351 if (diff < -60)
352 diff = -60;
353 if (diff > 60)
354 diff = 60;
355
356 lastsf += diff;
357 coder->sf[cnt] = lastsf;
358 }
359 }
360 return 1;
361 }
362 return 0;
363 }
364
CalcBW(unsigned * bw,int rate,SR_INFO * sr,AACQuantCfg * aacquantCfg)365 void CalcBW(unsigned *bw, int rate, SR_INFO *sr, AACQuantCfg *aacquantCfg)
366 {
367 // find max short frame band
368 int max = *bw * (BLOCK_LEN_SHORT << 1) / rate;
369 int cnt;
370 int l;
371
372 l = 0;
373 for (cnt = 0; cnt < sr->num_cb_short; cnt++)
374 {
375 if (l >= max)
376 break;
377 l += sr->cb_width_short[cnt];
378 }
379 aacquantCfg->max_cbs = cnt;
380 if (aacquantCfg->pnslevel)
381 *bw = (double)l * rate / (BLOCK_LEN_SHORT << 1);
382
383 // find max long frame band
384 max = *bw * (BLOCK_LEN_LONG << 1) / rate;
385 l = 0;
386 for (cnt = 0; cnt < sr->num_cb_long; cnt++)
387 {
388 if (l >= max)
389 break;
390 l += sr->cb_width_long[cnt];
391 }
392 aacquantCfg->max_cbl = cnt;
393 aacquantCfg->max_l = l;
394
395 *bw = (double)l * rate / (BLOCK_LEN_LONG << 1);
396 }
397
398 enum {MINSFB = 2};
399
calce(double * xr,int * bands,double e[NSFB_SHORT],int maxsfb,int maxl)400 static void calce(double *xr, int *bands, double e[NSFB_SHORT], int maxsfb,
401 int maxl)
402 {
403 int sfb;
404 int l;
405
406 // mute lines above cutoff freq
407 for (l = maxl; l < bands[maxsfb]; l++)
408 xr[l] = 0.0;
409
410 for (sfb = MINSFB; sfb < maxsfb; sfb++)
411 {
412 e[sfb] = 0;
413 for (l = bands[sfb]; l < bands[sfb + 1]; l++)
414 e[sfb] += xr[l] * xr[l];
415 }
416 }
417
resete(double min[NSFB_SHORT],double max[NSFB_SHORT],double e[NSFB_SHORT],int maxsfb)418 static void resete(double min[NSFB_SHORT], double max[NSFB_SHORT],
419 double e[NSFB_SHORT], int maxsfb)
420 {
421 int sfb;
422 for (sfb = MINSFB; sfb < maxsfb; sfb++)
423 min[sfb] = max[sfb] = e[sfb];
424 }
425
426 #define PRINTSTAT 0
427 #if PRINTSTAT
428 static int groups = 0;
429 static int frames = 0;
430 #endif
BlocGroup(double * xr,CoderInfo * coderInfo,AACQuantCfg * cfg)431 void BlocGroup(double *xr, CoderInfo *coderInfo, AACQuantCfg *cfg)
432 {
433 int win, sfb;
434 double e[NSFB_SHORT];
435 double min[NSFB_SHORT];
436 double max[NSFB_SHORT];
437 const double thr = 3.0;
438 int win0;
439 int fastmin;
440 int maxsfb, maxl;
441
442 if (coderInfo->block_type != ONLY_SHORT_WINDOW)
443 {
444 coderInfo->groups.n = 1;
445 coderInfo->groups.len[0] = 1;
446 return;
447 }
448
449 maxl = cfg->max_l / 8;
450 maxsfb = cfg->max_cbs;
451 fastmin = ((maxsfb - MINSFB) * 3) >> 2;
452
453 #ifdef DRM
454 coderInfo->groups.n = 1;
455 coderInfo->groups.len[0] = 8;
456 return;
457 #endif
458
459 #if PRINTSTAT
460 frames++;
461 #endif
462 calce(xr, coderInfo->sfb_offset, e, maxsfb, maxl);
463 resete(min, max, e, maxsfb);
464 win0 = 0;
465 coderInfo->groups.n = 0;
466 for (win = 1; win < MAX_SHORT_WINDOWS; win++)
467 {
468 int fast = 0;
469
470 calce(xr + win * BLOCK_LEN_SHORT, coderInfo->sfb_offset, e, maxsfb, maxl);
471 for (sfb = MINSFB; sfb < maxsfb; sfb++)
472 {
473 if (min[sfb] > e[sfb])
474 min[sfb] = e[sfb];
475 if (max[sfb] < e[sfb])
476 max[sfb] = e[sfb];
477
478 if (max[sfb] > thr * min[sfb])
479 fast++;
480 }
481 if (fast > fastmin)
482 {
483 coderInfo->groups.len[coderInfo->groups.n++] = win - win0;
484 win0 = win;
485 resete(min, max, e, maxsfb);
486 }
487 }
488 coderInfo->groups.len[coderInfo->groups.n++] = win - win0;
489 #if PRINTSTAT
490 groups += coderInfo->groups.n;
491 #endif
492 }
493
BlocStat(void)494 void BlocStat(void)
495 {
496 #if PRINTSTAT
497 printf("frames:%d; groups:%d; g/f:%f\n", frames, groups, (double)groups/frames);
498 #endif
499 }
500