[ncl-talk] Overlay profile of topography

Mary Haley haley at ucar.edu
Sun May 7 09:08:31 MDT 2017


Hi,

If by "optimized for this plot" you mean is the plotting being done
correctly, then I would say "yes".

However, I did notice that you have two slightly different scales on your
color bar. The top one goes from 30 to 90 in steps of 5, while the bottom
one goes from 20 to 90 in steps of 5.

If somebody is not paying attention when looking at these images, they may
not realize that the bottom purple color represents values < 30 in the top
plot, while it represents values < 20 in the bottom plot.

I recommend using the same color scale for both plots, and then using a
single panel labelbar so that it's more clear what the colors mean.

To make sure you use the same contour levels for both plots, set these
resources:

 res at cnLevelSelectionMode = "ManualLevels" ; manually set contour levels
 res at cnMinLevelValF       = 20.            ; minimum contour level
 res at cnMaxLevelValF       = 90.            ; maximum contour level
 res at cnLevelSpacingF        = 5.	

You then you need to turn off the individual labelbars for both plots,
so you can set a common labelbar in the panel:

 res at lbLabelBarOn = False

Now, in the panel resources, set:

pres at gsnPanelLabelBar = True
pres2 at lbOrientation   = "vertical" ; set this if you want a vertical label bar

If you need more customization of your panel plot, then we have some
examples on our panel page:

http://www.ncl.ucar.edu/Applications/panel.shtml

--Mary




On Fri, May 5, 2017 at 4:07 AM, Appo derbetini <appopson4 at gmail.com> wrote:

> Dear Mary,
> I applied your suggestions.
> Now, things are going like a charm.
> Here attached figure obtained and script.
> I guess that this script is optimized  for this plot?
>
> Thank you very much
>
> begin
>
>     in1  = addfile("uvwr.mean.JJAS.nc", "r")                         ;
> open netcdf file
>     input4  = addfile("topo.nc", "r")
>
>    u00     =  in1->u
>    v00     =  in1->v
>    w00     =  in1->w
>    rhum00  =  in1->r
>
>    topo1 = input4->HT
>
>
>    lat = in1->latitude
>    lon = in1->longitude
>
>
>    dv00    = uv2dvF (u00, v00)
>    uvd00   = dv2uvF (dv00)
>    ucomp_div00 = uvd00(0, :, : , :, :)
>    copy_VarCoords(u00, ucomp_div00)
>
>
>    sellat = 3.0
>    latmin = -5.0
>    latmax = 3.0
>    lonmin = min(lon)
>    lonmax = 20.0
>
>    sellat2 = -1.0
>    latmin2 = -4.0
>    latmax2 = -1.0
>    lonmin = min(lon)
>    lonmax = 20.0
>
>    udiv00 = dim_avg_n_Wrap(ucomp_div00(:, :, {latmin:latmax},
> {lonmin:lonmax}), 2)
>    omega00 = dim_avg_n_Wrap(w00(:, :, {latmin:latmax}, {lonmin:lonmax}), 2)
>    r00 = dim_avg_n_Wrap(rhum00(:, :, {latmin:latmax},{lonmin:lonmax}), 2)
>
>    udiv002 = dim_avg_n_Wrap(ucomp_div00(:, :, {latmin2:latmax2},
> {lonmin:lonmax}), 2)
>    omega002 = dim_avg_n_Wrap(w00(:, :, {latmin2:latmax2},
> {lonmin:lonmax}), 2)
>    r002 = dim_avg_n_Wrap(rhum00(:, :, {latmin2:latmax2},{lonmin:lonmax}),
> 2)
>
>
>
>    elev = 1013.25*(1 - topo1*0.0065/288.15)^5.25145
>
>    copy_VarCoords(topo1, elev)
>
> ;  topo = dim_avg_n_Wrap(elev({latmin:latmax},{lonmin:lonmax}), 0)
>
>    topo = elev({sellat}, {lonmin:lonmax})
>    topo2 = elev({sellat2}, {lonmin:lonmax})
>
> ;---create plot
>    wks   = gsn_open_wks ("eps", "Overlay_Rhum_udiv_JJAS_Timmean_5S3N")
> ; open workstation
>    plot = new(2, graphic)
>    res                     = True                  ; plot mods desired
>    res at gsnDraw             = False    ; turn off draw
>    res at gsnFrame            = False    ; turn off frame
>    res at gsnMaximize         = True
>
>    res at gsnRightString      = ""
>    res at gsnLeftString       = ""
>
>    cnres = res
>    cnres at cnFillOn          = True               ; turn on color fill
>    cnres at cnFillPalette     = "MPL_rainbow"
>    cnres at lbOrientation     = "Vertical"
>    cnres at pmLabelBarOrthogonalPosF = 0.08
>    cnres at tiYAxisString      = "Pressure (hPa)"
>
>    vcres = res
>    vcres at vcRefMagnitudeF = 10.0                ; define vector ref mag
>    vcres at vcRefLengthF    = 0.1              ; define length of vec ref
>    vcres at vcGlyphStyle    = "CurlyVector"      ; turn on curley vectors
>    vcres at vcMinDistanceF  = 0.04               ; thin out vectors
>    vcres at vcMapDirection  = False
>    vcres at vcLineArrowThicknessF   = 3.0
>    vcres at vcVectorDrawOrder  = "Draw"        ; draw vectors last
>
>
>
>    contour_plota1  =  gsn_csm_pres_hgt(wks, r00(0, :, {lonmin:lonmax}),
> cnres)
>    vector_plota1   = gsn_csm_vector(wks, udiv00(0, :, {lonmin:lonmax}),
> -200.0*omega00(0, :, {lonmin:lonmax}), vcres )
>
>    contour_plota2  =  gsn_csm_pres_hgt(wks, r002(0, :, {lonmin:lonmax}),
> cnres)
>    vector_plota2   = gsn_csm_vector(wks, udiv002(0, :, {lonmin:lonmax}),
> -200.0*omega002(0, :, {lonmin:lonmax}), vcres )
>
>
> ;---Add topo field using a filled polygon.
>    getvalues contour_plota1
>       "trYMinF" : ymin
>       "trYMaxF" : ymax
>    end getvalues
>
> ;---Create new X,Y arrays that form a closed polygon.
>    nlon  = dimsizes(topo&lon)
>    xtopo = new(nlon+3,typeof(topo&lon))
>    ytopo = new(nlon+3,typeof(topo))
>
>    xtopo(0:nlon-1) = topo&lon
>    ytopo(0:nlon-1) = topo
>    xtopo(nlon)     = topo&lon(nlon-1)
>    ytopo(nlon)     = ymax          ; Use actual Y max of contour plot
>    xtopo(nlon+1)   = topo&lon(0)
>    ytopo(nlon+1)   = ymax
>    xtopo(nlon+2)   = topo&lon(0)   ; This last point closes
>    ytopo(nlon+2)   = topo(0)       ; the polygon.
>
>
>
>
> ;---Add topo field using a filled polygon.
>    getvalues contour_plota2
>       "trYMinF" : ymin2
>       "trYMaxF" : ymax2
>    end getvalues
>
> ;---Create new X,Y arrays that form a closed polygon.
>    nlon  = dimsizes(topo2&lon)
>    xtopo2 = new(nlon+3,typeof(topo2&lon))
>    ytopo2 = new(nlon+3,typeof(topo2))
>
>    xtopo2(0:nlon-1) = topo2&lon
>    ytopo2(0:nlon-1) = topo2
>    xtopo2(nlon)     = topo2&lon(nlon-1)
>    ytopo2(nlon)     = ymax2          ; Use actual Y max of contour plot
>    xtopo2(nlon+1)   = topo2&lon(0)
>    ytopo2(nlon+1)   = ymax2
>    xtopo2(nlon+2)   = topo2&lon(0)   ; This last point closes
>    ytopo2(nlon+2)   = topo2(0)       ; the polygon.
>
> ;---Add the polygon to the contour plot.
>    gnres = True
>    gnres at gsFillColor = "gray20"
>
>    id = gsn_add_polygon(wks, contour_plota1, xtopo, ytopo, gnres)
>    i2 = gsn_add_polygon(wks, contour_plota2, xtopo2, ytopo2, gnres)
>
>     overlay(contour_plota1, vector_plota1)
>     overlay(contour_plota2, vector_plota2)
>
>      plot(0)  = contour_plota1
>      plot(1)  = contour_plota2
>
>    ;draw(contour_plot)     ; This draws everything
>
>
>
>
>    pres = True
>    pres at gsnMaximize = True
>    gsn_panel(wks, plot, (/2, 1/), pres)
>
>
>
> end
>
>
>
> 2017-05-05 0:39 GMT+01:00 Mary Haley <haley at ucar.edu>:
>
>> I'm so sorry I didn't get back to you sooner on this.  I've been swamped
>> and then was out-of-town for a few days.
>>
>> It's a little hard to follow what's going on in this script without being
>> able to run it.
>>
>> I think the issue is that you are using the "res" variable in your call
>> to gsn_panel, and "res" has both gsnDraw and gsnFrame set to False. This
>> will cause *no* plot to appear when you call gsn_panel.  I'm surprised you
>> are getting any plot at all.
>>
>> I recommend creating a new resource variable for gsn_panel, because you
>> really don't want to use the one that you used for the plots themselves:
>>
>> pres = True
>> pres at gsnMaximize = True
>> gsn_panel(wks, plot, (/2, 1/), pres)
>>
>> If you
>> ​ continue to have problems, then it would help if you send me the two
>> data files. If you can't send the files, then send me the image you're
>> getting, and any error messages ​that the script produces.
>>
>> ​Thanks,
>>
>> --Mary
>>>>
>>
>>
>> On Tue, Apr 25, 2017 at 12:56 AM, Appo derbetini <appopson4 at gmail.com>
>> wrote:
>>
>>> Good morning Mary,
>>> Sorry , i forgot to add ncl-users email to the reply.
>>>
>>> I am attaching the modified version of the script were I'm trying to
>>> create a panel plot for two transects at latitude sellat and sellat2.
>>> Unfortunately, only one plot appears.
>>> Really, I don't know what is wrong this script.
>>>
>>> Your help will be appreciated.
>>>
>>>
>>>
>>> begin
>>>
>>>     in1  = addfile("uvwr.mean.JJAS.nc", "r")                         ;
>>> open netcdf file
>>>     input4  = addfile("atanas_topo.nc", "r")
>>>
>>>    u00     =  in1->u
>>>    v00     =  in1->v
>>>    w00     =  in1->w
>>>    rhum00  =  in1->r
>>>
>>>    topo1 = input4->HT
>>>
>>>
>>>    lat = in1->latitude
>>>    lon = in1->longitude
>>>
>>>
>>>    dv00    = uv2dvF (u00, v00)
>>>    uvd00   = dv2uvF (dv00)
>>>    ucomp_div00 = uvd00(0, :, : , :, :)
>>>    copy_VarCoords(u00, ucomp_div00)
>>>
>>>
>>>    sellat = 3.0
>>>    sellat2 = -2.5
>>>    latmin = -5.0
>>>    latmax = 3.0
>>>    lonmin = min(lon)
>>>    lonmax = 20.0
>>>
>>>    udiv00 = dim_avg_n_Wrap(ucomp_div00(:, :, {latmin:latmax},
>>> {lonmin:lonmax}), 2)
>>>    omega00 = dim_avg_n_Wrap(w00(:, :, {latmin:latmax}, {lonmin:lonmax}),
>>> 2)
>>>    r00 = dim_avg_n_Wrap(rhum00(:, :, {latmin:latmax},{lonmin:lonmax}),
>>> 2)
>>>
>>>    ;udiv00  = ucomp_div00(:, :, {sellat}, {lonmin:lonmax})
>>>    ;omega00 = w00(:, :, {sellat}, {lonmin:lonmax})
>>>    ;r00     = rhum00(:, :, {sellat},{lonmin:lonmax})
>>>
>>>
>>>    elev = 1013.25*(1 - topo1*0.0065/288.15)^5.25145
>>>
>>>    copy_VarCoords(topo1, elev)
>>>
>>> ;  topo = dim_avg_n_Wrap(elev({latmin:latmax},{lonmin:lonmax}), 0)
>>>
>>>    topo  = elev({sellat}, {lonmin:lonmax})
>>>    topo2 = elev({sellat2}, {lonmin:lonmax})
>>>
>>> ;---create plot
>>>    wks   = gsn_open_wks ("eps", "Overlay_Rhum_udiv_topo")        ; open
>>> workstation
>>>     plot = new(2, graphic)
>>>    res                     = True                  ; plot mods desired
>>>    res at gsnDraw             = False    ; turn off draw
>>>    res at gsnFrame            = False    ; turn off frame
>>>    res at gsnMaximize         = True
>>>
>>>    res at gsnRightString      = ""
>>>    res at gsnLeftString       = ""
>>>
>>>    cnres = res
>>>    cnres at cnFillOn          = True               ; turn on color fill
>>>    cnres at cnFillPalette     = "MPL_rainbow"
>>>    cnres at lbOrientation     = "Vertical"
>>>    cnres at pmLabelBarOrthogonalPosF = 0.08
>>>    cnres at tiYAxisString      = "Pressure (hPa)"
>>>
>>>    vcres = res
>>>    vcres at vcRefMagnitudeF = 10.0                ; define vector ref mag
>>>    vcres at vcRefLengthF    = 0.1              ; define length of vec ref
>>>    vcres at vcGlyphStyle    = "CurlyVector"      ; turn on curley vectors
>>>    vcres at vcMinDistanceF  = 0.04               ; thin out vectors
>>>    vcres at vcMapDirection  = False
>>>    vcres at vcLineArrowThicknessF   = 3.0
>>>    vcres at vcVectorDrawOrder  = "Draw"        ; draw vectors last
>>>
>>>    contour_plot_1 = gsn_csm_pres_hgt(wks, r00(0, :, {lonmin:lonmax}),
>>> cnres)
>>>    vector_plot_1  = gsn_csm_vector(wks, udiv00(0, :, {lonmin:lonmax}),
>>> -200.0*omega00(0, :, {lonmin:lonmax}), vcres )
>>>
>>>    contour_plot_2 = gsn_csm_pres_hgt(wks, r00(0, :, {lonmin:lonmax}),
>>> cnres)
>>>    vector_plot_2  = gsn_csm_vector(wks, udiv00(0, :, {lonmin:lonmax}),
>>> -200.0*omega00(0, :, {lonmin:lonmax}), vcres )
>>>
>>>
>>> ;---Add topo field using a filled polygon.
>>>    getvalues contour_plot_1
>>>       "trYMinF" : ymin
>>>       "trYMaxF" : ymax
>>>    end getvalues
>>>
>>>
>>> ;---Add topo field using a filled polygon.
>>>    getvalues contour_plot_2
>>>       "trYMinF" : ymin2
>>>       "trYMaxF" : ymax2
>>>    end getvalues
>>>
>>>
>>>
>>> ;---Create new X,Y arrays that form a closed polygon.
>>>    nlon  = dimsizes(topo&lon)
>>>    xtopo = new(nlon+3,typeof(topo&lon))
>>>    ytopo = new(nlon+3,typeof(topo))
>>>
>>>    xtopo(0:nlon-1) = topo&lon
>>>    ytopo(0:nlon-1) = topo
>>>    xtopo(nlon)     = topo&lon(nlon-1)
>>>    ytopo(nlon)     = ymax          ; Use actual Y max of contour plot
>>>    xtopo(nlon+1)   = topo&lon(0)
>>>    ytopo(nlon+1)   = ymax
>>>    xtopo(nlon+2)   = topo&lon(0)   ; This last point closes
>>>    ytopo(nlon+2)   = topo(0)       ; the polygon.
>>>
>>>
>>>    nlon2  = dimsizes(topo2&lon)
>>>    xtopo2 = new(nlon2+3,typeof(topo2&lon))
>>>    ytopo2 = new(nlon2+3,typeof(topo2))
>>>
>>>    xtopo2(0:nlon-1) = topo2&lon
>>>    ytopo2(0:nlon-1) = topo2
>>>    xtopo2(nlon)     = topo2&lon(nlon-1)
>>>    ytopo2(nlon)     = ymax2         ; Use actual Y max of contour plot
>>>    xtopo2(nlon+1)   = topo2&lon(0)
>>>    ytopo2(nlon+1)   = ymax2
>>>    xtopo2(nlon+2)   = topo2&lon(0)   ; This last point closes
>>>    ytopo2(nlon+2)   = topo2(0)       ; the polygon.
>>>
>>>
>>> ;---Add the polygon to the contour plot.
>>>    gnres = True
>>>    gnres at gsFillColor = "gray20"
>>>
>>>    id = gsn_add_polygon(wks, contour_plot_1, xtopo, ytopo, gnres)
>>>    id2 = gsn_add_polygon(wks, contour_plot_2, xtopo2, ytopo2, gnres)
>>>    overlay(contour_plot_1, vector_plot_1)
>>>    overlay(contour_plot_2, vector_plot_2)
>>>
>>>
>>>      draw(contour_plot_1)     ; This draws everything
>>>      draw(contour_plot_2)     ; This draws everything
>>>
>>>    plot(0) = contour_plot_1
>>>
>>>    plot(1) = contour_plot_2
>>>
>>>     gsn_panel(wks,plot, (/2, 1/), res)               ; now draw as one
>>> plot
>>>
>>>
>>> end
>>>
>>>
>>> 2017-04-24 21:58 GMT+01:00 Mary Haley <haley at ucar.edu>:
>>>
>>>> Hi,
>>>>
>>>> In the future, please email ncl-talk with follow-up questions so
>>>> everybody can benefit from the answers.
>>>>
>>>> What is the exact error message you're getting?
>>>>
>>>> --Mary
>>>>
>>>>
>>>> On Mon, Apr 24, 2017 at 4:08 AM, Appo derbetini <appopson4 at gmail.com>
>>>> wrote:
>>>>
>>>>> Hi Mary,
>>>>>
>>>>> Since you send me the script that plot masked fields, I am working to
>>>>> do it for many latitude in a panel plot.
>>>>> Always i got many errors saying that i cannot overlay.
>>>>> How to solve it?
>>>>> Thank you in advance.
>>>>>
>>>>> Appo
>>>>>
>>>>> 2017-04-10 8:23 GMT+01:00 Appo derbetini <appopson4 at gmail.com>:
>>>>>
>>>>>> Thank you very much.
>>>>>> It's look like what I'm trying to do.
>>>>>> regards
>>>>>>
>>>>>> 2017-04-09 22:11 GMT+01:00 Mary Haley <haley at ucar.edu>:
>>>>>>
>>>>>>> I'm sorry I didn't get back to this sooner, but I was having trouble
>>>>>>> trying to understand what you are doing in the script.
>>>>>>>
>>>>>>> First, I'm not sure it makes sense to take an average of
>>>>>>> topographical data across a set of latitudes and plot that as a topo layer.
>>>>>>>
>>>>>>> But, that aside, I think you just want to plot the topo line as a
>>>>>>> filled polygon instead of a separate XY plot.
>>>>>>>
>>>>>>> I'm not sure what I've attached is correct, but hopefully it gives
>>>>>>> you an idea of what to do.
>>>>>>>
>>>>>>> Note that I'm using "sellat" to select a single latitude value,
>>>>>>> rather than averaging across a range of latitudes.
>>>>>>>
>>>>>>> --Mary
>>>>>>>
>>>>>>>
>>>>>>> On Sun, Apr 9, 2017 at 2:10 AM, Appo derbetini <appopson4 at gmail.com>
>>>>>>> wrote:
>>>>>>>
>>>>>>>> Dear Mary,
>>>>>>>>
>>>>>>>> Despite applying what you suggest the problem remain the same
>>>>>>>> Thank you very much
>>>>>>>>
>>>>>>>> 2017-04-06 8:53 GMT+01:00 Appo derbetini <appopson4 at gmail.com>:
>>>>>>>>
>>>>>>>>> Dear Mary,
>>>>>>>>>
>>>>>>>>> Datasets are uploaded as indicated.
>>>>>>>>>
>>>>>>>>> What I really trying to do is something like Mask Example 14.
>>>>>>>>>
>>>>>>>>> But I don't want to interpolate topography to grid of wind.
>>>>>>>>>
>>>>>>>>> Here attached figure produced by the script
>>>>>>>>>
>>>>>>>>> Thank you
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> 2017-04-05 15:36 GMT+01:00 Mary Haley <haley at ucar.edu>:
>>>>>>>>>
>>>>>>>>>> I don't think I've ever seen this error before.
>>>>>>>>>>
>>>>>>>>>> It's hard to debug this particular script without having the data
>>>>>>>>>> so I can run it. But, I'm wondering if this line might be part of the
>>>>>>>>>> problem:
>>>>>>>>>>
>>>>>>>>>>   plot  = gsn_csm_pres_hgt_vector(wks, r00(0, :,
>>>>>>>>>> {lonmin:lonmax}), udiv00(0, :, {lonmin:lonmax}), -200.0*omega00(0, :,
>>>>>>>>>> {lonmin:lonmax}), res )
>>>>>>>>>>
>>>>>>>>>> When you do an arithmetic operation on an array, like:
>>>>>>>>>>
>>>>>>>>>>   -200. * omega00
>>>>>>>>>>
>>>>>>>>>> you end up stripping all of the metadata off omega00.  What's
>>>>>>>>>> getting passed into gsn_csm_pres_hgt has no metadata. This might be okay,
>>>>>>>>>> however, because I think this routine might just use the metadata from
>>>>>>>>>> "udiv00" instead.
>>>>>>>>>>
>>>>>>>>>> Still, it's worth trying this:
>>>>>>>>>>
>>>>>>>>>>  omega00_scale = omega00(0, :, {lonmin:lonmax})   ; subsets
>>>>>>>>>> omega00 *and* copies all metadata
>>>>>>>>>>  omega00_scale = -200.0 * omega00_scale
>>>>>>>>>>
>>>>>>>>>>   plot  = gsn_csm_pres_hgt_vector(wks, r00(0, :,
>>>>>>>>>> {lonmin:lonmax}), udiv00(0, :, {lonmin:lonmax}), omega00_scale, res )
>>>>>>>>>>
>>>>>>>>>> If you continue to have problems with this script, could you
>>>>>>>>>> upload your data (if it's not too large) to our ftp:
>>>>>>>>>>
>>>>>>>>>> ftp ftp.cgd.ucar.edu
>>>>>>>>>> anonymous
>>>>>>>>>> <use your email address for the password>
>>>>>>>>>> cd incoming
>>>>>>>>>> put uvwr.mean.JJAS.nc
>>>>>>>>>> put atanas_topo.nc
>>>>>>>>>> quit
>>>>>>>>>>
>>>>>>>>>> Thanks,
>>>>>>>>>>
>>>>>>>>>> --Mary
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> On Tue, Apr 4, 2017 at 9:18 AM, Appo derbetini <
>>>>>>>>>> appopson4 at gmail.com> wrote:
>>>>>>>>>>
>>>>>>>>>>> Dear all,
>>>>>>>>>>> I am trying to overlay topography profil on a profil of wind and
>>>>>>>>>>> relative humidity.
>>>>>>>>>>> I don't want to interpolate topography to show some details.
>>>>>>>>>>>
>>>>>>>>>>> Atached are script and datasets used for the plot.
>>>>>>>>>>>
>>>>>>>>>>> But I having errors.
>>>>>>>>>>>
>>>>>>>>>>> warning:LlDataPolygon: point 6.150009,0.000000 outside data
>>>>>>>>>>> domain
>>>>>>>>>>> warning:LlDataPolygon: point 19.983337,0.000000 outside data
>>>>>>>>>>> domain
>>>>>>>>>>> warning:LlDataLineTo: point 6.150009,0.000000 outside data domain
>>>>>>>>>>> warning:LlDataLineTo: point 19.983337,0.000000 outside data
>>>>>>>>>>> domain
>>>>>>>>>>>
>>>>>>>>>>> Any help will be appreciated.
>>>>>>>>>>>
>>>>>>>>>>> Best regards.
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> ncl-talk mailing list
>>>>>>>>>>> ncl-talk at ucar.edu
>>>>>>>>>>> List instructions, subscriber options, unsubscribe:
>>>>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>
>>>>>
>>>>
>>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20170507/173186d4/attachment.html 


More information about the ncl-talk mailing list