home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 10.ddi / CONTROL.DI$ / DSTEP.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.9 KB  |  106 lines

  1. function [yout,x,n] = dstep(a,b,c,d,iu,n)
  2. %DSTEP    Step response of discrete-time linear systems.
  3. %    DSTEP(A,B,C,D,IU)  plots the response of the discrete system:
  4. %
  5. %        x[n+1] = Ax[n] + Bu[n]
  6. %        y[n]   = Cx[n] + Du[n]
  7. %
  8. %    to a step applied to the single input IU.  The number of points is
  9. %    determined automatically.
  10. %
  11. %    DSTEP(NUM,DEN)  plots the step response of the polynomial transfer
  12. %    function  G(z) = NUM(z)/DEN(z)  where NUM and DEN contain the 
  13. %    polynomial coefficients in descending powers of z.
  14. %
  15. %    DSTEP(A,B,C,D,IU,N) or DSTEP(NUM,DEN,N) uses the user-supplied 
  16. %    number of points, N.  When invoked with left hand arguments,
  17. %        [Y,X] = DSTEP(A,B,C,D,...)
  18. %        [Y,X] = DSTEP(NUM,DEN,...)
  19. %    returns the output and state time history in the matrices Y and X.
  20. %    No plot is drawn on the screen.  Y has as many columns as there 
  21. %    are outputs and X has as many columns as there are states.
  22. %
  23. %    See also: DIMPULSE,DINITIAL,DLSIM and STEP.
  24.  
  25. %    J.N. Little 4-21-85
  26. %    Revised JNL 7-18-88, CMT 7-31-90, ACWG 6-21-92
  27. %    Copyright (c) 1986-93 by the MathWorks, Inc.
  28.  
  29. if nargin==0, eval('exresp(''dstep'');'), return,end
  30.  
  31. error(nargchk(2,6,nargin));
  32.  
  33. if (nargin==2),        % Transfer function without number of points
  34.     [num,den] = tfchk(a,b);
  35.     iu = 1;
  36.     [a,b,c,d] = tf2ss(num,den);
  37.  
  38. elseif (nargin==3),    % Transfer function with number of points
  39.     [num,den] = tfchk(a,b);
  40.     n = c;
  41.     iu = 1;
  42.     [a,b,c,d] = tf2ss(num,den);
  43.  
  44. elseif (nargin>=4)
  45.     error(abcdchk(a,b,c,d));
  46. end
  47.  
  48. [ny,nu] = size(d);
  49. if nu*ny==0, x = []; n = []; if nargout~=0, yout=[]; end, return, end
  50.  
  51. if nargin>4
  52.     if ~isempty(b), b=b(:,iu); end
  53.     d=d(:,iu);
  54. end
  55.  
  56. % Workout time vector if not supplied.
  57. if nargin==4 | nargin==5 | nargin==2
  58.     if isempty(a),
  59.         n = 11;
  60.     else
  61.         % The next line controls the number of samples in the plot if N not specified
  62.         st=0.005; % Set settling time bound  = 0.5%
  63.         % Step response is effectively equal to placing initial conditions
  64.         % on the plant as follows:
  65.         [ns,mu]=size(b);
  66.         x0=(eye(ns,ns)-a)\(b*ones(mu,1));
  67.         % Cater for pure integrator case
  68.         infind=find(~finite(x0));
  69.         x0(infind)=ones(length(infind),1);
  70.         n=dtimvec(a,b,c,x0,st);
  71.  
  72.         %  Multivariable systems
  73.         if nargin==4
  74.             [iu,nargin,y]=mulresp('dstep',a,b,c,d,n,nargout,0);
  75.             if ~iu, if nargout, yout = y; end, return, end
  76.         end
  77.     end
  78. end
  79.  
  80. [ny,nu] = size(d);
  81. if (nargin <= 3)&(nargout <= 1),    % transfer function description 
  82.     y = dlsim(num,den,ones(n,1));        % More efficient: uses FILTER
  83. else
  84.     [y,x] = dlsim(a,b,c,d,ones(n,nu));
  85. end
  86.  
  87. if nargout==0,        % If not output arguments, plot graph
  88.     status = ishold;
  89.     if rcond(a-eye(size(a))) > eps
  90.         dcgain=-c/(a-eye(size(a)))*b+d;
  91.     else
  92.         dcgain = 0;
  93.     end
  94.  
  95.     stairs([0:n-1],y)
  96.     hold on
  97.     plot([0,n-1],[dcgain';dcgain'],'w:')
  98.  
  99.     xlabel('No. of Samples')
  100.     ylabel('Amplitude')
  101.  
  102.     if ~status, hold off, end    % Return hold to previous status
  103.     return % Suppress output
  104. end
  105. yout = y; 
  106.