wang_rtiqu2.m

用于心电分析的小波分析程序

%---------------提取心电信号特征程序---------------------------------
clear;clc;
%参数设置记录R个数
gg=0;
M=0;
%截取数据的长度
shu=1080;
%数据导入
load('D:\matlab\matlab7.21\work\E.mat');
s=E(1:2160);
%数据段的个数
for mm=1:fix(length(s)/shu)
    %每段数据的截取赋值
    for jj=1:shu
        ll=(mm-1)*shu+jj;
        ss(jj)=s(ll);
    end;
    %---------------小波分解-----------------
    [c l]=wavedec(ss,3,'db5');
    cd3=detcoef(c,l,3);
   
    %---------------求模极值对阈值-------------
    dd=0;xx=0;%累加参数
    for tt=2:length(cd3)-1
        %求极大值点
        if cd3(tt)>cd3(tt-1)&cd3(tt)>cd3(tt+1)
           dd=dd+1;
           cd3_Mod(dd)=cd3(tt);
        end;
        %求极小值点
        if cd3(tt)<cd3(tt-1)&cd3(tt)<cd3(tt+1)
            xx=xx+1;
            cd3_Mox(xx)=cd3(tt);
        end;
    end;
    %极大值集与极小值集分别排序（升序）
    cd3_Mod=sort(cd3_Mod);
    cd3_Mox=sort(cd3_Mox);
    max_cd3=0;
    min_cd3=0;
    %取最小的10个累加
    for ii=1:10
        min_cd3=min_cd3+cd3_Mox(ii);
    end;
    %取最大的10个累加
    for ii=length(cd3_Mod):-1:length(cd3_Mod)-10
       max_cd3=max_cd3+cd3_Mod(ii);
    end;
    %分别求每段的平均值
    max_yuzhi(mm)=max_cd3/10;
    min_yuzhi(mm)=min_cd3/10;
    %取自适应的阈值
    if mm==1
      Max_yuzhi=1/2*max_yuzhi(1);
      Min_yuzhi=1/3*min_yuzhi(1);
    else
      Max_yuzhi=1/4*(max_yuzhi(mm-1)+max_yuzhi(mm));
      Min_yuzhi=1/100*(min_yuzhi(mm-1)+min_yuzhi(mm));
    end;
    %-------------------确定R位置-----------
    %找大于最大阈值的极大值点
    jj=0;%循环计数值
    while jj<length(cd3)
        jj=jj+1;
        if jj>=length(cd3)
            break;
        end;
        if cd3(jj)>=Max_yuzhi&cd3(jj)>=cd3(jj-1)&cd3(jj)>=cd3(jj+1)
            %记录模极大值及位置
            Pmax=cd3(jj);
            P_max=jj;
            %判断是否为最后一点
            if jj==length(cd3)-1;
                break;
            end;
            %找极小值点
            for ii=1:length(cd3)-jj-1
                if cd3(jj+ii)<=cd3(jj+ii-1)&cd3(jj+ii)<=cd3(jj+ii+1)
                    break;
                end;
            end;
            %判断是否循环结束而跳出
            if ii+jj>length(cd3)
                break;
            end;
            %记录极小值及位置
            Pmin=cd3(jj+ii);
            P_min=jj+ii;
            kk=P_min;%设置循环跳转
            %模极值对计数
            M=M+1;
            %求斜率
            Rp(M)=(Pmax-Pmin)./(P_min-P_max);
            %------------求自适应斜率阈值-----------
            RR=0;%累加参数
            if M<=10
                for t=1:M
                   RR=RR+Rp(t);
                end;
                R_yuzhi=RR/(5*t);
            else
                for t=M-10:M-1
                   RR=RR+Rp(t);
                end;
                R_yuzhi=RR/100;
            end;
            %斜率阈值判断
            if Rp(M)>R_yuzhi
                %R波个数加1
               gg=gg+1;
               %找过零点
               for n=jj:jj+ii
                  if cd3(n)==0
                      break;
                  end;
               end;
             
               %记录过零点位置
               nn=n;
            else
                %消除极值对
                M=M-1;
                continue;
            end;
            Rw(gg)=nn;%cd3中R的位置
            
            %判断相邻R间的距离
            if gg==1|Rw(gg)-Rw(gg-1)<=90
                %------确定在原数据序列中的位置---------
                
                for k=1:length(cd3)
                   if k==Rw(gg)
                       cd_l(k)=1;
                   else
                       cd_l(k)=0;
                   end;
               end;
               
               %连续3次插值
               cd_l1=dyadup(cd_l);
               cd_l2=dyadup(cd_l1);
               cd_l3=dyadup(cd_l2);
               %找cd3中过零点相对应的位置
               for k=1:length(cd_l3)
                  if cd_l3(k)==1
                    break;
                  end;
               end;
              
               %计算原序列中的位置并修正 
               R_weizhi(gg)=(mm-1)*shu+k-25;
               R_fuzhi(gg)=E(R_weizhi(gg));
              
               while R_fuzhi(gg)<0.55
                 if E(R_weizhi(gg))>0&E(R_weizhi(gg))>=E(R_weizhi(gg)-1)&E(R_weizhi(gg))>=E(R_weizhi(gg)+1)
                       break;
                   end;
                   R_weizhi(gg)=R_weizhi(gg)-1;
                   R_fuzhi(gg)=E(R_weizhi(gg));
                   
                   
               end;
               jj=kk+1;
            else
                %------------重新设置参数继续返回确定过零点的位置--------
                M=M-1;
                gg=gg-1;
                if gg==0
                    jj=1;
                else
                    jj=Rw(gg);
                    if jj==0;
                        jj=jj+1;
                    end;
                   Max_yuzhi=1/2*Max_yuzhi;
                   Min_yuzhi=1/2*Min_yuzhi;
                   continue;
                   %—————————确定在原序列中的位置----------------
                   for k=1:length(cd3)
                      if k==Rw(gg)
                          cd_l(k)=1;
                      else
                          cd_l(k)=0;
                      end;
                  end;
                  %连续3次插值
                  cd_l1=dyadup(cd_l);
                  cd_l2=dyadup(cd_l1);
                  cd_l3=dyadup(cd_l2);
                  %找cd3中过零点相对应的位置
                  for k=1:length(cd_l3)
                    if cd_l3(k)==1
                        break;
                    end;
                  end;
                   R_weizhi(gg)=(mm-1)*shu+k-25;
                   R_fuzhi(gg)=E(R_weizhi(gg));
                   while R_fuzhi(gg)<0.55
                       if E(R_weizhi(gg))>0&E(R_weizhi(gg))>E(R_weizhi(gg)-1)&E(R_weizhi(gg))>E(R_weizhi(gg)+1)
                          break;
                       end;
                       R_weizhi(gg)=R_weizhi(gg)-1;
                       R_fuzhi(gg)=E(R_weizhi(gg));
                    end;
                   jj=kk;
                end;
            end;
        end;
    end;
end;
subplot(2,1,1);
plot(s);
subplot(2,1,2);
stem(R_weizhi,R_fuzhi);
M