home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 9.ddi / IDENT.DI$ / NYQPLOT.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.4 KB  |  89 lines
  1. function nyqplot(G,sd,mode)
  2. %NYQPLOT  plots the Nyquist diagram of a transfer function 
  3. %
  4. %    nyqplot(G)  or  nyqplot(G,SD)
  5. %
  6. %    G is entered as the transfer function(s) of the standard frequency 
  7. %    function format, produced by TH2FF, SPA or ETFE. 
  8. %    See HELP FREQFUNC.
  9. %    
  10. %    Several plots in the same diagram are obtained by
  11. %    nyqplot([G1 G2 .. Gn]) 
  12. %    The transfer functions Gk need not be given at the same frequencies.
  13. %    If SD is specified to a number larger than zero, also confidence
  14. %    intervals for the functions (corresponding to SD standard deviations)
  15. %    are shown as dash-dotted  curves.
  16. %
  17. %    The default mode is that ENTER will advance the plot from one
  18. %    input to the next (if any). To obtain all plots in the same diagram
  19. %    use nyqplot(G,SD,'same').
  20. %    Frequency functions for different outputs will always be given in
  21. %    separate plots.
  22.  
  23. %    L. Ljung 10-7-89
  24. %    Copyright (c) 1989-90 by the MathWorks, Inc.
  25. %    All Rights Reserved.
  26.  
  27. if (nargin==1), sd=0; mode='sep';end
  28. if (nargin==2),mode='sep';end
  29. clf
  30. mode=mode(1:3);
  31.  
  32. [m,n]=size(G);
  33. inputs1=G(1,:);mo=fix(inputs1/1000);mmo=max(mo);
  34. kko=[]; % The output indices
  35. for k=[0:mmo],if length(find(mo==k))>0,kko=[kko k];end,end
  36. G=G(2:m,:);
  37. for ko=kko+1
  38. inputs=inputs1-(ko-1)*1000;
  39. mi=max(inputs(find(inputs<19)));
  40.  
  41. kki=[]; % The input indices
  42. for k=[1:mi 0], if length(find(inputs==k))>0,kki=[kki k];end,end
  43.  
  44. if mode=='sep'
  45. for k=kki
  46. if k~=0,      
  47. clg
  48.  
  49.     gp=G(:,find(inputs==20+k))*pi/180;
  50.     polar(gp,G(:,find(inputs==k)))
  51.     ax=axis; hold on, polar([0 0],ax(1:2)),
  52.     polar([pi pi]/2,ax(1:2)),hold off
  53.     if sd>0,hold on,
  54.     ind=find(k+50==inputs);
  55.     polar(gp,G(:,ind-1)+sd*G(:,ind),'-.')
  56.     polar(gp, G(:,ind-1)-sd*G(:,ind),'-.')
  57.     hold off,end
  58.  
  59.     
  60.  
  61.     title(['NYQUIST PLOT, input # ',int2str(k),' output # ',int2str(ko)])
  62.     
  63. else disp('No Nyquist plot for spectra!'),end
  64.  
  65.  
  66. if k~=kki(length(kki)),pause,end
  67. end
  68. end
  69.  
  70. if mode=='sam' 
  71.         gp=G(:,find(inputs<39 & inputs > 19))*pi/180;
  72.           polar(gp,G(:,find(inputs<20)))
  73.     ax=axis; hold on, polar([0 0],ax(1:2)),
  74.     polar([pi pi]/2,ax(1:2)),hold off
  75.     if sd>0,hold on,
  76.         ind=find(inputs<69 & inputs>49);
  77.     polar(gp,G(:,ind-1)+sd*G(:,ind),'-.')
  78.     polar(gp, G(:,ind-1)-sd*G(:,ind),'-.')
  79.     
  80.     hold off,end
  81.     
  82.     title(['NYQUIST PLOT output # ',int2str(ko)])
  83.       end
  84. if ko~=kko(length(kko))+1,pause,end
  85. end
  86. hold off
  87. set(gcf,'NextPlot','replace');
  88.  
  89.