function RespObj = axeslims(RespObj)
%AXESLIMS  Sets the axes limits for Nyquist Response Objects.
%
%   RespObj = AXESLIMS(RespObj)

%   Authors: Karen Gondoly, 5-14-99
%   Revised: P. Gahinet, A. DiVergilio
%   Copyright 1986-2000 The MathWorks, Inc. 
%   $Revision: 1.17 $  $Date: 2000/06/14 20:09:11 $

AllRespProps = get(RespObj,'Undocumented');
DispAx = AllRespProps.PlotAxes;  
ResponseHandles = AllRespProps.ResponseHandles;
Xlims = AllRespProps.Xlim;
Ylims = AllRespProps.Ylim;
ZoomMode = AllRespProps.Zoom;

% Find axis lines and hide them (may interfer with axis limits)
% RE: Watch for axis grouping in MIMO case (handles get reparented)
SHH = get(0,'ShowHiddenHandles');
set(0,'ShowHiddenHandles','on');
XaxisLines = findobj(DispAx(:),'Tag','PZxaxisLine');
YaxisLines = findobj(DispAx(:),'Tag','PZyaxisLine');
set(0,'ShowHiddenHandles',SHH)
set(XaxisLines,'Xdata',[NaN NaN])
set(YaxisLines,'Ydata',[NaN NaN])

switch ZoomMode
case {'tight','critical'}
   % Compute desirable axis limits
   [Xlims,Ylims] = nylims(RespObj,ZoomMode);
   
otherwise
   % Default axis limit management
   AutoX = strcmp(AllRespProps.XlimMode,'auto');
   AutoY = strcmp(AllRespProps.YlimMode,'auto');
   %---Set axis lim mode
   if AutoX,
      set(DispAx(:),'XlimMode','auto')
   end
   if AutoY,
      set(DispAx(:),'YlimMode','auto')
   end  
   %---Adjust X and Y limits
   if AutoX,
      %---Enforce uniform X-axis limits across each column
      for ctax=1:size(DispAx,2),
         Xlims{ctax} = findlims(DispAx(:,ctax),'Xlim');
      end 
   end 
   if AutoY,
      %---Enforce uniform Y-axis limits across each row
      for ctax=1:size(DispAx,1),
         Ylims{ctax} = findlims(DispAx(ctax,:),'Ylim');
         if AutoY,
            % Enforce symmetry
            Ylims{ctax} = max(abs(Ylims{ctax})) * [-1 1];
         end
      end 
   end 
end

% Protect against axis grouping
switch AllRespProps.AxesGrouping
case 'all'  % limits of (1,1) axis apply to all plots
   XlimsG = Xlims(1);
   YlimsG = Ylims(1);
case 'inputs'  % (i,1) axis limits apply to all plots in row i
   XlimsG = Xlims(1);
   YlimsG = Ylims;
case 'outputs' % (1,j) axis limits apply to all plots in col j
   XlimsG = Xlims;
   YlimsG = Ylims(1);
otherwise
   XlimsG = Xlims;
   YlimsG = Ylims;
end

% (Re)draw Nyquist arrows
iosizes = size(DispAx);
for ksys=1:length(ResponseHandles),
    for kax=1:prod(iosizes),
        [i,j] = ind2sub(iosizes,kax);
        Xlim = XlimsG{min(j,end)};
        Ylim = YlimsG{min(i,end)};
        % All lines in (i,j) axis for system #ksys
        LineHandles = ResponseHandles{ksys}{kax}; % all lines in (i,j) axis
        for m=1:prod(size(LineHandles))
            % Reposition arrows
            LocalPositionArrow(LineHandles{m},Xlim,Ylim,1);
            LocalPositionArrow(LineHandles{m},Xlim,Ylim,-1);
        end
    end
end

%---Toggle grid
RespObj = plotgrid(RespObj);

% Update limits and ResponseHandles property
% RE: Leaves XlimMode and YlimMode untouched
RespObj = respfcn('setxylims',RespObj,Xlims,Ylims,ResponseHandles);

%---Adjust axis line extent
for ctAx=1:length(XaxisLines),
   ax = get(XaxisLines(ctAx),'Parent');
   set(XaxisLines(ctAx),'Xdata',get(ax,'Xlim'))
   set(YaxisLines(ctAx),'Ydata',get(ax,'Ylim'))
end   


%=============================
% subfunction GetArrowPosition
%=============================

function LocalPositionArrow(PlotHandles,Xlim,Ylim,Sign)
% RE: Assumes Ylim(1)=-Ylim(2)

% Find line handle
hline = findobj(PlotHandles,'Tag','LTIresponseLines');
ud = get(hline,{'UserData'});
hline = hline(1 + (ud{2}.Sign == Sign));

% Find arrow handle
harrow = findobj(PlotHandles,'Tag','NyquistArrow','UserData',Sign);
            
% Get line data      
data = get(hline,{'xdata','ydata'});
[xdata,ydata] = deal(data{:});
npts = length(xdata);

% Geometric parameters
xrange = Xlim(2)-Xlim(1);
yrange = Ylim(2)-Ylim(1);

% Shrink range to leave room for arrows
dX = xrange/40;  
dY = yrange/40;
Xlim = Xlim + [dX -dX];
Ylim = Ylim + [dY -dY];

% Select segment [A,B] where to plot the arrow 
% RE: A is always inside current axis limits
inscope = (xdata>Xlim(1) & xdata<Xlim(2) & ydata>Ylim(1) & ydata<Ylim(2));
[Ymax,iy] = max(abs(ydata));
if ~any(inscope),
   % No data in scope
   set(harrow,'Xdata',[NaN,NaN,NaN],'Ydata',[NaN,NaN,NaN])
   return
elseif all(inscope) | (inscope(iy) & Ymax>max(abs(Ylim))/100),
   % All data is in scope, or peak Y value is nonzero and in scope
   if ~any(ydata),
      % All data is real: pick mid point in x data range
      [junk,ia] = min(abs(xdata-(min(xdata)+max(xdata))/2));
   else
      % Pick point with largest |y| value
      ia = iy;
   end
   ib = ia-1+2*(ia==1);
elseif inscope(1),
   % Find first data point out of scope
   outs = find(~inscope);
   ib = outs(1);
   ia = ib-1;
else
   % Find first data point in scope   
   ins = find(inscope);
   ia = ins(1);
   ib = ia-1;
end
xa = xdata(ia);   
ya = ydata(ia);
xb = xdata(ib);  
yb = ydata(ib);

% Adjust position to place arrow inside and near axis limits
% (Find smallest t>0 s.t. A+t(B-A) intersect box [Xmin+dX Xmax-dX Ymin+dY Ymax-dY]
if xa==xb & ya==yb,
   % No arrow if (xa,ya)=(xb,yb) (e.g., sys = tf(1) )
   set(harrow,'Xdata',[NaN,NaN,NaN],'Ydata',[NaN,NaN,NaN])
   return
elseif prod(Xlim-xb)>=0 | prod(Ylim-yb)>=0,
   % (XB,YB) is outside of axis limits
   t = [];
   if xa~=xb,
      t = [t , (Xlim-xa)/(xb-xa)];
   end
   if ya~=yb, 
      t = [t , (Ylim-ya)/(yb-ya)];
   end
   t = min(t(t>0));   % smallest t>0 s.t. A+t(B-A) intersect box 
   if length(t),
      xa = xa + t * (xb-xa);
      ya = ya + t * (yb-ya);
   end
end

% Compute arrow extent
if all(inscope)
   sf = max((max(xdata)-min(xdata))/xrange,(max(ydata)-min(ydata))/yrange);
else
   sf = 0.9;
end
xsf = sf * dX;  % X arrow size
ysf = sf * dY;  % Y arrow size

% Determine (normalized) arrow direction DIR = B-A
dir = atan2((yb-ya)/yrange,(xb-xa)/xrange) + (Sign*(ib-ia)<0) * pi;

% Compute arrow position (tip at A, direction DIR)
Xarrow = xa + xsf * [cos(dir+5*pi/6);0;cos(dir-5*pi/6)];
Yarrow = ya + ysf * [sin(dir+5*pi/6);0;sin(dir-5*pi/6)];

% Compute offset to prevent arrow overlap when YA=0
offset = xsf/2 * (abs(ya)<ysf); 
if (xa-sum(Xlim)/2) * (xa-Xarrow(1))>0
    Xarrow = Xarrow + sign(Xarrow(1)-xa) * offset;
end

% Set position
set(harrow,'Xdata',Xarrow,'Ydata',Yarrow)
