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"

begin

; This program correlated observed monthly temperature and precipitation
; for 1982-2000.  The correlations are shaded over the United States
; using grey shadings/hatchings.  Shown only over ocean here.

 a = addfile("/ptmp/notaro/reanal/air.mon.mean_2.5x2.5_jan82dec00.nc","r")  ; open data files
 b = addfile("/ptmp/notaro/xiearkin/XieArkin-2.5x2.5-1982to2000.nc","r")
 lat = a->lat            ; read in lat, lon
 lon = a->lon
 times=ispan(1,228,1)    ; create a time loop
 sfct = a->air           ; read in temperature
 pcp = (b->PRC2)*30./10. ; read in precip and convert to cm/mon 

 sfct!0="time"           ; assign attributes
 sfct!1="lat"
 sfct!2="lon"
 sfct&time=times
 sfct&lat=lat
 sfct&lon=lon
 copy_VarCoords(sfct,pcp)   ; copy attributes from temp to pcp

;;;;;;;;;;;;;;;

  correl=new((/73,144/),float)
  correl=escorc(sfct(lat|:,lon|:,time|:),pcp(lat|:,lon|:,time|:))    ; correlate temp and pcp
  correl!0="lat"
  correl!1="lon"
  correl&lat=lat
  correl&lon=lon
  correl@_FillValue=1e+35

  wks = gsn_open_wks("ps","correl_sfct_pcp_bw")          ; open a ps file
  plot = new(1,graphic)
  setvalues wks
    "wkColorMap"        : "gsltod"             ; color map specified
    "wkForegroundColor" : (/0.,0.,0./)
    "wkBackgroundColor" : (/1.,1.,1./)
  end setvalues
  res                  = True                
  res@gsnDraw         = False
  res@gsnFrame        = False          
  res@cnFillOn         = True         
  res@cnLinesOn        =  True       
  res@cnLevelSelectionMode = "ExplicitLevels" 
  res@cnLevels             = (/-.3,-.1,.1,.3/)   ; contour levels
  res@cnMonoFillPattern    = False            
  res@cnFillPatterns       = (/17,17,-1,0,0/)    ; shading patterns
  res@cnMonoFillScale      = False            
  res@cnFillScales         = (/.45,.75,1,.3,.1/) 
  res@cnFillColors = (/30,30,-1,10,22/)          ; shading colors
  res@tiMainFontHeightF =  0.016
  res@tiMainString         = "R(Mn Reanal Sfc Air T, Mn Xie Arkin Pcp)"
  res@tmYROn = False
  res@tmXTOn = False
  res@cnFillDrawOrder      = "PreDraw"      ; this and the next 4 lines cover up the ocean
  res@cnLineDrawOrder = "PreDraw"
  res@mpOceanFillColor     = "white"
  res@mpLandFillColor = -1
  res@pmLabelBarDisplayMode = True
  res@lbLabelBarOn = True
  res@lbOrientation        = "vertical"
  res@mpMinLatF            = 25.                 ; set area to U.S.
  res@mpMaxLatF            = 49.5
  res@mpMinLonF = -126.
  res@mpMaxLonF = -66.
  res@mpOutlineBoundarySets = "AllBoundaries"    ; draw state lines
  res@mpGeophysicalLineThicknessF = 2.5
  res@mpOutlineOn = True
  plot(0) = gsn_csm_contour_map_ce(wks,correl,res)

  pres             = True                
  pres@gsnDraw=True
  pres@gsnFrame=True
  pres@gsnMaximize = True               
  pres@gsnPaperWidth=8.5
  pres@gsnPaperHeight=11.0
  pres@gsnPaperOrientation="portrait"
  pres@gsnPaperMargin=0.04
  pres@gsnPanelLeft = 0.02
  pres@gsnPanelRight = 0.98
  pres@gsnPanelBottom = 0.02
  pres@gsnPanelTop = 0.96
  pres@gsnPanelYWhiteSpacePercent = 3
  pres@gsnPanelXWhiteSpacePercent = 3
  pres@gsnPanelRowSpec=True
  pres@txFontHeightF = .014
  gsn_panel(wks,plot,(/1/),pres)        ; 1 panel

end