[Sel, Theta0, ModelVar, ItterVar] = MIMO_ML_DefaultValues(na, nb, nu, ny, PlantPlane, ModelStruct, Recip, Transient); Generates the default values and structures for the estimation of MIMO common denominator plant model structures: see the MIMO_ML function for the definition of the variables Box-Jenkins: Y = B/A U + Ig/A + C/D E + Jh/D Output parameters Sel = selection of the estimated and the known coefficients Sel = struct('A', [],'B', [], 'Ig', []) Sel.A = 1 x (OrderA+1) Sel.A(r) = 1 if coeff. a(r-1) is unknown Sel.A(r) = 0 if coeff. a(r-1) = 0 Sel.B = ny x nu x (OrderB+1) Sel.B(i,j,r) = 1 if coeff. b(i,j,r-1) is unknown Sel.B(i,j,r) = 0 if coeff. b(i,j,r-1) = 0 Sel.Ig = ny x (OrderIg+1) Sel.Ig(i,r) = 1 if coeff. ig(i,r-1) is unknown Sel.Ig(i,r) = 0 if coeff. ig(i,r-1) = 0 Sel.C = ny x ny x (OrderC+1) Theta0 = starting value plant, noise, and initial conditions parameters Theta0 = struct('A', [], 'B', [], 'Ig', []) Theta0.A = 1 x (OrderA+1) Theta0.A(r) = coefficient a(r-1) of Omega^(r-1) Theta0.B = ny x nu x (OrderB+1) Theta0.B(i,j,r) = coefficient b(i,j,r-1) of Omega^(r-1) Theta0.Ig = ny x (OrderIg+1) Theta0.Ig(i,r) = coefficient ig(i,r-1) of Omega^(r-1) ModelVar = contains the information about the model to be identified ModelVar = struct('Transient', [], 'PlantPlane', [], 'Struct', [], 'Reciprocal',[]) ModelVar.Transient = 1 then the initial conditions of the plant and/or noise are estimated ModelVar.PlantPlane = plane of the plant model 's': continuous-time; 'w': sqrt(s)-domain 'z': discrete-time; '': plane not defined ModelVar.Struct = model structure 'EIV': errors-in-variables (noisy input-output data) 'OE': generalised output error (known input, noisy output) ModelVar.RecipPlant = 1 if plant model is reciprocal: G(i,j) = G(j,i) IterVar = contains the information about the minimization procedure IterVar = struct('LM', 1, 'MaxIter', 100, 'TolParam', 1e-6, 'TolCost', 1e-10, ... 'TraceOn', 1, 'NormJacob', 1) IterVar.LM = 1 then then Levenberg-Marquardt minimization scheme is used IterVar.MaxIter = maximum number of itterations of the minimization procedure IterVar.TolParam = relative precision on parameter estimates IterVar.TolCost = relative precision on cost function IterVar.TraceOn = 1 then output iterations IterVar.NormJacob = 1 scales the columns of the Jacobian matrix such that the columns have norm one Input parameters na = order A polynomial nb = order B matrix polynomial nu = number of inputs ny = number of outputs PlantPlane = domain plant model: 'z' for z-domain, 's' for s-domain, 'w' for sqrt(s)-domain (optional, default 'z') ModelStruct = model structure: 'OE', 'EIV' (optional, default EIV) Recip = 1 if plant model reciprocal; otherwise zero (optional, default 0, non-reciprocal plant) Transient = 1 if the plant transient parameters are estimated (optional, default 0, transient parameters) Copyright (c) Rik Pintelon, Vrije Universiteit Brussel - dept. ELEC, November 2009 All rights reserved. Software can be used freely for non-commercial applications only. version 11 October 2011
0001 function [Sel, Theta0, ModelVar, IterVar] = MIMO_ML_DefaultValues(na, nb, nu, ny, PlantPlane, ModelStruct, Recip, Transient); 0002 % 0003 % [Sel, Theta0, ModelVar, ItterVar] = MIMO_ML_DefaultValues(na, nb, nu, ny, PlantPlane, ModelStruct, Recip, Transient); 0004 % 0005 % Generates the default values and structures for the estimation of MIMO common denominator plant model structures: 0006 % see the MIMO_ML function for the definition of the variables 0007 % 0008 % Box-Jenkins: Y = B/A U + Ig/A + C/D E + Jh/D 0009 % 0010 % Output parameters 0011 % 0012 % Sel = selection of the estimated and the known coefficients 0013 % Sel = struct('A', [],'B', [], 'Ig', []) 0014 % Sel.A = 1 x (OrderA+1) 0015 % Sel.A(r) = 1 if coeff. a(r-1) is unknown 0016 % Sel.A(r) = 0 if coeff. a(r-1) = 0 0017 % Sel.B = ny x nu x (OrderB+1) 0018 % Sel.B(i,j,r) = 1 if coeff. b(i,j,r-1) is unknown 0019 % Sel.B(i,j,r) = 0 if coeff. b(i,j,r-1) = 0 0020 % Sel.Ig = ny x (OrderIg+1) 0021 % Sel.Ig(i,r) = 1 if coeff. ig(i,r-1) is unknown 0022 % Sel.Ig(i,r) = 0 if coeff. ig(i,r-1) = 0 0023 % Sel.C = ny x ny x (OrderC+1) 0024 % 0025 % Theta0 = starting value plant, noise, and initial conditions parameters 0026 % Theta0 = struct('A', [], 'B', [], 'Ig', []) 0027 % Theta0.A = 1 x (OrderA+1) 0028 % Theta0.A(r) = coefficient a(r-1) of Omega^(r-1) 0029 % Theta0.B = ny x nu x (OrderB+1) 0030 % Theta0.B(i,j,r) = coefficient b(i,j,r-1) of Omega^(r-1) 0031 % Theta0.Ig = ny x (OrderIg+1) 0032 % Theta0.Ig(i,r) = coefficient ig(i,r-1) of Omega^(r-1) 0033 % 0034 % ModelVar = contains the information about the model to be identified 0035 % ModelVar = struct('Transient', [], 'PlantPlane', [], 'Struct', [], 'Reciprocal',[]) 0036 % ModelVar.Transient = 1 then the initial conditions of the plant and/or noise are estimated 0037 % ModelVar.PlantPlane = plane of the plant model 0038 % 's': continuous-time; 0039 % 'w': sqrt(s)-domain 0040 % 'z': discrete-time; 0041 % '': plane not defined 0042 % ModelVar.Struct = model structure 0043 % 'EIV': errors-in-variables (noisy input-output data) 0044 % 'OE': generalised output error (known input, noisy output) 0045 % ModelVar.RecipPlant = 1 if plant model is reciprocal: G(i,j) = G(j,i) 0046 % 0047 % IterVar = contains the information about the minimization procedure 0048 % IterVar = struct('LM', 1, 'MaxIter', 100, 'TolParam', 1e-6, 'TolCost', 1e-10, ... 0049 % 'TraceOn', 1, 'NormJacob', 1) 0050 % IterVar.LM = 1 then then Levenberg-Marquardt minimization scheme is used 0051 % IterVar.MaxIter = maximum number of itterations of the minimization procedure 0052 % IterVar.TolParam = relative precision on parameter estimates 0053 % IterVar.TolCost = relative precision on cost function 0054 % IterVar.TraceOn = 1 then output iterations 0055 % IterVar.NormJacob = 1 scales the columns of the Jacobian matrix such that the columns have norm one 0056 % 0057 % Input parameters 0058 % 0059 % na = order A polynomial 0060 % nb = order B matrix polynomial 0061 % nu = number of inputs 0062 % ny = number of outputs 0063 % PlantPlane = domain plant model: 'z' for z-domain, 's' for s-domain, 'w' for sqrt(s)-domain 0064 % (optional, default 'z') 0065 % ModelStruct = model structure: 'OE', 'EIV' 0066 % (optional, default EIV) 0067 % Recip = 1 if plant model reciprocal; otherwise zero 0068 % (optional, default 0, non-reciprocal plant) 0069 % Transient = 1 if the plant transient parameters are estimated 0070 % (optional, default 0, transient parameters) 0071 % 0072 % Copyright (c) Rik Pintelon, Vrije Universiteit Brussel - dept. ELEC, November 2009 0073 % All rights reserved. 0074 % Software can be used freely for non-commercial applications only. 0075 % version 11 October 2011 0076 % 0077 0078 0079 %%%%%%%%%%%%%%%%%%%%%%%%%%%% 0080 % initialisation variables % 0081 %%%%%%%%%%%%%%%%%%%%%%%%%%%% 0082 0083 try 0084 if isempty(PlantPlane) 0085 PlantPlane = 'z'; % default discrete-time 0086 end % if 0087 catch 0088 PlantPlane = 'z'; % default discrete-time 0089 end % try 0090 0091 try 0092 if isempty(ModelStruct) 0093 ModelStruct = 'EIV'; % EIV 0094 end % if 0095 catch 0096 ModelStruct = 'EIV'; % EIV 0097 end % try 0098 0099 try 0100 if isempty(Recip) 0101 Recip = 0; % default non-reciprocal 0102 end % if 0103 catch 0104 Recip = 0; % default non-reciprocal 0105 end % try 0106 0107 try 0108 if isempty(Transient) 0109 Transient = 0; % default no transient removal (has been removed non-parametrically 0110 end % if 0111 catch 0112 Transient = 0; % default no transient removal (has been removed non-parametrically 0113 end % try 0114 0115 nig = max(na, nb) - 1; % correct in z-domain; should be increased for s- and sqrt(s)-domains 0116 0117 Sel = struct('A',[], 'B',[], 'Ig', []); 0118 Sel.A = ones(1, na+1); 0119 Sel.B = ones(ny, nu, nb+1); 0120 if (nig < 0) | (Transient == 0) 0121 nig = 0; 0122 Sel.Ig = zeros(ny, nig+1); 0123 else 0124 Sel.Ig = ones(ny, nig+1); 0125 end 0126 0127 Theta0 = struct('A',[], 'B',[], 'Ig', []); 0128 Theta0.A = zeros(1, na+1); 0129 Theta0.A(1) = 1; 0130 Theta0.B = zeros(ny, nu, nb+1); 0131 Theta0.B(:,:,1) = eye(ny, nu); 0132 Theta0.Ig = zeros(ny, nig+1); 0133 0134 ModelVar = struct('Transient', [], 'PlantPlane', [], 'Struct', [], 'RecipPlant', [], ... 0135 'nu', [], 'ny', [], 'na', [], 'nb', [], 'nig', []); 0136 ModelVar.Transient = Transient; 0137 ModelVar.PlantPlane = PlantPlane; 0138 ModelVar.Struct = ModelStruct; 0139 ModelVar.RecipPlant = Recip; 0140 ModelVar.nu = nu; 0141 ModelVar.ny = ny; 0142 ModelVar.na = na; 0143 ModelVar.nb = nb; 0144 ModelVar.nig = nig; 0145 0146 IterVar = struct('LM', [], 'MaxIter', [], 'TolParam', [], 'TolCost', [], 'TraceOn', [], 'NormJacob', []); 0147 IterVar.LM = 1; 0148 IterVar.MaxIter = 100; 0149 IterVar.TolParam = 1e-6; 0150 IterVar.TolCost = 1e-10; 0151 IterVar.TraceOn = 1; 0152 IterVar.NormJacob = 1; 0153 0154 [Theta0, Sel, ModelVar] = MIMO_ML_ModelCompatibility(Theta0, Sel, ModelVar);