Skip to content

Commit 3adc610

Browse files
bnavigatorbnavigator
authored
Merge pull request #531 from sawyerbfuller/sisotooling
sisotool small visual cleanup, new feature to show step response of different input-output than loop
2 parents d75c395 + abca69d commit 3adc610

File tree

4 files changed

+153
-48
lines changed

4 files changed

+153
-48
lines changed

‎control/bdalg.py

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -334,7 +334,8 @@ def connect(sys, Q, inputv, outputv):
334334
interconnecting multiple systems.
335335
336336
"""
337-
inputv, outputv, Q = np.asarray(inputv), np.asarray(outputv), np.asarray(Q)
337+
inputv, outputv, Q = \
338+
np.atleast_1d(inputv), np.atleast_1d(outputv), np.atleast_1d(Q)
338339
# check indices
339340
index_errors = (inputv - 1 > sys.ninputs) | (inputv < 1)
340341
if np.any(index_errors):

‎control/rlocus.py

Lines changed: 16 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -137,8 +137,10 @@ def root_locus(sys, kvect=None, xlim=None, ylim=None,
137137
print_gain = config._get_param(
138138
'rlocus', 'print_gain', print_gain, _rlocus_defaults)
139139

140+
sys_loop = sys if sys.issiso() else sys[0,0]
141+
140142
# Convert numerator and denominator to polynomials if they aren't
141-
(nump, denp) = _systopoly1d(sys)
143+
(nump, denp) = _systopoly1d(sys_loop)
142144

143145
# if discrete-time system and if xlim and ylim are not given,
144146
# that we a view of the unit circle
@@ -241,12 +243,12 @@ def root_locus(sys, kvect=None, xlim=None, ylim=None,
241243
else:
242244
_sgrid_func()
243245
else:
244-
ax.axhline(0., linestyle=':', color='k', zorder=-20)
245-
ax.axvline(0., linestyle=':', color='k', zorder=-20)
246+
ax.axhline(0., linestyle=':', color='k', linewidth=.75, zorder=-20)
247+
ax.axvline(0., linestyle=':', color='k', linewidth=.75, zorder=-20)
246248
if isdtime(sys, strict=True):
247249
ax.add_patch(plt.Circle(
248250
(0, 0), radius=1.0, linestyle=':', edgecolor='k',
249-
linewidth=1.5, fill=False, zorder=-20))
251+
linewidth=0.75, fill=False, zorder=-20))
250252

251253
return mymat, kvect
252254

@@ -540,8 +542,9 @@ def _RLSortRoots(mymat):
540542

541543
def _RLZoomDispatcher(event, sys, ax_rlocus, plotstr):
542544
"""Rootlocus plot zoom dispatcher"""
545+
sys_loop = sys if sys.issiso() else sys[0,0]
543546

544-
nump, denp = _systopoly1d(sys)
547+
nump, denp = _systopoly1d(sys_loop)
545548
xlim, ylim = ax_rlocus.get_xlim(), ax_rlocus.get_ylim()
546549

547550
kvect, mymat, xlim, ylim = _default_gains(
@@ -573,21 +576,23 @@ def _RLClickDispatcher(event, sys, fig, ax_rlocus, plotstr, sisotool=False,
573576

574577
def _RLFeedbackClicksPoint(event, sys, fig, ax_rlocus, sisotool=False):
575578
"""Display root-locus gain feedback point for clicks on root-locus plot"""
576-
(nump, denp) = _systopoly1d(sys)
579+
sys_loop = sys if sys.issiso() else sys[0,0]
580+
581+
(nump, denp) = _systopoly1d(sys_loop)
577582

578583
xlim = ax_rlocus.get_xlim()
579584
ylim = ax_rlocus.get_ylim()
580-
x_tolerance = 0.05 * abs((xlim[1] - xlim[0]))
581-
y_tolerance = 0.05 * abs((ylim[1] - ylim[0]))
585+
x_tolerance = 0.1 * abs((xlim[1] - xlim[0]))
586+
y_tolerance = 0.1 * abs((ylim[1] - ylim[0]))
582587
gain_tolerance = np.mean([x_tolerance, y_tolerance])*0.1
583588

584589
# Catch type error when event click is in the figure but not in an axis
585590
try:
586591
s = complex(event.xdata, event.ydata)
587-
K = -1. / sys(s)
588-
K_xlim = -1. / sys(
592+
K = -1. / sys_loop(s)
593+
K_xlim = -1. / sys_loop(
589594
complex(event.xdata + 0.05 * abs(xlim[1] - xlim[0]), event.ydata))
590-
K_ylim = -1. / sys(
595+
K_ylim = -1. / sys_loop(
591596
complex(event.xdata, event.ydata + 0.05 * abs(ylim[1] - ylim[0])))
592597

593598
except TypeError:

‎control/sisotool.py

Lines changed: 59 additions & 33 deletions
Original file line numberDiff line numberDiff line change
@@ -1,17 +1,19 @@
11
__all__ = ['sisotool']
22

3+
from control.exception import ControlMIMONotImplemented
34
from .freqplot import bode_plot
45
from .timeresp import step_response
56
from .lti import issiso, isdtime
7+
from .xferfcn import TransferFunction
8+
from .bdalg import append, connect
69
import matplotlib
710
import matplotlib.pyplot as plt
811
import warnings
912

10-
def sisotool(sys, kvect = None, xlim_rlocus = None, ylim_rlocus = None,
11-
plotstr_rlocus = 'b' if int(matplotlib.__version__[0]) == 1 else 'C0',
12-
rlocus_grid = False, omega = None, dB = None, Hz = None,
13-
deg = None, omega_limits = None, omega_num = None,
14-
margins_bode = True, tvect=None):
13+
def sisotool(sys, kvect=None, xlim_rlocus=None, ylim_rlocus=None,
14+
plotstr_rlocus='C0', rlocus_grid=False, omega=None, dB=None,
15+
Hz=None, deg=None, omega_limits=None, omega_num=None,
16+
margins_bode=True, tvect=None):
1517
"""
1618
Sisotool style collection of plots inspired by MATLAB's sisotool.
1719
The left two plots contain the bode magnitude and phase diagrams.
@@ -22,7 +24,16 @@ def sisotool(sys, kvect = None, xlim_rlocus = None, ylim_rlocus = None,
2224
Parameters
2325
----------
2426
sys : LTI object
25-
Linear input/output systems (SISO only)
27+
Linear input/output systems. If sys is SISO, use the same
28+
system for the root locus and step response. If it is desired to
29+
see a different step response than feedback(K*loop,1), sys can be
30+
provided as a two-input, two-output system (e.g. by using
31+
:func:`bdgalg.connect' or :func:`iosys.interconnect`). Sisotool
32+
inserts the negative of the selected gain K between the first output
33+
and first input and uses the second input and output for computing
34+
the step response. This allows you to see the step responses of more
35+
complex systems, for example, systems with a feedforward path into the
36+
plant or in which the gain appears in the feedback path.
2637
kvect : list or ndarray, optional
2738
List of gains to use for plotting root locus
2839
xlim_rlocus : tuple or list, optional
@@ -32,21 +43,23 @@ def sisotool(sys, kvect = None, xlim_rlocus = None, ylim_rlocus = None,
3243
control of y-axis range
3344
plotstr_rlocus : :func:`matplotlib.pyplot.plot` format string, optional
3445
plotting style for the root locus plot(color, linestyle, etc)
35-
rlocus_grid: boolean (default = False)
36-
If True plot s-plane grid.
37-
omega : freq_range
38-
Range of frequencies in rad/sec for the bode plot
46+
rlocus_grid : boolean (default = False)
47+
If True plot s- or z-plane grid.
48+
omega : array_like
49+
List of frequencies in rad/sec to be used for bode plot
3950
dB : boolean
4051
If True, plot result in dB for the bode plot
4152
Hz : boolean
4253
If True, plot frequency in Hz for the bode plot (omega must be provided in rad/sec)
4354
deg : boolean
4455
If True, plot phase in degrees for the bode plot (else radians)
45-
omega_limits: tuple, list, ... of two values
56+
omega_limits : array_like of two values
4657
Limits of the to generate frequency vector.
47-
If Hz=True the limits are in Hz otherwise in rad/s.
48-
omega_num: int
49-
number of samples
58+
If Hz=True the limits are in Hz otherwise in rad/s. Ignored if omega
59+
is provided, and auto-generated if omitted.
60+
omega_num : int
61+
Number of samples to plot. Defaults to
62+
config.defaults['freqplot.number_of_samples'].
5063
margins_bode : boolean
5164
If True, plot gain and phase margin in the bode plot
5265
tvect : list or ndarray, optional
@@ -60,8 +73,11 @@ def sisotool(sys, kvect = None, xlim_rlocus = None, ylim_rlocus = None,
6073
"""
6174
from .rlocus import root_locus
6275

63-
# Check if it is a single SISO system
64-
issiso(sys,strict=True)
76+
# sys as loop transfer function if SISO
77+
if not sys.issiso():
78+
if not (sys.ninputs == 2 and sys.noutputs == 2):
79+
raise ControlMIMONotImplemented(
80+
'sys must be SISO or 2-input, 2-output')
6581

6682
# Setup sisotool figure or superimpose if one is already present
6783
fig = plt.gcf()
@@ -84,64 +100,74 @@ def sisotool(sys, kvect = None, xlim_rlocus = None, ylim_rlocus = None,
84100
}
85101

86102
# First time call to setup the bode and step response plots
87-
_SisotoolUpdate(sys, fig,1 if kvect is None else kvect[0],bode_plot_params)
103+
_SisotoolUpdate(sys, fig,
104+
1 if kvect is None else kvect[0], bode_plot_params)
88105

89106
# Setup the root-locus plot window
90-
root_locus(sys,kvect=kvect,xlim=xlim_rlocus,ylim = ylim_rlocus,plotstr=plotstr_rlocus,grid = rlocus_grid,fig=fig,bode_plot_params=bode_plot_params,tvect=tvect,sisotool=True)
107+
root_locus(sys, kvect=kvect, xlim=xlim_rlocus,
108+
ylim=ylim_rlocus, plotstr=plotstr_rlocus, grid=rlocus_grid,
109+
fig=fig, bode_plot_params=bode_plot_params, tvect=tvect, sisotool=True)
91110

92-
def _SisotoolUpdate(sys,fig,K,bode_plot_params,tvect=None):
111+
def _SisotoolUpdate(sys, fig, K, bode_plot_params, tvect=None):
93112

94-
if int(matplotlib.__version__[0]) == 1:
95-
title_font_size = 12
96-
label_font_size = 10
97-
else:
98-
title_font_size = 10
99-
label_font_size = 8
113+
title_font_size = 10
114+
label_font_size = 8
100115

101116
# Get the subaxes and clear them
102-
ax_mag,ax_rlocus,ax_phase,ax_step = fig.axes[0],fig.axes[1],fig.axes[2],fig.axes[3]
117+
ax_mag, ax_rlocus, ax_phase, ax_step = \
118+
fig.axes[0], fig.axes[1], fig.axes[2], fig.axes[3]
103119

104120
# Catch matplotlib 2.1.x and higher userwarnings when clearing a log axis
105121
with warnings.catch_warnings():
106122
warnings.simplefilter("ignore")
107123
ax_step.clear(), ax_mag.clear(), ax_phase.clear()
108124

125+
sys_loop = sys if sys.issiso() else sys[0,0]
126+
109127
# Update the bodeplot
110-
bode_plot_params['syslist'] = sys*K.real
128+
bode_plot_params['syslist'] = sys_loop*K.real
111129
bode_plot(**bode_plot_params)
112130

113131
# Set the titles and labels
114132
ax_mag.set_title('Bode magnitude',fontsize = title_font_size)
115133
ax_mag.set_ylabel(ax_mag.get_ylabel(), fontsize=label_font_size)
134+
ax_mag.tick_params(axis='both', which='major', labelsize=label_font_size)
116135

117136
ax_phase.set_title('Bode phase',fontsize=title_font_size)
118137
ax_phase.set_xlabel(ax_phase.get_xlabel(),fontsize=label_font_size)
119138
ax_phase.set_ylabel(ax_phase.get_ylabel(),fontsize=label_font_size)
120139
ax_phase.get_xaxis().set_label_coords(0.5, -0.15)
121140
ax_phase.get_shared_x_axes().join(ax_phase, ax_mag)
141+
ax_phase.tick_params(axis='both', which='major', labelsize=label_font_size)
122142

123143
ax_step.set_title('Step response',fontsize = title_font_size)
124144
ax_step.set_xlabel('Time (seconds)',fontsize=label_font_size)
125-
ax_step.set_ylabel('Amplitude',fontsize=label_font_size)
145+
ax_step.set_ylabel('Output',fontsize=label_font_size)
126146
ax_step.get_xaxis().set_label_coords(0.5, -0.15)
127147
ax_step.get_yaxis().set_label_coords(-0.15, 0.5)
148+
ax_step.tick_params(axis='both', which='major', labelsize=label_font_size)
128149

129150
ax_rlocus.set_title('Root locus',fontsize = title_font_size)
130151
ax_rlocus.set_ylabel('Imag', fontsize=label_font_size)
131152
ax_rlocus.set_xlabel('Real', fontsize=label_font_size)
132153
ax_rlocus.get_xaxis().set_label_coords(0.5, -0.15)
133154
ax_rlocus.get_yaxis().set_label_coords(-0.15, 0.5)
134-
135-
155+
ax_rlocus.tick_params(axis='both', which='major',labelsize=label_font_size)
136156

137157
# Generate the step response and plot it
138-
sys_closed = (K*sys).feedback(1)
158+
if sys.issiso():
159+
sys_closed = (K*sys).feedback(1)
160+
else:
161+
sys_closed = append(sys, -K)
162+
connects = [[1, 3],
163+
[3, 1]]
164+
sys_closed = connect(sys_closed, connects, 2, 2)
139165
if tvect is None:
140166
tvect, yout = step_response(sys_closed, T_num=100)
141167
else:
142-
tvect, yout = step_response(sys_closed,tvect)
168+
tvect, yout = step_response(sys_closed, tvect)
143169
if isdtime(sys_closed, strict=True):
144-
ax_step.plot(tvect, yout, 'o')
170+
ax_step.plot(tvect, yout, '.')
145171
else:
146172
ax_step.plot(tvect, yout)
147173
ax_step.axhline(1.,linestyle=':',color='k',zorder=-20)

‎control/tests/sisotool_test.py

Lines changed: 76 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,6 @@
11
"""sisotool_test.py"""
22

3+
from control.exception import ControlMIMONotImplemented
34
import matplotlib.pyplot as plt
45
import numpy as np
56
from numpy.testing import assert_array_almost_equal
@@ -8,6 +9,7 @@
89
from control.sisotool import sisotool
910
from control.rlocus import _RLClickDispatcher
1011
from control.xferfcn import TransferFunction
12+
from control.statesp import StateSpace
1113

1214

1315
@pytest.mark.usefixtures("mplcleanup")
@@ -19,6 +21,35 @@ def sys(self):
1921
"""Return a generic SISO transfer function"""
2022
return TransferFunction([1000], [1, 25, 100, 0])
2123

24+
@pytest.fixture
25+
def sysdt(self):
26+
"""Return a generic SISO transfer function"""
27+
return TransferFunction([1000], [1, 25, 100, 0], True)
28+
29+
@pytest.fixture
30+
def sys222(self):
31+
"""2-states square system (2 inputs x 2 outputs)"""
32+
A222 = [[4., 1.],
33+
[2., -3]]
34+
B222 = [[5., 2.],
35+
[-3., -3.]]
36+
C222 = [[2., -4],
37+
[0., 1.]]
38+
D222 = [[3., 2.],
39+
[1., -1.]]
40+
return StateSpace(A222, B222, C222, D222)
41+
42+
@pytest.fixture
43+
def sys221(self):
44+
"""2-states, 2 inputs x 1 output"""
45+
A222 = [[4., 1.],
46+
[2., -3]]
47+
B222 = [[5., 2.],
48+
[-3., -3.]]
49+
C221 = [[0., 1.]]
50+
D221 = [[1., -1.]]
51+
return StateSpace(A222, B222, C221, D221)
52+
2253
def test_sisotool(self, sys):
2354
sisotool(sys, Hz=False)
2455
fig = plt.gcf()
@@ -27,7 +58,7 @@ def test_sisotool(self, sys):
2758
# Check the initial root locus plot points
2859
initial_point_0 = (np.array([-22.53155977]), np.array([0.]))
2960
initial_point_1 = (np.array([-1.23422011]), np.array([-6.54667031]))
30-
initial_point_2 = (np.array([-1.23422011]), np.array([06.54667031]))
61+
initial_point_2 = (np.array([-1.23422011]), np.array([6.54667031]))
3162
assert_array_almost_equal(ax_rlocus.lines[0].get_data(),
3263
initial_point_0, 4)
3364
assert_array_almost_equal(ax_rlocus.lines[1].get_data(),
@@ -88,7 +119,49 @@ def test_sisotool(self, sys):
88119
step_response_moved = np.array(
89120
[0., 0.0072, 0.0516, 0.1554, 0.3281, 0.5681, 0.8646, 1.1987,
90121
1.5452, 1.875])
91-
# old: array([0., 0.0239, 0.161 , 0.4547, 0.8903, 1.407,
92-
# 1.9121, 2.2989, 2.4686, 2.353])
93122
assert_array_almost_equal(
94123
ax_step.lines[0].get_data()[1][:10], step_response_moved, 4)
124+
125+
def test_sisotool_tvect(self, sys):
126+
# test supply tvect
127+
tvect = np.linspace(0, 1, 10)
128+
sisotool(sys, tvect=tvect)
129+
fig = plt.gcf()
130+
ax_rlocus, ax_step = fig.axes[1], fig.axes[3]
131+
132+
# Move the rootlocus to another point and confirm same tvect
133+
event = type('test', (object,), {'xdata': 2.31206868287,
134+
'ydata': 15.5983051046,
135+
'inaxes': ax_rlocus.axes})()
136+
_RLClickDispatcher(event=event, sys=sys, fig=fig,
137+
ax_rlocus=ax_rlocus, sisotool=True, plotstr='-',
138+
bode_plot_params=dict(), tvect=tvect)
139+
assert_array_almost_equal(tvect, ax_step.lines[0].get_data()[0])
140+
141+
def test_sisotool_tvect_dt(self, sysdt):
142+
# test supply tvect
143+
tvect = np.linspace(0, 1, 10)
144+
sisotool(sysdt, tvect=tvect)
145+
fig = plt.gcf()
146+
ax_rlocus, ax_step = fig.axes[1], fig.axes[3]
147+
148+
# Move the rootlocus to another point and confirm same tvect
149+
event = type('test', (object,), {'xdata': 2.31206868287,
150+
'ydata': 15.5983051046,
151+
'inaxes': ax_rlocus.axes})()
152+
_RLClickDispatcher(event=event, sys=sysdt, fig=fig,
153+
ax_rlocus=ax_rlocus, sisotool=True, plotstr='-',
154+
bode_plot_params=dict(), tvect=tvect)
155+
assert_array_almost_equal(tvect, ax_step.lines[0].get_data()[0])
156+
157+
def test_sisotool_mimo(self, sys222, sys221):
158+
# a 2x2 should not raise an error:
159+
sisotool(sys222)
160+
161+
# but 2 input, 1 output should
162+
with pytest.raises(ControlMIMONotImplemented):
163+
sisotool(sys221)
164+
165+
166+
167+

0 commit comments

Comments
 (0)