Skip to content

changes to rlocus to be compatible with discrete-time systems #410

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Jul 14, 2020
Merged

changes to rlocus to be compatible with discrete-time systems #410

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
e2d8f9d
added link to lqe (linear quadratic estimator) function in docs
sawyerbfuller May 19, 2020
0854ee4
Merge branch 'master' of https://github.com/python-control/python-con…
murrayrm Jul 11, 2020
07e9c6e
changes to rlocus and sisotool to make them compatible with discrete-…
sawyerbfuller May 24, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
11 changes: 6 additions & 5 deletions control/grid.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -136,11 +136,12 @@ def nogrid():
return ax, f


def zgrid(zetas=None, wns=None):
def zgrid(zetas=None, wns=None, ax=None):
'''Draws discrete damping and frequency grid'''

fig = plt.gcf()
ax = fig.gca()
if ax is None:
ax = fig.gca()

# Constant damping lines
if zetas is None:
Expand All @@ -154,11 +155,11 @@ def zgrid(zetas=None, wns=None):
# Draw upper part in retangular coordinates
xret = mag*cos(ang)
yret = mag*sin(ang)
ax.plot(xret, yret, 'k:', lw=1)
ax.plot(xret, yret, ':', color='grey', lw=0.75)
# Draw lower part in retangular coordinates
xret = mag*cos(-ang)
yret = mag*sin(-ang)
ax.plot(xret, yret, 'k:', lw=1)
ax.plot(xret, yret, ':', color='grey', lw=0.75)
# Annotation
an_i = int(len(xret)/2.5)
an_x = xret[an_i]
Expand All @@ -177,7 +178,7 @@ def zgrid(zetas=None, wns=None):
# Draw in retangular coordinates
xret = mag*cos(ang)
yret = mag*sin(ang)
ax.plot(xret, yret, 'k:', lw=1)
ax.plot(xret, yret, ':', color='grey', lw=0.75)
# Annotation
an_i = -1
an_x = xret[an_i]
Expand Down
76 changes: 54 additions & 22 deletions control/rlocus.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,9 @@
# RMM, 2 April 2011: modified to work with new LTI structure (see ChangeLog)
# * Not tested: should still work on signal.ltisys objects
#
# Sawyer B. Fuller (minster@uw.edu) 21 May 2020:
# * added compatibility with discrete-time systems.
#
# $Id$

# Packages used by this module
Expand All @@ -52,9 +55,11 @@
import matplotlib.pyplot as plt
from numpy import array, poly1d, row_stack, zeros_like, real, imag
import scipy.signal # signal processing toolbox
from .lti import isdtime
from .xferfcn import _convert_to_transfer_function
from .exception import ControlMIMONotImplemented
from .sisotool import _SisotoolUpdate
from .grid import sgrid, zgrid
from . import config

__all__ = ['root_locus', 'rlocus']
Expand Down Expand Up @@ -131,6 +136,13 @@ def root_locus(sys, kvect=None, xlim=None, ylim=None,
# Convert numerator and denominator to polynomials if they aren't
(nump, denp) = _systopoly1d(sys)

# if discrete-time system and if xlim and ylim are not given,
# that we a view of the unit circle
if xlim is None and isdtime(sys, strict=True):
xlim = (-1.2, 1.2)
if ylim is None and isdtime(sys, strict=True):
xlim = (-1.3, 1.3)

if kvect is None:
start_mat = _RLFindRoots(nump, denp, [1])
kvect, mymat, xlim, ylim = _default_gains(nump, denp, xlim, ylim)
Expand Down Expand Up @@ -163,10 +175,14 @@ def root_locus(sys, kvect=None, xlim=None, ylim=None,
[root.real for root in start_mat],
[root.imag for root in start_mat],
'm.', marker='s', markersize=8, zorder=20, label='gain_point')
s = start_mat[0][0]
if isdtime(sys, strict=True):
zeta = -np.cos(np.angle(np.log(s)))
else:
zeta = -1 * s.real / abs(s)
fig.suptitle(
"Clicked at: %10.4g%+10.4gj gain: %10.4g damp: %10.4g" %
(start_mat[0][0].real, start_mat[0][0].imag,
1, -1 * start_mat[0][0].real / abs(start_mat[0][0])),
(s.real, s.imag, 1, zeta),
fontsize=12 if int(mpl.__version__[0]) == 1 else 10)
fig.canvas.mpl_connect(
'button_release_event',
Expand Down Expand Up @@ -199,20 +215,31 @@ def root_locus(sys, kvect=None, xlim=None, ylim=None,
ax.plot(real(col), imag(col), plotstr, label='rootlocus')

# Set up plot axes and labels
if xlim:
ax.set_xlim(xlim)
if ylim:
ax.set_ylim(ylim)

ax.set_xlabel('Real')
ax.set_ylabel('Imaginary')

if grid and sisotool:
_sgrid_func(f)
if isdtime(sys, strict=True):
zgrid(ax=ax)
else:
_sgrid_func(f)
elif grid:
_sgrid_func()
if isdtime(sys, strict=True):
zgrid(ax=ax)
else:
_sgrid_func()
else:
ax.axhline(0., linestyle=':', color='k', zorder=-20)
ax.axvline(0., linestyle=':', color='k')
ax.axvline(0., linestyle=':', color='k', zorder=-20)
if isdtime(sys, strict=True):
ax.add_patch(plt.Circle((0,0), radius=1.0,
linestyle=':', edgecolor='k', linewidth=1.5,
fill=False, zorder=-20))

if xlim:
ax.set_xlim(xlim)
if ylim:
ax.set_ylim(ylim)

return mymat, kvect

Expand Down Expand Up @@ -567,12 +594,17 @@ def _RLFeedbackClicksPoint(event, sys, fig, ax_rlocus, sisotool=False):
if abs(K.real) > 1e-8 and abs(K.imag / K.real) < gain_tolerance and \
event.inaxes == ax_rlocus.axes and K.real > 0.:

if isdtime(sys, strict=True):
zeta = -np.cos(np.angle(np.log(s)))
else:
zeta = -1 * s.real / abs(s)

# Display the parameters in the output window and figure
print("Clicked at %10.4g%+10.4gj gain %10.4g damp %10.4g" %
(s.real, s.imag, K.real, -1 * s.real / abs(s)))
(s.real, s.imag, K.real, zeta))
fig.suptitle(
"Clicked at: %10.4g%+10.4gj gain: %10.4g damp: %10.4g" %
(s.real, s.imag, K.real, -1 * s.real / abs(s)),
(s.real, s.imag, K.real, zeta),
fontsize=12 if int(mpl.__version__[0]) == 1 else 10)

# Remove the previous line
Expand Down Expand Up @@ -616,13 +648,13 @@ def _sgrid_func(fig=None, zeta=None, wn=None):
if zeta is None:
zeta = _default_zetas(xlim, ylim)

angules = []
angles = []
for z in zeta:
if (z >= 1e-4) and (z <= 1):
angules.append(np.pi/2 + np.arcsin(z))
angles.append(np.pi/2 + np.arcsin(z))
else:
zeta.remove(z)
y_over_x = np.tan(angules)
y_over_x = np.tan(angles)

# zeta-constant lines

Expand All @@ -647,30 +679,30 @@ def _sgrid_func(fig=None, zeta=None, wn=None):
ax.plot([0, 0], [ylim[0], ylim[1]],
color='gray', linestyle='dashed', linewidth=0.5)

angules = np.linspace(-90, 90, 20)*np.pi/180
angles = np.linspace(-90, 90, 20)*np.pi/180
if wn is None:
wn = _default_wn(xlocator(), ylim)

for om in wn:
if om < 0:
yp = np.sin(angules)*np.abs(om)
xp = -np.cos(angules)*np.abs(om)
yp = np.sin(angles)*np.abs(om)
xp = -np.cos(angles)*np.abs(om)
ax.plot(xp, yp, color='gray',
linestyle='dashed', linewidth=0.5)
an = "%.2f" % -om
ax.annotate(an, textcoords='data', xy=[om, 0], fontsize=8)


def _default_zetas(xlim, ylim):
"""Return default list of dumps coefficients"""
"""Return default list of damping coefficients"""
sep1 = -xlim[0]/4
ang1 = [np.arctan((sep1*i)/ylim[1]) for i in np.arange(1, 4, 1)]
sep2 = ylim[1] / 3
ang2 = [np.arctan(-xlim[0]/(ylim[1]-sep2*i)) for i in np.arange(1, 3, 1)]

angules = np.concatenate((ang1, ang2))
angules = np.insert(angules, len(angules), np.pi/2)
zeta = np.sin(angules)
angles = np.concatenate((ang1, ang2))
angles = np.insert(angles, len(angles), np.pi/2)
zeta = np.sin(angles)
return zeta.tolist()


Expand Down
7 changes: 5 additions & 2 deletions control/sisotool.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

from .freqplot import bode_plot
from .timeresp import step_response
from .lti import issiso
from .lti import issiso, isdtime
import matplotlib
import matplotlib.pyplot as plt
import warnings
Expand Down Expand Up @@ -139,7 +139,10 @@ def _SisotoolUpdate(sys,fig,K,bode_plot_params,tvect=None):
tvect, yout = step_response(sys_closed, T_num=100)
else:
tvect, yout = step_response(sys_closed,tvect)
ax_step.plot(tvect, yout)
if isdtime(sys_closed, strict=True):
ax_step.plot(tvect, yout, 'o')
else:
ax_step.plot(tvect, yout)
ax_step.axhline(1.,linestyle=':',color='k',zorder=-20)

# Manually adjust the spacing and draw the canvas
Expand Down