IIR数字滤波器的设计及软件实现.docx
《IIR数字滤波器的设计及软件实现.docx》由会员分享，可在线阅读，更多相关《IIR数字滤波器的设计及软件实现.docx（17页珍藏版）》请在冰豆网上搜索。
IIR数字滤波器的设计及软件实现
IIR数字滤波器的设计及软件实现
1.实验目的
（1）熟悉用双线性变换法设计IIR数字滤波器的原理与方法；
（2）学会用MATLAB信号处理工具箱中的滤波器设计函数（或滤波器设计分析工具FDAtool）设汁各种滤波器，学会根据滤波需求确定滤波器指标参数；
（3）掌握IIR数字滤波器的MATLAB实现方法；
（4）通过观察滤波器输入、输出信号的时域波形及其频谱，建立数字滤波的概念。
2.实验原理
设计IIR数字滤波器一般采用间接法（脉冲响应不变法和双线性不变法），应用最广泛的是双线性变换法。
基本的设计过程是：
将给定的数字滤波器指标转换成模拟滤波器的指标；
涉及模拟滤波器；
将模拟滤波器的系统函数转换成数字滤波器的系统函数。
MATLAB信号处理工具箱中的各种IIR数字滤波器设计函数都是采用双线性变换法。
本实验的数字滤波器的MATLAB实验是调用MATLAB信号处理工具箱的函数filter对给定的输入信号x（n）进行滤波，得到滤波后的输出信号y（n）・
3.实验内容及步骤
1.信号处产生函数mstg产生山三路抑制载波调幅信号相加构成的复合信号st,该函数还会自动回图显示st的时域波形和幅频特性曲线，山后图可见，三路信号时域混叠无法在时域分离。
但频域是分离的，所以可通过滤波的方法在频域分离。
2.将st中三路调幅信号分离，通过观察st的幅频特性曲线，分别确定可以分离st中三路抑制载波单频调幅信号的三个滤波器（低通滤波器、带通滤波器、高通滤波器）的通带截止频率和阻带截止频率。
且滤波器的通带最大衰减为0.ldB,阻带最小衰减为60bBo
提示：
抑制载波单频调幅信号的数学表示式为
=CO
d/八6（2兀fg[co^/-/（}）+co©"/。
+/J）]
其中，称为载波，f为载波频率，cos（2^/o/）称为单频调制信号,几为调制正弦波信号频率，且满足f>f^由上式可见，所谓抑制载波单频调幅信号，就是2个正弦信号相乘，它有2个频率成分：
和频/和差频
f-f，这2个频率成分关于载波频率f对称。
所以，1路抑制载波单频调幅"fJ0Jc
信号的频谱图是关于载波频率力对称的2根谱线，其中没有载频成分，故取名为抑制载波单频调幅信号。
图中三路调幅信号的载波频率分别为250Hz、500Hz.
0/13
lOOOHzo
13
3.编程调用MATLAB滤波器涉及函数ellipord和ellip分别设讣这三个椭圆滤波器，并绘图显示其损耗函数曲线；
4.调用滤波器实验函数f订ter,用三个滤波器分别对信号产生函数mstg产生的信号st进行滤波，分离出st中的三路不同载波频率的调幅信号yl（n）、y2（n）、y3（n）。
滤波器参数的选取：
•对于载波频率为250Hz的条幅信号，可以选用低通滤波器分离，其指标为：
通带截止频率f=280Hz,通带最大衰减=O.ldZ?
；
阻带截止频率/=450Hz,阻带最小衰减5S=60dB;
•对于载波频率为500Hz的条幅信号，可以选用带通滤波器分离，其指标为：
通带截止频率九厂440Hz,fh=560Hz,通带最大衰减Q=0.1dB：
阻带截止频率f=215Hz.,f=900/7z，阻带最小衰减a=60dB;
•对于载波频率为1000Hz的条幅信号，可以选用高通滤波器分离，其指标为:
通带截止频率刀=890*,通带最大衰减d，=0.1dB；
阻带截止频率/=550^，阻带最小衰减Q=60dB;
说明：
（1）为了使滤波器阶数尽可能低，每个滤波器的边界频率选择原则是尽可能使滤波器过渡带宽尽可能宽；
（2）与信号产生函数mstg相同釆用频率Fs=10kHz；
（3）为了滤波器阶数最低，选用椭圆滤波器。
4.试验程序框图
1、阅读信号产生函数mstg,确定三路调幅信号的载波频率和调制信号频率答：
第一路调幅信号的载波频率fcl二1000Hz；
第一路调幅信号的调制频率fml=100Hz；
笫二路调幅信号的载波频率fc2=500Hz；
第二路调幅信号的调制频率fm2=50Hz；
第三路调幅信号的载波频率fc3=250Hz；
第三路调幅信号的调制频率fm3=25Hz；
2、信号产生函数mstg中采样点数\二1600,对st进行'点FFT就可以得到6根理想谱线。
如果取N二1800,可否得到6根理想谱线？
为什么？
N二2000呢？
请改变釆样点数N的值，观察频谱图验证判断是否正确？
答：
因为信号St是周期序列，谱分析时要求观察时间为整数倍周期。
分析可知，st的每个频率成分都是25Hz的整数倍。
采样频率Fs=10kHz=25X400Hz,即在25Hz的正弦波的1个周期中采样400点。
所以，当N为400的整数倍时一定为st的整数个周期。
因此，釆样点数N二1600和N二2000时，对st进行N点FFT可以得到6根理想谱线。
如果取21800,不是400的整数倍，不能得到6根理想谱线。
（1）
N二1600时：
（a）s（t）的波形
•
At
A
A
11
I1A
I
1
5
I•I
V
00.0020.0040.0060.0080.010.0120.0140.0160.0180.02
（2）N=1800时:
⑻S（t）的波形
t/s
（b）s（t）的频谱
"HZ
（3）N二2000S时:
1\
i'\
11
FVW
i/“
二S
（a）s（t）的波形
00.0050.010.0150.020.025
（b）s（t）的频谱
f/HZ
3、修改信号产生函数mstg,给每路调幅信号加入载波成分，产生调幅（AM）信号，重复本实验，观察AM信号与抑制载波调幅信号的时域波形及其频谱的差别。
AM信号表示式:
呦=h+A”co〈2/r九脸耘门）人＞Am取值：
Ad=10，£=5,结果见（附录I）
六.实验结果（程序附录II）
原信号输出：
2
u1
To
0
-1
•
\/
9
（a）s⑴的波形
3
00.0010.0020.0030.0040.0050.0060.0070.0080.0090.01
t/s
（b）s（t）的频谱
t/HZ
（b）s（t）的波形
（c）s（t）的频谱
t/Hz
带通输出:
0-20-40-60-80
（a）s（t）损耗函数曲线
0.07
0
0.20.30.40.5
co/n
（b）s（t）的波形
00.010.020.030.040.05
t/s
（c）s（t）的频谱
rw
「—
»
■
\
•厂
---—…
0.7
0.9
1
0.5
0
-0.5
高通输出:
■1
（a）s（t）损耗函数曲线
1
0.5
0-0.5
/1:
*,
\-f\
\!
v
u\
V1•
<2
-
111/
（b）s（t）的波形
0.010.020.030.040.050.06
t/s
（c）s（t）的频谱
附录（I）:
原信号输出：
（a）s（t）的波形
t/s
（b）s（t）的频谱
f/HZ
低通输出:
（a）s（t）损耗函数曲线
（b）s（t）的波形
A
■a
厂、II
AA.
P\A•
./7V
00.010.020.030.040.050.060.07
t/s
带通输出:
高通输出:
（a）s（t）损耗函数曲线
（b）s（t）的波形
A
•门A
AA.
i\：
\
;Iji
P\A•
八;
./7V
V
■
00.010.020.030.040.050.060.07
t/s
（c）s（t）的频谱
f/Hz
附录＜n）:
%主函数
%HR数字滤波器设计及软件实现
clearall;
%调用信号产生函数mstg产生又三路抑制载波调幅信号相加构成的复合信号
symsst;
symst;
st=mstg;%低通滤波器设计与实现
Fs=10000;T=l/Fs;
n=800;Tp=n:
f:
T:
k=0:
n-l;f=k/Tp;
fp=280;fs=450;
wp=2*fp/Fs;ws=2*fs/Fs;
rp=0.1;rs=60;%DF指标;（低通滤波器的通阻带边界频率）
[N,wpO]=ellipord（wp,ws,rp,rs）;%调用ellipod讣算椭圆DF阶数N和通带截止频率wp[B.A]=ellip（N.rp,rs,wpO）;%调用ellip汁算椭圆带通DF系统函数系数向量B和Aylt=filter（B,A,st）;%滤波器的软件实现
fyt=fft（yltn）;
%下而为绘图部分
figure
（2）;
subplot（3J,l）:
myplot（B,A）;
yt-y」⑴：
subplot（3J,2）;
tplot（ylt.T,yt）；
subplot（3丄3）;
stem（f.abs（fyt）/max（abs（fyt））A,）;grid：
title（,（c）s⑴的频谱）；
axis（[0.Fs/5.0,l・2]）;xlabcl（T/Hz）;ylabclC幅度J;
%带通滤波器的实现与设计
fpl=440;fpu=560;fsl=275;fsu=900;wp=[2*fpl/Fs,2*fpu/Fs];
ws=[2*fsI/Fs,2*fsii/Fs];
rp=0.1:
rs=60;
[N.wpO]=ellipord（wp,wstrp.rs）;fB.A]=ellip（N.rp.rs,wpO）;
y2t=filter（B,A»st）;
fyt=fft（y2t,n）;
figure（3）;
subplot（3,lJ）;
myplot（B,A）;
yt='y_Kt）';
subplot（3,L2）;
tplot（ylt.T,yt）；
subplot（3,l,3）;steni（fabs（fyt）/max（abs（fyt））A,）;grid;
titleC（c）s（t）的频谱）;
axis（（0.Fs/5,0J.2J）;xlabelff/H^JjylabelC幅度）;
%髙通滤波器的实现与设讣fp=890;fs=600;
wp=2*fp/Fs;ws=2*fs/Fs;
rp=0.1;rs=60;
[N.wpO]=enipord（wp,wsjp.rs）;[B,A]=elIip（N,rp,rs,wpO；highf）;y3t=filter（B.A,st）;
fyt=fft（y3t.n）;
figure（4）;
subplot（3,l,l）;
myplot（B.A）;
yt=y.i（ty；
subplot（3,l,2）;
tplot（ylt.T,yt）；
subplot（3丄3）;
stem（f.abs（fyt）/max（abs（fyt））A,）;
grid;
titleCXc）s⑴的频谱）；
axis（[0.Fs/5.0,l・2]）;
xlabeI（,f/Hz,）;ylabel（,幅度）；
clc;clear
%子程序
%产生信号程序
functionst=mstg
N=800
FS=10000;T=1/FS:
TP=N*T;
t=0:
T:
（N-l）*T;K=0:
N-l;
f=K/TP:
fcl=FS/10;%第一路调幅信号的载波频率fcl=1000HZfml=fcl/10;%第一路调幅信号的调制信号频率为fml=100hz.
fc2=FS/20:
%第二路调幅信号的载波频率fc2=500HZfm2=fc2/10:
%第二路调幅信号的调制信号频率为fm2=50h乙fc3=FS/40;%第三路调幅信号的载波频率fc3=250HZfm3=fc3/10;%第三路调幅信号的调制信号频率为fm3=25h乙xtl=cos（2*pi*fml*t）.*cos（2*pi*fcl*t）;
Xt2=cos（2*pi*fm2*t）.*cos（2*pi*fc2*t）;xt3=cos（2*pi*fm3*t）.*cos（2*pi*fc3*t）;
st=xtl+xt2+xt3:
fxt=fft（st,N）;%计算信号st的频谱.％绘图
subplot（2JJ）
plot（t,st）;
gridon;
xlabcl（'t/s）;ylabcl（'s⑴）；
axis（[0,TP/8,min（st）,max（st）]）;
titleC（a）s（t）的波形）
subplot（2,1,2）stem（fabs（fxt）/max（abs（fxt））A,）;gridon;titleC（b）s⑴的频谱）axis（[0.FS/5,0,l・2]）;
xlabcl（'f/HZ'）;ylabclC幅度'）；
%损耗输出波形
functionmyplot（B.A）
[H.W]=freqz（B,AJ000）;
m=abs（H）;
plot（W/pi,20*logl0（m/max（m）））;gridon;
xlabeI（\omegaApi,）;ylabcl（'幅度（dB）'）;axis（|OJr8O,5]）;
titleC（a）s⑴损耗函数曲线J;
%滤波器输出波形
functiontpIot（xnTyn）n=O:
length（xn）-l;
t=n*T;
plot（t,xn）;
gridon;
xlabcK't/s）;
ylabcl（'y（n）‘）；
axis（[0,t（end）5min（xn）,1.2*max（xn）]）titleC（b）s⑴的波形J;
…精心整理，希望对您有所帮助
升级会员 