forrest@0: ; ***********************************************
forrest@0: ; xy_4.ncl
forrest@0: ; ***********************************************
forrest@0: load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
forrest@0: load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
forrest@0: load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
forrest@0: ;************************************************
forrest@0: begin
forrest@0: ;************************************************
forrest@0: ; read in observed data
forrest@0: ;************************************************
forrest@0:  f     = addfile ("Npp_T42_mean.nc","r")
forrest@0:  y     = f->NPP
forrest@0: ;************************************************
forrest@0: ; read in model data
forrest@0: ;************************************************
forrest@0:  g     = addfile ("i01.03cn_1545-1569_ANN_climo.nc","r")
forrest@0: ;g     = addfile ("i01.04casa_1605-1629_ANN_climo.nc","r")
forrest@0:  x     = g->NPP
forrest@0: 
forrest@0:  delta = 0.00000000001
forrest@0:  x0    = x(0,:,:)
forrest@0:  y = where(ismissing(y).and.(ismissing(x0).or.(x0.lt.delta)),0.,y)
forrest@0: 
forrest@0:  p     = zonalAve(y)        
forrest@0: 
forrest@0:  sec_to_year = 86400.*365.
forrest@0:  x0 = x0 * sec_to_year
forrest@0:  q = zonalAve(x0)
forrest@0: 
forrest@0:  good = ind(p .ne. 0.)
forrest@0:  u = p(good)
forrest@0:  v = q(good)
forrest@0:  
forrest@0:  ccr = esccr(u,v,0)
forrest@0:  print (ccr)
forrest@0:  bias = sum(((v-u)/u)^2)
forrest@0:  print (bias)
forrest@0:  M    = 1.- sqrt(bias/dimsizes(u))
forrest@0:  print (M)
forrest@0: 
forrest@0: end