function bathyout = bathyset(L,meshele,meshnod,bathyin);
% bathyout = bathyset(L,meshele,meshnod,bathy);
disp('bathyset.m  Do not click before the beep!')
clf
sou = cos(2*pi*[1:200]/3);
x = meshnod(:,2);
y = meshnod(:,3);

if nargin~=4
bathy = 20*ones(size(x));
else
bathy = bathyin(:,2);
end

x3 = meshele(:,2:4);
x3(:)=x(meshele(:,2:4));
y3 = meshele(:,2:4);
y3(:)=y(meshele(:,2:4));
field3 = meshele(:,2:4);
field3(:)=bathy(meshele(:,2:4));

lplot = find(x3(:,1)>L(1) & x3(:,1)<L(2));
lplot = lplot(find(y3(lplot,1)>L(3) & y3(lplot,1)>L(3) ));

lplot= 1:length(x3);

hps = patch(x3(lplot,:)',y3(lplot,:)',field3(lplot,:)');
axis('equal')
axis(L)
title('Click on nodes to give depth')
maxdepth = 20;
colormap(jet);
map = colormap;
map = map(length(map):-1:1,:);
colormap(map)
caxis([0 maxdepth]);

Hcolorbar = colorbar('horiz');
ax = get(gcf,'children');
HSLIDERom = uicontrol('Style','slider','Max',maxdepth,'Min',00,...
'Position',[200 10 380 20]);
xlim = get(Hcolorbar,'XLim');
set(HSLIDERom,'Max',xlim(2))
set(HSLIDERom,'Min',xlim(1))
set(HSLIDERom,'units',get(Hcolorbar,'units'))
  pos = get(Hcolorbar,'position');
  posc = get(Hcolorbar,'position');
  posc(2) = posc(2)-posc(4);
  pos(2) = pos(2);
  pos(1) = pos(1)-pos(4)*.85;
  pos(3) = pos(3)+2*pos(4)*.85;
  set(HSLIDERom,'position',pos)
  set(Hcolorbar,'position',posc)
depth = get(HSLIDERom,'Value');
olddepth = depth;
drawnow
pbquit = uicontrol(gcf,'Style','radio','String','QUIT' ...
 ,'units','normalized','position',[.01 .95 .107 .047] );

pbredraw = uicontrol(gcf,'Style','radio','String','DRAW','Value',0 ...
 ,'units','normalized','position',[.01 .17 .12 .047] );



while(1)
sound(sou)
[xin, yin] = ginput(1);
drawnow
if get(pbquit,'value')
  break
end

if get(pbredraw,'value')==1
  disp('Redraw graph')
  axes(ax(1))
  set(gca,'nextplot','replace')
  delete(hps)
  hps=  patch(x3(lplot,:)',y3(lplot,:)',field3(lplot,:)');
  axis('equal')
  axis(L)
  ylabel('REDRAWN')
  title('Click on nodes to give depth')
  set(pbredraw,'value',0)
elseif inaxis(xin,yin,L)
%  New point selected
  d = (x-xin).^2+(y-yin).^2;
  nod1 = find(d==min(d));
  depth = get(HSLIDERom,'Value');
  bathy(nod1) = depth;
  disp(['Node ',int2str(nod1),' at ',num2str(x3(nod1)),' ',num2str(y3(nod1)),...
     '  depth=  ',num2str(depth)])
  c = (bathy(meshele(:,2))+bathy(meshele(:,3))+bathy(meshele(:,4)) )/3;
  field3 = meshele(:,2:4);
  field3(:)=bathy(meshele(:,2:4));
end

%disp(['End of whileloop  depth=  ',num2str(depth)])

end

bathyout = [[1:length(bathy)]' bathy];

delete(HSLIDERom)
delete(pbredraw)
delete(pbquit)
clf

disp('End of bathyset')

