如何在matlab中使用全通滤波器去除声音中的噪音？

Question

大家好，我有一个有噪音的声音。我想消除这种噪音如何消除它？

原声：Zamfir-EinsamerHirte

嘈杂的声音：Zamfir-EinsamerHirte_noisy

[y4,Fs]=audioread('Zamfir-EinsamerHirte_noisy.ogg');
ffty4=fft(y4);

首先我分析了信号

shiftedffty4=fftshift(ffty4);
spectrumy41=abs(shiftedffty4);
phaseffty41 = angle(shiftedffty4);

N4=length(spectrumy41);
t4=-Fs/2:Fs/N4:Fs/2-Fs/N4;
spectrumy42=abs(fftshift(ffty4))/N4;
phaseffty42=angle(fftshift(ffty4));

其次，我制作了一个具有相同频谱长度的全通滤波器，并产生了 fft 的嘈杂声音，并制作了逆 fft 并去除了虚部并播放了声音。声音仍然有噪音

allpassfilter=ones([N4,2]);
allpassfilter(spectrumy42>1e+06)=0;
filteredy4=allpassfilter.*ffty4;
filteredyeni4=ifft(filteredy4);
filteredyy4=real(filteredyeni4);
sound(filteredyy4,Fs);

但我无法消除噪音。问题是我不知道如何将 allpassfilter 中的噪音值（有噪音的位置）设为零，如下所示：

allpassfilter(spectrumy42>1e+06)=0;

我怎样才能做到？ ！！！任何帮助将不胜感激！！！！提前致谢。

Answer 1

我下载了干净而嘈杂的音频文件。 首先让我们分析一小部分音频。

n=1024*8; % a small portion of data
w1=1e5;
w2=w1+n-1;

sig_noisy=data_n(w1:w2,1); % noisy audio
sig_clean=data_c(w1:w2,1); % clean audio

figure; hold all
plot(sig_noisy,'b')
plot(sig_clean,'r','LineWidth',2)
ylim([-1.5 1.5])
legend('Noisy','Clean')

正如这里所见，嘈杂的音频不知何故饱和了， 干净信号的截断版本。截断信号会导致谐波 在更大的频率。让我们看看功率谱 信号的密度。

n=1024*1; % a smaller portion of data
w1=1e5;
w2=w1+n-1;

sig_noisy=data_n(w1:w2,1); % noisy 
sig_clean=data_c(w1:w2,1); % clean  

[psd_noisy, f] = pwelch(sig_noisy);
[psd_clean, ~] = pwelch(sig_clean);

figure; hold all
plot(f/pi,db(psd_noisy),'b')
plot(f/pi,db(psd_clean),'r')
xlabel('Normalized Freq.')
legend('Noisy','Clean')

您会看到嘈杂的音频在高频处具有谐波和噪声。好吧，现在如果您假设噪声的特征在音频结束时没有改变，那么您可以设计一个过滤器来查看这小部分音频。既然你已经有了噪声和干净的信号，为什么不使用反卷积方法。

例如，如果您将干净的信号与嘈杂的信号去卷积，那么 您获得系统的逆响应（h_inv），这也是您可以用来过滤噪声信号的滤波器系数

(sig_noisy = sig_clean * h).

这里我使用维纳反卷积方法。另请注意，此函数并非仅用于图像，您还可以将 Matlab 中的反卷积方法用于一维信号。

h_inv=deconvwnr(sig_clean,sig_noisy);

figure,plot(h_inv)
legend('h^-^1')

正如我所说，这是您需要的滤波器系数。例如，如果我用 h_inv 过滤噪声信号：

sig_filtered=conv(sig_noisy,h_inv,'same');
[psd_filtered, ~] = pwelch(sig_filtered);

figure; hold all
plot(f/pi,db(psd_noisy),'b')
plot(f/pi,db(psd_clean),'r')
plot(f/pi,db(psd_filtered),'k')
xlabel('Normalized Freq.')
legend('Noisy','Clean','Filtered')

滤波后的信号频谱非常接近于干净的信号频谱。现在您有了滤波器系数，只需使用 h_inv 过滤整个嘈杂的音频并收听结果。

filtered_all=conv(data_n(:,1),h_inv,'same');
sound(filtered_all,48000)

您可以尝试其他反卷积方法并查看结果。您还可以在傅立叶域中将不需要的频谱归零，并为干净的信号取逆傅立叶。但是，由于信号太长，您将不得不在滑动窗口中进行。或者，您可以设计级联陷波滤波器来分别过滤每个谐波。

我看到有四个强谐波。因此，分别设计四个陷波滤波器和一个低通滤波器来滤除高频噪声。

% First notch
fc1=0.0001; bw1=0.05; N=4;
f = fdesign.notch('N,F0,BW',N,fc1,bw1); h = design(f);

% Second notch
fc2=0.21; bw2=0.2;
f = fdesign.notch('N,F0,BW',N,fc2,bw2); h2 = design(f);

% Third notch
fc3=0.41; bw3=0.2;
f = fdesign.notch('N,F0,BW',N,fc3,bw3); h3 = design(f);

% Fourth notch
fc4=0.58; bw4=0.2;
f = fdesign.notch('N,F0,BW',N,fc4,bw4); h4 = design(f);

% A Final lowpass filter
f = fdesign.lowpass('Fp,Fst,Ap,Ast',0.6,0.65,1,30);  h5 = design(f);

% Cascade the filters
hd = dfilt.cascade(h, h2, h3, h4, h5);

% See the filter characterisctic
ff=fvtool(hd,'Color','white');

% Now we can filter our 
sig_filtered2 = filter(hd,sig_noisy);
[psd_filtered2,f] = pwelch(sig_filtered2);

figure; hold all
plot(f/pi,db(psd_noisy),'b');
plot(f/pi,db(psd_clean),'r');
plot(f/pi,db(psd_filtered2),'k');
xlabel('Normalized Freq.')
legend('Noisy','Clean','Filtered')

现在您可以过滤整个音频

filtered_all2 = filter(hd,data_n(:,1));
sound(filtered_all2,48000)

希望我能帮上忙。