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/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"

begin

; for use of high resolution map database, execute
; 'setenv NCARG_RANGS /network/rit/lab/atm419lab/zhou/rangs' (on snow)
; 'setenv NCARG_RANGS /atm419lab/zhou/rangs' (on ash, reed)
; before running the script

c="./cmap_RH.txt"
cmap=asciiread(c,(/72,3/),"float")
cmap=cmap/255

data1=addfile("./RH_composite.nc","r")
data2=addfile("../../../NCL_FAILURE_START/RH/analysis/RH_composite_sub.nc","r")
RH1=data1->RH
RH2=data2->RH
RH1!0="case"
RH1!1="lat"
RH1!2="lon"
;printVarSummary(RH1)
;print(num(ismissing(RH)))
;print(RH&lev)
;print(RH&lat)
RH1&lat=fspan(-15,15,43)
RH1&lon=fspan(-15,15,43)

copy_VarMeta(RH1,RH2)

dsizes=dimsizes(RH1)
;print(dsizes)
dcases=dsizes(0)
dlat=dsizes(1)
dlon=dsizes(2)

RH1_avg=new((/dlat,dlon/),"float")
RH2_avg=new((/dlat,dlon/),"float")

do j=0,dlat-1,1
  do i=0,dlon-1,1
    RH1_avg(j,i)=avg(RH1(:,j,i))
    RH2_avg(j,i)=avg(RH2(:,j,i))
  end do
end do

copy_VarMeta(RH1(0,:,:),RH1_avg)
copy_VarMeta(RH2(0,:,:),RH2_avg)

RH_avg_diff=RH1_avg-RH2_avg

copy_VarMeta(RH1_avg,RH_avg_diff)

wks=gsn_open_wks("png","composite_diff")
gsn_define_colormap(wks,"MPL_BrBG")

res=True

res@gsnAddCyclic=False
res@gsnFrame=False
res@gsnMaximize=False
res@gsnStringFont=21
res@gsnLeftString="850-700hPa Average RH (%)"
res@gsnLeftStringOrthogonalPosF=0.01
res@gsnLeftStringFontHeightF=0.015
res@gsnRightString="difference"
;res@gsnRightString=name(case)+" @ "+sprinti("%10.0i",yyyymmddhh(case))
;res@gsnRightString="+"+sprinti("%0.2i",counter)+"h"
;res@gsnRightString="+"+sprintf("%4.1f",counter*0.5)+"h"
res@gsnRightStringOrthogonalPosF=0.01
res@gsnRightStringFontHeightF=0.015

res@vpXF=0.15
res@vpYF=0.95
res@vpHeightF=0.7
res@vpWidthF=0.7

res@trYMinF=-15
res@trYMaxF=15
res@trXMinF=-15
res@trXMaxF=15

res@cnFillOn=True
res@cnLevelSelectionMode="ExplicitLevels"
res@cnLevels=fspan(-15,15,61)
res@cnFillOpacityF=1
res@cnLinesOn=False

res@lbBoxSeparatorLinesOn=False
res@lbLabelStride=10
res@lbLabelFontHeightF=0.015

res@pmLabelBarHeightF=0.1

;res@mpDataResolution="FinestResolution"
;res@mpDataBaseVersion="HighRes"
;res@mpLimitMode="LatLon"
;res@mpCenterLatF=centlat
;res@mpCenterLonF=centlon
;res@mpMinLatF=centlat-10
;res@mpMaxLatF=centlat+10
;res@mpMinLonF=centlon-10
;res@mpMaxLonF=centlon+10
;res@mpMonoFillColor=False
;res@mpFillColors=(/-1,-1,-1/)
;res@mpOutlineOn=True

res@tmXTOn                = False
res@tmYROn                = False

res@tmXBMode              = "Explicit"
res@tmXBValues            = ispan(-12,12,4)
res@tmXBLabels            = res@tmXBValues
;res@tmXBLabels            = (/"100~S~o~N~E","110~S~o~N~","120~S~o~N~E",\
;                              "130~S~o~N~","140~S~o~N~E","150~S~o~N~",\
;                              "160~S~o~N~E","170~S~o~N~","180~S~o~N~E"/)
res@tmXBMajorLengthF      = 0.01
res@tmXBMajorOutwardLengthF=0.0
res@tmXBMajorThicknessF   = 3.0
res@tmXBMinorOn           = False
res@tmXBMinorPerMajor     = 4
res@tmXBMinorLengthF      = 0.004
res@tmXBMinorOutwardLengthF=0.0
res@tmXBMinorThicknessF   = 2.0
res@tmXBLabelFont         = 21
res@tmXBLabelFontHeightF  = 0.015
res@tmXBLabelDeltaF       = -0.6

res@tmYLMode              = "Explicit"
res@tmYLValues            = ispan(-12,12,4)
res@tmYLLabels            = res@tmYLValues
;res@tmYLLabels            = (/"0~S~o~N~","10~S~o~N~N","20~S~o~N~",\
;                              "30~S~o~N~N","40~S~o~N~","50~S~o~N~N"/)
res@tmYLMajorLengthF      = 0.01
res@tmYLMajorOutwardLengthF=0.0
res@tmYLMajorThicknessF   = 3.0
res@tmYLMinorOn           = False
res@tmYLMinorPerMajor     = 4
res@tmYLMinorLengthF      = 0.004
res@tmYLMinorOutwardLengthF=0.0
res@tmYLMinorThicknessF   = 2.0
res@tmYLLabelFont         = 21
res@tmYLLabelFontHeightF  = 0.015
res@tmYLLabelDeltaF       = -0.3

plot=gsn_csm_contour(wks,RH_avg_diff,res)

frame(wks)

end