BLT_DT02 = TTD_DTm02 - MLD_DReqDTm02
TTD_DTm02 = TTD_DTm02 = depth where (θ = θ10m - 0.2 °C)
http://www.ifremer.fr/cerweb/deboyer/mld/Surface_Warm_Layer_Depth.php#WLD_DTm02
MLD_DReqDTm02 = depth where ( σ0 = σ010m + Δσ0 )
with Δσ0 = σ0(θ10m - 0.2°C, S10m, P0) - σ0(θ10m, S10m, P0)
http://www.ifremer.fr/cerweb/deboyer/mld/Surface_Mixed_Layer_Depth.php#MLD_DReqDTm02
1.1. Global Scale, Annual Climatology (12 month)
::::::::::::::
pl.clim.MLD.deBoyer.run.sh
::::::::::::::
#!/bin/sh
export LANC=C
exe=./runncl.sh
ncl=$(basename $0 .run.sh).ncl
if [ ! -f $ncl ]; then
echo
echo Error in $0 : No such file, $ncl
echo
exit 1
fi
# Default value
dataset=dr080g
yyyy=1993
yearday=n10day_01_09
dthr=240
day=01
vlev=1
iys=1993 #2005
iye=2014 #2012
indir=/work4/data/ECCO.dr080g
if [ ! -d $indir ]; then
echo
echo Error in $0 : No such directory, $indir
echo
exit 1
fi
outdir=${HOME}/2015.Ibnu.Indian.Ocean/Fig/$(basename $(pwd))/\
$(basename $0 .sh)/
if [ ! -d $outdir ]; then mkdir -vp $outdir ; fi
yearday_list="\
n10day_01_09 \
n10day_10_18 \
n10day_19_27 \
n10day_28_37 \
"
for yearday in $yearday_list; do
if [ $yearday = "n10day_01_09" ]; then
datainfo="_08_08.00001_02160_"
elif [ $yearday = "n10day_10_18" ]; then
datainfo="_08_08.02160_04320_"
elif [ $yearday = "n10day_19_27" ]; then
datainfo="_08_08.04320_06480_"
elif [ $yearday = "n10day_28_37" ]; then
datainfo="_08_08.06480_08880_"
fi
$exe $ncl ${indir} ${yearday} ${datainfo} ${iys} ${iye} ${vlev} ${outdir}
done #yearday
echo
echo Done $(basename $0).
echo
exit 0
::::::::::::::
runncl.sh
::::::::::::::
#!/bin/bash
#
# Universal wrapper script for ncl.
# Pass arguments from the command line to environment variables
#
# version 0.1, Thierry Corti, C2SM ETH Zurich
#
E_BADARGS=65
if [ ! -n "$1" ]
then
echo "Usage: `basename $0` script.ncl argument1 argument2 etc."
exit $E_BADARGS
fi
# save number of arguments to environment variable NCL_N_ARG
export NCL_N_ARGS=$#
# save command line arguments to environment variable NCL_ARG_#
for ((index=1; index<=$#; index++))
do
eval export NCL_ARG_$index=\$$index
done
# run ncl
ncl -Q -n $1
exit 0
::::::::::::::
pl.clim.MLD.deBoyer.ncl
::::::::::::::
;
; Plot horizontal map of Monthly Climatology of Temperature and Salinity
;
; Refer to http://dx.doi.org/10.1002/2014JC010336 for details on dr080
;
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
begin
; Default values
dataset="dr080g"
dthr=240
iys=1993
iye=2014
ndy=9
nmo=3 ; Number of months in each file
figtype="png" ;"ps"
;
; Read command line arguments
;
scname = getenv("NCL_ARG_1")
indir = getenv("NCL_ARG_2")
yearday = getenv("NCL_ARG_3")
datainfo = getenv("NCL_ARG_4")
iys = stringtointeger(getenv("NCL_ARG_5"))
iye = stringtointeger(getenv("NCL_ARG_6"))
vlev = getenv("NCL_ARG_7")
outdir = getenv("NCL_ARG_8")
;iys=(iys_in)
;iye=floattointeger(iye_in)
;Set depth for plotting
klev=stringtointeger(vlev)-1
print("")
print("yearday = "+yearday)
;Set dimension
iy=iys
id=0
indir1=indir+"/"+dataset+"_"+iy+"/"+yearday+"/"
f = addfile(indir1 + "mld" + datainfo + dthr + ".cdf" ,"r")
mld = f->mld(id,:,:)
lat =f->lat
lon =f->lon
time_in=f->time
dt =f->mld@dt
zref =f->mld@zref
;print(dt + " " +zref)
delete(f)
dim=dimsizes(mld)
nyr=iye-iys+1
ttd_in=new((/dim(0),dim(1)/),float)
ttd_in@_FillValue=mld@_FillValue
mld_in=new((/dim(0),dim(1)/),float)
mld_in@_FillValue=mld@_FillValue
blt_in=new((/dim(0),dim(1)/),float)
blt_in@_FillValue=mld@_FillValue
ttd_temp=new((/dim(0),dim(1),nyr/),float)
mld_temp=new((/dim(0),dim(1),nyr/),float)
blt_temp=new((/dim(0),dim(1),nyr/),float)
ttdclim=new((/nmo,dim(0),dim(1)/),float)
ttdclim@_FillValue=mld@_FillValue
mldclim=new((/nmo,dim(0),dim(1)/),float)
mldclim@_FillValue=mld@_FillValue
bltclim=new((/nmo,dim(0),dim(1)/),float)
bltclim@_FillValue=mld@_FillValue
ttdclim@long_name = "Top of thermocilne depth (Monthly climatology)"
ttdclim@units = "m"
ttdclim@dt = dt
ttdclim@zref = zref
ttdclim!0="time"
ttdclim!1="lat"
ttdclim!2="lon"
ttdclim&lon=lon
ttdclim&lat=lat
mldclim@long_name = "Top of thermocilne depth (Monthly climatology)"
mldclim@units = "m"
mldclim@dt = dt
mldclim@zref = zref
mldclim!0="time"
mldclim!1="lat"
mldclim!2="lon"
mldclim&lon=lon
mldclim&lat=lat
bltclim@long_name = "Top of thermocilne depth (Monthly climatology)"
bltclim@units = "m"
bltclim@dt = dt
bltclim@zref = zref
bltclim!0="time"
bltclim!1="lat"
bltclim!2="lon"
bltclim&lon=lon
bltclim&lat=lat
; Month
mm=new( 3, integer)
mmm=new( 3, string )
if (yearday .eq. "n10day_01_09") then
ms=1
mmm(0)="Jan"
mmm(1)="Feb"
mmm(2)="Mar"
end if
if (yearday .eq. "n10day_10_18") then
ms=4
mmm(0)="Apr"
mmm(1)="May"
mmm(2)="Jun"
end if
if (yearday .eq. "n10day_19_27") then
ms=7
mmm(0)="Jul"
mmm(1)="Aug"
mmm(2)="Sep"
end if
if (yearday .eq. "n10day_28_37") then
ms=10
mmm(0)="Oct"
mmm(1)="Nov"
mmm(2)="Dec"
end if
mm(0)=ms
mm(1)=ms + 1
mm(2)=ms + 2
;
; Set NCL resources
;
res = True
res@gsnDraw = False ; don't draw
res@gsnFrame = False ; don't advance frame
;Figure size
res@vpHeightF = 0.2 ;0.35
res@vpWidthF = 0.25 ;0.55
res@vpYF = 0.6 ;0.2
res@vpXF = 0.1 ;0.2
;Map
res@mpFillOn = True
res@gsnPaperOrientation = "landscape"
res@mpShapeMode = "FreeAspect"
res@gsnAddCyclic = False
;Contour
res@cnFillOn = True ; turn on color fill
res@cnLineColor = "white"
res@lbOrientation = "horizontal"
res@cnLevelSelectionMode = "ManualLevels" ; manually set the contour levels with the following 3 resources
res@gsnContourLineThicknessesScale=1.5
;Choose subregion
res@mpMinLatF =-20
res@mpMaxLatF = 20
res@mpMinLonF = 40
res@mpMaxLonF =110
res@mpDataBaseVersion = "Ncarg4_1" ; use GMT coastline
do im=0,2 ; LOOP FOR MONTH
do iy=0,nyr-1 ; LOOP FOR YEAR
yyyy=iy+iys
print(mmm(im) + " " + yyyy)
;
; Open input files
;
indir1=indir+"/"+dataset+"_"+yyyy+"/"+yearday+"/"
fmld=addfile(indir1 + "mld" + datainfo + dthr + ".cdf" ,"r")
fttd=addfile(indir1 + "ttd" + datainfo + dthr + ".cdf" ,"r")
fblt=addfile(indir1 + "blt" + datainfo + dthr + ".cdf" ,"r")
do id=0,2 ; LOOP FOR BI-PENTAD DAYS
;
; Read input data
;
idx=im*3 + id
ttd_in = fttd->ttd(idx,:,:)
mld_in = fmld->mld(idx,:,:)
blt_in = fblt->blt(idx,:,:)
; For climatrogical average
ttd_temp(:,:,iy) = ttd_in(:,:)
mld_temp(:,:,iy) = mld_in(:,:)
blt_temp(:,:,iy) = blt_in(:,:)
end do ;id
end do ; iy
ttdclim(im,:,:)= dim_avg_n(ttd_temp, 2)
mldclim(im,:,:)= dim_avg_n(mld_temp, 2)
bltclim(im,:,:)= dim_avg_n(blt_temp, 2)
;
; Plotting
;
figdir1=outdir+"/"+dataset
foo=systemfunc("mkdir -vp "+figdir1)
figfle=figdir1+"/"+"ECCO.dr080_Clim.MLD.deBoyer" + "."+ mm(im) +"."+figtype
print("Figure: "+figfle)
wks = gsn_open_wks(figtype,figfle)
;
; Left panel
;
res1=res
;res1@vpYF = 0.9
res1@vpXF = 0.05
res1@cnMinLevelValF = 10 ; set min contour level
res1@cnMaxLevelValF = 100 ; set max contour level
res1@cnLevelSpacingF = 10 ; set contour spacing
res1@gsnLeftString = "" ; add the gsn titles
res1@gsnCenterString = "MLD [m]"
res1@gsnRightString = ""
plot1 = gsn_csm_contour_map(wks,mldclim(im,:,:),res1) ; Draw a contour plot.
; now assign plotA to array
draw(plot1)
delete(res1)
;
; Center panel
;
res2=res
;res2@vpYF = 0.45
res2@vpXF = 0.35
res2@cnMinLevelValF = 10 ; set min contour level
res2@cnMaxLevelValF = 100 ; set max contour level
res2@cnLevelSpacingF = 10 ; set contour spacing
res2@gsnLeftString = "" ; add the gsn titles
res2@gsnCenterString = "TTD [m]"
res2@gsnRightString = ""
plot2 = gsn_csm_contour_map(wks,ttdclim(im,:,:),res2) ; Draw a contour plot.
draw(plot2)
delete(res2)
;
; Right panel
;
gsn_define_colormap(wks,"precip3_16lev")
res3=res
;res3@vpYF = 0.45
res3@vpXF = 0.65
res3@cnMinLevelValF = 0 ; set min contour level
res3@cnMaxLevelValF = 40 ; set max contour level
res3@cnLevelSpacingF = 5 ; set contour spacing
res3@gsnLeftString = "" ; add the gsn titles
res3@gsnCenterString = "BLT [m]"
res3@gsnRightString = ""
plot3 = gsn_csm_contour_map(wks,bltclim(im,:,:),res3) ; Draw a contour plot.
draw(plot3)
delete(res3)
txres = True ; text mods desired
txres@txBackgroundFillColor = "white"
txres@txFontColor="black"
txres@txFontHeightF = 0.02
text=dataset
gsn_text_ndc(wks,text,0.475,0.8,txres)
txres@txFontHeightF = 0.015
text="Climatologies of TTD, MLD, and BLT using de Boyer's definition"
gsn_text_ndc(wks,text,0.475,0.75,txres)
text="Z~B~ref~N~="+zref+" ~F33~D~F21~T="+dt
gsn_text_ndc(wks,text,0.475,0.7,txres)
text= mmm(im) + " ("+iys+"-"+iye+")"
gsn_text_ndc(wks,text,0.475,0.65,txres)
;
; Header and footer
;
;drawNDCGrid(wks)
txres = True ; text mods desired
txres@txFontHeightF = 0.01 ; font smaller. default big
gsn_text_ndc(wks,outdir,0.5,0.99,txres)
text=systemfunc("export LANG=C; date")
gsn_text_ndc(wks,text,0.2,0.010,txres)
text=systemfunc("pwd")+"/"+scname
gsn_text_ndc(wks,text,0.5,0.03,txres)
frame(wks)
end do ; im
;
; Output netCDF files
;
time=new(1,float)
time!0="time"
time@units="Hours from Jan 1st"
do im=0,2
time=(mm(im)-0.5)*(30.0*24.0)
ofle1=indir+"/"+"ECCO.dr080g.ttdclm."+iys+"-"+iye+"."+sprinti("%0.2i",mm(im))+".nc"
ofle2=indir+"/"+"ECCO.dr080g.mldclm."+iys+"-"+iye+"."+sprinti("%0.2i",mm(im))+".nc"
ofle3=indir+"/"+"ECCO.dr080g.bltclm."+iys+"-"+iye+"."+sprinti("%0.2i",mm(im))+".nc"
system("/bin/rm -f " + ofle1)
fout1 =addfile (ofle1, "c")
system("/bin/rm -f " + ofle2)
fout2 =addfile (ofle2, "c")
system("/bin/rm -f " + ofle3)
fout3 =addfile (ofle3, "c")
setfileoption(fout1,"DefineMode",True)
setfileoption(fout1,"MissingToFillValue",True)
setfileoption(fout2,"DefineMode",True)
setfileoption(fout2,"MissingToFillValue",True)
setfileoption(fout3,"DefineMode",True)
setfileoption(fout3,"MissingToFillValue",True)
fAtt = True ; assign file attributes
fAtt@title = "ECCO dr080g (Monthly Climatology)"
fAtt@month = mmm(im)
fAtt@description = "Averaging period: "+iys+"-"+iye
fAtt@Conventions = "None"
fAtt@history= "Converted from MITGCM output on 07/24/2015\n"
fAtt@creation_date = systemfunc ("export LANG=C; date")
fAtt@hostname = systemfunc("hostname")
fAtt@cwd = systemfunc("pwd")
fAtt@script = scname
fAtt@input_dir = indir
fileattdef( fout1, fAtt ) ; copy file attributes
fileattdef( fout2, fAtt ) ; copy file attributes
fileattdef( fout3, fAtt ) ; copy file attributes
ntim = 1; dimsizes(time) ; get dimension sizes
nlat = dimsizes(lat)
nlon = dimsizes(lon)
dimNames = (/"time", "lat", "lon" /)
dimSizes = (/ ntim, nlat, nlon /)
dimUnlim = (/False, False, False /)
filedimdef(fout1,dimNames,dimSizes,dimUnlim)
filedimdef(fout2,dimNames,dimSizes,dimUnlim)
filedimdef(fout3,dimNames,dimSizes,dimUnlim)
filevardef(fout1, "time" ,typeof(time),getvardims(time))
filevardef(fout1, "lat" ,typeof(lat), getvardims(lat))
filevardef(fout1, "lon" ,typeof(lon), getvardims(lon))
filevardef(fout1, "ttdclim" ,typeof(ttdclim),getvardims(ttdclim(im,:,:)))
filevardef(fout2, "time" ,typeof(time),getvardims(time))
filevardef(fout2, "lat" ,typeof(lat), getvardims(lat))
filevardef(fout2, "lon" ,typeof(lon), getvardims(lon))
filevardef(fout2, "mldclim" ,typeof(mldclim),getvardims(mldclim(im,:,:)))
filevardef(fout3, "time" ,typeof(time),getvardims(time))
filevardef(fout3, "lat" ,typeof(lat), getvardims(lat))
filevardef(fout3, "lon" ,typeof(lon), getvardims(lon))
filevardef(fout3, "bltclim" ,typeof(bltclim),getvardims(bltclim(im,:,:)))
filevarattdef(fout1,"time" ,time)
filevarattdef(fout1,"lat" ,lat)
filevarattdef(fout1,"lon" ,lon)
filevarattdef(fout1,"ttdclim" ,ttdclim)
filevarattdef(fout2,"time" ,time)
filevarattdef(fout2,"lat" ,lat)
filevarattdef(fout2,"lon" ,lon)
filevarattdef(fout2,"mldclim" ,mldclim)
filevarattdef(fout3,"time" ,time)
filevarattdef(fout3,"lat" ,lat)
filevarattdef(fout3,"lon" ,lon)
filevarattdef(fout3,"bltclim" ,bltclim)
fout1->time=(/time/)
fout1->lat =(/lat/)
fout1->lon = (/lon/)
fout1->ttdclim = (/ttdclim(im,:,:)/)
fout2->time=(/time/)
fout2->lat =(/lat/)
fout2->lon = (/lon/)
fout2->mldclim = (/mldclim(im,:,:)/)
fout3->time=(/time/)
fout3->lat =(/lat/)
fout3->lon = (/lon/)
fout3->bltclim = (/bltclim(im,:,:)/)
print("Output file:"+ofle1)
print("Output file:"+ofle2)
print("Output file:"+ofle3)
print("")
end do ;im
end