home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 9.ddi / IDENT.DI$ / POLY2TH.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  3.3 KB  |  91 lines
  1. function th=poly2th(a,b,c,d,f,LAM,T)
  2. %POLY2TH constructs a "theta-matrix" from given polynomials
  3. %
  4. %    TH = poly2th(A,B,C,D,F,LAM,T)
  5. %
  6. %    TH: returned as a matrix of the standard format, describing the model
  7. %    A(q) y(t) = [B(q)/F(q)] u(t-nk) + [C(q)/D(q)] e(t)
  8. %    The exact format is given in HELP THETA.
  9. %
  10. %    A,B,C,D and F are entered as the polynomials. A,C,D and F start with
  11. %    1 , while B contains leading zeros to indicate the delay(s) for discr.
  12. %    time models. For multi-input systems B and F are matrices with the 
  13. %    number of rows equal to the number of inputs. For a time series, B and
  14. %    F are entered as []. POLY2TH is thus the inverse of TH2POLY.
  15. %
  16. %    LAM is the variance of the noise term e, and T is the sampling interval
  17. %    A negative value of T indicates a time-continuous model. Then the poly-
  18. %    nomials are entered in descending powers of s. Example: A=1,B=[1 2;0 3]
  19. %    C=1;D=1;F=[1 0;0 1]; T=-1 corresponds to the time-continuous system
  20. %    Y = (s+2)/s U1 + 3 U2.
  21. % 
  22. %    Trailing C,D,F, LAM, and T can be omitted, in which case they are
  23. %    taken as 1's (if B=[], then F=[]).
  24.  
  25. %    L. Ljung 10-1-86
  26. %       Revised 4-21-91
  27. %    Copyright (c) 1986-90 by the MathWorks, Inc.
  28. %    All Rights Reserved.
  29.  
  30. [nu,nb1]=size(b); nu=nu(1);
  31.  
  32. if nargin<7, T=1;end
  33. if nargin<6, LAM=1;end
  34. if nargin<5, f=ones(nu,1);end
  35. if nargin<4, d=1;end
  36. if nargin<3, c=1;end
  37. if a(1)~=1, error('The A-polynomial should be monic (start with "1")'),end 
  38. if T>0,for ku=1:nu,if f(ku,1)~=1,error('Each of the F-polynomials should be monic (start with "1") for discrete time models!'),end,end,end 
  39. na=length(a)-1;nc=length(c)-1;nd=length(d)-1;
  40. if nu>0,ib=b~=0; if nb1>1,nk=sum(cumsum(ib')==0);else nk=(cumsum(ib')==0);end
  41.     if T>0,    ib=(b(:,nb1:-1:1)~=0); 
  42.     if nb1>1,nb=-sum(cumsum(ib')==0)-nk+nb1;else nb=-(cumsum(ib')==0)-nk+nb1;end
  43.     nb=max(nb,zeros(1,nu));
  44.     else nb=nb1-nk;end,end
  45. if nu==0, nb=0;nk=0;end
  46. [nuf,nf1]=size(f);
  47. if nuf~=nu, error('f and b must have the same number of rows!'),return,end
  48. if nf1==1 | nf1==0
  49.    nf=zeros(1,nu);
  50. else
  51.    if T>0
  52.       ih=(f(:,nf1:-1:1)~=0); nf=-sum(cumsum(ih')==0)+nf1-1;
  53.    else 
  54.       ih=f~=0;
  55.       for ku=1:nu 
  56.         nf(ku)=nf1-min(find(ih(ku,:)~=0));
  57.         if f(ku,nf1-nf(ku))~=1
  58.          error('For continuous time systems the first non zero element of the rows of f must be a "1", corresponding to a monic F-polynomial')
  59.         end
  60.       end
  61.    end
  62. end
  63. n=na+sum(nb)+nc+nd+sum(nf);
  64. th=zeros(3,max([7 6+3*nu n]));
  65. if nu>0
  66.     th(1,1:6+3*nu)=[LAM T nu na nb nc nd nf nk];
  67. else
  68.     th(1,1:6)=[LAM T nu na nc nd];
  69. end
  70.  
  71. if na>0, th(3,1:na)=a(2:na+1);end
  72.  if nc>0, th(3,na+sum(nb)+1:na+nc+sum(nb))=c(2:nc+1);end
  73. if nd>0, th(3,na+nc+1+sum(nb):na+nc+nd+sum(nb))=d(2:nd+1);end
  74.  
  75. sb=na;sf=na+sum(nb)+nc+nd;
  76. for k=1:nu
  77.    if nb(k)>0,th(3,sb+1:sb+nb(k))=b(k,nk(k)+1:nk(k)+nb(k));end
  78.    if nf(k)>0,
  79.      if T>0,th(3,sf+1:sf+nf(k))=f(k,2:nf(k)+1);
  80.      else   th(3,sf+1:sf+nf(k))=f(k,nf1-nf(k)+1:nf1);
  81.      end
  82.    end
  83.    sb=sb+nb(k);
  84.    sf=sf+nf(k);
  85. end
  86. ti=fix(clock); ti(1)=ti(1)/100;
  87. th(2,2:6)=ti(1:5); th(2,7)=6;
  88. if T<0 & length(c)>length(conv(d,a)),disp('This model has differentiated white measurement noise. Expect problems!'),end
  89. if T<0 & length(c)<length(conv(d,a))
  90. disp('This model has no white  component in the measurement noise. Expect problems with sampling and state space representations!'),end
  91.