Merge pull request #438 from DHI/fix_elements_order

Fix problem with dfsu elements read out of order

modelskill/model/dfsu.py

@@ -4,6 +4,7 @@
 from typing import Optional, get_args
 
 import mikeio
+import numpy as np
 import pandas as pd
 
 from ._base import SpatialField, _validate_overlap_in_time, SelectedItems
@@ -169,7 +170,9 @@ def _extract_point(
 
         method = spatial_method or "inverse_distance"
         assert method in ["nearest", "contained", "inverse_distance"]
-        n_nearest = 5 if method == "inverse_distance" else 1
+        n_nearest = (
+            min(5, self.data.geometry.n_elements) if method == "inverse_distance" else 1
+        )
 
         x, y, z = observation.x, observation.y, observation.z
         if not self._in_domain(x, y):
@@ -196,6 +199,10 @@ def _extract_point(
                 elemids = self.data.geometry.find_nearest_elements(
                     x, y, n_nearest=n_nearest
                 )
+                # sort elemids, to ensure consistent results with all versions of mikeio
+                if isinstance(elemids, np.ndarray):
+                    elemids = np.sort(elemids)
+
                 ds = self.data.read(elements=elemids, items=self.sel_items.all)
                 ds_model = (
                     ds.interp(x=x, y=y, n_nearest=n_nearest) if n_nearest > 1 else ds

tests/test_config.py

@@ -59,5 +59,5 @@ def test_comparison_from_yml():
     assert cc["Klagshamn"].quantity.name == "Water Level"
     assert cc["Drogden"].quantity.name == "Water Level"
 
-    assert cc.score()["HD"] == pytest.approx(0.2331039)
+    assert cc.score()["HD"] == pytest.approx(0.233009409631)
     # spatial interp nearest: 0.23287298772

tests/test_multivariable_compare.py

@@ -89,9 +89,9 @@ def test_mv_skill(cc_1model):
     assert list(df.reset_index().observation) == cc_1model.obs_names
     assert df.index.names[0] == "observation"
     assert df.index.names[1] == "quantity"
-    assert pytest.approx(df.iloc[1].rmse) == 0.22792652
+    assert pytest.approx(df.iloc[1].rmse) == 0.22492342229
     idx = ("HKNA_wind", "Wind speed")
-    assert pytest.approx(df.loc[idx].rmse) == 1.5610323
+    assert pytest.approx(df.loc[idx].rmse) == 1.305358970478
     # spatial_interp nearest: 0.22359663 and 1.2761789
 
 
@@ -101,34 +101,34 @@ def test_mv_mm_skill(cc):
     assert df.index.names[1] == "observation"
     assert df.index.names[2] == "quantity"
     idx = ("SW_1", "HKNA_wind", "Wind speed")
-    assert pytest.approx(df.loc[idx].rmse) == 1.5610323
+    assert pytest.approx(df.loc[idx].rmse) == 1.30535897
     # spatial_interp nearest: 1.27617894455
 
     df = cc.sel(model="SW_1").skill().to_dataframe()
     assert df.index.names[0] == "observation"
     assert df.index.names[1] == "quantity"
-    assert pytest.approx(df.iloc[1].rmse) == 0.2279265
+    assert pytest.approx(df.iloc[1].rmse) == 0.2249234222997
     # spatial interp nearest: 0.22359663
     idx = ("HKNA_wind", "Wind speed")
-    assert pytest.approx(df.loc[idx].rmse) == 1.5610323
+    assert pytest.approx(df.loc[idx].rmse) == 1.30535897
 
     df = cc.sel(quantity="Wind speed").skill().to_dataframe()
     assert df.index.names[0] == "model"
     assert df.index.names[1] == "observation"
     idx = ("SW_1", "HKNA_wind")
-    assert pytest.approx(df.loc[idx].rmse) == 1.5610323
+    assert pytest.approx(df.loc[idx].rmse) == 1.30535897
 
 
 def test_mv_mm_mean_skill(cc):
     df = cc.mean_skill().to_dataframe()
     assert df.index.names[0] == "model"
     assert df.index.names[1] == "quantity"
     idx = ("SW_1", "Wind speed")
-    assert pytest.approx(df.loc[idx].r2) == 0.6246696
+    assert pytest.approx(df.loc[idx].r2) == 0.643293404624
     # spatial_interp nearest: 0.63344531
 
     df = cc.sel(quantity="Significant wave height").mean_skill().to_dataframe()
-    assert pytest.approx(df.loc["SW_1"].cc) == 0.9627803
+    assert pytest.approx(df.loc["SW_1"].cc) == 0.9640104274957
     # spatial_interp nearest: 0.963095
 
 

tests/test_pointcompare.py

@@ -100,16 +100,39 @@ def test_skill_from_observation_with_missing_values(modelresult_oresund_WL):
     assert not np.any(np.isnan(df))
 
 
+def test_score_two_elements():
+    mr = ms.model_result("tests/testdata/two_elements.dfsu", item=0)
+
+    obs_df = pd.DataFrame(
+        [2.0, 2.0], index=pd.date_range("2020-01-01", periods=2, freq="D")
+    )
+
+    # observation is in the center of the second element
+    obs = ms.PointObservation(obs_df, item=0, x=2.0, y=2.0, name="obs")
+
+    cc = ms.match(obs, mr, spatial_method="contained")
+
+    assert cc.score(metric=root_mean_squared_error)["two_elements"] == pytest.approx(
+        0.0
+    )
+
+    cc_default = ms.match(obs, mr)
+
+    assert cc_default.score(metric=root_mean_squared_error)[
+        "two_elements"
+    ] == pytest.approx(0.0)
+
+
 def test_score(modelresult_oresund_WL, klagshamn, drogden):
     mr = modelresult_oresund_WL
 
     cc = ms.match([klagshamn, drogden], mr)
 
     assert cc.score(metric=root_mean_squared_error)[
         "Oresund2D_subset"
-    ] == pytest.approx(0.19870695)
+    ] == pytest.approx(0.198637164895926)
     sk = cc.skill(metrics=[root_mean_squared_error, mean_absolute_error])
-    sk.root_mean_squared_error.data.mean() == pytest.approx(0.19870695)
+    sk.root_mean_squared_error.data.mean() == pytest.approx(0.198637164895926)
 
 
 def test_weighted_score(modelresult_oresund_WL):