%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%
%%% Institut fuer Regelungstechnik, TU Braunschweig
%%%
%%% Identifikation dynamischer Systeme
%%% Identifikation mittels Korelationsanalyse
%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%% Strecke G (System das identifiziert werden soll)
G = tf([1 -5], [1 0.5 1])

%%% Sprungantwort
figure(1)
clf
step(G)

%%% Wahl der Abtastzeit fuer Identifikation
T = 0.25

%%% diskrete Regelstrecke mit Halteglied
Gd = c2d(G, T)

%%% Sprungantwort der diskreten Regelstrecke mit Halteglied
figure(1)
hold on
step(Gd)


%%% PRBS fuer Identifikation

%%% Ordnung PRBS
n = 5;

%%% Taktzeit PRBS als Vielfaches von T
k = 4;

%%% Laenge einer Periode des PRBS
N = k * (2^n - 1);

%%% Periodendauer PRBS
Tp = N * T;

%%% Perioden PRBS
p = 10;

%%% Amplitude PRBS
V = 1;

%%% Generierung des PRBS
u = prbs(n, k, p, V, 0);

%%% Mittelwert PRBS entfernen
u = u - mean(u);

%%% Generierung eines Zeitvektors
t = (0:T:p*Tp-T)';

%%% Plotten des PRBS
figure(2)
clf
stairs(t, u, 'b-')
zoom

%%% Simulation Regelstrecke mit PRBS-Anregung
[y, tsim] = lsim(Gd, u);

%%% Plotten des Regelstrecken-Ausgangssignals
figure(3)
hold on
stairs(t, y, 'r-')


%%% Impulsantwort Regelstrecke
figure(4)
clf
impulse(G, 'b-')


%%% Identifikation mittels Korrelation
[cuy, cuylags] = xcorr(u, y, 'biased');

%%% Mitte
m = p * N;

%%% Impulsantwort durch Identifikation
figure(4)
hold on
plot(t(1:N), cuy(m:m-1+N)/(V^~2*k*T), 'r-')



%%% Addition von Stoerungen zum Ausgang der Regelstrecke
var = 1.0
yn = y + var * randn(length(y), 1);

%%% Plotten des gestoerten Ausgangssignals
figure(3)
stairs(t, yn, 'g-')

%%% Identifikation mittels Korrelation: gestoerter Fall
[cuyn, cuylagsn] = xcorr(u, yn, 'biased');

%%% Impulsantwort durch Identifikation: gestoerter Fall
figure(4)
plot(t(1:N), cuyn(m:m-1+N)/(V^~2*k*T), 'g-')
