home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 9.ddi / IDENT.DI$ / TH2SSAUX.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  1.6 KB  |  62 lines
  1. function [a,b,c,d,k,x0]=th2ssaux(th);
  2. %TH2SSAUX    An auxiliary routine to TH2SS
  3. %
  4. %    [a,b,c,d,k,x0]=th2ssaux(th);
  5.  
  6. %    L. Ljung 10-2-90
  7. %    Copyright (c) 1990 by the MathWorks, Inc.
  8. %    All Rights Reserved.
  9.  
  10. T=th(1,2);
  11. [ap,bp,cp,dp,fp]=th2poly(th);
  12. nu=th(1,3);na=th(1,4);nb=th(1,5:4+nu);nc=th(1,5+nu);nd=th(1,6+nu);
  13. nf=th(1,7+nu:6+2*nu);nk=th(1,7+2*nu:6+3*nu);nf1=max(nd,max(nf))+1;
  14.  
  15. if T>0
  16.   fp(nu+1,1:length(dp))=dp;
  17. else
  18.   fp(nu+1,:)=zeros(1,nf1);fp(nu+1,nf1-nd:nf1)=dp;
  19. end
  20. nf(nu+1)=nd;
  21.  
  22. for ku=1:nu+1
  23.   if T>0 findx=1:nf(ku)+1;else findx=nf1-nf(ku):nf1;end
  24.   ffp=fp(ku,findx);
  25.   ap=conv(ap,ffp);
  26.   if ku==nu+1
  27.     btemp=cp;
  28.   else
  29.     if T>0,btemp=bp(ku,1:nb(ku)+nk(ku));
  30.     else   btemp=bp(ku,nk(ku)+1:nb(ku)+nk(ku));
  31.     end
  32.   end
  33.   for kku=1:nu+1
  34.     if T>0 findx=1:nf(kku)+1;else findx=nf1-nf(kku):nf1;end
  35.     if kku~=ku,btemp=conv(btemp,fp(kku,findx));end
  36.   end
  37.   bpp(ku,1:length(btemp))=btemp;nbpp(ku)=length(btemp);
  38. end
  39. [dum,maxindb]=size(bpp);nap=length(ap);nab=max(maxindb,nap);
  40. bp=zeros(nu+1,nab);
  41. if T<0,
  42. if maxindb>nap,error('Pure differentiations of inputs. No state space model produced!'),end
  43.  
  44. for ku=1:nu+1,bp(ku,:)=[zeros(1,nab-nbpp(ku)) bpp(ku,1:nbpp(ku))];end
  45.  
  46. else
  47.  
  48. bp=[bpp,zeros(nu+1,nab-maxindb)];end
  49. app=zeros(1,nab);app(1:nap)=ap;
  50. d=bp(1:nu,1)';
  51. if bp(nu+1,1)~=1,disp('WARNING: Noise model indicates no white noise component,Resulting state space model has unreliable noise structure'),end
  52. bp=bp-bp(:,1)*app;
  53. nx=nab-1;
  54. bt=bp(:,2:nab);
  55. at=app(2:nab)';
  56.  
  57. a=[-at,eye(nx,nx-1)]; % Transform to state-space
  58. b=bt(1:nu,:)';
  59. c=[1 zeros(1,nx-1)];
  60. k=bt(nu+1,:)';
  61. x0=zeros(nx,1);
  62.