function argout = isolabel(cs, contours)
%       ISOLABEL Add labels to a contour plot. ISOLABEL is the Matlab
%       CLABEL function, slightly modified for the staplot package.
%
%	CLABEL(CS) adds height labels to the current contour plot
%	using the contour structure CS output from the CONTOUR routine.
%
%	The label positions are selected randomly.  For example:
%
%	   cs = contour(rand(10)); clabel(cs)
%
%	CLABEL(CS,V) labels just those contour levels given in
%	vector V.  The default action is to label all known contours.
%
%	See also ISOPYCS, CONTOUR, GINPUT.

%	Charles R. Denham, MathWorks, 1988, 1989.
%	Revised 4/12/90 by CRD: XYLIST now in screen coordinates.
%	Revised 4/16/90 by CRD: Removed an inappropriate test.
%	Copyright (c) 1984-94 by The MathWorks, Inc.
%       This version returns handles to all lines and text! (D.B. 11/6/95)
%

if nargin == 0
  error('Not enough input arguments.')
end
if min(size(cs)) > 2
  error('First input must be a valid contour description matrix.')
end
VERSION = version;
VERSION = str2num(VERSION(1));
cax = gca;
next = lower(get(cax,'NextPlot'));
manual = 0;
choice = 0;
if nargin > 1
  choice = 1;
  contours = sort(contours(:));
end

hold on;
[mcs, ncs] = size(cs);

% Find range of levels.
k = 1; i = 1;
while k <= ncs
   levels(i) = cs(1,k);
   i = i + 1;
   k = k + cs(2,k) + 1;
end
cmin = min(levels);
cmax = max(levels);
crange = max(abs(levels));
cdelta = abs(diff(levels)); 
cdelta = min(cdelta(cdelta > eps))/max(eps,crange); % Minimum significant change
if isempty(cdelta), cdelta = 0; end

% Decompose contour data structure if manual mode.

k = 0; n = 0; flip = 0;

ENT = 3;   % Ascii for "ENTER" key.
CR = 13;   % Ascii for "RETURN" key.
argout = [];
while (1)
  % Select a labeling point randomly if not manual.
  k = k + n + 1; if k > ncs, break, end
  c = cs(1, k); n = cs(2, k);
  if choice
    f = find(abs(c-contours)/max(eps+abs(contours)) < .00001);
    okay = length(f) > 0;
  else
    okay = 1;
  end
  if okay
    % put next line in in case rand('normal') is set - and to avoid
    % obsolescence message
    r = rand; while (r > 1 | r < 0) r = rand; end
    j = fix(r.* (n - 1)) + 1;
    if flip, j = n - j; end
    flip = ~flip;
    if n == 1 				% if there is only one point
      xx = cs(1, j+k); yy = cs(2, j+k);
    else
      x1 = cs(1, j+k); y1 = cs(2, j+k);
      x2 = cs(1, j+k+1); y2 = cs(2, j+k+1);
      xx = (x1 + x2) ./ 2; yy = (y1 + y2) ./ 2; % Test was here; removed.
    end
  end
  % Label the point.
  if okay
    % Set tiny labels to zero.
    if abs(c) <= 10*eps*crange, c = 0; end
    % Determine format string number of digits
    if cdelta > 0, 
      ndigits = max(3,ceil(-log10(cdelta)));
    else
      ndigits = 3;
    end
    s = [' ' sprintf(['%0.' num2str(ndigits) 'g'], c)];
    if VERSION < 5
      h = line('xdata',xx,'ydata',yy,'linesty','+','erasemode','none');
      argout = [argout; h];
    else
      h = line('xdata',xx,'ydata',yy,'marker','+', ...
	  'linestyle','none', 'erasemode','none');
      argout = [argout; h];
    end
    h = text(xx, yy, s, 'verticalalignment', 'bottom', ...
	'horizontalalignment', 'left','erasemode','none');
    argout = [argout; h];
  end
end
argout = argout(:);
return
