基于matlab的IIR数字滤波器设计.docx
《基于matlab的IIR数字滤波器设计.docx》由会员分享，可在线阅读，更多相关《基于matlab的IIR数字滤波器设计.docx（28页珍藏版）》请在冰点文库上搜索。
基于matlab的IIR数字滤波器设计
基于matlab的IIR数字滤波器设计
一．IIR数字滤波器介绍
所谓数字滤波器，是指输入，输出均为数字信号，通过一定运算关系改变输入信号所含频率成分的相比照例或者滤除某些频率成分的硬件。
实质上就是一个由有限精度算法实现的线性时不变离散系统。
它的根本工作原理是利用离散系统的特性对系统输入信号进行加工和变换，改变输入序列的频谱或信号波形，让有用的频率分量通过，抑制无用的信号分量输出，因此数字滤波与模拟滤波的概念相同，根据其频率特性同样可以分为低通，高通，带通，带阻，只是信号的形式和实现滤波方式有所不同。
如果要处理的信号是模拟信号，就可以通过A/D或者D/A转换，在信号形式上进行匹配转换，同样可以使用数字滤波器对模拟信号进行滤波。
数字滤波器滤波的数学表达式：
y〔n〕=x（n）*h（n）;如果滤波器的输入输出信号都是离散信号，那么该滤波器的脉冲响应也一定是离散信号，这样的滤波器就成为了数字滤波器。
上面的系统为时域离散系统时，其频域特性为：
其中
分别是数字滤波器的输出序列和输入序列的频域响应，
是数字滤波器的频域响应。
可以看见按照输入信号的频谱特点和处理信号的目的适中选择滤波器的频域响应
，使得滤波后的输出信号满足设计性能要求，就是滤波器的滤波原理。
2．IIR数字滤波器传输特性
IIR数字滤波器的系统函数可以表示为：
H（Z）=
，式中H（Z）称为N阶IIR滤波器函数。
3．.数字滤波器的技术要求.
我们通常设计的数字滤波器一般属于选频滤波器，
。
我们的目的是要设计一个因果可实现的滤波器，另外买也要考虑到本钱和复杂性问题，因此实用中通带和阻带都允许一定的误差容限，即通带不一定是完全水平的，阻带也不可能完全衰减到零。
而且，通带和阻带之间还要设置一定带宽的过渡带。
如以下图表示低通滤波器的技术要求：
图中，
分别表示通带截止频率和阻带截止频率，通带频率范围为0≤w≤
，通带中要求〔1-δ1〕≤|H
≤1，阻带截止频率范围
≤w≤Π，再阻带中要求
≤δ2，从
到
称为过渡带，在这个频带内，幅度响应从通带平滑的下落到阻带。
二．IIR数字滤波器的设计方法
1．常用的IIR滤波器设计方法
〔1〕.以模拟滤波器函数为根底的变换法；即先设计一满足指定条件的模拟滤波器H（s），再将该模拟滤波器转化为数字滤波器H（z）。
〔2〕.直接设计法：
在z平面内，根据零、极点对系统特性的影响，调整零极点位置得H（z）。
〔3〕.最优化设计法（计算机辅助设计），在某种最小化误差准那么下，建立差分方程系数ak、bi对理想特性的逼近方程，使用迭代方法解方程组得到最正确逼近系统。
由于此方法计算量大，需要借助于计算机进行设计。
下面分别介绍函数设计法和信号处理图形界面来设计IIR数字滤波器，得出最优化的设计方法。
2.基于matlab的函数设计IIR数字滤波器
IIR数字滤波器设计的一般方法是先设计低通模拟滤波器，进行频率变换，将其转换为相应的〔高通，带通等〕模拟滤波器，在转换为高通，带通或带阻数字滤波器，由模拟滤波器设计数字滤波器的方法。
这是因为模拟滤波器设计方法已经很成熟，它不仅有完整的设计公式，还有完善的图表供查阅，另外，还有一些典型的滤波器类型可供我们使用。
对设计的全过程的各个步骤，MATLAB都提供了了相应的工具箱函数，使IIR数字滤波器设计变得非常简单。
.基于matlab函数的IIR滤波器设计
.1.设计方法选择：
程序设计法是基于MATLAB中相应的工具箱函数来实现的，IIR数字滤波器的设计步骤可由图所示的流程图来表示。
这个图也清晰的表示了5类20个信号处理工具箱函数的作用，
在MATLAB中，模拟滤波器的系统函数
H（S）=
数字滤波器的系统函数
H（Z）=
在实际工程中，需要的设计结果是系数向量B和A，用B和A来综合滤波器的硬件实现结构或软件运算结构，为了直观的看出设计结果，本文的实例均以滤波器幅频响应曲线作为设计结果输出。
如果需要滤波器系数，在运行程序后，只要在MATLAB命令窗口键入系数向量名，那么相应的系数就显示出来了。
.2.程序设计实例分析
〔a〕设计高通和带通Butterworth数字滤波器
我们给出四阶归一化Butterworth模拟滤波器的系统函数
用双线性变换法从Ha〔s〕设计四阶带通butterworth数字滤波器
并图示|
|,
设计采样周期T=1s，指标如下
现在我们分步进行：
■建模
由于本例主要涉及三个问题：
（1）由数字滤波器指标求相应的模拟滤波器指标；
（2）模拟滤波器频率变换〔因为已给定阶数和模拟滤波器的归一化低通原型〕；
（3）由相应的模拟滤波器到数字滤波器〔双线性变换法〕。
由于调用bilinear函数将模拟滤波器转换成数字滤波器非常容易，并且有效抑制频率失真的问题，本例给定了数字滤波器指标，所以首先要设计处与该指标相应的四阶Butterworth模拟滤波器，然后调用bilinear函数将其转换为数字滤波器即可，应当特别注意的是，对于双线性变换法，由数字边界频率求相应的模拟边界频率时，一定要考虑预畸变矫正。
只有这样，最终设计结果才能满足所给指标，
〔ⅰ〕首先按照步骤一的要求
设计高通数字滤波器时，相应的模拟高通滤波器3dB截止频率为
设计带通数字滤波器时，相应的模拟滤波器的3dB截止频率为
〔ⅱ〕步骤二的的原理
可调用MATLAB频率变换函数lp2lp，lp2hp，lp2bp，分别实现从模拟低通到模拟低通，高通，带通，带阻的频率变换。
▲【Bt，At】=lp2hp〔B，A,wc〕，将系数向量为B和A的模拟滤波器归一化低通原型〔3Db〕截止频频为1rad/s〕，变换成3dB截止频率为wc的高通模拟滤波器，返回高通模拟滤波器系数向量Bt和At。
▲【Bt，At】=lp2bp〔B,A，wo，Bw〕将系数向量为B和A的模拟滤波器归一化低通原型变换成中心频率为wo，带宽为Bw的带通模拟滤波器，返回带通模拟滤波器的系数向量Bt和At。
其中，wo=
，
由以上原理我们来编写如下程序：
%用双线性变换法设计数字高通和带通滤波器
clear;closeall
T=1;wch=pi/2;%T为采样间隔，wch位数字高通3dB截止频率
wlc=0.35*pi;wuc=0.65*pi;%wlc，wuc；数字高通3dB截止频率
B=1;A=[1,2.6131,3.4142,2.6131,1];
[h,w]=freqs（B,A,512）;%求原归一化模拟滤波器的频率响应
subplot（3,2,1）;plot（w,20*log10（abs（h）））;%画模拟滤波器幅频特性
grid;axis（[0,10,-90,0]）
xlabel（'w/'）;ylabel（'模拟低通幅度〔dB〕'）
%
（1）设计高通
omegach=2*tan（wch/2）/T;%预畸变求模拟高通3dB截止频率
[Bhs,Ahs]=lp2hp（B,A,omegach）;%模拟域低通转换为高通系数
[Bhz,Ahz]=bilinear（Bhs,Ahs,1/T）;%模拟转换位数字高通系数变量
[h,w]=freqz（Bhz,Ahs,512）;%求画出数字滤波器幅频特性
Subplot（3,2,3）;plot（w/pi,20*log10（abs（h）））;
grid;axis（[0,1,-150,0]）
xlabel（'w/'）;ylabel（'数字滤波器幅度〔dB〕'）
%
（2）设计带通
omegalc=2*tan（wlc/2）/T;%预畸变求滤波器通带低端截止频率
omegauc=2*tan（wuc/2）/T;%预畸变求滤波器通带高端截止频率
wo=sqrt（omegalc*omegauc）;Bw=omegauc-omegalc;
[Bbs,Abs]=lp2bp（B,A,wo,Bw）;%模拟域低通转换为带通系数
[Bbz,Abz]=bilinear（Bbs,Abs,1/T）;%模拟转换为数字带通系数变量
[h,w]=freqz（Bbz,Abz,512）;%求并画出数字滤波器幅频特性
subplot（3,2,4）;plot（w/pi,20*log10（abs（h）））;
grid;axis（[0,1,-150,0]）
xlabel（'w/pi'）;ylabel（'数字滤波器幅度〔dB〕'）
程序运行结果如图，，所示分别表示模拟低通幅度，数字高通幅度和数字带通幅度特性。
模拟低通幅度数字高通幅度
图图
数字带通幅度
图
〔b〕设计滤波器进行图像去噪处理
以上是我们运用MATLAB工具箱函数来编程设计IIR数字滤波器并对其频响特性进行分析的实例，下面我们就来引入数字滤波器处理图像的实例具体分析。
我们用卷积定理来说明数字滤波器进行图像去噪处理的原理，设图像信号e（x,y）通过线性不变系统h（x,y）的结果是g（x,y），即r〔x,y〕=e（x,y）*h（x,y）;由时域卷积，频域乘积定理得R〔U，V〕=E〔U，V〕H〔U，V〕。
其中，R（u,v），E（u,v），H（u,v）分别是r（x,y），e〔x,y）和h（x,y）的傅里叶变换。
实际上，图像的能量大局部集中在幅度频的低频和中频段，而图像的边缘和噪声对应于高频局部。
因此，能降低高频成分幅度的滤波器那么能过滤噪声，减弱噪声的影响。
而Butterworth低通滤波器在物理上是可以实现的滤波器，它的转移函数模的平方为
为N阶，截止频率为
。
下面我们就来设计Butterwirth滤波器来对加高斯白噪声的图像进行去噪处理
I=imread（'K:
\祖国好.jpg'）;
figure,imshow（I）,
D=imnoise（I,'gaussian'）;
figure,imshow（D）
D=double（D）;
F=fft2（D）;%傅里叶变换
F=fftshift（F）;%转换数据矩阵
[N1,N2]=size（F）;
n=2;
d0=40;
n1=fix（N1/2）;
n2=fix（N2/2）;
fori=1:
N1
forj=1:
N2
d=sqrt（（i-n1）^2+（j-n2）^2）;
h=1/（1+（d/d0）^（2*n））;%计算低通转换函数
FD（i,j）=h*F（i,j）;%低通滤波
end
end
FD=ifftshift（FD）;
FD=ifft2（FD）;
FD=uint8（real（FD））;
figure,imshow（FD）
运行结果如下，我们可以从图与图看出Butterworth滤波器能有效地过滤图像中高频加性噪声，增强图像。
2.1.3基于函数设计法的总结
从以上一系列函数设计中，我们由matlab函数来设计滤波器是次优化的，它的设计步骤为：
1．先选择设计方法
2.猜想滤波器参数，后进行设计
3.观察滤波器的响应，判断其是否符合要求
4.反复这一尝试与失败过程直到符合要求。
这种设计方法，很显然在设计要求上进行权衡分析是不是很有效，它更多时候是凭借操作者来设计的。
因此我们下面将探讨更优化的设计方法。
3.基于信号处理图形用户界面设计IIR数字滤波器
基于matlab函数的滤波器设计完成后，需要对已设计的滤波器的频率响应要进行校核。
要得到幅频、相频响应特性，运算量也是很大的。
而利用MATLAB强大的信号处理界面工具进行计算机辅助设计，可以快速有效地设计数字滤波器，大大地简化了计算量。
3.1．FDATool设计IIR数字滤波器
3.FDATool工具包的介绍及使用
FDATool（FilterDesign&AnalysisTool）是MATLAB信号处理工具箱里专用的滤波器设计分析工具，以上的版本还专门增加了滤波器设计工具箱（FilterDesignToolbox）。
FDATool可以设计几乎所有的常规滤波器，包括FIR和IIR的各种设计方法。
它操作简单，方便灵活。
FDATool界面总共分两大局部，一局部是DesignFilter。
在界面的下半部，用来设置滤波器的设计参数；另一局部那么是特性区，在界面的上半局部，用来显示滤波器的各种特性。
DesignFilter局部主要分为：
FilterType（滤波器类型）选项，包括Lowpass（低通）、Highpass（高通）、Bandpass（带通）、Bandstop（带阻）和特殊的FIR滤波器。
DesignMethod（设计方法）选项，包括IIR滤波器的Butterwotth（巴特沃思）法、Chebyshev TypeI（切比雪夫I型）法、ChebyshevTypeII（切比雪夫II型）法、Elliptic（椭圆滤波器）法和FIR滤波器的Equiripple法、Least-Squares（最小乘方）法、Window（窗函数）法。
FilterOrder（滤波器阶数）选项，定义滤波器的阶数，包括SpecifyOrder（指定阶数）和MinimumOrder（最小阶数）。
在SpecifyOrder中填入所要设计的滤波器的阶数（N阶滤波器，SpecifyOrder=N-1）。
如果选择MinimumOrder，那么MATLAB根据所选择的滤波器类型自动使用最小阶数。
FrequencySpecifications选项，可以详细定义频带的各参数，包括采样频率和频带的截止频率。
它的具体选项由FilterType选项和DesignMethod选项决定。
例如Bandpass（带通）滤波器需要定义Fstop1（下阻带截止频率）、Fpass1（通带下限截止频率）、Fpass2（通带上限截止频率）、Fstop2（上阻带截止频率），而Lowpass（低通）滤波器只需要定义Fstop1、Fpass1。
采用窗函数设计滤波器时，由于过渡带是由窗函数的类型和阶数所决定，所以只需定义通带截止频率，而不必定义阻带参数。
MagnitudeSpecifications选项，可以定义幅值衰减的情况。
例如设计带通滤波器时，可以定义Wstop1（频率Fstop1处的幅值衰减）、Wpass（通带范围内的幅值衰减）、Wstop2（频率Fstop2处的幅值衰减）。
当采用窗函数设计时，通带截止频率处的幅值衰减固定为6db，所以不必定义。
3.FDATool设计IIR数字滤波器
我们将以一个IIR滤波器的设计实例来具体说明使用matlab工具箱的方便。
要求设计elliptic〔椭圆〕带通数字低通滤波器满足以下指标：
它的通带范围从100HZ到150HZ,采样频率Fs=1000HZ,通带最大衰减Rp=1dB，阻带最大衰减
Rs=60dB,阶数为10。
在该实例中，首先在ResponseType中选择Bandpass高通滤波器，然后在下面的DesighMethod中选择IIR类型，并且指定FilterOrder项中的阶数SpecifyOrder=10，由于是设计elliptic椭圆滤波器，其下面Option就不必选择。
然后在FrequencySpecifications中选择Unit为Hz，给出采样频率Fs=1000，通带Fpass1=100和Fpass2=150；最后在MagnitudeSpecifications中选择Unit为dB，Apass=1,Astop=60.设置完成后点击DesignFilter即可得到所设计的IIR滤波器。
通过菜单项选择项Analysis可以在特性区看到所设计的幅频响应、相频响应、冲击响应和零极点配置等特性，如图1。
下面即是运用FDATOOl对elliptic椭圆滤波器的设计界面：
其中幅频特性如图1.8所示
图1.8
相频特性
图1.9
冲击响应特性
图
从以上这些界面中我们可以清晰明了的看到设计的椭圆滤波器各种特性：
由以上图中我们能够很容易的分析，图1.8中椭圆滤波器具有等纹波的通频带、等纹波的抑止频带，而且过渡带宽非常狭窄，总之，使用FDATool工具包设计和分析滤波器，是非常方便易行的，而且交互性良好，不需要极其复杂的程序编制就可以实现。
在工程中也是广泛应用。
3.3．基于fdesign更加优化的设计方法
.1.fdesign设计方法概述
Fdesign是一种面向对象的滤波器设计方法。
这种设计方法的设计思路是：
2.因为MATLAB提供可符合这些要求的设计方法，例如fvtool，sptool等，使用这些工具箱进行设计
3.然后在各种方法中找到最优化的设计方法。
.2.使用滤波器对象的优点
1.设计的filter.dfilt为对象表示，对象中包含所有滤波器的特性及可供操作的函数
2.fvtool，sptool提供了滤波器分析及视觉化集成环境。
3.方便对参数和结构进行滤波功能上的权衡分析，包括：
延时，滤波器设计复杂度，阻带衰减权衡分析，支持多种滤波器结构包括：
direct-formFIR，transposed，Overlap-addFIR以及之间的转换用convert函数进行操作。
4.仿真与自动代码生成的途径：
它可以生成simulink模型
5.自动估算计算复杂度，使用cost函数：
df.cost。
.3.fvtool设计IIR数字滤波器
FVTool可用于查看设计或导入的滤波器的特性，包括其幅度响应、相位响应、群延迟、极点-零点图、冲激响应和阶跃响应等。
下面我就来用一个程序来分析如何用fvtool设计滤波器：
%用工具箱的画图工具进行画图
clear;
Wp=60/600;Ws=90/600;Rp=1;Rs=15;%滤波器参数设定
[N,Wn]=BUTTORD（Wp,Ws,Rp,Rs）;
[B,A]=BUTTER（N,Wn）;%巴特沃斯模拟滤波器
[num,den]=iirlp2bp（B,A,0.15,[0.1,0.2]）;%IIR低通向带通转换
fvtool（B,A,num,den）;%进入fvtool界面进行滤波器设计与分析
FT=dfilt.df1（B,A）;%将要设计的滤波器参数传递给变量FT
realizemdl（FT）;%生成simulink仿真模型界面
运行此程序那么会出现fvtool界面，如下图
Butterwoth滤波器幅频特性界面
通过这个界面我们不仅可以分析道所要设计的滤波器的幅频特性，还可以分析它的相频特性，以及零极点图，单位脉冲响应，单位阶跃响应，群延迟，相位延迟以及参加高斯噪声后的频谱等等。
在实际的语言编码通信中，解调后信号和原传递信号的差异是因幅度和时间的量化而产生的，而滤波器那么会引起这种差异的产生，而这种差异就是量化噪声，在fvtool界面中我们也能分析到滤波器
的量化噪声功率谱。
如以下图所示为Round-offNoisePowerSpectrum-该滤波器的量化噪声功率谱
图2.2
我们还可以在fvtool的信息栏看到看到滤波器的各项数据如下图，我们发现它将两项不同参数的滤波器进行比拟filter1为6阶的，而filter2为12阶的，因为是IIR滤波器的设计，所以两者都不是线性相位的，这是对同一个传递函数以不同参数进行自动比拟，从而选出最优化的设计。
图2.3
通过这个设计巴特沃斯滤波器的程序，我们用realizemdl〔FT〕命令可以得出该我们所需要设计的滤波器的仿真模型，进而出现simulink界面如下图;这样一个模块可以直接用于信号传输中滤波器模块的建立。
图2.4
用鼠标双击simulink界面中的该模型，我们可以得到滤波器的设计模型如下图：
图2.5
我可以看到该模型用到了12个延迟器，13个乘法器，12个加法器，这是一个典型的优化设计滤波器模型。
.4.Sptool设计IIR数字滤波器〔面向对象设计〕
SPTool是MATLAB信号处理工具箱中自带的交互式图形用户界面工具，它包含了信号处理工具箱中的大局部函数，可以方便快捷地完成对信号、滤波器及频谱的分析、设计和浏览，因此只需要操作界面就可以载入，观察，分析，和打印数字信号，分析和设计数字滤波器。
SPTool提供对信号、滤波器和频谱分析函数的访问入口。
借助其可以：
①设计和编辑各种长度和类型并具有标准配置的FIR和IIR滤波器
②查看设计或导入的滤波器的特性，包括其幅度响应、相位响应、群延迟、极点-零点图、冲激响应和阶跃响应等
③将滤波器应用于选定的信号使用不同频谱估计方法进行频域数据的图形化分析，其中包括Burg、
FFT、多正弦窗（MTM）、MUSIC、特征向量、Welch和Yule-WalkerAR。
■Sptool设计IIR滤波器实例分析：
首先在MATLAB命令窗口输入命令：
Fs＝500；t=（0:
500）/Fs;
f=sin（2*pi*t*40）+sin（3*pi*t*50）+sin（2*pi*t*100）;
此时，变量Fs、t、s将显示在workspace列表中。
在命令窗口键入Sptool，将弹出Sptool主界面，如下图；我们按照以下步骤操作：
〔1〕点击菜单File/Import将信号f导入并取名为f。
〔2〕单击Filters列表下的New，按照参数要求设计出滤波器filt1
〔3〕将滤波器filt1应用到f信号序列。
分别在Signals、Filters、Spectra列表中选择f、filt1、mtlbseauto单击Filters列表下的Apply按钮，在弹出的ApplyFilter对话框中将输出信号命名为信号3
〔4〕进行频谱分析。
在Signals中选滤波后的信号信号3，单击Spectra下的Create按钮，在弹出的SpectraViewer界面中选择Method为FFT，Nfft=512，单击Apply按钮生成滤波后信号的频谱。
Sptool主界面
图2.6
IIRChebyshev1型低通滤波器filt1设计界面
图
模拟信号源f
图
经滤波后信号3
图
滤波后经FFT处理后频谱
图3.0
分别选中原信号f、滤波后信号3，信号3的频谱，单击各自列表下方的View按钮，即可观察他们的波形，如下图。
低通滤波器filt1使输入信号f中频率为40hz的正弦波信号通过，而将频率为75hz和100hz的正弦波信号大大衰减。
滤出的信号3集中在40HZ的频率区，说明滤波的效果比拟理想。
这样滤波后的信号3波形非常清楚的展现在用户面前。
4．滤波器设计方法总结
在对滤波器实际设计时，运用函数设计法，整个过程的运算量是很大的。
设计阶数较高的IIR滤波器时，计算量更大，设计过程中要改变参数或滤波器类型时都要重新计算。
它需要反复的实验，而且需要设计者凭借经验设定参数，平时所要设计的数字滤波器，阶数和类型并不一定是完全给定的，很多时候要根据设计要求和滤波效果不断地调整，以到达设计的最优化。
在这种情况下，滤波器设计就要进行大量复杂的运算，单纯的靠公式计算和编制简单的程序很难在短时间内完成。
因此，基于对象的信号处理工具FDATool，Fvtool以及S
升级会员 