%LIN_ID Linear system identification.

% {ARS}, 10-20-99

clf;
figure(gcf)
echo on

%    LINEAR SYSTEM IDENTIFICATION:

%    Using NEWLIN/TRAIN or NEWLIND a linear neuron is designed
%    to model a linear system.

pause % Strike any key to continue...

%    DEFINING A WAVE FORM
%    ====================

%    TIME defines the time steps of this simulation.

time = 0:0.025:5;  % from 0 to 6 seconds

%    X defines a systems input signal:

X = sin(sin(time).*time*10);

%    The input to the neural network is the most recent
%    three values of the systems input signal:

%Q=size(X,2);
%P = zeros(3,Q);
%P(1,1:Q) = X(1,1:Q);
%P(2,2:Q) = X(1,1:(Q-1));
%P(3,3:Q) = X(1,1:(Q-2));

P=X;

%    The target T is the output of the original system.

T = filter([1 0.5 -1.5],1,X);

pause % Strike any key to see these signals...

%    PLOTTING THE INPUT SIGNAL
%    =========================

%    Here is a plot of the input signal to the system:

plot(time,X)
axis([0 5 -2 2]);
xlabel('Time');
ylabel('Input Signal');
title('Input Signal to the System');

pause % Strike any key to see these signals...

%    PLOTTING THE OUTPUT SIGNAL
%    ==========================

%    Here is a plot of the output signal of the system:

plot(time,T)
axis([0 5 -2 2]);
xlabel('Time');
ylabel('Output Signal');
title('Output Signal of the System');

pause % Strike any key to design the network...

%    DESIGNING THE NETWORK
%    =====================

%    NEWLIND solves for weights and biases which will let
%    the linear neuron model the system.

net = newlin(minmax(P),1,[0 1 2]);
%net = newlin(minmax(P),1);

net.IW{1}=[1.1 0.5 -1.5];
net.trainParam.epochs=100;
net.trainFcn='trainlm';
%net.trainFcn='trainwb';
net.trainParam.goal=1e-6;
net = train(net,con2seq(P),con2seq(T));

pause % Strike any key to test the neuron model...

%    TESTING THE MODEL
%    =================

%    SIM simulates the linear neuron which models the
%    original system.

pf=[0 0];
[a pf] = sim(net,con2seq(P),con2seq(pf));
a=seq2con(a);
a=a{1};

%    The output signal is plotted with the targets.

plot(time,a,time,T,'+')
axis([0 5 -2 2]);
xlabel('Time');
ylabel('Network Output - System Output +');
title('Network and System Output Signals');

%    The linear neuron does a good job.

pause % Strike any key to see the error signal...

%    THE ERROR SIGNAL
%    ================

%    Error is the difference between output and target signals.

e = T-a;

%    This error can be plotted.

plot(time,e,[min(time) max(time)],[0 0],':r')
axis([0 5 -2 2]);
xlabel('Time');
ylabel('Error');
title('Error Signal');

echo off
disp('End of LIN_ID')
