MITgcm_contrib/gmaze_pv/test_volbet2iso.m

% Test of the function volbet2iso
%

clear

% Theoritical fields:
eg = 1;

switch eg
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
 case 1 % The more simple:
  % Axis:
  lon = [200:1/8:300]; nlon = length(lon);
  lat = [0:1/8:20];   nlat = length(lat);
  dpt = [5:5:1000];    ndpt = length(dpt);
    
  % chp goes linearly from 10 at 30N to 0 at 40N
  % uniformely between the surface and the bottom:
  [a chp c] = meshgrid(lon,-lat+lat(nlat),dpt); clear a c
  chp = permute(chp,[3 1 2]);
  %chp(:,:,1:400) = chp(:,:,1:400).*NaN;
  
  % Define limits:
  LIMITS(1) = 18 ; % Between 1.75N and 2N
  LIMITS(2) = 18.2 ;
  LIMITS(3) = dpt(ndpt) ;
  LIMITS(4:5) = lat([1 nlat]) ;
  LIMITS(6:7) = lon([1 nlon]) ;
   
  % Expected volume: 
  dx = m_lldist([200 300],[1 1]*1.875)./1000;
  dy = m_lldist([1 1],[1.75 2])./1000;
  dz = dpt(ndpt)./1000;
  Vexp = dx*dy*dz; % Unit is km^3
  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
  
end %switch


% Get volume:
[V Vmat dV] = volbet2iso(chp,LIMITS,dpt,lat,lon);

disp('Computed:')
disp(num2str(V/1000^3))
disp('Approximatly expected:')
disp(num2str(Vexp))
1	gmaze	1.1	% Test of the function volbet2iso
2			%
3
4			clear
5
6			% Theoritical fields:
7			eg = 1;
8
9			switch eg
10			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
11			case 1 % The more simple:
12			% Axis:
13			lon = [200:1/8:300]; nlon = length(lon);
14			lat = [0:1/8:20]; nlat = length(lat);
15			dpt = [5:5:1000]; ndpt = length(dpt);
16
17			% chp goes linearly from 10 at 30N to 0 at 40N
18			% uniformely between the surface and the bottom:
19			[a chp c] = meshgrid(lon,-lat+lat(nlat),dpt); clear a c
20			chp = permute(chp,[3 1 2]);
21			%chp(:,:,1:400) = chp(:,:,1:400).*NaN;
22
23			% Define limits:
24			LIMITS(1) = 18 ; % Between 1.75N and 2N
25			LIMITS(2) = 18.2 ;
26			LIMITS(3) = dpt(ndpt) ;
27			LIMITS(4:5) = lat([1 nlat]) ;
28			LIMITS(6:7) = lon([1 nlon]) ;
29
30			% Expected volume:
31			dx = m_lldist([200 300],[1 1]*1.875)./1000;
32			dy = m_lldist([1 1],[1.75 2])./1000;
33			dz = dpt(ndpt)./1000;
34			Vexp = dxdydz; % Unit is km^3
35
36			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
37
38			end %switch
39
40
41
42			% Get volume:
43			[V Vmat dV] = volbet2iso(chp,LIMITS,dpt,lat,lon);
44
45			disp('Computed:')
46			disp(num2str(V/1000^3))
47			disp('Approximatly expected:')
48			disp(num2str(Vexp))