;*************************************************
; demonstrates the plotting of magnetic variables
;************************************************
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"   
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"    
load "/home/tgcm/ncl/web_ex/func_lib/hao_contrib.ncl"
load "/home/tgcm/ncl/web_ex/func_lib/axis.ncl"
;************************************************
begin
;************************************************
; before the data can be plotted, there is a 
; sequence of things that have to be done to it:
;************************************************
  in    = addfile("/chewy/d/murphys/files/mag_june95_migdisdi_wmc.nc","r")
  x     = mag2latlon(in->jmr)
  u     = mag2latlon(in->kqphi)
  v     = mag2latlon(in->kqlam)
; note the (/.../) syntax here means read in w/o any meta data. since 
; mtime has a back coordinate variable, reading it in directly produces
; some annoying error messages.
  mtime = (/in->mtime/)        
  t2plt = 14 
  ut    = mtime(t2plt,1)                      ; no ut on file
;************************************************
; create plot
;************************************************
  wks = gsn_open_wks("ps" ,"ce")              ; open a ps file
  gsn_define_colormap(wks,"BlWhRe") 

  res                          = True         ; plot mods desired
  res@cnFillOn                 = True         ; add color
  res@gsnSpreadColors          = True         ; use full range of color map
  res@lbLabelAutoStride        = True         ; nice lb labels
  res@mpFillOn                 = False        ; turn off map fill

  title = x@long_name+" ("+x@units+") day= "+ mtime(t2plt,0)+","+\
  mtime(t2plt,1)+","+ mtime(t2plt,2)          
  res@tiMainString             = title  

; these resources must exist when plotting on magnetic coordinates
  res@mpOutlineOn              = False
  res@mpFillOn                 = False
  res@gsnAddCyclic             = False

; these 5 resources zoom in on the map and increase the label size. The
; labels are set for a full globe, and sometime when we zoom in, the
; labels shrink and must be enlarged.
  res@mpMinLatF                = -60.
  res@mpMaxLatF                = 60.
  res@tmXBLabelFontHeightF     = 0.014        ; resize tick labels
  res@tmYLLabelFontHeightF     = 0.014
  res@tiMainFontHeightF        = 0.022        ; resize title

; set explict contour levels
  res@cnLevelSelectionMode     = "ExplicitLevels" 
  res@cnLevels                 =  0.000000001*ispan(-10,10,2)

; turn on local time axis and put mlt noon in the center
  res  = set_lt_axis_ceplot(ut,True,res)        ; turn on local time axis

  res@pmLabelBarOrthogonalPosF  =  0.25         ; mv label bar. room for title
  res@tiXAxisOffsetYF           = .11           ; move xaxis title 

  res@gsnScalarContour          = True          ; contours desired
  res@cnFillDrawOrder           = "Draw"        ; draw on top of map
  res@vcVectorDrawOrder         = "PostDraw"    ; draw vectors last
  res@vcRefMagnitudeF           = 4e-02         ; define vector ref mag
  res@vcRefLengthF              = 0.04          ; define length of vec ref
  res@vcMinDistanceF            = 0.02          ; thin vectors
  res@vcRefAnnoString2On        = True          ; turn on vec ref string
  res@vcRefAnnoString2          = u@units       ; string for ref vec
  res@vcRefAnnoArrowLineColor   = "black"       ; change ref vector color
  res@vcRefAnnoArrowUseVecColor = False         ; don't use vec color for ref
  res@vcLineArrowColor          = "white"       ; change vector color
  res@vcRefAnnoOrthogonalPosF   = -0.27         ; move reference vec up
  
  plot = gsn_csm_vector_scalar_map_ce(wks,u(t2plt,:,:),v(t2plt,:,:),\
  x(t2plt,:,:),res)
end