1 |
gforget |
1.1 |
|
2 |
|
|
if strcmp(setDiags,'SEAICE')&userStep==1;%diags to be computed |
3 |
|
|
listDiags='seaIceConc seaIceThick seaSnowThick seaIceStrFun seaIceDiv'; |
4 |
|
|
listDiags=[listDiags ' iceThickNorth iceThickSouth']; |
5 |
|
|
listDiags=[listDiags ' iceAreaNorth iceAreaSouth']; |
6 |
|
|
listDiags=[listDiags ' snowThickNorth snowThickSouth']; |
7 |
|
|
listDiags=[listDiags ' snowAreaNorth snowAreaSouth']; |
8 |
|
|
end; |
9 |
|
|
|
10 |
|
|
if strcmp(setDiags,'SEAICE')&userStep==2;%input files and variables |
11 |
|
|
listFlds={'SIarea ','SIheff ','SIhsnow ',... |
12 |
|
|
'ADVxHEFF','ADVyHEFF','ADVxSNOW','ADVySNOW',... |
13 |
|
|
'DFxEHEFF','DFyEHEFF','DFxESNOW','DFyESNOW'}; |
14 |
|
|
listFldsNames=deblank(listFlds); |
15 |
|
|
listFiles={'state_2d_set1'}; |
16 |
|
|
listSubdirs={[dirModel 'diags/STATE/' ]}; |
17 |
|
|
end; |
18 |
|
|
|
19 |
|
|
if strcmp(setDiags,'SEAICE')&userStep==3;%computational part; |
20 |
|
|
|
21 |
|
|
%ice masks |
22 |
|
|
mskC=mygrid.mskC(:,:,1); mskC(SIheff==0)=NaN; |
23 |
|
|
mskW=mygrid.mskW(:,:,1); mskW(ADVxHEFF==0)=NaN; |
24 |
|
|
mskS=mygrid.mskS(:,:,1); mskS(ADVyHEFF==0)=NaN; |
25 |
|
|
|
26 |
|
|
%ice and snow thickness |
27 |
|
|
seaIceConc=mskC.*SIarea; |
28 |
|
|
seaIceThick=SIheff./seaIceConc; |
29 |
|
|
seaSnowThick=SIhsnow./seaIceConc; |
30 |
|
|
|
31 |
|
|
%transports in kg/s |
32 |
|
|
tmpU=(myparms.rhoi*DFxEHEFF+myparms.rhosn*DFxESNOW+... |
33 |
|
|
myparms.rhoi*ADVxHEFF+myparms.rhosn*ADVxSNOW); |
34 |
|
|
tmpV=(myparms.rhoi*DFyEHEFF+myparms.rhosn*DFyESNOW+... |
35 |
|
|
myparms.rhoi*ADVyHEFF+myparms.rhosn*ADVySNOW); |
36 |
|
|
|
37 |
|
|
%divergence in kiloton/s |
38 |
|
|
seaIceDiv=1e-6*mskC.*calc_UV_div(tmpU,tmpV); %no dh needed here |
39 |
|
|
|
40 |
|
|
%streamfunction in kg/s |
41 |
|
|
%(factor 1e6 to compensare for calc_barostream internal conversion) |
42 |
|
|
seaIceStrFun=1e6*calc_barostream(tmpU.*mskW,tmpV.*mskS,1,{}); |
43 |
|
|
%convert to megaton/s |
44 |
|
|
seaIceStrFun=1e-9*mskC.*seaIceStrFun; |
45 |
|
|
|
46 |
|
|
%thickness distributions per hemisphere |
47 |
|
|
xx=convert2vector(seaIceThick); |
48 |
|
|
yy=convert2vector(1+0*seaIceThick); |
49 |
|
|
zz=convert2vector(seaIceConc.*mygrid.RAC); |
50 |
|
|
|
51 |
|
|
mm=convert2vector(mygrid.YC>0); jj=find(mm.*zz>0); |
52 |
|
|
[x,y,z,n]=MITprof_stats(xx(jj),[0:0.05:10],yy(jj),[0 1 2],'sum',zz(jj)); |
53 |
|
|
iceThickNorth=x(:,2); iceAreaNorth=z(:,2); |
54 |
|
|
|
55 |
|
|
mm=convert2vector(mygrid.YC<0); jj=find(mm.*zz>0); |
56 |
|
|
[x,y,z,n]=MITprof_stats(xx(jj),[0:0.05:10],yy(jj),[0 1 2],'sum',zz(jj)); |
57 |
|
|
iceThickSouth=x(:,2); iceAreaSouth=z(:,2); |
58 |
|
|
|
59 |
|
|
xx=convert2vector(seaSnowThick); |
60 |
|
|
|
61 |
|
|
mm=convert2vector(mygrid.YC>0); jj=find(mm.*zz>0); |
62 |
|
|
[x,y,z,n]=MITprof_stats(xx(jj),[0:0.05:10],yy(jj),[0 1 2],'sum',zz(jj)); |
63 |
|
|
snowThickNorth=x(:,2); snowAreaNorth=z(:,2); |
64 |
|
|
|
65 |
|
|
mm=convert2vector(mygrid.YC<0); jj=find(mm.*zz>0); |
66 |
|
|
[x,y,z,n]=MITprof_stats(xx(jj),[0:0.05:10],yy(jj),[0 1 2],'sum',zz(jj)); |
67 |
|
|
snowThickSouth=x(:,2); snowAreaSouth=z(:,2); |
68 |
|
|
|
69 |
|
|
end; |
70 |
|
|
|
71 |
|
|
if strcmp(setDiags,'SEAICE')&userStep==-1;%computational part; |
72 |
|
|
|
73 |
|
|
if ~doAnomalies; |
74 |
|
|
|
75 |
|
|
%compute seasonal cycle in thickness distributions |
76 |
|
|
|
77 |
|
|
y=alldiag.iceThickNorth(:,1:12); |
78 |
|
|
x=ones(size(y,1),1)*[1:12]; |
79 |
|
|
iceNorth=zeros(size(x)); iceSouth=zeros(size(x)); |
80 |
|
|
snowNorth=zeros(size(x)); snowSouth=zeros(size(x)); |
81 |
|
|
for mm=1:12; |
82 |
|
|
tmp1=alldiag.iceAreaNorth(:,mm:12:nt); |
83 |
|
|
tmp1(isnan(tmp1))=0; |
84 |
|
|
iceNorth(:,mm)=log10(median(tmp1,2)); |
85 |
|
|
% |
86 |
|
|
tmp1=alldiag.iceAreaSouth(:,mm:12:nt); |
87 |
|
|
tmp1(isnan(tmp1))=0; |
88 |
|
|
iceSouth(:,mm)=log10(median(tmp1,2)); |
89 |
|
|
% |
90 |
|
|
tmp1=alldiag.snowAreaNorth(:,mm:12:nt); |
91 |
|
|
tmp1(isnan(tmp1))=0; |
92 |
|
|
snowNorth(:,mm)=log10(median(tmp1,2)); |
93 |
|
|
% |
94 |
|
|
tmp1=alldiag.snowAreaSouth(:,mm:12:nt); |
95 |
|
|
tmp1(isnan(tmp1))=0; |
96 |
|
|
snowSouth(:,mm)=log10(median(tmp1,2)); |
97 |
|
|
end; |
98 |
|
|
iceNorth(~isfinite(iceNorth))=NaN; |
99 |
|
|
iceSouth(~isfinite(iceSouth))=NaN; |
100 |
|
|
snowNorth(~isfinite(snowNorth))=NaN; |
101 |
|
|
snowSouth(~isfinite(snowSouth))=NaN; |
102 |
|
|
|
103 |
|
|
%now display seasonal cycle |
104 |
|
|
x=[x x+12]-0.5; y=[y y]; z=[iceNorth iceNorth]; |
105 |
|
|
|
106 |
|
|
figureL; set(gcf,'Renderer','zbuffer'); |
107 |
|
|
subplot(2,1,1); |
108 |
|
|
z=[iceNorth iceNorth]; pcolor(x,y,z); shading interp; |
109 |
|
|
axis([1 13 0 6]); grid on; caxis([9.5 12]); colormap(jet(25)); |
110 |
|
|
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
111 |
|
|
title('log(ice area) -- North. Hemi.'); |
112 |
|
|
subplot(2,1,2); |
113 |
|
|
z=[snowNorth snowNorth]; pcolor(x,y,z); shading interp; |
114 |
|
|
axis([1 13 0 1.2]); grid on; caxis([9.5 12]+0.5); colormap(jet(25)); |
115 |
|
|
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
116 |
|
|
title('log(snow area) -- North. Hemi.'); |
117 |
|
|
|
118 |
|
|
myCaption={myYmeanTxt,'Northern Hemisphere :',... |
119 |
|
|
'monthly mean ice (top) and snow (bottom)',... |
120 |
|
|
' thickness distribution (in log(m2))'}; |
121 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
122 |
|
|
|
123 |
|
|
figureL; set(gcf,'Renderer','zbuffer'); |
124 |
|
|
subplot(2,1,1); |
125 |
|
|
z=[iceSouth iceSouth]; pcolor(x,y,z); shading interp; |
126 |
|
|
axis([1 13 0 4]); grid on; caxis([9.5 12]); colormap(jet(25)); |
127 |
|
|
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
128 |
|
|
title('log(ice area) -- South. Hemi.'); |
129 |
|
|
subplot(2,1,2); |
130 |
|
|
z=[snowSouth snowSouth]; pcolor(x,y,z); shading interp; |
131 |
|
|
axis([1 13 0 1.2]); grid on; caxis([9.5 12]+0.5); colormap(jet(30)); |
132 |
|
|
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
133 |
|
|
title('log(snow area) -- South. Hemi.'); |
134 |
|
|
|
135 |
|
|
myCaption={myYmeanTxt,'Southern Hemisphere :',... |
136 |
|
|
'monthly mean ice (top) and snow (bottom)',... |
137 |
|
|
' thickness distribution (in log(m2))'}; |
138 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
139 |
|
|
|
140 |
|
|
end; |
141 |
|
|
|
142 |
|
|
%maps of mean march and september ice fields |
143 |
|
|
|
144 |
|
|
list_var={'seaIceConc',... |
145 |
|
|
'seaIceThick',... |
146 |
|
|
'seaSnowThick',... |
147 |
|
|
'seaIceStrFun',... |
148 |
|
|
'seaIceDiv'}; |
149 |
|
|
|
150 |
|
|
list_cc={[0:0.1:1],... |
151 |
|
|
[0:0.1:0.5 0.75 1:0.5:2.5 3:5],... |
152 |
|
|
[0:0.1:0.5 0.75 1:0.5:2.5 3:5]/5,... |
153 |
|
|
1/50*[-5:-3 -2.5:0.5:-1 -0.5:0.2:0.5 1:0.5:2.5 3:5],... |
154 |
|
|
1/5*[-5:-3 -2.5:0.5:-1 -0.5:0.2:0.5 1:0.5:2.5 3:5]}; |
155 |
|
|
|
156 |
|
|
list_ccAno={2*[-0.2:0.04:0.2],... |
157 |
|
|
4*[-0.2:0.04:0.2],... |
158 |
|
|
[-0.2:0.04:0.2],... |
159 |
|
|
0.25*[-0.2:0.04:0.2],... |
160 |
|
|
[-0.2:0.04:0.2]}; |
161 |
|
|
|
162 |
|
|
list_tit={' ice concentration (unitless)',... |
163 |
|
|
' ice thickness (m)',... |
164 |
|
|
' snow thickness (m)',... |
165 |
|
|
' ice+snow streamfunction (megaton/s)',... |
166 |
|
|
' ice+snow divergence (kiloton/s)'}; |
167 |
|
|
|
168 |
|
|
for hemi=1:2; |
169 |
|
|
for seas=1:2; |
170 |
|
|
for vv=1:length(list_var); |
171 |
|
|
|
172 |
|
|
eval(['fld=alldiag.' list_var{vv} ';']); |
173 |
|
|
|
174 |
|
|
%compute mean march and september fields |
175 |
|
|
if seas==1; fld_seas=nanmedian(fld(:,:,3:12:240),3); mon='March'; |
176 |
|
|
else; fld_seas=nanmedian(fld(:,:,9:12:240),3); mon='September'; |
177 |
|
|
end; |
178 |
|
|
fld_seas(fld_seas==0)=NaN; |
179 |
|
|
|
180 |
|
|
%choose projection |
181 |
|
|
if hemi==1; pp=2; txt='Northern '; |
182 |
|
|
else; pp=3; txt='Southern '; |
183 |
|
|
end; |
184 |
|
|
|
185 |
|
|
%ice concentration |
186 |
|
|
figureL; set(gcf,'Renderer','zbuffer'); |
187 |
|
|
cc=list_cc{vv}; if doAnomalies; cc=list_ccAno{vv}; end; |
188 |
|
|
caxi={'myCaxis',cc}; docb={'doCbar',1}; |
189 |
|
|
m_map_gcmfaces(fld_seas,pp,caxi,docb); |
190 |
|
|
myCaption={myYmeanTxt,mon,' mean -- ',list_tit{vv}}; |
191 |
|
|
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
192 |
|
|
|
193 |
|
|
end;%for vv=1:length(list_var); |
194 |
|
|
end;%for seas=1:2; |
195 |
|
|
end;%for hemi=1:2; |
196 |
|
|
|
197 |
|
|
end; |
198 |
|
|
|