lib/lame/vbrquantize.c

   1 /*
   2  *      MP3 quantization
   3  *
   4  *      Copyright (c) 1999 Mark Taylor
   5  *
   6  * This library is free software; you can redistribute it and/or
   7  * modify it under the terms of the GNU Library General Public
   8  * License as published by the Free Software Foundation; either
   9  * version 2 of the License, or (at your option) any later version.
  10  *
  11  * This library is distributed in the hope that it will be useful,
  12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14  * Library General Public License for more details.
  15  *
  16  * You should have received a copy of the GNU Library General Public
  17  * License along with this library; if not, write to the
  18  * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
  19  * Boston, MA 02111-1307, USA.
  20  */
  21
  22 /* $Id: vbrquantize.c,v 1.1 2002/04/28 17:30:30 kramm Exp $ */
  23
  24 #include "config_static.h"
  25
  26 #include <assert.h>
  27 #include "util.h"
  28 #include "l3side.h"
  29 #include "reservoir.h"
  30 #include "quantize_pvt.h"
  31 #include "vbrquantize.h"
  32
  33 #ifdef WITH_DMALLOC
  34 #include <dmalloc.h>
  35 #endif
  36
  37
  38
  39 typedef union {
  40     float f;
  41     int   i;
  42 } fi_union;
  43
  44 #define MAGIC_FLOAT (65536*(128))
  45 #define MAGIC_INT    0x4b000000
  46
  47 #ifdef TAKEHIRO_IEEE754_HACK
  48
  49 #define DUFFBLOCKMQ() do { \
  50         xp = xr34[0] * sfpow34_p1; \
  51         xe = xr34[0] * sfpow34_eq; \
  52         xm = xr34[0] * sfpow34_m1; \
  53         if (xm > IXMAX_VAL)  \
  54             return -1; \
  55         xp += MAGIC_FLOAT; \
  56         xe += MAGIC_FLOAT; \
  57         xm += MAGIC_FLOAT; \
  58         fi[0].f = xp; \
  59         fi[1].f = xe; \
  60         fi[2].f = xm; \
  61         fi[0].f = xp + (adj43asm - MAGIC_INT)[fi[0].i]; \
  62         fi[1].f = xe + (adj43asm - MAGIC_INT)[fi[1].i]; \
  63         fi[2].f = xm + (adj43asm - MAGIC_INT)[fi[2].i]; \
  64         fi[0].i -= MAGIC_INT; \
  65         fi[1].i -= MAGIC_INT; \
  66         fi[2].i -= MAGIC_INT; \
  67         x0 = fabs(xr[0]); \
  68         xp = x0 - pow43[fi[0].i] * sfpow_p1; \
  69         xe = x0 - pow43[fi[1].i] * sfpow_eq; \
  70         xm = x0 - pow43[fi[2].i] * sfpow_m1; \
  71         xp *= xp; \
  72         xe *= xe; \
  73         xm *= xm; \
  74         xfsf_eq = Max(xfsf_eq, xe); \
  75         xfsf_p1 = Max(xfsf_p1, xp); \
  76         xfsf_m1 = Max(xfsf_m1, xm); \
  77         ++xr; \
  78         ++xr34; \
  79     } while(0)
  80
  81 #define DUFFBLOCK() do { \
  82         xp = xr34[0] * sfpow34_p1; \
  83         xe = xr34[0] * sfpow34_eq; \
  84         xm = xr34[0] * sfpow34_m1; \
  85         if (xm > IXMAX_VAL)  \
  86             return -1; \
  87         xp += MAGIC_FLOAT; \
  88         xe += MAGIC_FLOAT; \
  89         xm += MAGIC_FLOAT; \
  90         fi[0].f = xp; \
  91         fi[1].f = xe; \
  92         fi[2].f = xm; \
  93         fi[0].f = xp + (adj43asm - MAGIC_INT)[fi[0].i]; \
  94         fi[1].f = xe + (adj43asm - MAGIC_INT)[fi[1].i]; \
  95         fi[2].f = xm + (adj43asm - MAGIC_INT)[fi[2].i]; \
  96         fi[0].i -= MAGIC_INT; \
  97         fi[1].i -= MAGIC_INT; \
  98         fi[2].i -= MAGIC_INT; \
  99         x0 = fabs(xr[0]); \
 100         xp = x0 - pow43[fi[0].i] * sfpow_p1; \
 101         xe = x0 - pow43[fi[1].i] * sfpow_eq; \
 102         xm = x0 - pow43[fi[2].i] * sfpow_m1; \
 103         xfsf_p1 += xp * xp; \
 104         xfsf_eq += xe * xe; \
 105         xfsf_m1 += xm * xm; \
 106         ++xr; \
 107         ++xr34; \
 108     } while(0)
 109
 110 #else
 111
 112 /*********************************************************************
 113  * XRPOW_FTOI is a macro to convert floats to ints.
 114  * if XRPOW_FTOI(x) = nearest_int(x), then QUANTFAC(x)=adj43asm[x]
 115  *                                         ROUNDFAC= -0.0946
 116  *
 117  * if XRPOW_FTOI(x) = floor(x), then QUANTFAC(x)=asj43[x]
 118  *                                   ROUNDFAC=0.4054
 119  *********************************************************************/
 120 #  define QUANTFAC(rx)  adj43[rx]
 121 #  define ROUNDFAC 0.4054
 122 #  define XRPOW_FTOI(src,dest) ((dest) = (int)(src))
 123
 124
 125 #endif
 126
 127 /*  caution: a[] will be resorted!!
 128  */
 129 FLOAT8 select_kth(FLOAT8 a[], int N, int k)
 130 {
 131     int i, j, l, r;
 132     FLOAT8 v, w;
 133
 134     l = 0;
 135     r = N-1;
 136     while (r > l) {
 137         v = a[r];
 138         i = l-1;
 139         j = r;
 140         for (;;) {
 141             while (a[++i] < v) /*empty*/;
 142             while (a[--j] > v) /*empty*/;
 143             if (i >= j)
 144                 break;
 145             /* swap i and j */
 146             w = a[i];
 147             a[i] = a[j];
 148             a[j] = w;
 149         }
 150         /* swap i and r */
 151         w = a[i];
 152         a[i] = a[r];
 153         a[r] = w;
 154         if (i >= k)
 155             r = i-1;
 156         if (i <= k)
 157             l = i+1;
 158     }
 159     return a[k];
 160 }
 161
 162
 163
 164
 165 static FLOAT8
 166 calc_sfb_noise(const FLOAT8 *xr, const FLOAT8 *xr34, int bw, int sf)
 167 {
 168     int j;
 169     fi_union fi;
 170     FLOAT8 temp;
 171     FLOAT8 xfsf = 0;
 172     FLOAT8 sfpow, sfpow34;
 173
 174     sfpow   =  POW20(sf+210); /*pow(2.0,sf/4.0); */
 175     sfpow34 = IPOW20(sf+210); /*pow(sfpow,-3.0/4.0);*/
 176
 177     for ( j = 0; j < bw ; ++j ) {
 178         if ( xr34[j]*sfpow34 > IXMAX_VAL ) return -1;
 179
 180 #ifdef TAKEHIRO_IEEE754_HACK
 181         temp   = sfpow34*xr34[j];
 182         temp  += MAGIC_FLOAT;
 183         fi.f  = temp;
 184         fi.f  = temp + (adj43asm - MAGIC_INT)[fi.i];
 185         fi.i -= MAGIC_INT;
 186 #else
 187         temp = xr34[j]*sfpow34;
 188         XRPOW_FTOI(temp, fi.i);
 189         XRPOW_FTOI(temp + QUANTFAC(fi.i), fi.i);
 190 #endif
 191
 192         temp = fabs(xr[j])- pow43[fi.i]*sfpow;
 193         xfsf += temp*temp;
 194     }
 195     return xfsf;
 196 }
 197
 198
 199
 200
 201 static FLOAT8
 202 calc_sfb_noise_mq( const FLOAT8 *xr, const FLOAT8 *xr34, int bw, int sf,
 203                    int mq, FLOAT8 *scratch )
 204 {
 205     int j, k;
 206     fi_union fi;
 207     FLOAT8 temp;
 208     FLOAT8 sfpow, sfpow34, xfsfm = 0, xfsf = 0;
 209
 210     sfpow   =  POW20(sf+210); /*pow(2.0,sf/4.0); */
 211     sfpow34 = IPOW20(sf+210); /*pow(sfpow,-3.0/4.0);*/
 212
 213     for ( j = 0; j < bw; ++j ) {
 214         if ( xr34[j]*sfpow34 > IXMAX_VAL ) return -1;
 215
 216 #ifdef TAKEHIRO_IEEE754_HACK
 217         temp  = sfpow34*xr34[j];
 218         temp += MAGIC_FLOAT;
 219         fi.f  = temp;
 220         fi.f  = temp + (adj43asm - MAGIC_INT)[fi.i];
 221         fi.i -= MAGIC_INT;
 222 #else
 223         temp = xr34[j]*sfpow34;
 224         XRPOW_FTOI(temp, fi.i);
 225         XRPOW_FTOI(temp + QUANTFAC(fi.i), fi.i);
 226 #endif
 227
 228         temp = fabs(xr[j])- pow43[fi.i]*sfpow;
 229         temp *= temp;
 230
 231         scratch[j] = temp;
 232         if ( xfsfm < temp ) xfsfm = temp;
 233         xfsf += temp;
 234     }
 235     if ( mq == 1 ) return bw*select_kth(scratch,bw,bw*13/16);
 236
 237     xfsf /= bw;
 238     for ( k = 1, j = 0; j < bw; ++j ) {
 239         if ( scratch[j] > xfsf ) {
 240             xfsfm += scratch[j];
 241             ++k;
 242         }
 243     }
 244     return xfsfm/k * bw;
 245 }
 246
 247
 248
 249
 250 static FLOAT8
 251 calc_sfb_noise_ave(const FLOAT8 *xr, const FLOAT8 *xr34, int bw, int sf)
 252 {
 253     double xp;
 254     double xe;
 255     double xm;
 256 #ifdef TAKEHIRO_IEEE754_HACK
 257     double x0;
 258 #endif
 259     int xx[3], j;
 260     fi_union *fi = (fi_union *)xx;
 261     FLOAT8 sfpow34_eq, sfpow34_p1, sfpow34_m1;
 262     FLOAT8 sfpow_eq, sfpow_p1, sfpow_m1;
 263     FLOAT8 xfsf_eq = 0, xfsf_p1 = 0, xfsf_m1 = 0;
 264
 265     sfpow_eq = POW20(sf + 210); /*pow(2.0,sf/4.0); */
 266     sfpow_m1 = sfpow_eq * .8408964153;  /* pow(2,(sf-1)/4.0) */
 267     sfpow_p1 = sfpow_eq * 1.189207115;
 268
 269     sfpow34_eq = IPOW20(sf + 210); /*pow(sfpow,-3.0/4.0);*/
 270     sfpow34_m1 = sfpow34_eq * 1.13878863476;       /* .84089 ^ -3/4 */
 271     sfpow34_p1 = sfpow34_eq * 0.878126080187;
 272
 273 #ifdef TAKEHIRO_IEEE754_HACK
 274     /*
 275      *  loop unrolled into "Duff's Device".   Robert Hegemann
 276      */
 277     j = (bw+3) / 4;
 278     switch (bw % 4) {
 279         default:
 280         case 0: do{ DUFFBLOCK();
 281         case 3:     DUFFBLOCK();
 282         case 2:     DUFFBLOCK();
 283         case 1:     DUFFBLOCK(); } while (--j);
 284     }
 285 #else
 286     for (j = 0; j < bw; ++j) {
 287
 288         if (xr34[j]*sfpow34_m1 > IXMAX_VAL) return -1;
 289
 290         xe = xr34[j]*sfpow34_eq;
 291         XRPOW_FTOI(xe, fi[0].i);
 292         XRPOW_FTOI(xe + QUANTFAC(fi[0].i), fi[0].i);
 293         xe = fabs(xr[j])- pow43[fi[0].i]*sfpow_eq;
 294         xe *= xe;
 295
 296         xp = xr34[j]*sfpow34_p1;
 297         XRPOW_FTOI(xp, fi[0].i);
 298         XRPOW_FTOI(xp + QUANTFAC(fi[0].i), fi[0].i);
 299         xp = fabs(xr[j])- pow43[fi[0].i]*sfpow_p1;
 300         xp *= xp;
 301
 302         xm = xr34[j]*sfpow34_m1;
 303         XRPOW_FTOI(xm, fi[0].i);
 304         XRPOW_FTOI(xm + QUANTFAC(fi[0].i), fi[0].i);
 305         xm = fabs(xr[j])- pow43[fi[0].i]*sfpow_m1;
 306         xm *= xm;
 307
 308         xfsf_eq += xe;
 309         xfsf_p1 += xp;
 310         xfsf_m1 += xm;
 311     }
 312 #endif
 313
 314     if ( xfsf_eq < xfsf_p1 ) xfsf_eq = xfsf_p1;
 315     if ( xfsf_eq < xfsf_m1 ) xfsf_eq = xfsf_m1;
 316     return xfsf_eq;
 317 }
 318
 319
 320
 321 INLINE int
 322 find_scalefac( const FLOAT8 *xr, const FLOAT8 *xr34, FLOAT8 l3_xmin, int bw )
 323 {
 324     FLOAT8 xfsf;
 325     int i, sf, sf_ok, delsf;
 326
 327     /* search will range from sf:  -209 -> 45  */
 328     sf = -82;
 329     delsf = 128;
 330
 331     sf_ok = 10000;
 332     for (i = 0; i < 7; ++i) {
 333         delsf /= 2;
 334         xfsf = calc_sfb_noise( xr, xr34, bw, sf );
 335
 336         if (xfsf < 0) {
 337             /* scalefactors too small */
 338             sf += delsf;
 339         }
 340         else {
 341             if (xfsf > l3_xmin) {
 342                 /* distortion.  try a smaller scalefactor */
 343                 sf -= delsf;
 344             }
 345             else {
 346                 sf_ok = sf;
 347                 sf += delsf;
 348             }
 349         }
 350     }
 351
 352     /*  returning a scalefac without distortion, if possible
 353      */
 354     return sf_ok > 45 ? sf : sf_ok;
 355 }
 356
 357 INLINE int
 358 find_scalefac_mq( const FLOAT8 *xr, const FLOAT8 *xr34, FLOAT8 l3_xmin,
 359                   int bw, int mq, FLOAT8 *scratch )
 360 {
 361     FLOAT8 xfsf;
 362     int i, sf, sf_ok, delsf;
 363
 364     /* search will range from sf:  -209 -> 45  */
 365     sf = -82;
 366     delsf = 128;
 367
 368     sf_ok = 10000;
 369     for (i = 0; i < 7; ++i) {
 370         delsf /= 2;
 371         xfsf = calc_sfb_noise_mq( xr, xr34, bw, sf, mq, scratch );
 372
 373         if (xfsf < 0) {
 374             /* scalefactors too small */
 375             sf += delsf;
 376         }
 377         else {
 378             if (xfsf > l3_xmin) {
 379                 /* distortion.  try a smaller scalefactor */
 380                 sf -= delsf;
 381             }
 382             else {
 383                 sf_ok = sf;
 384                 sf += delsf;
 385             }
 386         }
 387     }
 388
 389     /*  returning a scalefac without distortion, if possible
 390      */
 391     return sf_ok > 45 ? sf : sf_ok;
 392 }
 393
 394 INLINE int
 395 find_scalefac_ave( const FLOAT8 *xr, const FLOAT8 *xr34, FLOAT8 l3_xmin, int bw )
 396 {
 397     FLOAT8 xfsf;
 398     int i, sf, sf_ok, delsf;
 399
 400     /* search will range from sf:  -209 -> 45  */
 401     sf = -82;
 402     delsf = 128;
 403
 404     sf_ok = 10000;
 405     for (i = 0; i < 7; ++i) {
 406         delsf /= 2;
 407         xfsf = calc_sfb_noise_ave( xr, xr34, bw, sf );
 408
 409         if (xfsf < 0) {
 410             /* scalefactors too small */
 411             sf += delsf;
 412         }
 413         else{
 414             if (xfsf > l3_xmin) {
 415                 /* distortion.  try a smaller scalefactor */
 416                 sf -= delsf;
 417             }
 418             else {
 419                 sf_ok = sf;
 420                 sf += delsf;
 421             }
 422         }
 423     }
 424
 425     /*  returning a scalefac without distortion, if possible
 426      */
 427     return sf_ok > 45 ? sf : sf_ok;
 428 }
 429
 430
 431 /**
 432  *  Robert Hegemann 2001-05-01
 433  *  calculates quantization step size determined by allowed masking
 434  */
 435 INLINE int
 436 calc_scalefac( FLOAT8 l3_xmin, int bw, FLOAT8 preset_tune )
 437 {
 438     FLOAT8 const c = (preset_tune > 0 ? preset_tune
 439                                               : 5.799142446);   // 10 * 10^(2/3) * log10(4/3)
 440     return (int)(c*log10(l3_xmin/bw)-.5);
 441 }
 442
 443
 444
 445 static const int max_range_short[SBMAX_s] =
 446 {15, 15, 15, 15, 15, 15, 7, 7, 7, 7, 7, 7, 0 };
 447
 448 static const int max_range_long[SBMAX_l] =
 449 {15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 0};
 450
 451 static const int max_range_long_lsf_pretab[SBMAX_l] =
 452 { 7,7,7,7,7,7, 3,3,3,3,3, 0,0,0,0, 0,0,0, 0,0,0, 0 };
 453
 454
 455
 456 /*
 457     sfb=0..5  scalefac < 16
 458     sfb>5     scalefac < 8
 459
 460     ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
 461     ol_sf =  (cod_info->global_gain-210.0);
 462     ol_sf -= 8*cod_info->subblock_gain[i];
 463     ol_sf -= ifqstep*scalefac[gr][ch].s[sfb][i];
 464
 465 */
 466 INLINE int
 467 compute_scalefacs_short(
 468     int sf[][3], const gr_info *cod_info, int scalefac[][3], int *sbg )
 469 {
 470     const int maxrange1 = 15, maxrange2 = 7;
 471     int maxrange, maxover = 0;
 472     int sfb, i;
 473     int ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
 474
 475     for (i = 0; i < 3; ++i) {
 476         int maxsf1 = 0, maxsf2 = 0, minsf = 1000;
 477         /* see if we should use subblock gain */
 478         for (sfb = 0; sfb < 6; ++sfb) {         /* part 1 */
 479             if (maxsf1 < -sf[sfb][i]) maxsf1 = -sf[sfb][i];
 480             if (minsf  > -sf[sfb][i]) minsf  = -sf[sfb][i];
 481         }
 482         for (; sfb < SBPSY_s; ++sfb) {          /* part 2 */
 483             if (maxsf2 < -sf[sfb][i]) maxsf2 = -sf[sfb][i];
 484             if (minsf  > -sf[sfb][i]) minsf  = -sf[sfb][i];
 485         }
 486
 487         /* boost subblock gain as little as possible so we can
 488          * reach maxsf1 with scalefactors
 489          * 8*sbg >= maxsf1
 490          */
 491         maxsf1 = Max (maxsf1-maxrange1*ifqstep, maxsf2-maxrange2*ifqstep);
 492         sbg[i] = 0;
 493         if (minsf  > 0) sbg[i] = floor (.125*minsf + .001);
 494         if (maxsf1 > 0) sbg[i] = Max (sbg[i], (maxsf1/8 + (maxsf1 % 8 != 0)));
 495         if (sbg[i] > 7) sbg[i] = 7;
 496
 497         for (sfb = 0; sfb < SBPSY_s; ++sfb) {
 498             sf[sfb][i] += 8*sbg[i];
 499
 500             if (sf[sfb][i] < 0) {
 501                 maxrange = sfb < 6 ? maxrange1 : maxrange2;
 502
 503                 scalefac[sfb][i]
 504                      = -sf[sfb][i]/ifqstep + (-sf[sfb][i]%ifqstep != 0);
 505
 506                 if (scalefac[sfb][i] > maxrange)
 507                     scalefac[sfb][i] = maxrange;
 508
 509                 if (maxover < -(sf[sfb][i] + scalefac[sfb][i]*ifqstep))
 510                     maxover = -(sf[sfb][i] + scalefac[sfb][i]*ifqstep);
 511             }
 512         }
 513         scalefac[sfb][i] = 0;
 514     }
 515
 516     return maxover;
 517 }
 518
 519
 520
 521 /*
 522           ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
 523           ol_sf =  (cod_info->global_gain-210.0);
 524           ol_sf -= ifqstep*scalefac[gr][ch].l[sfb];
 525           if (cod_info->preflag && sfb>=11)
 526           ol_sf -= ifqstep*pretab[sfb];
 527 */
 528 INLINE int
 529 compute_scalefacs_long_lsf( int *sf, const gr_info * cod_info, int *scalefac )
 530 {
 531     const int * max_range = max_range_long;
 532     int ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
 533     int sfb;
 534     int maxover;
 535
 536     if (cod_info->preflag) {
 537         max_range = max_range_long_lsf_pretab;
 538         for (sfb = 11; sfb < SBPSY_l; ++sfb)
 539             sf[sfb] += pretab[sfb] * ifqstep;
 540     }
 541
 542     maxover = 0;
 543     for (sfb = 0; sfb < SBPSY_l; ++sfb) {
 544
 545         if (sf[sfb] < 0) {
 546             /* ifqstep*scalefac >= -sf[sfb], so round UP */
 547             scalefac[sfb] = -sf[sfb]/ifqstep  + (-sf[sfb] % ifqstep != 0);
 548             if (scalefac[sfb] > max_range[sfb])
 549                 scalefac[sfb] = max_range[sfb];
 550
 551             /* sf[sfb] should now be positive: */
 552             if (-(sf[sfb] + scalefac[sfb]*ifqstep) > maxover) {
 553                 maxover = -(sf[sfb] + scalefac[sfb]*ifqstep);
 554             }
 555         }
 556     }
 557     scalefac[sfb] = 0;  /* sfb21 */
 558
 559     return maxover;
 560 }
 561
 562
 563
 564
 565
 566 /*
 567           ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
 568           ol_sf =  (cod_info->global_gain-210.0);
 569           ol_sf -= ifqstep*scalefac[gr][ch].l[sfb];
 570           if (cod_info->preflag && sfb>=11)
 571           ol_sf -= ifqstep*pretab[sfb];
 572 */
 573 INLINE int
 574 compute_scalefacs_long( int *sf, const gr_info * cod_info, int *scalefac )
 575 {
 576     int ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
 577     int sfb;
 578     int maxover;
 579
 580     if (cod_info->preflag) {
 581         for (sfb = 11; sfb < SBPSY_l; ++sfb)
 582             sf[sfb] += pretab[sfb] * ifqstep;
 583     }
 584
 585     maxover = 0;
 586     for (sfb = 0; sfb < SBPSY_l; ++sfb) {
 587
 588         if (sf[sfb] < 0) {
 589             /* ifqstep*scalefac >= -sf[sfb], so round UP */
 590             scalefac[sfb] = -sf[sfb]/ifqstep  + (-sf[sfb] % ifqstep != 0);
 591             if (scalefac[sfb] > max_range_long[sfb])
 592                 scalefac[sfb] = max_range_long[sfb];
 593
 594             /* sf[sfb] should now be positive: */
 595             if (-(sf[sfb] + scalefac[sfb]*ifqstep) > maxover) {
 596                 maxover = -(sf[sfb] + scalefac[sfb]*ifqstep);
 597             }
 598         }
 599     }
 600     scalefac[sfb] = 0;  /* sfb21 */
 601
 602     return maxover;
 603 }
 604
 605
 606
 607
 608
 609
 610
 611
 612 /************************************************************************
 613  *
 614  * quantize and encode with the given scalefacs and global gain
 615  *
 616  * compute scalefactors, l3_enc, and return number of bits needed to encode
 617  *
 618  *
 619  ************************************************************************/
 620
 621 static int
 622 VBR_quantize_granule(
 623           lame_internal_flags * gfc,
 624           FLOAT8              * xr34,
 625           int                 * l3_enc,
 626           III_scalefac_t      * scalefac,
 627           int gr, int ch )
 628 {
 629     int status;
 630     III_side_info_t * l3_side  = & gfc->l3_side;
 631     gr_info         * cod_info = & l3_side->gr[gr].ch[ch].tt;
 632
 633     /* encode scalefacs */
 634     if ( gfc->is_mpeg1 )
 635         status = scale_bitcount( scalefac, cod_info );
 636     else
 637         status = scale_bitcount_lsf( gfc, scalefac, cod_info );
 638
 639     if (status != 0) {
 640         return -1;
 641     }
 642
 643     /* quantize xr34 */
 644     cod_info->part2_3_length = count_bits(gfc,l3_enc,xr34,cod_info);
 645     if (cod_info->part2_3_length >= LARGE_BITS) return -2;
 646     cod_info->part2_3_length += cod_info->part2_length;
 647
 648
 649     if (gfc->use_best_huffman == 1) {
 650         best_huffman_divide(gfc, cod_info, l3_enc);
 651     }
 652     return 0;
 653 }
 654
 655
 656
 657 /***********************************************************************
 658  *
 659  *      calc_short_block_vbr_sf()
 660  *      calc_long_block_vbr_sf()
 661  *
 662  *  Mark Taylor 2000-??-??
 663  *  Robert Hegemann 2000-10-25 made functions of it
 664  *
 665  ***********************************************************************/
 666 static const int MAX_SF_DELTA = 4;
 667
 668 static void
 669 short_block_sf (
 670     const lame_internal_flags * gfc,
 671     const III_psy_xmin   * l3_xmin,
 672     const FLOAT8         * xr34_orig,
 673     const FLOAT8         * xr34,
 674           III_scalefac_t * vbrsf,
 675           int            * vbrmin,
 676           int            * vbrmax )
 677 {
 678     int j, sfb, b;
 679     int vbrmean, vbrmn;
 680     int sf_cache[SBMAX_s];
 681     int scalefac_criteria;
 682
 683     if (gfc->presetTune.use) {
 684             /* map experimentalX settings to internal selections */
 685         static char const map[] = {2,1,0,3,6};
 686         scalefac_criteria = map[gfc->presetTune.quantcomp_current];
 687     }
 688         else {
 689         scalefac_criteria = gfc->VBR->quality;
 690     }
 691
 692     for (j = 0, sfb = 0; sfb < SBMAX_s; ++sfb) {
 693         for (b = 0; b < 3; ++b) {
 694             const  int start = gfc->scalefac_band.s[ sfb ];
 695             const  int end   = gfc->scalefac_band.s[ sfb+1 ];
 696             const  int width = end - start;
 697
 698             switch( scalefac_criteria ) {
 699             default:
 700                 /*  the fastest
 701                  */
 702                 vbrsf->s[sfb][b] = calc_scalefac( l3_xmin->s[sfb][b], width,
 703                                                                       gfc->presetTune.quantcomp_adjust_mtrh );
 704                 break;
 705             case 5:
 706             case 4:
 707             case 3:
 708                 /*  the faster and sloppier mode to use at lower quality
 709                  */
 710                 vbrsf->s[sfb][b] = find_scalefac (&xr34[j], &xr34_orig[j],
 711                                               l3_xmin->s[sfb][b], width);
 712                 break;
 713             case 2:
 714                 /*  the slower and better mode to use at higher quality
 715                  */
 716                 vbrsf->s[sfb][b] = find_scalefac_ave (&xr34[j], &xr34_orig[j],
 717                                                     l3_xmin->s[sfb][b], width);
 718                 break;
 719             case 1:
 720                 /*  maxnoise mode to use at higher quality
 721                  */
 722                 vbrsf->s[sfb][b] = find_scalefac_mq (&xr34[j], &xr34_orig[j],
 723                                               l3_xmin->s[sfb][b], width,
 724                                               1, gfc->VBR->scratch);
 725                 break;
 726             case 0:
 727                 /*  maxnoise mode to use at higher quality
 728                  */
 729                 vbrsf->s[sfb][b] = find_scalefac_mq (&xr34[j], &xr34_orig[j],
 730                                               l3_xmin->s[sfb][b], width,
 731                                               0, gfc->VBR->scratch);
 732                 break;
 733             }
 734             j += width;
 735         }
 736     }
 737
 738     *vbrmax = -10000;
 739     *vbrmin = +10000;
 740
 741     for (b = 0; b < 3; ++b) {
 742
 743         /*  smoothing
 744          */
 745         switch( gfc->VBR->smooth ) {
 746         default:
 747         case  0:
 748             break;
 749
 750         case  1:
 751             /*  make working copy, get min value, select_kth_int will reorder!
 752              */
 753             for (vbrmn = +10000, sfb = 0; sfb < SBMAX_s; ++sfb) {
 754                 sf_cache[sfb] = vbrsf->s[sfb][b];
 755                 if (vbrmn > vbrsf->s[sfb][b])
 756                     vbrmn = vbrsf->s[sfb][b];
 757             }
 758
 759             /*  find median value, take it as mean
 760              */
 761             vbrmean = select_kth_int (sf_cache, SBMAX_s, (SBMAX_s+1)/2);
 762
 763             /*  cut peaks
 764              */
 765             for (sfb = 0; sfb < SBMAX_s; ++sfb) {
 766                 if (vbrsf->s[sfb][b] > vbrmean+(vbrmean-vbrmn))
 767                     vbrsf->s[sfb][b] = vbrmean+(vbrmean-vbrmn);
 768             }
 769             break;
 770
 771         case  2:
 772             for (sfb = 0; sfb < SBMAX_s; ++sfb) {
 773                 if (sfb > 0)
 774                     if (vbrsf->s[sfb][b] > vbrsf->s[sfb-1][b]+MAX_SF_DELTA)
 775                         vbrsf->s[sfb][b] = vbrsf->s[sfb-1][b]+MAX_SF_DELTA;
 776                 if (sfb < SBMAX_s-1)
 777                     if (vbrsf->s[sfb][b] > vbrsf->s[sfb+1][b]+MAX_SF_DELTA)
 778                         vbrsf->s[sfb][b] = vbrsf->s[sfb+1][b]+MAX_SF_DELTA;
 779             }
 780             break;
 781         }
 782
 783         /*  get max value
 784          */
 785         for (sfb = 0; sfb < SBMAX_s; ++sfb) {
 786             if (*vbrmax < vbrsf->s[sfb][b])
 787                 *vbrmax = vbrsf->s[sfb][b];
 788             if (*vbrmin > vbrsf->s[sfb][b])
 789                 *vbrmin = vbrsf->s[sfb][b];
 790         }
 791     }
 792 }
 793
 794
 795     /* a variation for vbr-mtrh */
 796 static void
 797 long_block_sf (
 798     const lame_internal_flags * gfc,
 799     const III_psy_xmin   * l3_xmin,
 800     const FLOAT8         * xr34_orig,
 801     const FLOAT8         * xr34,
 802           III_scalefac_t * vbrsf,
 803           int            * vbrmin,
 804           int            * vbrmax )
 805 {
 806     int sfb;
 807     int vbrmean, vbrmn;
 808     int sf_cache[SBMAX_l];
 809     int scalefac_criteria;
 810
 811     if (gfc->presetTune.use) {
 812             /* map experimentalX settings to internal selections */
 813         static char const map[] = {2,1,0,3,6};
 814         scalefac_criteria = map[gfc->presetTune.quantcomp_current];
 815     }
 816         else {
 817         scalefac_criteria = gfc->VBR->quality;
 818     }
 819
 820     for (sfb = 0; sfb < SBMAX_l; ++sfb) {
 821         const  int start = gfc->scalefac_band.l[ sfb ];
 822         const  int end   = gfc->scalefac_band.l[ sfb+1 ];
 823         const  int width = end - start;
 824
 825         switch( scalefac_criteria ) {
 826         default:
 827             /*  the fastest
 828              */
 829             vbrsf->l[sfb] = calc_scalefac( l3_xmin->l[sfb], width,
 830                                                            gfc->presetTune.quantcomp_adjust_mtrh );
 831             break;
 832         case 5:
 833         case 4:
 834         case 3:
 835             /*  the faster and sloppier mode to use at lower quality
 836              */
 837             vbrsf->l[sfb] = find_scalefac (&xr34[start], &xr34_orig[start],
 838                                            l3_xmin->l[sfb], width);
 839             break;
 840         case 2:
 841             /*  the slower and better mode to use at higher quality
 842              */
 843             vbrsf->l[sfb] = find_scalefac_ave (&xr34[start], &xr34_orig[start],
 844                                                l3_xmin->l[sfb], width);
 845             break;
 846         case 1:
 847             /*  maxnoise mode to use at higher quality
 848              */
 849             vbrsf->l[sfb] = find_scalefac_mq (&xr34[start], &xr34_orig[start],
 850                                            l3_xmin->l[sfb], width,
 851                                            1, gfc->VBR->scratch);
 852             break;
 853         case 0:
 854             /*  maxnoise mode to use at higher quality
 855              */
 856             vbrsf->l[sfb] = find_scalefac_mq (&xr34[start], &xr34_orig[start],
 857                                            l3_xmin->l[sfb], width,
 858                                            0, gfc->VBR->scratch);
 859             break;
 860         }
 861     }
 862
 863     switch( gfc->VBR->smooth ) {
 864     default:
 865     case  0:
 866         break;
 867
 868     case  1:
 869         /*  make working copy, get min value, select_kth_int will reorder!
 870          */
 871         for (vbrmn = +10000, sfb = 0; sfb < SBMAX_l; ++sfb) {
 872             sf_cache[sfb] = vbrsf->l[sfb];
 873             if (vbrmn > vbrsf->l[sfb])
 874                 vbrmn = vbrsf->l[sfb];
 875         }
 876         /*  find median value, take it as mean
 877          */
 878         vbrmean = select_kth_int (sf_cache, SBMAX_l, (SBMAX_l+1)/2);
 879
 880         /*  cut peaks
 881          */
 882         for (sfb = 0; sfb < SBMAX_l; ++sfb) {
 883             if (vbrsf->l[sfb] > vbrmean+(vbrmean-vbrmn))
 884                 vbrsf->l[sfb] = vbrmean+(vbrmean-vbrmn);
 885         }
 886         break;
 887
 888     case  2:
 889         for (sfb = 0; sfb < SBMAX_l; ++sfb) {
 890             if (sfb > 0)
 891                 if (vbrsf->l[sfb] > vbrsf->l[sfb-1]+MAX_SF_DELTA)
 892                     vbrsf->l[sfb] = vbrsf->l[sfb-1]+MAX_SF_DELTA;
 893             if (sfb < SBMAX_l-1)
 894                 if (vbrsf->l[sfb] > vbrsf->l[sfb+1]+MAX_SF_DELTA)
 895                     vbrsf->l[sfb] = vbrsf->l[sfb+1]+MAX_SF_DELTA;
 896         }
 897         break;
 898     }
 899
 900     /*  get max value
 901      */
 902     for (*vbrmin = +10000, *vbrmax = -10000, sfb = 0; sfb < SBMAX_l; ++sfb) {
 903         if (*vbrmax < vbrsf->l[sfb])
 904             *vbrmax = vbrsf->l[sfb];
 905         if (*vbrmin > vbrsf->l[sfb])
 906             *vbrmin = vbrsf->l[sfb];
 907     }
 908 }
 909
 910
 911
 912 /******************************************************************
 913  *
 914  *  short block scalefacs
 915  *
 916  ******************************************************************/
 917
 918 static void
 919 short_block_scalefacs (
 920     const lame_internal_flags * gfc,
 921           gr_info             * cod_info,
 922           III_scalefac_t      * scalefac,
 923           III_scalefac_t      * vbrsf,
 924           int                 * VBRmax )
 925 {
 926     int sfb, maxsfb, b;
 927     int maxover, maxover0, maxover1, mover;
 928     int v0, v1;
 929     int minsfb;
 930     int vbrmax = *VBRmax;
 931
 932     maxover0 = 0;
 933     maxover1 = 0;
 934     maxsfb = gfc->sfb21_extra ? SBMAX_s : SBPSY_s;
 935     for (sfb = 0; sfb < maxsfb; ++sfb) {
 936         for (b = 0; b < 3; ++b) {
 937             v0 = (vbrmax - vbrsf->s[sfb][b]) - (4*14 + 2*max_range_short[sfb]);
 938             v1 = (vbrmax - vbrsf->s[sfb][b]) - (4*14 + 4*max_range_short[sfb]);
 939             if (maxover0 < v0)
 940                 maxover0 = v0;
 941             if (maxover1 < v1)
 942                 maxover1 = v1;
 943         }
 944     }
 945
 946     if ((gfc->noise_shaping == 2) && (gfc->presetTune.use && !(gfc->presetTune.athadjust_safe_noiseshaping || gfc->ATH->adjust < 1.0)))
 947         /* allow scalefac_scale=1 */
 948         mover = Min (maxover0, maxover1);
 949     else
 950         mover = maxover0;
 951
 952     vbrmax   -= mover;
 953     maxover0 -= mover;
 954     maxover1 -= mover;
 955
 956     if (maxover0 == 0)
 957         cod_info->scalefac_scale = 0;
 958     else if (maxover1 == 0)
 959         cod_info->scalefac_scale = 1;
 960
 961     /* sf =  (cod_info->global_gain-210.0) */
 962     cod_info->global_gain = vbrmax + 210;
 963     assert(cod_info->global_gain < 256);
 964
 965     if (cod_info->global_gain < 0) {
 966         cod_info->global_gain = 0;
 967     }
 968     else
 969     if (cod_info->global_gain > 255) {
 970         cod_info->global_gain = 255;
 971     }
 972     for (sfb = 0; sfb < SBMAX_s; ++sfb) {
 973         for (b = 0; b < 3; ++b) {
 974             vbrsf->s[sfb][b] -= vbrmax;
 975         }
 976     }
 977     maxover = compute_scalefacs_short (vbrsf->s, cod_info, scalefac->s,
 978                                            cod_info->subblock_gain);
 979
 980     assert (maxover <= 0);
 981
 982     /* adjust global_gain so at least 1 subblock gain = 0 */
 983     minsfb = 999; /* prepare for minimum search */
 984     for (b = 0; b < 3; ++b)
 985         if (minsfb > cod_info->subblock_gain[b])
 986             minsfb = cod_info->subblock_gain[b];
 987
 988     if (minsfb > cod_info->global_gain/8)
 989         minsfb = cod_info->global_gain/8;
 990
 991     vbrmax                -= 8*minsfb;
 992     cod_info->global_gain -= 8*minsfb;
 993
 994     for (b = 0; b < 3; ++b)
 995         cod_info->subblock_gain[b] -= minsfb;
 996
 997     *VBRmax = vbrmax;
 998 }
 999
1000
1001
1002 /******************************************************************
1003  *
1004  *  long block scalefacs
1005  *
1006  ******************************************************************/
1007
1008 static void
1009 long_block_scalefacs (
1010     const lame_internal_flags * gfc,
1011           gr_info             * cod_info,
1012           III_scalefac_t      * scalefac,
1013           III_scalefac_t      * vbrsf,
1014           int                 * VBRmax )
1015 {
1016     const int * max_rangep;
1017     int sfb, maxsfb;
1018     int maxover, maxover0, maxover1, maxover0p, maxover1p, mover;
1019     int v0, v1, v0p, v1p;
1020     int vbrmax = *VBRmax;
1021
1022     max_rangep = gfc->is_mpeg1 ? max_range_long : max_range_long_lsf_pretab;
1023
1024     maxover0  = 0;
1025     maxover1  = 0;
1026     maxover0p = 0; /* pretab */
1027     maxover1p = 0; /* pretab */
1028
1029     maxsfb = gfc->sfb21_extra ? SBMAX_l : SBPSY_l;
1030     for ( sfb = 0; sfb < maxsfb; ++sfb ) {
1031         v0  = (vbrmax - vbrsf->l[sfb]) - 2*max_range_long[sfb];
1032         v1  = (vbrmax - vbrsf->l[sfb]) - 4*max_range_long[sfb];
1033         v0p = (vbrmax - vbrsf->l[sfb]) - 2*(max_rangep[sfb]+pretab[sfb]);
1034         v1p = (vbrmax - vbrsf->l[sfb]) - 4*(max_rangep[sfb]+pretab[sfb]);
1035         if (maxover0 < v0)
1036             maxover0 = v0;
1037         if (maxover1 < v1)
1038             maxover1 = v1;
1039         if (maxover0p < v0p)
1040             maxover0p = v0p;
1041         if (maxover1p < v1p)
1042             maxover1p = v1p;
1043     }
1044
1045     mover = Min (maxover0, maxover0p);
1046     if ((gfc->noise_shaping == 2) && (gfc->presetTune.use && !(gfc->presetTune.athadjust_safe_noiseshaping || gfc->ATH->adjust < 1.0))) {
1047         /* allow scalefac_scale=1 */
1048         mover = Min (mover, maxover1);
1049         mover = Min (mover, maxover1p);
1050     }
1051
1052     vbrmax    -= mover;
1053     maxover0  -= mover;
1054     maxover0p -= mover;
1055     maxover1  -= mover;
1056     maxover1p -= mover;
1057
1058     if (maxover0 <= 0) {
1059         cod_info->scalefac_scale = 0;
1060         cod_info->preflag = 0;
1061         vbrmax -= maxover0;
1062     } else if (maxover0p <= 0) {
1063         cod_info->scalefac_scale = 0;
1064         cod_info->preflag = 1;
1065         vbrmax -= maxover0p;
1066     } else if (maxover1 == 0) {
1067         cod_info->scalefac_scale = 1;
1068         cod_info->preflag = 0;
1069     } else if (maxover1p == 0) {
1070         cod_info->scalefac_scale = 1;
1071         cod_info->preflag = 1;
1072     } else {
1073         assert(0); /* this should not happen */
1074     }
1075
1076     /* sf =  (cod_info->global_gain-210.0) */
1077     cod_info->global_gain = vbrmax + 210;
1078     assert (cod_info->global_gain < 256);
1079
1080     if (cod_info->global_gain < 0) {
1081         cod_info->global_gain = 0;
1082     }
1083     else
1084     if (cod_info->global_gain > 255)
1085         cod_info->global_gain = 255;
1086
1087     for (sfb = 0; sfb < SBMAX_l; ++sfb)
1088         vbrsf->l[sfb] -= vbrmax;
1089
1090     if ( gfc->is_mpeg1 == 1 )
1091         maxover = compute_scalefacs_long (vbrsf->l, cod_info, scalefac->l);
1092     else
1093         maxover = compute_scalefacs_long_lsf (vbrsf->l, cod_info, scalefac->l);
1094
1095     assert (maxover <= 0);
1096
1097     *VBRmax = vbrmax;
1098 }
1099
1100
1101
1102 /***********************************************************************
1103  *
1104  *  quantize xr34 based on scalefactors
1105  *
1106  *  calc_short_block_xr34
1107  *  calc_long_block_xr34
1108  *
1109  *  Mark Taylor 2000-??-??
1110  *  Robert Hegemann 2000-10-20 made functions of them
1111  *
1112  ***********************************************************************/
1113
1114 static void
1115 short_block_xr34 (
1116     const lame_internal_flags * gfc,
1117     const gr_info             * cod_info,
1118     const III_scalefac_t      * scalefac,
1119     const FLOAT8              * xr34_orig,
1120           FLOAT8              * xr34 )
1121 {
1122     int    sfb, l, j, b;
1123     int    ifac, ifqstep, start, end;
1124     FLOAT8 fac;
1125
1126     /* even though there is no scalefactor for sfb12
1127      * subblock gain affects upper frequencies too, that's why
1128      * we have to go up to SBMAX_s
1129      */
1130     ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
1131     for ( j = 0, sfb = 0; sfb < SBMAX_s; ++sfb ) {
1132         start = gfc->scalefac_band.s[ sfb ];
1133         end   = gfc->scalefac_band.s[ sfb+1 ];
1134         for (b = 0; b < 3; ++b) {
1135             ifac = 8*cod_info->subblock_gain[b]+ifqstep*scalefac->s[sfb][b];
1136
1137             if ( ifac == 0 ) {  /* just copy */
1138                 l = (end-start+7) / 8;
1139                 switch ((end-start) % 8) {
1140                     default:
1141                     case 0: do{ xr34[j] = xr34_orig[j]; ++j;
1142                     case 7:     xr34[j] = xr34_orig[j]; ++j;
1143                     case 6:     xr34[j] = xr34_orig[j]; ++j;
1144                     case 5:     xr34[j] = xr34_orig[j]; ++j;
1145                     case 4:     xr34[j] = xr34_orig[j]; ++j;
1146                     case 3:     xr34[j] = xr34_orig[j]; ++j;
1147                     case 2:     xr34[j] = xr34_orig[j]; ++j;
1148                     case 1:     xr34[j] = xr34_orig[j]; ++j; } while (--l);
1149                 }
1150                 continue;
1151             }
1152             if (ifac < Q_MAX-210)
1153                 fac = IIPOW20_(ifac);
1154             else
1155                 fac = pow (2.0, 0.1875*ifac);
1156
1157             /*
1158              *  loop unrolled into "Duff's Device".  Robert Hegemann
1159              */
1160             l = (end-start+7) / 8;
1161             switch ((end-start) % 8) {
1162                 default:
1163                 case 0: do{ xr34[j] = xr34_orig[j]*fac; ++j;
1164                 case 7:     xr34[j] = xr34_orig[j]*fac; ++j;
1165                 case 6:     xr34[j] = xr34_orig[j]*fac; ++j;
1166                 case 5:     xr34[j] = xr34_orig[j]*fac; ++j;
1167                 case 4:     xr34[j] = xr34_orig[j]*fac; ++j;
1168                 case 3:     xr34[j] = xr34_orig[j]*fac; ++j;
1169                 case 2:     xr34[j] = xr34_orig[j]*fac; ++j;
1170                 case 1:     xr34[j] = xr34_orig[j]*fac; ++j; } while (--l);
1171             }
1172         }
1173     }
1174 }
1175
1176
1177
1178 static void
1179 long_block_xr34 (
1180     const lame_internal_flags * gfc,
1181     const gr_info             * cod_info,
1182     const III_scalefac_t      * scalefac,
1183     const FLOAT8              * xr34_orig,
1184           FLOAT8              * xr34 )
1185 {
1186     int    sfb, l, j;
1187     int    ifac, ifqstep, start, end;
1188     FLOAT8 fac;
1189
1190     ifqstep = ( cod_info->scalefac_scale == 0 ) ? 2 : 4;
1191     for ( sfb = 0; sfb < SBMAX_l; ++sfb ) {
1192
1193         ifac = ifqstep*scalefac->l[sfb];
1194         if (cod_info->preflag)
1195             ifac += ifqstep*pretab[sfb];
1196
1197         start = gfc->scalefac_band.l[ sfb ];
1198         end   = gfc->scalefac_band.l[ sfb+1 ];
1199
1200         if ( ifac == 0 ) {  /* just copy */
1201             j = start;
1202             l = (end-start+7) / 8;
1203             switch ((end-start) % 8) {
1204                 default:
1205                 case 0: do{ xr34[j] = xr34_orig[j]; ++j;
1206                 case 7:     xr34[j] = xr34_orig[j]; ++j;
1207                 case 6:     xr34[j] = xr34_orig[j]; ++j;
1208                 case 5:     xr34[j] = xr34_orig[j]; ++j;
1209                 case 4:     xr34[j] = xr34_orig[j]; ++j;
1210                 case 3:     xr34[j] = xr34_orig[j]; ++j;
1211                 case 2:     xr34[j] = xr34_orig[j]; ++j;
1212                 case 1:     xr34[j] = xr34_orig[j]; ++j; } while (--l);
1213             }
1214             continue;
1215         }
1216         if (ifac < Q_MAX-210)
1217             fac = IIPOW20_(ifac);
1218         else
1219             fac = pow (2.0, 0.1875*ifac);
1220
1221         /*
1222          *  loop unrolled into "Duff's Device".  Robert Hegemann
1223          */
1224         j = start;
1225         l = (end-start+7) / 8;
1226         switch ((end-start) % 8) {
1227             default:
1228             case 0: do{ xr34[j] = xr34_orig[j]*fac; ++j;
1229             case 7:     xr34[j] = xr34_orig[j]*fac; ++j;
1230             case 6:     xr34[j] = xr34_orig[j]*fac; ++j;
1231             case 5:     xr34[j] = xr34_orig[j]*fac; ++j;
1232             case 4:     xr34[j] = xr34_orig[j]*fac; ++j;
1233             case 3:     xr34[j] = xr34_orig[j]*fac; ++j;
1234             case 2:     xr34[j] = xr34_orig[j]*fac; ++j;
1235             case 1:     xr34[j] = xr34_orig[j]*fac; ++j; } while (--l);
1236         }
1237     }
1238 }
1239
1240
1241
1242
1243
1244
1245
1246 /************************************************************************
1247  *
1248  *  VBR_noise_shaping2()
1249  *
1250  *  may result in a need of too many bits, then do it CBR like
1251  *
1252  *  Robert Hegemann 2000-10-25
1253  *
1254  ***********************************************************************/
1255
1256 int
1257 VBR_noise_shaping2 (
1258     lame_global_flags * gfp,
1259     FLOAT8            * xr,
1260     FLOAT8            * xr34orig,
1261     int               * l3_enc,
1262     int                 minbits,
1263     int                 maxbits,
1264     III_scalefac_t    * scalefac,
1265     III_psy_xmin      * l3_xmin,
1266     int                 gr,
1267     int                 ch )
1268 {
1269     lame_internal_flags *gfc = gfp->internal_flags;
1270     III_scalefac_t vbrsf;
1271     III_scalefac_t vbrsf2;
1272     gr_info *cod_info;
1273     FLOAT8 xr34[576];
1274     int shortblock, ret, bits, huffbits;
1275     int vbrmin, vbrmax, vbrmin2, vbrmax2;
1276     int best_huffman = gfc->use_best_huffman;
1277     int count=6;
1278
1279     gfc->use_best_huffman = 0; /* we will do it later */
1280
1281     cod_info   = &gfc->l3_side.gr[gr].ch[ch].tt;
1282     shortblock = (cod_info->block_type == SHORT_TYPE);
1283
1284     if (shortblock) {
1285         short_block_sf (gfc, l3_xmin, xr34orig, xr, &vbrsf2, &vbrmin2, &vbrmax2);
1286     } else {
1287         long_block_sf (gfc, l3_xmin, xr34orig, xr, &vbrsf2, &vbrmin2, &vbrmax2);
1288     }
1289     vbrsf = vbrsf2;
1290     vbrmin = vbrmin2;
1291     vbrmax = vbrmax2;
1292
1293 do {
1294 --count;
1295
1296     if (shortblock) {
1297         short_block_scalefacs (gfc, cod_info, scalefac, &vbrsf, &vbrmax);
1298         short_block_xr34      (gfc, cod_info, scalefac, xr34orig, xr34);
1299     } else {
1300         long_block_scalefacs (gfc, cod_info, scalefac, &vbrsf, &vbrmax);
1301         long_block_xr34      (gfc, cod_info, scalefac, xr34orig, xr34);
1302     }
1303
1304     ret = VBR_quantize_granule (gfc, xr34, l3_enc, scalefac, gr, ch);
1305
1306     if (vbrmin == vbrmax) break;
1307     else if (cod_info->part2_3_length < minbits) {
1308         int i;
1309         vbrmax = vbrmin2 + (vbrmax2-vbrmin2) * count/6;
1310         vbrmin = vbrmin2;
1311         if (shortblock) {
1312             for (i = 0; i < SBMAX_s; ++i) {
1313                 //vbrsf.s[i][0] = vbrmin2 + (vbrsf2.s[i][0]-vbrmin2) * count/6;
1314                 //vbrsf.s[i][1] = vbrmin2 + (vbrsf2.s[i][1]-vbrmin2) * count/6;
1315                 //vbrsf.s[i][2] = vbrmin2 + (vbrsf2.s[i][2]-vbrmin2) * count/6;
1316                 vbrsf.s[i][0] = Min(vbrsf2.s[i][0], vbrmax);
1317                 vbrsf.s[i][1] = Min(vbrsf2.s[i][1], vbrmax);
1318                 vbrsf.s[i][2] = Min(vbrsf2.s[i][2], vbrmax);
1319             }
1320         }
1321         else {
1322             for (i = 0; i < SBMAX_l; ++i) {
1323                 //vbrsf.l[i] = vbrmin2 + (vbrsf2.l[i]-vbrmin2) * count/6;
1324                 vbrsf.l[i] = Min(vbrsf2.l[i], vbrmax);
1325             }
1326         }
1327     }
1328     else if (cod_info->part2_3_length > maxbits) {
1329         int i;
1330         vbrmax = vbrmax2;
1331         vbrmin = vbrmax2 + (vbrmin2-vbrmax2) * count/6;
1332         if (shortblock) {
1333             for (i = 0; i < SBMAX_s; ++i) {
1334                 //vbrsf.s[i][0] = vbrmax2 + (vbrsf2.s[i][0]-vbrmax2) * count/6;
1335                 //vbrsf.s[i][1] = vbrmax2 + (vbrsf2.s[i][1]-vbrmax2) * count/6;
1336                 //vbrsf.s[i][2] = vbrmax2 + (vbrsf2.s[i][2]-vbrmax2) * count/6;
1337                 vbrsf.s[i][0] = Max(vbrsf2.s[i][0], vbrmin);
1338                 vbrsf.s[i][1] = Max(vbrsf2.s[i][1], vbrmin);
1339                 vbrsf.s[i][2] = Max(vbrsf2.s[i][2], vbrmin);
1340             }
1341         }
1342         else {
1343             for (i = 0; i < SBMAX_l; ++i) {
1344                 //vbrsf.l[i] = vbrmax2 + (vbrsf2.l[i]-vbrmax2) * count/6;
1345                 vbrsf.l[i] = Max(vbrsf2.l[i], vbrmin);
1346             }
1347         }
1348     }
1349     else break;
1350 } while(1 && ret != -1);
1351
1352
1353     gfc->use_best_huffman = best_huffman;
1354
1355     if (ret == -1) /* Houston, we have a problem */
1356         return -1;
1357
1358     if (cod_info->part2_3_length < minbits) {
1359         huffbits = minbits - cod_info->part2_length;
1360         bits = bin_search_StepSize (gfc, cod_info, huffbits,
1361                                     gfc->OldValue[ch], xr34, l3_enc);
1362         gfc->OldValue[ch] = cod_info->global_gain;
1363         cod_info->part2_3_length  = bits + cod_info->part2_length;
1364     }
1365     if (cod_info->part2_3_length > maxbits) {
1366         huffbits = maxbits - cod_info->part2_length;
1367         if (huffbits < 0) huffbits = 0;
1368         bits = bin_search_StepSize (gfc, cod_info, huffbits,
1369                                     gfc->OldValue[ch], xr34, l3_enc);
1370         gfc->OldValue[ch] = cod_info->global_gain;
1371         cod_info->part2_3_length = bits;
1372         if (bits > huffbits) {
1373             bits = inner_loop (gfc, cod_info, huffbits, xr34, l3_enc);
1374             cod_info->part2_3_length  = bits;
1375         }
1376         if (bits >= LARGE_BITS) /* Houston, we have a problem */
1377             return -2;
1378         cod_info->part2_3_length += cod_info->part2_length;
1379     }
1380
1381     if (cod_info->part2_length >= LARGE_BITS) /* Houston, we have a problem */
1382         return -2;
1383
1384     assert (cod_info->global_gain < 256);
1385
1386     return 0;
1387 }
1388
1389
1390
1391
1392
1393
1394 /************************************************************************
1395  *
1396  * VBR_noise_shaping()
1397  *
1398  * compute scalefactors, l3_enc, and return number of bits needed to encode
1399  *
1400  * return code:    0   scalefactors were found with all noise < masking
1401  *
1402  *               n>0   scalefactors required too many bits.  global gain
1403  *                     was decreased by n
1404  *                     If n is large, we should probably recompute scalefacs
1405  *                     with a lower quality.
1406  *
1407  *               n<0   scalefactors used less than minbits.
1408  *                     global gain was increased by n.
1409  *                     If n is large, might want to recompute scalefacs
1410  *                     with a higher quality setting?
1411  *
1412  ************************************************************************/
1413 static int
1414 VBR_noise_shaping (
1415     lame_global_flags *gfp,
1416     FLOAT8             xr       [576],
1417     FLOAT8             xr34orig [576],
1418     int                l3_enc   [576],
1419     int                digital_silence,
1420     int                minbits,
1421     int                maxbits,
1422     III_scalefac_t    *scalefac,
1423     III_psy_xmin      *l3_xmin,
1424     int                gr,
1425     int                ch )
1426 {
1427     lame_internal_flags *gfc=gfp->internal_flags;
1428     III_scalefac_t save_sf;
1429     III_scalefac_t vbrsf;
1430     gr_info *cod_info;
1431     FLOAT8 xr34[576];
1432     int shortblock;
1433     int vbrmax, vbrmin;
1434     int global_gain_adjust = 0;
1435
1436     cod_info   = &gfc->l3_side.gr[gr].ch[ch].tt;
1437     shortblock = (cod_info->block_type == SHORT_TYPE);
1438
1439     if (shortblock)
1440         short_block_sf (gfc, l3_xmin, xr34orig, xr, &vbrsf, &vbrmin, &vbrmax);
1441     else
1442         long_block_sf (gfc, l3_xmin, xr34orig, xr, &vbrsf, &vbrmin, &vbrmax);
1443
1444     /* save a copy of vbrsf, incase we have to recomptue scalefacs */
1445     memcpy (&save_sf, &vbrsf, sizeof(III_scalefac_t));
1446
1447     do {
1448         memset (scalefac, 0, sizeof(III_scalefac_t));
1449
1450         if (shortblock) {
1451             short_block_scalefacs (gfc, cod_info, scalefac, &vbrsf, &vbrmax);
1452             short_block_xr34      (gfc, cod_info, scalefac, xr34orig, xr34);
1453         } else {
1454             long_block_scalefacs (gfc, cod_info, scalefac, &vbrsf, &vbrmax);
1455             long_block_xr34      (gfc, cod_info, scalefac, xr34orig, xr34);
1456         }
1457         VBR_quantize_granule (gfc, xr34, l3_enc, scalefac, gr, ch);
1458
1459
1460         /* decrease noise until we use at least minbits
1461          */
1462         if (cod_info->part2_3_length < minbits) {
1463             if (digital_silence) break;
1464             //if (cod_info->part2_3_length == cod_info->part2_length) break;
1465             if (vbrmax+210 == 0) break;
1466
1467             /* decrease global gain, recompute scale factors */
1468             --vbrmax;
1469             --global_gain_adjust;
1470             memcpy (&vbrsf, &save_sf, sizeof(III_scalefac_t));
1471         }
1472
1473     } while (cod_info->part2_3_length < minbits);
1474
1475     /* inject noise until we meet our bit limit
1476      */
1477     while (cod_info->part2_3_length > Min (maxbits, MAX_BITS)) {
1478         /* increase global gain, keep existing scale factors */
1479         ++cod_info->global_gain;
1480         if (cod_info->global_gain > 255)
1481             ERRORF (gfc,"%ld impossible to encode ??? frame! bits=%d\n",
1482                     //  gfp->frameNum, cod_info->part2_3_length);
1483                              -1,       cod_info->part2_3_length);
1484         VBR_quantize_granule (gfc, xr34, l3_enc, scalefac, gr, ch);
1485
1486         ++global_gain_adjust;
1487     }
1488
1489     return global_gain_adjust;
1490 }
1491
1492
1493
1494 void
1495 VBR_quantize(lame_global_flags *gfp,
1496                 FLOAT8 pe[2][2], FLOAT8 ms_ener_ratio[2],
1497                 FLOAT8 xr[2][2][576], III_psy_ratio ratio[2][2],
1498                 int l3_enc[2][2][576],
1499                 III_scalefac_t scalefac[2][2])
1500 {
1501   lame_internal_flags *gfc=gfp->internal_flags;
1502   III_psy_xmin l3_xmin[2][2];
1503   int minbits,maxbits,max_frame_bits,totbits,gr,ch,i,bits_ok;
1504   int bitsPerFrame,mean_bits;
1505   int analog_silence;
1506   FLOAT8 qadjust;
1507   III_side_info_t * l3_side;
1508   gr_info *cod_info;
1509   int digital_silence[2][2];
1510   FLOAT8 masking_lower_db=0;
1511   FLOAT8 xr34[2][2][576];
1512
1513   qadjust=0;   /* start with -1 db quality improvement over quantize.c VBR */
1514
1515   l3_side = &gfc->l3_side;
1516
1517   /* now find out: if the frame can be considered analog silent
1518    *               if each granule can be considered digital silent
1519    * and calculate l3_xmin and the fresh xr34 array
1520    */
1521
1522   assert( gfp->VBR_q <= 9 );
1523   assert( gfp->VBR_q >= 0 );
1524   analog_silence=1;
1525   for (gr = 0; gr < gfc->mode_gr; ++gr) {
1526     /* copy data to be quantized into xr */
1527     if (gfc->mode_ext==MPG_MD_MS_LR) {
1528       ms_convert(xr[gr],xr[gr]);
1529     }
1530     for (ch = 0; ch < gfc->channels_out; ++ch) {
1531       /* if in the following sections the quality would not be adjusted
1532        * then we would only have to call calc_xmin once here and
1533        * could drop subsequently calls (rh 2000/07/17)
1534        */
1535       int over_ath;
1536       cod_info = &l3_side->gr[gr].ch[ch].tt;
1537       cod_info->part2_3_length=LARGE_BITS;
1538
1539       if (cod_info->block_type == SHORT_TYPE) {
1540           cod_info->sfb_lmax = 0; /* No sb*/
1541           cod_info->sfb_smin = 0;
1542       } else {
1543           /* MPEG 1 doesnt use last scalefactor band */
1544           cod_info->sfb_lmax = SBPSY_l;
1545           cod_info->sfb_smin = SBPSY_s;    /* No sb */
1546           if (cod_info->mixed_block_flag) {
1547             cod_info->sfb_lmax        = gfc->is_mpeg1 ? 8 : 6;
1548             cod_info->sfb_smin        = 3;
1549           }
1550       }
1551
1552       /* quality setting */
1553       masking_lower_db = gfc->VBR->mask_adjust;
1554       if (pe[gr][ch]>750) {
1555         masking_lower_db -= Min(10,4*(pe[gr][ch]-750.)/750.);
1556       }
1557       gfc->masking_lower = pow(10.0,masking_lower_db/10);
1558
1559       /* masking thresholds */
1560       over_ath = calc_xmin(gfp,xr[gr][ch],&ratio[gr][ch],cod_info,&l3_xmin[gr][ch]);
1561
1562       /* if there are bands with more energy than the ATH
1563        * then we say the frame is not analog silent */
1564       if (over_ath) {
1565         analog_silence = 0;
1566       }
1567
1568       /* if there is no line with more energy than 1e-20
1569        * then this granule is considered to be digital silent
1570        * plus calculation of xr34 */
1571       digital_silence[gr][ch] = 1;
1572       for(i=0;i<576;++i) {
1573         FLOAT8 temp=fabs(xr[gr][ch][i]);
1574         xr34[gr][ch][i]=sqrt(sqrt(temp)*temp);
1575         digital_silence[gr][ch] &= temp < 1E-20;
1576       }
1577     } /* ch */
1578   }  /* gr */
1579
1580
1581   /* compute minimum allowed bits from minimum allowed bitrate */
1582   if (analog_silence) {
1583     gfc->bitrate_index=1;
1584   } else {
1585     gfc->bitrate_index=gfc->VBR_min_bitrate;
1586   }
1587   getframebits(gfp, &bitsPerFrame, &mean_bits);
1588   minbits = (mean_bits/gfc->channels_out);
1589
1590   /* compute maximum allowed bits from max allowed bitrate */
1591   gfc->bitrate_index=gfc->VBR_max_bitrate;
1592   getframebits(gfp, &bitsPerFrame, &mean_bits);
1593   max_frame_bits = ResvFrameBegin(gfp, l3_side, mean_bits, bitsPerFrame);
1594   maxbits=2.5*(mean_bits/gfc->channels_out);
1595
1596   {
1597   /* compute a target  mean_bits based on compression ratio
1598    * which was set based on VBR_q
1599    */
1600   int bit_rate = gfp->out_samplerate*16*gfc->channels_out/(1000.0*gfp->compression_ratio);
1601   bitsPerFrame = (bit_rate*gfp->framesize*1000)/gfp->out_samplerate;
1602   mean_bits = (bitsPerFrame - 8*gfc->sideinfo_len) / gfc->mode_gr;
1603   }
1604
1605
1606   minbits = Max(minbits,125);
1607   minbits=Max(minbits,.40*(mean_bits/gfc->channels_out));
1608   maxbits=Min(maxbits,2.5*(mean_bits/gfc->channels_out));
1609
1610
1611
1612
1613
1614
1615
1616   /*
1617    * loop over all ch,gr, encoding anything with bits > .5*(max_frame_bits/4)
1618    *
1619    * If a particular granule uses way too many bits, it will be re-encoded
1620    * on the next iteration of the loop (with a lower quality setting).
1621    * But granules which dont use
1622    * use too many bits will not be re-encoded.
1623    *
1624    * minbits:  minimum allowed bits for 1 granule 1 channel
1625    * maxbits:  maximum allowwed bits for 1 granule 1 channel
1626    * max_frame_bits:  maximum allowed bits for entire frame
1627    * (max_frame_bits/4)   estimate of average bits per granule per channel
1628    *
1629    */
1630
1631   do {
1632
1633     totbits=0;
1634     for (gr = 0; gr < gfc->mode_gr; ++gr) {
1635       int minbits_lr[2];
1636       minbits_lr[0]=minbits;
1637       minbits_lr[1]=minbits;
1638
1639 #if 0
1640       if (gfc->mode_ext==MPG_MD_MS_LR) {
1641         FLOAT8 fac;
1642         fac = .33*(.5-ms_ener_ratio[gr])/.5;
1643         if (fac<0) fac=0;
1644         if (fac>.5) fac=.5;
1645         minbits_lr[0] = (1+fac)*minbits;
1646         minbits_lr[1] = Max(125,(1-fac)*minbits);
1647       }
1648 #endif
1649
1650
1651       for (ch = 0; ch < gfc->channels_out; ++ch) {
1652         int adjusted,shortblock;
1653         cod_info = &l3_side->gr[gr].ch[ch].tt;
1654
1655         /* ENCODE this data first pass, and on future passes unless it uses
1656          * a very small percentage of the max_frame_bits  */
1657         if (cod_info->part2_3_length > (max_frame_bits/(2*gfc->channels_out*gfc->mode_gr))) {
1658
1659           shortblock = (cod_info->block_type == SHORT_TYPE);
1660
1661           /* Adjust allowed masking based on quality setting */
1662           if (qadjust!=0 /*|| shortblock*/) {
1663             masking_lower_db = gfc->VBR->mask_adjust + qadjust;
1664
1665             /*
1666             if (shortblock) masking_lower_db -= 4;
1667             */
1668
1669             if (pe[gr][ch]>750)
1670               masking_lower_db -= Min(10,4*(pe[gr][ch]-750.)/750.);
1671             gfc->masking_lower = pow(10.0,masking_lower_db/10);
1672             calc_xmin( gfp, xr[gr][ch], ratio[gr]+ch, cod_info, l3_xmin[gr]+ch);
1673           }
1674
1675           /* digital silent granules do not need the full round trip,
1676            * but this can be optimized later on
1677            */
1678           adjusted = VBR_noise_shaping (gfp,xr[gr][ch],xr34[gr][ch],
1679                                         l3_enc[gr][ch],
1680                                         digital_silence[gr][ch],
1681                                         minbits_lr[ch],
1682                                         maxbits,scalefac[gr]+ch,
1683                                         l3_xmin[gr]+ch,gr,ch);
1684           if (adjusted>10) {
1685             /* global_gain was changed by a large amount to get bits < maxbits */
1686             /* quality is set to high.  we could set bits = LARGE_BITS
1687              * to force re-encoding.  But most likely the other channels/granules
1688              * will also use too many bits, and the entire frame will
1689              * be > max_frame_bits, forcing re-encoding below.
1690              */
1691             // cod_info->part2_3_bits = LARGE_BITS;
1692           }
1693         }
1694         totbits += cod_info->part2_3_length;
1695       }
1696     }
1697     bits_ok=1;
1698     if (totbits>max_frame_bits) {
1699       /* lower quality */
1700       qadjust += Max(.125,Min(1,(totbits-max_frame_bits)/300.0));
1701       /* adjusting minbits and maxbits is necessary too
1702        * cos lowering quality is not enough in rare cases
1703        * when each granule still needs almost maxbits, it wont fit */
1704       minbits = Max(125,minbits*0.975);
1705       maxbits = Max(minbits,maxbits*0.975);
1706       //      DEBUGF("%i totbits>max_frame_bits   totbits=%i  maxbits=%i \n",gfp->frameNum,totbits,max_frame_bits);
1707       //      DEBUGF("next masking_lower_db = %f \n",masking_lower_db + qadjust);
1708       bits_ok=0;
1709     }
1710
1711   } while (!bits_ok);
1712
1713
1714
1715   /* find optimal scalefac storage.  Cant be done above because
1716    * might enable scfsi which breaks the interation loops */
1717   totbits=0;
1718   for (gr = 0; gr < gfc->mode_gr; ++gr) {
1719     for (ch = 0; ch < gfc->channels_out; ++ch) {
1720       best_scalefac_store(gfc, gr, ch, l3_enc, l3_side, scalefac);
1721       totbits += l3_side->gr[gr].ch[ch].tt.part2_3_length;
1722     }
1723   }
1724
1725
1726
1727
1728   if (analog_silence && !gfp->VBR_hard_min) {
1729     gfc->bitrate_index = 1;
1730   } else {
1731     gfc->bitrate_index = gfc->VBR_min_bitrate;
1732   }
1733   for( ; gfc->bitrate_index < gfc->VBR_max_bitrate; ++gfc->bitrate_index ) {
1734
1735     getframebits (gfp, &bitsPerFrame, &mean_bits);
1736     maxbits = ResvFrameBegin(gfp, l3_side, mean_bits, bitsPerFrame);
1737     if (totbits <= maxbits) break;
1738   }
1739   if (gfc->bitrate_index == gfc->VBR_max_bitrate) {
1740     getframebits (gfp, &bitsPerFrame, &mean_bits);
1741     maxbits = ResvFrameBegin(gfp, l3_side, mean_bits, bitsPerFrame);
1742   }
1743
1744   //  DEBUGF("%i total_bits=%i max_frame_bits=%i index=%i  \n",gfp->frameNum,totbits,max_frame_bits,gfc->bitrate_index);
1745
1746   for (gr = 0; gr < gfc->mode_gr; ++gr) {
1747     for (ch = 0; ch < gfc->channels_out; ++ch) {
1748       cod_info = &l3_side->gr[gr].ch[ch].tt;
1749
1750
1751       ResvAdjust (gfc, cod_info, l3_side, mean_bits);
1752
1753       /*******************************************************************
1754        * set the sign of l3_enc from the sign of xr
1755        *******************************************************************/
1756       for ( i = 0; i < 576; ++i) {
1757         if (xr[gr][ch][i] < 0) l3_enc[gr][ch][i] *= -1;
1758       }
1759     }
1760   }
1761   ResvFrameEnd (gfc, l3_side, mean_bits);
1762
1763
1764
1765 }
1766
1767
1768
1769