From 3716b03ff198681a2cc21b5995ed065e6236fc1b Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 13 Aug 2026 12:47:33 +0100 Subject: [PATCH 01/22] refactoring and changes to recipe # Conflicts: # src/CSET/operators/scoreswrappers.py --- src/CSET/operators/read.py | 5 + src/CSET/operators/scoreswrappers.py | 206 ++++++++++++++---- ...series_surface_difference_scores_RMSE.yaml | 4 +- 3 files changed, 171 insertions(+), 44 deletions(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index 96791d383..2fd39555b 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -167,7 +167,12 @@ def read_cubes( paths = iter_maybe(file_paths) model_names = iter_maybe(model_names) + # flattens if second element is a list + if isinstance(model_names[1], list): + model_names = [model_names[0]] + model_names[1] + # Check we have appropriate number of model names. + if model_names != (None,) and len(model_names) != len(paths): raise ValueError( f"The number of model names ({len(model_names)}) should equal " diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 681f7f5be..a0446c7c3 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -41,6 +41,29 @@ logger = logging.getLogger(__name__) +def _sort_cube_into_base_and_other(cubes): + base: Cube = cubes.extract_cube(iris.AttributeConstraint(cset_comparison_base=1)) + others: CubeList = cubes.extract( + iris.Constraint( + cube_func=lambda cube: "cset_comparison_base" not in cube.attributes + ) + ) + + return base, others + + +def _ensure_increasing_pressure_coordinates(cubes): + for cube in cubes: + try: + if len(cube.coord("pressure").points) > 2 and not is_increasing( + cube.coord("pressure").points + ): + reverse(cube, "pressure") + + except iris.exceptions.CoordinateNotFoundError: + pass + + def _sort_cubes_for_verification(cubes: CubeList): """Prepare cubes ready for verification in scores. @@ -149,6 +172,98 @@ def _sort_cubes_for_verification(cubes: CubeList): return base, other +def _process_cubes_for_verification(base: Cube, other: Cube): + """Prepare cubes ready for verification in scores. + + Parameters + ---------- + cubes: iris.cube.CubeList + A CubeList of exact 2 cubes, one from each model. + + Returns + ------- + base: iris.cube.Cube + The cube from the "analysis" in the same format as the other model. + other: iris.cube.Cube + The cube from the model in the same format as the base model. + + Raises + ------ + ValueError: "cubes should contain exactly 2 cubes." + If any other number of cubes are present. + + Notes + ----- + This operator is used for sorting the data into the correct format. It + is likely going to need to be refactored out of CSET and perhaps moved into + `CSET._utils` given common code between here and `misc.difference`. + """ + # Set cubes into correct format using code from difference operator + + # Extract just common time points. + other_model_name = other.attributes["model_name"] + base, other = _extract_common_time_points(base, other) + + # Get spatial coord names. + base_lat_name, base_lon_name = get_cube_yxcoordname(base) + other_lat_name, other_lon_name = get_cube_yxcoordname(other) + + # Ensure cubes to compare are on common differencing grid. + # This is triggered if either + # i) latitude and longitude shapes are not the same. Note grid points + # are not compared directly as these can differ through rounding + # errors. + # ii) or variables are known to often sit on different grid staggering + # in different models (e.g. cell center vs cell edge), as is the case + # for UM and LFRic comparisons. + # In future greater choice of regridding method might be applied depending + # on variable type. Linear regridding can in general be appropriate for smooth + # variables. Care should be taken with interpretation of differences + # given this dependency on regridding. + if ( + base.coord(base_lat_name).shape != other.coord(other_lat_name).shape + or base.coord(base_lon_name).shape != other.coord(other_lon_name).shape + ) or ( + base.long_name + in [ + "eastward_wind_at_10m", + "northward_wind_at_10m", + "northward_wind_at_cell_centres", + "eastward_wind_at_cell_centres", + "zonal_wind_at_pressure_levels", + "meridional_wind_at_pressure_levels", + "potential_vorticity_at_pressure_levels", + "vapour_specific_humidity_at_pressure_levels_for_climate_averaging", + ] + ): + logger.debug("Linear regridding base cube to other grid to compute differences") + base = regrid_onto_cube(base, other, method="Linear") + + # Figure out if we are comparing between UM and LFRic; flip array if so. + base_lat_direction = is_increasing(base.coord(base_lat_name).points) + other_lat_direction = is_increasing(other.coord(other_lat_name).points) + if base_lat_direction != other_lat_direction: + # Copy base cube for correct coordinate information. + other_tmp = base.copy() + # Flip the data and place in the copied cube. + other_tmp.data = np.flip( + other.data, other.coord(other_lat_name).cube_dims(other) + ) + # Use original name and units from the other cube. + other_tmp.rename(other.name()) + other_tmp.units = other.units + # Replace the cube. + other = other_tmp + + # Equalise attributes so we can merge. + fully_equalise_attributes(CubeList([base, other])) + + other.attributes["model_name"] = other_model_name + logger.debug("Base: %s\nOther: %s", base, other) + + return base, other + + def _resolve_preserve_dims( cube: Cube, data_array: xr.DataArray, @@ -274,6 +389,10 @@ def scores_rmse( scores_cube: iris.cube.Cube A cube containing the RMSE between the base and other cube. """ + base, others = _sort_cube_into_base_and_other(cubes) + rmse_cubelist = CubeList() + for other in others: + base, other = _process_cubes_for_verification(base, other) if obs_model_comparison: for cb in cubes: if "observed" in cb.long_name: @@ -283,53 +402,56 @@ def scores_rmse( else: base, other = _sort_cubes_for_verification(cubes) - # Copy the coordinates of the input cubes. - other_xr = xr.DataArray.from_iris(other) - base_xr = xr.DataArray.from_iris(base) - preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) + # Copy the coordinates of the input cubes. + other_xr = xr.DataArray.from_iris(other) + base_xr = xr.DataArray.from_iris(base) + preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) - # Scores operates on xarray data arrays, so we transform the iris cube into an array, - # apply scores, and then transform it back. - scores_cube = xr.DataArray.to_iris( - scores.continuous.rmse( - other_xr, - base_xr, - preserve_dims=preserve_dims, - ) - ) - - # If time is aggregated out, attach a scalar time coordinate with bounds - # so plotting can display the aggregated period in the title. - try: - if not scores_cube.coords("time"): - base_time = base.coord("time") - time_vals = ( - base_time.bounds.flatten() - if base_time.has_bounds() - else base_time.points + # Scores operates on xarray data arrays, so we transform the iris cube into an array, + # apply scores, and then transform it back. + scores_cube = xr.DataArray.to_iris( + scores.continuous.rmse( + other_xr, + base_xr, + preserve_dims=preserve_dims, ) - t_start = float(time_vals[0]) - t_end = float(time_vals[-1]) - t_mid = 0.5 * (t_start + t_end) + ) - scores_cube.add_aux_coord( - iris.coords.AuxCoord( - t_mid, - standard_name=base_time.standard_name, - long_name=base_time.long_name, - var_name=base_time.var_name, - units=base_time.units, - bounds=np.array([t_start, t_end]), - attributes=base_time.attributes.copy(), + # If time is aggregated out, attach a scalar time coordinate with bounds + # so plotting can display the aggregated period in the title. + try: + if not scores_cube.coords("time"): + base_time = base.coord("time") + time_vals = ( + base_time.bounds.flatten() + if base_time.has_bounds() + else base_time.points ) - ) - except iris.exceptions.CoordinateNotFoundError: - pass + t_start = float(time_vals[0]) + t_end = float(time_vals[-1]) + t_mid = 0.5 * (t_start + t_end) + + scores_cube.add_aux_coord( + iris.coords.AuxCoord( + t_mid, + standard_name=base_time.standard_name, + long_name=base_time.long_name, + var_name=base_time.var_name, + units=base_time.units, + bounds=np.array([t_start, t_end]), + attributes=base_time.attributes.copy(), + ) + ) + except iris.exceptions.CoordinateNotFoundError: + pass - scores_cube.rename(f"RMSE_of_{base.name()}") - # if preserved_coordinates == ["grid_latitude", "grid_longitude"]: - # scores_cube.add_aux_coord(time_coord) - return scores_cube + scores_cube.rename(f"RMSE_of_{base.name()}") + rmse_cubelist.append(scores_cube) + + model_name = other.attributes["model_name"] + scores_cube.attributes["model_name"] = model_name + + return rmse_cubelist def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = None): diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml index 7915c332d..fee2dd4f4 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores RMSE timeseries between $OTHER_MODEL and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores RMSE timeseries between $BASE_MODEL and $MODEL_NAMES and $SUBAREA_NAME" description: | Extracts and plots the Root Mean Square Error in $VARNAME computed over the domain for each timestep. The RMSE is calculated based on that used in the @@ -16,7 +16,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $MODEL_NAMES] constraint: operator: constraints.combine_constraints varname_constraint: From 67f6249cab826e1aee3b4c8207b31e187ad5d927 Mon Sep 17 00:00:00 2001 From: James Ross Date: Fri, 14 Aug 2026 12:46:25 +0100 Subject: [PATCH 02/22] changing more scores to account for a multi model comparison --- src/CSET/loaders/verification.py | 4 +- src/CSET/operators/scoreswrappers.py | 296 ++++++++++-------- .../surface_difference_scores_MAE.yaml | 2 +- .../surface_difference_scores_RMSE.yaml | 2 +- ...rface_difference_scores_additive_bias.yaml | 4 +- ...eseries_surface_difference_scores_MAE.yaml | 4 +- ...series_surface_difference_scores_RMSE.yaml | 4 +- ...rface_difference_scores_additive_bias.yaml | 4 +- ...ifference_scores_correlation_pearsonr.yaml | 4 +- .../timeseries_surface_rmse_scores.yaml | 43 --- 10 files changed, 172 insertions(+), 195 deletions(-) delete mode 100644 src/CSET/recipes/verification/timeseries_surface_rmse_scores.yaml diff --git a/src/CSET/loaders/verification.py b/src/CSET/loaders/verification.py index c0b136a30..f85e00b81 100644 --- a/src/CSET/loaders/verification.py +++ b/src/CSET/loaders/verification.py @@ -124,7 +124,7 @@ def load(conf: Config): variables={ "VARNAME": field, "BASE_MODEL": base_model["name"], - "OTHER_MODEL": model["name"], + "OTHER_MODELS": model["name"], "METHOD": method, "PRESERVED_COORDS": preserved_coords, "SUBAREA_NAME": conf.SUBAREA_NAME if conf.SELECT_SUBAREA else "", @@ -148,7 +148,7 @@ def load(conf: Config): variables={ "VARNAME": field, "BASE_MODEL": base_model["name"], - "OTHER_MODEL": model["name"], + "OTHER_MODELS": model["name"], "SUBAREA_NAME": conf.SUBAREA_NAME if conf.SELECT_SUBAREA else "", "SUBAREA_TYPE": conf.SUBAREA_TYPE if conf.SELECT_SUBAREA else None, "SUBAREA_EXTENT": conf.SUBAREA_EXTENT diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index a0446c7c3..380fc5160 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -386,11 +386,11 @@ def scores_rmse( Returns ------- - scores_cube: iris.cube.Cube - A cube containing the RMSE between the base and other cube. + scores_cubelist: iris.cube.CubeList + A cubelist containing the RMSE between the base and other cube. """ base, others = _sort_cube_into_base_and_other(cubes) - rmse_cubelist = CubeList() + scores_cubelist = CubeList() for other in others: base, other = _process_cubes_for_verification(base, other) if obs_model_comparison: @@ -446,12 +446,12 @@ def scores_rmse( pass scores_cube.rename(f"RMSE_of_{base.name()}") - rmse_cubelist.append(scores_cube) + scores_cubelist.append(scores_cube) model_name = other.attributes["model_name"] scores_cube.attributes["model_name"] = model_name - return rmse_cubelist + return scores_cubelist def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = None): @@ -473,56 +473,63 @@ def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = Returns ------- - scores_cube: iris.cube.Cube - A cube containing the MAE between the base and other cube. + scores_cubelist: iris.cube.CubeList + A cubelist containing the MAE between the base and other cubes. """ - base, other = _sort_cubes_for_verification(cubes) - - # Copy the coordinates of the input cubes. - other_xr = xr.DataArray.from_iris(other) - base_xr = xr.DataArray.from_iris(base) - preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) - - # Scores operates on xarray data arrays, so we transform the iris cube into an array, - # apply scores, and then transform it back. - scores_cube = xr.DataArray.to_iris( - scores.continuous.mae( - other_xr, - base_xr, - preserve_dims=preserve_dims, - ) - ) + base, others = _sort_cube_into_base_and_other(cubes) + scores_cubelist = CubeList() + for other in others: + base, other = _process_cubes_for_verification(base, other) - # If time is aggregated out, attach a scalar time coordinate with bounds - # so plotting can display the aggregated period in the title. - try: - if not scores_cube.coords("time"): - base_time = base.coord("time") - time_vals = ( - base_time.bounds.flatten() - if base_time.has_bounds() - else base_time.points + # Copy the coordinates of the input cubes. + other_xr = xr.DataArray.from_iris(other) + base_xr = xr.DataArray.from_iris(base) + preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) + + # Scores operates on xarray data arrays, so we transform the iris cube into an array, + # apply scores, and then transform it back. + scores_cube = xr.DataArray.to_iris( + scores.continuous.mae( + other_xr, + base_xr, + preserve_dims=preserve_dims, ) - t_start = float(time_vals[0]) - t_end = float(time_vals[-1]) - t_mid = 0.5 * (t_start + t_end) - - scores_cube.add_aux_coord( - iris.coords.AuxCoord( - t_mid, - standard_name=base_time.standard_name, - long_name=base_time.long_name, - var_name=base_time.var_name, - units=base_time.units, - bounds=np.array([t_start, t_end]), - attributes=base_time.attributes.copy(), + ) + + # If time is aggregated out, attach a scalar time coordinate with bounds + # so plotting can display the aggregated period in the title. + try: + if not scores_cube.coords("time"): + base_time = base.coord("time") + time_vals = ( + base_time.bounds.flatten() + if base_time.has_bounds() + else base_time.points ) - ) - except iris.exceptions.CoordinateNotFoundError: - pass + t_start = float(time_vals[0]) + t_end = float(time_vals[-1]) + t_mid = 0.5 * (t_start + t_end) - scores_cube.rename(f"MAE_of_{base.name()}") - return scores_cube + scores_cube.add_aux_coord( + iris.coords.AuxCoord( + t_mid, + standard_name=base_time.standard_name, + long_name=base_time.long_name, + var_name=base_time.var_name, + units=base_time.units, + bounds=np.array([t_start, t_end]), + attributes=base_time.attributes.copy(), + ) + ) + except iris.exceptions.CoordinateNotFoundError: + pass + + scores_cube.rename(f"MAE_of_{base.name()}") + scores_cubelist.append(scores_cube) + model_name = other.attributes["model_name"] + scores_cube.attributes["model_name"] = model_name + + return scores_cubelist def scores_additive_bias( @@ -546,55 +553,62 @@ def scores_additive_bias( Returns ------- - scores_cube: iris.cube.Cube - A cube containing the ME between the base and other cube. + scores_cubelist: iris.cube.CubeList + A cubelist containing the ME between the base and other cube. """ - base, other = _sort_cubes_for_verification(cubes) - - # Copy the coordinates of the input cubes. - other_xr = xr.DataArray.from_iris(other) - base_xr = xr.DataArray.from_iris(base) - preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) - - # Scores operates on xarray data arrays, so we transform the iris cube into an array, - # apply scores, and then transform it back. - scores_cube = xr.DataArray.to_iris( - scores.continuous.additive_bias( - other_xr, - base_xr, - preserve_dims=preserve_dims, - ) - ) + base, others = _sort_cube_into_base_and_other(cubes) + scores_cubelist = CubeList() + for other in others: + base, other = _process_cubes_for_verification(base, other) - # If time is aggregated out, attach a scalar time coordinate with bounds - # so plotting can display the aggregated period in the title. - try: - if not scores_cube.coords("time"): - base_time = base.coord("time") - time_vals = ( - base_time.bounds.flatten() - if base_time.has_bounds() - else base_time.points + # Copy the coordinates of the input cubes. + other_xr = xr.DataArray.from_iris(other) + base_xr = xr.DataArray.from_iris(base) + preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) + + # Scores operates on xarray data arrays, so we transform the iris cube into an array, + # apply scores, and then transform it back. + scores_cube = xr.DataArray.to_iris( + scores.continuous.additive_bias( + other_xr, + base_xr, + preserve_dims=preserve_dims, ) - t_start = float(time_vals[0]) - t_end = float(time_vals[-1]) - t_mid = 0.5 * (t_start + t_end) - - scores_cube.add_aux_coord( - iris.coords.AuxCoord( - t_mid, - standard_name=base_time.standard_name, - long_name=base_time.long_name, - var_name=base_time.var_name, - units=base_time.units, - bounds=np.array([t_start, t_end]), - attributes=base_time.attributes.copy(), + ) + + # If time is aggregated out, attach a scalar time coordinate with bounds + # so plotting can display the aggregated period in the title. + try: + if not scores_cube.coords("time"): + base_time = base.coord("time") + time_vals = ( + base_time.bounds.flatten() + if base_time.has_bounds() + else base_time.points ) - ) - except iris.exceptions.CoordinateNotFoundError: - pass - scores_cube.rename(f"Additive_Bias_of_{base.name()}") - return scores_cube + t_start = float(time_vals[0]) + t_end = float(time_vals[-1]) + t_mid = 0.5 * (t_start + t_end) + + scores_cube.add_aux_coord( + iris.coords.AuxCoord( + t_mid, + standard_name=base_time.standard_name, + long_name=base_time.long_name, + var_name=base_time.var_name, + units=base_time.units, + bounds=np.array([t_start, t_end]), + attributes=base_time.attributes.copy(), + ) + ) + except iris.exceptions.CoordinateNotFoundError: + pass + scores_cube.rename(f"Additive_Bias_of_{base.name()}") + scores_cubelist.append(scores_cube) + model_name = other.attributes["model_name"] + scores_cube.attributes["model_name"] = model_name + + return scores_cubelist def scores_correlation_pearsonr( @@ -621,53 +635,59 @@ def scores_correlation_pearsonr( scores_cube: iris.cube.Cube A cube containing the PC between the base and other cube. """ - base, other = _sort_cubes_for_verification(cubes) - - # Copy the coordinates of the input cubes. - other_xr = xr.DataArray.from_iris(other) - base_xr = xr.DataArray.from_iris(base) - preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) - - # Scores operates on xarray data arrays, so we transform the iris cube into an array, - # apply scores, and then transform it back. - scores_cube = xr.DataArray.to_iris( - scores.continuous.correlation.pearsonr( - other_xr, - base_xr, - preserve_dims=preserve_dims, - ) - ) + base, others = _sort_cube_into_base_and_other(cubes) + scores_cubelist = CubeList() + for other in others: + base, other = _process_cubes_for_verification(base, other) - # If time is aggregated out, attach a scalar time coordinate with bounds - # so plotting can display the aggregated period in the title. - try: - if not scores_cube.coords("time"): - base_time = base.coord("time") - time_vals = ( - base_time.bounds.flatten() - if base_time.has_bounds() - else base_time.points + # Copy the coordinates of the input cubes. + other_xr = xr.DataArray.from_iris(other) + base_xr = xr.DataArray.from_iris(base) + preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) + + # Scores operates on xarray data arrays, so we transform the iris cube into an array, + # apply scores, and then transform it back. + scores_cube = xr.DataArray.to_iris( + scores.continuous.correlation.pearsonr( + other_xr, + base_xr, + preserve_dims=preserve_dims, ) - t_start = float(time_vals[0]) - t_end = float(time_vals[-1]) - t_mid = 0.5 * (t_start + t_end) - - scores_cube.add_aux_coord( - iris.coords.AuxCoord( - t_mid, - standard_name=base_time.standard_name, - long_name=base_time.long_name, - var_name=base_time.var_name, - units=base_time.units, - bounds=np.array([t_start, t_end]), - attributes=base_time.attributes.copy(), + ) + + # If time is aggregated out, attach a scalar time coordinate with bounds + # so plotting can display the aggregated period in the title. + try: + if not scores_cube.coords("time"): + base_time = base.coord("time") + time_vals = ( + base_time.bounds.flatten() + if base_time.has_bounds() + else base_time.points ) - ) - except iris.exceptions.CoordinateNotFoundError: - pass + t_start = float(time_vals[0]) + t_end = float(time_vals[-1]) + t_mid = 0.5 * (t_start + t_end) + + scores_cube.add_aux_coord( + iris.coords.AuxCoord( + t_mid, + standard_name=base_time.standard_name, + long_name=base_time.long_name, + var_name=base_time.var_name, + units=base_time.units, + bounds=np.array([t_start, t_end]), + attributes=base_time.attributes.copy(), + ) + ) + except iris.exceptions.CoordinateNotFoundError: + pass - scores_cube.rename(f"Pearson_Correlation_of_{base.name()}") - return scores_cube + scores_cube.rename(f"Pearson_Correlation_of_{base.name()}") + scores_cubelist.append(scores_cube) + model_name = other.attributes["model_name"] + scores_cube.attributes["model_name"] = model_name + return scores_cubelist def scores_crps_for_ensemble( diff --git a/src/CSET/recipes/verification/surface_difference_scores_MAE.yaml b/src/CSET/recipes/verification/surface_difference_scores_MAE.yaml index cc77db011..a20e69095 100644 --- a/src/CSET/recipes/verification/surface_difference_scores_MAE.yaml +++ b/src/CSET/recipes/verification/surface_difference_scores_MAE.yaml @@ -29,7 +29,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/surface_difference_scores_RMSE.yaml b/src/CSET/recipes/verification/surface_difference_scores_RMSE.yaml index 9b4642cf4..65d7b5e42 100644 --- a/src/CSET/recipes/verification/surface_difference_scores_RMSE.yaml +++ b/src/CSET/recipes/verification/surface_difference_scores_RMSE.yaml @@ -26,7 +26,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/surface_difference_scores_additive_bias.yaml b/src/CSET/recipes/verification/surface_difference_scores_additive_bias.yaml index 9cc6b22ca..4b5d57bbb 100644 --- a/src/CSET/recipes/verification/surface_difference_scores_additive_bias.yaml +++ b/src/CSET/recipes/verification/surface_difference_scores_additive_bias.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$METHOD $VARNAME\nScores Mean Error between $OTHER_MODEL and $BASE_MODEL $SUBAREA_NAME" +title: "$METHOD $VARNAME\nScores Mean Error between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Error (also known as the Additive Bias) in $METHOD of $VARNAME. The ME is calculated based on that used in the @@ -24,7 +24,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml index d5973fbbb..f5d36359d 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Absolute Error timeseries between $OTHER_MODEL and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Absolute Error timeseries between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Absolute Error in $VARNAME computed over the domain for each timestep. The MAE is calculated based on that used in the @@ -17,7 +17,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml index fee2dd4f4..38eecfc6e 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores RMSE timeseries between $BASE_MODEL and $MODEL_NAMES and $SUBAREA_NAME" +title: "$VARNAME\nScores RMSE timeseries between $OTHER_MODELS and $BASE_MODEL and $SUBAREA_NAME" description: | Extracts and plots the Root Mean Square Error in $VARNAME computed over the domain for each timestep. The RMSE is calculated based on that used in the @@ -16,7 +16,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $MODEL_NAMES] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml index df9ee7c15..dbc60032c 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Error timeseries between $OTHER_MODEL and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Error timeseries between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Error (also known as the Additive Bias) in $VARNAME computed over the domain for each timestep. The ME is calculated based on that used in the @@ -15,7 +15,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml index 179a735f9..38597ec7f 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Pearson's Correlation timeseries between $OTHER_MODEL and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Pearson's Correlation timeseries between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Pearson's Correlation coefficient in $VARNAME computed over the domain for each timestep. The PC is calculated based on that used in the @@ -14,7 +14,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints varname_constraint: diff --git a/src/CSET/recipes/verification/timeseries_surface_rmse_scores.yaml b/src/CSET/recipes/verification/timeseries_surface_rmse_scores.yaml deleted file mode 100644 index db1ed972a..000000000 --- a/src/CSET/recipes/verification/timeseries_surface_rmse_scores.yaml +++ /dev/null @@ -1,43 +0,0 @@ -category: Scores -title: "$VARNAME\nRMSE timeseries between $OTHER_MODEL and $BASE_MODEL" -description: | - Extracts and plots the Root Mean Square Error in $VARNAME - for all times. The RMSE is calculated based on that used in the - package [`scores`](https://scores.readthedocs.io/en/stable/api.html#scores.continuous.rmse). - This recipe allows the preservation of the time coordinate to produce a timeseries. Therefore, the RMSE is - collapsed over all other coordinates in the cube and calculated for every timestep. - - A larger RMSE implies a greater error than a smaller RMSE. An - RMSE of zero implies the two fields match. The RMSE is calculated - on the grid point and thus a spatial view of the RMSE provides useful - information about where the differences are, or if placement errors - are domininating the score (usually indicated by dipoles). - -steps: - - operator: read.read_cubes - file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] - constraint: - operator: constraints.combine_constraints - varname_constraint: - operator: constraints.generate_var_constraint - varname: $VARNAME - cell_methods_constraint: - operator: constraints.generate_cell_methods_constraint - cell_methods: [] - varname: $VARNAME - pressure_level_constraint: - operator: constraints.generate_level_constraint - coordinate: "pressure" - levels: [] - subarea_type: $SUBAREA_TYPE - subarea_extent: $SUBAREA_EXTENT - - - operator: scoreswrappers.scores_rmse - preserved_coordinates: ["time"] - - - operator: plot.plot_line_series - series_coordinate: time - - - operator: write.write_cube_to_nc - overwrite: True From 3ffc1c6223dc26b3f4caa94c91450b467aae9fae Mon Sep 17 00:00:00 2001 From: James Ross Date: Fri, 14 Aug 2026 15:22:39 +0100 Subject: [PATCH 03/22] fix bug --- src/CSET/loaders/verification.py | 2 +- src/CSET/operators/read.py | 3 +- src/CSET/operators/scoreswrappers.py | 108 ------------------ .../generic_level_rmse_scores_profile.yaml | 4 +- 4 files changed, 5 insertions(+), 112 deletions(-) diff --git a/src/CSET/loaders/verification.py b/src/CSET/loaders/verification.py index f85e00b81..969d447da 100644 --- a/src/CSET/loaders/verification.py +++ b/src/CSET/loaders/verification.py @@ -251,7 +251,7 @@ def load(conf: Config): variables={ "VARNAME": field, "BASE_MODEL": base_model["name"], - "OTHER_MODEL": model["name"], + "OTHER_MODELS": model["name"], "PRESERVED_COORDS": ["pressure"], "AGGREGATION_MODE": "Case-study RMSE", "SUBAREA_TYPE": conf.SUBAREA_TYPE if conf.SELECT_SUBAREA else None, diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index 2fd39555b..d1e459509 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -32,6 +32,7 @@ import iris.util import numpy as np from iris.analysis.cartography import rotate_pole, rotate_winds +from shapely.measurement import length from CSET._common import iter_maybe from CSET.operators._stash_to_lfric import STASH_TO_LFRIC @@ -168,7 +169,7 @@ def read_cubes( model_names = iter_maybe(model_names) # flattens if second element is a list - if isinstance(model_names[1], list): + if length(model_names) == 2 and isinstance(model_names[1], list): model_names = [model_names[0]] + model_names[1] # Check we have appropriate number of model names. diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 380fc5160..703390c30 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -64,114 +64,6 @@ def _ensure_increasing_pressure_coordinates(cubes): pass -def _sort_cubes_for_verification(cubes: CubeList): - """Prepare cubes ready for verification in scores. - - Parameters - ---------- - cubes: iris.cube.CubeList - A CubeList of exact 2 cubes, one from each model. - - Returns - ------- - base: iris.cube.Cube - The cube from the "analysis" in the same format as the other model. - other: iris.cube.Cube - The cube from the model in the same format as the base model. - - Raises - ------ - ValueError: "cubes should contain exactly 2 cubes." - If any other number of cubes are present. - - Notes - ----- - This operator is used for sorting the data into the correct format. It - is likely going to need to be refactored out of CSET and perhaps moved into - `CSET._utils` given common code between here and `misc.difference`. - """ - # Set cubes into correct format using code from difference operator - if len(cubes) != 2: - raise ValueError("cubes should contain exactly 2 cubes.") - base: Cube = cubes.extract_cube(iris.AttributeConstraint(cset_comparison_base=1)) - other: Cube = cubes.extract_cube( - iris.Constraint( - cube_func=lambda cube: "cset_comparison_base" not in cube.attributes - ) - ) - - # If cubes contain a pressure coordinate, ensure it is increasing. - for cube in cubes: - try: - if len(cube.coord("pressure").points) > 2 and not is_increasing( - cube.coord("pressure").points - ): - reverse(cube, "pressure") - - except iris.exceptions.CoordinateNotFoundError: - pass - - # Extract just common time points. - base, other = _extract_common_time_points(base, other) - - # Get spatial coord names. - base_lat_name, base_lon_name = get_cube_yxcoordname(base) - other_lat_name, other_lon_name = get_cube_yxcoordname(other) - - # Ensure cubes to compare are on common differencing grid. - # This is triggered if either - # i) latitude and longitude shapes are not the same. Note grid points - # are not compared directly as these can differ through rounding - # errors. - # ii) or variables are known to often sit on different grid staggering - # in different models (e.g. cell center vs cell edge), as is the case - # for UM and LFRic comparisons. - # In future greater choice of regridding method might be applied depending - # on variable type. Linear regridding can in general be appropriate for smooth - # variables. Care should be taken with interpretation of differences - # given this dependency on regridding. - if ( - base.coord(base_lat_name).shape != other.coord(other_lat_name).shape - or base.coord(base_lon_name).shape != other.coord(other_lon_name).shape - ) or ( - base.long_name - in [ - "eastward_wind_at_10m", - "northward_wind_at_10m", - "northward_wind_at_cell_centres", - "eastward_wind_at_cell_centres", - "zonal_wind_at_pressure_levels", - "meridional_wind_at_pressure_levels", - "potential_vorticity_at_pressure_levels", - "vapour_specific_humidity_at_pressure_levels_for_climate_averaging", - ] - ): - logger.debug("Linear regridding base cube to other grid to compute differences") - base = regrid_onto_cube(base, other, method="Linear") - - # Figure out if we are comparing between UM and LFRic; flip array if so. - base_lat_direction = is_increasing(base.coord(base_lat_name).points) - other_lat_direction = is_increasing(other.coord(other_lat_name).points) - if base_lat_direction != other_lat_direction: - # Copy base cube for correct coordinate information. - other_tmp = base.copy() - # Flip the data and place in the copied cube. - other_tmp.data = np.flip( - other.data, other.coord(other_lat_name).cube_dims(other) - ) - # Use original name and units from the other cube. - other_tmp.rename(other.name()) - other_tmp.units = other.units - # Replace the cube. - other = other_tmp - - # Equalise attributes so we can merge. - fully_equalise_attributes(CubeList([base, other])) - logger.debug("Base: %s\nOther: %s", base, other) - - return base, other - - def _process_cubes_for_verification(base: Cube, other: Cube): """Prepare cubes ready for verification in scores. diff --git a/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml b/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml index d4b2788d7..db0a735bd 100644 --- a/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml +++ b/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nRMSE vertical profile between $OTHER_MODEL and $BASE_MODEL ($AGGREGATION_MODE)" +title: "$VARNAME\nRMSE vertical profile between $OTHER_MODELS and $BASE_MODEL ($AGGREGATION_MODE)" description: | Extracts and plots the Root Mean Square Error (RMSE) of $VARNAME as a vertical profile on pressure levels. The RMSE is calculated @@ -24,7 +24,7 @@ description: | steps: - operator: read.read_cubes file_paths: $INPUT_PATHS - model_names: [$BASE_MODEL, $OTHER_MODEL] + model_names: [$BASE_MODEL, $OTHER_MODELS] constraint: operator: constraints.combine_constraints variable_constraint: From 615066b1ee4d32ffd3af858e8bae6c1cc20a17af Mon Sep 17 00:00:00 2001 From: James Ross Date: Fri, 14 Aug 2026 15:26:01 +0100 Subject: [PATCH 04/22] make flatten simpler --- src/CSET/operators/read.py | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index d1e459509..4963ec3e2 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -32,7 +32,6 @@ import iris.util import numpy as np from iris.analysis.cartography import rotate_pole, rotate_winds -from shapely.measurement import length from CSET._common import iter_maybe from CSET.operators._stash_to_lfric import STASH_TO_LFRIC @@ -164,14 +163,14 @@ def read_cubes( FileNotFoundError If the provided path does not exist """ + # flattens if any element is a list + model_names = list(itertools.chain.from_iterable(model_names)) + # Get iterable of paths. Each path corresponds to 1 model. + paths = iter_maybe(file_paths) model_names = iter_maybe(model_names) - # flattens if second element is a list - if length(model_names) == 2 and isinstance(model_names[1], list): - model_names = [model_names[0]] + model_names[1] - # Check we have appropriate number of model names. if model_names != (None,) and len(model_names) != len(paths): From 9498fe9b14f712586d64ddecd26f763b841f1b9b Mon Sep 17 00:00:00 2001 From: James Ross Date: Fri, 14 Aug 2026 15:30:57 +0100 Subject: [PATCH 05/22] revert --- src/CSET/operators/read.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index 4963ec3e2..d1e459509 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -32,6 +32,7 @@ import iris.util import numpy as np from iris.analysis.cartography import rotate_pole, rotate_winds +from shapely.measurement import length from CSET._common import iter_maybe from CSET.operators._stash_to_lfric import STASH_TO_LFRIC @@ -163,14 +164,14 @@ def read_cubes( FileNotFoundError If the provided path does not exist """ - # flattens if any element is a list - model_names = list(itertools.chain.from_iterable(model_names)) - # Get iterable of paths. Each path corresponds to 1 model. - paths = iter_maybe(file_paths) model_names = iter_maybe(model_names) + # flattens if second element is a list + if length(model_names) == 2 and isinstance(model_names[1], list): + model_names = [model_names[0]] + model_names[1] + # Check we have appropriate number of model names. if model_names != (None,) and len(model_names) != len(paths): From 26e8f58915ab99af6a042c3776b7eb15fdcde626 Mon Sep 17 00:00:00 2001 From: James Ross Date: Fri, 14 Aug 2026 15:34:15 +0100 Subject: [PATCH 06/22] another way to flatten --- src/CSET/operators/read.py | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index d1e459509..db4c4f1cc 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -32,7 +32,6 @@ import iris.util import numpy as np from iris.analysis.cartography import rotate_pole, rotate_winds -from shapely.measurement import length from CSET._common import iter_maybe from CSET.operators._stash_to_lfric import STASH_TO_LFRIC @@ -168,9 +167,14 @@ def read_cubes( paths = iter_maybe(file_paths) model_names = iter_maybe(model_names) - # flattens if second element is a list - if length(model_names) == 2 and isinstance(model_names[1], list): - model_names = [model_names[0]] + model_names[1] + # flattens model_names if needed + flat = [] + for item in model_names: + if isinstance(item, list): + flat.extend(item) + else: + flat.append(item) + model_names = flat # Check we have appropriate number of model names. From 3c3eb1bed7d0b318653721dd5da5fcc8ecb2f2ff Mon Sep 17 00:00:00 2001 From: James Ross Date: Mon, 17 Aug 2026 08:54:54 +0100 Subject: [PATCH 07/22] trying to fix tests --- src/CSET/operators/read.py | 1 + src/CSET/operators/scoreswrappers.py | 14 +++++++++----- 2 files changed, 10 insertions(+), 5 deletions(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index db4c4f1cc..8dee699b0 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -168,6 +168,7 @@ def read_cubes( model_names = iter_maybe(model_names) # flattens model_names if needed + flat = [] for item in model_names: if isinstance(item, list): diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 703390c30..3e5d5242f 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -93,7 +93,11 @@ def _process_cubes_for_verification(base: Cube, other: Cube): # Set cubes into correct format using code from difference operator # Extract just common time points. - other_model_name = other.attributes["model_name"] + + if "model_name" in other.attributes: + other_model_name = other.attributes["model_name"] + else: + other_model_name = "modelname" base, other = _extract_common_time_points(base, other) # Get spatial coord names. @@ -343,7 +347,7 @@ def scores_rmse( model_name = other.attributes["model_name"] scores_cube.attributes["model_name"] = model_name - return scores_cubelist + return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = None): @@ -421,7 +425,7 @@ def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = model_name = other.attributes["model_name"] scores_cube.attributes["model_name"] = model_name - return scores_cubelist + return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist def scores_additive_bias( @@ -500,7 +504,7 @@ def scores_additive_bias( model_name = other.attributes["model_name"] scores_cube.attributes["model_name"] = model_name - return scores_cubelist + return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist def scores_correlation_pearsonr( @@ -579,7 +583,7 @@ def scores_correlation_pearsonr( scores_cubelist.append(scores_cube) model_name = other.attributes["model_name"] scores_cube.attributes["model_name"] = model_name - return scores_cubelist + return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist def scores_crps_for_ensemble( From 0c48ae2f336f9953efc26e221ccff7c1ad9ee51e Mon Sep 17 00:00:00 2001 From: James Ross Date: Mon, 17 Aug 2026 09:45:55 +0100 Subject: [PATCH 08/22] fixing tests --- src/CSET/operators/read.py | 17 ++++++++--------- 1 file changed, 8 insertions(+), 9 deletions(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index 8dee699b0..c6f68c545 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -168,17 +168,16 @@ def read_cubes( model_names = iter_maybe(model_names) # flattens model_names if needed - - flat = [] - for item in model_names: - if isinstance(item, list): - flat.extend(item) - else: - flat.append(item) - model_names = flat + if model_names != (None,): + flat = [] + for item in model_names: + if isinstance(item, list): + flat.extend(item) + else: + flat.append(item) + model_names = flat # Check we have appropriate number of model names. - if model_names != (None,) and len(model_names) != len(paths): raise ValueError( f"The number of model names ({len(model_names)}) should equal " From 4dbd70865bccb69575aee7288008e67fcd3ad958 Mon Sep 17 00:00:00 2001 From: James Ross Date: Mon, 17 Aug 2026 10:04:25 +0100 Subject: [PATCH 09/22] fixing --- src/CSET/operators/scoreswrappers.py | 97 ++++++++++++++-------------- 1 file changed, 49 insertions(+), 48 deletions(-) diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 3e5d5242f..bf63680c2 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -285,69 +285,70 @@ def scores_rmse( scores_cubelist: iris.cube.CubeList A cubelist containing the RMSE between the base and other cube. """ - base, others = _sort_cube_into_base_and_other(cubes) scores_cubelist = CubeList() - for other in others: - base, other = _process_cubes_for_verification(base, other) if obs_model_comparison: for cb in cubes: if "observed" in cb.long_name: base = cb else: - other = cb + others = cb else: - base, other = _sort_cubes_for_verification(cubes) - - # Copy the coordinates of the input cubes. - other_xr = xr.DataArray.from_iris(other) - base_xr = xr.DataArray.from_iris(base) - preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) - - # Scores operates on xarray data arrays, so we transform the iris cube into an array, - # apply scores, and then transform it back. - scores_cube = xr.DataArray.to_iris( - scores.continuous.rmse( - other_xr, - base_xr, - preserve_dims=preserve_dims, + base, others = _sort_cube_into_base_and_other(cubes) + for other in others: + base, other = _process_cubes_for_verification(base, other) + + # Copy the coordinates of the input cubes. + other_xr = xr.DataArray.from_iris(other) + base_xr = xr.DataArray.from_iris(base) + preserve_dims = _resolve_preserve_dims( + other, other_xr, preserved_coordinates ) - ) - # If time is aggregated out, attach a scalar time coordinate with bounds - # so plotting can display the aggregated period in the title. - try: - if not scores_cube.coords("time"): - base_time = base.coord("time") - time_vals = ( - base_time.bounds.flatten() - if base_time.has_bounds() - else base_time.points + # Scores operates on xarray data arrays, so we transform the iris cube into an array, + # apply scores, and then transform it back. + scores_cube = xr.DataArray.to_iris( + scores.continuous.rmse( + other_xr, + base_xr, + preserve_dims=preserve_dims, ) - t_start = float(time_vals[0]) - t_end = float(time_vals[-1]) - t_mid = 0.5 * (t_start + t_end) + ) - scores_cube.add_aux_coord( - iris.coords.AuxCoord( - t_mid, - standard_name=base_time.standard_name, - long_name=base_time.long_name, - var_name=base_time.var_name, - units=base_time.units, - bounds=np.array([t_start, t_end]), - attributes=base_time.attributes.copy(), + # If time is aggregated out, attach a scalar time coordinate with bounds + # so plotting can display the aggregated period in the title. + try: + if not scores_cube.coords("time"): + base_time = base.coord("time") + time_vals = ( + base_time.bounds.flatten() + if base_time.has_bounds() + else base_time.points ) - ) - except iris.exceptions.CoordinateNotFoundError: - pass + t_start = float(time_vals[0]) + t_end = float(time_vals[-1]) + t_mid = 0.5 * (t_start + t_end) + + scores_cube.add_aux_coord( + iris.coords.AuxCoord( + t_mid, + standard_name=base_time.standard_name, + long_name=base_time.long_name, + var_name=base_time.var_name, + units=base_time.units, + bounds=np.array([t_start, t_end]), + attributes=base_time.attributes.copy(), + ) + ) + except iris.exceptions.CoordinateNotFoundError: + pass - scores_cube.rename(f"RMSE_of_{base.name()}") - scores_cubelist.append(scores_cube) + scores_cube.rename(f"RMSE_of_{base.name()}") + scores_cubelist.append(scores_cube) - model_name = other.attributes["model_name"] - scores_cube.attributes["model_name"] = model_name + model_name = other.attributes["model_name"] + scores_cube.attributes["model_name"] = model_name - return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist + return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = None): From 900b3a1e4523de61c91d22f6d4356599f9684d68 Mon Sep 17 00:00:00 2001 From: James Ross Date: Mon, 17 Aug 2026 10:12:46 +0100 Subject: [PATCH 10/22] fix tests --- src/CSET/operators/scoreswrappers.py | 95 ++++++++++++++-------------- 1 file changed, 47 insertions(+), 48 deletions(-) diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index bf63680c2..eca7ae9e3 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -291,64 +291,63 @@ def scores_rmse( if "observed" in cb.long_name: base = cb else: - others = cb + others = [cb] else: base, others = _sort_cube_into_base_and_other(cubes) - for other in others: - base, other = _process_cubes_for_verification(base, other) - - # Copy the coordinates of the input cubes. - other_xr = xr.DataArray.from_iris(other) - base_xr = xr.DataArray.from_iris(base) - preserve_dims = _resolve_preserve_dims( - other, other_xr, preserved_coordinates + + for other in others: + base, other = _process_cubes_for_verification(base, other) + + # Copy the coordinates of the input cubes. + other_xr = xr.DataArray.from_iris(other) + base_xr = xr.DataArray.from_iris(base) + preserve_dims = _resolve_preserve_dims(other, other_xr, preserved_coordinates) + + # Scores operates on xarray data arrays, so we transform the iris cube into an array, + # apply scores, and then transform it back. + scores_cube = xr.DataArray.to_iris( + scores.continuous.rmse( + other_xr, + base_xr, + preserve_dims=preserve_dims, ) + ) - # Scores operates on xarray data arrays, so we transform the iris cube into an array, - # apply scores, and then transform it back. - scores_cube = xr.DataArray.to_iris( - scores.continuous.rmse( - other_xr, - base_xr, - preserve_dims=preserve_dims, + # If time is aggregated out, attach a scalar time coordinate with bounds + # so plotting can display the aggregated period in the title. + try: + if not scores_cube.coords("time"): + base_time = base.coord("time") + time_vals = ( + base_time.bounds.flatten() + if base_time.has_bounds() + else base_time.points ) - ) + t_start = float(time_vals[0]) + t_end = float(time_vals[-1]) + t_mid = 0.5 * (t_start + t_end) - # If time is aggregated out, attach a scalar time coordinate with bounds - # so plotting can display the aggregated period in the title. - try: - if not scores_cube.coords("time"): - base_time = base.coord("time") - time_vals = ( - base_time.bounds.flatten() - if base_time.has_bounds() - else base_time.points - ) - t_start = float(time_vals[0]) - t_end = float(time_vals[-1]) - t_mid = 0.5 * (t_start + t_end) - - scores_cube.add_aux_coord( - iris.coords.AuxCoord( - t_mid, - standard_name=base_time.standard_name, - long_name=base_time.long_name, - var_name=base_time.var_name, - units=base_time.units, - bounds=np.array([t_start, t_end]), - attributes=base_time.attributes.copy(), - ) + scores_cube.add_aux_coord( + iris.coords.AuxCoord( + t_mid, + standard_name=base_time.standard_name, + long_name=base_time.long_name, + var_name=base_time.var_name, + units=base_time.units, + bounds=np.array([t_start, t_end]), + attributes=base_time.attributes.copy(), ) - except iris.exceptions.CoordinateNotFoundError: - pass + ) + except iris.exceptions.CoordinateNotFoundError: + pass - scores_cube.rename(f"RMSE_of_{base.name()}") - scores_cubelist.append(scores_cube) + scores_cube.rename(f"RMSE_of_{base.name()}") + scores_cubelist.append(scores_cube) - model_name = other.attributes["model_name"] - scores_cube.attributes["model_name"] = model_name + model_name = other.attributes["model_name"] + scores_cube.attributes["model_name"] = model_name - return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist + return scores_cubelist[0] if len(scores_cubelist) == 1 else scores_cubelist def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = None): From 54a028c4a9f60c2dfa449d54456e787d24bff1cc Mon Sep 17 00:00:00 2001 From: James Ross Date: Mon, 17 Aug 2026 10:25:32 +0100 Subject: [PATCH 11/22] remove unneeded test --- tests/operators/test_scoreswrappers.py | 15 --------------- 1 file changed, 15 deletions(-) diff --git a/tests/operators/test_scoreswrappers.py b/tests/operators/test_scoreswrappers.py index a4430a27a..8f8863aa6 100644 --- a/tests/operators/test_scoreswrappers.py +++ b/tests/operators/test_scoreswrappers.py @@ -157,21 +157,6 @@ def test_scores_rmse_no_common_points(cube): scoreswrappers.scores_rmse(cubes) -def test_scores_rmse_incorrect_number_of_cubes(cube): - """Test exception when incorrect number of cubes provided.""" - no_cubes = CubeList([]) - with pytest.raises(ValueError, match="cubes should contain exactly 2 cubes."): - scoreswrappers.scores_rmse(no_cubes) - - one_cube = CubeList([cube]) - with pytest.raises(ValueError, match="cubes should contain exactly 2 cubes."): - scoreswrappers.scores_rmse(one_cube) - - three_cubes = CubeList([cube, cube, cube]) - with pytest.raises(ValueError, match="cubes should contain exactly 2 cubes."): - scoreswrappers.scores_rmse(three_cubes) - - def test_scores_rmse_different_data_shape_regrid(cube): """Test when data shape differs, but gets regridded. From ed3f396a11ee9e0d41ad130bf3dcfa3c21fcd56a Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 12:35:04 +0100 Subject: [PATCH 12/22] fixing loaders, fixing tests, applying suggestions --- src/CSET/loaders/verification.py | 8 ++++---- src/CSET/operators/plot.py | 1 + src/CSET/operators/scoreswrappers.py | 5 +---- tests/operators/test_scoreswrappers.py | 22 ++++++++++++++++++++++ 4 files changed, 28 insertions(+), 8 deletions(-) diff --git a/src/CSET/loaders/verification.py b/src/CSET/loaders/verification.py index 969d447da..627cc4b3c 100644 --- a/src/CSET/loaders/verification.py +++ b/src/CSET/loaders/verification.py @@ -140,22 +140,22 @@ def load(conf: Config): scores_timeseries_methods = _get_scores_timeseries_methods(conf) if scores_timeseries_methods: # Produce timeseries plots of scores metrics averaged over the domain for each case study. - for model, field, scores_method in itertools.product( - models[1:], conf.SURFACE_FIELDS, scores_timeseries_methods + for field, scores_method in itertools.product( + conf.SURFACE_FIELDS, scores_timeseries_methods ): yield RawRecipe( recipe=f"timeseries_surface_difference_scores_{scores_method}.yaml", variables={ "VARNAME": field, "BASE_MODEL": base_model["name"], - "OTHER_MODELS": model["name"], + "OTHER_MODELS": [model["name"] for model in models[1:]], "SUBAREA_NAME": conf.SUBAREA_NAME if conf.SELECT_SUBAREA else "", "SUBAREA_TYPE": conf.SUBAREA_TYPE if conf.SELECT_SUBAREA else None, "SUBAREA_EXTENT": conf.SUBAREA_EXTENT if conf.SELECT_SUBAREA else None, }, - model_ids=[base_model["id"], model["id"]], + model_ids=[base_model["id"]] + [model["id"] for model in models[1:]], aggregation=False, ) diff --git a/src/CSET/operators/plot.py b/src/CSET/operators/plot.py index 9b60ff8d9..c107df924 100644 --- a/src/CSET/operators/plot.py +++ b/src/CSET/operators/plot.py @@ -2360,6 +2360,7 @@ def plot_line_series( if nplot == 1 and seq_coord.has_bounds and np.size(seq_coord.bounds) > 1: title = f"{recipe_title}\n [{seq_coord.units.title(seq_coord.bounds[0][0])} to {seq_coord.units.title(seq_coord.bounds[0][1])}]" + print(plot_filename) # Do the actual plotting. plotting_func( cube_slice, diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index eca7ae9e3..2e2b1b80d 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -93,11 +93,8 @@ def _process_cubes_for_verification(base: Cube, other: Cube): # Set cubes into correct format using code from difference operator # Extract just common time points. + other_model_name = other.attributes["model_name"] - if "model_name" in other.attributes: - other_model_name = other.attributes["model_name"] - else: - other_model_name = "modelname" base, other = _extract_common_time_points(base, other) # Get spatial coord names. diff --git a/tests/operators/test_scoreswrappers.py b/tests/operators/test_scoreswrappers.py index 8f8863aa6..7c2661d67 100644 --- a/tests/operators/test_scoreswrappers.py +++ b/tests/operators/test_scoreswrappers.py @@ -38,6 +38,8 @@ def test_scores_correlation_pearsonr(cube: Cube): # Data preparation. other_cube = cube.copy() del other_cube.attributes["cset_comparison_base"] + cube.attributes["model_name"] = "model1" + other_cube.attributes["model_name"] = "model2" cubes = CubeList([cube, other_cube]) # Take difference. @@ -61,6 +63,8 @@ def test_scores_additive_bias(cube: Cube): # Data preparation. other_cube = cube.copy() del other_cube.attributes["cset_comparison_base"] + cube.attributes["model_name"] = "model1" + other_cube.attributes["model_name"] = "model2" cubes = CubeList([cube, other_cube]) # Take difference. @@ -80,6 +84,8 @@ def test_scores_mae(cube: Cube): # Data preparation. other_cube = cube.copy() del other_cube.attributes["cset_comparison_base"] + cube.attributes["model_name"] = "model1" + other_cube.attributes["model_name"] = "model2" cubes = CubeList([cube, other_cube]) # Take difference. @@ -97,6 +103,8 @@ def test_scores_rmse(cube: Cube): # Data preparation. other_cube = cube.copy() del other_cube.attributes["cset_comparison_base"] + cube.attributes["model_name"] = "model1" + other_cube.attributes["model_name"] = "model2" cubes = CubeList([cube, other_cube]) # Take difference. @@ -121,6 +129,8 @@ def test_scores_rmse_nonzero(): ) other_cube = cube.copy(data=np.ones((2, 2))) cube.attributes["cset_comparison_base"] = 1 + cube.attributes["model_name"] = "model1" + other_cube.attributes["model_name"] = "model2" different_cubes = CubeList((cube, other_cube)) # Take difference. rmse_cube = scoreswrappers.scores_rmse(different_cubes) @@ -138,6 +148,8 @@ def test_scores_rmse_no_time_coord(cube): c1.remove_coord("time") c2 = c1.copy() del c2.attributes["cset_comparison_base"] + c1.attributes["model_name"] = "model1" + c2.attributes["model_name"] = "model2" cubes = CubeList([c1, c2]) rmse_cube = scoreswrappers.scores_rmse(cubes) assert isinstance(rmse_cube, Cube) @@ -152,6 +164,8 @@ def test_scores_rmse_no_common_points(cube): new_times += 6 other_cube.coord("time").points = new_times del other_cube.attributes["cset_comparison_base"] + cube.attributes["model_name"] = "model1" + other_cube.attributes["model_name"] = "model2" cubes = CubeList([cube, other_cube]) with pytest.raises(ValueError, match="No common time points found!"): scoreswrappers.scores_rmse(cubes) @@ -165,6 +179,8 @@ def test_scores_rmse_different_data_shape_regrid(cube): rearranged_cube = cube.copy() rearranged_cube = rearranged_cube[:, :, 1:] del cube.attributes["cset_comparison_base"] + rearranged_cube.attributes["model_name"] = "model1" + cube.attributes["model_name"] = "model2" cubes = CubeList([rearranged_cube, cube]) # Need to preserve coordinates to test shape. rmse = scoreswrappers.scores_rmse( @@ -180,6 +196,8 @@ def test_rmse_grid_staggering_regrid(cube): rearranged_cube = cube.copy() rearranged_cube.rename("eastward_wind_at_10m") del cube.attributes["cset_comparison_base"] + rearranged_cube.attributes["model_name"] = "model1" + cube.attributes["model_name"] = "model2" cubes = CubeList([rearranged_cube, cube]) # Need to preserve coordinates to test shape. rmse = scoreswrappers.scores_rmse( @@ -194,6 +212,8 @@ def test_difference_different_model_types(cube): flipped = cube.copy() reverse(flipped, "grid_latitude") del flipped.attributes["cset_comparison_base"] + flipped.attributes["model_name"] = "model1" + cube.attributes["model_name"] = "model2" cubes = CubeList([cube, flipped]) # Take rmse. @@ -209,6 +229,8 @@ def test_difference_flip_pressure_order(transect_source_cube_readonly): flipped = transect_source_cube_readonly.copy() reverse(flipped, "pressure") del flipped.attributes["cset_comparison_base"] + flipped.attributes["model_name"] = "model1" + transect_source_cube_readonly.attributes["model_name"] = "model2" cubes = CubeList([transect_source_cube_readonly, flipped]) # Take rmse. From 7fc466b4601dddf469dd120ca195387696286f81 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 13:36:37 +0100 Subject: [PATCH 13/22] adding doc strings --- src/CSET/operators/scoreswrappers.py | 27 +++++++++++++++++++++++++++ 1 file changed, 27 insertions(+) diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 2e2b1b80d..1ec21686e 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -42,6 +42,21 @@ def _sort_cube_into_base_and_other(cubes): + """Sorts cube into base and other models. + + Parameters + ---------- + cubes: iris.cube.CubeList + A CubeList of multiple cubes. One base cube and other model cubes. + + Returns + ------- + base: iris.cube.Cube + The cube from the "analysis" in the same format as the other model. + others: iris.cube.CubeList + The cube list of containing the cube(s) from the model in the same format as the base model. + + """ base: Cube = cubes.extract_cube(iris.AttributeConstraint(cset_comparison_base=1)) others: CubeList = cubes.extract( iris.Constraint( @@ -53,6 +68,18 @@ def _sort_cube_into_base_and_other(cubes): def _ensure_increasing_pressure_coordinates(cubes): + """Ensure the pressure coordinate is increasing. + + Parameters + ---------- + cubes: iris.cube.CubeList + A CubeList of n cubes + + Returns + ------- + Cubes: iris.cube.CubeList + The original cube list but where each cube is ensured to have an increasing pressure coordinate. + """ for cube in cubes: try: if len(cube.coord("pressure").points) > 2 and not is_increasing( From dee04ff336b5c14bffd6bccb25865a3096a15e1a Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 13:52:41 +0100 Subject: [PATCH 14/22] change recipe title --- src/CSET/operators/scoreswrappers.py | 12 ++++++------ .../timeseries_surface_difference_scores_MAE.yaml | 2 +- .../timeseries_surface_difference_scores_RMSE.yaml | 2 +- ...ries_surface_difference_scores_additive_bias.yaml | 2 +- ...rface_difference_scores_correlation_pearsonr.yaml | 2 +- 5 files changed, 10 insertions(+), 10 deletions(-) diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 1ec21686e..6052e702a 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -254,8 +254,8 @@ def scores_rmse_model_obs( Returns ------- - scores_cube: iris.cube.Cube - A cube containing the RMSE between the models and observation cube. + scores_cubelist: iris.cube.CubeList + A cubelist containing the RMSE between the models and observation cube(s). """ rmse_cubes = CubeList() model_list = CubeList() @@ -394,7 +394,7 @@ def scores_mae(cubes: CubeList, preserved_coordinates: list[str] | str | None = Returns ------- scores_cubelist: iris.cube.CubeList - A cubelist containing the MAE between the base and other cubes. + A cubelist containing the MAE between the base and other cube(s). """ base, others = _sort_cube_into_base_and_other(cubes) scores_cubelist = CubeList() @@ -474,7 +474,7 @@ def scores_additive_bias( Returns ------- scores_cubelist: iris.cube.CubeList - A cubelist containing the ME between the base and other cube. + A cubelist containing the ME between the base and other cube(s). """ base, others = _sort_cube_into_base_and_other(cubes) scores_cubelist = CubeList() @@ -552,8 +552,8 @@ def scores_correlation_pearsonr( Returns ------- - scores_cube: iris.cube.Cube - A cube containing the PC between the base and other cube. + scores_cubelist: iris.cube.CubeList + A cubelist containing the PC between the base and other cube(s). """ base, others = _sort_cube_into_base_and_other(cubes) scores_cubelist = CubeList() diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml index f5d36359d..dc25e8176 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Absolute Error timeseries between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Absolute Error timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Absolute Error in $VARNAME computed over the domain for each timestep. The MAE is calculated based on that used in the diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml index 38eecfc6e..0306fee24 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores RMSE timeseries between $OTHER_MODELS and $BASE_MODEL and $SUBAREA_NAME" +title: "$VARNAME\nScores RMSE timeseries between models and $BASE_MODEL and $SUBAREA_NAME" description: | Extracts and plots the Root Mean Square Error in $VARNAME computed over the domain for each timestep. The RMSE is calculated based on that used in the diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml index dbc60032c..cf0d6a911 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Error timeseries between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Error timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Error (also known as the Additive Bias) in $VARNAME computed over the domain for each timestep. The ME is calculated based on that used in the diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml index 38597ec7f..539fe0cbe 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Pearson's Correlation timeseries between $OTHER_MODELS and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Pearson's Correlation timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Pearson's Correlation coefficient in $VARNAME computed over the domain for each timestep. The PC is calculated based on that used in the From a2aeb6ca936cd6f979c0b167dd5558922ba8bce3 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 13:53:53 +0100 Subject: [PATCH 15/22] Update src/CSET/operators/read.py Co-authored-by: James Warner <62252918+jwarner8@users.noreply.github.com> --- src/CSET/operators/read.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/operators/read.py b/src/CSET/operators/read.py index c6f68c545..1a9f49ea7 100644 --- a/src/CSET/operators/read.py +++ b/src/CSET/operators/read.py @@ -167,7 +167,7 @@ def read_cubes( paths = iter_maybe(file_paths) model_names = iter_maybe(model_names) - # flattens model_names if needed + # flattens model_names if needed into one dimensional list. if model_names != (None,): flat = [] for item in model_names: From 902700387810c9895b59ee692eb6f148224534c3 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 14:05:18 +0100 Subject: [PATCH 16/22] fixing doc string --- src/CSET/operators/scoreswrappers.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 6052e702a..499de1695 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -303,6 +303,8 @@ def scores_rmse( ["time","grid_latitude", "grid_longitude"] or if you want a time series you can preserve ["time"], if you want to collapse to a single value use `None`. The default is `None`. + obs_model_comparison: bool, default False + Set true if doing model-obs comparison. Returns ------- From f76811b4cd5ef898ad12e818f003d4eb86835d35 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 15:31:55 +0100 Subject: [PATCH 17/22] Update src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml Co-authored-by: James Warner <62252918+jwarner8@users.noreply.github.com> --- ...meseries_surface_difference_scores_correlation_pearsonr.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml index 539fe0cbe..2eb849e4f 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Pearson's Correlation timeseries between models and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Pearson's Correlation Difference timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Pearson's Correlation coefficient in $VARNAME computed over the domain for each timestep. The PC is calculated based on that used in the From 03cbacc013211294b2d84beca1795d5a52cfa886 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 15:32:37 +0100 Subject: [PATCH 18/22] Update src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml Co-authored-by: James Warner <62252918+jwarner8@users.noreply.github.com> --- .../timeseries_surface_difference_scores_additive_bias.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml index cf0d6a911..ef30e3941 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Error timeseries between models and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Error Difference timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Error (also known as the Additive Bias) in $VARNAME computed over the domain for each timestep. The ME is calculated based on that used in the From e487ff37967601fcc118100f80275d7fc8f69bd6 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 15:32:46 +0100 Subject: [PATCH 19/22] Update src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml Co-authored-by: James Warner <62252918+jwarner8@users.noreply.github.com> --- .../verification/timeseries_surface_difference_scores_RMSE.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml index 0306fee24..b3df030ad 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores RMSE timeseries between models and $BASE_MODEL and $SUBAREA_NAME" +title: "$VARNAME\nScores RMSE Difference timeseries between models and $BASE_MODEL and $SUBAREA_NAME" description: | Extracts and plots the Root Mean Square Error in $VARNAME computed over the domain for each timestep. The RMSE is calculated based on that used in the From 9b90969d1a5e71b25c9e981730d991f46473d992 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 15:33:59 +0100 Subject: [PATCH 20/22] Update src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml Co-authored-by: James Warner <62252918+jwarner8@users.noreply.github.com> --- .../verification/timeseries_surface_difference_scores_MAE.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml index dc25e8176..c6f8d4047 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Absolute Error timeseries between models and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Absolute Error Difference timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Absolute Error in $VARNAME computed over the domain for each timestep. The MAE is calculated based on that used in the From 58580a70bcb0e1d5b34d9ec9b97fb439dd9d0bf0 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 15:35:04 +0100 Subject: [PATCH 21/22] remove print --- src/CSET/operators/plot.py | 1 - 1 file changed, 1 deletion(-) diff --git a/src/CSET/operators/plot.py b/src/CSET/operators/plot.py index c107df924..9b60ff8d9 100644 --- a/src/CSET/operators/plot.py +++ b/src/CSET/operators/plot.py @@ -2360,7 +2360,6 @@ def plot_line_series( if nplot == 1 and seq_coord.has_bounds and np.size(seq_coord.bounds) > 1: title = f"{recipe_title}\n [{seq_coord.units.title(seq_coord.bounds[0][0])} to {seq_coord.units.title(seq_coord.bounds[0][1])}]" - print(plot_filename) # Do the actual plotting. plotting_func( cube_slice, From 7bbbebe7f25d01eb759512eef389aa25ae76d425 Mon Sep 17 00:00:00 2001 From: James Ross Date: Thu, 20 Aug 2026 15:41:44 +0100 Subject: [PATCH 22/22] revert title changes --- .../verification/timeseries_surface_difference_scores_MAE.yaml | 2 +- .../verification/timeseries_surface_difference_scores_RMSE.yaml | 2 +- .../timeseries_surface_difference_scores_additive_bias.yaml | 2 +- ...meseries_surface_difference_scores_correlation_pearsonr.yaml | 2 +- 4 files changed, 4 insertions(+), 4 deletions(-) diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml index c6f8d4047..dc25e8176 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_MAE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Absolute Error Difference timeseries between models and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Absolute Error timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Absolute Error in $VARNAME computed over the domain for each timestep. The MAE is calculated based on that used in the diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml index b3df030ad..0306fee24 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_RMSE.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores RMSE Difference timeseries between models and $BASE_MODEL and $SUBAREA_NAME" +title: "$VARNAME\nScores RMSE timeseries between models and $BASE_MODEL and $SUBAREA_NAME" description: | Extracts and plots the Root Mean Square Error in $VARNAME computed over the domain for each timestep. The RMSE is calculated based on that used in the diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml index ef30e3941..cf0d6a911 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_additive_bias.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Mean Error Difference timeseries between models and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Mean Error timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Mean Error (also known as the Additive Bias) in $VARNAME computed over the domain for each timestep. The ME is calculated based on that used in the diff --git a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml index 2eb849e4f..539fe0cbe 100644 --- a/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml +++ b/src/CSET/recipes/verification/timeseries_surface_difference_scores_correlation_pearsonr.yaml @@ -1,5 +1,5 @@ category: Scores -title: "$VARNAME\nScores Pearson's Correlation Difference timeseries between models and $BASE_MODEL $SUBAREA_NAME" +title: "$VARNAME\nScores Pearson's Correlation timeseries between models and $BASE_MODEL $SUBAREA_NAME" description: | Extracts and plots the Pearson's Correlation coefficient in $VARNAME computed over the domain for each timestep. The PC is calculated based on that used in the