标签：en end buf1 电信号 %% 算法 QRS AMP
参考链接：ECG ×AI: 机器/深度学习的ECG应用入门（3）_Aiwiscal的博客-CSDN博客_ecg 深度学习
function [QRS_amp,QRS_ind] = DS_detect( ecg_i,gr ) %% function [QRS_amp,QRS_ind]=DS_detect(ecg_i,gr) %% 输入 % ecg_i : 原信号，一维向量;1*N形式 % gr : 绘图与否，0：不绘图，1：绘图 %% 输出 % QRS_amp:QRS波振幅. % QRS_ind:QRS波索引. % 绘制图像. %% 作者：刘文涵（WHliu@whu.edu.cn） %% 版本：1.1 %% 04-24-2018 %% 代码: if nargin < 2 gr = 1; if nargin<1 error('The algorithm need a input:ecg_i.'); end end if ~isvector(ecg_i) error('ecg_i must be a row or column vector.'); end fs=360; size(ecg_i,2) if size(ecg_i,2)<round(1.5*fs)+1 error('The algorithm need a longer input.'); end tic, s=ecg_i; N=size(s,2); ECG=s; FIR_c1=[0.0041,0.0053,0.0068,0.0080,0.0081,0.0058,-0.0000,-0.0097,-0.0226,... -0.0370,-0.0498,-0.0577,-0.0576,-0.0477,-0.0278,0,0.0318,0.0625,0.0867,... 0.1000,0.1000,0.0867,0.0625,0.0318,0,-0.0278,-0.0477,-0.0576,-0.0577,... -0.0498,-0.0370,-0.0226,-0.0097,-0.0000,0.0058,0.0081,0.0080,0.0068,... 0.0053,0.0041]; % 使用fdatool设计并导出的滤波器系数,带通FIR,15~25Hz,详情使用fdatool打开DS1.fda查看 FIR_c2=[0.0070,0.0094,0.0162,0.0269,0.0405,0.0555,0.0703,0.0833,0.0928,... 0.0979,0.0979,0.0928,0.0833,0.0703,0.0555,0.0405,0.0269,0.0162,0.0094,... 0.0070]; % 使用fdatool设计并导出的滤波器系数,低通FIR,截止频率5Hz,详情使用fdatool打开DS2.fda查看 l1=size(FIR_c1,2); ECG_l=[ones(1,l1)*ECG(1) ECG ones(1,l1)*ECG(N)]; % 数据点延拓，防止滤波边缘效应； ECG=filter(FIR_c1,1,ECG_l); % 使用filter滤波； ECG=ECG((l1+1):(N+l1)); % 前面延拓了数据点，这里截取有用的部分； %% 双斜率处理 a=round(0.015*fs); % 两侧目标区间0.015~0.060s; b=round(0.060*fs); Ns=N-2*b; % 确保在不超过信号长度； S_l=zeros(1,b-a+1); S_r=zeros(1,b-a+1); S_dmax=zeros(1,Ns); for i=1:Ns % 对每个点双斜率处理； for k=a:b S_l(k-a+1)=(ECG(i+b)-ECG(i+b-k))./k; S_r(k-a+1)=(ECG(i+b)-ECG(i+b+k))./k; end S_lmax=max(S_l); S_lmin=min(S_l); S_rmax=max(S_r); S_rmin=min(S_r); C1=S_rmax-S_lmin; C2=S_lmax-S_rmin; S_dmax(i)=max([C1 C2]); end %% 再次进行低通滤波，思路与上述带通滤波一致 l2=size(FIR_c2,2); S_dmaxl=[ones(1,l2)*S_dmax(1) S_dmax ones(1,l2)*S_dmax(Ns)]; S_dmaxt=filter(FIR_c2,1,S_dmaxl); S_dmaxt=S_dmaxt((l2+1):(Ns+l2)); %% 滑动窗口积分 w=8;wd=7; d_l=[zeros(1,w) S_dmaxt zeros(1,w)]; % 零延拓，确保所有的点都可以进行窗口积分 m=zeros(1,Ns); for n=(w+1):(Ns+w) % 滑动窗口； m(n-w)=sum(d_l(n-w:n+w)); % 积分； end m_l=[ones(1,wd)*m(1) m ones(1,wd)*m(Ns)]; %% 双阈值检测与动态变化 QRS_buf1=[]; % 存储检测到的QRS波索引 AMP_buf1=[]; % 存储最近检测到的8个QRS波对应特征信号的波峰值 thr_init0=0.4;thr_lim0=0.23; thr_init1=0.6;thr_lim1=0.3; %% 阈值变化的初始值和下限设置 en=-1; % 标记波峰检出情况，高于高阈值--1，高低阈值之间--0，未检出-- -1 thr0=thr_init0; thr1=thr_init1; thr1_buf=[]; % 阈值缓存，记录阈值变化情况； thr0_buf=[]; for j=8:Ns t=1; cri=1; while t<=wd&&cri>0 % 检测候选波峰； cri=((m_l(j)-m_l(j-t))>0)&&(m_l(j)-m_l(j+t)>0); t=t+1; end if t==wd+1 N1=size(QRS_buf1,2); %N1:已经检测到的QRS波个数 if m_l(j)>thr1 % 高于高阈值时的处理 if N1<2 % N1小于2时直接存储； QRS_buf1=[QRS_buf1 (j-wd)]; % j-wd 减去了滑动窗口积分带来的延迟； AMP_buf1=[AMP_buf1 m_l(j)]; en=1; else dist=j-wd-QRS_buf1(N1); if dist>0.24*fs % 检测波峰距离； QRS_buf1=[QRS_buf1 (j-wd)]; AMP_buf1=[AMP_buf1 m_l(j)]; en=1; else if m_l(j)>AMP_buf1(end) % 不应期处理 QRS_buf1(end)=j-wd; AMP_buf1(end)=m_l(j); en=1; end end end else % 特征峰值低于高阈值 if N1<2&&m_l(j)>thr0 % 特征峰值在两阈值之间 QRS_buf1=[QRS_buf1 (j-wd)]; AMP_buf1=[AMP_buf1 m_l(j)]; en=0; else if m_l(j)>thr0 % 特征峰值在两阈值之间 dist_m=mean(diff(QRS_buf1)); dist=j-wd-QRS_buf1(N1); if dist>0.24*fs && dist>0.5*dist_m % 不应期检测，并且，波峰要距离足够远（> 平均距离的一半） QRS_buf1=[QRS_buf1 (j-wd)]; AMP_buf1=[AMP_buf1 m_l(j)]; en=0; else if m_l(j)>AMP_buf1(end) QRS_buf1(end)=j-wd; AMP_buf1(end)=m_l(j); en=0; end end else en=-1; end end end N2=size(AMP_buf1,2); if N2>8 AMP_buf1=AMP_buf1(2:9); % 确保只存储最近的8个特征波峰； end % 下面的if与博文中的公式对应 if en==1 thr1=0.7*mean(AMP_buf1); thr0=0.25*mean(AMP_buf1); else if en==0 thr1=thr1-(abs(m_l(j)-mean(AMP_buf1)))/2; thr0=0.4*m_l(j); end end end if thr1<=thr_lim1 % 确保阈值高于下限 thr1=thr_lim1; end if thr0<=thr_lim0 thr0=thr_lim0; end thr1_buf=[thr1_buf thr1]; thr0_buf=[thr0_buf thr0]; end delay=round(l1/2)-2*w+2; QRS_ind=QRS_buf1-delay; % 减去延迟，得到最终结果； QRS_amp=s(QRS_ind); toc if gr==1 %绘图 subplot(2,1,1);plot(m);axis([1 size(m,2) -0.3 1.6*max(m)]); hold on;title('Feature signal and thresholds');grid on; plot(QRS_buf1,m(QRS_buf1),'ro'); plot(thr1_buf,'r'); plot(thr0_buf,'k'); legend('Feature Signal','QRS Locations','Threshold1','Threshold0'); subplot(2,1,2);plot(s);%axis([1 size(s,2) min(s) 1.5*max(s)]); xlabel('n');ylabel('Voltage / mV'); hold on;title('The result on the raw ECG');grid on; plot(QRS_ind,QRS_amp,'ro'); legend('Raw ECG','QRS Locations'); end end
标签：en,end,buf1,电信号,%%,算法,QRS,AMP 来源： https://www.cnblogs.com/JiangYue-04/p/16070848.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。