contourfnu.m

function varargout = contourfnu(x,y,data,v,cmap,pos_colorbar,overticklabel,method,ninterp,nancolor)
% non-uniform contourf/imagesc/colorbar
%
%    Input variables
% x             : x-coordinates of grid, vector or 2d matrix
% y             : y-coordinates of grid, vector or 2d matrix
%                 If x and y are vectors, then length(x)==size(z,2) and length(y)==size(Z,1). 
%                 If x and y are 2d matrix, they are generated by meshgrid
% data          : 2d matrix to be ploted
%
%     Optional:
% v             : vector of contour levels (default:linspace(datamin,datamax,10))
% cmap          : color map array (default:jet)
% pos_colorbar  : 'none', or location with respect to the axes (default:'eastoutside')
% overticklabel : whether or not label the overrange ticks at colorbar (default:true)
% method        : imagesc, contourf, contour or pcolor (default:imagesc)
% ninterp       : repeatedly interp times in each dimension (default:0)
% nancolor      : axis backgroud color
%
%     Output variable
% hout          : structure with handles
%                 .h     plot handle
%                 .c     contour matrix (method='contourf')
%                 .hc    colorbar handle (pos_colorbar~='none')


datamax = max(data(:));
datamin = min(data(:));

if(~exist('v','var')||isempty(v)),                         v = linspace(datamin,datamax,10); end
if(~exist('cmap','var')||isempty(cmap)),                   cmap = jet;                       end
if(~exist('pos_colorbar','var')||isempty(pos_colorbar)),   pos_colorbar = 'eastoutside';     end
if(~exist('overticklabel','var')||isempty(overticklabel)), overticklabel = true;             end
if(~exist('method','var')||isempty(method)),               method = 'imagesc';               end
if(~exist('ninterp','var')||isempty(ninterp)),             ninterp = 0;                      end
if(ninterp>0)
    data=interp2(data,ninterp);
    if( strcmp(method,'contourf')||strcmp(method,'contour')||strcmp(method,'pcolor') )
        if( isvector(x)&&isvector(y) )  % vector to meshgrid
            [x,y] = meshgrid(x,y);
        end
        x = interp2(x,ninterp);
        y = interp2(y,ninterp);
    end
end
if( strcmp(method,'imagesc')&&~isvector(x) ) % meshgrid to vector
    x = x(1,:);
    y = y(:,1);
end

if(isrow(v)), v = v'; end
v = sort(v);
% extendmin = datamin<v(1);
% extendmax = datamax>v(end);
if(datamin<v(1))
    v = [-inf;v];
end
if(datamax>v(end))
    v = [v;inf];
end
nlev = length(v);

% piecewise linear mapping,  v -> 1:nlev,    data -> 1 to nlev
z = zeros(size(data));
for i=1:nlev-1
    if( i==1 && v(i)==-inf )
        index = data<v(2);
        z(index) = i+(data(index)-datamin)/(v(2)-datamin);
    elseif( i==(nlev-1) && v(i)==inf )
        index = data>=v(end-1);
        z(index) = i+(data(index)-v(end-1))/(datamax-v(end-1));
    else
        index = data>=v(i) & data<v(i+1);
        z(index) = i+(data(index)-v(i))/(v(i+1)-v(i));
    end
end
z(isnan(data)) = nan;

% draw
if(strcmp(method,'imagesc'))
%     hout.h = imagesc(x,y,z);axis xy
    hout.h = imagesc(x,y,z,'AlphaData',~isnan(z));axis xy
elseif(strcmp(method,'contourf'))
    [hout.c,hout.h] = contourf(x,y,z,1:nlev-1);
elseif(strcmp(method,'contour'))
    [hout.c,hout.h] = contour(x,y,z,1:nlev-1);
elseif(strcmp(method,'pcolor'))
    hout.h = pcolor(x,y,z);shading flat
end
if(exist('nancolor','var')), set(gca,'color',nancolor); end
    
% set colormap and colorbar
Ncmap=size(cmap,1);
Ndraw=nlev-1;
cmapdraw(1:Ndraw,:)=cmap( round(linspace(1,Ncmap,Ndraw)) ,:);
colormap(gca,cmapdraw)
caxis([1,nlev])
if(~strcmp(pos_colorbar,'none'))
    hc = colorbar('location',pos_colorbar);
    if(overticklabel)
        vlabel = v;
        if(v(1)==-inf), vlabel(1) = datamin; end
        if(v(end)==inf), vlabel(end) = datamax; end
    else
        vlabel = cellstr(num2str(v));
        if(v(1)==-inf), vlabel(1) = {''}; end
        if(v(end)==inf), vlabel(end) = {''}; end
    end
    if(any(strcmp(pos_colorbar,{'eastoutside','westoutside','east','west'})))
        ylimits = get(hc,'Ylim');
        ystep = (ylimits(2)-ylimits(1))/Ndraw;
        set(hc,'ytick',ylimits(1):ystep:ylimits(2))
        set(hc,'yticklabel',vlabel)
    elseif(any(strcmp(pos_colorbar,{'southoutside','northoutside','south','north'})))
        xlimits = get(hc,'Xlim');
        xstep = (xlimits(2)-xlimits(1))/Ndraw;
        set(hc,'xtick',xlimits(1):xstep:xlimits(2))
        set(hc,'xticklabel',vlabel)
    end
    hout.hc = hc;
end

if nargout > 0
    varargout{1} = hout;
end