1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/lai/03.lai_ob_to_T31.ncl.x Mon Jan 26 22:08:20 2009 -0500
1.3 @@ -0,0 +1,128 @@
1.4 +; ***********************************************
1.5 +; interpolate into model grids (T31)
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 +;************************************************
1.14 +; output data
1.15 +;************************************************
1.16 + diro = "/fis/cgd/cseg/people/jeff/clamp_data/lai/ob/"
1.17 + filo = "LAI_2000-2005_ensemble_T31.nc"
1.18 + c = addfile(diro+filo,"c")
1.19 + filedimdef(c,"time",-1,True)
1.20 +
1.21 +;************************************************
1.22 +; read in observed data
1.23 +;************************************************
1.24 + diri = "/fis/cgd/cseg/people/jeff/clamp_data/lai/ob/"
1.25 + fili = "LAI_2000-2005_ensemble_0.05deg.nc"
1.26 + g = addfile (diri+fili,"r")
1.27 + bi = g->LAI
1.28 + xi = g->lon
1.29 + yi = g->lat
1.30 +
1.31 +;************************************************
1.32 +; change into 0-360E, 90S-90N
1.33 +;************************************************
1.34 +
1.35 + yi = (/ yi(::-1) /)
1.36 + bi = (/ bi(:,::-1,:) /)
1.37 + printVarSummary(bi)
1.38 +
1.39 + b2 = bi
1.40 + x2 = xi
1.41 +
1.42 + nx = dimsizes(xi)
1.43 + do i= 0,nx-1
1.44 + if (i .lt. 3600) then
1.45 + p = i + 3600
1.46 + xi(p) = x2(i) + 360.
1.47 + else
1.48 + p = i - 3600
1.49 + xi(p) = x2(i)
1.50 + end if
1.51 + bi(:,:,p)= b2(:,:,i)
1.52 + end do
1.53 +
1.54 + bi&lat = yi
1.55 + bi&lon = xi
1.56 +
1.57 +;print (xi)
1.58 +;print (yi)
1.59 +;************************************************
1.60 +; read in model data
1.61 +;************************************************
1.62 + diri2 = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
1.63 + fili2 = "b30.061n_1995-2004_MONS_climo_lnd.nc"
1.64 + f = addfile (diri2+fili2,"r")
1.65 +
1.66 + lon = f->lon
1.67 + lat = f->lat
1.68 + nlon = dimsizes(lon)
1.69 + nlat = dimsizes(lat)
1.70 +
1.71 +; print (xi)
1.72 +; print (yi)
1.73 +; print (xo)
1.74 +; print (yo)
1.75 +
1.76 + bo = new((/12,nlat,nlon/),float)
1.77 +
1.78 + do m = 0,11
1.79 + do j=0,nlat-1
1.80 + if (j.eq.0 .or. j.eq.nlat-1) then
1.81 + if (j.eq.0) then
1.82 + LATS = -90.
1.83 + LATN = lat(j)+0.5*(lat(j+1)-lat(j))
1.84 + end if
1.85 + if (j.eq.nlat-1) then
1.86 + LATS = lat(j)-0.5*(lat(j)-lat(j-1))
1.87 + LATN = 90.
1.88 + end if
1.89 + else
1.90 + LATS = lat(j)-0.5*(lat(j)-lat(j-1))
1.91 + LATN = lat(j)+0.5*(lat(j+1)-lat(j))
1.92 + end if
1.93 +
1.94 + do i=0,nlon-1
1.95 + if (i.eq.0 .or. i.eq.nlon-1) then
1.96 + if (i.eq.0) then
1.97 + LONL = 0.
1.98 + LONR = lon(i)+0.5*(lon(i+1)-lon(i))
1.99 + end if
1.100 + if (i.eq.nlon-1) then
1.101 + LONL = lon(i)-0.5*(lon(i)-lon(i-1))
1.102 + LONR = 360.
1.103 + end if
1.104 + else
1.105 + LONL = lon(i)-0.5*(lon(i)-lon(i-1))
1.106 + LONR = lon(i)+0.5*(lon(i+1)-lon(i))
1.107 + end if
1.108 +
1.109 +;print (LATS)
1.110 +;print (LATN)
1.111 +;print (LONL)
1.112 +;print (LONR)
1.113 + bo(m,j,i) = avg(bi(m,{LATS:LATN},{LONL:LONR}))
1.114 + end do
1.115 + end do
1.116 + end do
1.117 +
1.118 + bo!0 = "time"
1.119 + bo!1 = "lat"
1.120 + bo!2 = "lon"
1.121 + bo&time= bi&time
1.122 + bo&lat = lat
1.123 + bo&lon = lon
1.124 + bo@units = bi@units
1.125 + bo@long_name = bi@long_name
1.126 + bo@_FillValue = bi@_FillValue
1.127 +; bo@units = "none"
1.128 +; bo@long_name = "Leaf Area Index"
1.129 +
1.130 + c->LAI = bo
1.131 +end
1.132 \ No newline at end of file