1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/lai/31.contour_ob_lai_max.ncl.1 Mon Jan 26 22:08:20 2009 -0500
1.3 @@ -0,0 +1,95 @@
1.4 +;*************************************************
1.5 +; ce_1.ncl
1.6 +;************************************************
1.7 +load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
1.8 +load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
1.9 +load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
1.10 +;************************************************
1.11 +begin
1.12 +;************************************************
1.13 +; read in observed data
1.14 +;************************************************
1.15 + diri = "/fis/cgd/cseg/people/jeff/clamp_data/"
1.16 +; fili = "LAI_2000-2005_ensemble_T42.nc"
1.17 + fili = "2000-2005_ensemble_T42.nc"
1.18 + f = addfile(diri+fili,"r")
1.19 +
1.20 + z = f->LAI
1.21 + y = z(0,:,:)
1.22 + y@long_name = "Leaf Area Index Max"
1.23 + p = z(0,:,:)
1.24 + p@long_name = "Leaf Area Index Max Month"
1.25 + s = z(:,0,0)
1.26 +; printVarSummary(s)
1.27 +
1.28 + dsizes_z = dimsizes(z)
1.29 + ntime = dsizes_z(0)
1.30 + nlat = dsizes_z(1)
1.31 + nlon = dsizes_z(2)
1.32 +
1.33 + do j = 0,nlat-1
1.34 + do i = 0,nlon-1
1.35 + s = z(:,j,i)
1.36 + y(j,i) = max(s)
1.37 + p(j,i) = maxind(s)+ 1
1.38 +
1.39 +; iMax = ind(s.eq.y(j,i)) + 1 ; may be multiple max
1.40 +; if (.not.ismissing(iMax(0))) then
1.41 +; p(j,i) = iMax(0) ; no meta data
1.42 +; end if
1.43 +; delete(iMax)
1.44 +
1.45 + end do
1.46 + end do
1.47 +
1.48 +; printVarSummary(p)
1.49 + print (min(y)+"/"+max(y))
1.50 + print (min(p)+"/"+max(p))
1.51 +exit
1.52 +;************************************************
1.53 +; read in model data
1.54 +;************************************************
1.55 + fili2 = "i01.03cn_1545-1569_ANN_climo.nc"
1.56 + g = addfile(diri+fili2,"r")
1.57 + x = g->TLAI
1.58 +
1.59 + delta = 0.00000000001
1.60 + x0 = x(0,:,:)
1.61 + y = where(ismissing(y).and.(ismissing(x0).or.(x0.lt.delta)),0.,y)
1.62 +;************************************************
1.63 +; create default plot
1.64 +;************************************************
1.65 +
1.66 +; setvalues NhlGetWorkspaceObjectId()
1.67 +; "wsMaximumSize" : 199999999
1.68 +; end setvalues
1.69 +
1.70 + wks = gsn_open_wks("ps","xy") ; open a ps file
1.71 +; gsn_define_colormap(wks,"wgne15") ; choose colormap
1.72 + gsn_define_colormap(wks,"gui_default") ; choose colormap
1.73 +
1.74 + res = True ; Use plot options
1.75 + res@cnFillOn = True ; Turn on color fill
1.76 + res@gsnSpreadColors = True ; use full colormap
1.77 +; res@cnFillMode = "RasterFill" ; Turn on raster color
1.78 +; res@lbLabelAutoStride = True
1.79 + res@cnLinesOn = False ; Turn off contourn lines
1.80 + res@mpFillOn = False ; Turn off map fill
1.81 + res@tiMainString = "MODIS MOD 15A2 2000-2005"
1.82 +
1.83 + res@gsnSpreadColors = True ; use full colormap
1.84 + res@cnLevelSelectionMode = "ManualLevels" ; Manual contour invtervals
1.85 + res@cnMinLevelValF = 0. ; Min level
1.86 + res@cnMaxLevelValF = 12. ; Max level
1.87 + res@cnLevelSpacingF = 1. ; interval
1.88 +
1.89 +; pres = True ; panel plot mods desired
1.90 +; pres@gsnMaximize = True ; fill the page
1.91 +
1.92 + plot = gsn_csm_contour_map_ce(wks,y,res) ; for observed
1.93 +
1.94 + plot = gsn_csm_contour_map_ce(wks,p,res) ; for observed
1.95 +
1.96 +; gsn_panel(wks,plot,(/1,1/),pres) ; create panel plot
1.97 + system("convert xy.ps xy.png")
1.98 +end
1.99 \ No newline at end of file