matlab_class/gcmfaces_calc/gcmfaces_bindata.m

function [varargout]=gcmfaces_bindata(varargin);

global mygrid mytri;

if nargin==0;
%generate delaunay triangulation:
XC=convert2array(mygrid.XC);
YC=convert2array(mygrid.YC);
%needed so that we do not loose data
XC(isnan(XC))=270;
YC(isnan(YC))=180;
%
TRI=delaunay(XC,YC); nxy = prod(size(XC));
Stri=sparse(TRI(:,[1 1 2 2 3 3]),TRI(:,[2 3 1 3 1 2]),1,nxy,nxy);
%
mytri.XC=XC; mytri.YC=YC; mytri.TRI=TRI; mytri.Stri=Stri;
%usage: ik=dsearch(mytri.XC,mytri.YC,mytri.TRI,lon,lat,mytri.Stri);
%       where lon and lat are vector of position

elseif nargin==2;
%compute grid point vector associated with lon/lat vectors
lon=varargin{1}; lat=varargin{2};
ik=dsearch(mytri.XC,mytri.YC,mytri.TRI,lon,lat,mytri.Stri);
if nargout==1; 
  varargout={ik};
elseif nargout==2;
  jj=ceil(ik/size(mytri.XC,1));
  ii=ik-(jj-1)*size(mytri.XC,1);
  varargout={ii}; varargout(2)={jj};
else;
  error('wrong output choice');
end;
elseif nargin==3;
%do the bin average (if nargout==1) or the bin sum+count (if nargout==2)
lon=varargin{1}; lat=varargin{2}; obs=varargin{3};
ii=find(~isnan(obs)); lon=lon(ii); lat=lat(ii); obs=obs(ii);
ik=dsearch(mytri.XC,mytri.YC,mytri.TRI,lon,lat,mytri.Stri);

OBS=0*mytri.XC; NOBS=OBS;
for k=1:length(ik)
        NOBS(ik(k))=NOBS(ik(k))+1;
        OBS(ik(k))=OBS(ik(k))+obs(k);
end  % k=1:length(ik)

if nargout==1;%output bin average
  in=find(NOBS); OBS(in)=OBS(in)./NOBS(in);
  in=find(~NOBS); OBS(in)=NaN; NOBS(in)=NaN;
  varargout={convert2array(OBS,mygrid.XC)};
elseif nargout==2;%output bin sum+count
  OBS=convert2array(OBS,mygrid.XC);
  NOBS=convert2array(NOBS,mygrid.XC);
  varargout={OBS}; varargout(2)={NOBS};
else;
  error('wrong output choice');
end;

else;
error('wrong input choice');
end;


1	gforget	1.1	function [varargout]=gcmfaces_bindata(varargin);
2
3			global mygrid mytri;
4
5			if nargin==0;
6			%generate delaunay triangulation:
7			XC=convert2array(mygrid.XC);
8			YC=convert2array(mygrid.YC);
9			%needed so that we do not loose data
10			XC(isnan(XC))=270;
11			YC(isnan(YC))=180;
12			%
13			TRI=delaunay(XC,YC); nxy = prod(size(XC));
14			Stri=sparse(TRI(:,[1 1 2 2 3 3]),TRI(:,[2 3 1 3 1 2]),1,nxy,nxy);
15			%
16			mytri.XC=XC; mytri.YC=YC; mytri.TRI=TRI; mytri.Stri=Stri;
17			%usage: ik=dsearch(mytri.XC,mytri.YC,mytri.TRI,lon,lat,mytri.Stri);
18			% where lon and lat are vector of position
19
20			elseif nargin==2;
21			%compute grid point vector associated with lon/lat vectors
22			lon=varargin{1}; lat=varargin{2};
23			ik=dsearch(mytri.XC,mytri.YC,mytri.TRI,lon,lat,mytri.Stri);
24			if nargout==1;
25			varargout={ik};
26			elseif nargout==2;
27			jj=ceil(ik/size(mytri.XC,1));
28			ii=ik-(jj-1)*size(mytri.XC,1);
29			varargout={ii}; varargout(2)={jj};
30			else;
31			error('wrong output choice');
32			end;
33			elseif nargin==3;
34			%do the bin average (if nargout==1) or the bin sum+count (if nargout==2)
35			lon=varargin{1}; lat=varargin{2}; obs=varargin{3};
36			ii=find(~isnan(obs)); lon=lon(ii); lat=lat(ii); obs=obs(ii);
37			ik=dsearch(mytri.XC,mytri.YC,mytri.TRI,lon,lat,mytri.Stri);
38
39			OBS=0*mytri.XC; NOBS=OBS;
40			for k=1:length(ik)
41			NOBS(ik(k))=NOBS(ik(k))+1;
42			OBS(ik(k))=OBS(ik(k))+obs(k);
43			end % k=1:length(ik)
44
45			if nargout==1;%output bin average
46			in=find(NOBS); OBS(in)=OBS(in)./NOBS(in);
47			in=find(~NOBS); OBS(in)=NaN; NOBS(in)=NaN;
48			varargout={convert2array(OBS,mygrid.XC)};
49			elseif nargout==2;%output bin sum+count
50			OBS=convert2array(OBS,mygrid.XC);
51			NOBS=convert2array(NOBS,mygrid.XC);
52			varargout={OBS}; varargout(2)={NOBS};
53			else;
54			error('wrong output choice');
55			end;
56
57			else;
58			error('wrong input choice');
59			end;
60
61
62
63