tns.c

jpeg and mpeg 编解码技术源代码

 
#include "tns.h"

/***********************************************/
/* TNS Profile/Frequency Dependent Parameters  */
/***********************************************/
static unsigned long tnsSupportedSamplingRates[13] = 
  {8000,11025,12000,16000,22050,24000,32000,44100,48000,64000,88200,96000,0};


       
/* 将带宽限制在2.0 kHz之内 */
static unsigned short tnsMinBandNumberLong[12] = 
  { 31, 30, 28, 24, 26, 25, 20, 17, 16, 15, 12, 11 };
static unsigned short tnsMinBandNumberShort[12] = 
  { 12, 10, 10, 8, 6, 6, 4, 3, 3, 2, 2, 2 };
       
/**************************************/
/* 采用基本配置和低复杂性配置时的TNS 参数 */
/**************************************/
static unsigned short tnsMaxBandsLongMainLow[12] = 
  { 39, 42, 42, 42, 46, 46, 51, 42, 40, 34, 31, 31 };

static unsigned short tnsMaxBandsShortMainLow[12] = 
  { 14, 14, 14, 14, 14, 14, 14, 14, 14, 10, 9, 9 };

static unsigned short tnsMaxOrderLongMain = 20;
static unsigned short tnsMaxOrderLongLow = 12;
static unsigned short tnsMaxOrderShortMainLow = 7;

/**************************************/
/*采用可变采样率配置时的TNS 参数*/
/**************************************/
static unsigned short tnsMaxBandsLongSSR[12] = 
  { 19, 23, 23, 23, 29, 29, 26, 26, 26, 27, 28, 28 };

static unsigned short tnsMaxBandsShortSSR[12] = 
  { 7, 8, 8, 8, 7, 7, 6, 6, 6, 7, 7, 7 };

static unsigned short tnsMaxOrderLongSSR = 12;
static unsigned short tnsMaxOrderShortSSR = 7;


/*****************************************************/
/* InitTns:                                          */
/*****************************************************/
void TnsInit(long samplingRate,enum AAC_PROFILE profile,TNS_INFO* tnsInfo) 
{
  int fsIndex=0;

  /* 判决采样率是否被支持*/
  while (samplingRate!=tnsSupportedSamplingRates[fsIndex]) {
    fsIndex++;
  }
  
  switch( profile ) {
    case MAIN :
      tnsInfo->tnsMaxBandsLong = tnsMaxBandsLongMainLow[fsIndex];
      tnsInfo->tnsMaxBandsShort = tnsMaxBandsShortMainLow[fsIndex];
      tnsInfo->tnsMaxOrderLong = tnsMaxOrderLongMain;
      tnsInfo->tnsMaxOrderShort = tnsMaxOrderShortMainLow;
      break;
    case LOW :
      tnsInfo->tnsMaxBandsLong = tnsMaxBandsLongMainLow[fsIndex];
      tnsInfo->tnsMaxBandsShort = tnsMaxBandsShortMainLow[fsIndex];
      tnsInfo->tnsMaxOrderLong = tnsMaxOrderLongLow;
      tnsInfo->tnsMaxOrderShort = tnsMaxOrderShortMainLow;
      break;
    case SSR :
      tnsInfo->tnsMaxBandsLong = tnsMaxBandsLongSSR[fsIndex];
      tnsInfo->tnsMaxBandsShort = tnsMaxBandsShortSSR[fsIndex];
      tnsInfo->tnsMaxOrderLong = tnsMaxOrderLongSSR;
      tnsInfo->tnsMaxOrderShort = tnsMaxOrderShortSSR;
      break;
  }
  tnsInfo->tnsMinBandNumberLong = tnsMinBandNumberLong[fsIndex];
  tnsInfo->tnsMinBandNumberShort = tnsMinBandNumberShort[fsIndex];
}


/*****************************************************/
/* TnsEncode:                                        */
/*****************************************************/
int TnsEncode(int numberOfBands,       /* 每一段的频带数 */
	       int maxSfb,              /* 最大sfb */
	       enum WINDOW_TYPE blockType,   /* 窗口块的类型（长段，短段） */
	       int* sfbOffsetTable,     /* 比例因子偏移表 */
	       double* spec,            /* 频谱数据数组 */
	       TNS_INFO* tnsInfo)       /* TNS信息 */

{
	int numberOfWindows,windowSize;
	int startBand,stopBand,order;    /* 要进行TNS的频带 */
	int lengthInBands;               /* 要进行滤波的长度，以频带为单位*/
	int w, error;
	int startIndex,length;
	double gain;

	switch( blockType ) {
	case ONLY_SHORT_WINDOW :
		numberOfWindows = NSHORT;
		windowSize = SN2;
		startBand = tnsInfo->tnsMinBandNumberShort;
		stopBand = numberOfBands; 
		lengthInBands = stopBand-startBand;
		order = tnsInfo->tnsMaxOrderShort;
		startBand = min(startBand,tnsInfo->tnsMaxBandsShort);
		stopBand = min(stopBand,tnsInfo->tnsMaxBandsShort);
		break;

	default:
		numberOfWindows = 1;
		windowSize = LN2;
		startBand = tnsInfo->tnsMinBandNumberLong;
		stopBand = numberOfBands;
		lengthInBands = stopBand - startBand;
		order = tnsInfo->tnsMaxOrderLong;
		startBand = min(startBand,tnsInfo->tnsMaxBandsLong);
		stopBand = min(stopBand,tnsInfo->tnsMaxBandsLong);
		break;
	}
	
	/* 确认起始和结束的频带数小于maxSfb 并且大于等于0*/
	startBand = min(startBand,maxSfb);
	stopBand = min(stopBand,maxSfb);
	startBand = max(startBand,0);
	stopBand = max(stopBand,0);

	tnsInfo->tnsDataPresent=0;     /* 缺省为没有使用 TNS */

#if 1
	if (blockType != ONLY_SHORT_WINDOW)
		/* 对每一段进行分析并滤波 */
	for (w=0;w<numberOfWindows;w++) {

		TNS_WINDOW_DATA* windowData = &tnsInfo->windowData[w];
		TNS_FILTER_DATA* tnsFilter = windowData->tnsFilter;
			double* k = tnsFilter->kCoeffs;    /*反射系数 */
		double* a = tnsFilter->aCoeffs;    /* 预测系数 */


		windowData->numFilters=0;
		windowData->coefResolution = DEF_TNS_COEFF_RES;
		startIndex = w * windowSize + sfbOffsetTable[startBand];
		length = sfbOffsetTable[stopBand] - sfbOffsetTable[startBand];
		gain = LevinsonDurbin(order,length,&spec[startIndex],k);

		if (gain>DEF_TNS_GAIN_THRESH) {  /* Use TNS */
			int truncatedOrder;
			windowData->numFilters++;
			tnsInfo->tnsDataPresent=1;
			tnsFilter->direction = 0;
			tnsFilter->coefCompress = 0;
			tnsFilter->length = lengthInBands;
			QuantizeReflectionCoeffs(order,DEF_TNS_COEFF_RES,k,tnsFilter->index);
			truncatedOrder = TruncateCoeffs(order,DEF_TNS_COEFF_THRESH,k);
			tnsFilter->order = truncatedOrder;
			StepUp(truncatedOrder,k,a);    /* Compute predictor coefficients */
			error = TnsInvFilter(length,&spec[startIndex],tnsFilter);      /* Filter */
			if (error == MBERROR)
				return MBERROR;
		}
	}
	return MBNO_ERROR;
#endif
}


/*****************************************************/
/* TnsFilter:                                        */
/*   对给定的频谱进行滤波          */
/*****************************************************/
void TnsFilter(int length,double* spec,TNS_FILTER_DATA* filter) {
	int i,j,k=0;
	int order=filter->order;
	double* a=filter->aCoeffs;

	/* Determine loop parameters for given direction */
	if (filter->direction) {

		/* Startup, initial state is zero */
		for (i=length-2;i>(length-1-order);i--) {
			k++;
			for (j=1;j<=k;j++) {
				spec[i]-=spec[i+j]*a[j];
			}
		}
		
		/* Now filter completely inplace */
		for	(i=length-1-order;i>=0;i--) {
			for (j=1;j<=order;j++) {
				spec[i]-=spec[i+j]*a[j];
			}
		}


	} else {

		/* Startup, initial state is zero */
		for (i=1;i<order;i++) {
			for (j=1;j<=i;j++) {
				spec[i]-=spec[i-j]*a[j];
			}
		}
		
		/* Now filter completely inplace */
		for	(i=order;i<length;i++) {
			for (j=1;j<=order;j++) {
				spec[i]-=spec[i-j]*a[j];
			}
		}
	}
}


/********************************************************/
/* TnsInvFilter:                                        */
/*   对指定的频谱进行反向滤波           */
/********************************************************/
int TnsInvFilter(int length,double* spec,TNS_FILTER_DATA* filter) {
	int i,j,k=0;
	int order=filter->order;
	double* a=filter->aCoeffs;
	double* temp;

    temp = (double *) malloc( length * sizeof (double));
    if (!temp) {
		return MBERROR;
//      CommonExit( 1, "TnsInvFilter: Could not allocate memory for TNS array\nExiting program...\n");
    }

	/* Determine loop parameters for given direction */
	if (filter->direction) {

		/* Startup, initial state is zero */
		temp[length-1]=spec[length-1];
		for (i=length-2;i>(length-1-order);i--) {
			temp[i]=spec[i];
			k++;
			for (j=1;j<=k;j++) {
				spec[i]+=temp[i+j]*a[j];
			}
		}
		
		/* Now filter the rest */
		for	(i=length-1-order;i>=0;i--) {
			temp[i]=spec[i];
			for (j=1;j<=order;j++) {
				spec[i]+=temp[i+j]*a[j];
			}
		}


	} else {

		/* Startup, initial state is zero */
		temp[0]=spec[0];
		for (i=1;i<order;i++) {
			temp[i]=spec[i];
			for (j=1;j<=i;j++) {
				spec[i]+=temp[i-j]*a[j];
			}
		}
		
		/* Now filter the rest */
		for	(i=order;i<length;i++) {
			temp[i]=spec[i];
			for (j=1;j<=order;j++) {
				spec[i]+=temp[i-j]*a[j];
			}
		}
	}
	free(temp);
	return MBNO_ERROR;
}





/*****************************************************/
/* TruncateCoeffs:                                   */
/*   对绝对值小于某一个阈值的系数进行截断，并返回截断
后的滤波器的阶数 */

/*****************************************************/
int TruncateCoeffs(int fOrder,double threshold,double* kArray) {
	int i;

	for (i=fOrder;i>=0;i--) {
		kArray[i] = (fabs(kArray[i])>threshold) ? kArray[i] : 0.0;
		if (kArray[i]!=0.0) return i;
	}
}

/*****************************************************/
/* QuantizeReflectionCoeffs:                         */
/*   将给定的反射系数量化为指定的阶数
/*****************************************************/
void QuantizeReflectionCoeffs(int fOrder,
			      int coeffRes,
			      double* kArray,
			      int* indexArray) {

	double iqfac,iqfac_m;
	int i;

	iqfac = ((1<<(coeffRes-1))-0.5)/(PI/2);
	iqfac_m = ((1<<(coeffRes-1))+0.5)/(PI/2);

	/* Quantize and inverse quantize */
	for (i=1;i<=fOrder;i++) {
		indexArray[i] = (int)(0.5+(asin(kArray[i])*((kArray[i]>=0)?iqfac:iqfac_m)));
		kArray[i] = sin((double)indexArray[i]/((indexArray[i]>=0)?iqfac:iqfac_m));
	}
}

/*****************************************************/
/* Autocorrelation,                                  */
/*   针对给定的数据计算自相关函数
/*****************************************************/
void Autocorrelation(int maxOrder,        /* Maximum autocorr order */
		     int dataSize,		  /* Size of the data array */
		     double* data,		  /* Data array */
		     double* rArray) {	  /* Autocorrelation array */

  int order,index;

  for (order=0;order<=maxOrder;order++) {
    rArray[order]=0.0;
    for (index=0;index<dataSize;index++) {
      rArray[order]+=data[index]*data[index+order];
    }
    dataSize--;
  }
}



/*****************************************************/
/* LevinsonDurbin:                                   */
/*   采用LevinsonDurbin迭代算法根据给定的数据计算反射
系数，并返回预测增益。         */
/*****************************************************/
double LevinsonDurbin(int fOrder,          /* 滤波器的阶数 */
					  int dataSize,		   /* 数据数组的大小 */
					  double* data,		   /* 数据数组*/
					  double* kArray) 	   /* 反射系数数组*/

{
	int order,i;
	double signal;
	double error, kTemp;				/* 预测误差 */
	double aArray1[TNS_MAX_ORDER+1];	/* 预测器系数数组1 */
	double aArray2[TNS_MAX_ORDER+1];	/* 预测器系数数组2 */
	double rArray[TNS_MAX_ORDER+1];		/* 自相关系数*/
	double* aPtr = aArray1;				/* 指向aArray1的指针 */
	double* aLastPtr = aArray2;			/*指向 aArray2 的指针*/

	double* aTemp;

	/* 计算自相关系数*/
	Autocorrelation(fOrder,dataSize,data,rArray);
	signal=rArray[0];	/* signal energy */

	/* 将指针指向当前和上一次迭代*/ 
	aPtr = aArray1;
	aLastPtr = aArray2;
	/* 如果没有信号能量，则返回*/
	if (!signal) {
		kArray[0]=1.0;
		for (order=1;order<=fOrder;order++) {
			kArray[order]=0.0;
		}
		return 0;
	} else {

		/* 建立第一个迭代*/
		kArray[0]=1.0;
		aPtr[0]=1.0;		/* 指向预测器系数, 当前循环*/
		aLastPtr[0]=1.0;	/* 指向预测器系数, 上一次循环*/

		error=rArray[0];

		for (order=1;order<=fOrder;order++) {
			kTemp = aLastPtr[0]*rArray[order-0];
			for (i=1;i<order;i++) {
				kTemp += aLastPtr[i]*rArray[order-i];
			}
			kTemp = -kTemp/error;
			kArray[order]=kTemp;
			aPtr[order]=kTemp;
			for (i=1;i<order;i++) {
				aPtr[i] = aLastPtr[i] + kTemp*aLastPtr[order-i];
			}
			error = error * (1 - kTemp*kTemp);

				/* 将当前循环设为最后一次循环*/
			aTemp=aLastPtr;
			aLastPtr=aPtr;		
			aPtr=aTemp;			
		}
		return signal/error;	/* 返回增益值 */
	}
}


/*****************************************************/
/* StepUp:                                           */
/*   将反射系数转化为预测器的系数                   */
/*****************************************************/
void StepUp(int fOrder,double* kArray,double* aArray) {

	double aTemp[TNS_MAX_ORDER+2];
	int i,order;

	aArray[0]=1.0;
	aTemp[0]=1.0;
	for (order=1;order<=fOrder;order++) {
		aArray[order]=0.0;
		for (i=1;i<=order;i++) {
			aTemp[i] = aArray[i] + kArray[order]*aArray[order-i];
		}
		for (i=1;i<=order;i++) {
			aArray[i]=aTemp[i];
		}
	}
}