noise_est.c

ears-0.32, linux下有用的语音信号处理工具包

				if (j != mmax_index)
					average1 += buffer->values[i][j];
			}
			average1 /= (float)(histo_step_framesize - 1);
			buffer->values[i][mmax_index] = average1;
			for(j = histo_step_framesize; j< FILT_LENGTH; j++)
				buffer->values[i][j] = average1;
        }
	}  /* end of if (i_call == histo_step_framesize) */
	
	if (i_call > histo_step_framesize)
	{
        for(i = pptr->first_good; i<lastfilt; i++)
        {
			
            if(det(buffer, i, filter, aspectrum) == 0)
            {
				/* spectral magnitude value is entered into the histogram */
				/* the oldest value in the histogram is stored so that it can be excluded 
				   from the histogram later */
				histo_index->values[i]++;
				m1 = histo_index->values[i]  % histo_length;
				waste_trajec->values[i] = histo_values->values[i][m1];
				histo_values->values[i][m1] = all_trajec->values[i][t];
            }
            else
				/* spectral magnitude value is not entered into the histogram */
				/* keep the oldest value in the histogram */ 
				waste_trajec->values[i] = -100;
        }
	}  /* end of if (i_call > histo_step_framesize) */  
	
	/* computation of the noise power for the first 100 ms */
	if ((i_call == 1) || (i_call == histo_step_framesize))
	{
        /* use parseval theorem in the beginning */
        for(i=pptr->first_good, sum = 0.0; i<lastfilt; i++)
        {
			for(j= 0, average = 0.0; j < i_call; j++)
			{
                average += all_trajec->values[i][j];
			}
			average /= i_call;
			sum += (average * average);
        }
        sum /= ( lastfilt-(pptr->first_good));
        noiseptr->values[0] = sum; 
	}  /* end of computation of the noise power for the first 100 ms */

	/* The x-axis of the histogram is scaled by the maximum spectral magnitude value. For e.g. 
	   given spectral magnitude value x and maximum value of y, the i th bin of the histogram 
	   is incremented by one where i = x/y * HISTO_RESOLUTION  */ 
	
	/* See if there is a new maximum in the trajectory */
	for (i = pptr->first_good; i< lastfilt; i++)
	{
        m1 = histo_index->values[i] % histo_length;
        /* if new_max_index->values[i] == m1, the maximum value is the oldest value of the histogram 
           which should be excluded from the histogram. Thus a new maximum is found */
        if ((new_max_index->values[i] == (short)m1)  && (i_call != 1) && (waste_trajec->values[i] != -100) )
        {
			new_max->values[i] = SET;
			for (j= 0, max->values[i] = 0; j < histo_length; j++)
			{
				if (histo_values->values[i][j] > max->values[i])
				{
					new_max_index->values[i] = (short)j;
					max->values[i] = histo_values->values[i][j];
				} 
			}
        }
        else
        {
			if( (waste_trajec->values[i] != -100) && (histo_values->values[i][m1] > max->values[i]))
			{
				new_max->values[i] = SET;
				max->values[i] = histo_values->values[i][m1]; 
				new_max_index->values[i] = (short)m1;
			}
        }
	}  /* end of for (i = pptr->first_good; i< lastfilt; i++)  */
	/* end of finding the new max */

	/* After 200ms, the histogram method is used for noise estimation */        
	if (i_call >= (2* histo_step_framesize))
	{
		for ( i= pptr->first_good; i< lastfilt; i++)
		{
			m1 = (histo_index->values[i]) % histo_length;
			/* decide the length of the histogram */
			if (histo_index->values[i] >= histo_length)
			{
                length = histo_length;
			}
			else
			{
                length = histo_index->values[i]+1; 
			}
			if (new_max->values[i] == SET)
			{
				/* completely new computation of the histogram */
				for (k=0; k < HISTO_RESOLUTION; k++) stat_lt[i][k] = 0;
				for (k=0; k < length; k++)
				{
					stat_lt[i][(int)( ((histo_values->values[i][k])/(max->values[i]))*(HISTO_RESOLUTION-EPSILON))]++;
				}
				new_max->values[i] = RESET;
			}
			else if (waste_trajec->values[i] != -100) 
			{
				stat_lt[i][(int)(((histo_values->values[i][m1])/(max->values[i]))*(HISTO_RESOLUTION-EPSILON))]++;
				if (histo_index->values[i] >=  histo_length)
				{
					/* the oldest value of the histogram is thrown away */
					stat_lt[i][(int)(((waste_trajec->values[i])/(max->values[i]))*(HISTO_RESOLUTION-EPSILON))]--;
				}
			}
            
			if ((i_call%histo_step_framesize) == 0) 
				/* Find the most frequently occurring value from the histogram every histo_step_framesize frames */ 
			{
				/* compute the actual average signal energy in this band */
				/* The possible estimated noise level is limited to the average spectral value. This is related
				   to the fact that no noise energy can occur which is higher than the total amount of energy
				   inside a band. */
				average = 0.0;
				if (i_call > histo_length)
                    length = histo_length;
				else
                    length = i_call;
				for(j=0; j<length; j++)
				{
					average += all_trajec->values[i][j];
				}
				average /= length;
				mmax_lp = 0;
				/* Find the most frequently occuring value from the histogram */
				for (k=0; k <  HISTO_RESOLUTION; k++)
				{
                    if ( stat_lt[i][k] > mmax_lp)
                    {
						mmax_lp = stat_lt[i][k];
						mindex = k;
                    }
				}   
				noise_level_lt->values[i] =   (mindex+0.5)/HISTO_RESOLUTION * (max->values[i]); 
				if (noise_level_lt->values[i] > average) noise_level_lt->values[i] = average;
				if (i_call == 2*histo_step_framesize)
					hptr->noiseOLD[i] = noise_level_lt->values[i];
				/* The noise level is smoothed by the preceding noise level */
				noise_level_lt->values[i] = 0.5*hptr->noiseOLD[i] + 0.5*noise_level_lt->values[i];
				hptr->noiseOLD[i] = noise_level_lt->values[i];
			}   /* end of if (i_call%histo_step_framesize == 0) */
		}   /* end loop over all frequencies */
		
		/* Noise levels for different subbands are squared and then the average of these 
		   squared values gives the noise power estimate. The noise power is computed 
		   every 100ms */
		if (i_call%histo_step_framesize == 0)
		{
            for(i = pptr->first_good, noisepower = 0.0; i<lastfilt; i++)
            {
				noisepower += pow(noise_level_lt->values[i],2.0);
            }
            noisepower /=  (lastfilt - (pptr->first_good)); 
		} /* end of if(i_call%histo_step_framesize == 0) */
		noiseptr->values[0] = noisepower; 
	} /* end of if i_call >= (2* histo_step_framesize)  */
	
	/* ------------------------------------------------------------------------------------------*/
	if ( (i_call == 1 || i_call%histo_step_framesize == 0) && (pptr->debug==TRUE) )
		fprintf(outdebug,"%d %f %e\n",i_call, noiseptr->values[0], 1/(C_CONSTANT*(noiseptr->values[0])));
	/*-------------------------------------------------------------------------------------------*/
	return(noiseptr);
}  /* end of subroutine comp_noisepower */            


/* create_filt is a subroutine for finding the coefficients of the filter for high energy 
   speech detection */
void create_filt(float over_est_factor, int time, float time_constant, struct fvec *filter)
{
	int i;
	float sum = 0;
	
	
	filter->values[0] = 1;
	for(i = 1; i< time; i++)
	{
		filter->values[i] = (float) (-pow((double)time_constant, (double)i));
		sum += fabs((double) filter->values[i]);
	}
	sum /= over_est_factor;
	for (i = 1; i< time; i++)
		filter->values[i] /= sum;
}


/* det is a subroutine for detecting high energy speech */
int det(struct fmat *buffer, int Crband,  struct fvec *filter, struct fvec *aspectrum)
{
	int j;
	float output;
	int flag;
	
	output = filter->values[0] * (float)sqrt((double)aspectrum->values[Crband]);
	for (j= 1; j<FILT_LENGTH; j++)
		output += filter->values[j]* buffer->values[Crband][j-1];
	/* When output >= 0 , high energy speech  is detected */
	for (j= FILT_LENGTH - 1; j>0; j--)
		buffer->values[Crband][j] = buffer->values[Crband][j-1];
	if (output < 0)
	{
		buffer->values[Crband][0] = (float)sqrt((double)aspectrum->values[Crband]);
		flag = 0;
	}
	else
	{
		buffer->values[Crband][0] = 1.05 * buffer->values[Crband][1];
		flag = 1;
	}
	return(flag);
}