
function [Data,TestData,Syst]=DataGeneration_Sect7_3
%[Data,TestData,Syst]=DataGeneration_Sect7_3
%
% Generating Estimation data and Test data as described in Section 7.3
% Syst: A randomly generated multi-input ARMAX system with polynomials
% of order 30. The number on inputs p is chosen randomly between 2 and 5.
% Certain constraints on the poles of the system and the sizes of the
% channel components are imposed.
%
% Data: Simulated data from this system with the inputs being white noises.
% 600 data points are generated and the last 300 are selected for estimation.
% Noise is added to the output according the the ARMAX noise model. 
%
% TestData: 3000 data points generated from the System with white noisea
% as inputs.  Noise added to the output according to the ARMAX model.
%
 
was=warning;
warning('off')
t=rand(1);
m=ceil(1+4*t);%# of inputs
SV=100;%upper bound on squared ir norms
snr=20;%snr check 
while snr<1 || snr>10
    P=2;%poles check
    V(1)=2*SV;%norm2 check
    while P>0.95 || V(1)>SV 
        [B(1,:),A]=drmodel(30);
        P=max(abs(roots(A)));
        V(1)=sum(dimpulse(B(1,:),A,300).^2);
    end  
    for k=2:m;
        V(k)=2*SV;
        while V(k)>SV;
            [B(k,:),~]=drmodel(30);
            V(k)=sum(dimpulse(B(k,:),A,300).^2);
        end
    end
    P=2;%poles check
    V(m+1)=2*SV;%norm2 check
    while P>0.95 || V(m+1)>SV 
        [~,C]=drmodel(30);
        P=max(abs(roots(C)));
        V(m+1)=sum(dimpulse(C,A,300).^2);
    end
    % update of snr
    snr=sum(V(1:m))/V(m+1);
end
Syst=idpoly(A,[zeros(size(B,1),1) B],C);
Syst.nk=Syst.nk./Syst.nk;%all the delays are set to one
 
%Identification and Test data
u = randn(600,m);
y = sim(Syst,u,'noise');
Data = iddata(y(301:end,1),u(301:end,:));
ut = randn(3300,m);
yt = sim(Syst,ut,'noise');
TestData = iddata(yt(301:end,1),ut(301:end,:));  
warning(was)
        

