1 |
gforget |
1.1 |
|
2 |
|
|
if userStep==1;%diags to be computed |
3 |
gforget |
1.2 |
%here _s stands for cumulative sum and _2s for cumulative sum of squares |
4 |
|
|
%the _s should be stated first |
5 |
|
|
listDiags=['THETA_s SALT_s SIatmQnt_s SIatmFW_s oceQnet_s oceFWflx_s ' ... |
6 |
|
|
'fldTZ_s fldTM_s curlTau_s fldETAN_s fldETANLEADS_s fldMLD_s ' ... |
7 |
|
|
'THETA_2s SALT_2s SIatmQnt_2s SIatmFW_2s oceQnet_2s oceFWflx_2s ' ... |
8 |
|
|
'fldTZ_2s fldTM_2s curlTau_2s fldETAN_2s fldETANLEADS_2s fldMLD_2s']; |
9 |
gforget |
1.4 |
|
10 |
gforget |
1.1 |
elseif userStep==2;%input files and variables |
11 |
|
|
listFlds={ 'THETA','SALT','SIatmQnt','SIatmFW ','oceQnet ','oceFWflx','oceTAUX','oceTAUY','ETAN','sIceLoad','MXLDEPTH'}; |
12 |
|
|
listFldsNames=deblank(listFlds); |
13 |
|
|
listFiles={'state_3d_set1','state_2d_set1','other_2d_set1'}; |
14 |
|
|
listSubdirs={[dirModel 'diags/OTHER/' ],[dirModel 'diags/STATE/' ]}; |
15 |
gforget |
1.4 |
|
16 |
|
|
elseif userStep==3;%computation |
17 |
gforget |
1.1 |
%mask fields: |
18 |
|
|
SIatmQnt=SIatmQnt.*mygrid.mskC(:,:,1); |
19 |
|
|
SIatmFW=SIatmFW.*mygrid.mskC(:,:,1); |
20 |
|
|
oceQnet=oceQnet.*mygrid.mskC(:,:,1); |
21 |
|
|
oceFWflx=oceFWflx.*mygrid.mskC(:,:,1); |
22 |
|
|
fldTX=oceTAUX.*mygrid.mskW(:,:,1); |
23 |
|
|
fldTY=oceTAUY.*mygrid.mskS(:,:,1); |
24 |
|
|
%compute Eastward/Northward wind stresses: |
25 |
|
|
[fldTZ,fldTM]=calc_UEVNfromUXVY(fldTX,fldTY); |
26 |
|
|
%compute wind stress curl: |
27 |
|
|
curlTau=calc_UV_curl(fldTX, fldTY,1 );%the doMask argument should not matter as msk was already applied |
28 |
|
|
%mask and re-arrange fields: |
29 |
|
|
fldETAN=ETAN.*mygrid.mskC(:,:,1); |
30 |
|
|
fldETANLEADS=(ETAN+sIceLoad/myparms.rhoconst).*mygrid.mskC(:,:,1); |
31 |
|
|
fldMLD=MXLDEPTH.*mygrid.mskC(:,:,1); |
32 |
|
|
% |
33 |
|
|
THETA=THETA.*mygrid.mskC; |
34 |
|
|
SALT=SALT.*mygrid.mskC; |
35 |
gforget |
1.2 |
% |
36 |
|
|
if ii>1; |
37 |
|
|
fileMatPrev=['diags_set_' tmp1 '_' num2str(listTimes(ii-1)) '.mat']; |
38 |
gforget |
1.4 |
listTimesBak=listTimes; |
39 |
gforget |
1.2 |
load([dirMat fileMatPrev]); |
40 |
gforget |
1.4 |
listTimes=listTimesBak; |
41 |
gforget |
1.2 |
end; |
42 |
|
|
for jj=1:length(listDiags)/2; |
43 |
|
|
myDiag=listDiags{jj}(1:end-2); |
44 |
|
|
if ii==1; |
45 |
|
|
eval([myDiag '_s=0*' myDiag ';']); |
46 |
|
|
eval([myDiag '_2s=0*' myDiag '.^2;']); |
47 |
|
|
end; |
48 |
|
|
eval([myDiag '_s=' myDiag '_s+' myDiag ';']); |
49 |
|
|
eval([myDiag '_2s=' myDiag '_2s+' myDiag '.^2;']); |
50 |
|
|
end; |
51 |
|
|
|
52 |
gforget |
1.4 |
%===================== COMPUTATIONAL SEQUENCE ENDS =========================% |
53 |
|
|
%===================== PLOTTING SEQUENCE BEGINS =========================% |
54 |
|
|
|
55 |
|
|
elseif userStep==0;%loading / post-processing of mat files |
56 |
|
|
|
57 |
|
|
%load last cumsum |
58 |
|
|
alldiag=alldiag_load(dirMat,[fileMat '_*.mat'],'',myparms.recInAve(2)); |
59 |
|
|
%load first cumsum |
60 |
|
|
if myparms.recInAve(1)>1; |
61 |
|
|
tmpdiag=alldiag_load(dirMat,[fileMat '_*.mat'],'',myparms.recInAve(1)-1); |
62 |
|
|
for ii=1:length(alldiag.listDiags); |
63 |
|
|
tmp0=alldiag.listDiags{ii}; |
64 |
|
|
if ~strcmp(tmp0,'listTimes')&~strcmp(tmp0,'listSteps'); |
65 |
|
|
tmp1=getfield(alldiag,alldiag.listDiags{ii}); |
66 |
|
|
tmp2=getfield(tmpdiag,alldiag.listDiags{ii}); |
67 |
|
|
alldiag=setfield(alldiag,tmp0,tmp1-tmp2); |
68 |
|
|
end; |
69 |
|
|
end; |
70 |
|
|
end; |
71 |
|
|
%accomodate missing fields (old MITgcm version) |
72 |
|
|
listFlds={'SIatmQnt_s','oceQnet_s','SIatmFW_s','oceFWflx_s',... |
73 |
|
|
'fldTZ_s','fldTM_s','curlTau_s','fldETAN_s','fldETANLEADS_s',... |
74 |
|
|
'SIatmQnt_2s','oceQnet_2s','SIatmFW_2s','oceFWflx_2s',... |
75 |
|
|
'fldTZ_2s','fldTM_2s','curlTau_2s','fldETAN_2s','fldETANLEADS_2s'}; |
76 |
|
|
for ii=1:length(listFlds); |
77 |
|
|
if ~sum(strcmp(fieldnames(alldiag),listFlds{ii})); |
78 |
|
|
eval(['alldiag.' listFlds{ii} '=NaN*alldiag.fldETAN_s;']); |
79 |
|
|
end; |
80 |
|
|
end; |
81 |
|
|
% |
82 |
|
|
n=diff(myparms.recInAve)+1; |
83 |
|
|
for ii=1:length(alldiag.listDiags)/2; |
84 |
|
|
tmp0=alldiag.listDiags{ii}; |
85 |
|
|
if ~strcmp(tmp0,'listTimes')&~strcmp(tmp0,'listSteps'); |
86 |
|
|
tmp1=1/n*getfield(alldiag,tmp0); |
87 |
|
|
tmp2=1/n*getfield(alldiag,[tmp0(1:end-1) '2s']); |
88 |
|
|
tmp2=(tmp2-tmp1.^2); |
89 |
|
|
tmp2=n/(n-1)*tmp2; |
90 |
|
|
%tmp2(tmp2<0)=0; |
91 |
|
|
tmp2=sqrt(tmp2); |
92 |
|
|
alldiag=setfield(alldiag,[tmp0(1:end-2) '_mean'],tmp1); |
93 |
|
|
alldiag=setfield(alldiag,[tmp0(1:end-2) '_std'],tmp2); |
94 |
|
|
end; |
95 |
|
|
end; |
96 |
|
|
|
97 |
|
|
diagsWereLoaded=1 |
98 |
|
|
|
99 |
|
|
elseif userStep==-1;%plotting |
100 |
|
|
|
101 |
|
|
if isempty(setDiagsParams); |
102 |
|
|
choicePlot={'all'}; |
103 |
|
|
else; |
104 |
|
|
choicePlot=setDiagsParams; |
105 |
|
|
end; |
106 |
|
|
|
107 |
|
|
if sum(strcmp(choicePlot,'all'))|sum(strcmp(choicePlot,'qnet')); |
108 |
|
|
|
109 |
|
|
%qnet from ocean+ice: |
110 |
|
|
fld=-alldiag.SIatmQnt_mean; |
111 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]; title0='QNET to ocean+ice'; |
112 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*10; end; |
113 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
114 |
|
|
myCaption={myYmeanTxt,'mean -- QNET to ocean+ice (W/m2)'}; |
115 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
116 |
|
|
|
117 |
|
|
%qnet to ocean: |
118 |
|
|
fld=alldiag.oceQnet_mean; |
119 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]; title0='QNET to ocean'; |
120 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*10; end; |
121 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
122 |
|
|
myCaption={myYmeanTxt,'mean -- QNET to ocean (W/m2)'}; |
123 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
124 |
|
|
|
125 |
|
|
if multiTimes; |
126 |
|
|
%qnet from ocean+ice: |
127 |
|
|
fld=alldiag.SIatmQnt_std; |
128 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]; title0='std(QNET to ocean+ice)'; |
129 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*5; end; |
130 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
131 |
|
|
myCaption={myYmeanTxt,' standard deviation -- QNET to ocean+ice (W/m2)'}; |
132 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
133 |
|
|
|
134 |
|
|
%qnet from ocean: |
135 |
|
|
fld=alldiag.oceQnet_std; |
136 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]; title0='std(QNET to ocean)'; |
137 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*5; end; |
138 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
139 |
|
|
myCaption={myYmeanTxt,' standard deviation -- QNET to ocean (W/m2)'}; |
140 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
141 |
|
|
end; |
142 |
|
|
|
143 |
|
|
end; |
144 |
|
|
|
145 |
|
|
|
146 |
|
|
if sum(strcmp(choicePlot,'all'))|sum(strcmp(choicePlot,'fwf')); |
147 |
|
|
|
148 |
|
|
%FW flux from ocean+ice: |
149 |
|
|
fld=-alldiag.SIatmFW_mean/1000;%conversion to m/s |
150 |
|
|
fld=fld*86400*1000;%conversion to mm/day |
151 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]*0.06; title0='EMPMR to ocean+ice'; |
152 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.5; end; |
153 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
154 |
|
|
myCaption={myYmeanTxt,'mean -- EMPMR from ocean+ice (mm/day)'}; |
155 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
156 |
|
|
|
157 |
|
|
%FW flux from ocean: |
158 |
|
|
fld=-alldiag.oceFWflx_mean/1000;%conversion to m/s |
159 |
|
|
fld=fld*86400*1000;%conversion to mm/day |
160 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]*0.06; title0='EMPMR to ocean'; |
161 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.5; end; |
162 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
163 |
|
|
myCaption={myYmeanTxt,'mean -- EMPMR from ocean (mm/day)'}; |
164 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
165 |
|
|
|
166 |
|
|
if multiTimes; |
167 |
|
|
%empmr from ocean+ice: |
168 |
|
|
fld=alldiag.SIatmFW_std*86400;%conversion to mm/day |
169 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]*0.04; title0='std(EMPMR to ocean+ice)'; |
170 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*0.5; end; |
171 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
172 |
|
|
myCaption={myYmeanTxt,' standard deviation -- EMPMR to ocean+ice (W/m2)'}; |
173 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
174 |
|
|
|
175 |
|
|
%empmr from ocean: |
176 |
|
|
fld=alldiag.oceFWflx_std*86400;%conversion to mm/day |
177 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]*0.04; title0='std(EMPMR to ocean)'; |
178 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*0.5; end; |
179 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
180 |
|
|
myCaption={myYmeanTxt,' standard deviation -- EMPMR to ocean (W/m2)'}; |
181 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
182 |
|
|
end; |
183 |
|
|
|
184 |
|
|
end; |
185 |
|
|
|
186 |
|
|
if sum(strcmp(choicePlot,'all'))|sum(strcmp(choicePlot,'tau')); |
187 |
|
|
|
188 |
|
|
%zonal wind stress: |
189 |
|
|
fld=alldiag.fldTZ_mean; |
190 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]/500; title0='zonal wind stress'; |
191 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.01; end; |
192 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
193 |
|
|
myCaption={myYmeanTxt,'mean -- zonal wind stress (N/m2)'}; |
194 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
195 |
|
|
|
196 |
|
|
%meridional wind stress: |
197 |
|
|
fld=alldiag.fldTM_mean; |
198 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]/500; title0='meridional wind stress'; |
199 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.01; end; |
200 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
201 |
|
|
myCaption={myYmeanTxt,'mean -- meridional wind stress (N/m2)'}; |
202 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
203 |
|
|
|
204 |
|
|
fld=alldiag.curlTau_mean; |
205 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]/5e8; title0='wind stress curl'; |
206 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.01; end; |
207 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
208 |
|
|
myCaption={myYmeanTxt,'mean -- wind stress curl (N/m3)'}; |
209 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
210 |
|
|
|
211 |
|
|
if multiTimes; |
212 |
|
|
%zonal wind stress: |
213 |
|
|
fld=alldiag.fldTZ_std; |
214 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]/2000; title0='std(tauZ)'; |
215 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*0.005; end; |
216 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
217 |
|
|
myCaption={myYmeanTxt,' standard deviation -- tauZ (W/m2)'}; |
218 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
219 |
|
|
|
220 |
|
|
%meridional wind stress: |
221 |
|
|
fld=alldiag.fldTM_std; |
222 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]/2000; title0='std(tauM)'; |
223 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*0.005; end; |
224 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
225 |
|
|
myCaption={myYmeanTxt,' standard deviation -- tauM (W/m2)'}; |
226 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
227 |
|
|
|
228 |
|
|
fld=alldiag.curlTau_std; |
229 |
|
|
cc=[[0:5:25] 35 [50:25:200] [250 300]]/1e9; title0='wind stress curl'; |
230 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.01; end; |
231 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
232 |
|
|
myCaption={myYmeanTxt,'standard deviation -- tauCurl (N/m3)'}; |
233 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
234 |
|
|
end; |
235 |
|
|
|
236 |
|
|
end; |
237 |
|
|
|
238 |
|
|
if sum(strcmp(choicePlot,'all'))|sum(strcmp(choicePlot,'ssh')); |
239 |
|
|
|
240 |
|
|
%ETAN: |
241 |
|
|
fld=alldiag.fldETAN_mean; |
242 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]/100; title0='sea surface height (EXCLUDING ice)'; |
243 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.05; end; |
244 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
245 |
|
|
myCaption={myYmeanTxt,'mean -- sea surface height (EXCLUDING ice, in m)'}; |
246 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
247 |
|
|
|
248 |
|
|
%ETANLEADS: |
249 |
|
|
fld=alldiag.fldETANLEADS_mean; |
250 |
|
|
cc=[[-250:50:-100] [-75 -50] [-35:10:35] [50 75] [100:50:250]]/100; title0='sea surface height (INCLUDING ice)'; |
251 |
|
|
if doAnomalies; cc=scaleAnom*[-1:0.1:1]*0.05; end; |
252 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
253 |
|
|
myCaption={myYmeanTxt,'mean -- sea surface height (INCLUDING ice, in m)'}; |
254 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
255 |
|
|
|
256 |
|
|
if multiTimes; |
257 |
|
|
%ETAN: |
258 |
|
|
fld=alldiag.fldETAN_std; |
259 |
|
|
cc=[0:25:500]/2500; title0='std(ETAN)'; |
260 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*0.02; end; |
261 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
262 |
|
|
myCaption={myYmeanTxt,' standard deviation -- sea surface height (EXCLUDING ice, in m)'}; |
263 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
264 |
|
|
|
265 |
|
|
%ETANLEADS: |
266 |
|
|
fld=alldiag.fldETANLEADS_std; |
267 |
|
|
cc=[0:25:500]/2500; title0='std(ETANLEADS)'; |
268 |
|
|
if doAnomalies; cc=scaleAnom*[0:0.1:1]*0.02; end; |
269 |
|
|
figureL; m_map_gcmfaces(fld,0,{'myCaxis',cc},{'do_m_coast',1},{'myTitle',title0}); |
270 |
|
|
myCaption={myYmeanTxt,' standard deviation -- sea surface height (INCLUDING ice, in m)'}; |
271 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
272 |
|
|
end; |
273 |
|
|
|
274 |
|
|
end; |
275 |
|
|
|
276 |
gforget |
1.1 |
end; |
277 |
gforget |
1.4 |
|
278 |
|
|
|