npp/32.contour_diff.ncl
author Forrest Hoffman <forrest@climatemodeling.org>
Mon, 26 Jan 2009 22:08:20 -0500
changeset 0 0c6405ab2ff4
permissions -rw-r--r--
Initial commit of C-LAMP Diagnostics from Jeff Lee
forrest@0
     1
;*************************************************
forrest@0
     2
; ce_1.ncl
forrest@0
     3
;************************************************
forrest@0
     4
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"  
forrest@0
     5
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"  
forrest@0
     6
;************************************************
forrest@0
     7
begin
forrest@0
     8
;************************************************
forrest@0
     9
; read in observed data
forrest@0
    10
;************************************************
forrest@0
    11
 f     = addfile ("Npp_T42_mean.nc","r")
forrest@0
    12
 y     = f->NPP
forrest@0
    13
;************************************************
forrest@0
    14
; read in model data
forrest@0
    15
;************************************************
forrest@0
    16
 g     = addfile ("i01.03cn_1545-1569_ANN_climo.nc","r")
forrest@0
    17
;g     = addfile ("i01.04casa_1605-1629_ANN_climo.nc","r")
forrest@0
    18
 x     = g->NPP
forrest@0
    19
forrest@0
    20
 nsec_per_year = 86400.*365.
forrest@0
    21
 x = x * nsec_per_year
forrest@0
    22
 
forrest@0
    23
 x@units = "gC/m^2/year" 
forrest@0
    24
;************************************************
forrest@0
    25
;fill ob missing grid the same as model
forrest@0
    26
;************************************************
forrest@0
    27
 delta = 0.00001
forrest@0
    28
 x0    = x(0,:,:)
forrest@0
    29
 y = where(ismissing(y).and.(ismissing(x0).or.(x0.lt.delta)),0.,y)
forrest@0
    30
forrest@0
    31
 y@units = "gC/m^2/year"
forrest@0
    32
;************************************************
forrest@0
    33
;model - observed
forrest@0
    34
;************************************************ 
forrest@0
    35
 z = x
forrest@0
    36
 z = x(0,:,:) - y(:,:)
forrest@0
    37
;************************************************
forrest@0
    38
; create 3 plots
forrest@0
    39
;************************************************
forrest@0
    40
  wks = gsn_open_wks("ps","xy")             ; open a ps file
forrest@0
    41
  gsn_define_colormap(wks,"gui_default")     ; choose colormap
forrest@0
    42
forrest@0
    43
  res                     = True             ; Use plot options
forrest@0
    44
  res@cnFillOn            = True             ; Turn on color fill
forrest@0
    45
  res@gsnSpreadColors      = True            ; use full colormap
forrest@0
    46
; res@cnFillMode          = "RasterFill"     ; Turn on raster color
forrest@0
    47
; res@lbLabelAutoStride   = True
forrest@0
    48
  res@cnLinesOn           = False            ; Turn off contourn lines
forrest@0
    49
  res@mpFillOn            = False            ; Turn off map fill
forrest@0
    50
  res@tiMainString        = "Observed MODIS MOD 17"
forrest@0
    51
forrest@0
    52
  res@gsnSpreadColors     = True             ; use full colormap
forrest@0
    53
  res@cnLevelSelectionMode = "ManualLevels"  ; Manual contour invtervals
forrest@0
    54
  res@cnMinLevelValF       = 0.              ; Min level
forrest@0
    55
  res@cnMaxLevelValF       = 2200.           ; Max level
forrest@0
    56
  res@cnLevelSpacingF      = 200.            ; interval
forrest@0
    57
forrest@0
    58
  res@gsnFrame             = False           ; Do not draw plot 
forrest@0
    59
  res@gsnDraw              = False           ; Do not advance frame
forrest@0
    60
forrest@0
    61
  plot=new(3,graphic)                        ; create graphic array
forrest@0
    62
forrest@0
    63
  plot(0) = gsn_csm_contour_map_ce(wks,y,res)        ; for observed
forrest@0
    64
forrest@0
    65
  res@tiMainString        = "Model i01.03cn"
forrest@0
    66
; res@tiMainString        = "Model i01.04casa"
forrest@0
    67
  plot(1) = gsn_csm_contour_map_ce(wks,x(0,:,:),res) ; for model
forrest@0
    68
forrest@0
    69
  res@cnMinLevelValF       = -500           ; Min level
forrest@0
    70
  res@cnMaxLevelValF       =  500.          ; Max level
forrest@0
    71
  res@cnLevelSpacingF      = 50.            ; interval
forrest@0
    72
  res@tiMainString        = "(Model i01.03cn) - (observed)"
forrest@0
    73
; res@tiMainString        = "(Model i01.04casa) - (observed)"
forrest@0
    74
  plot(2) = gsn_csm_contour_map_ce(wks,z(0,:,:),res) ; for model - ob
forrest@0
    75
forrest@0
    76
;***********************************************
forrest@0
    77
; create panel plot
forrest@0
    78
;***********************************************
forrest@0
    79
  pres                            = True        ; panel plot mods desired
forrest@0
    80
  pres@gsnPanelYWhiteSpacePercent = 5           ; increase white space around
forrest@0
    81
                                                ; indiv. plots in panel
forrest@0
    82
  pres@gsnMaximize                = True        ; fill the page
forrest@0
    83
forrest@0
    84
  gsn_panel(wks,plot,(/3,1/),pres)              ; create panel plot
forrest@0
    85
forrest@0
    86
  system("convert xy.ps xy.png")
forrest@0
    87
end