Skip to content

Commit 11a6264

Browse files
author
Jed Frey
authored
Merge pull request #1 from murrayrm/sstf_input_sanitization
Updates to input sanitization
2 parents d44489f + 68b5dd1 commit 11a6264

18 files changed

+442
-137
lines changed

‎.gitignore‎

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,7 @@ __conda_*.txt
99
record.txt
1010
build.log
1111
*.egg-info/
12+
.eggs/
1213
.coverage
1314
doc/_build
1415
doc/generated
@@ -18,4 +19,7 @@ examples/.ipynb_checkpoints/
1819
.project
1920
Untitled*.ipynb
2021
*.idea/
21-
.eggs
22+
23+
# Files created by or for emacs (RMM, 29 Dec 2017)
24+
*~
25+
TAGS

‎.travis.yml‎

Lines changed: 6 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -33,7 +33,7 @@ before_install:
3333
- conda install conda-build
3434
- conda config --add channels python-control
3535
- conda info -a
36-
- conda create -q -n test-environment python="$TRAVIS_PYTHON_VERSION" pip coverage slycot
36+
- conda create -q -n test-environment python="$TRAVIS_PYTHON_VERSION" pip coverage
3737
- source activate test-environment
3838
# coveralls not in conda repos
3939
- pip install coveralls
@@ -45,6 +45,11 @@ install:
4545

4646
# command to run tests
4747
script:
48+
# Before installing Slycot
49+
- python setup.py test
50+
51+
# Now, get and use Slycot
52+
- conda install slycot
4853
- coverage run setup.py test
4954

5055
after_success:

‎MANIFEST.in‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3,3 +3,4 @@ include tests/*.py
33
include README.rst
44
include ChangeLog
55
include Pending
6+
include LICENSE

‎control/bdalg.py‎

Lines changed: 6 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -54,6 +54,7 @@
5454
"""
5555

5656
import scipy as sp
57+
import numpy as np
5758
from . import xferfcn as tf
5859
from . import statesp as ss
5960
from . import frdata as frd
@@ -221,18 +222,18 @@ def feedback(sys1, sys2=1, sign=-1):
221222
"""
222223

223224
# Check for correct input types.
224-
if not isinstance(sys1, (int, float, complex, tf.TransferFunction,
225-
ss.StateSpace, frd.FRD)):
225+
if not isinstance(sys1, (int, float, complex, np.number,
226+
tf.TransferFunction, ss.StateSpace, frd.FRD)):
226227
raise TypeError("sys1 must be a TransferFunction, StateSpace " +
227228
"or FRD object, or a scalar.")
228-
if not isinstance(sys2, (int, float, complex, tf.TransferFunction,
229-
ss.StateSpace, frd.FRD)):
229+
if not isinstance(sys2, (int, float, complex, np.number,
230+
tf.TransferFunction, ss.StateSpace, frd.FRD)):
230231
raise TypeError("sys2 must be a TransferFunction, StateSpace " +
231232
"or FRD object, or a scalar.")
232233

233234
# If sys1 is a scalar, convert it to the appropriate LTI type so that we can
234235
# its feedback member function.
235-
if isinstance(sys1, (int, float, complex)):
236+
if isinstance(sys1, (int, float, complex, np.number)):
236237
if isinstance(sys2, tf.TransferFunction):
237238
sys1 = tf._convertToTransferFunction(sys1)
238239
elif isinstance(sys2, ss.StateSpace):

‎control/frdata.py‎

Lines changed: 6 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -49,6 +49,7 @@
4949
"""
5050

5151
# External function declarations
52+
import numpy as np
5253
from numpy import angle, array, empty, ones, \
5354
real, imag, matrix, absolute, eye, linalg, where, dot
5455
from scipy.interpolate import splprep, splev
@@ -219,7 +220,7 @@ def __mul__(self, other):
219220
"""Multiply two LTI objects (serial connection)."""
220221

221222
# Convert the second argument to a transfer function.
222-
if isinstance(other, (int, float, complex)):
223+
if isinstance(other, (int, float, complex, np.number)):
223224
return FRD(self.fresp * other, self.omega,
224225
smooth=(self.ifunc is not None))
225226
else:
@@ -245,7 +246,7 @@ def __rmul__(self, other):
245246
"""Right Multiply two LTI objects (serial connection)."""
246247

247248
# Convert the second argument to an frd function.
248-
if isinstance(other, (int, float, complex)):
249+
if isinstance(other, (int, float, complex, np.number)):
249250
return FRD(self.fresp * other, self.omega,
250251
smooth=(self.ifunc is not None))
251252
else:
@@ -272,7 +273,7 @@ def __rmul__(self, other):
272273
def __truediv__(self, other):
273274
"""Divide two LTI objects."""
274275

275-
if isinstance(other, (int, float, complex)):
276+
if isinstance(other, (int, float, complex, np.number)):
276277
return FRD(self.fresp * (1/other), self.omega,
277278
smooth=(self.ifunc is not None))
278279
else:
@@ -295,7 +296,7 @@ def __div__(self, other):
295296
# TODO: Division of MIMO transfer function objects is not written yet.
296297
def __rtruediv__(self, other):
297298
"""Right divide two LTI objects."""
298-
if isinstance(other, (int, float, complex)):
299+
if isinstance(other, (int, float, complex, np.number)):
299300
return FRD(other / self.fresp, self.omega,
300301
smooth=(self.ifunc is not None))
301302
else:
@@ -449,7 +450,7 @@ def _convertToFRD(sys, omega, inputs=1, outputs=1):
449450

450451
return FRD(fresp, omega, smooth=True)
451452

452-
elif isinstance(sys, (int, float, complex)):
453+
elif isinstance(sys, (int, float, complex, np.number)):
453454
fresp = ones((outputs, inputs, len(omega)), dtype=float)*sys
454455
return FRD(fresp, omega, smooth=True)
455456

‎control/lti.py‎

Lines changed: 5 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@
1212
timebaseEqual()
1313
"""
1414

15+
import numpy as np
1516
from numpy import absolute, real
1617

1718
__all__ = ['issiso', 'timebase', 'timebaseEqual', 'isdtime', 'isctime',
@@ -96,7 +97,7 @@ def dcgain(self):
9697

9798
# Test to see if a system is SISO
9899
def issiso(sys, strict=False):
99-
if isinstance(sys, (int, float, complex)) and not strict:
100+
if isinstance(sys, (int, float, complex, np.number)) and not strict:
100101
return True
101102
elif not isinstance(sys, LTI):
102103
raise ValueError("Object is not an LTI system")
@@ -114,7 +115,7 @@ def timebase(sys, strict=True):
114115
set to False, dt = True will be returned as 1.
115116
"""
116117
# System needs to be either a constant or an LTI system
117-
if isinstance(sys, (int, float, complex)):
118+
if isinstance(sys, (int, float, complex, np.number)):
118119
return None
119120
elif not isinstance(sys, LTI):
120121
raise ValueError("Timebase not defined")
@@ -162,7 +163,7 @@ def isdtime(sys, strict=False):
162163
"""
163164

164165
# Check to see if this is a constant
165-
if isinstance(sys, (int, float, complex)):
166+
if isinstance(sys, (int, float, complex, np.number)):
166167
# OK as long as strict checking is off
167168
return True if not strict else False
168169

@@ -187,7 +188,7 @@ def isctime(sys, strict=False):
187188
"""
188189

189190
# Check to see if this is a constant
190-
if isinstance(sys, (int, float, complex)):
191+
if isinstance(sys, (int, float, complex, np.number)):
191192
# OK as long as strict checking is off
192193
return True if not strict else False
193194

‎control/margins.py‎

Lines changed: 38 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -2,10 +2,11 @@
22
33
Functions for computing stability margins and related functions.
44
5-
Routeins in this module:
5+
Routines in this module:
66
77
margin.stability_margins
88
margin.phase_crossover_frequencies
9+
margin.margin
910
"""
1011

1112
# Python 3 compatibility (needs to go here)
@@ -87,7 +88,17 @@ def _polysqr(pol):
8788
# RvP, July 8, 2015, augmented to calculate all phase/gain crossings with
8889
# frd data. Correct to return smallest phase
8990
# margin, smallest gain margin and their frequencies
90-
def stability_margins(sysdata, returnall=False, epsw=1e-8):
91+
# RvP, Jun 10, 2017, modified the inclusion of roots found for phase
92+
# crossing to include all >= 0, made subsequent calc
93+
# insensitive to div by 0
94+
# also changed the selection of which crossings to
95+
# return on basis of "A note on the Gain and Phase
96+
# Margin Concepts" Journal of Control and Systems
97+
# Engineering, Yazdan Bavafi-Toosi, Dec 2015, vol 3
98+
# issue 1, pp 51-59, closer to Matlab behavior, but
99+
# not completely identical in edge cases, which don't
100+
# cross but touch gain=1
101+
def stability_margins(sysdata, returnall=False, epsw=0.0):
91102
"""Calculate stability margins and associated crossover frequencies.
92103
93104
Parameters
@@ -104,7 +115,7 @@ def stability_margins(sysdata, returnall=False, epsw=1e-8):
104115
minimum stability margins. For frequency data or FRD systems, only one
105116
margin is found and returned.
106117
epsw: float, optional
107-
Frequencies below this value (default 1e-8) are considered static gain,
118+
Frequencies below this value (default 0.0) are considered static gain,
108119
and not returned as margin.
109120
110121
Returns
@@ -161,12 +172,13 @@ def stability_margins(sysdata, returnall=False, epsw=1e-8):
161172
#print ('2:w_180', w_180)
162173

163174
# evaluate response at remaining frequencies, to test for phase 180 vs 0
164-
resp_w_180 = np.real(np.polyval(sys.num[0][0], 1.j*w_180) /
165-
np.polyval(sys.den[0][0], 1.j*w_180))
166-
#print ('resp_w_180', resp_w_180)
175+
with np.errstate(all='ignore'):
176+
resp_w_180 = np.real(
177+
np.polyval(sys.num[0][0], 1.j*w_180) /
178+
np.polyval(sys.den[0][0], 1.j*w_180))
167179

168180
# only keep frequencies where the negative real axis is crossed
169-
w_180 = w_180[np.real(resp_w_180) < 0.0]
181+
w_180 = w_180[np.real(resp_w_180) <= 0.0]
170182

171183
# and sort
172184
w_180.sort()
@@ -253,20 +265,30 @@ def dstab(w):
253265

254266
# margins, as iterables, converted frdata and xferfcn calculations to
255267
# vector for this
256-
GM = 1/np.abs(sys.evalfr(w_180)[0][0])
268+
with np.errstate(all='ignore'):
269+
gain_w_180 = np.abs(sys.evalfr(w_180)[0][0])
270+
GM = 1.0/gain_w_180
257271
SM = np.abs(sys.evalfr(wstab)[0][0]+1)
258-
PM = np.angle(sys.evalfr(wc)[0][0], deg=True) + 180
259-
272+
PM = np.remainder(np.angle(sys.evalfr(wc)[0][0], deg=True), 360.0) - 180.0
273+
260274
if returnall:
261275
return GM, PM, SM, w_180, wc, wstab
262276
else:
277+
if GM.shape[0] and not np.isinf(GM).all():
278+
with np.errstate(all='ignore'):
279+
gmidx = np.where(np.abs(np.log(GM)) ==
280+
np.min(np.abs(np.log(GM))))
281+
else:
282+
gmidx = -1
283+
if PM.shape[0]:
284+
pmidx = np.where(np.abs(PM) == np.amin(np.abs(PM)))[0]
263285
return (
264-
(GM.shape[0] or None) and np.amin(GM),
265-
(PM.shape[0] or None) and np.amin(PM),
266-
(SM.shape[0] or None) and np.amin(SM),
267-
(w_180.shape[0] or None) and w_180[GM==np.amin(GM)][0],
268-
(wc.shape[0] or None) and wc[PM==np.amin(PM)][0],
269-
(wstab.shape[0] or None) and wstab[SM==np.amin(SM)][0])
286+
(not gmidx != -1 and float('inf')) or GM[gmidx][0],
287+
(not PM.shape[0] and float('inf')) or PM[pmidx][0],
288+
(not SM.shape[0] and float('inf')) or np.amin(SM),
289+
(not gmidx != -1 and float('nan')) or w_180[gmidx][0],
290+
(not wc.shape[0] and float('nan')) or wc[pmidx][0],
291+
(not wstab.shape[0] and float('nan')) or wstab[SM==np.amin(SM)][0])
270292

271293

272294
# Contributed by Steffen Waldherr <waldherr@ist.uni-stuttgart.de>

‎control/statesp.py‎

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -58,6 +58,7 @@
5858
from numpy.random import rand, randn
5959
from numpy.linalg import solve, eigvals, matrix_rank
6060
from numpy.linalg.linalg import LinAlgError
61+
import scipy as sp
6162
from scipy.signal import lti, cont2discrete
6263
# from exceptions import Exception
6364
import warnings
@@ -218,7 +219,7 @@ def __add__(self, other):
218219
"""Add two LTI systems (parallel connection)."""
219220

220221
# Check for a couple of special cases
221-
if (isinstance(other, (int, float, complex))):
222+
if (isinstance(other, (int, float, complex, np.number))):
222223
# Just adding a scalar; put it in the D matrix
223224
A, B, C = self.A, self.B, self.C;
224225
D = self.D + other;
@@ -275,7 +276,7 @@ def __mul__(self, other):
275276
"""Multiply two LTI objects (serial connection)."""
276277

277278
# Check for a couple of special cases
278-
if isinstance(other, (int, float, complex)):
279+
if isinstance(other, (int, float, complex, np.number)):
279280
# Just multiplying by a scalar; change the output
280281
A, B = self.A, self.B
281282
C = self.C * other
@@ -316,7 +317,7 @@ def __rmul__(self, other):
316317
"""Right multiply two LTI objects (serial connection)."""
317318

318319
# Check for a couple of special cases
319-
if isinstance(other, (int, float, complex)):
320+
if isinstance(other, (int, float, complex, np.number)):
320321
# Just multiplying by a scalar; change the input
321322
A, C = self.A, self.C;
322323
B = self.B * other;
@@ -699,11 +700,10 @@ def _convertToStateSpace(sys, **kw):
699700
# TODO: do we want to squeeze first and check dimenations?
700701
# I think this will fail if num and den aren't 1-D after
701702
# the squeeze
702-
lti_sys = lti(squeeze(sys.num), squeeze(sys.den))
703-
return StateSpace(lti_sys.A, lti_sys.B, lti_sys.C, lti_sys.D,
704-
sys.dt)
703+
A, B, C, D = sp.signal.tf2ss(squeeze(sys.num), squeeze(sys.den))
704+
return StateSpace(A, B, C, D, sys.dt)
705705

706-
elif isinstance(sys, (int, float, complex)):
706+
elif isinstance(sys, (int, float, complex, np.number)):
707707
if "inputs" in kw:
708708
inputs = kw["inputs"]
709709
else:

‎control/tests/convert_test.py‎

Lines changed: 23 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -206,6 +206,29 @@ def testTf2ssStaticMimo(self):
206206
d = np.matrix([[0.5, 30, 0.0625], [-0.5, -1.25, 101.3]])
207207
np.testing.assert_array_equal(d, gmimo.D)
208208

209+
def testSs2tfStaticSiso(self):
210+
"""Regression: ss2tf for SISO static gain"""
211+
import control
212+
gsiso = control.ss2tf(control.ss([], [], [], 0.5))
213+
np.testing.assert_array_equal([[[0.5]]], gsiso.num)
214+
np.testing.assert_array_equal([[[1.]]], gsiso.den)
215+
216+
def testSs2tfStaticMimo(self):
217+
"""Regression: ss2tf for MIMO static gain"""
218+
import control
219+
# 2x3 TFM
220+
a = []
221+
b = []
222+
c = []
223+
d = np.matrix([[0.5, 30, 0.0625], [-0.5, -1.25, 101.3]])
224+
gtf = control.ss2tf(control.ss(a,b,c,d))
225+
226+
# we need a 3x2x1 array to compare with gtf.num
227+
# np.testing.assert_array_equal doesn't seem to like a matrices
228+
# with an extra dimension, so convert to ndarray
229+
numref = np.asarray(d)[...,np.newaxis]
230+
np.testing.assert_array_equal(numref, np.array(gtf.num) / np.array(gtf.den))
231+
209232

210233
def suite():
211234
return unittest.TestLoader().loadTestsFromTestCase(TestConvert)

0 commit comments

Comments
 (0)