npp/52.zonavg_to_T31.ncl
changeset 1 4be95183fbcd
equal deleted inserted replaced
-1:000000000000 0:6b86d090bd36
       
     1 ; ***********************************************
       
     2 ; zonal average plot
       
     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 ; read in observed data
       
    11 ;************************************************
       
    12  g     = addfile ("Npp_0.05deg_mean.nc","r")
       
    13  bi    = g->NPP   
       
    14  xi    = g->lon 
       
    15  yi    = g->lat
       
    16  yi    = (/ yi(::-1) /)
       
    17  bi    =  (/ bi(::-1,:) /)
       
    18  b2    = bi
       
    19  x2    = xi   
       
    20  
       
    21  nx = dimsizes(xi)
       
    22  do i= 0,nx-1
       
    23     if (i .lt. 3600) then
       
    24        p = i + 3600
       
    25        xi(p) = x2(i) + 360.      
       
    26     else
       
    27        p = i - 3600
       
    28        xi(p) = x2(i)
       
    29     end if
       
    30     bi(:,p)= b2(:,i) 
       
    31  end do
       
    32 ;print (xi)
       
    33 ;print (yi)
       
    34 ;exit
       
    35 
       
    36  f     = addfile ("i01.03cn_1545-1569_ANN_climo.nc","r")
       
    37 ;f     = addfile ("i01.04casa_1605-1629_ANN_climo.nc","r")
       
    38 
       
    39  xo    = f->lon     
       
    40  yo    = f->lat      
       
    41 
       
    42  print (xi)
       
    43  print (yi)
       
    44  print (xo)
       
    45  print (yo)
       
    46 
       
    47  bo = linint2_Wrap(xi,yi,bi,True,xo,yo,0) 
       
    48  
       
    49  v = zonalAve(bo)
       
    50 
       
    51  v@long_name = "NPP (gC/m2/year)"
       
    52  
       
    53 ;************************************************
       
    54 ; plotting parameters
       
    55 ;************************************************
       
    56  wks   = gsn_open_wks ("png","xy")                ; open workstation
       
    57 
       
    58  res                  = True                     ; plot mods desired
       
    59  res@tiMainString     = "Observed T31"          ; add title
       
    60 
       
    61  plot  = gsn_csm_xy (wks,v&lat,v,res) 
       
    62 end