Coverage for src/CSET/operators/plot.py: 80%
851 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-17 15:44 +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 ax.set_xlim(vmin, vmax)
1337 ax.tick_params(axis="both", labelsize=12)
1339 # Overlay grid-lines onto histogram plot.
1340 ax.grid(linestyle="--", color="grey", linewidth=1)
1341 if model_colors_map: 1341 ↛ 1342line 1341 didn't jump to line 1342 because the condition on line 1341 was never true
1342 ax.legend(loc="best", ncol=1, frameon=False, fontsize=16)
1344 # Save plot.
1345 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1346 logging.info("Saved histogram plot to %s", filename)
1347 plt.close(fig)
1350def _plot_and_save_postage_stamp_histogram_series(
1351 cube: iris.cube.Cube,
1352 filename: str,
1353 title: str,
1354 stamp_coordinate: str,
1355 vmin: float,
1356 vmax: float,
1357 **kwargs,
1358):
1359 """Plot and save postage (ensemble members) stamps for a histogram series.
1361 Parameters
1362 ----------
1363 cube: Cube
1364 2 dimensional Cube of the data to plot as histogram.
1365 filename: str
1366 Filename of the plot to write.
1367 title: str
1368 Plot title.
1369 stamp_coordinate: str
1370 Coordinate that becomes different plots.
1371 vmin: float
1372 minimum for pdf x-axis
1373 vmax: float
1374 maximum for pdf x-axis
1375 """
1376 # Use the smallest square grid that will fit the members.
1377 nmember = len(cube.coord(stamp_coordinate).points)
1378 grid_rows = int(math.sqrt(nmember))
1379 grid_size = math.ceil(nmember / grid_rows)
1381 fig = plt.figure(
1382 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
1383 )
1384 # Make a subplot for each member.
1385 for member, subplot in zip(
1386 cube.slices_over(stamp_coordinate),
1387 range(1, grid_size * grid_rows + 1),
1388 strict=False,
1389 ):
1390 # Implicit interface is much easier here, due to needing to have the
1391 # cartopy GeoAxes generated.
1392 plt.subplot(grid_rows, grid_size, subplot)
1393 # Reshape cube data into a single array to allow for a single histogram.
1394 # Otherwise we plot xdim histograms stacked.
1395 member_data_1d = (member.data).flatten()
1396 plt.hist(member_data_1d, density=True, stacked=True)
1397 axes = plt.gca()
1398 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
1399 axes.set_title(f"{mtitle}")
1400 axes.set_xlim(vmin, vmax)
1402 # Overall figure title.
1403 fig.suptitle(title, fontsize=16)
1405 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1406 logging.info("Saved histogram postage stamp plot to %s", filename)
1407 plt.close(fig)
1410def _plot_and_save_postage_stamps_in_single_plot_histogram_series(
1411 cube: iris.cube.Cube,
1412 filename: str,
1413 title: str,
1414 stamp_coordinate: str,
1415 vmin: float,
1416 vmax: float,
1417 **kwargs,
1418):
1419 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1420 ax.set_title(title, fontsize=16)
1421 ax.set_xlim(vmin, vmax)
1422 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1423 ax.set_ylabel("normalised probability density", fontsize=14)
1424 # Loop over all slices along the stamp_coordinate
1425 for member in cube.slices_over(stamp_coordinate):
1426 # Flatten the member data to 1D
1427 member_data_1d = member.data.flatten()
1428 # Plot the histogram using plt.hist
1429 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
1430 plt.hist(
1431 member_data_1d,
1432 density=True,
1433 stacked=True,
1434 label=f"{mtitle}",
1435 )
1437 # Add a legend
1438 ax.legend(fontsize=16)
1440 # Save the figure to a file
1441 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1442 logging.info("Saved histogram postage stamp plot to %s", filename)
1444 # Close the figure
1445 plt.close(fig)
1448def _plot_and_save_scattermap_plot(
1449 cube: iris.cube.Cube, filename: str, title: str, projection=None, **kwargs
1450):
1451 """Plot and save a geographical scatter plot.
1453 Parameters
1454 ----------
1455 cube: Cube
1456 1 dimensional Cube of the data points with auxiliary latitude and
1457 longitude coordinates,
1458 filename: str
1459 Filename of the plot to write.
1460 title: str
1461 Plot title.
1462 projection: str
1463 Mapping projection to be used by cartopy.
1464 """
1465 # Setup plot details, size, resolution, etc.
1466 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1467 if projection is not None:
1468 # Apart from the default, the only projection we currently support is
1469 # a stereographic projection over the North Pole.
1470 if projection == "NP_Stereo":
1471 axes = plt.axes(projection=ccrs.NorthPolarStereo(central_longitude=0.0))
1472 else:
1473 raise ValueError(f"Unknown projection: {projection}")
1474 else:
1475 axes = plt.axes(projection=ccrs.PlateCarree())
1477 # Scatter plot of the field. The marker size is chosen to give
1478 # symbols that decrease in size as the number of observations
1479 # increases, although the fraction of the figure covered by
1480 # symbols increases roughly as N^(1/2), disregarding overlaps,
1481 # and has been selected for the default figure size of (10, 10).
1482 # Should this be changed, the marker size should be adjusted in
1483 # proportion to the area of the figure.
1484 mrk_size = int(np.sqrt(2500000.0 / len(cube.data)))
1485 klon = None
1486 klat = None
1487 for kc in range(len(cube.aux_coords)):
1488 if cube.aux_coords[kc].standard_name == "latitude":
1489 klat = kc
1490 elif cube.aux_coords[kc].standard_name == "longitude":
1491 klon = kc
1492 scatter_map = iplt.scatter(
1493 cube.aux_coords[klon],
1494 cube.aux_coords[klat],
1495 c=cube.data[:],
1496 s=mrk_size,
1497 cmap="jet",
1498 edgecolors="k",
1499 )
1501 # Add coastlines and borderlines.
1502 try:
1503 axes.coastlines(resolution="10m")
1504 axes.add_feature(cfeature.BORDERS)
1505 except AttributeError:
1506 pass
1508 # Add title.
1509 axes.set_title(title, fontsize=16)
1511 # Add colour bar.
1512 cbar = fig.colorbar(scatter_map)
1513 cbar.set_label(label=f"{cube.name()} ({cube.units})", size=20)
1515 # Save plot.
1516 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1517 logging.info("Saved geographical scatter plot to %s", filename)
1518 plt.close(fig)
1521def _plot_and_save_power_spectrum_series(
1522 cubes: iris.cube.Cube | iris.cube.CubeList,
1523 filename: str,
1524 title: str,
1525 **kwargs,
1526):
1527 """Plot and save a power spectrum series.
1529 Parameters
1530 ----------
1531 cubes: Cube or CubeList
1532 2 dimensional Cube or CubeList of the data to plot as power spectrum.
1533 filename: str
1534 Filename of the plot to write.
1535 title: str
1536 Plot title.
1537 """
1538 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1539 ax = plt.gca()
1541 model_colors_map = get_model_colors_map(cubes)
1543 for cube in iter_maybe(cubes):
1544 # Calculate power spectrum
1546 # Extract time coordinate and convert to datetime
1547 time_coord = cube.coord("time")
1548 time_points = time_coord.units.num2date(time_coord.points)
1550 # Choose one time point (e.g., the first one)
1551 target_time = time_points[0]
1553 # Bind target_time inside the lambda using a default argument
1554 time_constraint = iris.Constraint(
1555 time=lambda cell, target_time=target_time: cell.point == target_time
1556 )
1558 cube = cube.extract(time_constraint)
1560 if cube.ndim == 2:
1561 cube_3d = cube.data[np.newaxis, :, :]
1562 logging.debug("Adding in new axis for a 2 dimensional cube.")
1563 elif cube.ndim == 3: 1563 ↛ 1564line 1563 didn't jump to line 1564 because the condition on line 1563 was never true
1564 cube_3d = cube.data
1565 else:
1566 raise ValueError("Cube dimensions unsuitable for power spectra code")
1567 raise ValueError(
1568 f"Cube is {cube.ndim} dimensional. Cube should be 2 or 3 dimensional."
1569 )
1571 # Calculate spectra
1572 ps_array = DCT_ps(cube_3d)
1574 ps_cube = iris.cube.Cube(
1575 ps_array,
1576 long_name="power_spectra",
1577 )
1579 ps_cube.attributes["model_name"] = cube.attributes.get("model_name")
1581 # Create a frequency/wavelength array for coordinate
1582 ps_len = ps_cube.data.shape[1]
1583 freqs = np.arange(1, ps_len + 1)
1584 freq_coord = iris.coords.DimCoord(freqs, long_name="frequency", units="1")
1586 # Convert datetime to numeric time using original units
1587 numeric_time = time_coord.units.date2num(time_points)
1588 # Create a new DimCoord with numeric time
1589 new_time_coord = iris.coords.DimCoord(
1590 numeric_time, standard_name="time", units=time_coord.units
1591 )
1593 # Add time and frequency coordinate to spectra cube.
1594 ps_cube.add_dim_coord(new_time_coord.copy(), 0)
1595 ps_cube.add_dim_coord(freq_coord.copy(), 1)
1597 # Extract data from the cube
1598 frequency = ps_cube.coord("frequency").points
1599 power_spectrum = ps_cube.data
1601 label = None
1602 color = "black"
1603 if model_colors_map: 1603 ↛ 1604line 1603 didn't jump to line 1604 because the condition on line 1603 was never true
1604 label = ps_cube.attributes.get("model_name")
1605 color = model_colors_map[label]
1606 ax.plot(frequency, power_spectrum[0], color=color, label=label)
1608 # Add some labels and tweak the style.
1609 ax.set_title(title, fontsize=16)
1610 ax.set_xlabel("Wavenumber", fontsize=14)
1611 ax.set_ylabel("Power", fontsize=14)
1612 ax.tick_params(axis="both", labelsize=12)
1614 # Set log-log scale
1615 ax.set_xscale("log")
1616 ax.set_yscale("log")
1618 # Overlay grid-lines onto power spectrum plot.
1619 ax.grid(linestyle="--", color="grey", linewidth=1)
1620 if model_colors_map: 1620 ↛ 1621line 1620 didn't jump to line 1621 because the condition on line 1620 was never true
1621 ax.legend(loc="best", ncol=1, frameon=False, fontsize=16)
1623 # Save plot.
1624 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1625 logging.info("Saved power spectrum plot to %s", filename)
1626 plt.close(fig)
1629def _plot_and_save_postage_stamp_power_spectrum_series(
1630 cube: iris.cube.Cube,
1631 filename: str,
1632 title: str,
1633 stamp_coordinate: str,
1634 **kwargs,
1635):
1636 """Plot and save postage (ensemble members) stamps for a power spectrum series.
1638 Parameters
1639 ----------
1640 cube: Cube
1641 2 dimensional Cube of the data to plot as power spectrum.
1642 filename: str
1643 Filename of the plot to write.
1644 title: str
1645 Plot title.
1646 stamp_coordinate: str
1647 Coordinate that becomes different plots.
1648 """
1649 # Use the smallest square grid that will fit the members.
1650 nmember = len(cube.coord(stamp_coordinate).points)
1651 grid_rows = int(math.sqrt(nmember))
1652 grid_size = math.ceil(nmember / grid_rows)
1654 fig = plt.figure(
1655 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
1656 )
1658 # Make a subplot for each member.
1659 for member, subplot in zip(
1660 cube.slices_over(stamp_coordinate),
1661 range(1, grid_size * grid_rows + 1),
1662 strict=False,
1663 ):
1664 # Implicit interface is much easier here, due to needing to have the
1665 # cartopy GeoAxes generated.
1666 plt.subplot(grid_rows, grid_size, subplot)
1668 frequency = member.coord("frequency").points
1670 axes = plt.gca()
1671 axes.plot(frequency, member.data)
1672 axes.set_title(f"Member #{member.coord(stamp_coordinate).points[0]}")
1674 # Overall figure title.
1675 fig.suptitle(title, fontsize=16)
1677 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1678 logging.info("Saved power spectra postage stamp plot to %s", filename)
1679 plt.close(fig)
1682def _plot_and_save_postage_stamps_in_single_plot_power_spectrum_series(
1683 cube: iris.cube.Cube,
1684 filename: str,
1685 title: str,
1686 stamp_coordinate: str,
1687 **kwargs,
1688):
1689 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1690 ax.set_title(title, fontsize=16)
1691 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1692 ax.set_ylabel("Power", fontsize=14)
1693 # Loop over all slices along the stamp_coordinate
1694 for member in cube.slices_over(stamp_coordinate):
1695 frequency = member.coord("frequency").points
1696 ax.plot(
1697 frequency,
1698 member.data,
1699 label=f"Member #{member.coord(stamp_coordinate).points[0]}",
1700 )
1702 # Add a legend
1703 ax.legend(fontsize=16)
1705 # Save the figure to a file
1706 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1707 logging.info("Saved power spectra plot to %s", filename)
1709 # Close the figure
1710 plt.close(fig)
1713def _spatial_plot(
1714 method: Literal["contourf", "pcolormesh"],
1715 cube: iris.cube.Cube,
1716 filename: str | None,
1717 sequence_coordinate: str,
1718 stamp_coordinate: str,
1719 overlay_cube: iris.cube.Cube | None = None,
1720 contour_cube: iris.cube.Cube | None = None,
1721 **kwargs,
1722):
1723 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1725 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1726 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1727 is present then postage stamp plots will be produced.
1729 If an overlay_cube and/or contour_cube are specified, multiple variables can
1730 be overplotted on the same figure.
1732 Parameters
1733 ----------
1734 method: "contourf" | "pcolormesh"
1735 The plotting method to use.
1736 cube: Cube
1737 Iris cube of the data to plot. It should have two spatial dimensions,
1738 such as lat and lon, and may also have a another two dimension to be
1739 plotted sequentially and/or as postage stamp plots.
1740 filename: str | None
1741 Name of the plot to write, used as a prefix for plot sequences. If None
1742 uses the recipe name.
1743 sequence_coordinate: str
1744 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1745 This coordinate must exist in the cube.
1746 stamp_coordinate: str
1747 Coordinate about which to plot postage stamp plots. Defaults to
1748 ``"realization"``.
1749 overlay_cube: Cube | None, optional
1750 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
1751 contour_cube: Cube | None, optional
1752 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
1754 Raises
1755 ------
1756 ValueError
1757 If the cube doesn't have the right dimensions.
1758 TypeError
1759 If the cube isn't a single cube.
1760 """
1761 recipe_title = get_recipe_metadata().get("title", "Untitled")
1763 # Ensure we've got a single cube.
1764 cube = check_single_cube(cube)
1766 # Check if there is a valid stamp coordinate in cube dimensions.
1767 if stamp_coordinate == "realization": 1767 ↛ 1772line 1767 didn't jump to line 1772 because the condition on line 1767 was always true
1768 stamp_coordinate = check_stamp_coordinate(cube)
1770 # Make postage stamp plots if stamp_coordinate exists and has more than a
1771 # single point.
1772 plotting_func = _plot_and_save_spatial_plot
1773 try:
1774 if cube.coord(stamp_coordinate).shape[0] > 1:
1775 plotting_func = _plot_and_save_postage_stamp_spatial_plot
1776 except iris.exceptions.CoordinateNotFoundError:
1777 pass
1779 # Produce a geographical scatter plot if the data have a
1780 # dimension called observation or model_obs_error
1781 if any( 1781 ↛ 1785line 1781 didn't jump to line 1785 because the condition on line 1781 was never true
1782 crd.var_name == "station" or crd.var_name == "model_obs_error"
1783 for crd in cube.coords()
1784 ):
1785 plotting_func = _plot_and_save_scattermap_plot
1787 # Must have a sequence coordinate.
1788 try:
1789 cube.coord(sequence_coordinate)
1790 except iris.exceptions.CoordinateNotFoundError as err:
1791 raise ValueError(f"Cube must have a {sequence_coordinate} coordinate.") from err
1793 # Create a plot for each value of the sequence coordinate.
1794 plot_index = []
1795 nplot = np.size(cube.coord(sequence_coordinate).points)
1797 for iseq, cube_slice in enumerate(cube.slices_over(sequence_coordinate)):
1798 # Set plot titles and filename
1799 seq_coord = cube_slice.coord(sequence_coordinate)
1800 plot_title, plot_filename = _set_title_and_filename(
1801 seq_coord, nplot, recipe_title, filename
1802 )
1804 # Extract sequence slice for overlay_cube and contour_cube if required.
1805 overlay_slice = slice_over_maybe(overlay_cube, sequence_coordinate, iseq)
1806 contour_slice = slice_over_maybe(contour_cube, sequence_coordinate, iseq)
1808 # Do the actual plotting.
1809 plotting_func(
1810 cube_slice,
1811 filename=plot_filename,
1812 stamp_coordinate=stamp_coordinate,
1813 title=plot_title,
1814 method=method,
1815 overlay_cube=overlay_slice,
1816 contour_cube=contour_slice,
1817 **kwargs,
1818 )
1819 plot_index.append(plot_filename)
1821 # Add list of plots to plot metadata.
1822 complete_plot_index = _append_to_plot_index(plot_index)
1824 # Make a page to display the plots.
1825 _make_plot_html_page(complete_plot_index)
1828####################
1829# Public functions #
1830####################
1833def spatial_contour_plot(
1834 cube: iris.cube.Cube,
1835 filename: str = None,
1836 sequence_coordinate: str = "time",
1837 stamp_coordinate: str = "realization",
1838 **kwargs,
1839) -> iris.cube.Cube:
1840 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1842 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1843 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1844 is present then postage stamp plots will be produced.
1846 Parameters
1847 ----------
1848 cube: Cube
1849 Iris cube of the data to plot. It should have two spatial dimensions,
1850 such as lat and lon, and may also have a another two dimension to be
1851 plotted sequentially and/or as postage stamp plots.
1852 filename: str, optional
1853 Name of the plot to write, used as a prefix for plot sequences. Defaults
1854 to the recipe name.
1855 sequence_coordinate: str, optional
1856 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1857 This coordinate must exist in the cube.
1858 stamp_coordinate: str, optional
1859 Coordinate about which to plot postage stamp plots. Defaults to
1860 ``"realization"``.
1862 Returns
1863 -------
1864 Cube
1865 The original cube (so further operations can be applied).
1867 Raises
1868 ------
1869 ValueError
1870 If the cube doesn't have the right dimensions.
1871 TypeError
1872 If the cube isn't a single cube.
1873 """
1874 _spatial_plot(
1875 "contourf", cube, filename, sequence_coordinate, stamp_coordinate, **kwargs
1876 )
1877 return cube
1880def spatial_pcolormesh_plot(
1881 cube: iris.cube.Cube,
1882 filename: str = None,
1883 sequence_coordinate: str = "time",
1884 stamp_coordinate: str = "realization",
1885 **kwargs,
1886) -> iris.cube.Cube:
1887 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1889 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1890 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1891 is present then postage stamp plots will be produced.
1893 This function is significantly faster than ``spatial_contour_plot``,
1894 especially at high resolutions, and should be preferred unless contiguous
1895 contour areas are important.
1897 Parameters
1898 ----------
1899 cube: Cube
1900 Iris cube of the data to plot. It should have two spatial dimensions,
1901 such as lat and lon, and may also have a another two dimension to be
1902 plotted sequentially and/or as postage stamp plots.
1903 filename: str, optional
1904 Name of the plot to write, used as a prefix for plot sequences. Defaults
1905 to the recipe name.
1906 sequence_coordinate: str, optional
1907 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1908 This coordinate must exist in the cube.
1909 stamp_coordinate: str, optional
1910 Coordinate about which to plot postage stamp plots. Defaults to
1911 ``"realization"``.
1913 Returns
1914 -------
1915 Cube
1916 The original cube (so further operations can be applied).
1918 Raises
1919 ------
1920 ValueError
1921 If the cube doesn't have the right dimensions.
1922 TypeError
1923 If the cube isn't a single cube.
1924 """
1925 _spatial_plot(
1926 "pcolormesh", cube, filename, sequence_coordinate, stamp_coordinate, **kwargs
1927 )
1928 return cube
1931def spatial_multi_pcolormesh_plot(
1932 cube: iris.cube.Cube,
1933 overlay_cube: iris.cube.Cube,
1934 contour_cube: iris.cube.Cube,
1935 filename: str = None,
1936 sequence_coordinate: str = "time",
1937 stamp_coordinate: str = "realization",
1938 **kwargs,
1939) -> iris.cube.Cube:
1940 """Plot a set of spatial variables onto a map from a 2D, 3D, or 4D cube.
1942 A 2D basis cube spatial field can be plotted, but if the sequence_coordinate is present
1943 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1944 is present then postage stamp plots will be produced.
1946 If specified, a masked overlay_cube can be overplotted on top of the base cube.
1948 If specified, contours of a contour_cube can be overplotted on top of those.
1950 For single-variable equivalent of this routine, use spatial_pcolormesh_plot.
1952 This function is significantly faster than ``spatial_contour_plot``,
1953 especially at high resolutions, and should be preferred unless contiguous
1954 contour areas are important.
1956 Parameters
1957 ----------
1958 cube: Cube
1959 Iris cube of the data to plot. It should have two spatial dimensions,
1960 such as lat and lon, and may also have a another two dimension to be
1961 plotted sequentially and/or as postage stamp plots.
1962 overlay_cube: Cube
1963 Iris cube of the data to plot as an overlay on top of basis cube. It should have two spatial dimensions,
1964 such as lat and lon, and may also have a another two dimension to be
1965 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.
1966 contour_cube: Cube
1967 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,
1968 such as lat and lon, and may also have a another two dimension to be
1969 plotted sequentially and/or as postage stamp plots.
1970 filename: str, optional
1971 Name of the plot to write, used as a prefix for plot sequences. Defaults
1972 to the recipe name.
1973 sequence_coordinate: str, optional
1974 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1975 This coordinate must exist in the cube.
1976 stamp_coordinate: str, optional
1977 Coordinate about which to plot postage stamp plots. Defaults to
1978 ``"realization"``.
1980 Returns
1981 -------
1982 Cube
1983 The original cube (so further operations can be applied).
1985 Raises
1986 ------
1987 ValueError
1988 If the cube doesn't have the right dimensions.
1989 TypeError
1990 If the cube isn't a single cube.
1991 """
1992 _spatial_plot(
1993 "pcolormesh",
1994 cube,
1995 filename,
1996 sequence_coordinate,
1997 stamp_coordinate,
1998 overlay_cube=overlay_cube,
1999 contour_cube=contour_cube,
2000 )
2001 return cube, overlay_cube, contour_cube
2004# TODO: Expand function to handle ensemble data.
2005# line_coordinate: str, optional
2006# Coordinate about which to plot multiple lines. Defaults to
2007# ``"realization"``.
2008def plot_line_series(
2009 cube: iris.cube.Cube | iris.cube.CubeList,
2010 filename: str = None,
2011 series_coordinate: str = "time",
2012 # line_coordinate: str = "realization",
2013 **kwargs,
2014) -> iris.cube.Cube | iris.cube.CubeList:
2015 """Plot a line plot for the specified coordinate.
2017 The Cube or CubeList must be 1D.
2019 Parameters
2020 ----------
2021 iris.cube | iris.cube.CubeList
2022 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
2023 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2024 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2025 filename: str, optional
2026 Name of the plot to write, used as a prefix for plot sequences. Defaults
2027 to the recipe name.
2028 series_coordinate: str, optional
2029 Coordinate about which to make a series. Defaults to ``"time"``. This
2030 coordinate must exist in the cube.
2032 Returns
2033 -------
2034 iris.cube.Cube | iris.cube.CubeList
2035 The original Cube or CubeList (so further operations can be applied).
2036 plotted data.
2038 Raises
2039 ------
2040 ValueError
2041 If the cubes don't have the right dimensions.
2042 TypeError
2043 If the cube isn't a Cube or CubeList.
2044 """
2045 # Ensure we have a name for the plot file.
2046 recipe_title = get_recipe_metadata().get("title", "Untitled")
2048 num_models = get_num_models(cube)
2050 validate_cube_shape(cube, num_models)
2052 # Iterate over all cubes and extract coordinate to plot.
2053 cubes = iter_maybe(cube)
2054 coords = []
2055 for cube in cubes:
2056 try:
2057 coords.append(cube.coord(series_coordinate))
2058 except iris.exceptions.CoordinateNotFoundError as err:
2059 raise ValueError(
2060 f"Cube must have a {series_coordinate} coordinate."
2061 ) from err
2062 if cube.ndim > 2 or not cube.coords("realization"):
2063 raise ValueError("Cube must be 1D or 2D with a realization coordinate.")
2065 # Format the title and filename using plotted series coordinate
2066 nplot = 1
2067 seq_coord = coords[0]
2068 plot_title, plot_filename = _set_title_and_filename(
2069 seq_coord, nplot, recipe_title, filename
2070 )
2072 # Do the actual plotting.
2073 _plot_and_save_line_series(cubes, coords, "realization", plot_filename, plot_title)
2075 # Add list of plots to plot metadata.
2076 plot_index = _append_to_plot_index([plot_filename])
2078 # Make a page to display the plots.
2079 _make_plot_html_page(plot_index)
2081 return cube
2084def plot_vertical_line_series(
2085 cubes: iris.cube.Cube | iris.cube.CubeList,
2086 filename: str = None,
2087 series_coordinate: str = "model_level_number",
2088 sequence_coordinate: str = "time",
2089 # line_coordinate: str = "realization",
2090 **kwargs,
2091) -> iris.cube.Cube | iris.cube.CubeList:
2092 """Plot a line plot against a type of vertical coordinate.
2094 The Cube or CubeList must be 1D.
2096 A 1D line plot with y-axis as pressure coordinate can be plotted, but if the sequence_coordinate is present
2097 then a sequence of plots will be produced.
2099 Parameters
2100 ----------
2101 iris.cube | iris.cube.CubeList
2102 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
2103 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2104 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2105 filename: str, optional
2106 Name of the plot to write, used as a prefix for plot sequences. Defaults
2107 to the recipe name.
2108 series_coordinate: str, optional
2109 Coordinate to plot on the y-axis. Can be ``pressure`` or
2110 ``model_level_number`` for UM, or ``full_levels`` or ``half_levels``
2111 for LFRic. Defaults to ``model_level_number``.
2112 This coordinate must exist in the cube.
2113 sequence_coordinate: str, optional
2114 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2115 This coordinate must exist in the cube.
2117 Returns
2118 -------
2119 iris.cube.Cube | iris.cube.CubeList
2120 The original Cube or CubeList (so further operations can be applied).
2121 Plotted data.
2123 Raises
2124 ------
2125 ValueError
2126 If the cubes doesn't have the right dimensions.
2127 TypeError
2128 If the cube isn't a Cube or CubeList.
2129 """
2130 # Ensure we have a name for the plot file.
2131 recipe_title = get_recipe_metadata().get("title", "Untitled")
2133 cubes = iter_maybe(cubes)
2134 # Initialise empty list to hold all data from all cubes in a CubeList
2135 all_data = []
2137 # Store min/max ranges for x range.
2138 x_levels = []
2140 num_models = get_num_models(cubes)
2142 validate_cube_shape(cubes, num_models)
2144 # Iterate over all cubes in cube or CubeList and plot.
2145 coords = []
2146 for cube in cubes:
2147 # Test if series coordinate i.e. pressure level exist for any cube with cube.ndim >=1.
2148 try:
2149 coords.append(cube.coord(series_coordinate))
2150 except iris.exceptions.CoordinateNotFoundError as err:
2151 raise ValueError(
2152 f"Cube must have a {series_coordinate} coordinate."
2153 ) from err
2155 try:
2156 if cube.ndim > 1 or not cube.coords("realization"): 2156 ↛ 2164line 2156 didn't jump to line 2164 because the condition on line 2156 was always true
2157 cube.coord(sequence_coordinate)
2158 except iris.exceptions.CoordinateNotFoundError as err:
2159 raise ValueError(
2160 f"Cube must have a {sequence_coordinate} coordinate or be 1D, or 2D with a realization coordinate."
2161 ) from err
2163 # Get minimum and maximum from levels information.
2164 _, levels, _ = colorbar_map_levels(cube, axis="x")
2165 if levels is not None: 2165 ↛ 2169line 2165 didn't jump to line 2169 because the condition on line 2165 was always true
2166 x_levels.append(min(levels))
2167 x_levels.append(max(levels))
2168 else:
2169 all_data.append(cube.data)
2171 if len(x_levels) == 0: 2171 ↛ 2173line 2171 didn't jump to line 2173 because the condition on line 2171 was never true
2172 # Combine all data into a single NumPy array
2173 combined_data = np.concatenate(all_data)
2175 # Set the lower and upper limit for the x-axis to ensure all plots have
2176 # same range. This needs to read the whole cube over the range of the
2177 # sequence and if applicable postage stamp coordinate.
2178 vmin = np.floor(combined_data.min())
2179 vmax = np.ceil(combined_data.max())
2180 else:
2181 vmin = min(x_levels)
2182 vmax = max(x_levels)
2184 # Matching the slices (matching by seq coord point; it may happen that
2185 # evaluated models do not cover the same seq coord range, hence matching
2186 # necessary)
2187 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2189 # Create a plot for each value of the sequence coordinate.
2190 # Allowing for multiple cubes in a CubeList to be plotted in the same plot for
2191 # similar sequence values. Passing a CubeList into the internal plotting function
2192 # for similar values of the sequence coordinate. cube_slice can be an iris.cube.Cube
2193 # or an iris.cube.CubeList.
2194 plot_index = []
2195 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2196 for cubes_slice in cube_iterables:
2197 # Format the coordinate value in a unit appropriate way.
2198 seq_coord = cubes_slice[0].coord(sequence_coordinate)
2199 plot_title, plot_filename = _set_title_and_filename(
2200 seq_coord, nplot, recipe_title, filename
2201 )
2203 # Do the actual plotting.
2204 _plot_and_save_vertical_line_series(
2205 cubes_slice,
2206 coords,
2207 "realization",
2208 plot_filename,
2209 series_coordinate,
2210 title=plot_title,
2211 vmin=vmin,
2212 vmax=vmax,
2213 )
2214 plot_index.append(plot_filename)
2216 # Add list of plots to plot metadata.
2217 complete_plot_index = _append_to_plot_index(plot_index)
2219 # Make a page to display the plots.
2220 _make_plot_html_page(complete_plot_index)
2222 return cubes
2225def qq_plot(
2226 cubes: iris.cube.CubeList,
2227 coordinates: list[str],
2228 percentiles: list[float],
2229 model_names: list[str],
2230 filename: str = None,
2231 one_to_one: bool = True,
2232 **kwargs,
2233) -> iris.cube.CubeList:
2234 """Plot a Quantile-Quantile plot between two models for common time points.
2236 The cubes will be normalised by collapsing each cube to its percentiles. Cubes are
2237 collapsed within the operator over all specified coordinates such as
2238 grid_latitude, grid_longitude, vertical levels, but also realisation representing
2239 ensemble members to ensure a 1D cube (array).
2241 Parameters
2242 ----------
2243 cubes: iris.cube.CubeList
2244 Two cubes of the same variable with different models.
2245 coordinate: list[str]
2246 The list of coordinates to collapse over. This list should be
2247 every coordinate within the cube to result in a 1D cube around
2248 the percentile coordinate.
2249 percent: list[float]
2250 A list of percentiles to appear in the plot.
2251 model_names: list[str]
2252 A list of model names to appear on the axis of the plot.
2253 filename: str, optional
2254 Filename of the plot to write.
2255 one_to_one: bool, optional
2256 If True a 1:1 line is plotted; if False it is not. Default is True.
2258 Raises
2259 ------
2260 ValueError
2261 When the cubes are not compatible.
2263 Notes
2264 -----
2265 The quantile-quantile plot is a variant on the scatter plot representing
2266 two datasets by their quantiles (percentiles) for common time points.
2267 This plot does not use a theoretical distribution to compare against, but
2268 compares percentiles of two datasets. This plot does
2269 not use all raw data points, but plots the selected percentiles (quantiles) of
2270 each variable instead for the two datasets, thereby normalising the data for a
2271 direct comparison between the selected percentiles of the two dataset distributions.
2273 Quantile-quantile plots are valuable for comparing against
2274 observations and other models. Identical percentiles between the variables
2275 will lie on the one-to-one line implying the values correspond well to each
2276 other. Where there is a deviation from the one-to-one line a range of
2277 possibilities exist depending on how and where the data is shifted (e.g.,
2278 Wilks 2011 [Wilks2011]_).
2280 For distributions above the one-to-one line the distribution is left-skewed;
2281 below is right-skewed. A distinct break implies a bimodal distribution, and
2282 closer values/values further apart at the tails imply poor representation of
2283 the extremes.
2285 References
2286 ----------
2287 .. [Wilks2011] Wilks, D.S., (2011) "Statistical Methods in the Atmospheric
2288 Sciences" Third Edition, vol. 100, Academic Press, Oxford, UK, 676 pp.
2289 """
2290 # Check cubes using same functionality as the difference operator.
2291 if len(cubes) != 2:
2292 raise ValueError("cubes should contain exactly 2 cubes.")
2293 base: Cube = cubes.extract_cube(iris.AttributeConstraint(cset_comparison_base=1))
2294 other: Cube = cubes.extract_cube(
2295 iris.Constraint(
2296 cube_func=lambda cube: "cset_comparison_base" not in cube.attributes
2297 )
2298 )
2300 # Get spatial coord names.
2301 base_lat_name, base_lon_name = get_cube_yxcoordname(base)
2302 other_lat_name, other_lon_name = get_cube_yxcoordname(other)
2304 # Ensure cubes to compare are on common differencing grid.
2305 # This is triggered if either
2306 # i) latitude and longitude shapes are not the same. Note grid points
2307 # are not compared directly as these can differ through rounding
2308 # errors.
2309 # ii) or variables are known to often sit on different grid staggering
2310 # in different models (e.g. cell center vs cell edge), as is the case
2311 # for UM and LFRic comparisons.
2312 # In future greater choice of regridding method might be applied depending
2313 # on variable type. Linear regridding can in general be appropriate for smooth
2314 # variables. Care should be taken with interpretation of differences
2315 # given this dependency on regridding.
2316 if (
2317 base.coord(base_lat_name).shape != other.coord(other_lat_name).shape
2318 or base.coord(base_lon_name).shape != other.coord(other_lon_name).shape
2319 ) or (
2320 base.long_name
2321 in [
2322 "eastward_wind_at_10m",
2323 "northward_wind_at_10m",
2324 "northward_wind_at_cell_centres",
2325 "eastward_wind_at_cell_centres",
2326 "zonal_wind_at_pressure_levels",
2327 "meridional_wind_at_pressure_levels",
2328 "potential_vorticity_at_pressure_levels",
2329 "vapour_specific_humidity_at_pressure_levels_for_climate_averaging",
2330 ]
2331 ):
2332 logging.debug(
2333 "Linear regridding base cube to other grid to compute differences"
2334 )
2335 base = regrid_onto_cube(base, other, method="Linear")
2337 # Extract just common time points.
2338 base, other = _extract_common_time_points(base, other)
2340 # Equalise attributes so we can merge.
2341 fully_equalise_attributes([base, other])
2342 logging.debug("Base: %s\nOther: %s", base, other)
2344 # Collapse cubes.
2345 base = collapse(
2346 base,
2347 coordinate=coordinates,
2348 method="PERCENTILE",
2349 additional_percent=percentiles,
2350 )
2351 other = collapse(
2352 other,
2353 coordinate=coordinates,
2354 method="PERCENTILE",
2355 additional_percent=percentiles,
2356 )
2358 # Ensure we have a name for the plot file.
2359 recipe_title = get_recipe_metadata().get("title", "Untitled")
2360 title = f"{recipe_title}"
2362 if filename is None:
2363 filename = slugify(recipe_title)
2365 # Add file extension.
2366 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2368 # Do the actual plotting on a scatter plot
2369 _plot_and_save_scatter_plot(
2370 base, other, plot_filename, title, one_to_one, model_names
2371 )
2373 # Add list of plots to plot metadata.
2374 plot_index = _append_to_plot_index([plot_filename])
2376 # Make a page to display the plots.
2377 _make_plot_html_page(plot_index)
2379 return iris.cube.CubeList([base, other])
2382def scatter_plot(
2383 cube_x: iris.cube.Cube | iris.cube.CubeList,
2384 cube_y: iris.cube.Cube | iris.cube.CubeList,
2385 filename: str = None,
2386 one_to_one: bool = True,
2387 **kwargs,
2388) -> iris.cube.CubeList:
2389 """Plot a scatter plot between two variables.
2391 Both cubes must be 1D.
2393 Parameters
2394 ----------
2395 cube_x: Cube | CubeList
2396 1 dimensional Cube of the data to plot on y-axis.
2397 cube_y: Cube | CubeList
2398 1 dimensional Cube of the data to plot on x-axis.
2399 filename: str, optional
2400 Filename of the plot to write.
2401 one_to_one: bool, optional
2402 If True a 1:1 line is plotted; if False it is not. Default is True.
2404 Returns
2405 -------
2406 cubes: CubeList
2407 CubeList of the original x and y cubes for further processing.
2409 Raises
2410 ------
2411 ValueError
2412 If the cube doesn't have the right dimensions and cubes not the same
2413 size.
2414 TypeError
2415 If the cube isn't a single cube.
2417 Notes
2418 -----
2419 Scatter plots are used for determining if there is a relationship between
2420 two variables. Positive relations have a slope going from bottom left to top
2421 right; Negative relations have a slope going from top left to bottom right.
2422 """
2423 # Iterate over all cubes in cube or CubeList and plot.
2424 for cube_iter in iter_maybe(cube_x):
2425 # Check cubes are correct shape.
2426 cube_iter = check_single_cube(cube_iter)
2427 if cube_iter.ndim > 1:
2428 raise ValueError("cube_x must be 1D.")
2430 # Iterate over all cubes in cube or CubeList and plot.
2431 for cube_iter in iter_maybe(cube_y):
2432 # Check cubes are correct shape.
2433 cube_iter = check_single_cube(cube_iter)
2434 if cube_iter.ndim > 1:
2435 raise ValueError("cube_y must be 1D.")
2437 # Ensure we have a name for the plot file.
2438 recipe_title = get_recipe_metadata().get("title", "Untitled")
2439 title = f"{recipe_title}"
2441 if filename is None:
2442 filename = slugify(recipe_title)
2444 # Add file extension.
2445 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2447 # Do the actual plotting.
2448 _plot_and_save_scatter_plot(cube_x, cube_y, plot_filename, title, one_to_one)
2450 # Add list of plots to plot metadata.
2451 plot_index = _append_to_plot_index([plot_filename])
2453 # Make a page to display the plots.
2454 _make_plot_html_page(plot_index)
2456 return iris.cube.CubeList([cube_x, cube_y])
2459def vector_plot(
2460 cube_u: iris.cube.Cube,
2461 cube_v: iris.cube.Cube,
2462 filename: str = None,
2463 sequence_coordinate: str = "time",
2464 **kwargs,
2465) -> iris.cube.CubeList:
2466 """Plot a vector plot based on the input u and v components."""
2467 recipe_title = get_recipe_metadata().get("title", "Untitled")
2469 # Cubes must have a matching sequence coordinate.
2470 try:
2471 # Check that the u and v cubes have the same sequence coordinate.
2472 if cube_u.coord(sequence_coordinate) != cube_v.coord(sequence_coordinate): 2472 ↛ anywhereline 2472 didn't jump anywhere: it always raised an exception.
2473 raise ValueError("Coordinates do not match.")
2474 except (iris.exceptions.CoordinateNotFoundError, ValueError) as err:
2475 raise ValueError(
2476 f"Cubes should have matching {sequence_coordinate} coordinate:\n{cube_u}\n{cube_v}"
2477 ) from err
2479 # Create a plot for each value of the sequence coordinate.
2480 plot_index = []
2481 nplot = np.size(cube_u[0].coord(sequence_coordinate).points)
2482 for cube_u_slice, cube_v_slice in zip(
2483 cube_u.slices_over(sequence_coordinate),
2484 cube_v.slices_over(sequence_coordinate),
2485 strict=True,
2486 ):
2487 # Format the coordinate value in a unit appropriate way.
2488 seq_coord = cube_u_slice.coord(sequence_coordinate)
2489 plot_title, plot_filename = _set_title_and_filename(
2490 seq_coord, nplot, recipe_title, filename
2491 )
2493 # Do the actual plotting.
2494 _plot_and_save_vector_plot(
2495 cube_u_slice,
2496 cube_v_slice,
2497 filename=plot_filename,
2498 title=plot_title,
2499 method="contourf",
2500 )
2501 plot_index.append(plot_filename)
2503 # Add list of plots to plot metadata.
2504 complete_plot_index = _append_to_plot_index(plot_index)
2506 # Make a page to display the plots.
2507 _make_plot_html_page(complete_plot_index)
2509 return iris.cube.CubeList([cube_u, cube_v])
2512def plot_histogram_series(
2513 cubes: iris.cube.Cube | iris.cube.CubeList,
2514 filename: str = None,
2515 sequence_coordinate: str = "time",
2516 stamp_coordinate: str = "realization",
2517 single_plot: bool = False,
2518 **kwargs,
2519) -> iris.cube.Cube | iris.cube.CubeList:
2520 """Plot a histogram plot for each vertical level provided.
2522 A histogram plot can be plotted, but if the sequence_coordinate (i.e. time)
2523 is present then a sequence of plots will be produced using the time slider
2524 functionality to scroll through histograms against time. If a
2525 stamp_coordinate is present then postage stamp plots will be produced. If
2526 stamp_coordinate and single_plot is True, all postage stamp plots will be
2527 plotted in a single plot instead of separate postage stamp plots.
2529 Parameters
2530 ----------
2531 cubes: Cube | iris.cube.CubeList
2532 Iris cube or CubeList of the data to plot. It should have a single dimension other
2533 than the stamp coordinate.
2534 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2535 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2536 filename: str, optional
2537 Name of the plot to write, used as a prefix for plot sequences. Defaults
2538 to the recipe name.
2539 sequence_coordinate: str, optional
2540 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2541 This coordinate must exist in the cube and will be used for the time
2542 slider.
2543 stamp_coordinate: str, optional
2544 Coordinate about which to plot postage stamp plots. Defaults to
2545 ``"realization"``.
2546 single_plot: bool, optional
2547 If True, all postage stamp plots will be plotted in a single plot. If
2548 False, each postage stamp plot will be plotted separately. Is only valid
2549 if stamp_coordinate exists and has more than a single point.
2551 Returns
2552 -------
2553 iris.cube.Cube | iris.cube.CubeList
2554 The original Cube or CubeList (so further operations can be applied).
2555 Plotted data.
2557 Raises
2558 ------
2559 ValueError
2560 If the cube doesn't have the right dimensions.
2561 TypeError
2562 If the cube isn't a Cube or CubeList.
2563 """
2564 recipe_title = get_recipe_metadata().get("title", "Untitled")
2566 cubes = iter_maybe(cubes)
2568 # Internal plotting function.
2569 plotting_func = _plot_and_save_histogram_series
2571 num_models = get_num_models(cubes)
2573 validate_cube_shape(cubes, num_models)
2575 # If several histograms are plotted, check sequence_coordinate
2576 check_sequence_coordinate(cubes, sequence_coordinate)
2578 # Get axis minimum and maximum from levels information.
2579 # If no levels set, derive minima and maxima from data in CubeList.
2580 vmin, vmax = _set_axis_range(cubes)
2582 # Make postage stamp plots if stamp_coordinate exists and has more than a
2583 # single point. If single_plot is True:
2584 # -- all postage stamp plots will be plotted in a single plot instead of
2585 # separate postage stamp plots.
2586 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
2587 # produced per single model only
2588 if num_models == 1: 2588 ↛ 2601line 2588 didn't jump to line 2601 because the condition on line 2588 was always true
2589 if ( 2589 ↛ 2593line 2589 didn't jump to line 2593 because the condition on line 2589 was never true
2590 stamp_coordinate in [c.name() for c in cubes[0].coords()]
2591 and cubes[0].coord(stamp_coordinate).shape[0] > 1
2592 ):
2593 if single_plot:
2594 plotting_func = (
2595 _plot_and_save_postage_stamps_in_single_plot_histogram_series
2596 )
2597 else:
2598 plotting_func = _plot_and_save_postage_stamp_histogram_series
2599 cube_iterables = cubes[0].slices_over(sequence_coordinate)
2600 else:
2601 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2603 plot_index = []
2604 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2605 # Create a plot for each value of the sequence coordinate. Allowing for
2606 # multiple cubes in a CubeList to be plotted in the same plot for similar
2607 # sequence values. Passing a CubeList into the internal plotting function
2608 # for similar values of the sequence coordinate. cube_slice can be an
2609 # iris.cube.Cube or an iris.cube.CubeList.
2610 for cube_slice in cube_iterables:
2611 single_cube = cube_slice
2612 if isinstance(cube_slice, iris.cube.CubeList): 2612 ↛ 2613line 2612 didn't jump to line 2613 because the condition on line 2612 was never true
2613 single_cube = cube_slice[0]
2615 # Ensure valid stamp coordinate in cube dimensions
2616 if stamp_coordinate == "realization": 2616 ↛ 2619line 2616 didn't jump to line 2619 because the condition on line 2616 was always true
2617 stamp_coordinate = check_stamp_coordinate(single_cube)
2618 # Set plot titles and filename, based on sequence coordinate
2619 seq_coord = single_cube.coord(sequence_coordinate)
2620 # Use time coordinate in title and filename if single histogram output.
2621 if sequence_coordinate == "realization" and nplot == 1: 2621 ↛ 2622line 2621 didn't jump to line 2622 because the condition on line 2621 was never true
2622 seq_coord = single_cube.coord("time")
2623 plot_title, plot_filename = _set_title_and_filename(
2624 seq_coord, nplot, recipe_title, filename
2625 )
2627 # Do the actual plotting.
2628 plotting_func(
2629 cube_slice,
2630 filename=plot_filename,
2631 stamp_coordinate=stamp_coordinate,
2632 title=plot_title,
2633 vmin=vmin,
2634 vmax=vmax,
2635 )
2636 plot_index.append(plot_filename)
2638 # Add list of plots to plot metadata.
2639 complete_plot_index = _append_to_plot_index(plot_index)
2641 # Make a page to display the plots.
2642 _make_plot_html_page(complete_plot_index)
2644 return cubes
2647def plot_power_spectrum_series(
2648 cubes: iris.cube.Cube | iris.cube.CubeList,
2649 filename: str = None,
2650 sequence_coordinate: str = "time",
2651 stamp_coordinate: str = "realization",
2652 single_plot: bool = False,
2653 **kwargs,
2654) -> iris.cube.Cube | iris.cube.CubeList:
2655 """Plot a power spectrum plot for each vertical level provided.
2657 A power spectrum plot can be plotted, but if the sequence_coordinate (i.e. time)
2658 is present then a sequence of plots will be produced using the time slider
2659 functionality to scroll through power spectrum against time. If a
2660 stamp_coordinate is present then postage stamp plots will be produced. If
2661 stamp_coordinate and single_plot is True, all postage stamp plots will be
2662 plotted in a single plot instead of separate postage stamp plots.
2664 Parameters
2665 ----------
2666 cubes: Cube | iris.cube.CubeList
2667 Iris cube or CubeList of the data to plot. It should have a single dimension other
2668 than the stamp coordinate.
2669 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2670 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2671 filename: str, optional
2672 Name of the plot to write, used as a prefix for plot sequences. Defaults
2673 to the recipe name.
2674 sequence_coordinate: str, optional
2675 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2676 This coordinate must exist in the cube and will be used for the time
2677 slider.
2678 stamp_coordinate: str, optional
2679 Coordinate about which to plot postage stamp plots. Defaults to
2680 ``"realization"``.
2681 single_plot: bool, optional
2682 If True, all postage stamp plots will be plotted in a single plot. If
2683 False, each postage stamp plot will be plotted separately. Is only valid
2684 if stamp_coordinate exists and has more than a single point.
2686 Returns
2687 -------
2688 iris.cube.Cube | iris.cube.CubeList
2689 The original Cube or CubeList (so further operations can be applied).
2690 Plotted data.
2692 Raises
2693 ------
2694 ValueError
2695 If the cube doesn't have the right dimensions.
2696 TypeError
2697 If the cube isn't a Cube or CubeList.
2698 """
2699 recipe_title = get_recipe_metadata().get("title", "Untitled")
2701 cubes = iter_maybe(cubes)
2703 # Internal plotting function.
2704 plotting_func = _plot_and_save_power_spectrum_series
2706 num_models = get_num_models(cubes)
2708 validate_cube_shape(cubes, num_models)
2710 # If several power spectra are plotted, check sequence_coordinate
2711 check_sequence_coordinate(cubes, sequence_coordinate)
2713 # Make postage stamp plots if stamp_coordinate exists and has more than a
2714 # single point. If single_plot is True:
2715 # -- all postage stamp plots will be plotted in a single plot instead of
2716 # separate postage stamp plots.
2717 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
2718 # produced per single model only
2719 if num_models == 1: 2719 ↛ 2732line 2719 didn't jump to line 2732 because the condition on line 2719 was always true
2720 if ( 2720 ↛ 2724line 2720 didn't jump to line 2724 because the condition on line 2720 was never true
2721 stamp_coordinate in [c.name() for c in cubes[0].coords()]
2722 and cubes[0].coord(stamp_coordinate).shape[0] > 1
2723 ):
2724 if single_plot:
2725 plotting_func = (
2726 _plot_and_save_postage_stamps_in_single_plot_power_spectrum_series
2727 )
2728 else:
2729 plotting_func = _plot_and_save_postage_stamp_power_spectrum_series
2730 cube_iterables = cubes[0].slices_over(sequence_coordinate)
2731 else:
2732 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2734 plot_index = []
2735 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2736 # Create a plot for each value of the sequence coordinate. Allowing for
2737 # multiple cubes in a CubeList to be plotted in the same plot for similar
2738 # sequence values. Passing a CubeList into the internal plotting function
2739 # for similar values of the sequence coordinate. cube_slice can be an
2740 # iris.cube.Cube or an iris.cube.CubeList.
2741 for cube_slice in cube_iterables:
2742 single_cube = cube_slice
2743 if isinstance(cube_slice, iris.cube.CubeList): 2743 ↛ 2744line 2743 didn't jump to line 2744 because the condition on line 2743 was never true
2744 single_cube = cube_slice[0]
2746 # Set stamp coordinate
2747 if stamp_coordinate == "realization": 2747 ↛ 2750line 2747 didn't jump to line 2750 because the condition on line 2747 was always true
2748 stamp_coordinate = check_stamp_coordinate(single_cube)
2749 # Set plot title and filenames based on sequence values
2750 seq_coord = single_cube.coord(sequence_coordinate)
2751 plot_title, plot_filename = _set_title_and_filename(
2752 seq_coord, nplot, recipe_title, filename
2753 )
2755 # Do the actual plotting.
2756 plotting_func(
2757 cube_slice,
2758 filename=plot_filename,
2759 stamp_coordinate=stamp_coordinate,
2760 title=plot_title,
2761 )
2762 plot_index.append(plot_filename)
2764 # Add list of plots to plot metadata.
2765 complete_plot_index = _append_to_plot_index(plot_index)
2767 # Make a page to display the plots.
2768 _make_plot_html_page(complete_plot_index)
2770 return cubes