;----------------------------------------------------------------------
; wrf_gsn_6.ncl
;----------------------------------------------------------------------
; Concepts illustrated:
;   - Using gsn_csm scripts to plot WRF-ARW data
;   - Overlaying filled contours, and vectors on a map
;   - Setting the correct WRF map projection using wrf_map_resources
;   - Setting lots of resources to customize a WRF plot
;----------------------------------------------------------------------
; This script is meant to show the difference between plotting WRF
; data using wrf_xxx scripts, and using gsn_csm_xxx scripts.
;----------------------------------------------------------------------
; You can safely ignore these warning messages:
;
; warning:start_lat is not a valid resource in wrf_gsn_contour at this time
; warning:start_lon is not a valid resource in wrf_gsn_contour at this time
; warning:end_lat is not a valid resource in wrf_gsn_contour at this time
; warning:end_lon is not a valid resource in wrf_gsn_contour at this time
; warning:mpNestTime is not a valid resource in map at this time
;----------------------------------------------------------------------
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"


begin

;---Parameters
  input_file         = "wrfout_d01_2015-02-05_00:00:00.nc"
  input_dir          = "/glade/u/home/mmasarik/WRF-NCL/input"
  input_path         = input_dir + "/" +  input_file

  plot_name          = "streamlines_station-points"



;---Open WRF output file
  a                  = addfile(input_path,"r")
  hgt = a->HGT(0,:,:)
  hgt@lat2d = a->XLAT(0,:,:)
  hgt@lon2d = a->XLONG(0,:,:)


;--Create Workstation
  wks                = gsn_open_wks("png","wrf_map_nonnative")


;--Start resources
  res                                = True
  res@gsnLeftString                  = ""
  res@gsnRightString                 = ""
  res@gsnAddCyclic                   = False

;---Map plot
  map_res                             = res
;  map_res                             = wrf_map_resources(a,map_res)
  map_res@gsnStringFontHeightF        = 0.015
  map_res@tiMainJust                  = "CenterLeft"
  map_res@tiMainPosition              = "Left"
  map_res@tiMainOffsetYF              = -0.01


;---Set map resources based on projection on WRF output file
  map_res@mpFillOn                    = False
  map_res@mpOutlineOn                 = True
  map_res@mpDataBaseVersion           = "MediumRes"
  map_res@mpOutlineDrawOrder          = "PostDraw"
  map_res@mpUSStateLineColor          = "Black"
  map_res@mpPerimLineColor            = "Black"
  map_res@mpNationalLineColor         = "Black"
  map_res@mpLimbLineColor             = "Black"
  map_res@mpGridLineColor             = "Black"
  map_res@mpGeophysicalLineColor      = "Black"
  map_res@mpUSStateLineThicknessF     = 3.0
  map_res@mpNationalLineThicknessF    = 3.0
  map_res@mpGeophysicalLineThicknessF = 3.0

  map_res@mpMinLatF = min(hgt@lat2d)-1
  map_res@mpMaxLatF = max(hgt@lat2d)+1
  map_res@mpMinLonF = min(hgt@lon2d)-1
  map_res@mpMaxLonF = max(hgt@lon2d)+1

  map_res@cnFillOn = True
  map_res@cnFillPalette = "OceanLakeLandSnow"
  map_res@cnLevelSelectionMode = "ExplicitLevels"
  map_res@cnLevels = (/5,20,50,100,150,200,250,300,350,400,450,500,550,600,650,700,750,800,850,900,950,1000,1025/)
  map_res@pmTickMarkDisplayMode = "Always"    ; nicer tickmarks
  map = gsn_csm_contour_map(wks,hgt,map_res)

  mkres = True
  mkres@gsMarkerIndex = 16
  mkres@gsMarkerSizeF = 2
  gsn_coordinates(wks,map,hgt,mkres)

end