Make sure a skill array with the right dimensions is returned when multiple models are used (#279)

* Make sure a skill array with the right dimensions is returned when multiple nwp models are used

* Add missing test

* add clim_kwargs test for blending/steps.py

Author: RubenImhoff

pysteps/blending/skill_scores.py

@@ -232,7 +232,7 @@ def clim_regr_values(n_cascade_levels, outdir_path, n_model=0, skill_kwargs=None
     """
 
     if skill_kwargs is None:
-        skill_kwargs = dict()
+        skill_kwargs = {"n_models": 1}
 
     # First, obtain climatological skill values
     try:
@@ -242,8 +242,8 @@ def clim_regr_values(n_cascade_levels, outdir_path, n_model=0, skill_kwargs=None
     except FileNotFoundError:
         # The climatological skill values file does not exist yet, so we'll
         # use the default values from BPS2004.
-        clim_cor_values = np.array(
-            [[0.848, 0.537, 0.237, 0.065, 0.020, 0.0044, 0.0052, 0.0040]]
+        clim_cor_values = clim.get_default_skill(
+            n_cascade_levels=n_cascade_levels, n_models=skill_kwargs["n_models"]
         )
 
     clim_cor_values = clim_cor_values[n_model, :]

pysteps/blending/steps.py

@@ -368,7 +368,8 @@ def forecast(
         vel_pert_kwargs = dict()
 
     if clim_kwargs is None:
-        clim_kwargs = dict()
+        # Make sure clim_kwargs at least contains the number of models
+        clim_kwargs = dict({"n_models": precip_models.shape[0]})
 
     if mask_kwargs is None:
         mask_kwargs = dict()

pysteps/tests/test_blending_skill_scores.py

@@ -69,6 +69,8 @@
     "expected_cor_t0",
     "expected_cor_nwp_lt",
     "expected_cor_nowcast_lt",
+    "n_model",
+    "number_of_models",
 )
 
 # Test function values
@@ -107,6 +109,8 @@
                 0.996475,
             ]
         ),
+        0,
+        None,
     ),
     (
         obs_6lev,
@@ -122,6 +126,8 @@
             [0.97455941, 0.9356775, 0.81972779, 0.55202975, 0.31534738, 0.02264599]
         ),
         np.array([0.996475, 0.996475, 0.996475, 0.996475, 0.996475, 0.996475]),
+        0,
+        1,
     ),
     (
         obs_9lev,
@@ -159,6 +165,34 @@
                 0.996475,
             ]
         ),
+        0,
+        1,
+    ),
+    (
+        obs_8lev,
+        mod_8lev,
+        0,
+        PHI_8lev,
+        None,
+        clim_cor_values_8lev,
+        regr_pars_8lev,
+        8,
+        np.repeat(1.0, 8),
+        np.repeat(1.0, 8),
+        np.array(
+            [
+                0.996475,
+                0.996475,
+                0.996475,
+                0.996475,
+                0.996475,
+                0.996475,
+                0.996475,
+                0.996475,
+            ]
+        ),
+        0,
+        1,
     ),
     (
         obs_8lev,
@@ -183,6 +217,8 @@
                 0.996475,
             ]
         ),
+        1,
+        2,
     ),
 ]
 
@@ -203,6 +239,8 @@ def test_blending_skill_scores(
     expected_cor_t0,
     expected_cor_nwp_lt,
     expected_cor_nowcast_lt,
+    n_model,
+    number_of_models,
 ):
     """Tests if the skill_score functions behave correctly. A dummy gridded
     model and observation field should be given for n_cascade_levels, which
@@ -212,6 +250,11 @@ def test_blending_skill_scores(
     extrapolation field.
 
     """
+    if number_of_models != None:
+        skill_kwargs = {"n_models": number_of_models}
+    else:
+        skill_kwargs = None
+
     domain_mask = np.full(obs[0, :, :].shape, False, dtype=bool)
 
     # Calculate the spatial correlation of the given model field
@@ -238,7 +281,10 @@ def test_blending_skill_scores(
     # a lead time in minutes
     # First, check if the climatological values are returned correctly
     correlations_clim, regr_clim = clim_regr_values(
-        n_cascade_levels=n_cascade_levels, outdir_path="./tmp/"
+        n_cascade_levels=n_cascade_levels,
+        outdir_path="./tmp/",
+        n_model=n_model,
+        skill_kwargs=skill_kwargs,
     )
     assert (
         correlations_clim.shape[0] == n_cascade_levels

pysteps/tests/test_blending_steps.py

@@ -114,7 +114,10 @@ def test_steps_blending(
 
     # Also set the outdir_path and clim_kwargs
     outdir_path_skill = "./tmp/"
-    clim_kwargs = dict({"n_models": n_models, "window_length": 30})
+    if n_models == 1:
+        clim_kwargs = None
+    else:
+        clim_kwargs = dict({"n_models": n_models, "window_length": 30})
 
     ###
     # First threshold the data and convert it to dBR