clc;
clear all;
close all;
 
%Prenos systému 
num = [1]; %čitateľ prenosu (numerator) 
den = [5 14 7 1]; %menovateľ prenosu (denominator)
Gs = tf(num,den); %prenos systému
Ts = 0.01; %vzorkovacia perióda 
t = (0:Ts:150); %časový interval
 
%Simulácia prenosu
u = [0 ones(1,length(t)-1)];
y = lsim(Gs,u,t);
title('Časová odozva dynamického systému');
plot(t,y,'r'); %vykreslenie simulácie systému
hold on;
 
%---------------------------------------------------------------------
%Aproximácia systému pomocou RMNŠ prenosom 2. rádu
 
N = 2;
%Inicializácia kovariančnej matice
value = 1e10; %inicializačná hodnota kovariančnej matice P
P = (zeros(4));
for i = 1:4 
    for j = 1:4
        if i == j
            P(i,j) = value; 
        end
    end
end
 
%Inicializácia
theta = [1;0;0;0]; %vektor parametrov
Mvst = [0 0]; %matica vstupných hodnôt
Mvyst = [0 0]; %matica výstupných hodnôt
 
%Výpočet neznámych parametrov pomocou rekurzívnej metódy najmenších štvorcov
for k = 3:length(u)
    Mvst = [Mvst u(k)];
    Mvyst = [Mvyst y(k)];
    [theta P]=RMNS_2_rad(Mvst,Mvyst,P,theta);
end
 
%Určenie prenosu systemu
GzIdent = tf([theta(3:4)'],[1 -theta(1:2)'],(Ts)) %diskrétný prenos
GsIdent = d2c(GzIdent) %prevod diskrétneho prenosu na spojitý
yIdent = lsim(GsIdent,u,t);
title('Časová odozva systému');
xlabel('Čas t [s]');
axis([0 50 -0.1 1.1]);
grid on;
plot(t,yIdent,'b'); %vykreslenie diskrétneho identifikovaného prenosu
legend('Prenos sústavy','Odhad prenosu sústavy modelom druhého rádu', 'Location', 'southoutside');