数字信号处理实验报告(实验三).doc

实验三 用双线性变换法设计IIR数字滤波器
1. 实验目的
(1) 熟悉用双线性变换法设计IIR数字滤波器的原理与方法。
(2) 掌握数字滤波器的计算机仿真方法。
(3) 通过观察对实际心电图信号的滤波作用， 获得数字滤波的感性知识。
2. 实验内容
(1) 用双线性变换法设计一个巴特沃斯低通IIR数字滤波器。设计指标参数为：在通带内频率低于0.2π时，最大衰减小于1dB；在阻带内[0.3π, π] 频率区间上，最小衰减大于15dB。
(2) 以 0.02π为采样间隔， 打印出数字滤波器在频率区间［0, π/2］上的幅频响应特性曲线。
(3) 用所设计的滤波器对实际心电图信号采样序列(在本实验后面给出)进行仿真滤波处理，并分别打印出滤波前后的心电图信号波形图， 观察总结滤波作用与效果。
3.实验原理
为了克服用脉冲响应不变法产生频谱混叠现象，可以采用非线性频率压缩方法（正切变换），从s平面映射到s1平面，再从s1平面映射到z平面,即实现了双线性变换。
4. 实验步骤
(1) 复习有关巴特沃斯模拟滤波器设计和用双线性变换法设计IIR数字滤波器的内容， 按照例 6.4.2， 用双线性变换法设计数字滤波器系统函数H(z)。
例 6.4.2 中已求出满足本实验要求的数字滤波器系统函数：

(2－1)

(2－2)

A=0.09036
B1=1.2686, C1=－0.7051
B2=1.0106, C2=－0.3583
B3=0.9044, C3=－0.2155
由(2－1)式和(2－2)式可见，滤波器H(z)由三个二阶滤波器H1(z)，H2(z)和H3(z)级联组成，如图 2－1 所示。

y3(n)=y(n)

图 2－1 滤波器H(z)的组成

(2) 编写滤波器仿真程序，计算H(z)对心电图信号采样序列x(n)的响应序列y(n)。设 为第k级二阶滤波器 的输出序列， 为输入序列，如图2-1 所示。 由(2-2)式可得到差分方程：

(2-3)
当k=1时， 所以 的总响应序列
可以用顺序迭代算法得到。即依次对k=1,2,3, 求解差分方程(2-3),
最后得到 仿真程序就是实现上述求解差分方程和顺序迭代算法的通用程序。也可以直接调用MATLAB filter函数实现仿真。
(3) 在通用计算机上运行仿真滤波程序，并调用通用绘图子程序， 完成实验内容(2)和(3)。

5. 思考题
用双线性变换法设计数字滤波器过程中， 变换公式 中T的取值， 对设计结果有无影响? 为什么?
答：的幅度与采样间隔成反比，当T很小时， 就会有太高的增益，会有一定的影响。
6. 实验报告要求
(1) 简述实验目的及原理。
(2) 由所打印的 特性曲线及设计过程简述双线性变换法的特点。
(3) 对比滤波前后的心电图信号波形， 说明数字滤波器的滤波过程与滤波作用。
(4) 简要回答思考题。
7. 心电图信号采样序列x(n)
人体心电图信号在测量过程中往往受到工业高频干扰， 所以必须经过低通滤波处理后， 才能作为判断心脏功能的有用信息。 下面给出一实际心电图信号采样序列样本x(n)， 其中存在高频干扰。 在实验中， 以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}
8. 实验结果
***以下为实验内容***
(1). 仿真程序：
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]; %向数组x中输入数据
k=1;
close all; %关闭所有打开的图形窗口
figure(1); %图形窗口命名
subplot(2,2,1); %对产生第一个图形划分坐标轴(确定方位)

n=0:55; % n从0开始以默认步长为1取到55
stem(n,x,'.'); %确定离散的序列图的参数
axis([0 55 -100 50]); %控制坐标值的大小
hold on; %保持当前的图形
n=0:60; % n从0开始以默认步长为1取到60（即滤波器可接受心电图信号采样序列x(n)的最大范围）
m=zeros(61); %m为61行61列的零矩阵
plot(n,m); %确定离散的序列图的参数
xlabel('n');ylabel('x(n)'); %将y轴命名为x(n), x轴命名为n
title('心电图信号采样序列x(n)'); %将该图形命名为“心电图信号采样序列x(n)”
B=[0.09036 2*0.09036 0.09036]; %向数组b中输入数据
A=[1.2686 -0.7051]; %向数组a中输入数据
A1=[1.0106 -0.3583]; %向数组a1中输入数据
A2=[0.9044 -0.2155]; %向数组a2中输入数据

第一种“循环”编程法：
while(k<=3)
y=filter(B,A,x);
x=y;
if k==2
A=A1;
end
if k==3
A=A2;
end
k=k+1;
end

第二种“顺序”编程法：
Y1=filter(B,A,x);
Y2=filter(B,A1,y1);
y=filter(B,A2,y2);

subplot(2,2,3); %对产生第三个图形划分坐标轴(确定方位)
n=0:55; % n从0开始以默认步长为1取到55
stem(n,y,'.'); %确定离散的序列图的参数
axis([0 55 -15 5]); %控制坐标值的大小
hold on; %保持当前的图形
n=0:60; % n从0开始以默认步长为1取到60（即滤波器可接受心电图信号采样序列x(n)的最大范围）
m=zeros(61); %m为61行61列的零矩阵
plot(n,m); %确定离散的序列图的参数
xlabel('n');ylabel('y(n)'); %将y轴命名为y(n), x轴命名为n
title('三级滤波后的心电图信号'); %将该图形命名为“三级滤波后的心电图信号”
A=[0.09036, 0.1872,0.09036]; %向数组a中输入数据
B1=[1,-1.2686,0.7051]; %向数组b1中输入数据
B2=[1,-1.0106,0.3583]; %向数组b2中输入数据
B3=[1,-0.9044,0.2155]; %向数组b3中输入数据

[H1,w]=freqz(A,B1,100); %数字滤波器（h1）100点幅频响应
[H2,w]=freqz(A,B2,100); %数字滤波器（h2）100点幅频响应
[H3,w]=freqz(A,B3,100); %数字滤波器（h3）100点幅频响应

H4=H1.*(H2);
H=H4.*(H3); %滤波器H(z)由三个二阶滤波器H1(z)，H2(z)和H3(z)级联组成

mag=abs(H); %求系统函数的绝对值
db=20*log10((mag+eps)/max(mag)); %转化为分贝数
subplot(2,2,2)； %对产生第二个图形划分坐标轴(确定方位)
plot(w/pi,db); %确定离散的序列图的参数
axis([0 0.5 -50 10]); %控制坐标值的大小
title('滤波器的幅频响应曲线'); %将该图形命名为“滤波器的幅频响应曲线”
(2). 经过matlab软件调试后的结果：

对实验内容的总结：滤波器H(z)由三个二阶滤波器H1(z)，H2(z)和H3(z)级联组成，且从其幅频响应曲线可以看出，其实质为一个低通滤波器,横轴为归一化后的频率，纵轴为转化后的分贝数。