aac_qc.c

AAC音频解码算法程序

#include "aacenc.h" 
#include "bitstream.h"
#include "tf_main.h"
#include "pulse.h"
#include "aac_qc.h"
#include "aac_se_enc.h"
#include <math.h>

#include "hufftab5.h"
#include "all.h"


double pow_quant[9000];
double adj_quant[9000];
int sign[1024];
int g_Count;
int old_startsf;
int pns_sfb_start = 1000;         /* lower border for Perceptual Noise Substitution 
                                      (off by default) */

//double ATH[SFB_NUM_MAX];

/*double ATHformula(double f)
{
	double ath;
	f  = max(0.02, f);
	/* from Painter & Spanias, 1997 */
	/* minimum: (i=77) 3.3kHz = -5db */
/*	ath=(3.640 * pow(f,-0.8)
		-  6.500 * exp(-0.6*pow(f-3.3,2.0))
		+  0.001 * pow(f,4.0));
	
	/* convert to energy */
/*	ath = pow( 10.0, ath/10.0 );
	return ath;
}
 */

/*void compute_ath(AACQuantInfo *quantInfo, double ATH[SFB_NUM_MAX])
{
	int sfb,i,start=0,end=0;
	double ATH_f;
	double samp_freq = 44.1;
	static int width[] = {0, 4,  4,  4,  4,  4,  8,  8,  8, 12, 12, 12, 16, 16, 16};
	if (quantInfo->block_type==ONLY_SHORT_WINDOW) {
		for ( sfb = 0; sfb < 14; sfb++ ) {
			start = start+(width[sfb]*8);
			end   = end+(width[sfb+1]*8);
			ATH[sfb]=1e99;
			for (i=start ; i < end; i++) {
				ATH_f = ATHformula(samp_freq*i/(128)); /* freq in kHz */
/*				ATH[sfb]=min(ATH[sfb],ATH_f);
			}
		}
	} else {
		for ( sfb = 0; sfb < quantInfo->nr_of_sfb; sfb++ ) {
			start = quantInfo->sfb_offset[sfb];
			end   = quantInfo->sfb_offset[sfb+1];
			ATH[sfb]=1e99;
			for (i=start ; i < end; i++) {
				ATH_f = ATHformula(samp_freq*i/(1024)); /* freq in kHz */
/*				ATH[sfb]=min(ATH[sfb],ATH_f);
			}
		}
	}
}*/


void tf_init_encode_spectrum_aac( int quality )
{
	int i;

	g_Count = quality;
	old_startsf = 0;

	for (i=0;i<9000;i++){
		pow_quant[i]=pow(i, ((double)4.0/(double)3.0));
	}
	for (i=0;i<8999;i++){
		adj_quant[i] = (i + 1) - pow(0.5 * (pow_quant[i] + pow_quant[i + 1]), 0.75);
	}
}

int quantize(AACQuantInfo *quantInfo,
			 double *p_spectrum,
			 double *pow_spectrum,
			 int quant[NUM_COEFF]
			 )
{
	int i, sb;
	double quantFac;
	double x;

	for (sb = 0; sb < quantInfo->nr_of_sfb; sb++) {

		quantFac = pow(2.0, 0.1875*(quantInfo->scale_factor[sb] -
			quantInfo->common_scalefac ));

		for (i = quantInfo->sfb_offset[sb]; i < quantInfo->sfb_offset[sb+1]; i++){
			x = pow_spectrum[i] * quantFac;
			if (x > MAX_QUANT) {
				return 1;
			}
			quant[i] = (int)(x + adj_quant[(int)x]);
			quant[i] = sgn(p_spectrum[i]) * quant[i];  /* restore the original sign */
		}
	}

	if (quantInfo->block_type==ONLY_SHORT_WINDOW)
		quantInfo->pulseInfo.pulse_data_present = 0;
	else
		PulseCoder(quantInfo, quant);

	return 0;
}

int calc_noise(AACQuantInfo *quantInfo,
				double *p_spectrum,
				int quant[NUM_COEFF],
				double requant[NUM_COEFF],
				double error_energy[SFB_NUM_MAX],
				double allowed_dist[SFB_NUM_MAX],
				double *over_noise,
				double *tot_noise,
				double *max_noise
				)
{
	int i, sb, sbw;
	int over = 0, count = 0;
	double invQuantFac;
	double linediff;

	*over_noise = 0.0;
	*tot_noise = 0.0;
	*max_noise = -999.0;

	if (quantInfo->block_type!=ONLY_SHORT_WINDOW)
		PulseDecoder(quantInfo, quant);

	for (sb = 0; sb < quantInfo->nr_of_sfb; sb++) {

		double max_sb_noise = 0.0;

		sbw = quantInfo->sfb_offset[sb+1] - quantInfo->sfb_offset[sb];

		invQuantFac = pow(2.0,-0.25 * (quantInfo->scale_factor[sb] - quantInfo->common_scalefac ));

		error_energy[sb] = 0.0;

		for (i = quantInfo->sfb_offset[sb]; i < quantInfo->sfb_offset[sb+1]; i++){
			requant[i] =  pow_quant[ABS(quant[i])] * invQuantFac; 

			/* measure the distortion in each scalefactor band */
			linediff = (double)(ABS(p_spectrum[i]) - ABS(requant[i]));
			linediff *= linediff;
			error_energy[sb] += linediff;
			max_sb_noise = max(max_sb_noise, linediff);
		}
		error_energy[sb] = error_energy[sb] / sbw;		

		if (allowed_dist[sb] != 0.0)
			error_energy[sb] = 10*log10(error_energy[sb] / allowed_dist[sb]);
		else error_energy[sb] = 0;
		if (error_energy[sb] > 0) {   //ERROR_ENERGY[sb] > allowed_dist[sb]
			over++;
			*over_noise += error_energy[sb];
		}
		*tot_noise += error_energy[sb];
		*max_noise = max(*max_noise, error_energy[sb]);
		count++;
  	}

	if (count>1) *tot_noise /= count;
	if (over>1) *over_noise /= over;
	return over;
}

int quant_compare(int best_over, double best_tot_noise, double best_over_noise,
				  double best_max_noise, int over, double tot_noise, double over_noise,
				  double max_noise)
{
	/*
	noise is given in decibals (db) relative to masking thesholds.

	over_noise:  sum of quantization noise > masking
	tot_noise:   sum of all quantization noise
	max_noise:   max quantization noise 

	*/
	int better=0;

#if 0
	better = ((over < best_over) ||
			((over==best_over) && (over_noise<best_over_noise)) ) ;
#else
#if 0
	better = max_noise < best_max_noise;
#else
#if 0
	better = tot_noise < best_tot_noise;
#else
#if 0
	better = (tot_noise < best_tot_noise) &&
		(max_noise < best_max_noise + 2);
#else
#if 0
	better = ( ( (0>=max_noise) && (best_max_noise>2)) ||
		( (0>=max_noise) && (best_max_noise<0) && ((best_max_noise+2)>max_noise) && (tot_noise<best_tot_noise) ) ||
		( (0>=max_noise) && (best_max_noise>0) && ((best_max_noise+2)>max_noise) && (tot_noise<(best_tot_noise+best_over_noise)) ) ||
		( (0<max_noise) && (best_max_noise>-0.5) && ((best_max_noise+1)>max_noise) && ((tot_noise+over_noise)<(best_tot_noise+best_over_noise)) ) ||
		( (0<max_noise) && (best_max_noise>-1) && ((best_max_noise+1.5)>max_noise) && ((tot_noise+over_noise+over_noise)<(best_tot_noise+best_over_noise+best_over_noise)) ) );
#else
#if 1
	better =   (over_noise <  best_over_noise)
		|| ((over_noise == best_over_noise)&&(tot_noise < best_tot_noise));
#if 0
	better = (over_noise < best_over_noise)
		||( (over_noise == best_over_noise)
		&&( (max_noise < best_max_noise)
		||( (max_noise == best_max_noise)
		&&(tot_noise <= best_tot_noise)
		)
		) 
		);
#endif
#endif
#endif
#endif
#endif
#endif
#endif

	return better;
}


int count_bits(AACQuantInfo* quantInfo,
			   int quant[NUM_COEFF],
			   int output_book_vector[SFB_NUM_MAX*2])
{
	int i, bits = 0;

	/* find a good method to section the scalefactor bands into huffman codebook sections */
	bit_search(quant,              /* Quantized spectral values */
		quantInfo);         /* Quantization information */

    /* Set special codebook for bands coded via PNS  */
    if (quantInfo->block_type != ONLY_SHORT_WINDOW) {     /* long blocks only */
		for(i=0;i<quantInfo->nr_of_sfb;i++) {
			if (quantInfo->pns_sfb_flag[i]) {
				quantInfo->book_vector[i] = PNS_HCB;
			}
		}
    }

	/* calculate the amount of bits needed for encoding the huffman codebook numbers */
	bits += sort_book_numbers(quantInfo,             /* Quantization information */
		output_book_vector,    /* Output codebook vector, formatted for bitstream */
		NULL,          /* Bitstream */
		0);                    /* Write flag: 0 count, 1 write */

	/* calculate the amount of bits needed for the spectral values */
	quantInfo -> spectralCount = 0;
	for(i=0;i< quantInfo -> nr_of_sfb;i++) {  
		bits += output_bits(
			quantInfo,
			quantInfo->book_vector[i],
			quant,
			quantInfo->sfb_offset[i], 
			quantInfo->sfb_offset[i+1]-quantInfo->sfb_offset[i],
			0);
	}

	/* the number of bits for the scalefactors */
	bits += write_scalefactor_bitstream(
		NULL,             /* Bitstream */  
		0,                        /* Write flag */
		quantInfo
		);

	/* the total amount of bits required */
	return bits;
}

int loop_count = 0;
int if_success = 0;

int tf_encode_spectrum_aac(
			   double      *p_spectrum[MAX_TIME_CHANNELS],   //输入1：MDCT的能量谱向量
			   double      *PsySigMaskRatio[MAX_TIME_CHANNELS], //输入2：心理声学模型输出:信掩比
			   double      allowed_dist[MAX_TIME_CHANNELS][MAX_SCFAC_BANDS], //输入2:允许的最大误差
			   double      energy[MAX_TIME_CHANNELS][MAX_SCFAC_BANDS],  //？
			   enum WINDOW_TYPE block_type[MAX_TIME_CHANNELS],	//块类型
			   int         sfb_width_table[MAX_TIME_CHANNELS][MAX_SCFAC_BANDS],  //输入5：SCALEFACTOR BAND WIDTH
			   int         nr_of_sfb[MAX_TIME_CHANNELS],	//sfb的数目
			   int         average_block_bits,           //输入3：平均可用比特数
			   int         available_bitreservoir_bits,  //输入4：可用比特数
			   int         padding_limit,
			   BsBitStream *fixed_stream,			//输出：
			   BsBitStream *var_stream,				//输出：
			   int         nr_of_chan,
			   double      *p_reconstructed_spectrum[MAX_TIME_CHANNELS],
			   int         useShortWindows,
			   int aacAllowScalefacs,
			   AACQuantInfo* quantInfo,      /* AAC quantization information */ 
			   Ch_Info* ch_info,
			   int varBitRate,
			   int bitRate)              //输入定义见标准P134
{
	int quant[NUM_COEFF];
	int i=0;
	int j=0;
	int k;
	double max_dct_line = 0;
	int global_gain;
	int store_common_scalefac;
	int best_scale_factor[SFB_NUM_MAX];
	double pow_spectrum[NUM_COEFF];
	double requant[NUM_COEFF];
	int sb;
	int extra_bits;
	int max_bits;
	int output_book_vector[SFB_NUM_MAX*2];
	int start_com_sf;
	double SigMaskRatio[SFB_NUM_MAX];
	IS_Info *is_info;
	int *ptr_book_vector;

	int inner_loop_count;
	char message[255];

	/* Set up local pointers to quantInfo elements for convenience */
	int* sfb_offset = quantInfo -> sfb_offset;
	int* scale_factor = quantInfo -> scale_factor;
	int* common_scalefac = &(quantInfo -> common_scalefac);

	int outer_loop_count;
	int quantizer_change;
	int over = 0, best_over = 100, better;
	int sfb_overflow, amp_over;
	int best_common_scalefac;
	int prev_scfac, prev_is_scfac;
	double noise_thresh;
	double over_noise, tot_noise, max_noise;
	double noise[SFB_NUM_MAX];
	double best_max_noise = 0;
	double best_over_noise = 0;
	double best_tot_noise = 0;
	static int init = -1;

	/* Set block type in quantization info */
	quantInfo -> block_type = block_type[MONO_CHAN];

#if 0
	if (init != quantInfo->block_type) {
		init = quantInfo->block_type;
//		compute_ath(quantInfo, ATH);
	}
#endif

	/** create the sfb_offset tables **/
	if (quantInfo -> block_type == ONLY_SHORT_WINDOW) {

		/* Now compute interleaved sf bands and spectrum */
		sort_for_grouping(
			quantInfo,                       /* ptr to quantization information */
			sfb_width_table[MONO_CHAN],      /* Widths of single window */
			p_spectrum,                      /* Spectral values, noninterleaved */
			SigMaskRatio,				/*output new interleaved sigmaskratio*/
			PsySigMaskRatio[MONO_CHAN]   /*input, from psycho calculate*/
			);

		extra_bits = 51;   //??????????????????????????
	} else{
		/* For long windows, band are not actually interleaved */
		if ((quantInfo -> block_type == ONLY_LONG_WINDOW) ||  
			(quantInfo -> block_type == LONG_SHORT_WINDOW) || 
			(quantInfo -> block_type == SHORT_LONG_WINDOW)) {
			quantInfo->nr_of_sfb = quantInfo->max_sfb;

			sfb_offset[0] = 0;
			k=0;
			for( i=0; i< quantInfo -> nr_of_sfb; i++ ){
				sfb_offset[i] = k;
				k +=sfb_width_table[MONO_CHAN][i];
				SigMaskRatio[i]=PsySigMaskRatio[MONO_CHAN][i];
			}
			sfb_offset[i] = k;
			extra_bits = 100; /* header bits and more ... *//*for what?*/

		} 
	}

	extra_bits += 1;  /*why?*/

    /* Take into account bits for TNS data */
    extra_bits += WriteTNSData(quantInfo,fixed_stream,0);    /* Count but don't write */

 // if(quantInfo->block_type!=ONLY_SHORT_WINDOW)
		/* Take into account bits for LTP data */
//		extra_bits += WriteLTP_PredictorData(quantInfo, fixed_stream, 0); /* Count but don't write */

    /* for short windows, compute interleaved energy here */
    if (quantInfo->block_type==ONLY_SHORT_WINDOW) {
		int numWindowGroups = quantInfo->num_window_groups;
		int maxBand = quantInfo->max_sfb;
		int windowOffset=0;
		int sfb_index=0;
		int g;
		for (g=0;g<numWindowGroups;g++) {
			int numWindowsThisGroup = quantInfo->window_group_length[g];
			int b;
			for (b=0;b<maxBand;b++) {
				double sum=0.0;
				int w;
				for (w=0;w<numWindowsThisGroup;w++) {
					int bandNum = (w+windowOffset)*maxBand + b;
					sum += energy[MONO_CHAN][bandNum];
				}
				energy[MONO_CHAN][sfb_index] = sum;
				sfb_index++;
			}
			windowOffset += numWindowsThisGroup;
		}
    } 

	/* initialize the scale_factors that aren't intensity stereo bands */
	is_info=&(ch_info->is_info);
	for(k=0; k< quantInfo -> nr_of_sfb ;k++) {
		scale_factor[k]=((is_info->is_present)&&(is_info->is_used[k])) ? scale_factor[k] : 0 /*min(15,(int)(1.0/SigMaskRatio[k]+0.5))*/;
	}

	/* Mark IS bands by setting book_vector to INTENSITY_HCB */
	ptr_book_vector=quantInfo->book_vector;
	for (k=0;k<quantInfo->nr_of_sfb;k++) {
		if ((is_info->is_present)&&(is_info->is_used[k])) {
			ptr_book_vector[k] = (is_info->sign[k]) ? INTENSITY_HCB2 : INTENSITY_HCB;
		} else {
			ptr_book_vector[k] = 0;
		}
	}

	/* PNS prepare */
    for(sb=0; sb < quantInfo->nr_of_sfb; sb++ )
		quantInfo->pns_sfb_flag[sb] = 0;

    if (block_type[MONO_CHAN] != ONLY_SHORT_WINDOW) {     /* long blocks only */
		for(sb = pns_sfb_start; sb < quantInfo->nr_of_sfb; sb++ ) {
			/* Calc. pseudo scalefactor */
			if (energy[0][sb] == 0.0) {
				quantInfo->pns_sfb_flag[sb] = 0;
				continue;
			}

			if ((10*log10(energy[MONO_CHAN][sb]*sfb_width_table[0][sb]+1e-60)<70)||(SigMaskRatio[sb] > 1.0)) {
				quantInfo->pns_sfb_flag[sb] = 1;
				quantInfo->pns_sfb_nrg[sb] = (int) (2.0 * log(energy[0][sb]*sfb_width_table[0][sb]+1e-60) / log(2.0) + 0.5) + PNS_SF_OFFSET;
				
				/* Erase spectral lines */
				for( i=sfb_offset[sb]; i<sfb_offset[sb+1]; i++ ) {
					p_spectrum[0][i] = 0.0;
				}
			}
		}
    }

	/* Compute allowed distortion */
	for(sb = 0; sb < quantInfo->nr_of_sfb; sb++) {