lai/31.contour_ob_lai_max.ncl.1
changeset 1 4be95183fbcd
equal deleted inserted replaced
-1:000000000000 0:188bb85b4497
       
     1 ;*************************************************
       
     2 ; ce_1.ncl
       
     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   diri  = "/fis/cgd/cseg/people/jeff/clamp_data/"
       
    13 ; fili  = "LAI_2000-2005_ensemble_T42.nc"
       
    14   fili  = "2000-2005_ensemble_T42.nc"
       
    15   f     = addfile(diri+fili,"r")
       
    16 
       
    17   z     = f->LAI
       
    18   y     = z(0,:,:)
       
    19   y@long_name = "Leaf Area Index Max"
       
    20   p     = z(0,:,:)
       
    21   p@long_name = "Leaf Area Index Max Month"
       
    22   s     = z(:,0,0)
       
    23 ; printVarSummary(s)
       
    24  
       
    25   dsizes_z = dimsizes(z)
       
    26   ntime    = dsizes_z(0)
       
    27   nlat     = dsizes_z(1)
       
    28   nlon     = dsizes_z(2)
       
    29   
       
    30   do j = 0,nlat-1
       
    31   do i = 0,nlon-1
       
    32      s      = z(:,j,i) 
       
    33      y(j,i) = max(s)
       
    34      p(j,i) = maxind(s)+ 1
       
    35 
       
    36 ;    iMax   = ind(s.eq.y(j,i)) + 1      ; may be multiple max
       
    37 ;    if (.not.ismissing(iMax(0))) then 
       
    38 ;        p(j,i) =  iMax(0)          ; no meta data
       
    39 ;    end if
       
    40 ;    delete(iMax)
       
    41 
       
    42   end do
       
    43   end do
       
    44 
       
    45 ; printVarSummary(p)
       
    46   print (min(y)+"/"+max(y))
       
    47   print (min(p)+"/"+max(p))
       
    48 exit
       
    49 ;************************************************
       
    50 ; read in model data
       
    51 ;************************************************
       
    52   fili2 = "i01.03cn_1545-1569_ANN_climo.nc"
       
    53   g     = addfile(diri+fili2,"r")
       
    54   x     = g->TLAI
       
    55 
       
    56   delta = 0.00000000001
       
    57   x0    = x(0,:,:)
       
    58   y = where(ismissing(y).and.(ismissing(x0).or.(x0.lt.delta)),0.,y)
       
    59 ;************************************************
       
    60 ; create default plot
       
    61 ;************************************************
       
    62   
       
    63 ; setvalues NhlGetWorkspaceObjectId()
       
    64 ;   "wsMaximumSize" : 199999999
       
    65 ; end setvalues
       
    66 
       
    67   wks = gsn_open_wks("ps","xy")         ; open a ps file
       
    68 ; gsn_define_colormap(wks,"wgne15")          ; choose colormap
       
    69   gsn_define_colormap(wks,"gui_default")     ; choose colormap
       
    70 
       
    71   res                     = True             ; Use plot options
       
    72   res@cnFillOn            = True             ; Turn on color fill
       
    73   res@gsnSpreadColors      = True            ; use full colormap
       
    74 ; res@cnFillMode          = "RasterFill"     ; Turn on raster color
       
    75 ; res@lbLabelAutoStride   = True
       
    76   res@cnLinesOn           = False            ; Turn off contourn lines
       
    77   res@mpFillOn            = False            ; Turn off map fill
       
    78   res@tiMainString        = "MODIS MOD 15A2 2000-2005"
       
    79 
       
    80   res@gsnSpreadColors     = True             ; use full colormap
       
    81   res@cnLevelSelectionMode = "ManualLevels"  ; Manual contour invtervals
       
    82   res@cnMinLevelValF       = 0.              ; Min level
       
    83   res@cnMaxLevelValF       = 12.             ; Max level
       
    84   res@cnLevelSpacingF      = 1.              ; interval
       
    85 
       
    86 ; pres                            = True     ; panel plot mods desired
       
    87 ; pres@gsnMaximize                = True     ; fill the page
       
    88 
       
    89   plot = gsn_csm_contour_map_ce(wks,y,res)   ; for observed
       
    90 
       
    91   plot = gsn_csm_contour_map_ce(wks,p,res)   ; for observed
       
    92 
       
    93 ; gsn_panel(wks,plot,(/1,1/),pres)           ; create panel plot
       
    94   system("convert xy.ps xy.png")
       
    95 end