Matlab实现IIR数字滤波器设计

时间：2024-01-22 12:18:59

低通滤波器的设计：

1、基本指标

　　低通滤波器顾名思义就是低频率成分通过，高频率成分截止，那么在设计一个低通滤波器时首先要明白想要截止多大的频率。比如想截掉4Hz以上的信号，理想状态下就是将4Hz以上的信号成分全部截止，4Hz以下的信号全部保留。然而事实上几乎不存在这样的滤波器，通常情况下总是在通过频率和截止频率之间存在一个过渡带。通过频率这部分称为通带，允许通过的最大频率为通带截止频率ωp，截至频率这部分称为阻带，阻带最小截至频率为ωs，通带和阻带之间的部分为过渡带，也即ωp~ωp。通带之间的波动称为通带波动δp，阻带之间的波动称为阻带波动δs如图所示：

来源：【北京交通大学陈后金教授数字信号处理课件】

2、Matlab函数介绍

声明：上面一节的基本指标均为模拟滤波器状态下的指标，Matlab不仅提供了模拟滤波器设计函数，同样提供了数字滤波器的设计函数，只是参数之间需要一些转换。

数字滤波器设计中，Wp和Ws为归一化角频率，如果工程要求设计低通数字滤波器满足通带截止频率为fp = 40，阻带截止频率fs = 50，则fp，fs和Wp，Ws之间的转换关系为：

Wp = 2*pi*fp/Fs;
Ws = 2*pi*fs/Fs;

截止频率的定义是当输出幅值响应下降到输入幅值的-3dB ( 20log(1/\sqrt{2}) )，也就是0.707(也就是 1/\sqrt{2} )时对应的频率，也称3db截频

巴特沃斯数字滤波器函数：
[n,Wc] = buttord(Wp/pi,Ws/pi,Ap,As);
[b,a] = butter(n,Wc,\'low\');%低通滤波器

Ap和As的计算用上图中的公式计算获得。

function FilteredSignal = filtered(signal,Fs,fp,fs,Ap,As,FilterType)
% fp=40; fs=50; Ap=1; As=20;
%=============================================
% Fs:信号的频率
% fp：滤波器的通带频率
% fs:滤波器的截止频率
% Ap：通带最大衰减
% As：阻带最小衰减
% FilterType:滤波器类型，\'low\',低通，\'high\',高通
% FilteredSignal：滤波后的信号
%==============================================
Wp = 2*pi*fp/Fs;
Ws = 2*pi*fs/Fs;
[n,Wc] = buttord(Wp/pi,Ws/pi,Ap,As);
[b,a] = butter(n,Wc,FilterType);
omega=[Wp Ws]; 
h = freqz(b,a,omega); %Compute Ap and As of AF
fprintf(\'Ap= %.4f\n\',-20*log10(abs(h(1))));
fprintf(\'As= %.4f\n\',-20*log10(abs(h(2))));
FilteredSignal = filter(b,a,signal);

　当 fp=40; fs=50; Ap=1; As=20;时，通过上面的代码设计的滤波器Ap= 0.9998，As= 20.0000，满足要求，设计的滤波器如下图：

　

切比雪夫I型数字滤波器函数：

[N,wc] = cheb1ord(Wp/pi,Ws/pi,Ap,As);
[b,a] = cheby1(N,Ap,wc,\'low\');
%% 切比雪夫1型低通滤波器
fp=40; fs=50; Ap=1; As=20;
Wp = 2*pi*fp/Fs;
Ws = 2*pi*fs/Fs;
[N,wc] = cheb1ord(Wp/pi,Ws/pi,Ap,As);
[b,a] = cheby1(N,Ap,wc,\'low\');
omega = [Wp,Ws];
h = freqz(b,a,omega);
fprintf(\'N = %.2f\n\',N);
fprintf(\'Ap=%.4f\n\',-20*log10(abs(h(1))));
fprintf(\'As=%.4f\n\',-20*log10(abs(h(2))));
figure(7)
freqz(b,a,512,200);title(\'CBI Lowpass Filter\')

　

N = 5.00；Ap=1.0000<=1；As=20.0093>20满足定义的Ap和As，通带刚好，阻带有余量。

%% 椭圆低通滤波器
[N,wc] = ellipord(Wp/pi,Ws/pi,Ap,As);
[b,a] = ellip(N,Ap,As,wc,\'low\');
omega = [Wp,Ws];
h = freqz(b,a,omega);
fprintf(\'N = %.2f\n\',N);
fprintf(\'Ap=%.4f\n\',-20*log10(abs(h(1))));
fprintf(\'As=%.4f\n\',-20*log10(abs(h(2))));
figure
freqz(b,a,512,200);title(\'elip Lowpass Filter\')

　　

　

　　

相关文章