%% Kittiphan Techakittiroj (27 Jan 1996)

%% Input
figure(1);
k = 0:1:1000;
x = [ zeros(11,1); 100*ones(990,1)];
plot(k*0.01,x);
title('an example input');
ylabel('Armature Voltage (V)');  xlabel('time (s)');  grid;
input = [ x ];


%% Model
figure(2);
B0 = [         0 ];   A0 = [      1];  %% delay 0
B1 = [ 4.9585e-5 ];   A1 = [-1.9752];  %% delay 1
B2 = [ 4.9174e-5 ];   A2 = [ 0.9753];  %% delay 2
AA = [ A0 A1 A2];
BB = [ B0 B1 B2];
Model = arx2th( AA, BB, 1, 1);
y = idsim(input, Model);
plot(k*0.01,y);
title('an example output from the simulation');
ylabel('Armature Voltage (V)');  xlabel('time (s)');  grid;


%% Add noise to the output
figure(3);
y_n = y + max(y)*0.1*(randn(size(y)) - 0.5);
plot(k*0.01,y_n);
title('an example output from the simulation');
ylabel('Armature Voltage (V)');  xlabel('time (s)');  grid;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% ARX Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Identification
ident_data = [y input];
na = [2];
nb = [2];
nc = [1];
[I_A I_B] = th2arx(arx(ident_data, [na nb nc]));

%% Compare
I_A
I_B

%% Identification
ident_data = [y_n input];
na = [2];
nb = [2];
nc = [1];
[I_A I_B] = th2arx(arx(ident_data, [na nb nc]));

%% Compare
I_A
I_B


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% RARX Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
lam = 0.9;


%% Identification
ident_data = [y input];
na = [2];
nb = [2];
nc = [1];
I_Model = rarx(ident_data, [na nb nc], 'ff', lam);

%% Compare
I_Model(length(ident_data),:)

%% Identification
ident_data = [y_n input];
na = [2];
nb = [2];
nc = [1];
I_Model = rarx(ident_data, [na nb nc], 'ff', lam);

%% Compare
I_Model(length(ident_data),:)
