From 38fc2e85fcbc3288db1b3825d8adec3c9aa451ed Mon Sep 17 00:00:00 2001 From: Erik van Sebille Date: Wed, 22 Jul 2026 11:41:00 +0200 Subject: [PATCH 1/4] Add _assert_no_nan_in_field_data function --- src/parcels/_core/model.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/src/parcels/_core/model.py b/src/parcels/_core/model.py index 79c9c6460..af208b3f1 100644 --- a/src/parcels/_core/model.py +++ b/src/parcels/_core/model.py @@ -149,6 +149,7 @@ def __init__(self, data: xr.Dataset, mesh: Mesh, vector_field_components: ptypin def assert_valid_field_data(self, field_data: xr.DataArray) -> None: # assert_all_field_dims_have_axis(field_data, self.grid.xgcm_grid) #! TODO v4: These checks should be revisited + _assert_no_nan_in_field_data(field_data) _assert_has_time_coordinate(field_data) @property @@ -507,3 +508,11 @@ def _assert_has_time_coordinate(da: xr.DataArray) -> None: if "time" not in da.coords: raise ValueError("Field data is missing a 'time' coordinate.") return + + +def _assert_no_nan_in_field_data(field_data: xr.DataArray) -> None: + arr = field_data.isel(time=0) if "time" in field_data.dims else field_data + if xr.ufuncs.isnan(arr).any(): + raise ValueError( + f"Field data {field_data.name!r} contains NaN values. Please fill or remove NaN values before using this field in a Parcels simulation." + ) From 99fcd822116b379daa3fd70460b205ef1d3729d6 Mon Sep 17 00:00:00 2001 From: Erik van Sebille Date: Wed, 22 Jul 2026 11:41:20 +0200 Subject: [PATCH 2/4] Use fillna in convert functions --- src/parcels/convert.py | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/parcels/convert.py b/src/parcels/convert.py index c8cea9ce3..5dfdc75a5 100644 --- a/src/parcels/convert.py +++ b/src/parcels/convert.py @@ -339,6 +339,8 @@ def nemo_to_sgrid(*, fields: dict[str, xr.Dataset | xr.DataArray], coords: xr.Da ds = _maybe_remove_depth_from_lonlat(ds) ds = _set_axis_attrs(ds, _NEMO_AXIS_VARNAMES) + ds = ds.fillna(0) # set NaN values to 0 + expected_axes = set("XYZT") # TODO: Update after we have support for 2D spatial fields if missing_axes := (expected_axes - set(ds.cf.axes)): raise ValueError( @@ -541,6 +543,9 @@ def copernicusmarine_to_sgrid( ds.attrs.clear() # Clear global attributes from the merging ds = _maybe_rename_coords(ds, _COPERNICUS_MARINE_AXIS_VARNAMES) + + ds = ds.fillna(0) # set NaN values to 0 + if "W" in ds.data_vars: # Negate W to convert from up positive to down positive (as that's the direction of positive z) ds["W"].data *= -1 From 2e43cb77ed670b73ad52d15ee0b421221c88d164 Mon Sep 17 00:00:00 2001 From: Erik van Sebille Date: Fri, 24 Jul 2026 15:05:43 +0200 Subject: [PATCH 3/4] Adding fill_nan_to_zero and assert_valid_fields options --- src/parcels/_core/fieldset.py | 10 ++++++++- src/parcels/_core/model.py | 38 +++++++++++++++++++++++------------ tests/test_fieldset.py | 24 +++++++++++++++++++--- 3 files changed, 55 insertions(+), 17 deletions(-) diff --git a/src/parcels/_core/fieldset.py b/src/parcels/_core/fieldset.py index 5864b9288..3532ee960 100644 --- a/src/parcels/_core/fieldset.py +++ b/src/parcels/_core/fieldset.py @@ -289,6 +289,8 @@ def from_sgrid_conventions( ds: xr.Dataset, mesh: ptyping.TMesh | None = None, vector_fields: ptyping.VectorFields | NotSetType = NOTSET, + fill_nan_to_zero: bool = True, + assert_valid_fields: bool = True, ): # TODO: Update mesh to be discovered from the dataset metadata """Create a FieldSet from a dataset using SGRID convention metadata. @@ -307,6 +309,10 @@ def from_sgrid_conventions( Mapping of vector field names to tuples of component variable names in the dataset. For example, ``{"UV": ("U", "V"), "UVW": ("U", "V", "W")}``. If omitted (default), vector fields are auto-discovered from standard variable names (``U``/``V``/``W``). + fill_nan_to_zero : bool, optional + If True (default), fill NaN values in the dataset with zeros, which avoids Interpolation errors. + assert_valid_fields : bool, optional + If True (default), asserts that the Fields are valid, including checks for dimensions names. In some cases, it can be useful to be able to override this check, but thread with caution, as invalid Fields can lead to unexpected behavior during execution. Returns ------- @@ -321,7 +327,9 @@ def from_sgrid_conventions( See https://sgrid.github.io/sgrid/ for more information on the SGRID conventions. """ - model = StructuredModelData.from_sgrid_conventions(ds, mesh, vector_fields) + model = StructuredModelData.from_sgrid_conventions( + ds, mesh, vector_fields, fill_nan_to_zero=fill_nan_to_zero, assert_valid_fields=assert_valid_fields + ) return cls([model]) def describe(self, buf: IO | None = None) -> None: diff --git a/src/parcels/_core/model.py b/src/parcels/_core/model.py index 3db1e784e..66cad282b 100644 --- a/src/parcels/_core/model.py +++ b/src/parcels/_core/model.py @@ -132,11 +132,20 @@ def preprocess_sgrid_model_data(ds: xr.Dataset) -> xr.Dataset: class StructuredModelData(ModelData): - def __init__(self, data: xr.Dataset, mesh: ptyping.TMesh, vector_field_components: ptyping.VectorFields): + def __init__( + self, + data: xr.Dataset, + mesh: ptyping.TMesh, + vector_field_components: ptyping.VectorFields, + fill_nan_to_zero: bool = True, + assert_valid_fields: bool = True, + ): if not isinstance(data, xr.Dataset): raise ValueError(f"Expected `data` to be an xarray.Dataset . Got {type(data)}") data = preprocess_sgrid_model_data(data) + if fill_nan_to_zero: + data = data.fillna(0) grid = XGrid(data, mesh) self.data = data @@ -144,11 +153,11 @@ def __init__(self, data: xr.Dataset, mesh: ptyping.TMesh, vector_field_component self.vector_field_components = vector_field_components self.field_to_interpolator = {} self._fields: list[Field | VectorField] | None = None - self.assert_valid_model_data() + if assert_valid_fields: + self.assert_valid_model_data() def assert_valid_field_data(self, field_data: xr.DataArray) -> None: # assert_all_field_dims_have_axis(field_data, self.grid.xgcm_grid) #! TODO v4: These checks should be revisited - _assert_no_nan_in_field_data(field_data) _assert_has_time_coordinate(field_data) @property @@ -182,7 +191,12 @@ def construct_fields(self) -> list[Field | VectorField]: @classmethod def from_sgrid_conventions( - cls, ds: xr.Dataset, mesh: ptyping.TMesh | None, vector_fields: ptyping.VectorFields | NotSetType + cls, + ds: xr.Dataset, + mesh: ptyping.TMesh | None, + vector_fields: ptyping.VectorFields | NotSetType, + fill_nan_to_zero: bool = True, + assert_valid_fields: bool = True, ) -> Self: ds = ds.copy() if mesh is None: @@ -213,7 +227,13 @@ def from_sgrid_conventions( vector_fields = resolve_vector_fields(ds, vector_fields) assert_valid_vector_fields(ds, vector_fields) - model = cls(ds, mesh=mesh, vector_field_components=vector_fields) + model = cls( + ds, + mesh=mesh, + vector_field_components=vector_fields, + fill_nan_to_zero=fill_nan_to_zero, + assert_valid_fields=assert_valid_fields, + ) model._fields = model.construct_fields() for f in model._fields: if isinstance(f, Field): @@ -509,11 +529,3 @@ def _assert_has_time_coordinate(da: xr.DataArray) -> None: if "time" not in da.coords: raise ValueError("Field data is missing a 'time' coordinate.") return - - -def _assert_no_nan_in_field_data(field_data: xr.DataArray) -> None: - arr = field_data.isel(time=0) if "time" in field_data.dims else field_data - if xr.ufuncs.isnan(arr).any(): - raise ValueError( - f"Field data {field_data.name!r} contains NaN values. Please fill or remove NaN values before using this field in a Parcels simulation." - ) diff --git a/tests/test_fieldset.py b/tests/test_fieldset.py index 015d4110e..ca492b4d6 100644 --- a/tests/test_fieldset.py +++ b/tests/test_fieldset.py @@ -330,6 +330,24 @@ def test_fieldset_from_sgrid_conventions(ds_name): assert len(fieldset.fields) > 0 +def test_fieldset_keep_nan_values(): + ds = datasets_structured["ds_2d_left"] + ds["U_A_grid"][:] = np.nan + + fieldset = FieldSet.from_sgrid_conventions(ds, mesh="flat") + assert np.isfinite(fieldset.U_A_grid.data).all() + fieldset = FieldSet.from_sgrid_conventions(ds, mesh="flat", fill_nan_to_zero=False) + assert np.isnan(fieldset.U_A_grid.data).all() + + +def test_fieldset_assert_valid_fields(): + ds = datasets_structured["ds_2d_left"].drop("time") + + FieldSet.from_sgrid_conventions(ds, mesh="flat", assert_valid_fields=False) + with pytest.raises(ValueError): + FieldSet.from_sgrid_conventions(ds, mesh="flat") + + def test_fieldset_add_error_on_duplicate_fields(): """Test that adding FieldSets with overlapping field names raises a ValueError.""" ds1 = datasets_structured["ds_2d_left"][["U_A_grid", "V_A_grid", "grid"]].rename({"U_A_grid": "U", "V_A_grid": "V"}) @@ -490,9 +508,9 @@ def test_fieldset_describe_backends(tmp_path): expected = """\ | Name | Type | Grid number | Interp method / value | Parcels backend | |:-------|:------------|--------------:|:------------------------|:------------------| -| U | Field | 0 | XLinear(...) | Zarr | -| V | Field | 0 | XLinear(...) | Zarr | -| W | Field | 0 | XLinear(...) | Zarr | +| U | Field | 0 | XLinear(...) | NumPy | +| V | Field | 0 | XLinear(...) | NumPy | +| W | Field | 0 | XLinear(...) | NumPy | | UV | VectorField | 0 | CGrid_Velocity(...) | - | | UVW | VectorField | 0 | CGrid_Velocity(...) | - | From 2bb3d10cb8b73bb484ba640e9931f33d916085f9 Mon Sep 17 00:00:00 2001 From: Erik van Sebille Date: Fri, 24 Jul 2026 15:22:22 +0200 Subject: [PATCH 4/4] Changing to one skip_field_data_validation option --- src/parcels/_core/fieldset.py | 11 ++++------- src/parcels/_core/model.py | 14 +++++--------- src/parcels/convert.py | 4 ---- tests/test_fieldset.py | 12 ++---------- 4 files changed, 11 insertions(+), 30 deletions(-) diff --git a/src/parcels/_core/fieldset.py b/src/parcels/_core/fieldset.py index 3532ee960..f0e48577f 100644 --- a/src/parcels/_core/fieldset.py +++ b/src/parcels/_core/fieldset.py @@ -289,8 +289,7 @@ def from_sgrid_conventions( ds: xr.Dataset, mesh: ptyping.TMesh | None = None, vector_fields: ptyping.VectorFields | NotSetType = NOTSET, - fill_nan_to_zero: bool = True, - assert_valid_fields: bool = True, + skip_field_data_validation: bool = False, ): # TODO: Update mesh to be discovered from the dataset metadata """Create a FieldSet from a dataset using SGRID convention metadata. @@ -309,10 +308,8 @@ def from_sgrid_conventions( Mapping of vector field names to tuples of component variable names in the dataset. For example, ``{"UV": ("U", "V"), "UVW": ("U", "V", "W")}``. If omitted (default), vector fields are auto-discovered from standard variable names (``U``/``V``/``W``). - fill_nan_to_zero : bool, optional - If True (default), fill NaN values in the dataset with zeros, which avoids Interpolation errors. - assert_valid_fields : bool, optional - If True (default), asserts that the Fields are valid, including checks for dimensions names. In some cases, it can be useful to be able to override this check, but thread with caution, as invalid Fields can lead to unexpected behavior during execution. + skip_field_data_validation : bool, optional + If True, skip validation of the field data. This can be useful for performance reasons, but may lead to unexpected behavior if the field data is invalid. Default is False. Returns ------- @@ -328,7 +325,7 @@ def from_sgrid_conventions( See https://sgrid.github.io/sgrid/ for more information on the SGRID conventions. """ model = StructuredModelData.from_sgrid_conventions( - ds, mesh, vector_fields, fill_nan_to_zero=fill_nan_to_zero, assert_valid_fields=assert_valid_fields + ds, mesh, vector_fields, skip_field_data_validation=skip_field_data_validation ) return cls([model]) diff --git a/src/parcels/_core/model.py b/src/parcels/_core/model.py index 66cad282b..9f31711d3 100644 --- a/src/parcels/_core/model.py +++ b/src/parcels/_core/model.py @@ -137,14 +137,13 @@ def __init__( data: xr.Dataset, mesh: ptyping.TMesh, vector_field_components: ptyping.VectorFields, - fill_nan_to_zero: bool = True, - assert_valid_fields: bool = True, + skip_field_data_validation: bool = False, ): if not isinstance(data, xr.Dataset): raise ValueError(f"Expected `data` to be an xarray.Dataset . Got {type(data)}") data = preprocess_sgrid_model_data(data) - if fill_nan_to_zero: + if not skip_field_data_validation: data = data.fillna(0) grid = XGrid(data, mesh) @@ -153,8 +152,7 @@ def __init__( self.vector_field_components = vector_field_components self.field_to_interpolator = {} self._fields: list[Field | VectorField] | None = None - if assert_valid_fields: - self.assert_valid_model_data() + self.assert_valid_model_data() def assert_valid_field_data(self, field_data: xr.DataArray) -> None: # assert_all_field_dims_have_axis(field_data, self.grid.xgcm_grid) #! TODO v4: These checks should be revisited @@ -195,8 +193,7 @@ def from_sgrid_conventions( ds: xr.Dataset, mesh: ptyping.TMesh | None, vector_fields: ptyping.VectorFields | NotSetType, - fill_nan_to_zero: bool = True, - assert_valid_fields: bool = True, + skip_field_data_validation: bool = False, ) -> Self: ds = ds.copy() if mesh is None: @@ -231,8 +228,7 @@ def from_sgrid_conventions( ds, mesh=mesh, vector_field_components=vector_fields, - fill_nan_to_zero=fill_nan_to_zero, - assert_valid_fields=assert_valid_fields, + skip_field_data_validation=skip_field_data_validation, ) model._fields = model.construct_fields() for f in model._fields: diff --git a/src/parcels/convert.py b/src/parcels/convert.py index 5dfdc75a5..be2a756bb 100644 --- a/src/parcels/convert.py +++ b/src/parcels/convert.py @@ -339,8 +339,6 @@ def nemo_to_sgrid(*, fields: dict[str, xr.Dataset | xr.DataArray], coords: xr.Da ds = _maybe_remove_depth_from_lonlat(ds) ds = _set_axis_attrs(ds, _NEMO_AXIS_VARNAMES) - ds = ds.fillna(0) # set NaN values to 0 - expected_axes = set("XYZT") # TODO: Update after we have support for 2D spatial fields if missing_axes := (expected_axes - set(ds.cf.axes)): raise ValueError( @@ -544,8 +542,6 @@ def copernicusmarine_to_sgrid( ds = _maybe_rename_coords(ds, _COPERNICUS_MARINE_AXIS_VARNAMES) - ds = ds.fillna(0) # set NaN values to 0 - if "W" in ds.data_vars: # Negate W to convert from up positive to down positive (as that's the direction of positive z) ds["W"].data *= -1 diff --git a/tests/test_fieldset.py b/tests/test_fieldset.py index ca492b4d6..26310f703 100644 --- a/tests/test_fieldset.py +++ b/tests/test_fieldset.py @@ -330,24 +330,16 @@ def test_fieldset_from_sgrid_conventions(ds_name): assert len(fieldset.fields) > 0 -def test_fieldset_keep_nan_values(): +def test_fieldset_skip_field_data_validation(): ds = datasets_structured["ds_2d_left"] ds["U_A_grid"][:] = np.nan fieldset = FieldSet.from_sgrid_conventions(ds, mesh="flat") assert np.isfinite(fieldset.U_A_grid.data).all() - fieldset = FieldSet.from_sgrid_conventions(ds, mesh="flat", fill_nan_to_zero=False) + fieldset = FieldSet.from_sgrid_conventions(ds, mesh="flat", skip_field_data_validation=True) assert np.isnan(fieldset.U_A_grid.data).all() -def test_fieldset_assert_valid_fields(): - ds = datasets_structured["ds_2d_left"].drop("time") - - FieldSet.from_sgrid_conventions(ds, mesh="flat", assert_valid_fields=False) - with pytest.raises(ValueError): - FieldSet.from_sgrid_conventions(ds, mesh="flat") - - def test_fieldset_add_error_on_duplicate_fields(): """Test that adding FieldSets with overlapping field names raises a ValueError.""" ds1 = datasets_structured["ds_2d_left"][["U_A_grid", "V_A_grid", "grid"]].rename({"U_A_grid": "U", "V_A_grid": "V"})