课题二基于MATLAB平台心电信号分析系统研发设计与仿真.docx
《课题二基于MATLAB平台心电信号分析系统研发设计与仿真.docx》由会员分享,可在线阅读,更多相关《课题二基于MATLAB平台心电信号分析系统研发设计与仿真.docx(22页珍藏版)》请在冰豆网上搜索。
课题二基于MATLAB平台心电信号分析系统研发设计与仿真
课题二基于MATLABDE的心电信号分析系统的设计与仿真
1、本课题的目的
本设计课题主要研究数字心电信号的初步分析及滤波器的应用。
通过完成本课题的设计,拟主要达到以下几个目的:
矚慫润厲钐瘗睞枥庑赖。
(1)了解MATLAB软件的特点和使用方法,熟悉基于Simulink的动态建模和仿真的步骤和过程;
(2)了解人体心电信号的时域特征和频谱特征;
(3)进一步了解数字信号的分析方法;
(4)通过应用具体的滤波器进一步加深对滤波器理解;
(5)通过本课题的设计,培养学生运用所学知识分析和解决实际问题的能力。
二、课题任务
设计一个简单的心电信号分析系统。
对输入的原始心电信号,进行一定的数字信号处理,进行频谱分析。
采用Matlab语言设计,要求分别采用两种方式进行仿真,即直接采用Matlab语言编程的静态仿真方式、采用Simulink进行动态建模和仿真的方式。
根据具体设计要求完成系统的程序编写、调试及功能测试。
聞創沟燴鐺險爱氇谴净。
(1)对原始数字心电信号进行读取,由数字信号数据绘制出其时域波形。
(2)对数字信号数据做一次线性插值,使其成为均匀数字信号,以便后面的信号分析。
(3)根据心电信号的频域特征(自己查阅相关资料),设计相应的低通和高通滤波器。
(4)编程绘制实现信号处理前后的频谱,做频谱分析,得出相关结论。
(5)对系统进行综合测试,整理数据,撰写设计报告。
三、主要设备和软件
(1)PC机一台。
(2)MATLAB6.5以上版本软件,一套。
四、设计内容、步骤和要求
4.1必做部分
4.1.1利用Matlab对MIT-BIH数据库提供的数字心电信号进行读取,并还原实际波形
美国麻省理工学院提供的MIT-BIH数据库是一个权威性的国际心电图检测标准库,近年来应用广泛,为我国的医学工程界所重视。
MIT-BIH数据库共有48个病例,每个病例数据长30min,总计约有116000多个心拍,包含有正常心拍和各种异常心拍,内容丰富完整。
残骛楼諍锩瀨濟溆塹籟。
为了读取简单方便,采用其txt格式的数据文件作为我们的原心电信号数据。
利用Matlab提供的文件textread或textscan函数,读取txt数据文件中的信号,并且还原实际波形。
酽锕极額閉镇桧猪訣锥。
4.1.2对原始心电信号做线性插值
由于原始心电信号数据不是通过等间隔采样得到的,也就是说原始的心电数据并不是均匀的,而用Matlab中提供的数字滤波器处理数据时,要求数据是等间隔的。
因此设计的系统首先应对原始心电信号做线性插值处理,使其变为等间隔的数字信号,否则直接处理后会出现偏差,根据心电信号的特点,把时间分隔成0.001s。
添加的幅值点采用一次线性插值。
对二维数据进行插值,相连幅值间数据的插值根据时间进行,运算公式如下:
彈贸摄尔霁毙攬砖卤庑。
,
,
,
,
其中
是第i个数据时间点,Ai是与之对应的数据,N是两数据之间需要的插值数,
是需要插值的两点数据差,
,
时数组
依次排列,即得到了插值后等间隔的新数据。
将插值后的幅值数据另存为txt文件,以供后续程序调用。
4.1.3根据心电信号的频域特征,设计相应的低通和带通滤波器
一般正常人的心电信号频率在0.7~100HZ范围内,幅度为
(胎儿)~5mV(成人)。
人体心电信号微弱,信噪比小,因此,在采集心电信号时,易受到仪器、人体活动等因素的影响,而且所采集的心电信号常伴有干扰。
采集心电数据时,由于人的说话呼吸,常常会混有约为0.1Hz到0.25Hz频段的干扰,对于这些低频干扰,可以让信号通过一个高频滤波器,低截止频率设置为0.25,来滤除低频信号,对于高频信号干扰,可以让信号再通过一个低频滤波器,其中截止频率设置为99Hz。
也可以直接应用带通滤波器设计。
謀荞抟箧飆鐸怼类蒋薔。
(1)根据以上指标,设计模拟巴特沃斯(切比雪夫)低通、高通或带通滤波器,画出幅频特性(模拟滤波器幅频特性freqs)。
厦礴恳蹒骈時盡继價骚。
(2)根据心电信号频谱范围设计一个3阶以上模拟滤波器对心电信号进行预滤波;
(3)采用直接、级联或并联方式,实现该系统,并画出系统的信号流图;
(4)分析系统的时域特性(阶跃响应、冲击响应等),并用Matlab绘出相关波形;
(5)用Matlab分析幅频特性,并绘出相关波形;
(6)分析系统函数零极点与幅频特性的关系。
4.1.4对处理前后的心电信号分别做频谱分析
利用Matlab软件对处理前后的心电信号编程显示其频谱,分析比较滤波前后的频谱,得出结论。
如果分析频谱,滤波效果不明显,则需变动滤波器参数指标,重新设计滤波器。
通过频谱分析,多次试验确定最合适的滤波器。
茕桢广鳓鯡选块网羈泪。
4.1.5Simulink仿真
根据前面的设计,进行基于Simulink的动态仿真设计。
实现心电信号的分析和处理。
给出系统的基于Simulink的动态建模和仿真的系统方框图,同时记录系统的各个输出点的波形和频谱图。
鹅娅尽損鹌惨歷茏鴛賴。
4.2选作部分
4.2.1减少分析数据的工作量试验
(1)只截取大约2.5s,三个周期左右,大约800个采样数据进行分析;
(2)利用Matlab编程实现对10s原始心电信号数据二次采样并进行结果分析。
4.2.2simulink仿真设计
(1)将Matlab语言编程进行线性插值后的心电信号结果,直接输入到simulink模块中进行滤波处理及系统仿真;籟丛妈羥为贍偾蛏练淨。
(2)在仿真环境下,自己设计simulink的线性插值模块,直接输入原始心电信号,进行插值处理再进行系统的功能仿真。
預頌圣鉉儐歲龈讶骅籴。
4.2.350Hz工频陷波器设计
由于电子设备采集到的信号经常会混有电源线干扰。
电源线干扰是以50Hz为中心的窄带噪声,带宽小于1Hz。
设计相应的带阻滤波器滤除电源线干扰,并对处理后的信号做频谱分析。
渗釤呛俨匀谔鱉调硯錦。
五、课程设计报告要求
(1)设计报告书包括内容:
课程设计题目,课程设计目的和意义,设计方案,详细设计步骤,设计结果(原理图等),测试和仿真结果(图形或数据)及其分析,其它有明确要求的设计内容,结论,参考文献等。
铙誅卧泻噦圣骋贶頂廡。
(2)提交课程设计报告时应同时提交相关设计和仿真分析材料(程序、结果等)的电子版。
六、参考文献
[1]北京迪阳正泰科技发展公司.综合通信实验系统——信号与系统指导书(第二版).2006,6
[2]丁玉美.数字信号处理(第二版).西安电子科技大学出版社,2001
[3]吴大正.信号与线性系统分析(第四版).高等教育出版社,2005,8
[4]谢嘉奎.电子线路--线性部分(第四版).高等教育出版社,2003,2
[5]陈后金.信号分析与处理实验.高等教育出版社,2006,8
七、附录——设计原理
1.心电信号的读取
txt格式的数据文件内容及格式如图1-1所示(以100.txt为例)。
图1-1txt格式的心电数据文件
其中文件的第一列为采样时间,第二列是在以MLII这种导联方式所得到的采样数据,第三列是以V5这种导联方式所得到的采样数据,全文件记录了约为10s的心电数据,3600个采样数据,每一行数据之间用Tab符分隔。
擁締凤袜备訊顎轮烂蔷。
由于数据文件中后两列数据是对同一种心电信号进行不同的导联方式所得到的采样数据,所以可以只采用其中的一种采样数据,摒弃另外一种,即可完成对此心电信号的分析。
全部的心电文件记录时间约为10s,共计12个左右周期的心电信号。
贓熱俣阃歲匱阊邺镓騷。
实际设计心电信号数据文件时应注意:
(1)数据文件的前两行为解释说明文字,不是真正的信号数据,读取信号程序要能够自动忽略前两行文字,只读取真正的数字信号数据(严禁自己手动删除原心电数据文件中的前两行数据,必须通过程序来实现忽略前两行文字的目的)。
坛摶乡囂忏蒌鍥铃氈淚。
(2)利用数组函数分别将文件的前两个列分别读入一个一维数组。
(3)最后利用已经转为数值的分别代表心电信号时间和幅值的两个一维数组,图形化还原原始心电信号波形,在此推荐利用plot(x,y)函数对数据做图形化显示。
蜡變黲癟報伥铉锚鈰赘。
2.心电信号的线性插值处理
根据上文中提到的插值公式,以此为原理,设计Matlab程序,对心电信号数据做线性插值处理。
插值完以后的数据应该是时间均匀的、以0.001秒为间隔的。
買鲷鴯譖昙膚遙闫撷凄。
此步骤的实现主要是基于Matlab中的数组操作函数来实现,建议一定首先熟悉并掌握Matlab中的所有数组操作函数的作用和操作方法。
綾镝鯛駕櫬鹕踪韦辚糴。
其中一种插值方法的思路是:
将第一步中读取的心电信号数据的时间数据和幅值数据分别存放在一个一维数组中。
然后利用for循环结构把所有数据依次读取进来。
判断时间数据数组中前后两个相邻的数据间隔是否为0.001s,如果是则判断下一对相邻两个数据;如果间隔大于0.001s则进行一维插值处理。
驅踬髏彦浃绥譎饴憂锦。
注意对时间数据做插值的同时一定不要忘记对幅值数据同样做插值处理,时间数据和幅值数据一定是相互对应的。
3.滤波器设计
3.1模拟滤波器设计原理
(1)模拟巴特沃思滤波器原理
巴特沃斯滤波器具有单调下降的幅频特性:
在小于截止频率
的范围内,具有最平幅度的响应,而在
后,幅频响应迅速下降。
猫虿驢绘燈鮒诛髅貺庑。
巴特沃思低通滤波器幅度平方函数为:
(3-1)
式中N为滤波器阶数,
为3dB截止角频率。
将幅度平方函数写成s的函数形式:
(3-2)
该幅度平方函数有2N个等间隔分布在半径为
的圆上的极点
,
为了形成稳定的滤波器,取左半平面的N个极点构成
,即:
锹籁饗迳琐筆襖鸥娅薔。
(3-3)
为使设计统一,将频率归一化,得到归一化极点
,相应的归一化系统函数为:
(3-4)
多项式形式为:
(3-5)
(2)模拟切比雪夫滤波器原理
切比雪夫滤波器的幅频特性具有等波纹特性,有两种形式,在通带内等波纹、阻带单调的是
型滤波器,在通带内单调、在阻带内等波纹的是Ⅱ滤波器。
以
型滤波器为例。
構氽頑黉碩饨荠龈话骛。
切比雪夫滤波器的幅度平方函数为:
(3-6)
ε为小于1的正数,表示通带内幅度波动的程度。
Ωp称为通带截止频率。
令λ=Ω/Ωp,称为对Ωp的归一化频率。
CN(x)为N阶切比雪夫多项式。
幅度平方函数的极点是分布在bΩp为长半轴,aΩp为短半轴的椭圆上的点。
同样取s平面左半平面的极点构成
:
輒峄陽檉簖疖網儂號泶。
(3-7)
进行归一化,得到:
其中
,
3.2模拟滤波器数字化原理
将模拟滤波器转化为数字滤波器在工程上常用的有脉冲响应不变法和双线性变换法。
脉冲响应不变法是一种时域上的转换方法,它是数字滤波器的单位取样响应在采样点上等于模拟滤波器的单位冲激响应,即:
尧侧閆繭絳闕绚勵蜆贅。
(3-8)
设模拟滤波器只有单阶极点,其系统函数为:
(3-9)
对
进行拉氏反变换得到
,对
进行等间隔采样,得到
,对
进行Z变换,得到数字滤波器系统函数:
(3-10)
这种方法s和z的关系是:
。
该方法的优点是频率坐标变换时线性的切数字滤波器的单位脉冲响应完全模仿模拟滤波器的单位冲激响应,时域特性逼近好;缺点是会产生频谱混叠现象,适合低通、带通滤波器的设计,不适合高通、带阻滤波器的设计。
识饒鎂錕缢灩筧嚌俨淒。
双线性变换法为了克服频谱混叠现象,采用非线性频率压缩方法,将整个频率轴上的频率范围压缩到
之间,再用
转换到Z平面上。
凍鈹鋨劳臘锴痫婦胫籴。
这种方法s和z的关系是:
。
该方法克服了频谱混叠现象,但带来了频率坐标变换的非线性:
,由模拟滤波器系统函数转换为数字滤波器系统函数公式为:
恥諤銪灭萦欢煬鞏鹜錦。
(3-11)
3.3数字高通、带通、带阻滤波器的设计
上述滤波器可以借助于模拟滤波器的频率变换设计一个所需类型的模拟滤波器,再通过双线性变换法将其转换成所需类型的数字滤波器。
鯊腎鑰诎褳鉀沩懼統庫。
首先确定所需类型数字滤波器的技术指标;然后将数字滤波器技术指标按照公式
转换成所需类型滤波器的模拟域技术指标;将所需类型滤波器的模拟域技术指标转换成低通滤波器技术指标;设计归一化模拟低通滤波器;去归一化得到模拟低通滤波器的系统函数;将模拟低通滤波器转换为所需类型的模拟滤波器;最后通过双线性变换法转换成所需类型的数字滤波器。
硕癘鄴颃诌攆檸攜驤蔹。
4.使用Simulink建模和仿真的过程
启动Matlab后,在命令窗口中输入命令“simulink”,打开simulink模块库窗口(使用命令‘simulink3’可以打开老版本的simulink模块库界面)。
典型的Simulink模块包括三个部分:
输入模块、状态模块、输出模块。
阌擻輳嬪諫迁择楨秘騖。
4.1Simulink模块库简介
(1)Continuous(连续模块)库
(2)Discrete(离散模块)库
(3)Math(数学模块)库
(4)Sinks(信号输出模块)库:
常用模块为Scope(示波器模块)、XYGraph(二维信号显示模块)、Display(显示模块)氬嚕躑竄贸恳彈瀘颔澩。
(5)Sources(信号源模块)库(如图4-1所示),常见模块有:
Constant(输入常数模块)、SignalGenerator(信号源发生器模块)。
釷鹆資贏車贖孙滅獅赘。
(6)SignalGenerator用于产生不同的信号波形,其中包括:
正弦波、方波、锯齿波信号。
Sources(信号源模块)还包括其它常用模块:
Ramp(斜坡输入信号)、SineWave(正弦波输入信号)、Step(阶跃输入信号)、Clock(时间信号)、Pulse(脉冲信号)等。
怂阐譜鯪迳導嘯畫長凉。
图4-1Sources(信号源模块)库
4.2利用Simulink建立滤波器仿真实例
设某一六阶Butterworth数字低通滤波器的系统函数为:
(1)采用级联型结构实现该数字低通滤波器;
(2)观察对实际心电图信号的滤波效果。
人体心电图信号在测量过程中往往受到工业高频干扰,所以必须经过低通滤波处理。
已知某一实际心电图信号的采样序列如下:
谚辞調担鈧谄动禪泻類。
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],共56点。
嘰觐詿缧铴嗫偽純铪锩。
解:
滤波器的系统函数H(z)可以写成以下形式:
其中,
H(z)由三个二阶子网络H1(z)、H2(z)和H3(z)级联组成,如图4-2所示。
图4-2滤波器H(z)的级联组成
仿真过程如下:
(1)建立模型
首先打开simulink模块库窗口,在simulink模块库窗口中单击菜单项“File/New/Model”,即可以建立一个新的simulink模型文件。
如图4-3所示。
熒绐譏钲鏌觶鷹緇機库。
利用鼠标单击Simulink模块库窗口中的Continuous子库,选取传递函数模块TransferFcn,将它拖动到新建模型文件窗口的合适的位置。
然后对模型模块进行参数设置和修改,单击右键从快捷菜单中选取“TransferFcnparameters…”修改传递函数参数,在弹出对话框中的传递函数分子系数“Numerator:
”栏填入[121],在传递函数分母系数“Denominator:
”栏填入[11.2686-0.7051],其余参数使用默认值。
若需要进一步了解该模块的参数设置说明,可以单击该对话框下方的“help”按钮。
然后确认,就得到了需要仿真的传递函数。
如图4-4所示。
鶼渍螻偉阅劍鲰腎邏蘞。
采用同样的方法,在Simulink模块库中的Sources子库中选取激励信号源,本例选取simin(FromWorkspace),并将它拖动到建模文件窗口的合适的位置。
然后在Sinks子库中选取示波器(Scope)模块作为系统输出波形显示。
接下来利用鼠标将三个模块连接起来。
纣忧蔣氳頑莶驅藥悯骛。
连接方法如下:
模块外部的大于符号“>”分别表示信号的输入输出节点,为了连接两个模块的输入输出节点,可以将鼠标置于节点处,这时鼠标显示为“十”字形状,拖动鼠标到另一个模块的端口,然后释放鼠标按钮,则可以形成带箭头的连线,箭头方向表示信号的流向。
颖刍莖蛺饽亿顿裊赔泷。
图4-3Simulink模块库窗口和新建模型文件窗口
图4-4修改仿真模型的参数
建立好的系统模型可以存盘为模型文件,扩展名为“mdl”。
如图4-5所示。
图4-5完成的建模方框图
(2)参数设置
建立好系统模型后,接着要对输入信号源进行参数设置。
输入信号源的参数设置界面如图4-6所示。
最后双击示波器模型图标,打开示波器显示窗口。
在快捷键设置菜单中设置为自动刻度。
濫驂膽閉驟羥闈詔寢賻。
(3)启动仿真
以上工作完成后,可通过建模窗口菜单项“Simulation/Start”启动仿真,也可以单击工具栏上的小三角按钮启动仿真。
銚銻縵哜鳗鸿锓謎諏涼。
图4-6输入信号源的参数设置界面
4.3SignalProcessingBlockset模块仿真过程
SignalProcessingBlockset提供了包括变换、矩阵运算、FIR、IIR、自适应和多速率滤波、谱分析和实时数据I/O等模块。
利用这个工具箱对数字滤波器进行仿真。
图4-7为SignalProcessingBlockset模块库,图4-8为DSPSources模块库,图4-9为FilterDesigns模块库。
挤貼綬电麥结鈺贖哓类。
图4-7SignalProcessingBlockset模块库
图4-8DSPSources模块库
图4-9FilterDesigns模块库
通过调用Simulink中的功能模块可以构成数字滤波器的仿真,在仿真过程中,可以双击各功能模块,随时改变参数,获得不同状态下的仿真结果。
赔荊紳谘侖驟辽輩袜錈。
(1)Matlab自带的滤波器设计与分析工具(FDAT001)介绍
在Matlabcommand窗口下,执行FDATool命令,将弹出一个FDATool窗口,Matlab自带的滤波器设计与分析工具(FDAT001)设计滤波器,FDATool是交互式信号分析处理环境,通过这一窗口用户能够方便地进行滤波器设计和分析。
FDATool(FilterDesign&AnalysisT001)是Matlab信号处理工具箱里专用的滤波器设计分析工具,Matlab6.0以上的版本还专门增加了滤波器设计工具箱(FilterDesignToolbox)。
FDATool可以设计几乎所有的基本的常规滤波器,包括FIR和IIR的各种设计方法。
它操作简单,方便灵活。
FDATool界面总共分两大部分,一部分是DesignFilter,在界面的下半部,用来设置滤波器的设计参数,另一部分则是特性区,在界面的上半部分,用来显示滤波器的各种特性。
DesignFilter部分主要分为:
FilterType(滤波器类型)选项,包括Lowpass(低通)、Highpass(高通)、Bandpass(带通)、Bandstop(带阻)和特殊的FIR滤波器。
DesignMethod(设计方法)选项,包括IIR滤波器的Butterworth(巴特沃思)法、ChebyshevTypeI(切比雪夫I型)法、ChebyshevTypeII(切比雪夫II型)法、Elliptic(椭圆滤波器)法和FIR滤波器的Equiripple法、Least—Squares(最小乘方)法、Window(窗函数)法。
FilterOrder(滤波器阶数)选项,定义滤波器的阶数,包括SpecifyOrder(指定阶数)和MinimumOrder(最小阶数)。
在SpecifyOrder中填入所要设计的滤波器的阶数(N阶滤波器,SpecifyOrder=N.1),如果选择MinimumOrder则Matlab根据所选择的滤波器类型自动使用最小阶数。
塤礙籟馐决穩賽釙冊庫。
心电滤波需要设计低通滤波器、高通滤波器和50Hz陷波器。
(2)低通滤波器设计
首先设计低通滤波器,在FilterType中选择Lowpass(低通滤波器);在DesignMethod选项中选择FIR(FIR滤波器),接着在WindowSpecifications选项中选取Equiripple;指定FilterOrder项中的Minimum;由于采用Equiripple法设计,只要给出Fstop(阻带截止频率)、Fpass(通带截止频率),选取Fpass=100Hz、Fstop=102Hz。
裊樣祕廬廂颤谚鍘羋蔺。
设置完以后点击DesignFilter即可得到所设计的FIR滤波器。
通过菜单选项Analysis可以在特性区看到所设计滤波器的幅频响应、相频响应、零极点配置和滤波器系数等各种特性。
设计完成后将结果保存为1.fda文件。
在设计过程中,可以对比滤波器幅频相频特性和设计要求,随时调整参数和滤波器类型,低通滤波器的幅频如图4-10所示。
仓嫗盤紲嘱珑詁鍬齊驁。
图4-10低通滤波器的参数设定及幅频特性
(3)高通滤波器设计
设计高通滤波器,在FilterType中选择highpass(高通滤波器);在DesignMethod选项中选择IIR(IIR滤波器),接着在WindowSpecifications选项中选取Butterworth;指定FilterOrder项中的Minimu;由于采用Butterworth法设计,只要给出Fstop(阻带截止频率)、Fpass(通带截止频率)。
选取Fstop=0.25Hz、Fpass=1Hz。
设置完以后点击DesignFilter即可得到所设计的IIR滤波器。
通过菜单选项Analysis可以在特性区看到所设计滤波器的幅频响应、相频响应、零极点配置和滤波器系数等各种特胜。
设计完成后将结果保存为2.fda文件。
在设计过程中,可以对比滤波器幅频相频特性和设计要求,随时调整参数和滤波器类型,高通滤波器的幅频和相频特性如图4-11示。
绽萬璉轆娛閬蛏鬮绾瀧。
图4-11高通滤波器的参数设定及幅频特性
(4)50Hz陷波器设计
设计50Hz陷波器,在FilterType中选择Bandstop(带阻滤波器);在DesignMethod选项中选择IIRWindow(HR滤波器巴特沃思法),接着在WindowSpecifications选项中选取Butterworth;指定FilterOrder项中的Minimu;由于采用Butterworth法设计,只要给出Fstopl(下阻带截止频率)、Fpassl(通带下限截止频率)、Fpass2(通带上限截止频率)、Fstop2(上阻带截止频率)。
选取Fstopl=49Hz、Fpassl=48Hz、Fstop2=51Hz、Fpass2=53Hz。
设置完以后点击DesignFilter即可得到所设计的IIR滤波器。
通过菜单选项Analysis可以在特性区看到