loop.c

MPEG 2的音频编码软件。喜欢多媒体的开发人员可以看看。

/**********************************************************************
 * 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
 **********************************************************************/
/**********************************************************************
 *   date   programmers                comment                        *
 *25. 6.92  Toshiyuki Ishino           Ver 1.0                        *
 *29.10.92  Masahiro Iwadare           Ver 2.0                        *
 *17. 4.93  Masahiro Iwadare           Updated for IS Modification    *
 *04.11.93  Seymour Shlien             Speed up inner loop            *
 *09.09.95  [e-mail]            Several changes for updated IS,*
 *                                     and some MPEG2-LSF support     *
 *20.12.96  [e-mail]  Fixed some bugs and improved   *
 *                                     the appearance                 *
 *********************************************************************/ 

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <assert.h>
#include "l3side.h"
#include "loop.h"
#include "huffman.h"
#include "l3bitstream.h"
#include "reservoir.h"
#include "loop-pvt.h"

/* #define DEBUG */
/* #define DEBUGSC */
/* #define PERFORM 3 */
/* If PERFORM is defined to some number, then a file encode.log
   is preduced showing the intermediate results produced by the
   outer_loop code for the frame number = PERFORM.
*/
#define BIN_SEARCH 
/* for speeding up the iteration_loop algorithm */

#ifdef PERFORM
FILE *log_output;
extern int frameNum;
float worst_xfsf_to_xmin_ratio(III_psy_xmin *l3_xmin, double xfsf[4][CBLIMIT]
                             ,int block_type,int gr,int ch);
#endif

/* New SS 20-12-96 */
#ifdef BIN_SEARCH 
int bin_search_StepSize(int desired_rate, double start, int ix[576],
           double xrs[576], gr_info * cod_info);
int count_bits();
float worst_xfsf_to_xmin_ratio();
#endif


/*
  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,45,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];



/************************************************************************/
/*  iteration_loop()                                                    */
/************************************************************************/
void
iteration_loop( double pe[][2], double xr_org[2][2][576], III_psy_ratio *ratio,
		III_side_info_t *l3_side, int l3_enc[2][2][576],
		int mean_bits, int stereo, double xr_dec[2][2][576],
		III_scalefac_t *scalefac, frame_params *fr_ps,
		int ancillary_pad, int bitsPerFrame ) 
{
    static int firstcall = 1;
    III_psy_xmin l3_xmin;
    gr_info *cod_info;
    layer *info;
    int *main_data_begin;

    int max_bits;
    int ch, gr, sfb, i, mode_gr;
    static int huffman_read_flag = 0; 

    double xr[2][2][576];
    I576  *ix;

    main_data_begin = &l3_side->main_data_begin;
    l3_side->resvDrain = 0;

    if ( firstcall )
    {
	*main_data_begin = 0;
	firstcall = 0;
#ifdef PERFORM
	log_output = fopen("encode.log","w");
#endif
    }

    info = fr_ps->header;
    mode_gr = (info->version == 1) ? 2 : 1;

    scalefac_band_long  = &sfBandIndex[info->sampling_frequency + (info->version * 3)].l[0];
    scalefac_band_short = &sfBandIndex[info->sampling_frequency + (info->version * 3)].s[0];

    /* reading huffman code table */
    if (huffman_read_flag == 0) {
        FILE *fi = OpenTableFile( "huffcode" );
        if ( fi == NULL )
            exit( EXIT_FAILURE );
        read_huffcodetab( fi );
        huffman_read_flag++;
        fclose( fi );
    }


    for ( gr = 0; gr < mode_gr; gr++ )
    {
        for ( ch = 0; ch < stereo; ch++ )
	{
            for ( i = 0; i < 576; i++ ) 
                xr[gr][ch][i] = xr_org[gr][ch][i];
	}
    }

    ResvFrameBegin( fr_ps, l3_side, mean_bits, bitsPerFrame );

    for ( gr = 0; gr < mode_gr; gr++ )
    {
        for ( ch = 0; ch < stereo; ch++ )
        {
            ix = (I576 *) l3_enc[gr][ch];
            cod_info = (gr_info *) &(l3_side->gr[gr].ch[ch]);
            gr_deco(cod_info);
            calc_xmin( xr, ratio, cod_info, &l3_xmin, gr, ch );
	    
#ifdef DEBUG
            printf( "----- start gr[%d] ch[%1d] : block_type=%1d, window_switching_flag=%1d (loop)\n",
                    gr,ch, cod_info->block_type, cod_info->window_switching_flag ); 
#endif
	    
            if ( info->version == 1 )
                calc_scfsi( xr[gr][ch], l3_side, &l3_xmin, ch, gr );
	    
            /* calculation of number of available bit( per granule ) */
	    max_bits = ResvMaxBits( fr_ps, l3_side, &pe[gr][ch], mean_bits );
#ifdef DEBUG
            printf( " max_bits = %d, mean_bits = %d (iteration_loop)\n", max_bits, mean_bits ); 
#endif
	    
            /* reset of iteration variables */
	    
            for ( sfb = 0; sfb < 21; sfb++ )
                scalefac->l[gr][ch][sfb] = 0;
            for ( sfb = 0; sfb < 13; sfb++ )
                for ( i = 0; i < 3; i++ )
                    scalefac->s[gr][ch][sfb][i] = 0;

	    for ( i = 0; i < 4; i++ )
		cod_info->slen[i] = 0;
	    cod_info->sfb_partition_table = &nr_of_sfb_block[0][0][0];

            cod_info->part2_3_length    = 0;
            cod_info->big_values        = 0;
            cod_info->count1            = 0;
            cod_info->scalefac_compress = 0;
            cod_info->table_select[0]   = 0;
            cod_info->table_select[1]   = 0;
            cod_info->table_select[2]   = 0;
            cod_info->subblock_gain[0]  = 0;
            cod_info->subblock_gain[1]  = 0;
            cod_info->subblock_gain[2]  = 0;
            cod_info->region0_count     = 0;
            cod_info->region1_count     = 0;
            cod_info->part2_length      = 0;
            cod_info->preflag           = 0;
            cod_info->scalefac_scale    = 0;
            cod_info->quantizerStepSize = 0.0;
            cod_info->count1table_select= 0;
            
            /* all spectral values zero ? */
            if ( fabs(xr_max(xr[gr][ch], 0, 576)) != 0.0 )
            {
                cod_info->quantizerStepSize =
			 (double) quantanf_init( xr[gr][ch] );
                cod_info->part2_3_length = outer_loop( xr, max_bits, &l3_xmin,
                                                       l3_enc, fr_ps, scalefac,
                                                       gr, ch, l3_side );
            }
	    ResvAdjust( fr_ps, cod_info, l3_side, mean_bits );

	    cod_info->global_gain = nint( cod_info->quantizerStepSize + 210.0 );
	    assert( cod_info->global_gain < 256 );
        } /* for ch */
    } /* for gr */
    ResvFrameEnd( fr_ps, l3_side, mean_bits );
}



/************************************************************************/
/*  quantanf_init                                                       */
/************************************************************************/
int quantanf_init( double xr[576] )
/* Function: Calculate the first quantization step quantanf.       */
{
    int i, tp = 0;
    double system_const, minlimit;
    double sfm = 0.0, sum1 = 0.0, sum2 = 0.0;
    
    system_const = 8.0;
    minlimit = -100.0;

    for ( i = 0; i < 576; i++ )
    {
        if ( xr[i] != 0 )
	{
            double tpd = xr[i] * xr[i];
            sum1 += log( tpd );
            sum2 += tpd;
        }
    }
    if ( sum2 != 0.0 )
    {
        sfm = exp( sum1 / 576.0 ) / (sum2 / 576.0);
        tp = nint( system_const * log(sfm) );
	if ( tp < minlimit )
	    tp = minlimit;
#ifdef DEBUG
        printf(" quantanf = %d (quantanf_init)\n",tp );
#endif
    }
      return(tp-70.0); /* SS 19-12-96. Starting value of
                          global_gain or quantizerStepSize 
                          has to be reduced for iteration_loop
                       */
}





/************************************************************************/
/*  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             */
/************************************************************************/
int outer_loop(
    double xr[2][2][576],     /*  magnitudes of the spectral values */
    int max_bits,
    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 ch, III_side_info_t *l3_side )
{
    int status ;
    int scalesave_l[CBLIMIT], scalesave_s[CBLIMIT][3];
    int sfb, bits, huff_bits, save_preflag, save_compress;
    double xfsf[4][CBLIMIT];
    int i, over, iteration;
    float max_ratio;


/* D576 *xrs; */ /* to eliminate warning messages from gcc compiler */
/* I576 *ix; */  /* replace this code with below. S. Shlien 15-1-97 */

  double *xrs; 
  int *ix;  
  gr_info *cod_info = &l3_side->gr[gr].ch[ch].tt;

/* xrs = (D576 *) &xr[gr][ch][0]; */ 
/* ix  = (I576 *) l3_enc[gr][ch]; */
 

xrs = (double *) &(xr[gr][ch][0]); 
ix  = (int *) &(l3_enc[gr][ch][0]);


    iteration = 0;
#ifdef PERFORM
    if(frameNum == PERFORM)
    fprintf(log_output,"\n\nframe = %d ch = %d gr= %d\n",frameNum,ch,gr);
#endif
    do 
    {
	iteration += 1;
	cod_info->part2_length = part2_length( scalefac, fr_ps, gr, ch, l3_side );
        huff_bits = max_bits - cod_info->part2_length;

#ifdef BIN_SEARCH
	if(iteration == 1)
        {
       bin_search_StepSize(max_bits,cod_info->quantizerStepSize,
          ix,xrs,cod_info); /* speeds things up a bit */
        }
#endif
#ifdef PERFORM
	if(frameNum==PERFORM)
          fprintf(log_output,"\n    Interim Results %d\n\n",iteration);
        bits = test_inner_loop( xr, l3_enc, huff_bits, cod_info, gr, ch,
         xfsf,l3_xmin);
#else
        bits = inner_loop( xr, l3_enc, huff_bits, cod_info, gr, ch );
#endif

        calc_noise( &xr[gr][ch][0], &l3_enc[gr][ch][0], cod_info, xfsf ); /* distortion calculation */

        for ( sfb = 0; sfb < CBLIMIT; sfb++ ) /* save scaling factors */
            scalesave_l[sfb] = scalefac->l[gr][ch][sfb];

        for ( sfb = 0; sfb < SFB_SMAX; sfb++ )
            for ( i = 0; i < 3; i++ )
                scalesave_s[sfb][i] = scalefac->s[gr][ch][sfb][i];
        
        save_preflag  = cod_info->preflag;
        save_compress = cod_info->scalefac_compress;

        preemphasis( &xr[gr][ch][0], xfsf, l3_xmin, gr, ch, l3_side );



#ifdef PERFORM 
    if(frameNum == PERFORM)  
    {
    fprintf(log_output,"\nbits = %d  huff_bits= %d ", bits,huff_bits);
    fprintf(log_output," max_bits = %d\n",max_bits);
    fprintf(log_output,"Stepsize = %f ",cod_info->quantizerStepSize);
    fprintf(log_output," scale_bits    = %d \n", cod_info->part2_length );  
    print_scalefacs(log_output,scalefac,cod_info->block_type,gr,ch);
    /*if (gr==0 && ch==0)
      print_quantized_values(log_output,l3_enc[gr][ch] ,cod_info);*/
    max_ratio = worst_xfsf_to_xmin_ratio(l3_xmin,xfsf,cod_info->block_type,gr,ch);
    fprintf(log_output,"max_ratio = %6.2f\n",max_ratio );
    print_ratios(log_output,l3_xmin,xfsf,cod_info->block_type,gr,ch);
    fprintf(log_output,"\n\n");