基于MATLAB的心电信号分析心电信号分析(自己做的)带程序 带图片

河北工业大学信息工程学院信号与线性系统课程设计1基于MATLAB的心电信号分析摘要：本课题设计了一个简单的心电信号分析系统。直接采用Matlab语言编程的静态仿真方式、采用Simulink进行动态建模和仿真的方式，对输入的原始心电信号，进行线性插值处理，并通过matlab语言编程设计对其进行时域和频域的波形频谱分析，根据具体设计要求完成系统的程序编写、调试及功能测试。得出一定的结论。关键字：matlab、心电信号提取、线性插值、滤波、simulink仿真。一、课题目的及意义心电信号是人类最早研究并应用于医学临床的生物信号之一，它比其它生物电信号更易于检测，并且具有较直观的规律性，因而心电图分析技术促进了医学的发展。然而，心电图自动诊断还未广泛应用于临床，从国内外的心电图机检测分析来看，自动分析精度还达不到可以替代医生的水平，仅可以为临床医生提供辅助信息。其主要原因是心电波形的识别不准，并且心电图诊断标准不统一。因此，探索新的方法以提高波形识别的准确率，寻找适合计算机实现又具诊断价值的诊断标准，是改进心电图自动诊断效果，扩大其应用范围的根本途径。如何把心电信号的特征更加精确的提取出来进行自动分析，判断出其异常的类型成了亟待解决的焦点问题。本课题通过matlab语言编程，对原始心电信号进行一定的分析处理。二、课题任务及要求1、必做部分（1）利用Matlab对MIT-BIH数据库提供的数字心电信号进行读取，并还原实际波形。（2）对原始心电信号做线性插值（3）对处理前后的心电信号分别做频谱分析利用Matlab软件对处理前后的心电信号编程显示其频谱，分析比较滤波前后的频谱，得出结论。（4）Simulink仿真根据前面的设计，进行基于Simulink的动态仿真设计。实现心电信号的分析和处理。2、选作部分河北工业大学信息工程学院信号与线性系统课程设计2（1）只截取大约2.5s，三个周期左右，大约800个采样数据进行分析。（2）60Hz工频陷波器设计三、设计技术指标四、设计方案论证1、必做部分2、选作部分河北工业大学信息工程学院信号与线性系统课程设计3五、设计内容及结果分析1、基于matlab编写的程序如下：%读取心电信号并转化成数组形式function[t,Xn]=duquexinhao1(w)fid=fopen(w);C=textscan(fid,'%8c%f%*f','headerlines',2);%去除前两行fclose(fid);a=C{2};b=C{1};k=length(b);fori=1:kc(i)=strread(b(i,:),'%*s%f','delimiter',':');endc=c';d=[c,a];t=d(:,1);%时间Xn=d(:,2);%幅度%线性插值function[t3,Xn3]=xianxingchazhi(t,Xn)m=max(t);t3=0:0.001:m;t3=t3';Xn3=interp1(t,Xn,t3);%保存插值前的信号functionbaocun1(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);%保存插值后的信号functionbaocun2(t1,Xn1)fid=fopen('t1.txt','wt');fprintf(fid,'%g\n',t1);fclose(fid);fid=fopen('Xn1.txt','wt');fprintf(fid,'%g\n',Xn1);fclose(fid);%画初始信号和即插值后信号频谱河北工业大学信息工程学院信号与线性系统课程设计4functionkeshehuatu(t,Xn,t1,Xn1)f=1000;T=1/f;m=1:length(Xn);k1=length(Xn1);m1=1:k1;q=f*m/length(Xn);q1=f*m1/k1;subplot(2,2,1)plot(t,Xn)title('初始信号时域波形')subplot(2,2,2)Y=fft(Xn);plot(q,abs(Y))title('初始信号频谱')subplot(2,2,3)axis([0,1000,0,1000])plot(t1,Xn1)title('插值信号时域波形')Y1=fft(Xn1);subplot(2,2,4)axis([0,1000,0,5000])plot(q1,abs(Y1))title('插值信号频谱')%低通滤波器function[H,f]=kesheditonglvboqi(wp,ws,Rp,As,Xn1)T=0.001;f=1/T;[N,Wc]=buttord(wp,ws,Rp,As,'s');[b,a]=butter(N,Wc,'s');f=(0:length(Xn1)-1)*f/length(Xn1);w=f*2*pi;H=freqs(b,a,w);%高通滤波器function[H,f]=keshegaotonglvboqi(wp,ws,Rp,As,Xn1)T=0.001;fs=1/T;[N,Wc]=buttord(wp,ws,Rp,As,'s');[b,a]=butter(N,Wc,'high','s');f=(0:length(Xn1)-1)*fs/length(Xn1);w=f*2*pi;H=freqs(b,a,w);%带阻滤波器function[H,f]=keshedaizulvboqi(wp,ws,p,s,Xn1)T=0.001;河北工业大学信息工程学院信号与线性系统课程设计5f=1/T;[N,Wc]=buttord(wp,ws,p,s,'s');[b,a]=butter(N,Wc,'stop','s');f=(0:length(Xn1)-1)*f/length(Xn1);w=f*2*pi;H=freqs(b,a,w);主函数如下（1）、将信号通过低通、高通、带阻滤波器程序[t,Xn]=duquexinhao1('117.txt');baocun1(t,Xn)%保存读取信号[t1,Xn1]=xianxingchazhi(t,Xn);baocun2(t1,Xn1)%保存插值后信号xy=[t1,Xn1];%仿真输入二维数组figure(1)keshehuatu(t,Xn,t1,Xn1)%画原始信号和插值后信号波形和频谱wp=90*2*pi;%低通滤波器滤波ws=99*2*pi;p=1;s=35;[H1,f]=kesheditonglvboqi(wp,ws,p,s,Xn1);wp=4*2*pi;%高通滤波器滤波ws=0.25*2*pi;p=1;s=35;[H2,f]=keshegaotonglvboqi(wp,ws,p,s,Xn1);wp=[58,62]*2*pi;%带阻滤波器ws=[59.9,60.1]*2*pi;[H3,f]=keshedaizulvboqi(wp,ws,p,s,Xn1);H=abs(H1).*abs(H2).*abs(H3);%低通和高通和带阻组合的滤波器Y=H'.*abs(fft(Xn1));%经过滤波后心电信号频谱y=ifft(Y);%滤波后心电信号时域波形figure(2)subplot(2,2,1)plot(f,abs(H1))axis([0,150,0,1.5])title('低通滤波器')subplot(2,2,2)plot(f,abs(H2))axis([0,50,0,1.5])title('高通滤波器')subplot(2,2,3)plot(f,abs(H3))axis([0,150,0,1.5])title('带阻滤波器')河北工业大学信息工程学院信号与线性系统课程设计6subplot(2,2,4)plot(f,abs(H))axis([0,100,0,1.5])title('组合后滤波器')figure(3)plot(f,abs(Y))axis([0,100,0,80])title('滤波后心电信号频谱')figure(4)subplot(2,1,1)plot(t1,Xn1)title('滤波前信号')subplot(2,1,2)plot(t1,y)title('滤波后信号')所出图形如下河北工业大学信息工程学院信号与线性系统课程设计7河北工业大学信息工程学院信号与线性系统课程设计8结果分析：（2）、直接通过带通滤波器程序[t,Xn]=duquexinhao1('117.txt');baocun1(t,Xn)%保存读取信号[t1,Xn1]=xianxingchazhi(t,Xn);baocun2(t1,Xn1)%保存插值后信号figure(1)keshehuatu(t,Xn,t1,Xn1)%画原始信号和插值后信号波形和频谱wp=[2,80]*2*pi;ws=[0.25,99]*2*pi;p=1;s=35;[H1,f]=kesheditonglvboqi(wp,ws,p,s,Xn1);H=abs(H1);%带通Y=H'.*abs(fft(Xn1));%经过滤波后心电信号频谱y=ifft(Y);%滤波后心电信号时域波形figure(2)subplot(1,2,1)plot(f,abs(H1))axis([0,200,0,1.5])title('带通滤波器')河北工业大学信息工程学院信号与线性系统课程设计9subplot(1,2,2)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,-1.5,1.5])title('滤波前截取一部分信号')subplot(2,2,4)plot(t1,y)axis([0,1.5,-1.5,1.5])title('滤波后截取一部分信号')所出图形如下河北工业大学信息工程学院信号与线性系统课程设计10结果分析：河北工业大学信息工程学院信号与线性系统课程设计11（3）、将信号通过低通、高通组合成的带通滤波器程序[t,Xn]=duquexinhao1('117.txt');baocun1(t,Xn)%保存读取信号[t1,Xn1]=xianxingchazhi(t,Xn);baocun2(t1,Xn1)%保存插值后信号figure(1)keshehuatu(t,Xn,t1,Xn1)%画原始信号和插值后信号波形和频谱xy=[t1,Xn1];wp=0.52*2*pi;%低通滤波器滤波ws=0.62*2*pi;p=1;s=35;[H1,f]=kesheditonglvboqi(wp,ws,p,s,Xn1);wp=0.10*2*pi;%高通滤波器滤波ws=0.25*2*pi;p=1;s=35;[H2,f]=keshegaotonglvboqi(wp,ws,p,s,Xn1);H=abs(H1).*abs(H2);%低通和高通组合的带通Y=H'.*abs(fft(Xn1));%经过滤波后心电信号频谱y=ifft(Y);%滤波后心电信号时域波形figure(2)subplot(2,2,1)plot(f,abs(H1))axis([0,1,0,1.5])title('低通滤波器')subplot(2,2,2)plot(f,abs(H2))axis([0,1,0,1.5])title('高通滤波器')subplot(2,2,3)plot(f,abs(H))axis([0,1,0,1.5])title('组合带通滤波器')subplot(2,2,4)plot(f,abs(Y))axis([0,1,0,260])title('滤波后心电信号频谱')figure(3)subplot(