home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 7.ddi / ROBUST.DI$ / AUGTF.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  5.1 KB  |  166 lines
  1. function [a,b1,b2,c1,c2,d11,d12,d21,d22] = augtf(Z1,Z2,Z3,Z4,Z5,Z6,Z7)
  2. % [TSS_] = AUGTF(SS_G,W1,W2,W3) or
  3. % [A,B1,B2,C1,C2,D11,D12,D21,D22] = AUGTF(AG,BG,CG,DG,W1,W2,W3) 
  4. % produces an augmented state-space realization for H2/Hinf control 
  5. % system design, where
  6. %         ss_g = mksys(ag,bg,cg,dg); (original plant state-space)
  7. %
  8. %         w1   = diagonal weighting matrix on "e"
  9. %              = [N11;D11;N22;D22;....]
  10. %
  11. %         w2   = diagonal weighting matrix on "u"
  12. %              = [N11;D11;N22;D22;....]
  13. %
  14. %         w3   = diagonal weighting matrix on "y", but allow polynomial 
  15. %         matrix at frequency infinity (i.e., w3(s) = sysw3 + w3poly(s))
  16. %              = [N11;D11;N22;D22;....]
  17. %
  18. %  The resulting state-space sysw1, sysw2, and W3(s) = sysw3 + w3poly(s) 
  19. %  are passed on to AUGSS.M to compute the final state-space.
  20. %
  21. %  If one of the weightings is missing, simply assign an empty bracket: [].
  22.  
  23. % R. Y. Chiang & M. G. Safonov 1/87
  24. % Copyright (c) 1988 by the MathWorks, Inc.
  25. % All Rights Reserved.
  26. % --------------------------------------------------------------------
  27.  
  28. inargs = '(ag,bg,cg,dg,w1,w2,w3)';
  29. eval(mkargs(inargs,nargin,'ss'))
  30.  
  31. if w1 == 0, w1 = []; end
  32. if w2 == 0, w2 = []; end
  33. if w3 == 0, w3 = []; end
  34. %
  35. % ------ Assemble the state-space quadruple of W1:
  36. %
  37. [mw1,nw1] = size(w1);
  38. %
  39. for i = 1:mw1/2
  40.     numw1 = w1(2*i-1,:);
  41.     denw1 = w1(2*i,:);
  42.     if norm(denw1,'fro') == 0 & norm(numw1,'fro') ~= 0
  43.        error('    W1 Denominator is zero !!')
  44.     end
  45.     numw1p = numw1; denw1p = denw1;
  46.     if norm(numw1,'fro') ~= 0
  47.        p = 1;                          % strip out leading zeros of 'numw1'
  48.        while numw1(1,p) == 0
  49.           numw1p = numw1(1,p+1:nw1);
  50.           p = p+1;
  51.        end
  52.        p = 1;                            
  53.        while denw1(1,p) == 0           % strip out leading zeros of 'denw1'
  54.           denw1p = denw1(1,p+1:nw1);
  55.           p = p+1;
  56.        end
  57.     else
  58.        numw1p = 0;
  59.        denw1p = 1;
  60.     end
  61.     if ((size(numw1p)-size(denw1p))*[0;1])>0
  62.         error('     Error: W1(s) is not proper !!')
  63.     end
  64.     [aw1p,bw1p,cw1p,dw1p] = tf2ss(numw1p,denw1p);  
  65.     if i == 1
  66.        aw1 = aw1p; bw1 = bw1p; cw1 = cw1p; dw1 = dw1p;
  67.     else
  68.        [aw1,bw1,cw1,dw1] = append(aw1,bw1,cw1,dw1,aw1p,bw1p,cw1p,dw1p);
  69.     end
  70. end
  71. %
  72. % ------ Assemble the state-space quadriple of W2:
  73. %
  74. [mw2,nw2] = size(w2);
  75. %
  76. for i = 1:mw2/2
  77.     numw2 = w2(2*i-1,:);
  78.     denw2 = w2(2*i,:);
  79.     if norm(denw2,'fro') == 0 & norm(numw2,'fro') ~= 0
  80.        error('    W2 Denominator is zero !!')
  81.     end
  82.     numw2p = numw2; denw2p = denw2;
  83.     if norm(numw2,'fro') ~= 0
  84.        p = 1;                          % strip out leading zeros of 'numw2'
  85.        while numw2(1,p) == 0
  86.           numw2p = numw2(1,p+1:nw2);
  87.           p = p+1;
  88.        end
  89.        p = 1;                            
  90.        while denw2(1,p) == 0           % strip out leading zeros of 'denw2'
  91.           denw2p = denw2(1,p+1:nw2);
  92.           p = p+1;
  93.        end
  94.     else
  95.        numw2p = 0;
  96.        denw2p = 1;
  97.     end
  98.     if ((size(numw2p)-size(denw2p))*[0;1])>0
  99.         error('     Error: W2(s) is not proper !!')
  100.     end
  101.     [aw2p,bw2p,cw2p,dw2p] = tf2ss(numw2p,denw2p);  
  102.     if i == 1
  103.        aw2 = aw2p; bw2 = bw2p; cw2 = cw2p; dw2 = dw2p;
  104.     else
  105.        [aw2,bw2,cw2,dw2] = append(aw2,bw2,cw2,dw2,aw2p,bw2p,cw2p,dw2p);
  106.     end
  107. end
  108. %
  109. % ------ Assemble the state-space quadriple of W3:
  110. %
  111. [mw3,nw3] = size(w3);
  112. n = mw3/2;
  113. N = nw3;
  114. w3poly = zeros(n,N*n);
  115. %
  116. for i = 1:n
  117.     numw3 = w3(2*i-1,:);
  118.     denw3 = w3(2*i,:);
  119.     if norm(denw3,'fro') == 0 & norm(numw3,'fro') ~= 0
  120.        error('    W3 Denominator is zero !!')
  121.     end
  122.     numw3p = numw3; denw3p = denw3;
  123.     if norm(numw3,'fro') ~= 0
  124.        p = 1;                            
  125.        while denw3(1,p) == 0           % strip out leading zeros of 'denw3'
  126.           denw3p = denw3(1,p+1:nw3);
  127.           p = p+1;
  128.        end
  129.     else
  130.        denw3p = 1;
  131.     end
  132.     [qw3,rw3] = deconv(numw3p,denw3p); % deconvolution of (numw3,denw3)
  133.     [mrw3,nrw3] = size(rw3);    
  134.     [mqw3,nqw3] = size(qw3);
  135.     if norm(rw3,'fro') < 1.e-14     % fix the numerical error from deconv
  136.        rw3 = zeros(mrw3,nrw3);
  137.     end                            % end of the fix
  138.     if nqw3 < N
  139.        qw3 = [zeros(1,N-nqw3) qw3];
  140.     end
  141.     aw3p = []; bw3p = []; cw3p = []; dw3p = 0;
  142.     if norm(rw3,'fro') ~= 0
  143.        p = 1;                   % strip out leading zeros of 'rw3'
  144.        while abs(rw3(1,p)) < 1.e-12
  145.              rw3p = rw3(1,p+1:nrw3);
  146.              p = p+1;
  147.        end
  148.        [aw3p,bw3p,cw3p,dw3p] = tf2ss(rw3p,denw3p);  % strictly proper tf 
  149.     end
  150.     [aw3,bw3,cw3,dw3] = append(aw3,bw3,cw3,dw3,aw3p,bw3p,cw3p,dw3p);
  151.     for j = 0:N-1                                % non-proper part
  152.         w3poly(i,j*n+i)=qw3(1,j+1);
  153.     end
  154. end
  155. %
  156. % ------ State-space augmentation with W(3) = sysw3 + w3poly:
  157. %
  158. [a,b1,b2,c1,c2,d11,d12,d21,d22] = augss(ag,bg,cg,dg,aw1,bw1,cw1,dw1,...
  159.         aw2,bw2,cw2,dw2,aw3,bw3,cw3,dw3,w3poly);
  160. %
  161. if xsflag
  162.    a = mksys(a,b1,b2,c1,c2,d11,d12,d21,d22,'tss');
  163. end   
  164. %
  165. % ------------- End of AUGTF.M  % RYC/MGS 10/19/1988 %
  166.