北京信息科技大学
数字信号处理实验报告
题 目：
数字信号处理课程设计实验
学 院： 信息与通信工程学院 专 业： 通信工程专业 姓 名： 班 级：
学 号： 指导老师：
实验目的
1、熟悉IIR数字滤波器的设计原理与方法。
2、掌握数字滤波器的计算机软件实现方法。
3、通过观察对实际心电图信号的滤波作用，学习数字滤波器在实际中的应用。
实验仪器及材料
计算机，MATLAB软件
实验内容及要求
1．设计巴特沃斯低通数字滤波器对人体心电信号进行滤波
(1) 人体心电图信号在测量过程中会受到工业高频干扰，所以必须经过低通滤波处理，才能作为判断心脏功能的有用信息。以下为一个实际心电图信号采样序列x(n)，其中存在高频干扰，抽样周期Ts=1秒。在实验中，以x(n)作为输入序列，滤除其中干扰成分。
x(n)=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0] 对序列x(n)用FFT做频谱分析，生成x(n)的频谱图。
(2)设计一个巴特沃斯低通IIR数字滤波器H(z)。
设计指标参数为：在通带内频率低于0.2π时，最大衰减小于1dB； 在阻带内 [0.3π, π]频率区间上，最小衰减大于15dB。
j|H(e)|。 写出数字滤波器H(z)的表达式，画出滤波器的幅频响应曲线
(3) 用所设计的滤波器对实际心电图信号采样序列x(n)进行滤波处理，编写程序，求滤波后的序列y(n)，并分别画出滤波前后的心电图信号波形图和频谱图。
y(n)= [0,0,0,0, 0,0,0,0,-0.14025,0.40279,-0.56085 ,0.33328,0.023981,-0.18809,0.11843,
-0.1038,0.11576,-0.1225,0.099815 ,-0.13769 ,0.095249,-0.0070273,0.018867,0.090543,-0.11257,-0.070884 ,0.17676,-0.55407,0.24813,-0.34732,-0.30428,0.59426,-0.29574,-0.063869,0.34018,-0.73334,1.0293,-0.57107,-0.2461,0.83605,-0.83026,0.45459,0.011551,-0.25667,0.23896,-0.17361,0.20829,-0.28417,0.28765 ,-0.2035,0.02865,0.066164,0.077916,-0.36052, 0.53517,-0.5571]
源程序
clear all，clc
x=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];%未经滤波的心电图信号 L=length(x); l=0:L-1;
y=fft(x,L); Wp=0.2*pi; Ws=0.3*pi; Rp=1; Rs=15;
[N,Wn] = buttord(Wp,Ws,Rp,Rs,\'s\'); [b,a] = butter(N,Wn,\'s\'); [numa,dena]=impinvar(b,a,1); w=linspace(0,pi,1024); h=freqz(numa,dena,w); norm=max(abs(h)); numa=numa/norm; [z,p]=tf2zp(b,a); figure(1)
plot(w,20*log10(abs(h)/norm)); grid;
xlabel(\'数字频率\'); ylabel(\'幅度响应dB\'); figure(2) plot(w,abs(h)); grid;
xlabel(\'数字频率\');
ylabel(\'幅度响应|H(e^(jw))|\'); figure(3) zplane(z,p); xx=filter(b,a,x); yy=fft(xx,L); figure(4) subplot(2,1,1) stem(l,x);
title(\'未经滤波的心电图信号\'); xlabel(\'n\'); subplot(2,1,2) stem(l,xx);
title(\'经滤波之后的心电图信号\'); xlabel(\'n\'); figure(5) subplot(2,1,1) plot(l,abs(y));
title(\'未经滤波的心电图信号的频谱\'); subplot(2,1,2) plot(l,abs(yy));
title(\'经滤波处理的心电图信号的频谱\');
2．用help查看内部函数cheb1ord.m及cheby1.m，了解调用格式。
编程设计教材习题6-2，求模拟滤波器Ha(s)的表达式。
源程序
close all clear all clc
Wp=2*pi*3000; Rp=2; Ws=2*pi*12000; Rs=50; Fs=24000;
w=linspace(0,pi,1024);
[N,Wn]=cheb1ord(Wp,Ws,Rp,Rs,\'s\'); e=sqrt(10^(Rp/10)-1); [b,a]=cheby1(N,e,Wn,\'s\') [numa,dena]=impinvar(b,a,Fs); h=freqz(numa,dena,w); norm=max(abs(h));
plot(w*Fs/pi,20*log10(abs(h)/norm)) title(\'幅度响应\') xlabel(\'频率(Hz)\') ylabel(\'幅度(dB)\') grid
3．模拟滤波器的数字化
用内部函数impinvar及bilinear实现教材习题6-5，求数字滤波器H(z)的表达式。
源程序
close all clear all clc b1=[0 0 1]; a1=[1 1 1]; b2=[0 0 1]; a2=[2 3 1];
[numa1,dena1]=impinvar(b1,a1,0.5) [numa2,dena2]=bilinear(b1,a1,0.5) [numa3,dena3]=impinvar(b2,a2,0.5)
[numa4,dena4]=bilinear(b2,a2,0.5)
本实验所用的部分MATLAB函数
L=length(x) ：求序列x长度。
y=fft(x,L) ：将序列x(n)做L点快速傅立叶变换，结果赋给序列y(n)。 [n,Wn] = buttord(Wp,Ws,Rp,Rs,\'s\') ：计算模拟Butterworth滤波器的最小阶次n和截止频率为Wn。
[b,a] = butter(n,Wn,\'s\')：设计模拟截止频率为Wn（rad/s）的n阶 Butterworth低通滤波器，返回值为模拟滤波器的系数。
y=filter(b,a,x): 将序列x(n) 通过滤波器滤波后生成序列y(n)，滤波器的分母多项式系数构成a向量，分子多项式系数构成b向量。
[BZ,AZ] = impinvar(B,A,Fs)：冲激响应不变法，返回值为数字滤波器的系数。 [BZ,AZ] = bilinear(B,A,fs) ：双线性变换，返回值为数字滤波器的系数。 [H w]=freqz(b,a)：由滤波器分母多项式系数构成的a向量和分子多项式系数构成的b向量求系统频响。
截图
实验体会
这次的实验让这让我看到了理论与实践相结合的优势与用处，让我受益匪浅。我认识到了自己理论知识的不足，也认识到了我们学习的基础知识究竟能运用于什么领域，如何运用。我们在老师的耐心指导下调试电路，直到得到要求的效果。让我们在学习电路、信号等理论知识的同时，明白如何把这些应用于实际。
