diff -r 000000000000 -r 0c6405ab2ff4 lai/32.contour_diff.ncl
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/lai/32.contour_diff.ncl	Mon Jan 26 22:08:20 2009 -0500
@@ -0,0 +1,87 @@
+;*************************************************
+; ce_1.ncl
+;************************************************
+load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"  
+load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"  
+;************************************************
+begin
+;************************************************
+; read in observed data
+;************************************************
+ f     = addfile ("Npp_T42_mean.nc","r")
+ y     = f->NPP
+;************************************************
+; read in model data
+;************************************************
+ g     = addfile ("i01.03cn_1545-1569_ANN_climo.nc","r")
+;g     = addfile ("i01.04casa_1605-1629_ANN_climo.nc","r")
+ x     = g->NPP
+
+ nsec_per_year = 86400.*365.
+ x = x * nsec_per_year
+ 
+ x@units = "gC/m^2/year" 
+;************************************************
+;fill ob missing grid the same as model
+;************************************************
+ delta = 0.00001
+ x0    = x(0,:,:)
+ y = where(ismissing(y).and.(ismissing(x0).or.(x0.lt.delta)),0.,y)
+
+ y@units = "gC/m^2/year"
+;************************************************
+;model - observed
+;************************************************ 
+ z = x
+ z = x(0,:,:) - y(:,:)
+;************************************************
+; create 3 plots
+;************************************************
+  wks = gsn_open_wks("ps","xy")             ; open a ps file
+  gsn_define_colormap(wks,"gui_default")     ; choose colormap
+
+  res                     = True             ; Use plot options
+  res@cnFillOn            = True             ; Turn on color fill
+  res@gsnSpreadColors      = True            ; use full colormap
+; res@cnFillMode          = "RasterFill"     ; Turn on raster color
+; res@lbLabelAutoStride   = True
+  res@cnLinesOn           = False            ; Turn off contourn lines
+  res@mpFillOn            = False            ; Turn off map fill
+  res@tiMainString        = "Observed MODIS MOD 17"
+
+  res@gsnSpreadColors     = True             ; use full colormap
+  res@cnLevelSelectionMode = "ManualLevels"  ; Manual contour invtervals
+  res@cnMinLevelValF       = 0.              ; Min level
+  res@cnMaxLevelValF       = 2200.           ; Max level
+  res@cnLevelSpacingF      = 200.            ; interval
+
+  res@gsnFrame             = False           ; Do not draw plot 
+  res@gsnDraw              = False           ; Do not advance frame
+
+  plot=new(3,graphic)                        ; create graphic array
+
+  plot(0) = gsn_csm_contour_map_ce(wks,y,res)        ; for observed
+
+  res@tiMainString        = "Model i01.03cn"
+; res@tiMainString        = "Model i01.04casa"
+  plot(1) = gsn_csm_contour_map_ce(wks,x(0,:,:),res) ; for model
+
+  res@cnMinLevelValF       = -500           ; Min level
+  res@cnMaxLevelValF       =  500.          ; Max level
+  res@cnLevelSpacingF      = 50.            ; interval
+  res@tiMainString        = "(Model i01.03cn) - (observed)"
+; res@tiMainString        = "(Model i01.04casa) - (observed)"
+  plot(2) = gsn_csm_contour_map_ce(wks,z(0,:,:),res) ; for model - ob
+
+;***********************************************
+; create panel plot
+;***********************************************
+  pres                            = True        ; panel plot mods desired
+  pres@gsnPanelYWhiteSpacePercent = 5           ; increase white space around
+                                                ; indiv. plots in panel
+  pres@gsnMaximize                = True        ; fill the page
+
+  gsn_panel(wks,plot,(/3,1/),pres)              ; create panel plot
+
+  system("convert xy.ps xy.png")
+end
\ No newline at end of file