eddy_flux/c20a/ll_regrid.m

function [v] = ll_regrid(loni,lati,vi,lon,lat)

% Function [v] = ll_regrid(loni,lati,vi,lon,lat)
%
% INPUTS
% loni   matrix of input lon coords on the interval [0,360)
% lati   matrix of input lat coords on the interval [0,180)
% vi     matrix of values on the lati,loni points
% lon    matrix of output lon coords on the interval [0,360)
% lat    matrix of output lat coords on the interval [0,180)
%
% OUTPUTS
% v      matrix of values on the lat,lon points
%
% This function implements a very simple regridding method for
% projecting scalars defined on general high-res grids
% (eg. cubesphere) onto lower-res Lat--Lon grids.

v = [];

if nargin ~= 5
  disp('Error: there must be exactly 5 arguments!');
  return
end

%  Mandatory arguments
if (prod(size(lati)) ~= prod(size(loni))) ...
    || (prod(size(lati)) ~= prod(size(vi))) 
  disp(['Error: the lati, loni, and vi must all ' ...
        'be the same size']);
  return
end

%  lat = 0:10:350
%  lon = -80:10:80

% slon = sort(lon(:));
% dlon = diff(slon);
% lons = slon(1)-dlon(1)/2;
% lone = slon(length(slon))+dlon(length(dlon))/2;
% mlon = [ lons ; lons+cumsum(dlon) ; lone ];
% klon = [ 1 ; [1:length(dlon)]' ; length(lon) ];

% slat = sort(lat(:));
% dlat = diff(slat);
% lats = slat(1)-dlat(1)/2;
% late = slat(length(slat))+dlat(length(dlat))/2;
% mlat = [ lats ; lats+cumsum(dlat) ; late ];
% klat = [ 1 ; [1:length(dlat)]' ; length(lat) ];

%  slon = sort(lon);
%  slat = sort(lat);
%  klon = [1:length(slon)];
%  klat = [1:length(slat)];

alon = min(lon);
blon = max(lon);
dlon = blon - alon;
nlon = length(lon);

alat = min(lat);
blat = max(lat);
dlat = blat - alat;
nlat = length(lat);

vsum = zeros(length(lon),length(lat));
vcnt = zeros(length(lon),length(lat));

for kk = 1:(prod(size(vi)))

  %  ii = interp1(slon,klon,loni(kk),'nearest');
  %  jj = interp1(slat,klat,lati(kk),'nearest');
  
  flon = (loni(kk) - alon) / dlon;
  flat = (lati(kk) - alat) / dlat;

  if  (flon >= 0) && (flon <= 1) && (flat >= 0) && (flat <= 1)
    ii = 1 + round(flon*(dlon-1));
    jj = 1 + round(flat*(dlat-1));
    vcnt(ii,jj) = vcnt(ii,jj) + 1;
    vsum(ii,jj) = vsum(ii,jj) + vi(kk);
  end
  
end

inz = find(vcnt ~= 0);
vsum(inz) = vsum(inz) ./ vcnt(inz);

iz = find(vcnt == 0);
vsum(iz) = NaN;

v = vsum;
1	function [v] = ll_regrid(loni,lati,vi,lon,lat)
2
3	% Function [v] = ll_regrid(loni,lati,vi,lon,lat)
4	%
5	% INPUTS
6	% loni matrix of input lon coords on the interval [0,360)
7	% lati matrix of input lat coords on the interval [0,180)
8	% vi matrix of values on the lati,loni points
9	% lon matrix of output lon coords on the interval [0,360)
10	% lat matrix of output lat coords on the interval [0,180)
11	%
12	% OUTPUTS
13	% v matrix of values on the lat,lon points
14	%
15	% This function implements a very simple regridding method for
16	% projecting scalars defined on general high-res grids
17	% (eg. cubesphere) onto lower-res Lat--Lon grids.
18
19	v = [];
20
21	if nargin ~= 5
22	disp('Error: there must be exactly 5 arguments!');
23	return
24	end
25
26	% Mandatory arguments
27	if (prod(size(lati)) ~= prod(size(loni))) ...
28	\|\| (prod(size(lati)) ~= prod(size(vi)))
29	disp(['Error: the lati, loni, and vi must all ' ...
30	'be the same size']);
31	return
32	end
33
34	% lat = 0:10:350
35	% lon = -80:10:80
36
37	% slon = sort(lon(:));
38	% dlon = diff(slon);
39	% lons = slon(1)-dlon(1)/2;
40	% lone = slon(length(slon))+dlon(length(dlon))/2;
41	% mlon = [ lons ; lons+cumsum(dlon) ; lone ];
42	% klon = [ 1 ; [1:length(dlon)]' ; length(lon) ];
43
44	% slat = sort(lat(:));
45	% dlat = diff(slat);
46	% lats = slat(1)-dlat(1)/2;
47	% late = slat(length(slat))+dlat(length(dlat))/2;
48	% mlat = [ lats ; lats+cumsum(dlat) ; late ];
49	% klat = [ 1 ; [1:length(dlat)]' ; length(lat) ];
50
51	% slon = sort(lon);
52	% slat = sort(lat);
53	% klon = [1:length(slon)];
54	% klat = [1:length(slat)];
55
56	alon = min(lon);
57	blon = max(lon);
58	dlon = blon - alon;
59	nlon = length(lon);
60
61	alat = min(lat);
62	blat = max(lat);
63	dlat = blat - alat;
64	nlat = length(lat);
65
66	vsum = zeros(length(lon),length(lat));
67	vcnt = zeros(length(lon),length(lat));
68
69	for kk = 1:(prod(size(vi)))
70
71	% ii = interp1(slon,klon,loni(kk),'nearest');
72	% jj = interp1(slat,klat,lati(kk),'nearest');
73
74	flon = (loni(kk) - alon) / dlon;
75	flat = (lati(kk) - alat) / dlat;
76
77	if (flon >= 0) && (flon <= 1) && (flat >= 0) && (flat <= 1)
78	ii = 1 + round(flon*(dlon-1));
79	jj = 1 + round(flat*(dlat-1));
80	vcnt(ii,jj) = vcnt(ii,jj) + 1;
81	vsum(ii,jj) = vsum(ii,jj) + vi(kk);
82	end
83
84	end
85
86	inz = find(vcnt ~= 0);
87	vsum(inz) = vsum(inz) ./ vcnt(inz);
88
89	iz = find(vcnt == 0);
90	vsum(iz) = NaN;
91
92	v = vsum;