/********************************************************************** * ISO MPEG Audio Subgroup Software Simulation Group (1996) * ISO 13818-3 MPEG-2 Audio Encoder - Lower Sampling Frequency Extension * * $Id: loop.c,v 1.2 1997/01/19 22:28:29 rowlands Exp $ * * $Log: loop.c,v $ * Revision 1.2 1997/01/19 22:28:29 rowlands * Layer 3 bug fixes from Seymour Shlien * * Revision 1.1 1996/02/14 04:04:23 rowlands * Initial revision * * Received from Mike Coleman **********************************************************************/ #include #include #include #include #include "globalflags.h" #include "common.h" #include "l3side.h" #include "loop.h" #include "huffman.h" #include "l3bitstream.h" #include "reservoir.h" #include "loop-pvt.h" #include "gtkanal.h" /* #define DEBUG */ /* #define DEBUGSC */ static int convert_mdct, convert_psy, reduce_sidechannel; /* mt 5/99. These global flags denote 4 possibilities: 1 MDCT input L/R, quantize L/R, psy-model thresholds: L/R -m s 2 MDCT input L/R, quantize M/S, psy-model thresholds: L/R -m j 3 MDCT input M/S, quantize M/S, psy-model thresholds: M/S -m f 4 MDCT input L/R, quantize M/S, psy-model thresholds: M/S coming soon 1: convert_mdct = 0, convert_psy=0, reduce_sidechannel=0 2: convert_mdct = 1, convert_psy=1, reduce_sidechannel=1 3: convert_mdct = 0, convert_psy=0, reduce_sidechannel=1 4: convert_mdct = 1, convert_psy=0, reduce_sidechannel=1 if (convert_mdct), then iteration_loop will quantize M/S data from the L/R input MDCT coefficients. if (convert_psy), then calc_noise will compute the noise for the L/R channels from M/S MDCT data and L/R psy-model threshold information. Distortion in ether L or R channel will be marked as distortion in both Mid and Side channels. if (reduce_sidechannel) then outer_loop will allocate less bits to the side channel and more bits to the mid channel based on relative energies. */ /* Here are MPEG1 Table B.8 and MPEG2 Table B.1 -- Layer III scalefactor bands. Index into this using a method such as: idx = fr_ps->header->sampling_frequency + (fr_ps->header->version * 3) */ struct scalefac_struct sfBandIndex[6] = { { /* Table B.2.b: 22.05 kHz */ {0,6,12,18,24,30,36,44,54,66,80,96,116,140,168,200,238,284,336,396,464,522,576}, {0,4,8,12,18,24,32,42,56,74,100,132,174,192} }, { /* Table B.2.c: 24 kHz */ {0,6,12,18,24,30,36,44,54,66,80,96,114,136,162,194,232,278,330,394,464,540,576}, {0,4,8,12,18,26,36,48,62,80,104,136,180,192} }, { /* Table B.2.a: 16 kHz */ {0,6,12,18,24,30,36,44,54,66,80,96,116,140,168,200,238,248,336,396,464,522,576}, {0,4,8,12,18,26,36,48,62,80,104,134,174,192} }, { /* Table B.8.b: 44.1 kHz */ {0,4,8,12,16,20,24,30,36,44,52,62,74,90,110,134,162,196,238,288,342,418,576}, {0,4,8,12,16,22,30,40,52,66,84,106,136,192} }, { /* Table B.8.c: 48 kHz */ {0,4,8,12,16,20,24,30,36,42,50,60,72,88,106,128,156,190,230,276,330,384,576}, {0,4,8,12,16,22,28,38,50,64,80,100,126,192} }, { /* Table B.8.a: 32 kHz */ {0,4,8,12,16,20,24,30,36,44,54,66,82,102,126,156,194,240,296,364,448,550,576}, {0,4,8,12,16,22,30,42,58,78,104,138,180,192} } }; /* The following table is used to implement the scalefactor partitioning for MPEG2 as described in section 2.4.3.2 of the IS. The indexing corresponds to the way the tables are presented in the IS: [table_number][row_in_table][column of nr_of_sfb] */ static unsigned nr_of_sfb_block[6][3][4] = { { {6, 5, 5, 5}, {9, 9, 9, 9}, {6, 9, 9, 9} }, { {6, 5, 7, 3}, {9, 9, 12, 6}, {6, 9, 12, 6} }, { {11, 10, 0, 0}, {18, 18, 0, 0}, {15,18,0,0} }, { {7, 7, 7, 0}, {12, 12, 12, 0}, {6, 15, 12, 0} }, { {6, 6, 6, 3}, {12, 9, 9, 6}, {6, 12, 9, 6} }, { {8, 8, 5, 0}, {15,12,9,0}, {6,18,9,0} } }; /* table of largest scalefactors for MPEG2 */ static unsigned max_sfac_tab[6][4] = { {4, 4, 3, 3}, {4, 4, 3, 0}, {3, 2, 0, 0}, {4, 5, 5, 0}, {3, 3, 3, 0}, {2, 2, 0, 0} }; /* Table B.6: layer3 preemphasis */ int pretab[21] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 3, 3, 3, 2 }; /* This is the scfsi_band table from 2.4.2.7 of the IS */ int scfsi_band_long[5] = { 0, 6, 11, 16, 21 }; int *scalefac_band_long = &sfBandIndex[3].l[0]; int *scalefac_band_short = &sfBandIndex[3].s[0]; /* convert from L/R <-> Mid/Side */ void ms_convert(double xr[2][576],double xr_org[2][576]) { int i; for ( i = 0; i < 576; i++ ) { xr[0][i] = (xr_org[0][i]+xr_org[1][i])/SQRT2; xr[1][i] = (xr_org[0][i]-xr_org[1][i])/SQRT2; } } /************************************************************************ * allocate bits among the 4 granules based on PE * mt 6/99 ************************************************************************/ int VBR_on_pe(frame_params *fr_ps, layer *info, III_side_info_t *l3_side, int VBRbits[2][2], double pe[][2], int *mean_bits) { int nshort=0; int index,gr,ch; int bitsPerFrame; int stereo = fr_ps->stereo; int fullframebits; int mode_gr = (info->version == 1) ? 2 : 1; getframebits(info,stereo,&bitsPerFrame,mean_bits); fullframebits= ResvFrameBegin( fr_ps, l3_side, *mean_bits, bitsPerFrame ); for (gr=0;grgr[gr].ch[ch].tt.block_type==2); if (shortblock) nshort ++; } index = nshort; if (index>3) index=3; /* impose a minimum bit rate based on number of short blocks */ if (index>0) { /* nshort=1: min 160kbs. (1 shortblock granule) * nshort=2: min 190kbs. (2 shortblock granules) * nshort=3: min 220kbs (3 shortblock granules) * nshort=4: min 256kbs (4 shortblock granules) */ if (info->bitrate_index < 9 + index) info->bitrate_index = 9+index; if (info->bitrate_index > 14) info->bitrate_index=14; getframebits(info,stereo,&bitsPerFrame,mean_bits); fullframebits= ResvFrameBegin( fr_ps, l3_side, *mean_bits, bitsPerFrame ); } /* allocate a minimum 125 bits per channel (approx 32kbs) */ for (gr=0;gr 0 ) { double pe2[2][2],pe_tot=.01; for (gr=0;grgr[gr].ch[ch].tt.block_type==2); pe2[gr][ch] = Min(pe[gr][ch],2800.); if (shortblock) pe2[gr][ch] = Max(pe2[gr][ch],1800.0); pe_tot += pe2[gr][ch]; } } for (gr=0;gr do not decrease bit rate further */ if (index) return (info->bitrate_index <= 9 + index); else return 0; } /************************************************************************ * increase or decrease frame bitrate based on max_over = number of scale * factor bands with audible distortion. If bitrate in increased, allocate * the additional bits among the granules based on "over" * mt 6/99 ************************************************************************/ void VBR_on_over(frame_params *fr_ps, layer *info, III_side_info_t *l3_side, int VBRbits[2][2], double pe[2][2], int max_over, int over[2][2], int *mean_bits, int *VBR_no_decrease) { int extrabits=0; int gr,ch; int bitsPerFrame; int stereo = fr_ps->stereo; int fullframebits; int mode_gr = (info->version == 1) ? 2 : 1; if (max_over <= VBR_q) { info->bitrate_index --; if (info->bitrate_index < VBR_min_bitrate) info->bitrate_index=VBR_min_bitrate; *VBR_no_decrease = VBR_on_pe(fr_ps,info,l3_side,VBRbits,pe,mean_bits); return; } getframebits(info,stereo,&bitsPerFrame,mean_bits); fullframebits= ResvFrameBegin( fr_ps, l3_side, *mean_bits, bitsPerFrame ); info->bitrate_index ++; if (info->bitrate_index > VBR_max_bitrate) info->bitrate_index=VBR_max_bitrate; getframebits(info,stereo,&bitsPerFrame,mean_bits); extrabits= ResvFrameBegin( fr_ps, l3_side, *mean_bits, bitsPerFrame ); extrabits -= fullframebits; /* divide extra bits based on OVER */ if (extrabits > 0) { double over2[2][2],over_tot=0; for (gr=0;gr4095) VBRbits[gr][ch]=4095; } } } /************************************************************************/ /* iteration_loop() */ /************************************************************************/ void iteration_loop( double pe[][2], double ms_ener_ratio[2], double xr_org[2][2][576], III_psy_ratio *ratio, III_side_info_t *l3_side, int l3_enc[2][2][576], III_scalefac_t *scalefac, frame_params *fr_ps) { static int firstcall = 1; III_psy_xmin l3_xmin; gr_info *cod_info; layer *info; double * pxr; int *pi; int VBRbits[2][2]; int best_over[2][2]; int bitsPerFrame; int mean_bits; int stereo = fr_ps->stereo; int ch, gr, i, mode_gr, bit_rate; int VBR_no_decrease; int samp_freq; double xr[2][2][576]; if ( firstcall ) { l3_side->main_data_begin = 0; memset((char *) &l3_xmin, 0, sizeof(l3_xmin)); firstcall = 0; } l3_side->resvDrain = 0; info = fr_ps->header; mode_gr = (info->version == 1) ? 2 : 1; bit_rate = bitrate[info->version][info->lay-1][info->bitrate_index]; samp_freq = s_freq[info->version][info->sampling_frequency]*1000 +.5; scalefac_band_long = &sfBandIndex[info->sampling_frequency + (info->version * 3)].l[0]; scalefac_band_short = &sfBandIndex[info->sampling_frequency + (info->version * 3)].s[0]; convert_mdct=0; convert_psy=0; reduce_sidechannel=0; if (info->mode_ext==2) { convert_mdct = 1; convert_psy = 1; reduce_sidechannel=1; } if (force_ms) { convert_mdct = 0; convert_psy = 0; reduce_sidechannel=1; } /* remove coefficients in scalefactor band 21 (12 for short blocks) * FhG does this for bps <= 128kbs, so we will too. * This amounds to a 16kHz low-pass filter. If that offends you, you * probably should not be encoding at 128kbs! * There is no ratio[21] or xfsf[21], so when these coefficients are * included they are just quantized as is. mt 5/99 */ if (sfb21 && (samp_freq>33000)) { if (bit_rate <= 64*stereo) { for (ch =0 ; ch < 2 ; ch++) for ( gr = 0; gr < mode_gr; gr++ ) { int shortblock = (l3_side->gr[gr].ch[0].tt.block_type==2); if (shortblock) { int j; for (j=0; j<3; j++) { int start = scalefac_band_short[ SFB_SMAX-1 ]; for ( i = start; i < 192; i++ ) { int i0 = 3*i+j; xr_org[gr][ch][i0]=0; } } }else{ int start = scalefac_band_long[ SFB_LMAX-1 ]; memset (xr_org[gr][ch],0,576-start); } } } } /* some intializations. */ for ( gr = 0; gr < mode_gr; gr++ ){ for ( ch = 0; ch < stereo; ch++ ){ cod_info = (gr_info *) &(l3_side->gr[gr].ch[ch]); gr_deco(cod_info); } } /* dont bother with scfsi. */ memset(l3_side->scfsi,0,4); VBR_no_decrease=0; if (VBR) { VBR_no_decrease=VBR_on_pe(fr_ps,info,l3_side,VBRbits,pe,&mean_bits); }else{ getframebits(info,stereo,&bitsPerFrame,&mean_bits); ResvFrameBegin( fr_ps, l3_side, mean_bits, bitsPerFrame ); } /* quantize! */ do { int max_over; for ( gr = 0; gr < mode_gr; gr++ ) outer_loop( xr, xr_org, mean_bits, VBRbits, bit_rate, best_over[gr], &l3_xmin,l3_enc, fr_ps, scalefac,gr,stereo, l3_side, ratio, pe, ms_ener_ratio); if (!VBR) break; /* see if we should try a higher bitrate quantization */ max_over=0; for ( gr = 0; gr < mode_gr; gr++ ) for (ch=0 ; ch < stereo ; ch ++ ) { if (best_over[gr][ch]>max_over) max_over=best_over[gr][ch]; } if ((max_over <= VBR_q) && (info->bitrate_index <= VBR_min_bitrate)) break; if ((max_over > VBR_q) && (info->bitrate_index >= VBR_max_bitrate)) break; /* do not decrease bit rate if we have ever increased bit rate * otherwise infinite loop */ if ((max_over <= VBR_q) && (VBR_no_decrease)) break; /* compute a new bitrate, and allocate bits between granules */ if (max_over > VBR_q) VBR_no_decrease=1; VBR_on_over(fr_ps,info,l3_side,VBRbits,pe,max_over,best_over,&mean_bits,&VBR_no_decrease); } while (1); /* update reservoir status after FINAL quantization/bitrate */ for ( gr = 0; gr < mode_gr; gr++ ) { for (ch=0 ; ch < stereo ; ch ++ ) { gr_info *cod_info = &l3_side->gr[gr].ch[ch].tt; ResvAdjust( fr_ps, cod_info, l3_side, mean_bits ); } } /* set the sign of l3_enc */ for ( gr = 0; gr < mode_gr; gr++ ) for ( ch = 0; ch < stereo; ch++ ) { pi= &l3_enc[gr][ch][0]; pxr=&xr[gr][ch][0]; for ( i = 0; i < 576; i++, pxr++, pi++) if ( (*pxr < 0) && (*pi > 0) ) *pi = -*pi; } ResvFrameEnd( fr_ps, l3_side, mean_bits ); } /************************************************************************/ /* init_outer_loop mt 6/99 */ /************************************************************************/ void init_outer_loop( double xr[2][2][576], /* could be L/R OR MID/SIDE */ double xr_org[2][2][576], III_psy_xmin *l3_xmin, /* the allowed distortion of the scalefactor */ III_scalefac_t *scalefac, /* scalefactors */ int gr, int stereo, III_side_info_t *l3_side, III_psy_ratio *ratio) { int ch,sfb,i; gr_info *cod_info[2]; cod_info[0] = &l3_side->gr[gr].ch[0].tt; cod_info[1] = &l3_side->gr[gr].ch[1].tt; /* copy data to be quantized into xr */ if (convert_mdct) ms_convert(xr[gr],xr_org[gr]); else memcpy(xr[gr],xr_org[gr],sizeof(double)*2*576); /* compute max allowed distortion */ for ( ch = 0; ch < stereo; ch++ ){ calc_xmin( xr_org, ratio, cod_info[ch], l3_xmin, gr, ch ); } for (ch=0 ; ch < stereo ; ch ++ ){ /* if ( info->version == 1 ) calc_scfsi( xr[gr][ch], l3_side, &l3_xmin, ch, gr ); */ /* reset of iteration variables */ for ( sfb = 0; sfb < SFB_LMAX; sfb++ ) scalefac->l[gr][ch][sfb] = 0; for ( sfb = 0; sfb < SFB_SMAX; sfb++ ) for ( i = 0; i < 3; i++ ) scalefac->s[gr][ch][sfb][i] = 0; for ( i = 0; i < 4; i++ ) cod_info[ch]->slen[i] = 0; cod_info[ch]->sfb_partition_table = &nr_of_sfb_block[0][0][0]; cod_info[ch]->part2_3_length = 0; cod_info[ch]->big_values = 0; cod_info[ch]->count1 = 0; cod_info[ch]->scalefac_compress = 0; cod_info[ch]->table_select[0] = 0; cod_info[ch]->table_select[1] = 0; cod_info[ch]->table_select[2] = 0; cod_info[ch]->subblock_gain[0] = 0; cod_info[ch]->subblock_gain[1] = 0; cod_info[ch]->subblock_gain[2] = 0; cod_info[ch]->region0_count = 0; cod_info[ch]->region1_count = 0; cod_info[ch]->part2_length = 0; cod_info[ch]->preflag = 0; cod_info[ch]->scalefac_scale = 0; cod_info[ch]->quantizerStepSize = 0.0; cod_info[ch]->count1table_select= 0; cod_info[ch]->address1 = 0; cod_info[ch]->address2 = 0; cod_info[ch]->address3 = 0; } #define SBGAINXXX #ifdef SBGAIN /* compute subblock gains */ for (ch=0 ; ch < stereo ; ch ++) { int j,b; double en[3],mx; if ((cod_info[ch]->block_type ==2) && count[ch]) { /* estimate energy within each subblock */ for (b=0; b<3; b++) en[b]=0; for ( i=0,j = 0; j < 192; j++ ) { for (b=0; b<3; b++) { en[b]+=xr[gr][ch][i]*xr[gr][ch][i]; i++; } } mx = 1e-12; for (b=0; b<3; b++) mx=Max(mx,en[b]); for (b=0; b<3; b++) en[b] = Max(en[b],1e-12)/mx; printf("ener = %4.2f %4.2f %4.2f \n",en[0],en[1],en[2]); /* pick gain so that 2^(2gain)*en[0] = 1 */ /* gain = .5* log( 1/en[0] )/log(2) = -.5*log(en[])/log(2) */ for (b=0; b<3; b++) { cod_info[ch]->subblock_gain[b]=nint(-.5*log(en[b])/log(2.0)); if (cod_info[ch]->subblock_gain[b] > 2) cod_info[ch]->subblock_gain[b]=2; if (cod_info[ch]->subblock_gain[b] < 0) cod_info[ch]->subblock_gain[b]=0; } } } #endif } /************************************************************************/ /* outer_loop */ /************************************************************************/ /* Function: The outer iteration loop controls the masking conditions */ /* of all scalefactorbands. It computes the best scalefac and */ /* global gain. This module calls the inner iteration loop * * mt 5/99 completely rewritten to allow for bit reservoir control, * mid/side channels with L/R or mid/side masking thresholds, * and chooses best quantization instead of last quantization when * no distortion free quantization can be found. * * added VBR support mt 5/99 ************************************************************************/ void outer_loop( double xr[2][2][576], /* could be L/R OR MID/SIDE */ double xr_org[2][2][576], int mean_bits, int VBRbits[2][2], int bit_rate, int best_over[2], III_psy_xmin *l3_xmin, /* the allowed distortion of the scalefactor */ int l3_enc[2][2][576], /* vector of quantized values ix(0..575) */ frame_params *fr_ps, III_scalefac_t *scalefac, /* scalefactors */ int gr, int stereo, III_side_info_t *l3_side, III_psy_ratio *ratio, double pe[2][2], double ms_ener_ratio[2]) { int status[2],notdone[2]={0,0},count[2]={0,0},bits_found[2]; int targ_bits[2],real_bits[2],tbits,extra_bits; int scalesave_l[2][CBLIMIT], scalesave_s[2][CBLIMIT][3]; int sfb, bits, huff_bits, save_preflag[2], save_compress[2]; double xfsf[2][4][CBLIMIT]; double xrpow[2][2][576],temp; int distort[2][4][CBLIMIT]; int save_l3_enc[2][576]; int save_real_bits[2]; int i,over[2], iteration, ch, try_scale; int better[2]; double tot_noise[2]; double best_noise[2]; gr_info save_cod_info[2]; gr_info *cod_info[2]; cod_info[0] = &l3_side->gr[gr].ch[0].tt; cod_info[1] = &l3_side->gr[gr].ch[1].tt; for (ch=0 ; ch < stereo ; ch ++) best_over[ch] = 100; init_outer_loop(xr,xr_org,l3_xmin,scalefac,gr,stereo,l3_side,ratio); for (ch=0 ; ch < stereo ; ch ++) { count[ch]=0; for (i=0; i<576; i++) { if ( fabs(xr[gr][ch][i]) > 0 ) count[ch]++; } notdone[ch]=count[ch]; if (count[ch]==0) best_over[ch]=0; } /****************************************************************** * allocate bits for each channel ******************************************************************/ if (VBR) { for (ch=0 ; ch < stereo ; ch ++ ) targ_bits[ch]=VBRbits[gr][ch]; }else { int add_bits[2]; double bits_needed; /* allocate targ_bits for granule */ ResvMaxBits2( mean_bits, &tbits, &extra_bits, gr); for (ch=0 ; ch < stereo ; ch ++ ) targ_bits[ch]=tbits/stereo; // allocate extra bits from reservoir based on PE bits=0; for (ch=0; ch 700 */ add_bits[ch]=(pe[gr][ch]-700)/2.0; /* 3.0; */ if (convert_psy) bits_needed = (pe[gr][0]+pe[gr][1]-1400)/2.0; /* short blocks need extra, no matter what the pe */ if (cod_info[ch]->block_type==2) if (add_bits[ch]<500) add_bits[ch]=500; if (add_bits[ch] < 0) add_bits[ch]=0; bits += add_bits[ch]; } for (ch=0; ch extra_bits) add_bits[ch] = (extra_bits*add_bits[ch])/bits; targ_bits[ch] = targ_bits[ch] + add_bits[ch]; } for (ch=0; ch 125) { targ_bits[0] += targ_bits[1]*fac; targ_bits[1] -= targ_bits[1]*fac; } } /* dont allow to many bits per channel */ for (ch=0; ch max_bits) { extra_bits += (targ_bits[ch] - max_bits); targ_bits[ch] = max_bits; } } /* BEGIN MAIN LOOP */ iteration = 0; while ( (notdone[0] || notdone[1]) ) { int pre_just_turned_on[2]; iteration ++; if (iteration==1) { /* compute initial quantization step */ for (ch=0 ; ch < stereo ; ch ++ ) if (count[ch]) { for(i=0;i<576;i++) { temp=fabs(xr[gr][ch][i]); xrpow[gr][ch][i]=sqrt(sqrt(temp)*temp); } bits_found[ch]=bin_search_StepSize2(targ_bits[ch],-211.0,46, l3_enc[gr][ch],xr[gr][ch],xrpow[gr][ch],cod_info[ch]); } } for (ch=0 ; ch < stereo ; ch ++ ) { /* inner_loop starts with the initial quantization step computed above * and slowly increases until the bits < huff_bits. * Thus is it important not to start with too large of an inital * quantization step. Too small is ok, but inner_loop will take longer */ for (ch=0 ; ch < stereo ; ch ++ ) { if (notdone[ch]) { cod_info[ch]->part2_length = part2_length( scalefac, fr_ps, gr, ch, l3_side ); huff_bits = targ_bits[ch] - cod_info[ch]->part2_length; if (huff_bits < 0) { if (iteration==1) { fprintf(stderr,"ERROR: outer_loop(): huff_bits < 0. \n"); exit(-5); }else{ /* scale factors too large, not enough bits. use previous quantizaton */ notdone[ch]=0; over[ch]=999; } }else{ /* if this is the first iteration, see if we can reuse the quantization * computed in bin_search_StepSize above */ if (iteration==1) { if(bits_found[ch]>huff_bits) { cod_info[ch]->quantizerStepSize+=1.0; real_bits[ch] = inner_loop( xr, xrpow, l3_enc, huff_bits, cod_info[ch], gr, ch ); } else real_bits[ch]=bits_found[ch]; } else real_bits[ch]=inner_loop( xr, xrpow, l3_enc, huff_bits, cod_info[ch], gr, ch ); } } } } if (fast_mode) { for (ch=0; chl[gr][ch][sfb]; for ( sfb = 0; sfb < SFB_SMAX; sfb++ ) for ( i = 0; i < 3; i++ ) scalesave_s[ch][sfb][i] = scalefac->s[gr][ch][sfb][i]; save_preflag[ch] = cod_info[ch]->preflag; save_compress[ch] = cod_info[ch]->scalefac_compress; memcpy(save_l3_enc[ch],l3_enc[gr][ch],sizeof(l3_enc[gr][ch])); memcpy(&save_cod_info[ch],cod_info[ch],sizeof(save_cod_info[ch])); save_real_bits[ch]=real_bits[ch]; #ifdef HAVEGTK if (gtkflag) { for ( i = 0; i < 3; i++ ) { for ( sfb = cod_info[ch]->sfb_smax; sfb < 12; sfb++ ) { pinfo->xfsf_s[gr][ch][3*sfb+i] = pinfo->thr_s[gr][ch][3*sfb+i]*xfsf[ch][i+1][sfb]/ (1e-20+l3_xmin->s[gr][ch][sfb][i]); } } for ( sfb = 0; sfb < cod_info[ch]->sfb_lmax; sfb++ ) { pinfo->xfsf[gr][ch][sfb] = pinfo->thr[gr][ch][sfb]*xfsf[ch][0][sfb]/ (1e-20 + l3_xmin->l[gr][ch][sfb]); } pinfo->over[gr][ch]=over[ch]; pinfo->noise[gr][ch]=tot_noise[ch]; } #endif } } } /* if no bands with distortion, we are done */ for (ch=0 ; ch < stereo ; ch ++ ) if (notdone[ch]) { if (convert_psy) notdone[ch] = (over[0] || over[1]); else notdone[ch] = over[ch]; } /* see if we should apply preemphasis */ for (ch=0 ; ch < stereo ; ch ++ ) { pre_just_turned_on[ch]=0; if (notdone[ch]) pre_just_turned_on[ch]= preemphasis2( xr[gr][ch], xrpow[gr][ch], l3_xmin, gr, ch, l3_side, distort[ch]); } /* if we didn't just apply pre-emph, let us see if we should * amplify some scale factor bands */ for (ch=0 ; ch < stereo ; ch ++ ) if (notdone[ch] && (!pre_just_turned_on[ch]) ) { amp_scalefac_bands2( xr[gr][ch], xrpow[gr][ch], l3_xmin, l3_side, scalefac, gr, ch, iteration,distort[ch]); } /* check to make sure we have not amplified too much */ try_scale=0; for (ch=0 ; ch < stereo ; ch ++ ) { if (notdone[ch]) { if ( (status[ch] = loop_break(scalefac, cod_info[ch], gr, ch)) == 0 ) { if ( fr_ps->header->version == 1 ) { status[ch] = scale_bitcount( scalefac, cod_info[ch], gr, ch ); if (status[ch] && (cod_info[ch]->scalefac_scale==0)) try_scale=1; }else{ status[ch] = scale_bitcount_lsf( scalefac, cod_info[ch], gr, ch ); if (status[ch] && (cod_info[ch]->scalefac_scale==0)) try_scale=1; } } notdone[ch] = !status[ch]; } } if (try_scale && experimentalY) { init_outer_loop(xr,xr_org,l3_xmin,scalefac,gr,stereo,l3_side,ratio); for (ch=0 ; ch < stereo ; ch ++ ) { iteration=0; notdone[ch]=1; cod_info[ch]->scalefac_scale=1; } } } /* done with main iteration */ /* restore some data */ for (ch=0 ; ch < stereo ; ch ++ ) { if (count[ch] ) { cod_info[ch]->preflag = save_preflag[ch]; cod_info[ch]->scalefac_compress = save_compress[ch]; for ( sfb = 0; sfb < CBLIMIT; sfb++ ) { scalefac->l[gr][ch][sfb] = scalesave_l[ch][sfb]; } for ( i = 0; i < 3; i++ ) for ( sfb = 0; sfb < SFB_SMAX; sfb++ ) { scalefac->s[gr][ch][sfb][i] = scalesave_s[ch][sfb][i]; } { real_bits[ch]=save_real_bits[ch]; memcpy(l3_enc[gr][ch],save_l3_enc[ch],sizeof(l3_enc[gr][ch])); memcpy(cod_info[ch],&save_cod_info[ch],sizeof(save_cod_info[ch])); if ( fr_ps->header->version == 1 ) status[ch] = scale_bitcount( scalefac, cod_info[ch], gr, ch ); else status[ch] = scale_bitcount_lsf( scalefac, cod_info[ch], gr, ch ); if (status[ch]) { fprintf(stderr,"Error recomputing scalefac_compress...this should not happen"); exit(-10); } } cod_info[ch]->part2_length = part2_length( scalefac, fr_ps, gr, ch, l3_side ); cod_info[ch]->part2_3_length = cod_info[ch]->part2_length + real_bits[ch]; #ifdef HAVEGTK if (gtkflag) pinfo->LAMEmainbits[gr][ch]=cod_info[ch]->part2_3_length; #endif } } /* finish up */ for (ch=0 ; ch < stereo ; ch ++ ) { if (!VBR) ResvAdjust( fr_ps, cod_info[ch], l3_side, mean_bits ); cod_info[ch]->global_gain = nint( cod_info[ch]->quantizerStepSize + 210.0 ); assert( cod_info[ch]->global_gain < 256 ); } } /*************************************************************************** * inner_loop * *************************************************************************** * The code selects the best quantizerStepSize for a particular set * of scalefacs */ int inner_loop( double xr[2][2][576], double xrpow[2][2][576], int l3_enc[2][2][576], int max_bits, gr_info *cod_info, int gr, int ch ) { int bits; assert( max_bits >= 0 ); cod_info->quantizerStepSize -= 1.0;; do { cod_info->quantizerStepSize += 1.0; quantize_xrpow( xrpow[gr][ch], l3_enc[gr][ch], cod_info ); bits = count_bits(l3_enc[gr][ch],cod_info); } while ( bits > max_bits ); return bits; } /***************************************************************************/ /* part2_length */ /***************************************************************************/ /* calculates the number of bits needed to encode the scalefacs in the */ /* main data block */ int part2_length( III_scalefac_t *scalefac, frame_params *fr_ps, int gr, int ch, III_side_info_t *si ) { int slen1, slen2, bits, partition; gr_info *gi = &si->gr[gr].ch[ch].tt; bits = 0; if ( fr_ps->header->version == 1 ) { static int slen1_tab[16] = { 0, 0, 0, 0, 3, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4 }; static int slen2_tab[16] = { 0, 1, 2, 3, 0, 1, 2, 3, 1, 2, 3, 1, 2, 3, 2, 3 }; slen1 = slen1_tab[ gi->scalefac_compress ]; slen2 = slen2_tab[ gi->scalefac_compress ]; if ( (gi->window_switching_flag == 1) && (gi->block_type == 2) ) { bits += (18 * slen1) + (18 * slen2); } else { if ( (gr == 0) || (si->scfsi[ch][0] == 0) ) bits += (6 * slen1); if ( (gr == 0) || (si->scfsi[ch][1] == 0) ) /* bits += (6 * slen1); This is wrong SS 19-12-96 */ bits += (5 * slen1); if ( (gr == 0) || (si->scfsi[ch][2] == 0) ) /* bits += (6 * slen1); This is wrong SS 19-12-96 */ bits += (5 * slen2); if ( (gr == 0) || (si->scfsi[ch][3] == 0) ) /* bits += (6 * slen1); This is wrong SS 19-12-96 */ bits += (5 * slen2); } } else { /* MPEG 2 */ assert( gi->sfb_partition_table ); for ( partition = 0; partition < 4; partition++ ) bits += gi->slen[partition] * gi->sfb_partition_table[partition]; } return bits; } /*************************************************************************/ /* scale_bitcount */ /*************************************************************************/ /* Also calculates the number of bits necessary to code the scalefactors. */ int scale_bitcount( III_scalefac_t *scalefac, gr_info *cod_info, int gr, int ch ) { int i, k, sfb, max_slen1 = 0, max_slen2 = 0, /*a, b, */ ep = 2; static int slen1[16] = { 0, 0, 0, 0, 3, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4 }; static int slen2[16] = { 0, 1, 2, 3, 0, 1, 2, 3, 1, 2, 3, 1, 2, 3, 2, 3 }; if ( cod_info->window_switching_flag != 0 && cod_info->block_type == 2 ) { /* a = 18; b = 18; */ for ( i = 0; i < 3; i++ ) { for ( sfb = 0; sfb < 6; sfb++ ) if ( scalefac->s[gr][ch][sfb][i] > max_slen1 ) max_slen1 = scalefac->s[gr][ch][sfb][i]; for (sfb = 6; sfb < 12; sfb++ ) if ( scalefac->s[gr][ch][sfb][i] > max_slen2 ) max_slen2 = scalefac->s[gr][ch][sfb][i]; } } else { /* block_type == 1,2,or 3 */ /* a = 11; b = 10; */ for ( sfb = 0; sfb < 11; sfb++ ) if ( scalefac->l[gr][ch][sfb] > max_slen1 ) max_slen1 = scalefac->l[gr][ch][sfb]; for ( sfb = 11; sfb < 21; sfb++ ) if ( scalefac->l[gr][ch][sfb] > max_slen2 ) max_slen2 = scalefac->l[gr][ch][sfb]; } for ( k = 0; k < 16; k++ ) { if ( (max_slen1 < (1<scalefac_compress = k; return ep; } /*************************************************************************/ /* scale_bitcount_lsf */ /*************************************************************************/ /* Also counts the number of bits to encode the scalefacs but for MPEG 2 */ /* Lower sampling frequencies (24, 22.05 and 16 kHz.) */ /* This is reverse-engineered from section 2.4.3.2 of the MPEG2 IS, */ /* "Audio Decoding Layer III" */ int scale_bitcount_lsf( III_scalefac_t *scalefac, gr_info *cod_info, int gr, int ch ) { int table_number, row_in_table, partition, nr_sfb, window, over; int i, sfb, max_sfac[ 4 ]; unsigned *partition_table; /* Set partition table. Note that should try to use table one, but do not yet... */ if ( cod_info->preflag ) table_number = 2; else table_number = 0; for ( i = 0; i < 4; i++ ) max_sfac[i] = 0; if ( cod_info->window_switching_flag != 0 && cod_info->block_type == 2 ) { row_in_table = 1; partition_table = &nr_of_sfb_block[table_number][row_in_table][0]; for ( sfb = 0, partition = 0; partition < 4; partition++ ) { nr_sfb = partition_table[ partition ] / 3; for ( i = 0; i < nr_sfb; i++, sfb++ ) for ( window = 0; window < 3; window++ ) if ( scalefac->s[gr][ch][sfb][window] > max_sfac[partition] ) max_sfac[partition] = scalefac->s[gr][ch][sfb][window]; } } else { row_in_table = 0; partition_table = &nr_of_sfb_block[table_number][row_in_table][0]; partition = 0; for ( sfb = 0, partition = 0; partition < 4; partition++ ) { nr_sfb = partition_table[ partition ]; for ( i = 0; i < nr_sfb; i++, sfb++ ) if ( scalefac->l[gr][ch][sfb] > max_sfac[partition] ) max_sfac[partition] = scalefac->l[gr][ch][sfb]; } } for ( over = 0, partition = 0; partition < 4; partition++ ) { if ( max_sfac[partition] > max_sfac_tab[table_number][partition] ) over++; } if ( !over ) { /* Since no bands have been over-amplified, we can set scalefac_compress and slen[] for the formatter */ static int log2tab[] = { 0, 1, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4 }; unsigned slen1, slen2, slen3, slen4; cod_info->sfb_partition_table = &nr_of_sfb_block[table_number][row_in_table][0]; for ( partition = 0; partition < 4; partition++ ) cod_info->slen[partition] = log2tab[max_sfac[partition]]; /* set scalefac_compress */ slen1 = cod_info->slen[ 0 ]; slen2 = cod_info->slen[ 1 ]; slen3 = cod_info->slen[ 2 ]; slen4 = cod_info->slen[ 3 ]; switch ( table_number ) { case 0: cod_info->scalefac_compress = (((slen1 * 5) + slen2) << 4) + (slen3 << 2) + slen4; break; case 1: cod_info->scalefac_compress = 400 + (((slen1 * 5) + slen2) << 2) + slen3; break; case 2: cod_info->scalefac_compress = 500 + (slen1 * 3) + slen2; break; default: fprintf( stderr, "intensity stereo not implemented yet\n" ); exit( EXIT_FAILURE ); break; } } #ifdef DEBUG if ( over ) printf( "---WARNING !! Amplification of some bands over limits\n" ); #endif return over; } /*************************************************************************/ /* calc_noise */ /*************************************************************************/ /* mt 5/99: Function: Improved calc_noise for a single channel */ int calc_noise1( double xr[576], int ix[576], gr_info *cod_info, double xfsf[4][CBLIMIT], int distort[4][CBLIMIT], III_psy_xmin *l3_xmin,int gr, int ch, double noise) { int start, end, sfb, l, i, over=0; double sum,step,bw; D192_3 *xr_s; I192_3 *ix_s; #define PRECALC_SIZE 1024 /* WAS 256 !!! */ static double pow43[PRECALC_SIZE]; static int init=0; noise=0; if (init==0) { init++; for(i=0;iquantizerStepSize) * 0.25 ); for ( sfb = 0; sfb < cod_info->sfb_lmax; sfb++ ) { start = scalefac_band_long[ sfb ]; end = scalefac_band_long[ sfb+1 ]; bw = end - start; for ( sum = 0.0, l = start; l < end; l++ ) { double temp; if (ix[l] %d\n",ix[l]); .. table is not big enough */ } sum += temp * temp; } xfsf[0][sfb] = sum / bw; distort[0][sfb] = ( xfsf[0][sfb] > l3_xmin->l[gr][ch][sfb] ); if (distort[0][sfb]) over++; if (distort[0][sfb]) noise += xfsf[0][sfb] - l3_xmin->l[gr][ch][sfb]; } for ( i = 0; i < 3; i++ ) { step = pow( 2.0, (cod_info->quantizerStepSize) * 0.25 ); /* subblock_gain ? */ if (cod_info->subblock_gain[i] ) step *= pow(2.0,-2.0*cod_info->subblock_gain[i]); for ( sfb = cod_info->sfb_smax; sfb < 12; sfb++ ) { start = scalefac_band_short[ sfb ]; end = scalefac_band_short[ sfb+1 ]; bw = end - start; for ( sum = 0.0, l = start; l < end; l++ ) { double temp; if((*ix_s)[l][i] %d\n",(*ix_s)[l][i]); */ } sum += temp * temp; } xfsf[i+1][sfb] = sum / bw; distort[i+1][sfb] = ( xfsf[i+1][sfb] > l3_xmin->s[gr][ch][sfb][i] ); if (distort[i+1][sfb]) over++; if (distort[i+1][sfb]) noise += xfsf[i+1][sfb]-l3_xmin->s[gr][ch][sfb][i]; } } return over; } /*************************************************************************/ /* calc_noise2 */ /*************************************************************************/ /* Improved version of calc_noise for dual channel. This routine is */ /* used when you are quantizaing mid and side channels using masking */ /* thresholds from L and R channels. mt 5/99 */ void calc_noise2( double xr[2][576], int ix[2][576], gr_info *cod_info[2], double xfsf[2][4][CBLIMIT], int distort[2][4][CBLIMIT], III_psy_xmin *l3_xmin,int gr,int stereo, int over[2], double noise[2]) { int start, end, sfb, l, i; double sum[2],step_s[3][2],step[2],bw; D192_3 *xr_s[2]; I192_3 *ix_s[2]; #define PRECALC_SIZE 1024 /* WAS 256 !!! */ static double pow43[PRECALC_SIZE]; int ch; static int init=0; double diff[2], qcoeff; if (init==0) { init++; for(i=0;iblock_type != 2) { for (ch=0 ; ch < stereo ; ch ++ ) { over[ch]=0; noise[ch]=0; step[ch] = pow( 2.0, (cod_info[ch]->quantizerStepSize) * 0.25 ); } for ( sfb = 0; sfb < SFB_LMAX-1; sfb++ ) { start = scalefac_band_long[ sfb ]; end = scalefac_band_long[ sfb+1 ]; bw = end - start; for (ch=0 ; ch < stereo ; ch ++ ) sum[ch]=0; for ( l = start; l < end; l++ ) { for (ch=0 ; ch < stereo ; ch ++ ) { qcoeff= (ix[ch][l] < PRECALC_SIZE) ? pow43[ix[ch][l]] * step[ch]: pow((double)ix[ch][l],4.0/3.0)*step[ch]; if (xr[ch][l]<0) qcoeff=-qcoeff; diff[ch]=xr[ch][l]-qcoeff; } sum[0] += (diff[0]+diff[1])*(diff[0]+diff[1])/(2.0); sum[1] += (diff[0]-diff[1])*(diff[0]-diff[1])/(2.0); } for (ch=0 ; ch < stereo ; ch ++ ) { xfsf[ch][0][sfb] = sum[ch] / bw; distort[ch][0][sfb] = ( xfsf[ch][0][sfb] > l3_xmin->l[gr][ch][sfb] ); if (distort[ch][0][sfb]) over[ch]++; if (distort[ch][0][sfb]) noise[ch] += xfsf[ch][0][sfb] - l3_xmin->l[gr][ch][sfb]; } /* if there is audible distortion in left or right channel, set flags * to denote distortion in both mid and side channels */ for (ch=0 ; ch < stereo ; ch ++ ) { distort[ch][0][sfb] = (distort[0][0][sfb] || distort[1][0][sfb]); } } } /* calc_noise2: we can assume block types of both channels must be the same */ if (cod_info[0]->block_type == 2) { for (ch=0 ; ch < stereo ; ch ++ ) { for (i=0;i<3;i++){ step_s[i][ch] = pow( 2.0, (cod_info[ch]->quantizerStepSize) * 0.25 ); /* subblock_gain ? */ if (cod_info[ch]->subblock_gain[i] ) step_s[i][ch] *= pow(2.0,-2.0*cod_info[ch]->subblock_gain[i]); } over[ch] = 0; xr_s[ch] = (D192_3 *) xr[ch]; ix_s[ch] = (I192_3 *) ix[ch]; } for ( sfb = 0 ; sfb < SFB_SMAX-1; sfb++ ) { start = scalefac_band_short[ sfb ]; end = scalefac_band_short[ sfb+1 ]; bw = end - start; for ( i = 0; i < 3; i++ ) { for (ch=0 ; ch < stereo ; ch ++ ) sum[ch] = 0.0; for ( l = start; l < end; l++ ) { for (ch=0 ; ch < stereo ; ch ++ ) { qcoeff = ((*ix_s[ch])[l][i] l3_xmin->s[gr][ch][sfb][i] ); if (distort[ch][i+1][sfb]) over[ch]++; if (distort[ch][i+1][sfb]) noise[ch] += xfsf[ch][i+1][sfb]-l3_xmin->s[gr][ch][sfb][i]; } /* if there is audible distortion in left or right channel, set flags * to denote distortion in both mid and side channels */ for (ch=0 ; ch < stereo ; ch ++ ) distort[ch][i+1][sfb] = (distort[0][i+1][sfb] || distort[1][i+1][sfb] ); } } } } /*************************************************************************/ /* calc_xmin */ /*************************************************************************/ /* Calculate the allowed distortion for each scalefactor band, as determined by the psychoacoustic model. xmin(sb) = ratio(sb) * en(sb) / bw(sb) */ void calc_xmin( double xr[2][2][576], III_psy_ratio *ratio, gr_info *cod_info, III_psy_xmin *l3_xmin, int gr, int ch ) { int start, end, sfb, l, b; double en0, bw; D192_3 *xr_s; xr_s = (D192_3 *) xr[gr][ch] ; for ( sfb = cod_info->sfb_smax; sfb < SFB_SMAX - 1; sfb++ ) { start = scalefac_band_short[ sfb ]; end = scalefac_band_short[ sfb + 1 ]; bw = end - start; for ( b = 0; b < 3; b++ ) { for ( en0 = 0.0, l = start; l < end; l++ ) en0 += (*xr_s)[l][b] * (*xr_s)[l][b]; l3_xmin->s[gr][ch][sfb][b] = ratio->s[gr][ch][sfb][b] * en0 / bw; } } for ( sfb = 0; sfb < cod_info->sfb_lmax; sfb++ ) { start = scalefac_band_long[ sfb ]; end = scalefac_band_long[ sfb+1 ]; bw = end - start; for ( en0 = 0.0, l = start; l < end; l++ ) en0 += xr[gr][ch][l] * xr[gr][ch][l]; l3_xmin->l[gr][ch][sfb] = ratio->l[gr][ch][sfb] * en0 / bw; } } /*************************************************************************/ /* loop_break */ /*************************************************************************/ /* Function: Returns zero if there is a scalefac which has not been amplified. Otherwise it returns one. */ int loop_break( III_scalefac_t *scalefac, gr_info *cod_info, int gr, int ch ) { int i, sfb; for ( sfb = 0; sfb < cod_info->sfb_lmax; sfb++ ) if ( scalefac->l[gr][ch][sfb] == 0 ) return 0; for ( sfb = cod_info->sfb_smax; sfb < 12; sfb++ ) for ( i = 0; i < 3; i++ ) if ( scalefac->s[gr][ch][sfb][i] == 0 ) return 0; return 1; } /*************************************************************************/ /* preemphasis */ /*************************************************************************/ /* See ISO 11172-3 section C.1.5.4.3.4 */ int preemphasis2( double xr[576], double xrpow[576], III_psy_xmin *l3_xmin, int gr, int ch, III_side_info_t *l3_side, int distort[4][CBLIMIT] ) { int i, sfb, start, end, scfsi_band, over; double ifqstep; gr_info *cod_info = &l3_side->gr[gr].ch[ch].tt; if ( gr == 1 ) { /* If the second granule is being coded and scfsi is active in at least one scfsi_band, the preemphasis in the second granule is set equal to the setting in the first granule */ for ( scfsi_band = 0; scfsi_band < 4; scfsi_band++ ) if ( l3_side->scfsi[ch][scfsi_band] ) { cod_info->preflag = l3_side->gr[0].ch[ch].tt.preflag; return 0; } } /* Preemphasis is switched on if in all the upper four scalefactor bands the actual distortion exceeds the threshold after the first call of the inner loop */ over = 0; if ( cod_info->block_type != 2 && cod_info->preflag == 0 ) { for ( sfb = 17; sfb < 21; sfb++ ) if ( distort[0][sfb] ) over++; if (over == 4 ) { double t,t34; cod_info->preflag = 1; ifqstep = ( cod_info->scalefac_scale == 0 ) ? SQRT2 : 2.0; for ( sfb = 11; sfb < cod_info->sfb_lmax; sfb++ ) { t=pow( ifqstep, (double) pretab[sfb] ); t34=sqrt(sqrt(t)*t); l3_xmin->l[gr][ch][sfb] *= t*t; start = scalefac_band_long[ sfb ]; end = scalefac_band_long[ sfb+1 ]; for( i = start; i < end; i++ ) xr[i]*=t; for( i = start; i < end; i++ ) xrpow[i]*=t34; } } } return (over == 4); } /*************************************************************************/ /* amp_scalefac_bands */ /*************************************************************************/ /* Amplify the scalefactor bands that violate the masking threshold. See ISO 11172-3 Section C.1.5.4.3.5 */ int amp_scalefac_bands2( double xr[576], double xrpow[576], III_psy_xmin *l3_xmin, III_side_info_t *l3_side, III_scalefac_t *scalefac, int gr, int ch, int iteration, int distort[4][CBLIMIT]) { int start, end, l, sfb, i, scfsi_band, over = 0; double ifqstep, ifqstep2, ifqstep34; D192_3 *xr_s; D192_3 *xrpow_s; gr_info *cod_info, *gr0; int copySF, preventSF; cod_info = &l3_side->gr[gr].ch[ch].tt; gr0 = &l3_side->gr[0].ch[ch].tt; xr_s = (D192_3 *) xr; xrpow_s = (D192_3 *) xrpow; copySF = 0; preventSF = 0; if ( cod_info->scalefac_scale == 0 ) ifqstep = SQRT2; else ifqstep = 2.0; if ( gr == 1 ) { /* If the second granule is being coded and scfsi is active in at least one scfsi_band... */ for ( scfsi_band = 0; scfsi_band < 4; scfsi_band++ ) if ( l3_side->scfsi[ch][scfsi_band] ) { /* a) ifqstep has to be set similar to the first granule... */ if ( gr0->scalefac_scale == 0 ) ifqstep = SQRT2; else ifqstep = 2.0; if ( iteration == 1 ) { /* b) If it is the first iteration, the scalefactors of scalefactor bands in which scfsi is enabled must be taken from the first granule */ copySF = 1; } else { /* c) If it is not the first iteration, the amplification must be prevented for scalefactor bands in which scfsi is enabled */ preventSF = 1; } break; } } ifqstep2 = ifqstep * ifqstep; scfsi_band = 0; ifqstep34=sqrt(sqrt(ifqstep)*ifqstep); for ( sfb = 0; sfb < cod_info->sfb_lmax; sfb++ ) { if ( copySF || preventSF ) { if ( sfb == scfsi_band_long[scfsi_band + 1] ) scfsi_band += 1; if ( l3_side->scfsi[ch][scfsi_band] ) { if ( copySF ) scalefac->l[gr][ch][sfb] = scalefac->l[0][ch][sfb]; continue; } } if ( distort[0][sfb] ) { over++; l3_xmin->l[gr][ch][sfb] *= ifqstep2; scalefac->l[gr][ch][sfb]++; start = scalefac_band_long[sfb]; end = scalefac_band_long[sfb+1]; for ( l = start; l < end; l++ ) xr[l] *= ifqstep; for ( l = start; l < end; l++ ) xrpow[l] *= ifqstep34; } } /* Note that scfsi is not enabled for frames containing short blocks */ for ( i = 0; i < 3; i++ ) for ( sfb = cod_info->sfb_smax; sfb < 12; sfb++ ) { if ( distort[i+1][sfb]) { over++; l3_xmin->s[gr][ch][sfb][i] *= ifqstep2; scalefac->s[gr][ch][sfb][i]++; #ifdef DEBUGSC printf( "cod_info->scalefac[%d][%d] = %d (amp_scale)\n", i,sfb,scalefac->s[gr][ch][sfb][i] ); #endif start = scalefac_band_short[sfb]; end = scalefac_band_short[sfb+1]; for ( l = start; l < end; l++ ) (*xr_s)[l][i] *= ifqstep; for ( l = start; l < end; l++ ) (*xrpow_s)[l][i] *= ifqstep34; } } return over; } void quantize_xrpow( double xr[576], int ix[576], gr_info *cod_info ) { /* quantize on xr^(3/4) instead of xr */ register int j; double quantizerStepSize; double istep_l,istep0,istep1,istep2; quantizerStepSize = cod_info->quantizerStepSize; istep_l = pow ( 2.0, quantizerStepSize * -0.1875 ); if ((cod_info->block_type==2)) { istep0 = istep_l * pow(2.0,1.5* (double) cod_info->subblock_gain[0]); istep1 = istep_l * pow(2.0,1.5* (double) cod_info->subblock_gain[1]); istep2 = istep_l * pow(2.0,1.5* (double) cod_info->subblock_gain[2]); } else { istep0 = istep_l; istep1 = istep_l; istep2 = istep_l; } for (j=192;j>0;j--) { #if defined(__GNUC__) && defined(__i386__) asm ("fistpl %0 ": "=m"(*(ix++)): "t"(istep0*(*(xr++)) - 0.0946): "st"); asm ("fistpl %0 ": "=m"(*(ix++)): "t"(istep1*(*(xr++)) - 0.0946): "st"); asm ("fistpl %0 ": "=m"(*(ix++)): "t"(istep2*(*(xr++)) - 0.0946): "st"); #else *(ix++) = (int)( istep0*(*(xr++)) + 0.4054); *(ix++) = (int)( istep1*(*(xr++)) + 0.4054); *(ix++) = (int)( istep2*(*(xr++)) + 0.4054); #endif } } /*************************************************************************/ /* ix_max */ /*************************************************************************/ /* Function: Calculate the maximum of ix from 0 to 575 */ int ix_max( int ix[576], unsigned int begin, unsigned int end ) { int i, max = 0; for ( i = begin; i < end; i++ ) { int x = ix[i]; if ( x > max ) max = x; } return max; } /*************************************************************************/ /* xr_max */ /*************************************************************************/ /* Function: Calculate the maximum of xr[576] from 0 to 575 */ double xr_max( double xr[576], unsigned int begin, unsigned int end ) { int i; double max = 0.0, temp; for ( i = begin; i < end; i++ ) if( (temp = fabs(xr[i])) > max ) max = temp; return max; } /* Noiseless coding -- Huffman coding */ /*************************************************************************/ /* calc_runlen */ /*************************************************************************/ /* Function: Calculation of count1, big_values (Partitions ix into big values, quadruples and zeros). */ void calc_runlen( int ix[576], gr_info *cod_info ) { int i; if ( cod_info->window_switching_flag && (cod_info->block_type == 2) ) { /* short blocks */ cod_info->count1 = 0; cod_info->big_values = 288; } else { for ( i = 576; i > 1; i -= 2 ) if ( ix[i-1] | ix[i-2] ) break; cod_info->count1 = 0 ; for ( ; i > 3; i -= 4 ) if ( (ix[i-1] | ix[i-2] | ix[i-3] | ix[i-4]) <= 1 ) cod_info->count1++; else break; cod_info->big_values = i/2; } } /*************************************************************************/ /* count1_bitcount */ /*************************************************************************/ /* Determines the number of bits to encode the quadruples. */ int count1_bitcount( int ix[ 576 ], gr_info *cod_info ) { int abs_and_sign( int *x ); int p, i, k, bitsum_count1; int v, w, x, y, signbits; int sum0 = 0, sum1 = 0; for ( i = cod_info->big_values * 2, k = 0; k < cod_info->count1; i += 4, k++ ) { v = ix[i]; w = ix[i+1]; x = ix[i+2]; y = ix[i+3]; /* v = ix[ i ]; w = ix[ i + 1 ]; x = ix[ i + 2 ]; y = ix[ i + 3 ]; abs_and_sign( &v ); abs_and_sign( &w ); abs_and_sign( &x ); abs_and_sign( &y ); */ p = v + (w << 1) + (x << 2) + (y << 3); signbits = 0; if ( v != 0 ) signbits ++; if ( w != 0 ) signbits ++; if ( x != 0 ) signbits ++; if ( y != 0 ) signbits ++; sum0 += signbits; sum1 += signbits; sum0 += ht[ 32 ].hlen[ p ]; sum1 += ht[ 33 ].hlen[ p ]; } if ( sum0 < sum1 ) { bitsum_count1 = sum0; cod_info->count1table_select = 0; } else { bitsum_count1 = sum1; cod_info->count1table_select = 1; } return( bitsum_count1 ); } struct { unsigned region0_count; unsigned region1_count; } subdv_table[ 23 ] = { {0, 0}, /* 0 bands */ {0, 0}, /* 1 bands */ {0, 0}, /* 2 bands */ {0, 0}, /* 3 bands */ {0, 0}, /* 4 bands */ {0, 1}, /* 5 bands */ {1, 1}, /* 6 bands */ {1, 1}, /* 7 bands */ {1, 2}, /* 8 bands */ {2, 2}, /* 9 bands */ {2, 3}, /* 10 bands */ {2, 3}, /* 11 bands */ {3, 4}, /* 12 bands */ {3, 4}, /* 13 bands */ {3, 4}, /* 14 bands */ {4, 5}, /* 15 bands */ {4, 5}, /* 16 bands */ {4, 6}, /* 17 bands */ {5, 6}, /* 18 bands */ {5, 6}, /* 19 bands */ {5, 7}, /* 20 bands */ {6, 7}, /* 21 bands */ {6, 7}, /* 22 bands */ }; /*************************************************************************/ /* subdivide */ /*************************************************************************/ /* presumable subdivides the bigvalue region which will use separate Huffman tables. */ void subdivide( gr_info *cod_info ) { int scfb_anz = 0; int bigvalues_region; if ( cod_info->big_values == 0 ) { /* no big_values region */ cod_info->region0_count = 0; cod_info->region1_count = 0; } else { bigvalues_region = 2 * cod_info->big_values; if ( (cod_info->window_switching_flag == 0) ) { /* long blocks */ int thiscount, index; /* Calculate scfb_anz */ while ( scalefac_band_long[scfb_anz] < bigvalues_region ) scfb_anz++; assert( scfb_anz < 23 ); cod_info->region0_count = subdv_table[scfb_anz].region0_count; thiscount = cod_info->region0_count; index = thiscount + 1; while ( thiscount && (scalefac_band_long[index] > bigvalues_region) ) { thiscount -= 1; index -= 1; } cod_info->region0_count = thiscount; cod_info->region1_count = subdv_table[scfb_anz].region1_count; index = cod_info->region0_count + cod_info->region1_count + 2; thiscount = cod_info->region1_count; while ( thiscount && (scalefac_band_long[index] > bigvalues_region) ) { thiscount -= 1; index -= 1; } cod_info->region1_count = thiscount; cod_info->address1 = scalefac_band_long[cod_info->region0_count+1]; cod_info->address2 = scalefac_band_long[cod_info->region0_count + cod_info->region1_count + 2 ]; cod_info->address3 = bigvalues_region; } else { if ( (cod_info->block_type == 2) ) { cod_info->region0_count = 8; cod_info->region1_count = 36; cod_info->address1 = 36; cod_info->address2 = bigvalues_region; cod_info->address3 = 0; } else { cod_info->region0_count = 7; cod_info->region1_count = 13; cod_info->address1 = scalefac_band_long[cod_info->region0_count+1]; cod_info->address2 = bigvalues_region; cod_info->address3 = 0; } } } } /*************************************************************************/ /* new_choose table */ /*************************************************************************/ /* Choose the Huffman table that will encode ix[begin..end] with the fewest bits. Note: This code contains knowledge about the sizes and characteristics of the Huffman tables as defined in the IS (Table B.7), and will not work with any arbitrary tables. */ int new_choose_table( int ix[576], unsigned int begin, unsigned int end, int * s ) { int i, max; int choice[ 2 ]; int sum[ 2 ]; max = ix_max( ix, begin, end ); if ( max == 0 ) return 0; choice[ 0 ] = 0; choice[ 1 ] = 0; if ( max < 15 ) { /* try tables with no linbits */ for ( i = 0; i < 14; i++ ) if ( ht[i].xlen > max ) { choice[ 0 ] = i; break; } sum[ 0 ] = count_bit( ix, begin, end, choice[0] ); switch ( choice[0] ) { case 2: sum[ 1 ] = count_bit( ix, begin, end, 3 ); if ( sum[1] <= sum[0] ) { sum[0]=sum[1]; choice[ 0 ] = 3; } break; case 5: sum[ 1 ] = count_bit( ix, begin, end, 6 ); if ( sum[1] <= sum[0] ) { sum[0]=sum[1]; choice[ 0 ] = 6; } break; case 7: sum[ 1 ] = count_bit( ix, begin, end, 8 ); if ( sum[1] <= sum[0] ) { choice[ 0 ] = 8; sum[ 0 ] = sum[ 1 ]; } sum[ 1 ] = count_bit( ix, begin, end, 9 ); if ( sum[1] <= sum[0] ) { sum[0]=sum[1]; choice[ 0 ] = 9; } break; case 10: sum[ 1 ] = count_bit( ix, begin, end, 11 ); if ( sum[1] <= sum[0] ) { choice[ 0 ] = 11; sum[ 0 ] = sum[ 1 ]; } sum[ 1 ] = count_bit( ix, begin, end, 12 ); if ( sum[1] <= sum[0] ) { sum[0]=sum[1]; choice[ 0 ] = 12; } break; case 13: sum[ 1 ] = count_bit( ix, begin, end, 15 ); if ( sum[1] <= sum[0] ) { sum[0]=sum[1]; choice[ 0 ] = 15; } break; default: break; } } else { /* try tables with linbits */ max -= 15; for ( i = 15; i < 24; i++ ) { if ( ht[i].linmax >= max ) { choice[ 0 ] = i; break; } } for ( i = 24; i < 32; i++ ) { if ( ht[i].linmax >= max ) { choice[ 1 ] = i; break; } } sum[ 0 ] = count_bit( ix, begin, end, choice[0] ); sum[ 1 ] = count_bit( ix, begin, end, choice[1] ); if ( sum[1] < sum[0] ) { sum[0]=sum[1]; choice[ 0 ] = choice[ 1 ]; } } *s=sum[0]; return choice[ 0 ]; } /*************************************************************************/ /* choose table */ /*************************************************************************/ int choose_table( int max ) { int i, choice; if ( max == 0 ) return 0; choice = 0; if ( max < 15 ) { for ( i = 0; i < 15; i++ ) { if ( ht[i].xlen > max ) { choice = i; break; } } } else { max -= 15; for (i = 15; i < 32; i++ ) { if ( ht[i].linmax >= max ) { choice = i; break; } } } assert( choice ); return choice; } /*************************************************************************/ /* bigv_bitcount */ /*************************************************************************/ /* Function: Select huffman code tables for bigvalues regions Function: Count the number of bits necessary to code the bigvalues region. */ int bigv_bitcount( int ix[576], gr_info *cod_info ) { int bits = 0; cod_info->table_select[0] = 0; cod_info->table_select[1] = 0; cod_info->table_select[2] = 0; if ( cod_info->window_switching_flag && (cod_info->block_type == 2) ) { /* Within each scalefactor band, data is given for successive time windows, beginning with window 0 and ending with window 2. Within each window, the quantized values are then arranged in order of increasing frequency... */ int sfb, window, line, start, end, max1, max2, x, y; int region1Start; int *pmax; I192_3 *ix_s; region1Start = 12; max1 = max2 = 0; for ( sfb = 0; sfb < 13; sfb++ ) { start = scalefac_band_short[ sfb ]; end = scalefac_band_short[ sfb+1 ]; if ( start < region1Start ) pmax = &max1; else pmax = &max2; for ( window = 0; window < 3; window++ ) for ( line = start; line < end; line += 2 ) { x = ix[ (line * 3) + window ]; y = ix[ ((line + 1) * 3) + window ]; *pmax = *pmax > x ? *pmax : x; *pmax = *pmax > y ? *pmax : y; } } cod_info->table_select[0] = choose_table( max1 ); cod_info->table_select[1] = choose_table( max2 ); /* Within each scalefactor band, data is given for successive time windows, beginning with window 0 and ending with window 2. Within each window, the quantized values are then arranged in order of increasing frequency... */ sfb = 0; ix_s = (I192_3 *) &ix[0]; for ( ; sfb < 13; sfb++ ) { unsigned tableindex = 100; start = scalefac_band_short[ sfb ]; end = scalefac_band_short[ sfb+1 ]; if ( start < 12 ) tableindex = cod_info->table_select[ 0 ]; else tableindex = cod_info->table_select[ 1 ]; for ( window = 0; window < 3; window++ ) for ( line = start; line < end; line += 2 ) { unsigned int code, ext; int cbits, xbits; int x = (*ix_s)[line][window]; int y = (*ix_s)[line + 1][window]; bits += HuffmanCode( tableindex, x, y, &code, &ext, &cbits, &xbits ); } } } else { unsigned int table; int s; if ( cod_info->address1 > 0 ) { table = cod_info->table_select[0] = new_choose_table( ix, 0, cod_info->address1, &s ); if( table != 0 ) /* region0 */ bits += s; } if ( cod_info->address2 > cod_info->address1 ) { table = cod_info->table_select[1] = new_choose_table( ix, cod_info->address1, cod_info->address2, &s ); if( table != 0 ) /* region1 */ bits += s; } if ( cod_info->big_values * 2 > cod_info->address2 ) { table = cod_info->table_select[2] = new_choose_table( ix, cod_info->address2, cod_info->big_values * 2, &s ); if( table != 0 ) /* region2 */ bits += s; } } return bits; } /*************************************************************************/ /* count_bit */ /*************************************************************************/ /* Function: Count the number of bits necessary to code the subregion. */ int count_bit( int ix[576], unsigned int start, unsigned int end, unsigned int table ) { unsigned linbits; register int i, sum; register int x,y; struct huffcodetab *h; unsigned char *hlen; if(table==0) return 0; h = &(ht[table]); sum = 0; linbits = h->linbits; hlen=h->hlen; if(table>15) { /* ESC-table is used */ for(i=start;i14) { x = 15; sum += linbits; } if(y>14) { y = 15; sum += linbits; } sum += hlen[(x*16)+y]; if(x!=0) sum++; if(y!=0) sum++; } } else { /* No ESC-words */ for(i=start;iwindow_switching_flag != 0 && cod_info->block_type == 2 ) { cod_info->sfb_lmax = 0; /* No sb*/ cod_info->sfb_smax = 0; } else { /* MPEG 1 doesnt use last scalefactor band? */ cod_info->sfb_lmax = SFB_LMAX - 1; cod_info->sfb_smax = SFB_SMAX - 1; /* No sb */ } } int count_bits(ix,cod_info) int *ix; /* I576 *ix; */ gr_info *cod_info; { int bits,i; for ( i = 0; i < 576; i++ ) { if ( ix[i] > 8191 + 14) return 100000; /* report unsuitable quantizer */ } calc_runlen(ix,cod_info); /*count1,big_values*/ bits = count1_bitcount(ix, cod_info); /*count1_table selection*/ subdivide(cod_info); /* bigvalues sfb division */ bits += bigv_bitcount(ix,cod_info); /* bit count */ /* printf("\nglobal_gain = %f bits= %d ",cod_info->quantizerStepSize,bits);*/ return bits; } static int OldValue = -30; /* guess it or so. */ typedef enum { BINSEARCH_NONE, BINSEARCH_UP, BINSEARCH_DOWN } binsearchDirection_t; /*-------------------------------------------------------------------------*/ int bin_search_StepSize2(int desired_rate, double start, int bot, int *ix, double xrs[576], double xrspow[576], gr_info *cod_info) /*-------------------------------------------------------------------------*/ { int flag_GoneOver = 0; int CurrentStep = 4; int nBits, Counter = 0; int StepSize = OldValue; binsearchDirection_t Direction = BINSEARCH_NONE; do { cod_info->quantizerStepSize = StepSize; quantize_xrpow(xrspow, ix, cod_info); nBits = count_bits(ix,cod_info); if (CurrentStep == 1 || Counter > 50) /* the 50 is for ''lock-up'' - havent found any.*/ { break; /* nothing to adjust anymore */ } if (flag_GoneOver) { CurrentStep /= 2; } if (nBits > desired_rate) /* increase Quantize_StepSize */ { if (Direction == BINSEARCH_DOWN && !flag_GoneOver) { flag_GoneOver = 1; CurrentStep /= 2; /* late adjust */ } Direction = BINSEARCH_UP; StepSize += CurrentStep; } else if (nBits < desired_rate) { if (Direction == BINSEARCH_UP && !flag_GoneOver) { flag_GoneOver = 1; CurrentStep /= 2; /* late adjust */ } Direction = BINSEARCH_DOWN; StepSize -= CurrentStep; } else break; /* nBits == desired_rate;; most unlikely to happen. */ } while (1); /* For-ever, break is adjusted. */ OldValue = StepSize; return nBits; } int bin_search_StepSize(int desired_rate, double start, int bot, int *ix, double xrs[576], double xrspow[576], gr_info * cod_info) { int top,next,last; int bit; top = start; next = start; do { last = next; next = (top+bot)/2.0; cod_info->quantizerStepSize = next; quantize_xrpow(xrspow,ix,cod_info); bit = count_bits(ix,cod_info); if (bit>desired_rate) top = next; else bot = next; // printf("\n%f %f %f %f %d %d",start,next, top,bot,bit,desired_rate); } while ((bit != desired_rate) && abs(last - next) > 1); //printf("\n done %f %d %d",next,bit,desired_rate); return bit; }