发布时间 : 星期二 文章数字信号处理上机指导 - 1dl更新完毕开始阅读a620bd4b852458fb770b56fc
magxk=abs(xk)/a; %归一化处理 n1=2/N*n0; subplot(4,1,1); stem(n1,magxk);
%信号序列的插值序列yn, 插值因子I, yn=chazhi(xn,I); M=length(yn);
yk=dft(yn,M); %序列插值后的频谱yk b=max(abs(yk));
magyk=abs(yk)/b; %归一化处理 m=[0:1:M-1]; m1=2/M*m; subplot(4,1,2) stem(m1,magyk);
%生成滤除高频分量的滤波器hk np=wp*M/2;
hk=[1,zeros(1,M-1)]; for i=2:1:np hk(i)=1;
hk(M-i+2)=1; end
ykk=yk.*hk; %经过滤波后的频谱ykk magykk=abs(ykk)/b; subplot(4,1,3)
stem(m1,magykk);
%计算经过滤波后的序列ynn
ykkk=conj(ykk); %取复共轭函数:conj ynn=conj(dft(ykkk,M))/M;
%ynn序列的抽取序列zn，抽取因子D zn=chouqu (ynn,D); M=length(zn);
zk=dft(zn,M); %序列抽取后的频谱zk c=max(abs(zk));
magzk=abs(zk)/c; %归一化处理 m=[0:1:M-1]; m1=2/M*m; subplot(4,1,4) stem(m1,magzk);
function yn=chazhi(xn,I) N=length(xn); M=I*N;
yn=zeros(1,M); for n=1:1:N
m=(n-1)*I+1; yn(m)=xn(n); end
function yn=chouqu(xn,D) N=length(xn); m=0;
for n=1:D:N m=m+1;
yn(m)=xn(n); end
实验四 IIR滤波器的实现与应用
- 28 -
一、实验目的
1．掌握用模拟滤波器设计IIR数字滤波器的方法。 2．应用IIR数字滤波器实现对信号的处理。 二、实验内容
设信号x(t)=cos(2π×50t)+2×cos(2π×400t)，试将它的两个频率分量分离，并绘制它们的时域波形及频谱图。 三、实验原理
（一）用模拟滤波器设计IIR数字滤波器
由模拟滤波器设计数字滤波器是先设计一个合适的模拟滤波器，然后将其数字化，即将S平面映射到Z平面，得到所需要的数字滤波器。由于模拟滤波器的设计方法已很成熟，故此法设计数字滤波器准确、简便，是目前最普遍采用的方法。 实现上述设计思想有多种途径，其中之一的设计步骤是：
1．确定所要设计的数字滤波器的性能指标及各数字特征频率{ωk}；
2．由冲激响应不变法或双线性变换法的变换关系将{ωk}变换为模拟特征频率{Ωk}；
冲激响应不变法：Ωk=ωk/2π
双线性变换法：Ωk=c×tan（ωk/2）（预畸）
3．把要求的Lp、Hp、Bp或Bs特征频率参数转化为模拟低通滤波器的设计参数{ΩLk }；
Hp→Lp：Ωp=Ωph；Ωst=Ωph/Ωsth Bp→Lp：Ωp=Ωp2－Ωp1
Ωst=min{|(Ωs2－Ω0)/Ωs2|，|(Ωs1－Ω0)/Ωs1|} 2
Ω0=Ωp2×Ωp1
Bs→Lp：Ωp=Ωp2×Ωp1/（Ωp2－Ωp1）
Ωst=min{|Ω0×Ωs2/(Ω0－Ωs2) |，|Ω0×Ωs1/(Ω0－Ωs1) |} Ω02=Ωp2×Ωp1
4．按指标{ΩLk }设计一个模拟低通滤波器原型，设其系统函数为Ha(s)； 5．将模拟低通滤波器原型Ha(s)转化为所需类型的模拟滤波器Ha(p)；
Lp→Hp：s=Ωph2/p
Lp→Bp：s=(P2－Ω02)/p Lp→Bs：s= pΩ02/(Ω02－p2)
6．将模拟滤波器Ha(p)转化为同类型的数字滤波器H(z)。
N2
2
2
2
2
2
2
2
2
2
2
冲激不变法： Ha(p)??k?1Akp?pk
N H(z)??1?ek?1TAkpkTz
?1
双线性变换法： H(z)?Hq(p)通常取c=2/T；T为采样周期。
p?c1?z1?z?1?1（二）应用IIR数字滤波器实现对信号的处理
对信号的处理包括滤波、变换、识别、压缩等，此处主要是应用IIR数字滤波器实现对信号的分离，即利用IIR数字滤波器选出有用的信号。 四、实验步骤
- 29 -
1．分析实验内容，设计实验方案
例如：
2．编程实现每一个子功能并调试之
例如： 设定采样频率fS，将连续信号x(t)离散为长度为N的序列x(n)。
设：fS=1kHz，
则：T=1/ fS=0.001秒
x(n)= x(t)| t=nT
=cos(0.1πn)+cos(0.8πn)
由于1/50与1/400的最小公倍数为1/50=0.02，因而x(t)的周期T0为0.02秒，序列的长度N=0.02/0.001=20。 其Matlab程序为[xn,wk,N]=shiyan40(fs,T0)
3．选择合理的滤波器参数，编程完成信号的分离。 五、实验要求
1．针对实验内容，首先进行理论分析，给出预测结果； 2．熟悉实验教材§6-4的内容；
3．编写每一个子功能程序，上机运行；
4．研究滤波器各参数值、采样频率、记录长度的大小对实验结果的影响。 六、程序示例
function [yl,yh]=shiyan49
fs=1600; %采样频率 Tt=0.02; %信号周期 T0=4*Tt; %记录长度 [xn,wk,N]=shiyan40(fs,T0); M=length(wk); if M==2
rp=1;rs=80;
f1=wk(1)*fs/N; f2=wk(2)*fs/N; f0=f2-f1;
fpl=f1+f0/20; fstl=f2-f0/10; %模拟低通滤波器的特征参数 [bzl,azl]=shiyan42(fpl,fstl,rp,rs,fs);
fph=f2-f0/2;fsth=f1+f0/20; %模拟高通滤波器的特征参数 [bzh,azh]=shiyan43(fph,fsth,rp,rs,fs); end
ynl=filter(bzl,azl,xn); %序列xn通过数字低通滤波器，输出为ynl ynh=filter(bzh,azh,xn); %序列xn通过数字高通滤波器，输出为ynh knl=abs(fft(ynl)); % ynl的频谱knl
kl=knl/max(knl); % ynl的幅度归一化频谱kl knh=abs(fft(ynh)); kh=knh/max(knh); T=1/fs; figure(8)
t=0:T:(T0-T);
w=0:2*pi/N:(2*pi-2*pi/N); subplot(2,2,1); plot(t,ynl); subplot(2,2,2); stem(w,kl);
- 30 -
subplot(2,2,3); plot(t,ynh); subplot(2,2,4); stem(w,kh);
%去掉滤波输出序列的头一个周期 yl=ynl(N/4+1:N); yh=ynh(N/4+1:N);
knl=abs(fft(yl)); kl=knl/max(knl); knh=abs(fft(yh)); kh=knh/max(knh); N=length(kl); t=Tt:T:(T0-T);
w=0:2*pi/N:(2*pi-2*pi/N); figure(9)
subplot(2,2,1) plot(t,yl);
Xlabel('秒'); title('低频分量时域波形') subplot(2,2,2); stem(w,kl)
Xlabel('数字频率/弧度 ');title('低频分量归一化频谱图') subplot(2,2,3) plot(t,yh);
Xlabel('秒'); title('高频分量时域波形') subplot(2,2,4) stem(w,kh)
Xlabel('数字频率/弧度 '); title('高频分量归一化频谱图')
function [xn,wk,N]=shiyan40(fs,T0) T=1/fs;
t=0:T:(T0-T);
xn=cos(100*pi*t)+2*cos(800*pi*t); xk=fft(xn); N=length(xn); i=1; wk=0;
for m=1:1:(N+1)/2
if (abs(xk(m))>0.0001) wk(i)=m-1; i=i+1; end end
n=0:1:N-1; figure(1)
subplot(2,1,1) plot(t,xn); subplot(2,1,2)
stem(n,abs(xk),'r')
function [Ba,Aa,Wn]=shiyan41(mp,ms,rp,rs) %生成归一化频率的模拟低通滤波器
[N,Wn]=buttord(mp,ms,rp,rs,'s');%mp/ms:通带/阻带截止频率(弧度%/秒) [z,p,k]=buttap(N); [Ba,Aa]=zp2tf(z,p,k);
function [bz,az]=shiyan42(fp,fst,rp,rs,fs) %数字低通滤波器的生成
W0=[0,fp,fst,fs/2]; rr=[rp,rp,rs,rs];%设计指标 mp=fp*2*pi;
ms=fst*2*pi; % mp/ms:通带/阻带截止频率（弧度/秒）
- 31 -