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

Ryan Krauss ryanlists@gmail....
Mon Jul 14 08:41:43 CDT 2008


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
<williamhpurcell@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-user@scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-user
>
>


More information about the SciPy-user mailing list