function [sel, i] = see40(x, y, presel, iprev)

% SEE40 Browse and mark data graphically.
%  [SEL,I]=SEE40(X,Y,PRESEL,IPREV) invites
%   point selection and zooming on a plot of
%   Y(X).  The logical vector PRESEL marks
%   pre-selected points.  The result is the
%   logical vector SEL showing which points
%   remain selected when the routine exits.
%   PRESEL defaults to no preselections.
%   The vector IPREV is a list of indices
%   previously output as I by SEE3, that
%   are used to restore the previous state
%   of zoom.  IPREV defaults to no previous
%   zooming.
%  SEE40(Z,PRESEL) derives X and Y from the
%   first two columns of array Z.  PRESEL
%   behaves as above.
%  SEE40 (no arguments) performs a demonstration
%   of the routine with random data.

% Copyright (C) 1991 Charles R. Denham, Zydeco.
%	Updated by Marinna Martini 9/18/93
%	to work with MATLAB version 4.0

invisible = 'c0';
black = 'c1';
red = 'c2';
green = 'c3';
blue = 'c4';
cyan = 'c5';
magenta = 'c6';
yellow = 'c7';

if nargin == 0   % Demonstration.
   help see3
   m = 10;
   x = (1:m).';
   y = rand(m, 1);
   seltemp = see3(x, y);
   disp(' ');
   disp('   Index Selected');
   disp([x seltemp]);
   return
end

if nargin == 1
   y = x(:, 2); x = x(:, 1);
   presel = 0 .* x;
end

if nargin == 2
   [m, n] = size(x);
   if m > 1 & n > 1
      presel = y;
      y = x(:, 2); x = x(:, 1);
     else
      x = x(:); y = y(:);
      presel = 0 .* x;
   end
end

% Setup.

x = x(:); y = y(:); sel = presel(:);
[m, n] = size(x);
if nargin == 4
   i = iprev;
  else
   i = [1 m];
end

% Monotonic behaviour in the x-direction is
%  controlled by the latter of the next two
%  lines.  When "ismonotonic" is true (non-zero),
%  points are selected by their x-position only.
%  Otherwise, selection is based on the nearest
%  point, assuming the plot is square visually.

ismonotonic = 0;                   % Disabled.
ismonotonic = all(diff(x) >= 0);   % Enabled.

xsave = []; ysave = [];
newplot = 1;
inited = 0;

preselected = sum(presel ~= 0);
disp([' Number preselected: ' int2str(preselected)]);

while 1

% Refresh the plot.

if newplot
   hold off
   i(1:2) = sort(i(1:2));
   imin = i(1); imax = i(2);
   ind = imin:imax;
   xi = x(ind); yi = y(ind);
   if i(1) ~= i(2)
      plot(xi, yi, [black '-'])   % Many points.
      hold on
      if length(ind) < 50
         plot(xi, yi,[black '+'])
      end
      ii = (sel(ind) == 1);
      if any(ii), plot(xi(ii), yi(ii), [red 'o']), end
      hold off
     else
      i(1:2) = [];   % Never plot just one point.
   end
   title('[> Zoom <]            [> Unzoom <]')
   xlabel('[> Select <]'), ylabel('[> Return <]')
   newplot = 0;
   disp([' Range: ' int2str(i(1)) ' ' int2str(i(2))])
end

% Current axes.
% This section changed by M. Martini
%[status, ax] = inquire('axis');
%xmin = ax(1); xmax = ax(2);
ax = get(gca, 'Xlim');
xmin=ax(1); xmax=ax(2);
%ymin = ax(3); ymax = ax(4);
ax = get(gca, 'Ylim');
ymin=ax(1); ymax=ax(2);
xrange = xmax - xmin;
yrange = ymax - ymin;
xmean = mean([xmin xmax]);
ymean = mean([ymin ymax]);

% Invitation.

if ~inited
   disp(' Ready to click!'), inited = 1;
end

% Get the click.
%gxprev = gx;
%gyprev = gy;
if exist('gx') ~=1, 
	gxprev=0;
else,
	gxprev = gx;
end
if exist('gy') ~=1, 
	gyprev=0;
else,
	gyprev = gy;
end
[gx, gy] = ginput(1);
gxy = gx + sqrt(-1) .* gy;


% Locate the click.

isreturn = gx < xmin & gy > ymin & gy < ymax;
iszoom = gy > ymax & gx < xmean;
isunzoom = gy > ymax & gx > xmean;
isselect = gy < ymin & gx > xmin & gx < xmax;
isgraph = ~any([isreturn isselect iszoom isunzoom]);

% Return.

if isreturn
   selected = sum(sel ~= 0);
   changed = sum((presel ~= 0) ~= (sel ~= 0));
   disp([' Number preselected: ' int2str(preselected)]);
   disp([' Number selected:    ' int2str(selected)]);
   disp([' Number changed:     ' int2str(changed)]);
   return
end

% Zoom.

if iszoom & length(i) > 2
   if ~rem(length(i), 2)
      newplot = 1;
     else
      disp(' Choose another point.')
   end
end

% Unzoom.

if isunzoom
   if rem(length(i), 2)
      i(1) = [];
   end
   if length(i) > 2
      i(1:2) = [];
   end
   newplot = 1;
end

% Mark or unmark the selected point.
% the following three lines added by M. Martini
if exist('f') ~= 1,
	f=0;
end
if isselect & any(f)
   if i(1) == f, i(1) = []; end
   sel(f) = ~sel(f);
   hold on
   if sel(f)
      plot(x(f), y(f), [red 'o']);
      disp('        Selected')
     else
      plot(x(f), y(f), [cyan 'o']);
      disp('        Deselected')
   end
   hold off
end

% Find the nearest data point.

if isgraph
   if ismonotonic   % Use x-distance only.
      temp = abs(xi - gx);
     else           % Use Euclidean distance.
      xxi = xi .* yrange; yyi = yi .* xrange;
      xytemp = xxi + sqrt(-1) .* yyi;
      gxtemp = gx .* yrange; gytemp = gy .* xrange;
      gxytemp = gxtemp + sqrt(-1) .* gytemp;
      temp = abs(xytemp - gxytemp);
   end
   f = find(temp == min(temp));
   if any(f), f = f(1) + imin - 1; end
   i = [f i];
   hold on
   plot(x(f), y(f), [black 'o']);
   hold off
   disp([' Point: ' int2str(f)])
end

end