FIR数字滤波器设计_MATLAB语音滤波实验02

通过FIR滤波器也可以实现对一段带噪语音的选频滤波：

clc,clear,close all;
%%读取音频并进行频谱分析
[x0,fs]=audioread('test.wav');
[x1,fs1]=audioread('airplane.mp3');
%下一步添加噪声，并自动对齐长度，向短对齐
x = x0(:,1) + x1(1:length(x0(:,1)),1); %信号len_x0 < len_x1
player = audioplayer(x,fs);
play(player)
%FFT变换
y=fft(x);    %对信号进行快速Fourier变换
mag=abs(y);    %求得Fourier变换后的振幅
figure(1)
subplot(411);
plot(mag);   %绘出随频率变化的振幅
grid; 
xlabel('频率/Hz');
ylabel('振幅');
title('原始信号FFT变换频谱');

%FIR低通滤波器
fp=1000;   %通带截止频率
fs=1200;   %阻带截止频率
At=100;    %带外衰减系数(dB)
wave=1;    %带内纹波系数(dB)
mval=[1 0];  %边界处的幅值
fedge=[fp fs];
derta1=1-10^(-wave/20); %δ1的计算
derta2=10^(-At/20);     %δ2的计算
dev=[derta1 derta2];   %基本参数δ1与δ2
fs=8000;     %FIR采样频率
%参数计算
[N,fpts,mag,wt]=remezord(fedge,mval,dev,fs);
%雷米兹交替算法设计最优FIR滤波器
b=remez(N,fpts,mag,wt);
%FIR滤波器单位脉冲响应
[h,w]=freqz(b,1,1024);
figure(2)
subplot(311);
plot(w*(fs/2)/pi,20*log10(abs(h)));
grid; 
xlabel('频率/Hz')
ylabel('幅度/dB')
title('FIR低通滤波器单位脉冲响应');

%FIR低通滤波器滤波
sf1=filter(b,1,x);
figure(1)
subplot(412);
x=sf1;      %信号
y=fft(x); %对信号进行快速Fourier变换
mag=abs(y);   %求得Fourier变换后的振幅
plot(mag);   %绘出随频率变化的振幅
grid; 
xlabel('频率/Hz');
ylabel('振幅');
title('FIR低通滤波后FFT变换频谱');

%FIR高通滤波器单位脉冲响应
fp=2800;fs=3000;Rp=1;As=100;   % 滤波器指标
fb=[fp,fs];m=[0,1]; % 计算remezord函数所需参数f,m,dev
derta1=1-10^(-Rp/20);  %δ1的计算
derta2=10^(-As/20); %δ2的计算
dev=[derta2 derta1];
[n,fo,mo,W]=remezord(fb,m,dev,8000);   % 确定remez函数所需参数
hn=remez(n,fo,mo,W);  % 调用remez函数进行设计,用于滤除噪声nt中的低频成分
%yt=filter(hn,1,1024);      %滤除随机噪声中低频成分，生成高通噪声yt
[h,w]=freqz(hn,1,1024);
figure(2)
subplot(312);
plot(w*(4000)/pi,20*log10(abs(h)));
grid; 
xlabel('频率/Hz')
ylabel('幅度/dB')
title('FIR高通滤波器单位脉冲响应');

%FIR高通滤波器滤波
sf2=filter(hn,1,x);
figure(1)
subplot(413);
x=sf2; %信号
y=fft(x);  %对信号进行快速Fourier变换
mag=abs(y);  %求得Fourier变换后的振幅
plot(mag); %绘出随频率变化的振幅
grid; 
xlabel('频率/Hz');
ylabel('振幅');
title('FIR高通滤波后FFT变换频谱');

%FIR带通滤波器
fp1=1200;   %通带截止频率
fc1=1000;   %阻带截止频率
fp2=3000;
fc2=3200;
At=100;     %带外衰减系数(dB)
wave=1;     %带内纹波系数(dB)
mval=[0 1 0];  %边界处的幅值
fedge=[fc1 fp1 fp2 fc2];
derta1=1-10^(-wave/20);  %δ1的计算
derta2=10^(-At/20);      %δ2的计算
dev=[derta2 derta1 derta2];  %基本参数δ1与δ2
fs=8000;   %FIR采样频率
%参数计算
[N,fpts,mag,wt]=remezord(fedge,mval,dev,fs);
%雷米兹交替算法设计最优FIR滤波器
b=remez(N,fpts,mag,wt);

%FIR滤波器单位脉冲响应
[h,w]=freqz(b,1,1024);
figure(2)
subplot(313);
plot(w*(fs/2)/pi,20*log10(abs(h)));
grid; 
xlabel('频率/Hz')
ylabel('幅度/dB')
title('FIR带通滤波器单位脉冲响应');

%FIR带通滤波器滤波
sf3=filter(b,1,x);
figure(1)
subplot(414);
x=sf3; %信号
y=fft(x);    %对信号进行快速Fourier变换
mag=abs(y);    %求得Fourier变换后的振幅
plot(mag);   %绘出随频率变化的振幅
grid; 
xlabel('频率/Hz');
ylabel('振幅');
title('FIR带通滤波后FFT变换频谱');

效果：

参考来源：少茗——MATLAB语音FIR滤波器处理