function COAST(axes,maxlat);
% COAST makes map of the world.
%
%        COAST([LONG_MIN LONG_MAX LAT_MIN LAT_MAX]) draws a
%        coastline map of the world at about 1/2 degree
%        resolution. COAST (called without arguments) draws
%        the whole world with the Pacific in the center.
% 
%        COAST('north',MAXLAT) draws a polar projection from the North
%        pole out to latitude MAXLAT. MAP('south') does the same for 
%        the South pole (MAXLAT is optional in both cases, but is taken
%        to lie in the same hemisphere as the pole).
%        

%Notes: RP (WHOI) 6/Dec/91
%                 7/Nov/92  Changed for matlab4.0
%                 17/Oct/93 New outline file.
%                 14/Mar/94 Fixed polar projections, and changed name
%                           to 'worldmap'.
%              

% change this line (if map.m is not in your default directory) to
% "load <path>/coasts".
load coastln
global lon_var lat_var proplist

% make proplist vars global
s = size(proplist); str = []; for i = 1:s(1); str = [str ' ']; end
str = sprintf('%s',[proplist str']'); eval(['global ' str]); 

VERSION = version; VERSION = str2num(VERSION(1));
  
shiftx=0;
drawtwice=0;
proj='rec';
if (nargin<2), maxlat=40; else maxlat=abs(maxlat); end;

if (maxlat>90); error('MAXLAT greater than 90!'); end;

if (nargin>0),
   if (isstr(axes)),
      if (axes(1:3)=='nor'),
         proj='npl';
         axes=(90-maxlat)*[-1 1 -1 1];
      elseif (axes(1:3)=='sou'),
         proj='spl';
         axes=(90-maxlat)*[-1 1 -1 1];
      else
         error('map: Unrecognized projection!');
      end;
   else         
      axes=axes(:)';
      if (max(size(axes)) ~= 4), 
         error('map: wrong number of limit args!');
      end;
      if ( (axes(4)-axes(3))>180. ),
         error('map: Lat range greater than 180 degrees');
      end;
      if ( (axes(2)-axes(1))>360. ),
         error('map: Long range greater than 360 degrees');
      else
         if (axes(1)<=-180. ),
            shiftx=-360;
            drawtwice=1;
         elseif (axes(2)>180.), 
            shiftx=360;
            drawtwice=1;
         else
            shiftx=0;
            drawtwice=0;
         end;
      end;
   end;
else
   axes=[-334 25 -90 90];
   shiftx=-360;
   drawtwice=1;
end;

if (proj=='rec'),
   plot(coastlines(1,:)+shiftx,coastlines(2,:),'w-',...
       'Erasemode','none')
elseif (proj=='npl'),
   kk=find(coastlines(2,:)<maxlat);
   xx=(90-coastlines(2,:)).*cos(coastlines(1,:)*pi/180-pi/2);
   yy=(90-coastlines(2,:)).*sin(coastlines(1,:)*pi/180-pi/2);
   xx(kk)=NaN*kk;
   plot(xx,yy,'w-','Erasemode','none')
   set(gca,'aspect',[1 1]);

elseif (proj=='spl'),
   kk=find(coastlines(2,:)>-maxlat);
   xx=(90+coastlines(2,:)).*cos(-coastlines(1,:)*pi/180+pi/2);
   yy=(90+coastlines(2,:)).*sin(-coastlines(1,:)*pi/180+pi/2);
   xx(kk)=NaN*kk;
   plot(xx,yy,'w-','Erasemode','none')
   set(gca,'aspect',[1 1]);
end;

if (drawtwice),
   line(coastlines(1,:),coastlines(2,:),'Linestyle','-',...
   'Color','w','Erasemode','none');
end;


if (proj=='rec'),
  axis(axes);
  xl = get(gca,'xtick'); [Nn,Mm] = size(xl);
  zz=rem(xl+540,360)-180; zz = zz(:);
  for kk=1:Mm,
    ixl = num2str(abs(zz(kk)));
    s = size(ixl); s = s(2);
    for ii=1:3-s
      ixl = [' ' ixl];
    end
    if zz(kk) < 0, 
      nxl(kk,:) = [ixl 'W'];
    elseif zz(kk) > 0  
      nxl(kk,:) = [ixl 'E'];
    else
      nxl(kk,:) = [ixl ' '];
    end;
  end;
  if VERSION < 5
    set(gca,'Xticklabels',nxl)
  else
    set(gca,'Xticklabel',nxl)
  end
  yl=get(gca,'Ytick');
 [Mm,Nn]=size(yl);
 for kk=1:Nn,
    iyl = num2str(abs(yl(kk)));
    s = size(iyl); s = s(2);
    for ii=1:3-s
      iyl = [' ' iyl];
    end
   if (yl(kk) < 0)
     nyl(kk,:) = [iyl 'S'];
   elseif (yl(kk) > 0)
     nyl(kk,:) = [iyl 'N'];
   else
     nyl(kk,:) = [iyl ' '];
   end;
 end;
 if VERSION < 5
   set(gca,'Yticklabels',nyl)
 else
   set(gca,'Yticklabel',nyl)
 end
 set(gca,'Ytick',yl);  % Need this to reset position limits

else % polar projection
  axis(axes);
  set(gca,'visible','off');
  for kk=10:10:(90-maxlat),
    xx=kk*cos([0:10:360]*pi/180);
    yy=kk*sin([0:10:360]*pi/180);
    line(xx,yy,'color','w','linestyle',':');
  end;
  xx=(90-maxlat)*cos([0:5:360]*pi/180);
  yy=(90-maxlat)*sin([0:5:360]*pi/180);
  line(xx,yy,'color','w','linestyle','-');
  
  % set tick labels
  for kk=10:10:(90-maxlat-2),
    if (proj=='npl'), 
       text(0,kk,sprintf('%2.0fN',90-kk),'horizontal','center');
    else
       text(0,kk,sprintf('%2.0fS',90-kk),'horizontal','center');
    end;
  end;
  for kk=0:30:359,
    xx=[10 90-maxlat]*cos(kk*pi/180);
    yy=[10 90-maxlat]*sin(kk*pi/180);
    line(xx,yy,'color','w','linestyle',':');
    if (kk<=180),nyl=sprintf('%3.0fE',kk);
    else nyl=sprintf('%3.0fW',abs(360-kk)); end;
    if (proj=='npl'), 
      text(yy(2),-xx(2),nyl,'rotation',kk,'horizontal','center','vertical','top');
    else
      text(yy(2),xx(2),nyl,'rotation',-kk,'horizontal','center','vertical','bottom');
    end;
  end;
end;

xout = []; yout = [];
if ~isempty(lon_var)
  xout=minmax(eval(lon_var));
end
if ~isempty(lat_var)
  yout=minmax(eval(lat_var));
end
if ~isempty(xout) & ~isempty(yout)
  if (proj=='rec'),
    line(xout([1 2 2 1 1]),yout([1 1 2 2 1]),'color','r');
  elseif proj=='npl',
    lx=[xout(1)+diff(xout)*[0:20 20:-1:0 0]/20 ];
    ly=[yout(ones(1,21)) yout(2*ones(1,21)) yout(1)];
    xx=(90-ly).*cos(lx*pi/180-pi/2);
    yy=(90-ly).*sin(lx*pi/180-pi/2);
    line(xx,yy,'color','r');
  else 					% proj is spl
    lx=[xout(1)+diff(xout)*[0:20 20:-1:0 0]/20 ];
    ly=[yout(ones(1,21)) yout(2*ones(1,21)) yout(1)];
    xx=(90+ly).*cos(-lx*pi/180+pi/2);
    yy=(90+ly).*sin(-lx*pi/180+pi/2);
    line(xx,yy,'color','r');
  end;
end

if (proj=='rec'), grid; end;
