;; -*- mode: emacs-lisp; coding: emacs-mule; -*-;  -*- mode: 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/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl" 

ana=addfile("nh4_w_w_co3.nc","r")
no3_1_w=ana->nh4_1(:,:,61:148,36:168)
no3_2_w=ana->nh4_2(:,:,61:148,36:168)
no3_3_w=ana->nh4_3(:,:,61:148,36:168)
no3_4_w=ana->nh4_4(:,:,61:148,36:168)
delete(ana)


ana=addfile("nh4_w_wo_co3.nc","r")
 no3_1_wo=ana->nh4_1(:,:,61:148,36:168)
 no3_2_wo=ana->nh4_2(:,:,61:148,36:168)
 no3_3_wo=ana->nh4_3(:,:,61:148,36:168)
 no3_4_wo=ana->nh4_4(:,:,61:148,36:168)
delete(ana)

nox1=(no3_2_w-no3_2_wo)
copy_VarCoords(no3_1_w,nox1)  
nox=dim_avg_n_Wrap(nox1,0)

ana=addfile("wrfout_d01_2015-11-12_21:00:00","r")
lat = wrf_user_getvar(ana,"lat",0)
lon = wrf_user_getvar(ana,"lon",0)
;delete(ana)

lat1=lat(61:148,:)
lon1=lon(:,36:168)


wks = gsn_open_wks("eps","nh4_2_new_flux")
gsn_define_colormap(wks,"BlWhRe")

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 
  res                     = True               ; plot mods desired
  res@cnFillOn            = True              ; turn on color   
  res@gsnSpreadColors     = True               ; use full color map
  res@cnLinesOn           = False              ; no contour lines
  res@cnLineLabelsOn      = False              ; no line labels

  res@gsnMaximize        = True
  res@mpMinLatF          = min(lat1)-1
  res@mpMaxLatF          = max(lat1)+1
  res@mpMinLonF          = min(lon1)-1
  res@mpMaxLonF          = max(lon1)+1

  res@gsnDraw             = False           ; don't draw
  res@gsnFrame            = False           ; don't advance frame
 
  res@gsnAddCyclic        = False 
  res@lbLabelBarOn        = True           ; turn off individual cb's 


;res@lbTitleOn        =  True                ; turn on title
;res@lbTitleString    = "ug/m^3"                ; title string
;res@lbTitlePosition  = "Right"              ; title position
;res@lbTitleFontHeightF= .02                ; make title smaller
;res@lbTitleDirection = "Across"             ; title direction


; res@cnLevelSelectionMode =  "AutomaticLevels"   ; AutomaticLevels

 res@cnLevelSelectionMode =  "ManualLevels"   ; AutomaticLevels
  res@cnMinLevelValF       = -0.064              ; min level
  res@cnMaxLevelValF       =  0.064              ; max level
  res@cnLevelSpacingF      = 0.008              ; interval



  res@mpFillOn                    = False        ; turn off gray fill
  res@mpOutlineBoundarySets       = "National"   ; turn on country boundaries
  res@mpGeophysicalLineColor      = "Black"       ; color of cont. outlines
  res@mpGeophysicalLineThicknessF = 1.5          ; thickness of outlines
  res@mpDataBaseVersion     = "MediumRes"        ; better map resolution

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;res@tiMainString        = "Without het reactions Nh4_4 NEI+Lake - NEI"
res@tiMainString        = "new_flux NH4_2 NEI+Lake - NEI "
res = wrf_map_resources(ana,res)
res@tfDoNDCOverlay   = True  
res@tmXBLabelFontHeightF = 0.02                 ; resize tick labels
res@tmYLLabelFontHeightF = 0.02
res@lbLabelFontHeightF = 0.02

plot= gsn_csm_contour_map(wks,nox(0,:,:),res)


draw (plot)
frame(wks)
delete(wks)