[ncl-talk] panel plots with overlays

mberdahl at envsci.rutgers.edu mberdahl at envsci.rutgers.edu
Mon Aug 3 13:19:47 MDT 2015


Actually, I take that back.  I don't get a blank plot (that was a stupid
error), but I do just get the same thing as before... Overlaid contours
and vectors, with no base map.

mira

> Hi Alex,
>
> Yes! The plot you attached is the map I'd like to have under my wind
> vectors and filled contours.  However, when I change the max and min lon
> in my map resources to what you did (0 and -90), I do not receive any map
> - it is a blank page with just a title.
>
> I wonder if it has to do with the lines where I extract the longitudes in
> my variables:
>
> u = a->uwnd(0:803,{500},{45:90},{270:357.5})
> v = b->vwnd(0:803,{500},{45:90},{270:357.5}) ; Get u, v, time (1),level
> (1000hpa),latitude(-90:90) and longitude(0:360) data.
> z = c->hgt(0:803,{500},{45:90},{270:357.5})  ; get geopotential heights...
>
> Do you think it is something to do with this?
>
> Mira
>
>
>
>> Hi Mira,
>>
>> I copied these lines from your script and ran them locally:
>>
>> wks = gsn_open_wks("ps","Panel_NAO_z_500") 		; open a ps file
>> gsn_define_colormap(wks,"temp1")
>>
>> plot = new(3,graphic)				; create a plot array
>>
>> ;---- set common resources for all plots
>> res 			= True
>> res at gsnDraw		= False 		; dont draw
>> res at gsnFrame		= False			; dont advance frame
>> res at cnInfoLabelOn	= False			; trn off cn info label
>> res at gsnAddCyclic = False			; has to do with wrapping the longitude at
>> 0/360
>> res at cnFillPalette = "matlab_jet"
>> ;************************************************
>> ; Choose a subregion
>> ;************************************************
>> res at mpMaxLatF = 90 				;maximum latitude
>> res at mpMinLatF = 45 				;minimum latitude
>> res at mpMaxLonF = 357.5 			;	;maximum longitude
>> res at mpMinLonF = 270 				;minimum longitude
>> ;res at mpFillBoundarySets = AllBoundaries
>> res at mpOutlineBoundarySets = "National"
>> res at mpOutlineOn = True
>> res at mpOutlineDrawOrder = "PostDraw"
>>
>> mpid = gsn_csm_map(wks,res)
>>
>>
>> I got an error and made a change in your lons.  I set them to 0 and -90
>> respectively and got this base map. Is this what you were looking for?
>>
>> -Alex
>>
>>
>>
>>
>>> On Jul 30, 2015, at 9:00 PM, mberdahl at envsci.rutgers.edu wrote:
>>>
>>> Hi Alex,
>>> Thanks again for the help.  Setting the workstation colormap worked
>>> once
>>> I
>>> realized I have an older version (6.0.0) of NCL so there is a limited
>>> set
>>> of colormaps available.
>>> I still, however, can't get a base map plotted.  I've tried a couple
>>> things, including what you suggested below, but this doesn't work.
>>> I've
>>> tried creating a map and overlaying it with the contour and vectors,
>>> but
>>> it doesn't work.  I receive the following error and warnings:
>>>
>>>
>>> ^Mfatal:NhlAddOverlay: tranform is already an annotation or overlay: 74
>>> ^Mwarning:ContourPlotSetValues: attempt to set overlay member plot view
>>> ignored
>>> ^Mwarning:ContourPlotSetValues: attempt to set overlay member plot view
>>> ignored
>>> ^Mwarning:ContourPlotSetValues: attempt to set overlay member plot view
>>> ignored
>>> ^Mwarning:NhlDraw: cannot draw Plot Member, ID 262, independently
>>> ^Mwarning:NhlDraw: cannot draw Plot Member, ID 298, independently
>>> ^Mwarning:NhlDraw: cannot draw Plot Member, ID 334, independently
>>>
>>>
>>>
>>> My code is copied below.
>>>
>>>
>>> Any thoughts?
>>> Thanks again,
>>> Mira
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>> 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/contributed.ncl"
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
>>> ;************************************************
>>> begin
>>> ;************************************************
>>> ; read in netCDF file s
>>> ;************************************************
>>> a = addfile("uwnd.mon.mean.alllevels.nc","r") ; u winds
>>> b = addfile("vwnd.mon.mean.alllevels.nc","r") ; v winds
>>>
>>> c = addfile("../Geopotential/hgt.mon.mean.nc","r") ; geopotential
>>> heights.
>>>
>>> ;************************************************
>>> ; read in zonal [u] and meridional [v] winds (July)
>>> ;************************************************
>>>
>>> u = a->uwnd(0:803,{500},{45:90},{270:357.5})
>>> v = b->vwnd(0:803,{500},{45:90},{270:357.5}) ; Get u, v, time (1),level
>>> (1000hpa),latitude(-90:90) and longitude(0:360) data.
>>> z = c->hgt(0:803,{500},{45:90},{270:357.5})  ; get geopotenial
>>> heights...
>>>
>>> printVarSummary(u)
>>> printVarSummary(v)
>>> printVarSummary(z)
>>>
>>> ; Calculate the seasonal averages.
>>> uDJF = month_to_season(u, "DJF")
>>> vDJF = month_to_season(v, "DJF")
>>> zDJF = month_to_season(z, "DJF")
>>>
>>> printVarSummary(uDJF)
>>> printVarSummary(vDJF)
>>> printVarSummary(zDJF)
>>>
>>> ; from the matlab script i wrote: findExtremeYrs, i pulled out the
>>> extreme
>>> years (> or < 1std) that i want to average and plot here.
>>>
>>> ; for ans =   4 (NAO)
>>> ; yearList_hi = 1973        1975        1983        1989        1995
>>> 2000        2007        2012
>>> ; yearList_lo = 1963        1964        1965        1969        1977
>>> 1979        1996        1997        2010        2011
>>>
>>>
>>> ; this data starts at 1948 (this is index 0), so 1953=5, 1963=15 etc.
>>>
>>> uDJF_NAO_hi = uDJF((/25,27,35,41,47,52,59,64/),:,:)
>>> uDJF_NAO_lo = uDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>
>>> vDJF_NAO_hi = vDJF((/25,27,35,41,47,52,59,64/),:,:)
>>> vDJF_NAO_lo = vDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>
>>> zDJF_NAO_hi = zDJF((/25,27,35,41,47,52,59,64/),:,:)
>>> zDJF_NAO_lo = zDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>
>>>
>>> uAvgTime_hi = dim_avg_n_Wrap(uDJF_NAO_hi,0)
>>> uAvgTime_lo = dim_avg_n_Wrap(uDJF_NAO_lo,0)
>>>
>>> printVarSummary(uAvgTime_hi)
>>> printVarSummary(uAvgTime_lo)
>>>
>>> vAvgTime_hi = dim_avg_n_Wrap(vDJF_NAO_hi,0)
>>> vAvgTime_lo = dim_avg_n_Wrap(vDJF_NAO_lo,0)
>>>
>>> printVarSummary(vAvgTime_hi)
>>> printVarSummary(vAvgTime_lo)
>>>
>>> zAvgTime_hi = dim_avg_n_Wrap(zDJF_NAO_hi,0)
>>> zAvgTime_lo = dim_avg_n_Wrap(zDJF_NAO_lo,0)
>>>
>>> printVarSummary(zAvgTime_hi)
>>> printVarSummary(zAvgTime_lo)
>>>
>>> ; dirty way to copy metadata over first.
>>> diff_u = uAvgTime_hi;
>>> diff_v = vAvgTime_hi;
>>> diff_z = zAvgTime_hi;
>>>
>>> diff_u = uAvgTime_hi - uAvgTime_lo
>>> diff_v = vAvgTime_hi - vAvgTime_lo
>>> diff_z = zAvgTime_hi - zAvgTime_lo
>>>
>>> printVarSummary(diff_u)
>>> printVarSummary(diff_v)
>>> printVarSummary(diff_z)
>>>
>>> ;************************************************
>>> ; create plot
>>> ;************************************************
>>> wks = gsn_open_wks("ps","Panel_NAO_z_500") 		; open a ps file
>>> gsn_define_colormap(wks,"temp1")
>>>
>>> plot = new(3,graphic)				; create a plot array
>>>
>>> ;---- set common resources for all plots
>>> res 			= True
>>> res at gsnDraw		= False 		; dont draw
>>> res at gsnFrame		= False			; dont advance frame
>>> res at cnInfoLabelOn	= False			; trn off cn info label
>>> res at gsnAddCyclic = False			; has to do with wrapping the longitude at
>>> 0/360
>>> res at cnFillPalette = "matlab_jet"
>>> ;************************************************
>>> ; Choose a subregion
>>> ;************************************************
>>> res at mpMaxLatF = 90 				;maximum latitude
>>> res at mpMinLatF = 45 				;minimum latitude
>>> res at mpMaxLonF = 357.5 			;	;maximum longitude
>>> res at mpMinLonF = 270 				;minimum longitude
>>> ;res at mpFillBoundarySets = AllBoundaries
>>> res at mpOutlineBoundarySets = "National"
>>> res at mpOutlineOn = True
>>> res at mpOutlineDrawOrder = "PostDraw"
>>>
>>> mpid = gsn_csm_map(wks,res)
>>>
>>> ;***********************************************
>>> ; ----wind  vector plot
>>> ;***********************************************
>>> vcres = res
>>> vcres at vcRefAnnoOrthogonalPosF = -1.0 		; move ref vector up
>>> vcres at vcRefMagnitudeF = 10.0 			; define vector ref mag
>>> vcres at vcRefLengthF = 0.045 			; define length of vec ref
>>> vcres at vcGlyphStyle = "CurlyVector" 		; turn on curly vectors
>>> vcres at vcMinDistanceF = 0.017
>>> vcres at mpFillOn = False 				; turn off gray fill
>>> vcres at mpOutlineBoundarySets = "National" 	; turn on country boundaries
>>> ;vcres at mpFillBoundarySets = AllBoundaries
>>> vcres at mpGeophysicalLineColor = "Navy" 		; color of cont. outlines
>>> vcres at mpGeophysicalLineThicknessF = 1.5 	; thickness of outlines
>>>
>>>
>>> ;vcres at gsnLeftString = "DJF High NAO"
>>> ; was previously winds_hi =
>>> gsn_csm_vector_map_ce(wks,uAvgTime_hi,vAvgTime_hi,vcres)
>>> winds_hi = gsn_csm_vector(wks,uAvgTime_hi,vAvgTime_hi,vcres)
>>> ;vcres at gsnLeftString = "DJF Low NAO"
>>> winds_lo = gsn_csm_vector(wks,uAvgTime_lo,vAvgTime_lo,vcres)
>>> ;vcres at gsnLeftString = "Difference of High - Low"
>>> winds_diff = gsn_csm_vector(wks, diff_u, diff_v,vcres)
>>> ;************************************************
>>> ;---- geopotential height filled contour plot
>>> ;***********************************************
>>> zfres                      = res
>>> zfres at cnFillOn             = True
>>> ;zfres at cnLevelSelectionMode = "ExplicitLevels
>>> ;zfres at cnLevels             = ispan(-20,90,5)
>>> zfres at lbLabelFontHeightF   = 0.015
>>> zfres at lbOrientation        = "Vertical"
>>> zfres at pmLabelBarOrthogonalPosF = -0.005
>>> zfres at cnFillPalette = "BlWhRe"
>>>
>>> contour_zf_hi = gsn_csm_contour(wks,zAvgTime_hi,zfres)
>>> contour_zf_lo = gsn_csm_contour(wks,zAvgTime_lo,zfres)
>>> contour_zf_diff = gsn_csm_contour(wks,diff_z,zfres)
>>>
>>> plot(0) = gsn_csm_contour(wks,zAvgTime_hi,zfres)
>>> plot(1) = gsn_csm_contour(wks,zAvgTime_lo,zfres)
>>> plot(2) = gsn_csm_contour(wks,diff_z,zfres)
>>>
>>> ;overlay(plot(0),winds_hi)
>>> overlay(mpid,plot(0))
>>> overlay(mpid,winds_hi)
>>>
>>> ;overlay(plot(1),winds_lo)
>>> overlay(mpid,plot(1))
>>> overlay(mpid,winds_lo)
>>>
>>> ;overlay(plot(2),winds_diff)
>>> overlay(mpid,plot(2))
>>> overlay(mpid,winds_lo)
>>>
>>>
>>> ;************************************************
>>> ; create panel
>>> ;************************************************
>>> resP = True					; modify the panel plot
>>> resP at txString = "NAO 500mb"
>>> gsn_panel(wks,plot,(/3,1/),resP)		; now draw as one plot;
>>>
>>>
>>> end
>>>
>>>
>>>
>>>
>>>
>>>
>>>> Hi Mira,
>>>>
>>>> As for the color bar, I would suggest either using cnFillColors or
>>>> setting
>>>> the color table for the workstation, I usually just set the color
>>>> table
>>>> for workstation.
>>>>
>>>> After you open the workstation:
>>>>
>>>> wks = gsn_open_wks("ps","Panel_NAO_z_500")
>>>>
>>>> Do:
>>>> gsn_define_colormap(wks, “matlab_jet”)
>>>>
>>>> You can choose the color palette that you want and put the name in
>>>> quotes.
>>>>
>>>> As for the coastlines, maybe set:
>>>>
>>>> mpDataSetName = “Earth..4”
>>>>
>>>> You can also try setting mpFillBoundarySets = AllBoundaries and see if
>>>> that helps.
>>>>
>>>> Hope that helps,
>>>> -Alex
>>>>
>>>>
>>>>
>>>>
>>>>> On Jul 29, 2015, at 8:28 PM, mberdahl at envsci.rutgers.edu wrote:
>>>>>
>>>>> Hi Alex,
>>>>>
>>>>> Great!   Now I have three panels, in the correct region, with filled
>>>>> contours and overlaid wind vectors.  Thanks!  I have two outstanding
>>>>> problems still:  I'm trying to change the color bar for the filled
>>>>> contours, as the default in the older version of ncl is not very
>>>>> appealing.  When I try to do res at cnFillPalette = "matlab_jet", or any
>>>>> other choice for that matter, I receive the following warning:
>>>>>
>>>>> warning:cnFillPalette is not a valid resource in
>>>>> Panel_NAO_z_500_vector
>>>>> at
>>>>> this time
>>>>>
>>>>> I'm not exactly sure why I can't seem to change the contour fill
>>>>> palette.
>>>>>
>>>>> Secondly, I still don't see a map of coastlines on my figures.  Do
>>>>> you
>>>>> know how to add one?
>>>>>
>>>>> Thank you very much again,
>>>>> Mira
>>>>>
>>>>>
>>>>>> Hi Mira,
>>>>>>
>>>>>> You arenÂ’t seeing vectors or map because they a buried under the
>>>>>> filled
>>>>>> contours.  You may want to flip the plots around and do the contours
>>>>>> as
>>>>>> a
>>>>>> map and overlay the vectors on top of the contours.
>>>>>>
>>>>>> Upon further analysis of your script, you were doing the overlay on
>>>>>> contour plots but then using different contour plots for the panel
>>>>>> plot.
>>>>>> Below I have altered your script, try it this way and see how it
>>>>>> goes.
>>>>>>
>>>>>> I am not sure if you have made any more recent changes to this
>>>>>> script
>>>>>> but
>>>>>> this is using the last one you sent.
>>>>>>
>>>>>> Hope that helps,
>>>>>> -Alex
>>>>>>
>>>>>>
>>>>>> ;*****************************************
>>>>>> ; ; plot average winds overlayed on geopotental filled contours for
>>>>>> the
>>>>>> years that are extreme
>>>>>> ; the years are found with my matlab script findExtremeYrs.m
>>>>>> ; this particlar script does the years that correspond to high
>>>>>> correlation
>>>>>> with DJF, SE precip and NAO
>>>>>> ;*****************************************
>>>>>> ; the original data goes from 1948 January to April 2015.
>>>>>> ; I will cut it out to 2014 December, so we have something divisible
>>>>>> by
>>>>>> 12
>>>>>> so we can do seasonal averages...
>>>>>>
>>>>>>
>>>>>> 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/contributed.ncl"
>>>>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
>>>>>> ;************************************************
>>>>>> begin
>>>>>> ;************************************************
>>>>>> ; read in netCDF file s
>>>>>> ;************************************************
>>>>>> a = addfile("uwnd.mon.mean.alllevels.nc","r") ; u winds
>>>>>> b = addfile("vwnd.mon.mean.alllevels.nc","r") ; v winds
>>>>>>
>>>>>> c = addfile("../Geopotential/hgt.mon.mean.nc","r") ; geopotential
>>>>>> heights.
>>>>>>
>>>>>> ;************************************************
>>>>>> ; read in zonal [u] and meridional [v] winds (July)
>>>>>> ;************************************************
>>>>>>
>>>>>> u = a->uwnd(0:803,{500},{45:90},{270:357.5})
>>>>>> v = b->vwnd(0:803,{500},{45:90},{270:357.5}) ; Get u, v, time
>>>>>> (1),level
>>>>>> (1000hpa),latitude(-90:90) and longitude(0:360) data.
>>>>>> z = c->hgt(0:803,{500},{45:90},{270:357.5})  ; get geopotenial
>>>>>> heights...
>>>>>>
>>>>>> printVarSummary(u)
>>>>>> printVarSummary(v)
>>>>>> printVarSummary(z)
>>>>>>
>>>>>> ; Calculate the seasonal averages.
>>>>>> uDJF = month_to_season(u, "DJF")
>>>>>> vDJF = month_to_season(v, "DJF")
>>>>>> zDJF = month_to_season(z, "DJF")
>>>>>>
>>>>>> printVarSummary(uDJF)
>>>>>> printVarSummary(vDJF)
>>>>>> printVarSummary(zDJF)
>>>>>>
>>>>>> ; from the matlab script i wrote: findExtremeYrs, i pulled out the
>>>>>> extreme
>>>>>> years (> or < 1std) that i want to average and plot here.
>>>>>>
>>>>>> ; for ans =   4 (NAO)
>>>>>> ; yearList_hi = 1973        1975        1983        1989        1995
>>>>>> 2000        2007        2012
>>>>>> ; yearList_lo = 1963        1964        1965        1969        1977
>>>>>> 1979        1996        1997        2010        2011
>>>>>>
>>>>>>
>>>>>> ; this data starts at 1948 (this is index 0), so 1953=5, 1963=15
>>>>>> etc.
>>>>>>
>>>>>> uDJF_NAO_hi = uDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>> uDJF_NAO_lo = uDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>
>>>>>> vDJF_NAO_hi = vDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>> vDJF_NAO_lo = vDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>
>>>>>> zDJF_NAO_hi = zDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>> zDJF_NAO_lo = zDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>
>>>>>>
>>>>>> uAvgTime_hi = dim_avg_n_Wrap(uDJF_NAO_hi,0)
>>>>>> uAvgTime_lo = dim_avg_n_Wrap(uDJF_NAO_lo,0)
>>>>>>
>>>>>> printVarSummary(uAvgTime_hi)
>>>>>> printVarSummary(uAvgTime_lo)
>>>>>>
>>>>>> vAvgTime_hi = dim_avg_n_Wrap(vDJF_NAO_hi,0)
>>>>>> vAvgTime_lo = dim_avg_n_Wrap(vDJF_NAO_lo,0)
>>>>>>
>>>>>> printVarSummary(vAvgTime_hi)
>>>>>> printVarSummary(vAvgTime_lo)
>>>>>>
>>>>>> zAvgTime_hi = dim_avg_n_Wrap(zDJF_NAO_hi,0)
>>>>>> zAvgTime_lo = dim_avg_n_Wrap(zDJF_NAO_lo,0)
>>>>>>
>>>>>> printVarSummary(zAvgTime_hi)
>>>>>> printVarSummary(zAvgTime_lo)
>>>>>>
>>>>>> ; dirty way to copy metadata over first.
>>>>>> diff_u = uAvgTime_hi;
>>>>>> diff_v = vAvgTime_hi;
>>>>>> diff_z = zAvgTime_hi;
>>>>>>
>>>>>> diff_u = uAvgTime_hi - uAvgTime_lo
>>>>>> diff_v = vAvgTime_hi - vAvgTime_lo
>>>>>> diff_z = zAvgTime_hi - zAvgTime_lo
>>>>>>
>>>>>> printVarSummary(diff_u)
>>>>>> printVarSummary(diff_v)
>>>>>> printVarSummary(diff_z)
>>>>>>
>>>>>> ;************************************************
>>>>>> ; create plot
>>>>>> ;************************************************
>>>>>> wks = gsn_open_wks("ps","Panel_NAO_z_500") 		; open a ps file
>>>>>> plot = new(3,graphic)				; create a plot array
>>>>>>
>>>>>> ;---- set common resources for all plots
>>>>>> res 			= True
>>>>>> res at gsnDraw		= False 		; dont draw
>>>>>> res at gsnFrame		= False			; dont advance frame
>>>>>> res at cnInfoLabelOn	= False			; trn off cn info label
>>>>>> res at gsnAddCyclic = False			; has to do with wrapping the longitude
>>>>>> at
>>>>>> 0/360
>>>>>> ;************************************************
>>>>>> ; Choose a subregion
>>>>>> ;************************************************
>>>>>> res at mpMaxLatF = 90 				;maximum latitude
>>>>>> res at mpMinLatF = 45 				;minimum latitude
>>>>>> res at mpMaxLonF = 357.5 			;	;maximum longitude
>>>>>> res at mpMinLonF = 270 				;minimum longitude
>>>>>>
>>>>>> ;************************************************
>>>>>> ; Map Options
>>>>>> ;************************************************
>>>>>> res at mpOutlineBoundarySets = "National" 	; turn on country boundaries
>>>>>> res at mpGeophysicalLineColor = "Navy" 		; color of cont. outlines
>>>>>> res at mpGeophysicalLineThicknessF = 1.5 	; thickness of outlines
>>>>>>
>>>>>> ;***********************************************
>>>>>> ; ----wind  vector plot
>>>>>> ;***********************************************
>>>>>> vcres = res
>>>>>> vcres at vcRefAnnoOrthogonalPosF = -1.0 		; move ref vector up
>>>>>> vcres at vcRefMagnitudeF = 10.0 			; define vector ref mag
>>>>>> vcres at vcRefLengthF = 0.045 			; define length of vec ref
>>>>>> vcres at vcGlyphStyle = "CurlyVector" 		; turn on curly vectors
>>>>>> vcres at vcMinDistanceF = 0.017
>>>>>> vcres at mpFillOn = False 				; turn off gray fill
>>>>>>
>>>>>> ;vcres at gsnLeftString = "DJF High NAO"
>>>>>> winds_hi = gsn_csm_vector(wks,uAvgTime_hi,vAvgTime_hi,vcres)
>>>>>> ;vcres at gsnLeftString = "DJF Low NAO"
>>>>>> winds_lo = gsn_csm_vector(wks,uAvgTime_lo,vAvgTime_lo,vcres)
>>>>>> ;vcres at gsnLeftString = "Difference of High - Low"
>>>>>> winds_diff = gsn_csm_vector(wks, diff_u, diff_v,vcres)
>>>>>> ;************************************************
>>>>>> ;---- geopotential height filled contour plot
>>>>>> ;***********************************************
>>>>>> zfres                      = res
>>>>>> zfres at cnFillOn             = True
>>>>>> ;zfres at cnLevelSelectionMode = "ExplicitLevels
>>>>>> ;zfres at cnLevels             = ispan(-20,90,5)
>>>>>> zfres at lbLabelFontHeightF   = 0.015
>>>>>> zfres at lbOrientation        = "Vertical"
>>>>>> zfres at pmLabelBarOrthogonalPosF = -0.005
>>>>>>
>>>>>> plot(0) = gsn_csm_contour_map_ce(wks,zAvgTime_hi,zfres)
>>>>>> plot(1) = gsn_csm_contour_map_ce(wks,zAvgTime_lo,zfres)
>>>>>> plot(2) = gsn_csm_contour_map_ce(wks,diff_z,zfres)
>>>>>>
>>>>>> overlay(plot(0),winds_hi)
>>>>>> overlay(plot(1),winds_lo)
>>>>>> overlay(plot(2),winds_diff)
>>>>>>
>>>>>> ;************************************************
>>>>>> ; create panel
>>>>>> ;************************************************
>>>>>> resP = True					; modify the panel plot
>>>>>> resP at txString = "NAO 500mb"
>>>>>> gsn_panel(wks,plot,(/3,1/),resP)		; now draw as one plot;
>>>>>>
>>>>>>
>>>>>> end
>>>>>>
>>>>>>
>>>>>>> On Jul 29, 2015, at 1:57 PM, mberdahl at envsci.rutgers.edu wrote:
>>>>>>>
>>>>>>> Hi Alex,
>>>>>>> Great - thanks for the suggestions!  That does help, but does not
>>>>>>> totally
>>>>>>> solve my problem.  Now I get three panel plots of the correct
>>>>>>> region
>>>>>>> and
>>>>>>> all on one page which is great.  But only the filled contours show.
>>>>>>> There
>>>>>>> are no overlaid vectors or map of coastal boundaries.  Any
>>>>>>> suggestions
>>>>>>> would be greatly appreciated.
>>>>>>> Thanks so much!
>>>>>>> Mira
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>> Hi Mira,
>>>>>>>>
>>>>>>>> I think the issue here is that you set a bunch of resources as res
>>>>>>>> but
>>>>>>>> then set vres to True instead of res.  I think switching it to
>>>>>>>> vres
>>>>>>>> =
>>>>>>>> res
>>>>>>>> will give you the maps you were expecting.  I see you also set res
>>>>>>>> =
>>>>>>>> True
>>>>>>>> twice, might be helpful to remove the double.
>>>>>>>>
>>>>>>>> Hope that helps,
>>>>>>>> -Alex
>>>>>>>>
>>>>>>>>
>>>>>>>>> On Jul 27, 2015, at 9:39 PM, mberdahl at envsci.rutgers.edu wrote:
>>>>>>>>>
>>>>>>>>> Hi Alex,
>>>>>>>>> Thanks for the suggestion.  I tried this and now I do not have
>>>>>>>>> the
>>>>>>>>> same
>>>>>>>>> error message any more, so that is good.  However, I am now just
>>>>>>>>> getting
>>>>>>>>> 3
>>>>>>>>> blank plots of the whole world, even though I'd like plots with
>>>>>>>>> overlaid
>>>>>>>>> geopotential height and winds, for just a region over Greenland.
>>>>>>>>> These
>>>>>>>>> are
>>>>>>>>> followed by a fourth figure that has filled contours of the
>>>>>>>>> correct
>>>>>>>>> region, but no map or wind vectors.  I'll attach the figures as
>>>>>>>>> an
>>>>>>>>> attachment here.  They are not panel plots as I want, but each
>>>>>>>>> figure
>>>>>>>>> is
>>>>>>>>> on its own page.  I'll copy the updated code below.
>>>>>>>>> Any thoughts?
>>>>>>>>> Thanks!
>>>>>>>>> Mira
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> ;*****************************************
>>>>>>>>> ; ; plot average winds overlayed on geopotental filled contours
>>>>>>>>> for
>>>>>>>>> the
>>>>>>>>> years that are extreme
>>>>>>>>> ; the years are found with my matlab script findExtremeYrs.m
>>>>>>>>> ; this particlar script does the years that correspond to high
>>>>>>>>> correlation
>>>>>>>>> with DJF, SE precip and NAO
>>>>>>>>> ;*****************************************
>>>>>>>>> ; the original data goes from 1948 January to April 2015.
>>>>>>>>> ; I will cut it out to 2014 December, so we have something
>>>>>>>>> divisible
>>>>>>>>> by
>>>>>>>>> 12
>>>>>>>>> so we can do seasonal averages...
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> 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/contributed.ncl"
>>>>>>>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
>>>>>>>>> ;************************************************
>>>>>>>>> begin
>>>>>>>>> ;************************************************
>>>>>>>>> ; read in netCDF file s
>>>>>>>>> ;************************************************
>>>>>>>>> a = addfile("uwnd.mon.mean.alllevels.nc","r") ; u winds
>>>>>>>>> b = addfile("vwnd.mon.mean.alllevels.nc","r") ; v winds
>>>>>>>>>
>>>>>>>>> c = addfile("../Geopotential/hgt.mon.mean.nc","r") ; geopotential
>>>>>>>>> heights.
>>>>>>>>>
>>>>>>>>> ;************************************************
>>>>>>>>> ; read in zonal [u] and meridional [v] winds (July)
>>>>>>>>> ;************************************************
>>>>>>>>>
>>>>>>>>> u = a->uwnd(0:803,{500},{45:90},{270:357.5})
>>>>>>>>> v = b->vwnd(0:803,{500},{45:90},{270:357.5}) ; Get u, v, time
>>>>>>>>> (1),level
>>>>>>>>> (1000hpa),latitude(-90:90) and longitude(0:360) data.
>>>>>>>>> z = c->hgt(0:803,{500},{45:90},{270:357.5})  ; get geopotenial
>>>>>>>>> heights...
>>>>>>>>>
>>>>>>>>> printVarSummary(u)
>>>>>>>>> printVarSummary(v)
>>>>>>>>> printVarSummary(z)
>>>>>>>>>
>>>>>>>>> ; Calculate the seasonal averages.
>>>>>>>>> uDJF = month_to_season(u, "DJF")
>>>>>>>>> vDJF = month_to_season(v, "DJF")
>>>>>>>>> zDJF = month_to_season(z, "DJF")
>>>>>>>>>
>>>>>>>>> printVarSummary(uDJF)
>>>>>>>>> printVarSummary(vDJF)
>>>>>>>>> printVarSummary(zDJF)
>>>>>>>>>
>>>>>>>>> ; from the matlab script i wrote: findExtremeYrs, i pulled out
>>>>>>>>> the
>>>>>>>>> extreme
>>>>>>>>> years (> or < 1std) that i want to average and plot here.
>>>>>>>>>
>>>>>>>>> ; for ans =   4 (NAO)
>>>>>>>>> ; yearList_hi = 1973        1975        1983        1989
>>>>>>>>> 1995
>>>>>>>>> 2000        2007        2012
>>>>>>>>> ; yearList_lo = 1963        1964        1965        1969
>>>>>>>>> 1977
>>>>>>>>> 1979        1996        1997        2010        2011
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> ; this data starts at 1948 (this is index 0), so 1953=5, 1963=15
>>>>>>>>> etc.
>>>>>>>>>
>>>>>>>>> uDJF_NAO_hi = uDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>>>>> uDJF_NAO_lo = uDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>>>>
>>>>>>>>> vDJF_NAO_hi = vDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>>>>> vDJF_NAO_lo = vDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>>>>
>>>>>>>>> zDJF_NAO_hi = zDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>>>>> zDJF_NAO_lo = zDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> uAvgTime_hi = dim_avg_n_Wrap(uDJF_NAO_hi,0)
>>>>>>>>> uAvgTime_lo = dim_avg_n_Wrap(uDJF_NAO_lo,0)
>>>>>>>>>
>>>>>>>>> printVarSummary(uAvgTime_hi)
>>>>>>>>> printVarSummary(uAvgTime_lo)
>>>>>>>>>
>>>>>>>>> vAvgTime_hi = dim_avg_n_Wrap(vDJF_NAO_hi,0)
>>>>>>>>> vAvgTime_lo = dim_avg_n_Wrap(vDJF_NAO_lo,0)
>>>>>>>>>
>>>>>>>>> printVarSummary(vAvgTime_hi)
>>>>>>>>> printVarSummary(vAvgTime_lo)
>>>>>>>>>
>>>>>>>>> zAvgTime_hi = dim_avg_n_Wrap(zDJF_NAO_hi,0)
>>>>>>>>> zAvgTime_lo = dim_avg_n_Wrap(zDJF_NAO_lo,0)
>>>>>>>>>
>>>>>>>>> printVarSummary(zAvgTime_hi)
>>>>>>>>> printVarSummary(zAvgTime_lo)
>>>>>>>>>
>>>>>>>>> ; dirty way to copy metadata over first.
>>>>>>>>> diff_u = uAvgTime_hi;
>>>>>>>>> diff_v = vAvgTime_hi;
>>>>>>>>> diff_z = zAvgTime_hi;
>>>>>>>>>
>>>>>>>>> diff_u = uAvgTime_hi - uAvgTime_lo
>>>>>>>>> diff_v = vAvgTime_hi - vAvgTime_lo
>>>>>>>>> diff_z = zAvgTime_hi - zAvgTime_lo
>>>>>>>>>
>>>>>>>>> printVarSummary(diff_u)
>>>>>>>>> printVarSummary(diff_v)
>>>>>>>>> printVarSummary(diff_z)
>>>>>>>>>
>>>>>>>>> ;************************************************
>>>>>>>>> ; create plot
>>>>>>>>> ;************************************************
>>>>>>>>> wks = gsn_open_wks("ps","Panel_NAO_z_500") 		; open a ps file
>>>>>>>>> plot = new(3,graphic)				; create a plot array
>>>>>>>>>
>>>>>>>>> ;---- set common resources for all plots
>>>>>>>>> res 			= True
>>>>>>>>> res at gsnDraw		= False 		; dont draw
>>>>>>>>> res at gsnFrame		= False			; dont advance frame
>>>>>>>>> res at cnInfoLabelOn	= False			; trn off cn info label
>>>>>>>>> res = True 					; plot mods desired
>>>>>>>>> res at gsnAddCyclic = False			; has to do with wrapping the
>>>>>>>>> longitude
>>>>>>>>> at
>>>>>>>>> 0/360
>>>>>>>>> ;************************************************
>>>>>>>>> ; Choose a subregion
>>>>>>>>> ;************************************************
>>>>>>>>> res at mpMaxLatF = 90 				;maximum latitude
>>>>>>>>> res at mpMinLatF = 45 				;minimum latitude
>>>>>>>>> res at mpMaxLonF = 357.5 			;	;maximum longitude
>>>>>>>>> res at mpMinLonF = 270 				;minimum longitude
>>>>>>>>>
>>>>>>>>> ;***********************************************
>>>>>>>>> ; ----wind  vector plot
>>>>>>>>> ;***********************************************
>>>>>>>>> vcres = True
>>>>>>>>> vcres at vcRefAnnoOrthogonalPosF = -1.0 		; move ref vector up
>>>>>>>>> vcres at vcRefMagnitudeF = 10.0 			; define vector ref mag
>>>>>>>>> vcres at vcRefLengthF = 0.045 			; define length of vec ref
>>>>>>>>> vcres at vcGlyphStyle = "CurlyVector" 		; turn on curly vectors
>>>>>>>>> vcres at vcMinDistanceF = 0.017
>>>>>>>>> vcres at mpFillOn = False 				; turn off gray fill
>>>>>>>>> vcres at mpOutlineBoundarySets = "National" 	; turn on country
>>>>>>>>> boundaries
>>>>>>>>> vcres at mpGeophysicalLineColor = "Navy" 		; color of cont. outlines
>>>>>>>>> vcres at mpGeophysicalLineThicknessF = 1.5 	; thickness of outlines
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> ;vcres at gsnLeftString = "DJF High NAO"
>>>>>>>>> winds_hi =
>>>>>>>>> gsn_csm_vector_map_ce(wks,uAvgTime_hi,vAvgTime_hi,vcres)
>>>>>>>>> ;vcres at gsnLeftString = "DJF Low NAO"
>>>>>>>>> winds_lo =
>>>>>>>>> gsn_csm_vector_map_ce(wks,uAvgTime_lo,vAvgTime_lo,vcres)
>>>>>>>>> ;vcres at gsnLeftString = "Difference of High - Low"
>>>>>>>>> winds_diff = gsn_csm_vector_map_ce(wks, diff_u, diff_v,vcres)
>>>>>>>>> ;************************************************
>>>>>>>>> ;---- geopotential height filled contour plot
>>>>>>>>> ;***********************************************
>>>>>>>>> zfres                      = res
>>>>>>>>> zfres at cnFillOn             = True
>>>>>>>>> ;zfres at cnLevelSelectionMode = "ExplicitLevels
>>>>>>>>> ;zfres at cnLevels             = ispan(-20,90,5)
>>>>>>>>> zfres at lbLabelFontHeightF   = 0.015
>>>>>>>>> zfres at lbOrientation        = "Vertical"
>>>>>>>>> zfres at pmLabelBarOrthogonalPosF = -0.005
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> contour_zf_hi = gsn_csm_contour(wks,zAvgTime_hi,zfres)
>>>>>>>>> contour_zf_lo = gsn_csm_contour(wks,zAvgTime_lo,zfres)
>>>>>>>>> contour_zf_diff = gsn_csm_contour(wks,diff_z,zfres)
>>>>>>>>>
>>>>>>>>> plot(0) = gsn_csm_contour(wks,zAvgTime_hi,zfres)
>>>>>>>>> plot(1) = gsn_csm_contour(wks,zAvgTime_lo,zfres)
>>>>>>>>> plot(2) = gsn_csm_contour(wks,diff_z,zfres)
>>>>>>>>>
>>>>>>>>> overlay(contour_zf_hi,winds_hi)
>>>>>>>>> overlay(contour_zf_lo,winds_lo)
>>>>>>>>> overlay(contour_zf_diff,winds_diff)
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> ;************************************************
>>>>>>>>> ; create panel
>>>>>>>>> ;************************************************
>>>>>>>>> resP = True					; modify the panel plot
>>>>>>>>> resP at txString = "NAO 500mb"
>>>>>>>>> gsn_panel(wks,plot,(/3,1/),resP)		; now draw as one plot;
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> end
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>> Mira,
>>>>>>>>>>
>>>>>>>>>> Try this:
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> plot(0) = gsn_csm_contour(wks,zAvgTime_hi,zfres)
>>>>>>>>>> plot(1) = gsn_csm_contour(wks,zAvgTime_lo,zfres)
>>>>>>>>>> plot(2) = gsn_csm_contour(wks,diff_z,zfres)
>>>>>>>>>>
>>>>>>>>>> overlay(contour_zf_hi,winds_hi)
>>>>>>>>>> overlay(contour_zf_lo,winds_lo)
>>>>>>>>>> overlay(contour_diff,winds_diff)
>>>>>>>>>>
>>>>>>>>>> I don't believe overlay returns anything and you are trying to
>>>>>>>>>> save
>>>>>>>>>> it
>>>>>>>>>> as
>>>>>>>>>> a
>>>>>>>>>> graphic. The code above would replace the lines in your script.
>>>>>>>>>>
>>>>>>>>>> Hope that helps,
>>>>>>>>>>
>>>>>>>>>> -Alex
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> On Sunday, July 26, 2015, <mberdahl at envsci.rutgers.edu
>>>>>>>>>> <javascript:_e(%7B%7D,'cvml','mberdahl at envsci.rutgers.edu
>>>>>>>>>> <mailto:mberdahl at envsci.rutgers.edu
>>>>>>>>>> <mailto:mberdahl at envsci.rutgers.edu
>>>>>>>>>> <mailto:mberdahl at envsci.rutgers.edu>>>');>> wrote:
>>>>>>>>>>
>>>>>>>>>>> Hi all,
>>>>>>>>>>>
>>>>>>>>>>> I'm trying to adapt a working code I had that was plotting 3
>>>>>>>>>>> panel
>>>>>>>>>>> plots
>>>>>>>>>>> of wind vectors in a certain region.  Now, all I want to do is
>>>>>>>>>>> add
>>>>>>>>>>> an
>>>>>>>>>>> overlay of geopotential height (filled contours) underneath
>>>>>>>>>>> each
>>>>>>>>>>> of
>>>>>>>>>>> these
>>>>>>>>>>> plots, but the way I'm approaching it is wrong.  I've tried
>>>>>>>>>>> simplifying
>>>>>>>>>>> it
>>>>>>>>>>> to do this without panels, just a single plot, but I still
>>>>>>>>>>> don't
>>>>>>>>>>> get
>>>>>>>>>>> what
>>>>>>>>>>> I need.  The error message I am getting now is:
>>>>>>>>>>>
>>>>>>>>>>> fatal:syntax error: line 153 in file plotWinds_z_NAO_level.ncl
>>>>>>>>>>> before
>>>>>>>>>>> or
>>>>>>>>>>> near overlay
>>>>>>>>>>> plot(0)= overlay
>>>>>>>>>>> ---------------^
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> My code is copied below.  Any suggestions are greatly
>>>>>>>>>>> appreciated.
>>>>>>>>>>> Thanks,
>>>>>>>>>>> Mira
>>>>>>>>>>>
>>>>>>>>>>> ;*****************************************
>>>>>>>>>>> ; ; plot average winds overlayed on geopotental filled contours
>>>>>>>>>>> for
>>>>>>>>>>> the
>>>>>>>>>>> years that are extreme
>>>>>>>>>>> ; the years are found with my matlab script findExtremeYrs.m
>>>>>>>>>>> ; this particlar script does the years that correspond to high
>>>>>>>>>>> correlation
>>>>>>>>>>> with DJF, SE precip and NAO
>>>>>>>>>>> ;*****************************************
>>>>>>>>>>> ; the original data goes from 1948 January to April 2015.
>>>>>>>>>>> ; I will cut it out to 2014 December, so we have something
>>>>>>>>>>> divisible
>>>>>>>>>>> by
>>>>>>>>>>> 12
>>>>>>>>>>> so we can do seasonal averages...
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> 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/contributed.ncl"
>>>>>>>>>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> begin
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> ; read in netCDF file s
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> a = addfile("uwnd.mon.mean.alllevels.nc","r") ; u winds
>>>>>>>>>>> b = addfile("vwnd.mon.mean.alllevels.nc","r") ; v winds
>>>>>>>>>>>
>>>>>>>>>>> c = addfile("../Geopotential/hgt.mon.mean.nc","r") ;
>>>>>>>>>>> geopotential
>>>>>>>>>>> heights.
>>>>>>>>>>>
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> ; read in zonal [u] and meridional [v] winds (July)
>>>>>>>>>>> ;************************************************
>>>>>>>>>>>
>>>>>>>>>>> u = a->uwnd(0:803,{500},{45:90},{270:357.5})
>>>>>>>>>>> v = b->vwnd(0:803,{500},{45:90},{270:357.5}) ; Get u, v, time
>>>>>>>>>>> (1),level
>>>>>>>>>>> (1000hpa),latitude(-90:90) and longitude(0:360) data.
>>>>>>>>>>> z = c->hgt(0:803,{500},{45:90},{270:357.5})  ; get geopotenial
>>>>>>>>>>> heights...
>>>>>>>>>>>
>>>>>>>>>>> printVarSummary(u)
>>>>>>>>>>> printVarSummary(v)
>>>>>>>>>>> printVarSummary(z)
>>>>>>>>>>>
>>>>>>>>>>> ; Calculate the seasonal averages.
>>>>>>>>>>> uDJF = month_to_season(u, "DJF")
>>>>>>>>>>> vDJF = month_to_season(v, "DJF")
>>>>>>>>>>> zDJF = month_to_season(z, "DJF")
>>>>>>>>>>>
>>>>>>>>>>> printVarSummary(uDJF)
>>>>>>>>>>> printVarSummary(vDJF)
>>>>>>>>>>> printVarSummary(zDJF)
>>>>>>>>>>>
>>>>>>>>>>> ; from the matlab script i wrote: findExtremeYrs, i pulled out
>>>>>>>>>>> the
>>>>>>>>>>> extreme
>>>>>>>>>>> years (> or < 1std) that i want to average and plot here.
>>>>>>>>>>>
>>>>>>>>>>> ; for ans =   4 (NAO)
>>>>>>>>>>> ; yearList_hi = 1973        1975        1983        1989
>>>>>>>>>>> 1995
>>>>>>>>>>> 2000        2007        2012
>>>>>>>>>>> ; yearList_lo = 1963        1964        1965        1969
>>>>>>>>>>> 1977
>>>>>>>>>>> 1979        1996        1997        2010        2011
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> ; this data starts at 1948 (this is index 0), so 1953=5,
>>>>>>>>>>> 1963=15
>>>>>>>>>>> etc.
>>>>>>>>>>>
>>>>>>>>>>> uDJF_NAO_hi = uDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>>>>>>> uDJF_NAO_lo = uDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>>>>>>
>>>>>>>>>>> vDJF_NAO_hi = vDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>>>>>>> vDJF_NAO_lo = vDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>>>>>>
>>>>>>>>>>> zDJF_NAO_hi = zDJF((/25,27,35,41,47,52,59,64/),:,:)
>>>>>>>>>>> zDJF_NAO_lo = zDJF((/15,16,17,21,29,31,48,49,62,63/),:,:)
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> uAvgTime_hi = dim_avg_n_Wrap(uDJF_NAO_hi,0)
>>>>>>>>>>> uAvgTime_lo = dim_avg_n_Wrap(uDJF_NAO_lo,0)
>>>>>>>>>>>
>>>>>>>>>>> printVarSummary(uAvgTime_hi)
>>>>>>>>>>> printVarSummary(uAvgTime_lo)
>>>>>>>>>>>
>>>>>>>>>>> vAvgTime_hi = dim_avg_n_Wrap(vDJF_NAO_hi,0)
>>>>>>>>>>> vAvgTime_lo = dim_avg_n_Wrap(vDJF_NAO_lo,0)
>>>>>>>>>>>
>>>>>>>>>>> printVarSummary(vAvgTime_hi)
>>>>>>>>>>> printVarSummary(vAvgTime_lo)
>>>>>>>>>>>
>>>>>>>>>>> zAvgTime_hi = dim_avg_n_Wrap(zDJF_NAO_hi,0)
>>>>>>>>>>> zAvgTime_lo = dim_avg_n_Wrap(zDJF_NAO_lo,0)
>>>>>>>>>>>
>>>>>>>>>>> printVarSummary(zAvgTime_hi)
>>>>>>>>>>> printVarSummary(zAvgTime_lo)
>>>>>>>>>>>
>>>>>>>>>>> ; dirty way to copy metadata over first.
>>>>>>>>>>> diff_u = uAvgTime_hi;
>>>>>>>>>>> diff_v = vAvgTime_hi;
>>>>>>>>>>> diff_z = zAvgTime_hi;
>>>>>>>>>>>
>>>>>>>>>>> diff_u = uAvgTime_hi - uAvgTime_lo
>>>>>>>>>>> diff_v = vAvgTime_hi - vAvgTime_lo
>>>>>>>>>>> diff_z = zAvgTime_hi - zAvgTime_lo
>>>>>>>>>>>
>>>>>>>>>>> printVarSummary(diff_u)
>>>>>>>>>>> printVarSummary(diff_v)
>>>>>>>>>>> printVarSummary(diff_z)
>>>>>>>>>>>
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> ; create plot
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> wks = gsn_open_wks("ps","Panel_NAO_z_500")              ; open
>>>>>>>>>>> a
>>>>>>>>>>> ps
>>>>>>>>>>> file
>>>>>>>>>>> plot = new(3,graphic)                           ; create a plot
>>>>>>>>>>> array
>>>>>>>>>>>
>>>>>>>>>>> ;---- set common resources for all plots
>>>>>>>>>>> res                     = True
>>>>>>>>>>> res at gsnDraw             = False                 ; dont draw
>>>>>>>>>>> res at gsnFrame            = False                 ; dont advance
>>>>>>>>>>> frame
>>>>>>>>>>> res at cnInfoLabelOn       = False                 ; trn off cn
>>>>>>>>>>> info
>>>>>>>>>>> label
>>>>>>>>>>> res = True                                      ; plot mods
>>>>>>>>>>> desired
>>>>>>>>>>> res at gsnAddCyclic = False                        ; has to do
>>>>>>>>>>> with
>>>>>>>>>>> wrapping
>>>>>>>>>>> the longitude at 0/360
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> ; Choose a subregion
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> res at mpMaxLatF = 90                              ;maximum
>>>>>>>>>>> latitude
>>>>>>>>>>> res at mpMinLatF = 45                              ;minimum
>>>>>>>>>>> latitude
>>>>>>>>>>> res at mpMaxLonF = 357.5                   ;       ;maximum
>>>>>>>>>>> longitude
>>>>>>>>>>> res at mpMinLonF = 270                             ;minimum
>>>>>>>>>>> longitude
>>>>>>>>>>>
>>>>>>>>>>> ;***********************************************
>>>>>>>>>>> ; ----wind  vector plot
>>>>>>>>>>> ;***********************************************
>>>>>>>>>>> vcres = True
>>>>>>>>>>> vcres at vcRefAnnoOrthogonalPosF = -1.0            ; move ref
>>>>>>>>>>> vector
>>>>>>>>>>> up
>>>>>>>>>>> vcres at vcRefMagnitudeF = 10.0                    ; define vector
>>>>>>>>>>> ref
>>>>>>>>>>> mag
>>>>>>>>>>> vcres at vcRefLengthF = 0.045                      ; define length
>>>>>>>>>>> of
>>>>>>>>>>> vec
>>>>>>>>>>> ref
>>>>>>>>>>> vcres at vcGlyphStyle = "CurlyVector"              ; turn on curly
>>>>>>>>>>> vectors
>>>>>>>>>>> vcres at vcMinDistanceF = 0.017
>>>>>>>>>>> vcres at mpFillOn = False                          ; turn off gray
>>>>>>>>>>> fill
>>>>>>>>>>> vcres at mpOutlineBoundarySets = "National"        ; turn on
>>>>>>>>>>> country
>>>>>>>>>>> boundaries
>>>>>>>>>>> vcres at mpGeophysicalLineColor = "Navy"           ; color of
>>>>>>>>>>> cont.
>>>>>>>>>>> outlines
>>>>>>>>>>> vcres at mpGeophysicalLineThicknessF = 1.5         ; thickness of
>>>>>>>>>>> outlines
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> ;vcres at gsnLeftString = "DJF High NAO"
>>>>>>>>>>> winds_hi =
>>>>>>>>>>> gsn_csm_vector_map_ce(wks,uAvgTime_hi,vAvgTime_hi,vcres)
>>>>>>>>>>> ;vcres at gsnLeftString = "DJF Low NAO"
>>>>>>>>>>> winds_lo =
>>>>>>>>>>> gsn_csm_vector_map_ce(wks,uAvgTime_lo,vAvgTime_lo,vcres)
>>>>>>>>>>> ;vcres at gsnLeftString = "Difference of High - Low"
>>>>>>>>>>> winds_diff = gsn_csm_vector_map_ce(wks, diff_u, diff_v,vcres)
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> ;---- geopotential height filled contour plot
>>>>>>>>>>> ;***********************************************
>>>>>>>>>>> zfres                      = res
>>>>>>>>>>> zfres at cnFillOn             = True
>>>>>>>>>>> ;zfres at cnLevelSelectionMode = "ExplicitLevels
>>>>>>>>>>> ;zfres at cnLevels             = ispan(-20,90,5)
>>>>>>>>>>> zfres at lbLabelFontHeightF   = 0.015
>>>>>>>>>>> zfres at lbOrientation        = "Vertical"
>>>>>>>>>>> zfres at pmLabelBarOrthogonalPosF = -0.005
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> contour_zf_hi = gsn_csm_contour(wks,zAvgTime_hi,zfres)
>>>>>>>>>>> contour_zf_lo = gsn_csm_contour(wks,zAvgTime_lo,zfres)
>>>>>>>>>>> contour_zf_diff = gsn_csm_contour(wks,diff_z,zfres)
>>>>>>>>>>>
>>>>>>>>>>> plot(0) = overlay(contour_zf_hi,winds_hi)
>>>>>>>>>>> plot(1) = overlay(contour_zf_lo,winds_lo)
>>>>>>>>>>> plot(2) = overlay(contour_diff,winds_diff)
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> ; create panel
>>>>>>>>>>> ;************************************************
>>>>>>>>>>> resP = True                                     ; modify the
>>>>>>>>>>> panel
>>>>>>>>>>> plot
>>>>>>>>>>> resP at txString = "NAO 500mb"
>>>>>>>>>>> gsn_panel(wks,plot,(/3,1/),resP)                ; now draw as
>>>>>>>>>>> one
>>>>>>>>>>> plot;
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> end
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> ncl-talk mailing list
>>>>>>>>>>> ncl-talk at ucar.edu
>>>>>>>>>>> List instructions, subscriber options, unsubscribe:
>>>>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>> <Panel_NAO_z_500.pdf>
>>>>
>>>>
>>>
>>
>>
>




More information about the ncl-talk mailing list