function mcolor(data, cint, mstyle);
%
%  mcolor(data, cint, mstyle);
%
%  This function uses the mapping toolbox to shade the
%  input 'data' matrix.  It assumes that the x and y axes
%  are defined in the global variables XAX and YAX, 
%  respectively.
%
%  Note that this function works best for cylindrical, or
%  pseudocylindrical projections.  (Try typing 'maps' to 
%  see the style of the projection you want to use).  It
%  may work for others too, but I haven't tested it.
%
%  The default mstyle is Global Isometric ('giso'), 
%  centered at the center of the variable FRAME.  This should
%  all be pretty easy to change, though.
%

%  Dan Vimont, 2 September, 1999

data = squeeze(data);

%  First, determine if hold is 'on', or 'off'.

next_ax = lower(get(gca, 'NextPlot'));
if strcmp(next_ax, 'replace');  cla;  end;

%  Get the global varaibles XAX YAX and FRAME

global XAX YAX FRAME
if (size(XAX, 1) == 1); XAX = XAX'; end;
if (size(YAX, 1) == 1); YAX = YAX'; end;

[m, n] = size(data);
if or((m~=length(YAX)), (n~=length(XAX)));
  error('[length(YAX) length(XAX)] must equal size(data)');
end

%  Set defaults

if nargin < 2; cint = (max(max(data))-min(min(data)))/10; end;
if nargin < 3; mstyle = 'giso'; end;
origin = [];
lab = 0;
line_style = 'k';

%  Define contour interval;

if isscalar(cint);
  clev = sort([-cint:-cint:min(min(data)) 0 cint:cint:max(max(data))]);
else
  clev = cint;
end

%  Redefine XAX, YAX and data for global contouring -- allowing 
%  for wrapping around the globe.

if (FRAME(2) - FRAME(1) == 360);
  data2 = [data data(:,1)];
  XAX2 = [XAX; XAX(1)+360];
else
  data2 = data;
  XAX2 = XAX;
end;

%  Define map axis:

axesm('mapprojection', mstyle, ...
      'origin', origin, ...
      'maplatlimit', FRAME(3:4), ...
      'maplonlimit', FRAME(1:2));

%  Plot the data

contourfm(YAX, XAX2, data2, clev);
axis_limits

set(gca, 'NextPlot', next_ax);
