| 1 |
|
| 2 |
if 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 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/' ],[dirModel 'diags/']}; |
| 17 |
end; |
| 18 |
|
| 19 |
if 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 |
%convergence in kiloton/s |
| 38 |
seaIceDiv=1e-6*mskC.*calc_UV_conv(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 userStep==-1&myparms.diagsAreMonthly==1;%plotting |
| 72 |
|
| 73 |
if ~doAnomalies; |
| 74 |
|
| 75 |
if addToTex; write2tex(fileTex,1,'Monthly Thickness Distribution',2); end; |
| 76 |
|
| 77 |
%compute seasonal cycle in thickness distributions |
| 78 |
|
| 79 |
y=alldiag.iceThickNorth(:,1:12); |
| 80 |
x=ones(size(y,1),1)*[1:12]; |
| 81 |
iceNorth=zeros(size(x)); iceSouth=zeros(size(x)); |
| 82 |
snowNorth=zeros(size(x)); snowSouth=zeros(size(x)); |
| 83 |
for mm=1:12; |
| 84 |
tmp1=alldiag.iceAreaNorth(:,mm:12:nt); |
| 85 |
tmp1(isnan(tmp1))=0; |
| 86 |
iceNorth(:,mm)=log10(mean(tmp1,2)); |
| 87 |
% |
| 88 |
tmp1=alldiag.iceAreaSouth(:,mm:12:nt); |
| 89 |
tmp1(isnan(tmp1))=0; |
| 90 |
iceSouth(:,mm)=log10(mean(tmp1,2)); |
| 91 |
% |
| 92 |
tmp1=alldiag.snowAreaNorth(:,mm:12:nt); |
| 93 |
tmp1(isnan(tmp1))=0; |
| 94 |
snowNorth(:,mm)=log10(mean(tmp1,2)); |
| 95 |
% |
| 96 |
tmp1=alldiag.snowAreaSouth(:,mm:12:nt); |
| 97 |
tmp1(isnan(tmp1))=0; |
| 98 |
snowSouth(:,mm)=log10(mean(tmp1,2)); |
| 99 |
end; |
| 100 |
iceNorth(~isfinite(iceNorth))=NaN; |
| 101 |
iceSouth(~isfinite(iceSouth))=NaN; |
| 102 |
snowNorth(~isfinite(snowNorth))=NaN; |
| 103 |
snowSouth(~isfinite(snowSouth))=NaN; |
| 104 |
|
| 105 |
%now display seasonal cycle |
| 106 |
x=[x x+12]-0.5; y=[y y]; z=[iceNorth iceNorth]; |
| 107 |
|
| 108 |
figureL; set(gcf,'Renderer','zbuffer'); |
| 109 |
subplot(2,1,1); |
| 110 |
z=[iceNorth iceNorth]; pcolor(x,y,z); shading interp; |
| 111 |
axis([1 13 0 6]); grid on; caxis([9.5 12]); colormap(jet(25)); |
| 112 |
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
| 113 |
title('log(ice area) -- North. Hemi.'); |
| 114 |
subplot(2,1,2); |
| 115 |
z=[snowNorth snowNorth]; pcolor(x,y,z); shading interp; |
| 116 |
axis([1 13 0 1.2]); grid on; caxis([9.5 12]+0.5); colormap(jet(25)); |
| 117 |
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
| 118 |
title('log(snow area) -- North. Hemi.'); |
| 119 |
|
| 120 |
myCaption={myYmeanTxt,'Northern Hemisphere :',... |
| 121 |
'monthly mean ice (top) and snow (bottom)',... |
| 122 |
' thickness distribution (in log(m$^2$))'}; |
| 123 |
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
| 124 |
|
| 125 |
figureL; set(gcf,'Renderer','zbuffer'); |
| 126 |
subplot(2,1,1); |
| 127 |
z=[iceSouth iceSouth]; pcolor(x,y,z); shading interp; |
| 128 |
axis([1 13 0 4]); grid on; caxis([9.5 12]); colormap(jet(25)); |
| 129 |
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
| 130 |
title('log(ice area) -- South. Hemi.'); |
| 131 |
subplot(2,1,2); |
| 132 |
z=[snowSouth snowSouth]; pcolor(x,y,z); shading interp; |
| 133 |
axis([1 13 0 1.2]); grid on; caxis([9.5 12]+0.5); colormap(jet(30)); |
| 134 |
colorbar; xlabel('month of year'); ylabel('log(thickness)'); |
| 135 |
title('log(snow area) -- South. Hemi.'); |
| 136 |
|
| 137 |
myCaption={myYmeanTxt,'Southern Hemisphere :',... |
| 138 |
'monthly mean ice (top) and snow (bottom)',... |
| 139 |
' thickness distribution (in log(m$^2$))'}; |
| 140 |
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
| 141 |
|
| 142 |
end; |
| 143 |
|
| 144 |
%maps of mean march and september ice fields |
| 145 |
|
| 146 |
list_var={'seaIceConc',... |
| 147 |
'seaIceThick',... |
| 148 |
'seaSnowThick',... |
| 149 |
'seaIceStrFun',... |
| 150 |
'seaIceDiv'}; |
| 151 |
|
| 152 |
list_cc={[0:0.1:1],... |
| 153 |
[0:0.1:0.5 0.75 1:0.5:2.5 3:5],... |
| 154 |
[0:0.1:0.5 0.75 1:0.5:2.5 3:5]/5,... |
| 155 |
1/50*[-5:-3 -2.5:0.5:-1 -0.5:0.2:0.5 1:0.5:2.5 3:5],... |
| 156 |
1/5*[-5:-3 -2.5:0.5:-1 -0.5:0.2:0.5 1:0.5:2.5 3:5]}; |
| 157 |
|
| 158 |
list_ccAno={2*[-0.2:0.04:0.2],... |
| 159 |
4*[-0.2:0.04:0.2],... |
| 160 |
[-0.2:0.04:0.2],... |
| 161 |
0.25*[-0.2:0.04:0.2],... |
| 162 |
[-0.2:0.04:0.2]}; |
| 163 |
|
| 164 |
list_tit={' ice concentration (unitless)',... |
| 165 |
' ice thickness (m)',... |
| 166 |
' snow thickness (m)',... |
| 167 |
' ice+snow streamfunction (megaton/s)',... |
| 168 |
' ice+snow convergence (kiloton/s)'}; |
| 169 |
|
| 170 |
for hemi=1:2; |
| 171 |
for seas=1:2; |
| 172 |
|
| 173 |
%select month |
| 174 |
if seas==1; mon='March'; else; mon='September'; end; |
| 175 |
|
| 176 |
%select projection |
| 177 |
if hemi==1; pp=2; txt='Northern '; else; pp=3; txt='Southern '; end; |
| 178 |
|
| 179 |
if addToTex; write2tex(fileTex,1,[txt 'Hem. in ' mon],2); end; |
| 180 |
|
| 181 |
for vv=1:length(list_var); |
| 182 |
|
| 183 |
eval(['fld=alldiag.' list_var{vv} ';']); |
| 184 |
|
| 185 |
%mask out negligible amounts (e.g. due to monthly averaging) |
| 186 |
msk=1+0*alldiag.seaIceConc; |
| 187 |
msk(alldiag.seaIceConc<1e-2)=NaN; |
| 188 |
msk(alldiag.seaIceThick<1e-2)=NaN; |
| 189 |
fld=fld.*msk; |
| 190 |
|
| 191 |
%compute mean march and september fields |
| 192 |
if seas==1; |
| 193 |
tmp1=fld(:,:,3:12:nt); |
| 194 |
tmp1(isnan(tmp1))=0; |
| 195 |
fld_seas=mean(tmp1,3); |
| 196 |
else; |
| 197 |
tmp1=fld(:,:,9:12:nt); |
| 198 |
tmp1(isnan(tmp1))=0; |
| 199 |
fld_seas=mean(tmp1,3); |
| 200 |
end; |
| 201 |
fld_seas(fld_seas==0)=NaN; |
| 202 |
|
| 203 |
%ice concentration |
| 204 |
figureL; set(gcf,'Renderer','zbuffer'); |
| 205 |
cc=list_cc{vv}; if doAnomalies; cc=list_ccAno{vv}; end; |
| 206 |
caxi={'myCaxis',cc}; docb={'doCbar',1}; |
| 207 |
m_map_gcmfaces(fld_seas,pp,caxi,docb); |
| 208 |
myCaption={myYmeanTxt,mon,' mean -- ',list_tit{vv}}; |
| 209 |
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
| 210 |
|
| 211 |
end;%for vv=1:length(list_var); |
| 212 |
end;%for seas=1:2; |
| 213 |
end;%for hemi=1:2; |
| 214 |
|
| 215 |
end; |
| 216 |
|
Parent Directory
| 
Revision Log
| 
Revision Graph