home *** CD-ROM | disk | FTP | other *** search
/ OS/2 Shareware BBS: 10 Tools / 10-Tools.zip / octa21eb.zip / octave / SCRIPTS.ZIP / scripts / control / zp2ss.m < prev    next >
Text File  |  1999-03-05  |  4KB  |  144 lines
  1. # Copyright (C) 1996 A. Scottedward Hodel 
  2. #
  3. # This file is part of Octave. 
  4. #
  5. # Octave is free software; you can redistribute it and/or modify it 
  6. # under the terms of the GNU General Public License as published by the 
  7. # Free Software Foundation; either version 2, or (at your option) any 
  8. # later version. 
  9. # 
  10. # Octave is distributed in the hope that it will be useful, but WITHOUT 
  11. # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 
  12. # FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License 
  13. # for more details.
  14. # 
  15. # You should have received a copy of the GNU General Public License 
  16. # along with Octave; see the file COPYING.  If not, write to the Free 
  17. # Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. 
  18.  
  19. function [a,b,c,d] = zp2ss(zer,pol,k)
  20. # [A,B,C,D] = zp2ss(zer,pol,k)
  21. # Conversion from zero / pole to state space.
  22. # Inputs: 
  23. #   zer,  pol: vectors of (possibly) complex poles and zeros of a transfer
  24. #              function.  Complex values must come in conjugate pairs
  25. #              (i.e., x+jy in zer means that x-jy is also in zer)
  26. #   k:  real scalar (leading coefficient)
  27. # Outputs:
  28. #  A, B, C, D:
  29. # The state space system
  30. #      .
  31. #      x = Ax + Bu
  32. #      y = Cx + Du
  33. #
  34. # is obtained from a vector of zeros and a vector of poles via the
  35. # function call [a,b,c,d] = zp2ss(zer,pol,k).  The vectors 'zer' and 
  36. # 'pol' may either be row or column vectors.  Each zero and pole that
  37. # has an imaginary part must have a conjugate in the list.
  38. # The number of poles must at least equal the number of zeros.
  39. # k is a gain that is associated with the zero vector.
  40.  
  41. # Written by David Clem August 15, 1994
  42.  
  43.   sav_val = empty_list_elements_ok;
  44.   empty_list_elements_ok = 1;
  45.  
  46.   if(nargin != 3)
  47.     error("Incorrect number of input arguments");
  48.   endif
  49.  
  50.   if(! (is_vector(zer) | isempty(zer)) )
  51.     error(["zer(",num2str(rows(zer)),",",num2str(columns(zer)), ...
  52.     ") should be a vector"]);
  53.   elseif(! (is_vector(pol) | isempty(pol) ) )
  54.     error(["pol(",num2str(rows(pol)),",",num2str(columns(pol)), ...
  55.     ") should be a vector"]);
  56.   elseif(! is_scalar(k))
  57.     error(["k(",num2str(rows(k)),",",num2str(columns(k)), ...
  58.     ") should be a scalar"]);
  59.   elseif( k != real(k))
  60.     warning("zp2ss: k is complex")
  61.   endif
  62.  
  63.   zpsys = ss2sys([],[],[],k);
  64.  
  65.   # Find the number of zeros and the number of poles
  66.   nzer=length(zer);
  67.   npol =length(pol);
  68.  
  69.   if(nzer > npol)
  70.     error([num2str(nzer)," zeros, exceeds number of poles=",num2str(npol)]);
  71.   endif
  72.  
  73.   # Sort to place complex conjugate pairs together
  74.   zer=sortcom(zer);
  75.   pol=sortcom(pol);
  76.  
  77.   # construct the system as a series connection of poles and zeros
  78.   # problem: poles and zeros may come in conjugate pairs, and not
  79.   # matched up!
  80.  
  81.   # approach: remove poles/zeros from the list as they are included in
  82.   # the ss system
  83.  
  84.   while(length(pol))
  85.  
  86.     # search for complex poles, zeros
  87.     cpol=[];    czer = [];
  88.     if(!isempty(pol))
  89.       cpol = find(imag(pol) != 0);
  90.     endif
  91.     if(!isempty(zer))
  92.       czer = find(imag(zer) != 0);
  93.     endif
  94.  
  95.     if(isempty(cpol) & isempty(czer))
  96.       pcnt = 1;
  97.     else 
  98.       pcnt = 2;
  99.     endif
  100.  
  101.     num=1;    # assume no zeros left.
  102.     switch(pcnt)
  103.     case(1)
  104.       # real pole/zero combination
  105.       if(length(zer))
  106.         num = [1 -zer(1)];  
  107.         zer = zer(2:length(zer));
  108.       endif
  109.       den = [1 -pol(1)];
  110.       pol = pol(2:length(pol));
  111.     case(2)
  112.       # got a complex pole or zero, need two roots (if available)
  113.       if(length(zer) > 1)
  114.         [num,zer] = zp2ssg2(zer);    # get two zeros
  115.       elseif(length(zer) == 1)
  116.         num = [1 -zer];            # use last zero (better be real!)
  117.         zer = [];
  118.       endif
  119.       [den,pol] = zp2ssg2(pol);        # get two poles
  120.     otherwise
  121.       error(["pcnt = ",num2str(pcnt)])
  122.     endswitch
  123.  
  124.     # pack tf into system form and put in series with earlier realization
  125.     zpsys1 = tf2sys(num,den,0,"u","yy");
  126.  
  127.     # change names to avoid warning messages from sysgroup
  128.     zpsys  = syssetsignals(zpsys,"in","u1",1);
  129.     zpsys1 = sysupdate(zpsys1,"ss");
  130.     nn     = sysdimensions(zpsys);        # working with continuous system
  131.     zpsys  = syssetsignals(zpsys,"st", sysdefioname(nn,"x"));
  132.     nn1    = sysdimensions(zpsys1);
  133.     zpsys1 = syssetsignals(zpsys1,"st",sysdefioname(nn1,"xx"));
  134.  
  135.     zpsys = sysmult(zpsys,zpsys1);
  136.  
  137.   endwhile 
  138.  
  139.   [a,b,c,d] = sys2ss(zpsys);
  140.  
  141.   empty_list_elements_ok = sav_val;
  142. endfunction
  143.  
  144.