From 287af3e21d3fc1049e9a753dae5ef171e5e3d4f7 Mon Sep 17 00:00:00 2001 From: Vecko <36369090+VeckoTheGecko@users.noreply.github.com> Date: Wed, 7 May 2025 11:07:53 +0200 Subject: [PATCH 1/5] Add fieldset tests for incompatible timebases --- tests/v4/test_fieldset.py | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/tests/v4/test_fieldset.py b/tests/v4/test_fieldset.py index 815b4da065..a50aa50fa5 100644 --- a/tests/v4/test_fieldset.py +++ b/tests/v4/test_fieldset.py @@ -2,7 +2,9 @@ import numpy as np import pytest +import xarray as xr +from parcels._datasets.structured.generic import T as T_structured from parcels._datasets.structured.generic import datasets as datasets_structured from parcels.field import Field, VectorField from parcels.fieldset import FieldSet @@ -87,3 +89,31 @@ def test_fieldset_time_interval(): assert fieldset.time_interval.left == np.datetime64("2000-01-02") assert fieldset.time_interval.right == np.datetime64("2001-01-01") + + +def test_fieldset_init_incompatible_timebases(): + ds1 = ds.copy() + ds1["time"] = xr.date_range("2000", "2001", T_structured, calendar="365_day", use_cftime=True) + + grid = Grid(ds1) + U = Field("U", ds1["U (A grid)"], grid, mesh_type="flat") + V = Field("V", ds1["V (A grid)"], grid, mesh_type="flat") + UV = VectorField("UV", U, V) + + ds2 = ds.copy() + ds2["time"] = xr.date_range("2000", "2001", T_structured, calendar="360_day", use_cftime=True) + grid2 = Grid(ds2) + incompatible_timebase = Field("test", ds2["data_g"], grid2, mesh_type="flat") + + with pytest.raises(ValueError): + FieldSet([U, V, UV, incompatible_timebase]) + + +def test_fieldset_add_field_incompatible_timebases(fieldset): + ds_test = ds.copy() + ds_test["time"] = xr.date_range("2000", "2001", T_structured, calendar="360_day", use_cftime=True) + grid = Grid(ds_test) + field = Field("test_field", ds_test["data_g"], grid, mesh_type="flat") + + with pytest.raises(ValueError): + fieldset.add_field(field, "test_field") From 543328ddf09cb380206b3b0fcd36fe6df0d8f3cd Mon Sep 17 00:00:00 2001 From: Vecko <36369090+VeckoTheGecko@users.noreply.github.com> Date: Wed, 7 May 2025 11:21:09 +0200 Subject: [PATCH 2/5] Add test_field_init_fail_on_bad_timebase and better error message --- tests/v4/test_field.py | 22 +++++++++++++++++++++- 1 file changed, 21 insertions(+), 1 deletion(-) diff --git a/tests/v4/test_field.py b/tests/v4/test_field.py index 6a77e088c9..91330e447f 100644 --- a/tests/v4/test_field.py +++ b/tests/v4/test_field.py @@ -4,6 +4,7 @@ import xarray as xr from parcels import Field +from parcels._datasets.structured.generic import T as T_structured from parcels._datasets.structured.generic import datasets as structured_datasets from parcels._datasets.unstructured.generic import datasets as unstructured_datasets from parcels.v4.grid import Grid @@ -58,7 +59,7 @@ def test_field_incompatible_combination(data, grid): ), # TODO: Perhaps this test should be expanded to cover more datasets? ], ) -def test_field_structured_grid_creation(data, grid): +def test_field_init_structured_grid(data, grid): """Test creating a field.""" field = Field( name="test_field", @@ -70,6 +71,25 @@ def test_field_structured_grid_creation(data, grid): assert field.grid == grid +@pytest.mark.parametrize("numpy_dtype", ["timedelta64[s]", "float64"]) +def test_field_init_fail_on_bad_timebase(numpy_dtype): + """Tests that field initialisation fails when the timebase isn't given as datetime object (i.e., is float or timedelta).""" + ds = structured_datasets["ds_2d_left"].copy() + ds["time"] = np.arange(0, T_structured, dtype=numpy_dtype) + + data = ds["data_g"] + grid = Grid(ds) + with pytest.raises( + ValueError, + match="Error getting time interval.*. Are you sure that the time dimension on the xarray dataset is stored as datetime or cftime datetime objects\?", + ): + Field( + name="test_field", + data=data, + grid=grid, + ) + + @pytest.mark.parametrize( "data,grid", [ From 6d70f7c50140c5298f420e1501dfd952d4ae4495 Mon Sep 17 00:00:00 2001 From: Vecko <36369090+VeckoTheGecko@users.noreply.github.com> Date: Wed, 7 May 2025 11:23:14 +0200 Subject: [PATCH 3/5] Update import alias --- tests/v4/test_field.py | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/tests/v4/test_field.py b/tests/v4/test_field.py index 91330e447f..847770c791 100644 --- a/tests/v4/test_field.py +++ b/tests/v4/test_field.py @@ -5,8 +5,8 @@ from parcels import Field from parcels._datasets.structured.generic import T as T_structured -from parcels._datasets.structured.generic import datasets as structured_datasets -from parcels._datasets.unstructured.generic import datasets as unstructured_datasets +from parcels._datasets.structured.generic import datasets as datasets_structured +from parcels._datasets.unstructured.generic import datasets as datasets_unstructured from parcels.v4.grid import Grid @@ -37,7 +37,7 @@ def test_field_init_param_types(): pytest.param(ux.UxDataArray(), Grid(xr.Dataset()), id="uxdata-grid"), pytest.param( xr.DataArray(), - unstructured_datasets["stommel_gyre_delaunay"].uxgrid, + datasets_unstructured["stommel_gyre_delaunay"].uxgrid, id="xarray-uxgrid", ), ], @@ -55,7 +55,7 @@ def test_field_incompatible_combination(data, grid): "data,grid", [ pytest.param( - structured_datasets["ds_2d_left"]["data_g"], Grid(structured_datasets["ds_2d_left"]), id="ds_2d_left" + datasets_structured["ds_2d_left"]["data_g"], Grid(datasets_structured["ds_2d_left"]), id="ds_2d_left" ), # TODO: Perhaps this test should be expanded to cover more datasets? ], ) @@ -74,7 +74,7 @@ def test_field_init_structured_grid(data, grid): @pytest.mark.parametrize("numpy_dtype", ["timedelta64[s]", "float64"]) def test_field_init_fail_on_bad_timebase(numpy_dtype): """Tests that field initialisation fails when the timebase isn't given as datetime object (i.e., is float or timedelta).""" - ds = structured_datasets["ds_2d_left"].copy() + ds = datasets_structured["ds_2d_left"].copy() ds["time"] = np.arange(0, T_structured, dtype=numpy_dtype) data = ds["data_g"] @@ -94,7 +94,7 @@ def test_field_init_fail_on_bad_timebase(numpy_dtype): "data,grid", [ pytest.param( - structured_datasets["ds_2d_left"]["data_g"], Grid(structured_datasets["ds_2d_left"]), id="ds_2d_left" + datasets_structured["ds_2d_left"]["data_g"], Grid(datasets_structured["ds_2d_left"]), id="ds_2d_left" ), ], ) From 53583e91c0bf39c7dd289fba908dbd64ae4e385b Mon Sep 17 00:00:00 2001 From: Vecko <36369090+VeckoTheGecko@users.noreply.github.com> Date: Wed, 7 May 2025 11:46:01 +0200 Subject: [PATCH 4/5] Add assert_compatible_calendars() to FieldSet() and FieldSet.add_field() This ensures that all fields in a FieldSet always have calendars that are compatible with each other This commit also adds CalendarError, DatetimeLike, and time helpers --- parcels/_core/utils/time.py | 30 ++++++++++++++++++++++++++++- parcels/_typing.py | 6 +++++- parcels/fieldset.py | 38 +++++++++++++++++++++++++++++++++++++ tests/v4/test_field.py | 4 ++-- tests/v4/test_fieldset.py | 8 ++++---- 5 files changed, 78 insertions(+), 8 deletions(-) diff --git a/parcels/_core/utils/time.py b/parcels/_core/utils/time.py index de53953097..caef07b122 100644 --- a/parcels/_core/utils/time.py +++ b/parcels/_core/utils/time.py @@ -1,13 +1,16 @@ from __future__ import annotations from datetime import datetime -from typing import TypeVar +from typing import TYPE_CHECKING, TypeVar import cftime import numpy as np T = TypeVar("T", datetime, cftime.datetime) +if TYPE_CHECKING: + from parcels._typing import DatetimeLike + class TimeInterval: """A class representing a time interval between two datetime objects. @@ -70,3 +73,28 @@ def is_compatible(t1: datetime | cftime.datetime, t2: datetime | cftime.datetime return False else: return True + + +def get_datetime_type_calendar( + example_datetime: DatetimeLike, +) -> tuple[type, str | None]: + """Get the type and calendar of a datetime object. + + Parameters + ---------- + example_datetime : datetime, cftime.datetime, or np.datetime64 + The datetime object to check. + + Returns + ------- + tuple[type, str | None] + A tuple containing the type of the datetime object and its calendar. + The calendar will be None if the datetime object is not a cftime datetime object. + """ + calendar = None + try: + calendar = example_datetime.calendar + except AttributeError: + # datetime isn't a cftime datetime object + pass + return type(example_datetime), calendar diff --git a/parcels/_typing.py b/parcels/_typing.py index 731969dbb0..cfa96dd376 100644 --- a/parcels/_typing.py +++ b/parcels/_typing.py @@ -8,8 +8,12 @@ import os from collections.abc import Callable +from datetime import datetime from typing import Any, Literal, get_args +import numpy as np +from cftime import datetime as cftime_datetime + InterpMethodOption = Literal[ "linear", "nearest", @@ -30,7 +34,7 @@ VectorType = Literal["3D", "3DSigma", "2D"] | None # corresponds with `vector_type` GridIndexingType = Literal["pop", "mom5", "mitgcm", "nemo", "croco"] # corresponds with `gridindexingtype` NetcdfEngine = Literal["netcdf4", "xarray", "scipy"] - +DatetimeLike = datetime | cftime_datetime | np.datetime64 KernelFunction = Callable[..., None] diff --git a/parcels/fieldset.py b/parcels/fieldset.py index ba76526794..71affd6cf2 100644 --- a/parcels/fieldset.py +++ b/parcels/fieldset.py @@ -1,13 +1,21 @@ +from __future__ import annotations + import functools +from collections.abc import Iterable +from typing import TYPE_CHECKING import numpy as np import xarray as xr +from parcels._core.utils.time import get_datetime_type_calendar +from parcels._core.utils.time import is_compatible as datetime_is_compatible from parcels._reprs import fieldset_repr from parcels._typing import Mesh from parcels.field import Field, VectorField from parcels.v4.grid import Grid +if TYPE_CHECKING: + from parcels._typing import DatetimeLike __all__ = ["FieldSet"] @@ -45,6 +53,7 @@ def __init__(self, fields: list[Field | VectorField]): for field in fields: if not isinstance(field, (Field, VectorField)): raise ValueError(f"Expected `field` to be a Field or VectorField object. Got {field}") + assert_compatible_calendars(fields) self.fields = {f.name: f for f in fields} self.constants = {} @@ -126,6 +135,7 @@ def add_field(self, field: Field, name: str | None = None): """ if not isinstance(field, (Field, VectorField)): raise ValueError(f"Expected `field` to be a Field or VectorField object. Got {type(field)}") + assert_compatible_calendars((*self.fields.values(), field)) name = field.name if name is None else name @@ -235,3 +245,31 @@ def add_constant(self, name, value): # return nextTime # else: # return time + nSteps * dt + + +class CalendarError(Exception): # TODO: Move to a parcels errors module + """Exception raised when the calendar of a field is not compatible with the rest of the Fields. The user should ensure that they only add fields to a FieldSet that have compatible CFtime calendars.""" + + +def assert_compatible_calendars(fields: Iterable[Field]): + time_intervals = [f.time_interval for f in fields if f.time_interval is not None] + reference_datetime_object = time_intervals[0].left + + for field in fields: + if field.time_interval is None: + continue + + if not datetime_is_compatible(reference_datetime_object, field.time_interval.left): + msg = format_calendar_error_message(field, reference_datetime_object) + raise CalendarError(msg) + + +def format_calendar_error_message(field: Field, reference_datetime: DatetimeLike) -> str: + def datetime_to_msg(example_datetime: DatetimeLike) -> str: + datetime_type, calendar = get_datetime_type_calendar(example_datetime) + msg = str(datetime_type) + if calendar is not None: + msg += f" with cftime calendar {calendar}'" + return msg + + return f"Expected field {field.name!r} to have calendar compatible with datetime object {datetime_to_msg(reference_datetime)}. Got field with calendar {datetime_to_msg(field.time_interval.left)}. Have you considered using xarray to update the time dimension of the dataset to have a compatible calendar?" diff --git a/tests/v4/test_field.py b/tests/v4/test_field.py index 847770c791..74daccfa70 100644 --- a/tests/v4/test_field.py +++ b/tests/v4/test_field.py @@ -72,8 +72,8 @@ def test_field_init_structured_grid(data, grid): @pytest.mark.parametrize("numpy_dtype", ["timedelta64[s]", "float64"]) -def test_field_init_fail_on_bad_timebase(numpy_dtype): - """Tests that field initialisation fails when the timebase isn't given as datetime object (i.e., is float or timedelta).""" +def test_field_init_fail_on_bad_time_type(numpy_dtype): + """Tests that field initialisation fails when the time isn't given as datetime object (i.e., is float or timedelta).""" ds = datasets_structured["ds_2d_left"].copy() ds["time"] = np.arange(0, T_structured, dtype=numpy_dtype) diff --git a/tests/v4/test_fieldset.py b/tests/v4/test_fieldset.py index a50aa50fa5..8ea5885dd2 100644 --- a/tests/v4/test_fieldset.py +++ b/tests/v4/test_fieldset.py @@ -91,7 +91,7 @@ def test_fieldset_time_interval(): assert fieldset.time_interval.right == np.datetime64("2001-01-01") -def test_fieldset_init_incompatible_timebases(): +def test_fieldset_init_incompatible_calendars(): ds1 = ds.copy() ds1["time"] = xr.date_range("2000", "2001", T_structured, calendar="365_day", use_cftime=True) @@ -103,13 +103,13 @@ def test_fieldset_init_incompatible_timebases(): ds2 = ds.copy() ds2["time"] = xr.date_range("2000", "2001", T_structured, calendar="360_day", use_cftime=True) grid2 = Grid(ds2) - incompatible_timebase = Field("test", ds2["data_g"], grid2, mesh_type="flat") + incompatible_calendar = Field("test", ds2["data_g"], grid2, mesh_type="flat") with pytest.raises(ValueError): - FieldSet([U, V, UV, incompatible_timebase]) + FieldSet([U, V, UV, incompatible_calendar]) -def test_fieldset_add_field_incompatible_timebases(fieldset): +def test_fieldset_add_field_incompatible_calendars(fieldset): ds_test = ds.copy() ds_test["time"] = xr.date_range("2000", "2001", T_structured, calendar="360_day", use_cftime=True) grid = Grid(ds_test) From 77395bcb779598b574bb40feaa18f9e9b021ac20 Mon Sep 17 00:00:00 2001 From: Vecko <36369090+VeckoTheGecko@users.noreply.github.com> Date: Wed, 7 May 2025 11:57:21 +0200 Subject: [PATCH 5/5] Add VectorField.time_interval --- parcels/field.py | 28 +++++++++++++++++++++++++++- tests/v4/test_field.py | 5 +++++ 2 files changed, 32 insertions(+), 1 deletion(-) diff --git a/parcels/field.py b/parcels/field.py index 66511da3d0..c0af36f4c9 100644 --- a/parcels/field.py +++ b/parcels/field.py @@ -166,7 +166,13 @@ def __init__( self.name = name self.data = data self.grid = grid - self.time_interval = get_time_interval(data) + try: + self.time_interval = get_time_interval(data) + except ValueError as e: + e.add_note( + f"Error getting time interval for field {name!r}. Are you sure that the time dimension on the xarray dataset is stored as datetime or cftime datetime objects?" + ) + raise e # For compatibility with parts of the codebase that rely on v3 definition of Grid. # Should be worked to be removed in v4 @@ -531,6 +537,13 @@ def __init__( self.V = V self.W = W + if W is None: + assert_same_time_interval((U, V)) + else: + assert_same_time_interval((U, V, W)) + + self.time_interval = U.time_interval + if self.W: self.vector_type = "3D" else: @@ -670,3 +683,16 @@ def get_time_interval(data: xr.DataArray | ux.UxDataArray) -> TimeInterval | Non return None return TimeInterval(data.time.values[0], data.time.values[-1]) + + +def assert_same_time_interval(fields: list[Field]) -> None: + if len(fields) == 0: + return + + reference_time_interval = fields[0].time_interval + + for field in fields[1:]: + if field.time_interval != reference_time_interval: + raise ValueError( + f"Fields must have the same time domain. {fields[0].name}: {reference_time_interval}, {field.name}: {field.time_interval}" + ) diff --git a/tests/v4/test_field.py b/tests/v4/test_field.py index 74daccfa70..74c255730d 100644 --- a/tests/v4/test_field.py +++ b/tests/v4/test_field.py @@ -105,6 +105,11 @@ def test_field_time_interval(data, grid): assert field.time_interval.right == np.datetime64("2001-01-01") +def test_vectorfield_init_different_time_intervals(): + # Tests that a VectorField raises a ValueError if the component fields have different time domains. + ... + + def test_field_unstructured_grid_creation(): ...