clc; % vymaže příkazový řádek
clear all; % vymaže proměnné ve workspace
close all; % zavře všechna otevřená okna (obrázky)

load EFV1.mat % výběr dat (1 min ablačního signálu o délce 120000 vzorků)
  
figure(1);

subplot(311)
time = linspace(0,60,120000); % vytvoření časové osy x (30s)
plot(time,signal); % vykreslení ablačního signálu (30s) v průběhu RFA
xlabel('t [s] \rightarrow') % popis osy x
ylabel('U [\muV] \rightarrow') % popis osy y

subplot(312)
hdp = fir1(18,0.03, 'low'); %  DP s mezní frekvencí 40 Hz
x = filter(hdp, 1, signal); % filtrace vstupního signálu
plot(time,x) % vykreslení signálu po filtraci
xlabel('t [s] \rightarrow') % popis osy x
ylabel('U [\muV] \rightarrow') % popis osy y


subplot(313)
time = linspace(0,30,60001); % vytvoření  30s časové osy x
plot(time,signal(end-60000:end)); % vykreslení 30s vstupního signálu
hold on
plot(time,x(end-60000:end),'r') % přikreslení 30s signálu po filtraci
xlabel('t [s] \rightarrow') % popis osy x
ylabel('U [\muV] \rightarrow')% popis osy y
hold off

% prvni_cast = x(60000:70000); % část ablačního signálu - RFA start
% prvni_cast = x(80000:90000); % RFA 10s
prvni_cast = x(105000:115000); % RFA 20s
% prvni_cast = x(105000:110000); % RFA 30s

figure(2);

% time2 = linspace(0,2,4001); % vytvoření  30s časové osy
plot(prvni_cast) % vykreslení úseku
axis([0 10000 -2e3 2e3])
xlabel(' pocet vzorku \rightarrow') % popis osy x
ylabel('U [\muV] \rightarrow') % popis osy y
hold on

% axis ([0 10000 min(prvni_cast)-1 max(prvni_cast)+1]) % rozsah osy y
[x1,y1] = ginput(1); % souřednice po prvním označení myší
[x2,y2] = ginput(1); % souřadnice po druhém označení myší

minimum = round(min(prvni_cast(x1:x2)));
maximum = round(max(prvni_cast(x1:x2)));
if maximum<0
    maximum=0;
else
    maximum;
end

if minimum>0
    stem(x1,maximum,'r') % vytvoření kladného kaliperu pro první kliknutí
    stem(x1,0,'r') % vytvoření záporného kaliperu pro první kliknutí
else
    stem(x1,maximum,'r')
    stem(x1,minimum,'r')
end    
hold on

if minimum>0
    stem(x2,maximum,'r') % vytvoření kladného kaliperu pro druhé kliknutí
    stem(x2,0,'r') % vytvoření záporného kaliperu pro druhé kliknutí
else
    stem(x2,maximum,'r')
    stem(x2,minimum,'r')
end
hold on

x1 = round(x1); % zaokrouhlení souřadnice prvního kliku na celé hodnoty
x2 = round(x2); % zaokrouhlení souřadnice druhého kliku na celé hodnoty
kladne = 0; % čítač pro kladné hodnoty vzorků
zaporne = 0; % čítač pro záporné hodnoty vzorků

if minimum<0
    for i = x1:x2 % for cyklus realizující výpočet plochy pod křivkou
        if prvni_cast(i)<0
        zaporne = (zaporne-(minimum-prvni_cast(i)));% sčítání záporných hodnot vzorků
        else
        zaporne = zaporne-(minimum);
        end
    end
end

okno = ((x2-x1)*(-minimum));
pomer = 1-(zaporne/okno);

% FDHM (šířka záporné výchylky signálu)
[minimum,pozice_minima] = min(prvni_cast(x1:x2)); % nalezní minima v rámci úseku označeném kalipery a uložení jeho hodnoty a pozice
[maximum,pozice_maxima] = max(prvni_cast(x1:x2)); % nalezení maxima a uložení jeho hodnoty a pozice
pozice_maxima = floor(pozice_maxima + x1); % zaokrouhlení pozice maxima
pozice_minima = floor(pozice_minima + x1); % zaokrouhlení pozice minima

hold on
if minimum<0
    plot(x1:x2,minimum,'r'); % vytvoření záporného vodorovného kaliperu
else
    plot(x1:x2,0,'r');
end
hold on
if maximum>0
    plot(x1:x2,maximum,'r'); % vytvoření kladného vodorovného kaliperu
else
    plot(x1:x2,0,'r');
end

polovina_min = minimum/2; % výpočet poloviny minima
polovina_max = maximum/2; % výpočet poloviny maxima

minimum2=-minimum/1000;
maximum2=maximum/1000;
tmi = minimum2/(maximum2+minimum2);

sirka1 = floor(x1); % zaokrouhlená pozice prvního kliku
while prvni_cast(sirka1)>polovina_min % nalezení poloviny záporné výchylky zleva
    sirka1 = sirka1+1;
end

sirka2 = floor(x2); % zaokrouhlená pozice druhého kliku
while prvni_cast(sirka2)>polovina_min % nalezení poloviny záporné výchylky zprava
    sirka2 = sirka2-1;
end

sirka1;
sirka2;
FDHM = sirka2-sirka1; % výpočet šířky záporné výchylky v polovině

if pomer==1
    vysl = (['Zadna zaporna plocha!']);
else
    vysl =(['Index procenta negativni plochy ', num2str(pomer), ',' ' FDHM odpovídá šířce ', num2str(FDHM),' a TMI je ', num2str(tmi), '.']); %
end
hold on
title((vysl),'FontSize',12) % nadpis grafu
hold off
      