lai/31.contour_ob_grow.ncl
author Forrest Hoffman <forrest@climatemodeling.org>
Thu, 26 Mar 2009 14:02:21 -0400
changeset 1 4be95183fbcd
permissions -rw-r--r--
Modifications to scoring and graphics production for the final version of code for the C-LAMP paper in GCB.
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
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"  
forrest@0
     7
;************************************************
forrest@0
     8
begin
forrest@0
     9
;************************************************
forrest@0
    10
; read in observed data
forrest@0
    11
;************************************************
forrest@0
    12
  diri  = "/fis/cgd/cseg/people/jeff/clamp_data/lai/"
forrest@0
    13
  fili  = "LAI_2000-2005_ensemble_T42.nc"
forrest@0
    14
  f     = addfile(diri+fili,"r")
forrest@0
    15
forrest@0
    16
  z     = f->LAI
forrest@0
    17
  y     = z(0,:,:)
forrest@0
    18
  y@long_name = "Days of Growing Season"
forrest@0
    19
 
forrest@0
    20
  dsizes_z = dimsizes(z)
forrest@0
    21
  ntime    = dsizes_z(0)
forrest@0
    22
  nlat     = dsizes_z(1)
forrest@0
    23
  nlon     = dsizes_z(2)
forrest@0
    24
forrest@0
    25
 day_of_data = (/31,28,31,30,31,30,31,31,30,31,30,31/)
forrest@0
    26
forrest@0
    27
  do j = 0,nlat-1
forrest@0
    28
  do i = 0,nlon-1
forrest@0
    29
     nday = 0.
forrest@0
    30
     do k = 0,ntime-1
forrest@0
    31
        if (.not. ismissing(z(k,j,i)) .and. z(k,j,i) .gt. 1.0) then
forrest@0
    32
           nday = nday + day_of_data(k)
forrest@0
    33
        end if
forrest@0
    34
     end do
forrest@0
    35
     y(j,i) = nday
forrest@0
    36
  end do
forrest@0
    37
  end do
forrest@0
    38
forrest@0
    39
  print (min(y)+"/"+max(y))
forrest@0
    40
forrest@0
    41
;************************************************
forrest@0
    42
; read in model data
forrest@0
    43
;************************************************
forrest@0
    44
  diri2  = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
forrest@0
    45
  fili2 = "i01.03cn_1545-1569_ANN_climo.nc"
forrest@0
    46
  g     = addfile(diri2+fili2,"r")
forrest@0
    47
  x     = g->TLAI
forrest@0
    48
forrest@0
    49
  delta = 0.000001
forrest@0
    50
  x0    = x(0,:,:)
forrest@0
    51
  y = where(ismissing(y).and.(ismissing(x0).or.(x0.lt.delta)),y@_FillValue,y)
forrest@0
    52
;************************************************
forrest@0
    53
; create default plot
forrest@0
    54
;************************************************
forrest@0
    55
forrest@0
    56
  wks = gsn_open_wks("ps","xy")         ; open a ps file
forrest@0
    57
  gsn_define_colormap(wks,"gui_default")     ; choose colormap
forrest@0
    58
forrest@0
    59
  res                     = True             ; Use plot options
forrest@0
    60
  res@cnFillOn            = True             ; Turn on color fill
forrest@0
    61
  res@gsnSpreadColors      = True            ; use full colormap
forrest@0
    62
; res@cnFillMode          = "RasterFill"     ; Turn on raster color
forrest@0
    63
; res@lbLabelAutoStride   = True
forrest@0
    64
  res@cnLinesOn           = False            ; Turn off contourn lines
forrest@0
    65
  res@mpFillOn            = False            ; Turn off map fill
forrest@0
    66
  res@tiMainString        = "MODIS MOD 15A2 2000-2005"
forrest@0
    67
forrest@0
    68
  res@gsnSpreadColors     = True             ; use full colormap
forrest@0
    69
  res@cnLevelSelectionMode = "ManualLevels"  ; Manual contour invtervals
forrest@0
    70
  res@cnMinLevelValF       = 0.              ; Min level
forrest@0
    71
  res@cnMaxLevelValF       = 390.            ; Max level
forrest@0
    72
  res@cnLevelSpacingF      = 30.             ; interval
forrest@0
    73
forrest@0
    74
; pres                            = True     ; panel plot mods desired
forrest@0
    75
; pres@gsnMaximize                = True     ; fill the page
forrest@0
    76
forrest@0
    77
  plot = gsn_csm_contour_map_ce(wks,y,res)   ; for observed
forrest@0
    78
forrest@0
    79
; gsn_panel(wks,plot,(/1,1/),pres)           ; create panel plot
forrest@0
    80
  system("convert xy.ps xy.png")
forrest@0
    81
end