% program ex_15_3.m to demo model selection
% programmed by B.Joseph Dec 27, 2000
%
% first generate data
%
% we will use a prbs signal on a known system to generate data
% a first order system with g(s)=1/(s+1) is considered
%
load ch15ex3.mat; % first load the data



% now fit some models to this data
%
% create z vector for identification
% the following requires Identification Toolbox in Matlab
%
z=[y1 u1];% display input output data
cla;
hold off;
idplot(z,[],.1);% this program is from identification toolbox
pause;
% now fit the data to a first order arx model
th1=arx(z,[1 1 1]); % arx is an id_toolbox function that implements 
% the least squares fit for the parameters the results are
% stored in a structure variable called th1
% the parameters are na,nb and nk. 

present(th1); % this prsents the model
pause;
% now fit to a [2,2] model
th2=arx(z,[2 2 1]);
% 2resent(th);
% now fit a [3 3] model
th3=arx(z,[3 3 1]);

th4=arx(z,[4 4 1]);

% generate some new data for testing the fit

y1sim=idsim(ut,th1);
y2sim=idsim(ut,th2);
y3sim=idsim(ut,th3);
y4sim=idsim(ut,th4);
pe1=yt-y1sim;pe1norm=sqrt(sum(pe1 .*pe1))
pe2=yt-y2sim;pe2norm=sqrt(sum(pe2 .*pe2))
pe3=yt-y3sim;pe3norm=sqrt(sum(pe3 .*pe3))
pe4=yt-y3sim;pe4norm=sqrt(sum(pe4 .*pe4))
% compare step responses of all models
close;step(sysd);hold on;
[num,den]=th2tf(th1); sysd1=tf(num,den,.1)
step(sysd1);
pause
[num,den]=th2tf(th2);sysd1=tf(num,den,.1)
step(sysd1);
pause
[num,den]=th2tf(th3);sysd1=tf(num,den,.1)
step(sysd1);
pause
[num,den]=th2tf(th4);sysd1=tf(num,den,.1)
step(sysd1);
hold off