matlab.py
MATLAB emulation functions.
This file contains a number of functions that emulate some of the functionality of MATLAB. The intent of these functions is to provide a simple interface to the python control systems library (pythoncontrol) for people who are familiar with the MATLAB Control Systems Toolbox (tm). Most of the functions are just calls to pythoncontrol functions defined elsewhere. Use ‘from control.matlab import *’ in python to include all of the functions defined here. Functions that are defined in other libraries that have the same names as their MATLAB equivalents are automatically imported here.
The following tables give an overview of the module control.matlab. They also show the implementation progress and the planned features of the module.
The symbols in the first column show the current state of a feature:
*  tf()  create transfer function (TF) models 
zpk  create zero/pole/gain (ZPK) models.  
*  ss()  create statespace (SS) models 
dss  create descriptor statespace models  
delayss  create statespace models with delayed terms  
*  frd()  create frequency response data (FRD) models 
lti/exp  create pure continuoustime delays (TF and ZPK only)  
filt  specify digital filters  
  lti/set  set/modify properties of LTI models 
  setdelaymodel  specify internal delay model (state space only) 
lti/tfdata  extract numerators and denominators  
lti/zpkdata  extract zero/pole/gain data  
lti/ssdata  extract statespace matrices  
lti/dssdata  descriptor version of SSDATA  
frd/frdata  extract frequency response data  
lti/get  access values of LTI model properties  
ss/getDelayModel  access internal delay model (state space) 
*  tf()  conversion to transfer function 
zpk  conversion to zero/pole/gain  
*  ss()  conversion to state space 
*  frd()  conversion to frequency data 
c2d  continuous to discrete conversion  
d2c  discrete to continuous conversion  
d2d  resample discretetime model  
upsample  upsample discretetime LTI systems  
*  ss2tf()  state space to transfer function 
s  ss2zpk  transfer function to zeropolegain 
*  tf2ss()  transfer function to state space 
s  tf2zpk  transfer function to zeropolegain 
s  zpk2ss  zeropolegain to state space 
s  zpk2tf  zeropolegain to transfer function 
*  append()  group LTI models by appending inputs/outputs 
*  parallel()  connect LTI models in parallel (see also overloaded +) 
*  series()  connect LTI models in series (see also overloaded *) 
*  feedback()  connect lti models with a feedback loop 
lti/lft  generalized feedback interconnection  
lti/connect  arbitrary interconnection of lti models  
sumblk  summing junction (for use with connect)  
strseq  builds sequence of indexed strings (for I/O naming) 
*  dcgain()  steadystate (D.C.) gain 
lti/bandwidth  system bandwidth  
lti/norm  h2 and Hinfinity norms of LTI models  
*  pole()  system poles 
*  zero()  system (transmission) zeros 
lti/order  model order (number of states)  
*  pzmap()  polezero map (TF only) 
lti/iopzmap  input/output polezero map  
*  damp()  natural frequency, damping of system poles 
esort  sort continuous poles by real part  
dsort  sort discrete poles by magnitude  
lti/stabsep  stable/unstable decomposition  
lti/modsep  regionbased modal decomposition 
*  step()  step response 
stepinfo  step response characteristics  
*  impulse()  impulse response 
*  initial()  free response with initial conditions 
*  lsim()  response to userdefined input signal 
lsiminfo  linear response characteristics  
gensig  generate input signal for LSIM  
covar  covariance of response to white noise 
*  bode()  Bode plot of the frequency response 
lti/bodemag  Bode magnitude diagram only  
sigma  singular value frequency plot  
*  nyquist()  Nyquist plot 
*  nichols()  Nichols plot 
*  margin()  gain and phase margins 
lti/allmargin  all crossover frequencies and margins  
*  freqresp()  frequency response over a frequency grid 
*  evalfr()  frequency response at single frequency 
*  minreal()  minimal realization; pole/zero cancellation 
ss/sminreal  structurally minimal realization  
*  hsvd()  hankel singular values (state contributions) 
*  balred()  reducedorder approximations of LTI models 
*  modred()  model order reduction 
*  rlocus()  evans root locus 
*  place()  pole placement 
estim  form estimator given estimator gain  
reg  form regulator given statefeedback and estimator gains 
ss/lqg  singlestep LQG design  
*  lqr()  linear quadratic (LQ) statefbk regulator 
dlqr  discretetime LQ statefeedback regulator  
lqry  LQ regulator with output weighting  
lqrd  discrete LQ regulator for continuous plant  
ss/lqi  LinearQuadraticIntegral (LQI) controller  
ss/kalman  Kalman state estimator  
ss/kalmd  discrete Kalman estimator for cts plant  
ss/lqgreg  build LQG regulator from LQ gain and Kalman estimator  
ss/lqgtrack  build LQG servocontroller  
augstate  augment output by appending states 
*  rss()  random stable ctstime statespace models 
*  drss()  random stable disctime statespace models 
ss2ss  state coordinate transformation  
canon  canonical forms of statespace models  
*  ctrb()  controllability matrix 
*  obsv()  observability matrix 
*  gram()  controllability and observability gramians 
ss/prescale  optimal scaling of statespace models.  
balreal  gramianbased input/output balancing  
ss/xperm  reorder states. 
frd/chgunits  change frequency vector units  
frd/fcat  merge frequency responses  
frd/fselect  select frequency range or subgrid  
frd/fnorm  peak gain as a function of frequency  
frd/abs  entrywise magnitude of frequency response  
frd/real  real part of the frequency response  
frd/imag  imaginary part of the frequency response  
frd/interp  interpolate frequency response data  
mag2db  convert magnitude to decibels (dB)  
db2mag  convert decibels (dB) to magnitude 
lti/hasdelay  true for models with time delays  
lti/totaldelay  total delay between each input/output pair  
lti/delay2z  replace delays by poles at z=0 or FRD phase shift  
*  pade()  pade approximation of time delays 
class  model type (‘tf’, ‘zpk’, ‘ss’, or ‘frd’)  
isa  test if model is of given type  
tf/size  model sizes  
lti/ndims  number of dimensions  
lti/isempty  true for empty models  
lti/isct  true for continuoustime models  
lti/isdt  true for discretetime models  
lti/isproper  true for proper models  
lti/issiso  true for singleinput/singleoutput models  
lti/isstable  true for models with stable dynamics  
lti/reshape  reshape array of linear models 
*  + and   add, subtract systems (parallel connection) 
*  *  multiply systems (series connection) 
/  right divide – sys1*inv(sys2)  
  \  left divide – inv(sys1)*sys2 
^  powers of a given system  
‘  pertransposition  
.’  transposition of input/output map  
.*  elementbyelement multiplication  
[..]  concatenate models along inputs or outputs  
lti/stack  stack models/arrays along some dimension  
lti/inv  inverse of an LTI system  
lti/conj  complex conjugation of model coefficients 
*  lyap()  solve continuoustime Lyapunov equations 
*  dlyap()  solve discretetime Lyapunov equations 
lyapchol, dlyapchol  squareroot Lyapunov solvers  
*  care()  solve continuoustime algebraic Riccati equations 
*  dare()  solve disctime algebraic Riccati equations 
gcare, gdare  generalized Riccati solvers  
bdschur  block diagonalization of a square matrix 
*  gangof4()  generate the Gang of 4 sensitivity plots 
*  linspace()  generate a set of numbers that are linearly spaced 
*  logspace()  generate a set of numbers that are logarithmically spaced 
*  unwrap()  unwrap phase angle to give continuous curve 
Bode plot of the frequency response
Plots a bode gain and phase diagram
Parameters:  sys : Lti, or list of Lti
omega: freq_range
dB : boolean
Hz : boolean
deg : boolean
Plot : boolean


Examples
>>> sys = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> mag, phase, omega = bode(sys)
Todo
Document these use cases
>>> bode(sys, w)
>>> bode(sys1, sys2, ..., sysN)
>>> bode(sys1, sys2, ..., sysN, w)
>>> bode(sys1, 'plotstyle1', ..., sysN, 'plotstyleN')
Compute natural frequency, damping and poles of a system
The function takes 1 or 2 parameters
Parameters:  sys: Lti (StateSpace or TransferFunction)
doprint:


Returns:  wn: array
damping: array
poles: array

See also
Compute the gain of the system in steady state.
The function takes either 1, 2, 3, or 4 parameters:
Parameters:  A, B, C, D: arraylike
Z, P, k: arraylike, arraylike, number
num, den: arraylike
sys: Lti (StateSpace or TransferFunction)


Returns:  gain: matrix

Notes
This function is only useful for systems with invertible system matrix A.
All systems are first converted to state space form. The function then computes:
Create a stable discrete random state space object.
Parameters:  states: integer
inputs: integer
outputs: integer


Returns:  sys: StateSpace

Raises:  ValueError

See also
Notes
If the number of states, inputs, or outputs is not specified, then the missing numbers are assumed to be 1. The poles of the returned system will always have a magnitude less than 1.
Evaluate the transfer function of an LTI system for a single complex number x.
To evaluate at a frequency, enter x = omega*j, where omega is the frequency in radians
Parameters:  sys: StateSpace or TransferFunction
x: scalar


Returns:  fresp: ndarray 
Notes
This function is a wrapper for StateSpace.evalfr and TransferFunction.evalfr.
Examples
>>> sys = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> evalfr(sys, 1j)
array([[ 44.821.4j]])
>>> # This is the transfer function matrix evaluated at s = i.
Todo
Add example with MIMO system
Construct a Frequency Response Data model, or convert a system
frd models store the (measured) frequency response of a system.
This function can be called in different ways:
Parameters:  response: array_like, or list
freq: array_lik or lis
sys: Lti (StateSpace or TransferFunction)


Returns:  sys: FRD

Frequency response of an LTI system at multiple angular frequencies.
Parameters:  sys: StateSpace or TransferFunction
omega: array_like


Returns:  mag: ndarray phase: ndarray omega: list, tuple, or ndarray 
Notes
This function is a wrapper for StateSpace.freqresp and TransferFunction.freqresp. The output omega is a sorted version of the input omega.
Examples
>>> sys = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> mag, phase, omega = freqresp(sys, [0.1, 1., 10.])
>>> mag
array([[[ 58.8576682 , 49.64876635, 13.40825927]]])
>>> phase
array([[[0.05408304, 0.44563154, 0.66837155]]])
Todo
Add example with MIMO system
#>>> sys = rss(3, 2, 2) #>>> mag, phase, omega = freqresp(sys, [0.1, 1., 10.]) #>>> mag[0, 1, :] #array([ 55.43747231, 42.47766549, 1.97225895]) #>>> phase[1, 0, :] #array([0.12611087, 1.14294316, 2.5764547 ]) #>>> # This is the magnitude of the frequency response from the 2nd #>>> # input to the 1st output, and the phase (in radians) of the #>>> # frequency response from the 1st input to the 2nd output, for #>>> # s = 0.1i, i, 10i.
Impulse response of a linear system
If the system has multiple inputs or outputs (MIMO), one input and one output must be selected for the simulation. The parameters input and output do this. All other inputs are set to 0, all other outputs are ignored.
Parameters:  sys: StateSpace, TransferFunction
T: arraylike object, optional
input: int
output: int
**keywords:


Returns:  yout: array
T: array

Examples
>>> T, yout = impulse(sys, T)
Initial condition response of a linear system
If the system has multiple inputs or outputs (MIMO), one input and one output have to be selected for the simulation. The parameters input and output do this. All other inputs are set to 0, all other outputs are ignored.
Parameters:  sys: StateSpace, or TransferFunction
T: arraylike object, optional
X0: arraylike object or number, optional
input: int
output: int
**keywords:


Returns:  yout: array
T: array

Examples
>>> T, yout = initial(sys, T, X0)
Simulate the output of a linear system.
As a convenience for parameters U, X0: Numbers (scalars) are converted to constant arrays with the correct shape. The correct shape is inferred from arguments sys and T.
Parameters:  sys: Lti (StateSpace, or TransferFunction)
U: arraylike or number, optional
T: arraylike
X0: arraylike or number, optional
**keywords:


Returns:  yout: array
T: array
xout: array

Examples
>>> T, yout, xout = lsim(sys, U, T, X0)
Calculate gain and phase margins and associated crossover frequencies
Function margin takes either 1 or 3 parameters.
Parameters:  sys : StateSpace or TransferFunction
mag, phase, w : array_like


Returns:  gm, pm, Wcg, Wcp : float

Examples
>>> sys = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> gm, pm, wg, wp = margin(sys)
margin: no magnitude crossings found
Todo
better ecample system!
#>>> gm, pm, wg, wp = margin(mag, phase, w)
Nichols chart grid
Parameters:  cl_mags : arraylike (dB), optional
cl_phases : arraylike (degrees), optional


Compute system poles.
Parameters:  sys: StateSpace or TransferFunction


Returns:  poles: ndarray

Raises:  NotImplementedError

See also
Notes
This function is a wrapper for StateSpace.pole and TransferFunction.pole.
Root locus plot
The rootlocus plot has a callback function that prints pole location, gain and damping to the Python consol on mouseclicks on the rootlocus graph.
Parameters:  sys: StateSpace or TransferFunction
klist:


Returns:  rlist:
klist:

Create a stable continuous random state space object.
Parameters:  states: integer
inputs: integer
outputs: integer


Returns:  sys: StateSpace

Raises:  ValueError

See also
Notes
If the number of states, inputs, or outputs is not specified, then the missing numbers are assumed to be 1. The poles of the returned system will always have a negative real part.
Create a state space system.
The function accepts either 1, 4 or 5 parameters:
Create a state space system from the matrices of its state and output equations:
Create a discretetime state space system from the matrices of its state and output equations:
The matrices can be given as array like data types or strings. Everything that the constructor of numpy.matrix accepts is permissible here too.
Parameters:  sys: Lti (StateSpace or TransferFunction)
A: array_like or string
B: array_like or string
C: array_like or string
D: array_like or string
dt: If present, specifies the sampling period and a discrete time


Returns:  out: StateSpace

Raises:  ValueError

Examples
>>> # Create a StateSpace object from four "matrices".
>>> sys1 = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> # Convert a TransferFunction to a StateSpace object.
>>> sys_tf = tf([2.], [1., 3])
>>> sys2 = ss(sys_tf)
Transform a state space system to a transfer function.
The function accepts either 1 or 4 parameters:
Create a state space system from the matrices of its state and output equations.
For details see: ss()
Parameters:  sys: StateSpace
A: array_like or string
B: array_like or string
C: array_like or string
D: array_like or string


Returns:  out: TransferFunction

Raises:  ValueError
TypeError

Examples
>>> A = [[1., 2], [3, 4]]
>>> B = [[5.], [7]]
>>> C = [[6., 8]]
>>> D = [[9.]]
>>> sys1 = ss2tf(A, B, C, D)
>>> sys_ss = ss(A, B, C, D)
>>> sys2 = ss2tf(sys_ss)
Return state space data objects for a system
Parameters:  sys: Lti (StateSpace, or TransferFunction)


Returns:  (A, B, C, D): list of matrices

Step response of a linear system
If the system has multiple inputs or outputs (MIMO), one input and one output have to be selected for the simulation. The parameters input and output do this. All other inputs are set to 0, all other outputs are ignored.
Parameters:  sys: StateSpace, or TransferFunction
T: arraylike object, optional
X0: arraylike or number, optional
input: int
output: int
**keywords:


Returns:  yout: array
T: array

Examples
>>> yout, T = step(sys, T, X0)
Create a transfer function system. Can create MIMO systems.
The function accepts either 1 or 2 parameters:
Create a transfer function system from its numerator and denominator polynomial coefficients.
If num and den are 1D array_like objects, the function creates a SISO system.
To create a MIMO system, num and den need to be 2D nested lists of array_like objects. (A 3 dimensional data structure in total.) (For details see note below.)
Parameters:  sys: Lti (StateSpace or TransferFunction)
num: array_like, or list of list of array_like
den: array_like, or list of list of array_like


Returns:  out: TransferFunction

Raises:  ValueError
TypeError

Notes
Todo
The next paragraph contradicts the comment in the example! Also “input” should come before “output” in the sentence:
“from the (j+1)st output to the (i+1)st input”
num[i][j] contains the polynomial coefficients of the numerator for the transfer function from the (j+1)st output to the (i+1)st input. den[i][j] works the same way.
The coefficients [2, 3, 4] denote the polynomial .
Examples
>>> # Create a MIMO transfer function object
>>> # The transfer function from the 2nd input to the 1st output is
>>> # (3s + 4) / (6s^2 + 5s + 4).
>>> num = [[[1., 2.], [3., 4.]], [[5., 6.], [7., 8.]]]
>>> den = [[[9., 8., 7.], [6., 5., 4.]], [[3., 2., 1.], [1., 2., 3.]]]
>>> sys1 = tf(num, den)
>>> # Convert a StateSpace to a TransferFunction object.
>>> sys_ss = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> sys2 = tf(sys1)
Transform a transfer function to a state space system.
The function accepts either 1 or 2 parameters:
Create a transfer function system from its numerator and denominator polynomial coefficients.
For details see: tf()
Parameters:  sys: Lti (StateSpace or TransferFunction)
num: array_like, or list of list of array_like
den: array_like, or list of list of array_like


Returns:  out: StateSpace

Raises:  ValueError
TypeError

Examples
>>> num = [[[1., 2.], [3., 4.]], [[5., 6.], [7., 8.]]]
>>> den = [[[9., 8., 7.], [6., 5., 4.]], [[3., 2., 1.], [1., 2., 3.]]]
>>> sys1 = tf2ss(num, den)
>>> sys_tf = tf(num, den)
>>> sys2 = tf2ss(sys_tf)
Return transfer function data objects for a system
Parameters:  sys: Lti (StateSpace, or TransferFunction)


Returns:  (num, den): numerator and denominator arrays

Compute system zeros.
Parameters:  sys: StateSpace or TransferFunction


Returns:  zeros: ndarray

Raises:  NotImplementedError

See also
Notes
This function is a wrapper for StateSpace.zero and TransferFunction.zero.
Todo
The following functions should be documented in their own modules! This is only a temporary solution.
Plot a pole/zero map for a linear system.
Parameters:  sys: Lti (StateSpace or TransferFunction)
Plot: bool


Returns:  pole: array
zeros: array

Nyquist plot for a system
Plots a Nyquist plot for the system over a (optional) frequency range.
Parameters:  syslist : list of Lti
omega : freq_range
Plot : boolean
labelFreq : int
*args, **kwargs:


Returns:  real : array
imag : array
freq : array

Examples
>>> sys = ss("1. 2; 3. 4", "5.; 7", "6. 8", "9.")
>>> real, imag, freq = nyquist_plot(sys)
Nichols plot for a system
Plots a Nichols plot for the system over a (optional) frequency range.
Parameters:  syslist : list of Lti, or Lti
omega : array_like
grid : boolean, optional


Returns:  None 
Place closed loop eigenvalues
Parameters:  A : 2d array
B : 2d array
p : 1d list


Returns:  K : 2d array

Examples
>>> A = [[1, 1], [0, 1]]
>>> B = [[0], [1]]
>>> K = place(A, B, [2, 5])
Linear quadratic regulator design
The lqr() function computes the optimal state feedback controller that minimizes the quadratic cost
The function can be called with either 3, 4, or 5 arguments:
Parameters:  A, B: 2d array
sys: Lti (StateSpace or TransferFunction)
Q, R: 2d array
N: 2d array, optional


Returns:  K: 2d array
S: 2d array
E: 1d array

Examples
>>> K, S, E = lqr(sys, Q, R, [N])
>>> K, S, E = lqr(A, B, Q, R, [N])
Controllabilty matrix
Parameters:  A, B: array_like or string


Returns:  C: matrix

Examples
>>> C = ctrb(A, B)
Observability matrix
Parameters:  A, C: array_like or string


Returns:  O: matrix

Examples
>>> O = obsv(A, C)
Gramian (controllability or observability)
Parameters:  sys: StateSpace
type: String


Returns:  gram: array

Raises:  ValueError
ImportError

Examples
>>> Wc = gram(sys,'c')
>>> Wo = gram(sys,'o')
Create a linear system that approximates a delay.
Return the numerator and denominator coefficients of the Pade approximation.
Parameters:  T : number
n : integer


Returns:  num, den : array

Notes
Based on an algorithm in Golub and van Loan, “Matrix Computation” 3rd. Ed. pp. 572574.
Plot the “Gang of 4” transfer functions for a system
Generates a 2x2 plot showing the “Gang of 4” sensitivity functions [T, PS; CS, S]
Parameters:  P, C : Lti
omega : array


Returns:  None 
Unwrap a phase angle to give a continuous curve
Parameters:  X : array_like
period : number


Returns:  Y : array_like

Examples
>>> import numpy as np
>>> X = [5.74, 5.97, 6.19, 0.13, 0.35, 0.57]
>>> unwrap(X, period=2 * np.pi)
[5.74, 5.97, 6.19, 6.413185307179586, 6.633185307179586, 6.8531853071795865]
X = lyap(A,Q) solves the continuoustime Lyapunov equation
A X + X A^T + Q = 0
where A and Q are square matrices of the same dimension. Further, Q must be symmetric.
X = lyap(A,Q,C) solves the Sylvester equation
A X + X Q + C = 0
where A and Q are square matrices.
X = lyap(A,Q,None,E) solves the generalized continuoustime Lyapunov equation
A X E^T + E X A^T + Q = 0
where Q is a symmetric matrix and A, Q and E are square matrices of the same dimension.
dlyap(A,Q) solves the discretetime Lyapunov equation
A X A^T  X + Q = 0
where A and Q are square matrices of the same dimension. Further Q must be symmetric.
dlyap(A,Q,C) solves the Sylvester equation
A X Q^T  X + C = 0
where A and Q are square matrices.
dlyap(A,Q,None,E) solves the generalized discretetime Lyapunov equation
A X A^T  E X E^T + Q = 0
where Q is a symmetric matrix and A, Q and E are square matrices of the same dimension.
(X,L,G) = care(A,B,Q) solves the continuoustime algebraic Riccati equation
A^T X + X A  X B B^T X + Q = 0
where A and Q are square matrices of the same dimension. Further, Q is a symmetric matrix. The function returns the solution X, the gain matrix G = B^T X and the closed loop eigenvalues L, i.e., the eigenvalues of A  B G.
(X,L,G) = care(A,B,Q,R,S,E) solves the generalized continuoustime algebraic Riccati equation
A^T X E + E^T X A  (E^T X B + S) R^1 (B^T X E + S^T) + Q = 0
where A, Q and E are square matrices of the same dimension. Further, Q and R are symmetric matrices. The function returns the solution X, the gain matrix G = R^1 (B^T X E + S^T) and the closed loop eigenvalues L, i.e., the eigenvalues of A  B G , E.
(X,L,G) = dare(A,B,Q,R) solves the discretetime algebraic Riccati equation
A^T X A  X  A^T X B (B^T X B + R)^1 B^T X A + Q = 0
where A and Q are square matrices of the same dimension. Further, Q is a symmetric matrix. The function returns the solution X, the gain matrix G = (B^T X B + R)^1 B^T X A and the closed loop eigenvalues L, i.e., the eigenvalues of A  B G.
(X,L,G) = dare(A,B,Q,R,S,E) solves the generalized discretetime algebraic Riccati equation
 A^T X A  E^T X E  (A^T X B + S) (B^T X B + R)^1 (B^T X A + S^T) +
 Q = 0
where A, Q and E are square matrices of the same dimension. Further, Q and R are symmetric matrices. The function returns the solution X, the gain matrix G = (B^T X B + R)^1 (B^T X A + S^T) and the closed loop eigenvalues L, i.e., the eigenvalues of A  B G , E.