X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?p=CreaPhase.git;a=blobdiff_plain;f=octave_packages%2Fcontrol-2.3.52%2F%40tf%2F__sys2ss__.m;fp=octave_packages%2Fcontrol-2.3.52%2F%40tf%2F__sys2ss__.m;h=d22dea43806220f82821de3dd2a8487dd371fa36;hp=0000000000000000000000000000000000000000;hb=f5f7a74bd8a4900f0b797da6783be80e11a68d86;hpb=1705066eceaaea976f010f669ce8e972f3734b05 diff --git a/octave_packages/control-2.3.52/@tf/__sys2ss__.m b/octave_packages/control-2.3.52/@tf/__sys2ss__.m new file mode 100644 index 0000000..d22dea4 --- /dev/null +++ b/octave_packages/control-2.3.52/@tf/__sys2ss__.m @@ -0,0 +1,177 @@ +## Copyright (C) 2009, 2011, 2012 Lukas F. Reichlin +## +## This file is part of LTI Syncope. +## +## LTI Syncope is free software: you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation, either version 3 of the License, or +## (at your option) any later version. +## +## LTI Syncope is distributed in the hope that it will be useful, +## but WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +## GNU General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with LTI Syncope. If not, see . + +## -*- texinfo -*- +## TF to SS conversion. +## Reference: +## Varga, A.: Computation of irreducible generalized state-space realizations. +## Kybernetika, 26:89-106, 1990 + +## Special thanks to Vasile Sima and Andras Varga for their advice. +## Author: Lukas Reichlin +## Created: October 2009 +## Version: 0.4.1 + +function [retsys, retlti] = __sys2ss__ (sys) + + ## TODO: determine appropriate tolerance from number of inputs + ## (since we multiply all denominators in a row), index, ... + ## default tolerance from TB01UD is TOLDEF = N*N*EPS + + ## SECRET WISH: a routine which accepts individual denominators for + ## each channel and which supports descriptor systems + + [p, m] = size (sys); + [num, den] = tfdata (sys); + + len_num = cellfun (@length, num); + len_den = cellfun (@length, den); + + ## check for properness + ## tfpoly ensures that there are no leading zeros + tmp = len_num > len_den; + if (any (tmp(:))) # non-proper transfer function + ## separation into strictly proper and polynomial part + [numq, numr] = cellfun (@deconv, num, den, "uniformoutput", false); + numq = cellfun (@__remove_leading_zeros__, numq, "uniformoutput", false); + numr = cellfun (@__remove_leading_zeros__, numr, "uniformoutput", false); + + ## minimal state-space realization for the proper part + [a1, b1, c1] = __proper_tf2ss__ (numr, den, p, m); + e1 = eye (size (a1)); + + ## minimal realization for the polynomial part + [e2, a2, b2, c2] = __polynomial_tf2ss__ (numq, p, m); + + ## assemble irreducible descriptor realization + e = blkdiag (e1, e2); + a = blkdiag (a1, a2); + b = vertcat (b1, b2); + c = horzcat (c1, c2); + retsys = dss (a, b, c, [], e); + else # proper transfer function + [a, b, c, d] = __proper_tf2ss__ (num, den, p, m); + retsys = ss (a, b, c, d); + endif + + retlti = sys.lti; # preserve lti properties such as tsam + +endfunction + + +## transfer function to state-space conversion for proper models +function [a, b, c, d] = __proper_tf2ss__ (num, den, p, m) + + ## new cells for the TF of same row denominators + numc = cell (p, m); + denc = cell (p, 1); + + ## multiply all denominators in a row and + ## update each numerator accordingly + ## except for single-input models and those + ## with equal denominators in a row + for i = 1 : p + if (m == 1 || isequal (den{i,:})) + denc(i) = den{i,1}; + numc(i,:) = num(i,:); + else + denc(i) = __conv__ (den{i,:}); + for j = 1 : m + idx = setdiff (1:m, j); + numc(i,j) = __conv__ (num{i,j}, den{i,idx}); + endfor + endif + endfor + + len_numc = cellfun (@length, numc); + len_denc = cellfun (@length, denc); + + ## check for properness + ## tfpoly ensures that there are no leading zeros + ## tmp = len_numc > repmat (len_denc, 1, m); + ## if (any (tmp(:))) + ## error ("tf: tf2ss: system must be proper"); + ## endif + + ## create arrays and fill in the data + ## in a way that Slicot TD04AD can use + max_len_denc = max (len_denc(:)); + ucoeff = zeros (p, m, max_len_denc); + dcoeff = zeros (p, max_len_denc); + index = len_denc-1; + + for i = 1 : p + len = len_denc(i); + dcoeff(i, 1:len) = denc{i}; + for j = 1 : m + ucoeff(i, j, len-len_numc(i,j)+1 : len) = numc{i,j}; + endfor + endfor + + tol = min (sqrt (eps), eps*prod (index)); + [a, b, c, d] = sltd04ad (ucoeff, dcoeff, index, tol); + +endfunction + + +## realization of the polynomial part according to Andras' paper +function [e2, a2, b2, c2] = __polynomial_tf2ss__ (numq, p, m) + + len_numq = cellfun (@length, numq); + max_len_numq = max (len_numq(:)); + numq = cellfun (@(x) prepad (x, max_len_numq, 0, 2), numq, "uniformoutput", false); + f = @(y) cellfun (@(x) x(y), numq); + s = 1 : max_len_numq; + D = arrayfun (f, s, "uniformoutput", false); + + e2 = diag (ones (p*(max_len_numq-1), 1), -p); + a2 = eye (p*max_len_numq); + b2 = vertcat (D{:}); + c2 = horzcat (zeros (p, p*(max_len_numq-1)), -eye (p)); + + ## remove uncontrollable part + [a2, e2, b2, c2] = sltg01jd (a2, e2, b2, c2, 0.0, true, 1, 2); + +endfunction + + +## convolution for more than two arguments +function vec = __conv__ (vec, varargin) + + if (nargin == 1) + return; + else + for k = 1 : nargin-1 + vec = conv (vec, varargin{k}); + endfor + endif + +endfunction + + +## remove leading zeros from polynomial vector +function p = __remove_leading_zeros__ (p) + + idx = find (p != 0); + + if (isempty (idx)) + p = 0; + else + p = p(idx(1) : end); # p(idx) would remove all zeros + endif + +endfunction \ No newline at end of file