Skip to content

Add error message for >4D field input#2605

Open
VeckoTheGecko wants to merge 1 commit intoParcels-code:mainfrom
VeckoTheGecko:push-stllrnoswvuw
Open

Add error message for >4D field input#2605
VeckoTheGecko wants to merge 1 commit intoParcels-code:mainfrom
VeckoTheGecko:push-stllrnoswvuw

Conversation

@VeckoTheGecko
Copy link
Copy Markdown
Contributor

Description

I have added a check in assert_all_field_dims_have_axis to assert that field input is <=4D.

I have manually tested this with the following:

code

import numpy as np
import parcels
from parcels.interpolators import XLinear

import intake
runtime = np.timedelta64(2, "D")
dt = np.timedelta64(15, "m")


def _load_ds(chunk):
    """Helper function to load xarray dataset from catalog with or without chunking"""
    cat = intake.open_catalog(
        "/Users/Hodgs004/coding/repos/parcels-benchmarks/data/surf-data/parcels-benchmarks/catalog.yml"
    )
    chunks = {"time_counter": 1, "depth": 2, "y": chunk, "x": chunk} if chunk else None

    ds_u = (
        cat.moi_u(chunks=chunks).to_dask()[["vozocrtx"]].rename_vars({"vozocrtx": "U"})
    )
    ds_v = (
        cat.moi_v(chunks=chunks).to_dask()[["vomecrty"]].rename_vars({"vomecrty": "V"})
    )
    da_depth = cat.moi_w(chunks=chunks).to_dask()["depthw"]
    ds_mesh = cat.moi_mesh(chunks=None).read()[["glamf", "gphif"]].isel(t=0)
    ds_mesh["depthw"] = da_depth
    ds = parcels.convert.nemo_to_sgrid(fields=dict(U=ds_u, V=ds_v), coords=ds_mesh)

    return ds

if __name__=="__main__":
    ds = _load_ds(256)
    fieldset = parcels.FieldSet.from_sgrid_conventions(ds)

Error

/Users/Hodgs004/coding/repos/parcels/t.py:2: UserWarning: This is an alpha version of Parcels v4. The API is not stable and may change without deprecation warnings.
  import parcels
/Users/Hodgs004/coding/repos/parcels/src/parcels/convert.py:126: UserWarning: No depth dimension found in your dataset. Assuming no depth (i.e., surface data).
  warnings.warn("No depth dimension found in your dataset. Assuming no depth (i.e., surface data).", stacklevel=1)
Traceback (most recent call last):
  File "/Users/Hodgs004/coding/repos/parcels/t.py", line 32, in <module>
    fieldset = parcels.FieldSet.from_sgrid_conventions(ds)
  File "/Users/Hodgs004/coding/repos/parcels/src/parcels/_core/fieldset.py", line 311, in from_sgrid_conventions
    fields["U"] = Field("U", ds["U"], grid, XLinear)
                  ~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/Hodgs004/coding/repos/parcels/src/parcels/_core/field.py", line 106, in __init__
    assert_all_field_dims_have_axis(data, grid.xgcm_grid)
    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^
  File "/Users/Hodgs004/coding/repos/parcels/src/parcels/_core/xgrid.py", line 60, in assert_all_field_dims_have_axis
    raise ValueError(
    ...<2 lines>...
    )
ValueError: Two dimensions ('depth_center' and 'depth') provide values in the axis direction 'Z'. This is not possible, a field cannot have two dimensions on a single axis.

Checklist

"This is not possible, a field cannot have two dimensions on a single axis."
)
seen_axes[ax] = dim_name
assert len(ax_dims) <= 4
Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also raise a more informative error on this assert?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

Status: Backlog

Development

Successfully merging this pull request may close these issues.

Add checking that Fields have expected number of dimensions when calling an interpolation method

2 participants