%% Demo1
% Slouží k vyzkoušení algoritmu proložení bi-exponenciálního modelu.
% Nejdříve jsou načtena vstupní data.
% Následně s pomocí funkce automatic_detection jsou získány společné
% pozice prahovaných parametrů.
% Na závěr jsou v pozicích vykresleny průběhy AIF ve všech voxelech řezu
% a průměrná AIF.

%% Příprava workspace
clear all;
clc;
close all;

%% Načtení datasetu
% Kvůli úspoře místa do přílohy jsou data načtena až po převodu do 4D
% matice. Pouze osmý řez z datasetu AMCAO9. Počet snímků je pouze 250. 

%% Převedení buňky do 4D matice (x,y,z,čas)
% nSlices=size(data2,1);
% nSamples=size(data2,2);
% [nRows,nCols]=size(data2{1,1});
% data4D=zeros(nRows,nCols,nSlices,nSamples);
% for sl=1:nSlices
%     for t=1:nSamples
%         data4D(:,:,sl,t)=imrotate(data2{sl,t},90); % Rotace o 90 stupňů
%     end
% end
% 
% outputFolder=[];
%nSamples=size(data4D,4);

nSamples=250;
timeAxis=(0:1:nSamples-1)*0.9075/60; % Časová osa
Ts = 0.9075; % Vzorkovací perioda

% data = data4D;

load('data_vst.mat');

%% Normalizace dat

AUC_data = sum(abs(data),4)/size(data,4); % Plocha pod křivkou
AUC_data = AUC_data(:);
AUC_data(AUC_data==0)=[];
med_data = median(AUC_data,'all'); % Výběr mediánu
data = data./med_data; % Data dělená nalezeným mediánem

% data(isinf(data))=0;
% data(isnan(data))=0;

%% Automatická detekce

[output_common] = automatic_detection(data,Ts);

figure;
imshow(output_common(:,:,:))

aifs_mean = zeros(size(data,3),size(data,4));

for j = 1
[rows, columns]=find(output_common(:,:,j)==1); % Nalezení indexů společných pozic

%% Získání AIF
aif = zeros(size(data,4)*size(data,3),length(rows));
aif = [];
aifs = [];
aif_s = [];

for i = 1:length(rows)
        aifs= squeeze(data(rows(i,1),columns(i,1),j,:)); % Vykreslení AIF ve voxelech podle nalezených společných pozic
        aif_s(i,:) = aifs;
end

figure;
plot(timeAxis,aif_s'); % Vykreslení všech AIF v řezu
xlabel('Časová osa [min.]')
ylabel('AIF [a.u.]')

aif_mean = mean(aif_s,1); % Průměrná AIF v řezu
figure;
plot(timeAxis,aif_mean); % Vykreslení průměrné AIF v řezu
xlabel('Časová osa [min.]')
ylabel('AIF [a.u.]')

end