0001 function varargout = ltpda_ss_simulate(varargin)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031 ALGONAME = mfilename;
0032 VERSION = '$Id: ltpda_ss_simulate.m,v 1.1 2008/03/11 16:52:56 adrien Exp $';
0033 CATEGORY = 'STATESPACE';
0034 display(['starting ' ALGONAME]);
0035
0036 if not(isempty(varargin))
0037 if isequal( varargin{1}, 'Version')
0038 varargout = VERSION;
0039 return;
0040 elseif isequal(varargin{1}, 'Params')
0041 varargout = plist();
0042 return;
0043 elseif isequal(varargin{1}, 'Category')
0044 varargout = CATEGORY;
0045 return;
0046 end
0047 end
0048 systs = varargin{2};
0049 ao_in_lists = varargin{1};
0050 ao_out = ao();
0051
0052 for i_syst = 1:length(systs)
0053 for i_ao_in = 1:length(ao_in_lists)
0054 syst =systs(i_syst);
0055 ao_in_list = ao_in_lists(i_ao_in);
0056
0057
0058 A = find(syst, 'AMAT');
0059 B = find(syst, 'BMATS');
0060 C = find(syst, 'CMAT');
0061 D = find(syst, 'DMATS');
0062 A = A{1};
0063 C = C{1};
0064 nss = size(A,1);
0065
0066 Xini = find(syst, 'XINI');
0067 timestep = find(syst, 'TIMESTEP');
0068 NBINPUTS = find(syst, 'NBINPUTS');
0069 INPUTISUSED = find(syst, 'INPUTISUSED');
0070 INPUTNAMES = find(syst, 'INPUTNAMES');
0071 INPUTSIZES = find(syst, 'INPUTSIZES');
0072 XISOUTPUT = find(syst, 'XISOUTPUT');
0073 YISOUTPUT = find(syst, 'YISOUTPUT');
0074
0075
0076
0077
0078 ao_in = cell(NBINPUTS,1);
0079 for i=1:NBINPUTS
0080 if INPUTISUSED(i)
0081 ao_in{i} = find(ao_in_list, INPUTNAMES{i});
0082 if isequal(class(ao_in{i}), 'ao')
0083 if not(isequal(ao_in{i}.data.x,{}))
0084 t = ao_in{i}.data.x;
0085 dt = diff(t);
0086 if sum(not(isequal(dt,timestep)))>1
0087 msg = ['error in input number ', i , ' called ', INPUTNAMES{i}, 'because wrong timestep for simulation'];
0088 error(msg);
0089 end
0090 if not( isequal( size(B{i},1), size(D{i},1), size(ao_in{i}.data.y,1) ) )
0091 msg = ['error in input number ', i , ' called ', INPUTNAMES{i}, 'because wrong input size'];
0092 error(msg);
0093 end
0094 end
0095 end
0096 end
0097 end
0098
0099 nsamples = inf;
0100 for i=1:NBINPUTS
0101 if isequal(class(ao_in{i}), 'ao')
0102 nsamples = min(nsamples,size(ao_in{i}.data.y,1));
0103 end
0104 end
0105 if nsamples == inf
0106 nsamples = 0;
0107 end
0108
0109
0110 for i=1:NBINPUTS
0111 if isequal(class(ao_in{i}), 'plist')
0112 sigma = find(ao_in{i}, 'SIGMA');
0113 mean = find(ao_in{i}, 'MEAN');
0114 if sigma==0
0115 ao_in{i} = ao(ones(nsamples,INPUTSIZES(i))*diag(mean));
0116 else
0117 ao_in{i} = ao(ones(nsamples,INPUTSIZES(i))*diag(mean) + randn(nsamples,INPUTSIZES(i))*sigma);
0118 end
0119 end
0120 end
0121
0122
0123 x = zeros(nss, nsamples);
0124 y = zeros(size(C,1), nsamples);
0125 x(:,1) = Xini;
0126 IncrementY = zeros(size(C,1),1);
0127 for i=1:NBINPUTS
0128 if INPUTISUSED(i)
0129 IncrementY = IncrementY + D{i}* ao_in{i}.data.y(1,:)' ;
0130 end
0131 end
0132 y(:,1) = C* x(:,1) +IncrementY;
0133
0134
0135 for k = 2:nsamples
0136 IncrementX = zeros(nss,1);
0137 IncrementY = zeros(size(C,1),1);
0138 for i=1:NBINPUTS
0139 if INPUTISUSED(i)
0140 IncrementX = IncrementX + B{i}*ao_in{i}.data.y(k,:)' ;
0141 IncrementY = IncrementY + D{i}*ao_in{i}.data.y(k,:)' ;
0142 end
0143 end
0144 x(:,k) = A* x(:,k-1) + IncrementX ;
0145 y(:,k) = C* x(:,k) + IncrementY ;
0146 end
0147
0148
0149 out_dat = zeros(0, nsamples);
0150 if XISOUTPUT
0151 out_dat = [out_dat; x];
0152 end
0153 if YISOUTPUT
0154 out_dat = [out_dat; y];
0155 end
0156
0157 if exist('t','var')
0158 ts = tsdata(t,out_dat) ;
0159 else
0160 ts = tsdata(1:nsamples,out_dat) ;
0161 end
0162
0163
0164
0165 hi_new = history('simulate_discrete','version', syst);
0166
0167
0168 h = history(hi_new);
0169 ao_out_loc = ao(ts,h);
0170 ao_out_loc.name = ALGONAME;
0171 ao_out(i_syst , i_ao_in) = ao_out_loc;
0172 end
0173 end
0174
0175 varargout = {ao_out} ;
0176 end