clc; clear all; close all;
%--  R    - kardiogramos R danteliu adresai
%--  APB  - arterines pulsines weivai reiksmiu vektorius
%--  Td   - diskretizavimo periodas =125Hz
%--  N    - is laikmenos nuskaitytu reiksmiu skaicius

Td = 0.008;

%--  Nuskaitomi duomenys
fid=    fopen('C:\G2V02\RR3_039.dat','r');
   [R_adr, N]= fread(fid,'uint32');
fclose(fid);
fid=    fopen('C:\G2V02\APB3_039.dat','r');
   [APB, M]= fread(fid,'int32');
fclose(fid);

%-- atvaizduojame signalus
RR           = zeros(1,N);
RR(R_adr) = 1;
asis         = 0:M-1;
figure(1)
plot(asis,APB,asis,230+50*RR,'r')
title('Duoti duomenys');box off;
axis tight; set(gca,'FontSize',12)

% * * * 1 Užduotis * * *
% Pulsinių bangų (PB) seką išskaidome į N atskirų pulsinių bangų.
for i = 2 : N
    weivai{i - 1} = APB(R_adr(i - 1) : R_adr(i));
end

% Kiekvienai pulsinei bangai (PB) apskaičiuojame diskrečiąją Furjė transformaciją
% tam, kad galėtume nustatyti harmonikų amplitudes.
for i = 1 : N - 1
    furje{i} = fft(weivai{i});
    for j = 1 : 12
        hrmAmpl{j}(i) = abs(furje{i}(j));
     end
end

% Sudarome 12 sekų po N reikšmių - pirmoji seka sudaryta iš pirmųjų
% harmonikų amplitudžių reikšmių, antroji iš antrųjų ir t.t.
for i = 1 : 12
    figure(1+i)
    stem(hrmAmpl{i})
    title([num2str(i),' seka']);box off;
    axis tight; set(gca,'FontSize',12)
end

% Atvaizduojame 4 pulsines bangas (PB).
for i = 1:4
    figure(14)
    subplot(4,1,i)
    plot(weivai{i})
    title([num2str(i),' pulsine banga']);box off;
    axis tight; set(gca,'FontSize',12)
end

% 2 užduotis
%Rasti PTT kiekvienai pulsinei bangai ir
%sudaryti N ilgio seką iš PTT reikšmių.
PTT = zeros(1,N-1);
for i = 1 : N - 1
    [y, x] = max(weivai{i});
    for j = x : -1 : 0
        if weivai{i}(j-1) >= weivai{i}(j)
            break;
        end
    end
    PTT(i) = j;
end
%Rezultatus atvaizduoti naudojant funkciją stem().
figure(15)
stem(PTT(1 : N - 1))
title('PPT reiksmes');box off;
axis tight; set(gca,'FontSize',12)

% 3 užduotis
% Paskaičiuoti koreliacijos koeficientus tarp PTT sekos ir
% kiekvienos 1-oje užduotyje sudarytos PB harmonikų sekos.

korr_koef = zeros(1, 12);
for i = 1 : 12
    temp = corrcoef(PTT, hrmAmpl{i});
    korr_koef(i) = temp(1, 2);
end
% Gautą koreliacijos koeficientų seką atvaizduoti naudojant funkciją stem().
figure(16)
stem(korr_koef(1 : 12))
title('Koreliacijos koeficientai');box off;
axis tight; set(gca,'FontSize',12)

% 4 užduotis
% Kiekvienai pirmos užduoties PB harmonikų sekai sudaryti ir
% atvaizduoti histogramas, reikšmių intervalą skirstant į 20 dalių.
for i = 1:12
    figure(16+i)
    histogram(hrmAmpl{i},20)
    title([num2str(i),' seka']);box off;
    axis tight; set(gca,'FontSize',12)
end