global_ocean.128x64x15/diags_matlab/mit_oceanmasks.m

function grid = mit_oceanmasks(grid,plotit)
% function grid = mit_oceanmasks(grid)
% create masks for the 4x4 degree run for the individual oceans
% hfacc/s/w have to be available for this
  
  if grid.nx ~= 90 | grid.ny ~= 40 
    error(sprintf( ...
        '%s\n%s\n%s\n%s', ...
        'Your horizontal grid has dimensions different from 90x40.', ...
        'mit_oceanmasks creates masks for the individual ocean basins,', ...
        'therefore it needs to be customized if you change the', ...
        'horizontal dimensions.'))
  end
    
  % masks for the zonal velocities:
  % hfacw
  
  % Atlantic mask
  grid.atlantic_hfacw = grid.hfacw;
  ixsow = 6:73; % southern ocean and most of pacific
  iysow = 1:24; 
  grid.atlantic_hfacw(ixsow,iysow,:) = NaN;
  ixnpw = 15:67; % northern pacific
  iynpw = 25:37;
  grid.atlantic_hfacw(ixnpw,iynpw,:) = NaN;
  % Pacific mask
  grid.pacific_hfacw = grid.hfacw;
  ixaw = 74:90; % most of the atlantic
  iyaw = 1:40;
  grid.pacific_hfacw(ixaw,iyaw,:) = NaN;
  ixnaw = 68:73; % north atlantic and gulf of Mexico
  iynaw = 27:40;
  grid.pacific_hfacw(ixnaw,iynaw,:) = NaN;
  ixcaw = 70:73; % carribean
  iycaw = 24:26;
  grid.pacific_hfacw(ixcaw,iycaw,:) = NaN;
  ixiw = 1:27; % most of the indic
  iyiw = 1:40;
  grid.pacific_hfacw(ixiw,iyiw,:) = NaN;
  ixauw = 28:35; % water south and west of Australia
  iyauw = 1:20;
  grid.pacific_hfacw(ixauw,iyauw,:) = NaN;
  
  
  % masks for the meridional overturning stream functions:
  % hfacs
  
  % Atlantic mask
  grid.atlantic_hfacs = grid.hfacs;
  ixso = 5:73; % southern ocean
  iyso = 1:24;
  grid.atlantic_hfacs(ixso,iyso,:) = NaN;
  ixnp = 14:65; % northern pacific
  iynp = 25:36; 
  grid.atlantic_hfacs(ixnp,iynp,:) = NaN;
  ixsoc = 1:90; % complete southern ocean
  iysoc = 1:12; % south of latg(12) = -34 degN;
  grid.atlantic_hfacs(ixsoc,iysoc,:) = NaN;
  
  % Pacific mask
  grid.pacific_hfacs = grid.hfacs;
  ixatl = [1:26 74:90]; % most of altantic and indian ocean
  iyatl = 1:grid.ny;
  grid.pacific_hfacs(ixatl,iyatl,:) = NaN;
  ixatlf = 73; % small fix
  iyatlf = 31;
  grid.pacific_hfacs(ixatlf,iyatlf,:) = NaN;
  ixgm = 66:73; % Gulf of Mexico
  iygm = 25:30;
  grid.pacific_hfacs(ixgm,iygm,:) = NaN;
  ixhb = 68:73; % Hudson Bay and Baffin Bay
  iyhb = 36:40;
  grid.pacific_hfacs(ixhb,iyhb,:) = NaN;
  ixsop = 1:grid.nx; % southern ocean
  iysop = 1:12;
  grid.pacific_hfacs(ixsop,iysop,:) = NaN;
  ixio = 27:32; % rest of indian ocean
  iyio = 13:19;
  grid.pacific_hfacs(ixio,iyio,:) = NaN;
  %pacific_hfacs(33:34,12,:) = NaN; % two points south of Australia
  
  % Indian ocean
  grid.indic_hfacs = grid.hfacs;
  ixsoi = 1:90; % southern ocean
  iysoi = 1:12;
  grid.indic_hfacs(ixsoi,iysoi,:) = NaN;
  ixap = [1:4 33:90]; % atlantic and pacific
  iyap = 1:40;
  grid.indic_hfacs(ixap,iyap,:) = NaN;
  iynps = 35:40; % north polar sea
  grid.indic_hfacs(:,iynps,:) = NaN;
  ixmeds = 1:10;  % mediterraenan
  iymeds = 29:31;
  grid.indic_hfacs(ixmeds,iymeds,:) = NaN;
  ixscs = 27:32; % south china sea
  iyscs = 21:29;
  grid.indic_hfacs(ixscs,iyscs,:) = NaN;
  grid.indic_hfacs(5,39,:) = NaN; % singular point near Murmansk
  
  % Southern Ocean
  grid.so_hfacs = change( ...
      change(grid.atlantic_hfacs,'==',NaN,-1) ...
      + change(grid.pacific_hfacs,'==',NaN,-1) ...
      + change(grid.indic_hfacs,'==',NaN,-1) ...
      + grid.hfacs-grid.hfacs,'==',-1,NaN);
  iso = find(~isnan(grid.so_hfacs));
  grid.so_hfacs(iso)= grid.hfacs(iso);

  if exist('plotit','var')
    figure;
    spy(isnan(grid.hfacs(:,:,1)'),'kx');
    hold on; 
    spy(~isnan(grid.atlantic_hfacs(:,:,1)')); 
    spy(~isnan(grid.indic_hfacs(:,:,1)'),'g.'); 
    spy(~isnan(grid.so_hfacs(:,:,1)'),'r.'); 
    spy(~isnan(grid.pacific_hfacs(:,:,1)'),'y.');
    title('hfacs')
    axis xy
  end
  
  % masks for zonal averages of C-properties:
  % hfacc
  
  % Atlantic mask
  grid.atlantic_hfacc = grid.hfacc;
  ixe = 6:65; % indian ocean
  iye = 1:28;
  grid.atlantic_hfacc(ixe,iye,:) = NaN;
  ixep = 26:65; % eastern pacific
  iyep = 1:grid.ny; 
  grid.atlantic_hfacc(ixep,iyep,:) = NaN;
  ixw = 6:73; % western pacific
  iyw = 1:23;
  grid.atlantic_hfacc(ixw,iyw,:) = NaN;
  ixgl = 66:68; % galapagos
  iygl = 24;
  grid.atlantic_hfacc(ixgl,iygl,:) = NaN;
  ixch = 73; % coast of Chile
  iych = 13:15;
  grid.atlantic_hfacc(ixch,iych,:) = NaN;
  
  % Pacific mask
  grid.pacific_hfacc = grid.hfacc;
  ixatl = [1:25 74:90]; % most of altantic and indian ocean
  iyatl = 1:grid.ny;
  grid.pacific_hfacc(ixatl,iyatl,:) = NaN;
  ixatlf = 73; % small fix
  iyatlf = 31;
  grid.pacific_hfacc(ixatlf,iyatlf,:) = NaN;
  ixgme = 69:73; % Gulf of Mexico
  iygme = 24:30;
  ixgmw = 66:69;
  iygmw = 26:28;
  grid.pacific_hfacc(ixgme,iygme,:) = NaN;
  grid.pacific_hfacc(ixgmw,iygmw,:) = NaN;
  ixhbb = 68:73; % Hudson and Baffin Bay
  iyhbb = 35:40;
  grid.pacific_hfacc(ixhbb,iyhbb,:) = NaN;
  ixio = 26:34; % rest of indian ocean
  iyio =  1:19;
  grid.pacific_hfacc(ixio,iyio,:) = NaN;
%  grid.pacific_hfacc(33:34,12,:) = NaN; % two points south of Australia
  
  % Indian ocean
  grid.indic_hfacc = ...
      change(change(grid.atlantic_hfacc,'==',NaN,0),'>',0,NaN) ...
      + change(change(grid.pacific_hfacc,'==',NaN,0),'>',0,NaN) ...
      + grid.hfacc;
  
  if exist('plotit','var')
    figure;
    spy(isnan(grid.hfacc(:,:,1)'),'kx'); axis xy
    hold on; 
    spy(~isnan(grid.atlantic_hfacc(:,:,1)')); 
    spy(~isnan(grid.pacific_hfacc(:,:,1)'),'r.');
    spy(~isnan(grid.indic_hfacc(:,:,1)'),'g.'); 
    axis xy
    title('hfacc');

  end
  
  return
1	cnh	1.1	function grid = mit_oceanmasks(grid,plotit)
2			% function grid = mit_oceanmasks(grid)
3			% create masks for the 4x4 degree run for the individual oceans
4			% hfacc/s/w have to be available for this
5
6			if grid.nx ~= 90 \| grid.ny ~= 40
7			error(sprintf( ...
8			'%s\n%s\n%s\n%s', ...
9			'Your horizontal grid has dimensions different from 90x40.', ...
10			'mit_oceanmasks creates masks for the individual ocean basins,', ...
11			'therefore it needs to be customized if you change the', ...
12			'horizontal dimensions.'))
13			end
14
15			% masks for the zonal velocities:
16			% hfacw
17
18			% Atlantic mask
19			grid.atlantic_hfacw = grid.hfacw;
20			ixsow = 6:73; % southern ocean and most of pacific
21			iysow = 1:24;
22			grid.atlantic_hfacw(ixsow,iysow,:) = NaN;
23			ixnpw = 15:67; % northern pacific
24			iynpw = 25:37;
25			grid.atlantic_hfacw(ixnpw,iynpw,:) = NaN;
26			% Pacific mask
27			grid.pacific_hfacw = grid.hfacw;
28			ixaw = 74:90; % most of the atlantic
29			iyaw = 1:40;
30			grid.pacific_hfacw(ixaw,iyaw,:) = NaN;
31			ixnaw = 68:73; % north atlantic and gulf of Mexico
32			iynaw = 27:40;
33			grid.pacific_hfacw(ixnaw,iynaw,:) = NaN;
34			ixcaw = 70:73; % carribean
35			iycaw = 24:26;
36			grid.pacific_hfacw(ixcaw,iycaw,:) = NaN;
37			ixiw = 1:27; % most of the indic
38			iyiw = 1:40;
39			grid.pacific_hfacw(ixiw,iyiw,:) = NaN;
40			ixauw = 28:35; % water south and west of Australia
41			iyauw = 1:20;
42			grid.pacific_hfacw(ixauw,iyauw,:) = NaN;
43
44
45			% masks for the meridional overturning stream functions:
46			% hfacs
47
48			% Atlantic mask
49			grid.atlantic_hfacs = grid.hfacs;
50			ixso = 5:73; % southern ocean
51			iyso = 1:24;
52			grid.atlantic_hfacs(ixso,iyso,:) = NaN;
53			ixnp = 14:65; % northern pacific
54			iynp = 25:36;
55			grid.atlantic_hfacs(ixnp,iynp,:) = NaN;
56			ixsoc = 1:90; % complete southern ocean
57			iysoc = 1:12; % south of latg(12) = -34 degN;
58			grid.atlantic_hfacs(ixsoc,iysoc,:) = NaN;
59
60			% Pacific mask
61			grid.pacific_hfacs = grid.hfacs;
62			ixatl = [1:26 74:90]; % most of altantic and indian ocean
63			iyatl = 1:grid.ny;
64			grid.pacific_hfacs(ixatl,iyatl,:) = NaN;
65			ixatlf = 73; % small fix
66			iyatlf = 31;
67			grid.pacific_hfacs(ixatlf,iyatlf,:) = NaN;
68			ixgm = 66:73; % Gulf of Mexico
69			iygm = 25:30;
70			grid.pacific_hfacs(ixgm,iygm,:) = NaN;
71			ixhb = 68:73; % Hudson Bay and Baffin Bay
72			iyhb = 36:40;
73			grid.pacific_hfacs(ixhb,iyhb,:) = NaN;
74			ixsop = 1:grid.nx; % southern ocean
75			iysop = 1:12;
76			grid.pacific_hfacs(ixsop,iysop,:) = NaN;
77			ixio = 27:32; % rest of indian ocean
78			iyio = 13:19;
79			grid.pacific_hfacs(ixio,iyio,:) = NaN;
80			%pacific_hfacs(33:34,12,:) = NaN; % two points south of Australia
81
82			% Indian ocean
83			grid.indic_hfacs = grid.hfacs;
84			ixsoi = 1:90; % southern ocean
85			iysoi = 1:12;
86			grid.indic_hfacs(ixsoi,iysoi,:) = NaN;
87			ixap = [1:4 33:90]; % atlantic and pacific
88			iyap = 1:40;
89			grid.indic_hfacs(ixap,iyap,:) = NaN;
90			iynps = 35:40; % north polar sea
91			grid.indic_hfacs(:,iynps,:) = NaN;
92			ixmeds = 1:10; % mediterraenan
93			iymeds = 29:31;
94			grid.indic_hfacs(ixmeds,iymeds,:) = NaN;
95			ixscs = 27:32; % south china sea
96			iyscs = 21:29;
97			grid.indic_hfacs(ixscs,iyscs,:) = NaN;
98			grid.indic_hfacs(5,39,:) = NaN; % singular point near Murmansk
99
100			% Southern Ocean
101			grid.so_hfacs = change( ...
102			change(grid.atlantic_hfacs,'==',NaN,-1) ...
103			+ change(grid.pacific_hfacs,'==',NaN,-1) ...
104			+ change(grid.indic_hfacs,'==',NaN,-1) ...
105			+ grid.hfacs-grid.hfacs,'==',-1,NaN);
106			iso = find(~isnan(grid.so_hfacs));
107			grid.so_hfacs(iso)= grid.hfacs(iso);
108
109			if exist('plotit','var')
110			figure;
111			spy(isnan(grid.hfacs(:,:,1)'),'kx');
112			hold on;
113			spy(~isnan(grid.atlantic_hfacs(:,:,1)'));
114			spy(~isnan(grid.indic_hfacs(:,:,1)'),'g.');
115			spy(~isnan(grid.so_hfacs(:,:,1)'),'r.');
116			spy(~isnan(grid.pacific_hfacs(:,:,1)'),'y.');
117			title('hfacs')
118			axis xy
119			end
120
121			% masks for zonal averages of C-properties:
122			% hfacc
123
124			% Atlantic mask
125			grid.atlantic_hfacc = grid.hfacc;
126			ixe = 6:65; % indian ocean
127			iye = 1:28;
128			grid.atlantic_hfacc(ixe,iye,:) = NaN;
129			ixep = 26:65; % eastern pacific
130			iyep = 1:grid.ny;
131			grid.atlantic_hfacc(ixep,iyep,:) = NaN;
132			ixw = 6:73; % western pacific
133			iyw = 1:23;
134			grid.atlantic_hfacc(ixw,iyw,:) = NaN;
135			ixgl = 66:68; % galapagos
136			iygl = 24;
137			grid.atlantic_hfacc(ixgl,iygl,:) = NaN;
138			ixch = 73; % coast of Chile
139			iych = 13:15;
140			grid.atlantic_hfacc(ixch,iych,:) = NaN;
141
142			% Pacific mask
143			grid.pacific_hfacc = grid.hfacc;
144			ixatl = [1:25 74:90]; % most of altantic and indian ocean
145			iyatl = 1:grid.ny;
146			grid.pacific_hfacc(ixatl,iyatl,:) = NaN;
147			ixatlf = 73; % small fix
148			iyatlf = 31;
149			grid.pacific_hfacc(ixatlf,iyatlf,:) = NaN;
150			ixgme = 69:73; % Gulf of Mexico
151			iygme = 24:30;
152			ixgmw = 66:69;
153			iygmw = 26:28;
154			grid.pacific_hfacc(ixgme,iygme,:) = NaN;
155			grid.pacific_hfacc(ixgmw,iygmw,:) = NaN;
156			ixhbb = 68:73; % Hudson and Baffin Bay
157			iyhbb = 35:40;
158			grid.pacific_hfacc(ixhbb,iyhbb,:) = NaN;
159			ixio = 26:34; % rest of indian ocean
160			iyio = 1:19;
161			grid.pacific_hfacc(ixio,iyio,:) = NaN;
162			% grid.pacific_hfacc(33:34,12,:) = NaN; % two points south of Australia
163
164			% Indian ocean
165			grid.indic_hfacc = ...
166			change(change(grid.atlantic_hfacc,'==',NaN,0),'>',0,NaN) ...
167			+ change(change(grid.pacific_hfacc,'==',NaN,0),'>',0,NaN) ...
168			+ grid.hfacc;
169
170			if exist('plotit','var')
171			figure;
172			spy(isnan(grid.hfacc(:,:,1)'),'kx'); axis xy
173			hold on;
174			spy(~isnan(grid.atlantic_hfacc(:,:,1)'));
175			spy(~isnan(grid.pacific_hfacc(:,:,1)'),'r.');
176			spy(~isnan(grid.indic_hfacc(:,:,1)'),'g.');
177			axis xy
178			title('hfacc');
179
180			end
181
182			return