Vectorfield#
import hvplot.xarray # noqa
vectorfield accepts 2d arrays of magnitude and angle on a grid and produces an array of vectors. x and y can be 2d or 1d coordinates.
import numpy as np
import xarray as xr
import holoviews as hv
import cartopy.crs as ccrs
def sample_data(shape=(20, 30)):
"""
Return ``(x, y, u, v, crs)`` of some vector data
computed mathematically. The returned crs will be a rotated
pole CRS, meaning that the vectors will be unevenly spaced in
regular PlateCarree space.
"""
crs = ccrs.RotatedPole(pole_longitude=177.5, pole_latitude=37.5)
x = np.linspace(311.9, 391.1, shape[1])
y = np.linspace(-23.6, 24.8, shape[0])
x2d, y2d = np.meshgrid(x, y)
u = 10 * (2 * np.cos(2 * np.deg2rad(x2d) + 3 * np.deg2rad(y2d + 30)) ** 2)
v = 20 * np.cos(6 * np.deg2rad(x2d))
return x, y, u, v, crs
xs, ys, U, V, crs = sample_data()
mag = np.sqrt(U**2 + V**2)
angle = (np.pi/2.) - np.arctan2(U/mag, V/mag)
ds = xr.Dataset({'mag': xr.DataArray(mag, dims=('y', 'x'), coords={'y': ys, 'x': xs}),
'angle': xr.DataArray(angle, dims=('y', 'x'), coords={'y': ys, 'x': xs})},
attrs={'crs': crs})
ds
<xarray.Dataset> Size: 10kB
Dimensions: (y: 20, x: 30)
Coordinates:
* y (y) float64 160B -23.6 -21.05 -18.51 -15.96 ... 19.71 22.25 24.8
* x (x) float64 240B 311.9 314.6 317.4 320.1 ... 385.6 388.4 391.1
Data variables:
mag (y, x) float64 5kB 6.459 2.149 5.899 11.25 ... 22.28 22.74 22.0
angle (y, x) float64 5kB 1.413 0.3677 -0.9793 ... -0.9368 -1.049 -1.127
Attributes:
crs: +proj=ob_tran +ellps=WGS84 +a=6378137.0 +o_proj=latlon +o_lon_p...ds.hvplot.vectorfield(x='x', y='y', angle='angle', mag='mag', hover=False).opts(magnitude='mag')
Geographic Data#
If a dataset has an attr called crs which is a cartopy object or a proj4 string, then just by setting the option geo=True will use the correct crs.
ds.hvplot.vectorfield(x='x', y='y', angle='angle', mag='mag',
hover=False, geo=True, tiles="CartoLight")
If you set coastline or features it will keep the original crs and transform the features to the data crs.
ds.hvplot.vectorfield(x='x', y='y', angle='angle', mag='mag',
hover=False, geo=True, coastline=True)
Large Data#
The visualization of vector fields from large datasets often presents a challenge. Direct plotting methods can quickly consume excessive memory, leading to crashes or unresponsive applications. To address this issue, we introduce a dynamic downsampling technique that enables interactive exploration of vector fields without sacrificing performance. We first create a sample_data that contains 4,000,000 points.
xs, ys, U, V, crs = sample_data(shape=(2000, 2000))
mag = np.sqrt(U**2 + V**2)
angle = (np.pi/2.) - np.arctan2(U/mag, V/mag)
ds = xr.Dataset({'mag': xr.DataArray(mag, dims=('y', 'x'), coords={'y': ys, 'x': xs}),
'angle': xr.DataArray(angle, dims=('y', 'x'), coords={'y': ys, 'x': xs})},
attrs={'crs': crs})
ds
<xarray.Dataset> Size: 64MB
Dimensions: (y: 2000, x: 2000)
Coordinates:
* y (y) float64 16kB -23.6 -23.58 -23.55 -23.53 ... 24.75 24.78 24.8
* x (x) float64 16kB 311.9 311.9 312.0 312.0 ... 391.0 391.1 391.1
Data variables:
mag (y, x) float64 32MB 6.459 6.383 6.307 6.232 ... 22.04 22.02 22.0
angle (y, x) float64 32MB 1.413 1.41 1.406 1.402 ... -1.125 -1.126 -1.127
Attributes:
crs: +proj=ob_tran +ellps=WGS84 +a=6378137.0 +o_proj=latlon +o_lon_p...If we just try to call ds.hvplot.vectorfield this is probably returning MemoryError. The alternative is to dynamically downsample the data based on the visible range. This helps manage memory consumption when dealing with large datasets, especially when plotting vector fields. We are going to use HoloViews to create a view (hv.DynamicMap) whose content is dynamically updated based on the data range displayed (tracked by the hv.streams.RangeXY stream), find out more about these concepts.
def downsample_quiver(x_range=None, y_range=None, nmax=10):
"""
Creates a HoloViews vector field plot from a dataset, dynamically downsampling
data based on the visible range to optimize memory usage.
Args:
x_range (tuple, optional): Range of x values to include. Defaults to None (full range).
y_range (tuple, optional): Range of y values to include. Defaults to None (full range).
nmax (int, optional): Maximum number of points along each axis after coarsening.
Defaults to 10.
Returns:
HoloViews DynamicMap: A dynamic vector field plot that updates based on the visible range.
"""
if x_range is None or y_range is None:
# No range provided, downsample the entire dataset for initial display
xs, ys = ds.x.size, ds.y.size # Get dataset dimensions
ix, iy = xs // nmax, ys // nmax # Calculate downsampling intervals
ix = max(1, ix) # Ensure interval is at least 1
iy = max(1, iy)
sub = ds.coarsen(x=ix, y=iy, side="center", boundary="trim").mean() # Downsample
else:
# Select data within the specified range
sub = ds.sel(x=slice(*x_range), y=slice(*y_range))
# Downsample the selected data
xs, ys = sub.x.size, sub.y.size
ix, iy = xs // nmax, ys // nmax
ix = max(1, ix)
iy = max(1, iy)
sub = sub.coarsen(x=ix, y=iy, side="center", boundary="trim").mean()
# Create the vector field plot
quiver = sub.hvplot.vectorfield(
x="x",
y="y",
mag="mag",
angle="angle",
hover=False,
).opts(magnitude="mag")
return quiver
# Create interactive plot components
range_xy = hv.streams.RangeXY() # Stream to capture range changes
filtered = hv.DynamicMap(downsample_quiver, streams=[range_xy]) # Dynamic plot
filtered # Display the plot