lai/03.lai_ob_to_T31.ncl
changeset 0 0c6405ab2ff4
equal deleted inserted replaced
-1:000000000000 0:77a0005c01b8
       
     1 ; ***********************************************
       
     2 ; interpolate into model grids (T31)
       
     3 ; ***********************************************
       
     4 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
       
     5 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
       
     6 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
       
     7 ;************************************************
       
     8 begin
       
     9 
       
    10   year = 2000
       
    11 
       
    12 ;************************************************
       
    13 ; output data
       
    14 ;************************************************
       
    15   diro  = "/fis/cgd/cseg/people/jeff/clamp_data/lai/ob/"
       
    16   filo  = "LAI_" + year + "_monthly_T31.nc"
       
    17 
       
    18   c = addfile(diro+filo,"c")
       
    19   filedimdef(c,"time",-1,True)
       
    20 
       
    21 ;************************************************
       
    22 ; read in observed data
       
    23 ;************************************************
       
    24   diri  = "/fis/cgd/cseg/people/jeff/clamp_data/lai/ob/"
       
    25   fili  = "LAI_" + year + "_monthly.nc"
       
    26 
       
    27   g     = addfile (diri+fili,"r")
       
    28   bi    = g->LAI   
       
    29   xi    = g->lon 
       
    30   yi    = g->lat
       
    31 
       
    32 ;************************************************
       
    33 ; change into 0-360E, 90S-90N
       
    34 ; Observed NPP*scale_factor
       
    35 ;************************************************
       
    36  
       
    37  yi    = (/ yi(::-1) /)
       
    38  bi    = (/ bi(:,::-1,:) /)
       
    39  printVarSummary(bi)
       
    40 
       
    41  b2    = bi
       
    42  x2    = xi   
       
    43  
       
    44  nx = dimsizes(xi)
       
    45  do i= 0,nx-1
       
    46     if (i .lt. 3600) then
       
    47        p = i + 3600
       
    48        xi(p) = x2(i) + 360.      
       
    49     else
       
    50        p = i - 3600
       
    51        xi(p) = x2(i)
       
    52     end if
       
    53     bi(:,:,p)= b2(:,:,i) 
       
    54  end do
       
    55 
       
    56  bi&lat =  yi
       
    57  bi&lon =  xi
       
    58 
       
    59 ;print (xi)
       
    60 ;print (yi)
       
    61 
       
    62 ;************************************************
       
    63 ; read in model data
       
    64 ;************************************************
       
    65  diri2  = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
       
    66  fili2  = "b30.061n_1995-2004_MONS_climo_lnd.nc"
       
    67  f      = addfile (diri2+fili2,"r")
       
    68 
       
    69  lon    = f->lon     
       
    70  lat    = f->lat
       
    71  nlon   = dimsizes(lon)
       
    72  nlat   = dimsizes(lat)      
       
    73 
       
    74 ; print (xi)
       
    75 ; print (yi)
       
    76 ; print (xo)
       
    77 ; print (yo)
       
    78 
       
    79  bo = new((/12,nlat,nlon/),float)
       
    80 
       
    81  do m = 0,11
       
    82     do j=0,nlat-1
       
    83        if (j.eq.0 .or. j.eq.nlat-1) then
       
    84           if (j.eq.0) then
       
    85              LATS = -90.          
       
    86              LATN = lat(j)+0.5*(lat(j+1)-lat(j))
       
    87           end if
       
    88           if (j.eq.nlat-1) then
       
    89              LATS = lat(j)-0.5*(lat(j)-lat(j-1))
       
    90              LATN = 90.                  
       
    91           end if
       
    92        else
       
    93           LATS = lat(j)-0.5*(lat(j)-lat(j-1))
       
    94           LATN = lat(j)+0.5*(lat(j+1)-lat(j))
       
    95        end if
       
    96  
       
    97       do i=0,nlon-1
       
    98        if (i.eq.0 .or. i.eq.nlon-1) then
       
    99           if (i.eq.0) then
       
   100              LONL = 0.          
       
   101              LONR = lon(i)+0.5*(lon(i+1)-lon(i))
       
   102           end if
       
   103           if (i.eq.nlon-1) then
       
   104              LONL = lon(i)-0.5*(lon(i)-lon(i-1))
       
   105              LONR = 360.                 
       
   106           end if
       
   107        else
       
   108           LONL = lon(i)-0.5*(lon(i)-lon(i-1))
       
   109           LONR = lon(i)+0.5*(lon(i+1)-lon(i))
       
   110        end if
       
   111 
       
   112 ;print (LATS)
       
   113 ;print (LATN)
       
   114 ;print (LONL)
       
   115 ;print (LONR)
       
   116 
       
   117          bo(m,j,i) = avg(bi(m,{LATS:LATN},{LONL:LONR}))  
       
   118       end do 
       
   119     end do
       
   120   end do
       
   121 
       
   122   bo!0   = "time"
       
   123   bo!1   = "lat"
       
   124   bo!2   = "lon"
       
   125   bo&time= bi&time
       
   126   bo&lat = lat
       
   127   bo&lon = lon
       
   128   bo@units      = bi@units
       
   129   bo@long_name  = bi@long_name
       
   130   bo@_FillValue = bi@_FillValue
       
   131 
       
   132   c->LAI  = bo
       
   133 end