npp/53.zonavg_3line.ncl
changeset 0 0c6405ab2ff4
equal deleted inserted replaced
-1:000000000000 0:ef07843fb511
       
     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 ;************************************************
       
    11 ; read in observed data
       
    12 ;************************************************
       
    13  
       
    14  a     = addfile ("Npp_T42_mean.nc","r")
       
    15  b     = addfile ("i01.03cn_1545-1569_ANN_climo.nc","r")
       
    16  c     = addfile ("i01.04casa_1605-1629_ANN_climo.nc","r")
       
    17 
       
    18  ta    = a->NPP
       
    19  tb    = b->NPP
       
    20  tc    = c->NPP
       
    21    
       
    22  lon   = a->lon 
       
    23  lat   = a->lat
       
    24  nlon  = dimsizes(lon)
       
    25  nlat  = dimsizes(lat)
       
    26 
       
    27  s  = new ((/3,nlat/), float)
       
    28 
       
    29 ;************************************************
       
    30 ; set value to 0. at missing point of observed
       
    31 ;************************************************
       
    32  delta = 0.00000000001
       
    33  x0    = tb(0,:,:)
       
    34  ta = where(ismissing(ta).and.(ismissing(x0).or.(x0.lt.delta)),0.,ta)
       
    35 
       
    36  s(0,:) = zonalAve(ta(:,:))
       
    37  delete (ta)
       
    38 
       
    39 ;************************************************
       
    40 ; for model data
       
    41 ;************************************************
       
    42  scale_factor = 86400.*365.
       
    43 
       
    44  s(1,:) = zonalAve(tb(0,:,:)) * scale_factor
       
    45  delete (tb)     
       
    46 
       
    47  s(2,:) = zonalAve(tc(0,:,:)) * scale_factor     
       
    48  delete (tc)     
       
    49 
       
    50  s@long_name = "NPP (gC/m2/year)"
       
    51  
       
    52 ;***************************************************** 
       
    53 ; create plot
       
    54 ;***************************************************** 
       
    55   wks = gsn_open_wks("png","xy")              ; create plot
       
    56   i   = NhlNewColor(wks,1.0,0.71,0.76)       ; add color to colormap
       
    57   j   = NhlNewColor(wks,0.64,0.71,0.8)       ; ditto
       
    58   
       
    59   res                    = True              ; plot mods desired
       
    60   res@gsnDraw            = False             ; don't draw yet
       
    61   res@gsnFrame           = False             ; don't advance frame yet
       
    62 
       
    63 ; res@vpHeightF 	 = 0.4               ; change aspect ratio of plot
       
    64 ; res@vpWidthF 	         = 0.7
       
    65 
       
    66 ; res@trXMinF	         = 1890              ; set x-axis minimum
       
    67 
       
    68   res@xyMonoLineColor    = "False"           ; want colored lines
       
    69   res@xyLineColors       = (/"Red","Blue","Black"/) ; colors chosen
       
    70 ; res@xyLineThicknesses	 = (/3.,3.,4./)      ; line thicknesses
       
    71   res@xyLineThicknesses	 = (/2.,2.,2./)      ; line thicknesses
       
    72   res@xyDashPatterns	 = (/0.,0.,0./)      ; make all lines solid
       
    73 
       
    74   res@tiMainString        = "Zoanl Average"
       
    75   res@tiYAxisString	 = "NPP (gC/m2/year)"      ; add a axis title    
       
    76   res@txFontHeightF	 = 0.0195            ; change title font heights
       
    77   
       
    78   plot  = gsn_csm_xy (wks,lat,s,res)       ; create plot
       
    79 
       
    80 ; plot  = fill_xy2(wks,plot(0),time,mnmx(2,:),mnmx(3,:),(/0.64,0.71,0.8/),\
       
    81 ; (/0.64,0.71,0.8/))
       
    82 ; plot  = fill_xy2(wks,plot(0),time,mnmx(0,:),mnmx(1,:),(/1.0,0.71,0.76/),\
       
    83 ; (/1.0,0.71,0.76/))
       
    84 ;*****************************************************   
       
    85 ; Manually create legend
       
    86 ;***************************************************** 
       
    87   res_text                    = True                  ; text mods desired
       
    88   res_text@txFontHeightF      = 0.015                 ; change text size
       
    89   res_text@txJust             = "CenterLeft"          ; text justification
       
    90 
       
    91   res_lines                   = True                  ; polyline mods desired
       
    92   res_lines@gsLineDashPattern = 0.                    ; solid line
       
    93   res_lines@gsLineThicknessF  = 5.                    ; line thicker
       
    94   res_lines@gsLineColor       = "red"                 ; line color
       
    95   xx = (/-85.,-75./)
       
    96   yy = (/1400.,1400./)
       
    97   gsn_polyline(wks,plot,xx,yy,res_lines)              ; add polyline
       
    98   gsn_text(wks,plot,"Observed",-70.,1400.,res_text); add text
       
    99   
       
   100   yy = (/1500.,1500./)
       
   101   res_lines@gsLineColor       = "blue"                ; change to blue
       
   102   gsn_polyline(wks,plot,xx,yy,res_lines)              ; add polyline
       
   103   gsn_text(wks,plot,"Model i01.03cn",-70.,1500.,res_text)     ; add text
       
   104   
       
   105   yy = (/1600.,1600./)
       
   106   res_lines@gsLineColor       = "black"               ; change to black
       
   107   gsn_polyline(wks,plot,xx,yy,res_lines)              ; add poly line
       
   108   gsn_text(wks,plot,"Model i01.04casa",-70.,1600.,res_text) ; add text
       
   109 ;*****************************************************   
       
   110 ; Manually create titles
       
   111 ;*****************************************************   
       
   112 ; res_text@txJust        = "CenterCenter"             ; change justification
       
   113 ; res_text@txFontHeightF = 0.03                       ; change font size
       
   114 ; gsn_text_ndc(wks,"Parallel Climate Model Ensembles",0.55,0.90,res_text)
       
   115 
       
   116 ; res_text@txFontHeightF = 0.02                       ; change font size
       
   117 ; gsn_text_ndc(wks,"Global Temperature Anomalies",0.55,0.86,res_text)
       
   118 
       
   119 ; res_text@txFontHeightF = 0.015                      ; change font size
       
   120 ; gsn_text_ndc(wks,"from 1890-1919 average",0.55,0.83,res_text)
       
   121   
       
   122   draw(plot)
       
   123   frame(wks)                                            ; advance frame
       
   124   
       
   125 end