$ srcdump.sh jofuro-2.monthly.trend.ncl
------------------------------
List of the following files:
------------------------------
jofuro-2.monthly.trend.ncl
------------------------------
Machine info
------------------------------
aofd165.bio.mie-u.ac.jp
/work2/am/2017.Heavy.Rain/J-OFURO2
Sat Jul 15 08:01:15 JST 2017
======================
jofuro-2.monthly.trend.ncl
======================
begin
indir="/misc/raid112MANDA/KakenA.H28/J-OFURO2"
prefix="J-OFURO2_LHF_HF004_MONTHLY_"
;mmm="JUN"
mmm="JUL"
lonw=122.
lone=130.
lats=25.
latn=33.
script_name = get_script_name()
script_no_suffix = systemfunc("basename "+ script_name +" .ncl")
;print("script_no_suffix =" +script_no_suffix)
;print("")
figdir="Fig."+script_no_suffix
system("mkdir -vp "+figdir)
if (mmm .eq. "JUN" )then
mm=6-1
end if
if (mmm .eq. "JUL" )then
mm=7-1
end if
yyyys=1988
yyyye=2008
yyyy=yyyys
infle=prefix+yyyy+".nc"
a = addfile(indir+"/"+infle, "r")
lon = a->XLONG
lat = a->YLAT
tmp = a->LHF
dim=dimsizes(tmp)
nt=yyyye-yyyys+1
ny=dim(1)
nx=dim(2)
lhflt=new((/nt/),typeof(tmp))
lat2d=new((/ny,nx/),float)
lon2d=new((/ny,nx/),float)
do j=0,ny-1
do i=0,nx-1
lon2d(j,i)=tofloat(lon(i))
lat2d(j,i)=tofloat(lat(j))
end do
end do
rad = 4.0*atan(1.0)/180.0
clat = cos(lat*rad)
it=-1
do yyyy=yyyys,yyyye
infle=prefix+yyyy+".nc"
a = addfile(indir+"/"+infle, "r")
; print(a)
lhf = a->LHF(mm,:,:)
it=it+1
lhf_sub=where(\
(lon2d .ge. lonw .and. lon2d .le. lone .and. \
lat2d .ge. lats .and. lat2d .le. latn),\
lhf, lhf@_FillValue)
lhflt(it)=wgt_areaave(lhf_sub, clat, 1.0, 0)
end do ;yyyy
nyr=yyyye-yyyys+1
year=fspan(yyyys,yyyye,nyr)
xa=year-yyyys
rc = regline(xa, lhflt)
rc@units = "Wm2/year"
print(rc)
ltrend = new ( (/nyr/), float, lhflt@_FillValue)
ltrend = rc*xa + rc@yintercept
p=trend_manken(lhflt,False,0)
ttrend = new ( (/nyr/), float, lhflt@_FillValue)
ttrend = p(1)*xa + rc@yintercept
print(rc)
; print(p)
x=year
y=new( (/3,nyr/), float, lhflt@_FillValue)
y(0,:)=lhflt
y(1,:)=ltrend
y(2,:)=ttrend
region=toint(lonw)+"-"+toint(lone)+"_"+toint(lats)+"-"+toint(latn)
figfile=figdir+"/"+script_no_suffix+"_"+prefix+mmm+"_"+region
figtype="png"
wks = gsn_open_wks (figtype,figfile)
res = True
opt=res
opt@gsnLeftString = "LHF"
opt@gsnCenterString = mmm +"("+yyyys+"-"+yyyye+")"
opt@gsnRightString = "W/m~S~2~N~"
opt@gsnFrame = False
opt@gsnDraw = False
plot=gsn_csm_xy(wks,x,y,opt)
; HEADER
system("export LANG=C")
txres = True
txres@txFontHeightF = 0.01
txres@txJust="CenterLeft"
text=indir
gsn_text_ndc(wks,text,0.05,0.97,txres)
text=systemfunc("ls -lh "+ indir +"/"+ infle)
gsn_text_ndc(wks,text,0.05,0.95,txres)
txres@txFontHeightF = 0.015
text=lonw+"E-"+lone+"E"+" "+lats+"N-"+latn+"N"
gsn_text_ndc(wks,text,0.05,0.92,txres)
text="Trend (Linear regression) = "+rc
gsn_text_ndc(wks,text,0.05,0.89,txres)
text="Trend (Theil-Sen) = "+p(1)+" (p="+p(0)+")"
gsn_text_ndc(wks,text,0.05,0.86,txres)
;FOOTER
txres@txFontHeightF = 0.01
txres@txJust="CenterLeft"
text = systemfunc("date")
gsn_text_ndc(wks,text, 0.05,0.010,txres)
cwd =systemfunc("pwd")
gsn_text_ndc(wks,cwd, 0.05,0.025,txres)
gsn_text_ndc(wks,figfile+"."+figtype,0.05,0.040,txres)
draw(plot)
frame(wks)
print("")
print("Input: "+infle)
print("Output: "+figfile)
print("")
end
----------------------
End of jofuro-2.monthly.trend.ncl
----------------------