home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 9.ddi / IDENT.DI$ / CONTIN.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.3 KB  |  72 lines
  1. function [bc,fc]=contin(th,ku,delays)
  2. %CONTIN   converts a model to continuous-time form. USE RATHER THD2THC!
  3. %
  4. %    [BC,FC] = contin(TH)
  5. %
  6. %    TH: A discrete time model defined in the format described by HELP THETA
  7. %
  8. %    BC and FC are returned as the numerator and denominator polynomials,
  9. %    respectively, of a continuous-time counterpart. Rows k contain
  10. %    the model associated with the k:th input.
  11. %
  12. %    With [BC,FC] = contin(TH,KU) the continuous-time models associated with
  13. %    the input numbers given by the entries of the row vector KU are
  14. %    computed (default is all the inputs). To obtain a continuous-time
  15. %    version of the noise model, let KU include the value zero. The noise 
  16. %    model will then be given in the corresponding row of [BC,FC].
  17. %
  18. %    When the model TH contains extra delays (nk>1) these are approximated 
  19. %    when forming the continuous time model. With a third, positive argument
  20. %    [BC,FC]=contin(TH,KU,1), the extra delays are removed before forming 
  21. %    the model. They should then be appended to the continuous model as a
  22. %    dead-time.
  23.  
  24. %    L. Ljung 10-1-86,4-20-87
  25. %    Copyright (c) 1986 by the MathWorks, Inc.
  26. %    All Rights Reserved.
  27.  
  28. T=th(1,2); nu=th(1,3);nk=th(1,7+2*nu:6+3*nu);
  29. nb=th(1,5:4+nu);nf=th(1,7+nu:6+2*nu);nc=th(1,5+nu);nd=th(1,6+nu);
  30. % *** Set up default values ***
  31.  
  32. if nargin<3,delays=0;end
  33. if nargin<2, ku=1:nu;end
  34. if length(ku)==1,if ku<0,ku=1:nu;end,end
  35. if delays<0,delays=0;end
  36.  
  37.  
  38. [a,b,c,d,f]=polyform(th);
  39.  
  40. if length(find(ku==0))>0
  41.      b(nu+1,1:length(c))=c;
  42.      f(nu+1,1:length(d))=d;
  43.      nk(nu+1)=0;
  44.      nb(nu+1)=nc+1;
  45.      nf(nu+1)=nd;    
  46. end
  47.  
  48. ss=1;
  49. for k=ku
  50.    if k==0, k=nu+1;end
  51.    den=conv(a,f(k,1:nf(k)+1)); 
  52.    if delays~=0,nkmod=max(nk(k),1);else nkmod=1;end %Corr 89-07-27
  53.    num=b(k,nkmod:nb(k)+nk(k));
  54.     nd=length(den);, nn=length(num);
  55.     n=max(nd,nn);
  56.     de=zeros(1,n);, de(1:nd)=den;
  57.     nume=zeros(1,n);, nume(1:nn)=num;
  58.     A=[-de(2:n);eye(n-2,n-1)]; % Transform to state-space
  59.     B=eye(n-1,1);
  60.     C=nume(2:n)-nume(1)*de(2:n);
  61.     D=nume(1);
  62.     s=real(logm([[A B];zeros(1,n-1) eye(1)]))/T;  % Sample
  63.     Ac=s(1:n-1,1:n-1);
  64.     Bcc=s(1:n-1,n);
  65.     ff=poly(Ac);                % Transform to i/o form
  66.     bb=poly(Ac-Bcc*C)+(D-1)*ff;
  67.     bc(ss,1:length(bb))=bb;
  68.     fc(ss,1:length(ff))=ff;
  69.     ss=ss+1;
  70. end
  71.  
  72.