matlab_class/sample_processing/example_griddata.m

function []=example_griddata(choiceV3orV4);
%object:    a griddata example within gcmfaces
%inputs:    choiceV3orV4 ('v3' or 'v4') selects the sample GRID

gcmfaces_global;
if myenv.verbose>0;
    gcmfaces_msg('===============================================');
    gcmfaces_msg(['*** entering example_griddata that will define ' ...
        'a 2x2 gridded variable (an SSH field here), and then ' ...
        'use griddata to map it to the chosen grid. '],'');
end;

%%%%%%%%%%%%%%%%%
%load parameters:
%%%%%%%%%%%%%%%%%

gcmfaces_global;
dir0=[myenv.gcmfaces_dir '/sample_input/'];
if strcmp(choiceV3orV4,'v4');
    nF=5;
    fileFormat='compact';
    dirGrid=[dir0 '../../GRID/'];
else;
    error('files are missing');
    nF=1;
    fileFormat='straight';
    dirGrid=[dir0 '/GRID' choiceV3orV4  '/'];
end;

mygrid=[]; grid_load(dirGrid,nF,fileFormat,1);

%%%%%%%%%%%%%%%%%%%%%%%
%get sample data: subsampled mygrid.Depth
if myenv.verbose>0;
    gcmfaces_msg('* define test case : subsampled bathymetry');
end;
original_Depth=mygrid.Depth;
original_Depth(original_Depth==0)=NaN;
lon=mygrid.XC; lat=mygrid.YC;
for ff=1:mygrid.nFaces;
    lon{ff}(1:2:end,:)=NaN;
    lat{ff}(1:2:end,:)=NaN;
    original_Depth{ff}(1:2:end,:)=NaN;
    lon{ff}(:,1:2:end)=NaN;
    lat{ff}(:,1:2:end)=NaN;
    original_Depth{ff}(:,1:2:end)=NaN;
end;
lon=convert2vector(lon);
lat=convert2vector(lat);
vec=convert2vector(original_Depth);

%%%%%%%%%%%%%%%%%%%%%%%
%make sure to accomodate the two longitude conventions:
if myenv.verbose>0;
    gcmfaces_msg('* accomodate both lon. conventions, then convert to data vector');
end;
x=[lon-360;lon;lon+360]; y=[lat;lat;lat]; z=[vec;vec;vec];
jj=find(~isnan(z)); x=x(jj); y=y(jj); z=z(jj);

%%%%%%%%%%%%%%%%%%%%%%%
%do the interpolation:
if myenv.verbose>0;
    gcmfaces_msg('* map/interpolate data using griddata');
end;
griddata_Depth=gcmfaces(5);
for ii=1:mygrid.nFaces;
    xi=mygrid.XC{ii}; yi=mygrid.YC{ii};
    zi = griddata(x',y',z',xi,yi);
    griddata_Depth{ii}=zi;
end;

%msk=mygrid.hFacC(:,:,1); griddata_Depth(find(msk==0))=NaN;

%%%%%%%%%%%%%%%%%%%%%%%
%illustrate the result:

figureL;
subplot(2,1,1); set(gca,'FontSize',14);
[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,original_Depth);
pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
title('input to griddata');
subplot(2,1,2); set(gca,'FontSize',14);
[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,griddata_Depth);
pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
title('output of griddata');

if myenv.verbose>0;
    gcmfaces_msg('* note of caution : griddata puts values on land. ');
    gcmfaces_msg('*** leaving example_griddata');
    gcmfaces_msg('===============================================');
end;
1	gforget	1.8	function []=example_griddata(choiceV3orV4);
2			%object: a griddata example within gcmfaces
3			%inputs: choiceV3orV4 ('v3' or 'v4') selects the sample GRID
4	gforget	1.1
5	gforget	1.11	gcmfaces_global;
6			if myenv.verbose>0;
7			gcmfaces_msg('===============================================');
8			gcmfaces_msg(['*** entering example_griddata that will define ' ...
9			'a 2x2 gridded variable (an SSH field here), and then ' ...
10			'use griddata to map it to the chosen grid. '],'');
11			end;
12
13	gforget	1.1	%%%%%%%%%%%%%%%%%
14			%load parameters:
15			%%%%%%%%%%%%%%%%%
16	gforget	1.12
17			gcmfaces_global;
18	gforget	1.9	dir0=[myenv.gcmfaces_dir '/sample_input/'];
19	gforget	1.12	if strcmp(choiceV3orV4,'v4');
20			nF=5;
21			fileFormat='compact';
22			dirGrid=[dir0 '../../GRID/'];
23			else;
24			error('files are missing');
25			nF=1;
26			fileFormat='straight';
27			dirGrid=[dir0 '/GRID' choiceV3orV4 '/'];
28			end;
29
30			mygrid=[]; grid_load(dirGrid,nF,fileFormat,1);
31	gforget	1.1
32			%%%%%%%%%%%%%%%%%%%%%%%
33	gforget	1.12	%get sample data: subsampled mygrid.Depth
34	gforget	1.11	if myenv.verbose>0;
35	gforget	1.12	gcmfaces_msg('* define test case : subsampled bathymetry');
36	gforget	1.11	end;
37	gforget	1.12	original_Depth=mygrid.Depth;
38			original_Depth(original_Depth==0)=NaN;
39			lon=mygrid.XC; lat=mygrid.YC;
40			for ff=1:mygrid.nFaces;
41			lon{ff}(1:2:end,:)=NaN;
42			lat{ff}(1:2:end,:)=NaN;
43			original_Depth{ff}(1:2:end,:)=NaN;
44			lon{ff}(:,1:2:end)=NaN;
45			lat{ff}(:,1:2:end)=NaN;
46			original_Depth{ff}(:,1:2:end)=NaN;
47			end;
48			lon=convert2vector(lon);
49			lat=convert2vector(lat);
50			vec=convert2vector(original_Depth);
51	gforget	1.1
52			%%%%%%%%%%%%%%%%%%%%%%%
53	gforget	1.11	%make sure to accomodate the two longitude conventions:
54			if myenv.verbose>0;
55			gcmfaces_msg('* accomodate both lon. conventions, then convert to data vector');
56			end;
57	gforget	1.12	x=[lon-360;lon;lon+360]; y=[lat;lat;lat]; z=[vec;vec;vec];
58	gforget	1.1	jj=find(~isnan(z)); x=x(jj); y=y(jj); z=z(jj);
59
60	gforget	1.11	%%%%%%%%%%%%%%%%%%%%%%%
61			%do the interpolation:
62			if myenv.verbose>0;
63			gcmfaces_msg('* map/interpolate data using griddata');
64			end;
65	gforget	1.12	griddata_Depth=gcmfaces(5);
66	gforget	1.8	for ii=1:mygrid.nFaces;
67	gforget	1.12	xi=mygrid.XC{ii}; yi=mygrid.YC{ii};
68			zi = griddata(x',y',z',xi,yi);
69			griddata_Depth{ii}=zi;
70	gforget	1.1	end;
71
72	gforget	1.12	%msk=mygrid.hFacC(:,:,1); griddata_Depth(find(msk==0))=NaN;
73	gforget	1.1
74			%%%%%%%%%%%%%%%%%%%%%%%
75			%illustrate the result:
76
77	gforget	1.12	figureL;
78			subplot(2,1,1); set(gca,'FontSize',14);
79			[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,original_Depth);
80			pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
81			title('input to griddata');
82			subplot(2,1,2); set(gca,'FontSize',14);
83			[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,griddata_Depth);
84			pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
85			title('output of griddata');
86	gforget	1.1
87	gforget	1.11	if myenv.verbose>0;
88			gcmfaces_msg('* note of caution : griddata puts values on land. ');
89			gcmfaces_msg('*** leaving example_griddata');
90			gcmfaces_msg('===============================================');
91			end;