course-idma2016/matlab/idma_interp_2d.m

function [val]=idma_interp_2d(lon,lat,fld,method);
%IDMA_INTERP_2D linearly interpolates gcmfaces field (fld) to
%   a set of locations (lon, lat) using one of several methods:
%   'natural', 'linear', 'nearest', or 'mix' (default).
%   The 'mix' is `natural' extended with 'nearest' when
%   the input field has been land-masked with NaNs.
%
%Example:
%     lon=[-179.9:0.2:179.9]; lat=[-89.9:0.2:89.9];
%     [lat,lon] = meshgrid(lat,lon);
%     fld=mygrid.Depth.*mygrid.mskC(:,:,1);
%     [val]=idma_interp_2d(lon,lat,fld);
%     figureL; pcolor(lon,lat,val); shading flat;

gcmfaces_global;
if isempty(whos('method'));
    method='mix';
end;

if isempty(which('DelaunayTri'));
    error('this code needs DelaunayTri that is not found');
end;

XC=convert2array(mygrid.XC);
YC=convert2array(mygrid.YC);
VEC=convert2array(fld);
val=NaN*lon;
for ii=1:3;
    if ii==1;
        myXC=XC; myXC(myXC<0)=myXC(myXC<0)+360;
        jj=find(lon>90);
    elseif ii==2;
        myXC=XC;
        jj=find(lon>=-90&lon<=90);
    else;
        myXC=XC; myXC(myXC>0)=myXC(myXC>0)-360;
        jj=find(lon<-90);
    end;
    kk=find(~isnan(myXC));
    TRI=DelaunayTri(myXC(kk),YC(kk));
    if strcmp(method,'mix');
        F = TriScatteredInterp(TRI, VEC(kk),'natural');
        tmp1=F(lon(jj),lat(jj));
        F = TriScatteredInterp(TRI, VEC(kk),'nearest');
        tmp2=F(lon(jj),lat(jj));
        tmp1(isnan(tmp1))=tmp2(isnan(tmp1));
        val(jj)=tmp1;
    else;
        F = TriScatteredInterp(TRI, VEC(kk),method);
        val(jj)=F(lon(jj),lat(jj));
    end;
end;


1	gforget	1.1	function [val]=idma_interp_2d(lon,lat,fld,method);
2			%IDMA_INTERP_2D linearly interpolates gcmfaces field (fld) to
3			% a set of locations (lon, lat) using one of several methods:
4			% 'natural', 'linear', 'nearest', or 'mix' (default).
5			% The 'mix' is `natural' extended with 'nearest' when
6			% the input field has been land-masked with NaNs.
7			%
8			%Example:
9			% lon=[-179.9:0.2:179.9]; lat=[-89.9:0.2:89.9];
10			% [lat,lon] = meshgrid(lat,lon);
11			% fld=mygrid.Depth.*mygrid.mskC(:,:,1);
12			% [val]=idma_interp_2d(lon,lat,fld);
13			% figureL; pcolor(lon,lat,val); shading flat;
14
15			gcmfaces_global;
16			if isempty(whos('method'));
17			method='mix';
18			end;
19
20			if isempty(which('DelaunayTri'));
21			error('this code needs DelaunayTri that is not found');
22			end;
23
24			XC=convert2array(mygrid.XC);
25			YC=convert2array(mygrid.YC);
26			VEC=convert2array(fld);
27			val=NaN*lon;
28			for ii=1:3;
29			if ii==1;
30			myXC=XC; myXC(myXC<0)=myXC(myXC<0)+360;
31			jj=find(lon>90);
32			elseif ii==2;
33			myXC=XC;
34			jj=find(lon>=-90&lon<=90);
35			else;
36			myXC=XC; myXC(myXC>0)=myXC(myXC>0)-360;
37			jj=find(lon<-90);
38			end;
39			kk=find(~isnan(myXC));
40			TRI=DelaunayTri(myXC(kk),YC(kk));
41			if strcmp(method,'mix');
42			F = TriScatteredInterp(TRI, VEC(kk),'natural');
43			tmp1=F(lon(jj),lat(jj));
44			F = TriScatteredInterp(TRI, VEC(kk),'nearest');
45			tmp2=F(lon(jj),lat(jj));
46			tmp1(isnan(tmp1))=tmp2(isnan(tmp1));
47			val(jj)=tmp1;
48			else;
49			F = TriScatteredInterp(TRI, VEC(kk),method);
50			val(jj)=F(lon(jj),lat(jj));
51			end;
52			end;
53
54