Skip to content

Feature/from icon #2482

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
fluidnumerics-joe merged 3 commits into from
Jan 17, 2026
Merged

Feature/from icon #2482

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
1b9f442
Cherry pick changes from https://github.com/Parcels-code/Parcels/pull…
fluidnumerics-joe Jan 17, 2026
80165b1
Patch interpolation tutorial notebook
fluidnumerics-joe Jan 17, 2026
8de06c8
Update nestedgrids and stommel_uxarray examples to use new interpolat…
fluidnumerics-joe Jan 17, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
26 changes: 7 additions & 19 deletions docs/user_guide/examples/tutorial_interpolation.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -259,8 +259,10 @@
"source": [
"## Interpolators on unstructured grids\n",
"Parcels v4 supports the use of general circulation model output that is defined on unstructured grids. We include basic interpolators to help you get started, including\n",
"- `UxPiecewiseConstantFace` - this interpolator implements piecewise constant interpolation and is appropriate for data that is registered to the face centers of the unstructured grid\n",
"- `UxPiecewiseLinearNode` - this interpolator implements barycentric interpolation and is appropriate for data that is registered to the corner vertices of the unstructured grid faces\n",
"- `UxConstantFaceConstantZC` - this interpolator implements piecewise constant interpolation in both the lateral and vertical directions. It is appropriate for data that is registered to the face centers of the unstructured grid and centered on vertical layers.\n",
"- `UxLinearNodeConstantZC` - this interpolator implements barycentric interpolation in the lateral direction and piecewise constant in the vertical direction. It is appropriate for data that is registered to the corner vertices of the unstructured grid faces and centered on vertical layers.\n",
"- `UxLinearNodeLinearZF` - this interpolator implements barycentric interpolation in the lateral direction and piecewise linear interpolation in the vertical direction. It is appropriate for data that is registered to the corner vertices of the unstructured grid faces and on vertical layer interfaces.\n",
"- `UxConstantFaceLinearZF` - this interpolator implements piecewise constant interpolation in the lateral direction and piecewise linear interpolation in the vertical direction. It is appropriate for data that is registered to the face centers of the unstructured grid and centered on vertical layers\n",
"\n",
"To get started, we use a very simple generated `UxArray.UxDataset` that is included with Parcels."
]
Expand All @@ -281,7 +283,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, we create the `Field` and `Fieldset` objects that will be used later in advancing particles. When creating a `Field` or `VectorField` object for unstructured grid data, we attach a {py:obj}`parcels.UxGrid` object and attach an `interp_method` to each object. For data that is defined on face centers, we use the `UxPiecewiseConstantFace` interpolator and for data that is defined on the face vertices, we use the `UxPiecewiseLinearNode` interpolator. In this example, we will look specifically at interpolating a tracer field that is defined by the same underlying analytical function, but is defined on both faces and vertices as separate fields."
"Next, we create the`Fieldset` object that will be used later in advancing particles. When using the `FieldSet.from_ugrid_conventions()` method, Parcels takes care of automatically assigning the appropriate interpolator based on the coordinates of each data variable."
]
},
{
Expand All @@ -294,21 +296,7 @@
"\n",
"import parcels\n",
"\n",
"grid = parcels.UxGrid(ds.uxgrid, ds.nz, mesh=\"flat\")\n",
"\n",
"Tnode = parcels.Field(\n",
" \"T_node\",\n",
" ds[\"T_node\"],\n",
" grid,\n",
" interp_method=parcels.interpolators.UxPiecewiseLinearNode,\n",
")\n",
"Tface = parcels.Field(\n",
" \"T_face\",\n",
" ds[\"T_face\"],\n",
" grid,\n",
" interp_method=parcels.interpolators.UxPiecewiseConstantFace,\n",
")\n",
"fieldset = parcels.FieldSet([Tnode, Tface])"
"fieldset = parcels.FieldSet.from_ugrid_conventions(ds, mesh=\"flat\")"
]
},
{
Expand Down Expand Up @@ -341,7 +329,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"Now, we can create a `ParticleSet` for each interpolation experiment. In one of the `ParticleSet` objects we use the `SampleTracer_Node` kernel to showcase interpolating with the `UxPiecewiseLinearNode` interpolator for node-registered data. In the other, we use the `SampleTracer_Face` to showcase interpolating with the `UxPiecewiseConstantFace` interpolator for face-registered data."
"Now, we can create a `ParticleSet` for each interpolation experiment. In one of the `ParticleSet` objects we use the `SampleTracer_Node` kernel to showcase interpolating with the `UxLinearNodeLinearZF` interpolator for node-registered data. In the other, we use the `SampleTracer_Face` to showcase interpolating with the `UxConstantFaceConstantZC` interpolator for face-registered data."
]
},
{
Expand Down
10 changes: 5 additions & 5 deletions docs/user_guide/examples/tutorial_nestedgrids.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -312,7 +312,7 @@
" \"face_polygon\": (\n",
" (\n",
" \"time\",\n",
" \"nz\",\n",
" \"zf\",\n",
" \"n_face\",\n",
" ),\n",
" face_poly[np.newaxis, np.newaxis, :],\n",
Expand All @@ -325,7 +325,7 @@
" },\n",
" coords={\n",
" \"time\": np.array([np.timedelta64(0, \"ns\")]),\n",
" \"nz\": np.array([0]),\n",
" \"zf\": np.array([0]),\n",
" \"n_node\": np.arange(n_node),\n",
" \"n_face\": np.arange(n_face),\n",
" },\n",
Expand All @@ -338,8 +338,8 @@
" \"GridID\",\n",
" uxda,\n",
" # TODO note that here we need to use mesh=\"flat\" otherwise the hashing doesn't work. See https://github.com/Parcels-code/Parcels/pull/2439#discussion_r2627664010\n",
" parcels.UxGrid(uxda.uxgrid, z=uxda[\"nz\"], mesh=\"flat\"),\n",
" interp_method=parcels.interpolators.UxPiecewiseConstantFace,\n",
" parcels.UxGrid(uxda.uxgrid, z=uxda[\"zf\"], mesh=\"flat\"),\n",
" interp_method=parcels.interpolators.UxConstantFaceConstantZC,\n",
")\n",
"fieldset = parcels.FieldSet([GridID])"
]
Expand Down Expand Up @@ -541,7 +541,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "docs",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand Down
10 changes: 5 additions & 5 deletions docs/user_guide/examples_v3/tutorial_stommel_uxarray.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -112,7 +112,7 @@
"\n",
"In Parcels, grid searching is conducted with respect to the faces. In other words, when a grid index `ei` is provided to an interpolation method, this refers the face index `fi` at vertical layer `zi` (when unraveled). Within the interpolation method, the `field.grid.uxgrid.face_node_connectivity` attribute can be used to obtain the node indices that surround the face. Using these connectivity tables is necessary for properly indexing node registered data.\n",
"\n",
"For the example Stommel gyre dataset in this tutorial, the `u` and `v` velocity components are face registered (similar to FESOM). Parcels includes a nearest neighbor interpolator for face registered unstructured grid data through `Parcels.application_kernels.interpolation.UxPiecewiseConstantFace`. Below, we create the `Field`s `U` and `V` and associate them with the `UxGrid` we created previously and this interpolation method."
"For the example Stommel gyre dataset in this tutorial, the `u` and `v` velocity components are face registered (similar to FESOM). Parcels includes a nearest neighbor interpolator for face registered unstructured grid data through `Parcels.application_kernels.interpolation.UxConstantFaceConstantZC`. Below, we create the `Field`s `U` and `V` and associate them with the `UxGrid` we created previously and this interpolation method."
]
},
{
Expand All @@ -121,26 +121,26 @@
"metadata": {},
"outputs": [],
"source": [
"from parcels.application_kernels.interpolation import UxPiecewiseConstantFace\n",
"from parcels.application_kernels.interpolation import UxConstantFaceConstantZC\n",
"from parcels.field import Field\n",
"\n",
"U = Field(\n",
" name=\"U\",\n",
" data=ds.U,\n",
" grid=grid,\n",
" interp_method=UxPiecewiseConstantFace,\n",
" interp_method=UxConstantFaceConstantZC,\n",
")\n",
"V = Field(\n",
" name=\"V\",\n",
" data=ds.V,\n",
" grid=grid,\n",
" interp_method=UxPiecewiseConstantFace,\n",
" interp_method=UxConstantFaceConstantZC,\n",
")\n",
"P = Field(\n",
" name=\"P\",\n",
" data=ds.p,\n",
" grid=grid,\n",
" interp_method=UxPiecewiseConstantFace,\n",
" interp_method=UxConstantFaceConstantZC,\n",
")"
]
},
Expand Down
21 changes: 2 additions & 19 deletions src/parcels/_core/field.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -406,9 +406,9 @@ def _assert_valid_uxdataarray(data: ux.UxDataArray):
uxarray.UxDataArray object.
"""
# Validate dimensions
if not ("nz1" in data.dims or "nz" in data.dims):
if not ("zf" in data.dims or "zc" in data.dims):
raise ValueError(
"Field is missing a 'nz1' or 'nz' dimension in the field's metadata. "
"Field is missing a 'zf' or 'zc' dimension in the field's metadata. "
"This attribute is required for xarray.DataArray objects."
)

Expand All @@ -427,23 +427,6 @@ def _assert_valid_uxdataarray(data: ux.UxDataArray):
"This attribute is required for xarray.DataArray objects."
)

_assert_valid_uxgrid(data.uxgrid)


def _assert_valid_uxgrid(grid):
"""Verifies that all the required attributes are present in the uxarray.UxDataArray.UxGrid object."""
if "Conventions" not in grid.attrs.keys():
raise ValueError(
"Field is missing a 'Conventions' attribute in the field's metadata. "
"This attribute is required for uxarray.UxDataArray objects."
)
if grid.attrs["Conventions"] != "UGRID-1.0":
raise ValueError(
"Field has a 'Conventions' attribute that is not 'UGRID-1.0'. "
"This attribute is required for uxarray.UxDataArray objects."
"See https://ugrid-conventions.github.io/ugrid-conventions/ for more information."
)


def _assert_compatible_combination(data: xr.DataArray | ux.UxDataArray, grid: ux.Grid | XGrid):
if isinstance(data, ux.UxDataArray):
Expand Down
113 changes: 90 additions & 23 deletions src/parcels/_core/fieldset.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -24,8 +24,10 @@
from parcels.interpolators import (
CGrid_Velocity,
Ux_Velocity,
UxPiecewiseConstantFace,
UxPiecewiseLinearNode,
UxConstantFaceConstantZC,
UxConstantFaceLinearZF,
UxLinearNodeConstantZC,
UxLinearNodeLinearZF,
XConstantField,
XLinear,
XLinear_Velocity,
Expand Down Expand Up @@ -187,35 +189,39 @@ def gridset(self) -> list[BaseGrid]:
return grids

@classmethod
def from_fesom2(cls, ds: ux.UxDataset):
"""Create a FieldSet from a FESOM2 uxarray.UxDataset.
def from_ugrid_conventions(cls, ds: ux.UxDataset, mesh: str = "spherical"):
"""Create a FieldSet from a Parcels compliant uxarray.UxDataset.
The main requirements for a uxDataset are naming conventions for vertical grid dimensions & coordinates

zf - Name for coordinate and dimension for vertical positions at layer interfaces
zc - Name for coordinate and dimension for vertical positions at layer centers

Parameters
----------
ds : uxarray.UxDataset
uxarray.UxDataset as obtained from the uxarray package.
uxarray.UxDataset as obtained from the uxarray package but with appropriate named vertical dimensions

Returns
-------
FieldSet
FieldSet object containing the fields from the dataset that can be used for a Parcels simulation.
"""
ds = ds.copy()
ds_dims = list(ds.dims)
if not all(dim in ds_dims for dim in ["time", "nz", "nz1"]):
if not all(dim in ds_dims for dim in ["time", "zf", "zc"]):
raise ValueError(
f"Dataset missing one of the required dimensions 'time', 'nz', or 'nz1'. Found dimensions {ds_dims}"
f"Dataset missing one of the required dimensions 'time', 'zf', or 'zc' for uxDataset. Found dimensions {ds_dims}"
)
grid = UxGrid(ds.uxgrid, z=ds.coords["nz"], mesh="spherical")
ds = _discover_fesom2_U_and_V(ds)

grid = UxGrid(ds.uxgrid, z=ds.coords["zf"], mesh=mesh)
ds = _discover_ux_U_and_V(ds)

fields = {}
if "U" in ds.data_vars and "V" in ds.data_vars:
fields["U"] = Field("U", ds["U"], grid, _select_uxinterpolator(ds["U"]))
fields["V"] = Field("V", ds["V"], grid, _select_uxinterpolator(ds["U"]))
fields["V"] = Field("V", ds["V"], grid, _select_uxinterpolator(ds["V"]))

if "W" in ds.data_vars:
fields["W"] = Field("W", ds["W"], grid, _select_uxinterpolator(ds["U"]))
fields["W"] = Field("W", ds["W"], grid, _select_uxinterpolator(ds["W"]))
fields["UVW"] = VectorField(
"UVW", fields["U"], fields["V"], fields["W"], vector_interp_method=Ux_Velocity
)
Expand All @@ -227,6 +233,63 @@ def from_fesom2(cls, ds: ux.UxDataset):

return cls(list(fields.values()))

@classmethod
def from_fesom2(cls, ds: ux.UxDataset, mesh: str = "spherical"):
"""Create a FieldSet from a FESOM2 uxarray.UxDataset.

Parameters
----------
ds : uxarray.UxDataset
uxarray.UxDataset as obtained from the uxarray package.

Returns
-------
FieldSet
FieldSet object containing the fields from the dataset that can be used for a Parcels simulation.
"""
ds = ds.copy()
ds_dims = list(ds.dims)
if not all(dim in ds_dims for dim in ["time", "nz", "nz1"]):
raise ValueError(
f"Dataset missing one of the required dimensions 'time', 'nz', or 'nz1' for FESOM data. Found dimensions {ds_dims}"
)
ds = ds.rename(
{
"nz": "zf", # Vertical Interface
"nz1": "zc", # Vertical Center
}
).set_index(zf="zf", zc="zc")

return FieldSet.from_ugrid_conventions(ds, mesh=mesh)

@classmethod
def from_icon(cls, ds: ux.UxDataset, mesh: str = "spherical"):
"""Create a FieldSet from a ICON uxarray.UxDataset.

Parameters
----------
ds : uxarray.UxDataset
uxarray.UxDataset as obtained from the uxarray package.

Returns
-------
FieldSet
FieldSet object containing the fields from the dataset that can be used for a Parcels simulation.
"""
ds = ds.copy()
ds_dims = list(ds.dims)
if not all(dim in ds_dims for dim in ["time", "depth", "depth_2"]):
raise ValueError(
f"Dataset missing one of the required dimensions 'time', 'depth', or 'depth_2' for ICON data. Found dimensions {ds_dims}"
)
ds = ds.rename(
{
"depth_2": "zf", # Vertical Interface
"depth": "zc", # Vertical Center
}
).set_index(zf="zf", zc="zc")
return FieldSet.from_ugrid_conventions(ds, mesh=mesh)

@classmethod
def from_sgrid_conventions(
cls, ds: xr.Dataset, mesh: Mesh | None = None
Expand Down Expand Up @@ -393,13 +456,13 @@ def _format_calendar_error_message(field: Field, reference_datetime: TimeLike) -
}


def _discover_fesom2_U_and_V(ds: ux.UxDataset) -> ux.UxDataset:
def _discover_ux_U_and_V(ds: ux.UxDataset) -> ux.UxDataset:
# Common variable names for U and V found in UxDatasets
common_fesom_UV = [("unod", "vnod"), ("u", "v")]
common_fesom_W = ["w"]
common_ux_UV = [("unod", "vnod"), ("u", "v")]
common_ux_W = ["w"]

if "W" not in ds:
for common_W in common_fesom_W:
for common_W in common_ux_W:
if common_W in ds:
ds = _ds_rename_using_standard_names(ds, {common_W: "W"})
break
Expand All @@ -411,7 +474,7 @@ def _discover_fesom2_U_and_V(ds: ux.UxDataset) -> ux.UxDataset:
"Dataset has only one of the two variables 'U' and 'V'. Please rename the appropriate variable in your dataset to have both 'U' and 'V' for Parcels simulation."
)

for common_U, common_V in common_fesom_UV:
for common_U, common_V in common_ux_UV:
if common_U in ds:
if common_V not in ds:
raise ValueError(
Expand All @@ -438,24 +501,28 @@ def _discover_fesom2_U_and_V(ds: ux.UxDataset) -> ux.UxDataset:
def _select_uxinterpolator(da: ux.UxDataArray):
"""Selects the appropriate uxarray interpolator for a given uxarray UxDataArray"""
supported_uxinterp_mapping = {
# (nz1,n_face): face-center laterally, layer centers vertically — piecewise constant
"nz1,n_face": UxPiecewiseConstantFace,
# (nz,n_node): node/corner laterally, layer interfaces vertically — barycentric lateral & linear vertical
"nz,n_node": UxPiecewiseLinearNode,
# (zc,n_face): face-center laterally, layer centers vertically — piecewise constant
"zc,n_face": UxConstantFaceConstantZC,
# (zc,n_node): node/corner laterally, layer centers vertically — barycentric lateral & piecewise constant vertical
"zc,n_node": UxLinearNodeConstantZC,
# (zf,n_node): node/corner laterally, layer interfaces vertically — barycentric lateral & linear vertical
"zf,n_node": UxLinearNodeLinearZF,
# (zf,n_face): face-center laterally, layer interfaces vertically — piecewise constant lateral & linear vertical
"zf,n_face": UxConstantFaceLinearZF,
}
# Extract only spatial dimensions, neglecting time
da_spatial_dims = tuple(d for d in da.dims if d not in ("time",))
if len(da_spatial_dims) != 2:
raise ValueError(
"Fields on unstructured grids must have two spatial dimensions, one vertical (nz or nz1) and one lateral (n_face, n_edge, or n_node)"
"Fields on unstructured grids must have two spatial dimensions, one vertical (zf or zc) and one lateral (n_face, n_edge, or n_node)"
)

# Construct key (string) for mapping to interpolator
# Find vertical and lateral tokens
vdim = None
ldim = None
for d in da_spatial_dims:
if d in ("nz", "nz1"):
if d in ("zf", "zc"):
vdim = d
if d in ("n_face", "n_node"):
ldim = d
Expand Down
Loading
Loading