matlab - 如何在Matlab中使用allpassfilter消除声音中的噪音？

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去除了虚部，然后播放了声音。声音仍然有噪音

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

但是我无法消除噪声，问题在于我不知道如何将allpassfilter中的noise(具有噪声的位置)的值设为零，如下所示:

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)

您可以尝试其他反卷积方法并查看结果。您也可以将傅立叶域中的无用频谱归零，并对原始信号进行逆傅立叶变换。但是，由于信号太长，您将不得不在滑动窗口中进行操作。另外，您可以设计级联陷波滤波器来分别滤波每个谐波。

我看到有四个Stron谐波。因此，为每个滤波器设计四个陷波滤波器和一个低通滤波器以滤除高频噪声。

% 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)

希望我能帮上忙。