Coverage for src/CSET/operators/plot.py: 80%
854 statements
« prev ^ index » next coverage.py v7.14.2, created at 2026-06-22 10:06 +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 importlib.resources
19import itertools
20import json
21import logging
22import math
23import os
24from typing import Literal
26import cartopy.crs as ccrs
27import cartopy.feature as cfeature
28import iris
29import iris.coords
30import iris.cube
31import iris.exceptions
32import iris.plot as iplt
33import matplotlib as mpl
34import matplotlib.pyplot as plt
35import numpy as np
36from cartopy.mpl.geoaxes import GeoAxes
37from iris.cube import Cube
38from markdown_it import MarkdownIt
39from mpl_toolkits.axes_grid1.inset_locator import inset_axes
41from CSET._common import (
42 filename_slugify,
43 get_recipe_metadata,
44 iter_maybe,
45 render_file,
46 slugify,
47)
48from CSET.operators._colormaps import (
49 colorbar_map_levels,
50 get_model_colors_map,
51)
52from CSET.operators._spectra import DCT_ps
53from CSET.operators._utils import (
54 check_sequence_coordinate,
55 check_single_cube,
56 check_stamp_coordinate,
57 fully_equalise_attributes,
58 get_cube_yxcoordname,
59 get_num_models,
60 is_transect,
61 slice_over_maybe,
62 validate_cube_shape,
63 validate_cubes_coords,
64)
65from CSET.operators.collapse import collapse
66from CSET.operators.misc import _extract_common_time_points
67from CSET.operators.regrid import regrid_onto_cube
69# Use a non-interactive plotting backend.
70mpl.use("agg")
73############################
74# Private helper functions #
75############################
78def _append_to_plot_index(plot_index: list) -> list:
79 """Add plots into the plot index, returning the complete plot index."""
80 with open("meta.json", "r+t", encoding="UTF-8") as fp:
81 fcntl.flock(fp, fcntl.LOCK_EX)
82 fp.seek(0)
83 meta = json.load(fp)
84 complete_plot_index = meta.get("plots", [])
85 complete_plot_index = complete_plot_index + plot_index
86 meta["plots"] = complete_plot_index
87 if os.getenv("CYLC_TASK_CYCLE_POINT") and not bool(
88 os.getenv("DO_CASE_AGGREGATION")
89 ):
90 meta["case_date"] = os.getenv("CYLC_TASK_CYCLE_POINT", "")
91 fp.seek(0)
92 fp.truncate()
93 json.dump(meta, fp, indent=2)
94 return complete_plot_index
97def _make_plot_html_page(plots: list):
98 """Create a HTML page to display a plot image."""
99 # Debug check that plots actually contains some strings.
100 assert isinstance(plots[0], str)
102 # Load HTML template file.
103 operator_files = importlib.resources.files()
104 template_file = operator_files.joinpath("_plot_page_template.html")
106 # Get some metadata.
107 meta = get_recipe_metadata()
108 title = meta.get("title", "Untitled")
109 description = MarkdownIt().render(meta.get("description", "*No description.*"))
111 # Prepare template variables.
112 variables = {
113 "title": title,
114 "description": description,
115 "initial_plot": plots[0],
116 "plots": plots,
117 "title_slug": slugify(title),
118 }
120 # Render template.
121 html = render_file(template_file, **variables)
123 # Save completed HTML.
124 with open("index.html", "wt", encoding="UTF-8") as fp:
125 fp.write(html)
128def _setup_spatial_map(
129 cube: iris.cube.Cube,
130 figure,
131 cmap,
132 grid_size: tuple[int, int] | None = None,
133 subplot: int | None = None,
134):
135 """Define map projections, extent and add coastlines and borderlines for spatial plots.
137 For spatial map plots, a relevant map projection for rotated or non-rotated inputs
138 is specified, and map extent defined based on the input data.
140 Parameters
141 ----------
142 cube: Cube
143 2 dimensional (lat and lon) Cube of the data to plot.
144 figure:
145 Matplotlib Figure object holding all plot elements.
146 cmap:
147 Matplotlib colormap.
148 grid_size: (int, int), optional
149 Size of grid (rows, cols) for subplots if multiple spatial subplots in figure.
150 subplot: int, optional
151 Subplot index if multiple spatial subplots in figure.
153 Returns
154 -------
155 axes:
156 Matplotlib GeoAxes definition.
157 """
158 # Identify min/max plot bounds.
159 try:
160 lat_axis, lon_axis = get_cube_yxcoordname(cube)
161 x1 = np.min(cube.coord(lon_axis).points)
162 x2 = np.max(cube.coord(lon_axis).points)
163 y1 = np.min(cube.coord(lat_axis).points)
164 y2 = np.max(cube.coord(lat_axis).points)
166 # Adjust bounds within +/- 180.0 if x dimension extends beyond half-globe.
167 if np.abs(x2 - x1) > 180.0:
168 x1 = x1 - 180.0
169 x2 = x2 - 180.0
170 logging.debug("Adjusting plot bounds to fit global extent.")
172 # Consider map projection orientation.
173 # Adapting orientation enables plotting across international dateline.
174 # Users can adapt the default central_longitude if alternative projections views.
175 if x2 > 180.0 or x1 < -180.0:
176 central_longitude = 180.0
177 else:
178 central_longitude = 0.0
180 # Define spatial map projection.
181 coord_system = cube.coord(lat_axis).coord_system
182 if isinstance(coord_system, iris.coord_systems.RotatedGeogCS):
183 # Define rotated pole map projection for rotated pole inputs.
184 projection = ccrs.RotatedPole(
185 pole_longitude=coord_system.grid_north_pole_longitude,
186 pole_latitude=coord_system.grid_north_pole_latitude,
187 central_rotated_longitude=central_longitude,
188 )
189 crs = projection
190 elif isinstance(coord_system, iris.coord_systems.TransverseMercator): 190 ↛ 192line 190 didn't jump to line 192 because the condition on line 190 was never true
191 # Define Transverse Mercator projection for TM inputs.
192 projection = ccrs.TransverseMercator(
193 central_longitude=coord_system.longitude_of_central_meridian,
194 central_latitude=coord_system.latitude_of_projection_origin,
195 false_easting=coord_system.false_easting,
196 false_northing=coord_system.false_northing,
197 scale_factor=coord_system.scale_factor_at_central_meridian,
198 )
199 crs = projection
200 else:
201 # Define regular map projection for non-rotated pole inputs.
202 # Alternatives might include e.g. for global model outputs:
203 # projection=ccrs.Robinson(central_longitude=X.y, globe=None)
204 # See also https://scitools.org.uk/cartopy/docs/v0.15/crs/projections.html.
205 projection = ccrs.PlateCarree(central_longitude=central_longitude)
206 crs = ccrs.PlateCarree()
208 # Define axes for plot (or subplot) with required map projection.
209 if subplot is not None:
210 axes = figure.add_subplot(
211 grid_size[0], grid_size[1], subplot, projection=projection
212 )
213 else:
214 axes = figure.add_subplot(projection=projection)
216 # Add coastlines and borderlines if cube contains x and y map coordinates.
217 # Avoid adding lines for masked data or specific fixed ancillary spatial plots.
218 if iris.util.is_masked(cube.data) or any( 218 ↛ 221line 218 didn't jump to line 221 because the condition on line 218 was never true
219 name in cube.name() for name in ["land_", "orography", "altitude"]
220 ):
221 pass
222 else:
223 if cmap.name in ["viridis", "Greys"]:
224 coastcol = "magenta"
225 else:
226 coastcol = "black"
227 logging.debug("Plotting coastlines and borderlines in colour %s.", coastcol)
228 axes.coastlines(resolution="10m", color=coastcol)
229 axes.add_feature(cfeature.BORDERS, edgecolor=coastcol)
231 # Add gridlines.
232 gl = axes.gridlines(
233 alpha=0.3,
234 draw_labels=True,
235 dms=False,
236 x_inline=False,
237 y_inline=False,
238 )
239 gl.top_labels = False
240 gl.right_labels = False
241 if subplot:
242 gl.bottom_labels = False
243 gl.left_labels = False
244 if subplot % grid_size[1] == 1:
245 gl.left_labels = True
246 if subplot > ((grid_size[0] - 1) * grid_size[1]): 246 ↛ 251line 246 didn't jump to line 251 because the condition on line 246 was always true
247 gl.bottom_labels = True
249 # If is lat/lon spatial map, fix extent to keep plot tight.
250 # Specifying crs within set_extent helps ensure only data region is shown.
251 if isinstance(coord_system, iris.coord_systems.GeogCS):
252 axes.set_extent([x1, x2, y1, y2], crs=crs)
254 except ValueError:
255 # Skip if not both x and y map coordinates.
256 axes = figure.gca()
257 pass
259 return axes
262def _get_plot_resolution() -> int:
263 """Get resolution of rasterised plots in pixels per inch."""
264 return get_recipe_metadata().get("plot_resolution", 100)
267def _get_start_end_strings(seq_coord: iris.coords.Coord, use_bounds: bool):
268 """Return title and filename based on start and end points or bounds."""
269 if use_bounds and seq_coord.has_bounds():
270 vals = seq_coord.bounds.flatten()
271 else:
272 vals = seq_coord.points
273 start = seq_coord.units.title(vals[0])
274 end = seq_coord.units.title(vals[-1])
276 if start == end:
277 sequence_title = f"\n [{start}]"
278 sequence_fname = f"_{filename_slugify(start)}"
279 else:
280 sequence_title = f"\n [{start} to {end}]"
281 sequence_fname = f"_{filename_slugify(start)}_{filename_slugify(end)}"
283 # Do not include time if coord set to zero.
284 if (
285 seq_coord.units == "hours since 0001-01-01 00:00:00"
286 and vals[0] == 0
287 and vals[-1] == 0
288 ):
289 sequence_title = ""
290 sequence_fname = ""
292 return sequence_title, sequence_fname
295def _set_title_and_filename(
296 seq_coord: iris.coords.Coord,
297 nplot: int,
298 recipe_title: str,
299 filename: str,
300):
301 """Set plot title and filename based on cube coordinate.
303 Parameters
304 ----------
305 sequence_coordinate: iris.coords.Coord
306 Coordinate about which to make a plot sequence.
307 nplot: int
308 Number of output plots to generate - controls title/naming.
309 recipe_title: str
310 Default plot title, potentially to update.
311 filename: str
312 Input plot filename, potentially to update.
314 Returns
315 -------
316 plot_title: str
317 Output formatted plot title string, based on plotted data.
318 plot_filename: str
319 Output formatted plot filename string.
320 """
321 ndim = seq_coord.ndim
322 npoints = np.size(seq_coord.points)
323 sequence_title = ""
324 sequence_fname = ""
326 # Case 1: Multiple dimension sequence input - list number of aggregated cases
327 # (e.g. aggregation histogram plots)
328 if ndim > 1:
329 ncase = np.shape(seq_coord)[0]
330 sequence_title = f"\n [{ncase} cases]"
331 sequence_fname = f"_{ncase}cases"
333 # Case 2: Single dimension input
334 else:
335 # Single sequence point
336 if npoints == 1:
337 if nplot > 1:
338 # Default labels for sequence inputs
339 sequence_value = seq_coord.units.title(seq_coord.points[0])
340 sequence_title = f"\n [{sequence_value}]"
341 sequence_fname = f"_{filename_slugify(sequence_value)}"
342 else:
343 # Aggregated attribute available where input collapsed over aggregation
344 try:
345 ncase = seq_coord.attributes["number_reference_times"]
346 sequence_title = f"\n [{ncase} cases]"
347 sequence_fname = f"_{ncase}cases"
348 except KeyError:
349 sequence_title, sequence_fname = _get_start_end_strings(
350 seq_coord, use_bounds=seq_coord.has_bounds()
351 )
352 # Multiple sequence (e.g. time) points
353 else:
354 sequence_title, sequence_fname = _get_start_end_strings(
355 seq_coord, use_bounds=False
356 )
358 # Set plot title and filename
359 plot_title = f"{recipe_title}{sequence_title}"
361 # Set plot filename, defaulting to user input if provided.
362 if filename is None:
363 filename = slugify(recipe_title)
364 plot_filename = f"{filename.rsplit('.', 1)[0]}{sequence_fname}.png"
365 else:
366 if nplot > 1:
367 plot_filename = f"{filename.rsplit('.', 1)[0]}{sequence_fname}.png"
368 else:
369 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
371 return plot_title, plot_filename
374def _set_postage_stamp_title(stamp_coord: iris.coords.Coord) -> str:
375 """Control postage stamp plot output titles based on stamp coordinate."""
376 if stamp_coord.name() == "realization":
377 mtitle = "Member"
378 else:
379 mtitle = stamp_coord.name().capitalize()
381 if stamp_coord.name() == "time":
382 mtitle = f"{stamp_coord.units.title(stamp_coord.points[0])}"
383 else:
384 mtitle = f"{mtitle} #{stamp_coord.points[0]}"
386 return mtitle
389def _set_axis_range(cubes):
390 """Get minimum and maximum from levels information."""
391 levels = None
392 for cube in cubes: 392 ↛ 408line 392 didn't jump to line 408 because the loop on line 392 didn't complete
393 # First check if user-specified "auto" range variable.
394 # This maintains the value of levels as None, so proceed.
395 _, levels, _ = colorbar_map_levels(cube, axis="y")
396 if levels is None:
397 break
398 # If levels is changed, recheck to use the vmin,vmax or
399 # levels-based ranges for histogram plots.
400 _, levels, _ = colorbar_map_levels(cube)
401 logging.debug("levels: %s", levels)
402 if levels is not None: 402 ↛ 392line 402 didn't jump to line 392 because the condition on line 402 was always true
403 vmin = min(levels)
404 vmax = max(levels)
405 logging.debug("Updated vmin, vmax: %s, %s", vmin, vmax)
406 break
408 if levels is None:
409 vmin = min(cb.data.min() for cb in cubes)
410 vmax = max(cb.data.max() for cb in cubes)
412 return vmin, vmax
415def _find_matched_slices(cubes, sequence_coordinate):
416 """Identify matched cubes in CubeList by sequence_coordinate values.
418 Ensures common points are compared for multiple cube inputs.
419 """
420 all_points = sorted(
421 set(
422 itertools.chain.from_iterable(
423 cb.coord(sequence_coordinate).points for cb in cubes
424 )
425 )
426 )
427 all_slices = list(
428 itertools.chain.from_iterable(
429 cb.slices_over(sequence_coordinate) for cb in cubes
430 )
431 )
432 # Matched slices (matched by seq coord point; it may happen that
433 # evaluated models do not cover the same seq coord range, hence matching
434 # necessary)
435 cube_iterables = [
436 iris.cube.CubeList(
437 s for s in all_slices if s.coord(sequence_coordinate).points[0] == point
438 )
439 for point in all_points
440 ]
442 return cube_iterables
445def _plot_and_save_spatial_plot(
446 cube: iris.cube.Cube,
447 filename: str,
448 title: str,
449 method: Literal["contourf", "pcolormesh"],
450 overlay_cube: iris.cube.Cube | None = None,
451 contour_cube: iris.cube.Cube | None = None,
452 **kwargs,
453):
454 """Plot and save a spatial plot.
456 Parameters
457 ----------
458 cube: Cube
459 2 dimensional (lat and lon) Cube of the data to plot.
460 filename: str
461 Filename of the plot to write.
462 title: str
463 Plot title.
464 method: "contourf" | "pcolormesh"
465 The plotting method to use.
466 overlay_cube: Cube, optional
467 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
468 contour_cube: Cube, optional
469 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
470 """
471 # Setup plot details, size, resolution, etc.
472 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
474 # Specify the color bar
475 cmap, levels, norm = colorbar_map_levels(cube)
477 # If overplotting, set required colorbars
478 if overlay_cube:
479 over_cmap, over_levels, over_norm = colorbar_map_levels(overlay_cube)
480 if contour_cube:
481 cntr_cmap, cntr_levels, cntr_norm = colorbar_map_levels(contour_cube)
483 # Setup plot map projection, extent and coastlines and borderlines.
484 axes = _setup_spatial_map(cube, fig, cmap)
486 # Set colorscale bounds
487 try:
488 vmin = min(levels)
489 vmax = max(levels)
490 except TypeError:
491 vmin, vmax = None, None
492 # Ensure to use norm and not vmin/vmax if levels are defined.
493 if norm is not None:
494 vmin = None
495 vmax = None
496 logging.debug("Plotting using defined levels.")
498 # Plot the field.
499 if method == "contourf":
500 plot = iplt.contourf(cube, cmap=cmap, levels=levels, norm=norm)
501 elif method == "pcolormesh":
502 plot = iplt.pcolormesh(cube, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
503 else:
504 raise ValueError(f"Unknown plotting method: {method}")
506 # Overplot overlay field, if required
507 if overlay_cube:
508 try:
509 over_vmin = min(over_levels)
510 over_vmax = max(over_levels)
511 except TypeError:
512 over_vmin, over_vmax = None, None
513 if over_norm is not None: 513 ↛ 514line 513 didn't jump to line 514 because the condition on line 513 was never true
514 over_vmin = None
515 over_vmax = None
516 overlay = iplt.pcolormesh(
517 overlay_cube,
518 cmap=over_cmap,
519 norm=over_norm,
520 alpha=0.8,
521 vmin=over_vmin,
522 vmax=over_vmax,
523 )
524 # Overplot contour field, if required, with contour labelling.
525 if contour_cube:
526 contour = iplt.contour(
527 contour_cube,
528 colors="darkgray",
529 levels=cntr_levels,
530 norm=cntr_norm,
531 alpha=0.5,
532 linestyles="--",
533 linewidths=1,
534 )
535 plt.clabel(contour)
537 # Check to see if transect, and if so, adjust y axis.
538 if is_transect(cube):
539 if "pressure" in [coord.name() for coord in cube.coords()]:
540 axes.invert_yaxis()
541 axes.set_yscale("log")
542 axes.set_ylim(1100, 100)
543 # If both model_level_number and level_height exists, iplt can construct
544 # plot as a function of height above orography (NOT sea level).
545 elif {"model_level_number", "level_height"}.issubset( 545 ↛ 550line 545 didn't jump to line 550 because the condition on line 545 was always true
546 {coord.name() for coord in cube.coords()}
547 ):
548 axes.set_yscale("log")
550 axes.set_title(
551 f"{title}\n"
552 f"Start Lat: {cube.attributes['transect_coords'].split('_')[0]}"
553 f" Start Lon: {cube.attributes['transect_coords'].split('_')[1]}"
554 f" End Lat: {cube.attributes['transect_coords'].split('_')[2]}"
555 f" End Lon: {cube.attributes['transect_coords'].split('_')[3]}",
556 fontsize=16,
557 )
559 # Inset code
560 axins = inset_axes(
561 axes,
562 width="20%",
563 height="20%",
564 loc="upper right",
565 axes_class=GeoAxes,
566 axes_kwargs={"map_projection": ccrs.PlateCarree()},
567 )
569 # Slightly transparent to reduce plot blocking.
570 axins.patch.set_alpha(0.4)
572 axins.coastlines(resolution="50m")
573 axins.add_feature(cfeature.BORDERS, linewidth=0.3)
575 SLat, SLon, ELat, ELon = (
576 float(coord) for coord in cube.attributes["transect_coords"].split("_")
577 )
579 # Draw line between them
580 axins.plot(
581 [SLon, ELon], [SLat, ELat], color="black", transform=ccrs.PlateCarree()
582 )
584 # Plot points (note: lon, lat order for Cartopy)
585 axins.plot(SLon, SLat, marker="x", color="green", transform=ccrs.PlateCarree())
586 axins.plot(ELon, ELat, marker="x", color="red", transform=ccrs.PlateCarree())
588 lon_min, lon_max = sorted([SLon, ELon])
589 lat_min, lat_max = sorted([SLat, ELat])
591 # Midpoints
592 lon_mid = (lon_min + lon_max) / 2
593 lat_mid = (lat_min + lat_max) / 2
595 # Maximum half-range
596 half_range = max(lon_max - lon_min, lat_max - lat_min) / 2
597 if half_range == 0: # points identical → provide small default 597 ↛ 601line 597 didn't jump to line 601 because the condition on line 597 was always true
598 half_range = 1
600 # Set square extent
601 axins.set_extent(
602 [
603 lon_mid - half_range,
604 lon_mid + half_range,
605 lat_mid - half_range,
606 lat_mid + half_range,
607 ],
608 crs=ccrs.PlateCarree(),
609 )
611 # Ensure square aspect
612 axins.set_aspect("equal")
614 else:
615 # Add title.
616 axes.set_title(title, fontsize=16)
618 # Adjust padding if spatial plot or transect
619 if is_transect(cube):
620 yinfopad = -0.1
621 ycbarpad = 0.1
622 else:
623 yinfopad = 0.01
624 ycbarpad = 0.042
626 # Add watermark with min/max/mean. Currently not user togglable.
627 # In the bbox dictionary, fc and ec are hex colour codes for grey shade.
628 axes.annotate(
629 f"Min: {np.min(cube.data):.3g} Max: {np.max(cube.data):.3g} Mean: {np.mean(cube.data):.3g}",
630 xy=(0.025, yinfopad),
631 xycoords="axes fraction",
632 xytext=(-5, 5),
633 textcoords="offset points",
634 ha="left",
635 va="bottom",
636 size=11,
637 bbox=dict(boxstyle="round", fc="#cccccc", ec="#808080", alpha=0.9),
638 )
640 # Add secondary colour bar for overlay_cube field if required.
641 if overlay_cube:
642 cbarB = fig.colorbar(
643 overlay, orientation="horizontal", location="bottom", pad=0.0, shrink=0.7
644 )
645 cbarB.set_label(label=f"{overlay_cube.name()} ({overlay_cube.units})", size=14)
646 # add ticks and tick_labels for every levels if less than 20 levels exist
647 if over_levels is not None and len(over_levels) < 20: 647 ↛ 648line 647 didn't jump to line 648 because the condition on line 647 was never true
648 cbarB.set_ticks(over_levels)
649 cbarB.set_ticklabels([f"{level:.2f}" for level in over_levels])
650 if "rainfall" or "snowfall" or "visibility" in overlay_cube.name():
651 cbarB.set_ticklabels([f"{level:.3g}" for level in over_levels])
652 logging.debug("Set secondary colorbar ticks and labels.")
654 # Add main colour bar.
655 cbar = fig.colorbar(
656 plot, orientation="horizontal", location="bottom", pad=ycbarpad, shrink=0.7
657 )
659 cbar.set_label(label=f"{cube.name()} ({cube.units})", size=14)
660 # add ticks and tick_labels for every levels if less than 20 levels exist
661 if levels is not None and len(levels) < 20:
662 cbar.set_ticks(levels)
663 cbar.set_ticklabels([f"{level:.2f}" for level in levels])
664 if "rainfall" or "snowfall" or "visibility" in cube.name(): 664 ↛ 666line 664 didn't jump to line 666 because the condition on line 664 was always true
665 cbar.set_ticklabels([f"{level:.3g}" for level in levels])
666 logging.debug("Set colorbar ticks and labels.")
668 # Save plot.
669 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
670 logging.info("Saved spatial plot to %s", filename)
671 plt.close(fig)
674def _plot_and_save_postage_stamp_spatial_plot(
675 cube: iris.cube.Cube,
676 filename: str,
677 stamp_coordinate: str,
678 title: str,
679 method: Literal["contourf", "pcolormesh"],
680 overlay_cube: iris.cube.Cube | None = None,
681 contour_cube: iris.cube.Cube | None = None,
682 **kwargs,
683):
684 """Plot postage stamp spatial plots from an ensemble.
686 Parameters
687 ----------
688 cube: Cube
689 Iris cube of data to be plotted. It must have the stamp coordinate.
690 filename: str
691 Filename of the plot to write.
692 stamp_coordinate: str
693 Coordinate that becomes different plots.
694 method: "contourf" | "pcolormesh"
695 The plotting method to use.
696 overlay_cube: Cube, optional
697 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
698 contour_cube: Cube, optional
699 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
701 Raises
702 ------
703 ValueError
704 If the cube doesn't have the right dimensions.
705 """
706 # Use the smallest square grid that will fit the members.
707 nmember = len(cube.coord(stamp_coordinate).points)
708 grid_rows = int(math.sqrt(nmember))
709 grid_size = math.ceil(nmember / grid_rows)
711 fig = plt.figure(
712 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
713 )
715 # Specify the color bar
716 cmap, levels, norm = colorbar_map_levels(cube)
717 # If overplotting, set required colorbars
718 if overlay_cube: 718 ↛ 719line 718 didn't jump to line 719 because the condition on line 718 was never true
719 over_cmap, over_levels, over_norm = colorbar_map_levels(overlay_cube)
720 if contour_cube: 720 ↛ 721line 720 didn't jump to line 721 because the condition on line 720 was never true
721 cntr_cmap, cntr_levels, cntr_norm = colorbar_map_levels(contour_cube)
723 # Make a subplot for each member.
724 for member, subplot in zip(
725 cube.slices_over(stamp_coordinate),
726 range(1, grid_size * grid_rows + 1),
727 strict=False,
728 ):
729 # Setup subplot map projection, extent and coastlines and borderlines.
730 axes = _setup_spatial_map(
731 member, fig, cmap, grid_size=(grid_rows, grid_size), subplot=subplot
732 )
733 if method == "contourf":
734 # Filled contour plot of the field.
735 plot = iplt.contourf(member, cmap=cmap, levels=levels, norm=norm)
736 elif method == "pcolormesh":
737 if levels is not None:
738 vmin = min(levels)
739 vmax = max(levels)
740 else:
741 raise TypeError("Unknown vmin and vmax range.")
742 vmin, vmax = None, None
743 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
744 # if levels are defined.
745 if norm is not None: 745 ↛ 746line 745 didn't jump to line 746 because the condition on line 745 was never true
746 vmin = None
747 vmax = None
748 # pcolormesh plot of the field.
749 plot = iplt.pcolormesh(member, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
750 else:
751 raise ValueError(f"Unknown plotting method: {method}")
753 # Overplot overlay field, if required
754 if overlay_cube: 754 ↛ 755line 754 didn't jump to line 755 because the condition on line 754 was never true
755 try:
756 over_vmin = min(over_levels)
757 over_vmax = max(over_levels)
758 except TypeError:
759 over_vmin, over_vmax = None, None
760 if over_norm is not None:
761 over_vmin = None
762 over_vmax = None
763 iplt.pcolormesh(
764 overlay_cube[member.coord(stamp_coordinate).points[0]],
765 cmap=over_cmap,
766 norm=over_norm,
767 alpha=0.6,
768 vmin=over_vmin,
769 vmax=over_vmax,
770 )
771 # Overplot contour field, if required
772 if contour_cube: 772 ↛ 773line 772 didn't jump to line 773 because the condition on line 772 was never true
773 iplt.contour(
774 contour_cube[member.coord(stamp_coordinate).points[0]],
775 colors="darkgray",
776 levels=cntr_levels,
777 norm=cntr_norm,
778 alpha=0.6,
779 linestyles="--",
780 linewidths=1,
781 )
782 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
783 axes.set_title(f"{mtitle}")
785 # Put the shared colorbar in its own axes.
786 colorbar_axes = fig.add_axes([0.15, 0.05, 0.7, 0.03])
787 colorbar = fig.colorbar(
788 plot, colorbar_axes, orientation="horizontal", pad=0.042, shrink=0.7
789 )
790 colorbar.set_label(f"{cube.name()} ({cube.units})", size=14)
792 # Overall figure title.
793 fig.suptitle(title, fontsize=16)
795 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
796 logging.info("Saved contour postage stamp plot to %s", filename)
797 plt.close(fig)
800def _plot_and_save_line_series(
801 cubes: iris.cube.CubeList,
802 coords: list[iris.coords.Coord],
803 ensemble_coord: str,
804 filename: str,
805 title: str,
806 **kwargs,
807):
808 """Plot and save a 1D line series.
810 Parameters
811 ----------
812 cubes: Cube or CubeList
813 Cube or CubeList containing the cubes to plot on the y-axis.
814 coords: list[Coord]
815 Coordinates to plot on the x-axis, one per cube.
816 ensemble_coord: str
817 Ensemble coordinate in the cube.
818 filename: str
819 Filename of the plot to write.
820 title: str
821 Plot title.
822 """
823 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
825 model_colors_map = get_model_colors_map(cubes)
827 # Store min/max ranges.
828 y_levels = []
830 # Check match-up across sequence coords gives consistent sizes
831 validate_cubes_coords(cubes, coords)
833 for cube, coord in zip(cubes, coords, strict=True):
834 label = None
835 color = "black"
836 if model_colors_map:
837 label = cube.attributes.get("model_name")
838 color = model_colors_map.get(label)
839 for cube_slice in cube.slices_over(ensemble_coord):
840 # Label with (control) if part of an ensemble or not otherwise.
841 if cube_slice.coord(ensemble_coord).points == [0]:
842 iplt.plot(
843 coord,
844 cube_slice,
845 color=color,
846 marker="o",
847 ls="-",
848 lw=3,
849 label=f"{label} (control)"
850 if len(cube.coord(ensemble_coord).points) > 1
851 else label,
852 )
853 # Label with (perturbed) if part of an ensemble and not the control.
854 else:
855 iplt.plot(
856 coord,
857 cube_slice,
858 color=color,
859 ls="-",
860 lw=1.5,
861 alpha=0.75,
862 label=f"{label} (member)",
863 )
865 # Calculate the global min/max if multiple cubes are given.
866 _, levels, _ = colorbar_map_levels(cube, axis="y")
867 if levels is not None: 867 ↛ 868line 867 didn't jump to line 868 because the condition on line 867 was never true
868 y_levels.append(min(levels))
869 y_levels.append(max(levels))
871 # Get the current axes.
872 ax = plt.gca()
874 # Add some labels and tweak the style.
875 # check if cubes[0] works for single cube if not CubeList
876 if coords[0].name() == "time":
877 ax.set_xlabel(f"{coords[0].name()}", fontsize=14)
878 else:
879 ax.set_xlabel(f"{coords[0].name()} / {coords[0].units}", fontsize=14)
880 ax.set_ylabel(f"{cubes[0].name()} / {cubes[0].units}", fontsize=14)
881 ax.set_title(title, fontsize=16)
883 ax.ticklabel_format(axis="y", useOffset=False)
884 ax.tick_params(axis="x", labelrotation=15)
885 ax.tick_params(axis="both", labelsize=12)
887 # Set y limits to global min and max, autoscale if colorbar doesn't exist.
888 if y_levels: 888 ↛ 889line 888 didn't jump to line 889 because the condition on line 888 was never true
889 ax.set_ylim(min(y_levels), max(y_levels))
890 # Add zero line.
891 if min(y_levels) < 0.0 and max(y_levels) > 0.0:
892 ax.axhline(y=0, xmin=0, xmax=1, ls="-", color="grey", lw=2)
893 logging.debug(
894 "Line plot with y-axis limits %s-%s", min(y_levels), max(y_levels)
895 )
896 else:
897 ax.autoscale()
899 # Add gridlines
900 ax.grid(linestyle="--", color="grey", linewidth=1)
901 # Ientify unique labels for legend
902 handles = list(
903 {
904 label: handle
905 for (handle, label) in zip(*ax.get_legend_handles_labels(), strict=True)
906 }.values()
907 )
908 ax.legend(handles=handles, loc="best", ncol=1, frameon=False, fontsize=16)
910 # Save plot.
911 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
912 logging.info("Saved line plot to %s", filename)
913 plt.close(fig)
916def _plot_and_save_vertical_line_series(
917 cubes: iris.cube.CubeList,
918 coords: list[iris.coords.Coord],
919 ensemble_coord: str,
920 filename: str,
921 series_coordinate: str,
922 title: str,
923 vmin: float,
924 vmax: float,
925 **kwargs,
926):
927 """Plot and save a 1D line series in vertical.
929 Parameters
930 ----------
931 cubes: CubeList
932 1 dimensional Cube or CubeList of the data to plot on x-axis.
933 coord: list[Coord]
934 Coordinates to plot on the y-axis, one per cube.
935 ensemble_coord: str
936 Ensemble coordinate in the cube.
937 filename: str
938 Filename of the plot to write.
939 series_coordinate: str
940 Coordinate to use as vertical axis.
941 title: str
942 Plot title.
943 vmin: float
944 Minimum value for the x-axis.
945 vmax: float
946 Maximum value for the x-axis.
947 """
948 # plot the vertical pressure axis using log scale
949 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
951 model_colors_map = get_model_colors_map(cubes)
953 # Check match-up across sequence coords gives consistent sizes
954 validate_cubes_coords(cubes, coords)
956 for cube, coord in zip(cubes, coords, strict=True):
957 label = None
958 color = "black"
959 if model_colors_map: 959 ↛ 960line 959 didn't jump to line 960 because the condition on line 959 was never true
960 label = cube.attributes.get("model_name")
961 color = model_colors_map.get(label)
963 for cube_slice in cube.slices_over(ensemble_coord):
964 # If ensemble data given plot control member with (control)
965 # unless single forecast.
966 if cube_slice.coord(ensemble_coord).points == [0]:
967 iplt.plot(
968 cube_slice,
969 coord,
970 color=color,
971 marker="o",
972 ls="-",
973 lw=3,
974 label=f"{label} (control)"
975 if len(cube.coord(ensemble_coord).points) > 1
976 else label,
977 )
978 # If ensemble data given plot perturbed members with (perturbed).
979 else:
980 iplt.plot(
981 cube_slice,
982 coord,
983 color=color,
984 ls="-",
985 lw=1.5,
986 alpha=0.75,
987 label=f"{label} (member)",
988 )
990 # Get the current axis
991 ax = plt.gca()
993 # Special handling for pressure level data.
994 if series_coordinate == "pressure": 994 ↛ 1016line 994 didn't jump to line 1016 because the condition on line 994 was always true
995 # Invert y-axis and set to log scale.
996 ax.invert_yaxis()
997 ax.set_yscale("log")
999 # Define y-ticks and labels for pressure log axis.
1000 y_tick_labels = [
1001 "1000",
1002 "850",
1003 "700",
1004 "500",
1005 "300",
1006 "200",
1007 "100",
1008 ]
1009 y_ticks = [1000, 850, 700, 500, 300, 200, 100]
1011 # Set y-axis limits and ticks.
1012 ax.set_ylim(1100, 100)
1014 # Test if series_coordinate is model level data. The UM data uses
1015 # model_level_number and lfric uses full_levels as coordinate.
1016 elif series_coordinate in ("model_level_number", "full_levels", "half_levels"):
1017 # Define y-ticks and labels for vertical axis.
1018 y_ticks = iter_maybe(cubes)[0].coord(series_coordinate).points
1019 y_tick_labels = [str(int(i)) for i in y_ticks]
1020 ax.set_ylim(min(y_ticks), max(y_ticks))
1022 ax.set_yticks(y_ticks)
1023 ax.set_yticklabels(y_tick_labels)
1025 # Set x-axis limits.
1026 ax.set_xlim(vmin, vmax)
1027 # Mark y=0 if present in plot.
1028 if vmin < 0.0 and vmax > 0.0: 1028 ↛ 1029line 1028 didn't jump to line 1029 because the condition on line 1028 was never true
1029 ax.axvline(x=0, ymin=0, ymax=1, ls="-", color="grey", lw=2)
1031 # Add some labels and tweak the style.
1032 ax.set_ylabel(f"{coord.name()} / {coord.units}", fontsize=14)
1033 ax.set_xlabel(
1034 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
1035 )
1036 ax.set_title(title, fontsize=16)
1037 ax.ticklabel_format(axis="x")
1038 ax.tick_params(axis="y")
1039 ax.tick_params(axis="both", labelsize=12)
1041 # Add gridlines
1042 ax.grid(linestyle="--", color="grey", linewidth=1)
1043 # Ientify unique labels for legend
1044 handles = list(
1045 {
1046 label: handle
1047 for (handle, label) in zip(*ax.get_legend_handles_labels(), strict=True)
1048 }.values()
1049 )
1050 ax.legend(handles=handles, loc="best", ncol=1, frameon=False, fontsize=16)
1052 # Save plot.
1053 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1054 logging.info("Saved line plot to %s", filename)
1055 plt.close(fig)
1058def _plot_and_save_scatter_plot(
1059 cube_x: iris.cube.Cube | iris.cube.CubeList,
1060 cube_y: iris.cube.Cube | iris.cube.CubeList,
1061 filename: str,
1062 title: str,
1063 one_to_one: bool,
1064 model_names: list[str] = None,
1065 **kwargs,
1066):
1067 """Plot and save a 2D scatter plot.
1069 Parameters
1070 ----------
1071 cube_x: Cube | CubeList
1072 1 dimensional Cube or CubeList of the data to plot on x-axis.
1073 cube_y: Cube | CubeList
1074 1 dimensional Cube or CubeList of the data to plot on y-axis.
1075 filename: str
1076 Filename of the plot to write.
1077 title: str
1078 Plot title.
1079 one_to_one: bool
1080 Whether a 1:1 line is plotted.
1081 """
1082 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1083 # plot the cube_x and cube_y 1D fields as a scatter plot. If they are CubeLists this ensures
1084 # to pair each cube from cube_x with the corresponding cube from cube_y, allowing to iterate
1085 # over the pairs simultaneously.
1087 # Ensure cube_x and cube_y are iterable
1088 cube_x_iterable = iter_maybe(cube_x)
1089 cube_y_iterable = iter_maybe(cube_y)
1091 for cube_x_iter, cube_y_iter in zip(cube_x_iterable, cube_y_iterable, strict=True):
1092 iplt.scatter(cube_x_iter, cube_y_iter)
1093 if one_to_one is True:
1094 plt.plot(
1095 [
1096 np.nanmin([np.nanmin(cube_y.data), np.nanmin(cube_x.data)]),
1097 np.nanmax([np.nanmax(cube_y.data), np.nanmax(cube_x.data)]),
1098 ],
1099 [
1100 np.nanmin([np.nanmin(cube_y.data), np.nanmin(cube_x.data)]),
1101 np.nanmax([np.nanmax(cube_y.data), np.nanmax(cube_x.data)]),
1102 ],
1103 "k",
1104 linestyle="--",
1105 )
1106 ax = plt.gca()
1108 # Add some labels and tweak the style.
1109 if model_names is None:
1110 ax.set_xlabel(f"{cube_x[0].name()} / {cube_x[0].units}", fontsize=14)
1111 ax.set_ylabel(f"{cube_y[0].name()} / {cube_y[0].units}", fontsize=14)
1112 else:
1113 # Add the model names, these should be order of base (x) and other (y).
1114 ax.set_xlabel(
1115 f"{model_names[0]}_{cube_x[0].name()} / {cube_x[0].units}", fontsize=14
1116 )
1117 ax.set_ylabel(
1118 f"{model_names[1]}_{cube_y[0].name()} / {cube_y[0].units}", fontsize=14
1119 )
1120 ax.set_title(title, fontsize=16)
1121 ax.ticklabel_format(axis="y", useOffset=False)
1122 ax.tick_params(axis="x", labelrotation=15)
1123 ax.tick_params(axis="both", labelsize=12)
1124 ax.autoscale()
1126 # Save plot.
1127 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1128 logging.info("Saved scatter plot to %s", filename)
1129 plt.close(fig)
1132def _plot_and_save_vector_plot(
1133 cube_u: iris.cube.Cube,
1134 cube_v: iris.cube.Cube,
1135 filename: str,
1136 title: str,
1137 method: Literal["contourf", "pcolormesh"],
1138 **kwargs,
1139):
1140 """Plot and save a 2D vector plot.
1142 Parameters
1143 ----------
1144 cube_u: Cube
1145 2 dimensional Cube of u component of the data.
1146 cube_v: Cube
1147 2 dimensional Cube of v component of the data.
1148 filename: str
1149 Filename of the plot to write.
1150 title: str
1151 Plot title.
1152 """
1153 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1155 # Create a cube containing the magnitude of the vector field.
1156 cube_vec_mag = (cube_u**2 + cube_v**2) ** 0.5
1157 cube_vec_mag.rename(f"{cube_u.name()}_{cube_v.name()}_magnitude")
1159 # Specify the color bar
1160 cmap, levels, norm = colorbar_map_levels(cube_vec_mag)
1162 # Setup plot map projection, extent and coastlines and borderlines.
1163 axes = _setup_spatial_map(cube_vec_mag, fig, cmap)
1165 if method == "contourf":
1166 # Filled contour plot of the field.
1167 plot = iplt.contourf(cube_vec_mag, cmap=cmap, levels=levels, norm=norm)
1168 elif method == "pcolormesh":
1169 try:
1170 vmin = min(levels)
1171 vmax = max(levels)
1172 except TypeError:
1173 vmin, vmax = None, None
1174 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
1175 # if levels are defined.
1176 if norm is not None:
1177 vmin = None
1178 vmax = None
1179 plot = iplt.pcolormesh(cube_vec_mag, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
1180 else:
1181 raise ValueError(f"Unknown plotting method: {method}")
1183 # Check to see if transect, and if so, adjust y axis.
1184 if is_transect(cube_vec_mag):
1185 if "pressure" in [coord.name() for coord in cube_vec_mag.coords()]:
1186 axes.invert_yaxis()
1187 axes.set_yscale("log")
1188 axes.set_ylim(1100, 100)
1189 # If both model_level_number and level_height exists, iplt can construct
1190 # plot as a function of height above orography (NOT sea level).
1191 elif {"model_level_number", "level_height"}.issubset(
1192 {coord.name() for coord in cube_vec_mag.coords()}
1193 ):
1194 axes.set_yscale("log")
1196 axes.set_title(
1197 f"{title}\n"
1198 f"Start Lat: {cube_vec_mag.attributes['transect_coords'].split('_')[0]}"
1199 f" Start Lon: {cube_vec_mag.attributes['transect_coords'].split('_')[1]}"
1200 f" End Lat: {cube_vec_mag.attributes['transect_coords'].split('_')[2]}"
1201 f" End Lon: {cube_vec_mag.attributes['transect_coords'].split('_')[3]}",
1202 fontsize=16,
1203 )
1205 else:
1206 # Add title.
1207 axes.set_title(title, fontsize=16)
1209 # Add watermark with min/max/mean. Currently not user togglable.
1210 # In the bbox dictionary, fc and ec are hex colour codes for grey shade.
1211 axes.annotate(
1212 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}",
1213 xy=(0.05, -0.05),
1214 xycoords="axes fraction",
1215 xytext=(-5, 5),
1216 textcoords="offset points",
1217 ha="right",
1218 va="bottom",
1219 size=11,
1220 bbox=dict(boxstyle="round", fc="#cccccc", ec="#808080", alpha=0.9),
1221 )
1223 # Add colour bar.
1224 cbar = fig.colorbar(plot, orientation="horizontal", pad=0.042, shrink=0.7)
1225 cbar.set_label(label=f"{cube_vec_mag.name()} ({cube_vec_mag.units})", size=14)
1226 # add ticks and tick_labels for every levels if less than 20 levels exist
1227 if levels is not None and len(levels) < 20:
1228 cbar.set_ticks(levels)
1229 cbar.set_ticklabels([f"{level:.1f}" for level in levels])
1231 # 30 barbs along the longest axis of the plot, or a barb per point for data
1232 # with less than 30 points.
1233 step = max(max(cube_u.shape) // 30, 1)
1234 iplt.quiver(cube_u[::step, ::step], cube_v[::step, ::step], pivot="middle")
1236 # Save plot.
1237 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1238 logging.info("Saved vector plot to %s", filename)
1239 plt.close(fig)
1242def _plot_and_save_histogram_series(
1243 cubes: iris.cube.Cube | iris.cube.CubeList,
1244 filename: str,
1245 title: str,
1246 vmin: float,
1247 vmax: float,
1248 **kwargs,
1249):
1250 """Plot and save a histogram series.
1252 Parameters
1253 ----------
1254 cubes: Cube or CubeList
1255 2 dimensional Cube or CubeList of the data to plot as histogram.
1256 filename: str
1257 Filename of the plot to write.
1258 title: str
1259 Plot title.
1260 vmin: float
1261 minimum for colorbar
1262 vmax: float
1263 maximum for colorbar
1264 """
1265 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1266 ax = plt.gca()
1268 model_colors_map = get_model_colors_map(cubes)
1270 # Set default that histograms will produce probability density function
1271 # at each bin (integral over range sums to 1).
1272 density = True
1274 for cube in iter_maybe(cubes):
1275 # Easier to check title (where var name originates)
1276 # than seeing if long names exist etc.
1277 # Exception case, where distribution better fits log scales/bins.
1278 if "surface_microphysical" in title:
1279 if "amount" in title: 1279 ↛ 1281line 1279 didn't jump to line 1281 because the condition on line 1279 was never true
1280 # Compute histogram following Klingaman et al. (2017): ASoP
1281 bin2 = np.exp(np.log(0.02) + 0.1 * np.linspace(0, 99, 100))
1282 bins = np.pad(bin2, (1, 0), "constant", constant_values=0)
1283 density = False
1284 else:
1285 bins = 10.0 ** (
1286 np.arange(-10, 27, 1) / 10.0
1287 ) # Suggestion from RMED toolbox.
1288 bins = np.insert(bins, 0, 0)
1289 ax.set_yscale("log")
1290 vmin = bins[1]
1291 vmax = bins[-1] # Manually set vmin/vmax to override json derived value.
1292 ax.set_xscale("log")
1293 elif "lightning" in title:
1294 bins = [0, 1, 2, 3, 4, 5]
1295 else:
1296 bins = np.linspace(vmin, vmax, 51)
1297 logging.debug(
1298 "Plotting histogram with %s bins %s - %s.",
1299 np.size(bins),
1300 np.min(bins),
1301 np.max(bins),
1302 )
1304 # Reshape cube data into a single array to allow for a single histogram.
1305 # Otherwise we plot xdim histograms stacked.
1306 cube_data_1d = (cube.data).flatten()
1308 label = None
1309 color = "black"
1310 if model_colors_map: 1310 ↛ 1311line 1310 didn't jump to line 1311 because the condition on line 1310 was never true
1311 label = cube.attributes.get("model_name")
1312 color = model_colors_map[label]
1313 x, y = np.histogram(cube_data_1d, bins=bins, density=density)
1315 # Compute area under curve.
1316 if "surface_microphysical" in title and "amount" in title: 1316 ↛ 1317line 1316 didn't jump to line 1317 because the condition on line 1316 was never true
1317 bin_mean = (bins[:-1] + bins[1:]) / 2.0
1318 x = x * bin_mean / x.sum()
1319 x = x[1:]
1320 y = y[1:]
1322 ax.plot(
1323 y[:-1], x, color=color, linewidth=3, marker="o", markersize=6, label=label
1324 )
1326 # Add some labels and tweak the style.
1327 ax.set_title(title, fontsize=16)
1328 ax.set_xlabel(
1329 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
1330 )
1331 ax.set_ylabel("Normalised probability density", fontsize=14)
1332 if "surface_microphysical" in title and "amount" in title: 1332 ↛ 1333line 1332 didn't jump to line 1333 because the condition on line 1332 was never true
1333 ax.set_ylabel(
1334 f"Contribution to mean ({iter_maybe(cubes)[0].units})", fontsize=14
1335 )
1336 try:
1337 ax.set_xlim(vmin, vmax)
1338 except ValueError:
1339 pass
1340 ax.tick_params(axis="both", labelsize=12)
1342 # Overlay grid-lines onto histogram plot.
1343 ax.grid(linestyle="--", color="grey", linewidth=1)
1344 if model_colors_map: 1344 ↛ 1345line 1344 didn't jump to line 1345 because the condition on line 1344 was never true
1345 ax.legend(loc="best", ncol=1, frameon=False, fontsize=16)
1347 # Save plot.
1348 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1349 logging.info("Saved histogram plot to %s", filename)
1350 plt.close(fig)
1353def _plot_and_save_postage_stamp_histogram_series(
1354 cube: iris.cube.Cube,
1355 filename: str,
1356 title: str,
1357 stamp_coordinate: str,
1358 vmin: float,
1359 vmax: float,
1360 **kwargs,
1361):
1362 """Plot and save postage (ensemble members) stamps for a histogram series.
1364 Parameters
1365 ----------
1366 cube: Cube
1367 2 dimensional Cube of the data to plot as histogram.
1368 filename: str
1369 Filename of the plot to write.
1370 title: str
1371 Plot title.
1372 stamp_coordinate: str
1373 Coordinate that becomes different plots.
1374 vmin: float
1375 minimum for pdf x-axis
1376 vmax: float
1377 maximum for pdf x-axis
1378 """
1379 # Use the smallest square grid that will fit the members.
1380 nmember = len(cube.coord(stamp_coordinate).points)
1381 grid_rows = int(math.sqrt(nmember))
1382 grid_size = math.ceil(nmember / grid_rows)
1384 fig = plt.figure(
1385 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
1386 )
1387 # Make a subplot for each member.
1388 for member, subplot in zip(
1389 cube.slices_over(stamp_coordinate),
1390 range(1, grid_size * grid_rows + 1),
1391 strict=False,
1392 ):
1393 # Implicit interface is much easier here, due to needing to have the
1394 # cartopy GeoAxes generated.
1395 plt.subplot(grid_rows, grid_size, subplot)
1396 # Reshape cube data into a single array to allow for a single histogram.
1397 # Otherwise we plot xdim histograms stacked.
1398 member_data_1d = (member.data).flatten()
1399 plt.hist(member_data_1d, density=True, stacked=True)
1400 axes = plt.gca()
1401 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
1402 axes.set_title(f"{mtitle}")
1403 axes.set_xlim(vmin, vmax)
1405 # Overall figure title.
1406 fig.suptitle(title, fontsize=16)
1408 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1409 logging.info("Saved histogram postage stamp plot to %s", filename)
1410 plt.close(fig)
1413def _plot_and_save_postage_stamps_in_single_plot_histogram_series(
1414 cube: iris.cube.Cube,
1415 filename: str,
1416 title: str,
1417 stamp_coordinate: str,
1418 vmin: float,
1419 vmax: float,
1420 **kwargs,
1421):
1422 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1423 ax.set_title(title, fontsize=16)
1424 ax.set_xlim(vmin, vmax)
1425 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1426 ax.set_ylabel("normalised probability density", fontsize=14)
1427 # Loop over all slices along the stamp_coordinate
1428 for member in cube.slices_over(stamp_coordinate):
1429 # Flatten the member data to 1D
1430 member_data_1d = member.data.flatten()
1431 # Plot the histogram using plt.hist
1432 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
1433 plt.hist(
1434 member_data_1d,
1435 density=True,
1436 stacked=True,
1437 label=f"{mtitle}",
1438 )
1440 # Add a legend
1441 ax.legend(fontsize=16)
1443 # Save the figure to a file
1444 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1445 logging.info("Saved histogram postage stamp plot to %s", filename)
1447 # Close the figure
1448 plt.close(fig)
1451def _plot_and_save_scattermap_plot(
1452 cube: iris.cube.Cube, filename: str, title: str, projection=None, **kwargs
1453):
1454 """Plot and save a geographical scatter plot.
1456 Parameters
1457 ----------
1458 cube: Cube
1459 1 dimensional Cube of the data points with auxiliary latitude and
1460 longitude coordinates,
1461 filename: str
1462 Filename of the plot to write.
1463 title: str
1464 Plot title.
1465 projection: str
1466 Mapping projection to be used by cartopy.
1467 """
1468 # Setup plot details, size, resolution, etc.
1469 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1470 if projection is not None:
1471 # Apart from the default, the only projection we currently support is
1472 # a stereographic projection over the North Pole.
1473 if projection == "NP_Stereo":
1474 axes = plt.axes(projection=ccrs.NorthPolarStereo(central_longitude=0.0))
1475 else:
1476 raise ValueError(f"Unknown projection: {projection}")
1477 else:
1478 axes = plt.axes(projection=ccrs.PlateCarree())
1480 # Scatter plot of the field. The marker size is chosen to give
1481 # symbols that decrease in size as the number of observations
1482 # increases, although the fraction of the figure covered by
1483 # symbols increases roughly as N^(1/2), disregarding overlaps,
1484 # and has been selected for the default figure size of (10, 10).
1485 # Should this be changed, the marker size should be adjusted in
1486 # proportion to the area of the figure.
1487 mrk_size = int(np.sqrt(2500000.0 / len(cube.data)))
1488 klon = None
1489 klat = None
1490 for kc in range(len(cube.aux_coords)):
1491 if cube.aux_coords[kc].standard_name == "latitude":
1492 klat = kc
1493 elif cube.aux_coords[kc].standard_name == "longitude":
1494 klon = kc
1495 scatter_map = iplt.scatter(
1496 cube.aux_coords[klon],
1497 cube.aux_coords[klat],
1498 c=cube.data[:],
1499 s=mrk_size,
1500 cmap="jet",
1501 edgecolors="k",
1502 )
1504 # Add coastlines and borderlines.
1505 try:
1506 axes.coastlines(resolution="10m")
1507 axes.add_feature(cfeature.BORDERS)
1508 except AttributeError:
1509 pass
1511 # Add title.
1512 axes.set_title(title, fontsize=16)
1514 # Add colour bar.
1515 cbar = fig.colorbar(scatter_map)
1516 cbar.set_label(label=f"{cube.name()} ({cube.units})", size=20)
1518 # Save plot.
1519 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1520 logging.info("Saved geographical scatter plot to %s", filename)
1521 plt.close(fig)
1524def _plot_and_save_power_spectrum_series(
1525 cubes: iris.cube.Cube | iris.cube.CubeList,
1526 filename: str,
1527 title: str,
1528 **kwargs,
1529):
1530 """Plot and save a power spectrum series.
1532 Parameters
1533 ----------
1534 cubes: Cube or CubeList
1535 2 dimensional Cube or CubeList of the data to plot as power spectrum.
1536 filename: str
1537 Filename of the plot to write.
1538 title: str
1539 Plot title.
1540 """
1541 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1542 ax = plt.gca()
1544 model_colors_map = get_model_colors_map(cubes)
1546 for cube in iter_maybe(cubes):
1547 # Calculate power spectrum
1549 # Extract time coordinate and convert to datetime
1550 time_coord = cube.coord("time")
1551 time_points = time_coord.units.num2date(time_coord.points)
1553 # Choose one time point (e.g., the first one)
1554 target_time = time_points[0]
1556 # Bind target_time inside the lambda using a default argument
1557 time_constraint = iris.Constraint(
1558 time=lambda cell, target_time=target_time: cell.point == target_time
1559 )
1561 cube = cube.extract(time_constraint)
1563 if cube.ndim == 2:
1564 cube_3d = cube.data[np.newaxis, :, :]
1565 logging.debug("Adding in new axis for a 2 dimensional cube.")
1566 elif cube.ndim == 3: 1566 ↛ 1567line 1566 didn't jump to line 1567 because the condition on line 1566 was never true
1567 cube_3d = cube.data
1568 else:
1569 raise ValueError("Cube dimensions unsuitable for power spectra code")
1570 raise ValueError(
1571 f"Cube is {cube.ndim} dimensional. Cube should be 2 or 3 dimensional."
1572 )
1574 # Calculate spectra
1575 ps_array = DCT_ps(cube_3d)
1577 ps_cube = iris.cube.Cube(
1578 ps_array,
1579 long_name="power_spectra",
1580 )
1582 ps_cube.attributes["model_name"] = cube.attributes.get("model_name")
1584 # Create a frequency/wavelength array for coordinate
1585 ps_len = ps_cube.data.shape[1]
1586 freqs = np.arange(1, ps_len + 1)
1587 freq_coord = iris.coords.DimCoord(freqs, long_name="frequency", units="1")
1589 # Convert datetime to numeric time using original units
1590 numeric_time = time_coord.units.date2num(time_points)
1591 # Create a new DimCoord with numeric time
1592 new_time_coord = iris.coords.DimCoord(
1593 numeric_time, standard_name="time", units=time_coord.units
1594 )
1596 # Add time and frequency coordinate to spectra cube.
1597 ps_cube.add_dim_coord(new_time_coord.copy(), 0)
1598 ps_cube.add_dim_coord(freq_coord.copy(), 1)
1600 # Extract data from the cube
1601 frequency = ps_cube.coord("frequency").points
1602 power_spectrum = ps_cube.data
1604 label = None
1605 color = "black"
1606 if model_colors_map: 1606 ↛ 1607line 1606 didn't jump to line 1607 because the condition on line 1606 was never true
1607 label = ps_cube.attributes.get("model_name")
1608 color = model_colors_map[label]
1609 ax.plot(frequency, power_spectrum[0], color=color, label=label)
1611 # Add some labels and tweak the style.
1612 ax.set_title(title, fontsize=16)
1613 ax.set_xlabel("Wavenumber", fontsize=14)
1614 ax.set_ylabel("Power", fontsize=14)
1615 ax.tick_params(axis="both", labelsize=12)
1617 # Set log-log scale
1618 ax.set_xscale("log")
1619 ax.set_yscale("log")
1621 # Overlay grid-lines onto power spectrum plot.
1622 ax.grid(linestyle="--", color="grey", linewidth=1)
1623 if model_colors_map: 1623 ↛ 1624line 1623 didn't jump to line 1624 because the condition on line 1623 was never true
1624 ax.legend(loc="best", ncol=1, frameon=False, fontsize=16)
1626 # Save plot.
1627 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1628 logging.info("Saved power spectrum plot to %s", filename)
1629 plt.close(fig)
1632def _plot_and_save_postage_stamp_power_spectrum_series(
1633 cube: iris.cube.Cube,
1634 filename: str,
1635 title: str,
1636 stamp_coordinate: str,
1637 **kwargs,
1638):
1639 """Plot and save postage (ensemble members) stamps for a power spectrum series.
1641 Parameters
1642 ----------
1643 cube: Cube
1644 2 dimensional Cube of the data to plot as power spectrum.
1645 filename: str
1646 Filename of the plot to write.
1647 title: str
1648 Plot title.
1649 stamp_coordinate: str
1650 Coordinate that becomes different plots.
1651 """
1652 # Use the smallest square grid that will fit the members.
1653 nmember = len(cube.coord(stamp_coordinate).points)
1654 grid_rows = int(math.sqrt(nmember))
1655 grid_size = math.ceil(nmember / grid_rows)
1657 fig = plt.figure(
1658 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
1659 )
1661 # Make a subplot for each member.
1662 for member, subplot in zip(
1663 cube.slices_over(stamp_coordinate),
1664 range(1, grid_size * grid_rows + 1),
1665 strict=False,
1666 ):
1667 # Implicit interface is much easier here, due to needing to have the
1668 # cartopy GeoAxes generated.
1669 plt.subplot(grid_rows, grid_size, subplot)
1671 frequency = member.coord("frequency").points
1673 axes = plt.gca()
1674 axes.plot(frequency, member.data)
1675 axes.set_title(f"Member #{member.coord(stamp_coordinate).points[0]}")
1677 # Overall figure title.
1678 fig.suptitle(title, fontsize=16)
1680 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1681 logging.info("Saved power spectra postage stamp plot to %s", filename)
1682 plt.close(fig)
1685def _plot_and_save_postage_stamps_in_single_plot_power_spectrum_series(
1686 cube: iris.cube.Cube,
1687 filename: str,
1688 title: str,
1689 stamp_coordinate: str,
1690 **kwargs,
1691):
1692 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1693 ax.set_title(title, fontsize=16)
1694 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1695 ax.set_ylabel("Power", fontsize=14)
1696 # Loop over all slices along the stamp_coordinate
1697 for member in cube.slices_over(stamp_coordinate):
1698 frequency = member.coord("frequency").points
1699 ax.plot(
1700 frequency,
1701 member.data,
1702 label=f"Member #{member.coord(stamp_coordinate).points[0]}",
1703 )
1705 # Add a legend
1706 ax.legend(fontsize=16)
1708 # Save the figure to a file
1709 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1710 logging.info("Saved power spectra plot to %s", filename)
1712 # Close the figure
1713 plt.close(fig)
1716def _spatial_plot(
1717 method: Literal["contourf", "pcolormesh"],
1718 cube: iris.cube.Cube,
1719 filename: str | None,
1720 sequence_coordinate: str,
1721 stamp_coordinate: str,
1722 overlay_cube: iris.cube.Cube | None = None,
1723 contour_cube: iris.cube.Cube | None = None,
1724 **kwargs,
1725):
1726 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1728 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1729 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1730 is present then postage stamp plots will be produced.
1732 If an overlay_cube and/or contour_cube are specified, multiple variables can
1733 be overplotted on the same figure.
1735 Parameters
1736 ----------
1737 method: "contourf" | "pcolormesh"
1738 The plotting method to use.
1739 cube: Cube
1740 Iris cube of the data to plot. It should have two spatial dimensions,
1741 such as lat and lon, and may also have a another two dimension to be
1742 plotted sequentially and/or as postage stamp plots.
1743 filename: str | None
1744 Name of the plot to write, used as a prefix for plot sequences. If None
1745 uses the recipe name.
1746 sequence_coordinate: str
1747 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1748 This coordinate must exist in the cube.
1749 stamp_coordinate: str
1750 Coordinate about which to plot postage stamp plots. Defaults to
1751 ``"realization"``.
1752 overlay_cube: Cube | None, optional
1753 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
1754 contour_cube: Cube | None, optional
1755 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
1757 Raises
1758 ------
1759 ValueError
1760 If the cube doesn't have the right dimensions.
1761 TypeError
1762 If the cube isn't a single cube.
1763 """
1764 recipe_title = get_recipe_metadata().get("title", "Untitled")
1766 # Ensure we've got a single cube.
1767 cube = check_single_cube(cube)
1769 # Check if there is a valid stamp coordinate in cube dimensions.
1770 if stamp_coordinate == "realization": 1770 ↛ 1775line 1770 didn't jump to line 1775 because the condition on line 1770 was always true
1771 stamp_coordinate = check_stamp_coordinate(cube)
1773 # Make postage stamp plots if stamp_coordinate exists and has more than a
1774 # single point.
1775 plotting_func = _plot_and_save_spatial_plot
1776 try:
1777 if cube.coord(stamp_coordinate).shape[0] > 1:
1778 plotting_func = _plot_and_save_postage_stamp_spatial_plot
1779 except iris.exceptions.CoordinateNotFoundError:
1780 pass
1782 # Produce a geographical scatter plot if the data have a
1783 # dimension called observation or model_obs_error
1784 if any( 1784 ↛ 1788line 1784 didn't jump to line 1788 because the condition on line 1784 was never true
1785 crd.var_name == "station" or crd.var_name == "model_obs_error"
1786 for crd in cube.coords()
1787 ):
1788 plotting_func = _plot_and_save_scattermap_plot
1790 # Must have a sequence coordinate.
1791 try:
1792 cube.coord(sequence_coordinate)
1793 except iris.exceptions.CoordinateNotFoundError as err:
1794 raise ValueError(f"Cube must have a {sequence_coordinate} coordinate.") from err
1796 # Create a plot for each value of the sequence coordinate.
1797 plot_index = []
1798 nplot = np.size(cube.coord(sequence_coordinate).points)
1800 for iseq, cube_slice in enumerate(cube.slices_over(sequence_coordinate)):
1801 # Set plot titles and filename
1802 seq_coord = cube_slice.coord(sequence_coordinate)
1803 plot_title, plot_filename = _set_title_and_filename(
1804 seq_coord, nplot, recipe_title, filename
1805 )
1807 # Extract sequence slice for overlay_cube and contour_cube if required.
1808 overlay_slice = slice_over_maybe(overlay_cube, sequence_coordinate, iseq)
1809 contour_slice = slice_over_maybe(contour_cube, sequence_coordinate, iseq)
1811 # Do the actual plotting.
1812 plotting_func(
1813 cube_slice,
1814 filename=plot_filename,
1815 stamp_coordinate=stamp_coordinate,
1816 title=plot_title,
1817 method=method,
1818 overlay_cube=overlay_slice,
1819 contour_cube=contour_slice,
1820 **kwargs,
1821 )
1822 plot_index.append(plot_filename)
1824 # Add list of plots to plot metadata.
1825 complete_plot_index = _append_to_plot_index(plot_index)
1827 # Make a page to display the plots.
1828 _make_plot_html_page(complete_plot_index)
1831####################
1832# Public functions #
1833####################
1836def spatial_contour_plot(
1837 cube: iris.cube.Cube,
1838 filename: str = None,
1839 sequence_coordinate: str = "time",
1840 stamp_coordinate: str = "realization",
1841 **kwargs,
1842) -> iris.cube.Cube:
1843 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1845 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1846 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1847 is present then postage stamp plots will be produced.
1849 Parameters
1850 ----------
1851 cube: Cube
1852 Iris cube of the data to plot. It should have two spatial dimensions,
1853 such as lat and lon, and may also have a another two dimension to be
1854 plotted sequentially and/or as postage stamp plots.
1855 filename: str, optional
1856 Name of the plot to write, used as a prefix for plot sequences. Defaults
1857 to the recipe name.
1858 sequence_coordinate: str, optional
1859 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1860 This coordinate must exist in the cube.
1861 stamp_coordinate: str, optional
1862 Coordinate about which to plot postage stamp plots. Defaults to
1863 ``"realization"``.
1865 Returns
1866 -------
1867 Cube
1868 The original cube (so further operations can be applied).
1870 Raises
1871 ------
1872 ValueError
1873 If the cube doesn't have the right dimensions.
1874 TypeError
1875 If the cube isn't a single cube.
1876 """
1877 _spatial_plot(
1878 "contourf", cube, filename, sequence_coordinate, stamp_coordinate, **kwargs
1879 )
1880 return cube
1883def spatial_pcolormesh_plot(
1884 cube: iris.cube.Cube,
1885 filename: str = None,
1886 sequence_coordinate: str = "time",
1887 stamp_coordinate: str = "realization",
1888 **kwargs,
1889) -> iris.cube.Cube:
1890 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1892 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1893 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1894 is present then postage stamp plots will be produced.
1896 This function is significantly faster than ``spatial_contour_plot``,
1897 especially at high resolutions, and should be preferred unless contiguous
1898 contour areas are important.
1900 Parameters
1901 ----------
1902 cube: Cube
1903 Iris cube of the data to plot. It should have two spatial dimensions,
1904 such as lat and lon, and may also have a another two dimension to be
1905 plotted sequentially and/or as postage stamp plots.
1906 filename: str, optional
1907 Name of the plot to write, used as a prefix for plot sequences. Defaults
1908 to the recipe name.
1909 sequence_coordinate: str, optional
1910 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1911 This coordinate must exist in the cube.
1912 stamp_coordinate: str, optional
1913 Coordinate about which to plot postage stamp plots. Defaults to
1914 ``"realization"``.
1916 Returns
1917 -------
1918 Cube
1919 The original cube (so further operations can be applied).
1921 Raises
1922 ------
1923 ValueError
1924 If the cube doesn't have the right dimensions.
1925 TypeError
1926 If the cube isn't a single cube.
1927 """
1928 _spatial_plot(
1929 "pcolormesh", cube, filename, sequence_coordinate, stamp_coordinate, **kwargs
1930 )
1931 return cube
1934def spatial_multi_pcolormesh_plot(
1935 cube: iris.cube.Cube,
1936 overlay_cube: iris.cube.Cube,
1937 contour_cube: iris.cube.Cube,
1938 filename: str = None,
1939 sequence_coordinate: str = "time",
1940 stamp_coordinate: str = "realization",
1941 **kwargs,
1942) -> iris.cube.Cube:
1943 """Plot a set of spatial variables onto a map from a 2D, 3D, or 4D cube.
1945 A 2D basis cube spatial field can be plotted, but if the sequence_coordinate is present
1946 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1947 is present then postage stamp plots will be produced.
1949 If specified, a masked overlay_cube can be overplotted on top of the base cube.
1951 If specified, contours of a contour_cube can be overplotted on top of those.
1953 For single-variable equivalent of this routine, use spatial_pcolormesh_plot.
1955 This function is significantly faster than ``spatial_contour_plot``,
1956 especially at high resolutions, and should be preferred unless contiguous
1957 contour areas are important.
1959 Parameters
1960 ----------
1961 cube: Cube
1962 Iris cube of the data to plot. It should have two spatial dimensions,
1963 such as lat and lon, and may also have a another two dimension to be
1964 plotted sequentially and/or as postage stamp plots.
1965 overlay_cube: Cube
1966 Iris cube of the data to plot as an overlay on top of basis cube. It should have two spatial dimensions,
1967 such as lat and lon, and may also have a another two dimension to be
1968 plotted sequentially and/or as postage stamp plots. This is likely to be a masked cube in order not to hide the underlying basis cube.
1969 contour_cube: Cube
1970 Iris cube of the data to plot as a contour overlay on top of basis cube and overlay_cube. It should have two spatial dimensions,
1971 such as lat and lon, and may also have a another two dimension to be
1972 plotted sequentially and/or as postage stamp plots.
1973 filename: str, optional
1974 Name of the plot to write, used as a prefix for plot sequences. Defaults
1975 to the recipe name.
1976 sequence_coordinate: str, optional
1977 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1978 This coordinate must exist in the cube.
1979 stamp_coordinate: str, optional
1980 Coordinate about which to plot postage stamp plots. Defaults to
1981 ``"realization"``.
1983 Returns
1984 -------
1985 Cube
1986 The original cube (so further operations can be applied).
1988 Raises
1989 ------
1990 ValueError
1991 If the cube doesn't have the right dimensions.
1992 TypeError
1993 If the cube isn't a single cube.
1994 """
1995 _spatial_plot(
1996 "pcolormesh",
1997 cube,
1998 filename,
1999 sequence_coordinate,
2000 stamp_coordinate,
2001 overlay_cube=overlay_cube,
2002 contour_cube=contour_cube,
2003 )
2004 return cube, overlay_cube, contour_cube
2007# TODO: Expand function to handle ensemble data.
2008# line_coordinate: str, optional
2009# Coordinate about which to plot multiple lines. Defaults to
2010# ``"realization"``.
2011def plot_line_series(
2012 cube: iris.cube.Cube | iris.cube.CubeList,
2013 filename: str = None,
2014 series_coordinate: str = "time",
2015 # line_coordinate: str = "realization",
2016 **kwargs,
2017) -> iris.cube.Cube | iris.cube.CubeList:
2018 """Plot a line plot for the specified coordinate.
2020 The Cube or CubeList must be 1D.
2022 Parameters
2023 ----------
2024 iris.cube | iris.cube.CubeList
2025 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
2026 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2027 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2028 filename: str, optional
2029 Name of the plot to write, used as a prefix for plot sequences. Defaults
2030 to the recipe name.
2031 series_coordinate: str, optional
2032 Coordinate about which to make a series. Defaults to ``"time"``. This
2033 coordinate must exist in the cube.
2035 Returns
2036 -------
2037 iris.cube.Cube | iris.cube.CubeList
2038 The original Cube or CubeList (so further operations can be applied).
2039 plotted data.
2041 Raises
2042 ------
2043 ValueError
2044 If the cubes don't have the right dimensions.
2045 TypeError
2046 If the cube isn't a Cube or CubeList.
2047 """
2048 # Ensure we have a name for the plot file.
2049 recipe_title = get_recipe_metadata().get("title", "Untitled")
2051 num_models = get_num_models(cube)
2053 validate_cube_shape(cube, num_models)
2055 # Iterate over all cubes and extract coordinate to plot.
2056 cubes = iter_maybe(cube)
2057 coords = []
2058 for cube in cubes:
2059 try:
2060 coords.append(cube.coord(series_coordinate))
2061 except iris.exceptions.CoordinateNotFoundError as err:
2062 raise ValueError(
2063 f"Cube must have a {series_coordinate} coordinate."
2064 ) from err
2065 if cube.ndim > 2 or not cube.coords("realization"):
2066 raise ValueError("Cube must be 1D or 2D with a realization coordinate.")
2068 # Format the title and filename using plotted series coordinate
2069 nplot = 1
2070 seq_coord = coords[0]
2071 plot_title, plot_filename = _set_title_and_filename(
2072 seq_coord, nplot, recipe_title, filename
2073 )
2075 # Do the actual plotting.
2076 _plot_and_save_line_series(cubes, coords, "realization", plot_filename, plot_title)
2078 # Add list of plots to plot metadata.
2079 plot_index = _append_to_plot_index([plot_filename])
2081 # Make a page to display the plots.
2082 _make_plot_html_page(plot_index)
2084 return cube
2087def plot_vertical_line_series(
2088 cubes: iris.cube.Cube | iris.cube.CubeList,
2089 filename: str = None,
2090 series_coordinate: str = "model_level_number",
2091 sequence_coordinate: str = "time",
2092 # line_coordinate: str = "realization",
2093 **kwargs,
2094) -> iris.cube.Cube | iris.cube.CubeList:
2095 """Plot a line plot against a type of vertical coordinate.
2097 The Cube or CubeList must be 1D.
2099 A 1D line plot with y-axis as pressure coordinate can be plotted, but if the sequence_coordinate is present
2100 then a sequence of plots will be produced.
2102 Parameters
2103 ----------
2104 iris.cube | iris.cube.CubeList
2105 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
2106 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2107 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2108 filename: str, optional
2109 Name of the plot to write, used as a prefix for plot sequences. Defaults
2110 to the recipe name.
2111 series_coordinate: str, optional
2112 Coordinate to plot on the y-axis. Can be ``pressure`` or
2113 ``model_level_number`` for UM, or ``full_levels`` or ``half_levels``
2114 for LFRic. Defaults to ``model_level_number``.
2115 This coordinate must exist in the cube.
2116 sequence_coordinate: str, optional
2117 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2118 This coordinate must exist in the cube.
2120 Returns
2121 -------
2122 iris.cube.Cube | iris.cube.CubeList
2123 The original Cube or CubeList (so further operations can be applied).
2124 Plotted data.
2126 Raises
2127 ------
2128 ValueError
2129 If the cubes doesn't have the right dimensions.
2130 TypeError
2131 If the cube isn't a Cube or CubeList.
2132 """
2133 # Ensure we have a name for the plot file.
2134 recipe_title = get_recipe_metadata().get("title", "Untitled")
2136 cubes = iter_maybe(cubes)
2137 # Initialise empty list to hold all data from all cubes in a CubeList
2138 all_data = []
2140 # Store min/max ranges for x range.
2141 x_levels = []
2143 num_models = get_num_models(cubes)
2145 validate_cube_shape(cubes, num_models)
2147 # Iterate over all cubes in cube or CubeList and plot.
2148 coords = []
2149 for cube in cubes:
2150 # Test if series coordinate i.e. pressure level exist for any cube with cube.ndim >=1.
2151 try:
2152 coords.append(cube.coord(series_coordinate))
2153 except iris.exceptions.CoordinateNotFoundError as err:
2154 raise ValueError(
2155 f"Cube must have a {series_coordinate} coordinate."
2156 ) from err
2158 try:
2159 if cube.ndim > 1 or not cube.coords("realization"): 2159 ↛ 2167line 2159 didn't jump to line 2167 because the condition on line 2159 was always true
2160 cube.coord(sequence_coordinate)
2161 except iris.exceptions.CoordinateNotFoundError as err:
2162 raise ValueError(
2163 f"Cube must have a {sequence_coordinate} coordinate or be 1D, or 2D with a realization coordinate."
2164 ) from err
2166 # Get minimum and maximum from levels information.
2167 _, levels, _ = colorbar_map_levels(cube, axis="x")
2168 if levels is not None: 2168 ↛ 2172line 2168 didn't jump to line 2172 because the condition on line 2168 was always true
2169 x_levels.append(min(levels))
2170 x_levels.append(max(levels))
2171 else:
2172 all_data.append(cube.data)
2174 if len(x_levels) == 0: 2174 ↛ 2176line 2174 didn't jump to line 2176 because the condition on line 2174 was never true
2175 # Combine all data into a single NumPy array
2176 combined_data = np.concatenate(all_data)
2178 # Set the lower and upper limit for the x-axis to ensure all plots have
2179 # same range. This needs to read the whole cube over the range of the
2180 # sequence and if applicable postage stamp coordinate.
2181 vmin = np.floor(combined_data.min())
2182 vmax = np.ceil(combined_data.max())
2183 else:
2184 vmin = min(x_levels)
2185 vmax = max(x_levels)
2187 # Matching the slices (matching by seq coord point; it may happen that
2188 # evaluated models do not cover the same seq coord range, hence matching
2189 # necessary)
2190 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2192 # Create a plot for each value of the sequence coordinate.
2193 # Allowing for multiple cubes in a CubeList to be plotted in the same plot for
2194 # similar sequence values. Passing a CubeList into the internal plotting function
2195 # for similar values of the sequence coordinate. cube_slice can be an iris.cube.Cube
2196 # or an iris.cube.CubeList.
2197 plot_index = []
2198 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2199 for cubes_slice in cube_iterables:
2200 # Format the coordinate value in a unit appropriate way.
2201 seq_coord = cubes_slice[0].coord(sequence_coordinate)
2202 plot_title, plot_filename = _set_title_and_filename(
2203 seq_coord, nplot, recipe_title, filename
2204 )
2206 # Do the actual plotting.
2207 _plot_and_save_vertical_line_series(
2208 cubes_slice,
2209 coords,
2210 "realization",
2211 plot_filename,
2212 series_coordinate,
2213 title=plot_title,
2214 vmin=vmin,
2215 vmax=vmax,
2216 )
2217 plot_index.append(plot_filename)
2219 # Add list of plots to plot metadata.
2220 complete_plot_index = _append_to_plot_index(plot_index)
2222 # Make a page to display the plots.
2223 _make_plot_html_page(complete_plot_index)
2225 return cubes
2228def qq_plot(
2229 cubes: iris.cube.CubeList,
2230 coordinates: list[str],
2231 percentiles: list[float],
2232 model_names: list[str],
2233 filename: str = None,
2234 one_to_one: bool = True,
2235 **kwargs,
2236) -> iris.cube.CubeList:
2237 """Plot a Quantile-Quantile plot between two models for common time points.
2239 The cubes will be normalised by collapsing each cube to its percentiles. Cubes are
2240 collapsed within the operator over all specified coordinates such as
2241 grid_latitude, grid_longitude, vertical levels, but also realisation representing
2242 ensemble members to ensure a 1D cube (array).
2244 Parameters
2245 ----------
2246 cubes: iris.cube.CubeList
2247 Two cubes of the same variable with different models.
2248 coordinate: list[str]
2249 The list of coordinates to collapse over. This list should be
2250 every coordinate within the cube to result in a 1D cube around
2251 the percentile coordinate.
2252 percent: list[float]
2253 A list of percentiles to appear in the plot.
2254 model_names: list[str]
2255 A list of model names to appear on the axis of the plot.
2256 filename: str, optional
2257 Filename of the plot to write.
2258 one_to_one: bool, optional
2259 If True a 1:1 line is plotted; if False it is not. Default is True.
2261 Raises
2262 ------
2263 ValueError
2264 When the cubes are not compatible.
2266 Notes
2267 -----
2268 The quantile-quantile plot is a variant on the scatter plot representing
2269 two datasets by their quantiles (percentiles) for common time points.
2270 This plot does not use a theoretical distribution to compare against, but
2271 compares percentiles of two datasets. This plot does
2272 not use all raw data points, but plots the selected percentiles (quantiles) of
2273 each variable instead for the two datasets, thereby normalising the data for a
2274 direct comparison between the selected percentiles of the two dataset distributions.
2276 Quantile-quantile plots are valuable for comparing against
2277 observations and other models. Identical percentiles between the variables
2278 will lie on the one-to-one line implying the values correspond well to each
2279 other. Where there is a deviation from the one-to-one line a range of
2280 possibilities exist depending on how and where the data is shifted (e.g.,
2281 Wilks 2011 [Wilks2011]_).
2283 For distributions above the one-to-one line the distribution is left-skewed;
2284 below is right-skewed. A distinct break implies a bimodal distribution, and
2285 closer values/values further apart at the tails imply poor representation of
2286 the extremes.
2288 References
2289 ----------
2290 .. [Wilks2011] Wilks, D.S., (2011) "Statistical Methods in the Atmospheric
2291 Sciences" Third Edition, vol. 100, Academic Press, Oxford, UK, 676 pp.
2292 """
2293 # Check cubes using same functionality as the difference operator.
2294 if len(cubes) != 2:
2295 raise ValueError("cubes should contain exactly 2 cubes.")
2296 base: Cube = cubes.extract_cube(iris.AttributeConstraint(cset_comparison_base=1))
2297 other: Cube = cubes.extract_cube(
2298 iris.Constraint(
2299 cube_func=lambda cube: "cset_comparison_base" not in cube.attributes
2300 )
2301 )
2303 # Get spatial coord names.
2304 base_lat_name, base_lon_name = get_cube_yxcoordname(base)
2305 other_lat_name, other_lon_name = get_cube_yxcoordname(other)
2307 # Ensure cubes to compare are on common differencing grid.
2308 # This is triggered if either
2309 # i) latitude and longitude shapes are not the same. Note grid points
2310 # are not compared directly as these can differ through rounding
2311 # errors.
2312 # ii) or variables are known to often sit on different grid staggering
2313 # in different models (e.g. cell center vs cell edge), as is the case
2314 # for UM and LFRic comparisons.
2315 # In future greater choice of regridding method might be applied depending
2316 # on variable type. Linear regridding can in general be appropriate for smooth
2317 # variables. Care should be taken with interpretation of differences
2318 # given this dependency on regridding.
2319 if (
2320 base.coord(base_lat_name).shape != other.coord(other_lat_name).shape
2321 or base.coord(base_lon_name).shape != other.coord(other_lon_name).shape
2322 ) or (
2323 base.long_name
2324 in [
2325 "eastward_wind_at_10m",
2326 "northward_wind_at_10m",
2327 "northward_wind_at_cell_centres",
2328 "eastward_wind_at_cell_centres",
2329 "zonal_wind_at_pressure_levels",
2330 "meridional_wind_at_pressure_levels",
2331 "potential_vorticity_at_pressure_levels",
2332 "vapour_specific_humidity_at_pressure_levels_for_climate_averaging",
2333 ]
2334 ):
2335 logging.debug(
2336 "Linear regridding base cube to other grid to compute differences"
2337 )
2338 base = regrid_onto_cube(base, other, method="Linear")
2340 # Extract just common time points.
2341 base, other = _extract_common_time_points(base, other)
2343 # Equalise attributes so we can merge.
2344 fully_equalise_attributes([base, other])
2345 logging.debug("Base: %s\nOther: %s", base, other)
2347 # Collapse cubes.
2348 base = collapse(
2349 base,
2350 coordinate=coordinates,
2351 method="PERCENTILE",
2352 additional_percent=percentiles,
2353 )
2354 other = collapse(
2355 other,
2356 coordinate=coordinates,
2357 method="PERCENTILE",
2358 additional_percent=percentiles,
2359 )
2361 # Ensure we have a name for the plot file.
2362 recipe_title = get_recipe_metadata().get("title", "Untitled")
2363 title = f"{recipe_title}"
2365 if filename is None:
2366 filename = slugify(recipe_title)
2368 # Add file extension.
2369 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2371 # Do the actual plotting on a scatter plot
2372 _plot_and_save_scatter_plot(
2373 base, other, plot_filename, title, one_to_one, model_names
2374 )
2376 # Add list of plots to plot metadata.
2377 plot_index = _append_to_plot_index([plot_filename])
2379 # Make a page to display the plots.
2380 _make_plot_html_page(plot_index)
2382 return iris.cube.CubeList([base, other])
2385def scatter_plot(
2386 cube_x: iris.cube.Cube | iris.cube.CubeList,
2387 cube_y: iris.cube.Cube | iris.cube.CubeList,
2388 filename: str = None,
2389 one_to_one: bool = True,
2390 **kwargs,
2391) -> iris.cube.CubeList:
2392 """Plot a scatter plot between two variables.
2394 Both cubes must be 1D.
2396 Parameters
2397 ----------
2398 cube_x: Cube | CubeList
2399 1 dimensional Cube of the data to plot on y-axis.
2400 cube_y: Cube | CubeList
2401 1 dimensional Cube of the data to plot on x-axis.
2402 filename: str, optional
2403 Filename of the plot to write.
2404 one_to_one: bool, optional
2405 If True a 1:1 line is plotted; if False it is not. Default is True.
2407 Returns
2408 -------
2409 cubes: CubeList
2410 CubeList of the original x and y cubes for further processing.
2412 Raises
2413 ------
2414 ValueError
2415 If the cube doesn't have the right dimensions and cubes not the same
2416 size.
2417 TypeError
2418 If the cube isn't a single cube.
2420 Notes
2421 -----
2422 Scatter plots are used for determining if there is a relationship between
2423 two variables. Positive relations have a slope going from bottom left to top
2424 right; Negative relations have a slope going from top left to bottom right.
2425 """
2426 # Iterate over all cubes in cube or CubeList and plot.
2427 for cube_iter in iter_maybe(cube_x):
2428 # Check cubes are correct shape.
2429 cube_iter = check_single_cube(cube_iter)
2430 if cube_iter.ndim > 1:
2431 raise ValueError("cube_x must be 1D.")
2433 # Iterate over all cubes in cube or CubeList and plot.
2434 for cube_iter in iter_maybe(cube_y):
2435 # Check cubes are correct shape.
2436 cube_iter = check_single_cube(cube_iter)
2437 if cube_iter.ndim > 1:
2438 raise ValueError("cube_y must be 1D.")
2440 # Ensure we have a name for the plot file.
2441 recipe_title = get_recipe_metadata().get("title", "Untitled")
2442 title = f"{recipe_title}"
2444 if filename is None:
2445 filename = slugify(recipe_title)
2447 # Add file extension.
2448 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2450 # Do the actual plotting.
2451 _plot_and_save_scatter_plot(cube_x, cube_y, plot_filename, title, one_to_one)
2453 # Add list of plots to plot metadata.
2454 plot_index = _append_to_plot_index([plot_filename])
2456 # Make a page to display the plots.
2457 _make_plot_html_page(plot_index)
2459 return iris.cube.CubeList([cube_x, cube_y])
2462def vector_plot(
2463 cube_u: iris.cube.Cube,
2464 cube_v: iris.cube.Cube,
2465 filename: str = None,
2466 sequence_coordinate: str = "time",
2467 **kwargs,
2468) -> iris.cube.CubeList:
2469 """Plot a vector plot based on the input u and v components."""
2470 recipe_title = get_recipe_metadata().get("title", "Untitled")
2472 # Cubes must have a matching sequence coordinate.
2473 try:
2474 # Check that the u and v cubes have the same sequence coordinate.
2475 if cube_u.coord(sequence_coordinate) != cube_v.coord(sequence_coordinate): 2475 ↛ anywhereline 2475 didn't jump anywhere: it always raised an exception.
2476 raise ValueError("Coordinates do not match.")
2477 except (iris.exceptions.CoordinateNotFoundError, ValueError) as err:
2478 raise ValueError(
2479 f"Cubes should have matching {sequence_coordinate} coordinate:\n{cube_u}\n{cube_v}"
2480 ) from err
2482 # Create a plot for each value of the sequence coordinate.
2483 plot_index = []
2484 nplot = np.size(cube_u[0].coord(sequence_coordinate).points)
2485 for cube_u_slice, cube_v_slice in zip(
2486 cube_u.slices_over(sequence_coordinate),
2487 cube_v.slices_over(sequence_coordinate),
2488 strict=True,
2489 ):
2490 # Format the coordinate value in a unit appropriate way.
2491 seq_coord = cube_u_slice.coord(sequence_coordinate)
2492 plot_title, plot_filename = _set_title_and_filename(
2493 seq_coord, nplot, recipe_title, filename
2494 )
2496 # Do the actual plotting.
2497 _plot_and_save_vector_plot(
2498 cube_u_slice,
2499 cube_v_slice,
2500 filename=plot_filename,
2501 title=plot_title,
2502 method="contourf",
2503 )
2504 plot_index.append(plot_filename)
2506 # Add list of plots to plot metadata.
2507 complete_plot_index = _append_to_plot_index(plot_index)
2509 # Make a page to display the plots.
2510 _make_plot_html_page(complete_plot_index)
2512 return iris.cube.CubeList([cube_u, cube_v])
2515def plot_histogram_series(
2516 cubes: iris.cube.Cube | iris.cube.CubeList,
2517 filename: str = None,
2518 sequence_coordinate: str = "time",
2519 stamp_coordinate: str = "realization",
2520 single_plot: bool = False,
2521 **kwargs,
2522) -> iris.cube.Cube | iris.cube.CubeList:
2523 """Plot a histogram plot for each vertical level provided.
2525 A histogram plot can be plotted, but if the sequence_coordinate (i.e. time)
2526 is present then a sequence of plots will be produced using the time slider
2527 functionality to scroll through histograms against time. If a
2528 stamp_coordinate is present then postage stamp plots will be produced. If
2529 stamp_coordinate and single_plot is True, all postage stamp plots will be
2530 plotted in a single plot instead of separate postage stamp plots.
2532 Parameters
2533 ----------
2534 cubes: Cube | iris.cube.CubeList
2535 Iris cube or CubeList of the data to plot. It should have a single dimension other
2536 than the stamp coordinate.
2537 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2538 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2539 filename: str, optional
2540 Name of the plot to write, used as a prefix for plot sequences. Defaults
2541 to the recipe name.
2542 sequence_coordinate: str, optional
2543 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2544 This coordinate must exist in the cube and will be used for the time
2545 slider.
2546 stamp_coordinate: str, optional
2547 Coordinate about which to plot postage stamp plots. Defaults to
2548 ``"realization"``.
2549 single_plot: bool, optional
2550 If True, all postage stamp plots will be plotted in a single plot. If
2551 False, each postage stamp plot will be plotted separately. Is only valid
2552 if stamp_coordinate exists and has more than a single point.
2554 Returns
2555 -------
2556 iris.cube.Cube | iris.cube.CubeList
2557 The original Cube or CubeList (so further operations can be applied).
2558 Plotted data.
2560 Raises
2561 ------
2562 ValueError
2563 If the cube doesn't have the right dimensions.
2564 TypeError
2565 If the cube isn't a Cube or CubeList.
2566 """
2567 recipe_title = get_recipe_metadata().get("title", "Untitled")
2569 cubes = iter_maybe(cubes)
2571 # Internal plotting function.
2572 plotting_func = _plot_and_save_histogram_series
2574 num_models = get_num_models(cubes)
2576 validate_cube_shape(cubes, num_models)
2578 # If several histograms are plotted, check sequence_coordinate
2579 check_sequence_coordinate(cubes, sequence_coordinate)
2581 # Get axis minimum and maximum from levels information.
2582 # If no levels set, derive minima and maxima from data in CubeList.
2583 vmin, vmax = _set_axis_range(cubes)
2585 # Make postage stamp plots if stamp_coordinate exists and has more than a
2586 # single point. If single_plot is True:
2587 # -- all postage stamp plots will be plotted in a single plot instead of
2588 # separate postage stamp plots.
2589 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
2590 # produced per single model only
2591 if num_models == 1: 2591 ↛ 2604line 2591 didn't jump to line 2604 because the condition on line 2591 was always true
2592 if ( 2592 ↛ 2596line 2592 didn't jump to line 2596 because the condition on line 2592 was never true
2593 stamp_coordinate in [c.name() for c in cubes[0].coords()]
2594 and cubes[0].coord(stamp_coordinate).shape[0] > 1
2595 ):
2596 if single_plot:
2597 plotting_func = (
2598 _plot_and_save_postage_stamps_in_single_plot_histogram_series
2599 )
2600 else:
2601 plotting_func = _plot_and_save_postage_stamp_histogram_series
2602 cube_iterables = cubes[0].slices_over(sequence_coordinate)
2603 else:
2604 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2606 plot_index = []
2607 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2608 # Create a plot for each value of the sequence coordinate. Allowing for
2609 # multiple cubes in a CubeList to be plotted in the same plot for similar
2610 # sequence values. Passing a CubeList into the internal plotting function
2611 # for similar values of the sequence coordinate. cube_slice can be an
2612 # iris.cube.Cube or an iris.cube.CubeList.
2613 for cube_slice in cube_iterables:
2614 single_cube = cube_slice
2615 if isinstance(cube_slice, iris.cube.CubeList): 2615 ↛ 2616line 2615 didn't jump to line 2616 because the condition on line 2615 was never true
2616 single_cube = cube_slice[0]
2618 # Ensure valid stamp coordinate in cube dimensions
2619 if stamp_coordinate == "realization": 2619 ↛ 2622line 2619 didn't jump to line 2622 because the condition on line 2619 was always true
2620 stamp_coordinate = check_stamp_coordinate(single_cube)
2621 # Set plot titles and filename, based on sequence coordinate
2622 seq_coord = single_cube.coord(sequence_coordinate)
2623 # Use time coordinate in title and filename if single histogram output.
2624 if sequence_coordinate == "realization" and nplot == 1: 2624 ↛ 2625line 2624 didn't jump to line 2625 because the condition on line 2624 was never true
2625 seq_coord = single_cube.coord("time")
2626 plot_title, plot_filename = _set_title_and_filename(
2627 seq_coord, nplot, recipe_title, filename
2628 )
2630 # Do the actual plotting.
2631 plotting_func(
2632 cube_slice,
2633 filename=plot_filename,
2634 stamp_coordinate=stamp_coordinate,
2635 title=plot_title,
2636 vmin=vmin,
2637 vmax=vmax,
2638 )
2639 plot_index.append(plot_filename)
2641 # Add list of plots to plot metadata.
2642 complete_plot_index = _append_to_plot_index(plot_index)
2644 # Make a page to display the plots.
2645 _make_plot_html_page(complete_plot_index)
2647 return cubes
2650def plot_power_spectrum_series(
2651 cubes: iris.cube.Cube | iris.cube.CubeList,
2652 filename: str = None,
2653 sequence_coordinate: str = "time",
2654 stamp_coordinate: str = "realization",
2655 single_plot: bool = False,
2656 **kwargs,
2657) -> iris.cube.Cube | iris.cube.CubeList:
2658 """Plot a power spectrum plot for each vertical level provided.
2660 A power spectrum plot can be plotted, but if the sequence_coordinate (i.e. time)
2661 is present then a sequence of plots will be produced using the time slider
2662 functionality to scroll through power spectrum against time. If a
2663 stamp_coordinate is present then postage stamp plots will be produced. If
2664 stamp_coordinate and single_plot is True, all postage stamp plots will be
2665 plotted in a single plot instead of separate postage stamp plots.
2667 Parameters
2668 ----------
2669 cubes: Cube | iris.cube.CubeList
2670 Iris cube or CubeList of the data to plot. It should have a single dimension other
2671 than the stamp coordinate.
2672 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2673 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2674 filename: str, optional
2675 Name of the plot to write, used as a prefix for plot sequences. Defaults
2676 to the recipe name.
2677 sequence_coordinate: str, optional
2678 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2679 This coordinate must exist in the cube and will be used for the time
2680 slider.
2681 stamp_coordinate: str, optional
2682 Coordinate about which to plot postage stamp plots. Defaults to
2683 ``"realization"``.
2684 single_plot: bool, optional
2685 If True, all postage stamp plots will be plotted in a single plot. If
2686 False, each postage stamp plot will be plotted separately. Is only valid
2687 if stamp_coordinate exists and has more than a single point.
2689 Returns
2690 -------
2691 iris.cube.Cube | iris.cube.CubeList
2692 The original Cube or CubeList (so further operations can be applied).
2693 Plotted data.
2695 Raises
2696 ------
2697 ValueError
2698 If the cube doesn't have the right dimensions.
2699 TypeError
2700 If the cube isn't a Cube or CubeList.
2701 """
2702 recipe_title = get_recipe_metadata().get("title", "Untitled")
2704 cubes = iter_maybe(cubes)
2706 # Internal plotting function.
2707 plotting_func = _plot_and_save_power_spectrum_series
2709 num_models = get_num_models(cubes)
2711 validate_cube_shape(cubes, num_models)
2713 # If several power spectra are plotted, check sequence_coordinate
2714 check_sequence_coordinate(cubes, sequence_coordinate)
2716 # Make postage stamp plots if stamp_coordinate exists and has more than a
2717 # single point. If single_plot is True:
2718 # -- all postage stamp plots will be plotted in a single plot instead of
2719 # separate postage stamp plots.
2720 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
2721 # produced per single model only
2722 if num_models == 1: 2722 ↛ 2735line 2722 didn't jump to line 2735 because the condition on line 2722 was always true
2723 if ( 2723 ↛ 2727line 2723 didn't jump to line 2727 because the condition on line 2723 was never true
2724 stamp_coordinate in [c.name() for c in cubes[0].coords()]
2725 and cubes[0].coord(stamp_coordinate).shape[0] > 1
2726 ):
2727 if single_plot:
2728 plotting_func = (
2729 _plot_and_save_postage_stamps_in_single_plot_power_spectrum_series
2730 )
2731 else:
2732 plotting_func = _plot_and_save_postage_stamp_power_spectrum_series
2733 cube_iterables = cubes[0].slices_over(sequence_coordinate)
2734 else:
2735 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2737 plot_index = []
2738 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2739 # Create a plot for each value of the sequence coordinate. Allowing for
2740 # multiple cubes in a CubeList to be plotted in the same plot for similar
2741 # sequence values. Passing a CubeList into the internal plotting function
2742 # for similar values of the sequence coordinate. cube_slice can be an
2743 # iris.cube.Cube or an iris.cube.CubeList.
2744 for cube_slice in cube_iterables:
2745 single_cube = cube_slice
2746 if isinstance(cube_slice, iris.cube.CubeList): 2746 ↛ 2747line 2746 didn't jump to line 2747 because the condition on line 2746 was never true
2747 single_cube = cube_slice[0]
2749 # Set stamp coordinate
2750 if stamp_coordinate == "realization": 2750 ↛ 2753line 2750 didn't jump to line 2753 because the condition on line 2750 was always true
2751 stamp_coordinate = check_stamp_coordinate(single_cube)
2752 # Set plot title and filenames based on sequence values
2753 seq_coord = single_cube.coord(sequence_coordinate)
2754 plot_title, plot_filename = _set_title_and_filename(
2755 seq_coord, nplot, recipe_title, filename
2756 )
2758 # Do the actual plotting.
2759 plotting_func(
2760 cube_slice,
2761 filename=plot_filename,
2762 stamp_coordinate=stamp_coordinate,
2763 title=plot_title,
2764 )
2765 plot_index.append(plot_filename)
2767 # Add list of plots to plot metadata.
2768 complete_plot_index = _append_to_plot_index(plot_index)
2770 # Make a page to display the plots.
2771 _make_plot_html_page(complete_plot_index)
2773 return cubes