endpointdetection.c

这是一个语音去除噪声的程序

/* ---------------------------------------------
The correlating functions in SpeechRecognition 
EndPoint Detection 
------------------------------------------------
*/

#include <math.h>
#include <stdio.h>
#include "EndPointDetection.h"

void fft(float *fAr,float *fAi);
void FmelFilter(pFilter pFilter1);
void FreqEner(OrignalData *pOrignaldata1, EnergyBand *pEnergyband1, Filter *pFilter);

//void SequBEner(EnergyBand *pEnergyband1);

float SequBEner(EnergyBand *pEnergyband1);//add on Mar 5th,2006
int NusefulBand(float fEmin, EnergyBand *pEnergyband1);//add on Mar 5th,2006
//int NusefulBand(float fNoiseFreq);

//float FrameEntropy(int iNub, EnergyBand *pEnergyband1);

float FrameEntropy(int iNub, float fEmin, EnergyBand *pEnergyband1);//add on Mar 5 th,2006

//float Mean(float *fHFrame);

float MeanVariance(float *fHFrame, float *mean);


/*-----------------  FFT ------------------
	            短时傅立叶变换
-------------------------------------------
*/
void fft(float *fAr,float *fAi)              
{
	int i, j, k, iB;
	float tr, ti, p;
	int iLh, iN1;
	int iL;

	iLh = HALFOFFRAME;
	j  = iLh;
	iN1 = FRAMELENGTH-2;

	for(i=1; i<=iN1; i++)
	{
		if(i >= j)
			k = iLh;
		else 
		{
			tr = fAr[i];
			ti = fAi[i];
			fAr[i] = fAr[j];
			fAi[i] = fAi[j];
			fAr[j] = tr;
			fAi[j] = ti;
			k = iLh;
		}

		while ( j >= k )
		{
			j = j - k;
			k = k/2;
		}
		j = j + k;
	}

	for (iL=0; iL<MLOG; iL++)
	{
		iB = (int)pow(2, iL);

		for (j=0; j<iB; j++)
		{
			p = (float)(pow(2, MLOG-iL-1) * j*2.0*PI/ (float)FRAMELENGTH);
			for (k=j; k<FRAMELENGTH; k=k+2*iB )
			{
				tr = fAr[k+iB] * (float)cos(p) + fAi[k+iB] * (float)sin(p);
				ti = fAi[k+iB] * (float)cos(p) - fAr[k+iB] * (float)sin(p);
				fAi[k+iB] = fAi[k] - ti;
				fAr[k+iB] = fAr[k] - tr;
				fAr[k] = fAr[k] + tr;
				fAi[k] = fAi[k] + ti;			
			}
		}
	}
}

/*-----------------------FmelFilter----------------------------------
                   构造mel刻度的滤波器组
---------------------------------------------------------------------
*/
void FmelFilter(pFilter pFilter1)
{
	float fFmelLow,fFmelHigh,fFmelBandWid;
	float fFlow,fHigh;
	float fFmel[BANDSUM + 2],fFrequency[BANDSUM + 2];
	float fBoundaryPoint[BANDSUM + 2];
	float fFrameTime;
	int i,j,k;

	fFlow = 1 + (float)FLOW/(float)700;
	fHigh = 1 + (float)FHIGH/(float)700;
	fFmelLow = (float)(1127 * log(fFlow));
	fFmelHigh = (float)(1127 * log(fHigh));
	fFmelBandWid = (fFmelHigh-fFmelLow)/(float)(BANDSUM + 1);

  	for(i = 0; i <= BANDSUM + 1; i++)
	{
		fFmel[i] = fFmelLow + i * fFmelBandWid;
		fFrequency[i] = (float)(700 * (exp(fFmel[i]/1127) - 1));
		fFrameTime =(float)FRAMELENGTH/(float)SAMPLING;
		fBoundaryPoint[i] = fFrequency[i] * fFrameTime;
	}
	
	j = 1;
	for(i = 0; i < BANDSUM; i ++)
	{
		for(k = 0; k <HALFOFFRAME; k ++)
		{
			if(k >= fBoundaryPoint[j-1] && k <= fBoundaryPoint[j])
			{	pFilter1->fFilter[i][k] = \
				((float)k - fBoundaryPoint[j-1])/(fBoundaryPoint[j] - \
				fBoundaryPoint[j-1]);
			}
			else if(k >= fBoundaryPoint[i] && k <= fBoundaryPoint[i + 1])
			{
				pFilter1->fFilter[i][k] = (fBoundaryPoint[i + 1] - \
				(float) k)/(fBoundaryPoint[i + 1] - fBoundaryPoint[i]);
			}
			else pFilter1->fFilter[i][k] = 0;

		}

		j++;
	}
	
}

/* ---------------------FreqEner-------------------------------
     计算给定帧的能量，每个带的能量 fEb，每帧的总能量fEsum	   
---------------------------------------------------------------------
*/
void FreqEner(OrignalData *pOrignaldata1, EnergyBand *pEnergyband1, Filter *pFilter)
{

	float fEpoint[HALFOFFRAME],fE[HALFOFFRAME];
	int i, k;
	
//	fEpoint[0] = 0; // debug on Mar 6th,2006
	for( k=0; k<HALFOFFRAME; k++)
	{
		fE[k] = (float)pow((double)(pOrignaldata1->fAr[k]), 2) + (float)pow((double)(pOrignaldata1->fAi[k]), 2);
		fEpoint[k] = (float)sqrt((double)fE[k]);
	}

	for(i=0; i<BANDSUM; i++)
	{
		pEnergyband1->fEb[i] = 0;
		for(k = 0; k<HALFOFFRAME; k++)
		{
			pEnergyband1->fEb[i] = pEnergyband1->fEb[i] + fEpoint[k] * pFilter->fFilter[i][k];
		}

	}

	pEnergyband1->fEsum = 0;
	for(i=0; i<BANDSUM; i++)
	{
		pEnergyband1->fEsum = pEnergyband1->fEsum + pEnergyband1->fEb[i];
	}

}

/* ----------------------SequBEner---------------
	            给定帧中带能量按从大到小排序
 ----------------------------------------------------
*/
//void SequBEner(EnergyBand *pEnergyband1)
float SequBEner(EnergyBand *pEnergyband1)//add on Mar 5th,2006
{

	int i,j,k;
	float t;
	float fEmin;
	
	for (i=0; i<BANDSUM; i++)
	{
		k = i;
		for(j = i + 1;j<BANDSUM;j ++)
		{
			if(pEnergyband1->fEb[j] > pEnergyband1->fEb[k])  k = j;
		}
		
		if (k != i)
		{
			t = pEnergyband1->fEb[i];
			pEnergyband1->fEb[i] = pEnergyband1->fEb[k];
			pEnergyband1->fEb[k] = t;
		}

	}

	fEmin = pEnergyband1->fEb[BANDSUM - 1];
	return fEmin;

}

/*---------------------NusefulBand----------------
	          找出给定帧中有用的带数iNub
  ------------------------------------------------
*/

//add on Mar 5th,2006
int NusefulBand(float fEmin, EnergyBand *pEnergyband1)
{

	float fI, fP;
	int iNub;

	fP = fEmin / pEnergyband1->fEsum;
	fI = - (float)(log10(fP)) / (float)(log10(2));

	if (fI <= 5) 
		iNub = 30;
	else if(fI>5 && fI<25) 
		iNub = (int)(36.5 - 1.3 * fI);
	else 
		iNub = 4;

	return iNub;
}

/*-----------------------------------------------------
int NusefulBand(float fBNoiseFreq)
{

	int iNub;

	iNub = (int)(-1.5 * fBNoiseFreq + 142.5);
	if( iNub<3 ) iNub = 3;
	if( iNub>18 ) iNub = 18;
	
	return iNub;
}
--------------------------------------------------------*/
/*--------------------FrameEntropy----------------------
	               计算该帧的熵dHFrame
--------------------------------------------------------
*/
//float FrameEntropy(int iNub, EnergyBand *pEnergyband1)
float FrameEntropy(int iNub, float fEmin, EnergyBand *pEnergyband1)//add on Mar 5 th,2006
{

	float fHFrame = 0, fPb[BANDSUM];
	float fPboffset[BANDSUM];//add on Mar 5th,2006
	float fW[BANDSUM], fMean[BANDSUM];//add on Mar 5th,2006
	int i;

/*-------------------debug on Mar 5th, 2006---------	
	for (i=0; i<iNub; i++)
	{
		fPb[i] = pEnergyband1->fEb[i] / pEnergyband1->fEsum;

		fHFrame = fHFrame - fPb[i] * (float)(log10(fPb[i])) / (float)(log10(2));
	}
---------------------------------------------------*/

//add on Mar 5th, 2006
/*---------------------------------------------------
	          每个子带能量所占该帧总能量的概率
----------------------------------------------------*/
	for ( i = 0; i<BANDSUM; i ++ )
	{
		fPb[i] = pEnergyband1->fEb[i] / pEnergyband1->fEsum;
		fPboffset[i] = fEmin / pEnergyband1->fEb[i];
	}

/*----------------------------------------------------
                   各子带加权系数
-----------------------------------------------------*/
	fW[0] = 1;
	fW[BANDSUM - 1] = 1;
	for( i = 1; i<BANDSUM - 1; i ++ )
	{
		fMean[i] = (fPboffset[i - 1] +  fPboffset[i] + fPboffset[i + 1]) / 3;
		fW[i] = ((float)(pow((fPboffset[i - 1] - fMean[i]),2)) + (float)(pow((fPboffset[i] - fMean[i]),2)) + (float)(pow((fPboffset[i + 1] - fMean[i]),2))) / 3;
		fW[i] = (float)sqrt((double)fW[i]);
	}
/*-----------------------------------------------------
                       各帧的熵
------------------------------------------------------*/
	for(i = 0; i<iNub; i++)
	{
		fHFrame = fHFrame - fPb[i] * fW[i] * (float)(log10(fPb[i])) / (float)(log10(2));
	}
	return fHFrame;
}
/*------------------------Mean----------------------
	                  计算均值fMean
--------------------------------------------------------
*/
/*float Mean(float *fHFrame)
{
	float fMean;	
	float fHb = 0;
	float *fHFrame1 = fHFrame;
	int  i;

	if ( fHFrame == NULL )			// add by ling
	{
		printf("There is some input pointer error!\r\n");
		return 0;
	}
	for (i = 0; i<FRAMESUM; i++)
	{
		fHb = fHb + *fHFrame1;
		fHFrame1--;
	}


	fMean = fHb / FRAMESUM;

	return fMean;
}
*/
/*--------------------MeanVariance----------------------
	       该函数计算均值fMean和方差fVariance
--------------------------------------------------------
*/
float MeanVariance(float *fHFrame,float *mean)
{
	float fVariance;	
	float fHb = 0;
	float fAmid, fBmid, fSum = 0;
	int  i;
	float *fHFrame1 = fHFrame;

	for (i=0; i<FRAMESUM; i++)
	{
		fHb = fHb + *fHFrame1;
//		fHFrame1++;
		fHFrame1--;		// debug on Mar 4th, 2006
	}


//	*mean = fHb / (FRAMESUM + 1);
	*mean = fHb / FRAMESUM ;	// debug on Mar 4th, 2006

	for (i=0; i<FRAMESUM; i++)
	{
/*		fAmid = fHFrame[i] - *mean;
		fBmid = (float)(pow(fAmid, 2));
		fSum = fSum + fBmid;             */

		fAmid = *fHFrame1 - *mean;// debug on Mar 4th, 2006
		fBmid = (float)(pow(fAmid, 2));// debug on Mar 4th, 2006
		fSum = fSum + fBmid; // debug on Mar 4th, 2006
		fHFrame1 ++ ;		// debug on Mar 4th, 2006
	}

//	fVariance = fSum / (FRAMESUM + 1);
	fVariance = fSum / FRAMESUM ; //	 debug on Mar 4th, 2006
	fVariance = (float)sqrt((float)fVariance);//	 debug on Mar 4th, 2006

	return fVariance;
}