function cfg = aplotit(data, config, varargin);
% cfg = aplotit(data, config, varargin);
% plots data vs/ dday as a line or pcolor (w/ colorbar) in current axes.
% if just one argument, it assumes abcissa, ordinate, and field are x,y,z.
% optional: 'config' holds plotting nformation, varargin will override
%
%% ==========================================================================
% In this context, structure "data" comes from agetmat:
%
%
% ---------- structure "data":
%          
% dday         %1xN or Nx1                                                   
% depth        %NBINSxNPRFS                                                  
% [fieldname]  %pg, sw, amp, u_meas, (see below for a detailed list)
%
% ---------- structure "config":
% fn           % plot this field name  (defaults to '')              
% abc_name     % fieldname to use as abcissa         
%              %      (if one argument, defaults to 'x')
%              %      (if two arguments, defaults to 'dday')        
%              %    can also set 'abcissa_title' for xlabel  (default: empty)
% ord_name     % fieldname to use as ordinate (defaults to data.bins)       
%              %    can also set 'ordinate_title' for ylabel
%              %    (line plot (1xN), or pcolor (nbinsxN))      
% yearbase     % should have been in config structure used in agetmat
% cur_ax       % handle to figure axes (defaults to new axes in gca)
% reverse_ydir % reverse  direction of ordinate                             
%              %       (defaults to 1 for bins or depth); otherwise 0
% index        % how to subsample abcissa :                        
% ilist        % (two names for the same thing; index wins)
%              %      - empty (use it all)
%              %      - scalar (subsample by this many) 
%              %      - length 2 vector (range in dday)
%              %      - length > 2 vector (indices of abcissa to use)
% xlim         % limits of abcissa  (default is autoscale)
% ylim         % limits of ordinate (default is autoscale)
% var_range    % variable range to use for field: [minval maxval] 
%              %      defaults to substructures data.minclipval and
%              %           data.maxclipval (values from fieldname)
% clevs        % specific ranges to use for colors and labels: 
%              %     - default depends on field
%              %     - for finer resolution in color, use config.cmapres
%              %     - if choosing manually, clevs should look like this:
%              %        [-60:20:60]   (for example, for velocity)
%              %     DO NOT include the outer bounds for the color map
% cmapname     % specific colormap to be used (default is 'jet')
% cmapres      % colormap resolution: if using direct color mapping,
%              %      subdivide coarse cmap (from clevs or best_ticks of 
%              %      var_range) by this number. default is 4. 
%              %      (use 1 to get crudest colors divisions)
%              %      Only used if 'cdirect' true (BEWARE: appends color maps)
% flagged      % matrix of flags for excluding data in plots. same scheme    
%              %     as PFLAGS, i.e. 0 is OK (not flagged), otherwise 
%              %     type of flagging is indicated by value in bin
%              %     (default is all zeros, i.e. none flagged)
%              %     (see aflagit_setup or below for parameters to set)
% lcolor       % color of text label  (defaults to black)                  
% markindex    % index into abcissa to draw vertical marks (defaults to [])
% markbinrange % first and last bin to draw the vertical line 
%              %     defaults to bins corresponding to 70% and 80% of bin range
% markcolor    % color to draw the vertical line (defaults to black)
% comment2     % additional comment; follows config.comment, specifically 
%              %     final_comment is [comment, (comment)]
% astitle      % use the comment as a title (put comment in the title field)
%              %     (default is false: comment is lower left corner)
% badrect      % Nx1 cellaray of structures with fields (abcissa ordinate)
%              %     If only one bad rectangle, badrect is a structure with
%              %     those fields, not a cellarray of sructures.
% badpoly      % Nx1 cellaray of structures with fields (abcissa ordinate)
%              %     delimiting the bad polygon vertices. 
%              %     If only one bad rectangle, badrect is a structure with
%              %     those fields, not a cellarray of sructures.
% badbins      % structure with Nx1 fields of bad abcissa and ordinate
% badprof      % Nx1 array of bad abcissa and ordinate values
% mabs         % structure with Nx1 fields of max amp bin
%
%
% examples:
% aplotit(data, config, 'flagged', flagged);
% %
% dd.x = data.dday; dd.y = data.bins; dd.z = data.pg; aplotit(dd);
% %
% data.z = something;aplotit(data, config, 'fn', 'z');

% JH 2000/12
% JH 2001/12/6 added cdata_direct.  Still need to make it subdivide colors 
%           better (assumes color ticks are regular)
%           -- also added set(gcf, 'renderer', 'painters')

% cdirect      % 0, 1: use direct color indexing (default is 0 (false)).
% historical: can take it out???


program_name = 'aplotit';%this program; for 'help'

if (nargin < 1)
  help(program_name)
end

if (nargin == 1)                 %the most basic tool. autoscale everything
   config.abc_name = 'x';        %
   config.ord_name = 'y';         
   config.fn       = 'z';        
   
   if (~isstruct(data))          %ok, well, bomb out if it isn't a structure
      help(program_name)
      error('first argunemt must be a structure')
   end

end

%===== set up the data: make sure config is up-to-date, but don't use 
%      anything yet.

if (~isfield(config, 'cur_ax')),           config.cur_ax = [];       end
if (~isfield(config, 'xlim')),             config.xlim = [];         end
if (~isfield(config, 'ylim')),             config.ylim = [];         end
if (~isfield(config, 'var_range')),        config.var_range = [];    end
if (~isfield(config, 'fn')),               config.fn = 'pg';         end
if (~isfield(config, 'flagged')),          config.flagged = 0;       end
if (~isfield(config, 'lcolor')),           config.lcolor = 'k';      end
if (~isfield(config, 'comment')),          config.comment = '';      end
if (~isfield(config, 'comment2')),         config.comment2 = '';     end
if (~isfield(config, 'abc_name')),         config.abc_name = 'dday'; end  
if (~isfield(config, 'ord_name')),         config.ord_name = 'bins'; end
if (~isfield(config, 'markbinrange')),     config.markbinrange = []; end
if (~isfield(config, 'markcolor')),        config.markcolor = 'k';   end
if (~isfield(config, 'markindex')),        config.markindex = [];    end
if (~isfield(config, 'index')),            config.index = [];        end
if (~isfield(config, 'ilist')),            config.ilist = [];        end
if (~isfield(config, 'badprof')),          config.badprof = [];    end
if (~isfield(config, 'badrect')),          config.badrect = [];    end
if (~isfield(config, 'badpoly')),          config.badpoly = [];    end
if (~isfield(config, 'badbins')),          config.badbins = [];    end
if (~isfield(config, 'mabs')),             config.mabs = [];    end
if (~isfield(config, 'clevs')),            config.clevs = [];         end
if (~isfield(config, 'cmapname')),         config.cmapname = '';   end
if (~isfield(config, 'cmapres')),          config.cmapres = 4;   end
if (~isfield(config, 'cdirect')),          config.cdirect = 1;   end
if (~isfield(config, 'astitle')),          config.astitle = 0;       end
if (~isfield(config, 'yearbase')),         config.yearbase = -9999;  end
if (~isfield(config, 'prefix')),           config.prefix = 'a_';     end



%  <----- (config comes in)    ....  ----> cfg goes out
%  set defaults and finalize config values
%        i.e. defaults from (1) agetmat, (2) config (3) arg list

if (nargin > 2)%   we'll assume the rest is a set of name,value pairs
   %override input structure with name,value pairs
   cfg = fillstruct(config, varargin,0);
   cfg.box_ax = cfg.cur_ax;  %better make them match.
else
   cfg = config;
end

%now we have the correctly updated field 
%get the data and then set the clipping


if (~isfield(data, cfg.fn))     
   warning(sprintf('%s: field %s does not exist\n', program_name, cfg.fn))
   return
else
   fieldname  = cfg.fn;               %contour this name
   field      = getfield(data, fieldname);
   if length(unique(field(:))) == 1;
      return
   end
   cfg.var_range =  getclips(field, fieldname, cfg);
end


if (isemptys(field))   %cfg will return unchanged
   warning(sprintf('%s: field %s is empty\n', program_name, cfg.fn))
   return
end

%get abcissa  ------------------------------
abcissa = getfield(data, getfield(cfg, 'abc_name'));
abcissa_title = '';
if (~isfield(cfg, 'abcissa_title'))  %no title yet
   if (isfield(cfg, 'title'))
      if (isfield(cfg.title, 'abc_name')) %there is a default title
         abcissa_title = getfield(cfg.title, 'abc_name');
      end
   end
end

if (length(abcissa) ~= length(abcissa(:))) %then it is not MxN
  abcissa = abcissa(1,:);                  %assume it is shaped like depth
end

% get abcissa indices (either name is OK; use 'ilist' values if both exist)
if ~isemptys(cfg.index)  
   cfg.ilist = [];
else
   if ~isemptys(cfg.ilist)
      cfg.index = cfg.ilist;
   end
end

% now deal with cfg.index only.
if (isemptys(cfg.index)) 
   abci = 1:length(abcissa);     % use all indices      
elseif (length(cfg.index) == 1)     % subsample by this number
   abci = [1:cfg.index:length(abcissa)];
elseif (length(cfg.index) == 2)     %it's a range
   abci = find (...
      (cfg.index(1) <= abcissa) & ...
      (cfg.index(2) >= abcissa) );
else                                      %it's already set of indices
   abci = cfg.index;
   abci(abci < 1) = [];
   abci(abci > length(abcissa)) = [];
end
cfg.index = abci;

% get ordinate ------------------------------
ordinate = getfield(data, getfield(cfg, 'ord_name'));
ordinate_title = '';
if (~isfield(cfg, 'ordinate_title'))  %no title yet
  if (isfield(cfg, 'title'))
    if (isfield(cfg.title, 'ord_name')) %there is a default title
      ordinate_title = getfield(cfg.title, 'ord_name');
    end
  end
end

if (length(ordinate) ~= length(ordinate(:))) %then it is not MxN
  ordinate = ordinate(:,1);                  %assume it is shaped like depth
end

if (~isfield(cfg, 'reverse_ydir') | ...
    (strcmp(cfg.ord_name, 'z')) | ...
    (strcmp(cfg.ord_name, 'bins')) | ...
    (strcmp(cfg.ord_name, 'depth')))
  cfg.reverse_ydir = 1;
else
  cfg.reverse_ydir = 0;
end

% get the bin range for vertical marks
if isemptys(cfg.markbinrange)
   cfg.markbinrange = [ floor(.7*length(ordinate)) ceil(.8*length(ordinate))];
end

%set the title ------------------------------
titlestr = '';
if (isfield(cfg, 'title'))
  if (isfield(cfg.title, fieldname))
    titlestr   =  getfield(cfg.title,fieldname);
  end
end

if (isemptys(cfg.comment)),  cfg.comment = '';  end
if (isemptys(cfg.comment2)), cfg.comment2 = '';  end

numcomments =  2 - sum( [strcmp(cfg.comment, '') strcmp(cfg.comment2, '')]);
if (numcomments == 0)
  comment = '';
elseif (numcomments  == 1  )
   comment = sprintf(' (%s%s)', cfg.comment, cfg.comment2);
else %numcomments = 2
   comment = sprintf(' (%s,%s)', cfg.comment, cfg.comment2);
end
titlestr = [titlestr,comment];


%multiplicative mask (from flagged) dot-multiply -------
if (isemptys(cfg.flagged))
  cfg.flagged = zeros(size(field));
end

okmask = (cfg.flagged == 0).*ones(size(field));
okmask(find(okmask==0)) = NaN;  

% axis  ------------------------------
if (isemptys(cfg.xlim))
   cfg.xlim =  [min(abcissa(cfg.index))  max(abcissa(cfg.index))];
end

if (isemptys(cfg.ylim))
   cfg.ylim =  [min(ordinate(:)) max(ordinate(:))];
end

%set the axes ------------------------------
if (isemptys(cfg.cur_ax))
  figure(gcf);
  clf
  cfg.cur_ax = axes;
end


%get the colormap
if (isemptys(cfg.cmapname))
   cfg.cmapname = 'jet';
   %override if possible
   if (isfield(cfg, 'cmaps'))
      if (isfield(cfg.cmaps, cfg.fn)) %there is a default title
         cfg.cmapname = getfield(cfg.cmaps, cfg.fn);
      end
   end
end
   
%set up the contour levels
if (isemptys(cfg.clevs))
   %set up colors
   cfg.clevs = best_ticks(cfg.var_range);
end
%this algorithm should be made better to deal with oddly spaced clev
if (~cfg.cdirect)
   cfg.cmapres = 2;  %override the use of it if using caxis
end
color_dtick = mean(diff(cfg.clevs))/cfg.cmapres;
colorticks = cfg.clevs(1):color_dtick:cfg.clevs(end);
if (cfg.cdirect) %pad by alot
   cfg.clevs = [-10000*max(abs(cfg.var_range))  ...
             colorticks 10000*max(abs(cfg.var_range))];
else %pad by a little bit
   cfg.clevs = [colorticks colorticks(end)+color_dtick];
end   
   
   
fprintf('about to plot %s\n',fieldname')

% set up plotting ==================================================

% it's a line: do a line plot
if (length(field) == length(field(:))) 
   h=plot(abcissa(cfg.index),field(cfg.index),'k');
   set(h,'linewidth',.1)
   %here abcissa = field
   if cfg.astitle
      title(titlestr)
   else
      h=text( cfg.xlim(1) + .05*diff(cfg.xlim),...    
              cfg.ylim(1) + .05*diff(cfg.ylim),titlestr);
      set(h,'fontweight','bold','fontsize',10)
      axis([cfg.xlim cfg.ylim])
   end
   xlabel(abcissa_title)
   %return control to original axes
   return  %get out if it is just a line plot 
end

%----------------------------------------------------------------------
%...  only going past here if it's a color plot

%set up the colorbar
nclev = length(cfg.clevs);
if (cfg.cdirect)
   cmap = mk_cmap(nclev - 1, cfg.cmapname); 
   orig_cmap =  get(gcf, 'colormap');
   allcmap = [orig_cmap; cmap];
   colormap(allcmap);
else
   cmap = mk_cmap(nclev + 1, cfg.cmapname); 
   colormap(cmap);
end

%make the plot
axes(cfg.cur_ax);
abcissa = abcissa(cfg.index);
cfg.h_ax = pcolor(abcissa, ordinate, field(:,cfg.index).*okmask(:,cfg.index));
shading flat
hold on
if (~isemptys(cfg.markindex))
   ajunk = abcissa(:);
   ojunk = ordinate(:);
   plot([ajunk(cfg.markindex) ajunk(cfg.markindex)]',...
        [ojunk(cfg.markbinrange(1)) ojunk(cfg.markbinrange(2))]',...
        'color', cfg.markcolor);
   clear ajunk
   clear ojunk
end

if (~isemptys(cfg.badprof))
   ojunk = ordinate(:);
   for ii=1:length(cfg.badprof)
      meanda = mean(diff(abcissa));  %meanda = mean delta abcissa
      plot((meanda/2 + cfg.badprof(ii))*[1 1]', [ojunk(1) ojunk(end)]','k');
   end
end


if (~isemptys(cfg.badbins))
   meanda = mean(diff(abcissa));  %meanda = mean delta abcissa
   meando = mean(diff(ordinate));  %meanda = mean delta ordinate
   hh= plot(cfg.badbins.abcissa + meanda/2, ...
            cfg.badbins.ordinate + meando/2, 'm.');
   set(hh,'markersize', 2)
end


if (~isemptys(cfg.mabs))
   meanda = mean(diff(abcissa));  %meanda = mean delta abcissa
   meando = mean(diff(ordinate));  %meanda = mean delta ordinate
   hh=plot(cfg.mabs.abcissa + meanda/2, cfg.mabs.ordinate + meando/2, 'w');
   set(hh,'linewidth',3)
   hh=plot(cfg.mabs.abcissa + meanda/2, cfg.mabs.ordinate + meando/2, 'k')
   set(hh,'linewidth',2)
end


if (~isemptys(cfg.badrect))
   if (isstruct(cfg.badrect))
      if (length(cfg.badrect.abcissa) > 1)
         badrect_da = diff(cfg.badrect.abcissa(:)');
         badrect_da = badrect_da([1:end end]);
      else
         badrect_da = mean(diff(abcissa));
      end
      badrect_abcissa = cfg.badrect.abcissa(:)' + badrect_da/2;      
      
         
      if (length(badrect.ordinate) > 1)
         badrect_do =  [diff(cfg.badrect.ordinate)];
         badrect_do = badrect_do([1:end end]);
      else
         badrect_do = mean(diff(ordinate));
      end
      badrect_ordinate = cfg.badrect.ordinate(:)' + badrect_do/2;
         

      [abM, ordM] = meshgrid(badrect_abcissa, badrect_ordinate);
      hh=plot(abM, ordM, 'm.');
      set(hh,'markersize', 2)
   else
      for ii=1:length(cfg.badrect)
         if (length(cfg.badrect{ii}.abcissa) > 1)
            badrect_da = diff(cfg.badrect{ii}.abcissa(:)');
            badrect_da = badrect_da([1:end end]);
         else
            badrect_da = mean(diff(abcissa));
         end
         badrect_abcissa = cfg.badrect{ii}.abcissa(:)' + badrect_da/2;
         
         
         if (length(cfg.badrect{ii}.ordinate) > 1)
            badrect_do =  [diff(cfg.badrect{ii}.ordinate)];
            badrect_do = badrect_do([1:end end]);
         else
            badrect_do = mean(diff(ordinate));
         end
         badrect_ordinate = cfg.badrect{ii}.ordinate(:)' + badrect_do/2;
      
         [abM, ordM] = meshgrid(badrect_abcissa, badrect_ordinate);
         hh=plot(abM, ordM, 'm.');
         set(hh,'markersize', 2)         
      end
   end
end

%% vertices stored: must use inpolygon
if (~isemptys(cfg.badpoly))
   [abcM, ordM] = meshgrid(abcissa, ordinate);
   allabc = abcM(:);
   allord = ordM(:);
   da = mstdgap(diff(abcissa));
   do = mstdgap(diff(ordinate));
   
   if (isstruct(cfg.badpoly))
      ini = inpolygon(allabc, allord, ...
                      cfg.badpoly.abcissa,...
                      cfg.badpoly.ordinate);
      if (length(find(ini)) > 0)
         polyabc = allabc(ini);
         polyord = allord(ini);
         hh=plot(polyabc + da/1, polyord + do/2, 'm.');
         set(hh,'markersize', 2)
      end
   else
      for ii=1:length(cfg.badpoly)
         ini = inpolygon(allabc, allord, ...
                         cfg.badpoly{ii}.abcissa,...
                         cfg.badpoly{ii}.ordinate);
         if (length(find(ini)) > 0)
            polyabc = allabc(ini);
            polyord = allord(ini);
            hh=plot(polyabc + da/1, polyord + do/2, 'm.');
            set(hh,'markersize', 2)
         end
      end
   end
end

set(cfg.cur_ax,'box','off')
set(cfg.cur_ax, 'TickDir', 'out', 'TickLength', [0.01 0.01])
set(cfg.cur_ax,'xlim', cfg.xlim);
set(cfg.cur_ax,'ylim', cfg.ylim);
if (cfg.reverse_ydir) 
  set(cfg.cur_ax,'ydir', 'reverse')
end

if cfg.astitle
   title(titlestr)
else
   h=text( cfg.xlim(1) + .05*diff(cfg.xlim),...    
           cfg.ylim(end) - .05*diff(cfg.ylim),titlestr,...
           'fontweight','bold','fontsize',10,'color',cfg.lcolor);
end
axis([cfg.xlim cfg.ylim])
xlabel(abcissa_title)
ylabel(ordinate_title)
if ( (strcmp(lower(cfg.abc_name), 'lon')  | ...
      strcmp(lower(cfg.abc_name), 'lat') ))
  lltick(cfg.cur_ax, 'x', lower(cfg.abc_name))
end
if ( (strcmp(lower(cfg.ord_name), 'lon')  | ...
      strcmp(lower(cfg.ord_name), 'lat') ))
  lltick(cfg.cur_ax, 'y', lower(cfg.abc_name))
end

if (strcmp(lower(cfg.abc_name), 'heading')) 
  lltick(cfg.cur_ax, 'x', 'none');
end
if (strcmp(lower(cfg.ord_name), 'heading')) 
  lltick(cfg.cur_ax, 'y', 'none');
end

%set up the color bar
if (cfg.cdirect)
   cmapi = length(orig_cmap(:,1)) + 1;
   cdata_direct(cfg.h_ax, cfg.clevs, cmapi);
   cfg.cbar_ax = cbar_direct(cfg.clevs, cmapi, 'vert', '');
   ticks = get(cfg.cbar_ax, 'YTick');
   set(cfg.cbar_ax, 'YTick', ticks(1:cfg.cmapres:end));
   set(findobj(cfg.cbar_ax, 'Type', 'patch'), 'EdgeColor', 'none');
   
   set(gcf, 'renderer', 'painters')
else
   caxis(cfg.var_range);
   cfg.cbar_ax = dcolbar(cfg.clevs, 'vert', ' ');
   h_dcbp = findobj(cfg.cbar_ax,'Type','patch');
%   cf_index(h_dcbp, cfg.clevs, [1:length(cfg.clevs)]);
   set(cfg.cbar_ax, 'XTick', [], 'YTick', colorticks);
   set(cfg.cbar_ax, 'YTicklabel', num2cell(colorticks));
   set(cfg.cbar_ax, 'TickDir', 'out', 'TickLength', [0.02 0.02])
   set(cfg.cbar_ax, 'yaxislocation','right')
end

cfg.box_ax = axes('pos',get(cfg.cur_ax,'pos'),'color','none','box','on',...
                     'xtick',[],'ytick',[]);

%======================================================================
function [var_range] = getclips(field, fieldname, config);

% get the specified clipping and plot ranges
%first try the data
var_range = [ min(field(:)) max(field(:))];
%now if exists in defaults, override
if (isfield(config, 'minclipval'))
  if (isfield(config.minclipval, fieldname))
    var_range(1) = getfield(config.minclipval, fieldname);
  end
end

if (isfield(config, 'maxclipval'))
  if (isfield(config.maxclipval, fieldname))
    var_range(2) = getfield(config.maxclipval, fieldname);
  end
end

%now test config (maybe it was overridden here?)
if (~isemptys(config.var_range))
      var_range = getfield(config, 'var_range');
end

%======================================================================

