Subsetting#
When working with large unstructured grids, it’s often more efficient to focus on a smaller, region-specific subset rather than the entire global grid. UXarray offers a suite of grid-aware subsetting operations to help you isolate exactly the area you need for your analysis.
Nearest Neighbor
Bounding Box
Bounding Circle
Constant Latitude/Longitude
Interval Latitude/Longitude
import warnings
import cartopy.crs as ccrs
import geocat.datafiles as geodf
import geoviews.feature as gf
import holoviews as hv
import uxarray as ux
plot_opts = {"width": 700, "height": 350}
hv.extension("bokeh")
warnings.filterwarnings("ignore")
/home/docs/checkouts/readthedocs.org/user_builds/uxarray/conda/1454/lib/python3.14/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
from .autonotebook import tqdm as notebook_tqdm
Data#
In this example, we will be using the geocat-datafiles package to obtain our grid and data files. The dataset used in this example is a 30km global MPAS meshes. We will be investigating the relative humidity vertically interpolated to 200hPa (relhum200hPa) data variable.
datafiles = (
geodf.get(
"netcdf_files/MPAS/FalkoJudt/dyamond_1/30km/diag.2016-08-20_00.00.00_subset.nc"
),
geodf.get("netcdf_files/MPAS/FalkoJudt/dyamond_1/30km/x1.655362.grid_subset.nc"),
)
uxds = ux.open_dataset(datafiles[1], datafiles[0])
clim = (uxds["relhum_200hPa"][0].values.min(), uxds["relhum_200hPa"][0].values.max())
features = gf.coastline(
projection=ccrs.PlateCarree(), line_width=1, scale="50m"
) * gf.states(projection=ccrs.PlateCarree(), line_width=1, scale="50m")
Global Grid#
Many unstructured grids, such as those from global climate models, span the entire surface of a sphere. UXarray supports working with these global grids, handling cases that arise with the spherical geometry of the earth (wrap around at the antimeridian, pole points, etc.)
uxds["relhum_200hPa"][0].plot(
rasterize=True, periodic_elements="exclude", title="Global Grid", **plot_opts
) * features
In addition to plotting global grids, we can perform analysis operations on the entire grid.
uxds["relhum_200hPa"][0].values.mean()
np.float32(46.81902)
Regional Subsets#
UXarray supports taking subsets of a grid, which allows us to select a region and perform analysis directly on that area, as opposed to the global grid.
There are currently three supported subsetting methods, both for the Grid and UxDataArray data structures.
uxds["relhum_200hPa"].subset
<uxarray.UxDataArray.subset>
Supported Methods:
* nearest_neighbor(center_coord, k, element)
* bounding_circle(center_coord, r, element)
* bounding_box(lon_bounds, lat_bounds)
* constant_latitude(lat, lon_range)
* constant_longitude(lon, lat_range)
* constant_latitude_interval(lats)
* constant_longitude_interval(lons)
uxds["relhum_200hPa"].uxgrid.subset
<uxarray.Grid.subset>
Supported Methods:
* nearest_neighbor(center_coord, k, element)
* bounding_circle(center_coord, r, element)
* bounding_box(lon_bounds, lat_bounds)
* constant_latitude(lat, lon_range)
* constant_longitude(lon, lat_range)
* constant_latitude_interval(lats)
* constant_longitude_interval(lons)
Bounding Box#
We can declare a bounding box centered about the Chicago area by specifying the minimum and maximum longitude and latitude bounds.
lon_bounds = (-87.6298 - 2, -87.6298 + 2)
lat_bounds = (41.8781 - 2, 41.8781 + 2)
bbox_subset_nodes = uxds["relhum_200hPa"][0].subset.bounding_box(
lon_bounds,
lat_bounds,
)
bbox_subset_nodes.plot(
rasterize=True,
periodic_elements="exclude",
clim=clim,
title="Bounding Box Subset (Corner Node Query)",
**plot_opts,
) * features
Bounding Circle#
A bounding circle is defined using a center coordinate (lon, lat) and a radius (in degrees). The resulting subset will contain all elements within the radius of that circle.
center_coord = [-87.6298, 41.8781]
r = 2
bcircle_subset = uxds["relhum_200hPa"][0].subset.bounding_circle(center_coord, r)
bcircle_subset.plot(
rasterize=True,
periodic_elements="exclude",
clim=clim,
title="Bounding Circle Subset",
**plot_opts,
) * features
Nearest Neighbor#
Similar to the bounding circle, we can perform a nearest neighbor subset at some center coordinate (lon, lat) and query for some number of elements k
center_coord = [-87.6298, 41.8781]
nn_subset = uxds["relhum_200hPa"][0].subset.nearest_neighbor(
center_coord, k=30, element="nodes"
)
nn_subset.plot(
rasterize=True,
periodic_elements="exclude",
clim=clim,
title="Nearest Neighbor Subset (Query 30 closest nodes)",
**plot_opts,
) * features
We can increase the number of neighbors k to make the region larger.
nn_subset_120 = uxds["relhum_200hPa"][0].subset.nearest_neighbor(
center_coord, k=120, element="face centers"
)
nn_subset_120.plot(
rasterize=True,
periodic_elements="exclude",
clim=clim,
title="Nearest Neighbor Subset (Query 120 Closest Faces)",
**plot_opts,
) * features
When we set k=1, it selects the closest neighbor.
nn_subset_1 = uxds["relhum_200hPa"][0].subset.nearest_neighbor(
center_coord, k=1, element="face centers"
)
nn_subset_1.plot(
rasterize=True,
periodic_elements="exclude",
clim=clim,
title="Nearest Neighbor Subset (Query Closest Face)",
**plot_opts,
) * features
Latitude & Longitude Slices#
Constant Longitude#
lon = 0.0
clon_subset = uxds["relhum_200hPa"][0].subset.constant_longitude(lon)
clon_subset.plot(
rasterize=True,
clim=clim,
title="Constant Longitude Subset",
global_extent=True,
**plot_opts,
) * features
Constant Latitude#
lat = 0.0
clat_subset = uxds["relhum_200hPa"][0].subset.constant_latitude(lat)
clat_subset.plot(
rasterize=True,
clim=clim,
title="Constant Latitude Subset",
global_extent=True,
**plot_opts,
) * features
Longitude Interval#
lons = (-50, 50)
clon_int_subset = uxds["relhum_200hPa"][0].subset.constant_longitude_interval(lons)
clon_int_subset.plot(
rasterize=True,
clim=clim,
title="Longitude Interval Subset",
global_extent=True,
**plot_opts,
) * features
Latitude Interval#
lats = (-50, 50)
clat_int_subset = uxds["relhum_200hPa"][0].subset.constant_latitude_interval(lats)
clat_int_subset.plot(
rasterize=True,
clim=clim,
title="Latitude Interval Subset",
global_extent=True,
**plot_opts,
) * features
Determining if a Grid is a Subset#
To check if a Grid (or dataset using .uxgrid) is a subset, we can use Grid.is_subset, which will return either True or False, depending on whether the Grid is a subset. Since nn_subset_120 is a subset, using this feature we will return True:
nn_subset_120.uxgrid.is_subset
True
The file we have been using to create these subsets, uxds, is not a subset, so using the same call we will return False:
uxds.uxgrid.is_subset
False