MITgcm_contrib/gmaze_pv/intbet2outcrops.m

%
% I  = intbet2outcrops(TRACER,LIMITS,LAT,LONG)
%
% This function computes the horizontal surface integral between two 
% outcrops of the TRACER field, given fixed limits eastward, westward 
% and southward.
%
% TRACER = TRACER(LAT,LONG) : surface tracer variable in 2D
% LIMITS = [OUTCROP1 OUTCROP2 MAX_LAT1 MAX_LAT2 MAX_LONG1 MAX_LONG2]
%          : limit's values (MAX_LAT2 is used only if
%            the outcrop's surfaces reach them).
% LAT                             : latitude axis (1D), degrees northward
% LONG                            : longitude axis (1D), degrees east
% I                               : single surface integral value
%
% 06/15/2006
% gmaze@mit.edu
% 


function varargout = intbet2outcrops(TRACER,LIMITS,LAT,LONG)


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PRE-PROCESS and ERROR CHECK %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Check number of input:
if nargin ~= 4
  help intbet2outcrops.m
  error('intbet2outcrops.m : Wrong number of parameters')
  return
end %if
  
% Check dimensions:
n = size(TRACER);
if length(n)==2
  [ny nx] = size(TRACER);
  if ny~=length(LAT) | nx~=length(LONG)
     help intbet2outcrops.m
     error('intbet2outcrops.m : Axis must have same dimensions than TRACER field');
     return
  end %if
else
  help intbet2outcrops.m
  error('intbet2outcrops.m : TRACER must be a 2D field')
  return
end %if

% Ensure that axis are of dim: (1,N) and well sorted (increasing values):
a=size(LAT);
if a(1) ~= 1,  LAT=LAT'; end 
S = sort(LAT);
if S ~= LAT
  help intbet2outcrops.m
  error('intbet2outcrops.m : LAT must be increasing values')
  return
end %if
a=size(LONG);
if a(1) ~= 1,  LONG=LONG'; end 
S = sort(LONG);
if S ~= LONG
  help intbet2outcrops.m
  error('intbet2outcrops.m : LONG must be increasing values')
  return
end %if

% LIMITS definition:
if length(LIMITS) ~= 6
  help intbet2outcrops.m
  error('intbet2outcrops.m : LIMITS must contains 6 values')
  return
end %if
OUTCROPS = sort( LIMITS(1:2) );
LAT_MAX  = sort( LIMITS(3:4) );
LONG_MAX = sort( LIMITS(5:6) ); 


%%%%%%%%%%%%%%%%%%%%
% COMPUTE INTEGRAL %
%%%%%%%%%%%%%%%%%%%%
% We first determine the element surface matrix and points to integrate:
[I1 I1mat dI1] = subfct_getsurf(TRACER,LAT,LONG,[OUTCROPS(1) LAT_MAX LONG_MAX]);
[I2 I2mat dI2] = subfct_getsurf(TRACER,LAT,LONG,[OUTCROPS(2) LAT_MAX LONG_MAX]);

% Then we determine the outcrop surface limits:
I1mat(find(I1mat==0))=NaN;
I2mat(find(I2mat==0))=NaN;
I1mat(find(I1mat==1))=0;
I2mat(find(I2mat==1))=0;
Imat = I1mat + I2mat;
Imat(find(Imat==0))=1;
Imat(isnan(Imat))=0;
Imat = logical(Imat);

% And the integral of the TRACER on it:
I = sum(TRACER(Imat).*dI1(Imat));


%%%%%%%%%%%
% OUTPUTS %
%%%%%%%%%%%
switch nargout
 case 1
  varargout(1) = {I};
 case 2
  varargout(1) = {I};
  varargout(2) = {Imat};
 case 3
  varargout(1) = {I};
  varargout(2) = {Imat};
  varargout(3) = {dI1};
end %switch nargout


1	gmaze	1.1	%
2			% I = intbet2outcrops(TRACER,LIMITS,LAT,LONG)
3			%
4			% This function computes the horizontal surface integral between two
5			% outcrops of the TRACER field, given fixed limits eastward, westward
6			% and southward.
7			%
8			% TRACER = TRACER(LAT,LONG) : surface tracer variable in 2D
9			% LIMITS = [OUTCROP1 OUTCROP2 MAX_LAT1 MAX_LAT2 MAX_LONG1 MAX_LONG2]
10			% : limit's values (MAX_LAT2 is used only if
11			% the outcrop's surfaces reach them).
12			% LAT : latitude axis (1D), degrees northward
13			% LONG : longitude axis (1D), degrees east
14			% I : single surface integral value
15			%
16			% 06/15/2006
17			% gmaze@mit.edu
18			%
19
20
21			function varargout = intbet2outcrops(TRACER,LIMITS,LAT,LONG)
22
23
24			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
25			% PRE-PROCESS and ERROR CHECK %
26			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
27
28			% Check number of input:
29			if nargin ~= 4
30			help intbet2outcrops.m
31			error('intbet2outcrops.m : Wrong number of parameters')
32			return
33			end %if
34
35			% Check dimensions:
36			n = size(TRACER);
37			if length(n)==2
38			[ny nx] = size(TRACER);
39			if ny~=length(LAT) \| nx~=length(LONG)
40			help intbet2outcrops.m
41			error('intbet2outcrops.m : Axis must have same dimensions than TRACER field');
42			return
43			end %if
44			else
45			help intbet2outcrops.m
46			error('intbet2outcrops.m : TRACER must be a 2D field')
47			return
48			end %if
49
50			% Ensure that axis are of dim: (1,N) and well sorted (increasing values):
51			a=size(LAT);
52			if a(1) ~= 1, LAT=LAT'; end
53			S = sort(LAT);
54			if S ~= LAT
55			help intbet2outcrops.m
56			error('intbet2outcrops.m : LAT must be increasing values')
57			return
58			end %if
59			a=size(LONG);
60			if a(1) ~= 1, LONG=LONG'; end
61			S = sort(LONG);
62			if S ~= LONG
63			help intbet2outcrops.m
64			error('intbet2outcrops.m : LONG must be increasing values')
65			return
66			end %if
67
68			% LIMITS definition:
69			if length(LIMITS) ~= 6
70			help intbet2outcrops.m
71			error('intbet2outcrops.m : LIMITS must contains 6 values')
72			return
73			end %if
74			OUTCROPS = sort( LIMITS(1:2) );
75			LAT_MAX = sort( LIMITS(3:4) );
76			LONG_MAX = sort( LIMITS(5:6) );
77
78
79
80			%%%%%%%%%%%%%%%%%%%%
81			% COMPUTE INTEGRAL %
82			%%%%%%%%%%%%%%%%%%%%
83			% We first determine the element surface matrix and points to integrate:
84			[I1 I1mat dI1] = subfct_getsurf(TRACER,LAT,LONG,[OUTCROPS(1) LAT_MAX LONG_MAX]);
85			[I2 I2mat dI2] = subfct_getsurf(TRACER,LAT,LONG,[OUTCROPS(2) LAT_MAX LONG_MAX]);
86
87			% Then we determine the outcrop surface limits:
88			I1mat(find(I1mat==0))=NaN;
89			I2mat(find(I2mat==0))=NaN;
90			I1mat(find(I1mat==1))=0;
91			I2mat(find(I2mat==1))=0;
92			Imat = I1mat + I2mat;
93			Imat(find(Imat==0))=1;
94			Imat(isnan(Imat))=0;
95			Imat = logical(Imat);
96
97			% And the integral of the TRACER on it:
98			I = sum(TRACER(Imat).*dI1(Imat));
99
100
101
102
103			%%%%%%%%%%%
104			% OUTPUTS %
105			%%%%%%%%%%%
106			switch nargout
107			case 1
108			varargout(1) = {I};
109			case 2
110			varargout(1) = {I};
111			varargout(2) = {Imat};
112			case 3
113			varargout(1) = {I};
114			varargout(2) = {Imat};
115			varargout(3) = {dI1};
116			end %switch nargout
117
118
119