之前的文章讲了使用经典滤波器来抑制噪声。
Review：噪声抑制之经典滤波器篇
里面提到，“用经典滤波器抑制噪声，非常简单。如果噪声的功率谱PSD和有用信号功率谱PSD没有重叠的话，那可以实现非常好的效果。
但是，如果有重叠，去噪的效果就不是特别理想了。因为在复指数信号空间里面，没办法把有用信号和噪声信号分离啊。”
正是由于“噪声谱和有用信号谱可能重叠”，所以发展了维纳滤波器。
前面的文章对维纳滤波器的设计也讲过了。
Review：维纳滤波器的设计
这篇文章，就是真实地来操作一下，设计一个维纳滤波器来抑制噪声。
因为没有去录音，所以噪声源还是matlab里的randn产生的高斯过程的数据。
再假设高斯过程并不是直接加入有用数据的，而是经过了一个“信道”，发生了一些变化，比如，AR过程。
然后就是分两步了。
1. 训练过程，用已知的训练信号和已知的接收信号，通过解方程，求得滤波器系数；
2. 滤波器系数不变，用这组系数对此后接收到的信号进行滤波，实现噪声抑制。
代码比较简单，所以直接出来了。
主函数。
%% 维纳去噪，基于训练的% 通过解方程求滤波器系数% 作者：qcy% 版本说明：补充滤波阶段% 版本：v1.1 % 时间：10月29日18:18:21% 版本说明：训练阶段% 版本：v1.0% 时间：10月29日16:52:37clear;close all;clc%% 读入音频filedir=[]; % 设置路径filename='bluesky1.wav';% 设置文件名fle=[filedir filename]; % 构成完整的路径和文件名[s, fs] = audioread(fle);% 读入数据文件s=s-mean(s); % 消除直流分量s=s/max(abs(s)); % 幅值归一N=length(s); % 语音长度time=(0:N-1)/fs; % 设置时间刻度%% 产生噪声并加入SNR = 0; % 设置信噪比r2=randn(size(s)); % 产生随机噪声b=fir1(31,0.5);% 设计FIR滤波器,代替Hr21=filter(b,1,r2); % FIR滤波% [signal_with_noise,noise]=add_noisedata(s,data,fs,fs1,snr);[r1,r22]=add_noisedata(s,r21,fs,fs,SNR);% 产生带噪语音，信噪比为SNR%% 解方程% 【训练阶段】h_length = 100; % 这个东西我目前所知道的，就只有凭感觉去设置了。desired_signal = s;observed_signal = r1;h = my_weiner_filter1( h_length,desired_signal,observed_signal );% 滤波% 用维纳滤波器，作用在接收到的训练信号上。看看效果-->这还是属于训练阶段y = filter(h,1,r1);output = y;snr1=SNR_singlech(s,r1);% 计算初始信噪比snr2=SNR_singlech(s,output); % 计算滤波后的信噪比snr=snr2-snr1;%% 画图figure;subplot 311; plot(time,s,'k'); ylabel('幅值');ylim([-1 1 ]); title('原始语音信号');subplot 312; plot(time,r1,'k'); ylabel('幅值') ylim([-1 1 ]); title('带噪语音信号');subplot 313; plot(time,output,'k'); ylim([-1 1 ]); title('滤波输出语音信号');xlabel('时间/s'); ylabel('幅值')%% 打印SNRfprintf('[训练]\n',snr1);fprintf('滤波前 SNR = %f [dB] \n',snr1);fprintf('滤波后 SNR = %f [dB] \n',snr2);fprintf('提升 %f [dB] \n',snr);fprintf('\n');%% 听效果% sound(s,fs); % 干净的语音% sound(r1,fs); % 含噪的语音% sound(output,fs); % 滤波后的语音%% 【降噪阶段】% 假设这是后来录的一段音，混入了性质类似的噪声。% 现在想用刚刚得到的滤波器系数，来去除掉现在这段含噪信号中的噪声。filedir=[]; % 设置路径filename='risesun.wav';% 设置文件名fle=[filedir filename]; % 构成完整的路径和文件名[s, fs] = audioread(fle);% 读入数据文件s=s-mean(s); % 消除直流分量s=s/max(abs(s)); % 幅值归一N=length(s); % 语音长度time=(0:N-1)/fs; % 设置时间刻度%% 产生噪声并加入% SNR = 0; % 设置信噪比r2=randn(size(s)); % 产生随机噪声b=fir1(31,0.5);% 设计FIR滤波器,代替Hr21=filter(b,1,r2); % FIR滤波% [signal_with_noise,noise]=add_noisedata(s,data,fs,fs1,snr);[r1,r22]=add_noisedata(s,r21,fs,fs,SNR);% 产生带噪语音，信噪比为SNR%% 用上一阶段求得的h，降噪% 滤波y = filter(h,1,r1);output = y;snr1=SNR_singlech(s,r1);% 计算初始信噪比snr2=SNR_singlech(s,output); % 计算滤波后的信噪比snr=snr2-snr1;%% 画图figure;subplot 311; plot(time,s,'k'); ylabel('幅值');ylim([-1 1 ]); title('原始语音信号');subplot 312; plot(time,r1,'k'); ylabel('幅值') ylim([-1 1 ]); title('带噪语音信号');subplot 313; plot(time,output,'k'); ylim([-1 1 ]); title('滤波输出语音信号');xlabel('时间/s'); ylabel('幅值')%% 打印SNRfprintf('去噪 \n',snr);fprintf('滤波前 SNR = %f [dB] \n',snr1);fprintf('滤波后 SNR = %f [dB] \n',snr2);fprintf('提升 %f [dB] \n',snr);%% 听效果% sound(s,fs); % 干净的语音% sound(r1,fs); % 含噪的语音sound(output,fs); % 滤波后的语音
维纳滤波器设计的函数。
function h = my_weiner_filter1( h_length,desired_signal,observed_signal )%function h = my_weiner_filter1( h_length,desired_signal,observed_signal )% 维纳滤波器的实现% 输入参数% h_length: 返回的FIR滤波器的长度% desired_signal: 所期望的信号,训练信号,干净信号% observed_signal: 观测到的信号% 返回参数% h: FIR滤波器系数%% 作者：qcy% 版本：v1.0% 时间：10月29日18:50:56% 0. 定义线性方程组的大小row_number = h_length;col_number = row_number;% 1. Rx --> observed_signalM = col_number;% lags = -(N-1):(N-1);Rx_c_full = xcorr(observed_signal);[~,k] = max(Rx_c_full);Rx_c = Rx_c_full(k:k+M-1);Rx_c = Rx_c.';% 2. RdxRdx_c_full = xcorr(desired_signal,observed_signal);Rdx_c = Rdx_c_full(k:k+M-1);% 3. 求h, Ax = b% (1) 生成自相关A阵A = toeplitz(Rx_c,Rx_c);b = Rdx_c;h = A\b;end
其中，加入噪声的代码如下。
function [signal,noise]=add_noisedata(s,data,fs,fs1,snr)s=s(:);% 把信号转换成列数据s=s-mean(s); % 消除直流分量sL=length(s);% 求出信号的长度if fs~=fs1 % 若纯语音信号的采样频率与噪声的采样频率不相等x=resample(data,fs,fs1);% 对噪声重采样，使噪声采样频率与纯语音信号的采样频率相同elsex=data;endx=x(:);% 把噪声数据转换成列数据x=x-mean(x); % 消除直流分量xL=length(x);% 求噪声数据长度if xL>=sL % 如果噪声数据长度与信号数据长度不等，把噪声数据截断或补足x=x(1:sL);elsedisp('Warning: 噪声数据短于信号数据，以补0来补足！')x=[x; zeros(sL-xL,1)];endSr=snr;Es=sum(x.*x);% 求出信号的能量Ev=sum(s.*s);% 求出噪声的能量a=sqrt(Ev/Es/(10^(Sr/10)));% 计算出噪声的比例因子noise=a*x; % 调整噪声的幅值signal=s+noise; % 构成带噪语音
主要是要注意采样率，还有要根据信噪比重新调整噪声幅度。
效果示意图。
训练阶段。
滤波阶段。
求解方程得到的滤波器的频率响应。
看得出来，因为是语音嘛，所以，设计出来的滤波器的性质可能更多的还是低通…
总之，这是一种在整个频率范围内的，均方意义下的最优……
最后，因为维纳滤波器的前提是，信号与噪声不相关，信号与噪声宽平稳，……
反正约束条件还是挺多的，所以，虽然对信号的噪声抑制效果，比经典滤波器的要好，
但是，为了得到更好的效果，还需要利用更先进的技术。比如，之后要讲的，自适应滤波器（adaptive filter）。
谢谢观赏。
敬请期待。
