From: <par...@us...> - 2012-09-15 15:16:36
|
Revision: 11028 http://octave.svn.sourceforge.net/octave/?rev=11028&view=rev Author: paramaniac Date: 2012-09-15 15:16:30 +0000 (Sat, 15 Sep 2012) Log Message: ----------- control: quicksave draft code Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-15 14:36:20 UTC (rev 11027) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-15 15:16:30 UTC (rev 11028) @@ -35,6 +35,7 @@ tmp = cellfun (@is_real_matrix, args); vec_idx = find (tmp); n_vec = length (vec_idx) + n_sys = length (sys_cell) if (n_vec >= 1) arg = args{vec_idx(1)}; @@ -57,8 +58,6 @@ tmp = (@c2d, sys_cell(ct_idx), dt, {"zoh"}, "uniformoutput", false) sys_dt_cell(ct_idx) = tmp; - -endfunction %{ if (! isa (sys, "ss")) sys = ss (sys); # sys must be proper @@ -99,6 +98,141 @@ l_t = length (t); %} + if (plotflag) # display plot + switch (resptype) + case "initial" + str = "Response to Initial Conditions"; + cols = 1; + case "step" + str = "Step Response"; + cols = m; + case "impulse" + str = "Impulse Response"; + cols = m; + otherwise + error ("time_response: invalid response type"); + endswitch + + 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)); + 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})); + endif + endfor + endfor + else # continuous system + for k = 1 : p + for j = 1 : cols + subplot (p, cols, (k-1)*cols+j); + plot (t, y(:, 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})); + endif + endfor + endfor + endif + endfor + xlabel ("Time [s]"); + 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, sys_dt, t, x0) [A, B, C, D] = ssdata (sys); @@ -207,82 +341,11 @@ endfunction -%{ - if (plotflag) # display plot - ## TODO: Set correct titles, especially for multi-input systems +function - 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 -%} - - % function [tfinal, dt] = __sim_horizon__ (A, discrete, tfinal, Ts) function [tfinal, dt] = __sim_horizon__ (sys, tfinal, Ts) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-16 07:35:12
|
Revision: 11034 http://octave.svn.sourceforge.net/octave/?rev=11034&view=rev Author: paramaniac Date: 2012-09-16 07:35:06 +0000 (Sun, 16 Sep 2012) Log Message: ----------- control: minor touch up in draft code Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-16 07:31:08 UTC (rev 11033) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-16 07:35:06 UTC (rev 11034) @@ -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_arr] = __time_response_2__ (resptype, args) +function [y, t, x] = __time_response_2__ (resptype, args) sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models sys_cell = cellfun (@ss, args(sys_idx)); # system must be proper @@ -96,8 +96,26 @@ ## time vector t = reshape (0 : dt : tfinal, [], 1); l_t = length (t); -%} +%} +%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_cell, sys_dt_cell, t, {x0} or x0); + case "step" + [y, x] = cellfun (@__step_response__, sys_dt_cell, t); + case "impulse" + [y, x] = cellfun (@__impulse_response__, sys_cell, sys_dt_cell, t); + otherwise + error ("time_response: invalid response type"); + endswitch + + + + if (plotflag) # display plot switch (resptype) case "initial" This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-17 05:51:51
|
Revision: 11039 http://octave.svn.sourceforge.net/octave/?rev=11039&view=rev Author: paramaniac Date: 2012-09-17 05:51:45 +0000 (Mon, 17 Sep 2012) Log Message: ----------- control: touch up time response draft code Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-16 22:39:23 UTC (rev 11038) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-17 05:51:45 UTC (rev 11039) @@ -104,11 +104,11 @@ switch (resptype) case "initial" - [y, x] = cellfun (@__initial_response__, sys_cell, sys_dt_cell, t, {x0} or x0); + [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); + [y, x] = cellfun (@__step_response__, sys_dt_cell, t, "uniformoutput", false); case "impulse" - [y, x] = cellfun (@__impulse_response__, sys_cell, sys_dt_cell, t); + [y, x] = cellfun (@__impulse_response__, sys_cell, sys_dt_cell, t, "uniformoutput", false); otherwise error ("time_response: invalid response type"); endswitch @@ -251,10 +251,9 @@ -function [y, x_arr] = __initial_response__ (sys, sys_dt, t, x0) +function [y, x_arr] = __initial_response__ (sys_dt, t, x0) - [A, B, C, D] = ssdata (sys); - [F, G] = ssdata (sys_dt); + [F, G, C, D] = ssdata (sys_dt); n = rows (F); # number of states m = columns (G); # number of inputs @@ -314,8 +313,8 @@ function [y, x_arr] = __impulse_response__ (sys, sys_dt, t) - [A, B, C, D, dt] = ssdata (sys); - [F, G] = ssdata (sys_dt); + [~, B, ~, ~, dt] = ssdata (sys); + [F, G, C, D] = ssdata (sys_dt); discrete = ! isct (sys); n = rows (F); # number of states This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-18 04:09:46
|
Revision: 11043 http://octave.svn.sourceforge.net/octave/?rev=11043&view=rev Author: paramaniac Date: 2012-09-18 04:09:40 +0000 (Tue, 18 Sep 2012) Log Message: ----------- control: style fixes Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.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 00:05:50 UTC (rev 11042) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-18 04:09:40 UTC (rev 11043) @@ -25,8 +25,8 @@ % 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) - sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models - sys_cell = cellfun (@ss, args(sys_idx), "uniformoutput", false); # convert to state-space + sys_idx = cellfun (@isa, args, {"lti"}); # look for LTI models + 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); @@ -63,18 +63,18 @@ %{ if (! isa (sys, "ss")) - sys = ss (sys); # sys must be proper + sys = ss (sys); # sys must be proper endif - if (is_real_vector (tfinal) && length (tfinal) > 1) # time vector t passed - dt = tfinal(2) - tfinal(1); # assume that t is regularly spaced + if (is_real_vector (tfinal) && length (tfinal) > 1) # time vector t passed + dt = tfinal(2) - tfinal(1); # assume that t is regularly spaced tfinal = tfinal(end); endif [A, B, C, D, tsam] = ssdata (sys); - discrete = ! isct (sys); # static gains are treated as analog systems - tsam = abs (tsam); # use 1 second if tsam is unspecified (-1) + discrete = ! isct (sys); # static gains are treated as analog systems + tsam = abs (tsam); # use 1 second if tsam is unspecified (-1) if (discrete) if (! isempty (dt)) @@ -90,11 +90,11 @@ sys = c2d (sys, dt, "zoh"); endif - [F, G] = ssdata (sys); # matrices C and D don't change + [F, G] = ssdata (sys); # matrices C and D don't change - n = rows (F); # number of states - m = columns (G); # number of inputs - p = rows (C); # number of outputs + n = rows (F); # number of states + m = columns (G); # number of inputs + p = rows (C); # number of outputs ## time vector t = reshape (0 : dt : tfinal, [], 1); @@ -122,7 +122,7 @@ endswitch - if (plotflag) # display plot + if (plotflag) # display plot [p, m] = size (sys_cell{1}); switch (resptype) case "initial" @@ -141,40 +141,40 @@ outname = get (sys_cell{end}, "outname"); outname = __labels__ (outname, "y"); - for i = 1 : n_sys - discrete = ! ct_idx(i); - 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)); + for k = 1 : n_sys # for every system + 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)); hold on; grid on; - if (i == n_sys) + if (k == n_sys) axis tight; ylim (__axis_margin__ (ylim)) if (j == 1) - ylabel (outname{k}); - if (k == 1) + ylabel (outname{i}); + if (i == 1) title (str); endif endif endif endfor endfor - 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)); + 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); + plot (t{k}, y{k}(:, i, j)); hold on; grid on; - if (i == n_sys) + if (k == n_sys) axis tight ylim (__axis_margin__ (ylim)) if (j == 1) - ylabel (outname{k}); - if (k == 1) + ylabel (outname{i}); + if (i == 1) title (str); endif endif @@ -194,11 +194,11 @@ function [y, x_arr] = __initial_response__ (sys_dt, t, x0) - [F, G, C, D] = ssdata (sys_dt); # system must be proper + [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 + n = rows (F); # number of states + m = columns (G); # number of inputs + p = rows (C); # number of outputs l_t = length (t); ## preallocate memory @@ -206,7 +206,7 @@ x_arr = zeros (l_t, n); ## initial conditions - x = reshape (x0, [], 1); # make sure that x is a column vector + 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); @@ -226,16 +226,16 @@ [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 + 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 + for j = 1 : m # for every input channel ## initial conditions x = zeros (n, 1); u = zeros (m, 1); @@ -255,30 +255,30 @@ 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 + [F, G, C, D, dt] = ssdata (sys_dt); # system must be proper discrete = ! isct (sys); - n = rows (F); # number of states - m = columns (G); # number of inputs - p = rows (C); # number of outputs + 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 + 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 + 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! + x = B * u; # B, not G! y(1, :, j) = C * x; x_arr(1, :, j) = x; x = F * x; @@ -306,14 +306,14 @@ ## 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 + 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 + n = length (ev); # number of states/poles continuous = isct (sys); discrete = ! continuous; Ts = get (sys, "tsam"); 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:41:53
|
Revision: 11044 http://octave.svn.sourceforge.net/octave/?rev=11044&view=rev Author: paramaniac Date: 2012-09-18 05:41:47 +0000 (Tue, 18 Sep 2012) Log Message: ----------- control: use colororder Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.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 04:09:40 UTC (rev 11043) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-18 05:41:47 UTC (rev 11044) @@ -140,14 +140,18 @@ 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 + color = colororder(1+rem (k-1, rc), :); + style = {"-", "color", color}; 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)); + stairs (t{k}, y{k}(:, i, j), style{:}); hold on; grid on; if (k == n_sys) @@ -166,7 +170,7 @@ 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); - plot (t{k}, y{k}(:, i, j)); + plot (t{k}, y{k}(:, i, j), style{:}); hold on; grid on; if (k == n_sys) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-18 06:25:36
|
Revision: 11046 http://octave.svn.sourceforge.net/octave/?rev=11046&view=rev Author: paramaniac Date: 2012-09-18 06:25:30 +0000 (Tue, 18 Sep 2012) Log Message: ----------- control: quicksave time response draft code Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.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:58:35 UTC (rev 11045) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-18 06:25:30 UTC (rev 11046) @@ -128,12 +128,15 @@ 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 @@ -170,7 +173,12 @@ 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); - plot (t{k}, y{k}(:, i, j), style{:}); + 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) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-09-21 05:46:27
|
Revision: 11059 http://octave.svn.sourceforge.net/octave/?rev=11059&view=rev Author: paramaniac Date: 2012-09-21 05:46:21 +0000 (Fri, 21 Sep 2012) Log Message: ----------- control: save minor changes Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-20 14:11:11 UTC (rev 11058) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-21 05:46:21 UTC (rev 11059) @@ -105,14 +105,13 @@ ## 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) -%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}, "uniformoutput", false); case "step" [y, x] = cellfun (@__step_response__, sys_dt_cell, t, "uniformoutput", false); case "impulse" @@ -128,15 +127,15 @@ case "initial" str = "Response to Initial Conditions"; cols = 1; - yfinal = zeros (p, 1); + ## yfinal = zeros (p, 1); case "step" str = "Step Response"; cols = m; - yfinal = dcgain (sys_cell{1}); + ## yfinal = dcgain (sys_cell{1}); case "impulse" str = "Impulse Response"; cols = m; - yfinal = zeros (p, m); + ## yfinal = zeros (p, m); otherwise error ("time_response: invalid response type"); endswitch @@ -173,12 +172,12 @@ 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 + ##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) 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:10:59
|
Revision: 11072 http://octave.svn.sourceforge.net/octave/?rev=11072&view=rev Author: paramaniac Date: 2012-09-22 10:10:53 +0000 (Sat, 22 Sep 2012) Log Message: ----------- control: handle vector arguments in multiplot time responses Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 10:03:17 UTC (rev 11071) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 10:10:53 UTC (rev 11072) @@ -22,93 +22,104 @@ ## Created: October 2009 ## 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, sysname, plotflag) +% 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", resptype); + 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 - - %if (n_vec >= 1) - % arg = args{vec_idx(1)}; - % - %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 - ## extract tfinal/t, dt, x0 - tfinal = []; - dt = {[]}; + dt = []; + x0 = []; - %[tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, dt, "uniformoutput", false); + ## 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}, "uniformoutput", false); - + [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, {dt}, "uniformoutput", false); tfinal = max ([tfinal{:}]); - -dt - - +dt ct_idx = cellfun (@isct, sys_cell); sys_dt_cell = sys_cell; - tmp = cellfun (@c2d, sys_cell(ct_idx), dt, {"zoh"}, "uniformoutput", false); + tmp = cellfun (@c2d, sys_cell(ct_idx), dt(ct_idx), {"zoh"}, "uniformoutput", false); sys_dt_cell(ct_idx) = tmp; -%{ - if (! isa (sys, "ss")) - sys = ss (sys); # sys must be proper - endif - - if (is_real_vector (tfinal) && length (tfinal) > 1) # time vector t passed - dt = tfinal(2) - tfinal(1); # assume that t is regularly spaced - tfinal = tfinal(end); - endif - - [A, B, C, D, tsam] = ssdata (sys); - - discrete = ! isct (sys); # static gains are treated as analog systems - tsam = abs (tsam); # use 1 second if tsam is unspecified (-1) - - if (discrete) - if (! isempty (dt)) - warning ("time_response: argument dt has no effect on sampling time of discrete system"); - endif - - dt = tsam; - endif - - [tfinal, dt] = __sim_horizon__ (A, discrete, tfinal, dt); - - if (! discrete) - sys = c2d (sys, dt, "zoh"); - endif - - [F, G] = ssdata (sys); # matrices C and D don't change - - n = rows (F); # number of states - m = columns (G); # number of inputs - p = rows (C); # number of outputs - ## time vector - t = reshape (0 : dt : tfinal, [], 1); - 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) + switch (response) case "initial" [y, x] = cellfun (@__initial_response__, sys_dt_cell, t, {x0}, "uniformoutput", false); case "step" @@ -122,7 +133,7 @@ if (plotflag) # display plot [p, m] = size (sys_cell{1}); - switch (resptype) + switch (response) case "initial" str = "Response to Initial Conditions"; cols = 1; @@ -146,7 +157,6 @@ rc = rows (colororder); for k = 1 : n_sys # for every system - color = colororder(1+rem (k-1, rc), :); if (k == n_sys) lim = numel (args); else @@ -154,9 +164,9 @@ 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 - style discrete = ! ct_idx(k); if (discrete) # discrete-time system for i = 1 : p # for every output @@ -213,8 +223,6 @@ endfunction - - function [y, x_arr] = __initial_response__ (sys_dt, t, x0) [F, G, C, D] = ssdata (sys_dt); # system must be proper @@ -279,6 +287,7 @@ [~, 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 @@ -323,8 +332,6 @@ endfunction - -% function [tfinal, dt] = __sim_horizon__ (A, discrete, tfinal, Ts) function [tfinal, dt] = __sim_horizon__ (sys, tfinal, Ts) ## code based on __stepimp__.m of Kai P. Mueller and A. Scottedward Hodel @@ -339,9 +346,9 @@ n = length (ev); # number of states/poles continuous = isct (sys); discrete = ! continuous; - Ts = get (sys, "tsam"); if (discrete) + dt = Ts = abs (get (sys, "tsam")); ## perform bilinear transformation on poles in z for k = 1 : n pol = ev(k); @@ -400,8 +407,8 @@ endif endif - %if (! isempty (Ts)) # catch case cont. system with dt specified - % dt = Ts; - %endif + if (continuous && ! isempty (Ts)) # catch case cont. system with dt specified + dt = Ts; + 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-22 10:14:18
|
Revision: 11073 http://octave.svn.sourceforge.net/octave/?rev=11073&view=rev Author: paramaniac Date: 2012-09-22 10:14:12 +0000 (Sat, 22 Sep 2012) Log Message: ----------- control: minor bug fix (stairs doesn't like the "-" style argument) Modified Paths: -------------- trunk/octave-forge/main/control/devel/__time_response_2__.m Modified: trunk/octave-forge/main/control/devel/__time_response_2__.m =================================================================== --- trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 10:10:53 UTC (rev 11072) +++ trunk/octave-forge/main/control/devel/__time_response_2__.m 2012-09-22 10:14:12 UTC (rev 11073) @@ -106,7 +106,7 @@ [tfinal, dt] = cellfun (@__sim_horizon__, sys_cell, {tfinal}, {dt}, "uniformoutput", false); tfinal = max ([tfinal{:}]); -dt + ct_idx = cellfun (@isct, sys_cell); sys_dt_cell = sys_cell; tmp = cellfun (@c2d, sys_cell(ct_idx), dt(ct_idx), {"zoh"}, "uniformoutput", false); @@ -165,7 +165,7 @@ 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}; + style = {"color", color}; endif discrete = ! ct_idx(k); if (discrete) # discrete-time system This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |