diff --git a/docs/examples/Metocean_track_comparison.qmd b/docs/examples/Metocean_track_comparison.qmd index 58219a79d..f518e30c5 100644 --- a/docs/examples/Metocean_track_comparison.qmd +++ b/docs/examples/Metocean_track_comparison.qmd @@ -36,7 +36,7 @@ Csv files have no metadata on which quantity it contains, we add this manually, ```{python} o1 = ms.TrackObservation(df, item="surface_elevation", name='alti', - quantity=ms.Quantity(name="Surface Elevation", unit="meter")) + quantity=ms.Quantity(name="Surface Elevation", unit="m")) o1 ``` @@ -67,7 +67,7 @@ mr df = pd.read_csv('../data/altimetry_NorthSea_20171027.csv', index_col=0, parse_dates=True) o1 = ms.TrackObservation(df, item=2, name='alti', - quantity=ms.Quantity(name="Surface Elevation", unit="meter")) + quantity=ms.Quantity(name="Surface Elevation", unit="m")) o1 ``` @@ -91,7 +91,7 @@ mr = ms.model_result('../data/NorthSeaHD_and_windspeed.dfsu', df = pd.read_csv('../data/altimetry_NorthSea_20171027.csv', index_col=0, parse_dates=True) o1 = ms.TrackObservation(df, item=2, name='alti', - quantity=ms.Quantity(name="Surface Elevation", unit="meter")) + quantity=ms.Quantity(name="Surface Elevation", unit="m")) cmp = ms.match(o1, mr) cmp ``` diff --git a/src/modelskill/matching.py b/src/modelskill/matching.py index b7ac4b983..41b514d23 100644 --- a/src/modelskill/matching.py +++ b/src/modelskill/matching.py @@ -204,6 +204,7 @@ def match( spatial_method: str | None = None, spatial_tolerance: float = 1e-3, obs_no_overlap: Literal["ignore", "error", "warn"] = "error", + check_quantity: Literal["error", "ignore"] = "ignore", ) -> Comparer: ... @@ -216,6 +217,7 @@ def match( spatial_method: str | None = None, spatial_tolerance: float = 1e-3, obs_no_overlap: Literal["ignore", "error", "warn"] = "error", + check_quantity: Literal["error", "ignore"] = "ignore", ) -> ComparerCollection: ... @@ -227,6 +229,7 @@ def match( spatial_method: str | None = None, spatial_tolerance: float = 1e-3, obs_no_overlap: Literal["ignore", "error", "warn"] = "error", + check_quantity: Literal["error", "ignore"] = "ignore", ): """Match observation and model result data in space and time @@ -260,6 +263,10 @@ def match( and TrackObservation, by default 1e-3. obs_no_overlap: str, optional How to handle observations with no overlap with model results. One of: 'ignore', 'error', 'warn', by default 'error'. + check_quantity: str, optional + How to handle a model result whose unit differs from the observation's, + see [](`~modelskill.Quantity.is_compatible`). One of: 'error', 'ignore', + by default 'ignore'. Returns ------- @@ -273,6 +280,11 @@ def match( from_matched - Create a Comparer from observation and model results that are already matched """ + if check_quantity not in ("error", "ignore"): + raise ValueError( + f"check_quantity must be 'error' or 'ignore', not {check_quantity!r}" + ) + if isinstance(obs, get_args(ObsInputType)): return _match_single_obs( obs, @@ -281,6 +293,7 @@ def match( spatial_method=spatial_method, spatial_tolerance=spatial_tolerance, obs_no_overlap=obs_no_overlap, + check_quantity=check_quantity, ) if isinstance(obs, Collection): @@ -322,6 +335,7 @@ def match( spatial_method=spatial_method, spatial_tolerance=spatial_tolerance, obs_no_overlap=obs_no_overlap, + check_quantity=check_quantity, ) for o in obs ] @@ -339,6 +353,7 @@ def _match_single_obs( spatial_method: str | None, spatial_tolerance: float, obs_no_overlap: Literal["ignore", "error", "warn"], + check_quantity: Literal["error", "ignore"], ) -> Comparer | None: if isinstance(mod, get_args(MRInputType)): models: list = [mod] @@ -364,6 +379,14 @@ def _match_single_obs( else: matching_obs = m + if check_quantity == "error" and not obs.quantity.is_compatible( + matching_obs.quantity + ): + raise ValueError( + f"Model result '{m.name}' has quantity {matching_obs.quantity}, " + f"which is not compatible with {obs.quantity} of observation " + f"'{obs.name}'." + ) raw_mod_data[m.name] = matching_obs matched_data = _match_space_time( diff --git a/src/modelskill/quantity.py b/src/modelskill/quantity.py index 99d2eb4c2..ed779d477 100644 --- a/src/modelskill/quantity.py +++ b/src/modelskill/quantity.py @@ -58,29 +58,33 @@ def __repr__(self): # hide is_directional if False to avoid clutter return f"Quantity(name='{self.name}', unit='{self.unit}')" - def is_compatible(self, other) -> bool: + def is_compatible(self, other: "Quantity") -> bool: """Check if the quantity is compatible with another quantity + Two quantities are compatible when their units agree. Names are not + compared, since the same physical quantity is often named differently + in observation and model, e.g. "Water Level" and "Surface Elevation". + Units are compared as written, so "meter" and "m" differ. An undefined + quantity, or one without a unit, is compatible with any other. + Examples -------- ```{python} - wl = Quantity(name="Water Level", unit="meter") - ws = Quantity(name="Wind Speed", unit="meter per second") + wl = Quantity(name="Water Level", unit="m") + ws = Quantity(name="Wind Speed", unit="m/s") wl.is_compatible(ws) ``` ```{python} - uq = Quantity(name="Undefined", unit="Undefined") - wl.is_compatible(uq) + wl.is_compatible(Quantity(name="Surface Elevation", unit="m")) + ``` + ```{python} + wl.is_compatible(Quantity.undefined()) ``` """ - - if self == other: + if _is_undefined_unit(self.unit) or _is_undefined_unit(other.unit): return True - if (self.name == "Undefined") or (other.name == "Undefined"): - return True - - return False + return self.unit == other.unit @staticmethod def undefined() -> "Quantity": @@ -185,7 +189,14 @@ def from_mikeio_eum_name(type_name: str) -> "Quantity": raise ValueError( f"{type_name=} is not recognized as a known type. Please create a Quantity(name='{type_name}' unit='')" ) - unit = etype.units[0].name + unit = etype.units[0].short_name is_directional = unit == "degree" warnings.warn(f"{unit=} was automatically set for {type_name=}") return Quantity(name=type_name, unit=unit, is_directional=is_directional) + + +def _is_undefined_unit(unit: str) -> bool: + # "" from Quantity.undefined() and from res1d/EPANET results, which carry no + # unit; "undefined" from mikeio items of EUM type Undefined; "Undefined" from + # earlier modelskill versions + return unit in ("", "undefined", "Undefined") diff --git a/tests/test_match.py b/tests/test_match.py index 3324f4788..c3ee757ea 100644 --- a/tests/test_match.py +++ b/tests/test_match.py @@ -7,6 +7,7 @@ import modelskill as ms from modelskill.comparison._comparison import ItemSelection from modelskill.model.dfsu import DfsuModelResult + try: from modelskill.network import _make_basic_network except ImportError: @@ -1064,7 +1065,8 @@ def test_network_match_multi_obs_multi_model_comprehensive( def test_network_match_error_non_node_observation(network_mr, point_obs_error): """Test that non-NodeObservation raises appropriate error""" with pytest.raises( - TypeError, match="NetworkModelResult supports NodeObservation and ReachObservation" + TypeError, + match="NetworkModelResult supports NodeObservation and ReachObservation", ): ms.match(point_obs_error, network_mr) @@ -1076,3 +1078,39 @@ def test_match_nodeobs_with_other_result(node_obs1, mr1): match="Extraction from .* to is not implemented.", ): ms.match(node_obs1, mr1) + + +@pytest.fixture +def wl_obs_discharge_model(): + time = pd.date_range("2020-01-01", periods=5, freq="h") + obs = ms.PointObservation( + pd.Series(np.arange(5.0), index=time, name="obs"), + x=0.0, + y=0.0, + quantity=ms.Quantity(name="Water Level", unit="m"), + ) + mod = ms.PointModelResult( + pd.Series(np.arange(5.0), index=time, name="mod"), + x=0.0, + y=0.0, + quantity=ms.Quantity(name="Discharge", unit="m^3/s"), + ) + return obs, mod + + +def test_match_raises_on_incompatible_quantity(wl_obs_discharge_model): + obs, mod = wl_obs_discharge_model + with pytest.raises(ValueError, match="not compatible"): + ms.match(obs, mod, check_quantity="error") + + +def test_match_check_quantity_ignore(wl_obs_discharge_model): + obs, mod = wl_obs_discharge_model + cmp = ms.match(obs, mod, check_quantity="ignore") + assert cmp.n_points == 5 + + +def test_match_check_quantity_invalid_value(wl_obs_discharge_model): + obs, mod = wl_obs_discharge_model + with pytest.raises(ValueError, match="check_quantity"): + ms.match(obs, mod, check_quantity="warn") diff --git a/tests/test_pointcompare.py b/tests/test_pointcompare.py index a75206574..e7219622c 100644 --- a/tests/test_pointcompare.py +++ b/tests/test_pointcompare.py @@ -22,7 +22,7 @@ def drogden(): item=0, x=355568.0, y=6156863.0, - quantity=ms.Quantity("Water Level", unit="meter"), + quantity=ms.Quantity("Water Level", unit="m"), ) @@ -149,7 +149,7 @@ def test_weighted_score(modelresult_oresund_WL): x=355568.0, y=6156863.0, quantity=ms.Quantity( - "Water Level", unit="meter" + "Water Level", unit="m" ), # not sure if this is relevant in this test ) @@ -176,7 +176,7 @@ def test_weighted_score(modelresult_oresund_WL): x=355568.0, y=6156863.0, quantity=ms.Quantity( - "Water Level", unit="meter" + "Water Level", unit="m" ), # not sure if this is relevant in this test weight=0.1, ) diff --git a/tests/test_quantity.py b/tests/test_quantity.py index efb6c6ca2..963ddedfc 100644 --- a/tests/test_quantity.py +++ b/tests/test_quantity.py @@ -12,14 +12,13 @@ def test_from_EUMType_string(): # mikeio.EUMType.Significant_wave_height q = ms.Quantity.from_mikeio_eum_name("Significant_wave_height") - assert q.unit == "meter" + assert q.unit == "m" with pytest.warns(match="unit"): # mikeio.EUMType.Discharge q = ms.Quantity.from_mikeio_eum_name("Discharge") - # TODO should this be "meter^3 per second"? or m3/s? - assert q.unit == "meter_pow_3_per_sec" + assert q.unit == "m^3/s" def test_unknown_quantity_raises_error(): @@ -45,3 +44,33 @@ def test_from_cf_attrs_incomplete(): q = ms.Quantity.from_cf_attrs({"long_name": "Wind speed"}) assert q.name == "" assert q.unit == "" + + +def test_is_compatible_ignores_name(): + wl = ms.Quantity(name="Water Level", unit="m") + assert wl.is_compatible(ms.Quantity(name="Surface Elevation", unit="m")) + + +def test_is_compatible_different_units(): + wl = ms.Quantity(name="Water Level", unit="m") + assert not wl.is_compatible(ms.Quantity(name="Discharge", unit="m^3/s")) + + +def test_is_compatible_compares_units_as_written(): + wl = ms.Quantity(name="Water Level", unit="meter") + assert not wl.is_compatible(ms.Quantity(name="Water Level", unit="m")) + + +@pytest.mark.parametrize( + "other", + [ + ms.Quantity.undefined(), + ms.Quantity(name="Undefined", unit="Undefined"), + ms.Quantity(name="Undefined", unit="undefined"), + ms.Quantity(name="Pressure", unit=""), + ], +) +def test_is_compatible_without_unit(other): + q = ms.Quantity(name="Pressure", unit="MetresWater") + assert q.is_compatible(other) + assert other.is_compatible(q)