# [SciPy-user] Fwd: signal.lti(A,B,C,D) with D=0

flyeng4 williamhpurcell@gmail....
Tue Feb 17 14:25:54 CST 2009

```I am getting back to using signal.lti for state-space representation.
was no follow up. Any thoughts on a solution for representing MIMO or
SIMO systems with signal.lti, with the issue described below with
ss2zpk?
-Bill

---------- Forwarded message ----------
From: "Ryan Krauss" <ryanli...@gmail.com>
Date: Jul 14 2008, 7:41 am
Subject: signal.lti(A,B,C,D) with D=0
To: SciPy-user

This seems like a bug in ltisys.  I think this line
149     type_test = A[:,0] + B[:,0] + C[0,:] + D

should be
149     type_test = A[:,0] + B[:,0] + C[0,:] + D[0,:]

for a multipe-input/multiple-output system with i inputs, n states,
and m outputs, A should be n by n, B should be n by i, C should be m
by n, and D should be m by i.

(actually, because of this code a few lines above:
# make MOSI from possibly MOMI system.
if B.shape[-1] != 0:
B = B[:,input]
B.shape = (B.shape[0],1)
if D.shape[-1] != 0:
D = D[:,input]

I guess line 149 should be
149     type_test = A[:,0] + B[:,0] + C[0,:] + D[0]
)

But making this change exposes another problem:

/usr/lib/python2.5/site-packages/scipy/signal/ltisys.py in
__init__(self, *args, **kwords)
224                                 self.__dict__['D'] =
abcd_normalize(*args)
225             self.__dict__['zeros'], self.__dict__['poles'], \
--> 226                                     self.__dict__['gain'] =
ss2zpk(*args)
227             self.__dict__['num'], self.__dict__['den'] = ss2tf
(*args)
228             self.inputs = self.B.shape[-1]

/usr/lib/python2.5/site-packages/scipy/signal/ltisys.py in ss2zpk(A,
B, C, D, input)
183     """
184     Pdb().set_trace()
--> 185     return tf2zpk(*ss2tf(A,B,C,D,input=input))
186
187 class lti(object):

/usr/lib/python2.5/site-packages/scipy/signal/filter_design.py in
tf2zpk(b, a)
128     k = b[0]
129     b /= b[0]
--> 130     z = roots(b)
131     p = roots(a)
132     return z, p, k

/usr/lib/python2.5/site-packages/numpy/lib/polynomial.py in roots(p)
98     p = atleast_1d(p)
99     if len(p.shape) != 1:
--> 100         raise ValueError,"Input must be a rank-1 array."
101
102     # find non-zero array entries

<type 'exceptions.ValueError'>: Input must be a rank-1 array.
WARNING: Failure executing file: <proj2.py>

In [2]: %debug> /usr/lib/python2.5/site-packages/numpy/lib/
polynomial.py(100)roots()

99     if len(p.shape) != 1:
--> 100         raise ValueError,"Input must be a rank-1 array."
101

ipdb> print p
[[  1.00000000e+00   1.00000000e+00   1.00000000e+00   1.00000000e+00
1.00000000e+00   1.00000000e+00]
[  1.81898940e-11  -3.06222400e+00  -2.79196150e+02  -5.87518033e+03
-1.59843721e+04   1.91386789e-07]
[  1.75910240e+01   1.60479206e+03   3.37698682e+04   9.12618857e+04
-1.23219975e+04  -3.29828641e+04]]

ss2tf seems to correctly handle MIMO (or at least SIMO systems)
correctly and returns one denominator polynomial with several (m)
numerator polynomials.  But tf2zpk in filter_design.py does not seem
able to handle more than siso systems (which makes sense, it is
expecting a transfer function which is just a SISO system).

How should this be fixed?

I understand why ss2tf converts a MIMO system to SIMO - trying to
represent a mulitple input, multiple output system with a transfer
function has some limitations.  I think the real offender is in the
__init__ method of signal.lti:

elif N == 4:
...
self.__dict__['zeros'], self.__dict__['poles'], \
self.__dict__['gain'] = ss2zpk
(*args)

For a MIMO system, what should the zpk representation be?  For a SIMO
system, I would expect a vector of poles and a matrix of zeros.  This
seems to be in line with what ss2tf does.

It seems like the init method could find the pole by passing C[0,:]
and D[0,:] to ss2zpk.  The zeros and gains could then be found by
calling ss2zpk in some vectorized way or simply in a for loop for each
row of C and D.  But there may well be a more elegant solution.

Any thoughts?

Ryan

On Sun, Jul 13, 2008 at 2:59 PM, William Purcell

<williamhpurc...@gmail.com> wrote:
> I am trying to set up state space representation of a system using
> signal.lti. The system has no feedforward so D = 0. I've tried the three
> options listed in the code below to represent D.

> The only way I can get it to work is option 2 if C has 1 row. If C has more
> than 1 row it won't work.

> Any thoughts?
> -Bill

> Code
> ---------------------------------------------------------------
> from numpy import ones,matrix
> from scipy import signal
> r = ones(3)
> A = matrix([r,r,r])
> B = matrix([r]).T
> C = matrix([r,r])

> #three options of D to make it 0
> #1) D=0
> #2) D = matrix([0,0]).T
> #3) D = None

> D = 0
> #D = None
> #D = matrix([0,0]).T

> Gss = signal.lti(A,B,C,D)
> -----------------------------------------------------------

> Tracebacks
> -----------------------------------------------------------
> Option 1
> /usr/lib/python2.5/site-packages/scipy/signal/ltisys.py in abcd_normalize(A,
> B, C, D)
>     101         raise ValueError, "A and C must have the same number of
> columns."
>     102     if MD != MC:
> --> 103         raise ValueError, "C and D must have the same number of
> rows."
>     104     if ND != NB:
>     105         raise ValueError, "B and D must have the same number of
> columns."

> <type 'exceptions.ValueError'>: C and D must have the same number of rows.
> WARNING: Failure executing file: <testlti.py>

> Option 2 (with C as two rows...if C is a single row I do not get this
> traceback)
> /usr/lib/python2.5/site-packages/scipy/signal/ltisys.py in ss2tf(A, B, C, D,
> input)
>     147
>     148     num_states = A.shape[0]
> --> 149     type_test = A[:,0] + B[:,0] + C[0,:] + D
>     150     num = numpy.zeros((nout, num_states+1), type_test.dtype)
>     151     for k in range(nout):

> <type 'exceptions.ValueError'>: shape mismatch: objects cannot be broadcast
> to a single shape
> WARNING: Failure executing file: <testlti.py>

> Option 3
> same as 1

> _______________________________________________
> SciPy-user mailing list
> SciPy-u...@scipy.org
>http://projects.scipy.org/mailman/listinfo/scipy-user

_______________________________________________
SciPy-user mailing list
SciPy-u...@scipy.orghttp://projects.scipy.org/mailman/listinfo/scipy-
user
```