1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/npp/53.zonavg_3line.ncl Mon Jan 26 22:08:20 2009 -0500
1.3 @@ -0,0 +1,125 @@
1.4 +; ***********************************************
1.5 +; zonal average plot
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 +; read in observed data
1.15 +;************************************************
1.16 +
1.17 + a = addfile ("Npp_T42_mean.nc","r")
1.18 + b = addfile ("i01.03cn_1545-1569_ANN_climo.nc","r")
1.19 + c = addfile ("i01.04casa_1605-1629_ANN_climo.nc","r")
1.20 +
1.21 + ta = a->NPP
1.22 + tb = b->NPP
1.23 + tc = c->NPP
1.24 +
1.25 + lon = a->lon
1.26 + lat = a->lat
1.27 + nlon = dimsizes(lon)
1.28 + nlat = dimsizes(lat)
1.29 +
1.30 + s = new ((/3,nlat/), float)
1.31 +
1.32 +;************************************************
1.33 +; set value to 0. at missing point of observed
1.34 +;************************************************
1.35 + delta = 0.00000000001
1.36 + x0 = tb(0,:,:)
1.37 + ta = where(ismissing(ta).and.(ismissing(x0).or.(x0.lt.delta)),0.,ta)
1.38 +
1.39 + s(0,:) = zonalAve(ta(:,:))
1.40 + delete (ta)
1.41 +
1.42 +;************************************************
1.43 +; for model data
1.44 +;************************************************
1.45 + scale_factor = 86400.*365.
1.46 +
1.47 + s(1,:) = zonalAve(tb(0,:,:)) * scale_factor
1.48 + delete (tb)
1.49 +
1.50 + s(2,:) = zonalAve(tc(0,:,:)) * scale_factor
1.51 + delete (tc)
1.52 +
1.53 + s@long_name = "NPP (gC/m2/year)"
1.54 +
1.55 +;*****************************************************
1.56 +; create plot
1.57 +;*****************************************************
1.58 + wks = gsn_open_wks("png","xy") ; create plot
1.59 + i = NhlNewColor(wks,1.0,0.71,0.76) ; add color to colormap
1.60 + j = NhlNewColor(wks,0.64,0.71,0.8) ; ditto
1.61 +
1.62 + res = True ; plot mods desired
1.63 + res@gsnDraw = False ; don't draw yet
1.64 + res@gsnFrame = False ; don't advance frame yet
1.65 +
1.66 +; res@vpHeightF = 0.4 ; change aspect ratio of plot
1.67 +; res@vpWidthF = 0.7
1.68 +
1.69 +; res@trXMinF = 1890 ; set x-axis minimum
1.70 +
1.71 + res@xyMonoLineColor = "False" ; want colored lines
1.72 + res@xyLineColors = (/"Red","Blue","Black"/) ; colors chosen
1.73 +; res@xyLineThicknesses = (/3.,3.,4./) ; line thicknesses
1.74 + res@xyLineThicknesses = (/2.,2.,2./) ; line thicknesses
1.75 + res@xyDashPatterns = (/0.,0.,0./) ; make all lines solid
1.76 +
1.77 + res@tiMainString = "Zoanl Average"
1.78 + res@tiYAxisString = "NPP (gC/m2/year)" ; add a axis title
1.79 + res@txFontHeightF = 0.0195 ; change title font heights
1.80 +
1.81 + plot = gsn_csm_xy (wks,lat,s,res) ; create plot
1.82 +
1.83 +; plot = fill_xy2(wks,plot(0),time,mnmx(2,:),mnmx(3,:),(/0.64,0.71,0.8/),\
1.84 +; (/0.64,0.71,0.8/))
1.85 +; plot = fill_xy2(wks,plot(0),time,mnmx(0,:),mnmx(1,:),(/1.0,0.71,0.76/),\
1.86 +; (/1.0,0.71,0.76/))
1.87 +;*****************************************************
1.88 +; Manually create legend
1.89 +;*****************************************************
1.90 + res_text = True ; text mods desired
1.91 + res_text@txFontHeightF = 0.015 ; change text size
1.92 + res_text@txJust = "CenterLeft" ; text justification
1.93 +
1.94 + res_lines = True ; polyline mods desired
1.95 + res_lines@gsLineDashPattern = 0. ; solid line
1.96 + res_lines@gsLineThicknessF = 5. ; line thicker
1.97 + res_lines@gsLineColor = "red" ; line color
1.98 + xx = (/-85.,-75./)
1.99 + yy = (/1400.,1400./)
1.100 + gsn_polyline(wks,plot,xx,yy,res_lines) ; add polyline
1.101 + gsn_text(wks,plot,"Observed",-70.,1400.,res_text); add text
1.102 +
1.103 + yy = (/1500.,1500./)
1.104 + res_lines@gsLineColor = "blue" ; change to blue
1.105 + gsn_polyline(wks,plot,xx,yy,res_lines) ; add polyline
1.106 + gsn_text(wks,plot,"Model i01.03cn",-70.,1500.,res_text) ; add text
1.107 +
1.108 + yy = (/1600.,1600./)
1.109 + res_lines@gsLineColor = "black" ; change to black
1.110 + gsn_polyline(wks,plot,xx,yy,res_lines) ; add poly line
1.111 + gsn_text(wks,plot,"Model i01.04casa",-70.,1600.,res_text) ; add text
1.112 +;*****************************************************
1.113 +; Manually create titles
1.114 +;*****************************************************
1.115 +; res_text@txJust = "CenterCenter" ; change justification
1.116 +; res_text@txFontHeightF = 0.03 ; change font size
1.117 +; gsn_text_ndc(wks,"Parallel Climate Model Ensembles",0.55,0.90,res_text)
1.118 +
1.119 +; res_text@txFontHeightF = 0.02 ; change font size
1.120 +; gsn_text_ndc(wks,"Global Temperature Anomalies",0.55,0.86,res_text)
1.121 +
1.122 +; res_text@txFontHeightF = 0.015 ; change font size
1.123 +; gsn_text_ndc(wks,"from 1890-1919 average",0.55,0.83,res_text)
1.124 +
1.125 + draw(plot)
1.126 + frame(wks) ; advance frame
1.127 +
1.128 +end