Skip to content

Commit 89dc937

Browse files
author
victor
committed
Refactored sgrid based on matplotlib projections and moved sgrid and zgrid to new python module.
1 parent 3d0bf03 commit 89dc937

File tree

3 files changed

+212
-136
lines changed

3 files changed

+212
-136
lines changed

‎control/grid.py‎

Lines changed: 182 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,182 @@
1+
import numpy as np
2+
from numpy import cos, sin, sqrt, linspace, pi, exp
3+
import matplotlib.pyplot as plt
4+
from mpl_toolkits.axisartist import SubplotHost
5+
from mpl_toolkits.axisartist.grid_helper_curvelinear import GridHelperCurveLinear
6+
import mpl_toolkits.axisartist.angle_helper as angle_helper
7+
from matplotlib.projections import PolarAxes
8+
from matplotlib.transforms import Affine2D
9+
10+
class FormatterDMS(object):
11+
'''Transforms angle ticks to damping ratios'''
12+
def __call__(self,direction,factor,values):
13+
angles_deg = values/factor
14+
damping_ratios = np.cos((180-angles_deg)*np.pi/180)
15+
ret = ["%.2f"%val for val in damping_ratios]
16+
return ret
17+
18+
class ModifiedExtremeFinderCycle(angle_helper.ExtremeFinderCycle):
19+
'''Changed to allow only left hand-side polar grid'''
20+
def __call__(self, transform_xy, x1, y1, x2, y2):
21+
x_, y_ = np.linspace(x1, x2, self.nx), np.linspace(y1, y2, self.ny)
22+
x, y = np.meshgrid(x_, y_)
23+
lon, lat = transform_xy(np.ravel(x), np.ravel(y))
24+
25+
with np.errstate(invalid='ignore'):
26+
if self.lon_cycle is not None:
27+
lon0 = np.nanmin(lon)
28+
lon -= 360. * ((lon - lon0) > 360.) # Changed from 180 to 360 to be able to span only 90-270 (left hand side)
29+
if self.lat_cycle is not None:
30+
lat0 = np.nanmin(lat)
31+
lat -= 360. * ((lat - lat0) > 360.) # Changed from 180 to 360 to be able to span only 90-270 (left hand side)
32+
33+
lon_min, lon_max = np.nanmin(lon), np.nanmax(lon)
34+
lat_min, lat_max = np.nanmin(lat), np.nanmax(lat)
35+
36+
lon_min, lon_max, lat_min, lat_max = \
37+
self._adjust_extremes(lon_min, lon_max, lat_min, lat_max)
38+
39+
return lon_min, lon_max, lat_min, lat_max
40+
41+
def sgrid():
42+
# From matplotlib demos:
43+
# https://matplotlib.org/gallery/axisartist/demo_curvelinear_grid.html
44+
# https://matplotlib.org/gallery/axisartist/demo_floating_axis.html
45+
46+
# PolarAxes.PolarTransform takes radian. However, we want our coordinate
47+
# system in degree
48+
tr = Affine2D().scale(np.pi/180., 1.) + PolarAxes.PolarTransform()
49+
# polar projection, which involves cycle, and also has limits in
50+
# its coordinates, needs a special method to find the extremes
51+
# (min, max of the coordinate within the view).
52+
53+
# 20, 20 : number of sampling points along x, y direction
54+
sampling_points = 20
55+
extreme_finder = ModifiedExtremeFinderCycle(sampling_points, sampling_points,
56+
lon_cycle=360,
57+
lat_cycle=None,
58+
lon_minmax=(90,270),
59+
lat_minmax=(0, np.inf),)
60+
61+
grid_locator1 = angle_helper.LocatorDMS(15)
62+
tick_formatter1 = FormatterDMS()
63+
grid_helper = GridHelperCurveLinear(tr,
64+
extreme_finder=extreme_finder,
65+
grid_locator1=grid_locator1,
66+
tick_formatter1=tick_formatter1
67+
)
68+
69+
fig = plt.figure()
70+
ax = SubplotHost(fig, 1, 1, 1, grid_helper=grid_helper)
71+
72+
# make ticklabels of right invisible, and top axis visible.
73+
visible = True
74+
ax.axis[:].major_ticklabels.set_visible(visible)
75+
ax.axis[:].major_ticks.set_visible(False)
76+
ax.axis[:].invert_ticklabel_direction()
77+
78+
ax.axis["wnxneg"] = axis = ax.new_floating_axis(0, 180)
79+
axis.set_ticklabel_direction("-")
80+
axis.label.set_visible(False)
81+
ax.axis["wnxpos"] = axis = ax.new_floating_axis(0, 0)
82+
axis.label.set_visible(False)
83+
ax.axis["wnypos"] = axis = ax.new_floating_axis(0, 90)
84+
axis.label.set_visible(False)
85+
axis.set_axis_direction("left")
86+
ax.axis["wnyneg"] = axis = ax.new_floating_axis(0, 270)
87+
axis.label.set_visible(False)
88+
axis.set_axis_direction("left")
89+
axis.invert_ticklabel_direction()
90+
axis.set_ticklabel_direction("-")
91+
92+
# let left axis shows ticklabels for 1st coordinate (angle)
93+
ax.axis["left"].get_helper().nth_coord_ticks = 0
94+
ax.axis["right"].get_helper().nth_coord_ticks = 0
95+
ax.axis["left"].get_helper().nth_coord_ticks = 0
96+
ax.axis["bottom"].get_helper().nth_coord_ticks = 0
97+
98+
fig.add_subplot(ax)
99+
100+
### RECTANGULAR X Y AXES WITH SCALE
101+
#par2 = ax.twiny()
102+
#par2.axis["top"].toggle(all=False)
103+
#par2.axis["right"].toggle(all=False)
104+
#new_fixed_axis = par2.get_grid_helper().new_fixed_axis
105+
#par2.axis["left"] = new_fixed_axis(loc="left",
106+
# axes=par2,
107+
# offset=(0, 0))
108+
#par2.axis["bottom"] = new_fixed_axis(loc="bottom",
109+
# axes=par2,
110+
# offset=(0, 0))
111+
### FINISH RECTANGULAR
112+
113+
ax.grid(True, zorder=0,linestyle='dotted')
114+
115+
_final_setup(ax)
116+
return ax, fig
117+
118+
def _final_setup(ax):
119+
ax.set_xlabel('Real')
120+
ax.set_ylabel('Imaginary')
121+
ax.axhline(y=0, color='black', lw=1)
122+
ax.axvline(x=0, color='black', lw=1)
123+
plt.axis('equal')
124+
125+
def nogrid():
126+
f = plt.figure()
127+
ax = plt.axes()
128+
129+
_final_setup(ax)
130+
return ax, f
131+
132+
def zgrid(zetas=None, wns=None):
133+
'''Draws discrete damping and frequency grid'''
134+
135+
fig = plt.figure()
136+
ax = fig.gca()
137+
138+
# Constant damping lines
139+
if zetas is None:
140+
zetas = linspace(0, 0.9, 10)
141+
for zeta in zetas:
142+
# Calculate in polar coordinates
143+
factor = zeta/sqrt(1-zeta**2)
144+
x = linspace(0, sqrt(1-zeta**2),200)
145+
ang = pi*x
146+
mag = exp(-pi*factor*x)
147+
# Draw upper part in retangular coordinates
148+
xret = mag*cos(ang)
149+
yret = mag*sin(ang)
150+
ax.plot(xret,yret, 'k:', lw=1)
151+
# Draw lower part in retangular coordinates
152+
xret = mag*cos(-ang)
153+
yret = mag*sin(-ang)
154+
ax.plot(xret,yret,'k:', lw=1)
155+
# Annotation
156+
an_i = int(len(xret)/2.5)
157+
an_x = xret[an_i]
158+
an_y = yret[an_i]
159+
ax.annotate(str(round(zeta,2)), xy=(an_x, an_y), xytext=(an_x, an_y), size=7)
160+
161+
# Constant natural frequency lines
162+
if wns is None:
163+
wns = linspace(0, 1, 10)
164+
for a in wns:
165+
# Calculate in polar coordinates
166+
x = linspace(-pi/2,pi/2,200)
167+
ang = pi*a*sin(x)
168+
mag = exp(-pi*a*cos(x))
169+
# Draw in retangular coordinates
170+
xret = mag*cos(ang)
171+
yret = mag*sin(ang)
172+
ax.plot(xret,yret,'k:', lw=1)
173+
# Annotation
174+
an_i = -1
175+
an_x = xret[an_i]
176+
an_y = yret[an_i]
177+
num = '{:1.1f}'.format(a)
178+
ax.annotate("$\\frac{"+num+"\pi}{T}$", xy=(an_x, an_y), xytext=(an_x, an_y), size=9)
179+
180+
_final_setup(ax)
181+
return ax, fig
182+

‎control/pzmap.py‎

Lines changed: 16 additions & 43 deletions
Original file line numberDiff line numberDiff line change
@@ -40,46 +40,17 @@
4040
#
4141
# $Id:pzmap.py 819 2009-05-29 21:28:07Z murray $
4242

43-
from numpy import real, imag, linspace, pi, exp, cos, sin, sqrt
43+
from numpy import real, imag, linspace, exp, cos, sin, sqrt
44+
from math import pi
4445
from .lti import LTI, isdtime, isctime
46+
from .grid import sgrid, zgrid, nogrid
4547

4648
__all__ = ['pzmap']
4749

48-
def zgrid(plt):
49-
for zeta in linspace(0, 0.9,10):
50-
factor = zeta/sqrt(1-zeta**2)
51-
x = linspace(0, sqrt(1-zeta**2),200)
52-
ang = pi*x
53-
mag = exp(-pi*factor*x)
54-
xret = mag*cos(ang)
55-
yret = mag*sin(ang)
56-
plt.plot(xret,yret, 'k:', lw=1)
57-
xret = mag*cos(-ang)
58-
yret = mag*sin(-ang)
59-
plt.plot(xret,yret,'k:', lw=1)
60-
an_i = int(len(xret)/2.5)
61-
an_x = xret[an_i]
62-
an_y = yret[an_i]
63-
plt.annotate(str(zeta), xy=(an_x, an_y), xytext=(an_x, an_y), size=7)
64-
65-
for a in linspace(0,1,10):
66-
x = linspace(-pi/2,pi/2,200)
67-
ang = pi*a*sin(x)
68-
mag = exp(-pi*a*cos(x))
69-
xret = mag*cos(ang)
70-
yret = mag*sin(ang)
71-
plt.plot(xret,yret,'k:', lw=1)
72-
an_i = -1 #int(len(xret)/2.5)
73-
an_x = xret[an_i]
74-
an_y = yret[an_i]
75-
num = '{:1.1f}'.format(a)
76-
plt.annotate("$\\frac{"+num+"\pi}{T}$", xy=(an_x, an_y), xytext=(an_x, an_y), size=8)
77-
78-
7950
# TODO: Implement more elegant cross-style axes. See:
8051
# http://matplotlib.sourceforge.net/examples/axes_grid/demo_axisline_style.html
8152
# http://matplotlib.sourceforge.net/examples/axes_grid/demo_curvelinear_grid.html
82-
def pzmap(sys, Plot=True, title='Pole Zero Map'):
53+
def pzmap(sys, Plot=True, grid=False, title='Pole Zero Map'):
8354
"""
8455
Plot a pole/zero map for a linear system.
8556
@@ -90,6 +61,8 @@ def pzmap(sys, Plot=True, title='Pole Zero Map'):
9061
Plot: bool
9162
If ``True`` a graph is generated with Matplotlib,
9263
otherwise the poles and zeros are only computed and returned.
64+
grid: boolean (default = False)
65+
If True plot omega-damping grid.
9366
9467
Returns
9568
-------
@@ -107,21 +80,21 @@ def pzmap(sys, Plot=True, title='Pole Zero Map'):
10780
if (Plot):
10881
import matplotlib.pyplot as plt
10982

110-
if isdtime(sys):
111-
zgrid(plt)
83+
if grid:
84+
if isdtime(sys, strict=True):
85+
ax, fig = zgrid()
86+
else:
87+
ax, fig = sgrid()
88+
else:
89+
ax, fig = nogrid()
11290

11391
# Plot the locations of the poles and zeros
11492
if len(poles) > 0:
115-
plt.scatter(real(poles), imag(poles), s=50, marker='x', facecolors='k')
93+
ax.scatter(real(poles), imag(poles), s=50, marker='x', facecolors='k')
11694
if len(zeros) > 0:
117-
plt.scatter(real(zeros), imag(zeros), s=50, marker='o',
95+
ax.scatter(real(zeros), imag(zeros), s=50, marker='o',
11896
facecolors='none', edgecolors='k')
119-
# Add axes
120-
#Somewhat silly workaround
121-
plt.axhline(y=0, color='black', lw=1)
122-
plt.axvline(x=0, color='black', lw=1)
123-
plt.xlabel('Re')
124-
plt.ylabel('Im')
97+
12598

12699
plt.title(title)
127100

0 commit comments

Comments
 (0)