From: <par...@us...> - 2012-06-21 07:08:39
|
Revision: 10648 http://octave.svn.sourceforge.net/octave/?rev=10648&view=rev Author: paramaniac Date: 2012-06-21 07:08:29 +0000 (Thu, 21 Jun 2012) Log Message: ----------- control: apply changes by Megan Zagrobelny Modified Paths: -------------- trunk/octave-forge/main/control/devel/dlqe.m trunk/octave-forge/main/control/devel/lqe.m Modified: trunk/octave-forge/main/control/devel/dlqe.m =================================================================== --- trunk/octave-forge/main/control/devel/dlqe.m 2012-06-20 21:48:29 UTC (rev 10647) +++ trunk/octave-forge/main/control/devel/dlqe.m 2012-06-21 07:08:29 UTC (rev 10648) @@ -23,6 +23,14 @@ ## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{z}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}, @var{s}) ## Kalman filter for discrete-time systems. ## +## @example +## @group +## x[k] = Ax[k] + Bu[k] + Gw[k] (State equation) +## y[k] = Cx[k] + Du[k] + v[k] (Measurement Equation) +## E(w) = 0, E(v) = 0, cov(w) = Q, cov(v) = R, cov(w,v) = S +## @end group +## @end example +## ## @strong{Inputs} ## @table @var ## @item sys @@ -30,15 +38,15 @@ ## @item a ## State transition matrix of discrete-time system (n-by-n). ## @item g -## Process noise matrix of discrete-time system (n-by-.). +## Process noise matrix of discrete-time system (n-by-g). ## @item c ## Measurement matrix of discrete-time system (p-by-n). ## @item q -## Process noise covariance matrix (.-by-.). +## Process noise covariance matrix (g-by-g). ## @item r ## Measurement noise covariance matrix (p-by-p). ## @item s -## Optional cross term covariance matrix (n-by-p), s = cov(w,v) If @var{s} is not specified, a zero matrix is assumed. +## Optional cross term covariance matrix (g-by-p), s = cov(w,v) If @var{s} is not specified, a zero matrix is assumed. ## @end table ## ## @strong{Outputs} @@ -89,14 +97,9 @@ else [~, p, e] = dlqr (a.', c.', g*q*g.', r, g*s); endif - - ## k computed by dlqr would be - ## k = (r + c*p*c.') \ (c*p*a.' + s.') - ## such that l = a \ k.' - ## what about the s term? l = p*c.' / (c*p*c.' + r); - ## z = p - p*c.' / (c*p*c.' + r) * c*p; + z = p - l*c*p; endfunction Modified: trunk/octave-forge/main/control/devel/lqe.m =================================================================== --- trunk/octave-forge/main/control/devel/lqe.m 2012-06-20 21:48:29 UTC (rev 10647) +++ trunk/octave-forge/main/control/devel/lqe.m 2012-06-21 07:08:29 UTC (rev 10648) @@ -1,4 +1,5 @@ ## Copyright (C) 2012 Lukas F. Reichlin +## Copyright (C) 2012 Megan Zagrobelny ## ## This file is part of LTI Syncope. ## @@ -22,50 +23,56 @@ ## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}, @var{s}) ## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}) ## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}, @var{s}) -## Linear-quadratic estimator. +## Kalman filter for continuous-time systems. ## +## @example +## @group +## . +## x = Ax + Bu + Gw (State equation) +## y = Cx + Du + v (Measurement Equation) +## E(w) = 0, E(v) = 0, cov(w) = Q, cov(v) = R, cov(w,v) = S +## @end group +## @end example +## ## @strong{Inputs} ## @table @var ## @item sys -## Continuous or discrete-time LTI model. +## Continuous or discrete-time LTI model (p-by-m, n states). ## @item a -## State transition matrix of continuous-time system. -## @item b -## Input matrix of continuous-time system. +## State transition matrix of continuous-time system (n-by-n). +## @item g +## Process noise matrix of continuous-time system (n-by-g). +## @item c +## Measurement matrix of continuous-time system (p-by-n). ## @item q -## State weighting matrix. +## Process noise covariance matrix (g-by-g). ## @item r -## Input weighting matrix. +## Measurement noise covariance matrix (p-by-p). ## @item s -## Optional cross term matrix. If @var{s} is not specified, a zero matrix is assumed. -## @item e -## Optional descriptor matrix. If @var{e} is not specified, an identity matrix is assumed. +## Optional cross term covariance matrix (g-by-p), s = cov(w,v) If @var{s} is not specified, a zero matrix is assumed. ## @end table ## ## @strong{Outputs} ## @table @var ## @item l -## Observer gain matrix. +## Kalman filter gain matrix (n-by-p). ## @item p -## Unique stabilizing solution of the continuous-time Riccati equation. +## Unique stabilizing solution of the continuous-time Riccati equation (n-by-n). ## @item e -## Closed-loop poles. +## Closed-loop poles (n-by-1). ## @end table ## ## @strong{Equations} ## @example ## @group ## . -## x = A x + B u, x(0) = x0 +## x = Ax + Bu + L(y - Cx -Du) ## -## inf -## J(x0) = INT (x' Q x + u' R u + 2 x' S u) dt -## 0 -## -## L = eig (A - B*G) +## E = eig(A - L*C) +## ## @end group ## @end example -## @seealso{lqr, care, dare, dlqe} +## @seealso{dare, care, dlqr, lqr, dlqe} ## @end deftypefn ## Author: Lukas Reichlin <luk...@gm...> This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-07-03 12:45:09
|
Revision: 10717 http://octave.svn.sourceforge.net/octave/?rev=10717&view=rev Author: paramaniac Date: 2012-07-03 12:45:00 +0000 (Tue, 03 Jul 2012) Log Message: ----------- control: remove custom lapack Removed Paths: ------------- trunk/octave-forge/main/control/devel/Makefile.ReferenceLapack trunk/octave-forge/main/control/devel/lapack-3.4.1.tgz Deleted: trunk/octave-forge/main/control/devel/Makefile.ReferenceLapack =================================================================== --- trunk/octave-forge/main/control/devel/Makefile.ReferenceLapack 2012-07-03 12:39:16 UTC (rev 10716) +++ trunk/octave-forge/main/control/devel/Makefile.ReferenceLapack 2012-07-03 12:45:00 UTC (rev 10717) @@ -1,77 +0,0 @@ -MKOCTFILE ?= mkoctfile - -ifndef LAPACK_LIBS -LAPACK_LIBS := $(shell $(MKOCTFILE) -p LAPACK_LIBS) -endif -ifndef BLAS_LIBS -BLAS_LIBS := $(shell $(MKOCTFILE) -p BLAS_LIBS) -endif -ifndef FLIBS -FLIBS := $(shell $(MKOCTFILE) -p FLIBS) -endif -MKOCTFILE ?= mkoctfile - -# LAPACK_LIBS := $(shell $(MKOCTFILE) -p LAPACK_LIBS) -# BLAS_LIBS := $(shell $(MKOCTFILE) -p BLAS_LIBS) -FLIBS := $(shell $(MKOCTFILE) -p FLIBS) - -all: control_slicot_functions.oct \ - is_real_scalar.oct \ - is_real_vector.oct \ - is_real_matrix.oct \ - is_real_square_matrix.oct - -# TODO: Private oct-files for control package. - -# unpack and compile SLICOT library -# Note that TG04BX is a custom routine. -# It has the extension .fortran such that -# it is not deleted by rm *.f when using -# the developer makefile makefile_control.m -slicotlibrary.a: slicot.tar.gz - tar -xzf slicot.tar.gz - mv slicot/src/*.f . - mv slicot/src_aux/*.f . - cp TG04BX.fortran TG04BX.f - $(MKOCTFILE) -c *.f - ar -rc slicotlibrary.a *.o - rm -rf *.o *.f slicot - -# reference lapack is just included for debugging purposes -# it will be removed again before an official release of the control package -lapacklibrary.a: lapack-3.4.1.tgz - tar -xzf lapack-3.4.1.tgz - mv lapack-3.4.1/BLAS/SRC/*.f . - mv lapack-3.4.1/INSTALL/*.f . - mv lapack-3.4.1/SRC/*.f . - $(MKOCTFILE) -c *.f - ar -rc lapacklibrary.a *.o - rm -rf *.o *.f lapack-3.4.1 - -# slicot functions -control_slicot_functions.oct: control_slicot_functions.cc slicotlibrary.a lapacklibrary.a - $(MKOCTFILE) control_slicot_functions.cc common.cc slicotlibrary.a lapacklibrary.a \ - ${FLIBS} -# slicot functions -#control_slicot_functions.oct: control_slicot_functions.cc slicotlibrary.a -# $(MKOCTFILE) control_slicot_functions.cc common.cc slicotlibrary.a \ -# ${LAPACK_LIBS} ${BLAS_LIBS} ${FLIBS} - -# helpers -is_real_scalar.oct: is_real_scalar.cc - $(MKOCTFILE) is_real_scalar.cc - -is_real_vector.oct: is_real_vector.cc - $(MKOCTFILE) is_real_vector.cc - -is_real_matrix.oct: is_real_matrix.cc - $(MKOCTFILE) is_real_matrix.cc - -is_real_square_matrix.oct: is_real_square_matrix.cc - $(MKOCTFILE) is_real_square_matrix.cc - -clean: - rm -rf *.o core octave-core *.oct *~ *.f slicot lapack-3.4.1 - -realclean: clean - rm -rf *.a \ No newline at end of file Deleted: trunk/octave-forge/main/control/devel/lapack-3.4.1.tgz =================================================================== (Binary files differ) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-12 11:20:25
|
Revision: 10998 http://octave.svn.sourceforge.net/octave/?rev=10998&view=rev Author: paramaniac Date: 2012-09-12 11:20:19 +0000 (Wed, 12 Sep 2012) Log Message: ----------- control: sync draft code Modified Paths: -------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/bode2.m Modified: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-12 10:39:10 UTC (rev 10997) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-12 11:20:19 UTC (rev 10998) @@ -1,4 +1,4 @@ -## Copyright (C) 2009, 2010 Lukas F. Reichlin +## Copyright (C) 2009, 2010, 2012 Lukas F. Reichlin ## ## This file is part of LTI Syncope. ## @@ -21,9 +21,9 @@ ## Author: Lukas Reichlin <luk...@gm...> ## Created: November 2009 -## Version: 0.2 +## Version: 0.4 -% function [H, w] = __frequency_response_2__ (sys, w = [], mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) +% function [H, w] = __frequency_response__ (sys, w = [], mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) function [H, w] = __frequency_response_2__ (mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false, varargin) sys_vec = cellfun (@(x) isa (x, "lti"), varargin) # true or false @@ -40,6 +40,7 @@ varargin{sys_idx} + ## check arguments if(! isa (sys, "lti")) error ("frequency_response: first argument sys must be an LTI system"); @@ -50,18 +51,21 @@ endif if (isa (sys, "frd")) - if (isempty (w)) - w = get (sys, "w"); - else + if (! isempty (w)) warning ("frequency_response: second argument w is ignored"); endif + w = get (sys, "w"); + H = __freqresp__ (sys, [], resptype, cellflag); elseif (isempty (w)) # find interesting frequency range w if not specified w = __frequency_vector__ (sys, wbounds); + H = __freqresp__ (sys, w, resptype, cellflag); + elseif (iscell (w) && numel (w) == 2 && issample (w{1}) && issample (w{2})) + w = __frequency_vector__ (sys, wbounds, w{1}, w{2}); + H = __freqresp__ (sys, w, resptype, cellflag); elseif (! is_real_vector (w)) error ("frequency_response: second argument w must be a vector of frequencies"); + else + H = __freqresp__ (sys, w, resptype, cellflag); endif - ## frequency response - H = __freqresp__ (sys, w, resptype, cellflag); - -endfunction \ No newline at end of file +endfunction Modified: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-12 10:39:10 UTC (rev 10997) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-12 11:20:19 UTC (rev 10998) @@ -1,4 +1,4 @@ -## Copyright (C) 2009, 2010 Lukas F. Reichlin +## Copyright (C) 2009, 2010, 2011, 2012 Lukas F. Reichlin ## ## This file is part of LTI Syncope. ## @@ -28,6 +28,9 @@ ## @item w ## Optional vector of frequency values. If @var{w} is not specified, ## it is calculated by the zeros and poles of the system. +## Alternatively, the cell @code{@{wmin, wmax@}} specifies a frequency range, +## where @var{wmin} and @var{wmax} denote minimum and maximum frequencies +## in rad/s. ## @end table ## ## @strong{Outputs} @@ -45,7 +48,7 @@ ## Author: Lukas Reichlin <luk...@gm...> ## Created: November 2009 -## Version: 0.2 +## Version: 0.4 function [mag_r, pha_r, w_r] = bode2 (varargin) @@ -64,12 +67,6 @@ if (! nargout) mag_db = 20 * log10 (mag); - wv = [min(w), max(w)]; - ax_vec_mag = __axis_limits__ ([reshape(w, [], 1), reshape(mag_db, [], 1)]); - ax_vec_mag(1:2) = wv; - ax_vec_pha = __axis_limits__ ([reshape(w, [], 1), reshape(pha, [], 1)]); - ax_vec_pha(1:2) = wv; - if (isct (sys)) xl_str = "Frequency [rad/s]"; else @@ -78,14 +75,16 @@ subplot (2, 1, 1) semilogx (w, mag_db) - axis (ax_vec_mag) + axis ("tight") + ylim (__axis_margin__ (ylim)) grid ("on") title (["Bode Diagram of ", inputname(1)]) ylabel ("Magnitude [dB]") subplot (2, 1, 2) semilogx (w, pha) - axis (ax_vec_pha) + axis ("tight") + ylim (__axis_margin__ (ylim)) grid ("on") xlabel (xl_str) ylabel ("Phase [deg]") This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-12 14:42:28
|
Revision: 10999 http://octave.svn.sourceforge.net/octave/?rev=10999&view=rev Author: paramaniac Date: 2012-09-12 14:42:16 +0000 (Wed, 12 Sep 2012) Log Message: ----------- control: get multiplot draft code working \ o / Modified Paths: -------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/__frequency_vector_2__.m trunk/octave-forge/main/control/devel/bode2.m Added Paths: ----------- trunk/octave-forge/main/control/devel/multiplot.m trunk/octave-forge/main/control/devel/tfs.dat Modified: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-12 11:20:19 UTC (rev 10998) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-12 14:42:16 UTC (rev 10999) @@ -26,21 +26,24 @@ % function [H, w] = __frequency_response__ (sys, w = [], mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) function [H, w] = __frequency_response_2__ (mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false, varargin) - sys_vec = cellfun (@(x) isa (x, "lti"), varargin) # true or false - sys_idx = find (sys_vec) + sys_idx = cellfun (@isa, varargin, {"lti"}); # true or false + w_idx = cellfun (@is_real_vector, varargin); # look for frequency vectors + % ? = cellfun (@iscell, varargin) + % varargin(?) (end) - w_vec = cellfun (@is_real_vector, varargin) - w_idx = find (w_vec) % w_idx(end) - if (isempty (w_idx)) - w = __frequency_vector_2__ (wbounds, varargin{sys_idx}) + if (! any (w_idx)) + w = __frequency_vector_2__ (varargin(sys_idx), wbounds); endif -varargin{sys_idx} +%varargin{sys_idx} + H = cellfun (@__freqresp__, varargin(sys_idx), {w}, {resptype}, {cellflag}, "uniformoutput", false); + +%{ ## check arguments if(! isa (sys, "lti")) error ("frequency_response: first argument sys must be an LTI system"); @@ -67,5 +70,7 @@ else H = __freqresp__ (sys, w, resptype, cellflag); endif + +%} endfunction Modified: trunk/octave-forge/main/control/devel/__frequency_vector_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-12 11:20:19 UTC (rev 10998) +++ trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-12 14:42:16 UTC (rev 10999) @@ -37,7 +37,7 @@ ## Date: October 2009 ## Version: 0.2 -function w = __frequency_vector__ (sys_cell, wbounds = "std") +function w = __frequency_vector_2__ (sys_cell, wbounds = "std") if (! iscell (sys_cell)) sys_cell = {sys_cell} Modified: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-12 11:20:19 UTC (rev 10998) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-12 14:42:16 UTC (rev 10999) @@ -59,22 +59,30 @@ % endif [H, w] = __frequency_response_2__ (false, 0, "std", false, varargin{:}); + + H = cellfun (@reshape, H, {[]}, {1}, "uniformoutput", false); + mag = cellfun (@abs, H, "uniformoutput", false); + pha = cellfun (@(H) unwrap (arg (H)) * 180 / pi, H, "uniformoutput", false); - H = reshape (H, [], 1); - mag = abs (H); - pha = unwrap (arg (H)) * 180 / pi; - if (! nargout) - mag_db = 20 * log10 (mag); + mag_db = cellfun (@(mag) 20 * log10 (mag), mag, "uniformoutput", false); - if (isct (sys)) + w_cell = repmat ({w}, 1, numel (H)); + + mag_args = vertcat (w_cell, mag_db)(:); + pha_args = vertcat (w_cell, pha)(:); + + %mag_args = cellfun (@horzcat, {w}, mag_db, "uniformoutput", false); + %pha_args = cellfun (@horzcat, {w}, pha, "uniformoutput", false); + + %if (isct (sys)) xl_str = "Frequency [rad/s]"; - else - xl_str = sprintf ("Frequency [rad/s] w_N = %g", pi / get (sys, "tsam")); - endif + %else + % xl_str = sprintf ("Frequency [rad/s] w_N = %g", pi / get (sys, "tsam")); + %endif subplot (2, 1, 1) - semilogx (w, mag_db) + semilogx (mag_args{:}) axis ("tight") ylim (__axis_margin__ (ylim)) grid ("on") @@ -82,15 +90,15 @@ ylabel ("Magnitude [dB]") subplot (2, 1, 2) - semilogx (w, pha) + semilogx (pha_args{:}) axis ("tight") ylim (__axis_margin__ (ylim)) grid ("on") xlabel (xl_str) ylabel ("Phase [deg]") else - mag_r = mag; - pha_r = pha; + mag_r = mag{1}; + pha_r = pha{1}; w_r = w; endif Added: trunk/octave-forge/main/control/devel/multiplot.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot.m (rev 0) +++ trunk/octave-forge/main/control/devel/multiplot.m 2012-09-12 14:42:16 UTC (rev 10999) @@ -0,0 +1,3 @@ +load tfs.dat + +bode2 (C_AH, C_opt) \ No newline at end of file Added: trunk/octave-forge/main/control/devel/tfs.dat =================================================================== --- trunk/octave-forge/main/control/devel/tfs.dat (rev 0) +++ trunk/octave-forge/main/control/devel/tfs.dat 2012-09-12 14:42:16 UTC (rev 10999) @@ -0,0 +1,3990 @@ +# Created by Octave 3.6.3, Wed Sep 12 16:18:22 2012 CEST <lukas@ra.local> +# name: C_AH +# type: class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 4 + 0.008979867577393098 0.2861205263052757 2.279124799660292 0 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 3 + 1.211107551841654 1.929442317265359 0.846571595181162 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: C_opt +# type: class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 4 + 0.0277727961257141 0.5755081966283047 2.981421846103795 0 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 3 + 1.742425436966594 1.805327984349914 0.605525845565654 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: C_par_0 +# type: matrix +# rows: 3 +# columns: 1 + 0.846571595181162 + 2.279124799660292 + 0.6276982426498158 + + +# name: C_par_opt +# type: matrix +# rows: 3 +# columns: 1 + 0.605525845565654 + 2.981421846103795 + 0.9651572745071195 + + +# name: L_AH +# type: class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 10 + 0.008979867577393098 0.3310198641922412 3.808505974537995 14.668667933743 29.17483870788839 35.09997232248909 26.50995529536698 11.68174452460673 2.279124799660292 0 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 3 + 1.211107551841654 1.929442317265359 0.846571595181162 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: L_opt +# type: class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 10 + 0.0277727961257141 0.7143721772568752 6.164463586628174 21.62651853919032 41.15825581732086 48.20935998899305 35.70095408640898 15.48261742714728 2.981421846103795 0 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 3 + 1.742425436966594 1.805327984349914 0.605525845565654 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: P +# type: global class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 7 + 1 5 11 14 11 5 1 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +1 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: T_AH +# type: class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 10 + 0.008979867577393098 0.3310198641922412 3.808505974537995 14.668667933743 29.17483870788839 35.09997232248909 26.50995529536698 12.89285207644839 4.208567116925651 0.846571595181162 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 3 + 1.211107551841654 1.929442317265359 0.846571595181162 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: T_opt +# type: class +# classname: tf +# length: 5 +# name: den +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 10 + 0.0277727961257141 0.7143721772568752 6.164463586628174 21.62651853919032 41.15825581732086 48.20935998899305 35.70095408640898 17.22504286411387 4.786749830453709 0.605525845565654 + + + + + + + + + + + + + +# name: inv +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: bool +0 + + + + + +# name: lti +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: lti +# length: 6 +# name: inname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: name +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + +# name: notes +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 0 +# columns: 0 + + + + + +# name: outname +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_string +# elements: 0 + + + + + + + + +# name: tsam +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: scalar +0 + + + + + +# name: userdata +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: null_matrix +# rows: 0 +# columns: 0 + + + + + + + + + + +# name: num +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: class +# classname: tfpoly +# length: 1 +# name: poly +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: matrix +# rows: 1 +# columns: 3 + 1.742425436966594 1.805327984349914 0.605525845565654 + + + + + + + + + + + + + +# name: tfvar +# type: cell +# rows: 1 +# columns: 1 +# name: <cell-element> +# type: string +# elements: 1 +# length: 1 +s + + + + + + + +# name: Td_AH +# type: scalar +0.6276982426498158 + + +# name: Td_opt +# type: scalar +0.9651572745071195 + + +# name: Ti_AH +# type: scalar +2.279124799660292 + + +# name: Ti_opt +# type: scalar +2.981421846103795 + + +# name: Tu +# type: scalar +9.693196537717052 + + +# name: denP +# type: matrix +# rows: 1 +# columns: 7 + 1 5 11 14 11 5 1 + + +# name: dt +# type: global scalar +0.05 + + +# name: gamma +# type: scalar +1.75407835279559 + + +# name: gamma_AH +# type: scalar +1.587864585936343 + + +# name: gamma_opt +# type: scalar +2.867548815608383 + + +# name: kappa +# type: scalar +1 + + +# name: kp_AH +# type: scalar +0.846571595181162 + + +# name: kp_opt +# type: scalar +0.605525845565654 + + +# name: ku +# type: scalar +1.75407835279559 + + +# name: mu_1 +# type: global scalar +1 + + +# name: mu_2 +# type: global scalar +10 + + +# name: mu_3 +# type: global scalar +20 + + +# name: numP +# type: scalar +1 + + +# name: phi +# type: scalar +180 + + +# name: phi_AH +# type: scalar +31.39693841431557 + + +# name: phi_opt +# type: scalar +64.20796755745674 + + +# name: st_AH +# type: bool +1 + + +# name: st_opt +# type: bool +1 + + +# name: t +# type: global range +# base, limit, increment +0 30 0.05 + + +# name: t_AH +# type: matrix +# rows: 601 +# columns: 1 + 0 + 0.05 + 0.1 + 0.15 + 0.2 + 0.25 + 0.3 + 0.35 + 0.4 + 0.45 + 0.5 + 0.55 + 0.6000000000000001 + 0.65 + 0.7000000000000001 + 0.75 + 0.8 + 0.8500000000000001 + 0.9 + 0.9500000000000001 + 1 + 1.05 + 1.1 + 1.15 + 1.2 + 1.25 + 1.3 + 1.35 + 1.4 + 1.45 + 1.5 + 1.55 + 1.6 + 1.65 + 1.7 + 1.75 + 1.8 + 1.85 + 1.9 + 1.95 + 2 + 2.05 + 2.1 + 2.15 + 2.2 + 2.25 + 2.3 + 2.35 + 2.4 + 2.45 + 2.5 + 2.55 + 2.6 + 2.65 + 2.7 + 2.75 + 2.8 + 2.85 + 2.9 + 2.95 + 3 + 3.05 + 3.1 + 3.15 + 3.2 + 3.25 + 3.3 + 3.35 + 3.4 + 3.45 + 3.5 + 3.55 + 3.6 + 3.65 + 3.7 + 3.75 + 3.8 + 3.85 + 3.9 + 3.95 + 4 + 4.05 + 4.100000000000001 + 4.15 + 4.2 + 4.25 + 4.3 + 4.350000000000001 + 4.4 + 4.45 + 4.5 + 4.55 + 4.600000000000001 + 4.65 + 4.7 + 4.75 + 4.800000000000001 + 4.850000000000001 + 4.9 + 4.95 + 5 + 5.050000000000001 + 5.100000000000001 + 5.15 + 5.2 + 5.25 + 5.300000000000001 + 5.350000000000001 + 5.4 + 5.45 + 5.5 + 5.550000000000001 + 5.600000000000001 + 5.65 + 5.7 + 5.75 + 5.800000000000001 + 5.850000000000001 + 5.9 + 5.95 + 6 + 6.050000000000001 + 6.100000000000001 + 6.15 + 6.2 + 6.25 + 6.300000000000001 + 6.350000000000001 + 6.4 + 6.45 + 6.5 + 6.550000000000001 + 6.600000000000001 + 6.65 + 6.7 + 6.75 + 6.800000000000001 + 6.850000000000001 + 6.9 + 6.95 + 7 + 7.050000000000001 + 7.100000000000001 + 7.15 + 7.2 + 7.25 + 7.300000000000001 + 7.350000000000001 + 7.4 + 7.45 + 7.5 + 7.550000000000001 + 7.600000000000001 + 7.65 + 7.7 + 7.75 + 7.800000000000001 + 7.850000000000001 + 7.9 + 7.95 + 8 + 8.050000000000001 + 8.1 + 8.15 + 8.200000000000001 + 8.25 + 8.300000000000001 + 8.35 + 8.4 + 8.450000000000001 + 8.5 + 8.550000000000001 + 8.6 + 8.65 + 8.700000000000001 + 8.75 + 8.800000000000001 + 8.85 + 8.9 + 8.950000000000001 + 9 + 9.050000000000001 + 9.1 + 9.15 + 9.200000000000001 + 9.25 + 9.300000000000001 + 9.35 + 9.4 + 9.450000000000001 + 9.5 + 9.550000000000001 + 9.600000000000001 + 9.65 + 9.700000000000001 + 9.75 + 9.800000000000001 + 9.850000000000001 + 9.9 + 9.950000000000001 + 10 + 10.05 + 10.1 + 10.15 + 10.2 + 10.25 + 10.3 + 10.35 + 10.4 + 10.45 + 10.5 + 10.55 + 10.6 + 10.65 + 10.7 + 10.75 + 10.8 + 10.85 + 10.9 + 10.95 + 11 + 11.05 + 11.1 + 11.15 + 11.2 + 11.25 + 11.3 + 11.35 + 11.4 + 11.45 + 11.5 + 11.55 + 11.6 + 11.65 + 11.7 + 11.75 + 11.8 + 11.85 + 11.9 + 11.95 + 12 + 12.05 + 12.1 + 12.15 + 12.2 + 12.25 + 12.3 + 12.35 + 12.4 + 12.45 + 12.5 + 12.55 + 12.6 + 12.65 + 12.7 + 12.75 + 12.8 + 12.85 + 12.9 + 12.95 + 13 + 13.05 + 13.1 + 13.15 + 13.2 + 13.25 + 13.3 + 13.35 + 13.4 + 13.45 + 13.5 + 13.55 + 13.6 + 13.65 + 13.7 + 13.75 + 13.8 + 13.85 + 13.9 + 13.95 + 14 + 14.05 + 14.1 + 14.15 + 14.2 + 14.25 + 14.3 + 14.35 + 14.4 + 14.45 + 14.5 + 14.55 + 14.6 + 14.65 + 14.7 + 14.75 + 14.8 + 14.85 + 14.9 + 14.95 + 15 + 15.05 + 15.1 + 15.15 + 15.2 + 15.25 + 15.3 + 15.35 + 15.4 + 15.45 + 15.5 + 15.55 + 15.6 + 15.65 + 15.7 + 15.75 + 15.8 + 15.85 + 15.9 + 15.95 + 16 + 16.05 + 16.1 + 16.15 + 16.2 + 16.25 + 16.3 + 16.35 + 16.4 + 16.45 + 16.5 + 16.55 + 16.6 + 16.65 + 16.7 + 16.75 + 16.8 + 16.85 + 16.9 + 16.95 + 17 + 17.05 + 17.1 + 17.15 + 17.2 + 17.25 + 17.3 + 17.35 + 17.4 + 17.45 + 17.5 + 17.55 + 17.6 + 17.65 + 17.7 + 17.75 + 17.8 + 17.85 + 17.9 + 17.95 + 18 + 18.05 + 18.1 + 18.15 + 18.2 + 18.25 + 18.3 + 18.35 + 18.4 + 18.45 + 18.5 + 18.55 + 18.6 + 18.65 + 18.7 + 18.75 + 18.8 + 18.85 + 18.9 + 18.95 + 19 + 19.05 + 19.1 + 19.15 + 19.2 + 19.25 + 19.3 + 19.35 + 19.4 + 19.45 + 19.5 + 19.55 + 19.6 + 19.65 + 19.7 + 19.75 + 19.8 + 19.85 + 19.9 + 19.95 + 20 + 20.05 + 20.1 + 20.15 + 20.2 + 20.25 + 20.3 + 20.35 + 20.4 + 20.45 + 20.5 + 20.55 + 20.6 + 20.65 + 20.7 + 20.75 + 20.8 + 20.85 + 20.9 + 20.95 + 21 + 21.05 + 21.1 + 21.15 + 21.2 + 21.25 + 21.3 + 21.35 + 21.4 + 21.45 + 21.5 + 21.55 + 21.6 + 21.65 + 21.7 + 21.75 + 21.8 + 21.85 + 21.9 + 21.95 + 22 + 22.05 + 22.1 + 22.15 + 22.2 + 22.25 + 22.3 + 22.35 + 22.4 + 22.45 + 22.5 + 22.55 + 22.6 + 22.65 + 22.7 + 22.75 + 22.8 + 22.85 + 22.9 + 22.95 + 23 + 23.05 + 23.1 + 23.15 + 23.2 + 23.25 + 23.3 + 23.35 + 23.4 + 23.45 + 23.5 + 23.55 + 23.6 + 23.65 + 23.7 + 23.75 + 23.8 + 23.85 + 23.9 + 23.95 + 24 + 24.05 + 24.1 + 24.15 + 24.2 + 24.25 + 24.3 + 24.35 + 24.4 + 24.45 + 24.5 + 24.55 + 24.6 + 24.65 + 24.7 + 24.75 + 24.8 + 24.85 + 24.9 + 24.95 + 25 + 25.05 + 25.1 + 25.15 + 25.2 + 25.25 + 25.3 + 25.35 + 25.4 + 25.45 + 25.5 + 25.55 + 25.6 + 25.65 + 25.7 + 25.75 + 25.8 + 25.85 + 25.9 + 25.95 + 26 + 26.05 + 26.1 + 26.15 + 26.2 + 26.25 + 26.3 + 26.35 + 26.4 + 26.45 + 26.5 + 26.55 + 26.6 + 26.65 + 26.7 + 26.75 + 26.8 + 26.85 + 26.9 + 26.95 + 27 + 27.05 + 27.1 + 27.15 + 27.2 + 27.25 + 27.3 + 27.35 + 27.4 + 27.45 + 27.5 + 27.55 + 27.6 + 27.65 + 27.7 + 27.75 + 27.8 + 27.85 + 27.9 + 27.95 + 28 + 28.05 + 28.1 + 28.15 + 28.2 + 28.25 + 28.3 + 28.35 + 28.4 + 28.45 + 28.5 + 28.55 + 28.6 + 28.65 + 28.7 + 28.75 + 28.8 + 28.85 + 28.9 + 28.95 + 29 + 29.05 + 29.1 + 29.15 + 29.2 + 29.25 + 29.3 + 29.35 + 29.4 + 29.45 + 29.5 + 29.55 + 29.6 + 29.65 + 29.7 + 29.75 + 29.8 + 29.85 + 29.9 + 29.95 + 30 + + +# name: t_opt +# type: matrix +# rows: 601 +# columns: 1 + 0 + 0.05 + 0.1 + 0.15 + 0.2 + 0.25 + 0.3 + 0.35 + 0.4 + 0.45 + 0.5 + 0.55 + 0.6000000000000001 + 0.65 + 0.7000000000000001 + 0.75 + 0.8 + 0.8500000000000001 + 0.9 + 0.9500000000000001 + 1 + 1.05 + 1.1 + 1.15 + 1.2 + 1.25 + 1.3 + 1.35 + 1.4 + 1.45 + 1.5 + 1.55 + 1.6 + 1.65 + 1.7 + 1.75 + 1.8 + 1.85 + 1.9 + 1.95 + 2 + 2.05 + 2.1 + 2.15 + 2.2 + 2.25 + 2.3 + 2.35 + 2.4 + 2.45 + 2.5 + 2.55 + 2.6 + 2.65 + 2.7 + 2.75 + 2.8 + 2.85 + 2.9 + 2.95 + 3 + 3.05 + 3.1 + 3.15 + 3.2 + 3.25 + 3.3 + 3.35 + 3.4 + 3.45 + 3.5 + 3.55 + 3.6 + 3.65 + 3.7 + 3.75 + 3.8 + 3.85 + 3.9 + 3.95 + 4 + 4.05 + 4.100000000000001 + 4.15 + 4.2 + 4.25 + 4.3 + 4.350000000000001 + 4.4 + 4.45 + 4.5 + 4.55 + 4.600000000000001 + 4.65 + 4.7 + 4.75 + 4.800000000000001 + 4.850000000000001 + 4.9 + 4.95 + 5 + 5.050000000000001 + 5.100000000000001 + 5.15 + 5.2 + 5.25 + 5.300000000000001 + 5.350000000000001 + 5.4 + 5.45 + 5.5 + 5.550000000000001 + 5.600000000000001 + 5.65 + 5.7 + 5.75 + 5.800000000000001 + 5.850000000000001 + 5.9 + 5.95 + 6 + 6.050000000000001 + 6.100000000000001 + 6.15 + 6.2 + 6.25 + 6.300000000000001 + 6.350000000000001 + 6.4 + 6.45 + 6.5 + 6.550000000000001 + 6.600000000000001 + 6.65 + 6.7 + 6.75 + 6.800000000000001 + 6.850000000000001 + 6.9 + 6.95 + 7 + 7.050000000000001 + 7.100000000000001 + 7.15 + 7.2 + 7.25 + 7.300000000000001 + 7.350000000000001 + 7.4 + 7.45 + 7.5 + 7.550000000000001 + 7.600000000000001 + 7.65 + 7.7 + 7.75 + 7.800000000000001 + 7.850000000000001 + 7.9 + 7.95 + 8 + 8.050000000000001 + 8.1 + 8.15 + 8.200000000000001 + 8.25 + 8.300000000000001 + 8.35 + 8.4 + 8.450000000000001 + 8.5 + 8.550000000000001 + 8.6 + 8.65 + 8.700000000000001 + 8.75 + 8.800000000000001 + 8.85 + 8.9 + 8.950000000000001 + 9 + 9.050000000000001 + 9.1 + 9.15 + 9.200000000000001 + 9.25 + 9.300000000000001 + 9.35 + 9.4 + 9.450000000000001 + 9.5 + 9.550000000000001 + 9.600000000000001 + 9.65 + 9.700000000000001 + 9.75 + 9.800000000000001 + 9.850000000000001 + 9.9 + 9.950000000000001 + 10 + 10.05 + 10.1 + 10.15 + 10.2 + 10.25 + 10.3 + 10.35 + 10.4 + 10.45 + 10.5 + 10.55 + 10.6 + 10.65 + 10.7 + 10.75 + 10.8 + 10.85 + 10.9 + 10.95 + 11 + 11.05 + 11.1 + 11.15 + 11.2 + 11.25 + 11.3 + 11.35 + 11.4 + 11.45 + 11.5 + 11.55 + 11.6 + 11.65 + 11.7 + 11.75 + 11.8 + 11.85 + 11.9 + 11.95 + 12 + 12.05 + 12.1 + 12.15 + 12.2 + 12.25 + 12.3 + 12.35 + 12.4 + 12.45 + 12.5 + 12.55 + 12.6 + 12.65 + 12.7 + 12.75 + 12.8 + 12.85 + 12.9 + 12.95 + 13 + 13.05 + 13.1 + 13.15 + 13.2 + 13.25 + 13.3 + 13.35 + 13.4 + 13.45 + 13.5 + 13.55 + 13.6 + 13.65 + 13.7 + 13.75 + 13.8 + 13.85 + 13.9 + 13.95 + 14 + 14.05 + 14.1 + 14.15 + 14.2 + 14.25 + 14.3 + 14.35 + 14.4 + 14.45 + 14.5 + 14.55 + 14.6 + 14.65 + 14.7 + 14.75 + 14.8 + 14.85 + 14.9 + 14.95 + 15 + 15.05 + 15.1 + 15.15 + 15.2 + 15.25 + 15.3 + 15.35 + 15.4 + 15.45 + 15.5 + 15.55 + 15.6 + 15.65 + 15.7 + 15.75 + 15.8 + 15.85 + 15.9 + 15.95 + 16 + 16.05 + 16.1 + 16.15 + 16.2 + 16.25 + 16.3 + 16.35 + 16.4 + 16.45 + 16.5 + 16.55 + 16.6 + 16.65 + 16.7 + 16.75 + 16.8 + 16.85 + 16.9 + 16.95 + 17 + 17.05 + 17.1 + 17.15 + 17.2 + 17.25 + 17.3 + 17.35 + 17.4 + 17.45 + 17.5 + 17.55 + 17.6 + 17.65 + 17.7 + 17.75 + 17.8 + 17.85 + 17.9 + 17.95 + 18 + 18.05 + 18.1 + 18.15 + 18.2 + 18.25 + 18.3 + 18.35 + 18.4 + 18.45 + 18.5 + 18.55 + 18.6 + 18.65 + 18.7 + 18.75 + 18.8 + 18.85 + 18.9 + 18.95 + 19 + 19.05 + 19.1 + 19.15 + 19.2 + 19.25 + 19.3 + 19.35 + 19.4 + 19.45 + 19.5 + 19.55 + 19.6 + 19.65 + 19.7 + 19.75 + 19.8 + 19.85 + 19.9 + 19.95 + 20 + 20.05 + 20.1 + 20.15 + 20.2 + 20.25 + 20.3 + 20.35 + 20.4 + 20.45 + 20.5 + 20.55 + 20.6 + 20.65 + 20.7 + 20.75 + 20.8 + 20.85 + 20.9 + 20.95 + 21 + 21.05 + 21.1 + 21.15 + 21.2 + 21.25 + 21.3 + 21.35 + 21.4 + 21.45 + 21.5 + 21.55 + 21.6 + 21.65 + 21.7 + 21.75 + 21.8 + 21.85 + 21.9 + 21.95 + 22 + 22.05 + 22.1 + 22.15 + 22.2 + 22.25 + 22.3 + 22.35 + 22.4 + 22.45 + 22.5 + 22.55 + 22.6 + 22.65 + 22.7 + 22.75 + 22.8 + 22.85 + 22.9 + 22.95 + 23 + 23.05 + 23.1 + 23.15 + 23.2 + 23.25 + 23.3 + 23.35 + 23.4 + 23.45 + 23.5 + 23.55 + 23.6 + 23.65 + 23.7 + 23.75 + 23.8 + 23.85 + 23.9 + 23.95 + 24 + 24.05 + 24.1 + 24.15 + 24.2 + 24.25 + 24.3 + 24.35 + 24.4 + 24.45 + 24.5 + 24.55 + 24.6 + 24.65 + 24.7 + 24.75 + 24.8 + 24.85 + 24.9 + 24.95 + 25 + 25.05 + 25.1 + 25.15 + 25.2 + 25.25 + 25.3 + 25.35 + 25.4 + 25.45 + 25.5 + 25.55 + 25.6 + 25.65 + 25.7 + 25.75 + 25.8 + 25.85 + 25.9 + 25.95 + 26 + 26.05 + 26.1 + 26.15 + 26.2 + 26.25 + 26.3 + 26.35 + 26.4 + 26.45 + 26.5 + 26.55 + 26.6 + 26.65 + 26.7 + 26.75 + 26.8 + 26.85 + 26.9 + 26.95 + 27 + 27.05 + 27.1 + 27.15 + 27.2 + 27.25 + 27.3 + 27.35 + 27.4 + 27.45 + 27.5 + 27.55 + 27.6 + 27.65 + 27.7 + 27.75 + 27.8 + 27.85 + 27.9 + 27.95 + 28 + 28.05 + 28.1 + 28.15 + 28.2 + 28.25 + 28.3 + 28.35 + 28.4 + 28.45 + 28.5 + 28.55 + 28.6 + 28.65 + 28.7 + 28.75 + 28.8 + 28.85 + 28.9 + 28.95 + 29 + 29.05 + 29.1 + 29.15 + 29.2 + 29.25 + 29.3 + 29.35 + 29.4 + 29.45 + 29.5 + 29.55 + 29.6 + 29.65 + 29.7 + 29.75 + 29.8 + 29.85 + 29.9 + 29.95 + 30 + + +# name: t_sim +# type: scalar +30 + + +# name: w_gamma +# type: scalar +0.6482057062116896 + + +# name: w_phi +# type: scalar +NaN + + +# name: y_AH +# type: matrix +# rows: 601 +# columns: 1 + 0 + 1.687098435041176e-11 + 1.762592811836889e-09 + 2.483913383974933e-08 + 1.549627947623688e-07 + 6.20791364277213e-07 + 1.883954924375331e-06 + 4.728984393313541e-06 + 1.034490277755657e-05 + 2.039434810762801e-05 + 3.706872164232781e-05 + 6.312934499169241e-05 + 0.0001019353803053033 + 0.0001574596225259848 + 0.000234293401465818 + 0.0003376418431104489 + 0.0004733106937070915 + 0.0006476858383064071 + 0.0008677065636424558 + 0.001140833531586371 + 0.001475012347760123 + 0.001878633531995059 + 0.002360489623863582 + 0.002929730087687233 + 0.003595814617173589 + 0.004368465379991284 + 0.005257618686945202 + 0.006273376518736743 + 0.007425958295361934 + 0.008725653228794414 + 0.01018277355851297 + 0.01180760893146569 + 0.01361038215302795 + 0.01560120650323242 + 0.01779004478285701 + 0.02018667022669558 + 0.02280062939635355 + 0.02564120714206698 + 0.02871739370220445 + 0.03203785399014898 + 0.03561089910105374 + 0.03944446005540498 + 0.04354606378230296 + 0.04792281133278319 + 0.0525813583022509 + 0.05752789743109955 + 0.06276814334374274 + 0.0683073193785274 + 0.07415014645423569 + 0.08030083391305265 + 0.0867630722749043 + 0.09354002783389394 + 0.1006343390241196 + 0.1080481144793844 + 0.1157829327091595 + 0.1238398433115748 + 0.1322193696431459 + 0.14092151286435 + 0.1499457572799974 + 0.1592910768935681 + 0.1689559430952517 + 0.1789383334043138 + 0.1892357411875784 + 0.1998451862772256 + 0.2107632264127363 + 0.2219859694336389 + 0.2335090861517019 + 0.2453278238333429 + 0.2574370202252816 + 0.2698311180588043 + 0.2825041799704461 + 0.2954499037793794 + 0.3086616380643445 + 0.3221323979855202 + 0.3358548812993249 + 0.3498214845167291 + 0.3640243191582495 + 0.3784552280613643 + 0.3931058016986362 + 0.4079673944673421 + 0.42303114091388 + 0.4382879718586457 + 0.4537286303894431 + 0.4693436876938045 + 0.4851235587028432 + 0.5010585175214434 + 0.5171387126217069 + 0.5333541817786119 + 0.5496948667288045 + 0.5661506275353336 + 0.5827112566429513 + 0.5993664926103306 + 0.6161060335072094 + 0.6329195499660468 + 0.6497966978792692 + 0.6667271307346126 + 0.6837005115824019 + 0.7007065246298773 + 0.7177348864588765 + 0.7347753568642889 + 0.7518177493117548 + 0.7688519410140591 + 0.7858678826265719 + 0.8028556075629429 + 0.8198052409330251 + 0.8367070081057341 + 0.8535512429001951 + 0.8703283954091398 + 0.8870290394590603 + 0.903643879712121 + 0.9201637584152702 + 0.9365796618023924 + 0.9528827261556939 + 0.9690642435328161 + 0.9851156671664417 + 1.001028616543385 + 1.016794882170343 + 1.032406430033653 + 1.047855405760514 + 1.063134138489228 + 1.078235144456086 + 1.093151130306565 + 1.107874996138512 + 1.122399838284982 + 1.136718951844389 + 1.150825832965572 + 1.164714180895293 + 1.178377899795658 + 1.19181110033881 + 1.20500810108617 + 1.217963429659361 + 1.230671823709864 + 1.243128231694268 + 1.255327813461882 + 1.26726594066131 + 1.278938196972423 + 1.290340378170039 + 1.30146849202543 + 1.31231875805163 + 1.322887607098342 + 1.333171680802104 + 1.343167830897153 + 1.352873118392328 + 1.362284812619146 + 1.371400390156022 + 1.380217533633467 + 1.388734130424931 + 1.396948271227785 + 1.404858248538809 + 1.412462555028383 + 1.419759881817461 + 1.426749116661232 + 1.433429342043253 + 1.439799833183726 + 1.445860055965419 + 1.451609664780644 + 1.457048500302573 + 1.462176587184064 + 1.466994131687043 + 1.471501519245407 + 1.475699311964305 + 1.479588246058548 + 1.483169229232812 + 1.486443338006241 + 1.489411814983919 + 1.492076066077661 + 1.494437657678463 + 1.496498313782898 + 1.498259913075671 + 1.499724485970504 + 1.500894211611441 + 1.501771414836618 + 1.502358563106524 + 1.502658263398668 + 1.502673259070596 + 1.502406426693104 + 1.5018607728555 + 1.501039430944698 + 1.49994565789993 + 1.49858283094479 + 1.496954444298351 + 1.495064105867014 + 1.492915533918777 + 1.490512553741553 + 1.48785909428716 + 1.484959184802603 + 1.481816951450222 + 1.478436613918288 + 1.474822482023605 + 1.47097895230767 + 1.466910504627921 + 1.462621698745595 + 1.458117170911714 + 1.453401630452694 + 1.448479856357084 + 1.44335669386489 + 1.438037051060989 + 1.432525895474072 + 1.426828250682576 + 1.420949192929068 + 1.414893847744482 + 1.408667386583665 + 1.402275023473631 + 1.395722011675936 + 1.389013640364563 + 1.382155231320701 + 1.375152135645797 + 1.368009730494228 + 1.360733415826955 + 1.353328611187479 + 1.34580075250143 + 1.338155288901083 + 1.330397679576106 + 1.322533390651798 + 1.314567892096092 + 1.306506654656546 + 1.298355146828562 + 1.290118831856022 + 1.281803164765535 + 1.273413589435458 + 1.264955535700831 + 1.256434416495348 + 1.24785562503147 + 1.239224532019752 + 1.230546482928431 + 1.22182679528432 + 1.213070756015989 + 1.204283618840234 + 1.195470601692768 + 1.186636884204069 + 1.177787605221274 + 1.168927860376997 + 1.160062699705907 + 1.151197125309872 + 1.142336089072464 + 1.133484490423565 + 1.1246471741548 + 1.115828928286492 + 1.107034481986787 + 1.09826850354359 + 1.0895355983899 + 1.080840307183105 + 1.072187103938783 + 1.06358039421949 + 1.055024513379008 + 1.046523724862492 + 1.038082218562907 + 1.029704109234125 + 1.021393434961015 + 1.013154155686826 + 1.004990151798124 + 0.9969052227675322 + 0.9889030858544485 + 0.9809873748639296 + 0.9731616389638602 + 0.9654293415605164 + 0.9577938592325803 + 0.9502584807236509 + 0.9428264059932429 + 0.9355007453262481 + 0.9282845185007935 + 0.9211806540144024 + 0.9141919883683307 + 0.9073212654099204 + 0.900571135732782 + 0.8939441561345822 + 0.887442789132195 + 0.8810694025339247 + 0.8748262690685025 + 0.8687155660705089 + 0.8627393752218617 + 0.8568996823489688 + 0.851198377275125 + 0.8456372537277042 + 0.8402180092996634 + 0.8349422454648647 + 0.8298114676466771 + 0.8248270853393093 + 0.8199904122812932 + 0.8153026666805145 + 0.8107649714901654 + 0.8063783547349704 + 0.8021437498870115 + 0.7980619962904635 + 0.79413383963452 + 0.7903599324737832 + 0.7867408347953537 + 0.7832770146318552 + 0.779968848719596 + 0.7768166232010623 + 0.7738205343709108 + 0.7709806894646216 + 0.7682971074889473 + 0.7657697200932875 + 0.7633983724810952 + 0.7611828243604137 + 0.7591227509326265 + 0.7572177439184892 + 0.755467312620505 + 0.7538708850206876 + 0.7524278089127489 + 0.7511373530677399 + 0.74999870843216 + 0.749010989357542 + 0.7481732348605172 + 0.7474844099123515 + 0.7469434067569375 + 0.7465490462562266 + 0.7463000792620738 + 0.7461951880134658 + 0.7462329875580996 + 0.7464120271972744 + 0.7467307919530557 + 0.7471877040566702 + 0.7477811244570883 + 0.7485093543487468 + 0.7493706367173724 + 0.7503631579028528 + 0.7514850491781196 + 0.7527343883429959 + 0.7541092013319656 + 0.7556074638348348 + 0.7572271029292409 + 0.7589659987239833 + 0.7608219860121493 + 0.7627928559330091 + 0.7648763576416656 + 0.7670701999854497 + 0.7693720531860512 + 0.7717795505263927 + 0.7742902900412512 + 0.7769018362106491 + 0.7796117216550386 + 0.7824174488313104 + 0.7853164917286765 + 0.7883062975634753 + 0.7913842884719643 + 0.7945478632001757 + 0.797794398789915 + 0.801121252260006 + 0.8045257622818833 + 0.8080052508486558 + 0.8115570249367741 + 0.8151783781594442 + 0.8188665924109508 + 0.8226189395010584 + 0.8264326827786773 + 0.8303050787439961 + 0.8342333786482932 + 0.8382148300806591 + 0.8422466785408737 + 0.8463261689976947 + 0.8504505474318415 + 0.8546170623629534 + 0.858822966359843 + 0.8630655175333574 + 0.8673419810111948 + 0.8716496303940309 + 0.8759857491923281 + 0.8803476322432198 + 0.8847325871068776 + 0.88913793544179 + 0.8935610143583889 + 0.8979991777504931 + 0.9024497976040429 + 0.9069102652826219 + 0.9113779927892829 + 0.9158504140042131 + 0.9203249858977813 + 0.924799189718547 + 0.9292705321558112 + 0.9337365464763203 + 0.9381947936347401 + 0.9426428633575532 + 0.9470783752000316 + 0.9514989795759714 + 0.9559023587598798 + 0.960286227861341 + 0.9646483357712866 + 0.9689864660799303 + 0.9732984379661339 + 0.9775821070580001 + 0.9818353662644923 + 0.9860561465779148 + 0.9902424178470948 + 0.9943921895211261 + 0.9985035113635583 + 1.002574474136927 + 1.006603210257539 + 1.010587894420451 + 1.01452674419458 + 1.018418020587929 + 1.022260028582893 + 1.026051117641668 + 1.029789682181761 + 1.033474162021641 + 1.037103042796596 + 1.04067485634484 + 1.044188181063962 + 1.04764164223783 + 1.051033912334033 + 1.054363711272019 + 1.057629806662061 + 1.060831014015204 + 1.063966196924378 + 1.067034267216863 + 1.070034185078301 + 1.072964959148474 + 1.075825646589085 + 1.078615353123774 + 1.081333233050629 + 1.083978489227468 + 1.086550373030158 + 1.089048184284277 + 1.09147127117042 + 1.093819030103471 + 1.096090905586155 + 1.098286390037225 + 1.100405023594628 + 1.102446393894014 + 1.104410135822946 + 1.106295931251213 + 1.108103508737616 + 1.109832643213636 + 1.111483155644398 + 1.11305491266733 + 1.114547826208962 + 1.115961853080283 + 1.117296994551106 + 1.118553295903886 + 1.11973084596744 + 1.120829776631048 + 1.121850262339371 + 1.122792519568704 + 1.123656806285002 + 1.124443421384185 + 1.125152704115214 + 1.125785033486424 + 1.126340827655613 + 1.126820543304399 + 1.127224674997337 + 1.127553754526332 + 1.127808350240824 + 1.127989066364302 + 1.12809654229763 + 1.128131451909729 + 1.12809450281613 + 1.127986435645925 + 1.127808023297626 + 1.127560070184488 + 1.127243411469795 + 1.12685891229265 + 1.126407466984794 + 1.125889998278978 + 1.125307456509423 + 1.124660818804878 + 1.123951088274819 + 1.123179293189296 + 1.122346486152959 + 1.121453743273782 + 1.120502163326992 + 1.119492866914732 + 1.118426995621961 + 1.1173057111691 + 1.116130194561933 + 1.114901645239259 + 1.113621280218804 + 1.112290333241876 + 1.110910053917254 + 1.10948170686481 + 1.108006570859322 + 1.106485937974976 + 1.10492111273101 + 1.103313411238983 + 1.101664160352107 + 1.099974696817115 + 1.098246366429108 + 1.096480523189809 + 1.094678528469689 + 1.092841750174358 + 1.090971561915673 + 1.089069342187961 + 1.087136473549777 + 1.085174341811587 + 1.083184335229776 + 1.081167843707371 + 1.079126258001848 + 1.077060968940407 + 1.074973366643063 + 1.07286483975392 + 1.070736774680973 + 1.068590554844772 + 1.066427559936282 + 1.064249165184262 + 1.062056740632474 + 1.059851650427027 + 1.057635252114155 + 1.055408895948715 + 1.053173924213678 + 1.050931670550894 + 1.048683459303383 + 1.046430604869401 + 1.044174411068534 + 1.041916170520048 + 1.039657164033714 + 1.037398660013329 + 1.03514191387314 + 1.032888167467358 + 1.030638648532969 + 1.028394570145988 + 1.026157130191361 + 1.023927510846645 + 1.021706878079633 + 1.019496381160056 + 1.0172971521855 + 1.015110305621653 + 1.012936937857002 + 1.010778126772079 + 1.008634931323353 + 1.006508391141851 + 1.004399526146587 + 1.002309336172869 + 1.000238800615529 + 0.9981888780871464 + 0.9961605060912904 + 0.9941546007108117 + 0.9921720563112212 + 0.9902137452591534 + 0.9882805176559332 + 0.9863732010862361 + 0.9844926003818322 + 0.9826394974003981 + 0.9808146508193644 + 0.9790187959447668 + 0.977252644535055 + 0.9755168846398091 + 0.9738121804533041 + 0.9721391721828536 + 0.9704984759318616 + 0.9688906835974953 + 0.9673163627828971 + 0.9657760567238292 + 0.9642702842296555 + 0.9627995396385473 + 0.9613642927867923 + 0.9599649889920884 + 0.9586020490506858 + 0.9572758692482461 + 0.9559868213842675 + 0.954735252809935 + 0.9535214864792328 + 0.9523458210131615 + 0.9512085307768907 + 0.9501098659696787 + 0.9490500527273741 + 0.9480292932373214 + 0.9470477658654821 + 0.9461056252955748 + 0.9452030026800378 + 0.9443400058026123 + 0.9435167192523407 + 0.9427332046087603 + 0.9419895006380891 + 0.9412856235001745 + 0.9406215669659845 + 0.9399973026454167 + 0.9394127802251908 + 0.9388679277165927 + 0.9383626517128353 + 0.9378968376557898 + 0.9374703501118536 + 0.9370830330567003 + 0.9367347101686696 + 0.9364251851305467 + 0.9361542419394752 + 0.9359216452247511 + 0.9357271405732431 + 0.9355704548621792 + 0.9354512965990393 + 0.9353693562682955 + 0.9353243066847362 + 0.9353158033531107 + 0.9353434848338315 + 0.9354069731144681 + 0.93550587398677 + 0.935639777428946 + 0.9358082579929439 + 0.9360108751964545 + 0.9362471739193825 + 0.9365166848045094 + 0.9368189246620917 + 0.9371533968781218 + 0.9375195918259925 + 0.9379169872812975 + 0.9383450488395062 + 0.9388032303362529 + 0.939290974269975 + 0.9398077122266471 + 0.9403528653063457 + 0.9409258445513955 + 0.941526051375836 + + +# name: y_opt +# type: matrix +# rows: 601 +# columns: 1 + 0 + 8.354783498726138e-12 + 9.230226583754954e-10 + 1.367237371013527e-08 + 8.917343257190273e-08 + 3.716712185197994e-07 + 1.168484152022217e-06 + 3.026913971063492e-06 + 6.810276799104703e-06 + 1.376720350363142e-05 + 2.559050225541717e-05 + 4.446361595196161e-05 + 7.309398386925997e-05 + 0.0001147334883461923 + 0.0001731867227412064 + 0.0002528081419857187 + 0.0003584893220036596 + 0.0004956376121735757 + 0.0006701474540875041 + 0.000888365586656117 + 0.001157051280068725 + 0.00148333265108637 + 0.001874660017265783 + 0.002338757152779788 + 0.002883571216429221 + 0.00351722203491431 + 0.004247951342339291 + 0.005084072500650072 + 0.006033921155327065 + 0.00710580721606093 + 0.008307968493105482 + 0.009648526266246949 + 0.01113544301454131 + 0.01277648249083131 + 0.01457917228524432 + 0.01655076898608275 + 0.01869822601445498 + 0.02102816418037576 + 0.0235468449826264 + 0.02626014665215722 + 0.02917354291900438 + 0.03229208446536171 + 0.03562038301239067 + 0.03916259797537757 + 0.04292242561077839 + 0.04690309056936139 + 0.05110733976191238 + 0.05553743843766176 + 0.06019516837059315 + 0.06508182804497535 + 0.0701982347287068 + 0.07554472832126834 + 0.08112117686214315 + 0.08692698358539533 + 0.09296109540660806 + 0.09922201272949839 + 0.1057078004611661 + 0.1124161001270413 + 0.1193441429790987 + 0.1264887639937556 + 0.1338464166590125 + 0.1414131884537744 + 0.1491848169258795 + 0.1571567062791031 + 0.1653239443832721 + 0.1736813201255887 + 0.1822233410252775 + 0.190944251037732 + 0.1998380484783996 + 0.2088985040007027 + 0.2181191785663198 + 0.2274934413501296 + 0.2370144875260407 + 0.2466753558837732 + 0.2564689462304136 + 0.2663880365342285 + 0.276425299771779 + 0.2865733204428215 + 0.2968246107208143 + 0.3071716262100505 + 0.3176067812835275 + 0.3281224639786114 + 0.3387110504303831 + 0.3493649188252488 + 0.3600764628599579 + 0.3708381046936077 + 0.3816423073825184 + 0.3924815867900338 + 0.4033485229653582 + 0.4142357709874612 + 0.4251360712718876 + 0.4360422593399935 + 0.4469472750517 + 0.4578441713043138 + 0.4687261222013105 + 0.4795864306962213 + 0.4904185357179028 + 0.5012160187845198 + 0.511972610114513 + 0.522682194243692 + 0.5333388151583646 + 0.5439366809551066 + 0.5544701680383948 + 0.5649338248678623 + 0.5753223752674093 + 0.5856307213088067 + 0.5958539457827696 + 0.6059873142707646 + 0.6160262768310338 + 0.6259664693125073 + 0.6358037143103856 + 0.6455340217772687 + 0.6551535893037382 + 0.6646588020823028 + 0.674046232568578 + 0.6833126398535024 + 0.6924549687602903 + 0.7014703486796899 + 0.7103560921569632 + 0.719109693243822 + 0.7277288256283587 + 0.7362113405557844 + 0.7445552645525576 + 0.7527587969662299 + 0.7608203073330748 + 0.7687383325852826 + 0.7765115741092264 + 0.7841388946660017 + 0.7916193151851404 + 0.7989520114420955 + 0.8061363106297803 + 0.8131716878341183 + 0.8200577624232581 + 0.8267942943597724 + 0.8333811804448492 + 0.8398184505031537 + 0.8461062635167312 + 0.8522449037159894 + 0.8582347766354952 + 0.864076405142002 + 0.8697704254418169 + 0.8753175830743094 + 0.880718728898065 + 0.8859748150758943 + 0.8910868910646067 + 0.8960560996151888 + 0.9008836727887291 + 0.9055709279931734 + 0.9101192640457221 + 0.9145301572654093 + 0.9188051576001679 + 0.9229458847924187 + 0.9269540245869909 + 0.9308313249849394 + 0.9345795925465992 + 0.9382006887469964 + 0.9416965263865129 + 0.9450690660595038 + 0.9483203126833546 + 0.9514523120902767 + 0.9544671476839444 + 0.9573669371629047 + 0.9601538293124982 + 0.9628300008668784 + 0.9653976534425389 + 0.9678590105446043 + 0.9702163146469946 + 0.9724718243474183 + 0.9746278115980145 + 0.9766865590123304 + 0.9786503572491869 + 0.9805215024738684 + 0.9823022938969481 + 0.9839950313909509 + 0.9856020131849487 + 0.9871255336370724 + 0.98856788108484 + 0.9899313357730909 + 0.9912181678592407 + 0.9924306354954753 + 0.9935709829874314 + 0.9946414390288271 + 0.9956442150114385 + 0.9965815034097476 + 0.9974554762395216 + 0.9982682835895277 + 0.9990220522255204 + 0.999718884265596 + 1.000360855925949 + 1.000950016336024 + 1.001488386422008 + 1.001977957857573 + 1.002420692080731 + 1.002818519375643 + 1.003173338018173 + 1.003487013483962 + 1.003761377717771 + 1.003998228462802 + 1.004199328648698 + 1.004366405836906 + 1.004501151722045 + 1.004605221687944 + 1.004680234416964 + 1.004727771551231 + 1.004749377404389 + 1.00474655872247 + 1.004720784492489 + 1.004673485797347 + 1.004606055715634 + 1.004519849264927 + 1.004416183387165 + 1.004296336974712 + 1.004161550935682 + 1.004013028297144 + 1.003851934344813 + 1.003679396797849 + 1.003496506017372 + 1.003304315247352 + 1.00310384088651 + 1.002896062789883 + 1.00268192459875 + 1.002462334097573 + 1.002238163596686 + 1.002010250339432 + 1.001779396932496 + 1.00154637179817 + 1.001311909647344 + 1.001076711971998 + 1.000841447556001 + 1.000606753003072 + 1.000373233280718 + 1.000141462279061 + 0.999911983383409 + 0.999685310059522 + 0.9994619264504879 + 0.9992422879841836 + 0.9990268219902994 + 0.998815928325941 + 0.9986099800088383 + 0.9984093238572179 + 0.9982142811354255 + 0.9980251482043975 + 0.9978421971761136 + 0.9976656765711909 + 0.9974958119787916 + 0.9973328067180539 + 0.9971768425002749 + 0.9970280800910979 + 0.9968866599719867 + 0.9967527030002874 + 0.9966263110672099 + 0.9965075677530797 + 0.9963965389792399 + 0.996293273656006 + 0.9961978043260968 + 0.9961101478029969 + 0.996030305803721 + 0.9959582655754806 + 0.9958940005157728 + 0.9958374707854379 + 0.9957886239142506 + 0.9957473953986373 + 0.9957137092911286 + 0.9956874787811855 + 0.9956686067670524 + 0.9956569864183129 + 0.9956525017288507 + 0.9956550280599293 + 0.9956644326731349 + 0.9956805752529357 + 0.9957033084186416 + 0.9957324782255553 + 0.9957679246551384 + 0.9958094820940181 + 0.9958569798016959 + 0.9959102423668194 + 0.9959690901519086 + 0.9960333397264391 + 0.9961028042882009 + 0.9961772940728648 + 0.9962566167517137 + 0.9963405778174947 + 0.9964289809583791 + 0.9965216284200161 + 0.9966183213556943 + 0.9967188601646243 + 0.9968230448183769 + 0.9969306751755205 + 0.9970415512845162 + 0.9971554736749321 + 0.9972722436370598 + 0.9973916634900178 + 0.9975135368384414 + 0.9976376688178659 + 0.997763866328921 + 0.9978919382604572 + 0.9980216957017424 + 0.9981529521438676 + 0.9982855236705072 + 0.9984192291381961 + 0.9985538903462763 + 0.9986893321966877 + 0.9988253828437754 + 0.9989618738342894 + 0.9990986402377643 + 0.9992355207674652 + 0.9993723578920951 + 0.9995089979384558 + 0.9996452911852681 + 0.9997810919483529 + 0.999916258657377 + 1.000050653924375 + 1.000184144604257 + 1.000316601847521 + 1.00044790114537 + 1.000577922367471 + 1.00070654979255 + 1.00083367213206 + 1.000959182547128 + 1.001082978659008 + 1.001204962553245 + 1.00132504077779 + 1.00144312433526 + 1.001559128669574 + 1.001672973647182 + 1.001784583533089 + 1.001893886961905 + 1.002000816904113 + 1.002105310627787 + 1.002207309655946 + 1.002306759719771 + 1.002403610707865 + 1.002497816611786 + 1.002589335468026 + 1.002678129296641 + 1.002764164036738 + 1.002847409478985 + 1.00292783919535 + 1.003005430466247 + 1.003080164205267 + 1.003152024881676 + 1.003221000440847 + 1.003287082222808 + 1.003350264879054 + 1.003410546287811 + 1.003467927467882 + 1.003522412491257 + 1.00357400839462 + 1.003622725089909 + 1.003668575274071 + 1.003711574338152 + 1.00375174027585 + 1.00378909359168 + 1.003823657208861 + 1.003855456377054 + 1.003884518580073 + 1.003910873443687 + 1.003934552643611 + 1.00395558981381 + 1.003974020455207 + 1.003989881844902 + 1.004003212945986 + 1.00401405431806 + 1.004022448028526 + 1.004028437564743 + 1.004032067747134 + 1.004033384643303 + 1.004032435483249 + 1.004029268575739 + 1.004023933225902 + 1.00401647965411 + 1.004006958916197 + 1.003995422825081 + 1.00398192387382 + 1.003966515160172 + 1.003949250312676 + 1.003930183418327 + 1.003909368951844 + 1.003886861706598 + 1.003862716727206 + 1.003836989243839 + 1.003809734608248 + 1.003781008231549 + 1.003750865523767 + 1.003719361835176 + 1.00368655239943 + 1.003652492278502 + 1.003617236309445 + 1.003580839052974 + 1.003543354743876 + 1.003504837243242 + 1.003465339992539 + 1.003424915969484 + 1.003383617645763 + 1.003341496946537 + 1.003298605211764 + 1.00325499315931 + 1.003210710849831 + 1.003165807653429 + 1.003120332218044 + 1.003074332439589 + 1.003027855433784 + 1.002980947509686 + 1.002933654144889 + 1.002886019962363 + 1.002838088708918 + 1.002789903235268 + 1.002741505477656 + 1.002692936441023 + 1.002644236183694 + 1.002595443803545 + 1.002546597425623 + 1.002497734191193 + 1.002448890248171 + 1.002400100742922 + 1.002351399813383 + 1.002302820583479 + 1.002254395158801 + 1.002206154623509 + 1.002158129038426 + 1.002110347440294 + 1.002062837842145 + 1.002015627234769 + 1.001968741589226 + 1.001922205860378 + 1.001876043991406 + 1.00183027891927 + 1.001784932581077 + 1.001740025921337 + 1.001695578900044 + 1.001651610501574 + 1.001608138744353 + 1.001565180691255 + 1.001522752460712 + 1.001480869238481 + 1.00143954529006 + 1.001398793973689 + 1.001358627753932 + 1.001319058215789 + 1.001280096079312 + 1.001241751214695 + 1.0012040326578 + 1.001166948626099 + 1.00113050653499 + 1.001094713014466 + 1.001059573926102 + 1.00102509438034 + 1.00099127875403 + 1.000958130708211 + 1.000925653206105 + 1.000893848531287 + 1.000862718306024 + 1.000832263509731 + 1.000802484497552 + 1.000773381019014 + 1.000744952236749 + 1.000717196745254 + 1.000690112589663 + 1.000663697284528 + 1.000637947832565 + 1.000612860743359 + 1.00058843205201 + 1.000564657337697 + 1.00054153174214 + 1.00051904998795 + 1.000497206396845 + 1.000475994907724 + 1.000455409094565 + 1.000435442184165 + 1.000416087073673 + 1.000397336347923 + 1.000379182296563 + 1.000361616930934 + 1.000344632000734 + 1.000328219010412 + 1.000312369235314 + 1.00029707373756 + 1.000282323381635 + 1.000268108849706 + 1.000254420656634 + 1.000241249164695 + 1.00022858459799 + 1.000216417056545 + 1.000204736530084 + 1.000193532911496 + 1.000182796009959 + 1.00017251556374 + 1.000162681252664 + 1.000153282710237 + 1.000144309535437 + 1.000135751304164 + 1.000127597580337 + 1.000119837926661 + 1.00011246191504 + 1.000105459136647 + 1.000098819211653 + 1.000092531798603 + 1.000086586603461 + 1.000080973388303 + 1.000075681979671 + 1.000070702276594 + 1.000066024258267 + 1.000061637991395 + 1.000057533637213 + 1.000053701458171 + 1.000050131824295 + 1.000046815219236 + 1.000043742245983 + 1.000040903632285 + 1.000038290235747 + 1.000035893048627 + 1.000033703202339 + 1.00003171197165 + 1.000029910778597 + 1.000028291196113 + 1.000026844951376 + 1.000025563928884 + 1.000024440173262 + 1.000023465891801 + 1.00002263345675 + 1.00002193540735 + 1.000021364451623 + 1.000020913467923 + 1.000020575506258 + 1.000020343789378 + 1.00002021171365 + 1.000020172849714 + 1.000020220942931 + 1.000020349913634 + 1.000020553857181 + 1.000020827043817 + 1.000021163918358 + 1.000021559099696 + 1.000022007380133 + 1.000022503724558 + 1.000023043269462 + 1.000023621321798 + 1.000024233357707 + 1.000024875021098 + 1.000025542122096 + 1.000026230635369 + 1.000026936698325 + 1.000027656609209 + 1.000028386825078 + 1.000029123959684 + 1.000029864781248 + 1.000030606210159 + 1.000031345316574 + 1.000032079317943 + 1.000032805576461 + 1.000033521596446 + 1.000034225021661 + 1.000034913632565 + 1.00003558534352 + 1.000036238199947 + 1.000036870375432 + 1.000037480168799 + 1.000038066001145 + 1.000038626412841 + 1.000039160060517 + 1.000039665714011 + 1.000040142253312 + 1.00004058866548 + 1.000041004041564 + 1.000041387573505 + 1.000041738551045 + 1.000042056358631 + 1.000042340472321 + 1.000042590456711 + 1.000042805961851 + 1.000042986720196 + 1.00004313254356 + 1.000043243320095 + 1.000043319011291 + 1.000043359649002 + 1.000043365332497 + 1.00004333622554 + 1.000043272553508 + 1.000043174600537 + 1.000043042706704 + 1.000042877265249 + 1.000042678719842 + 1.00004244756188 + 1.000042184327834 + 1.000041889596643 + 1.000041563987149 + 1.000041208155583 + 1.000040822793095 + 1.000040408623339 + 1.000039966400105 + 1.000039496905005 + 1.000039000945208 + 1.000038479351236 + 1.000037932974803 + 1.000037362686714 + 1.000036769374825 + 1.000036153942044 + 1.000035517304404 + 1.000034860389179 + 1.000034184133062 + 1.000033489480401 + + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-13 08:09:18
|
Revision: 11001 http://octave.svn.sourceforge.net/octave/?rev=11001&view=rev Author: paramaniac Date: 2012-09-13 08:09:07 +0000 (Thu, 13 Sep 2012) Log Message: ----------- control: work on multiplot feature Modified Paths: -------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/__frequency_vector_2__.m trunk/octave-forge/main/control/devel/bode2.m Added Paths: ----------- trunk/octave-forge/main/control/devel/multiplot2.m Modified: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-12 19:03:09 UTC (rev 11000) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-13 08:09:07 UTC (rev 11001) @@ -24,24 +24,33 @@ ## Version: 0.4 % function [H, w] = __frequency_response__ (sys, w = [], mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) -function [H, w] = __frequency_response_2__ (mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false, varargin) +function [H, w] = __frequency_response_2__ (args, mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) - sys_idx = cellfun (@isa, varargin, {"lti"}); # true or false - w_idx = cellfun (@is_real_vector, varargin); # look for frequency vectors - % ? = cellfun (@iscell, varargin) - % varargin(?) (end) + sys_idx = cellfun (@isa, args, {"lti"}); # true or false + w_idx = cellfun (@is_real_vector, args); # look for frequency vectors + c_idx = cellfun (@iscell, args); + % args(?) (end) % w_idx(end) - if (! any (w_idx)) - w = __frequency_vector_2__ (varargin(sys_idx), wbounds); + if (any (c_idx)) + w = args(c_idx){end}; + if (numel (w) == 2 && issample (w{1}) && issample (w{2})) + w = __frequency_vector_2__ (args(sys_idx), wbounds, w{1}, w{2}); + else + error ("frequency_response: invalid cell"); + endif + elseif (any (w_idx)) + w = args(w_idx){end}; + else + w = __frequency_vector_2__ (args(sys_idx), wbounds); endif -%varargin{sys_idx} +%args{sys_idx} - H = cellfun (@__freqresp__, varargin(sys_idx), {w}, {resptype}, {cellflag}, "uniformoutput", false); + H = cellfun (@__freqresp__, args(sys_idx), {w}, {resptype}, {cellflag}, "uniformoutput", false); %{ ## check arguments Modified: trunk/octave-forge/main/control/devel/__frequency_vector_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-12 19:03:09 UTC (rev 11000) +++ trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-13 08:09:07 UTC (rev 11001) @@ -37,7 +37,7 @@ ## Date: October 2009 ## Version: 0.2 -function w = __frequency_vector_2__ (sys_cell, wbounds = "std") +function w = __frequency_vector_2__ (sys_cell, wbounds = "std", wmin, wmax) if (! iscell (sys_cell)) sys_cell = {sys_cell} @@ -48,8 +48,16 @@ [dec_min, dec_max, zp] = cellfun (@(x) __frequency_range__ (x, wbounds), sys_cell, "uniformoutput", false); - dec_min = min (cell2mat (dec_min)); - dec_max = max (cell2mat (dec_max)); + if (nargin == 2) + dec_min = min (cell2mat (dec_min)); + dec_max = max (cell2mat (dec_max)); + elseif (nargin == 4) # w = {wmin, wmax} + dec_min = log10 (wmin); + dec_max = log10 (wmax); + else + print_usage (); + endif + zp = horzcat (zp{:}); ## include zeros and poles for nice peaks in plots Modified: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-12 19:03:09 UTC (rev 11000) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-13 08:09:07 UTC (rev 11001) @@ -58,7 +58,7 @@ print_usage (); endif - [H, w] = __frequency_response_2__ (false, 0, "std", false, varargin{:}); + [H, w] = __frequency_response_2__ (varargin, false, 0, "std", false); H = cellfun (@reshape, H, {[]}, {1}, "uniformoutput", false); mag = cellfun (@abs, H, "uniformoutput", false); Added: trunk/octave-forge/main/control/devel/multiplot2.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot2.m (rev 0) +++ trunk/octave-forge/main/control/devel/multiplot2.m 2012-09-13 08:09:07 UTC (rev 11001) @@ -0,0 +1,110 @@ +%% -*- texinfo -*- +%% Frequency-weighted controller reduction. + +% =============================================================================== +% Frequency Weighted Controller Reduction Lukas Reichlin December 2011 +% =============================================================================== +% Reference: Madievski, A.G. and Anderson, B.D.O. +% Sampled-Data Controller Reduction Procedure +% IEEE Transactions of Automatic Control +% Vol. 40, No. 11, November 1995 +% =============================================================================== + +% Tabula Rasa +clear all, close all, clc + +% Plant +Ap1 = [ 0.0 1.0 + 0.0 0.0 ]; + +Ap2 = [ -0.015 0.765 + -0.765 -0.015 ]; + +Ap3 = [ -0.028 1.410 + -1.410 -0.028 ]; + +Ap4 = [ -0.04 1.85 + -1.85 -0.04 ]; + +Ap = blkdiag (Ap1, Ap2, Ap3, Ap4); + +Bp = [ 0.026 + -0.251 + 0.033 + -0.886 + -4.017 + 0.145 + 3.604 + 0.280 ]; + +Cp = [ -0.996 -0.105 0.261 0.009 -0.001 -0.043 0.002 -0.026 ]; + +Dp = [ 0.0 ]; + +P = ss (Ap, Bp, Cp, Dp); + +% Controller +Ac = [ -0.4077 0.9741 0.1073 0.0131 0.0023 -0.0186 -0.0003 -0.0098 + -0.0977 -0.1750 0.0215 -0.0896 -0.0260 0.0057 0.0109 -0.0105 + 0.0011 0.0218 -0.0148 0.7769 0.0034 -0.0013 -0.0014 0.0011 + -0.0361 -0.5853 -0.7701 -0.3341 -0.0915 0.0334 0.0378 -0.0290 + -0.1716 -2.6546 -0.0210 -1.4467 -0.4428 1.5611 0.1715 -0.1318 + -0.0020 0.0950 0.0029 0.0523 -1.3950 -0.0338 -0.0062 0.0045 + 0.1607 2.3824 0.0170 1.2979 0.3721 -0.1353 -0.1938 1.9685 + -0.0006 0.1837 0.0048 0.1010 0.0289 -0.0111 -1.8619 -0.0311 ]; + +Bc = [ -0.4105 + -0.0868 + -0.0004 + 0.0036 + 0.0081 + -0.0085 + -0.0004 + -0.0132 ]; + +Cc = [ -0.0447 -0.6611 -0.0047 -0.3601 -0.1033 0.0375 0.0427 -0.0329 ]; + +Dc = [ 0.0 ]; + +K = ss (Ac, Bc, Cc, Dc); + +% Controller Reduction +Kr4 = spaconred (P, K, 4, 'feedback', '-') +Kr2 = spaconred (P, K, 2, 'feedback', '-') + +%{ +% Open Loop +L = P * K; +Lr4 = P * Kr4; +Lr2 = P * Kr2; + +% Closed Loop +T = feedback (L); +Tr4 = feedback (Lr4); +Tr2 = feedback (Lr2); +%} + +% Frequency Range +w = {1e-2, 1e1}; +% w = logspace (-2, 1, 500); + +% Bode Plot of Controller +figure (1) +bode2 (K, Kr4, Kr2, w) +title ('Bode Diagrams of K and Kr') +legend ('K (8 states)', 'Kr (4 states)', 'Kr (2 states)', 'location', 'southwest') + +%{ +% Step Response of Closed Loop +[y, t] = step (T, 100); +[yr4, tr4] = step (Tr4, 100); +[yr2, tr2] = step (Tr2, 100); + +figure (2) +plot (t, y, tr4, yr4, tr2, yr2) +grid ('on') +title ('Step Response of Closed Loop') +xlabel ('Time [s]') +ylabel ('Output [-]') +legend ('K (8 states)', 'Kr (4 states)', 'Kr (2 states)', 'Location', 'SouthEast') +%} \ No newline at end of file This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-13 14:34:37
|
Revision: 11002 http://octave.svn.sourceforge.net/octave/?rev=11002&view=rev Author: paramaniac Date: 2012-09-13 14:34:26 +0000 (Thu, 13 Sep 2012) Log Message: ----------- control: handle FRD models in multiplot function Modified Paths: -------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/__frequency_vector_2__.m trunk/octave-forge/main/control/devel/bode2.m trunk/octave-forge/main/control/devel/multiplot.m Modified: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-13 08:09:07 UTC (rev 11001) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-13 14:34:26 UTC (rev 11002) @@ -26,32 +26,40 @@ % function [H, w] = __frequency_response__ (sys, w = [], mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) function [H, w] = __frequency_response_2__ (args, mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) - sys_idx = cellfun (@isa, args, {"lti"}); # true or false - w_idx = cellfun (@is_real_vector, args); # look for frequency vectors - c_idx = cellfun (@iscell, args); - % args(?) (end) + sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models + w_idx = cellfun (@is_real_vector, args); # look for frequency vectors + r_idx = cellfun (@iscell, args); # look for frequency ranges {wmin, wmax} + + sys_cell = args(sys_idx); # extract LTI models + frd_idx = cellfun (@isa, sys_cell, {"frd"}); # look for FRD models - -% w_idx(end) - - if (any (c_idx)) - w = args(c_idx){end}; - if (numel (w) == 2 && issample (w{1}) && issample (w{2})) - w = __frequency_vector_2__ (args(sys_idx), wbounds, w{1}, w{2}); + if (any (r_idx)) # if there are frequency ranges + r = args(r_idx){end}; # take the last one + if (numel (r) == 2 && issample (r{1}) && issample (r{2})) + w = __frequency_vector_2__ (sys_cell, wbounds, r{1}, r{2}); else error ("frequency_response: invalid cell"); endif - elseif (any (w_idx)) + elseif (any (w_idx)) # are there any frequency vectors? w = args(w_idx){end}; - else - w = __frequency_vector_2__ (args(sys_idx), wbounds); + else # there are neither frequency ranges nor vectors + w = __frequency_vector_2__ (sys_cell, wbounds); endif -%args{sys_idx} + w = repmat ({w}, 1, numel (sys_cell)); # return cell of frequency vectors + w(frd_idx) = {[]}; # freqresp returns all frequencies of FRD models for w=[] + ## compute frequency response H for all LTI models + H = cellfun (@__freqresp__, sys_cell, w, {resptype}, {cellflag}, "uniformoutput", false); - H = cellfun (@__freqresp__, args(sys_idx), {w}, {resptype}, {cellflag}, "uniformoutput", false); + ## save frequency vectors of FRD models in w + w_frd = cellfun (@get, sys_cell(frd_idx), {"w"}, "uniformoutput", false); + w(frd_idx) = w_frd; +%w +%w = get (sys, "w"); +%args{sys_idx} + %{ ## check arguments if(! isa (sys, "lti")) Modified: trunk/octave-forge/main/control/devel/__frequency_vector_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-13 08:09:07 UTC (rev 11001) +++ trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-13 14:34:26 UTC (rev 11002) @@ -46,7 +46,7 @@ idx = cellfun (@(x) isa (x, "lti"), sys_cell); sys_cell = sys_cell(idx); - [dec_min, dec_max, zp] = cellfun (@(x) __frequency_range__ (x, wbounds), sys_cell, "uniformoutput", false); + [dec_min, dec_max, zp] = cellfun (@__frequency_range__, sys_cell, {wbounds}, "uniformoutput", false); if (nargin == 2) dec_min = min (cell2mat (dec_min)); @@ -78,6 +78,14 @@ function [dec_min, dec_max, zp] = __frequency_range__ (sys, wbounds = "std") + if (isa (sys, "frd")) + w = get (sys, "w"); + dec_min = log10 (w(1)); + dec_max = log10 (w(end)); + zp = []; + return; + endif + zer = zero (sys); pol = pole (sys); tsam = abs (get (sys, "tsam")); # tsam could be -1 Modified: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-13 08:09:07 UTC (rev 11001) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-13 14:34:26 UTC (rev 11002) @@ -67,10 +67,10 @@ if (! nargout) mag_db = cellfun (@(mag) 20 * log10 (mag), mag, "uniformoutput", false); - w_cell = repmat ({w}, 1, numel (H)); + %w_cell = repmat ({w}, 1, numel (H)); - mag_args = vertcat (w_cell, mag_db)(:); - pha_args = vertcat (w_cell, pha)(:); + mag_args = vertcat (w, mag_db)(:); + pha_args = vertcat (w, pha)(:); %if (isct (sys)) xl_str = "Frequency [rad/s]"; Modified: trunk/octave-forge/main/control/devel/multiplot.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot.m 2012-09-13 08:09:07 UTC (rev 11001) +++ trunk/octave-forge/main/control/devel/multiplot.m 2012-09-13 14:34:26 UTC (rev 11002) @@ -1,3 +1,9 @@ load tfs.dat -bode2 (C_AH, C_opt) \ No newline at end of file +figure (1) +bode2 (C_AH, C_opt) + +figure (2) +bode2 (5*C_AH, frd (C_AH), frd (C_opt)) + +% bode2 (5*C_AH, frd (C_AH, 1:10), frd (C_opt, 11:20)) \ No newline at end of file This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-13 16:47:35
|
Revision: 11003 http://octave.svn.sourceforge.net/octave/?rev=11003&view=rev Author: paramaniac Date: 2012-09-13 16:47:24 +0000 (Thu, 13 Sep 2012) Log Message: ----------- control: support plot styles (very basic) Modified Paths: -------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/bode2.m trunk/octave-forge/main/control/devel/multiplot.m Modified: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-13 14:34:26 UTC (rev 11002) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-13 16:47:24 UTC (rev 11003) @@ -21,11 +21,14 @@ ## Author: Lukas Reichlin <luk...@gm...> ## Created: November 2009 -## Version: 0.4 +## Version: 0.5 -% function [H, w] = __frequency_response__ (sys, w = [], mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) function [H, w] = __frequency_response_2__ (args, mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) + %if (! iscell (args)) + % args = {args}; + %endif + sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models w_idx = cellfun (@is_real_vector, args); # look for frequency vectors r_idx = cellfun (@iscell, args); # look for frequency ranges {wmin, wmax} @@ -33,6 +36,12 @@ sys_cell = args(sys_idx); # extract LTI models frd_idx = cellfun (@isa, sys_cell, {"frd"}); # look for FRD models + ## check arguments + if (! mimoflag && ! all (cellfun (@issiso, sys_cell))) + error ("frequency_response: require SISO systems"); + endif + + ## determine frequencies if (any (r_idx)) # if there are frequency ranges r = args(r_idx){end}; # take the last one if (numel (r) == 2 && issample (r{1}) && issample (r{2})) @@ -55,11 +64,7 @@ ## save frequency vectors of FRD models in w w_frd = cellfun (@get, sys_cell(frd_idx), {"w"}, "uniformoutput", false); w(frd_idx) = w_frd; -%w -%w = get (sys, "w"); -%args{sys_idx} - %{ ## check arguments if(! isa (sys, "lti")) Modified: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-13 14:34:26 UTC (rev 11002) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-13 16:47:24 UTC (rev 11003) @@ -48,12 +48,10 @@ ## Author: Lukas Reichlin <luk...@gm...> ## Created: November 2009 -## Version: 0.4 +## Version: 0.5 function [mag_r, pha_r, w_r] = bode2 (varargin) - ## TODO: multiplot feature: bode (sys1, "b", sys2, "r", ...) - if (nargin == 0) print_usage (); endif @@ -67,23 +65,23 @@ if (! nargout) mag_db = cellfun (@(mag) 20 * log10 (mag), mag, "uniformoutput", false); - %w_cell = repmat ({w}, 1, numel (H)); + style = repmat ({""}, 1, numel (H)); - mag_args = vertcat (w, mag_db)(:); - pha_args = vertcat (w, pha)(:); - - %if (isct (sys)) - xl_str = "Frequency [rad/s]"; - %else - % xl_str = sprintf ("Frequency [rad/s] w_N = %g", pi / get (sys, "tsam")); - %endif + tmp = cellfun (@ischar, varargin); + char_idx = find (tmp); + char_idx = char_idx(1:min (numel (H), end)); + style(1:length (char_idx)) = varargin(char_idx); + + mag_args = vertcat (w, mag_db, style)(:); + pha_args = vertcat (w, pha, style)(:); + subplot (2, 1, 1) semilogx (mag_args{:}) axis ("tight") ylim (__axis_margin__ (ylim)) grid ("on") - title (["Bode Diagram of ", inputname(1)]) + title ("Bode Diagram") ylabel ("Magnitude [dB]") subplot (2, 1, 2) @@ -91,7 +89,7 @@ axis ("tight") ylim (__axis_margin__ (ylim)) grid ("on") - xlabel (xl_str) + xlabel ("Frequency [rad/s]") ylabel ("Phase [deg]") else mag_r = mag{1}; Modified: trunk/octave-forge/main/control/devel/multiplot.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot.m 2012-09-13 14:34:26 UTC (rev 11002) +++ trunk/octave-forge/main/control/devel/multiplot.m 2012-09-13 16:47:24 UTC (rev 11003) @@ -6,4 +6,7 @@ figure (2) bode2 (5*C_AH, frd (C_AH), frd (C_opt)) -% bode2 (5*C_AH, frd (C_AH, 1:10), frd (C_opt, 11:20)) \ No newline at end of file +% bode2 (5*C_AH, frd (C_AH, 1:10), frd (C_opt, 11:20)) + +figure (3) +bode2 (5*C_AH, '*r', C_AH, ':b', C_opt, '-.k') \ No newline at end of file This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-14 08:24:42
|
Revision: 11006 http://octave.svn.sourceforge.net/octave/?rev=11006&view=rev Author: paramaniac Date: 2012-09-14 08:24:31 +0000 (Fri, 14 Sep 2012) Log Message: ----------- control: handle plots with implicit frequencies in backend of multiplot feature Modified Paths: -------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/__frequency_vector_2__.m Modified: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-13 18:01:43 UTC (rev 11005) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-14 08:24:31 UTC (rev 11006) @@ -55,44 +55,13 @@ w = __frequency_vector_2__ (sys_cell, wbounds); endif - w = repmat ({w}, 1, numel (sys_cell)); # return cell of frequency vectors + w_frd = w(frd_idx); # temporarily save frequency vectors of FRD models w(frd_idx) = {[]}; # freqresp returns all frequencies of FRD models for w=[] ## compute frequency response H for all LTI models H = cellfun (@__freqresp__, sys_cell, w, {resptype}, {cellflag}, "uniformoutput", false); - ## save frequency vectors of FRD models in w - w_frd = cellfun (@get, sys_cell(frd_idx), {"w"}, "uniformoutput", false); + ## restore frequency vectors of FRD models in w w(frd_idx) = w_frd; -%{ - ## check arguments - if(! isa (sys, "lti")) - error ("frequency_response: first argument sys must be an LTI system"); - endif - - if (! mimoflag && ! issiso (sys)) - error ("frequency_response: require SISO system"); - endif - - if (isa (sys, "frd")) - if (! isempty (w)) - warning ("frequency_response: second argument w is ignored"); - endif - w = get (sys, "w"); - H = __freqresp__ (sys, [], resptype, cellflag); - elseif (isempty (w)) # find interesting frequency range w if not specified - w = __frequency_vector__ (sys, wbounds); - H = __freqresp__ (sys, w, resptype, cellflag); - elseif (iscell (w) && numel (w) == 2 && issample (w{1}) && issample (w{2})) - w = __frequency_vector__ (sys, wbounds, w{1}, w{2}); - H = __freqresp__ (sys, w, resptype, cellflag); - elseif (! is_real_vector (w)) - error ("frequency_response: second argument w must be a vector of frequencies"); - else - H = __freqresp__ (sys, w, resptype, cellflag); - endif - -%} - endfunction Modified: trunk/octave-forge/main/control/devel/__frequency_vector_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-13 18:01:43 UTC (rev 11005) +++ trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-14 08:24:31 UTC (rev 11006) @@ -35,44 +35,68 @@ ## Adapted-By: Lukas Reichlin <luk...@gm...> ## Date: October 2009 -## Version: 0.2 +## Version: 0.3 function w = __frequency_vector_2__ (sys_cell, wbounds = "std", wmin, wmax) - if (! iscell (sys_cell)) + isc = iscell (sys_cell); + + if (! isc) # __sys2frd__ methods pass LTI models not in cells sys_cell = {sys_cell} endif idx = cellfun (@(x) isa (x, "lti"), sys_cell); sys_cell = sys_cell(idx); + len = numel (sys_cell); [dec_min, dec_max, zp] = cellfun (@__frequency_range__, sys_cell, {wbounds}, "uniformoutput", false); - if (nargin == 2) - dec_min = min (cell2mat (dec_min)); - dec_max = max (cell2mat (dec_max)); - elseif (nargin == 4) # w = {wmin, wmax} - dec_min = log10 (wmin); - dec_max = log10 (wmax); + if (strcmpi (wbounds, "std")) # plots with explicit frequencies + + if (nargin == 2) + dec_min = min (cell2mat (dec_min)); + dec_max = max (cell2mat (dec_max)); + elseif (nargin == 4) # w = {wmin, wmax} + dec_min = log10 (wmin); + dec_max = log10 (wmax); + else + print_usage (); + endif + + zp = horzcat (zp{:}); + + ## include zeros and poles for nice peaks in plots + idx = find (zp > 10^dec_min & zp < 10^dec_max); + zp = zp(idx); + + w = logspace (dec_min, dec_max, 500); + w = unique ([w, zp]); # unique also sorts frequency vector + + w = repmat ({w}, 1, len); # return cell of frequency vectors + + elseif (strcmpi (wbounds, "ext")) # plots with implicit frequencies + + if (nargin == 4) + dec_min = repmat ({log10 (wmin)}, 1, len); + dec_max = repmat ({log10 (wmax)}, 1, len); + endif + + idx = cellfun (@(zp, dec_min, dec_max) find (zp > 10^dec_min & zp < 10^dec_max), \ + zp, dec_min, dec_max, "uniformoutput", false); + zp = cellfun (@(zp, idx) zp(idx), zp, idx, "uniformoutput", false); + + w = cellfun (@logspace, dec_min, dec_max, {500}, "uniformoutput", false); + w = cellfun (@(w, zp) unique ([w, zp]), w, zp, "uniformoutput", false); + ## unique also sorts frequency vector + else - print_usage (); + error ("frequency_vector: invalid argument 'wbounds'"); endif - zp = horzcat (zp{:}); - - ## include zeros and poles for nice peaks in plots - idx = find (zp > 10^dec_min & zp < 10^dec_max); - zp = zp(idx); + if (! isc) # for __sys2frd__ methods + w = w{1}; + endif - w = logspace (dec_min, dec_max, 500); - w = unique ([w, zp]); # unique also sorts frequency vector - - ## TODO: nyquist diagrams may need individual dec_min and dec_max - ## if curve goes to infinity - - ## TODO: handle FRD models (they have fixed w), maybe frequency_response - ## is a better place for this? - endfunction @@ -88,8 +112,8 @@ zer = zero (sys); pol = pole (sys); - tsam = abs (get (sys, "tsam")); # tsam could be -1 - discrete = ! isct (sys); # static gains (tsam = -2) are assumed continuous + tsam = abs (get (sys, "tsam")); # tsam could be -1 + discrete = ! isct (sys); # static gains (tsam = -2) are assumed continuous ## make sure zer, pol are row vectors pol = reshape (pol, 1, []); @@ -132,8 +156,8 @@ if (isempty (iip) && isempty (iiz)) ## no poles/zeros away from omega = 0; pick defaults - dec_min = 0; # -1 - dec_max = 2; # 3 + dec_min = 0; # -1 + dec_max = 2; # 3 else dec_min = floor (log10 (min (abs ([cpol, czer])))); dec_max = ceil (log10 (max (abs ([cpol, czer])))); @@ -141,7 +165,7 @@ ## expand to show the entirety of the "interesting" portion of the plot switch (wbounds) - case "std" # standard + case "std" # standard if (dec_min == dec_max) dec_min -= 2; dec_max += 2; @@ -149,8 +173,8 @@ dec_min--; dec_max++; endif - case "ext" # extended (for nyquist) - if (any (abs (pol) < sqrt (eps))) # look for integrators + case "ext" # extended (for nyquist) + if (any (abs (pol) < sqrt (eps))) # look for integrators ## dec_min -= 0.5; dec_max += 2; else This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-14 09:45:26
|
Revision: 11008 http://octave.svn.sourceforge.net/octave/?rev=11008&view=rev Author: paramaniac Date: 2012-09-14 09:45:19 +0000 (Fri, 14 Sep 2012) Log Message: ----------- control: quicksave multiplot nyquist (2) Modified Paths: -------------- trunk/octave-forge/main/control/devel/multiplot.m trunk/octave-forge/main/control/devel/nyquist2.m Modified: trunk/octave-forge/main/control/devel/multiplot.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot.m 2012-09-14 08:57:52 UTC (rev 11007) +++ trunk/octave-forge/main/control/devel/multiplot.m 2012-09-14 09:45:19 UTC (rev 11008) @@ -13,4 +13,8 @@ figure (4) -nyquist2 (C_AH, C_opt) \ No newline at end of file +nyquist2 (C_AH, C_opt) + + +figure (5) +nyquist2 (C_AH, "r-", C_opt, "b:") Modified: trunk/octave-forge/main/control/devel/nyquist2.m =================================================================== --- trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 08:57:52 UTC (rev 11007) +++ trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 09:45:19 UTC (rev 11008) @@ -64,9 +64,9 @@ if (! nargout) tmp = cellfun (@isa, varargin, {"lti"}); - sys_idx = find (tmp); + sys_idx = find (tmp) tmp = cellfun (@ischar, varargin); - style_idx = find (tmp); + style_idx = find (tmp) len = numel (H); plot_args_pos = {}; @@ -77,11 +77,17 @@ for k = 1:len col = colororder(1+rem (k-1, rc), :); - plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, {"-", "color", col}); - plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, {"-.", "color", col}); + if (k == len) + lim = nargin; + else + lim = sys_idx(k+1); + endif + style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)) + plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, {"-", "color", col}, style); + plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, {"-.", "color", col}, style); legend_args{k} = inputname(sys_idx(k)); endfor - +%plot_args_pos h = plot (plot_args_pos{:}, plot_args_neg{:}); axis ("tight") xlim (__axis_margin__ (xlim)) @@ -98,4 +104,4 @@ w_r = w{1}; endif -endfunction \ No newline at end of file +endfunction This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-14 10:07:40
|
Revision: 11009 http://octave.svn.sourceforge.net/octave/?rev=11009&view=rev Author: paramaniac Date: 2012-09-14 10:07:29 +0000 (Fri, 14 Sep 2012) Log Message: ----------- control: finish multiplot nyquist Modified Paths: -------------- trunk/octave-forge/main/control/devel/multiplot.m trunk/octave-forge/main/control/devel/nyquist2.m Modified: trunk/octave-forge/main/control/devel/multiplot.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot.m 2012-09-14 09:45:19 UTC (rev 11008) +++ trunk/octave-forge/main/control/devel/multiplot.m 2012-09-14 10:07:29 UTC (rev 11009) @@ -17,4 +17,5 @@ figure (5) -nyquist2 (C_AH, "r-", C_opt, "b:") +nyquist2 (C_AH, "xr", C_opt, "ob") +legend ("Test C_AH", "Test C_opt") \ No newline at end of file Modified: trunk/octave-forge/main/control/devel/nyquist2.m =================================================================== --- trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 09:45:19 UTC (rev 11008) +++ trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 10:07:29 UTC (rev 11009) @@ -64,9 +64,9 @@ if (! nargout) tmp = cellfun (@isa, varargin, {"lti"}); - sys_idx = find (tmp) + sys_idx = find (tmp); tmp = cellfun (@ischar, varargin); - style_idx = find (tmp) + style_idx = find (tmp); len = numel (H); plot_args_pos = {}; @@ -82,12 +82,17 @@ else lim = sys_idx(k+1); endif - style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)) - plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, {"-", "color", col}, style); - plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, {"-.", "color", col}, style); + style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); + if (isempty (style)) + plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, {"-", "color", col}); + plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, {"-.", "color", col}); + else + plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, style); + plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, style); + endif legend_args{k} = inputname(sys_idx(k)); endfor -%plot_args_pos + h = plot (plot_args_pos{:}, plot_args_neg{:}); axis ("tight") xlim (__axis_margin__ (xlim)) @@ -97,7 +102,6 @@ xlabel ("Real Axis") ylabel ("Imaginary Axis") legend (h(1:len), legend_args) - % legend (h(1:2:2*len), legend_args) else re_r = re{1}; im_r = im{1}; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-14 11:33:50
|
Revision: 11010 http://octave.svn.sourceforge.net/octave/?rev=11010&view=rev Author: paramaniac Date: 2012-09-14 11:33:41 +0000 (Fri, 14 Sep 2012) Log Message: ----------- control: finish multiplot bode Modified Paths: -------------- trunk/octave-forge/main/control/devel/bode2.m trunk/octave-forge/main/control/devel/nyquist2.m Modified: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-14 10:07:29 UTC (rev 11009) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-14 11:33:41 UTC (rev 11010) @@ -65,16 +65,27 @@ if (! nargout) mag_db = cellfun (@(mag) 20 * log10 (mag), mag, "uniformoutput", false); - style = repmat ({""}, 1, numel (H)); - + tmp = cellfun (@isa, varargin, {"lti"}); + sys_idx = find (tmp); tmp = cellfun (@ischar, varargin); - char_idx = find (tmp); - char_idx = char_idx(1:min (numel (H), end)); + style_idx = find (tmp); - style(1:length (char_idx)) = varargin(char_idx); + len = numel (H); + mag_args = {}; + pha_args = {}; + legend_args = cell (len, 1); - mag_args = vertcat (w, mag_db, style)(:); - pha_args = vertcat (w, pha, style)(:); + for k = 1:len + if (k == len) + lim = nargin; + else + lim = sys_idx(k+1); + endif + style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); + mag_args = cat (2, mag_args, w(k), mag_db(k), style); + pha_args = cat (2, pha_args, w(k), pha(k), style); + legend_args{k} = inputname(sys_idx(k)); + endfor subplot (2, 1, 1) semilogx (mag_args{:}) @@ -91,6 +102,7 @@ grid ("on") xlabel ("Frequency [rad/s]") ylabel ("Phase [deg]") + legend (legend_args) else mag_r = mag{1}; pha_r = pha{1}; Modified: trunk/octave-forge/main/control/devel/nyquist2.m =================================================================== --- trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 10:07:29 UTC (rev 11009) +++ trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 11:33:41 UTC (rev 11010) @@ -69,8 +69,8 @@ style_idx = find (tmp); len = numel (H); - plot_args_pos = {}; - plot_args_neg = {}; + pos_args = {}; + neg_args = {}; legend_args = cell (len, 1); colororder = get (gca, "colororder"); rc = rows (colororder); @@ -84,16 +84,20 @@ endif style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); if (isempty (style)) - plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, {"-", "color", col}); - plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, {"-.", "color", col}); + pos_args = cat (2, pos_args, re{k}, im{k}, {"-", "color", col}); + neg_args = cat (2, neg_args, re{k}, -im{k}, {"-.", "color", col}); else - plot_args_pos = cat (2, plot_args_pos, re{k}, im{k}, style); - plot_args_neg = cat (2, plot_args_neg, re{k}, -im{k}, style); + pos_args = cat (2, pos_args, re{k}, im{k}, style); + neg_args = cat (2, neg_args, re{k}, -im{k}, style); endif legend_args{k} = inputname(sys_idx(k)); endfor + + ## FIXME: pos_args = cat (2, pos_args, re{k}, im{k}, {"-", "color", col}, style); + ## doesn't work! it would be nice to have default arguments that can be + ## (partially) overwritten by user-specified plot styles. - h = plot (plot_args_pos{:}, plot_args_neg{:}); + h = plot (pos_args{:}, neg_args{:}); axis ("tight") xlim (__axis_margin__ (xlim)) ylim (__axis_margin__ (ylim)) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-14 11:44:46
|
Revision: 11012 http://octave.svn.sourceforge.net/octave/?rev=11012&view=rev Author: paramaniac Date: 2012-09-14 11:44:35 +0000 (Fri, 14 Sep 2012) Log Message: ----------- control: cleaning up Modified Paths: -------------- trunk/octave-forge/main/control/devel/multiplot.m trunk/octave-forge/main/control/devel/multiplot2.m Removed Paths: ------------- trunk/octave-forge/main/control/devel/__frequency_response_2__.m trunk/octave-forge/main/control/devel/__frequency_vector_2__.m trunk/octave-forge/main/control/devel/bode2.m trunk/octave-forge/main/control/devel/nyquist2.m Deleted: trunk/octave-forge/main/control/devel/__frequency_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-14 11:40:39 UTC (rev 11011) +++ trunk/octave-forge/main/control/devel/__frequency_response_2__.m 2012-09-14 11:44:35 UTC (rev 11012) @@ -1,67 +0,0 @@ -## Copyright (C) 2009, 2010, 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 <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## Return frequency response H and frequency vector w. -## If w is empty, it will be calculated by __frequency_vector__. - -## Author: Lukas Reichlin <luk...@gm...> -## Created: November 2009 -## Version: 0.5 - -function [H, w] = __frequency_response_2__ (args, mimoflag = 0, resptype = 0, wbounds = "std", cellflag = false) - - %if (! iscell (args)) - % args = {args}; - %endif - - sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models - w_idx = cellfun (@is_real_vector, args); # look for frequency vectors - r_idx = cellfun (@iscell, args); # look for frequency ranges {wmin, wmax} - - sys_cell = args(sys_idx); # extract LTI models - frd_idx = cellfun (@isa, sys_cell, {"frd"}); # look for FRD models - - ## check arguments - if (! mimoflag && ! all (cellfun (@issiso, sys_cell))) - error ("frequency_response: require SISO systems"); - endif - - ## determine frequencies - if (any (r_idx)) # if there are frequency ranges - r = args(r_idx){end}; # take the last one - if (numel (r) == 2 && issample (r{1}) && issample (r{2})) - w = __frequency_vector_2__ (sys_cell, wbounds, r{1}, r{2}); - else - error ("frequency_response: invalid cell"); - endif - elseif (any (w_idx)) # are there any frequency vectors? - w = args(w_idx){end}; - else # there are neither frequency ranges nor vectors - w = __frequency_vector_2__ (sys_cell, wbounds); - endif - - w_frd = w(frd_idx); # temporarily save frequency vectors of FRD models - w(frd_idx) = {[]}; # freqresp returns all frequencies of FRD models for w=[] - - ## compute frequency response H for all LTI models - H = cellfun (@__freqresp__, sys_cell, w, {resptype}, {cellflag}, "uniformoutput", false); - - ## restore frequency vectors of FRD models in w - w(frd_idx) = w_frd; - -endfunction Deleted: trunk/octave-forge/main/control/devel/__frequency_vector_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-14 11:40:39 UTC (rev 11011) +++ trunk/octave-forge/main/control/devel/__frequency_vector_2__.m 2012-09-14 11:44:35 UTC (rev 11012) @@ -1,196 +0,0 @@ -## Copyright (C) 1996, 2000, 2004, 2005, 2006, 2007 -## Auburn University. All rights reserved. -## Copyright (C) 2009 - 2012 Lukas F. Reichlin -## -## This program 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. -## -## This program 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 this program; see the file COPYING. If not, see -## <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## @deftypefn {Function File} {@var{w} =} __frequency_vector__ (@var{sys}) -## Get default range of frequencies based on cutoff frequencies of system -## poles and zeros. -## Frequency range is the interval -## @iftex -## @tex -## $ [ 10^{w_{min}}, 10^{w_{max}} ] $ -## @end tex -## @end iftex -## @ifnottex -## [10^@var{wmin}, 10^@var{wmax}] -## @end ifnottex -## -## Used by @command{__frequency_response__} -## @end deftypefn - -## Adapted-By: Lukas Reichlin <luk...@gm...> -## Date: October 2009 -## Version: 0.3 - -function w = __frequency_vector_2__ (sys_cell, wbounds = "std", wmin, wmax) - - isc = iscell (sys_cell); - - if (! isc) # __sys2frd__ methods pass LTI models not in cells - sys_cell = {sys_cell} - endif - - idx = cellfun (@(x) isa (x, "lti"), sys_cell); - sys_cell = sys_cell(idx); - len = numel (sys_cell); - - [dec_min, dec_max, zp] = cellfun (@__frequency_range__, sys_cell, {wbounds}, "uniformoutput", false); - - if (strcmpi (wbounds, "std")) # plots with explicit frequencies - - if (nargin == 2) - dec_min = min (cell2mat (dec_min)); - dec_max = max (cell2mat (dec_max)); - elseif (nargin == 4) # w = {wmin, wmax} - dec_min = log10 (wmin); - dec_max = log10 (wmax); - else - print_usage (); - endif - - zp = horzcat (zp{:}); - - ## include zeros and poles for nice peaks in plots - idx = find (zp > 10^dec_min & zp < 10^dec_max); - zp = zp(idx); - - w = logspace (dec_min, dec_max, 500); - w = unique ([w, zp]); # unique also sorts frequency vector - - w = repmat ({w}, 1, len); # return cell of frequency vectors - - elseif (strcmpi (wbounds, "ext")) # plots with implicit frequencies - - if (nargin == 4) - dec_min = repmat ({log10 (wmin)}, 1, len); - dec_max = repmat ({log10 (wmax)}, 1, len); - endif - - idx = cellfun (@(zp, dec_min, dec_max) find (zp > 10^dec_min & zp < 10^dec_max), \ - zp, dec_min, dec_max, "uniformoutput", false); - zp = cellfun (@(zp, idx) zp(idx), zp, idx, "uniformoutput", false); - - w = cellfun (@logspace, dec_min, dec_max, {500}, "uniformoutput", false); - w = cellfun (@(w, zp) unique ([w, zp]), w, zp, "uniformoutput", false); - ## unique also sorts frequency vector - - else - error ("frequency_vector: invalid argument 'wbounds'"); - endif - - if (! isc) # for __sys2frd__ methods - w = w{1}; - endif - -endfunction - - -function [dec_min, dec_max, zp] = __frequency_range__ (sys, wbounds = "std") - - if (isa (sys, "frd")) - w = get (sys, "w"); - dec_min = log10 (w(1)); - dec_max = log10 (w(end)); - zp = []; - return; - endif - - zer = zero (sys); - pol = pole (sys); - tsam = abs (get (sys, "tsam")); # tsam could be -1 - discrete = ! isct (sys); # static gains (tsam = -2) are assumed continuous - - ## make sure zer, pol are row vectors - pol = reshape (pol, 1, []); - zer = reshape (zer, 1, []); - - ## check for natural frequencies away from omega = 0 - if (discrete) - ## The 2nd conditions prevents log(0) in the next log command - iiz = find (abs(zer-1) > norm(zer)*eps && abs(zer) > norm(zer)*eps); - iip = find (abs(pol-1) > norm(pol)*eps && abs(pol) > norm(pol)*eps); - - ## avoid dividing empty matrices, it would work but looks nasty - if (! isempty (iiz)) - czer = log (zer(iiz))/tsam; - else - czer = []; - endif - - if (! isempty (iip)) - cpol = log (pol(iip))/tsam; - else - cpol = []; - endif - else - ## continuous - iip = find (abs(pol) > norm(pol)*eps); - iiz = find (abs(zer) > norm(zer)*eps); - - if (! isempty (zer)) - czer = zer(iiz); - else - czer = []; - endif - if (! isempty (pol)) - cpol = pol(iip); - else - cpol = []; - endif - endif - - if (isempty (iip) && isempty (iiz)) - ## no poles/zeros away from omega = 0; pick defaults - dec_min = 0; # -1 - dec_max = 2; # 3 - else - dec_min = floor (log10 (min (abs ([cpol, czer])))); - dec_max = ceil (log10 (max (abs ([cpol, czer])))); - endif - - ## expand to show the entirety of the "interesting" portion of the plot - switch (wbounds) - case "std" # standard - if (dec_min == dec_max) - dec_min -= 2; - dec_max += 2; - else - dec_min--; - dec_max++; - endif - case "ext" # extended (for nyquist) - if (any (abs (pol) < sqrt (eps))) # look for integrators - ## dec_min -= 0.5; - dec_max += 2; - else - dec_min -= 2; - dec_max += 2; - endif - otherwise - error ("frequency_range: second argument invalid"); - endswitch - - ## run discrete frequency all the way to pi - if (discrete) - dec_max = log10 (pi/tsam); - endif - - ## include zeros and poles for nice peaks in plots - zp = [abs(zer), abs(pol)]; - -endfunction Deleted: trunk/octave-forge/main/control/devel/bode2.m =================================================================== --- trunk/octave-forge/main/control/devel/bode2.m 2012-09-14 11:40:39 UTC (rev 11011) +++ trunk/octave-forge/main/control/devel/bode2.m 2012-09-14 11:44:35 UTC (rev 11012) @@ -1,112 +0,0 @@ -## Copyright (C) 2009, 2010, 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 <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## @deftypefn {Function File} {[@var{mag}, @var{pha}, @var{w}] =} bode (@var{sys}) -## @deftypefnx {Function File} {[@var{mag}, @var{pha}, @var{w}] =} bode (@var{sys}, @var{w}) -## Bode diagram of frequency response. If no output arguments are given, -## the response is printed on the screen. -## -## @strong{Inputs} -## @table @var -## @item sys -## LTI system. Must be a single-input and single-output (SISO) system. -## @item w -## Optional vector of frequency values. If @var{w} is not specified, -## it is calculated by the zeros and poles of the system. -## Alternatively, the cell @code{@{wmin, wmax@}} specifies a frequency range, -## where @var{wmin} and @var{wmax} denote minimum and maximum frequencies -## in rad/s. -## @end table -## -## @strong{Outputs} -## @table @var -## @item mag -## Vector of magnitude. Has length of frequency vector @var{w}. -## @item pha -## Vector of phase. Has length of frequency vector @var{w}. -## @item w -## Vector of frequency values used. -## @end table -## -## @seealso{nichols, nyquist, sigma} -## @end deftypefn - -## Author: Lukas Reichlin <luk...@gm...> -## Created: November 2009 -## Version: 0.5 - -function [mag_r, pha_r, w_r] = bode2 (varargin) - - if (nargin == 0) - print_usage (); - endif - - [H, w] = __frequency_response_2__ (varargin, false, 0, "std", false); - - H = cellfun (@reshape, H, {[]}, {1}, "uniformoutput", false); - mag = cellfun (@abs, H, "uniformoutput", false); - pha = cellfun (@(H) unwrap (arg (H)) * 180 / pi, H, "uniformoutput", false); - - if (! nargout) - mag_db = cellfun (@(mag) 20 * log10 (mag), mag, "uniformoutput", false); - - tmp = cellfun (@isa, varargin, {"lti"}); - sys_idx = find (tmp); - tmp = cellfun (@ischar, varargin); - style_idx = find (tmp); - - len = numel (H); - mag_args = {}; - pha_args = {}; - legend_args = cell (len, 1); - - for k = 1:len - if (k == len) - lim = nargin; - else - lim = sys_idx(k+1); - endif - style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); - mag_args = cat (2, mag_args, w(k), mag_db(k), style); - pha_args = cat (2, pha_args, w(k), pha(k), style); - legend_args{k} = inputname(sys_idx(k)); - endfor - - subplot (2, 1, 1) - semilogx (mag_args{:}) - axis ("tight") - ylim (__axis_margin__ (ylim)) - grid ("on") - title ("Bode Diagram") - ylabel ("Magnitude [dB]") - - subplot (2, 1, 2) - semilogx (pha_args{:}) - axis ("tight") - ylim (__axis_margin__ (ylim)) - grid ("on") - xlabel ("Frequency [rad/s]") - ylabel ("Phase [deg]") - legend (legend_args) - else - mag_r = mag{1}; - pha_r = pha{1}; - w_r = w{1}; - endif - -endfunction Modified: trunk/octave-forge/main/control/devel/multiplot.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot.m 2012-09-14 11:40:39 UTC (rev 11011) +++ trunk/octave-forge/main/control/devel/multiplot.m 2012-09-14 11:44:35 UTC (rev 11012) @@ -1,21 +1,21 @@ load tfs.dat figure (1) -bode2 (C_AH, C_opt) +bode (C_AH, C_opt) figure (2) -bode2 (5*C_AH, frd (C_AH), frd (C_opt)) +bode (5*C_AH, frd (C_AH), frd (C_opt)) % bode2 (5*C_AH, frd (C_AH, 1:10), frd (C_opt, 11:20)) figure (3) -bode2 (5*C_AH, '*r', C_AH, 'xb', C_opt, 'ok') +bode (5*C_AH, '*r', C_AH, 'xb', C_opt, 'ok') figure (4) -nyquist2 (C_AH, C_opt) +nyquist (C_AH, C_opt) figure (5) -nyquist2 (C_AH, "xr", C_opt, "ob") +nyquist (C_AH, "xr", C_opt, "ob") legend ("Test C_AH", "Test C_opt") \ No newline at end of file Modified: trunk/octave-forge/main/control/devel/multiplot2.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot2.m 2012-09-14 11:40:39 UTC (rev 11011) +++ trunk/octave-forge/main/control/devel/multiplot2.m 2012-09-14 11:44:35 UTC (rev 11012) @@ -90,8 +90,8 @@ % Bode Plot of Controller figure (1) -bode2 (K, Kr4, Kr2, w) -title ('Bode Diagrams of K and Kr') +bode (K, Kr4, Kr2, w) +% title ('Bode Diagrams of K and Kr') legend ('K (8 states)', 'Kr (4 states)', 'Kr (2 states)', 'location', 'southwest') %{ Deleted: trunk/octave-forge/main/control/devel/nyquist2.m =================================================================== --- trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 11:40:39 UTC (rev 11011) +++ trunk/octave-forge/main/control/devel/nyquist2.m 2012-09-14 11:44:35 UTC (rev 11012) @@ -1,115 +0,0 @@ -## Copyright (C) 2009, 2010, 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 <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## @deftypefn {Function File} {[@var{re}, @var{im}, @var{w}] =} nyquist (@var{sys}) -## @deftypefnx {Function File} {[@var{re}, @var{im}, @var{w}] =} nyquist (@var{sys}, @var{w}) -## Nyquist diagram of frequency response. If no output arguments are given, -## the response is printed on the screen. -## -## @strong{Inputs} -## @table @var -## @item sys -## LTI system. Must be a single-input and single-output (SISO) system. -## @item w -## Optional vector of frequency values. If @var{w} is not specified, -## it is calculated by the zeros and poles of the system. -## Alternatively, the cell @code{@{wmin, wmax@}} specifies a frequency range, -## where @var{wmin} and @var{wmax} denote minimum and maximum frequencies -## in rad/s. -## @end table -## -## @strong{Outputs} -## @table @var -## @item re -## Vector of real parts. Has length of frequency vector @var{w}. -## @item im -## Vector of imaginary parts. Has length of frequency vector @var{w}. -## @item w -## Vector of frequency values used. -## @end table -## -## @seealso{bode, nichols, sigma} -## @end deftypefn - -## Author: Lukas Reichlin <luk...@gm...> -## Created: November 2009 -## Version: 0.4 - -function [re_r, im_r, w_r] = nyquist2 (varargin) - - if (nargin == 0) - print_usage (); - endif - - [H, w] = __frequency_response_2__ (varargin, false, 0, "ext"); - - H = cellfun (@reshape, H, {[]}, {1}, "uniformoutput", false); - re = cellfun (@real, H, "uniformoutput", false); - im = cellfun (@imag, H, "uniformoutput", false); - - if (! nargout) - tmp = cellfun (@isa, varargin, {"lti"}); - sys_idx = find (tmp); - tmp = cellfun (@ischar, varargin); - style_idx = find (tmp); - - len = numel (H); - pos_args = {}; - neg_args = {}; - legend_args = cell (len, 1); - colororder = get (gca, "colororder"); - rc = rows (colororder); - - for k = 1:len - col = colororder(1+rem (k-1, rc), :); - if (k == len) - lim = nargin; - else - lim = sys_idx(k+1); - endif - style = varargin(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); - if (isempty (style)) - pos_args = cat (2, pos_args, re{k}, im{k}, {"-", "color", col}); - neg_args = cat (2, neg_args, re{k}, -im{k}, {"-.", "color", col}); - else - pos_args = cat (2, pos_args, re{k}, im{k}, style); - neg_args = cat (2, neg_args, re{k}, -im{k}, style); - endif - legend_args{k} = inputname(sys_idx(k)); - endfor - - ## FIXME: pos_args = cat (2, pos_args, re{k}, im{k}, {"-", "color", col}, style); - ## doesn't work! it would be nice to have default arguments that can be - ## (partially) overwritten by user-specified plot styles. - - h = plot (pos_args{:}, neg_args{:}); - axis ("tight") - xlim (__axis_margin__ (xlim)) - ylim (__axis_margin__ (ylim)) - grid ("on") - title ("Nyquist Diagram") - xlabel ("Real Axis") - ylabel ("Imaginary Axis") - legend (h(1:len), legend_args) - else - re_r = re{1}; - im_r = im{1}; - w_r = w{1}; - endif - -endfunction This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-17 15:10:05
|
Revision: 11040 http://octave.svn.sourceforge.net/octave/?rev=11040&view=rev Author: paramaniac Date: 2012-09-17 15:09:54 +0000 (Mon, 17 Sep 2012) Log Message: ----------- control: get draft code for multiplot step response working Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m trunk/octave-forge/main/control/devel/step2.m Added Paths: ----------- trunk/octave-forge/main/control/devel/multiplot3.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-17 05:51:45 UTC (rev 11039) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-17 15:09:54 UTC (rev 11040) @@ -23,10 +23,10 @@ ## Version: 0.3 % function [y, t, x_arr] = __time_response_2__ (sys, resptype, plotflag, tfinal, dt, x0, sysname) -function [y, t, x] = __time_response_2__ (resptype, args) +function [y, t, x] = __time_response_2__ (resptype, args, plotflag) sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models - sys_cell = cellfun (@ss, args(sys_idx)); # system must be proper + sys_cell = cellfun (@ss, args(sys_idx), "uniformoutput", false); # convert to state-space if (! size_equal (sys_cell{:})) error ("%s: models must have equal sizes", resptype); @@ -37,25 +37,28 @@ n_vec = length (vec_idx) n_sys = length (sys_cell) - if (n_vec >= 1) - arg = args{vec_idx(1)}; - - endif + %if (n_vec >= 1) + % arg = args{vec_idx(1)}; + % + %endif ## extract tfinal/t, dt, x0 - [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, dt); + tfinal = []; + dt = {[]}; + + [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, dt, "uniformoutput", false); - tfinal = max (tfinal); + tfinal = max ([tfinal{:}]); % __sim_horizon__ (sys, tfinal, dt); - hier sim_horizon + %hier sim_horizon - ct_idx = cellfun (@isct, sys_cell) + ct_idx = cellfun (@isct, sys_cell); sys_dt_cell = sys_cell; - tmp = (@c2d, sys_cell(ct_idx), dt, {"zoh"}, "uniformoutput", false) + tmp = cellfun (@c2d, sys_cell(ct_idx), dt, {"zoh"}, "uniformoutput", false); sys_dt_cell(ct_idx) = tmp; %{ @@ -98,13 +101,18 @@ l_t = length (t); %} + + ## time vector + t = @cellfun (@(dt) reshape (0 : dt : tfinal, [], 1), dt, "uniformoutput", false); + + %function [y, x_arr] = __initial_response__ (sys, sys_dt, t, x0) %function [y, x_arr] = __step_response__ (sys_dt, t) %function [y, x_arr] = __impulse_response__ (sys, sys_dt, t) switch (resptype) case "initial" - [y, x] = cellfun (@__initial_response__, sys_dt_cell, t, {x0} or x0, "uniformoutput", false); + %[y, x] = cellfun (@__initial_response__, sys_dt_cell, t, {x0} or x0, "uniformoutput", false); case "step" [y, x] = cellfun (@__step_response__, sys_dt_cell, t, "uniformoutput", false); case "impulse" @@ -114,9 +122,8 @@ endswitch - - if (plotflag) # display plot + [p, m] = size (sys_cell{1}); switch (resptype) case "initial" str = "Response to Initial Conditions"; @@ -130,23 +137,24 @@ otherwise error ("time_response: invalid response type"); endswitch + + outname = get (sys_cell{end}, "outname"); + outname = __labels__ (outname, "y"); for i = 1 : n_sys - t = t{i}; - y = y{i}; discrete = ! ct_idx(i); if (discrete) # discrete system for k = 1 : p for j = 1 : cols subplot (p, cols, (k-1)*cols+j); - stairs (t, y(:, k, j)); + stairs (t{i}, y{i}(:, k, j)); hold on; grid ("on"); if (i == n_sys && k == 1 && j == 1) title (str); endif if (i == n_sys && j == 1) - ylabel (sprintf ("Amplitude %s", outname{k})); + ylabel (outname{k}); endif endfor endfor @@ -154,14 +162,14 @@ for k = 1 : p for j = 1 : cols subplot (p, cols, (k-1)*cols+j); - plot (t, y(:, k, j)); + plot (t{i}, y{i}(:, k, j)); hold on; grid ("on"); if (i == n_sys && k == 1 && j == 1) title (str); endif if (i == n_sys && j == 1) - ylabel (sprintf ("Amplitude %s", outname{k})); + ylabel (outname{k}); endif endfor endfor @@ -171,89 +179,14 @@ hold off; endif -%{ - if (plotflag) # display plot - - ## TODO: Set correct titles, especially for multi-input systems - - stable = isstable (sys); - outname = get (sys, "outname"); - outname = __labels__ (outname, "y_"); - - if (strcmp (resptype, "initial")) - cols = 1; - else - cols = m; - endif - - if (discrete) # discrete system - for k = 1 : p - for j = 1 : cols - - subplot (p, cols, (k-1)*cols+j); - - if (stable) - stairs (t, [y(:, k, j), yfinal(k, j) * ones(l_t, 1)]); - else - stairs (t, y(:, k, j)); - endif - - grid ("on"); - - if (k == 1 && j == 1) - title (str); - endif - - if (j == 1) - ylabel (sprintf ("Amplitude %s", outname{k})); - endif - - endfor - endfor - - xlabel ("Time [s]"); - - else # continuous system - for k = 1 : p - for j = 1 : cols - - subplot (p, cols, (k-1)*cols+j); - - if (stable) - plot (t, [y(:, k, j), yfinal(k, j) * ones(l_t, 1)]); - else - plot (t, y(:, k, j)); - endif - - grid ("on"); - - if (k == 1 && j == 1) - title (str); - endif - - if (j == 1) - ylabel (sprintf ("Amplitude %s", outname{k})); - endif - - endfor - endfor - - xlabel ("Time [s]"); - - endif - endif - endfunction -%} -endfunction - function [y, x_arr] = __initial_response__ (sys_dt, t, x0) - [F, G, C, D] = ssdata (sys_dt); + [F, G, C, D] = ssdata (sys_dt); # system must be proper n = rows (F); # number of states m = columns (G); # number of inputs @@ -283,7 +216,7 @@ function [y, x_arr] = __step_response__ (sys_dt, t) - [F, G, C, D] = ssdata (sys_dt); + [F, G, C, D] = ssdata (sys_dt); # system must be proper n = rows (F); # number of states m = columns (G); # number of inputs @@ -313,8 +246,8 @@ function [y, x_arr] = __impulse_response__ (sys, sys_dt, t) - [~, B, ~, ~, dt] = ssdata (sys); - [F, G, C, D] = ssdata (sys_dt); + [~, B] = ssdata (sys); + [F, G, C, D, dt] = ssdata (sys_dt); # system must be proper discrete = ! isct (sys); n = rows (F); # number of states @@ -359,10 +292,7 @@ endfunction -function - - % function [tfinal, dt] = __sim_horizon__ (A, discrete, tfinal, Ts) function [tfinal, dt] = __sim_horizon__ (sys, tfinal, Ts) @@ -376,6 +306,7 @@ ev = pole (sys); n = length (ev); + discrete = ! isct (sys); if (discrete) ## perform bilinear transformation on poles in z Added: trunk/octave-forge/main/control/devel/multiplot3.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot3.m (rev 0) +++ trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-17 15:09:54 UTC (rev 11040) @@ -0,0 +1,4 @@ +load tfs.dat + +figure (1) +step2 (T_AH, T_opt) Modified: trunk/octave-forge/main/control/devel/step2.m =================================================================== --- trunk/octave-forge/main/control/devel/step2.m 2012-09-17 05:51:45 UTC (rev 11039) +++ trunk/octave-forge/main/control/devel/step2.m 2012-09-17 15:09:54 UTC (rev 11040) @@ -69,7 +69,7 @@ %sys_names = arrayfun (@inputname, sys_idx); - [y, t, x] = __time_response__ ("step", varargin, ! nargout); + [y, t, x] = __time_response_2__ ("step", varargin, ! nargout); if (nargout) y_r = y{1}; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-17 20:17:24
|
Revision: 11041 http://octave.svn.sourceforge.net/octave/?rev=11041&view=rev Author: paramaniac Date: 2012-09-17 20:17:18 +0000 (Mon, 17 Sep 2012) Log Message: ----------- control: tighten axes Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m trunk/octave-forge/main/control/devel/multiplot2.m trunk/octave-forge/main/control/devel/multiplot3.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-17 15:09:54 UTC (rev 11040) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-17 20:17:18 UTC (rev 11041) @@ -34,8 +34,8 @@ tmp = cellfun (@is_real_matrix, args); vec_idx = find (tmp); - n_vec = length (vec_idx) - n_sys = length (sys_cell) + n_vec = length (vec_idx); + n_sys = length (sys_cell); %if (n_vec >= 1) % arg = args{vec_idx(1)}; @@ -47,14 +47,14 @@ tfinal = []; dt = {[]}; - [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, dt, "uniformoutput", false); + %[tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, dt, "uniformoutput", false); + + [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, "uniformoutput", false); tfinal = max ([tfinal{:}]); - % __sim_horizon__ (sys, tfinal, dt); - - - %hier sim_horizon +dt + ct_idx = cellfun (@isct, sys_cell); sys_dt_cell = sys_cell; @@ -143,34 +143,42 @@ for i = 1 : n_sys discrete = ! ct_idx(i); - if (discrete) # discrete system + if (discrete) # discrete-time system for k = 1 : p for j = 1 : cols subplot (p, cols, (k-1)*cols+j); stairs (t{i}, y{i}(:, k, j)); hold on; - grid ("on"); - if (i == n_sys && k == 1 && j == 1) - title (str); + grid on; + if (i == n_sys) + axis tight; + ylim (__axis_margin__ (ylim)) + if (j == 1) + ylabel (outname{k}); + if (k == 1) + title (str); + endif + endif endif - if (i == n_sys && j == 1) - ylabel (outname{k}); - endif endfor endfor - else # continuous system + else # continuous-time system for k = 1 : p for j = 1 : cols subplot (p, cols, (k-1)*cols+j); plot (t{i}, y{i}(:, k, j)); hold on; - grid ("on"); - if (i == n_sys && k == 1 && j == 1) - title (str); + grid on; + if (i == n_sys) + axis tight + ylim (__axis_margin__ (ylim)) + if (j == 1) + ylabel (outname{k}); + if (k == 1) + title (str); + endif + endif endif - if (i == n_sys && j == 1) - ylabel (outname{k}); - endif endfor endfor endif @@ -305,8 +313,10 @@ T_DEF = 10; # default simulation time ev = pole (sys); - n = length (ev); - discrete = ! isct (sys); + n = length (ev); # number of states/poles + continuous = isct (sys); + discrete = ! continuous; + Ts = get (sys, "tsam"); if (discrete) ## perform bilinear transformation on poles in z @@ -334,14 +344,14 @@ tfinal = T_DEF; endif - if (! discrete) + if (continuous) dt = tfinal / N_DEF; endif else ev = ev(find (ev)); ev_max = max (abs (ev)); - if (! discrete) + if (continuous) dt = 0.2 * pi / ev_max; endif @@ -354,7 +364,7 @@ tfinal = yy * ceil (tfinal / yy); endif - if (! discrete) + if (continuous) N = tfinal / dt; if (N < N_MIN) @@ -367,8 +377,8 @@ endif endif - if (! isempty (Ts)) # catch case cont. system with dt specified - dt = Ts; - endif + %if (! isempty (Ts)) # catch case cont. system with dt specified + % dt = Ts; + %endif endfunction Modified: trunk/octave-forge/main/control/devel/multiplot2.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot2.m 2012-09-17 15:09:54 UTC (rev 11040) +++ trunk/octave-forge/main/control/devel/multiplot2.m 2012-09-17 20:17:18 UTC (rev 11041) @@ -72,7 +72,6 @@ Kr4 = spaconred (P, K, 4, 'feedback', '-') Kr2 = spaconred (P, K, 2, 'feedback', '-') -%{ % Open Loop L = P * K; Lr4 = P * Kr4; @@ -82,7 +81,6 @@ T = feedback (L); Tr4 = feedback (Lr4); Tr2 = feedback (Lr2); -%} % Frequency Range w = {1e-2, 1e1}; @@ -94,6 +92,10 @@ % title ('Bode Diagrams of K and Kr') legend ('K (8 states)', 'Kr (4 states)', 'Kr (2 states)', 'location', 'southwest') +% Step Response of Closed Loop +figure (2) +step2 (T, Tr4, Tr2) + %{ % Step Response of Closed Loop [y, t] = step (T, 100); Modified: trunk/octave-forge/main/control/devel/multiplot3.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-17 15:09:54 UTC (rev 11040) +++ trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-17 20:17:18 UTC (rev 11041) @@ -2,3 +2,7 @@ figure (1) step2 (T_AH, T_opt) + + +figure (2) +step2 (T_AH, c2d (T_opt, 1)) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-18 05:58:41
|
Revision: 11045 http://octave.svn.sourceforge.net/octave/?rev=11045&view=rev Author: paramaniac Date: 2012-09-18 05:58:35 +0000 (Tue, 18 Sep 2012) Log Message: ----------- control: print legend with system names Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m trunk/octave-forge/main/control/devel/step2.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-18 05:41:47 UTC (rev 11044) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-18 05:58:35 UTC (rev 11045) @@ -23,7 +23,7 @@ ## Version: 0.3 % function [y, t, x_arr] = __time_response_2__ (sys, resptype, plotflag, tfinal, dt, x0, sysname) -function [y, t, x] = __time_response_2__ (resptype, args, plotflag) +function [y, t, x] = __time_response_2__ (resptype, args, sysname, plotflag) sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models sys_cell = cellfun (@ss, args(sys_idx), "uniformoutput", false); # convert to state-space @@ -188,6 +188,9 @@ endif endfor xlabel ("Time [s]"); + if (p == 1 && m == 1) + legend (sysname) + endif hold off; endif Modified: trunk/octave-forge/main/control/devel/step2.m =================================================================== --- trunk/octave-forge/main/control/devel/step2.m 2012-09-18 05:41:47 UTC (rev 11044) +++ trunk/octave-forge/main/control/devel/step2.m 2012-09-18 05:58:35 UTC (rev 11045) @@ -62,14 +62,23 @@ if (nargin == 0) print_usage (); endif - - %tmp = cellfun (@isa, varargin, {"lti"}); - %sys_idx = find (tmp); - - %sys_names = arrayfun (@inputname, sys_idx); + if (nargout) + sysname = {}; + else + sys_idx = find (cellfun (@isa, varargin, {"lti"})); + len = length (sys_idx); + sysname = cell (len, 1); + for k = 1 : len + try + sysname{k} = inputname(sys_idx(k)); + catch + sysname{k} = ""; + end_try_catch + endfor + endif - [y, t, x] = __time_response_2__ ("step", varargin, ! nargout); + [y, t, x] = __time_response_2__ ("step", varargin, sysname, ! nargout); if (nargout) y_r = y{1}; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-22 08:15:12
|
Revision: 11069 http://octave.svn.sourceforge.net/octave/?rev=11069&view=rev Author: paramaniac Date: 2012-09-22 08:15:06 +0000 (Sat, 22 Sep 2012) Log Message: ----------- control: support custom plot styles Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m trunk/octave-forge/main/control/devel/multiplot3.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-21 18:32:16 UTC (rev 11068) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 08:15:06 UTC (rev 11069) @@ -25,17 +25,16 @@ % function [y, t, x_arr] = __time_response_2__ (sys, resptype, plotflag, tfinal, dt, x0, sysname) function [y, t, x] = __time_response_2__ (resptype, args, sysname, plotflag) - sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models + sys_idx = find (cellfun (@isa, args, {"lti"})); # look for LTI models, 'find' needed for plot styles sys_cell = cellfun (@ss, args(sys_idx), "uniformoutput", false); # convert to state-space if (! size_equal (sys_cell{:})) error ("%s: models must have equal sizes", resptype); endif - tmp = cellfun (@is_real_matrix, args); - vec_idx = find (tmp); - n_vec = length (vec_idx); - n_sys = length (sys_cell); + vec_idx = find (cellfun (@is_real_matrix, args)); # indices of vector arguments + n_vec = length (vec_idx); # number of vector arguments + n_sys = length (sys_cell); # number of LTI systems %if (n_vec >= 1) % arg = args{vec_idx(1)}; @@ -140,6 +139,7 @@ error ("time_response: invalid response type"); endswitch + style_idx = find (cellfun (@ischar, args)); outname = get (sys_cell{end}, "outname"); outname = __labels__ (outname, "y"); colororder = get (gca, "colororder"); @@ -147,7 +147,16 @@ for k = 1 : n_sys # for every system color = colororder(1+rem (k-1, rc), :); - style = {"-", "color", color}; + if (k == n_sys) + lim = numel (args); + else + lim = sys_idx(k+1); + endif + style = args(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); + if (isempty (style)) + style = {"-", "color", color}; + endif + style discrete = ! ct_idx(k); if (discrete) # discrete-time system for i = 1 : p # for every output Modified: trunk/octave-forge/main/control/devel/multiplot3.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-21 18:32:16 UTC (rev 11068) +++ trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-22 08:15:06 UTC (rev 11069) @@ -3,6 +3,8 @@ figure (1) step2 (T_AH, T_opt) +figure (2) +step2 (T_AH, '-.r', T_opt, '-b') -figure (2) +figure (3) step2 (T_AH, c2d (T_opt, 1)) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-22 10:25:14
|
Revision: 11075 http://octave.svn.sourceforge.net/octave/?rev=11075&view=rev Author: paramaniac Date: 2012-09-22 10:25:04 +0000 (Sat, 22 Sep 2012) Log Message: ----------- control: remove cruft Modified Paths: -------------- trunk/octave-forge/main/control/devel/multiplot3.m Removed Paths: ------------- trunk/octave-forge/main/control/devel/__time_response_2__.m trunk/octave-forge/main/control/devel/step2.m Deleted: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 10:23:23 UTC (rev 11074) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 10:25:04 UTC (rev 11075) @@ -1,414 +0,0 @@ -## Copyright (C) 2009, 2010, 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 <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## Common code for the time response functions step, impulse and initial. - -## Author: Lukas Reichlin <luk...@gm...> -## Created: October 2009 -## Version: 0.3 - -% function [y, t, x_arr] = __time_response_2__ (sys, response, plotflag, tfinal, dt, x0, sysname) -function [y, t, x] = __time_response_2__ (response, args, sysname, plotflag) - - sys_idx = find (cellfun (@isa, args, {"lti"})); # look for LTI models, 'find' needed for plot styles - sys_cell = cellfun (@ss, args(sys_idx), "uniformoutput", false); # convert to state-space - - if (! size_equal (sys_cell{:})) - error ("%s: models must have equal sizes", response); - endif - - vec_idx = find (cellfun (@is_real_matrix, args)); # indices of vector arguments - n_vec = length (vec_idx); # number of vector arguments - n_sys = length (sys_cell); # number of LTI systems - - tfinal = []; - dt = []; - x0 = []; - - ## extract tfinal/t, dt, x0 from args - if (strcmpi (response, "initial")) - if (n_vec < 1) - error ("initial: require initial state vector 'x0'"); - else # initial state vector x0 specified - arg = args{vec_idx(1)}; - if (is_real_vector (arg)) - x0 = arg; - else - error ("initial: initial state vector 'x0' must be a vector of real values"); - endif - if (n_vec > 1) # tfinal or time vector t specified - arg = args{vec_idx(2)}; - if (issample (arg)) - tfinal = arg; - elseif (isempty (arg)) - ## tfinal = []; # nothing to do here - elseif (is_real_vector (arg)) - dt = abs (arg(2) - arg(1)); # assume that t is regularly spaced - tfinal = arg(end); - else - warning ("initial: argument number %d ignored", vec_idx(2)); - endif - if (n_vec > 2) # sampling time dt specified - arg = args{vec_idx(3)}; - if (issample (arg)) - dt = arg; - else - warning ("initial: argument number %d ignored", vec_idx(3)); - endif - if (n_vec > 3) - warning ("initial: ignored"); - endif - endif - endif - endif - else # step or impulse response - if (n_vec > 0) # tfinal or time vector t specified - arg = args{vec_idx(1)}; - if (issample (arg)) - tfinal = arg; - elseif (isempty (arg)) - ## tfinal = []; # nothing to do here - elseif (is_real_vector (arg)) - dt = abs (arg(2) - arg(1)); # assume that t is regularly spaced - tfinal = arg(end); - else - warning ("%s: argument number %d ignored", response, vec_idx(1)); - endif - if (n_vec > 1) # sampling time dt specified - arg = args{vec_idx(2)}; - if (issample (arg)) - dt = arg; - else - warning ("%s: argument number %d ignored", response, vec_idx(2)); - endif - if (n_vec > 2) - warning ("%s: ignored", response); - endif - endif - endif - endif - ## TODO: share common code between initial and step/impulse - - [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, {dt}, "uniformoutput", false); - tfinal = max ([tfinal{:}]); - - ct_idx = cellfun (@isct, sys_cell); - sys_dt_cell = sys_cell; - tmp = cellfun (@c2d, sys_cell(ct_idx), dt(ct_idx), {"zoh"}, "uniformoutput", false); - sys_dt_cell(ct_idx) = tmp; - - ## time vector - t = @cellfun (@(dt) reshape (0 : dt : tfinal, [], 1), dt, "uniformoutput", false); - - ## function [y, x_arr] = __initial_response__ (sys, sys_dt, t, x0) - ## function [y, x_arr] = __step_response__ (sys_dt, t) - ## function [y, x_arr] = __impulse_response__ (sys, sys_dt, t) - - switch (response) - case "initial" - [y, x] = cellfun (@__initial_response__, sys_dt_cell, t, {x0}, "uniformoutput", false); - case "step" - [y, x] = cellfun (@__step_response__, sys_dt_cell, t, "uniformoutput", false); - case "impulse" - [y, x] = cellfun (@__impulse_response__, sys_cell, sys_dt_cell, t, "uniformoutput", false); - otherwise - error ("time_response: invalid response type"); - endswitch - - - if (plotflag) # display plot - [p, m] = size (sys_cell{1}); - switch (response) - case "initial" - str = "Response to Initial Conditions"; - cols = 1; - ## yfinal = zeros (p, 1); - case "step" - str = "Step Response"; - cols = m; - ## yfinal = dcgain (sys_cell{1}); - case "impulse" - str = "Impulse Response"; - cols = m; - ## yfinal = zeros (p, m); - otherwise - error ("time_response: invalid response type"); - endswitch - - style_idx = find (cellfun (@ischar, args)); - outname = get (sys_cell{end}, "outname"); - outname = __labels__ (outname, "y"); - colororder = get (gca, "colororder"); - rc = rows (colororder); - - for k = 1 : n_sys # for every system - if (k == n_sys) - lim = numel (args); - else - lim = sys_idx(k+1); - endif - style = args(style_idx(style_idx > sys_idx(k) & style_idx <= lim)); - if (isempty (style)) - color = colororder(1+rem (k-1, rc), :); - style = {"color", color}; - endif - discrete = ! ct_idx(k); - if (discrete) # discrete-time system - for i = 1 : p # for every output - for j = 1 : cols # for every input (except for initial where cols=1) - subplot (p, cols, (i-1)*cols+j); - stairs (t{k}, y{k}(:, i, j), style{:}); - hold on; - grid on; - if (k == n_sys) - axis tight; - ylim (__axis_margin__ (ylim)) - if (j == 1) - ylabel (outname{i}); - if (i == 1) - title (str); - endif - endif - endif - endfor - endfor - else # continuous-time system - for i = 1 : p # for every output - for j = 1 : cols # for every input (except for initial where cols=1) - subplot (p, cols, (i-1)*cols+j); - ##if (n_sys == 1 && isstable (sys_cell{1})) - ## plot (t{k}, y{k}(:, i, j), style{:}, [t{k}(1), t{k}(end)], repmat (yfinal(i,j), 1, 2)); - ## ## TODO: plot final value first such that its line doesn't overprint the response - ##else - plot (t{k}, y{k}(:, i, j), style{:}); - ##endif - hold on; - grid on; - if (k == n_sys) - axis tight - ylim (__axis_margin__ (ylim)) - if (j == 1) - ylabel (outname{i}); - if (i == 1) - title (str); - endif - endif - endif - endfor - endfor - endif - endfor - xlabel ("Time [s]"); - if (p == 1 && m == 1) - legend (sysname) - endif - hold off; - endif - -endfunction - - -function [y, x_arr] = __initial_response__ (sys_dt, t, x0) - - [F, G, C, D] = ssdata (sys_dt); # system must be proper - - n = rows (F); # number of states - m = columns (G); # number of inputs - p = rows (C); # number of outputs - l_t = length (t); - - ## preallocate memory - y = zeros (l_t, p); - x_arr = zeros (l_t, n); - - ## initial conditions - x = reshape (x0, [], 1); # make sure that x is a column vector - - if (n != length (x0) || ! is_real_vector (x0)) - error ("initial: x0 must be a real vector with %d elements", n); - endif - - ## simulation - for k = 1 : l_t - y(k, :) = C * x; - x_arr(k, :) = x; - x = F * x; - endfor - -endfunction - - -function [y, x_arr] = __step_response__ (sys_dt, t) - - [F, G, C, D] = ssdata (sys_dt); # system must be proper - - n = rows (F); # number of states - m = columns (G); # number of inputs - p = rows (C); # number of outputs - l_t = length (t); - - ## preallocate memory - y = zeros (l_t, p, m); - x_arr = zeros (l_t, n, m); - - for j = 1 : m # for every input channel - ## initial conditions - x = zeros (n, 1); - u = zeros (m, 1); - u(j) = 1; - - ## simulation - for k = 1 : l_t - y(k, :, j) = C * x + D * u; - x_arr(k, :, j) = x; - x = F * x + G * u; - endfor - endfor - -endfunction - - -function [y, x_arr] = __impulse_response__ (sys, sys_dt, t) - - [~, B] = ssdata (sys); - [F, G, C, D, dt] = ssdata (sys_dt); # system must be proper - dt = abs (dt); # use 1 second if tsam is unspecified (-1) - discrete = ! isct (sys); - - n = rows (F); # number of states - m = columns (G); # number of inputs - p = rows (C); # number of outputs - l_t = length (t); - - ## preallocate memory - y = zeros (l_t, p, m); - x_arr = zeros (l_t, n, m); - - for j = 1 : m # for every input channel - ## initial conditions - u = zeros (m, 1); - u(j) = 1; - - if (discrete) - x = zeros (n, 1); # zero by definition - y(1, :, j) = D * u / dt; - x_arr(1, :, j) = x; - x = G * u / dt; - else - x = B * u; # B, not G! - y(1, :, j) = C * x; - x_arr(1, :, j) = x; - x = F * x; - endif - - ## simulation - for k = 2 : l_t - y (k, :, j) = C * x; - x_arr(k, :, j) = x; - x = F * x; - endfor - endfor - - if (discrete) - y *= dt; - x_arr *= dt; - endif - -endfunction - - -function [tfinal, dt] = __sim_horizon__ (sys, tfinal, Ts) - - ## code based on __stepimp__.m of Kai P. Mueller and A. Scottedward Hodel - - TOL = 1.0e-10; # values below TOL are assumed to be zero - N_MIN = 50; # min number of points - N_MAX = 2000; # max number of points - N_DEF = 1000; # default number of points - T_DEF = 10; # default simulation time - - ev = pole (sys); - n = length (ev); # number of states/poles - continuous = isct (sys); - discrete = ! continuous; - - if (discrete) - dt = Ts = abs (get (sys, "tsam")); - ## perform bilinear transformation on poles in z - for k = 1 : n - pol = ev(k); - if (abs (pol + 1) < TOL) - ev(k) = 0; - else - ev(k) = 2 / Ts * (pol - 1) / (pol + 1); - endif - endfor - endif - - ## remove poles near zero from eigenvalue array ev - nk = n; - for k = 1 : n - if (abs (real (ev(k))) < TOL) - ev(k) = 0; - nk -= 1; - endif - endfor - - if (nk == 0) - if (isempty (tfinal)) - tfinal = T_DEF; - endif - - if (continuous) - dt = tfinal / N_DEF; - endif - else - ev = ev(find (ev)); - ev_max = max (abs (ev)); - - if (continuous) - dt = 0.2 * pi / ev_max; - endif - - if (isempty (tfinal)) - ev_min = min (abs (real (ev))); - tfinal = 5.0 / ev_min; - - ## round up - yy = 10^(ceil (log10 (tfinal)) - 1); - tfinal = yy * ceil (tfinal / yy); - endif - - if (continuous) - N = tfinal / dt; - - if (N < N_MIN) - dt = tfinal / N_MIN; - endif - - if (N > N_MAX) - dt = tfinal / N_MAX; - endif - endif - endif - - if (continuous && ! isempty (Ts)) # catch case cont. system with dt specified - dt = Ts; - endif - -endfunction Modified: trunk/octave-forge/main/control/devel/multiplot3.m =================================================================== --- trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-22 10:23:23 UTC (rev 11074) +++ trunk/octave-forge/main/control/devel/multiplot3.m 2012-09-22 10:25:04 UTC (rev 11075) @@ -1,10 +1,10 @@ load tfs.dat figure (1) -step2 (T_AH, T_opt) +step (T_AH, T_opt) figure (2) -step2 (T_AH, '-.r', T_opt, '-b') +step (T_AH, '-.r', T_opt, '-b') figure (3) -step2 (T_AH, c2d (T_opt, 1)) +step (T_AH, c2d (T_opt, 1)) Deleted: trunk/octave-forge/main/control/devel/step2.m =================================================================== --- trunk/octave-forge/main/control/devel/step2.m 2012-09-22 10:23:23 UTC (rev 11074) +++ trunk/octave-forge/main/control/devel/step2.m 2012-09-22 10:25:04 UTC (rev 11075) @@ -1,89 +0,0 @@ -## Copyright (C) 2009, 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 <http://www.gnu.org/licenses/>. - -## -*- texinfo -*- -## @deftypefn{Function File} {[@var{y}, @var{t}, @var{x}] =} step (@var{sys}) -## @deftypefnx{Function File} {[@var{y}, @var{t}, @var{x}] =} step (@var{sys}, @var{t}) -## @deftypefnx{Function File} {[@var{y}, @var{t}, @var{x}] =} step (@var{sys}, @var{tfinal}) -## @deftypefnx{Function File} {[@var{y}, @var{t}, @var{x}] =} step (@var{sys}, @var{tfinal}, @var{dt}) -## Step response of LTI system. -## If no output arguments are given, the response is printed on the screen. -## -## @strong{Inputs} -## @table @var -## @item sys -## LTI model. -## @item t -## Time vector. Should be evenly spaced. If not specified, it is calculated by -## the poles of the system to reflect adequately the response transients. -## @item tfinal -## Optional simulation horizon. If not specified, it is calculated by -## the poles of the system to reflect adequately the response transients. -## @item dt -## Optional sampling time. Be sure to choose it small enough to capture transient -## phenomena. If not specified, it is calculated by the poles of the system. -## @end table -## -## @strong{Outputs} -## @table @var -## @item y -## Output response array. Has as many rows as time samples (length of t) -## and as many columns as outputs. -## @item t -## Time row vector. -## @item x -## State trajectories array. Has @code{length (t)} rows and as many columns as states. -## @end table -## -## @seealso{impulse, initial, lsim} -## @end deftypefn - -## Author: Lukas Reichlin <luk...@gm...> -## Created: October 2009 -## Version: 0.2 - -% function [y_r, t_r, x_r] = step2 (sys, tfinal = [], dt = []) -function [y_r, t_r, x_r] = step2 (varargin) - - if (nargin == 0) - print_usage (); - endif - - if (nargout) - sysname = {}; - else - sys_idx = find (cellfun (@isa, varargin, {"lti"})); - len = length (sys_idx); - sysname = cell (len, 1); - for k = 1 : len - try - sysname{k} = inputname(sys_idx(k)); - catch - sysname{k} = ""; - end_try_catch - endfor - endif - - [y, t, x] = __time_response_2__ ("step", varargin, sysname, ! nargout); - - if (nargout) - y_r = y{1}; - t_r = t{1}; - x_r = x{1}; - endif - -endfunction This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |