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

  1. function [a,b,c,d]=drmodel(n,p,m)
  2. %DRMODEL Generates random stable discrete nth order test models.
  3. %
  4. %    [NUM,DEN]=DRMODEL(N) forms an Nth order SISO transfer function 
  5. %        model.
  6. %
  7. %    [NUM,DEN]=DRMODEL(N,P) generates an Nth order single input, P 
  8. %        output transfer function model.
  9. %
  10. %    [A,B,C,D]=DRMODEL(N) generates an Nth order SISO state space model
  11. %
  12. %    [A,B,C,D]=DRMODEL(N,P,M) generates an Nth order, P output, M 
  13. %        input, state space model.
  14. %        
  15. %    See also: RMODEL.
  16.  
  17. % Clay M. Thompson  12-27-90
  18. % Copyright (c) 1986-93 by the MathWorks, Inc.
  19.  
  20. error(nargchk(0,3,nargin));
  21.  
  22. if nargin==0,
  23.   n=max([1,round(abs(10*randn(1,1)))]);
  24.   if nargout==4,
  25.     p=max([1,round(4*randn(1,1))]);
  26.     m=max([1,round(4*randn(1,1))]);
  27.   else
  28.     m=1;
  29.     p=1;
  30.   end
  31. end
  32. if nargin<2, p=1; end
  33. if nargin<3, m=1; end
  34.  
  35. rand('uniform')
  36. % Prob of an integrator is 0.10 for the first and 0.01 for all others
  37.  nint = (rand(1,1)<0.10)+sum(rand(n-1,1)<0.01);
  38. % Prob of repeated roots is 0.05
  39.  nrepeated = floor(sum(rand(n-nint,1)<0.05)/2);
  40. % Prob of complex roots is 0.5
  41.  ncomplex = floor(sum(rand(n-nint-2*nrepeated,1)<0.5)/2);
  42.  nreal = n-nint-2*nrepeated-2*ncomplex;
  43.  
  44. % Determine random poles
  45. rep = 2*rand(nrepeated,1)-1;
  46. mag = rand(ncomplex,1);
  47. ang = pi*rand(ncomplex,1);
  48. jay = sqrt(-1);
  49. complex = mag.*exp(jay*ang);
  50. re = real(complex); im = imag(complex);
  51.  
  52. if nargout==2, % Random transfer function
  53.   poles = [re+jay*im;re-jay*im;ones(nint,1);rep;rep;2*rand(nreal,1)-1];
  54.   
  55.   % Determine random zeros
  56.   zer = inf*ones(n,p);
  57.   jay = sqrt(-1);
  58.   for i=1:p
  59.     rand('uniform')
  60.     % Prob of complex zeros is 0.35
  61.     ncomplex = floor(sum(rand(n,1)<0.35)/2);
  62.     % Prob of real zero is 0.35
  63.     nreal = sum(rand(n-2*ncomplex,1)<0.35);
  64.     nzeros = 2*ncomplex+nreal;
  65.     if nzeros>0,
  66.       re = randn(ncomplex,1); im = randn(ncomplex,1);
  67.       zer(1:nzeros,i) = [re+jay*im;re-jay*im;randn(nreal,1)];
  68.     end
  69.   end
  70.   [a,b] = zp2tf(zer,poles,randn(p,1));
  71.  
  72. else % Random state space model
  73.   a = zeros(n,n);
  74.   for i=1:ncomplex,
  75.     ndx = [2*i-1,2*i];
  76.     a(ndx,ndx) = [re(i),im(i);-im(i),re(i)];
  77.   end
  78.   ndx = [2*ncomplex+1:n];
  79.   if ~isempty(ndx),
  80.      a(ndx,ndx) = diag([ones(nint,1);rep;rep;2*rand(nreal,1)-1]);
  81.   end
  82.   T = orth(rand(n,n));
  83.   a = T\a*T;
  84.   b = rand(n,m);
  85.   c = rand(p,n);
  86.   d = rand(p,m);
  87.   rand('uniform')
  88.   b = b .* (rand(n,m)<0.75);
  89.   c = c .* (rand(p,n)<0.75);
  90.   d = d .* (rand(p,m)<0.5);
  91. end
  92.