数字信号处理课程设计二心电信号分析系统设计-matlab.doc

数字信号处理课程设计二心电信号分析系统设计-matlab

———————————————————————————————— 作者：
———————————————————————————————— 日期：

2

个人收集整理 勿做商业用途

数字信号处理
课程设计报告

课题名称：

心电信号分析系统的设计 –matlab

班 级: 通信103
姓 名： 滕冲
学 号： 102189
成 绩:
指导教师： 王宝珠
日 期： 2011年12月28日

目录
摘要 ………………………………………………………………………3
一、 课程设计的目的、意义…………………………………………………3
二、 设计任务及技术指标 ……………………………………………………3
三、 设计方案论证……………………………………………………………3
三、 相关函数及程……………………………………………………………3
（1） 读取txt形式的心电信号M文件………………………………………3
(2） 数据保存：分别用来保存插值前和插值后的数据………………………4
(3) 线性插值：………………………………………………………………4
（4)对以上步骤进行画图………………………………………………………4
（5) 两个滤波器………………………………………………………………5
(6）主程序……………………………………………………………………5
五、 实验结果与分析………………………………………………………10
六、 结论（只对设计课题内容进行总结并给出结论）……………………14
七、课程设计心得与体会……………………………………………………15
八、参考文献…………………………………………………………………16

摘要：设计一个心电信号分析系统。读取原始心电信号txt文件并将其前两行内容用程序的方法去除，对其做一次线性插值并保存在新的txt文件里面。由心电信号特点得知需要设计带通滤波器还有一个带阻滤波器，设置初步滤波技术指标画出初步滤波图形后看出的结果好不好，出图是否与设计指标相符，若计数指标不好应该怎样调节,对处理完的时域和频域波形进行分析,得出结论，做过程记录.该系统做完后,启用simulink功能,建立模型，进行连线，设置参数，最后启动仿真。观察示波器模块的输入和输出波形是否真确。

关键词：数字心电信号，线性插值,滤波器，功能仿真

一、 课程设计的目的、意义

通过本课题的设计学会用matlab软件实现对心电信号进行读取，设计滤波器并进行滤波,分析心电信号的特点并反馈给滤波器的设计指标,同时强化对matlab语句应用的训练更好地掌握该软件；提高对实际问题的分析解决能力，锻炼综合能力;通过实际的动手过程发现问题解决问题，加深对所学知识的理解和掌握，重在过程.
二、 设计任务及技术指标
（1）对原始数字心电信号进行读取，由数字信号数据绘制出其时域波形。
（2）对数字信号数据做一次线性插值，使其成为均匀数字信号，以便后面的信号分析。
（3）根据心电信号的频域特征（自己查阅相关资料)，设计相应的低通和高通滤波器。
(4)编程绘制实现信号处理前后的频谱，做频谱分析，得出相关结论。
（5)对系统进行综合测试，整理数据，撰写设计报告.
三、 设计方案论证
（1）采用线性插值函数interp1完成插值处理
（2）设计一个带通滤波器,经画图分析,需加上一个带阻滤波器，滤除60Hz处的干扰信号
(3）初始技术指标：
带通：wp=[1,90]＊2*pi（滤波器通带截止频率); ws=[0。25,99］＊2＊pi(阻带截止频率）； Rp=1db(通带衰减）； As=25db（阻带衰减)；
带阻：wp1=［58，62］＊2＊pi;（滤波器通带截止频率); ws1=[59.9，60.1]*2＊pi;（阻带截止频率）; Rp=1db（通带衰减）; As=25db（阻带衰减）；

四、 相关函数及程序
（1） 读取txt形式的心电信号M文件

%读取心电信号并转化成数组形式
function ［t,Xn］=duquxinhao(w)
fid=fopen(w）;
A=textscan(fid,'%8c ％f ％＊f'，'headerlines’，2); ％将其前两行内容去除
fclose（fid);
a=A{1｝;
b=A{2};
k=length(a)
for i=1：k
c（i)=strread（a(i,:），'％＊s ％f’，’delimiter'，’:');
end
c=c'
d=［c,b];
%为了减少分析数据的工作量截取2。5s(三个周期左右)的心电信号
for i=1：k
if c(i)〈=2.5
e（i,：）=d（i，：);
else break;
end
end
t=e(：，1)； %时间
Xn=e(：，2）; ％幅度
（2） 数据保存：分别用来保存插值前和插值后的数据

function shujubaocun(t,Xn)
fid = fopen（'t。txt'，'wt');
fprintf（fid,'%g\n’，t）;
fclose(fid)；
fid = fopen（'Xn。txt’，’wt')；
fprintf(fid，'％g\n'，Xn）;
fclose(fid）;
（3） 线性插值：

function [t3，Xn3]=xianxingchazhi(t，Xn）
m=max(t)；
t3=0：0.001：m；
t3=t3';
Xn3=interp1（t,Xn，t3）;
（4） 对以上步骤进行画图：

％画初始信号和即插值后信号频谱
function huizhiboxing（t，Xn，t1，Xn1)
fs=1000；
T=1/fs；
n1=length（Xn);
n=1:n1;
m1=length(Xn1）；
m=1：m1;
q=fs*n/n1；
q1=fs*m/m1;
subplot（2,2，1）
plot（t,Xn）
title('初始信号时域波形'）
subplot（2，2，2）
Y=fft（Xn）;
plot(q,abs（Y）)
axis（［0 100 0 200]）；
title（’初始信号频谱')
subplot（2,2,3)
plot(t1,Xn1)
title(’插值后信号时域波形'）
Y1=fft(Xn1）；
subplot(2,2，4）
plot（q1,abs（Y1））
title('均匀信号频谱’）
axis([0 100 0 200］);
（5） 两个滤波器：

%带通滤波器
function ［H,f］=daitonglvboqi（wp，ws,Rp,As,Xn2）
T=0。001;
f=1/T;
［N,Wc］=buttord（wp,ws,Rp，As，’s’）;
[b，a］=butter(N,Wc,'s’)；
%[b0,B,A］=dir2cas(b，a)
f=（0：length(Xn2)-1)＊f/length（Xn2）；w=f＊2＊pi；
H=freqs（b，a,w);
%plot（w/（2＊pi)，20*log10（abs（H））);
%xlabel(’频率Hz'）
％ylabel（’db')
%title([’巴特沃斯N=’,num2str（N）］）
%grid on
%带阻滤波器
function [H,f]=daizulvboqi（wp，ws,p,s，Xn2)
T=0。001;
fs=1/T;
[N，Wc］=buttord（wp,ws,p,s,'s');
[b，a］=butter（N,Wc,'stop','s'）；
f=(0：length（Xn2）-1）*fs/length（Xn2）；w=f*2*pi;
H=freqs(b，a,w）；
（6）主程序：
[t，Xn］=duquxinhao（'113.txt’)
shujubaocun（t，Xn) ％保存读取信号
［t1,Xn1]=xianxingchazhi(t，Xn);
xin=[t1，Xn1］； %用二维数组读取信号源
shujubaocun(t1，Xn1） ％保存插值后信号
figure(1)
huizhiboxing(t,Xn，t1,Xn1） ％画原始信号和插值后信号波形和频谱
wp=[5，80]*2*pi;
ws=[0。25,99］*2＊pi；

wp1=［58，62］＊2*pi；
ws1=[59.9，60.1]*2*pi；
Rp=1;
As=20；
[H1,f］=daitonglvboqi（wp,ws，Rp,As,Xn1);
[H2，f］=daizulvboqi(wp1，ws1,Rp，As,Xn1);
H=abs(H1） ；%带通
HH=abs（H2）%带阻
Y=(H.＊HH)'。＊abs（fft(Xn1）);%经过滤波后心电信号频谱
y=ifft(Y）; %滤波后心电信号时域波形
figure(2）
subplot（3，1，1)
plot(f,abs（H1）)
axis（［-20，140，0，1。5]）
title（'带通滤波器'）
subplot（3,1，2）
plot（f，abs(HH))
axis（［0，200，0，1.5]）
title('带阻滤波器’）
subplot(3,1,3)
plot（f，abs（Y))
axis([0，100，0,80])
title('滤波后心电信号频谱'）
figure(3）
subplot(2,2，1)
plot（t1，Xn1)
title（'滤波前信号'）
subplot（2，2,2）
plot（t1,y）
title('滤波后信号’）
subplot(2，2，3)
plot(t1，Xn1)
axis(［0,1.5,0,1.5］)
title('滤波前截取一部分信号')
subplot(2,2，4）
plot(t1，y）
axis(［0，1。5，0，1。5]）
title（'滤波后截取一部分信号’)
% 分析带通滤波器
figure(4)
wp=［5,80］*pi＊2 ％模拟通带角频率
ws=[0。25，99]＊pi＊2 ％模拟阻带角频率
Rp=1,As=5%通带波动、阻带衰减％巴特沃斯型模拟带通滤波器设计
［N，wc］=buttord（wp,ws，Rp,As，'s'）％计算滤波器阶数N和3db截止频率wc
［b,a]=butter(N，wc,’bandpass',’s’) ％计算模拟带通滤波器函数系数％巴特沃斯型模拟带通滤波器频率响应
［H，w］=freqs(b,a）
mag=abs（H）;
fs=1000 %采样频率
［bn，an］=impinvar(b，a,1000) %用冲击响应不变法将模拟滤波器转换成数字滤波器;数字带通滤波器系数％巴特沃斯型数字带通滤波器频率响应
ht=impz(bn,an) %数字滤波器冲击响应曲线
subplot（3，1，1)
plot（w/(2*pi），mag）％模拟带通滤波器的幅频响应
axis([—50,150，0,1。5］）
title（’模拟带通滤波器幅频曲线’）
grid ％模拟幅频曲线
subplot(3,1，2);plot（ht)；
title（'数字带通滤波器冲击响应曲线’)
axis（[0,100,—0.5，0。5]);
grid
subplot（3,1,3)；zplane（bn，an）;
axis（［—2,2,—1.5，1。5])
title（’数字带通滤波器的零极点图’）％zplane函数画出零极点图
grid
％ 分析带阻滤波器
figure（5）
wp1=[58,62］*2*pi；%模拟通带角频率
ws1=［59.9,60。1]＊2＊pi； %模拟阻带角频率
Rp=1，As=20%通带波动、阻带衰减%巴特沃斯型模拟带通滤波器设计
［N,wc］=buttord(wp1，ws1,Rp,As，'s')%计算滤波器阶数N和3db截止频率wc
［b，a］=butter(N，wc，’stop’，'s') ％计算模拟带通滤波器函数系数％巴特沃斯型模拟带通滤波器频率响应
［H，w］=freqs（b，a)
mag=abs(H）；
fs=1000 ％采样频率
[bn,an］=impinvar(b,a，1000） %用冲击响应不变法将模拟滤波器转换成数字滤波器；数字带阻滤波器系数
ht=impz（bn，an) %数字滤波器冲击响应曲线
subplot(3，1，1)
plot(w/（2＊pi)，mag）%模拟带阻滤波器的幅频响应
axis（[0，150，0,1.5]）
title('模拟带阻滤波器幅频曲线’）
grid ％模拟幅频曲线
subplot（3，1,2)；plot(ht）；
title（'数字带阻滤波器冲击响应曲线')
axis（［0,600,-0.05,0。05］）;
grid
subplot(3，1，3)；zplane(bn，an）;
title（’数字带阻滤波器的零极点图’)%zplane函数画出零极点图
grid

五、 实验结果与分析

分析：由图可以看出,经过线性插值后的时域信号和均匀信号频谱图，在均匀信号频谱中可以看出在以60Hz为中心的窄带处存电源线干扰，所以在后面添加限制60Hz频率通过的带阻滤波器。

分析：经过带通滤波器后再通过带阻滤波器将以60Hz为中心的窄带处存电源线干扰滤除后的心电信号频谱图。

分析：只截取大约2。5s，三个周期左右的时域信号，观察并对比滤波效果
带通滤波器指标为：wp=[0.5,90]*2*pi；ws=［0。25,99］*2*pi； Rp=1；As=20；
由图所示主频附近还有很高的干扰信号，开始以为是阻带衰减取值过大，以下是减小阻带衰减后的滤波
分析：带通滤波器指标为:wp=［0。5，90］*2＊pi;ws=[0.25，99]*2*pi； Rp=1；As=5；
由图所示主频附近仍然有很高的干扰信号，发现改变阻带衰减并无影响，觉得是滤波器通带过宽所致,以下是把通带变窄后的滤波

带通滤波器指标为：wp=[20,90］＊2＊pi；ws=［0。25，99]*2＊pi; Rp=1;As=25；
由图所示主频附近还有很高的干扰信号，应该是滤波器通带过窄以下是把通带稍加宽后的滤波

带通滤波器指标为：wp=［5，90］＊2*pi；ws=［0。25,99]*2*pi； Rp=1；As=25;
由图所示:在主频附近已经基本没有了干扰信号，基本满足要求。
5阶巴特沃斯带通滤波器的时域分析：

1阶巴特沃斯带通滤波器的时域分析：

六、 结论（只对设计课题内容进行总结并给出结论）
总结：课题设计最重要的部分是滤波器的设计，起承上启下的作用，同时也是最难和最费时间的内容。Matlab系统本身自带了多种滤波器供选择，通过调用和设置可以达到目的。

结论：所设计的滤波器基本满足滤波要求，所得波形较为合理。

七、课程设计心得与体会

历时两周的课程设计最终结束，在这段时间里我们一起朝着一个共同的目标努力.在这个过程里我付出了努力，取得了收获，学到了知识，同时也得到了深刻的教训。课程设计是一个完整的系统工程,需要做大量和充分的前期准备和工作……
课程设计后我最大的收获是学会了大量未知的matlab相关函数和语句语法，比如读取和保存txt文件，对非均匀信号用线性插值函数处理得到新的均匀信号;还有滤波器设计方面的有关方法；学会了采用simulin进行动态建模和仿真的步骤和过程等等.在集体的学习氛围下，通过和同学们的交流和相互学习，体会到了团队的重要性，它可以使大家相互发现自己的不足和优点，大大降低了卡壳的几率,防止因为一个小的羁绊而无法前进。同时也可以增加自己对胜利的信心。课程设计是一个理论联系实践的机会，通过这个学习的过程，可以将书本上所学的知识通过这个过程进行检验，从而发现知识的盲点并及时弥补，通过实践更好地将知识进行整合。回归到课题本身，在参考了同学的思路后，我翻了很多有关滤波器设计的书上例题还有实验的实现方法，并对自己的程序进行了仔细的检查，虽然最终的结果不是太过理想,不过还是有了很大的体会，磨刀不误砍柴工，工欲善其事必先利其器很好的表达了我对这段时期的感想。
既然是数字信号处理的课程设计，安排在课程学完之后，那么就应该掌握本书所要求的知识技能，涉及到的是有关滤波器的设计。在上第一节课时，有关滤波器的设计我还是掌握的不太熟练，于是我用课下时间抓紧复习了书上的五六七章，并通过习题来强化.有些东西虽然看起来很难,但其实思路很清晰，只要踏实做下去就会很好的完成它;有些东西虽然看起来很简单，但如果没有扎实的理论基础就会陷入迷茫之中；学而不思则罔,思而不学则殆.
通信专业的课程都是环环相扣的，一门课所学总会在以后的课程中用到，如果在前面的学习中没有好好付出努力就会在将来的课程中付出代价。所以说对于通信的孩子更不能放松学习，需要时刻绷紧这根弦。考试有重点，但是对于知识是不存在所谓的重点的,所以在以后的学习中要端正这个态度是非常重要的事情。此外，进行这场课程设计里一个目的就是学习新的东西,比如说我们组的simulink工具使用，还有兄弟组里的labview软件的使用，这些东西先前都没有接触过，不过需要我们花点时间先学会和掌握，然后很好地和课题融为一体,作为课设的一部分总有一些创新的东西在里面.这样做可以锻炼我们自主学习的主动性和应用能力的提高。所以说花费了这么长的时间进行它关键在于对综合创新能力的锻炼和提高。
大二第一个学期学习了matlab这门课程,只有简单的几节课，而matlab这门课的教科书却是那么的厚，更不用说它能实现的强大功能，而通过十来节课的学习就要求对这门将来会经常用到的计数全面掌握是非常可笑的事情。当然这是学校安排的问题，无法去具体评价,我只是想说，在后来的学科中我们一次又一次地与它见面了.类似的问题还有就是示波器，从一开始就在使用和学习的实验仪器，如果大一没有掌握它那么就会发现在做实验就会困难重重问题不断;有些学科比如说低频电子线路是高频的基础,而高频又是下学期通信原理的理论基础,都说通信难其实是基础没到位的缘故。所以作为通信专业的基本工具和要求是必须要掌握的，不能是简单的考试合格就行了的,而是要像小学时的乘法口诀那样做梦都能倒背如流的牢牢掌握.鉴于大学的情况，考完试就忘了所学的现实，最主要的还是我们学生自己的问题,如果每个学生都能主动地学习，时时刻刻都不放松自己,那么什么样的情况都不会对我们造成影响。作为学校的一方,安排了课程设计这门课就应该从实际出发，一方面完成既定的教学大纲，另一方面也应该贴近现实着手准备能够起到强烈督促学生自主学习的目的,而不应该是采用简单的布置任务后就基本不管的压力政策，在现代化的企业中是对的,而办学的首要目的是教书育人，最起码在学校这段时期里能有所提高。个人收集整理,勿做商业用途文档为个人收集整理，来源于网络

八、参考文献

（1） 张志涌。MATLAB教程（R2010a）。北京航空航天大学出版社，2010.6
(2） 丁玉美。数字信号处理（第三版）.西安电子科技大学出版社，2008
(3） 吴大正. 信号与线性系统分析(第四版）. 高等教育出版社，2005,8
(4） 谢嘉奎。 电子线路--线性部分(第四版）. 高等教育出版社，2007，12
（5） 陈后金. 信号分析与处理实验。 高等教育出版社，2006，8

第 16页 共 18页