sysout.m

来自「用得到的已知输入输出数据」· M 代码 · 共 51 行

function [Y]=SysOut(ParamA,ParamB,ParamC,U,Sigma)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 系统模型A(q)y(t)=B(q)u(t)+C(q)e(t) %
% Y(k)=-ParamA(1)*Y(k-1)-...         %
%      -ParamA(Na)*Y(k-Na)           %
%      +ParamB(1)*U(k-1)+...         %
%      +ParamB(Nb)*U(k-Nb)           %
%      +ParamC(1)*E(k-1)+...         %
%      +ParamC(Nc)*E(k-Nc)           %
%      +E(k)                         %
% U:输入数据                         %
% Sigma:零均值正态分布噪声方差开方   %
% Y:计算输出                         %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Na=size(ParamA,1);
Nb=size(ParamB,1);
Nc=size(ParamC,1);
WholeTime=size(U,1);
Y=[];E=[];
for t=1:WholeTime
    temp1=[];temp2=[];temp3=[];
    for i=1:Na
        if t-i>0
            temp1=[temp1;-Y(t-i)];
        else
            temp1=[temp1;0];
        end
    end
    for i=1:Nb
        if t-i>0
            temp2=[temp2;U(t-i)];
        else
            temp2=[temp2;0];
        end
    end
    for i=1:Nc
        if t-i>0
            temp3=[temp3;E(t-i)];
        else
            temp3=[temp3;0];
        end
    end
    y=ParamA'*temp1+ParamB'*temp2+ParamC'*temp3;
    e=normrnd(0,Sigma);
    y=y+e;
    Y=[Y;y];
    E=[E;e];
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Designed by Ling Qing  2002.01.16  %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%