燕山大学课程设计说明书
4.2 QRS 波检测方法与程序
Q 波和 S 波通常是低幅高频波,一般 Q 波位于 S 波之前,S 波位于 R
波之后 ,由于他们是一般向下的波,所以他们的峰值点和极值是对应的。因
次在检测到 R 波向左和向右分别搜寻到极值点,对应的就是 Q 波和 S 波。
具体程序如下:
clear all;
clc;
z=textread('ECG.txt');
ECG=z(:,1);
input=ECG(1:256);
rate=ECG(100);
sig=input;
lensig=length(sig);
wtsig1=cwt(sig,6,'mexh');
lenwtsig1=length(wtsig1);
wtsig1(1:20)=0;
wtsig1(lenwtsig1-20:lenwtsig1)=0;
y=wtsig1;
yabs=abs(y);
%?
sigtemp=y;
siglen=length(y);
sigmax=[];
for i=1:siglen-2
if (y(i+1)>y(i)&y(i+1)>y(i+2))|(y(i+1)
燕山大学课程设计说明书
subplot(2,1,1),plot(sig);
subplot(2,1,2),plot(wtsig1);
%取阈值,阈值为相对幅值的差的 60%
thrtemp=sort(sigmax);
thrlen=length(sigmax);
thr=0;
for i=(thrlen-7):thrlen
thr=thr+thrtemp(i);
end;
thrmax=thr/8;
zerotemp=sort(y);
zerovalue=0;
for i=1:100
%最大幅度平均值,8 个最大幅值点的平均值
zerovalue=zerovalue+zerotemp(i);
end;
zerovalue=zerovalue/100;
值点的平均值
%最小幅度平均值,对消幅度,100 个最小幅
thr=(thrmax-zerovalue)*0.3; %最大、最小幅度的差值的 30%为判别 R 波的
阈值
%定位 R 波
rvalue=[];
for i=1:thrlen
if sigmax(i,1)>thr
rvalue=[rvalue;sigmax(i,2)];
end;
end;
rvalue_1=rvalue;
%排除误检,如果相邻两个极大值间距小于 0.4,则去掉幅度较小的一个
lenvalue=length(rvalue);
i=2;
while i<=lenvalue
if (rvalue(i)-rvalue(i-1))*rate<0.4
- 2 -
燕山大学课程设计说明书
if yabs(rvalue(i))>yabs(rvalue(i-1))
rvalue(i-1)=[];
rvalue(i)=[];
else
end;
lenvalue=length(rvalue);
i=i-1;
end;
i=i+1;
end;
lenvalue=length(rvalue);
%在原信号上精确校准
for i=1:lenvalue
if (wtsig1(rvalue(i))>0)
k=(rvalue(i)-5):(rvalue(i)+5);
[a,b]=max(sig(k));
rvalue(i)=rvalue(i)-6+b;
k=(rvalue(i)-5):(rvalue(i)+5);
[a,b]=min(sig(k));
rvalue(i)=rvalue(i)-6+b;
else
end;
end;
%打印纠正及校准前后的 R 波信号
figure(2);
subplot(2,1,1),plot(1:lensig,wtsig1,rvalue_1,wtsig1(rvalue_1),'r.')
;
subplot(2,1,2),plot(1:lensig,sig,rvalue,sig(rvalue),'r.');
%检测 Q 波
wtsig2=cwt(sig,8,'mexh');
lenrvalue=length(rvalue);
qvalue=[];
for i=1:lenrvalue
for j=rvalue(i):-1:(rvalue(i)-30)
- 3 -
燕山大学课程设计说明书
if wtsig1(rvalue(i))>0
if wtsig2(j)wtsig2(j-1)&wtsig2(j)>wtsig2(j+1)
tempqvalue=j-10;
break;
end;
end;
end;
x1=tempqvalue;
y1=sig(tempqvalue);
x2=rvalue(i);
y2=sig(rvalue(i));
a0=(y2-y1)/(x2-x1);
b0=-1;
c0=-a0*x1+y1;
dist=[];
for k=tempqvalue:rvalue(i)
%确定检测窗的起点
%倒置 R 波,取第一个正极大值
%求直线公式参数 ax+by+c=0
tempdist=(abs(a0*k+b0*sig(k)+c0))/sqrt(a0^2+b0^2);
dist=[dist;tempdist];
end;
[a,b]=max(dist);
附近
tempqvalue=tempqvalue+b-1;
l=(tempqvalue-5):rvalue(i);
[c,d]=min(sig(l));
tempqvalue=tempqvalue-6+d;
%
%
%
得到结果
qvalue=[qvalue;tempqvalue];
end;
%检测 S 波
svalue=[];
for i=1:lenrvalue-1
for j=rvalue(i):1:(rvalue(i)+100)
%求点到直线距离
%找到距离最大值,Q 波就在
%在最大值附近修正 Q 波,
- 4 -
燕山大学课程设计说明书
if wtsig1(rvalue(i))>0
if (wtsig2(j)wtsig2(j-1))&(wtsig2(j)>wtsig2(j+1))
tempsvalue=j+10;
%在小波变换域从 R 波开始向后寻找第一个极大值
break;
end;
end;
end;
x1=tempsvalue;
y1=sig(tempsvalue);
x2=rvalue(i);
y2=sig(rvalue(i));
a0=(y2-y1)/(x2-x1);
b0=-1;
c0=-a0*x1+y1;
dist=[];
for k=rvalue(i):tempsvalue
%求直线公式参数 ax+by+c=0
tempdist=(abs(a0*k+b0*sig(k)+c0))/sqrt(a0^2+b0^2);
dist=[dist;tempdist];
%求点到直线距离
%找到距离最大值,S 波就在附近
%在最大值附近修正 S 波,得到结果
end;
[a,b]=max(dist);
tempsvalue=rvalue(i)+b-1;
l=rvalue(i):(tempsvalue+10);
[c,d]=min(sig(l));
tempsvalue=rvalue(i)+d-1;
%
%
%
svalue=[svalue;tempsvalue];
end;
%检测 QRS 起点
start=[];
for i=1:lenrvalue
for j=qvalue(i):-1:(qvalue(i)-100)
if wtsig1(j)>0
start=[start;j];
- 5 -
燕山大学课程设计说明书
break;
end;
end;
end;
%打印 Q,S 波信号
qrvalue=[qvalue;rvalue];
qrvalue=sort(qrvalue);
qrsvalue=[qvalue;rvalue;svalue;start];
qrsvalue=sort(qrsvalue);
figure(3);
subplot(2,1,1),plot(1:lensig,sig,qrvalue,sig(qrvalue),'r.');
subplot(2,1,2),plot(1:lensig,sig,qrsvalue,sig(qrsvalue),'r.');
运行结果图如下:
原信号及变换信号
- 6 -
燕山大学课程设计说明书
纠正及校准前后的 R 波信号
- 7 -
燕山大学课程设计说明书
- 8 -