%% fltrace vstupního ekg
function [Y,yy]=ekgfil(ekg,atr)
    delka=length(atr);
    cinzk=1; %cinitel zkraceni atr o nepotrebna data
    for i=1:delka
        if atr(i,1)>12 %
            atr(i:delka-1,:)=atr(i+1:delka,:);
            cinzk=cinzk+1;
        end
    end
    test=isempty(find(atr(:,1)>12,1));
    if test==0% kontrola predesleho
        for i=1:delka
            if atr(i,1)>12 
                atr(i:delka-1,:)=atr(i+1:delka,:);
                cinzk=cinzk+1;
            end
        end
    end
cinzk=cinzk-1;
atr=atr(1:delka-cinzk,:);
pozice=int32(atr(:,3));
pozice=double(pozice);
delka=find(pozice>length(ekg),1);
    if isempty(delka)==1
        delka=length(pozice);
        pozice=pozice(1:delka);
    else
        pozice=pozice(1:delka-1);
    end
r=1;    
pr=1000;
a=40;
if length(atr)>10   % detektor QRS s vyuzitim anotaci
    r=2;            %urcuje pozici ve vektorech
    pr=find(atr(:,1)==1,2)'; %najde prvna dva radne prubehy
    prp=[pozice(pr(:,1)) pozice(pr(:,2))];
    Rhod(:,1)=max(ekg(prp(1)-30:prp(1)+30));
    Rpoz(:,1)=(find(ekg(prp(1)-30:prp(1)+30)==Rhod(:,1),1))+prp(1)-31;
    Rhod(:,2)=max(ekg(prp(2)-30:prp(2)+30));
    Rpoz(:,2)=(find(ekg(prp(2)-30:prp(2)+30)==Rhod(:,2),1))+prp(2)-31;
    stredy(:,1)=fix(Rpoz(:,1)+(Rpoz(:,2)-Rpoz(:,1))/2);
    delky(:,1)=Rpoz(:,2)-Rpoz(:,1);
    for i=pr(2):length(pozice)-1
        if pozice(i+1)<(Rpoz(r)+delky(r-1)+30) && pozice(i+1)>(Rpoz(r)+delky(r-1)-30)%existuje atr v oblasti predpokladane R
            Rhod(:,r+1)=max(ekg(pozice(i+1)-30:pozice(i+1)+30));%najde pozici R vlny
            Rpoz(:,r+1)=(find(ekg(pozice(i+1)-30:pozice(i+1)+30)==Rhod(:,r+1),1))+pozice(i+1)-31;
            stredy(:,r)=fix(Rpoz(r)+(Rpoz(r+1)-Rpoz(r))/2);
            delky(:,r)=Rpoz(r+1)-Rpoz(r);
        else % nenalezena atr v predpokladane pozici, hledá se v celem EKG cyklu
            e=find(pozice>Rpoz(r)+100,1);%zjisti se, zda je v danem rozsahu atr pozice
            if pozice(e)<(Rpoz(r)+3*fix(sum(delky(r-5:r-1)/5)))% pokud projde nalezena atr v rozsahu EKG cyklu
                %if (max(ekg(pozice(e)-30:pozice(e)+30)))>(0.5*Rhod(r)) % je nebo není R vlna
                    Rhod(:,r+1)=max(ekg(pozice(e)-30:pozice(e)+30));%najde pozici R vlny
                    Rpoz(:,r+1)=find(ekg(pozice(e)-30:pozice(e)+30)==Rhod(:,r+1),1)+pozice(i+1)-31;
                    stredy(:,r)=fix(Rpoz(r)+(Rpoz(r+1)-Rpoz(r))/2);
                    delky(:,r)=Rpoz(r+1)-Rpoz(r);
            else % nenalezena atr, hledá se v predpokladane pozici
                r=r-1; 
            end
        end
        r=r+1;
    end
end
   pocetR=length(Rpoz);
   
   for i=1:(pocetR)   
   Qp=0;   % hodnoty pro zamezení opakovaného zápisu 
       for j=1:70
            if Qp==0
                    if ekg(Rpoz(i)-j-1)>ekg(Rpoz(i)-j)
                    Qhod(:,i)=ekg(Rpoz(i)-j);
                    Qpoz(:,i)=Rpoz(i)-j+20;%+20 hodnota pro splajn
                    Qp=1;
                    end
            end
        end
   end
%% vlnková transformace rozklad s podvzorkovanim
% prvni odebrani driftu
for i=1:pocetR
hod(i)=mean(ekg(Qpoz(i)-23:Qpoz(i)-19));
end
yy = spline(Qpoz,hod,1:length(ekg))';
ekg2=ekg-yy;
%wavedec - jednodimensionalni rozklad signalu, vyuziva DWT s podvzorkovanim
% rozklad do osme urovne
% vlnka db4
[C,L]=wavedec(ekg2,8,'db4');
C=C';
%sohr - typ prahovani - meke tvrde...
[thr,sorh,keepapp]=ddencmp('den','wv',ekg2);
thr2=wthrmngr('dw1ddenoLVL','penalhi',C,L,8);
%THR = wthrmngr('sw1ddenoLVL',METHOD,SWTDEC,'one')
Y1=wdencmp('gbl',C,L,'db4',8,thr2(1),sorh,keepapp);

for i=1:pocetR % hodnoty pro vypocet slinu
hod(i)=mean(ekg2(Qpoz(i)-23:Qpoz(i)-19));
end
%kubicky spline pro odecet driftu nulove izolinie
yy2 = spline(Qpoz,hod,1:length(ekg2));
Y=Y1-yy2;

%% fourierova transformace vstupniho signalu
%     Fs=250;
%     delka=length(ekg);
%     NFFT = 2^nextpow2(delka);
%     a=(fft(ekk));
%     f = Fs/2*linspace(0,1,NFFT/2);
%     a(1:2,1)=0;
%     plot(f,abs(a(1:NFFT/2))) 

%% vykresleni vysledku filtrace

% plot(1:1000,ekg(1:1000)-mean(ekg(1:1000)),1:1000,Y(1:1000)-mean(Y(1:1000)))
% legend('Nefiltrovany signal','Filtorvany signal')


