home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 9.ddi / IDENT.DI$ / ROTVAR.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  1.5 KB  |  55 lines
  1. function zv=rotvar(pol,covm)
  2. %ROTVAR    computes the standard deviations of roots to polynomials
  3. %
  4. %    ZV = rotvar(POL,COVM)
  5. %
  6. %    POL is the entered polynomial and
  7. %    COVM is the covariance matrix of its coefficients
  8. %    ZV is returned as the roots of POL and their standard deviations.
  9. %    The first column of ZV contains the roots and the corresponding
  10. %    elements in the second column give the standard deviations.
  11. %    For complex conjugate pairs the second entry is the correlation
  12. %    coefficient between the real and imaginary parts.
  13. %    
  14. %    The routine is primarily a help subroutine to th2zp.
  15.  
  16. %    L. Ljung 87-7-2
  17. %    Copyright (c) 1987 by the MathWorks, Inc.
  18. %    All Rights Reserved.
  19.  
  20. [dr,dc]=size(covm);lp=length(pol);
  21. r=roots(pol);
  22. zv(:,1)=r;lr=length(r);
  23. if isempty(covm),zv(:,2)=zeros(lr,1);return,end
  24. k=1;
  25. while k<lr+1
  26.     if imag(r(k))==0
  27.         D(k,:)=real(-poly(r([1:k-1,k+1:lr])));
  28.         k=k+1;
  29.     else
  30.     
  31.         D(k,:)=real(-2*poly([r([1:k-1,k+2:lr]);real(r(k))]));
  32.         D(k+1,:)=real([0, 2*imag(r(k))*poly(r([1:k-1,k+2:lr]))]);
  33.         k=k+2;
  34.     end
  35. end
  36. D=inv(D)';
  37. if lr~=dr,
  38.      DZ(:,2:lp)=D/pol(1);
  39.      DZ(:,1)=-D*pol(2:lp)'/(pol(1)^2);
  40. else
  41.      DZ=D;
  42. end
  43. PV=DZ*covm*DZ';
  44. k=1;l=1;i=sqrt(-1);
  45. while k<lp
  46.      %if r(k)==0,zv(k,2)=0;k=k+1;else 
  47.      if imag(r(k))==0,
  48.      zv(k,2)=sqrt(PV(l,l));k=k+1;l=l+1;
  49.           else zv(k,2)=sqrt(PV(l,l))+i*sqrt(PV(l+1,l+1));
  50.           zv(k+1,2)=PV(l,l+1)/sqrt(PV(l,l)*PV(l+1,l+1));
  51.       k=k+2;l=l+2;
  52.           end
  53.      end
  54. end
  55.