以前有類似主題的文章。

這一篇講的是，如果有輔助資料，用LMS演算法對訊號濾波去噪的一個應用示例。

注意，現在的信噪比已經只有-40dB了！噪聲功率是訊號功率的1萬倍！

模型還是這樣的。v2是在離噪聲源很近的地方，採集到的噪聲（因此可以假設它就是v）。

v1是因為噪聲源發出的噪聲v，經過了一段通道傳輸以後，發生了一些變化，因此變成了v1。

%% 基於LMS演算法，有輔助資料的噪聲消除
% 作者：qcy
% 版本：v1.0
% 時間：2016年10月30日15:39:25

%% 
close all;clear; 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 = -40;                              % 設定信噪比
r2=randn(size(s));                      % 產生隨機噪聲
b=fir1(31,0.5);                         % 設計FIR濾波器,代替H
r21=filter(b,1,r2);                     % FIR濾波
[r1,r22]=add_noisedata(s,r21,fs,fs,SNR);% 產生帶噪語音，信噪比為SNR　

%% LMS濾波

h_length = 35;
h = zeros(h_length,1); % 濾波器的初始化
miu = 8e-4;

v1_hat = zeros(size(s));
err2 = zeros(size(s));
y_output = r1;

% 開始濾波
for kk = 1:1 % 這裡設定濾波過程需要重複幾次。
    for k = h_length : N
        idx = k: -1 :(k-h_length+1);
        r2_in_sub = r2(idx); % 參考噪聲 --> 原始噪聲
        v1_hat_out = h.' * r2_in_sub; % 由 原始噪聲估計的 --> 訊號中的噪聲
        v1_hat(k) = v1_hat_out;
        dk = r1(k); % 觀察到的含噪語音，是參考訊號
        err = dk - v1_hat_out;
        err2(k) = err^2;
        h = h + miu * err * r2_in_sub;
        y_output(k) = y_output(k) - v1_hat_out;
    end
end

% 誤差平方
figure;
plot(err2);

%% 濾波器的幅頻特性
[H,w] = freqz(h,1,1e4);
figure;
plot(w/2/pi*fs,20*log10(abs(H)/max(abs(H))+eps));
xlabel('頻率(Hz)');

%% 計算信噪比的提升
snr1=SNR_singlech(s,r1);                % 計算初始信噪比
snr2=SNR_singlech(s,y_output);            % 計算濾波後的信噪比
snr=snr2-snr1;
SN1=snr1; SN2=snr2; SN3=snr;

%% 列印SNR
fprintf('[LMS]\n',snr1);
fprintf('濾波前 SNR = %f [dB] \n',snr1);
fprintf('濾波後 SNR = %f [dB] \n',snr2);
fprintf('提升 %f [dB] \n',snr);

%% 聽效果
% sound(s,fs); % 乾淨的語音
% sound(r1,fs); % 含噪的語音
sound(y_output,fs); % 濾波後的語音

%% 作圖
figure;
subplot 311; plot(time,s); ylabel('幅值') 
ylim([-1 1 ]); title('原始語音訊號');
subplot 312; plot(time,r1); ylabel('幅值') 
ylim([-100 100]); title('帶噪語音訊號');
subplot 313; plot(time,y_output); 
ylim([-1 1 ]); title('LMS濾波輸出語音訊號');
xlabel('時間/s'); ylabel('幅值');
 

可以看出，原始語音訊號已經完全淹沒在噪聲中了。在這麼惡劣的條件下，竟然能分辨出微弱的有用訊號，

在以前簡直是不可想象的。感謝偉大的科學家與工程師們。

LMS演算法對濾波器進行更新，當濾波器係數收斂（基本不變）後，

可以看看這個濾波器的幅頻特性。

由於有用訊號是語音訊號（300 Hz - 3400 Hz），我們的最優濾波器也很符合邏輯地顯示出其低通特性。

聽一聽效果，真是不錯。

但，還是以前說的，

1. 需要另外增加一個訊號採集裝置，增加成本；

2. 這個額外的訊號採集裝置還要能很好地接收到基本純淨的噪聲或通過一定線性通道（我們用的FIR濾波器）的噪聲；

這兩個因素，能不能實現，就要具體情況具體分析了。