Joe,
We had great results recently with almost the identical plot that
you requested.
We took advantage of the high level gsn routine
gsn_csm_vector_scalar_map_ce, which draws shaded contours and curly
vectors together in one function. This simplifies the code. The
only difference from yours is that the top overlay is the line
contours, rather than the curly vectors.
We thought this was a reasonable compromise, since those are both
line plots. However, you could instead draw the vectors as a
separate overlay and get the exact layering that you want.
The very simplified code sequence is:
res_at_gsnDraw = False
res_at_gsnFrame = False
res2_at_gsnDraw = False
res2_at_gsnFrame = False
plot1 = gsn_csm_vector_scalar_map_ce (wks, u, v, grid, res)
plot2 = gsn_csm_contour (wks, sfunc, res2)
overlay (plot1, plot2)
draw (plot1)
frame (wks)
Attached is the full working script. This one also adds a fourth
"overlay" on top, in the form of a circle symbol using wmlabs.
Apologies for the clutter, this version has been through a lot.
--Dave A.
CU/CIRES Climate Diagnostics Center (CDC)
NOAA/ESRL/PSD, Climate Analysis Branch (CAB)
Joe Metzger wrote:
> I've searched the website and archives for this capability and may have
> missed it, but I'm wondering if NCL can plot *three* overlays? I want a
> color-filled plot on the bottom, then line contours on top of that, and
> finally curly vectors on the very top. Is this possible? Thanks for any
> advice or sample scripts.
>
> Joe
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk_at_ucar.edu
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
;----------------------------------------------------------------------------
;
; plotcomp.ncl -- Make offset plot series of composite fields.
;
; 2006-jun-22 bl, convert to south america. Adapted from sstcor.ncl.
; 2006-jun-29 dra, add both binary and Netcdf readers for mulitple offsets.
; 2007-jan-23 bl, add vector capability (program copied from plotcomp.ncl)
; 2007-jun-05 dra, add stream function from composite U and V winds.
; 2007-jun-25 dra, add contour line overlay for stream function, on
; same plot as shaded contour data field and wind vectors.
; Add command output selectors "ps=1" and "pdf=1", default = X11.
; Fix overlay dots.
;
;----------------------------------------------------------------------------
print ("Load ncl scripts")
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
; load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
print ("Load complete")
begin
; Parameters.
olrdata = False ; select data field, set only one to True
raindata = True
if(olrdata)then
infile = "/mac_home/bl/transient/composite/composite.olr.nc"
varname = "compolr"
end if
if(raindata)then
infile = "/mac_home/bl/transient/composite/composite.rain.lanczos.nc"
varname = "comprain"
end if
idir = "/mac_home/bl/transient/composite/"
uwindfile = idir + "composite.uwnd.lanczos.nc"
varnameu = "compuwnd"
vwindfile = idir + "composite.vwnd.lanczos.nc"
varnamev = "compvwnd"
sf_scale = "1e6" ; scale factor for stream function
;; sf_units = fs_scale + " x unity" ; units label for stream function
;; output = "x" ; x=window, ps=postscript, pdf=PDF
;; ; not used any more; use "ps=1" or "pdf=1"
; ; in command line; default = X11 window
outfile = "kindexrain850-2.5N95W-3day10w"
; outfile = "kindexstream200-largest"
; name of output file; leave off .nc extension
title = "Composite from index"
; Plot borders. Independent from data grid borders.
lonw = -120 ;
lone = -5
lats = -50
latn = 30
; lat_values = (/ -40, -30, -20, -10, 0, 10 /) ; lat labels with "EQ"
; lat_labels = (/ "40S", "30S", "20S", "10S", " EQ", "10N" /)
; contour_levels = (/ -7,-6,-5,-4,-3,-2, 2, 3, 4, 5, 6, 7 /)
; colormap = "WhViBlGrYeOrRe"
colormap = "testcmap"
black = 1 ; colors are specific to each color map
blue = 35
green = 80
red = 200
;----------------------------------------------------------------------------
; Get command line parameters.
if (isdefined ("ps")) then ; get output selector
output = "ps" ; default = X11 window
else
if (isdefined ("pdf")) then
output = "pdf"
else
output = "x11"
end if
end if
;----------------------------------------------------------------------------
; Derived parameters.
chars = stringtochar (infile) ; determine file type
len = dimsizes (chars) - 1 ; shortcut: assume 3 char ".xx" exten.
e1 = max ( (/ 0, len-3 /) )
e2 = max ( (/ e1, len-1 /) )
suffix = chartostring (chars(e1:e2))
netcdf = (suffix .eq. ".nc")
binary = .not. netcdf
if (binary) then
nt = offset2 - offset_base + 1
end if
;; Replaced by command line output selector.
;; if (output .ne. "ps" .and. output .ne. "pdf") then
;; output = "x11"
;; end if
sf_scale_F = stringtofloat (sf_scale)
;----------------------------------------------------------------------------
; Read composite main input file.
; if (netcdf) then
print ("Read main input file.")
ifile = addfile (infile, "r")
in = ifile->$varname$
; end if
if (any (isnan_ieee (in))) then ; check input grids for NaN's
print ("*** Abort, input contains NaN's.")
exit
end if
if (.not. isdefined ("_FillValue")) then
in@_FillValue = in_at_missing_value
end if
; printVarSummary (in)
print ("Min data = " + sprintf ("%0.3f", min (in)))
print ("Max data = " + sprintf ("%0.3f", max (in)))
print ("")
;----------------------------------------------------------------------------
print ("Read composite wind files.")
wfu = addfile (uwindfile, "r")
inu = wfu->$varnameu$
wfv = addfile (vwindfile, "r")
inv = wfv->$varnamev$
;----------------------------------------------------------------------
; Common setup for all plots.
;----------------------------------------------------------------------
; print ("Start output file or window.")
wks = gsn_open_wks (output, outfile) ; open output "device"
gsn_define_colormap (wks, colormap)
res = True ; plot mods desired
res_at_gsnDraw = False ; hold draw and auto frame for
res_at_gsnFrame = False ; overlays (contours and dots)
res_at_cnFillOn = True ; turn on color fill
res_at_cnLinesOn = True ; turn on contour lines
res_at_gsnSpreadColors = True ; choose colors across entire range
; res_at_gsnSpreadColorStart = -1 ; reverse colors
; res_at_gsnSpreadColorEnd = 2
; res_at_gsnSpreadColorEnd = 10
; res_at_gsnSpreadColorEnd = 30
res_at_lbOrientation = "horizontal" ; horizontal label bars
; res_at_lbOrientation = "vertical" ; vertical label bars
res_at_pmLabelBarOrthogonalPosF = 0.3 ; label bar farther from plot
res_at_tiMainFontHeightF = 0.018 ; adjust main title font size
if(olrdata)then ; set contours for OLR data
res_at_cnLevelSelectionMode = "ExplicitLevels"
contour_levels = (/ -30,-25.,-20.,-15.,-10.,-5.,5.,10.,15.,20.,25.,30. /)
res_at_cnLevels = contour_levels
; res_at_cnLevelSelectionMode = "ManualLevels" ; manual levels
; res_at_cnMinLevelValF = -10 ; min level
; res_at_cnMaxLevelValF = 10 ; max level
; res_at_cnLevelSpacingF = 2 ; interval
; res_at_cnMinLevelValF = -7 ; min level ; good for rain
; res_at_cnMaxLevelValF = 7 ; max level
; res_at_cnLevelSpacingF = 2 ; interval
; res_at_cnMinLevelValF = -30 ; min level
; res_at_cnMaxLevelValF = 30 ; max level
; res_at_cnLevelSpacingF = 5 ; interval
end if
if(raindata)then ; set contours for rain data
res_at_cnLevelSelectionMode = "ExplicitLevels"
contour_levels = (/ -10.,-8.,-6.,-4.,-2.,2.,4.,6.,8.,10. /)
res_at_cnLevels = contour_levels
; res_at_cnLevelSelectionMode = "ManualLevels" ; manual levels
; res_at_cnMinLevelValF = -8 ; min level
; res_at_cnMaxLevelValF = 8 ; max level
; res_at_cnLevelSpacingF = 1 ; interval
end if
; Display color bar without end labels.
; res_at_cnLevels = contour_levels
; Display color bar with forced end labels.
; Warning: Data beyond end values will be plotted to appear as if inside
; the end values. No warning is printed unless code is added.
; This version: Forced bottom end, open top end.
; nclevs = dimsizes (contour_levels)
;; res_at_cnLevels = contour_levels(1:nclevs-2) ; hide end vals from NCL
; res_at_cnLevels = contour_levels(1:) ; hide only bottom value
; res_at_cnExplicitLabelBarLabelsOn = True
; res_at_cnLabelBarEndLabelsOn = True
; label_strings = new (nclevs+1, string, "")
; label_strings(0:nclevs-1) = sprintf ("%1.0f", contour_levels)
; res_at_lbLabelStrings = label_strings
; end of forced bottom and top labels
res_at_cnMissingValPerimOn = True
res_at_cnMissingValPerimDashPattern = 1
res_at_mpFillOn = False ; turn off gray continents
res_at_gsnMaximize = True ; makes figure fit in frame properly
res_at_gsnPaperOrientation = "portrait" ; always portrait
; Plot boundaries.
res_at_mpMinLonF = lonw
res_at_mpMaxLonF = lone
res_at_mpMinLatF = lats
res_at_mpMaxLatF = latn
res_at_mpCenterLonF = 0.5 * (lonw + lone) ; get accurate center longitude
res_at_gsnAddCyclic = False ; must be false if not wrapping
; res_at_gsnAddCyclic = True ; must be false if not wrapping
; data at Greenwich meridian
; res_at_mpOutlineBoundarySets = "AllBoundaries"
; res_at_mpOutlineBoundarySets = "USStates"
; res_at_mpOutlineBoundarySets = "geophysicalandusstates"; turn on states
res_at_mpOutlineBoundarySets = "national" ; countries, continents, islands
res_at_mpDataBaseVersion = "mediumres" ; select database
res_at_mpDataSetName = "Earth..2"
;------------------------------------------------------------------------------
; Common setup for stream function overlay plot.
res2 = True ; make resource object for plot2
res2_at_gsnDraw = False ; also hold draw and auto frame
res2_at_gsnFrame = False ; for overlay plot
res2_at_cnLevelSelectionMode = "ManualLevels" ; manual contour levels
res2_at_lbLabelAutoStride = True
res2_at_cnMinLevelValF = -6 ; min level
res2_at_cnMaxLevelValF = 6 ; max level
res2_at_cnLevelSpacingF = 0.5 ; interval
res2_at_gsnLeftString = "" ; no auto labels
res2_at_gsnRightString = ""
res2_at_gsnCenterString = ""
res2_at_cnInfoLabelOrthogonalPosF = 0.15 ; move info label down
; (to be fixed later)
;----------------------------------------------------------------------
; Loop on offset, for each frame.
;----------------------------------------------------------------------
offset1= in_at_first_offset
offset2= in_at_last_offset
dims = dimsizes(in)
noffsets = dims(0)
do ti = 0, noffsets -1
offset = offset1+ti
print("offset in loop = " + offset)
grid = in(ti,:,:) ; get composite grid for this offfset
u = inu(ti,::2,::2) ; get U and V grids for this offset
v = inv(ti,::2,::2)
; v = 0. ;NOTE SETTING V COMPONENT to 0
res_at_gsnCenterString = sprinti ("offset %0i", offset)
res_at_gsnStringFontHeightF = 0.018
printVarSummary(in)
degrees_north = in_at_key_lat
degrees_west = 360 - in_at_key_lon
comp_days = in_at_comp_count
mblevel = in_at_level
; res_at_gsnLeftString = in_at_title3 + " (" + comp_days + " days)"
title3_plus = in_at_title3 + " (" + comp_days + " days)"
res_at_gsnLeftString = ""
if (in_at_key_lat .eq. 0) then
print_lat = "Equator"
else
if (in_at_key_lat .gt. 0) then
print_lat = degrees_north + " deg N"
else
print_lat = -degrees_north + " deg S"
end if
end if
res_at_tiMainString = in_at_title1 + "~C~ " \
+ in_at_title2 + " " + degrees_west + " deg W, " + print_lat \
+ "~C~ " + in_at_title4 + " " + in_at_title5 \
+ "~C~" + title3_plus
res_at_tiMainOffsetYF = 0.005
print (sprinti("offset %3i",offset) + sprintf (" min %5.2f", min(grid)) \
+ sprintf (" max %5.2f", max(grid)))
res_at_gsnScalarContour = True ; contours desired
; res_at_vcRefMagnitudeF = 8.0 ; define vector ref mag ;200 mb
; res_at_vcRefMagnitudeF = 2.0 ; define vector ref mag 1000 mb
res_at_vcRefMagnitudeF = 4.0 ; define vector ref mag 1000 mb
res_at_vcRefLengthF = 0.05 ; define length of vec ref
res_at_vcRefAnnoOrthogonalPosF = -1.0 ; move ref vector
res_at_vcRefAnnoArrowLineColor = "black" ; change ref vector color
res_at_vcRefAnnoArrowUseVecColor = False ; don't use vec color for ref
res_at_vcGlyphStyle = "CurlyVector" ; turn on curley vectors
res_at_vcLineArrowColor = "black" ; change vector color
res_at_vcLineArrowThicknessF = 2.0 ; change vector thickness
res_at_vcVectorDrawOrder = "PostDraw" ; draw vectors last
res_at_vcMinDistanceF = 0.005 ; Thin the vectors
plot1 = gsn_csm_vector_scalar_map_ce (wks, u, v, grid, res)
; plot1 = gsn_csm_vector_scalar_map_ce (wks, u, v, u, res)
; make contour plot with wind vector overlay.
; plot1 = gsn_csm_contour_map_ce (wks, grid, res) ; make next plot frame
;----------------------------------------------------------------------
; Add overlay contour lines for stream function.
sf_vp = uv2sfvpF (inu(ti,:,:), inv(ti,:,:)) ; compute stream function
; for current offset
sfunc = sf_vp(0,:,:) / sf_scale_F ; get stream function alone
ugrid = inu(ti,:,:) ; copy over lat and lon coords
copy_VarCoords (ugrid, sfunc)
;; grid_at_units = sf_units ; SF not labeled in this version
plot2 = gsn_csm_contour (wks, sfunc, res2) ; add stream func overlay
overlay (plot1, plot2) ; on top of main plot
draw (plot1)
;----------------------------------------------------------------------
; Draw dot at specified location.
; Multiple dots may be drawn with multiple wmlabs commands.
; Color and size may be changed between dots.
dot_lat = degrees_north
dot_lon = -degrees_west ; must be -180 to +180 ???
wmsetp ("EZF", 1)
wmsetp ("DBC", blue) ; set color of dot outer ring
wmsetp ("DTC", red) ; set color of dot center
wmsetp ("DTS", 0.012) ; set dot size
wmlabs (wks, dot_lat, dot_lon, "DT") ; draw dot at given lat, lon
print("Dot lat = " + dot_lat + ", lon = " + dot_lon)
;----------------------------------------------------------------------
frame (wks) ; send final plot to output
delete (grid) ; must delete working grid
delete (plot1) ; and working plot objects
delete (plot2) ; between loops
end do
print ("Done.")
end
_______________________________________________
ncl-talk mailing list
ncl-talk_at_ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Wed Jul 11 2007 - 13:12:27 MDT
This archive was generated by hypermail 2.2.0 : Tue Jul 17 2007 - 06:52:00 MDT