Demo of trying Vorticity_Divergence_Inversion.ipynb via MyBinder.org-served session¶
How this was generated in response discouse conversation about running a notebook.
(NOTE THAT BECAUSE I RAN INTO THINGS I DIDN'T UNDERSTAND I ADDED A PLOT USING A DEMO AT INPUT CELL #5 and stopped. Maybe the OP will have domain knowledge to get past this apparent hurdle that is probably due to different versions of the software being used. Ideally, the deveoper of the source repo should have listed versions in the repo with requiremetns.txt or environment.yml or in notebook with output of %pip list or %conda list, %pip freeze or something along those lines.)
Went to here and clicked 'launch binder'.
In the session that came up, I opened a terminal and ran git clone https://github.com/winash12/vortdivinversion.git to clone the related notebook and content.
In the file navigation panel on the left, I double clicked on the directory that was listed as vortdivinversion.
In the top of that notebook, I made a new cell and pasted in the following code that will install the necessary packages to run the first few cells and ran it.
In [1]:
%pip install cartopy
%pip install boto3
%pip install metpy
Requirement already satisfied: cartopy in /srv/conda/envs/notebook/lib/python3.10/site-packages (0.23.0) Requirement already satisfied: numpy>=1.21 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from cartopy) (1.26.4) Requirement already satisfied: matplotlib>=3.5 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from cartopy) (3.8.4) Requirement already satisfied: shapely>=1.7 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from cartopy) (2.0.4) Requirement already satisfied: packaging>=20 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from cartopy) (24.0) Requirement already satisfied: pyshp>=2.3 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from cartopy) (2.3.1) Requirement already satisfied: pyproj>=3.3.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from cartopy) (3.6.1) Requirement already satisfied: contourpy>=1.0.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (1.2.1) Requirement already satisfied: cycler>=0.10 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (0.12.1) Requirement already satisfied: fonttools>=4.22.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (4.51.0) Requirement already satisfied: kiwisolver>=1.3.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (1.4.5) Requirement already satisfied: pillow>=8 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (10.3.0) Requirement already satisfied: pyparsing>=2.3.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (3.1.2) Requirement already satisfied: python-dateutil>=2.7 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5->cartopy) (2.8.2) Requirement already satisfied: certifi in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pyproj>=3.3.1->cartopy) (2024.2.2) Requirement already satisfied: six>=1.5 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from python-dateutil>=2.7->matplotlib>=3.5->cartopy) (1.16.0) Note: you may need to restart the kernel to use updated packages. Requirement already satisfied: boto3 in /srv/conda/envs/notebook/lib/python3.10/site-packages (1.34.96) Requirement already satisfied: botocore<1.35.0,>=1.34.96 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from boto3) (1.34.97) Requirement already satisfied: jmespath<2.0.0,>=0.7.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from boto3) (1.0.1) Requirement already satisfied: s3transfer<0.11.0,>=0.10.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from boto3) (0.10.1) Requirement already satisfied: python-dateutil<3.0.0,>=2.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from botocore<1.35.0,>=1.34.96->boto3) (2.8.2) Requirement already satisfied: urllib3!=2.2.0,<3,>=1.25.4 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from botocore<1.35.0,>=1.34.96->boto3) (2.0.6) Requirement already satisfied: six>=1.5 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from python-dateutil<3.0.0,>=2.1->botocore<1.35.0,>=1.34.96->boto3) (1.16.0) Note: you may need to restart the kernel to use updated packages. Requirement already satisfied: metpy in /srv/conda/envs/notebook/lib/python3.10/site-packages (1.6.2) Requirement already satisfied: matplotlib>=3.5.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (3.8.4) Requirement already satisfied: numpy>=1.20.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (1.26.4) Requirement already satisfied: pandas>=1.4.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (2.2.2) Requirement already satisfied: pint>=0.17 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (0.23) Requirement already satisfied: pooch>=1.2.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (1.8.1) Requirement already satisfied: pyproj>=3.0.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (3.6.1) Requirement already satisfied: scipy>=1.8.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (1.13.0) Requirement already satisfied: traitlets>=5.0.5 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (5.11.2) Requirement already satisfied: xarray>=0.21.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from metpy) (2024.3.0) Requirement already satisfied: contourpy>=1.0.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (1.2.1) Requirement already satisfied: cycler>=0.10 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (0.12.1) Requirement already satisfied: fonttools>=4.22.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (4.51.0) Requirement already satisfied: kiwisolver>=1.3.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (1.4.5) Requirement already satisfied: packaging>=20.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (24.0) Requirement already satisfied: pillow>=8 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (10.3.0) Requirement already satisfied: pyparsing>=2.3.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (3.1.2) Requirement already satisfied: python-dateutil>=2.7 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from matplotlib>=3.5.0->metpy) (2.8.2) Requirement already satisfied: pytz>=2020.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pandas>=1.4.0->metpy) (2023.3.post1) Requirement already satisfied: tzdata>=2022.7 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pandas>=1.4.0->metpy) (2024.1) Requirement already satisfied: typing-extensions in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pint>=0.17->metpy) (4.8.0) Requirement already satisfied: platformdirs>=2.5.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pooch>=1.2.0->metpy) (3.11.0) Requirement already satisfied: requests>=2.19.0 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pooch>=1.2.0->metpy) (2.31.0) Requirement already satisfied: certifi in /srv/conda/envs/notebook/lib/python3.10/site-packages (from pyproj>=3.0.0->metpy) (2024.2.2) Requirement already satisfied: six>=1.5 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from python-dateutil>=2.7->matplotlib>=3.5.0->metpy) (1.16.0) Requirement already satisfied: charset-normalizer<4,>=2 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from requests>=2.19.0->pooch>=1.2.0->metpy) (3.3.0) Requirement already satisfied: idna<4,>=2.5 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from requests>=2.19.0->pooch>=1.2.0->metpy) (3.4) Requirement already satisfied: urllib3<3,>=1.21.1 in /srv/conda/envs/notebook/lib/python3.10/site-packages (from requests>=2.19.0->pooch>=1.2.0->metpy) (2.0.6) Note: you may need to restart the kernel to use updated packages.
In [2]:
import os.path
import xarray as xr
import numpy as np
import boto3
import metpy.calc as mpcalc
from botocore import UNSIGNED
from botocore.config import Config
import matplotlib.pyplot as plt
from matplotlib.cm import get_cmap
#from matplotlib.colormaps import get_cmap
import matplotlib.ticker as mticker
import cartopy.crs as crs
from cartopy.feature import NaturalEarthFeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import sys
import time
if (not os.path.isfile('gfs.t12z.pgrb2.0p25.f000')):
client = boto3.client('s3', config=Config(signature_version=UNSIGNED))
client.download_file('noaa-gfs-bdp-pds', 'gfs.20230809/12/atmos/gfs.t12z.pgrb2.0p25.f000', 'gfs.t12z.pgrb2.0p25.f000')
In [3]:
u850 = xr.open_dataset('gfs.t12z.pgrb2.0p25.f000', engine='cfgrib',backend_kwargs={'filter_by_keys':{'typeOfLevel': 'isobaricInhPa', 'shortName': 'u', 'level': 850}})
u = u850.u
print(u.shape)
v850 = xr.open_dataset('gfs.t12z.pgrb2.0p25.f000', engine='cfgrib', backend_kwargs={'filter_by_keys':{'typeOfLevel': 'isobaricInhPa', 'shortName': 'v', 'level': 850}})
v = v850.v
(721, 1440)
Compute the 850 hPa relative vorticity.
In [4]:
vort850 = mpcalc.vorticity(u, v)
fig = plt.figure(figsize=(12,9), dpi=300.)
# Create a set of axes for the figure and set
# its map projection to that of the input data.
ax = plt.axes(projection=crs.PlateCarree())
/tmp/ipykernel_816/816657038.py:1: UserWarning: More than one time coordinate present for variable "u". vort850 = mpcalc.vorticity(u, v)
I don't know why that doesn't work to show any plot. (I do note that I see UserWarning: More than one time coordinate present for variable "u" and a lot time warnings are just warnings and something usually shows.
The vorticity dpcumentation here does say, "Changed in version 1.0: Changed signature from (u, v, dx, dy)". I cannot tell if maybe that is the issue and the notebook was using older version of metpy?
I note that the example here works fine and gives a plot, as can be seen by running the code below:
In [5]:
import matplotlib.pyplot as plt
import metpy.calc as mpcalc
from metpy.cbook import example_data
# load example data
ds = example_data()
# Calculate the vertical vorticity of the flow
vort = mpcalc.vorticity(ds.uwind, ds.vwind)
# start figure and set axis
fig, ax = plt.subplots(figsize=(5, 5))
# scale vorticity by 1e5 for plotting purposes
cf = ax.contourf(ds.lon, ds.lat, vort * 1e5, range(-80, 81, 1), cmap=plt.cm.PuOr_r)
plt.colorbar(cf, pad=0, aspect=50)
ax.barbs(ds.lon.values, ds.lat.values, ds.uwind, ds.vwind, color='black', length=5, alpha=0.5)
ax.set(xlim=(260, 270), ylim=(30, 40))
ax.set_title('Relative Vorticity Calculation')
plt.show()
Add country borders and coastlines.
In [ ]:
countries = NaturalEarthFeature(category="cultural", scale="50m",
facecolor="none",
name="admin_0_countries")
ax.add_feature(countries, linewidth=.5, edgecolor="black")
ax.coastlines('50m', linewidth=0.8)
In [ ]:
plot = vort850.plot(levels=np.arange(-1.e-4, 1.e-4, 0.2e-5), cmap=get_cmap('PRGn'), transform=crs.PlateCarree(), cbar_kwargs={'label':'relative vorticity (x$10^{-5} s^{-1}$)', 'shrink': 0.98})
Set the map's extent to cover just Hurricane Dora.
In [ ]:
ax.set_extent([-180.,-150.,0.,20.],crs=crs.PlateCarree())
Add latitude/longitude gridlines.
In [ ]:
gridlines = ax.gridlines(color="grey", linestyle="dotted", draw_labels=True)
gridlines.xlabels_top = False
gridlines.ylabels_right = False
gridlines.xlocator = mticker.FixedLocator(np.arange(-180.,149.,5.))
gridlines.ylocator = mticker.FixedLocator(np.arange(0.,21.,5.))
gridlines.xlabel_style = {'size':12, 'color':'black'}
gridlines.ylabel_style = {'size':12, 'color':'black'}
gridlines.xformatter = LONGITUDE_FORMATTER
gridlines.yformatter = LATITUDE_FORMATTER
Add a plot title, then show the image.
In [ ]:
plt.title("GFS 0-h 850 hPa relative vorticity (x$10^{-5} s^{-1}$) at 1200 UTC 9 August 2023")
plt.show()
Compute the 850 hPa divergence.
In [ ]:
div850 = mpcalc.divergence(u, v)
Create a figure instance.
In [ ]:
fig = plt.figure(figsize=(12,9), dpi=300.)
Create a set of axes for the figure and set
its map projection to that of the input data.
In [ ]:
ax = plt.axes(projection=crs.PlateCarree())
Add country borders and coastlines.
In [ ]:
countries = NaturalEarthFeature(category="cultural", scale="50m",
facecolor="none",
name="admin_0_countries")
ax.add_feature(countries, linewidth=.5, edgecolor="black")
ax.coastlines('50m', linewidth=0.8)
Plot the 850 hPa divergence using xarray's plot functionality.
In [ ]:
plot = div850.plot(levels=np.arange(-1.e-4, 1.e-4, 0.2e-5), cmap=get_cmap('PRGn'), transform=crs.PlateCarree(), cbar_kwargs={'label':'relative vorticity (x$10^{-5} s^{-1}$)', 'shrink': 0.98})
Set the map's extent to cover just Hurricane Dora.
In [ ]:
ax.set_extent([-180.,-150.,0.,20.],crs=crs.PlateCarree())
Add latitude/longitude gridlines.
In [ ]:
gridlines = ax.gridlines(color="grey", linestyle="dotted", draw_labels=True)
gridlines.xlabels_top = False
gridlines.ylabels_right = False
gridlines.xlocator = mticker.FixedLocator(np.arange(-180.,149.,5.))
gridlines.ylocator = mticker.FixedLocator(np.arange(0.,21.,5.))
gridlines.xlabel_style = {'size':12, 'color':'black'}
gridlines.ylabel_style = {'size':12, 'color':'black'}
gridlines.xformatter = LONGITUDE_FORMATTER
gridlines.yformatter = LATITUDE_FORMATTER
Add a plot title, then show the image.
In [ ]:
plt.title("GFS 0-h 850 hPa divergence (x$10^{-5} s^{-1}$) at 1200 UTC 9 August 2023")
plt.show()
In [ ]:
vortmask = mpcalc.bounding_box_mask(vort850,5.,13.5,191.,202.)
In [ ]:
divmask = mpcalc.bounding_box_mask(div850,5.,13.5,191.,202.)
In [ ]:
i_bb_indices = mpcalc.find_bounding_box_indices(vortmask,5.,13.5,191.,202.)
In [ ]:
o_bb_indices = mpcalc.find_bounding_box_indices(vortmask,0.,30,180.,220)
In [ ]:
dx, dy = mpcalc.lat_lon_grid_deltas(vortmask.longitude, vortmask.latitude)
In [ ]:
upsi,vpsi = mpcalc.rotational_wind_from_inversion(vortmask,dx,dy,o_bb_indices,i_bb_indices)
Create a figure instance.
In [ ]:
fig = plt.figure(figsize=(12,9), dpi=300.)
Create a set of axes for the figure and set
its map projection to that of the input data.
In [ ]:
ax = plt.axes(projection=crs.PlateCarree())
Add country borders and coastlines.
In [ ]:
countries = NaturalEarthFeature(category="cultural", scale="50m",
facecolor="none",
name="admin_0_countries")
ax.add_feature(countries, linewidth=.5, edgecolor="black")
ax.coastlines('50m', linewidth=0.8)
Compute the magnitude of the non-divergent component of the 850 hPa wind.
In [ ]:
nd_spd = np.sqrt(upsi**2 + vpsi**2)
Plot this using xarray's plot functionality.
In [ ]:
plot = nd_spd.plot(levels=np.arange(0., 13., 1.), cmap=get_cmap('YlGnBu'), transform=crs.PlateCarree(), cbar_kwargs={'label':'non-divergent wind ($m s^{-1}$)', 'shrink': 0.98})
Set the map's extent to match that over which we computed the non-divergent wind.
In [ ]:
ax.set_extent([-180.,-140.,0.,30.],crs=crs.PlateCarree())
Add latitude/longitude gridlines.
In [ ]:
gridlines = ax.gridlines(color="grey", linestyle="dotted", draw_labels=True)
gridlines.xlabels_top = False
gridlines.ylabels_right = False
gridlines.xlocator = mticker.FixedLocator(np.arange(-180.,139.,5.))
gridlines.ylocator = mticker.FixedLocator(np.arange(0.,31.,5.))
gridlines.xlabel_style = {'size':12, 'color':'black'}
gridlines.ylabel_style = {'size':12, 'color':'black'}
gridlines.xformatter = LONGITUDE_FORMATTER
gridlines.yformatter = LATITUDE_FORMATTER
Add a plot title, then show the image.
In [ ]:
plt.title("GFS 0-h 850 hPa non-divergent wind magnitude ($m s^{-1}$) due to Dora at 1200 UTC 9 August 2023")
plt.show()
In [ ]:
uchi,vchi = mpcalc.divergent_wind_from_inversion(divmask,dx,dy,o_bb_indices,i_bb_indices)
Create a set of axes for the figure and set
its map projection to that of the input data.
In [ ]:
ax = plt.axes(projection=crs.PlateCarree())
Add country borders and coastlines.
In [ ]:
countries = NaturalEarthFeature(category="cultural", scale="50m",
facecolor="none",
name="admin_0_countries")
ax.add_feature(countries, linewidth=.5, edgecolor="black")
ax.coastlines('50m', linewidth=0.8)
Compute the magnitude of the non-divergent component of the 850 hPa wind.
In [ ]:
nd_spd = np.sqrt(uchi**2 + vchi**2)
Plot this using xarray's plot functionality.
In [ ]:
plot = nd_spd.plot(levels=np.arange(0., 13., 1.), cmap=get_cmap('YlGnBu'), transform=crs.PlateCarree(), cbar_kwargs={'label':'non-divergent wind ($m s^{-1}$)', 'shrink': 0.98})
Set the map's extent to match that over which we computed the non-divergent wind.
In [ ]:
ax.set_extent([-180.,-140.,0.,30.],crs=crs.PlateCarree())
Add latitude/longitude gridlines.
In [ ]:
gridlines = ax.gridlines(color="grey", linestyle="dotted", draw_labels=True)
gridlines.top_labels = False
gridlines.right_labels = False
gridlines.xlocator = mticker.FixedLocator(np.arange(-180.,139.,5.))
gridlines.ylocator = mticker.FixedLocator(np.arange(0.,31.,5.))
gridlines.xlabel_style = {'size':12, 'color':'black'}
gridlines.ylabel_style = {'size':12, 'color':'black'}
gridlines.xformatter = LONGITUDE_FORMATTER
gridlines.yformatter = LATITUDE_FORMATTER
Add a plot title, then show the image.
In [ ]:
plt.title("GFS 0-h 850 hPa divergent wind magnitude ($m s^{-1}$) due to Dora at 1200 UTC 9 August 2023")
plt.savefig('vectorized_version')
plt.show()