2 * quantize_pvt source file
4 * Copyright (c) 1999-2002 Takehiro Tominaga
5 * Copyright (c) 2000-2011 Robert Hegemann
6 * Copyright (c) 2001 Naoki Shibata
7 * Copyright (c) 2002-2005 Gabriel Bouvigne
9 * This library is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Library General Public
11 * License as published by the Free Software Foundation; either
12 * version 2 of the License, or (at your option) any later version.
14 * This library is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Library General Public License for more details.
19 * You should have received a copy of the GNU Library General Public
20 * License along with this library; if not, write to the
21 * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
22 * Boston, MA 02111-1307, USA.
25 /* $Id: quantize_pvt.c,v 1.169 2011/05/24 20:45:55 robert Exp $ */
35 #include "quantize_pvt.h"
36 #include "reservoir.h"
37 #include "lame-analysis.h"
41 #define NSATHSCALE 100 /* Assuming dynamic range=96dB, this value should be 92 */
44 The following table is used to implement the scalefactor
45 partitioning for MPEG2 as described in section
46 2.4.3.2 of the IS. The indexing corresponds to the
47 way the tables are presented in the IS:
49 [table_number][row_in_table][column of nr_of_sfb]
51 const int nr_of_sfb_block[6][3][4] = {
85 /* Table B.6: layer3 preemphasis */
86 const int pretab[SBMAX_l] = {
87 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
88 1, 1, 1, 1, 2, 2, 3, 3, 3, 2, 0
92 Here are MPEG1 Table B.8 and MPEG2 Table B.1
93 -- Layer III scalefactor bands.
94 Index into this using a method such as:
95 idx = fr_ps->header->sampling_frequency
96 + (fr_ps->header->version * 3)
100 const scalefac_struct sfBandIndex[9] = {
101 { /* Table B.2.b: 22.05 kHz */
102 {0, 6, 12, 18, 24, 30, 36, 44, 54, 66, 80, 96, 116, 140, 168, 200, 238, 284, 336, 396, 464,
104 {0, 4, 8, 12, 18, 24, 32, 42, 56, 74, 100, 132, 174, 192}
105 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
106 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
108 { /* Table B.2.c: 24 kHz */ /* docs: 332. mpg123(broken): 330 */
109 {0, 6, 12, 18, 24, 30, 36, 44, 54, 66, 80, 96, 114, 136, 162, 194, 232, 278, 332, 394, 464,
111 {0, 4, 8, 12, 18, 26, 36, 48, 62, 80, 104, 136, 180, 192}
112 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
113 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
115 { /* Table B.2.a: 16 kHz */
116 {0, 6, 12, 18, 24, 30, 36, 44, 54, 66, 80, 96, 116, 140, 168, 200, 238, 284, 336, 396, 464,
118 {0, 4, 8, 12, 18, 26, 36, 48, 62, 80, 104, 134, 174, 192}
119 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
120 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
122 { /* Table B.8.b: 44.1 kHz */
123 {0, 4, 8, 12, 16, 20, 24, 30, 36, 44, 52, 62, 74, 90, 110, 134, 162, 196, 238, 288, 342, 418,
125 {0, 4, 8, 12, 16, 22, 30, 40, 52, 66, 84, 106, 136, 192}
126 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
127 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
129 { /* Table B.8.c: 48 kHz */
130 {0, 4, 8, 12, 16, 20, 24, 30, 36, 42, 50, 60, 72, 88, 106, 128, 156, 190, 230, 276, 330, 384,
132 {0, 4, 8, 12, 16, 22, 28, 38, 50, 64, 80, 100, 126, 192}
133 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
134 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
136 { /* Table B.8.a: 32 kHz */
137 {0, 4, 8, 12, 16, 20, 24, 30, 36, 44, 54, 66, 82, 102, 126, 156, 194, 240, 296, 364, 448, 550,
139 {0, 4, 8, 12, 16, 22, 30, 42, 58, 78, 104, 138, 180, 192}
140 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
141 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
143 { /* MPEG-2.5 11.025 kHz */
144 {0, 6, 12, 18, 24, 30, 36, 44, 54, 66, 80, 96, 116, 140, 168, 200, 238, 284, 336, 396, 464,
146 {0 / 3, 12 / 3, 24 / 3, 36 / 3, 54 / 3, 78 / 3, 108 / 3, 144 / 3, 186 / 3, 240 / 3, 312 / 3,
147 402 / 3, 522 / 3, 576 / 3}
148 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
149 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
151 { /* MPEG-2.5 12 kHz */
152 {0, 6, 12, 18, 24, 30, 36, 44, 54, 66, 80, 96, 116, 140, 168, 200, 238, 284, 336, 396, 464,
154 {0 / 3, 12 / 3, 24 / 3, 36 / 3, 54 / 3, 78 / 3, 108 / 3, 144 / 3, 186 / 3, 240 / 3, 312 / 3,
155 402 / 3, 522 / 3, 576 / 3}
156 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
157 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
159 { /* MPEG-2.5 8 kHz */
160 {0, 12, 24, 36, 48, 60, 72, 88, 108, 132, 160, 192, 232, 280, 336, 400, 476, 566, 568, 570,
162 {0 / 3, 24 / 3, 48 / 3, 72 / 3, 108 / 3, 156 / 3, 216 / 3, 288 / 3, 372 / 3, 480 / 3, 486 / 3,
163 492 / 3, 498 / 3, 576 / 3}
164 , {0, 0, 0, 0, 0, 0, 0} /* sfb21 pseudo sub bands */
165 , {0, 0, 0, 0, 0, 0, 0} /* sfb12 pseudo sub bands */
171 FLOAT pow20[Q_MAX + Q_MAX2 + 1];
173 FLOAT pow43[PRECALC_SIZE];
174 /* initialized in first call to iteration_init */
175 #ifdef TAKEHIRO_IEEE754_HACK
176 FLOAT adj43asm[PRECALC_SIZE];
178 FLOAT adj43[PRECALC_SIZE];
182 compute the ATH for each scalefactor band
185 Input: 3.3kHz signal 32767 amplitude (3.3kHz is where ATH is smallest = -5db)
186 longblocks: sfb=12 en0/bw=-11db max_en0 = 1.3db
187 shortblocks: sfb=5 -9db 0db
189 Input: 1 1 1 1 1 1 1 -1 -1 -1 -1 -1 -1 -1 (repeated)
190 longblocks: amp=1 sfb=12 en0/bw=-103 db max_en0 = -92db
191 amp=32767 sfb=12 -12 db -1.4db
193 Input: 1 1 1 1 1 1 1 -1 -1 -1 -1 -1 -1 -1 (repeated)
194 shortblocks: amp=1 sfb=5 en0/bw= -99 -86
195 amp=32767 sfb=5 -9 db 4db
198 MAX energy of largest wave at 3.3kHz = 1db
199 AVE energy of largest wave at 3.3kHz = -11db
200 Let's take AVE: -11db = maximum signal in sfb=12.
201 Dynamic range of CD: 96db. Therefor energy of smallest audible wave
202 in sfb=12 = -11 - 96 = -107db = ATH at 3.3kHz.
204 ATH formula for this wave: -5db. To adjust to LAME scaling, we need
205 ATH = ATH_formula - 103 (db)
206 ATH = ATH * 2.5e-10 (ener)
211 ATHmdct(SessionConfig_t const *cfg, FLOAT f)
215 ath = ATHformula(cfg, f);
217 if (cfg->ATHfixpoint > 0) {
218 ath -= cfg->ATHfixpoint;
223 ath += cfg->ATH_offset_db;
225 /* modify the MDCT scaling for the ATH and convert to energy */
226 ath = powf(10.0f, ath * 0.1f);
231 compute_ath(lame_internal_flags const* gfc)
233 SessionConfig_t const *const cfg = &gfc->cfg;
234 FLOAT *const ATH_l = gfc->ATH->l;
235 FLOAT *const ATH_psfb21 = gfc->ATH->psfb21;
236 FLOAT *const ATH_s = gfc->ATH->s;
237 FLOAT *const ATH_psfb12 = gfc->ATH->psfb12;
238 int sfb, i, start, end;
240 FLOAT const samp_freq = cfg->samplerate_out;
242 for (sfb = 0; sfb < SBMAX_l; sfb++) {
243 start = gfc->scalefac_band.l[sfb];
244 end = gfc->scalefac_band.l[sfb + 1];
245 ATH_l[sfb] = FLOAT_MAX;
246 for (i = start; i < end; i++) {
247 FLOAT const freq = i * samp_freq / (2 * 576);
248 ATH_f = ATHmdct(cfg, freq); /* freq in kHz */
249 ATH_l[sfb] = Min(ATH_l[sfb], ATH_f);
253 for (sfb = 0; sfb < PSFB21; sfb++) {
254 start = gfc->scalefac_band.psfb21[sfb];
255 end = gfc->scalefac_band.psfb21[sfb + 1];
256 ATH_psfb21[sfb] = FLOAT_MAX;
257 for (i = start; i < end; i++) {
258 FLOAT const freq = i * samp_freq / (2 * 576);
259 ATH_f = ATHmdct(cfg, freq); /* freq in kHz */
260 ATH_psfb21[sfb] = Min(ATH_psfb21[sfb], ATH_f);
264 for (sfb = 0; sfb < SBMAX_s; sfb++) {
265 start = gfc->scalefac_band.s[sfb];
266 end = gfc->scalefac_band.s[sfb + 1];
267 ATH_s[sfb] = FLOAT_MAX;
268 for (i = start; i < end; i++) {
269 FLOAT const freq = i * samp_freq / (2 * 192);
270 ATH_f = ATHmdct(cfg, freq); /* freq in kHz */
271 ATH_s[sfb] = Min(ATH_s[sfb], ATH_f);
273 ATH_s[sfb] *= (gfc->scalefac_band.s[sfb + 1] - gfc->scalefac_band.s[sfb]);
276 for (sfb = 0; sfb < PSFB12; sfb++) {
277 start = gfc->scalefac_band.psfb12[sfb];
278 end = gfc->scalefac_band.psfb12[sfb + 1];
279 ATH_psfb12[sfb] = FLOAT_MAX;
280 for (i = start; i < end; i++) {
281 FLOAT const freq = i * samp_freq / (2 * 192);
282 ATH_f = ATHmdct(cfg, freq); /* freq in kHz */
283 ATH_psfb12[sfb] = Min(ATH_psfb12[sfb], ATH_f);
285 /*not sure about the following */
286 ATH_psfb12[sfb] *= (gfc->scalefac_band.s[13] - gfc->scalefac_band.s[12]);
291 * reduce ATH to -200 dB
295 for (sfb = 0; sfb < SBMAX_l; sfb++) {
298 for (sfb = 0; sfb < PSFB21; sfb++) {
299 ATH_psfb21[sfb] = 1E-20;
301 for (sfb = 0; sfb < SBMAX_s; sfb++) {
304 for (sfb = 0; sfb < PSFB12; sfb++) {
305 ATH_psfb12[sfb] = 1E-20;
309 /* work in progress, don't rely on it too much
311 gfc->ATH->floor = 10. * log10(ATHmdct(cfg, -1.));
314 { FLOAT g=10000, t=1e30, x;
315 for ( f = 100; f < 10000; f++ ) {
316 x = ATHmdct( cfg, f );
317 if ( t > x ) t = x, g = f;
319 printf("min=%g\n", g);
324 static float const payload_long[2][4] =
325 { {-0.000f, -0.000f, -0.000f, +0.000f}
326 , {-0.500f, -0.250f, -0.025f, +0.500f}
328 static float const payload_short[2][4] =
329 { {-0.000f, -0.000f, -0.000f, +0.000f}
330 , {-2.000f, -1.000f, -0.050f, +0.500f}
333 /************************************************************************/
334 /* initialization for iteration_loop */
335 /************************************************************************/
337 iteration_init(lame_internal_flags * gfc)
339 SessionConfig_t const *const cfg = &gfc->cfg;
340 III_side_info_t *const l3_side = &gfc->l3_side;
344 if (gfc->iteration_init_init == 0) {
345 gfc->iteration_init_init = 1;
347 l3_side->main_data_begin = 0;
351 for (i = 1; i < PRECALC_SIZE; i++)
352 pow43[i] = pow((FLOAT) i, 4.0 / 3.0);
354 #ifdef TAKEHIRO_IEEE754_HACK
356 for (i = 1; i < PRECALC_SIZE; i++)
357 adj43asm[i] = i - 0.5 - pow(0.5 * (pow43[i - 1] + pow43[i]), 0.75);
359 for (i = 0; i < PRECALC_SIZE - 1; i++)
360 adj43[i] = (i + 1) - pow(0.5 * (pow43[i] + pow43[i + 1]), 0.75);
363 for (i = 0; i < Q_MAX; i++)
364 ipow20[i] = pow(2.0, (double) (i - 210) * -0.1875);
365 for (i = 0; i <= Q_MAX + Q_MAX2; i++)
366 pow20[i] = pow(2.0, (double) (i - 210 - Q_MAX2) * 0.25);
369 init_xrpow_core_init(gfc);
371 sel = 1;/* RH: all modes like vbr-new (cfg->vbr == vbr_mt || cfg->vbr == vbr_mtrh) ? 1 : 0;*/
374 db = cfg->adjust_bass_db + payload_long[sel][0];
375 adjust = powf(10.f, db * 0.1f);
376 for (i = 0; i <= 6; ++i) {
377 gfc->sv_qnt.longfact[i] = adjust;
379 db = cfg->adjust_alto_db + payload_long[sel][1];
380 adjust = powf(10.f, db * 0.1f);
381 for (; i <= 13; ++i) {
382 gfc->sv_qnt.longfact[i] = adjust;
384 db = cfg->adjust_treble_db + payload_long[sel][2];
385 adjust = powf(10.f, db * 0.1f);
386 for (; i <= 20; ++i) {
387 gfc->sv_qnt.longfact[i] = adjust;
389 db = cfg->adjust_sfb21_db + payload_long[sel][3];
390 adjust = powf(10.f, db * 0.1f);
391 for (; i < SBMAX_l; ++i) {
392 gfc->sv_qnt.longfact[i] = adjust;
396 db = cfg->adjust_bass_db + payload_short[sel][0];
397 adjust = powf(10.f, db * 0.1f);
398 for (i = 0; i <= 2; ++i) {
399 gfc->sv_qnt.shortfact[i] = adjust;
401 db = cfg->adjust_alto_db + payload_short[sel][1];
402 adjust = powf(10.f, db * 0.1f);
403 for (; i <= 6; ++i) {
404 gfc->sv_qnt.shortfact[i] = adjust;
406 db = cfg->adjust_treble_db + payload_short[sel][2];
407 adjust = powf(10.f, db * 0.1f);
408 for (; i <= 11; ++i) {
409 gfc->sv_qnt.shortfact[i] = adjust;
411 db = cfg->adjust_sfb21_db + payload_short[sel][3];
412 adjust = powf(10.f, db * 0.1f);
413 for (; i < SBMAX_s; ++i) {
414 gfc->sv_qnt.shortfact[i] = adjust;
423 /************************************************************************
424 * allocate bits among 2 channels based on PE
426 * bugfixes rh 8/01: often allocated more than the allowed 4095 bits
427 ************************************************************************/
429 on_pe(lame_internal_flags * gfc, const FLOAT pe[][2], int targ_bits[2], int mean_bits, int gr, int cbr)
431 SessionConfig_t const *const cfg = &gfc->cfg;
432 int extra_bits = 0, tbits, bits;
433 int add_bits[2] = {0, 0};
434 int max_bits; /* maximum allowed bits for this granule */
437 /* allocate targ_bits for granule */
438 ResvMaxBits(gfc, mean_bits, &tbits, &extra_bits, cbr);
439 max_bits = tbits + extra_bits;
440 if (max_bits > MAX_BITS_PER_GRANULE) /* hard limit per granule */
441 max_bits = MAX_BITS_PER_GRANULE;
443 for (bits = 0, ch = 0; ch < cfg->channels_out; ++ch) {
444 /******************************************************************
445 * allocate bits for each channel
446 ******************************************************************/
447 targ_bits[ch] = Min(MAX_BITS_PER_CHANNEL, tbits / cfg->channels_out);
449 add_bits[ch] = targ_bits[ch] * pe[gr][ch] / 700.0 - targ_bits[ch];
451 /* at most increase bits by 1.5*average */
452 if (add_bits[ch] > mean_bits * 3 / 4)
453 add_bits[ch] = mean_bits * 3 / 4;
454 if (add_bits[ch] < 0)
457 if (add_bits[ch] + targ_bits[ch] > MAX_BITS_PER_CHANNEL)
458 add_bits[ch] = Max(0, MAX_BITS_PER_CHANNEL - targ_bits[ch]);
460 bits += add_bits[ch];
462 if (bits > extra_bits && bits > 0) {
463 for (ch = 0; ch < cfg->channels_out; ++ch) {
464 add_bits[ch] = extra_bits * add_bits[ch] / bits;
468 for (ch = 0; ch < cfg->channels_out; ++ch) {
469 targ_bits[ch] += add_bits[ch];
470 extra_bits -= add_bits[ch];
473 for (bits = 0, ch = 0; ch < cfg->channels_out; ++ch) {
474 bits += targ_bits[ch];
476 if (bits > MAX_BITS_PER_GRANULE) {
478 for (ch = 0; ch < cfg->channels_out; ++ch) {
479 targ_bits[ch] *= MAX_BITS_PER_GRANULE;
480 targ_bits[ch] /= bits;
481 sum += targ_bits[ch];
483 assert(sum <= MAX_BITS_PER_GRANULE);
493 reduce_side(int targ_bits[2], FLOAT ms_ener_ratio, int mean_bits, int max_bits)
498 assert(max_bits <= MAX_BITS_PER_GRANULE);
499 assert(targ_bits[0] + targ_bits[1] <= MAX_BITS_PER_GRANULE);
501 /* ms_ener_ratio = 0: allocate 66/33 mid/side fac=.33
502 * ms_ener_ratio =.5: allocate 50/50 mid/side fac= 0 */
503 /* 75/25 split is fac=.5 */
504 /* float fac = .50*(.5-ms_ener_ratio[gr])/.5; */
505 fac = .33 * (.5 - ms_ener_ratio) / .5;
511 /* number of bits to move from side channel to mid channel */
512 /* move_bits = fac*targ_bits[1]; */
513 move_bits = fac * .5 * (targ_bits[0] + targ_bits[1]);
515 if (move_bits > MAX_BITS_PER_CHANNEL - targ_bits[0]) {
516 move_bits = MAX_BITS_PER_CHANNEL - targ_bits[0];
521 if (targ_bits[1] >= 125) {
522 /* dont reduce side channel below 125 bits */
523 if (targ_bits[1] - move_bits > 125) {
525 /* if mid channel already has 2x more than average, dont bother */
526 /* mean_bits = bits per granule (for both channels) */
527 if (targ_bits[0] < mean_bits)
528 targ_bits[0] += move_bits;
529 targ_bits[1] -= move_bits;
532 targ_bits[0] += targ_bits[1] - 125;
537 move_bits = targ_bits[0] + targ_bits[1];
538 if (move_bits > max_bits) {
539 targ_bits[0] = (max_bits * targ_bits[0]) / move_bits;
540 targ_bits[1] = (max_bits * targ_bits[1]) / move_bits;
542 assert(targ_bits[0] <= MAX_BITS_PER_CHANNEL);
543 assert(targ_bits[1] <= MAX_BITS_PER_CHANNEL);
544 assert(targ_bits[0] + targ_bits[1] <= MAX_BITS_PER_GRANULE);
549 * Robert Hegemann 2001-04-27:
550 * this adjusts the ATH, keeping the original noise floor
551 * affects the higher frequencies more than the lower ones
555 athAdjust(FLOAT a, FLOAT x, FLOAT athFloor, float ATHfixpoint)
559 FLOAT const o = 90.30873362f;
560 FLOAT const p = (ATHfixpoint < 1.f) ? 94.82444863f : ATHfixpoint;
561 FLOAT u = FAST_LOG10_X(x, 10.0f);
562 FLOAT const v = a * a;
564 u -= athFloor; /* undo scaling */
566 w = 1.f + FAST_LOG10_X(v, 10.0f / o);
570 u += athFloor + o - p; /* redo scaling */
572 return powf(10.f, 0.1f * u);
577 /*************************************************************************/
579 /*************************************************************************/
582 Calculate the allowed distortion for each scalefactor band,
583 as determined by the psychoacoustic model.
584 xmin(sb) = ratio(sb) * en(sb) / bw(sb)
586 returns number of sfb's with energy > ATH
590 calc_xmin(lame_internal_flags const *gfc,
591 III_psy_ratio const *const ratio, gr_info * const cod_info, FLOAT * pxmin)
593 SessionConfig_t const *const cfg = &gfc->cfg;
594 int sfb, gsfb, j = 0, ath_over = 0, k;
595 ATH_t const *const ATH = gfc->ATH;
596 const FLOAT *const xr = cod_info->xr;
599 for (gsfb = 0; gsfb < cod_info->psy_lmax; gsfb++) {
604 xmin = athAdjust(ATH->adjust_factor, ATH->l[gsfb], ATH->floor, cfg->ATHfixpoint);
605 xmin *= gfc->sv_qnt.longfact[gsfb];
607 width = cod_info->width[gsfb];
612 rh2 = 2.2204460492503131e-016;
615 for (l = 0; l < width; ++l) {
616 FLOAT const xa = xr[j++];
617 FLOAT const x2 = xa * xa;
619 rh2 += (x2 < rh1) ? x2 : rh1;
627 else if (rh2 < xmin) {
635 FLOAT const e = ratio->en.l[gsfb];
638 x = en0 * ratio->thm.l[gsfb] / e;
639 x *= gfc->sv_qnt.longfact[gsfb];
645 } /* end of long block loop */
650 /*use this function to determine the highest non-zero coeff */
652 if (cod_info->block_type != SHORT_TYPE) { /* NORM, START or STOP type, but not SHORT */
654 while (k-- && fabs(xr[k]) < 1e-12f) {
658 cod_info->max_nonzero_coeff = max_nonzero;
662 for (sfb = cod_info->sfb_smin; gsfb < cod_info->psymax; sfb++, gsfb += 3) {
666 tmpATH = athAdjust(ATH->adjust_factor, ATH->s[sfb], ATH->floor, cfg->ATHfixpoint);
667 tmpATH *= gfc->sv_qnt.shortfact[sfb];
669 width = cod_info->width[gsfb];
670 for (b = 0; b < 3; b++) {
671 FLOAT en0 = 0.0, xmin = tmpATH;
674 rh1 = tmpATH / width;
678 rh2 = 2.2204460492503131e-016;
680 for (l = 0; l < width; ++l) {
681 FLOAT const xa = xr[j++];
682 FLOAT const x2 = xa * xa;
684 rh2 += (x2 < rh1) ? x2 : rh1;
692 else if (rh2 < tmpATH) {
700 FLOAT const e = ratio->en.s[sfb][b];
703 x = en0 * ratio->thm.s[sfb][b] / e;
704 x *= gfc->sv_qnt.shortfact[sfb];
711 if (cfg->use_temporal_masking_effect) {
712 if (pxmin[-3] > pxmin[-3 + 1])
713 pxmin[-3 + 1] += (pxmin[-3] - pxmin[-3 + 1]) * gfc->cd_psy->decay;
714 if (pxmin[-3 + 1] > pxmin[-3 + 2])
715 pxmin[-3 + 2] += (pxmin[-3 + 1] - pxmin[-3 + 2]) * gfc->cd_psy->decay;
717 } /* end of short block sfb loop */
724 calc_noise_core_c(const gr_info * const cod_info, int *startline, int l, FLOAT step)
728 const int *const ix = cod_info->l3_enc;
730 if (j > cod_info->count1) {
733 temp = cod_info->xr[j];
735 noise += temp * temp;
736 temp = cod_info->xr[j];
738 noise += temp * temp;
741 else if (j > cod_info->big_values) {
747 temp = fabs(cod_info->xr[j]) - ix01[ix[j]];
749 noise += temp * temp;
750 temp = fabs(cod_info->xr[j]) - ix01[ix[j]];
752 noise += temp * temp;
758 temp = fabs(cod_info->xr[j]) - pow43[ix[j]] * step;
760 noise += temp * temp;
761 temp = fabs(cod_info->xr[j]) - pow43[ix[j]] * step;
763 noise += temp * temp;
772 /*************************************************************************/
774 /*************************************************************************/
776 /* -oo dB => -1.00 */
777 /* - 6 dB => -0.97 */
778 /* - 3 dB => -0.80 */
779 /* - 2 dB => -0.64 */
780 /* - 1 dB => -0.38 */
782 /* + 1 dB => +0.49 */
783 /* + 2 dB => +1.06 */
784 /* + 3 dB => +1.68 */
785 /* + 6 dB => +3.69 */
786 /* +10 dB => +6.45 */
789 calc_noise(gr_info const *const cod_info,
790 FLOAT const *l3_xmin,
791 FLOAT * distort, calc_noise_result * const res, calc_noise_data * prev_noise)
793 int sfb, l, over = 0;
794 FLOAT over_noise_db = 0;
795 FLOAT tot_noise_db = 0; /* 0 dB relative to masking */
796 FLOAT max_noise = -20.0; /* -200 dB relative to masking */
798 const int *scalefac = cod_info->scalefac;
803 for (sfb = 0; sfb < cod_info->psymax; sfb++) {
805 cod_info->global_gain - (((*scalefac++) + (cod_info->preflag ? pretab[sfb] : 0))
806 << (cod_info->scalefac_scale + 1))
807 - cod_info->subblock_gain[cod_info->window[sfb]] * 8;
810 if (prev_noise && (prev_noise->step[sfb] == s)) {
812 /* use previously computed values */
813 noise = prev_noise->noise[sfb];
814 j += cod_info->width[sfb];
815 *distort++ = noise / *l3_xmin++;
817 noise = prev_noise->noise_log[sfb];
821 FLOAT const step = POW20(s);
822 l = cod_info->width[sfb] >> 1;
824 if ((j + cod_info->width[sfb]) > cod_info->max_nonzero_coeff) {
826 usefullsize = cod_info->max_nonzero_coeff - j + 1;
829 l = usefullsize >> 1;
834 noise = calc_noise_core_c(cod_info, &j, l, step);
838 /* save noise values */
839 prev_noise->step[sfb] = s;
840 prev_noise->noise[sfb] = noise;
843 noise = *distort++ = noise / *l3_xmin++;
845 /* multiplying here is adding in dB, but can overflow */
846 noise = FAST_LOG10(Max(noise, 1E-20));
849 /* save noise values */
850 prev_noise->noise_log[sfb] = noise;
855 /* save noise values */
856 prev_noise->global_gain = cod_info->global_gain;;
860 /*tot_noise *= Max(noise, 1E-20); */
861 tot_noise_db += noise;
866 tmp = Max((int) (noise * 10 + .5), 1);
867 res->over_SSD += tmp * tmp;
870 /* multiplying here is adding in dB -but can overflow */
871 /*over_noise *= noise; */
872 over_noise_db += noise;
874 max_noise = Max(max_noise, noise);
878 res->over_count = over;
879 res->tot_noise = tot_noise_db;
880 res->over_noise = over_noise_db;
881 res->max_noise = max_noise;
893 /************************************************************************
897 * updates plotting data
899 * Mark Taylor 2000-??-??
901 * Robert Hegemann: moved noise/distortion calc into it
903 ************************************************************************/
906 set_pinfo(lame_internal_flags const *gfc,
907 gr_info * const cod_info, const III_psy_ratio * const ratio, const int gr, const int ch)
909 SessionConfig_t const *const cfg = &gfc->cfg;
911 int j, i, l, start, end, bw;
913 FLOAT const ifqstep = (cod_info->scalefac_scale == 0) ? .5 : 1.0;
914 int const *const scalefac = cod_info->scalefac;
916 FLOAT l3_xmin[SFBMAX], xfsf[SFBMAX];
917 calc_noise_result noise;
919 (void) calc_xmin(gfc, ratio, cod_info, l3_xmin);
920 (void) calc_noise(cod_info, l3_xmin, xfsf, &noise, 0);
923 sfb2 = cod_info->sfb_lmax;
924 if (cod_info->block_type != SHORT_TYPE && !cod_info->mixed_block_flag)
926 for (sfb = 0; sfb < sfb2; sfb++) {
927 start = gfc->scalefac_band.l[sfb];
928 end = gfc->scalefac_band.l[sfb + 1];
930 for (en0 = 0.0; j < end; j++)
931 en0 += cod_info->xr[j] * cod_info->xr[j];
933 /* convert to MDCT units */
934 en1 = 1e15; /* scaling so it shows up on FFT plot */
935 gfc->pinfo->en[gr][ch][sfb] = en1 * en0;
936 gfc->pinfo->xfsf[gr][ch][sfb] = en1 * l3_xmin[sfb] * xfsf[sfb] / bw;
938 if (ratio->en.l[sfb] > 0 && !cfg->ATHonly)
939 en0 = en0 / ratio->en.l[sfb];
943 gfc->pinfo->thr[gr][ch][sfb] = en1 * Max(en0 * ratio->thm.l[sfb], gfc->ATH->l[sfb]);
945 /* there is no scalefactor bands >= SBPSY_l */
946 gfc->pinfo->LAMEsfb[gr][ch][sfb] = 0;
947 if (cod_info->preflag && sfb >= 11)
948 gfc->pinfo->LAMEsfb[gr][ch][sfb] = -ifqstep * pretab[sfb];
951 assert(scalefac[sfb] >= 0); /* scfsi should be decoded by caller side */
952 gfc->pinfo->LAMEsfb[gr][ch][sfb] -= ifqstep * scalefac[sfb];
956 if (cod_info->block_type == SHORT_TYPE) {
958 for (sfb = cod_info->sfb_smin; sfb < SBMAX_s; sfb++) {
959 start = gfc->scalefac_band.s[sfb];
960 end = gfc->scalefac_band.s[sfb + 1];
962 for (i = 0; i < 3; i++) {
963 for (en0 = 0.0, l = start; l < end; l++) {
964 en0 += cod_info->xr[j] * cod_info->xr[j];
967 en0 = Max(en0 / bw, 1e-20);
968 /* convert to MDCT units */
969 en1 = 1e15; /* scaling so it shows up on FFT plot */
971 gfc->pinfo->en_s[gr][ch][3 * sfb + i] = en1 * en0;
972 gfc->pinfo->xfsf_s[gr][ch][3 * sfb + i] = en1 * l3_xmin[sfb2] * xfsf[sfb2] / bw;
973 if (ratio->en.s[sfb][i] > 0)
974 en0 = en0 / ratio->en.s[sfb][i];
977 if (cfg->ATHonly || cfg->ATHshort)
980 gfc->pinfo->thr_s[gr][ch][3 * sfb + i] =
981 en1 * Max(en0 * ratio->thm.s[sfb][i], gfc->ATH->s[sfb]);
983 /* there is no scalefactor bands >= SBPSY_s */
984 gfc->pinfo->LAMEsfb_s[gr][ch][3 * sfb + i]
985 = -2.0 * cod_info->subblock_gain[i];
987 gfc->pinfo->LAMEsfb_s[gr][ch][3 * sfb + i] -= ifqstep * scalefac[sfb2];
992 } /* block type short */
993 gfc->pinfo->LAMEqss[gr][ch] = cod_info->global_gain;
994 gfc->pinfo->LAMEmainbits[gr][ch] = cod_info->part2_3_length + cod_info->part2_length;
995 gfc->pinfo->LAMEsfbits[gr][ch] = cod_info->part2_length;
997 gfc->pinfo->over[gr][ch] = noise.over_count;
998 gfc->pinfo->max_noise[gr][ch] = noise.max_noise * 10.0;
999 gfc->pinfo->over_noise[gr][ch] = noise.over_noise * 10.0;
1000 gfc->pinfo->tot_noise[gr][ch] = noise.tot_noise * 10.0;
1001 gfc->pinfo->over_SSD[gr][ch] = noise.over_SSD;
1005 /************************************************************************
1009 * updates plotting data for a whole frame
1011 * Robert Hegemann 2000-10-21
1013 ************************************************************************/
1016 set_frame_pinfo(lame_internal_flags * gfc, const III_psy_ratio ratio[2][2])
1018 SessionConfig_t const *const cfg = &gfc->cfg;
1022 /* for every granule and channel patch l3_enc and set info
1024 for (gr = 0; gr < cfg->mode_gr; gr++) {
1025 for (ch = 0; ch < cfg->channels_out; ch++) {
1026 gr_info *const cod_info = &gfc->l3_side.tt[gr][ch];
1027 int scalefac_sav[SFBMAX];
1028 memcpy(scalefac_sav, cod_info->scalefac, sizeof(scalefac_sav));
1030 /* reconstruct the scalefactors in case SCFSI was used
1034 for (sfb = 0; sfb < cod_info->sfb_lmax; sfb++) {
1035 if (cod_info->scalefac[sfb] < 0) /* scfsi */
1036 cod_info->scalefac[sfb] = gfc->l3_side.tt[0][ch].scalefac[sfb];
1040 set_pinfo(gfc, cod_info, &ratio[gr][ch], gr, ch);
1041 memcpy(cod_info->scalefac, scalefac_sav, sizeof(scalefac_sav));