
function [Data,TestData_w,TestData_u,Syst]=DataGeneration_Sect7_2
%[Data,TestData_w,TestData_u,Syst]=DataGeneration_Sect7_2
%
% Generating Estimation data and Test data as described in Section 7.2
% Syst: A randomly generated system or order 30. Certain constraints on the
% poles of the system.
%
% Data: Simulated data from this system with an input being white noise
% filtered through a randomly generated 2nd order filter. 1000 data points
% are generated. White noise is added to the output so that the SNR gets a 
% value that is randomly chosen in the range from 1 to 10.
% In the estimation experiment the 200 first are used as
% well as the whole record. In both cases the data is prepended with 200
% zero inputs and outputs to mark that the initial state of the system is
% zero.
%
% TestData_w: 3000 data points generated from the System with a white noise
% input. No noise added to the output.
%
% TestData_u: 3000 data points generated from the System with an input
% being white noise filtered through the same 2nd order filter that
% generated the Data sequence.  No noise added to the output.


was=warning;
warning('off')

P=2;%poles check

while P>0.95  % Sufficiently stable system
    bw=inf;
    while isinf(bw)||isnan(bw) % finite bandwidth of generated system
        mc = rss(30,1,1);
        mc.d=0;
        bw = bandwidth(mc);
    end
    f = bw*3*2*pi;
    md = c2d(mc,1/f,'zoh');
    Syst = idpoly(md);
    %P = max(abs(eig(md.a)));
    P = max(abs(roots(Syst.f)));
    
end
P = 2;
% input generation
while P>0.95
     bw=inf;
    while isinf(bw)||isnan(bw) % finite bandwidth of generated system
        mc = rss(2,1,1);
        
        mc.d=0;
        bw = bandwidth(mc);
    end
    f = bw*3*2*pi;
    mdu = c2d(mc,1/f,'zoh');
    mdu = idpoly(mdu);
    %P = max(abs(eig(md.a)));
    %P = max(abs(roots(Syst.f)));
    %mdu = drss(2,1,1);
    
   
    P = max(abs(roots(mdu.f)));
    
end
ugen = idpoly(mdu);
%Syst = idpoly(md);
we = randn(1000,1);
ue = sim(ugen,we);
wv = randn(1000,1);
uv = sim(ugen,wv);
ut = randn(1000,1);
TestData_w = iddata(sim(Syst,ut),ut);
TestData_u = iddata(sim(Syst,uv),uv);
ynf = sim(Syst,ue); %Noise-free system output
ynorm = ynf'*ynf/1000;
snr = 10*(0.1+0.9*rand(1));
enorm = ynorm/snr;
e = sqrt(enorm)*randn(1000,1);
Data = iddata(ynf+e,ue);
warning(was)


