Coverage for src/CSET/operators/plot.py: 88%
765 statements
« prev ^ index » next coverage.py v7.11.0, created at 2025-10-30 15:17 +0000
1# © Crown copyright, Met Office (2022-2025) and CSET contributors.
2#
3# Licensed under the Apache License, Version 2.0 (the "License");
4# you may not use this file except in compliance with the License.
5# You may obtain a copy of the License at
6#
7# http://www.apache.org/licenses/LICENSE-2.0
8#
9# Unless required by applicable law or agreed to in writing, software
10# distributed under the License is distributed on an "AS IS" BASIS,
11# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12# See the License for the specific language governing permissions and
13# limitations under the License.
15"""Operators to produce various kinds of plots."""
17import fcntl
18import functools
19import importlib.resources
20import itertools
21import json
22import logging
23import math
24import os
25import sys
26from typing import Literal
28import cartopy.crs as ccrs
29import iris
30import iris.coords
31import iris.cube
32import iris.exceptions
33import iris.plot as iplt
34import matplotlib as mpl
35import matplotlib.colors as mcolors
36import matplotlib.pyplot as plt
37import numpy as np
38from markdown_it import MarkdownIt
40from CSET._common import (
41 combine_dicts,
42 get_recipe_metadata,
43 iter_maybe,
44 render_file,
45 slugify,
46)
47from CSET.operators._utils import get_cube_yxcoordname, is_transect
49# Use a non-interactive plotting backend.
50mpl.use("agg")
52DEFAULT_DISCRETE_COLORS = mpl.colormaps["tab10"].colors + mpl.colormaps["Accent"].colors
54############################
55# Private helper functions #
56############################
59def _append_to_plot_index(plot_index: list) -> list:
60 """Add plots into the plot index, returning the complete plot index."""
61 with open("meta.json", "r+t", encoding="UTF-8") as fp:
62 fcntl.flock(fp, fcntl.LOCK_EX)
63 fp.seek(0)
64 meta = json.load(fp)
65 complete_plot_index = meta.get("plots", [])
66 complete_plot_index = complete_plot_index + plot_index
67 meta["plots"] = complete_plot_index
68 if os.getenv("CYLC_TASK_CYCLE_POINT") and not bool(
69 os.getenv("DO_CASE_AGGREGATION")
70 ):
71 meta["case_date"] = os.getenv("CYLC_TASK_CYCLE_POINT", "")
72 fp.seek(0)
73 fp.truncate()
74 json.dump(meta, fp, indent=2)
75 return complete_plot_index
78def _check_single_cube(cube: iris.cube.Cube | iris.cube.CubeList) -> iris.cube.Cube:
79 """Ensure a single cube is given.
81 If a CubeList of length one is given that the contained cube is returned,
82 otherwise an error is raised.
84 Parameters
85 ----------
86 cube: Cube | CubeList
87 The cube to check.
89 Returns
90 -------
91 cube: Cube
92 The checked cube.
94 Raises
95 ------
96 TypeError
97 If the input cube is not a Cube or CubeList of a single Cube.
98 """
99 if isinstance(cube, iris.cube.Cube):
100 return cube
101 if isinstance(cube, iris.cube.CubeList):
102 if len(cube) == 1:
103 return cube[0]
104 raise TypeError("Must have a single cube", cube)
107def _py312_importlib_resources_files_shim():
108 """Importlib behaviour changed in 3.12 to avoid circular dependencies.
110 This shim is needed until python 3.12 is our oldest supported version, after
111 which it can just be replaced by directly using importlib.resources.files.
112 """
113 if sys.version_info.minor >= 12:
114 files = importlib.resources.files()
115 else:
116 import CSET.operators
118 files = importlib.resources.files(CSET.operators)
119 return files
122def _make_plot_html_page(plots: list):
123 """Create a HTML page to display a plot image."""
124 # Debug check that plots actually contains some strings.
125 assert isinstance(plots[0], str)
127 # Load HTML template file.
128 operator_files = _py312_importlib_resources_files_shim()
129 template_file = operator_files.joinpath("_plot_page_template.html")
131 # Get some metadata.
132 meta = get_recipe_metadata()
133 title = meta.get("title", "Untitled")
134 description = MarkdownIt().render(meta.get("description", "*No description.*"))
136 # Prepare template variables.
137 variables = {
138 "title": title,
139 "description": description,
140 "initial_plot": plots[0],
141 "plots": plots,
142 "title_slug": slugify(title),
143 }
145 # Render template.
146 html = render_file(template_file, **variables)
148 # Save completed HTML.
149 with open("index.html", "wt", encoding="UTF-8") as fp:
150 fp.write(html)
153@functools.cache
154def _load_colorbar_map(user_colorbar_file: str = None) -> dict:
155 """Load the colorbar definitions from a file.
157 This is a separate function to make it cacheable.
158 """
159 colorbar_file = _py312_importlib_resources_files_shim().joinpath(
160 "_colorbar_definition.json"
161 )
162 with open(colorbar_file, "rt", encoding="UTF-8") as fp:
163 colorbar = json.load(fp)
165 logging.debug("User colour bar file: %s", user_colorbar_file)
166 override_colorbar = {}
167 if user_colorbar_file:
168 try:
169 with open(user_colorbar_file, "rt", encoding="UTF-8") as fp:
170 override_colorbar = json.load(fp)
171 except FileNotFoundError:
172 logging.warning("Colorbar file does not exist. Using default values.")
174 # Overwrite values with the user supplied colorbar definition.
175 colorbar = combine_dicts(colorbar, override_colorbar)
176 return colorbar
179def _get_model_colors_map(cubes: iris.cube.CubeList | iris.cube.Cube) -> dict:
180 """Get an appropriate colors for model lines in line plots.
182 For each model in the list of cubes colors either from user provided
183 color definition file (so-called style file) or from default colors are mapped
184 to model_name attribute.
186 Parameters
187 ----------
188 cubes: CubeList or Cube
189 Cubes with model_name attribute
191 Returns
192 -------
193 model_colors_map:
194 Dictionary mapping model_name attribute to colors
195 """
196 user_colorbar_file = get_recipe_metadata().get("style_file_path", None)
197 colorbar = _load_colorbar_map(user_colorbar_file)
198 model_names = sorted(
199 filter(
200 lambda x: x is not None,
201 (cube.attributes.get("model_name", None) for cube in iter_maybe(cubes)),
202 )
203 )
204 if not model_names:
205 return {}
206 use_user_colors = all(mname in colorbar.keys() for mname in model_names)
207 if use_user_colors: 207 ↛ 208line 207 didn't jump to line 208 because the condition on line 207 was never true
208 return {mname: colorbar[mname] for mname in model_names}
210 color_list = itertools.cycle(DEFAULT_DISCRETE_COLORS)
211 return {mname: color for mname, color in zip(model_names, color_list, strict=False)}
214def _colorbar_map_levels(cube: iris.cube.Cube, axis: Literal["x", "y"] | None = None):
215 """Get an appropriate colorbar for the given cube.
217 For the given variable the appropriate colorbar is looked up from a
218 combination of the built-in CSET colorbar definitions, and any user supplied
219 definitions. As well as varying on variables, these definitions may also
220 exist for specific pressure levels to account for variables with
221 significantly different ranges at different heights. The colorbars also exist
222 for masks and mask differences for considering variable presence diagnostics.
223 Specific variable ranges can be separately set in user-supplied definition
224 for x- or y-axis limits, or indicate where automated range preferred.
226 Parameters
227 ----------
228 cube: Cube
229 Cube of variable for which the colorbar information is desired.
230 axis: "x", "y", optional
231 Select the levels for just this axis of a line plot. The min and max
232 can be set by xmin/xmax or ymin/ymax respectively. For variables where
233 setting a universal range is not desirable (e.g. temperature), users
234 can set ymin/ymax values to "auto" in the colorbar definitions file.
235 Where no additional xmin/xmax or ymin/ymax values are provided, the
236 axis bounds default to use the vmin/vmax values provided.
238 Returns
239 -------
240 cmap:
241 Matplotlib colormap.
242 levels:
243 List of levels to use for plotting. For continuous plots the min and max
244 should be taken as the range.
245 norm:
246 BoundaryNorm information.
247 """
248 # Grab the colorbar file from the recipe global metadata.
249 user_colorbar_file = get_recipe_metadata().get("style_file_path", None)
250 colorbar = _load_colorbar_map(user_colorbar_file)
251 cmap = None
253 try:
254 # We assume that pressure is a scalar coordinate here.
255 pressure_level_raw = cube.coord("pressure").points[0]
256 # Ensure pressure_level is a string, as it is used as a JSON key.
257 pressure_level = str(int(pressure_level_raw))
258 except iris.exceptions.CoordinateNotFoundError:
259 pressure_level = None
261 # First try long name, then standard name, then var name. This order is used
262 # as long name is the one we correct between models, so it most likely to be
263 # consistent.
264 varnames = list(filter(None, [cube.long_name, cube.standard_name, cube.var_name]))
265 for varname in varnames:
266 # Get the colormap for this variable.
267 try:
268 var_colorbar = colorbar[varname]
269 cmap = plt.get_cmap(colorbar[varname]["cmap"], 51)
270 varname_key = varname
271 break
272 except KeyError:
273 logging.debug("Cube name %s has no colorbar definition.", varname)
275 # Get colormap if it is a mask.
276 if any("mask_for_" in name for name in varnames):
277 cmap, levels, norm = _custom_colormap_mask(cube, axis=axis)
278 return cmap, levels, norm
279 # If winds on Beaufort Scale use custom colorbar and levels
280 if any("Beaufort_Scale" in name for name in varnames):
281 cmap, levels, norm = _custom_beaufort_scale(cube, axis=axis)
282 return cmap, levels, norm
283 # If probability is plotted use custom colorbar and levels
284 if any("probability_of_" in name for name in varnames):
285 cmap, levels, norm = _custom_colormap_probability(cube, axis=axis)
286 return cmap, levels, norm
287 # If aviation colour state use custom colorbar and levels
288 if any("aviation_colour_state" in name for name in varnames): 288 ↛ 289line 288 didn't jump to line 289 because the condition on line 288 was never true
289 cmap, levels, norm = _custom_colormap_aviation_colour_state(cube)
290 return cmap, levels, norm
292 # If no valid colormap has been defined, use defaults and return.
293 if not cmap:
294 logging.warning("No colorbar definition exists for %s.", cube.name())
295 cmap, levels, norm = mpl.colormaps["viridis"], None, None
296 return cmap, levels, norm
298 # Test if pressure-level specific settings are provided for cube.
299 if pressure_level:
300 try:
301 var_colorbar = colorbar[varname_key]["pressure_levels"][pressure_level]
302 except KeyError:
303 logging.debug(
304 "%s has no colorbar definition for pressure level %s.",
305 varname,
306 pressure_level,
307 )
309 # Check for availability of x-axis or y-axis user-specific overrides
310 # for setting level bounds for line plot types and return just levels.
311 # Line plots do not need a colormap, and just use the data range.
312 if axis:
313 if axis == "x":
314 try:
315 vmin, vmax = var_colorbar["xmin"], var_colorbar["xmax"]
316 except KeyError:
317 vmin, vmax = var_colorbar["min"], var_colorbar["max"]
318 if axis == "y":
319 try:
320 vmin, vmax = var_colorbar["ymin"], var_colorbar["ymax"]
321 except KeyError:
322 vmin, vmax = var_colorbar["min"], var_colorbar["max"]
323 # Check if user-specified auto-scaling for this variable
324 if vmin == "auto" or vmax == "auto":
325 levels = None
326 else:
327 levels = [vmin, vmax]
328 return None, levels, None
329 # Get and use the colorbar levels for this variable if spatial or histogram.
330 else:
331 try:
332 levels = var_colorbar["levels"]
333 # Use discrete bins when levels are specified, rather
334 # than a smooth range.
335 norm = mpl.colors.BoundaryNorm(levels, ncolors=cmap.N)
336 logging.debug("Using levels for %s colorbar.", varname)
337 logging.info("Using levels: %s", levels)
338 except KeyError:
339 # Get the range for this variable.
340 vmin, vmax = var_colorbar["min"], var_colorbar["max"]
341 logging.debug("Using min and max for %s colorbar.", varname)
342 # Calculate levels from range.
343 levels = np.linspace(vmin, vmax, 101)
344 norm = None
346 # Overwrite cmap, levels and norm for specific variables that
347 # require custom colorbar_map as these can not be defined in the
348 # JSON file.
349 cmap, levels, norm = _custom_colourmap_precipitation(cube, cmap, levels, norm)
350 cmap, levels, norm = _custom_colourmap_visibility_in_air(
351 cube, cmap, levels, norm
352 )
353 cmap, levels, norm = _custom_colormap_celsius(cube, cmap, levels, norm)
354 return cmap, levels, norm
357def _setup_spatial_map(
358 cube: iris.cube.Cube,
359 figure,
360 cmap,
361 grid_size: int | None = None,
362 subplot: int | None = None,
363):
364 """Define map projections, extent and add coastlines for spatial plots.
366 For spatial map plots, a relevant map projection for rotated or non-rotated inputs
367 is specified, and map extent defined based on the input data.
369 Parameters
370 ----------
371 cube: Cube
372 2 dimensional (lat and lon) Cube of the data to plot.
373 figure:
374 Matplotlib Figure object holding all plot elements.
375 cmap:
376 Matplotlib colormap.
377 grid_size: int, optional
378 Size of grid for subplots if multiple spatial subplots in figure.
379 subplot: int, optional
380 Subplot index if multiple spatial subplots in figure.
382 Returns
383 -------
384 axes:
385 Matplotlib GeoAxes definition.
386 """
387 # Identify min/max plot bounds.
388 try:
389 lat_axis, lon_axis = get_cube_yxcoordname(cube)
390 x1 = np.min(cube.coord(lon_axis).points)
391 x2 = np.max(cube.coord(lon_axis).points)
392 y1 = np.min(cube.coord(lat_axis).points)
393 y2 = np.max(cube.coord(lat_axis).points)
395 # Adjust bounds within +/- 180.0 if x dimension extends beyond half-globe.
396 if np.abs(x2 - x1) > 180.0:
397 x1 = x1 - 180.0
398 x2 = x2 - 180.0
399 logging.debug("Adjusting plot bounds to fit global extent.")
401 # Consider map projection orientation.
402 # Adapting orientation enables plotting across international dateline.
403 # Users can adapt the default central_longitude if alternative projections views.
404 if x2 > 180.0:
405 central_longitude = 180.0
406 else:
407 central_longitude = 0.0
409 # Define spatial map projection.
410 coord_system = cube.coord(lat_axis).coord_system
411 if isinstance(coord_system, iris.coord_systems.RotatedGeogCS):
412 # Define rotated pole map projection for rotated pole inputs.
413 projection = ccrs.RotatedPole(
414 pole_longitude=coord_system.grid_north_pole_longitude,
415 pole_latitude=coord_system.grid_north_pole_latitude,
416 central_rotated_longitude=0.0,
417 )
418 crs = projection
419 else:
420 # Define regular map projection for non-rotated pole inputs.
421 # Alternatives might include e.g. for global model outputs:
422 # projection=ccrs.Robinson(central_longitude=X.y, globe=None)
423 # See also https://scitools.org.uk/cartopy/docs/v0.15/crs/projections.html.
424 projection = ccrs.PlateCarree(central_longitude=central_longitude)
425 crs = ccrs.PlateCarree()
427 # Define axes for plot (or subplot) with required map projection.
428 if subplot is not None:
429 axes = figure.add_subplot(
430 grid_size, grid_size, subplot, projection=projection
431 )
432 else:
433 axes = figure.add_subplot(projection=projection)
435 # Add coastlines if cube contains x and y map coordinates.
436 if cmap.name in ["viridis", "Greys"]:
437 coastcol = "magenta"
438 else:
439 coastcol = "black"
440 logging.debug("Plotting coastlines in colour %s.", coastcol)
441 axes.coastlines(resolution="10m", color=coastcol)
443 # If is lat/lon spatial map, fix extent to keep plot tight.
444 # Specifying crs within set_extent helps ensure only data region is shown.
445 if isinstance(coord_system, iris.coord_systems.GeogCS):
446 axes.set_extent([x1, x2, y1, y2], crs=crs)
448 except ValueError:
449 # Skip if not both x and y map coordinates.
450 axes = figure.gca()
451 pass
453 return axes
456def _get_plot_resolution() -> int:
457 """Get resolution of rasterised plots in pixels per inch."""
458 return get_recipe_metadata().get("plot_resolution", 100)
461def _plot_and_save_spatial_plot(
462 cube: iris.cube.Cube,
463 filename: str,
464 title: str,
465 method: Literal["contourf", "pcolormesh"],
466 **kwargs,
467):
468 """Plot and save a spatial plot.
470 Parameters
471 ----------
472 cube: Cube
473 2 dimensional (lat and lon) Cube of the data to plot.
474 filename: str
475 Filename of the plot to write.
476 title: str
477 Plot title.
478 method: "contourf" | "pcolormesh"
479 The plotting method to use.
480 """
481 # Setup plot details, size, resolution, etc.
482 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
484 # Specify the color bar
485 cmap, levels, norm = _colorbar_map_levels(cube)
487 # Setup plot map projection, extent and coastlines.
488 axes = _setup_spatial_map(cube, fig, cmap)
490 # Plot the field.
491 if method == "contourf":
492 # Filled contour plot of the field.
493 plot = iplt.contourf(cube, cmap=cmap, levels=levels, norm=norm)
494 elif method == "pcolormesh":
495 try:
496 vmin = min(levels)
497 vmax = max(levels)
498 except TypeError:
499 vmin, vmax = None, None
500 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
501 # if levels are defined.
502 if norm is not None:
503 vmin = None
504 vmax = None
505 logging.debug("Plotting using defined levels.")
506 plot = iplt.pcolormesh(cube, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
507 else:
508 raise ValueError(f"Unknown plotting method: {method}")
510 # Check to see if transect, and if so, adjust y axis.
511 if is_transect(cube):
512 if "pressure" in [coord.name() for coord in cube.coords()]:
513 axes.invert_yaxis()
514 axes.set_yscale("log")
515 axes.set_ylim(1100, 100)
516 # If both model_level_number and level_height exists, iplt can construct
517 # plot as a function of height above orography (NOT sea level).
518 elif {"model_level_number", "level_height"}.issubset( 518 ↛ 523line 518 didn't jump to line 523 because the condition on line 518 was always true
519 {coord.name() for coord in cube.coords()}
520 ):
521 axes.set_yscale("log")
523 axes.set_title(
524 f"{title}\n"
525 f"Start Lat: {cube.attributes['transect_coords'].split('_')[0]}"
526 f" Start Lon: {cube.attributes['transect_coords'].split('_')[1]}"
527 f" End Lat: {cube.attributes['transect_coords'].split('_')[2]}"
528 f" End Lon: {cube.attributes['transect_coords'].split('_')[3]}",
529 fontsize=16,
530 )
532 else:
533 # Add title.
534 axes.set_title(title, fontsize=16)
536 # Add watermark with min/max/mean. Currently not user togglable.
537 # In the bbox dictionary, fc and ec are hex colour codes for grey shade.
538 axes.annotate(
539 f"Min: {np.min(cube.data):.3g} Max: {np.max(cube.data):.3g} Mean: {np.mean(cube.data):.3g}",
540 xy=(1, -0.05),
541 xycoords="axes fraction",
542 xytext=(-5, 5),
543 textcoords="offset points",
544 ha="right",
545 va="bottom",
546 size=11,
547 bbox=dict(boxstyle="round", fc="#cccccc", ec="#808080", alpha=0.9),
548 )
550 # Add colour bar.
551 cbar = fig.colorbar(plot, orientation="horizontal", pad=0.042, shrink=0.7)
552 cbar.set_label(label=f"{cube.name()} ({cube.units})", size=14)
553 # add ticks and tick_labels for every levels if less than 20 levels exist
554 if levels is not None and len(levels) < 20:
555 cbar.set_ticks(levels)
556 cbar.set_ticklabels([f"{level:.2f}" for level in levels])
557 if "visibility" in cube.name(): 557 ↛ 558line 557 didn't jump to line 558 because the condition on line 557 was never true
558 cbar.set_ticklabels([f"{level:.3g}" for level in levels])
559 logging.debug("Set colorbar ticks and labels.")
561 # Save plot.
562 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
563 logging.info("Saved spatial plot to %s", filename)
564 plt.close(fig)
567def _plot_and_save_postage_stamp_spatial_plot(
568 cube: iris.cube.Cube,
569 filename: str,
570 stamp_coordinate: str,
571 title: str,
572 method: Literal["contourf", "pcolormesh"],
573 **kwargs,
574):
575 """Plot postage stamp spatial plots from an ensemble.
577 Parameters
578 ----------
579 cube: Cube
580 Iris cube of data to be plotted. It must have the stamp coordinate.
581 filename: str
582 Filename of the plot to write.
583 stamp_coordinate: str
584 Coordinate that becomes different plots.
585 method: "contourf" | "pcolormesh"
586 The plotting method to use.
588 Raises
589 ------
590 ValueError
591 If the cube doesn't have the right dimensions.
592 """
593 # Use the smallest square grid that will fit the members.
594 grid_size = int(math.ceil(math.sqrt(len(cube.coord(stamp_coordinate).points))))
596 fig = plt.figure(figsize=(10, 10))
598 # Specify the color bar
599 cmap, levels, norm = _colorbar_map_levels(cube)
601 # Make a subplot for each member.
602 for member, subplot in zip(
603 cube.slices_over(stamp_coordinate), range(1, grid_size**2 + 1), strict=False
604 ):
605 # Setup subplot map projection, extent and coastlines.
606 axes = _setup_spatial_map(
607 member, fig, cmap, grid_size=grid_size, subplot=subplot
608 )
609 if method == "contourf":
610 # Filled contour plot of the field.
611 plot = iplt.contourf(member, cmap=cmap, levels=levels, norm=norm)
612 elif method == "pcolormesh":
613 if levels is not None:
614 vmin = min(levels)
615 vmax = max(levels)
616 else:
617 raise TypeError("Unknown vmin and vmax range.")
618 vmin, vmax = None, None
619 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
620 # if levels are defined.
621 if norm is not None: 621 ↛ 622line 621 didn't jump to line 622 because the condition on line 621 was never true
622 vmin = None
623 vmax = None
624 # pcolormesh plot of the field.
625 plot = iplt.pcolormesh(member, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
626 else:
627 raise ValueError(f"Unknown plotting method: {method}")
628 axes.set_title(f"Member #{member.coord(stamp_coordinate).points[0]}")
629 axes.set_axis_off()
631 # Put the shared colorbar in its own axes.
632 colorbar_axes = fig.add_axes([0.15, 0.07, 0.7, 0.03])
633 colorbar = fig.colorbar(
634 plot, colorbar_axes, orientation="horizontal", pad=0.042, shrink=0.7
635 )
636 colorbar.set_label(f"{cube.name()} ({cube.units})", size=14)
638 # Overall figure title.
639 fig.suptitle(title, fontsize=16)
641 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
642 logging.info("Saved contour postage stamp plot to %s", filename)
643 plt.close(fig)
646def _plot_and_save_line_series(
647 cubes: iris.cube.CubeList,
648 coords: list[iris.coords.Coord],
649 ensemble_coord: str,
650 filename: str,
651 title: str,
652 **kwargs,
653):
654 """Plot and save a 1D line series.
656 Parameters
657 ----------
658 cubes: Cube or CubeList
659 Cube or CubeList containing the cubes to plot on the y-axis.
660 coords: list[Coord]
661 Coordinates to plot on the x-axis, one per cube.
662 ensemble_coord: str
663 Ensemble coordinate in the cube.
664 filename: str
665 Filename of the plot to write.
666 title: str
667 Plot title.
668 """
669 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
671 model_colors_map = _get_model_colors_map(cubes)
673 # Store min/max ranges.
674 y_levels = []
676 # Check match-up across sequence coords gives consistent sizes
677 _validate_cubes_coords(cubes, coords)
679 for cube, coord in zip(cubes, coords, strict=True):
680 label = None
681 color = "black"
682 if model_colors_map:
683 label = cube.attributes.get("model_name")
684 color = model_colors_map.get(label)
685 for cube_slice in cube.slices_over(ensemble_coord):
686 # Label with (control) if part of an ensemble or not otherwise.
687 if cube_slice.coord(ensemble_coord).points == [0]:
688 iplt.plot(
689 coord,
690 cube_slice,
691 color=color,
692 marker="o",
693 ls="-",
694 lw=3,
695 label=f"{label} (control)"
696 if len(cube.coord(ensemble_coord).points) > 1
697 else label,
698 )
699 # Label with (perturbed) if part of an ensemble and not the control.
700 else:
701 iplt.plot(
702 coord,
703 cube_slice,
704 color=color,
705 ls="-",
706 lw=1.5,
707 alpha=0.75,
708 label=f"{label} (member)",
709 )
711 # Calculate the global min/max if multiple cubes are given.
712 _, levels, _ = _colorbar_map_levels(cube, axis="y")
713 if levels is not None: 713 ↛ 714line 713 didn't jump to line 714 because the condition on line 713 was never true
714 y_levels.append(min(levels))
715 y_levels.append(max(levels))
717 # Get the current axes.
718 ax = plt.gca()
720 # Add some labels and tweak the style.
721 # check if cubes[0] works for single cube if not CubeList
722 ax.set_xlabel(f"{coords[0].name()} / {coords[0].units}", fontsize=14)
723 ax.set_ylabel(f"{cubes[0].name()} / {cubes[0].units}", fontsize=14)
724 ax.set_title(title, fontsize=16)
726 ax.ticklabel_format(axis="y", useOffset=False)
727 ax.tick_params(axis="x", labelrotation=15)
728 ax.tick_params(axis="both", labelsize=12)
730 # Set y limits to global min and max, autoscale if colorbar doesn't exist.
731 if y_levels: 731 ↛ 732line 731 didn't jump to line 732 because the condition on line 731 was never true
732 ax.set_ylim(min(y_levels), max(y_levels))
733 # Add zero line.
734 if min(y_levels) < 0.0 and max(y_levels) > 0.0:
735 ax.axhline(y=0, xmin=0, xmax=1, ls="-", color="grey", lw=2)
736 logging.debug(
737 "Line plot with y-axis limits %s-%s", min(y_levels), max(y_levels)
738 )
739 else:
740 ax.autoscale()
742 # Add gridlines
743 ax.grid(linestyle="--", color="grey", linewidth=1)
744 # Ientify unique labels for legend
745 handles = list(
746 {
747 label: handle
748 for (handle, label) in zip(*ax.get_legend_handles_labels(), strict=True)
749 }.values()
750 )
751 ax.legend(handles=handles, loc="best", ncol=1, frameon=False, fontsize=16)
753 # Save plot.
754 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
755 logging.info("Saved line plot to %s", filename)
756 plt.close(fig)
759def _plot_and_save_vertical_line_series(
760 cubes: iris.cube.CubeList,
761 coords: list[iris.coords.Coord],
762 ensemble_coord: str,
763 filename: str,
764 series_coordinate: str,
765 title: str,
766 vmin: float,
767 vmax: float,
768 **kwargs,
769):
770 """Plot and save a 1D line series in vertical.
772 Parameters
773 ----------
774 cubes: CubeList
775 1 dimensional Cube or CubeList of the data to plot on x-axis.
776 coord: list[Coord]
777 Coordinates to plot on the y-axis, one per cube.
778 ensemble_coord: str
779 Ensemble coordinate in the cube.
780 filename: str
781 Filename of the plot to write.
782 series_coordinate: str
783 Coordinate to use as vertical axis.
784 title: str
785 Plot title.
786 vmin: float
787 Minimum value for the x-axis.
788 vmax: float
789 Maximum value for the x-axis.
790 """
791 # plot the vertical pressure axis using log scale
792 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
794 model_colors_map = _get_model_colors_map(cubes)
796 # Check match-up across sequence coords gives consistent sizes
797 _validate_cubes_coords(cubes, coords)
799 for cube, coord in zip(cubes, coords, strict=True):
800 label = None
801 color = "black"
802 if model_colors_map: 802 ↛ 803line 802 didn't jump to line 803 because the condition on line 802 was never true
803 label = cube.attributes.get("model_name")
804 color = model_colors_map.get(label)
806 for cube_slice in cube.slices_over(ensemble_coord):
807 # If ensemble data given plot control member with (control)
808 # unless single forecast.
809 if cube_slice.coord(ensemble_coord).points == [0]:
810 iplt.plot(
811 cube_slice,
812 coord,
813 color=color,
814 marker="o",
815 ls="-",
816 lw=3,
817 label=f"{label} (control)"
818 if len(cube.coord(ensemble_coord).points) > 1
819 else label,
820 )
821 # If ensemble data given plot perturbed members with (perturbed).
822 else:
823 iplt.plot(
824 cube_slice,
825 coord,
826 color=color,
827 ls="-",
828 lw=1.5,
829 alpha=0.75,
830 label=f"{label} (member)",
831 )
833 # Get the current axis
834 ax = plt.gca()
836 # Special handling for pressure level data.
837 if series_coordinate == "pressure": 837 ↛ 859line 837 didn't jump to line 859 because the condition on line 837 was always true
838 # Invert y-axis and set to log scale.
839 ax.invert_yaxis()
840 ax.set_yscale("log")
842 # Define y-ticks and labels for pressure log axis.
843 y_tick_labels = [
844 "1000",
845 "850",
846 "700",
847 "500",
848 "300",
849 "200",
850 "100",
851 ]
852 y_ticks = [1000, 850, 700, 500, 300, 200, 100]
854 # Set y-axis limits and ticks.
855 ax.set_ylim(1100, 100)
857 # Test if series_coordinate is model level data. The UM data uses
858 # model_level_number and lfric uses full_levels as coordinate.
859 elif series_coordinate in ("model_level_number", "full_levels", "half_levels"):
860 # Define y-ticks and labels for vertical axis.
861 y_ticks = iter_maybe(cubes)[0].coord(series_coordinate).points
862 y_tick_labels = [str(int(i)) for i in y_ticks]
863 ax.set_ylim(min(y_ticks), max(y_ticks))
865 ax.set_yticks(y_ticks)
866 ax.set_yticklabels(y_tick_labels)
868 # Set x-axis limits.
869 ax.set_xlim(vmin, vmax)
870 # Mark y=0 if present in plot.
871 if vmin < 0.0 and vmax > 0.0: 871 ↛ 872line 871 didn't jump to line 872 because the condition on line 871 was never true
872 ax.axvline(x=0, ymin=0, ymax=1, ls="-", color="grey", lw=2)
874 # Add some labels and tweak the style.
875 ax.set_ylabel(f"{coord.name()} / {coord.units}", fontsize=14)
876 ax.set_xlabel(
877 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
878 )
879 ax.set_title(title, fontsize=16)
880 ax.ticklabel_format(axis="x")
881 ax.tick_params(axis="y")
882 ax.tick_params(axis="both", labelsize=12)
884 # Add gridlines
885 ax.grid(linestyle="--", color="grey", linewidth=1)
886 # Ientify unique labels for legend
887 handles = list(
888 {
889 label: handle
890 for (handle, label) in zip(*ax.get_legend_handles_labels(), strict=True)
891 }.values()
892 )
893 ax.legend(handles=handles, loc="best", ncol=1, frameon=False, fontsize=16)
895 # Save plot.
896 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
897 logging.info("Saved line plot to %s", filename)
898 plt.close(fig)
901def _plot_and_save_scatter_plot(
902 cube_x: iris.cube.Cube | iris.cube.CubeList,
903 cube_y: iris.cube.Cube | iris.cube.CubeList,
904 filename: str,
905 title: str,
906 one_to_one: bool,
907 **kwargs,
908):
909 """Plot and save a 2D scatter plot.
911 Parameters
912 ----------
913 cube_x: Cube | CubeList
914 1 dimensional Cube or CubeList of the data to plot on x-axis.
915 cube_y: Cube | CubeList
916 1 dimensional Cube or CubeList of the data to plot on y-axis.
917 filename: str
918 Filename of the plot to write.
919 title: str
920 Plot title.
921 one_to_one: bool
922 Whether a 1:1 line is plotted.
923 """
924 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
925 # plot the cube_x and cube_y 1D fields as a scatter plot. If they are CubeLists this ensures
926 # to pair each cube from cube_x with the corresponding cube from cube_y, allowing to iterate
927 # over the pairs simultaneously.
929 # Ensure cube_x and cube_y are iterable
930 cube_x_iterable = iter_maybe(cube_x)
931 cube_y_iterable = iter_maybe(cube_y)
933 for cube_x_iter, cube_y_iter in zip(cube_x_iterable, cube_y_iterable, strict=True):
934 iplt.scatter(cube_x_iter, cube_y_iter)
935 if one_to_one is True:
936 plt.plot(
937 [
938 np.nanmin([np.nanmin(cube_y.data), np.nanmin(cube_x.data)]),
939 np.nanmax([np.nanmax(cube_y.data), np.nanmax(cube_x.data)]),
940 ],
941 [
942 np.nanmin([np.nanmin(cube_y.data), np.nanmin(cube_x.data)]),
943 np.nanmax([np.nanmax(cube_y.data), np.nanmax(cube_x.data)]),
944 ],
945 "k",
946 linestyle="--",
947 )
948 ax = plt.gca()
950 # Add some labels and tweak the style.
951 ax.set_xlabel(f"{cube_x[0].name()} / {cube_x[0].units}", fontsize=14)
952 ax.set_ylabel(f"{cube_y[0].name()} / {cube_y[0].units}", fontsize=14)
953 ax.set_title(title, fontsize=16)
954 ax.ticklabel_format(axis="y", useOffset=False)
955 ax.tick_params(axis="x", labelrotation=15)
956 ax.tick_params(axis="both", labelsize=12)
957 ax.autoscale()
959 # Save plot.
960 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
961 logging.info("Saved scatter plot to %s", filename)
962 plt.close(fig)
965def _plot_and_save_vector_plot(
966 cube_u: iris.cube.Cube,
967 cube_v: iris.cube.Cube,
968 filename: str,
969 title: str,
970 method: Literal["contourf", "pcolormesh"],
971 **kwargs,
972):
973 """Plot and save a 2D vector plot.
975 Parameters
976 ----------
977 cube_u: Cube
978 2 dimensional Cube of u component of the data.
979 cube_v: Cube
980 2 dimensional Cube of v component of the data.
981 filename: str
982 Filename of the plot to write.
983 title: str
984 Plot title.
985 """
986 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
988 # Create a cube containing the magnitude of the vector field.
989 cube_vec_mag = (cube_u**2 + cube_v**2) ** 0.5
990 cube_vec_mag.rename(f"{cube_u.name()}_{cube_v.name()}_magnitude")
992 # Specify the color bar
993 cmap, levels, norm = _colorbar_map_levels(cube_vec_mag)
995 # Setup plot map projection, extent and coastlines.
996 axes = _setup_spatial_map(cube_vec_mag, fig, cmap)
998 if method == "contourf": 998 ↛ 1001line 998 didn't jump to line 1001 because the condition on line 998 was always true
999 # Filled contour plot of the field.
1000 plot = iplt.contourf(cube_vec_mag, cmap=cmap, levels=levels, norm=norm)
1001 elif method == "pcolormesh":
1002 try:
1003 vmin = min(levels)
1004 vmax = max(levels)
1005 except TypeError:
1006 vmin, vmax = None, None
1007 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
1008 # if levels are defined.
1009 if norm is not None:
1010 vmin = None
1011 vmax = None
1012 plot = iplt.pcolormesh(cube_vec_mag, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
1013 else:
1014 raise ValueError(f"Unknown plotting method: {method}")
1016 # Check to see if transect, and if so, adjust y axis.
1017 if is_transect(cube_vec_mag): 1017 ↛ 1018line 1017 didn't jump to line 1018 because the condition on line 1017 was never true
1018 if "pressure" in [coord.name() for coord in cube_vec_mag.coords()]:
1019 axes.invert_yaxis()
1020 axes.set_yscale("log")
1021 axes.set_ylim(1100, 100)
1022 # If both model_level_number and level_height exists, iplt can construct
1023 # plot as a function of height above orography (NOT sea level).
1024 elif {"model_level_number", "level_height"}.issubset(
1025 {coord.name() for coord in cube_vec_mag.coords()}
1026 ):
1027 axes.set_yscale("log")
1029 axes.set_title(
1030 f"{title}\n"
1031 f"Start Lat: {cube_vec_mag.attributes['transect_coords'].split('_')[0]}"
1032 f" Start Lon: {cube_vec_mag.attributes['transect_coords'].split('_')[1]}"
1033 f" End Lat: {cube_vec_mag.attributes['transect_coords'].split('_')[2]}"
1034 f" End Lon: {cube_vec_mag.attributes['transect_coords'].split('_')[3]}",
1035 fontsize=16,
1036 )
1038 else:
1039 # Add title.
1040 axes.set_title(title, fontsize=16)
1042 # Add watermark with min/max/mean. Currently not user togglable.
1043 # In the bbox dictionary, fc and ec are hex colour codes for grey shade.
1044 axes.annotate(
1045 f"Min: {np.min(cube_vec_mag.data):.3g} Max: {np.max(cube_vec_mag.data):.3g} Mean: {np.mean(cube_vec_mag.data):.3g}",
1046 xy=(1, -0.05),
1047 xycoords="axes fraction",
1048 xytext=(-5, 5),
1049 textcoords="offset points",
1050 ha="right",
1051 va="bottom",
1052 size=11,
1053 bbox=dict(boxstyle="round", fc="#cccccc", ec="#808080", alpha=0.9),
1054 )
1056 # Add colour bar.
1057 cbar = fig.colorbar(plot, orientation="horizontal", pad=0.042, shrink=0.7)
1058 cbar.set_label(label=f"{cube_vec_mag.name()} ({cube_vec_mag.units})", size=14)
1059 # add ticks and tick_labels for every levels if less than 20 levels exist
1060 if levels is not None and len(levels) < 20: 1060 ↛ 1061line 1060 didn't jump to line 1061 because the condition on line 1060 was never true
1061 cbar.set_ticks(levels)
1062 cbar.set_ticklabels([f"{level:.1f}" for level in levels])
1064 # 30 barbs along the longest axis of the plot, or a barb per point for data
1065 # with less than 30 points.
1066 step = max(max(cube_u.shape) // 30, 1)
1067 iplt.quiver(cube_u[::step, ::step], cube_v[::step, ::step], pivot="middle")
1069 # Save plot.
1070 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1071 logging.info("Saved vector plot to %s", filename)
1072 plt.close(fig)
1075def _plot_and_save_histogram_series(
1076 cubes: iris.cube.Cube | iris.cube.CubeList,
1077 filename: str,
1078 title: str,
1079 vmin: float,
1080 vmax: float,
1081 **kwargs,
1082):
1083 """Plot and save a histogram series.
1085 Parameters
1086 ----------
1087 cubes: Cube or CubeList
1088 2 dimensional Cube or CubeList of the data to plot as histogram.
1089 filename: str
1090 Filename of the plot to write.
1091 title: str
1092 Plot title.
1093 vmin: float
1094 minimum for colorbar
1095 vmax: float
1096 maximum for colorbar
1097 """
1098 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1099 ax = plt.gca()
1101 model_colors_map = _get_model_colors_map(cubes)
1103 # Set default that histograms will produce probability density function
1104 # at each bin (integral over range sums to 1).
1105 density = True
1107 for cube in iter_maybe(cubes):
1108 # Easier to check title (where var name originates)
1109 # than seeing if long names exist etc.
1110 # Exception case, where distribution better fits log scales/bins.
1111 if "surface_microphysical" in title:
1112 if "amount" in title: 1112 ↛ 1114line 1112 didn't jump to line 1114 because the condition on line 1112 was never true
1113 # Compute histogram following Klingaman et al. (2017): ASoP
1114 bin2 = np.exp(np.log(0.02) + 0.1 * np.linspace(0, 99, 100))
1115 bins = np.pad(bin2, (1, 0), "constant", constant_values=0)
1116 density = False
1117 else:
1118 bins = 10.0 ** (
1119 np.arange(-10, 27, 1) / 10.0
1120 ) # Suggestion from RMED toolbox.
1121 bins = np.insert(bins, 0, 0)
1122 ax.set_yscale("log")
1123 vmin = bins[1]
1124 vmax = bins[-1] # Manually set vmin/vmax to override json derived value.
1125 ax.set_xscale("log")
1126 elif "lightning" in title:
1127 bins = [0, 1, 2, 3, 4, 5]
1128 else:
1129 bins = np.linspace(vmin, vmax, 51)
1130 logging.debug(
1131 "Plotting histogram with %s bins %s - %s.",
1132 np.size(bins),
1133 np.min(bins),
1134 np.max(bins),
1135 )
1137 # Reshape cube data into a single array to allow for a single histogram.
1138 # Otherwise we plot xdim histograms stacked.
1139 cube_data_1d = (cube.data).flatten()
1141 label = None
1142 color = "black"
1143 if model_colors_map: 1143 ↛ 1144line 1143 didn't jump to line 1144 because the condition on line 1143 was never true
1144 label = cube.attributes.get("model_name")
1145 color = model_colors_map[label]
1146 x, y = np.histogram(cube_data_1d, bins=bins, density=density)
1148 # Compute area under curve.
1149 if "surface_microphysical" in title and "amount" in title: 1149 ↛ 1150line 1149 didn't jump to line 1150 because the condition on line 1149 was never true
1150 bin_mean = (bins[:-1] + bins[1:]) / 2.0
1151 x = x * bin_mean / x.sum()
1152 x = x[1:]
1153 y = y[1:]
1155 ax.plot(
1156 y[:-1], x, color=color, linewidth=3, marker="o", markersize=6, label=label
1157 )
1159 # Add some labels and tweak the style.
1160 ax.set_title(title, fontsize=16)
1161 ax.set_xlabel(
1162 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
1163 )
1164 ax.set_ylabel("Normalised probability density", fontsize=14)
1165 if "surface_microphysical" in title and "amount" in title: 1165 ↛ 1166line 1165 didn't jump to line 1166 because the condition on line 1165 was never true
1166 ax.set_ylabel(
1167 f"Contribution to mean ({iter_maybe(cubes)[0].units})", fontsize=14
1168 )
1169 ax.set_xlim(vmin, vmax)
1170 ax.tick_params(axis="both", labelsize=12)
1172 # Overlay grid-lines onto histogram plot.
1173 ax.grid(linestyle="--", color="grey", linewidth=1)
1174 if model_colors_map: 1174 ↛ 1175line 1174 didn't jump to line 1175 because the condition on line 1174 was never true
1175 ax.legend(loc="best", ncol=1, frameon=False, fontsize=16)
1177 # Save plot.
1178 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1179 logging.info("Saved line plot to %s", filename)
1180 plt.close(fig)
1183def _plot_and_save_postage_stamp_histogram_series(
1184 cube: iris.cube.Cube,
1185 filename: str,
1186 title: str,
1187 stamp_coordinate: str,
1188 vmin: float,
1189 vmax: float,
1190 **kwargs,
1191):
1192 """Plot and save postage (ensemble members) stamps for a histogram series.
1194 Parameters
1195 ----------
1196 cube: Cube
1197 2 dimensional Cube of the data to plot as histogram.
1198 filename: str
1199 Filename of the plot to write.
1200 title: str
1201 Plot title.
1202 stamp_coordinate: str
1203 Coordinate that becomes different plots.
1204 vmin: float
1205 minimum for pdf x-axis
1206 vmax: float
1207 maximum for pdf x-axis
1208 """
1209 # Use the smallest square grid that will fit the members.
1210 grid_size = int(math.ceil(math.sqrt(len(cube.coord(stamp_coordinate).points))))
1212 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1213 # Make a subplot for each member.
1214 for member, subplot in zip(
1215 cube.slices_over(stamp_coordinate), range(1, grid_size**2 + 1), strict=False
1216 ):
1217 # Implicit interface is much easier here, due to needing to have the
1218 # cartopy GeoAxes generated.
1219 plt.subplot(grid_size, grid_size, subplot)
1220 # Reshape cube data into a single array to allow for a single histogram.
1221 # Otherwise we plot xdim histograms stacked.
1222 member_data_1d = (member.data).flatten()
1223 plt.hist(member_data_1d, density=True, stacked=True)
1224 ax = plt.gca()
1225 ax.set_title(f"Member #{member.coord(stamp_coordinate).points[0]}")
1226 ax.set_xlim(vmin, vmax)
1228 # Overall figure title.
1229 fig.suptitle(title, fontsize=16)
1231 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1232 logging.info("Saved histogram postage stamp plot to %s", filename)
1233 plt.close(fig)
1236def _plot_and_save_postage_stamps_in_single_plot_histogram_series(
1237 cube: iris.cube.Cube,
1238 filename: str,
1239 title: str,
1240 stamp_coordinate: str,
1241 vmin: float,
1242 vmax: float,
1243 **kwargs,
1244):
1245 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1246 ax.set_title(title, fontsize=16)
1247 ax.set_xlim(vmin, vmax)
1248 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1249 ax.set_ylabel("normalised probability density", fontsize=14)
1250 # Loop over all slices along the stamp_coordinate
1251 for member in cube.slices_over(stamp_coordinate):
1252 # Flatten the member data to 1D
1253 member_data_1d = member.data.flatten()
1254 # Plot the histogram using plt.hist
1255 plt.hist(
1256 member_data_1d,
1257 density=True,
1258 stacked=True,
1259 label=f"Member #{member.coord(stamp_coordinate).points[0]}",
1260 )
1262 # Add a legend
1263 ax.legend(fontsize=16)
1265 # Save the figure to a file
1266 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1268 # Close the figure
1269 plt.close(fig)
1272def _spatial_plot(
1273 method: Literal["contourf", "pcolormesh"],
1274 cube: iris.cube.Cube,
1275 filename: str | None,
1276 sequence_coordinate: str,
1277 stamp_coordinate: str,
1278):
1279 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1281 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1282 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1283 is present then postage stamp plots will be produced.
1285 Parameters
1286 ----------
1287 method: "contourf" | "pcolormesh"
1288 The plotting method to use.
1289 cube: Cube
1290 Iris cube of the data to plot. It should have two spatial dimensions,
1291 such as lat and lon, and may also have a another two dimension to be
1292 plotted sequentially and/or as postage stamp plots.
1293 filename: str | None
1294 Name of the plot to write, used as a prefix for plot sequences. If None
1295 uses the recipe name.
1296 sequence_coordinate: str
1297 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1298 This coordinate must exist in the cube.
1299 stamp_coordinate: str
1300 Coordinate about which to plot postage stamp plots. Defaults to
1301 ``"realization"``.
1303 Raises
1304 ------
1305 ValueError
1306 If the cube doesn't have the right dimensions.
1307 TypeError
1308 If the cube isn't a single cube.
1309 """
1310 recipe_title = get_recipe_metadata().get("title", "Untitled")
1312 # Ensure we have a name for the plot file.
1313 if filename is None:
1314 filename = slugify(recipe_title)
1316 # Ensure we've got a single cube.
1317 cube = _check_single_cube(cube)
1319 # Make postage stamp plots if stamp_coordinate exists and has more than a
1320 # single point.
1321 plotting_func = _plot_and_save_spatial_plot
1322 try:
1323 if cube.coord(stamp_coordinate).shape[0] > 1:
1324 plotting_func = _plot_and_save_postage_stamp_spatial_plot
1325 except iris.exceptions.CoordinateNotFoundError:
1326 pass
1328 # Must have a sequence coordinate.
1329 try:
1330 cube.coord(sequence_coordinate)
1331 except iris.exceptions.CoordinateNotFoundError as err:
1332 raise ValueError(f"Cube must have a {sequence_coordinate} coordinate.") from err
1334 # Create a plot for each value of the sequence coordinate.
1335 plot_index = []
1336 nplot = np.size(cube.coord(sequence_coordinate).points)
1337 for cube_slice in cube.slices_over(sequence_coordinate):
1338 # Use sequence value so multiple sequences can merge.
1339 sequence_value = cube_slice.coord(sequence_coordinate).points[0]
1340 plot_filename = f"{filename.rsplit('.', 1)[0]}_{sequence_value}.png"
1341 coord = cube_slice.coord(sequence_coordinate)
1342 # Format the coordinate value in a unit appropriate way.
1343 title = f"{recipe_title}\n [{coord.units.title(coord.points[0])}]"
1344 # Use sequence (e.g. time) bounds if plotting single non-sequence outputs
1345 if nplot == 1 and coord.has_bounds:
1346 if np.size(coord.bounds) > 1:
1347 title = f"{recipe_title}\n [{coord.units.title(coord.bounds[0][0])} to {coord.units.title(coord.bounds[0][1])}]"
1348 # Do the actual plotting.
1349 plotting_func(
1350 cube_slice,
1351 filename=plot_filename,
1352 stamp_coordinate=stamp_coordinate,
1353 title=title,
1354 method=method,
1355 )
1356 plot_index.append(plot_filename)
1358 # Add list of plots to plot metadata.
1359 complete_plot_index = _append_to_plot_index(plot_index)
1361 # Make a page to display the plots.
1362 _make_plot_html_page(complete_plot_index)
1365def _custom_colormap_mask(cube: iris.cube.Cube, axis: Literal["x", "y"] | None = None):
1366 """Get colourmap for mask.
1368 If "mask_for_" appears anywhere in the name of a cube this function will be called
1369 regardless of the name of the variable to ensure a consistent plot.
1371 Parameters
1372 ----------
1373 cube: Cube
1374 Cube of variable for which the colorbar information is desired.
1375 axis: "x", "y", optional
1376 Select the levels for just this axis of a line plot. The min and max
1377 can be set by xmin/xmax or ymin/ymax respectively. For variables where
1378 setting a universal range is not desirable (e.g. temperature), users
1379 can set ymin/ymax values to "auto" in the colorbar definitions file.
1380 Where no additional xmin/xmax or ymin/ymax values are provided, the
1381 axis bounds default to use the vmin/vmax values provided.
1383 Returns
1384 -------
1385 cmap:
1386 Matplotlib colormap.
1387 levels:
1388 List of levels to use for plotting. For continuous plots the min and max
1389 should be taken as the range.
1390 norm:
1391 BoundaryNorm information.
1392 """
1393 if "difference" not in cube.long_name:
1394 if axis:
1395 levels = [0, 1]
1396 # Complete settings based on levels.
1397 return None, levels, None
1398 else:
1399 # Define the levels and colors.
1400 levels = [0, 1, 2]
1401 colors = ["white", "dodgerblue"]
1402 # Create a custom color map.
1403 cmap = mcolors.ListedColormap(colors)
1404 # Normalize the levels.
1405 norm = mcolors.BoundaryNorm(levels, cmap.N)
1406 logging.debug("Colourmap for %s.", cube.long_name)
1407 return cmap, levels, norm
1408 else:
1409 if axis:
1410 levels = [-1, 1]
1411 return None, levels, None
1412 else:
1413 # Search for if mask difference, set to +/- 0.5 as values plotted <
1414 # not <=.
1415 levels = [-2, -0.5, 0.5, 2]
1416 colors = ["goldenrod", "white", "teal"]
1417 cmap = mcolors.ListedColormap(colors)
1418 norm = mcolors.BoundaryNorm(levels, cmap.N)
1419 logging.debug("Colourmap for %s.", cube.long_name)
1420 return cmap, levels, norm
1423def _custom_beaufort_scale(cube: iris.cube.Cube, axis: Literal["x", "y"] | None = None):
1424 """Get a custom colorbar for a cube in the Beaufort Scale.
1426 Specific variable ranges can be separately set in user-supplied definition
1427 for x- or y-axis limits, or indicate where automated range preferred.
1429 Parameters
1430 ----------
1431 cube: Cube
1432 Cube of variable with Beaufort Scale in name.
1433 axis: "x", "y", optional
1434 Select the levels for just this axis of a line plot. The min and max
1435 can be set by xmin/xmax or ymin/ymax respectively. For variables where
1436 setting a universal range is not desirable (e.g. temperature), users
1437 can set ymin/ymax values to "auto" in the colorbar definitions file.
1438 Where no additional xmin/xmax or ymin/ymax values are provided, the
1439 axis bounds default to use the vmin/vmax values provided.
1441 Returns
1442 -------
1443 cmap:
1444 Matplotlib colormap.
1445 levels:
1446 List of levels to use for plotting. For continuous plots the min and max
1447 should be taken as the range.
1448 norm:
1449 BoundaryNorm information.
1450 """
1451 if "difference" not in cube.long_name:
1452 if axis:
1453 levels = [0, 12]
1454 return None, levels, None
1455 else:
1456 levels = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13]
1457 colors = [
1458 "black",
1459 (0, 0, 0.6),
1460 "blue",
1461 "cyan",
1462 "green",
1463 "yellow",
1464 (1, 0.5, 0),
1465 "red",
1466 "pink",
1467 "magenta",
1468 "purple",
1469 "maroon",
1470 "white",
1471 ]
1472 cmap = mcolors.ListedColormap(colors)
1473 norm = mcolors.BoundaryNorm(levels, cmap.N)
1474 logging.info("change colormap for Beaufort Scale colorbar.")
1475 return cmap, levels, norm
1476 else:
1477 if axis:
1478 levels = [-4, 4]
1479 return None, levels, None
1480 else:
1481 levels = [
1482 -3.5,
1483 -2.5,
1484 -1.5,
1485 -0.5,
1486 0.5,
1487 1.5,
1488 2.5,
1489 3.5,
1490 ]
1491 cmap = plt.get_cmap("bwr", 8)
1492 norm = mcolors.BoundaryNorm(levels, cmap.N)
1493 return cmap, levels, norm
1496def _custom_colormap_celsius(cube: iris.cube.Cube, cmap, levels, norm):
1497 """Return altered colourmap for temperature with change in units to Celsius.
1499 If "Celsius" appears anywhere in the name of a cube this function will be called.
1501 Parameters
1502 ----------
1503 cube: Cube
1504 Cube of variable for which the colorbar information is desired.
1505 cmap: Matplotlib colormap.
1506 levels: List
1507 List of levels to use for plotting. For continuous plots the min and max
1508 should be taken as the range.
1509 norm: BoundaryNorm.
1511 Returns
1512 -------
1513 cmap: Matplotlib colormap.
1514 levels: List
1515 List of levels to use for plotting. For continuous plots the min and max
1516 should be taken as the range.
1517 norm: BoundaryNorm.
1518 """
1519 varnames = filter(None, [cube.long_name, cube.standard_name, cube.var_name])
1520 if any("temperature" in name for name in varnames) and "Celsius" == cube.units:
1521 levels = np.array(levels)
1522 levels -= 273
1523 levels = levels.tolist()
1524 else:
1525 # Do nothing keep the existing colourbar attributes
1526 levels = levels
1527 cmap = cmap
1528 norm = norm
1529 return cmap, levels, norm
1532def _custom_colormap_probability(
1533 cube: iris.cube.Cube, axis: Literal["x", "y"] | None = None
1534):
1535 """Get a custom colorbar for a probability cube.
1537 Specific variable ranges can be separately set in user-supplied definition
1538 for x- or y-axis limits, or indicate where automated range preferred.
1540 Parameters
1541 ----------
1542 cube: Cube
1543 Cube of variable with probability in name.
1544 axis: "x", "y", optional
1545 Select the levels for just this axis of a line plot. The min and max
1546 can be set by xmin/xmax or ymin/ymax respectively. For variables where
1547 setting a universal range is not desirable (e.g. temperature), users
1548 can set ymin/ymax values to "auto" in the colorbar definitions file.
1549 Where no additional xmin/xmax or ymin/ymax values are provided, the
1550 axis bounds default to use the vmin/vmax values provided.
1552 Returns
1553 -------
1554 cmap:
1555 Matplotlib colormap.
1556 levels:
1557 List of levels to use for plotting. For continuous plots the min and max
1558 should be taken as the range.
1559 norm:
1560 BoundaryNorm information.
1561 """
1562 if axis:
1563 levels = [0, 1]
1564 return None, levels, None
1565 else:
1566 cmap = mcolors.ListedColormap(
1567 [
1568 "#FFFFFF",
1569 "#636363",
1570 "#e1dada",
1571 "#B5CAFF",
1572 "#8FB3FF",
1573 "#7F97FF",
1574 "#ABCF63",
1575 "#E8F59E",
1576 "#FFFA14",
1577 "#FFD121",
1578 "#FFA30A",
1579 ]
1580 )
1581 levels = [0.0, 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0]
1582 norm = mcolors.BoundaryNorm(levels, cmap.N)
1583 return cmap, levels, norm
1586def _custom_colourmap_precipitation(cube: iris.cube.Cube, cmap, levels, norm):
1587 """Return a custom colourmap for the current recipe."""
1588 varnames = filter(None, [cube.long_name, cube.standard_name, cube.var_name])
1589 if (
1590 any("surface_microphysical" in name for name in varnames)
1591 and "difference" not in cube.long_name
1592 and "mask" not in cube.long_name
1593 ):
1594 # Define the levels and colors
1595 levels = [0, 0.125, 0.25, 0.5, 1, 2, 4, 8, 16, 32, 64, 128, 256]
1596 colors = [
1597 "w",
1598 (0, 0, 0.6),
1599 "b",
1600 "c",
1601 "g",
1602 "y",
1603 (1, 0.5, 0),
1604 "r",
1605 "pink",
1606 "m",
1607 "purple",
1608 "maroon",
1609 "gray",
1610 ]
1611 # Create a custom colormap
1612 cmap = mcolors.ListedColormap(colors)
1613 # Normalize the levels
1614 norm = mcolors.BoundaryNorm(levels, cmap.N)
1615 logging.info("change colormap for surface_microphysical variable colorbar.")
1616 else:
1617 # do nothing and keep existing colorbar attributes
1618 cmap = cmap
1619 levels = levels
1620 norm = norm
1621 return cmap, levels, norm
1624def _custom_colormap_aviation_colour_state(cube: iris.cube.Cube):
1625 """Return custom colourmap for aviation colour state.
1627 If "aviation_colour_state" appears anywhere in the name of a cube
1628 this function will be called.
1630 Parameters
1631 ----------
1632 cube: Cube
1633 Cube of variable for which the colorbar information is desired.
1635 Returns
1636 -------
1637 cmap: Matplotlib colormap.
1638 levels: List
1639 List of levels to use for plotting. For continuous plots the min and max
1640 should be taken as the range.
1641 norm: BoundaryNorm.
1642 """
1643 levels = [-0.5, 0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5]
1644 colors = [
1645 "#87ceeb",
1646 "#ffffff",
1647 "#8ced69",
1648 "#ffff00",
1649 "#ffd700",
1650 "#ffa500",
1651 "#fe3620",
1652 ]
1653 # Create a custom colormap
1654 cmap = mcolors.ListedColormap(colors)
1655 # Normalise the levels
1656 norm = mcolors.BoundaryNorm(levels, cmap.N)
1657 return cmap, levels, norm
1660def _custom_colourmap_visibility_in_air(cube: iris.cube.Cube, cmap, levels, norm):
1661 """Return a custom colourmap for the current recipe."""
1662 varnames = filter(None, [cube.long_name, cube.standard_name, cube.var_name])
1663 if (
1664 any("visibility_in_air" in name for name in varnames)
1665 and "difference" not in cube.long_name
1666 and "mask" not in cube.long_name
1667 ):
1668 # Define the levels and colors (in km)
1669 levels = [0, 0.05, 0.1, 0.2, 1.0, 2.0, 5.0, 10.0, 20.0, 30.0, 50.0, 70.0, 100.0]
1670 norm = mcolors.BoundaryNorm(levels, cmap.N)
1671 colours = [
1672 "#8f00d6",
1673 "#d10000",
1674 "#ff9700",
1675 "#ffff00",
1676 "#00007f",
1677 "#6c9ccd",
1678 "#aae8ff",
1679 "#37a648",
1680 "#8edc64",
1681 "#c5ffc5",
1682 "#dcdcdc",
1683 "#ffffff",
1684 ]
1685 # Create a custom colormap
1686 cmap = mcolors.ListedColormap(colours)
1687 # Normalize the levels
1688 norm = mcolors.BoundaryNorm(levels, cmap.N)
1689 logging.info("change colormap for visibility_in_air variable colorbar.")
1690 else:
1691 # do nothing and keep existing colorbar attributes
1692 cmap = cmap
1693 levels = levels
1694 norm = norm
1695 return cmap, levels, norm
1698def _get_num_models(cube: iris.cube.Cube | iris.cube.CubeList) -> int:
1699 """Return number of models based on cube attributes."""
1700 model_names = list(
1701 filter(
1702 lambda x: x is not None,
1703 {cb.attributes.get("model_name", None) for cb in iter_maybe(cube)},
1704 )
1705 )
1706 if not model_names:
1707 logging.debug("Missing model names. Will assume single model.")
1708 return 1
1709 else:
1710 return len(model_names)
1713def _validate_cube_shape(
1714 cube: iris.cube.Cube | iris.cube.CubeList, num_models: int
1715) -> None:
1716 """Check all cubes have a model name."""
1717 if isinstance(cube, iris.cube.CubeList) and len(cube) != num_models: 1717 ↛ 1718line 1717 didn't jump to line 1718 because the condition on line 1717 was never true
1718 raise ValueError(
1719 f"The number of model names ({num_models}) should equal the number "
1720 f"of cubes ({len(cube)})."
1721 )
1724def _validate_cubes_coords(
1725 cubes: iris.cube.CubeList, coords: list[iris.coords.Coord]
1726) -> None:
1727 """Check same number of cubes as sequence coordinate for zip functions."""
1728 if len(cubes) != len(coords): 1728 ↛ 1729line 1728 didn't jump to line 1729 because the condition on line 1728 was never true
1729 raise ValueError(
1730 f"The number of CubeList entries ({len(cubes)}) should equal the number "
1731 f"of sequence coordinates ({len(coords)})."
1732 f"Check that number of time entries in input data are consistent if "
1733 f"performing time-averaging steps prior to plotting outputs."
1734 )
1737####################
1738# Public functions #
1739####################
1742def spatial_contour_plot(
1743 cube: iris.cube.Cube,
1744 filename: str = None,
1745 sequence_coordinate: str = "time",
1746 stamp_coordinate: str = "realization",
1747 **kwargs,
1748) -> iris.cube.Cube:
1749 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1751 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1752 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1753 is present then postage stamp plots will be produced.
1755 Parameters
1756 ----------
1757 cube: Cube
1758 Iris cube of the data to plot. It should have two spatial dimensions,
1759 such as lat and lon, and may also have a another two dimension to be
1760 plotted sequentially and/or as postage stamp plots.
1761 filename: str, optional
1762 Name of the plot to write, used as a prefix for plot sequences. Defaults
1763 to the recipe name.
1764 sequence_coordinate: str, optional
1765 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1766 This coordinate must exist in the cube.
1767 stamp_coordinate: str, optional
1768 Coordinate about which to plot postage stamp plots. Defaults to
1769 ``"realization"``.
1771 Returns
1772 -------
1773 Cube
1774 The original cube (so further operations can be applied).
1776 Raises
1777 ------
1778 ValueError
1779 If the cube doesn't have the right dimensions.
1780 TypeError
1781 If the cube isn't a single cube.
1782 """
1783 _spatial_plot("contourf", cube, filename, sequence_coordinate, stamp_coordinate)
1784 return cube
1787def spatial_pcolormesh_plot(
1788 cube: iris.cube.Cube,
1789 filename: str = None,
1790 sequence_coordinate: str = "time",
1791 stamp_coordinate: str = "realization",
1792 **kwargs,
1793) -> iris.cube.Cube:
1794 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1796 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1797 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1798 is present then postage stamp plots will be produced.
1800 This function is significantly faster than ``spatial_contour_plot``,
1801 especially at high resolutions, and should be preferred unless contiguous
1802 contour areas are important.
1804 Parameters
1805 ----------
1806 cube: Cube
1807 Iris cube of the data to plot. It should have two spatial dimensions,
1808 such as lat and lon, and may also have a another two dimension to be
1809 plotted sequentially and/or as postage stamp plots.
1810 filename: str, optional
1811 Name of the plot to write, used as a prefix for plot sequences. Defaults
1812 to the recipe name.
1813 sequence_coordinate: str, optional
1814 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1815 This coordinate must exist in the cube.
1816 stamp_coordinate: str, optional
1817 Coordinate about which to plot postage stamp plots. Defaults to
1818 ``"realization"``.
1820 Returns
1821 -------
1822 Cube
1823 The original cube (so further operations can be applied).
1825 Raises
1826 ------
1827 ValueError
1828 If the cube doesn't have the right dimensions.
1829 TypeError
1830 If the cube isn't a single cube.
1831 """
1832 _spatial_plot("pcolormesh", cube, filename, sequence_coordinate, stamp_coordinate)
1833 return cube
1836# TODO: Expand function to handle ensemble data.
1837# line_coordinate: str, optional
1838# Coordinate about which to plot multiple lines. Defaults to
1839# ``"realization"``.
1840def plot_line_series(
1841 cube: iris.cube.Cube | iris.cube.CubeList,
1842 filename: str = None,
1843 series_coordinate: str = "time",
1844 # line_coordinate: str = "realization",
1845 **kwargs,
1846) -> iris.cube.Cube | iris.cube.CubeList:
1847 """Plot a line plot for the specified coordinate.
1849 The Cube or CubeList must be 1D.
1851 Parameters
1852 ----------
1853 iris.cube | iris.cube.CubeList
1854 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
1855 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
1856 We do not support different data such as temperature and humidity in the same CubeList for plotting.
1857 filename: str, optional
1858 Name of the plot to write, used as a prefix for plot sequences. Defaults
1859 to the recipe name.
1860 series_coordinate: str, optional
1861 Coordinate about which to make a series. Defaults to ``"time"``. This
1862 coordinate must exist in the cube.
1864 Returns
1865 -------
1866 iris.cube.Cube | iris.cube.CubeList
1867 The original Cube or CubeList (so further operations can be applied).
1868 plotted data.
1870 Raises
1871 ------
1872 ValueError
1873 If the cubes don't have the right dimensions.
1874 TypeError
1875 If the cube isn't a Cube or CubeList.
1876 """
1877 # Ensure we have a name for the plot file.
1878 title = get_recipe_metadata().get("title", "Untitled")
1880 if filename is None:
1881 filename = slugify(title)
1883 # Add file extension.
1884 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
1886 num_models = _get_num_models(cube)
1888 _validate_cube_shape(cube, num_models)
1890 # Iterate over all cubes and extract coordinate to plot.
1891 cubes = iter_maybe(cube)
1892 coords = []
1893 for cube in cubes:
1894 try:
1895 coords.append(cube.coord(series_coordinate))
1896 except iris.exceptions.CoordinateNotFoundError as err:
1897 raise ValueError(
1898 f"Cube must have a {series_coordinate} coordinate."
1899 ) from err
1900 if cube.ndim > 2 or not cube.coords("realization"):
1901 raise ValueError("Cube must be 1D or 2D with a realization coordinate.")
1903 # Do the actual plotting.
1904 _plot_and_save_line_series(cubes, coords, "realization", plot_filename, title)
1906 # Add list of plots to plot metadata.
1907 plot_index = _append_to_plot_index([plot_filename])
1909 # Make a page to display the plots.
1910 _make_plot_html_page(plot_index)
1912 return cube
1915def plot_vertical_line_series(
1916 cubes: iris.cube.Cube | iris.cube.CubeList,
1917 filename: str = None,
1918 series_coordinate: str = "model_level_number",
1919 sequence_coordinate: str = "time",
1920 # line_coordinate: str = "realization",
1921 **kwargs,
1922) -> iris.cube.Cube | iris.cube.CubeList:
1923 """Plot a line plot against a type of vertical coordinate.
1925 The Cube or CubeList must be 1D.
1927 A 1D line plot with y-axis as pressure coordinate can be plotted, but if the sequence_coordinate is present
1928 then a sequence of plots will be produced.
1930 Parameters
1931 ----------
1932 iris.cube | iris.cube.CubeList
1933 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
1934 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
1935 We do not support different data such as temperature and humidity in the same CubeList for plotting.
1936 filename: str, optional
1937 Name of the plot to write, used as a prefix for plot sequences. Defaults
1938 to the recipe name.
1939 series_coordinate: str, optional
1940 Coordinate to plot on the y-axis. Can be ``pressure`` or
1941 ``model_level_number`` for UM, or ``full_levels`` or ``half_levels``
1942 for LFRic. Defaults to ``model_level_number``.
1943 This coordinate must exist in the cube.
1944 sequence_coordinate: str, optional
1945 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1946 This coordinate must exist in the cube.
1948 Returns
1949 -------
1950 iris.cube.Cube | iris.cube.CubeList
1951 The original Cube or CubeList (so further operations can be applied).
1952 Plotted data.
1954 Raises
1955 ------
1956 ValueError
1957 If the cubes doesn't have the right dimensions.
1958 TypeError
1959 If the cube isn't a Cube or CubeList.
1960 """
1961 # Ensure we have a name for the plot file.
1962 recipe_title = get_recipe_metadata().get("title", "Untitled")
1964 if filename is None:
1965 filename = slugify(recipe_title)
1967 cubes = iter_maybe(cubes)
1968 # Initialise empty list to hold all data from all cubes in a CubeList
1969 all_data = []
1971 # Store min/max ranges for x range.
1972 x_levels = []
1974 num_models = _get_num_models(cubes)
1976 _validate_cube_shape(cubes, num_models)
1978 # Iterate over all cubes in cube or CubeList and plot.
1979 coords = []
1980 for cube in cubes:
1981 # Test if series coordinate i.e. pressure level exist for any cube with cube.ndim >=1.
1982 try:
1983 coords.append(cube.coord(series_coordinate))
1984 except iris.exceptions.CoordinateNotFoundError as err:
1985 raise ValueError(
1986 f"Cube must have a {series_coordinate} coordinate."
1987 ) from err
1989 try:
1990 if cube.ndim > 1 or not cube.coords("realization"): 1990 ↛ 1998line 1990 didn't jump to line 1998 because the condition on line 1990 was always true
1991 cube.coord(sequence_coordinate)
1992 except iris.exceptions.CoordinateNotFoundError as err:
1993 raise ValueError(
1994 f"Cube must have a {sequence_coordinate} coordinate or be 1D, or 2D with a realization coordinate."
1995 ) from err
1997 # Get minimum and maximum from levels information.
1998 _, levels, _ = _colorbar_map_levels(cube, axis="x")
1999 if levels is not None: 1999 ↛ 2003line 1999 didn't jump to line 2003 because the condition on line 1999 was always true
2000 x_levels.append(min(levels))
2001 x_levels.append(max(levels))
2002 else:
2003 all_data.append(cube.data)
2005 if len(x_levels) == 0: 2005 ↛ 2007line 2005 didn't jump to line 2007 because the condition on line 2005 was never true
2006 # Combine all data into a single NumPy array
2007 combined_data = np.concatenate(all_data)
2009 # Set the lower and upper limit for the x-axis to ensure all plots have
2010 # same range. This needs to read the whole cube over the range of the
2011 # sequence and if applicable postage stamp coordinate.
2012 vmin = np.floor(combined_data.min())
2013 vmax = np.ceil(combined_data.max())
2014 else:
2015 vmin = min(x_levels)
2016 vmax = max(x_levels)
2018 # Matching the slices (matching by seq coord point; it may happen that
2019 # evaluated models do not cover the same seq coord range, hence matching
2020 # necessary)
2021 def filter_cube_iterables(cube_iterables) -> bool:
2022 return len(cube_iterables) == len(coords)
2024 cube_iterables = filter(
2025 filter_cube_iterables,
2026 (
2027 iris.cube.CubeList(
2028 s
2029 for s in itertools.chain.from_iterable(
2030 cb.slices_over(sequence_coordinate) for cb in cubes
2031 )
2032 if s.coord(sequence_coordinate).points[0] == point
2033 )
2034 for point in sorted(
2035 set(
2036 itertools.chain.from_iterable(
2037 cb.coord(sequence_coordinate).points for cb in cubes
2038 )
2039 )
2040 )
2041 ),
2042 )
2044 # Create a plot for each value of the sequence coordinate.
2045 # Allowing for multiple cubes in a CubeList to be plotted in the same plot for
2046 # similar sequence values. Passing a CubeList into the internal plotting function
2047 # for similar values of the sequence coordinate. cube_slice can be an iris.cube.Cube
2048 # or an iris.cube.CubeList.
2049 plot_index = []
2050 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2051 for cubes_slice in cube_iterables:
2052 # Use sequence value so multiple sequences can merge.
2053 seq_coord = cubes_slice[0].coord(sequence_coordinate)
2054 sequence_value = seq_coord.points[0]
2055 plot_filename = f"{filename.rsplit('.', 1)[0]}_{sequence_value}.png"
2056 # Format the coordinate value in a unit appropriate way.
2057 title = f"{recipe_title}\n [{seq_coord.units.title(sequence_value)}]"
2058 # Use sequence (e.g. time) bounds if plotting single non-sequence outputs
2059 if nplot == 1 and seq_coord.has_bounds: 2059 ↛ 2060line 2059 didn't jump to line 2060 because the condition on line 2059 was never true
2060 if np.size(seq_coord.bounds) > 1:
2061 title = f"{recipe_title}\n [{seq_coord.units.title(seq_coord.bounds[0][0])} to {seq_coord.units.title(seq_coord.bounds[0][1])}]"
2062 # Do the actual plotting.
2063 _plot_and_save_vertical_line_series(
2064 cubes_slice,
2065 coords,
2066 "realization",
2067 plot_filename,
2068 series_coordinate,
2069 title=title,
2070 vmin=vmin,
2071 vmax=vmax,
2072 )
2073 plot_index.append(plot_filename)
2075 # Add list of plots to plot metadata.
2076 complete_plot_index = _append_to_plot_index(plot_index)
2078 # Make a page to display the plots.
2079 _make_plot_html_page(complete_plot_index)
2081 return cubes
2084def scatter_plot(
2085 cube_x: iris.cube.Cube | iris.cube.CubeList,
2086 cube_y: iris.cube.Cube | iris.cube.CubeList,
2087 filename: str = None,
2088 one_to_one: bool = True,
2089 **kwargs,
2090) -> iris.cube.CubeList:
2091 """Plot a scatter plot between two variables.
2093 Both cubes must be 1D.
2095 Parameters
2096 ----------
2097 cube_x: Cube | CubeList
2098 1 dimensional Cube of the data to plot on y-axis.
2099 cube_y: Cube | CubeList
2100 1 dimensional Cube of the data to plot on x-axis.
2101 filename: str, optional
2102 Filename of the plot to write.
2103 one_to_one: bool, optional
2104 If True a 1:1 line is plotted; if False it is not. Default is True.
2106 Returns
2107 -------
2108 cubes: CubeList
2109 CubeList of the original x and y cubes for further processing.
2111 Raises
2112 ------
2113 ValueError
2114 If the cube doesn't have the right dimensions and cubes not the same
2115 size.
2116 TypeError
2117 If the cube isn't a single cube.
2119 Notes
2120 -----
2121 Scatter plots are used for determining if there is a relationship between
2122 two variables. Positive relations have a slope going from bottom left to top
2123 right; Negative relations have a slope going from top left to bottom right.
2125 A variant of the scatter plot is the quantile-quantile plot. This plot does
2126 not use all data points, but the selected quantiles of each variable
2127 instead. Quantile-quantile plots are valuable for comparing against
2128 observations and other models. Identical percentiles between the variables
2129 will lie on the one-to-one line implying the values correspond well to each
2130 other. Where there is a deviation from the one-to-one line a range of
2131 possibilities exist depending on how and where the data is shifted (e.g.,
2132 Wilks 2011 [Wilks2011]_).
2134 For distributions above the one-to-one line the distribution is left-skewed;
2135 below is right-skewed. A distinct break implies a bimodal distribution, and
2136 closer values/values further apart at the tails imply poor representation of
2137 the extremes.
2139 References
2140 ----------
2141 .. [Wilks2011] Wilks, D.S., (2011) "Statistical Methods in the Atmospheric
2142 Sciences" Third Edition, vol. 100, Academic Press, Oxford, UK, 676 pp.
2143 """
2144 # Iterate over all cubes in cube or CubeList and plot.
2145 for cube_iter in iter_maybe(cube_x):
2146 # Check cubes are correct shape.
2147 cube_iter = _check_single_cube(cube_iter)
2148 if cube_iter.ndim > 1:
2149 raise ValueError("cube_x must be 1D.")
2151 # Iterate over all cubes in cube or CubeList and plot.
2152 for cube_iter in iter_maybe(cube_y):
2153 # Check cubes are correct shape.
2154 cube_iter = _check_single_cube(cube_iter)
2155 if cube_iter.ndim > 1:
2156 raise ValueError("cube_y must be 1D.")
2158 # Ensure we have a name for the plot file.
2159 title = get_recipe_metadata().get("title", "Untitled")
2161 if filename is None:
2162 filename = slugify(title)
2164 # Add file extension.
2165 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2167 # Do the actual plotting.
2168 _plot_and_save_scatter_plot(cube_x, cube_y, plot_filename, title, one_to_one)
2170 # Add list of plots to plot metadata.
2171 plot_index = _append_to_plot_index([plot_filename])
2173 # Make a page to display the plots.
2174 _make_plot_html_page(plot_index)
2176 return iris.cube.CubeList([cube_x, cube_y])
2179def vector_plot(
2180 cube_u: iris.cube.Cube,
2181 cube_v: iris.cube.Cube,
2182 filename: str = None,
2183 sequence_coordinate: str = "time",
2184 **kwargs,
2185) -> iris.cube.CubeList:
2186 """Plot a vector plot based on the input u and v components."""
2187 recipe_title = get_recipe_metadata().get("title", "Untitled")
2189 # Ensure we have a name for the plot file.
2190 if filename is None: 2190 ↛ 2191line 2190 didn't jump to line 2191 because the condition on line 2190 was never true
2191 filename = slugify(recipe_title)
2193 # Cubes must have a matching sequence coordinate.
2194 try:
2195 # Check that the u and v cubes have the same sequence coordinate.
2196 if cube_u.coord(sequence_coordinate) != cube_v.coord(sequence_coordinate): 2196 ↛ 2197line 2196 didn't jump to line 2197 because the condition on line 2196 was never true
2197 raise ValueError("Coordinates do not match.")
2198 except (iris.exceptions.CoordinateNotFoundError, ValueError) as err:
2199 raise ValueError(
2200 f"Cubes should have matching {sequence_coordinate} coordinate:\n{cube_u}\n{cube_v}"
2201 ) from err
2203 # Create a plot for each value of the sequence coordinate.
2204 plot_index = []
2205 for cube_u_slice, cube_v_slice in zip(
2206 cube_u.slices_over(sequence_coordinate),
2207 cube_v.slices_over(sequence_coordinate),
2208 strict=True,
2209 ):
2210 # Use sequence value so multiple sequences can merge.
2211 sequence_value = cube_u_slice.coord(sequence_coordinate).points[0]
2212 plot_filename = f"{filename.rsplit('.', 1)[0]}_{sequence_value}.png"
2213 coord = cube_u_slice.coord(sequence_coordinate)
2214 # Format the coordinate value in a unit appropriate way.
2215 title = f"{recipe_title}\n{coord.units.title(coord.points[0])}"
2216 # Do the actual plotting.
2217 _plot_and_save_vector_plot(
2218 cube_u_slice,
2219 cube_v_slice,
2220 filename=plot_filename,
2221 title=title,
2222 method="contourf",
2223 )
2224 plot_index.append(plot_filename)
2226 # Add list of plots to plot metadata.
2227 complete_plot_index = _append_to_plot_index(plot_index)
2229 # Make a page to display the plots.
2230 _make_plot_html_page(complete_plot_index)
2232 return iris.cube.CubeList([cube_u, cube_v])
2235def plot_histogram_series(
2236 cubes: iris.cube.Cube | iris.cube.CubeList,
2237 filename: str = None,
2238 sequence_coordinate: str = "time",
2239 stamp_coordinate: str = "realization",
2240 single_plot: bool = False,
2241 **kwargs,
2242) -> iris.cube.Cube | iris.cube.CubeList:
2243 """Plot a histogram plot for each vertical level provided.
2245 A histogram plot can be plotted, but if the sequence_coordinate (i.e. time)
2246 is present then a sequence of plots will be produced using the time slider
2247 functionality to scroll through histograms against time. If a
2248 stamp_coordinate is present then postage stamp plots will be produced. If
2249 stamp_coordinate and single_plot is True, all postage stamp plots will be
2250 plotted in a single plot instead of separate postage stamp plots.
2252 Parameters
2253 ----------
2254 cubes: Cube | iris.cube.CubeList
2255 Iris cube or CubeList of the data to plot. It should have a single dimension other
2256 than the stamp coordinate.
2257 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2258 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2259 filename: str, optional
2260 Name of the plot to write, used as a prefix for plot sequences. Defaults
2261 to the recipe name.
2262 sequence_coordinate: str, optional
2263 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2264 This coordinate must exist in the cube and will be used for the time
2265 slider.
2266 stamp_coordinate: str, optional
2267 Coordinate about which to plot postage stamp plots. Defaults to
2268 ``"realization"``.
2269 single_plot: bool, optional
2270 If True, all postage stamp plots will be plotted in a single plot. If
2271 False, each postage stamp plot will be plotted separately. Is only valid
2272 if stamp_coordinate exists and has more than a single point.
2274 Returns
2275 -------
2276 iris.cube.Cube | iris.cube.CubeList
2277 The original Cube or CubeList (so further operations can be applied).
2278 Plotted data.
2280 Raises
2281 ------
2282 ValueError
2283 If the cube doesn't have the right dimensions.
2284 TypeError
2285 If the cube isn't a Cube or CubeList.
2286 """
2287 recipe_title = get_recipe_metadata().get("title", "Untitled")
2289 cubes = iter_maybe(cubes)
2291 # Ensure we have a name for the plot file.
2292 if filename is None:
2293 filename = slugify(recipe_title)
2295 # Internal plotting function.
2296 plotting_func = _plot_and_save_histogram_series
2298 num_models = _get_num_models(cubes)
2300 _validate_cube_shape(cubes, num_models)
2302 # If several histograms are plotted with time as sequence_coordinate for the
2303 # time slider option.
2304 for cube in cubes:
2305 try:
2306 cube.coord(sequence_coordinate)
2307 except iris.exceptions.CoordinateNotFoundError as err:
2308 raise ValueError(
2309 f"Cube must have a {sequence_coordinate} coordinate."
2310 ) from err
2312 # Get minimum and maximum from levels information.
2313 levels = None
2314 for cube in cubes: 2314 ↛ 2330line 2314 didn't jump to line 2330 because the loop on line 2314 didn't complete
2315 # First check if user-specified "auto" range variable.
2316 # This maintains the value of levels as None, so proceed.
2317 _, levels, _ = _colorbar_map_levels(cube, axis="y")
2318 if levels is None:
2319 break
2320 # If levels is changed, recheck to use the vmin,vmax or
2321 # levels-based ranges for histogram plots.
2322 _, levels, _ = _colorbar_map_levels(cube)
2323 logging.debug("levels: %s", levels)
2324 if levels is not None: 2324 ↛ 2314line 2324 didn't jump to line 2314 because the condition on line 2324 was always true
2325 vmin = min(levels)
2326 vmax = max(levels)
2327 logging.debug("Updated vmin, vmax: %s, %s", vmin, vmax)
2328 break
2330 if levels is None:
2331 vmin = min(cb.data.min() for cb in cubes)
2332 vmax = max(cb.data.max() for cb in cubes)
2334 # Make postage stamp plots if stamp_coordinate exists and has more than a
2335 # single point. If single_plot is True:
2336 # -- all postage stamp plots will be plotted in a single plot instead of
2337 # separate postage stamp plots.
2338 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
2339 # produced per single model only
2340 if num_models == 1: 2340 ↛ 2353line 2340 didn't jump to line 2353 because the condition on line 2340 was always true
2341 if ( 2341 ↛ 2345line 2341 didn't jump to line 2345 because the condition on line 2341 was never true
2342 stamp_coordinate in [c.name() for c in cubes[0].coords()]
2343 and cubes[0].coord(stamp_coordinate).shape[0] > 1
2344 ):
2345 if single_plot:
2346 plotting_func = (
2347 _plot_and_save_postage_stamps_in_single_plot_histogram_series
2348 )
2349 else:
2350 plotting_func = _plot_and_save_postage_stamp_histogram_series
2351 cube_iterables = cubes[0].slices_over(sequence_coordinate)
2352 else:
2353 all_points = sorted(
2354 set(
2355 itertools.chain.from_iterable(
2356 cb.coord(sequence_coordinate).points for cb in cubes
2357 )
2358 )
2359 )
2360 all_slices = list(
2361 itertools.chain.from_iterable(
2362 cb.slices_over(sequence_coordinate) for cb in cubes
2363 )
2364 )
2365 # Matched slices (matched by seq coord point; it may happen that
2366 # evaluated models do not cover the same seq coord range, hence matching
2367 # necessary)
2368 cube_iterables = [
2369 iris.cube.CubeList(
2370 s for s in all_slices if s.coord(sequence_coordinate).points[0] == point
2371 )
2372 for point in all_points
2373 ]
2375 plot_index = []
2376 nplot = np.size(cube.coord(sequence_coordinate).points)
2377 # Create a plot for each value of the sequence coordinate. Allowing for
2378 # multiple cubes in a CubeList to be plotted in the same plot for similar
2379 # sequence values. Passing a CubeList into the internal plotting function
2380 # for similar values of the sequence coordinate. cube_slice can be an
2381 # iris.cube.Cube or an iris.cube.CubeList.
2382 for cube_slice in cube_iterables:
2383 single_cube = cube_slice
2384 if isinstance(cube_slice, iris.cube.CubeList): 2384 ↛ 2385line 2384 didn't jump to line 2385 because the condition on line 2384 was never true
2385 single_cube = cube_slice[0]
2387 # Use sequence value so multiple sequences can merge.
2388 sequence_value = single_cube.coord(sequence_coordinate).points[0]
2389 plot_filename = f"{filename.rsplit('.', 1)[0]}_{sequence_value}.png"
2390 coord = single_cube.coord(sequence_coordinate)
2391 # Format the coordinate value in a unit appropriate way.
2392 title = f"{recipe_title}\n [{coord.units.title(coord.points[0])}]"
2393 # Use sequence (e.g. time) bounds if plotting single non-sequence outputs
2394 if nplot == 1 and coord.has_bounds: 2394 ↛ 2395line 2394 didn't jump to line 2395 because the condition on line 2394 was never true
2395 if np.size(coord.bounds) > 1:
2396 title = f"{recipe_title}\n [{coord.units.title(coord.bounds[0][0])} to {coord.units.title(coord.bounds[0][1])}]"
2397 # Do the actual plotting.
2398 plotting_func(
2399 cube_slice,
2400 filename=plot_filename,
2401 stamp_coordinate=stamp_coordinate,
2402 title=title,
2403 vmin=vmin,
2404 vmax=vmax,
2405 )
2406 plot_index.append(plot_filename)
2408 # Add list of plots to plot metadata.
2409 complete_plot_index = _append_to_plot_index(plot_index)
2411 # Make a page to display the plots.
2412 _make_plot_html_page(complete_plot_index)
2414 return cubes