clear all; close all; clc; g=9.81;
set(groot,'defaulttextinterpreter','latex');  
set(groot, 'defaultAxesTickLabelInterpreter','latex');  
set(groot, 'defaultLegendInterpreter','latex');

% vars = importdata('Odozva1.csv');
vars = importdata('Vzorka2dlhsia0.csv');
t=vars.data(:,1);
V=vars.data(:,2);
t=t-t(1);
grid on;
Ts = mean(diff(t));                                                 
Fs = 1/Ts;
Fn=Fs/2;
L = length(V);
Y = fft(V);
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
f = Fs*(0:(L/2))/L;
[Vsort,beginidx]=sort(abs(V),'descend');

toffset = t(beginidx(1))-t(1);
figure
subplot(2,1,1)
    plot(f,P1) 
    grid on
    title('A-f charakteristika')
    [~,idx]=max(P1);
    xlabel('f [Hz]')
    ylabel('Mag(f)')
    xlim([0 30])
    
subplot(2,1,2)
    pk = [7.6707, 6.8369, 6.0032,5.4195, 5.3361,5.068, 4.6691, 4.2522,4.0855,3.752, 3.5081, 3.4185, 3.085,  2.8348, 2.7514,2.6681, 2.3346, 2.2512,...%minusove peaky
        2.1678,2.0011, 1.9177, 1.8343, 1.6675, 1.5842, 1.4174,1.334];
    ksi = zeros(1,numel(pk)-1);
    for i =1:numel(pk)-1
        delta = (1/(i))*log( pk(1)/pk(i+1) ) ;
        ksi(i) = delta/sqrt((4*pi^2)+(delta^2));
    end
    envdown = exp(-mean(ksi)*f(idx)*2*pi.*t)*max(abs(Vsort(1)));
    plot(t-toffset,V)
    hold on
    grid on
    plot(t,+envdown)
    plot(t,-envdown)
    xlabel('t [s]');
    ylabel('V [V]');
    xi=mean(ksi);

figure('units','normalized','outerposition',[0 0 1 1])
plot(t-toffset,V)
hold on
grid on
plot(t,+envdown)
plot(t,-envdown)
xlabel('t [s]');
ylabel('V [V]');
title('Tlmené harmonické kmity')
t_mer=t-toffset;
V_mer = V;
save('KSI_vz2','xi')