;*************************************************
;plot annual surface temperature stsu
;Qian Zhang @ 11/11/2011
;************************************************
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"

begin

; READ land areas
  fm = addfile("csiro_landsea.nc","r")                 ; read the mask
  f0 = addfile("csiro_agcm_area.nc","r")               ; read grid area
  oro  = fm->mask
  area = f0->area                                      ; unit in m^2
  asum = sum(area)                                      

; READ out_cable*.nc 

  f = addfile("/short/p66/yxw599/coal_1/run1/stsu_testc.nc","r")
  Tmn = f->tsu
  yr  = f->year

  startyr = 1
  endyr   = dimsizes(yr)

  Tann = dim_avg_n(Tmn,1)
  T    = dim_sum_n(Tann*conform(Tann,area,(/1,2/)),(/1,2/))/asum

  T@units = Tmn@units
;  printVarSummary (T)

  y = ispan(startyr,endyr,1)
  
;---------------------------------------------------------------
; graph plotting

  wks      = gsn_open_wks("pdf","tsu_coal_20111213")

  res                  = True
  res@gsnMaximize      = True
;  res@tiMainString     = "Surface Air Temperature"  ; add title
  res@vpWidthF         = 0.8                  ; set width of plot
  res@vpHeightF        = 0.6                  ; set height of plot
  res@xyLineThicknessF = 3
;  res@xyDashPatterns   = (/0,0,0,1,1,1/)
;  res@xyLineColors     = (/"black","darkgreen","red","black","darkgreen","red"/) 
  res@txFontHeightF    = 0.15                  ; change title font heights
  res@trXMinF          = startyr               ; y axis min
  res@trXMaxF          = endyr                 ; y axis max
  res@trYMinF          = 285.8                 ; y axis min
  res@trYMaxF          = 286.4                 ; y axis max
  res@tiXAxisString    = "year"
  res@tiYAxisString    = "K"

; Manually create and attach legend
  res@pmLegendDisplayMode    = "Always"

  res@pmLegendSide           = "Top"               ; Change location of
  res@pmLegendParallelPosF   = .2                  ; move units right
  res@pmLegendOrthogonalPosF = -0.25                ; move units down

  res@pmLegendWidthF         = 0.12                ; Change width and
  res@pmLegendHeightF        = 0.14                ; height of legend.
  res@lgLabelFontHeightF     = .03                ; change font height
;  res@lgOrientation         = "horizontal"
  res@lgPerimOn              = False               ; no box around
  res@xyExplicitLegendLabels = (/"  Ctl-tsu"/)

  plot = gsn_csm_xy(wks,y,T,res)

end