diff --git a/src/CSET/cset_workflow/includes/metplus_point_stat.cylc b/src/CSET/cset_workflow/includes/metplus_point_stat.cylc index 00a01ab97..8880908e4 100644 --- a/src/CSET/cset_workflow/includes/metplus_point_stat.cylc +++ b/src/CSET/cset_workflow/includes/metplus_point_stat.cylc @@ -80,7 +80,7 @@ # Runs VerPy metloader utility to create the VerPy databases inherit = VERPY_METLOADER [[[environment]]] - VERPY_DIR = {{VERPY_DIR}} + VERPY_DIR = {{VERPY_DIR|default("")}} VER_METHOD = area STAT_TYPE = cnt STREAM = point_stat diff --git a/src/CSET/cset_workflow/meta/verification/rose-meta.conf b/src/CSET/cset_workflow/meta/verification/rose-meta.conf index ec57b3ad6..0b5b037d3 100644 --- a/src/CSET/cset_workflow/meta/verification/rose-meta.conf +++ b/src/CSET/cset_workflow/meta/verification/rose-meta.conf @@ -77,7 +77,7 @@ trigger=template variables=SCORES_SPATIAL_RMSE: False; template variables=SCORES_TIMESERIES_AB: False; template variables=SCORES_TIMESERIES_PC: False; compulsory=true -sort-key=scores1 +sort-key=scoresspatial1 # Scores spatial plots @@ -226,3 +226,31 @@ help=Usually this should be set to 0, but in some case it may be 1. type=integer compulsory=true sort-key=scores4b + +[Verification/Scores vertical profiles] +ns=Verification/Scores vertical profiles +sort-key=scoresvertical +title=Scores vertical profiles + +[Verification/Scores vertical profiles] +[template variables=SCORES_RMSE_VERTICAL_PROFILES] +ns=Verification/Scores vertical profiles +description=Create vertical profile RMSE plots for the specified level fields. + RMSE profiles are calculated at each grid point for all timesteps within a case study at each pressure level. +type=python_boolean +compulsory=true +sort-key=scoresprofile1 + +[template variables=SCORES_RMSE_VERTICAL_PROFILES_SEQUENCE] +ns=Verification/Scores vertical profiles +description=Also produce one RMSE vertical profile per time step. +help=When True, the RMSE is calculated at each grid point and level for every + time step, collapsing only the spatial (horizontal) dimensions. This produces + a sequence of vertical RMSE profile plots, one per time step, in addition to + the profile that is aggregated across the whole case study. + + When False (default), only the aggregated vertical RMSE profile is produced. +type=python_boolean +compulsory=true +#trigger=template variables=SCORES_RMSE_VERTICAL_PROFILES: True; +sort-key=scoresprofile2 diff --git a/src/CSET/cset_workflow/rose-suite.conf.example b/src/CSET/cset_workflow/rose-suite.conf.example index a43f7a5f6..636747ce3 100644 --- a/src/CSET/cset_workflow/rose-suite.conf.example +++ b/src/CSET/cset_workflow/rose-suite.conf.example @@ -150,6 +150,8 @@ SCORES_ALL=False SCORES_CRPS_FOR_ENSEMBLE=False SCORES_RMSE_SPATIAL=False SCORES_RMSE_TIMESERIES=False +SCORES_RMSE_VERTICAL_PROFILES_SEQUENCE=False +SCORES_RMSE_VERTICAL_PROFILES=False SCORES_SPATIAL_AB=False SCORES_SPATIAL_MAE=False SCORES_SPATIAL_RMSE=False diff --git a/src/CSET/cset_workflow/site/nci-gadi.cylc b/src/CSET/cset_workflow/site/nci-gadi.cylc index 543ef4ea5..3e83a92f0 100644 --- a/src/CSET/cset_workflow/site/nci-gadi.cylc +++ b/src/CSET/cset_workflow/site/nci-gadi.cylc @@ -23,7 +23,7 @@ """ [[[ environment ]]] PROJECT = {{ PROJECT }} - PYTHONPATH = "{{VERPY_DIR}}" + PYTHONPATH = "{{VERPY_DIR|default("")}}" {% if RUN_METPLUS_GRID_STAT|default(False) or RUN_METPLUS_POINT_STAT|default(False) %} [[METPLUS]] diff --git a/src/CSET/loaders/verification.py b/src/CSET/loaders/verification.py index 18623c243..7b99e75ec 100644 --- a/src/CSET/loaders/verification.py +++ b/src/CSET/loaders/verification.py @@ -49,12 +49,14 @@ def load(conf: Config): """Yield recipes from the given workflow configuration.""" # Load a list of model detail dictionaries. models = get_models(conf.asdict()) + if not models: + return # Models are listed in order, so model 1 is the first element. + base_model = models[0] scores_spatial_methods = _get_scores_spatial_methods(conf) if scores_spatial_methods: # Produce 2D spatial plots of scores metrics. - base_model = models[0] for model, field, method, scores_method in itertools.product( models[1:], conf.SURFACE_FIELDS, @@ -101,7 +103,6 @@ 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. - base_model = models[0] for model, field, scores_method in itertools.product( models[1:], conf.SURFACE_FIELDS, scores_timeseries_methods ): @@ -137,3 +138,43 @@ def load(conf: Config): model_ids=[model["id"]], aggregation=False, ) + + # including AGGREGATION_MODE in the variables dictionary to allow for clearer labeling of plots. + if conf.SCORES_RMSE_VERTICAL_PROFILES: + for model, field in itertools.product(models[1:], conf.PRESSURE_LEVEL_FIELDS): + yield RawRecipe( + recipe="generic_level_rmse_scores_profile.yaml", + variables={ + "VARNAME": field, + "BASE_MODEL": base_model["name"], + "OTHER_MODEL": model["name"], + "PRESERVED_COORDS": ["pressure"], + "AGGREGATION_MODE": "Case-study RMSE", + "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"]], + aggregation=False, + ) + + # including AGGREGATION_MODE in the variables dictionary to allow for clearer labeling of plots. + if conf.SCORES_RMSE_VERTICAL_PROFILES_SEQUENCE: + for model, field in itertools.product(models[1:], conf.PRESSURE_LEVEL_FIELDS): + yield RawRecipe( + recipe="generic_level_rmse_scores_profile.yaml", + variables={ + "VARNAME": field, + "BASE_MODEL": base_model["name"], + "OTHER_MODEL": model["name"], + "PRESERVED_COORDS": ["time", "pressure"], + "AGGREGATION_MODE": "Time-step RMSE", + "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"]], + aggregation=False, + ) diff --git a/src/CSET/operators/plot.py b/src/CSET/operators/plot.py index a4740a4c1..d294f3ab5 100644 --- a/src/CSET/operators/plot.py +++ b/src/CSET/operators/plot.py @@ -2469,28 +2469,74 @@ def plot_vertical_line_series( vmin = min(x_levels) vmax = max(x_levels) - # Matching the slices (matching by seq coord point; it may happen that - # evaluated models do not cover the same seq coord range, hence matching - # necessary) - cube_iterables = _find_matched_slices(cubes, sequence_coordinate) + # Check if the cube has a sequence coordinate (e.g. time). If not, plot + # a single profile directly without iterating over a sequence. + sequence_coords = [ + cube.coord(sequence_coordinate) + for cube in cubes + if cube.coords(sequence_coordinate) + ] + has_sequence_coord = len(sequence_coords) == len(cubes) and all( + np.size(coord.points) > 1 for coord in sequence_coords + ) + has_scalar_sequence_coord = len(sequence_coords) == len(cubes) and all( + np.size(coord.points) == 1 for coord in sequence_coords + ) - # Create a plot for each value of the sequence coordinate. - # Allowing for multiple cubes in a CubeList to be plotted in the same plot for - # similar sequence values. Passing a CubeList into the internal plotting function - # for similar values of the sequence coordinate. cube_slice can be an iris.cube.Cube - # or an iris.cube.CubeList. plot_index = [] - nplot = np.size(cubes[0].coord(sequence_coordinate).points) - for cubes_slice in cube_iterables: - # Format the coordinate value in a unit appropriate way. - seq_coord = cubes_slice[0].coord(sequence_coordinate) + if has_sequence_coord: + # Matching the slices (matching by seq coord point; it may happen that + # evaluated models do not cover the same seq coord range, hence matching + # necessary) + cube_iterables = _find_matched_slices(cubes, sequence_coordinate) + nplot = np.size(cubes[0].coord(sequence_coordinate).points) + for cubes_slice in cube_iterables: + # Format the coordinate value in a unit appropriate way. + seq_coord = cubes_slice[0].coord(sequence_coordinate) + plot_title, plot_filename = _set_title_and_filename( + seq_coord, nplot, recipe_title, filename + ) + + # Do the actual plotting. + _plot_and_save_vertical_line_series( + cubes_slice, + coords, + "realization", + plot_filename, + series_coordinate, + title=plot_title, + vmin=vmin, + vmax=vmax, + ) + plot_index.append(plot_filename) + elif has_scalar_sequence_coord: + # Scalar sequence coordinate (typically aggregated time bounds): + # make one plot and include sequence period in title/filename. plot_title, plot_filename = _set_title_and_filename( - seq_coord, nplot, recipe_title, filename + sequence_coords[0], 1, recipe_title, filename ) - # Do the actual plotting. _plot_and_save_vertical_line_series( - cubes_slice, + cubes, + coords, + "realization", + plot_filename, + series_coordinate, + title=plot_title, + vmin=vmin, + vmax=vmax, + ) + plot_index.append(plot_filename) + else: + # 1D case: no sequence coordinate, plot a single profile. + plot_title = recipe_title + if filename: + plot_filename = filename + else: + plot_filename = f"{slugify(plot_title)}.png" + + _plot_and_save_vertical_line_series( + cubes, coords, "realization", plot_filename, diff --git a/src/CSET/operators/scoreswrappers.py b/src/CSET/operators/scoreswrappers.py index 977380250..bee585144 100644 --- a/src/CSET/operators/scoreswrappers.py +++ b/src/CSET/operators/scoreswrappers.py @@ -208,7 +208,7 @@ def scores_rmse(cubes: CubeList, preserved_coordinates: list[str] | str | None = A CubeList containing exactly two cubes: a base and an "other" model, this can be an analysis and the model. preserved_coordinates: list[str] | str | None, default is None. - The coordinates that you wish to preserve in the calculaiton of the + The coordinates (or xarray dimension names) that you wish to preserve in the calculaiton of the RMSE. For example if you want a map of each time you can preserve ["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 diff --git a/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml b/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml new file mode 100644 index 000000000..d4b2788d7 --- /dev/null +++ b/src/CSET/recipes/verification/generic_level_rmse_scores_profile.yaml @@ -0,0 +1,48 @@ +category: Scores +title: "$VARNAME\nRMSE vertical profile between $OTHER_MODEL 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 + pairwise between $BASE_MODEL and $OTHER_MODEL for each pressure + level field utilising the preserved coordinates in the scores package. + + By adding a preserved coordinate to the recipe, the RMSE can be calculated for each + time step in a case study or for all time steps combined. + + References: + https://scores.readthedocs.io/en/stable/ + https://scores.readthedocs.io/en/stable/api.html#scores.continuous.rmse + + + Aggregation mode: + - Case-study RMSE: Single RMSE profile calculated using all time points in the case study. + - Time-step RMSE: RMSE profiles calculated at each time point in a case study. + + A larger RMSE implies a greater error than a smaller RMSE. An + RMSE of zero implies the two fields match exactly. + +steps: + - operator: read.read_cubes + file_paths: $INPUT_PATHS + model_names: [$BASE_MODEL, $OTHER_MODEL] + constraint: + operator: constraints.combine_constraints + variable_constraint: + operator: constraints.generate_var_constraint + varname: $VARNAME + level_constraint: + operator: constraints.generate_level_constraint + coordinate: pressure + levels: "*" + subarea_type: $SUBAREA_TYPE + subarea_extent: $SUBAREA_EXTENT + + - operator: scoreswrappers.scores_rmse + preserved_coordinates: $PRESERVED_COORDS + + - operator: plot.plot_vertical_line_series + series_coordinate: pressure + sequence_coordinate: time + + - operator: write.write_cube_to_nc + overwrite: True