Compatibility with HoloViz Tools#

This user guide section showcases how UXarray data structures (Grid, UxDataArray) can be converted in a Spatialpandas GeoDataFrame, which can be used directly by packages from the HoloViz stack, such as hvPlot, Datashader, Holoviews, and Geoviews. In addition to showing how to convert to a GeoDataFrame, a series of visualizations using hvPlot and GeoViews is showcased based around the converted data structure.

import uxarray as ux
import numpy as np
import warnings

import matplotlib
import matplotlib.pyplot as plt

import hvplot.pandas
import holoviews as hv

warnings.filterwarnings("ignore")

Data#

For this notebook, we will be using E3SM model output, which uses a cubed-sphere grid. However, this notebook can be adapted to datasets in any of our supported formats.

base_path = "../../test/meshfiles/ugrid/outCSne30/"
grid_path = base_path + "outCSne30.ug"
data_path = base_path + "outCSne30_vortex.nc"

uxds = ux.open_dataset(grid_path, data_path)

uxds
<xarray.UxDataset> Size: 43kB
Dimensions:  (n_face: 5400)
Dimensions without coordinates: n_face
Data variables:
    psi      (n_face) float64 43kB ...

Conversion to spatialpandas.GeoDataFrame for Visualization#

In order to support visualization with the popular HoloViz stack of libraries (hvPlot, HoloViews, Datashader, etc.), UXarray provides methods for converting Grid and UxDataArray objects into a SpatialPandas GeoDataFrame, which can be used for visualizing the nodes, edges, and polygons that represent each grid, in addition to data variables.

Grid Conversion#

If you wish to plot only the grid’s geometrical elements such as nodes and edges (i.e. without any data variables mapped to them), you can directly convert a Grid instance into a GeoDataFrame.

Each UxDataset and UxDataArray has a Grid instance paired to it, which can be accessed using the .uxgrid attribute.

You can use the Grid.to_geodataframe() method to construct a GeoDataFrame.

gdf_grid = uxds.uxgrid.to_geodataframe()
gdf_grid
geometry
0 Polygon([[-45.0, -35.26438903808594, -42.0, -3...
1 Polygon([[-42.0, -36.61769485473633, -39.0, -3...
2 Polygon([[-39.0, -37.852420806884766, -36.0, -...
3 Polygon([[-36.0, -38.973445892333984, -33.0, -...
4 Polygon([[-33.0, -39.98556900024414, -30.0, -4...
... ...
5333 Polygon([[147.3314971923828, 43.07368087768555...
5334 Polygon([[144.1993408203125, 42.01158142089844...
5335 Polygon([[141.0996856689453, 40.83760452270508...
5336 Polygon([[138.03317260742188, 39.5490798950195...
5337 Polygon([[135.0, 38.143314361572266, 132.0, 36...

5338 rows × 1 columns

UxDataArray & UxDataset Conversion#

If you are interested in mapping data to each face, you can index the UxDataset with the variable of interest (in this case “psi”) to return the same GeoDataFrame as above, but now with data mapped to each face.

Warning

UXarray currently only supports mapping data variables that are mapped to each face. Variables that reside on the corner nodes or edges are currently not supported for visualization.

gdf_data = uxds["psi"].to_geodataframe()
gdf_data
geometry psi
0 Polygon([[-45.0, -35.26438903808594, -42.0, -3... 1.351317
1 Polygon([[-42.0, -36.61769485473633, -39.0, -3... 1.330915
2 Polygon([[-39.0, -37.852420806884766, -36.0, -... 1.310140
3 Polygon([[-36.0, -38.973445892333984, -33.0, -... 1.289056
4 Polygon([[-33.0, -39.98556900024414, -30.0, -4... 1.267717
... ... ...
5333 Polygon([[147.3314971923828, 43.07368087768555... 0.733539
5334 Polygon([[144.1993408203125, 42.01158142089844... 0.712085
5335 Polygon([[141.0996856689453, 40.83760452270508... 0.690883
5336 Polygon([[138.03317260742188, 39.5490798950195... 0.669989
5337 Polygon([[135.0, 38.143314361572266, 132.0, 36... 0.649467

5338 rows × 2 columns

Challenges with Representing Geoscience Data as Geometries#

When we convert to a GeoDataFrame, we internally represent the surface of a sphere as a collection of polygons over a 2D projection. However, this leads to issues around the Antimeridian (180 degrees east or west) where polygons are incorrectly constructed and have incorrect geometries. When constructing the GeoDataFrame, UXarray detects and corrects any polygon that touches or crosses the antimeridian. An array of indices of these faces can be accessed as part of the Grid object.

Antimeridian Visual
uxds.uxgrid.antimeridian_face_indices
array([1814, 1844, 1874, 1904, 1934, 1964, 1994, 2024, 2054, 2084, 2114,
       2144, 2174, 2204, 2234, 2264, 2294, 2324, 2354, 2384, 2414, 2444,
       2474, 2504, 2534, 2564, 2594, 2624, 2654, 2684, 3615, 3645, 3675,
       3705, 3735, 3765, 3795, 3825, 3855, 3885, 3915, 3945, 3975, 4005,
       4034, 4035, 4964, 4965, 4995, 5025, 5055, 5085, 5115, 5145, 5175,
       5205, 5235, 5265, 5295, 5325, 5355, 5385])

Taking a look at one of these faces that crosses or touches the antimeridian, we can see that it’s split across the antimeridian and represented as a MultiPolygon, which allows us to properly render this two dimensional grid.

gdf_data.geometry[uxds.uxgrid.antimeridian_face_indices[0]]
Polygon([[-180.0, -45.0, -177.0, -44.96071243286133, -177.0, -41.96092987060547, -180.0, -42.0, -180.0, -45.0]])

For more details about the algorithm used for splitting these polygons, see the Antimeridian Python Package, which we use internally for correcting polygons on the antimeridian.

Visualizing Geometries#

The next sections will show how the GeoDataFrame representing our unstructured grid can be used with hvPlot and GeoViews to generate visualizations. Examples using both the matplotlib and bokeh backends are shown, with bokeh examples allowing for interactive plots.

Nodes#

# use matplotlib backend for rendering
hv.extension("matplotlib")

plot_kwargs = {
    "size": 6.0,
    "xlabel": "Longitude",
    "ylabel": "Latitude",
    "coastline": True,
    "width": 1600,
    "title": "Node Plot (Matplotlib Backend)",
}


gdf_grid.hvplot.points(**plot_kwargs)
# use bokeh backend for rendering
hv.extension("bokeh")

plot_kwargs = {
    "s": 1.0,
    "xlabel": "Longitude",
    "ylabel": "Latitude",
    "coastline": True,
    "frame_width": 700,
    "title": "Node Plot (Bokeh Backend)",
}

gdf_grid.hvplot.points(**plot_kwargs)