matlab_class/sample_processing/example_interp.m

function []=example_interp();
% EXAMPLE_INTERP illustrates interpolation capabilities
%   by going back and forth between gcmfaces grid and
%   longitude-latitude arrays

gcmfaces_global;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if myenv.verbose>0;
    gcmfaces_msg('===============================================');
    gcmfaces_msg(['*** entering example_interp: illustrate ' ...
        'interpolation capabilities to get in and out of gcmfaces'],'');
end;

%%%%%%%%%%%%%%%%%
%load grid and setup test case:
%%%%%%%%%%%%%%%%%

if isempty(mygrid);
   grid_load;
end;

%target locations:
lon=[-179.75:0.5:179.75]; lat=[-89.75:0.5:89.75];
[lat,lon] = meshgrid(lat,lon);

%original field:
m=mygrid.mskC(:,:,1);
fld=m.*mygrid.Depth;

%%%%%%%%%%%%%%%%%%%%%%%
%interpolate mygrid.Depth to lon-lat grid:

if myenv.verbose>0;
    gcmfaces_msg('* interpolate mygrid.Depth to lon-lat grid');
end;

fld_interp=gcmfaces_interp_2d(fld,lon,lat);

if 0;
%use sparse matrix method:
interp=gcmfaces_interp_coeffs(lon(:),lat(:));
%
tmp1=convert2vector(fld);
tmp0=1*~isnan(tmp1);
tmp1(isnan(tmp1))=0;
%
tmp0=interp.SPM*tmp0;
tmp1=interp.SPM*tmp1;
%
fld_SPM=reshape(tmp1./tmp0,size(lon));
end;

%%%%%%%%%%%%%%%%%%%%%%%
%interpolate back to mygrid.XC, mygrid.YC:

if myenv.verbose>0;
    gcmfaces_msg('* interpolate back to gcmfaces grid');
end;

fld_reinterp=gcmfaces_interp_2d(fld_interp,lon,lat,'linear');

%%%%%%%%%%%%%%%%%%%%%%%
%remap to gcmfaces grid using extrapolation:

if myenv.verbose>0;
    gcmfaces_msg('* remap to gcmfaces grid using extrapolation');
end;

fld_remap=gcmfaces_remap_2d(lon,lat,fld_interp,0,m);

%%%%%%%%%%%%%%%%%%%%%%%
%illustrate results:

figureL;
subplot(2,1,1); set(gca,'FontSize',14);
[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,fld);
pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
title('original field');
subplot(2,1,2); set(gca,'FontSize',14);
pcolor(lon,lat,fld_interp); axis([-180 180 -90 90]); shading flat;
title('interpolated field');

figureL;
subplot(2,1,1); set(gca,'FontSize',14);
[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,fld_reinterp);
pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
title('reinterpolated field');
subplot(2,1,2); set(gca,'FontSize',14);
[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,fld_remap);
pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
title('remapped field');


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if myenv.verbose>0;
    gcmfaces_msg('*** leaving example_interp');
    gcmfaces_msg('===============================================');
end;


1	function []=example_interp();
2	% EXAMPLE_INTERP illustrates interpolation capabilities
3	% by going back and forth between gcmfaces grid and
4	% longitude-latitude arrays
5
6	gcmfaces_global;
7
8	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
9	if myenv.verbose>0;
10	gcmfaces_msg('===============================================');
11	gcmfaces_msg(['*** entering example_interp: illustrate ' ...
12	'interpolation capabilities to get in and out of gcmfaces'],'');
13	end;
14
15	%%%%%%%%%%%%%%%%%
16	%load grid and setup test case:
17	%%%%%%%%%%%%%%%%%
18
19	if isempty(mygrid);
20	grid_load;
21	end;
22
23	%target locations:
24	lon=[-179.75:0.5:179.75]; lat=[-89.75:0.5:89.75];
25	[lat,lon] = meshgrid(lat,lon);
26
27	%original field:
28	m=mygrid.mskC(:,:,1);
29	fld=m.*mygrid.Depth;
30
31	%%%%%%%%%%%%%%%%%%%%%%%
32	%interpolate mygrid.Depth to lon-lat grid:
33
34	if myenv.verbose>0;
35	gcmfaces_msg('* interpolate mygrid.Depth to lon-lat grid');
36	end;
37
38	fld_interp=gcmfaces_interp_2d(fld,lon,lat);
39
40	if 0;
41	%use sparse matrix method:
42	interp=gcmfaces_interp_coeffs(lon(:),lat(:));
43	%
44	tmp1=convert2vector(fld);
45	tmp0=1*~isnan(tmp1);
46	tmp1(isnan(tmp1))=0;
47	%
48	tmp0=interp.SPM*tmp0;
49	tmp1=interp.SPM*tmp1;
50	%
51	fld_SPM=reshape(tmp1./tmp0,size(lon));
52	end;
53
54	%%%%%%%%%%%%%%%%%%%%%%%
55	%interpolate back to mygrid.XC, mygrid.YC:
56
57	if myenv.verbose>0;
58	gcmfaces_msg('* interpolate back to gcmfaces grid');
59	end;
60
61	fld_reinterp=gcmfaces_interp_2d(fld_interp,lon,lat,'linear');
62
63	%%%%%%%%%%%%%%%%%%%%%%%
64	%remap to gcmfaces grid using extrapolation:
65
66	if myenv.verbose>0;
67	gcmfaces_msg('* remap to gcmfaces grid using extrapolation');
68	end;
69
70	fld_remap=gcmfaces_remap_2d(lon,lat,fld_interp,0,m);
71
72	%%%%%%%%%%%%%%%%%%%%%%%
73	%illustrate results:
74
75	figureL;
76	subplot(2,1,1); set(gca,'FontSize',14);
77	[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,fld);
78	pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
79	title('original field');
80	subplot(2,1,2); set(gca,'FontSize',14);
81	pcolor(lon,lat,fld_interp); axis([-180 180 -90 90]); shading flat;
82	title('interpolated field');
83
84	figureL;
85	subplot(2,1,1); set(gca,'FontSize',14);
86	[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,fld_reinterp);
87	pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
88	title('reinterpolated field');
89	subplot(2,1,2); set(gca,'FontSize',14);
90	[X,Y,FLD]=convert2pcol(mygrid.XC,mygrid.YC,fld_remap);
91	pcolor(X,Y,FLD); axis([-180 180 -90 90]); shading flat;
92	title('remapped field');
93
94
95	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
96	if myenv.verbose>0;
97	gcmfaces_msg('*** leaving example_interp');
98	gcmfaces_msg('===============================================');
99	end;
100
101