噪声估计对语音滤波/语音增强效果好坏很重要。若噪声估计过小，滤波后有烦人的残留噪声；若噪声估计过高，滤波结果会导致语音失真。我这实验，根据不同噪声环境的说话场合，采用可选的其中一种噪声估计算法配合谱减滤波，滤波后语音清晰度较好。 7种可选噪声估计算法是：
1.martin 最小跟踪算法
2.mcra 最小控制递归平均算法
3.mcra2 最小控制递归平均变体算法
4.imcra 改进的最小控制递归平均算法
5.doblinger 连续谱最小跟踪算法
6.hirsch 加权频谱平均值算法
7.conn_freq 连接时间_频率区域算法
参考书：《语音增强--理论与实践》 作者：罗艾洲
% -----------------------------------------------
% 我编辑的主程序： ns_specsubs.m
close; clear ; clc;
nbits =16;
alpha=0.9; %(幂指数)过减因子,值小则去噪效果增强但失真度增大
FLOOR=0.002;
load methods.mat;
disp('以下是7种噪声估计算法：');
disp('1:martin, 2:mcra, 3:mcra2, 4:imcra,');
disp('5:doblinger, 6:hirsch, 7:conn_freq');
while 1
i = input('请选择噪声估计算法的编号 1～7: (输入0 则退出) ') ;
if i == 0
return;
end
while (i<0 | i>7 )
i = input('请选择噪声估计算法的编号 1～7: (输入0 则退出) ') ;
end
method =char(meth(i));
[x0,fs] =audioread('signal.mp3'); % 读入带人群噪声的语音文件
x0 = x0(:,1); % 取1声道
x0 = x0-mean(x0); % 消除直流分量
x0 = x0/max(abs(x0)); % 幅值归一化
x = x0 ;
head_time =0.9; % 裁掉起始段无语音的秒数，依具体情况调整
head = fs*head_time;
x =x(head+1: length(x));
% 初始化变量 ===============
len=floor(40*fs/1000); % 样本的大小
if rem(len,2)==1
len=len+1;
end
PERC=50; % 窗口重叠以帧大小的百分比表示
len1=floor(len*PERC/100);
len2=len-len1;
win=hamming(len); % 汉明窗
k=1;
nFFT=2*len;
img=sqrt(-1);
x_old=zeros(len1,1);
Nframes=floor(length(x)/len2)-1;
xfinal=zeros(Nframes*len2,1);
for n=1:Nframes %从第1帧至最末帧，逐幅处理
insign=win.*x(k:k+len-1); % Windowing 加窗
spec=fft(insign,nFFT); % 计算一帧的傅立叶变换
sig=abs(spec);
ns_ps=sig.^2;
% 估计/更新噪声的psd --------
if n == 1
parameters = initialise_parameters(ns_ps,fs,method);
else
parameters = noise_estimation(ns_ps,method,parameters);
end
noise_ps = parameters.noise_ps;
noise_mu=sqrt(noise_ps); % 噪声幅值谱
theta=angle(spec); %保存每帧的相位信息
% 计算信噪比：
SNRseg=10*log10(norm(sig,2)^2/norm(noise_mu,2)^2);
% 依据过减因子alpha及信噪比选择下限系数beta:
if alpha==1.0
beta=berouti1(SNRseg);
%beta = 4 ;
else
beta=berouti2(SNRseg);
end
% 功率谱减去噪：
sub_speech=sig.^alpha - beta*noise_mu.^alpha;
diffw = sub_speech-FLOOR*noise_mu.^alpha;
% 将值朝负无穷方向舍入到最邻近的整数组件：
z=find(diffw <0);
if~isempty(z)
sub_speech(z)=FLOOR*noise_mu(z).^alpha;
end
% 复数的实部重构确保共轭对称
sub_speech(nFFT/2+2:nFFT)=flipud(sub_speech(2:nFFT/2));
% 将全帧fft与相位信息相乘：
x_phase=(sub_speech.^(1/alpha)).*(cos(theta)+img*(sin(theta)));
xi=real(ifft(x_phase)); % 逆傅立叶变换
% 重叠相加：
xfinal(k:k+len2-1)=x_old+xi(1:len1);
x_old=xi(1+len1:len);
k=k+len2;
end % EOF for n=1:Nframes
out = xfinal;
N =length(x); %含噪信号长度
Nout=length(out); %滤波后信号长度
if Nout>N
out=out(1:N);
elseif Nout<N
%尾部补零，使输出信号与与原带噪信号等长
out=[out; zeros(N-Nout,1)];
end
out=out/max(abs(out)); %输出信号幅值归一化
B = 30; % debug的输入参数，B值依具体情况调整
out2 = debug(out,B); %清除残留噪声
% 绘处理前后对比波形图：
figure('position',[50,50,600,700] );
subplot(4,1,1);
plot(x0(1:length(x0)));
title('原带噪语音波形');
axis([0 N, -1 1]); grid on;
subplot(4,1,2);
plot(x(1:length(x)));
axis([0 N, -1 1]); grid on;
subplot(4,1,3);
plot(out(1:length(out)));
title([ num2str(method) '估计噪声算法>谱减降噪后的波形' ]);
axis([0 N, -1 1]); grid on;
subplot(4,1,4);
plot(out2(1:length(out2)));
title('清除残留噪声后的波形' );
axis([0 N, -1 1]); grid on;
% sound(x0,fs) ;%播放原带噪语音
% sound(out,fs) ;%播放去噪后的语音
sound(out2,fs) ;%播放清除残留噪声后的语音
end % EOF while 1
% ----------------------------------------
处理前后的语音波形对比图:
|
楼主
标题置顶
标题高亮
