Stepzwn.m
上传用户:eighthdate
上传日期:2014-05-24
资源大小:270k
文件大小:3k
源码类别:

其他行业

开发平台:

Matlab

  1. function c = stepres(Z, wn, R, A, T, t)
  2. % Hadi Saadat,  1998
  3. clc
  4. discr=[
  5. '                                                                      '
  6. 'This program obtains the step response for the s-domain function given'
  7. 'below:  The graphical time response c(t) is plotted on the screen and '
  8. 'the time-domain specifications are obtained.                          '
  9. '                                                                      '
  10. '                    wn^2(1 + as)                                      '
  11. '   C(s) = -------------------------------- R(s) ,    where R(s) = R/s '
  12. '           (1 + Ts)(s^2 + 2 z wn s + wn^2)                            '];
  13. disp(discr)
  14. l=length(t);
  15. K=R;
  16. if Z < 1
  17.   Z1=1.-Z^2;
  18.   CN1=1-2.*Z*A*wn + A^2*wn^2;
  19.   CD11=1.-2.*T*Z*wn+T^2*wn^2;
  20.   CD1=Z1*CD11;
  21.   C1=sqrt(CN1/CD1);
  22.   CN2=wn^2*T*(A-T);
  23.   C2=CN2/CD11;
  24.   B=sqrt(Z1);
  25.   AFR=A*wn*B;
  26.   AFI=1.-A*Z*wn;
  27.   F1=atan2(AFR,AFI);
  28.   BFR=T*wn*B;
  29.   BFI=1.-T*Z*wn;
  30.   F2=atan2(BFR,BFI);
  31.   CFR=B;
  32.   CFI=-Z;
  33.   F3=atan2(CFR,CFI);
  34.   F= F1-F2-F3;
  35.   ALF=Z*wn ;
  36.     if T ~=  0
  37.     ALF2=1./T;
  38.     else
  39.     end
  40.   FT=3;
  41.   if T ~= 0
  42.    c=K*(1.+C1*exp(-ALF*t).*sin(wn*B*t+F)+C2*exp(-ALF2*t));
  43.    else
  44.    c=K*(1.+C1*exp(-ALF*t).*sin(wn*B*t+F));
  45.    end
  46. else
  47. end
  48. if Z > 1
  49.   bta=sqrt(Z^2 - 1);  p1 = Z*wn -wn*bta; p2 = Z*wn+wn*bta;
  50.   if T ~= 0
  51.    k1 = 1/(p1*p2); k2 = -(1-A/T)/( (-1/T+p1)*(-1/T+p2) );
  52.    k3 = -(1-A*p1)/( T*p1*(-p1+1/T)*(-p1+p2) );
  53.    k4 = -(1-A*p2)/( T*p2*(-p2+1/T)*(-p2+p1) );
  54.    c= K*wn^2*( k1 + k2*exp(-1/T*t) + k3*exp(-p1*t) + k4*exp(-p2*t) );
  55.    else
  56.    k1 = 1/(p1*p2);
  57.    k3 = -(1-A*p1)/( p1*(-p1+p2));
  58.    k4 = -(1-A*p2)/( p2*(-p2+p1)) ;
  59.    c= K*wn^2*( k1 +  k3*exp(-p1*t) + k4*exp(-p2*t) );
  60.    end
  61. else
  62. end
  63. if Z ==1
  64.    p=Z*wn;
  65.    if T ~= 0
  66.    k1 = 1/(p^2); k2 = -(1-A/T)/((-1/T+p)^2);
  67.    k3 = -(1-A*p)/( T*p*(-p+1/T) );
  68.    k4 = (-A*p*(-p+1/T)-(-2*p+1/T)*(-A*p+1))/( T*(-p)^2*(-p+1/T)^2 );
  69.    c= K*wn^2*( k1 + k2*exp(-1/T*t) + k3*t.*exp(-p*t) + k4*exp(-p*t) );
  70.    else
  71.    k1 = 1/(p^2); k3 =-(1-A*p)/p;
  72.    k4 =(-A*p - (1-A*p) )/( (-p)^2 );
  73.    c= K*wn^2*( k1 +  k3*t.*exp(-p*t) + k4*exp(-p*t) );
  74.    end
  75. end
  76. j=0 ;m =0;
  77. Cmax=max(c);
  78. if Cmax > K
  79.    overshoot= (Cmax - K)/K * 100;
  80.    for i=1:l
  81.        if c(i) == Cmax , tmax=t(i);
  82.        end
  83.        if A ==0
  84.           if T==0
  85.              sz=sqrt(1-Z^2);
  86.              overshoot = 100*exp(-pi*Z/sz);
  87.              tmax = pi/(wn*sz);
  88.              else, end
  89.         else,end
  90.    end
  91. fprintf('nn'),
  92. fprintf('Peak time = %g',tmax),fprintf('         Percent overshoot = %g',overshoot),
  93. %else
  94. else,end
  95.      for i =1:l
  96.           if  j == 0
  97.               if c(i) >= 0.1*K, t1=t(i); j=1;
  98.               end
  99.           else, end
  100.           if  m ==0
  101.               if c(i) >= 0.9*K t9 = t(i); m =1;
  102.               end
  103.           else, end
  104.       end
  105. end
  106. if t9 ~= 0
  107.   trise = t9 -t1;
  108. fprintf('n')
  109. fprintf('Rise time = %g',trise),fprintf('nn')
  110. end
  111. clg
  112. plot(t,c)
  113. xlabel('  Time    Sec.  '), ylabel('  c(t)  ')
  114. title('Step Response'),grid