[ncl-talk] Wind barb orientation in python plots

[Zilore Mumba]

Ok, I will try. Thanks Dennis.

On Thu, Aug 26, 2021 at 6:01 PM Dennis Shea <shea at ucar.edu> wrote:

> I believe you will have to ask  'python support  group' not NCL.
>
> On Thu, Aug 26, 2021 at 2:06 AM Zilore Mumba via ncl-talk <
> ncl-talk at mailman.ucar.edu> wrote:
>
>> I am wondering if python does not recognise that wind barbs are directed
>> to the right in the southern hemisphere and to the left in the northern
>> hemisphere. This is the same problem in ncl which seems to be carried over
>> to python.
>> If it is possible to change this I would appreciate assistance. I
>> reproduce below a script to plot 850 winds, if it is possible to indicate
>> where to make any change.
>>
>> from netCDF4 import Dataset
>> import cartopy.io.shapereader as shpreader
>> import numpy as np
>> import matplotlib.pyplot as plt
>> from mpl_toolkits.axes_grid1 import make_axes_locatable
>> from matplotlib.cm import get_cmap
>> import cartopy.crs as crs
>> from cartopy.feature import NaturalEarthFeature
>>
>> from wrf import (getvar, interplevel, to_np, latlon_coords, get_cartopy,
>>  cartopy_xlim, cartopy_ylim)
>>
>> #Open the NetCDF file
>> #--------------------
>> ncfile = Dataset("wrfout_d01_2021-08-19_06:00:00")
>>
>> #Extract the pressure, geopotential height, and wind variables
>> #-------------------------------------------------------------
>> p = getvar(ncfile, "pressure")
>> z = getvar(ncfile, "z", units="dm")
>> ua = getvar(ncfile, "ua", units="kt")
>> va = getvar(ncfile, "va", units="kt")
>> wspd = getvar(ncfile, "wspd_wdir", units="kts")[0,:]
>>
>> #Interpolate geopotential height, u, and v winds to 850 hPa
>> #-----------------------------------------------------------
>> ht_850 = interplevel(z, p, 850)
>> u_850 = interplevel(ua, p, 850)
>> v_850 = interplevel(va, p, 850)
>> wspd_850 = interplevel(wspd, p, 850)
>>
>> #Get the lat/lon coordinates
>> #---------------------------
>> lats, lons = latlon_coords(ht_850)
>>
>> #Get the map projection information
>> #----------------------------------
>> cart_proj = get_cartopy(ht_850)
>>
>> #Create the figure
>> #-----------------
>> fig = plt.figure(figsize=(12,9))
>> ax = plt.axes(projection=cart_proj)
>>
>> #Download and add the states and coastlines
>> #-----------------------------------------
>> #states = NaturalEarthFeature(category="cultural", scale="50m",
>> # facecolor="none",
>> # name="admin_1_states_provinces_shp")
>> #ax.add_feature(states, linewidth=0.5, edgecolor="black")
>> #ax.coastlines('50m', linewidth=0.8)
>>
>> fname = './ZMB_ADM1.shp'
>> adm1_shapes = list(shpreader.Reader(fname).geometries())
>>
>> #Add the 850 hPa geopotential height contours
>> #--------------------------------------------
>> levels = np.arange(150., 155., 1.)
>> contours = plt.contour(to_np(lons), to_np(lats), to_np(ht_850),
>>    levels=levels, colors="blue",
>>    transform=crs.PlateCarree())
>> plt.clabel(contours, inline=1, fontsize=10, fmt="%i")
>>
>> #Add the wind speed contours
>> #---------------------------
>> levels = [0, 5, 10, 15, 20, 25, 30]
>> wspd_contours = plt.contourf(to_np(lons), to_np(lats), to_np(wspd_850),
>>  levels=levels,
>>  cmap=get_cmap("rainbow"),
>>  transform=crs.PlateCarree())
>> plt.colorbar(wspd_contours, ax=ax, orientation="vertical", pad=.05)
>>
>> #divider = make_axes_locatable(ax)
>> #ax_cb = divider.new_horizontal(size="5%", pad=0.1, axes_class=plt.Axes)
>>
>> #plt.add_axes(ax_cb)
>> #plt.colorbar(wspd_contours, cax=ax_cb, pad=.05)
>>
>> #Add the 850 hPa wind barbs, only plotting every 5th data point.
>> #-----------------------------------------------------------------
>> plt.barbs(to_np(lons[::10,::10]), to_np(lats[::10,::10]),
>>   to_np(u_850[::10, ::10]), to_np(v_850[::10, ::10]),
>>   transform=crs.PlateCarree(), length=6)
>>
>>
>> ax.add_geometries(adm1_shapes, crs.PlateCarree(),
>>                   edgecolor='black', facecolor='gray', alpha=0.5)
>>
>> #Set the map bounds
>> #------------------
>> ax.set_xlim(cartopy_xlim(ht_850))
>> ax.set_ylim(cartopy_ylim(ht_850))
>>
>> ax.gridlines()
>>
>> plt.title("850 MB Height (dm), Wind Speed (kt), Barbs (kt)")
>>
>> plt.show()
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at mailman.ucar.edu
>> List instructions, subscriber options, unsubscribe:
>> https://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20210826/471e8aba/attachment.html>