; load libraries, twice for good measure

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/wrf/WRFUserARW.ncl"
;load "/home/torn/ncl/wrf_map_info.ncl"
;load "/home/torn/ncl/array_stat.ncl"
;load "/home/torn/ncl/image_conv.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/wrf/WRFUserARW.ncl"
;load "/home/torn/ncl/wrf_map_info.ncl"
;load "/home/torn/ncl/array_stat.ncl"
;load "/home/torn/ncl/image_conv.ncl"
;load "wrf_dart_util.ncl"

begin
quote = inttochar(34) 
; ------ Load Data
delim = ","

lonmin = -81 ;round(min(gustlons),3)-2
lonmax = -71 ;round(max(gustlons),3)+2
latmin = 38  ;round(min(gustlats),3)-2
latmax = 48  ;round(max(gustlats),3)+2

wks_type = "png"
wks_type@wkWidth = 4000
wks_type@wkHeight = 3200 
wks = gsn_open_wks(wks_type,"plot1")

res = True
res@gsnFrame      = False
res@mpProjection = "CylindricalEquidistant"
;res@mpOutlineBoundarySets = "GeophysicalAndUSStates"
res@gsnMaximize = True
res@mpMinLonF = lonmin ;-- min longitude
res@mpMaxLonF = lonmax ;-- max longitude
res@mpMinLatF = latmin ;-- min latitude
res@mpMaxLatF = latmax ;-- max latitude
res@tfDoNDCOverlay = True ;this is a resource for overlaying multiple plots

res@mpDataBaseVersion           = "Ncarg4_1"
res@mpDataSetName               = "Earth..4"   ; For counties
res@mpFillAreaSpecifiers        = (/"Land", "United States:States","Water"/) ; for states only
res@mpFillAreaSpecifiers        = (/"Land", "Land:North America:United States:Conterminous US:New York","Water"/)
res@mpSpecifiedFillColors       = (/"transparent","transparent",       "transparent"/)
res@mpInlandWaterFillColor      = 238
res@mpLandFillColor             = -1     ; transparent
res@mpOutlineOn                 = True
res@mpMaskOutlineSpecifiers =   (/"Land"/)
res@mpOutlineSpecifiers =   (/ "Land:North America:United States:Conterminous US:New York:counties"/)

res@mpGeophysicalLineThicknessF = 20
res@mpUSStateLineThicknessF = 20
res@mpNationalLineThicknessF = 20
res@mpCountyLineThicknessF = 5
res@mpGridLonSpacingF = 2
res@mpDataBaseVersion = "MediumRes" ;-- better map resolution


  res@tmXBBorderOn      = False
  res@tmXBMinorPerMajor = 1
  res@tmXBOn            = False
  res@tmXTBorderOn      = False
  res@tmXTOn            = False
  res@tmYLBorderOn      = False
  res@tmYLMinorPerMajor = 1
  res@tmYLOn            = False
  res@tmYRBorderOn      = False
  res@tmYROn            = False 
res@gsnDraw             = False

res@mpOutlineDrawOrder   = "PreDraw"       ;      = "PostDraw" 


; make feelslike plot below
wks_type = "png"
wks_type@wkWidth = 2000
wks_type@wkHeight = 1600 
wks = gsn_open_wks(wks_type,"highlights")
gsn_define_colormap(wks,"GMT_wysiwygcont")

plot = gsn_csm_map(wks,res)

sitefile  = "maxwind1day.dat"
sitedata  = asciiread(sitefile,-1,"string")
gustlats  = tofloat(str_get_field(sitedata,1,delim))
gustlons  = tofloat(str_get_field(sitedata,2,delim))
gustdata  = tofloat(str_get_field(sitedata,3,delim))

gust = .1*round(max(10*gustdata),3)
gustlat = gustlats(maxind(gustdata))
gustlon = gustlons(maxind(gustdata))

sitefile  = "maxdewp1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
dplats  = tofloat(str_get_field(sitedata,1,delim))
dplons  = tofloat(str_get_field(sitedata,2,delim))
dpdata  = tofloat(str_get_field(sitedata,3,delim))

dp = .1*round(max(10*dpdata),3)
dplat = dplats(maxind(dpdata))
dplon = dplons(maxind(dpdata))

sitefile  = "precip1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
preclats  = tofloat(str_get_field(sitedata,1,delim))
preclons  = tofloat(str_get_field(sitedata,2,delim))
precdata  = tofloat(str_get_field(sitedata,3,delim))

prec = .01*round(max(100*precdata),3)
preclat = preclats(maxind(precdata))
preclon = preclons(maxind(precdata))

sitefile  = "mintemp1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
mintlats  := tofloat(str_get_field(sitedata,1,delim))
mintlons  := tofloat(str_get_field(sitedata,2,delim))
mintdata  := tofloat(str_get_field(sitedata,3,delim))

mint = .1*round(min(10*mintdata),3)
mintlat = mintlats(minind(mintdata))
mintlon = mintlons(minind(mintdata))

sitefile  = "maxtemp1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
maxtlats  := tofloat(str_get_field(sitedata,1,delim))
maxtlons  := tofloat(str_get_field(sitedata,2,delim))
maxtdata  := tofloat(str_get_field(sitedata,3,delim))

maxt = .1*round(max(10*maxtdata),3)
maxtlat = maxtlats(maxind(maxtdata))
maxtlon = maxtlons(maxind(maxtdata))

sitefile  = "maxhi1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
maxhilats  := tofloat(str_get_field(sitedata,1,delim))
maxhilons  := tofloat(str_get_field(sitedata,2,delim))
maxhidata  := tofloat(str_get_field(sitedata,3,delim))

maxhi = .1*round(max(10*maxhidata),3)
maxhilat = maxhilats(maxind(maxhidata))
maxhilon = maxhilons(maxind(maxhidata))

sitefile  = "maxwc1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
minwclats  := tofloat(str_get_field(sitedata,1,delim))
minwclons  := tofloat(str_get_field(sitedata,2,delim))
minwcdata  := tofloat(str_get_field(sitedata,3,delim))

minwc = .1*round(min(10*minwcdata),3)
minwclat = minwclats(minind(minwcdata))
minwclon = minwclons(minind(minwcdata))

sitefile  = "minpres1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
minslplats  := tofloat(str_get_field(sitedata,1,delim))
minslplons  := tofloat(str_get_field(sitedata,2,delim))
minslpdata  := tofloat(str_get_field(sitedata,3,delim))

minslp = round(min(minslpdata),3)
minslplat = minslplats(minind(minslpdata))
minslplon = minslplons(minind(minslpdata))

sitefile  = "maxpres1day.dat"
sitedata  := asciiread(sitefile,-1,"string")
maxslplats  := tofloat(str_get_field(sitedata,1,delim))
maxslplons  := tofloat(str_get_field(sitedata,2,delim))
maxslpdata  := tofloat(str_get_field(sitedata,3,delim))

maxslp = round(max(maxslpdata),3)
maxslplat = maxslplats(maxind(maxslpdata))
maxslplon = maxslplons(maxind(maxslpdata))


;if fileexists("current3hsnowchange.dat")
;     sitefile  = "current3hsnowchange.dat"
;else 
;     sitefile  = "currentmissingsnow.dat"
;end if

;sitedata  := asciiread(sitefile,-1,"string")
;snowlats  = tofloat(str_get_field(sitedata,1,delim))
;snowlons  = tofloat(str_get_field(sitedata,2,delim))
;snowdata  = tofloat(str_get_field(sitedata,3,delim))*1000

mres                = True         ; marker resources for best track
mres@gsMarkerIndex  = 16           ; marker style (filled circle)
mres@gsMarkerSizeF  = 11.0          ; marker size
mres@gsMarkerColor  := "gray50"      ; marker color
mres@tfPolyDrawOrder= "PostDraw"

tres                      = True                ; text mods desired
tres@tfPolyDrawOrder  = "PostDraw"
tres@txFont = "helvetica-bold"
tres@txFontHeightF = .01 
;gustmarker = gsn_add_polymarker(wks,plot,gustlon,gustlat,mres)
;precmarker = gsn_add_polymarker(wks,plot,preclon,preclat,mres)
mres@gsMarkerSizeF  = 32.0          ; marker size

  ;  mres@gsMarkerOpacityF =.7 ;
    mres@gsMarkerColor  := "red"
    maxtmarker = gsn_add_polymarker(wks,plot,maxtlon+.0,maxtlat-.0,mres)
;    gusttext  =     gsn_add_text(wks,plot,round(gustdata(n),3),gustlons(n)+.0,gustlats(n)-.07,tres)

 ;   mres@gsMarkerOpacityF =.7 ;
    mres@gsMarkerColor  := "hot pink"
    maxhimarker = gsn_add_polymarker(wks,plot,maxhilon+.0,maxhilat-.0,mres)
;    gusttext  =     gsn_add_text(wks,plot,round(gustdata(n),3),gustlons(n)+.0,gustlats(n)-.07,tres)
if prec.gt.(.005)
 ;   mres@gsMarkerOpacityF =  .7 ; 
    mres@gsMarkerColor  := "green"
    precmarker = gsn_add_polymarker(wks,plot,preclon-.0,preclat+.0,mres)
;    prectext  =     gsn_add_text(wks,plot,.1*round(10*precdata(n),3),preclons(n)-.07,preclats(n)+.07,tres)
end if
 ;   mres@gsMarkerOpacityF =  .7 ; 
    mres@gsMarkerColor  := "orange"
    gustmarker = gsn_add_polymarker(wks,plot,gustlon-.0,gustlat+.0,mres)
;    prectext  =     gsn_add_text(wks,plot,.1*round(10*precdata(n),3),preclons(n)-.07,preclats(n)+.07,tres)

  ;  mres@gsMarkerOpacityF = .7
    mres@gsMarkerColor  := "blue"
    mintmarker = gsn_add_polymarker(wks,plot,mintlon+.0,mintlat+.0,mres)
 ;   snowtext  =     gsn_add_text(wks,plot,round(snowdata(n),3),snowlons(n)+.07,snowlats(n)+.07,tres)

if minwc.lt.50
  ;  mres@gsMarkerOpacityF = .7
    mres@gsMarkerColor  := "skyblue"
    wcmarker = gsn_add_polymarker(wks,plot,minwclon+.0,minwclat+.0,mres)
 ;   snowtext  =     gsn_add_text(wks,plot,round(snowdata(n),3),snowlons(n)+.07,snowlats(n)+.07,tres)
end if

  ;  mres@gsMarkerOpacityF = .7
    mres@gsMarkerColor  := "seagreen"
    dpmarker = gsn_add_polymarker(wks,plot,dplon+.0,dplat+.0,mres)
 ;   snowtext  =     gsn_add_text(wks,plot,round(snowdata(n),3),snowlons(n)+.07,snowlats(n)+.07,tres)

  ;  mres@gsMarkerOpacityF = .7
    mres@gsMarkerColor  := "brown"
    maxslpmarker = gsn_add_polymarker(wks,plot,maxslplon+.0,maxslplat+.0,mres)
 ;   snowtext  =     gsn_add_text(wks,plot,round(snowdata(n),3),snowlons(n)+.07,snowlats(n)+.07,tres)

  ;  mres@gsMarkerOpacityF = .7
    mres@gsMarkerColor  := "purple"
    minslpmarker = gsn_add_polymarker(wks,plot,minslplon+.0,minslplat+.0,mres)
 ;   snowtext  =     gsn_add_text(wks,plot,round(snowdata(n),3),snowlons(n)+.07,snowlats(n)+.07,tres)

    maxttext  =     gsn_add_text(wks,plot,maxt,maxtlon+.0,maxtlat-.0,tres)
    maxhitext =     gsn_add_text(wks,plot,maxhi,maxhilon+.0,maxhilat-.0,tres)
if prec.gt.(.005)
    prectext  =     gsn_add_text(wks,plot,prec,preclon-.0,preclat+.0,tres)
end if

    gusttext  =     gsn_add_text(wks,plot,gust,gustlon-.0,gustlat+.0,tres)
if minwc.lt.50
    minwctext  =     gsn_add_text(wks,plot,minwc,minwclon+.0,minwclat+.0,tres)
end if
    dptext  =     gsn_add_text(wks,plot,dp,dplon+.0,dplat+.0,tres)
    tres@txFontColor := "White"
    minttext  =     gsn_add_text(wks,plot,mint,mintlon+.0,mintlat+.0,tres)
    maxslptext  =     gsn_add_text(wks,plot,maxslp,maxslplon+.0,maxslplat-.0,tres)
    minslptext  =     gsn_add_text(wks,plot,minslp,minslplon+.0,minslplat-.0,tres)


;mres@gsMarkerSizeF  = 25.0          ; marker size
;mres@gsMarkerColor  := "gray50"
;bigwstationmarker = gsn_add_polymarker(wks,plot,-76,41.22,mres)
;mres@gsMarkerOpacityF = .75
;mres@gsMarkerSizeF  = 80.0          ; marker size
;mres@gsMarkerColor  := "green"
;bigprecmarker = gsn_add_polymarker(wks,plot,-76-.41,41.2+.28,mres)
;mres@gsMarkerColor  := "blue"
;bigsnowmarker = gsn_add_polymarker(wks,plot,-76+.41,41.2+.28,mres)
;mres@gsMarkerColor  := "red"
;biggustmarker = gsn_add_polymarker(wks,plot,-76,41.2-.45,mres)

;tres@txFontHeightF = .008
;tres@txFontColor := "Black"
;ptext1  =     gsn_add_text(wks,plot,"3 Hour",-76-.41,41.2+.43,tres)
;ptext2  =     gsn_add_text(wks,plot,"Precipitation",-76-.41,41.2+.27,tres)
;ptext3  =     gsn_add_text(wks,plot,">1"+quote,-76-.41,41.2+.13,tres)
;tres@txFontColor := "White"
;stext1  =     gsn_add_text(wks,plot,"3 Hour",-76+.41,41.2+.43,tres)
;stext2  =     gsn_add_text(wks,plot,"Snow Depth",-76+.41,41.2+.27,tres)
;stext3  =     gsn_add_text(wks,plot,"Change >2"+quote,-76+.41,41.2+.13,tres)
;tres@txFontColor := "Black"
;wtext1  =     gsn_add_text(wks,plot,"3 Hour",-76,41.2-.3,tres)
;wtext2  =     gsn_add_text(wks,plot,"Wind Gust",-76,41.2-.45,tres)
;wtext3  =     gsn_add_text(wks,plot,">30 mph",-76,41.2-.6,tres)

txres               = True                            
txres@txFontHeightF = .023          ; Set the font height
txres@txJust        = "CenterCenter"
txres@txFont = "helvetica-bold"
label = "Daily NYS Mesonet Extremes"
gsn_text_ndc(wks,label,.36,.36 ,txres)

txres@txFontHeightF = .02 
label = "On " + systemfunc("date -u -d -2days +%A") + ", " + systemfunc("date -u -d -2days +%F")
gsn_text_ndc(wks,label,.36,.32 ,txres)

txres@txFontHeightF = .018 
txres@txJust        = "CenterLeft"
gsn_text_ndc(wks,"State Maximums",.77,.69 ,txres)
txres@txFontHeightF = .014 
txres@txFontColor := "Red"
gsn_text_ndc(wks,"Temperature: "+maxt+"~S~o~N~F",.77,.658 ,txres)
txres@txFontColor := "seagreen"
gsn_text_ndc(wks,"Dew Point: "+dp+"~S~o~N~F",.77,.63 ,txres)
txres@txFontColor := "brown"
gsn_text_ndc(wks,"Pressure: "+maxslp+" mb",.77 ,.60 ,txres)
txres@txFontColor := "orange"
gsn_text_ndc(wks,"Wind Gust: "+gust+" mph",.77,.568 ,txres)

if maxhi.gt.80.and.prec.gt.(.005)
txres@txFontColor := "hot pink"
gsn_text_ndc(wks,"Heat Index: "+maxhi+"~S~o~N~F",.77,.54 ,txres)
txres@txFontColor := "green"
gsn_text_ndc(wks,"Precipitation: "+prec+" in.",.77,.508 ,txres)
elseif  maxhi.lt.80.and.prec.gt.(.005)
txres@txFontColor := "green"
gsn_text_ndc(wks,"Precipitation: "+prec+" in",.77,.538 ,txres)
elseif  maxhi.gt.80.and.prec.lt.(.005)
txres@txFontColor := "hot pink"
gsn_text_ndc(wks,"Heat Index: "+maxhi+"~S~o~N~F",.77,.54 ,txres)

elseif  maxhi.lt.80.and.prec.lt.(.005)
   here = "yes"
end if

txres@txFontHeightF = .018 
txres@txJust        = "CenterLeft"
txres@txFontColor := "black"
gsn_text_ndc(wks,"State Minimums",.77,.46 ,txres)
txres@txFontHeightF = .014 
txres@txFontColor := "Blue"
gsn_text_ndc(wks,"Temperature: "+mint+"~S~o~N~F",.77,.428 ,txres)
txres@txFontColor := "purple"
gsn_text_ndc(wks,"Pressure: "+minslp+" mb",.77 ,.40 ,txres)
if minwc.lt.50
txres@txFontColor := "skyblue"
gsn_text_ndc(wks,"Wind Chill: "+minwc+"~S~o~N~F",.77 ,.37 ,txres)
end if

txres@txFontHeightF = .007
txres@txFontColor := "Blue"
label = "Snow Depth Change Can Be Either Snowfall Or Snow Drifting"; - Strong Winds Can Affect These Estimates"
;gsn_text_ndc(wks,label,.10,.25 ,txres)

draw(plot)
frame(wks)

end


