function [annt,annAmp,banot,QT,QTC,QW,RW,SW,QRS]=wavelimits(dirhea,dirsig,dirann,ecgnr,anot,typerec,res,nl,ti,tf,nbo_flag,Kq,Kr,Ks,Krr,Kpb,Kpe,Ktb,Kte,pco)
%----- Principal program ----------------------
%ECG wave limits detection.
%
%Input parameters:
% dirhea: header directory
% dirsig: signal directory
% dirann: annotator directory
% ecgnr: record name
% anot: annotator name; 0 if not DB-annotation file
% typerec: type of recording (0:MIT-DB, 1:Lund-Siemens)
% res: results format (0: struct, 1: text file)
% nl: leads (1:12)
% ti: begin time of processing in s.
% tf: end time of processing in s.
% nbo_flag: rejection of non-normal beats and neighbours (0:not applied, 1:applied)
% Kq: threshold Q-wave begin.
% Kr: threshold R-wave begin and end.
% Ks: threshold S-wave end.
% Krr: threshold R'-wave end.
% Kpb: threshold P-wave begin.
% Kpe: threshold P-wave end.
% Ktb: threshold T-wave begin.
% Kte: threshold T-wave end.
% pco: threshold T-wave morphology classification.
%
%Output parameter:
% banot: struct with the annotator.
% annt: struct with wave limits: Ponset, P, Poffset,
% QRSonset, Q, R, fiducial point, S, R2, QRSoffset, Tonset, T, T2,
% Toffset.
% annAmp: struct with wave amplitudes: P, Q, R, S, R2, T, T2.
% QT, QTC, QW, RW, SW, QRS: structs with position and values of the
% different waves and intervals.
%Default parameters.
if nargin<7
res=0;
end
if nargin<8
nl=1;
end
if nargin<9
ti=0;
end
if nargin<10
tf=Inf;
end
if nargin<11
nbo_flag=0;
end
if nargin<12 Kq=1.5; end
if nargin<13 Kr=5; end
if nargin<14 Ks=3; end
if nargin<15 Krr=5; end
if nargin<16 Kpb=1.35; end
if nargin<17 Kpe=2; end
if nargin<18 Ktb=2; end
if nargin<19 Kte=3.5; end
if nargin<20 pco=8; end
fidan=0;
t=[ti tf];
% ---- Getting header of the signal ----
heasig=readheader([dirhea ecgnr '.hea']);
if ~isempty(t) t=round(t*heasig.freq);
end
if isempty(t) t=[1 round(heasig.nsamp-0.1*heasig.freq)]; end
if t(1) < 1, t(1) = 1; end
if t(2) == Inf, t(2) = round(heasig.nsamp-0.1*heasig.freq); end
if anot~=0,
anname=[dirann ecgnr '.' anot];
annot=readannot([dirann ecgnr '.' anot],t);
else
disp('DB-annotation file needed')
return;
end
no_leads=heasig.nsig;
if nl>no_leads
disp('There are not nl signals at the record');
return;
end
Fs=round(heasig.freq);
nsamp=heasig.nsamp;
timeoffset=0;
if tf>nsamp tf=nsamp; end
if typerec==0,
if (heasig.fmt(nl)==16) | (heasig.fmt(nl)==212) | (heasig.fmt(nl)==61),
formato = num2str(heasig.fmt(nl));
else
error('This format is not supported by the program')
end
end
if strcmp(formato,'16')|strcmp(formato,'61'),
fid = fopen([dirsig heasig.fname(nl,:)],'rb');
fseek(fid,0,-1); % Rewind the file
if strcmp(heasig.fname(nl,1),'_'), % Siemens card recordings with MIT-type header
timeoffset=512;
fseek (fid, timeoffset,-1); % offset
end
end
if typerec==0,
fid = fopen([dirsig heasig.fname(nl,:)],'rb');
gain=heasig.gain(nl);
ioffset=heasig.adczero(nl);
if gain>0 gain=round(gain);
else gain=200;
end
% if ioffset==0
% v=getvec(fid,heasig,1,2*Fs);
% ioffset=mean(v(:,nl));
% end
end
rmax=0;
index=0;
ilastqrs=0;
fin=0;
ns=0;
prewindt=0;
a1=1;
qt1=1;
q1=1;
r1=1;
s1=1;
qrs1=1; nt1=1;
QT.pos=[]; QT.val=[];
QTC.pos=[]; QTC.val=[];
QW.pos=[]; QW.val=[];
RW.pos=[]; RW.val=[];
SW.pos=[]; SW.val=[];
QRS.pos=[]; QRS.val=[];
annt.Ponset=[]; annt.P=[]; annt.Poffset=[]; annt.QRSonset=[]; annt.Q=[];
annt.R=[]; annt.fiducial=[]; annt.S=[]; annt.R2=[]; annt.QRSoffset=[];
annt.Tonset=[]; annt.T=[]; annt.T2=[]; annt.Toffset=[];
annAmp.P=[]; annAmp.Q=[]; annAmp.R=[]; annAmp.S=[]; annAmp.R2=[]; annAmp.T=[]; annAmp.T2=[];
banot.time=[]; banot.anntyp=[]; banot.subtyp=[]; banot.num=[]; banot.chan=[];
% ---- Removing of non-QRS annotations ----
annot=isqrs(annot,heasig,t);
if isempty(annot.time), disp('no events detected'), return, end
no_beats=length(annot.time);
% ---- Removing of non-NORMAL labelled beats and adjacents ----
if nbo_flag==1
[annot,rr]=isnormal(annot,heasig,t,rr);
no_beats=length(annot.time);
end
s=sprintf('%s.b%d',ecgnr,nl-1);
if res==1
fidan=fopen(s,'wt+');
end
% ---- Reading signal segments ----
ti=t(1);
nlat=100;
iqrs=zeros(1,nlat);
atyp=char(zeros(1,nlat));
ilat=1;
nlat=min(nlat,length(annot.time));
while ilat<no_beats-3
flat=ilat+nlat+3;
if nlat+3>=length(annot.time)
tf=annot.time(nlat)+round(2*heasig.freq);
flat=length(annot.time);
nlat=nlat-3;
else tf=annot.time(flat)+round(2*heasig.freq);
end
ta=annot.time(ilat:flat);
antyp=annot.anntyp(ilat:flat);
if typerec==0 % MIT format
if strcmp(formato,'212'), % !!!!! Different formats (heasig)
X = rdsign212([dirsig ecgnr '.dat'],heasig.nsig,ti+1,tf);
elseif strcmp(formato,'16')|strcmp(formato,'61'),
fseek(fid, timeoffset+2*(ti)*heasig.nsig, -1); % Locate the pointer
X = fread(fid,[heasig.nsig tf-ti] ,'int16')';
if (strcmp(computer,'SOL2')&strcmp(formato,'16'))|(~strcmp(computer,'SOL2')&strcmp(formato,'61')),
X=swap16(sig);
end
end
X = X(:,nl); % Only selected lead
X = (X - heasig.adczero(nl))/heasig.gain(nl); % mV
elseif typerec==1 % LUND format
X=getsig(dirsig,ecgnr,[ti+1,tf],nl)';
X=X.'/1000; % mV
end
% ---- Get signals for processing ----
[Xpa,Xpb,D,F,Der]=lynfilt2(index,X,Fs,ns);
n=1;
nqrs=0;
% ---- Determine peak position of QRS complex ----
for i=1:length(ta)
if nbo_flag==0|(nbo_flag==1&antyp(i)=='N')
ibe=max(1,ta(i)-ti-round(0.2*Fs));
ien=ta(i)-ti+round(0.17*Fs);
if ien<length(Xpa)
[ymax,imax]=max(Xpa(ibe:ien)); imax=ibe+imax-1;
[ymin,imin]=min(Xpa(ibe:ien)); imin=ibe+imin-1;
%if abs(ymin)>abs(ymax)
if abs(ymin)>1.3*abs(ymax) %JGM
iqrs(n)=imin;
else iqrs(n)=imax;
end
atyp(n)=antyp(i);
nqrs=nqrs+1;
n=n+1;
end
end
end
index=1;
if iqrs(1)+ti>ilastqrs+0.4*Fs iprimerqrs=1;
elseif iqrs(2)+ti>ilastqrs+0.4*Fs iprimerqrs=2;
elseif iqrs(3)>ilastqrs+0.4*Fs iprimerqrs=3;
elseif iqrs(3)<=ilastqrs+0.4*Fs iprimerqrs=4;
end
T=linspace(ti,tf,ns);
% ---- Detection of limits and peaks ----
[POS,AMP,ANNOT,POS_ANNOT,NUMFIELD,SUBTYPEFIELD,CHANFIELD,POS_QT,VAL_QT,VAL_QTC,AMP_Q,POS_Q,AMP_R,POS_R,AMP_S,POS_S,VAL_QRS,POS_QRS,prewindt]=proces...
(fidan,X,Xpa,Xpb,D,F,Der,ti,tf,iprimerqrs,nqrs,iqrs,atyp,ns,Fs,nl,res,prewindt,Kq,Kr,Ks,Krr,Kpb,Kpe,Ktb,Kte,pco);
% ---- Recording of results ----
if (~isempty(POS_ANNOT))
a2=a1+length(POS_ANNOT)-1;
banot.time(a1:a2,1)=POS_ANNOT;
banot.anntyp(a1:a2,1)=setstr(ANNOT);
banot.subtyp(a1:a2,1)=num2str(SUBTYPEFIELD);
banot.chan(a1:a2,1)=num2str(CHANFIELD);
banot.num(a1:a2,1)=num2str(NUMFIELD);
a1=a2+1;
if ~isempty(POS_QT)
qt2=qt1+length(POS_QT)-1;
QT.pos(qt1:qt2,1)=POS_QT;
QT.val(qt1:qt2,1)=VAL_QT;
QTC.pos(qt1:qt2,1)=POS_QT;
QTC.val(qt1:qt2,1)=VAL_QTC;
qt1=qt2+1;
end
if ~isempty(POS_Q)
q2=q1+length(POS_Q)-1;
QW.pos(q1:q2,1)=POS_Q;
QW.val(q1:q2,1)=AMP_Q;
q1=q2+1;
end
if ~isempty(POS_R)
r2=r1+length(POS_R)-1;
RW.pos(r1:r2,1)=POS_R;
RW.val(r1:r2,1)=AMP_R;
r1=r2+1;
end
if ~isempty(POS_S)
s2=s1+length(POS_S)-1;
SW.pos(s1:s2,1)=POS_S;
SW.val(s1:s2,1)=AMP_S;
s1=s2+1;
end
if ~isempty(POS_QRS)
qrs2=qrs1+length(POS_QRS)-1;
QRS.pos(qrs1:qrs2,1)=POS_QRS;
QRS.val(qrs1:qrs2,1)=VAL_QRS;
qrs1=qrs2+1;
end
nt2=nt1+length(POS.Ponset)-1;
annt.Ponset(nt1:nt2)=POS.Ponset;
annt.P(nt1:nt2)=POS.P;
annt.Poffset(nt1:nt2)=POS.Poffset;
annt.QRSonset(nt1:nt2)=POS.QRSonset;
annt.Q(nt1:nt2)=POS.Q;
annt.R(nt1:nt2)=POS.R;
annt.fiducial(nt1:nt2)=POS.fiducial;
annt.S(nt1:nt2)=POS.S;
annt.R2(nt1:nt2)=POS.R2;
annt.QRSoffset(nt1:nt2)=POS.QRSoffset;
annt.Tonset(nt1:nt2)=POS.Tonset;
annt.T(nt1:nt2)=POS.T;
annt.T2(nt1:nt2)=POS.T2;
annt.Toffset(nt1:nt2)=POS.Toffset;
annAmp.P(nt1:nt2)=AMP.P;
annAmp.Q(nt1:nt2)=AMP.Q;
annAmp.R(nt1:nt2)=AMP.R;
annAmp.S(nt1:nt2)=AMP.S;
annAmp.R2(nt1:nt2)=AMP.R2;
annAmp.T(nt1:nt2)=AMP.T;
annAmp.T2(nt1:nt2)=AMP.T2;
nt1=nt2+1;
end
if n>=1
ilastqrs=iqrs(nqrs-2)+ti;
end
ilat=ilat+nlat;
nlat=min(nlat,no_beats-ilat-3);
ti=annot.time(ilat)-round(2*heasig.freq);
end
if res==1 fclose(fidan);
end
%if res==0 save('fanot','banot'); end
fclose('all');
keyboard