| 1 |
|
| 2 |
if userStep==1;%diags to be computed |
| 3 |
listDiags='fldMldBoyer fldMldSuga fldMldKara'; |
| 4 |
elseif userStep==2;%input files and variables |
| 5 |
listFlds={ 'THETA','SALT'}; |
| 6 |
listFldsNames=deblank(listFlds); |
| 7 |
listFiles={'monthly_2d_set1','monthly_3d_set1','state_2d_set1','other_2d_set1','state_3d_set1'}; |
| 8 |
elseif userStep==3;%computational part; |
| 9 |
fldT=THETA.*mygrid.mskC; fldS=SALT.*mygrid.mskC; |
| 10 |
% |
| 11 |
%prepare to compute potential density: |
| 12 |
fldP=0*mygrid.mskC; for kk=1:length(mygrid.RC); fldP(:,:,kk)=-mygrid.RC(kk); end; |
| 13 |
T=convert2vector(fldT); |
| 14 |
S=convert2vector(fldS); |
| 15 |
msk=convert2vector(mygrid.mskC); |
| 16 |
P=convert2vector(fldP); |
| 17 |
%compute potential density: |
| 18 |
RHO=0*msk; alpha=0*msk; |
| 19 |
tmp1=find(~isnan(msk)); |
| 20 |
RHO(tmp1) = density(T(tmp1),S(tmp1),P(tmp1)); |
| 21 |
fldRhoPot=convert2vector(RHO); |
| 22 |
alpha(tmp1) = density(T(tmp1)+1e-4,S(tmp1),P(tmp1)); |
| 23 |
fldAlpha=(convert2vector(alpha)-fldRhoPot)/1e-4; |
| 24 |
|
| 25 |
clear T S P msk RHO RHOis tmp1; |
| 26 |
|
| 27 |
%compute mld: |
| 28 |
tmp1=NaN*mygrid.mskC(:,:,1); |
| 29 |
for kk=1:50; |
| 30 |
tmp2=fldRhoPot(:,:,kk)-fldRhoPot(:,:,1); |
| 31 |
%if we pass RHO(1)+0.03 for the first time (or we reach the bottom) |
| 32 |
%then mld is the velocity point above RC(kk), which is RF(kk) |
| 33 |
jj=find((tmp2>0.03|isnan(tmp2))&isnan(tmp1)); |
| 34 |
tmp1(jj)=-mygrid.RF(kk); |
| 35 |
end; |
| 36 |
fldMldBoyer=tmp1; |
| 37 |
|
| 38 |
%compute mld: |
| 39 |
tmp1=NaN*mygrid.mskC(:,:,1); |
| 40 |
for kk=1:50; |
| 41 |
tmp2=fldRhoPot(:,:,kk)-fldRhoPot(:,:,1); |
| 42 |
%if we pass RHO(1)+0.125 for the first time (or we reach the bottom) |
| 43 |
%then mld is the velocity point above RC(kk), which is RF(kk) |
| 44 |
jj=find((tmp2>0.125|isnan(tmp2))&isnan(tmp1)); |
| 45 |
tmp1(jj)=-mygrid.RF(kk); |
| 46 |
end; |
| 47 |
fldMldSuga=tmp1; |
| 48 |
|
| 49 |
%compute mld: |
| 50 |
tmp1=NaN*mygrid.mskC(:,:,1); |
| 51 |
fldRhoPotMax=fldRhoPot(:,:,1)-0.8*fldAlpha(:,:,1); |
| 52 |
for kk=1:50; |
| 53 |
tmp2=fldRhoPot(:,:,kk)-fldRhoPotMax; |
| 54 |
%if we pass RHO(1)+0.8*alpha(1) for the first time (or we reach the bottom) |
| 55 |
%then mld is the velocity point above RC(kk), which is RF(kk) |
| 56 |
jj=find((tmp2>0|isnan(tmp2))&isnan(tmp1)); |
| 57 |
tmp1(jj)=-mygrid.RF(kk); |
| 58 |
end; |
| 59 |
fldMldKara=tmp1; |
| 60 |
|
| 61 |
elseif userStep==-1;%plotting; |
| 62 |
|
| 63 |
list_var={'fldMldKara','fldMldSuga','fldMldBoyer'}; |
| 64 |
|
| 65 |
list_tit={' mixed layer depth per Kara formula (m)',... |
| 66 |
' mixed layer depth per Suga formula (m)',... |
| 67 |
' mixed layer depth per Boyer M. formula (m)'}; |
| 68 |
|
| 69 |
if myparms.diagsAreMonthly==1; |
| 70 |
for seas=1:2; |
| 71 |
for vv=1:length(list_var); |
| 72 |
|
| 73 |
eval(['fld=alldiag.' list_var{vv} ';']); |
| 74 |
|
| 75 |
%compute mean march and september fields |
| 76 |
if seas==1; fld_seas=nanmedian(fld(:,:,3:12:nt),3); mon='March'; |
| 77 |
else; fld_seas=nanmedian(fld(:,:,9:12:nt),3); mon='September'; |
| 78 |
end; |
| 79 |
fld_seas(fld_seas==0)=NaN; |
| 80 |
|
| 81 |
%plot |
| 82 |
cc=[[0:20:100] [150:50:300] 400 [500:200:1100] [1500:500:2000]]; |
| 83 |
if doAnomalies; cc=scaleAnom*[-5:0.5:5]; end; |
| 84 |
figureL; set(gcf,'Renderer','zbuffer'); |
| 85 |
m_map_gcmfaces(fld_seas,0,{'myCaxis',cc}); |
| 86 |
myCaption={myYmeanTxt,mon,' mean -- ',list_tit{vv}}; |
| 87 |
if addToTex; write2tex(fileTex,2,myCaption,gcf); end; |
| 88 |
|
| 89 |
end;%for vv=1:length(list_var); |
| 90 |
end;%for seas=1:2; |
| 91 |
end;%if myparms.diagsAreMonthly==1 |
| 92 |
|
| 93 |
end; |
| 94 |
|
| 95 |
|
Parent Directory
| 
Revision Log
| 
Revision Graph