Skip to content

Commit cd87eda

Browse files
murrayrmmurrayrm
committed
Merge pull request #62 from cwrowley/matlab
Move matlab compatibility routines to subpackage. Reviewed all changes and confirmed that these should give the same functionality, modulo the changes @cwrowley noted in the PR.
2 parents 58c7f85 + 06e5c34 commit cd87eda

31 files changed

+1829
-1657
lines changed

‎control/__init__.py

Lines changed: 19 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -45,28 +45,26 @@
4545
"""
4646

4747
# Import functions from within the control system library
48-
# Should probably only import the exact functions we use...
49-
from .bdalg import series, parallel, negate, feedback
50-
from .delay import pade
51-
from .dtime import sample_system
52-
from .freqplot import bode_plot, nyquist_plot, gangof4_plot
53-
from .freqplot import bode, nyquist, gangof4
54-
from .lti import issiso, timebase, timebaseEqual, isdtime, isctime
55-
from .margins import stability_margins, phase_crossover_frequencies
56-
from .mateqn import lyap, dlyap, care, dare
57-
from .modelsimp import hsvd, modred, balred, era, markov, minreal
58-
from .nichols import nichols_plot, nichols
59-
from .phaseplot import phase_plot, box_grid
60-
from .pzmap import pzmap
61-
from .rlocus import root_locus
62-
from .statefbk import place, lqr, ctrb, obsv, gram, acker
63-
from .statesp import StateSpace
64-
from .timeresp import forced_response, initial_response, step_response, \
65-
impulse_response
66-
from .xferfcn import TransferFunction
48+
# Note: the functions we use are specified as __all__ variables in the modules
49+
from .bdalg import *
50+
from .delay import *
51+
from .dtime import *
52+
from .freqplot import *
53+
from .lti import *
54+
from .margins import *
55+
from .mateqn import *
56+
from .modelsimp import *
57+
from .nichols import *
58+
from .phaseplot import *
59+
from .pzmap import *
60+
from .rlocus import *
61+
from .statefbk import *
62+
from .statesp import *
63+
from .timeresp import *
64+
from .xferfcn import *
6765
from .ctrlutil import *
68-
from .frdata import FRD
69-
from .canonical import canonical_form, reachable_form
66+
from .frdata import *
67+
from .canonical import *
7068

7169
# Exceptions
7270
from .exception import *
@@ -77,22 +75,6 @@
7775
except ImportError:
7876
__version__ = "dev"
7977

80-
# Import some of the more common (and benign) MATLAB shortcuts
81-
# By default, don't import conflicting commands here
82-
#! TODO (RMM, 4 Nov 2012): remove MATLAB dependencies from __init__.py
83-
#!
84-
#! Eventually, all functionality should be in modules *other* than matlab.
85-
#! This will allow inclusion of the matlab module to set up a different set
86-
#! of defaults from the main package. At that point, the matlab module will
87-
#! allow provide compatibility with MATLAB but no package functionality.
88-
#!
89-
from .matlab import ss, tf, ss2tf, tf2ss, drss
90-
from .matlab import pole, zero, evalfr, freqresp, dcgain
91-
from .matlab import nichols, rlocus, margin
92-
# bode and nyquist come directly from freqplot.py
93-
from .matlab import step, impulse, initial, lsim
94-
from .matlab import ssdata, tfdata
95-
9678
# The following is to use Numpy's testing framework
9779
# Tests go under directory tests/, benchmarks under directory benchmarks/
9880
from numpy.testing import Tester

‎control/bdalg.py

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -58,6 +58,8 @@
5858
from . import statesp as ss
5959
from . import frdata as frd
6060

61+
__all__ = ['series', 'parallel', 'negate', 'feedback', 'append', 'connect']
62+
6163
def series(sys1, sys2):
6264
"""Return the series connection sys2 * sys1 for --> sys1 --> sys2 -->.
6365

‎control/canonical.py

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,7 @@
99
from numpy import zeros, shape, poly
1010
from numpy.linalg import inv
1111

12+
__all__ = ['canonical_form', 'reachable_form']
1213

1314
def canonical_form(xsys, form='reachable'):
1415
"""Convert a system into canonical form

‎control/delay.py

Lines changed: 15 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -3,7 +3,7 @@
33
#
44
# Author: Sawyer Fuller
55
# Date: 26 Aug 2010
6-
#
6+
#
77
# This file contains functions for implementing time delays (currently
88
# only the pade() function).
99
#
@@ -16,16 +16,16 @@
1616
#
1717
# 1. Redistributions of source code must retain the above copyright
1818
# notice, this list of conditions and the following disclaimer.
19-
#
19+
#
2020
# 2. Redistributions in binary form must reproduce the above copyright
2121
# notice, this list of conditions and the following disclaimer in the
2222
# documentation and/or other materials provided with the distribution.
23-
#
23+
#
2424
# 3. Neither the name of the California Institute of Technology nor
2525
# the names of its contributors may be used to endorse or promote
2626
# products derived from this software without specific prior
2727
# written permission.
28-
#
28+
#
2929
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
3030
# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
3131
# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
@@ -38,30 +38,32 @@
3838
# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
3939
# OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
4040
# SUCH DAMAGE.
41-
#
41+
#
4242
# $Id$
4343

4444
# Python 3 compatability (needs to go here)
4545
from __future__ import print_function
4646

47+
__all__ = ['pade']
48+
4749
def pade(T, n=1):
48-
"""
50+
"""
4951
Create a linear system that approximates a delay.
50-
52+
5153
Return the numerator and denominator coefficients of the Pade approximation.
52-
54+
5355
Parameters
5456
----------
5557
T : number
5658
time delay
5759
n : integer
5860
order of approximation
59-
61+
6062
Returns
61-
-------
63+
-------
6264
num, den : array
6365
Polynomial coefficients of the delay model, in descending powers of s.
64-
66+
6567
Notes
6668
-----
6769
Based on an algorithm in Golub and van Loan, "Matrix Computation" 3rd.
@@ -79,7 +81,7 @@ def pade(T, n=1):
7981
for k in range(1, n+1):
8082
c = T * c * (n - k + 1)/(2 * n - k + 1)/k
8183
num[n - k] = c * (-1)**k
82-
den[n - k] = c
84+
den[n - k] = c
8385
num = [coeff/den[0] for coeff in num]
8486
den = [coeff/den[0] for coeff in den]
85-
return num, den
87+
return num, den

‎control/dtime.py

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -47,6 +47,9 @@
4747
"""
4848

4949
from .lti import isctime
50+
from .statesp import StateSpace, _convertToStateSpace
51+
52+
__all__ = ['sample_system', 'c2d']
5053

5154
# Sample a continuous time system
5255
def sample_system(sysc, Ts, method='zoh', alpha=None):
@@ -85,3 +88,30 @@ def sample_system(sysc, Ts, method='zoh', alpha=None):
8588
raise ValueError("First argument must be continuous time system")
8689

8790
return sysc.sample(Ts, method, alpha)
91+
92+
93+
def c2d(sysc, Ts, method='zoh'):
94+
'''
95+
Return a discrete-time system
96+
97+
Parameters
98+
----------
99+
sysc: LTI (StateSpace or TransferFunction), continuous
100+
System to be converted
101+
102+
Ts: number
103+
Sample time for the conversion
104+
105+
method: string, optional
106+
Method to be applied,
107+
'zoh' Zero-order hold on the inputs (default)
108+
'foh' First-order hold, currently not implemented
109+
'impulse' Impulse-invariant discretization, currently not implemented
110+
'tustin' Bilinear (Tustin) approximation, only SISO
111+
'matched' Matched pole-zero method, only SISO
112+
'''
113+
# Call the sample_system() function to do the work
114+
sysd = sample_system(sysc, Ts, method)
115+
if isinstance(sysc, StateSpace) and not isinstance(sysd, StateSpace):
116+
return _convertToStateSpace(sysd)
117+
return sysd

‎control/frdata.py

Lines changed: 40 additions & 28 deletions
Original file line numberDiff line numberDiff line change
@@ -1,35 +1,9 @@
11
from __future__ import division
2-
"""frdata.py
3-
2+
"""
43
Frequency response data representation and functions.
54
6-
This file contains the FRD class and also functions that operate on
5+
This module contains the FRD class and also functions that operate on
76
FRD data.
8-
9-
Routines in this module:
10-
11-
FRD.__init__
12-
FRD.copy
13-
FRD.__str__
14-
FRD.__neg__
15-
FRD.__add__
16-
FRD.__radd__
17-
FRD.__sub__
18-
FRD.__rsub__
19-
FRD.__mul__
20-
FRD.__rmul__
21-
FRD.__div__
22-
FRD.__rdiv__
23-
FRD.__truediv__
24-
FRD.__rtruediv__
25-
FRD.evalfr
26-
FRD.freqresp
27-
FRD.pole
28-
FRD.zero
29-
FRD.feedback
30-
FRD._common_den
31-
_convertToFRD
32-
337
"""
348

359
"""Copyright (c) 2010 by California Institute of Technology
@@ -80,6 +54,8 @@
8054
from scipy.interpolate import splprep, splev
8155
from .lti import LTI
8256

57+
__all__ = ['FRD', 'frd']
58+
8359
class FRD(LTI):
8460
"""A class for models defined by Frequency Response Data (FRD)
8561
@@ -480,3 +456,39 @@ def _convertToFRD(sys, omega, inputs=1, outputs=1):
480456

481457
raise TypeError('''Can't convert given type "%s" to FRD system.''' %
482458
sys.__class__)
459+
460+
def frd(*args):
461+
'''
462+
Construct a Frequency Response Data model, or convert a system
463+
464+
frd models store the (measured) frequency response of a system.
465+
466+
This function can be called in different ways:
467+
468+
``frd(response, freqs)``
469+
Create an frd model with the given response data, in the form of
470+
complex response vector, at matching frequency freqs [in rad/s]
471+
472+
``frd(sys, freqs)``
473+
Convert an LTI system into an frd model with data at frequencies
474+
freqs.
475+
476+
Parameters
477+
----------
478+
response: array_like, or list
479+
complex vector with the system response
480+
freq: array_lik or lis
481+
vector with frequencies
482+
sys: LTI (StateSpace or TransferFunction)
483+
A linear system
484+
485+
Returns
486+
-------
487+
sys: FRD
488+
New frequency response system
489+
490+
See Also
491+
--------
492+
ss, tf
493+
'''
494+
return FRD(*args)

‎control/freqplot.py

Lines changed: 14 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -49,6 +49,9 @@
4949
from .bdalg import feedback
5050
from .lti import isdtime, timebaseEqual
5151

52+
__all__ = ['bode_plot', 'nyquist_plot', 'gangof4_plot',
53+
'bode', 'nyquist', 'gangof4']
54+
5255
#
5356
# Main plotting functions
5457
#
@@ -74,7 +77,7 @@ def bode_plot(syslist, omega=None, dB=None, Hz=None, deg=None,
7477
Hz : boolean
7578
If True, plot frequency in Hz (omega must be provided in rad/sec)
7679
deg : boolean
77-
If True, return phase in degrees (else radians)
80+
If True, plot phase in degrees (else radians)
7881
Plot : boolean
7982
If True, plot magnitude and phase
8083
*args, **kwargs:
@@ -85,9 +88,9 @@ def bode_plot(syslist, omega=None, dB=None, Hz=None, deg=None,
8588
mag : array (list if len(syslist) > 1)
8689
magnitude
8790
phase : array (list if len(syslist) > 1)
88-
phase
91+
phase in radians
8992
omega : array (list if len(syslist) > 1)
90-
frequency
93+
frequency in rad/sec
9194
9295
Notes
9396
-----
@@ -126,16 +129,15 @@ def bode_plot(syslist, omega=None, dB=None, Hz=None, deg=None,
126129
omega = default_frequency_range(syslist)
127130

128131
# Get the magnitude and phase of the system
132+
omega = np.array(omega)
129133
mag_tmp, phase_tmp, omega_sys = sys.freqresp(omega)
130134
mag = np.atleast_1d(np.squeeze(mag_tmp))
131135
phase = np.atleast_1d(np.squeeze(phase_tmp))
132136
phase = unwrap(phase)
133137
if Hz:
134-
omega_plot = omega_sys/(2*sp.pi)
138+
omega_plot = omega_sys / (2 * np.pi)
135139
else:
136140
omega_plot = omega_sys
137-
if dB: mag = 20*sp.log10(mag)
138-
if deg: phase = phase * 180 / sp.pi
139141

140142
mags.append(mag)
141143
phases.append(phase)
@@ -147,7 +149,7 @@ def bode_plot(syslist, omega=None, dB=None, Hz=None, deg=None,
147149
# Magnitude plot
148150
plt.subplot(211);
149151
if dB:
150-
plt.semilogx(omega_plot, mag, *args, **kwargs)
152+
plt.semilogx(omega_plot, 20 * np.log10(mag), *args, **kwargs)
151153
else:
152154
plt.loglog(omega_plot, mag, *args, **kwargs)
153155
plt.hold(True);
@@ -159,7 +161,11 @@ def bode_plot(syslist, omega=None, dB=None, Hz=None, deg=None,
159161

160162
# Phase plot
161163
plt.subplot(212);
162-
plt.semilogx(omega_plot, phase, *args, **kwargs)
164+
if deg:
165+
phase_plot = phase * 180 / np.pi
166+
else:
167+
phase_plot = phase
168+
plt.semilogx(omega_plot, phase_plot, *args, **kwargs)
163169
plt.hold(True);
164170

165171
# Add a grid to the plot + labeling

0 commit comments

Comments
 (0)