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_uv200.txt"
cmap=asciiread(c,(/83,3/),"float")
cmap=cmap/255

center=asciiread("../../DATA_SUCCESS_START/latlon_success_start.txt",(/92,2/),"float")
;print(center(:,1))
tinfo=asciiread("../../DATA_SUCCESS_START/time_success_start.txt",(/92,2/),"integer")
yyyymmddhh=tinfo(:,1)
;print(yyyymmddhh)
name=asciiread("../../DATA_SUCCESS_START/name_success.txt",(/92/),"string")
;print(name)

dcase=dimsizes(yyyymmddhh)
U_comp=new((/dcase,57,57/),"float")
V_comp=new((/dcase,57,57/),"float")
Z_comp=new((/dcase,57,57/),"float")

do case=0,dcase-1,1

print("working on "+yyyymmddhh(case)+" -- "+name(case))

data1=addfile("../../DATA_SUCCESS_START/ei.oper.an.pl.regn128uv."+yyyymmddhh(case),"r")
data2=addfile("../../DATA_SUCCESS_START/ei.oper.an.pl.regn128sc."+yyyymmddhh(case),"r")
;names=getfilevarnames(data)
;print(names)

U:=data1->U_GDS4_ISBL
U!0="lev"
U!1="lat"
U!2="lon"

V:=data1->V_GDS4_ISBL
V!0="lev"
V!1="lat"
V!2="lon"

Z:=data2->Z_GDS4_ISBL
;printVarSummary(Z)
Z!0="lev"
Z!1="lat"
Z!2="lon"

;print(U&lev) 
;print(U&lat)

lat=U&lat
lon=U&lon
centlat_ind=ind_nearest_coord(center(case,0),lat,0)
centlon_ind=ind_nearest_coord(center(case,1),lon,0)
centlat=lat(centlat_ind)
centlon=lon(centlon_ind)

U200=U(14,:,:)
V200=V(14,:,:)
MAG=sqrt(U200^2+V200^2)
Z200=Z(14,:,:)/100

copy_VarMeta(U(0,:,:),U200)
copy_VarMeta(U(0,:,:),V200)
copy_VarMeta(U(0,:,:),MAG)
copy_VarMeta(Z(0,:,:),Z200)

;print(" "+centlat+" "+lat(centlat_ind+28)+" "+lat(centlat_ind-28))
;print(" "+centlon+" "+lon(centlon_ind+28)+" "+lon(centlon_ind-28))

U_comp(case,:,:)=U200(centlat_ind+28:centlat_ind-28,centlon_ind-28:centlon_ind+28)
V_comp(case,:,:)=V200(centlat_ind+28:centlat_ind-28,centlon_ind-28:centlon_ind+28)
Z_comp(case,:,:)=Z200(centlat_ind+28:centlat_ind-28,centlon_ind-28:centlon_ind+28)

wks=gsn_open_wks("png","success_start_case_"+case)
gsn_define_colormap(wks,cmap)

res=True

res@gsnAddCyclic=False
res@gsnFrame=False
res@gsnDraw=False
res@gsnMaximize=False
res@gsnStringFont=21
res@gsnLeftString="200hPa Wind (m/s), Geopotential Height (dam)"
res@gsnLeftStringOrthogonalPosF=0.01
res@gsnLeftStringFontHeightF=0.015
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@gsnScalarContour=True

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

res@vcRefMagnitudeF=50
res@vcRefLengthF=0.035
res@vcGlyphStyle="CurlyVector"
res@vcLineArrowThicknessF=3.0
res@vcRefAnnoOrthogonalPosF=-0.079
res@vcRefAnnoParallelPosF=1.0
res@vcMinDistanceF=0.03
res@vcMonoLineArrowColor=True
res@vcLineArrowColor=1
res@vcRefAnnoArrowUseVecColor=False
res@vcRefAnnoArrowLineColor=1

;res@stLineColor=0
;res@stArrowLengthF=0.004
;res@stMinArrowSpacingF=0.004
;res@stMinDistanceF=0.1
;res@stLineThicknessF=2.0
;res@stStepSizeF=0.01

res@cnFillOn=True
res@cnLevelSelectionMode="ExplicitLevels"
res@cnLevels:=fspan(0,80,81)
res@cnFillOpacityF=1
res@cnLinesOn=False

res@lbBoxSeparatorLinesOn=False
res@lbLabelStride=20

res@pmLabelBarHeightF=0.1

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

res@tmXTOn                = False
res@tmYROn                = False

res@tmXBMode              = "Explicit"
res@tmXBValues            = ispan(100,180,10)
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(0,50,10)
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

plot1=gsn_csm_contour_map(wks,MAG,res)

res@cnFillOn=False
res@cnLinesOn=True
res@cnInfoLabelOn=False

zmin=min(Z200)
zmax=max(Z200)

res@cnLevels:=ispan(floattoint(zmin),floattoint(zmax),5)
res@cnLineLabelPlacementMode="constant"
res@cnLineLabelInterval=2
res@cnLineThicknessF=3.0

plot2=gsn_csm_contour_map(wks,smth9_Wrap(Z200,0.5,0.25,False),res)

plot3=gsn_csm_vector_map(wks,U200,V200,res)

draw(plot1)
draw(plot2)
draw(plot3)

frame(wks)

end do

rm=systemfunc("rm -rf UVZ_composite.nc")

fout=addfile("./UVZ_composite.nc","c")
fout->U=U_comp
fout->V=V_comp
fout->Z=Z_comp

end