home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 10.ddi / CONTROL.DI$ / FREQRESP.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  1.1 KB  |  48 lines
  1. function g=freqresp(a,b,c,d,iu,s)
  2. %FREQRESP Low level frequency response function.
  3. %
  4. %    G=FREQRESP(A,B,C,D,IU,S)
  5. %    G=FREQRESP(NUM,DEN,S)
  6.  
  7. %    Clay M. Thompson 7-10-90
  8. %    Copyright (c) 1986-93 by the MathWorks, Inc.
  9.  
  10. if (nargin==3)
  11.   % It is in transfer function form.  Do directly, using Horner's method
  12.   % of polynomial evaluation at the frequency points, for each row in
  13.   % the numerator.  Then divide by the denominator.
  14.   [ny,nx] = size(a);
  15.   s=c(:);
  16.   for i=1:ny
  17.     g(:,i) = polyval(a(i,:),s);
  18.   end
  19.   g = polyval(b,s)*ones(1,ny).\g;
  20. else
  21.   % It is in state space form.  Reduce to Hessenberg form then directly 
  22.   % evaluate frequency response.
  23.   [ny,nu] = size(d);
  24.   [nx,na] = size(a);
  25.   nw = max(size(s));
  26.  
  27.   % Balance A
  28.   [t,a] = balance(a);
  29.   b = t \ b;
  30.   c = c * t;
  31.   
  32.   % Reduce A to Hessenburg form
  33.   [p,a] = hess(a);
  34.  
  35.   % Apply similarity transformations from Hessenberg
  36.   % reduction to B and C:
  37.   if nx>0,
  38.     b = p' * b(:,iu);
  39.     c = c * p;
  40.     d = d(:,iu);
  41.     g = ltifr(a,b,s(:));
  42.     g = (c * g + diag(d) * ones(ny,nw)).';
  43.   else
  44.     d = d(:,iu);
  45.     g = (diag(d) * ones(ny,nw)).';
  46.   end
  47. end
  48.