Coverage for src/CSET/operators/plot.py: 84%
943 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-19 11:18 +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.nanmin(cube.coord(lon_axis).points)
162 x2 = np.nanmax(cube.coord(lon_axis).points)
163 y1 = np.nanmin(cube.coord(lat_axis).points)
164 y2 = np.nanmax(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 # Assume polar projection for regional grids encompassing N. Pole
202 if y1 > 20.0 and y2 > 80.0:
203 projection = ccrs.NorthPolarStereo(central_longitude=0.0)
204 crs = ccrs.PlateCarree()
205 elif y1 < -80.0 and y2 < -20.0:
206 projection = ccrs.SouthPolarStereo(central_longitude=0.0)
207 crs = ccrs.PlateCarree()
208 # Define regular map projection for non-rotated pole inputs.
209 # Alternatives might include e.g. for global model outputs:
210 else:
211 projection = ccrs.Robinson(central_longitude=0.0, globe=None)
212 projection = ccrs.NearsidePerspective(
213 central_longitude=180.0,
214 central_latitude=0,
215 satellite_height=35785831,
216 )
217 # See also https://scitools.org.uk/cartopy/docs/v0.15/crs/projections.html.
218 # else:
219 projection = ccrs.PlateCarree(central_longitude=central_longitude)
220 crs = ccrs.PlateCarree()
222 # Define axes for plot (or subplot) with required map projection.
223 if subplot is not None:
224 axes = figure.add_subplot(
225 grid_size[0], grid_size[1], subplot, projection=projection
226 )
227 else:
228 axes = figure.add_subplot(projection=projection)
230 # Add coastlines and borderlines if cube contains x and y map coordinates.
231 # Avoid adding lines for 2D masked data or specific fixed ancillary spatial plots.
232 if (cube.ndim > 1 and iris.util.is_masked(cube.data)) or any(
233 name in cube.name() for name in ["land_", "orography", "altitude"]
234 ):
235 pass
236 else:
237 if cmap.name in ["viridis", "Greys"]:
238 coastcol = "magenta"
239 else:
240 coastcol = "black"
241 logging.debug("Plotting coastlines and borderlines in colour %s.", coastcol)
242 axes.coastlines(resolution="10m", color=coastcol, alpha=0.8)
243 axes.add_feature(cfeature.BORDERS, edgecolor=coastcol, alpha=0.3)
245 # Add gridlines.
246 gl = axes.gridlines(
247 alpha=0.3,
248 draw_labels=True,
249 dms=False,
250 x_inline=False,
251 y_inline=False,
252 )
253 gl.top_labels = False
254 gl.right_labels = False
255 if subplot:
256 gl.bottom_labels = False
257 gl.left_labels = False
258 if subplot % grid_size[1] == 1:
259 gl.left_labels = True
260 if subplot > ((grid_size[0] - 1) * grid_size[1]): 260 ↛ 265line 260 didn't jump to line 265 because the condition on line 260 was always true
261 gl.bottom_labels = True
263 # If is lat/lon spatial map, fix extent to keep plot tight.
264 # Specifying crs within set_extent helps ensure only data region is shown.
265 if isinstance(coord_system, iris.coord_systems.GeogCS):
266 axes.set_extent([x1, x2, y1, y2], crs=crs)
268 except ValueError:
269 # Skip if not both x and y map coordinates.
270 axes = figure.gca()
271 pass
273 return axes
276def _get_plot_resolution() -> int:
277 """Get resolution of rasterised plots in pixels per inch."""
278 return get_recipe_metadata().get("plot_resolution", 100)
281def _get_start_end_strings(seq_coord: iris.coords.Coord, use_bounds: bool):
282 """Return title and filename based on start and end points or bounds."""
283 if use_bounds and seq_coord.has_bounds():
284 vals = seq_coord.bounds.flatten()
285 else:
286 vals = seq_coord.points
287 start = seq_coord.units.title(vals[0])
288 end = seq_coord.units.title(vals[-1])
290 if start == end:
291 sequence_title = f"\n [{start}]"
292 sequence_fname = f"_{filename_slugify(start)}"
293 else:
294 sequence_title = f"\n [{start} to {end}]"
295 sequence_fname = f"_{filename_slugify(start)}_{filename_slugify(end)}"
297 # Do not include time if coord set to zero.
298 if (
299 seq_coord.units == "hours since 0001-01-01 00:00:00"
300 and vals[0] == 0
301 and vals[-1] == 0
302 ):
303 sequence_title = ""
304 sequence_fname = ""
306 return sequence_title, sequence_fname
309def _set_title_and_filename(
310 seq_coord: iris.coords.Coord,
311 nplot: int,
312 recipe_title: str,
313 filename: str,
314):
315 """Set plot title and filename based on cube coordinate.
317 Parameters
318 ----------
319 sequence_coordinate: iris.coords.Coord
320 Coordinate about which to make a plot sequence.
321 nplot: int
322 Number of output plots to generate - controls title/naming.
323 recipe_title: str
324 Default plot title, potentially to update.
325 filename: str
326 Input plot filename, potentially to update.
328 Returns
329 -------
330 plot_title: str
331 Output formatted plot title string, based on plotted data.
332 plot_filename: str
333 Output formatted plot filename string.
334 """
335 ndim = seq_coord.ndim
336 npoints = np.size(seq_coord.points)
337 sequence_title = ""
338 sequence_fname = ""
340 # Case 1: Multiple dimension sequence input - list number of aggregated cases
341 # (e.g. aggregation histogram plots)
342 if ndim > 1:
343 ncase = np.shape(seq_coord)[0]
344 sequence_title = f"\n [{ncase} cases]"
345 sequence_fname = f"_{ncase}cases"
347 # Case 2: Single dimension input
348 else:
349 # Single sequence point
350 if npoints == 1:
351 if nplot > 1:
352 # Default labels for sequence inputs
353 sequence_value = seq_coord.units.title(seq_coord.points[0])
354 sequence_value = sequence_value.replace(" unknown", "")
355 sequence_title = f"\n [{sequence_value}]"
356 sequence_fname = f"_{filename_slugify(sequence_value)}"
357 else:
358 # Aggregated attribute available where input collapsed over aggregation
359 try:
360 ncase = seq_coord.attributes["number_reference_times"]
361 sequence_title = f"\n [{ncase} cases]"
362 sequence_fname = f"_{ncase}cases"
363 except KeyError:
364 sequence_title, sequence_fname = _get_start_end_strings(
365 seq_coord, use_bounds=seq_coord.has_bounds()
366 )
367 # Multiple sequence (e.g. time) points
368 else:
369 sequence_title, sequence_fname = _get_start_end_strings(
370 seq_coord, use_bounds=False
371 )
373 # Set plot title and filename
374 plot_title = f"{recipe_title}{sequence_title}"
376 # Set plot filename, defaulting to user input if provided.
377 if filename is None:
378 filename = slugify(recipe_title)
379 plot_filename = f"{filename.rsplit('.', 1)[0]}{sequence_fname}.png"
380 else:
381 if nplot > 1:
382 plot_filename = f"{filename.rsplit('.', 1)[0]}{sequence_fname}.png"
383 else:
384 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
386 return plot_title, plot_filename
389def _set_postage_stamp_title(stamp_coord: iris.coords.Coord) -> str:
390 """Control postage stamp plot output titles based on stamp coordinate."""
391 if stamp_coord.name() == "realization":
392 mtitle = "Member"
393 else:
394 mtitle = stamp_coord.name().capitalize()
396 if stamp_coord.name() == "time":
397 mtitle = f"{stamp_coord.units.title(stamp_coord.points[0])}"
398 else:
399 mtitle = f"{mtitle} #{stamp_coord.points[0]}"
401 return mtitle
404def _set_axis_range(cubes):
405 """Get minimum and maximum from levels information."""
406 levels = None
407 for cube in cubes: 407 ↛ 423line 407 didn't jump to line 423 because the loop on line 407 didn't complete
408 # First check if user-specified "auto" range variable.
409 # This maintains the value of levels as None, so proceed.
410 _, levels, _ = colorbar_map_levels(cube, axis="y")
411 if levels is None:
412 break
413 # If levels is changed, recheck to use the vmin,vmax or
414 # levels-based ranges for histogram plots.
415 _, levels, _ = colorbar_map_levels(cube)
416 logging.debug("levels: %s", levels)
417 if levels is not None: 417 ↛ 407line 417 didn't jump to line 407 because the condition on line 417 was always true
418 vmin = min(levels)
419 vmax = max(levels)
420 logging.debug("Updated vmin, vmax: %s, %s", vmin, vmax)
421 break
423 if levels is None:
424 vmin = min(cb.data.min() for cb in cubes)
425 vmax = max(cb.data.max() for cb in cubes)
427 return vmin, vmax
430def _find_matched_slices(cubes, sequence_coordinate):
431 """Identify matched cubes in CubeList by sequence_coordinate values.
433 Ensures common points are compared for multiple cube inputs.
434 """
435 all_points = sorted(
436 set(
437 itertools.chain.from_iterable(
438 cb.coord(sequence_coordinate).points for cb in cubes
439 )
440 )
441 )
442 all_slices = list(
443 itertools.chain.from_iterable(
444 cb.slices_over(sequence_coordinate) for cb in cubes
445 )
446 )
447 # Matched slices (matched by seq coord point; it may happen that
448 # evaluated models do not cover the same seq coord range, hence matching
449 # necessary)
450 cube_iterables = [
451 iris.cube.CubeList(
452 s for s in all_slices if s.coord(sequence_coordinate).points[0] == point
453 )
454 for point in all_points
455 ]
457 return cube_iterables
460def _plot_and_save_spatial_plot(
461 cube: iris.cube.Cube,
462 filename: str,
463 title: str,
464 method: Literal["contourf", "pcolormesh", "scatter"],
465 overlay_cube: iris.cube.Cube | None = None,
466 contour_cube: iris.cube.Cube | None = None,
467 scatter_cube: iris.cube.Cube | None = None,
468 **kwargs,
469):
470 """Plot and save a spatial plot.
472 Parameters
473 ----------
474 cube: Cube
475 2 dimensional (lat and lon) Cube of the data to plot.
476 filename: str
477 Filename of the plot to write.
478 title: str
479 Plot title.
480 method: "contourf" | "pcolormesh" | "scatter"
481 The plotting method to use.
482 overlay_cube: Cube, optional
483 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
484 contour_cube: Cube, optional
485 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
486 scatter_cube: Cube, optional
487 Optional 1 (station list) or 2 dimensional (lat and lon) Cube of data to overplot as map of scatter points over base cube
488 """
489 # Setup plot details, size, resolution, etc.
490 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
492 # Specify the color bar
493 cmap, levels, norm = colorbar_map_levels(cube)
495 # If overplotting, set required colorbars
496 if overlay_cube:
497 over_cmap, over_levels, over_norm = colorbar_map_levels(overlay_cube)
498 if contour_cube:
499 cntr_cmap, cntr_levels, cntr_norm = colorbar_map_levels(contour_cube)
501 # Setup plot map projection, extent and coastlines and borderlines.
502 axes = _setup_spatial_map(cube, fig, cmap)
504 # Set colorscale bounds
505 try:
506 vmin = min(levels)
507 vmax = max(levels)
508 except TypeError:
509 vmin, vmax = None, None
510 # Ensure to use norm and not vmin/vmax if levels are defined.
511 if norm is not None:
512 vmin = None
513 vmax = None
514 logging.debug("Plotting using defined levels.")
516 # Plot the field.
517 try:
518 vmin = min(levels)
519 vmax = max(levels)
520 except TypeError:
521 vmin, vmax = None, None
522 # Ensure to use norm and not vmin/vmax if levels are defined.
523 if norm is not None:
524 vmin = None
525 vmax = None
526 logging.debug("Plotting using defined levels.")
527 if method == "contourf":
528 plot = iplt.contourf(cube, cmap=cmap, levels=levels, norm=norm)
529 elif method == "pcolormesh":
530 plot = iplt.pcolormesh(cube, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
531 elif method == "scatter": 531 ↛ 539line 531 didn't jump to line 539 because the condition on line 531 was never true
532 # Scatter plot of the field. The marker size is chosen to give
533 # symbols that decrease in size as the number of observations
534 # increases, although the fraction of the figure covered by
535 # symbols increases roughly as N^(1/2), disregarding overlaps,
536 # and has been selected for the default figure size of (10, 10).
537 # Should this be changed, the marker size should be adjusted in
538 # proportion to the area of the figure.
539 mrk_size = int(np.sqrt(2500000.0 / len(cube.data)))
540 lat_axis, lon_axis = get_cube_yxcoordname(cube)
541 plot = iplt.scatter(
542 cube.coord(lon_axis),
543 cube.coord(lat_axis),
544 c=cube.data[:],
545 s=mrk_size,
546 cmap=cmap,
547 edgecolors="k",
548 norm=norm,
549 vmin=vmin,
550 vmax=vmax,
551 )
552 else:
553 raise ValueError(f"Unknown plotting method: {method}")
555 # Overplot overlay field, if required
556 if overlay_cube:
557 try:
558 over_vmin = min(over_levels)
559 over_vmax = max(over_levels)
560 except TypeError:
561 over_vmin, over_vmax = None, None
562 if over_norm is not None: 562 ↛ 563line 562 didn't jump to line 563 because the condition on line 562 was never true
563 over_vmin = None
564 over_vmax = None
565 overlay = iplt.pcolormesh(
566 overlay_cube,
567 cmap=over_cmap,
568 norm=over_norm,
569 alpha=0.8,
570 vmin=over_vmin,
571 vmax=over_vmax,
572 )
573 # Overplot contour field, if required, with contour labelling.
574 if contour_cube:
575 contour = iplt.contour(
576 contour_cube,
577 colors="darkgray",
578 levels=cntr_levels,
579 norm=cntr_norm,
580 alpha=0.5,
581 linestyles="--",
582 linewidths=1,
583 )
584 plt.clabel(contour)
585 # Overplot valid elements of scatter field, if required
586 if scatter_cube: 586 ↛ 587line 586 didn't jump to line 587 because the condition on line 586 was never true
587 mrk_size = int(np.sqrt(2500000.0 / len(scatter_cube.data)))
588 lat_axis, lon_axis = get_cube_yxcoordname(scatter_cube)
589 lon_coord = scatter_cube.coord(lon_axis)
590 lat_coord = scatter_cube.coord(lat_axis)
591 valid = ~scatter_cube.data.mask
592 valid_lon = iris.coords.AuxCoord(
593 lon_coord.points[valid],
594 standard_name=lon_coord.standard_name,
595 units=lon_coord.units,
596 coord_system=lon_coord.coord_system,
597 )
598 valid_lat = iris.coords.AuxCoord(
599 lat_coord.points[valid],
600 standard_name=lat_coord.standard_name,
601 units=lat_coord.units,
602 coord_system=lat_coord.coord_system,
603 )
604 iplt.scatter(
605 valid_lon,
606 valid_lat,
607 c=scatter_cube.data[valid],
608 s=mrk_size,
609 cmap=cmap,
610 edgecolors="k",
611 norm=norm,
612 vmin=vmin,
613 vmax=vmax,
614 )
616 # Check to see if transect, and if so, adjust y axis.
617 if is_transect(cube):
618 if "pressure" in [coord.name() for coord in cube.coords()]:
619 axes.invert_yaxis()
620 axes.set_yscale("log")
621 axes.set_ylim(1100, 100)
622 # If both model_level_number and level_height exists, iplt can construct
623 # plot as a function of height above orography (NOT sea level).
624 elif {"model_level_number", "level_height"}.issubset( 624 ↛ 629line 624 didn't jump to line 629 because the condition on line 624 was always true
625 {coord.name() for coord in cube.coords()}
626 ):
627 axes.set_yscale("log")
629 axes.set_title(
630 f"{title}\n"
631 f"Start Lat: {cube.attributes['transect_coords'].split('_')[0]}"
632 f" Start Lon: {cube.attributes['transect_coords'].split('_')[1]}"
633 f" End Lat: {cube.attributes['transect_coords'].split('_')[2]}"
634 f" End Lon: {cube.attributes['transect_coords'].split('_')[3]}",
635 fontsize=16,
636 )
638 # Inset code
639 axins = inset_axes(
640 axes,
641 width="20%",
642 height="20%",
643 loc="upper right",
644 axes_class=GeoAxes,
645 axes_kwargs={"map_projection": ccrs.PlateCarree()},
646 )
648 # Slightly transparent to reduce plot blocking.
649 axins.patch.set_alpha(0.4)
651 axins.coastlines(resolution="50m")
652 axins.add_feature(cfeature.BORDERS, linewidth=0.3)
654 SLat, SLon, ELat, ELon = (
655 float(coord) for coord in cube.attributes["transect_coords"].split("_")
656 )
658 # Draw line between them
659 axins.plot(
660 [SLon, ELon], [SLat, ELat], color="black", transform=ccrs.PlateCarree()
661 )
663 # Plot points (note: lon, lat order for Cartopy)
664 axins.plot(SLon, SLat, marker="x", color="green", transform=ccrs.PlateCarree())
665 axins.plot(ELon, ELat, marker="x", color="red", transform=ccrs.PlateCarree())
667 lon_min, lon_max = sorted([SLon, ELon])
668 lat_min, lat_max = sorted([SLat, ELat])
670 # Midpoints
671 lon_mid = (lon_min + lon_max) / 2
672 lat_mid = (lat_min + lat_max) / 2
674 # Maximum half-range
675 half_range = max(lon_max - lon_min, lat_max - lat_min) / 2
676 if half_range == 0: # points identical → provide small default 676 ↛ 680line 676 didn't jump to line 680 because the condition on line 676 was always true
677 half_range = 1
679 # Set square extent
680 axins.set_extent(
681 [
682 lon_mid - half_range,
683 lon_mid + half_range,
684 lat_mid - half_range,
685 lat_mid + half_range,
686 ],
687 crs=ccrs.PlateCarree(),
688 )
690 # Ensure square aspect
691 axins.set_aspect("equal")
693 else:
694 # Add title.
695 axes.set_title(title, fontsize=16)
697 # Adjust padding if spatial plot or transect
698 if is_transect(cube):
699 yinfopad = -0.1
700 ycbarpad = 0.1
701 else:
702 yinfopad = 0.01
703 ycbarpad = 0.042
705 # Add watermark with min/max/mean. Currently not user togglable.
706 # In the bbox dictionary, fc and ec are hex colour codes for grey shade.
707 axes.annotate(
708 f"Min: {np.min(cube.data):.3g} Max: {np.max(cube.data):.3g} Mean: {np.mean(cube.data):.3g}",
709 xy=(0.025, yinfopad),
710 xycoords="axes fraction",
711 xytext=(-5, 5),
712 textcoords="offset points",
713 ha="left",
714 va="bottom",
715 size=11,
716 bbox=dict(boxstyle="round", fc="#cccccc", ec="#808080", alpha=0.9),
717 )
719 # Add secondary colour bar for overlay_cube field if required.
720 if overlay_cube:
721 cbarB = fig.colorbar(
722 overlay, orientation="horizontal", location="bottom", pad=0.0, shrink=0.7
723 )
724 cbarB.set_label(label=f"{overlay_cube.name()} ({overlay_cube.units})", size=14)
725 # add ticks and tick_labels for every levels if less than 20 levels exist
726 if over_levels is not None and len(over_levels) < 20: 726 ↛ 727line 726 didn't jump to line 727 because the condition on line 726 was never true
727 cbarB.set_ticks(over_levels)
728 cbarB.set_ticklabels([f"{level:.2f}" for level in over_levels])
729 if "rainfall" or "snowfall" or "visibility" in overlay_cube.name():
730 cbarB.set_ticklabels([f"{level:.3g}" for level in over_levels])
731 logging.debug("Set secondary colorbar ticks and labels.")
733 # Add main colour bar.
734 cbar = fig.colorbar(
735 plot, orientation="horizontal", location="bottom", pad=ycbarpad, shrink=0.7
736 )
738 cbar.set_label(label=f"{cube.name()} ({cube.units})", size=14)
739 # add ticks and tick_labels for every levels if less than 20 levels exist
740 if levels is not None and len(levels) < 20:
741 cbar.set_ticks(levels)
742 cbar.set_ticklabels([f"{level:.2f}" for level in levels])
743 if "rainfall" or "snowfall" or "visibility" in cube.name(): 743 ↛ 745line 743 didn't jump to line 745 because the condition on line 743 was always true
744 cbar.set_ticklabels([f"{level:.3g}" for level in levels])
745 logging.debug("Set colorbar ticks and labels.")
747 # Save plot.
748 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
749 logging.info("Saved spatial plot to %s", filename)
750 plt.close(fig)
753def _plot_and_save_postage_stamp_spatial_plot(
754 cube: iris.cube.Cube,
755 filename: str,
756 stamp_coordinate: str,
757 title: str,
758 method: Literal["contourf", "pcolormesh"],
759 overlay_cube: iris.cube.Cube | None = None,
760 contour_cube: iris.cube.Cube | None = None,
761 **kwargs,
762):
763 """Plot postage stamp spatial plots from an ensemble.
765 Parameters
766 ----------
767 cube: Cube
768 Iris cube of data to be plotted. It must have the stamp coordinate.
769 filename: str
770 Filename of the plot to write.
771 stamp_coordinate: str
772 Coordinate that becomes different plots.
773 method: "contourf" | "pcolormesh"
774 The plotting method to use.
775 overlay_cube: Cube, optional
776 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
777 contour_cube: Cube, optional
778 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
780 Raises
781 ------
782 ValueError
783 If the cube doesn't have the right dimensions.
784 """
785 # Use the smallest square grid that will fit the members.
786 nmember = len(cube.coord(stamp_coordinate).points)
787 grid_rows = int(math.sqrt(nmember))
788 grid_size = math.ceil(nmember / grid_rows)
790 fig = plt.figure(
791 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
792 )
794 # Specify the color bar
795 cmap, levels, norm = colorbar_map_levels(cube)
796 # If overplotting, set required colorbars
797 if overlay_cube: 797 ↛ 798line 797 didn't jump to line 798 because the condition on line 797 was never true
798 over_cmap, over_levels, over_norm = colorbar_map_levels(overlay_cube)
799 if contour_cube: 799 ↛ 800line 799 didn't jump to line 800 because the condition on line 799 was never true
800 cntr_cmap, cntr_levels, cntr_norm = colorbar_map_levels(contour_cube)
802 # Make a subplot for each member.
803 for member, subplot in zip(
804 cube.slices_over(stamp_coordinate),
805 range(1, grid_size * grid_rows + 1),
806 strict=False,
807 ):
808 # Setup subplot map projection, extent and coastlines and borderlines.
809 axes = _setup_spatial_map(
810 member, fig, cmap, grid_size=(grid_rows, grid_size), subplot=subplot
811 )
812 if method == "contourf":
813 # Filled contour plot of the field.
814 plot = iplt.contourf(member, cmap=cmap, levels=levels, norm=norm)
815 elif method == "pcolormesh":
816 if levels is not None:
817 vmin = min(levels)
818 vmax = max(levels)
819 else:
820 raise TypeError("Unknown vmin and vmax range.")
821 vmin, vmax = None, None
822 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
823 # if levels are defined.
824 if norm is not None: 824 ↛ 825line 824 didn't jump to line 825 because the condition on line 824 was never true
825 vmin = None
826 vmax = None
827 # pcolormesh plot of the field.
828 plot = iplt.pcolormesh(member, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
829 else:
830 raise ValueError(f"Unknown plotting method: {method}")
832 # Overplot overlay field, if required
833 if overlay_cube: 833 ↛ 834line 833 didn't jump to line 834 because the condition on line 833 was never true
834 try:
835 over_vmin = min(over_levels)
836 over_vmax = max(over_levels)
837 except TypeError:
838 over_vmin, over_vmax = None, None
839 if over_norm is not None:
840 over_vmin = None
841 over_vmax = None
842 iplt.pcolormesh(
843 overlay_cube[member.coord(stamp_coordinate).points[0]],
844 cmap=over_cmap,
845 norm=over_norm,
846 alpha=0.6,
847 vmin=over_vmin,
848 vmax=over_vmax,
849 )
850 # Overplot contour field, if required
851 if contour_cube: 851 ↛ 852line 851 didn't jump to line 852 because the condition on line 851 was never true
852 iplt.contour(
853 contour_cube[member.coord(stamp_coordinate).points[0]],
854 colors="darkgray",
855 levels=cntr_levels,
856 norm=cntr_norm,
857 alpha=0.6,
858 linestyles="--",
859 linewidths=1,
860 )
861 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
862 axes.set_title(f"{mtitle}")
864 # Put the shared colorbar in its own axes.
865 colorbar_axes = fig.add_axes([0.15, 0.05, 0.7, 0.03])
866 colorbar = fig.colorbar(
867 plot, colorbar_axes, orientation="horizontal", pad=0.042, shrink=0.7
868 )
869 colorbar.set_label(f"{cube.name()} ({cube.units})", size=14)
871 # Overall figure title.
872 fig.suptitle(title, fontsize=16)
874 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
875 logging.info("Saved contour postage stamp plot to %s", filename)
876 plt.close(fig)
879def _plot_and_save_line_series(
880 cubes: iris.cube.CubeList,
881 coords: list[iris.coords.Coord],
882 ensemble_coord: str,
883 filename: str,
884 title: str,
885 **kwargs,
886):
887 """Plot and save a 1D line series.
889 Parameters
890 ----------
891 cubes: Cube or CubeList
892 Cube or CubeList containing the cubes to plot on the y-axis.
893 coords: list[Coord]
894 Coordinates to plot on the x-axis, one per cube.
895 ensemble_coord: str
896 Ensemble coordinate in the cube.
897 filename: str
898 Filename of the plot to write.
899 title: str
900 Plot title.
901 """
902 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
904 model_colors_map = get_model_colors_map(cubes)
906 # Store min/max ranges.
907 y_levels = []
909 # Check match-up across sequence coords gives consistent sizes
910 validate_cubes_coords(cubes, coords)
912 for cube, coord in zip(cubes, coords, strict=True):
913 label = None
914 color = "black"
915 if model_colors_map:
916 label = cube.attributes.get("model_name")
917 color = model_colors_map.get(label)
918 for cube_slice in cube.slices_over(ensemble_coord):
919 # Label with (control) if part of an ensemble or not otherwise.
920 if cube_slice.coord(ensemble_coord).points == [0]:
921 iplt.plot(
922 coord,
923 cube_slice,
924 color=color,
925 marker="o",
926 ls="-",
927 lw=3,
928 label=f"{label} (control)"
929 if len(cube.coord(ensemble_coord).points) > 1
930 else label,
931 )
932 # Label with (perturbed) if part of an ensemble and not the control.
933 else:
934 iplt.plot(
935 coord,
936 cube_slice,
937 color=color,
938 ls="-",
939 lw=1.5,
940 alpha=0.75,
941 label=f"{label} (member)",
942 )
944 # Calculate the global min/max if multiple cubes are given.
945 _, levels, _ = colorbar_map_levels(cube, axis="y")
946 if levels is not None: 946 ↛ 947line 946 didn't jump to line 947 because the condition on line 946 was never true
947 y_levels.append(min(levels))
948 y_levels.append(max(levels))
950 # Get the current axes.
951 ax = plt.gca()
953 # Add some labels and tweak the style.
954 # check if cubes[0] works for single cube if not CubeList
955 if coords[0].name() == "time":
956 ax.set_xlabel(f"{coords[0].name()}", fontsize=14)
957 else:
958 ax.set_xlabel(f"{coords[0].name()} / {coords[0].units}", fontsize=14)
959 ax.set_ylabel(f"{cubes[0].name()} / {cubes[0].units}", fontsize=14)
960 ax.set_title(title, fontsize=16)
962 ax.ticklabel_format(axis="y", useOffset=False)
963 ax.tick_params(axis="x", labelrotation=15)
964 ax.tick_params(axis="both", labelsize=12)
966 # Set y limits to global min and max, autoscale if colorbar doesn't exist.
967 if y_levels: 967 ↛ 968line 967 didn't jump to line 968 because the condition on line 967 was never true
968 ax.set_ylim(min(y_levels), max(y_levels))
969 # Add zero line.
970 if min(y_levels) < 0.0 and max(y_levels) > 0.0:
971 ax.axhline(y=0, xmin=0, xmax=1, ls="-", color="grey", lw=2)
972 logging.debug(
973 "Line plot with y-axis limits %s-%s", min(y_levels), max(y_levels)
974 )
975 else:
976 ax.autoscale()
978 # Add gridlines
979 ax.grid(linestyle="--", color="grey", linewidth=1)
980 # Ientify unique labels for legend
981 handles = list(
982 {
983 label: handle
984 for (handle, label) in zip(*ax.get_legend_handles_labels(), strict=True)
985 }.values()
986 )
987 ax.legend(handles=handles, loc="best", ncol=1, frameon=True, fontsize=16)
989 # Save plot.
990 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
991 logging.info("Saved line plot to %s", filename)
992 plt.close(fig)
995def _plot_and_save_vertical_line_series(
996 cubes: iris.cube.CubeList,
997 coords: list[iris.coords.Coord],
998 ensemble_coord: str,
999 filename: str,
1000 series_coordinate: str,
1001 title: str,
1002 vmin: float,
1003 vmax: float,
1004 **kwargs,
1005):
1006 """Plot and save a 1D line series in vertical.
1008 Parameters
1009 ----------
1010 cubes: CubeList
1011 1 dimensional Cube or CubeList of the data to plot on x-axis.
1012 coord: list[Coord]
1013 Coordinates to plot on the y-axis, one per cube.
1014 ensemble_coord: str
1015 Ensemble coordinate in the cube.
1016 filename: str
1017 Filename of the plot to write.
1018 series_coordinate: str
1019 Coordinate to use as vertical axis.
1020 title: str
1021 Plot title.
1022 vmin: float
1023 Minimum value for the x-axis.
1024 vmax: float
1025 Maximum value for the x-axis.
1026 """
1027 # plot the vertical pressure axis using log scale
1028 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1030 model_colors_map = get_model_colors_map(cubes)
1032 # Check match-up across sequence coords gives consistent sizes
1033 validate_cubes_coords(cubes, coords)
1035 for cube, coord in zip(cubes, coords, strict=True):
1036 label = None
1037 color = "black"
1038 if model_colors_map: 1038 ↛ 1039line 1038 didn't jump to line 1039 because the condition on line 1038 was never true
1039 label = cube.attributes.get("model_name")
1040 color = model_colors_map.get(label)
1042 for cube_slice in cube.slices_over(ensemble_coord):
1043 # If ensemble data given plot control member with (control)
1044 # unless single forecast.
1045 if cube_slice.coord(ensemble_coord).points == [0]:
1046 iplt.plot(
1047 cube_slice,
1048 coord,
1049 color=color,
1050 marker="o",
1051 ls="-",
1052 lw=3,
1053 label=f"{label} (control)"
1054 if len(cube.coord(ensemble_coord).points) > 1
1055 else label,
1056 )
1057 # If ensemble data given plot perturbed members with (perturbed).
1058 else:
1059 iplt.plot(
1060 cube_slice,
1061 coord,
1062 color=color,
1063 ls="-",
1064 lw=1.5,
1065 alpha=0.75,
1066 label=f"{label} (member)",
1067 )
1069 # Get the current axis
1070 ax = plt.gca()
1072 # Special handling for pressure level data.
1073 if series_coordinate == "pressure": 1073 ↛ 1095line 1073 didn't jump to line 1095 because the condition on line 1073 was always true
1074 # Invert y-axis and set to log scale.
1075 ax.invert_yaxis()
1076 ax.set_yscale("log")
1078 # Define y-ticks and labels for pressure log axis.
1079 y_tick_labels = [
1080 "1000",
1081 "850",
1082 "700",
1083 "500",
1084 "300",
1085 "200",
1086 "100",
1087 ]
1088 y_ticks = [1000, 850, 700, 500, 300, 200, 100]
1090 # Set y-axis limits and ticks.
1091 ax.set_ylim(1100, 100)
1093 # Test if series_coordinate is model level data. The UM data uses
1094 # model_level_number and lfric uses full_levels as coordinate.
1095 elif series_coordinate in ("model_level_number", "full_levels", "half_levels"):
1096 # Define y-ticks and labels for vertical axis.
1097 y_ticks = iter_maybe(cubes)[0].coord(series_coordinate).points
1098 y_tick_labels = [str(int(i)) for i in y_ticks]
1099 ax.set_ylim(min(y_ticks), max(y_ticks))
1101 ax.set_yticks(y_ticks)
1102 ax.set_yticklabels(y_tick_labels)
1104 # Set x-axis limits.
1105 ax.set_xlim(vmin, vmax)
1106 # Mark y=0 if present in plot.
1107 if vmin < 0.0 and vmax > 0.0: 1107 ↛ 1108line 1107 didn't jump to line 1108 because the condition on line 1107 was never true
1108 ax.axvline(x=0, ymin=0, ymax=1, ls="-", color="grey", lw=2)
1110 # Add some labels and tweak the style.
1111 ax.set_ylabel(f"{coord.name()} / {coord.units}", fontsize=14)
1112 ax.set_xlabel(
1113 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
1114 )
1115 ax.set_title(title, fontsize=16)
1116 ax.ticklabel_format(axis="x")
1117 ax.tick_params(axis="y")
1118 ax.tick_params(axis="both", labelsize=12)
1120 # Add gridlines
1121 ax.grid(linestyle="--", color="grey", linewidth=1)
1122 # Ientify unique labels for legend
1123 handles = list(
1124 {
1125 label: handle
1126 for (handle, label) in zip(*ax.get_legend_handles_labels(), strict=True)
1127 }.values()
1128 )
1129 ax.legend(handles=handles, loc="best", ncol=1, frameon=True, fontsize=16)
1131 # Save plot.
1132 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1133 logging.info("Saved line plot to %s", filename)
1134 plt.close(fig)
1137def _plot_and_save_scatter_plot(
1138 cube_x: iris.cube.Cube | iris.cube.CubeList,
1139 cube_y: iris.cube.Cube | iris.cube.CubeList,
1140 filename: str,
1141 title: str,
1142 one_to_one: bool,
1143 model_names: list[str] = None,
1144 **kwargs,
1145):
1146 """Plot and save a 2D scatter plot.
1148 Parameters
1149 ----------
1150 cube_x: Cube | CubeList
1151 1 dimensional Cube or CubeList of the data to plot on x-axis.
1152 cube_y: Cube | CubeList
1153 1 dimensional Cube or CubeList of the data to plot on y-axis.
1154 filename: str
1155 Filename of the plot to write.
1156 title: str
1157 Plot title.
1158 one_to_one: bool
1159 Whether a 1:1 line is plotted.
1160 """
1161 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1162 # plot the cube_x and cube_y 1D fields as a scatter plot. If they are CubeLists this ensures
1163 # to pair each cube from cube_x with the corresponding cube from cube_y, allowing to iterate
1164 # over the pairs simultaneously.
1166 # Ensure cube_x and cube_y are iterable
1167 cube_x_iterable = iter_maybe(cube_x)
1168 cube_y_iterable = iter_maybe(cube_y)
1170 for cube_x_iter, cube_y_iter in zip(cube_x_iterable, cube_y_iterable, strict=True):
1171 iplt.scatter(cube_x_iter, cube_y_iter)
1172 if one_to_one is True:
1173 plt.plot(
1174 [
1175 np.nanmin([np.nanmin(cube_y.data), np.nanmin(cube_x.data)]),
1176 np.nanmax([np.nanmax(cube_y.data), np.nanmax(cube_x.data)]),
1177 ],
1178 [
1179 np.nanmin([np.nanmin(cube_y.data), np.nanmin(cube_x.data)]),
1180 np.nanmax([np.nanmax(cube_y.data), np.nanmax(cube_x.data)]),
1181 ],
1182 "k",
1183 linestyle="--",
1184 )
1185 ax = plt.gca()
1187 # Add some labels and tweak the style.
1188 if model_names is None:
1189 ax.set_xlabel(f"{cube_x[0].name()} / {cube_x[0].units}", fontsize=14)
1190 ax.set_ylabel(f"{cube_y[0].name()} / {cube_y[0].units}", fontsize=14)
1191 else:
1192 # Add the model names, these should be order of base (x) and other (y).
1193 ax.set_xlabel(
1194 f"{model_names[0]}_{cube_x[0].name()} / {cube_x[0].units}", fontsize=14
1195 )
1196 ax.set_ylabel(
1197 f"{model_names[1]}_{cube_y[0].name()} / {cube_y[0].units}", fontsize=14
1198 )
1199 ax.set_title(title, fontsize=16)
1200 ax.ticklabel_format(axis="y", useOffset=False)
1201 ax.tick_params(axis="x", labelrotation=15)
1202 ax.tick_params(axis="both", labelsize=12)
1203 ax.autoscale()
1205 # Save plot.
1206 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1207 logging.info("Saved scatter plot to %s", filename)
1208 plt.close(fig)
1211def _plot_and_save_vector_plot(
1212 cube_u: iris.cube.Cube,
1213 cube_v: iris.cube.Cube,
1214 filename: str,
1215 title: str,
1216 method: Literal["contourf", "pcolormesh"],
1217 **kwargs,
1218):
1219 """Plot and save a 2D vector plot.
1221 Parameters
1222 ----------
1223 cube_u: Cube
1224 2 dimensional Cube of u component of the data.
1225 cube_v: Cube
1226 2 dimensional Cube of v component of the data.
1227 filename: str
1228 Filename of the plot to write.
1229 title: str
1230 Plot title.
1231 """
1232 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1234 # Create a cube containing the magnitude of the vector field.
1235 cube_vec_mag = (cube_u**2 + cube_v**2) ** 0.5
1236 cube_vec_mag.rename(f"{cube_u.name()}_{cube_v.name()}_magnitude")
1238 # Specify the color bar
1239 cmap, levels, norm = colorbar_map_levels(cube_vec_mag)
1241 # Setup plot map projection, extent and coastlines and borderlines.
1242 axes = _setup_spatial_map(cube_vec_mag, fig, cmap)
1244 if method == "contourf":
1245 # Filled contour plot of the field.
1246 plot = iplt.contourf(cube_vec_mag, cmap=cmap, levels=levels, norm=norm)
1247 elif method == "pcolormesh":
1248 try:
1249 vmin = min(levels)
1250 vmax = max(levels)
1251 except TypeError:
1252 vmin, vmax = None, None
1253 # pcolormesh plot of the field and ensure to use norm and not vmin/vmax
1254 # if levels are defined.
1255 if norm is not None:
1256 vmin = None
1257 vmax = None
1258 plot = iplt.pcolormesh(cube_vec_mag, cmap=cmap, norm=norm, vmin=vmin, vmax=vmax)
1259 else:
1260 raise ValueError(f"Unknown plotting method: {method}")
1262 # Check to see if transect, and if so, adjust y axis.
1263 if is_transect(cube_vec_mag):
1264 if "pressure" in [coord.name() for coord in cube_vec_mag.coords()]:
1265 axes.invert_yaxis()
1266 axes.set_yscale("log")
1267 axes.set_ylim(1100, 100)
1268 # If both model_level_number and level_height exists, iplt can construct
1269 # plot as a function of height above orography (NOT sea level).
1270 elif {"model_level_number", "level_height"}.issubset(
1271 {coord.name() for coord in cube_vec_mag.coords()}
1272 ):
1273 axes.set_yscale("log")
1275 axes.set_title(
1276 f"{title}\n"
1277 f"Start Lat: {cube_vec_mag.attributes['transect_coords'].split('_')[0]}"
1278 f" Start Lon: {cube_vec_mag.attributes['transect_coords'].split('_')[1]}"
1279 f" End Lat: {cube_vec_mag.attributes['transect_coords'].split('_')[2]}"
1280 f" End Lon: {cube_vec_mag.attributes['transect_coords'].split('_')[3]}",
1281 fontsize=16,
1282 )
1284 else:
1285 # Add title.
1286 axes.set_title(title, fontsize=16)
1288 # Add watermark with min/max/mean. Currently not user togglable.
1289 # In the bbox dictionary, fc and ec are hex colour codes for grey shade.
1290 axes.annotate(
1291 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}",
1292 xy=(0.05, -0.05),
1293 xycoords="axes fraction",
1294 xytext=(-5, 5),
1295 textcoords="offset points",
1296 ha="right",
1297 va="bottom",
1298 size=11,
1299 bbox=dict(boxstyle="round", fc="#cccccc", ec="#808080", alpha=0.9),
1300 )
1302 # Add colour bar.
1303 cbar = fig.colorbar(plot, orientation="horizontal", pad=0.042, shrink=0.7)
1304 cbar.set_label(label=f"{cube_vec_mag.name()} ({cube_vec_mag.units})", size=14)
1305 # add ticks and tick_labels for every levels if less than 20 levels exist
1306 if levels is not None and len(levels) < 20:
1307 cbar.set_ticks(levels)
1308 cbar.set_ticklabels([f"{level:.1f}" for level in levels])
1310 # 30 barbs along the longest axis of the plot, or a barb per point for data
1311 # with less than 30 points.
1312 step = max(max(cube_u.shape) // 30, 1)
1313 iplt.quiver(cube_u[::step, ::step], cube_v[::step, ::step], pivot="middle")
1315 # Save plot.
1316 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1317 logging.info("Saved vector plot to %s", filename)
1318 plt.close(fig)
1321def _plot_and_save_histogram_series(
1322 cubes: iris.cube.Cube | iris.cube.CubeList,
1323 filename: str,
1324 title: str,
1325 vmin: float,
1326 vmax: float,
1327 **kwargs,
1328):
1329 """Plot and save a histogram series.
1331 Parameters
1332 ----------
1333 cubes: Cube or CubeList
1334 2 dimensional Cube or CubeList of the data to plot as histogram.
1335 filename: str
1336 Filename of the plot to write.
1337 title: str
1338 Plot title.
1339 vmin: float
1340 minimum for colorbar
1341 vmax: float
1342 maximum for colorbar
1343 """
1344 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1345 ax = plt.gca()
1347 model_colors_map = get_model_colors_map(cubes)
1349 # Set default that histograms will produce probability density function
1350 # at each bin (integral over range sums to 1).
1351 density = True
1353 for cube in iter_maybe(cubes):
1354 # Easier to check title (where var name originates)
1355 # than seeing if long names exist etc.
1356 # Exception case, where distribution better fits log scales/bins.
1357 if "surface_microphysical" in title:
1358 if "amount" in title: 1358 ↛ 1360line 1358 didn't jump to line 1360 because the condition on line 1358 was never true
1359 # Compute histogram following Klingaman et al. (2017): ASoP
1360 bin2 = np.exp(np.log(0.02) + 0.1 * np.linspace(0, 99, 100))
1361 bins = np.pad(bin2, (1, 0), "constant", constant_values=0)
1362 density = False
1363 else:
1364 bins = 10.0 ** (
1365 np.arange(-10, 27, 1) / 10.0
1366 ) # Suggestion from RMED toolbox.
1367 bins = np.insert(bins, 0, 0)
1368 ax.set_yscale("log")
1369 vmin = bins[1]
1370 vmax = bins[-1] # Manually set vmin/vmax to override json derived value.
1371 ax.set_xscale("log")
1372 elif "lightning" in title:
1373 bins = [0, 1, 2, 3, 4, 5]
1374 else:
1375 bins = np.linspace(vmin, vmax, 51)
1376 logging.debug(
1377 "Plotting histogram with %s bins %s - %s.",
1378 np.size(bins),
1379 np.min(bins),
1380 np.max(bins),
1381 )
1383 # Reshape cube data into a single array to allow for a single histogram.
1384 # Otherwise we plot xdim histograms stacked.
1385 cube_data_1d = (cube.data).flatten()
1387 label = None
1388 color = "black"
1389 if model_colors_map:
1390 label = cube.attributes.get("model_name")
1391 color = model_colors_map[label]
1392 x, y = np.histogram(cube_data_1d, bins=bins, density=density)
1394 # Compute area under curve.
1395 if "surface_microphysical" in title and "amount" in title: 1395 ↛ 1396line 1395 didn't jump to line 1396 because the condition on line 1395 was never true
1396 bin_mean = (bins[:-1] + bins[1:]) / 2.0
1397 x = x * bin_mean / x.sum()
1398 x = x[1:]
1399 y = y[1:]
1401 ax.plot(
1402 y[:-1], x, color=color, linewidth=3, marker="o", markersize=6, label=label
1403 )
1405 # Add some labels and tweak the style.
1406 ax.set_title(title, fontsize=16)
1407 ax.set_xlabel(
1408 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
1409 )
1410 ax.set_ylabel("Normalised probability density", fontsize=14)
1411 if "surface_microphysical" in title and "amount" in title: 1411 ↛ 1412line 1411 didn't jump to line 1412 because the condition on line 1411 was never true
1412 ax.set_ylabel(
1413 f"Contribution to mean ({iter_maybe(cubes)[0].units})", fontsize=14
1414 )
1415 ax.set_xlim(vmin, vmax)
1416 ax.tick_params(axis="both", labelsize=12)
1418 # Overlay grid-lines onto histogram plot.
1419 ax.grid(linestyle="--", color="grey", linewidth=1)
1420 if model_colors_map:
1421 ax.legend(loc="best", ncol=1, frameon=True, fontsize=16)
1423 # Save plot.
1424 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1425 logging.info("Saved histogram plot to %s", filename)
1426 plt.close(fig)
1429def _plot_and_save_postage_stamp_histogram_series(
1430 cube: iris.cube.Cube,
1431 filename: str,
1432 title: str,
1433 stamp_coordinate: str,
1434 vmin: float,
1435 vmax: float,
1436 **kwargs,
1437):
1438 """Plot and save postage (ensemble members) stamps for a histogram series.
1440 Parameters
1441 ----------
1442 cube: Cube
1443 2 dimensional Cube of the data to plot as histogram.
1444 filename: str
1445 Filename of the plot to write.
1446 title: str
1447 Plot title.
1448 stamp_coordinate: str
1449 Coordinate that becomes different plots.
1450 vmin: float
1451 minimum for pdf x-axis
1452 vmax: float
1453 maximum for pdf x-axis
1454 """
1455 # Use the smallest square grid that will fit the members.
1456 nmember = len(cube.coord(stamp_coordinate).points)
1457 grid_rows = int(math.sqrt(nmember))
1458 grid_size = math.ceil(nmember / grid_rows)
1460 fig = plt.figure(
1461 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
1462 )
1463 # Make a subplot for each member.
1464 for member, subplot in zip(
1465 cube.slices_over(stamp_coordinate),
1466 range(1, grid_size * grid_rows + 1),
1467 strict=False,
1468 ):
1469 # Implicit interface is much easier here, due to needing to have the
1470 # cartopy GeoAxes generated.
1471 plt.subplot(grid_rows, grid_size, subplot)
1472 # Reshape cube data into a single array to allow for a single histogram.
1473 # Otherwise we plot xdim histograms stacked.
1474 member_data_1d = (member.data).flatten()
1475 plt.hist(member_data_1d, density=True, stacked=True)
1476 axes = plt.gca()
1477 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
1478 axes.set_title(f"{mtitle}")
1479 axes.set_xlim(vmin, vmax)
1481 # Overall figure title.
1482 fig.suptitle(title, fontsize=16)
1484 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1485 logging.info("Saved histogram postage stamp plot to %s", filename)
1486 plt.close(fig)
1489def _plot_and_save_postage_stamps_in_single_plot_histogram_series(
1490 cube: iris.cube.Cube,
1491 filename: str,
1492 title: str,
1493 stamp_coordinate: str,
1494 vmin: float,
1495 vmax: float,
1496 **kwargs,
1497):
1498 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1499 ax.set_title(title, fontsize=16)
1500 ax.set_xlim(vmin, vmax)
1501 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1502 ax.set_ylabel("normalised probability density", fontsize=14)
1503 # Loop over all slices along the stamp_coordinate
1504 for member in cube.slices_over(stamp_coordinate):
1505 # Flatten the member data to 1D
1506 member_data_1d = member.data.flatten()
1507 # Plot the histogram using plt.hist
1508 mtitle = _set_postage_stamp_title(member.coord(stamp_coordinate))
1509 plt.hist(
1510 member_data_1d,
1511 density=True,
1512 stacked=True,
1513 label=f"{mtitle}",
1514 )
1516 # Add a legend
1517 ax.legend(fontsize=16)
1519 # Save the figure to a file
1520 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1521 logging.info("Saved histogram postage stamp plot to %s", filename)
1523 # Close the figure
1524 plt.close(fig)
1527def _plot_and_save_scatter_series(
1528 cubes: iris.cube.Cube | iris.cube.CubeList,
1529 filename: str,
1530 title: str,
1531 vmin: float,
1532 vmax: float,
1533 hexbin: bool,
1534 **kwargs,
1535):
1536 """Plot and save a scatter plot series.
1538 Parameters
1539 ----------
1540 cubes: Cube or CubeList
1541 2 dimensional Cube or CubeList of the data to plot as scatter.
1542 filename: str
1543 Filename of the plot to write.
1544 title: str
1545 Plot title.
1546 vmin: float
1547 minimum for colorbar
1548 vmax: float
1549 maximum for colorbar
1550 """
1551 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1552 ax = plt.gca()
1554 model_colors_map = get_model_colors_map(cubes)
1556 percentiles = np.arange(0, 100, 5)
1557 percentiles[0] = 1
1558 percentiles[-1] = 99
1559 quantiles = iris.cube.CubeList()
1561 for plottype in ["points", "quantiles"]:
1562 nplot = 0
1563 for cube in iter_maybe(cubes):
1564 label = None
1565 color = "black"
1566 if model_colors_map: 1566 ↛ 1571line 1566 didn't jump to line 1571 because the condition on line 1566 was always true
1567 label = cube.attributes.get("model_name")
1568 color = model_colors_map[label]
1570 # Plot all data points
1571 if plottype == "points":
1572 if nplot > 0:
1573 if hexbin:
1574 hb = plt.hexbin(
1575 cubes[0].data.flatten(),
1576 cube.data.flatten(),
1577 alpha=0.3,
1578 gridsize=100,
1579 mincnt=1,
1580 )
1581 else:
1582 plt.scatter(
1583 cubes[0].data.flatten(),
1584 cube.data.flatten(),
1585 color=color,
1586 marker="+",
1587 label=None,
1588 alpha=0.3,
1589 )
1591 elif plottype == "quantiles": 1591 ↛ 1610line 1591 didn't jump to line 1610 because the condition on line 1591 was always true
1592 # Construct Q-Q plot
1593 quantiles.append(
1594 cube.collapsed(
1595 cube.coords(dim_coords=True),
1596 iris.analysis.PERCENTILE,
1597 percent=percentiles,
1598 )
1599 )
1600 if nplot > 0:
1601 iplt.scatter(
1602 quantiles[0],
1603 quantiles[-1],
1604 color=color,
1605 marker="o",
1606 label=label,
1607 edgecolors="black",
1608 )
1610 nplot = nplot + 1
1612 # Add some labels and tweak the style.
1613 ax.set_title(title, fontsize=16)
1614 ax.set_xlabel(
1615 f"{iter_maybe(cubes)[0].name()} / {iter_maybe(cubes)[0].units}", fontsize=14
1616 )
1617 ax.set_ylabel(
1618 f"{iter_maybe(cubes)[1].name()} / {iter_maybe(cubes)[1].units}", fontsize=14
1619 )
1620 ax.tick_params(axis="both", labelsize=12)
1621 ax.autoscale()
1623 lims = [
1624 np.min([ax.get_xlim(), ax.get_ylim()]), # min of both axes
1625 np.max([ax.get_xlim(), ax.get_ylim()]), # max of both axes
1626 ]
1627 ax.plot(lims, lims, "k-", alpha=0.75, zorder=0)
1628 ax.set_aspect("equal")
1629 ax.set_xlim(lims)
1630 ax.set_ylim(lims)
1632 # Overlay grid-lines onto scatter plot.
1633 ax.grid(linestyle="--", color="grey", linewidth=1)
1634 if model_colors_map: 1634 ↛ 1638line 1634 didn't jump to line 1638 because the condition on line 1634 was always true
1635 ax.legend(loc="upper left", ncol=1, frameon=True, fontsize=16)
1637 # Add colorbar if hexbin output
1638 if hexbin:
1639 cb = plt.colorbar(
1640 hb, orientation="horizontal", location="bottom", pad=0.08, shrink=0.7
1641 )
1642 cb.set_label("Number of data points", size=12)
1644 # Save plot.
1645 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1646 logging.info("Saved scatter plot to %s", filename)
1647 plt.close(fig)
1650def _plot_and_save_power_spectrum_series(
1651 cubes: iris.cube.Cube | iris.cube.CubeList,
1652 filename: str,
1653 title: str,
1654 **kwargs,
1655):
1656 """Plot and save a power spectrum series.
1658 Parameters
1659 ----------
1660 cubes: Cube or CubeList
1661 2 dimensional Cube or CubeList of the data to plot as power spectrum.
1662 filename: str
1663 Filename of the plot to write.
1664 title: str
1665 Plot title.
1666 """
1667 fig = plt.figure(figsize=(10, 10), facecolor="w", edgecolor="k")
1668 ax = plt.gca()
1670 model_colors_map = get_model_colors_map(cubes)
1672 for cube in iter_maybe(cubes):
1673 # Calculate power spectrum
1675 # Extract time coordinate and convert to datetime
1676 time_coord = cube.coord("time")
1677 time_points = time_coord.units.num2date(time_coord.points)
1679 # Choose one time point (e.g., the first one)
1680 target_time = time_points[0]
1682 # Bind target_time inside the lambda using a default argument
1683 time_constraint = iris.Constraint(
1684 time=lambda cell, target_time=target_time: cell.point == target_time
1685 )
1687 cube = cube.extract(time_constraint)
1689 if cube.ndim == 2:
1690 cube_3d = cube.data[np.newaxis, :, :]
1691 logging.debug("Adding in new axis for a 2 dimensional cube.")
1692 elif cube.ndim == 3: 1692 ↛ 1693line 1692 didn't jump to line 1693 because the condition on line 1692 was never true
1693 cube_3d = cube.data
1694 else:
1695 raise ValueError("Cube dimensions unsuitable for power spectra code")
1696 raise ValueError(
1697 f"Cube is {cube.ndim} dimensional. Cube should be 2 or 3 dimensional."
1698 )
1700 # Calculate spectra
1701 ps_array = DCT_ps(cube_3d)
1703 ps_cube = iris.cube.Cube(
1704 ps_array,
1705 long_name="power_spectra",
1706 )
1708 ps_cube.attributes["model_name"] = cube.attributes.get("model_name")
1710 # Create a frequency/wavelength array for coordinate
1711 ps_len = ps_cube.data.shape[1]
1712 freqs = np.arange(1, ps_len + 1)
1713 freq_coord = iris.coords.DimCoord(freqs, long_name="frequency", units="1")
1715 # Convert datetime to numeric time using original units
1716 numeric_time = time_coord.units.date2num(time_points)
1717 # Create a new DimCoord with numeric time
1718 new_time_coord = iris.coords.DimCoord(
1719 numeric_time, standard_name="time", units=time_coord.units
1720 )
1722 # Add time and frequency coordinate to spectra cube.
1723 ps_cube.add_dim_coord(new_time_coord.copy(), 0)
1724 ps_cube.add_dim_coord(freq_coord.copy(), 1)
1726 # Extract data from the cube
1727 frequency = ps_cube.coord("frequency").points
1728 power_spectrum = ps_cube.data
1730 label = None
1731 color = "black"
1732 if model_colors_map: 1732 ↛ 1733line 1732 didn't jump to line 1733 because the condition on line 1732 was never true
1733 label = ps_cube.attributes.get("model_name")
1734 color = model_colors_map[label]
1735 ax.plot(frequency, power_spectrum[0], color=color, label=label)
1737 # Add some labels and tweak the style.
1738 ax.set_title(title, fontsize=16)
1739 ax.set_xlabel("Wavenumber", fontsize=14)
1740 ax.set_ylabel("Power", fontsize=14)
1741 ax.tick_params(axis="both", labelsize=12)
1743 # Set log-log scale
1744 ax.set_xscale("log")
1745 ax.set_yscale("log")
1747 # Overlay grid-lines onto power spectrum plot.
1748 ax.grid(linestyle="--", color="grey", linewidth=1)
1749 if model_colors_map: 1749 ↛ 1750line 1749 didn't jump to line 1750 because the condition on line 1749 was never true
1750 ax.legend(loc="best", ncol=1, frameon=True, fontsize=16)
1752 # Save plot.
1753 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1754 logging.info("Saved power spectrum plot to %s", filename)
1755 plt.close(fig)
1758def _plot_and_save_postage_stamp_power_spectrum_series(
1759 cube: iris.cube.Cube,
1760 filename: str,
1761 title: str,
1762 stamp_coordinate: str,
1763 **kwargs,
1764):
1765 """Plot and save postage (ensemble members) stamps for a power spectrum series.
1767 Parameters
1768 ----------
1769 cube: Cube
1770 2 dimensional Cube of the data to plot as power spectrum.
1771 filename: str
1772 Filename of the plot to write.
1773 title: str
1774 Plot title.
1775 stamp_coordinate: str
1776 Coordinate that becomes different plots.
1777 """
1778 # Use the smallest square grid that will fit the members.
1779 nmember = len(cube.coord(stamp_coordinate).points)
1780 grid_rows = int(math.sqrt(nmember))
1781 grid_size = math.ceil(nmember / grid_rows)
1783 fig = plt.figure(
1784 figsize=(10, 10 * max(grid_rows / grid_size, 0.5)), facecolor="w", edgecolor="k"
1785 )
1787 # Make a subplot for each member.
1788 for member, subplot in zip(
1789 cube.slices_over(stamp_coordinate),
1790 range(1, grid_size * grid_rows + 1),
1791 strict=False,
1792 ):
1793 # Implicit interface is much easier here, due to needing to have the
1794 # cartopy GeoAxes generated.
1795 plt.subplot(grid_rows, grid_size, subplot)
1797 frequency = member.coord("frequency").points
1799 axes = plt.gca()
1800 axes.plot(frequency, member.data)
1801 axes.set_title(f"Member #{member.coord(stamp_coordinate).points[0]}")
1803 # Overall figure title.
1804 fig.suptitle(title, fontsize=16)
1806 fig.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1807 logging.info("Saved power spectra postage stamp plot to %s", filename)
1808 plt.close(fig)
1811def _plot_and_save_postage_stamps_in_single_plot_power_spectrum_series(
1812 cube: iris.cube.Cube,
1813 filename: str,
1814 title: str,
1815 stamp_coordinate: str,
1816 **kwargs,
1817):
1818 fig, ax = plt.subplots(figsize=(10, 10), facecolor="w", edgecolor="k")
1819 ax.set_title(title, fontsize=16)
1820 ax.set_xlabel(f"{cube.name()} / {cube.units}", fontsize=14)
1821 ax.set_ylabel("Power", fontsize=14)
1822 # Loop over all slices along the stamp_coordinate
1823 for member in cube.slices_over(stamp_coordinate):
1824 frequency = member.coord("frequency").points
1825 ax.plot(
1826 frequency,
1827 member.data,
1828 label=f"Member #{member.coord(stamp_coordinate).points[0]}",
1829 )
1831 # Add a legend
1832 ax.legend(fontsize=16)
1834 # Save the figure to a file
1835 plt.savefig(filename, bbox_inches="tight", dpi=_get_plot_resolution())
1836 logging.info("Saved power spectra plot to %s", filename)
1838 # Close the figure
1839 plt.close(fig)
1842def _spatial_plot(
1843 method: Literal["contourf", "pcolormesh"],
1844 cube: iris.cube.Cube,
1845 filename: str | None,
1846 sequence_coordinate: str,
1847 stamp_coordinate: str,
1848 overlay_cube: iris.cube.Cube | None = None,
1849 contour_cube: iris.cube.Cube | None = None,
1850 scatter_cube: iris.cube.Cube | None = None,
1851 **kwargs,
1852):
1853 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1855 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1856 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1857 is present then postage stamp plots will be produced.
1859 If an overlay_cube and/or contour_cube and/or scatter_cube are specified,
1860 multiple variables or outputs can be overplotted on the same figure.
1862 Parameters
1863 ----------
1864 method: "contourf" | "pcolormesh"
1865 The plotting method to use.
1866 cube: Cube
1867 Iris cube of the data to plot. It should have two spatial dimensions,
1868 such as lat and lon, and may also have a another two dimension to be
1869 plotted sequentially and/or as postage stamp plots.
1870 filename: str | None
1871 Name of the plot to write, used as a prefix for plot sequences. If None
1872 uses the recipe name.
1873 sequence_coordinate: str
1874 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1875 This coordinate must exist in the cube.
1876 stamp_coordinate: str
1877 Coordinate about which to plot postage stamp plots. Defaults to
1878 ``"realization"``.
1879 overlay_cube: Cube | None, optional
1880 Optional 2 dimensional (lat and lon) Cube of data to overplot on top of base cube
1881 contour_cube: Cube | None, optional
1882 Optional 2 dimensional (lat and lon) Cube of data to overplot as contours over base cube
1883 scatter_cube: Cube | None, optional
1884 Optional 1 dimensional (station) or 2 dimensional (lat and lon) Cube of data to overplot as points on scatter map over base cube
1886 Raises
1887 ------
1888 ValueError
1889 If the cube doesn't have the right dimensions.
1890 TypeError
1891 If the cube isn't a single cube.
1892 """
1893 recipe_title = get_recipe_metadata().get("title", "Untitled")
1895 # Ensure we've got a single cube.
1896 cube = check_single_cube(cube)
1898 # Check if there is a valid stamp coordinate in cube dimensions.
1899 if stamp_coordinate == "realization": 1899 ↛ 1904line 1899 didn't jump to line 1904 because the condition on line 1899 was always true
1900 stamp_coordinate = check_stamp_coordinate(cube)
1902 # Make postage stamp plots if stamp_coordinate exists and has more than a
1903 # single point.
1904 plotting_func = _plot_and_save_spatial_plot
1905 try:
1906 if cube.coord(stamp_coordinate).shape[0] > 1:
1907 plotting_func = _plot_and_save_postage_stamp_spatial_plot
1908 except iris.exceptions.CoordinateNotFoundError:
1909 pass
1911 # Produce a geographical scatter plot if the data have a
1912 # dimension called observation or model_obs_error
1913 if any( 1913 ↛ 1917line 1913 didn't jump to line 1917 because the condition on line 1913 was never true
1914 crd.var_name == "station" or crd.var_name == "model_obs_error"
1915 for crd in cube.coords()
1916 ):
1917 plotting_func = _plot_and_save_spatial_plot
1918 method = "scatter"
1920 # Must have a sequence coordinate.
1921 try:
1922 cube.coord(sequence_coordinate)
1923 except iris.exceptions.CoordinateNotFoundError as err:
1924 raise ValueError(f"Cube must have a {sequence_coordinate} coordinate.") from err
1926 # Create a plot for each value of the sequence coordinate.
1927 plot_index = []
1928 nplot = np.size(cube.coord(sequence_coordinate).points)
1930 for iseq, cube_slice in enumerate(cube.slices_over(sequence_coordinate)):
1931 # Set plot titles and filename
1932 seq_coord = cube_slice.coord(sequence_coordinate)
1933 plot_title, plot_filename = _set_title_and_filename(
1934 seq_coord, nplot, recipe_title, filename
1935 )
1937 # Extract sequence slice for overlay cubes if required.
1938 overlay_slice = slice_over_maybe(overlay_cube, sequence_coordinate, iseq)
1939 contour_slice = slice_over_maybe(contour_cube, sequence_coordinate, iseq)
1940 scatter_slice = slice_over_maybe(scatter_cube, sequence_coordinate, iseq)
1942 # Do the actual plotting.
1943 plotting_func(
1944 cube_slice,
1945 filename=plot_filename,
1946 stamp_coordinate=stamp_coordinate,
1947 title=plot_title,
1948 method=method,
1949 overlay_cube=overlay_slice,
1950 contour_cube=contour_slice,
1951 scatter_cube=scatter_slice,
1952 **kwargs,
1953 )
1954 plot_index.append(plot_filename)
1956 # Add list of plots to plot metadata.
1957 complete_plot_index = _append_to_plot_index(plot_index)
1959 # Make a page to display the plots.
1960 _make_plot_html_page(complete_plot_index)
1963####################
1964# Public functions #
1965####################
1968def spatial_contour_plot(
1969 cube: iris.cube.Cube,
1970 filename: str = None,
1971 sequence_coordinate: str = "time",
1972 stamp_coordinate: str = "realization",
1973 **kwargs,
1974) -> iris.cube.Cube:
1975 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
1977 A 2D spatial field can be plotted, but if the sequence_coordinate is present
1978 then a sequence of plots will be produced. Similarly if the stamp_coordinate
1979 is present then postage stamp plots will be produced.
1981 Parameters
1982 ----------
1983 cube: Cube
1984 Iris cube of the data to plot. It should have two spatial dimensions,
1985 such as lat and lon, and may also have a another two dimension to be
1986 plotted sequentially and/or as postage stamp plots.
1987 filename: str, optional
1988 Name of the plot to write, used as a prefix for plot sequences. Defaults
1989 to the recipe name.
1990 sequence_coordinate: str, optional
1991 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
1992 This coordinate must exist in the cube.
1993 stamp_coordinate: str, optional
1994 Coordinate about which to plot postage stamp plots. Defaults to
1995 ``"realization"``.
1997 Returns
1998 -------
1999 Cube
2000 The original cube (so further operations can be applied).
2002 Raises
2003 ------
2004 ValueError
2005 If the cube doesn't have the right dimensions.
2006 TypeError
2007 If the cube isn't a single cube.
2008 """
2009 _spatial_plot(
2010 "contourf", cube, filename, sequence_coordinate, stamp_coordinate, **kwargs
2011 )
2012 return cube
2015def spatial_pcolormesh_plot(
2016 cube: iris.cube.Cube,
2017 filename: str = None,
2018 sequence_coordinate: str = "time",
2019 stamp_coordinate: str = "realization",
2020 **kwargs,
2021) -> iris.cube.Cube:
2022 """Plot a spatial variable onto a map from a 2D, 3D, or 4D cube.
2024 A 2D spatial field can be plotted, but if the sequence_coordinate is present
2025 then a sequence of plots will be produced. Similarly if the stamp_coordinate
2026 is present then postage stamp plots will be produced.
2028 This function is significantly faster than ``spatial_contour_plot``,
2029 especially at high resolutions, and should be preferred unless contiguous
2030 contour areas are important.
2032 Parameters
2033 ----------
2034 cube: Cube
2035 Iris cube of the data to plot. It should have two spatial dimensions,
2036 such as lat and lon, and may also have a another two dimension to be
2037 plotted sequentially and/or as postage stamp plots.
2038 filename: str, optional
2039 Name of the plot to write, used as a prefix for plot sequences. Defaults
2040 to the recipe name.
2041 sequence_coordinate: str, optional
2042 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2043 This coordinate must exist in the cube.
2044 stamp_coordinate: str, optional
2045 Coordinate about which to plot postage stamp plots. Defaults to
2046 ``"realization"``.
2048 Returns
2049 -------
2050 Cube
2051 The original cube (so further operations can be applied).
2053 Raises
2054 ------
2055 ValueError
2056 If the cube doesn't have the right dimensions.
2057 TypeError
2058 If the cube isn't a single cube.
2059 """
2060 _spatial_plot(
2061 "pcolormesh", cube, filename, sequence_coordinate, stamp_coordinate, **kwargs
2062 )
2063 return cube
2066def spatial_multi_pcolormesh_plot(
2067 cube: iris.cube.Cube,
2068 overlay_cube: iris.cube.Cube | None = None,
2069 contour_cube: iris.cube.Cube | None = None,
2070 scatter_cube: iris.cube.Cube | None = None,
2071 filename: str = None,
2072 sequence_coordinate: str = "time",
2073 stamp_coordinate: str = "realization",
2074 **kwargs,
2075) -> iris.cube.Cube:
2076 """Plot a set of spatial variables onto a map from a 2D, 3D, or 4D cube.
2078 A 2D basis cube spatial field can be plotted, but if the sequence_coordinate is present
2079 then a sequence of plots will be produced. Similarly if the stamp_coordinate
2080 is present then postage stamp plots will be produced.
2082 If specified, a masked overlay_cube can be overplotted on top of the base cube.
2084 If specified, contours of a contour_cube can be overplotted on top of those.
2086 If specified, a scattermap of scatter_cube can be overplotted.
2088 For single-variable equivalent of this routine, use spatial_pcolormesh_plot.
2090 This function is significantly faster than ``spatial_contour_plot``,
2091 especially at high resolutions, and should be preferred unless contiguous
2092 contour areas are important.
2094 Parameters
2095 ----------
2096 cube: Cube
2097 Iris cube of the data to plot. It should have two spatial dimensions,
2098 such as lat and lon, and may also have a another two dimension to be
2099 plotted sequentially and/or as postage stamp plots.
2100 overlay_cube: Cube
2101 Iris cube of the data to plot as an overlay on top of basis cube. It should have two spatial dimensions,
2102 such as lat and lon, and may also have a another two dimension to be
2103 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.
2104 contour_cube: Cube
2105 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,
2106 such as lat and lon, and may also have a another two dimension to be
2107 plotted sequentially and/or as postage stamp plots.
2108 scatter_cube: Cube
2109 Iris cube of the data to plot as a scatter map overlay on top of basis cube, overlay_cube and contour_cube. It should have two spatial dimensions,
2110 such as lat and lon, but these can describe a 1-D cube (e.g. list of
2111 observation stations with lat/lon coordinates) and may also have a another
2112 two dimension to be plotted sequentially and/or as postage stamp plots.
2113 filename: str, optional
2114 Name of the plot to write, used as a prefix for plot sequences. Defaults
2115 to the recipe name.
2116 sequence_coordinate: str, optional
2117 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2118 This coordinate must exist in the cube.
2119 stamp_coordinate: str, optional
2120 Coordinate about which to plot postage stamp plots. Defaults to
2121 ``"realization"``.
2123 Returns
2124 -------
2125 Cube
2126 The original cube (so further operations can be applied).
2128 Raises
2129 ------
2130 ValueError
2131 If the cube doesn't have the right dimensions.
2132 TypeError
2133 If the cube isn't a single cube.
2134 """
2135 _spatial_plot(
2136 "pcolormesh",
2137 cube,
2138 filename,
2139 sequence_coordinate,
2140 stamp_coordinate,
2141 overlay_cube=overlay_cube,
2142 contour_cube=contour_cube,
2143 scatter_cube=scatter_cube,
2144 )
2145 return cube, overlay_cube, contour_cube, scatter_cube
2148# TODO: Expand function to handle ensemble data.
2149# line_coordinate: str, optional
2150# Coordinate about which to plot multiple lines. Defaults to
2151# ``"realization"``.
2152def plot_line_series(
2153 cube: iris.cube.Cube | iris.cube.CubeList,
2154 filename: str = None,
2155 series_coordinate: str = "time",
2156 # line_coordinate: str = "realization",
2157 **kwargs,
2158) -> iris.cube.Cube | iris.cube.CubeList:
2159 """Plot a line plot for the specified coordinate.
2161 The Cube or CubeList must be 1D.
2163 Parameters
2164 ----------
2165 iris.cube | iris.cube.CubeList
2166 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
2167 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2168 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2169 filename: str, optional
2170 Name of the plot to write, used as a prefix for plot sequences. Defaults
2171 to the recipe name.
2172 series_coordinate: str, optional
2173 Coordinate about which to make a series. Defaults to ``"time"``. This
2174 coordinate must exist in the cube.
2176 Returns
2177 -------
2178 iris.cube.Cube | iris.cube.CubeList
2179 The original Cube or CubeList (so further operations can be applied).
2180 plotted data.
2182 Raises
2183 ------
2184 ValueError
2185 If the cubes don't have the right dimensions.
2186 TypeError
2187 If the cube isn't a Cube or CubeList.
2188 """
2189 # Ensure we have a name for the plot file.
2190 recipe_title = get_recipe_metadata().get("title", "Untitled")
2192 num_models = get_num_models(cube)
2194 validate_cube_shape(cube, num_models)
2196 # Iterate over all cubes and extract coordinate to plot.
2197 cubes = iris.cube.CubeList(iter_maybe(cube))
2198 coords = []
2199 for cube in cubes:
2200 try:
2201 coords.append(cube.coord(series_coordinate))
2202 except iris.exceptions.CoordinateNotFoundError as err:
2203 raise ValueError(
2204 f"Cube must have a {series_coordinate} coordinate."
2205 ) from err
2206 if cube.ndim > 2 or not cube.coords("realization"):
2207 raise ValueError("Cube must be 1D or 2D with a realization coordinate.")
2209 # Format the title and filename using plotted series coordinate
2210 plot_index = []
2211 nplot = 1
2212 seq_coord = coords[0]
2213 plot_title, plot_filename = _set_title_and_filename(
2214 seq_coord, nplot, recipe_title, filename
2215 )
2217 # Do the actual plotting.
2219 # Treat cubes with station coordinate as point observation timeseries, looping over available points
2220 if (
2221 "station" in [c.name() for c in cubes[0].coords()]
2222 and len(cubes[0].coord("station").points) > 1
2223 ):
2224 for station in cubes[0].coord("station").points:
2225 station_cubes = cubes.extract(iris.Constraint(station=station))
2226 station_name = station_cubes[0].coord("Station_Name").points[0]
2227 station_plotname = plot_filename.replace(
2228 ".png", "_" + station_name + ".png"
2229 )
2230 _plot_and_save_line_series(
2231 station_cubes,
2232 coords,
2233 "realization",
2234 station_plotname,
2235 f"{plot_title} {station_name}",
2236 )
2237 plot_index.append(station_plotname)
2239 # Default line plotting
2240 else:
2241 _plot_and_save_line_series(
2242 cubes, coords, "realization", plot_filename, plot_title
2243 )
2244 plot_index.append(plot_filename)
2246 # Add list of plots to plot metadata.
2247 complete_plot_index = _append_to_plot_index(plot_index)
2249 # Make a page to display the plots.
2250 _make_plot_html_page(complete_plot_index)
2252 return cube
2255def plot_vertical_line_series(
2256 cubes: iris.cube.Cube | iris.cube.CubeList,
2257 filename: str = None,
2258 series_coordinate: str = "model_level_number",
2259 sequence_coordinate: str = "time",
2260 # line_coordinate: str = "realization",
2261 **kwargs,
2262) -> iris.cube.Cube | iris.cube.CubeList:
2263 """Plot a line plot against a type of vertical coordinate.
2265 The Cube or CubeList must be 1D.
2267 A 1D line plot with y-axis as pressure coordinate can be plotted, but if the sequence_coordinate is present
2268 then a sequence of plots will be produced.
2270 Parameters
2271 ----------
2272 iris.cube | iris.cube.CubeList
2273 Cube or CubeList of the data to plot. The individual cubes should have a single dimension.
2274 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2275 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2276 filename: str, optional
2277 Name of the plot to write, used as a prefix for plot sequences. Defaults
2278 to the recipe name.
2279 series_coordinate: str, optional
2280 Coordinate to plot on the y-axis. Can be ``pressure`` or
2281 ``model_level_number`` for UM, or ``full_levels`` or ``half_levels``
2282 for LFRic. Defaults to ``model_level_number``.
2283 This coordinate must exist in the cube.
2284 sequence_coordinate: str, optional
2285 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2286 This coordinate must exist in the cube.
2288 Returns
2289 -------
2290 iris.cube.Cube | iris.cube.CubeList
2291 The original Cube or CubeList (so further operations can be applied).
2292 Plotted data.
2294 Raises
2295 ------
2296 ValueError
2297 If the cubes doesn't have the right dimensions.
2298 TypeError
2299 If the cube isn't a Cube or CubeList.
2300 """
2301 # Ensure we have a name for the plot file.
2302 recipe_title = get_recipe_metadata().get("title", "Untitled")
2304 cubes = iter_maybe(cubes)
2305 # Initialise empty list to hold all data from all cubes in a CubeList
2306 all_data = []
2308 # Store min/max ranges for x range.
2309 x_levels = []
2311 num_models = get_num_models(cubes)
2313 validate_cube_shape(cubes, num_models)
2315 # Iterate over all cubes in cube or CubeList and plot.
2316 coords = []
2317 for cube in cubes:
2318 # Test if series coordinate i.e. pressure level exist for any cube with cube.ndim >=1.
2319 try:
2320 coords.append(cube.coord(series_coordinate))
2321 except iris.exceptions.CoordinateNotFoundError as err:
2322 raise ValueError(
2323 f"Cube must have a {series_coordinate} coordinate."
2324 ) from err
2326 try:
2327 if cube.ndim > 1 or not cube.coords("realization"): 2327 ↛ 2335line 2327 didn't jump to line 2335 because the condition on line 2327 was always true
2328 cube.coord(sequence_coordinate)
2329 except iris.exceptions.CoordinateNotFoundError as err:
2330 raise ValueError(
2331 f"Cube must have a {sequence_coordinate} coordinate or be 1D, or 2D with a realization coordinate."
2332 ) from err
2334 # Get minimum and maximum from levels information.
2335 _, levels, _ = colorbar_map_levels(cube, axis="x")
2336 if levels is not None: 2336 ↛ 2340line 2336 didn't jump to line 2340 because the condition on line 2336 was always true
2337 x_levels.append(min(levels))
2338 x_levels.append(max(levels))
2339 else:
2340 all_data.append(cube.data)
2342 if len(x_levels) == 0: 2342 ↛ 2344line 2342 didn't jump to line 2344 because the condition on line 2342 was never true
2343 # Combine all data into a single NumPy array
2344 combined_data = np.concatenate(all_data)
2346 # Set the lower and upper limit for the x-axis to ensure all plots have
2347 # same range. This needs to read the whole cube over the range of the
2348 # sequence and if applicable postage stamp coordinate.
2349 vmin = np.floor(combined_data.min())
2350 vmax = np.ceil(combined_data.max())
2351 else:
2352 vmin = min(x_levels)
2353 vmax = max(x_levels)
2355 # Matching the slices (matching by seq coord point; it may happen that
2356 # evaluated models do not cover the same seq coord range, hence matching
2357 # necessary)
2358 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2360 # Create a plot for each value of the sequence coordinate.
2361 # Allowing for multiple cubes in a CubeList to be plotted in the same plot for
2362 # similar sequence values. Passing a CubeList into the internal plotting function
2363 # for similar values of the sequence coordinate. cube_slice can be an iris.cube.Cube
2364 # or an iris.cube.CubeList.
2365 plot_index = []
2366 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2367 for cubes_slice in cube_iterables:
2368 # Format the coordinate value in a unit appropriate way.
2369 seq_coord = cubes_slice[0].coord(sequence_coordinate)
2370 plot_title, plot_filename = _set_title_and_filename(
2371 seq_coord, nplot, recipe_title, filename
2372 )
2374 # Do the actual plotting.
2375 _plot_and_save_vertical_line_series(
2376 cubes_slice,
2377 coords,
2378 "realization",
2379 plot_filename,
2380 series_coordinate,
2381 title=plot_title,
2382 vmin=vmin,
2383 vmax=vmax,
2384 )
2385 plot_index.append(plot_filename)
2387 # Add list of plots to plot metadata.
2388 complete_plot_index = _append_to_plot_index(plot_index)
2390 # Make a page to display the plots.
2391 _make_plot_html_page(complete_plot_index)
2393 return cubes
2396def qq_plot(
2397 cubes: iris.cube.CubeList,
2398 coordinates: list[str],
2399 percentiles: list[float],
2400 model_names: list[str],
2401 filename: str = None,
2402 one_to_one: bool = True,
2403 **kwargs,
2404) -> iris.cube.CubeList:
2405 """Plot a Quantile-Quantile plot between two models for common time points.
2407 The cubes will be normalised by collapsing each cube to its percentiles. Cubes are
2408 collapsed within the operator over all specified coordinates such as
2409 grid_latitude, grid_longitude, vertical levels, but also realisation representing
2410 ensemble members to ensure a 1D cube (array).
2412 Parameters
2413 ----------
2414 cubes: iris.cube.CubeList
2415 Two cubes of the same variable with different models.
2416 coordinate: list[str]
2417 The list of coordinates to collapse over. This list should be
2418 every coordinate within the cube to result in a 1D cube around
2419 the percentile coordinate.
2420 percent: list[float]
2421 A list of percentiles to appear in the plot.
2422 model_names: list[str]
2423 A list of model names to appear on the axis of the plot.
2424 filename: str, optional
2425 Filename of the plot to write.
2426 one_to_one: bool, optional
2427 If True a 1:1 line is plotted; if False it is not. Default is True.
2429 Raises
2430 ------
2431 ValueError
2432 When the cubes are not compatible.
2434 Notes
2435 -----
2436 The quantile-quantile plot is a variant on the scatter plot representing
2437 two datasets by their quantiles (percentiles) for common time points.
2438 This plot does not use a theoretical distribution to compare against, but
2439 compares percentiles of two datasets. This plot does
2440 not use all raw data points, but plots the selected percentiles (quantiles) of
2441 each variable instead for the two datasets, thereby normalising the data for a
2442 direct comparison between the selected percentiles of the two dataset distributions.
2444 Quantile-quantile plots are valuable for comparing against
2445 observations and other models. Identical percentiles between the variables
2446 will lie on the one-to-one line implying the values correspond well to each
2447 other. Where there is a deviation from the one-to-one line a range of
2448 possibilities exist depending on how and where the data is shifted (e.g.,
2449 Wilks 2011 [Wilks2011]_).
2451 For distributions above the one-to-one line the distribution is left-skewed;
2452 below is right-skewed. A distinct break implies a bimodal distribution, and
2453 closer values/values further apart at the tails imply poor representation of
2454 the extremes.
2456 References
2457 ----------
2458 .. [Wilks2011] Wilks, D.S., (2011) "Statistical Methods in the Atmospheric
2459 Sciences" Third Edition, vol. 100, Academic Press, Oxford, UK, 676 pp.
2460 """
2461 # Check cubes using same functionality as the difference operator.
2462 if len(cubes) != 2:
2463 raise ValueError("cubes should contain exactly 2 cubes.")
2464 base: Cube = cubes.extract_cube(iris.AttributeConstraint(cset_comparison_base=1))
2465 other: Cube = cubes.extract_cube(
2466 iris.Constraint(
2467 cube_func=lambda cube: "cset_comparison_base" not in cube.attributes
2468 )
2469 )
2471 # Get spatial coord names.
2472 base_lat_name, base_lon_name = get_cube_yxcoordname(base)
2473 other_lat_name, other_lon_name = get_cube_yxcoordname(other)
2475 # Ensure cubes to compare are on common differencing grid.
2476 # This is triggered if either
2477 # i) latitude and longitude shapes are not the same. Note grid points
2478 # are not compared directly as these can differ through rounding
2479 # errors.
2480 # ii) or variables are known to often sit on different grid staggering
2481 # in different models (e.g. cell center vs cell edge), as is the case
2482 # for UM and LFRic comparisons.
2483 # In future greater choice of regridding method might be applied depending
2484 # on variable type. Linear regridding can in general be appropriate for smooth
2485 # variables. Care should be taken with interpretation of differences
2486 # given this dependency on regridding.
2487 if (
2488 base.coord(base_lat_name).shape != other.coord(other_lat_name).shape
2489 or base.coord(base_lon_name).shape != other.coord(other_lon_name).shape
2490 ) or (
2491 base.long_name
2492 in [
2493 "eastward_wind_at_10m",
2494 "northward_wind_at_10m",
2495 "northward_wind_at_cell_centres",
2496 "eastward_wind_at_cell_centres",
2497 "zonal_wind_at_pressure_levels",
2498 "meridional_wind_at_pressure_levels",
2499 "potential_vorticity_at_pressure_levels",
2500 "vapour_specific_humidity_at_pressure_levels_for_climate_averaging",
2501 ]
2502 ):
2503 logging.debug(
2504 "Linear regridding base cube to other grid to compute differences"
2505 )
2506 base = regrid_onto_cube(base, other, method="Linear")
2508 # Extract just common time points.
2509 base, other = _extract_common_time_points(base, other)
2511 # Equalise attributes so we can merge.
2512 fully_equalise_attributes([base, other])
2513 logging.debug("Base: %s\nOther: %s", base, other)
2515 # Collapse cubes.
2516 base = collapse(
2517 base,
2518 coordinate=coordinates,
2519 method="PERCENTILE",
2520 additional_percent=percentiles,
2521 )
2522 other = collapse(
2523 other,
2524 coordinate=coordinates,
2525 method="PERCENTILE",
2526 additional_percent=percentiles,
2527 )
2529 # Ensure we have a name for the plot file.
2530 recipe_title = get_recipe_metadata().get("title", "Untitled")
2531 title = f"{recipe_title}"
2533 if filename is None:
2534 filename = slugify(recipe_title)
2536 # Add file extension.
2537 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2539 # Do the actual plotting on a scatter plot
2540 _plot_and_save_scatter_plot(
2541 base, other, plot_filename, title, one_to_one, model_names
2542 )
2544 # Add list of plots to plot metadata.
2545 plot_index = _append_to_plot_index([plot_filename])
2547 # Make a page to display the plots.
2548 _make_plot_html_page(plot_index)
2550 return iris.cube.CubeList([base, other])
2553def scatter_plot(
2554 cube_x: iris.cube.Cube | iris.cube.CubeList,
2555 cube_y: iris.cube.Cube | iris.cube.CubeList,
2556 filename: str = None,
2557 one_to_one: bool = True,
2558 **kwargs,
2559) -> iris.cube.CubeList:
2560 """Plot a scatter plot between two variables.
2562 Both cubes must be 1D.
2564 Parameters
2565 ----------
2566 cube_x: Cube | CubeList
2567 1 dimensional Cube of the data to plot on y-axis.
2568 cube_y: Cube | CubeList
2569 1 dimensional Cube of the data to plot on x-axis.
2570 filename: str, optional
2571 Filename of the plot to write.
2572 one_to_one: bool, optional
2573 If True a 1:1 line is plotted; if False it is not. Default is True.
2575 Returns
2576 -------
2577 cubes: CubeList
2578 CubeList of the original x and y cubes for further processing.
2580 Raises
2581 ------
2582 ValueError
2583 If the cube doesn't have the right dimensions and cubes not the same
2584 size.
2585 TypeError
2586 If the cube isn't a single cube.
2588 Notes
2589 -----
2590 Scatter plots are used for determining if there is a relationship between
2591 two variables. Positive relations have a slope going from bottom left to top
2592 right; Negative relations have a slope going from top left to bottom right.
2593 """
2594 # Iterate over all cubes in cube or CubeList and plot.
2595 for cube_iter in iter_maybe(cube_x):
2596 # Check cubes are correct shape.
2597 cube_iter = check_single_cube(cube_iter)
2598 if cube_iter.ndim > 1:
2599 raise ValueError("cube_x must be 1D.")
2601 # Iterate over all cubes in cube or CubeList and plot.
2602 for cube_iter in iter_maybe(cube_y):
2603 # Check cubes are correct shape.
2604 cube_iter = check_single_cube(cube_iter)
2605 if cube_iter.ndim > 1:
2606 raise ValueError("cube_y must be 1D.")
2608 # Ensure we have a name for the plot file.
2609 recipe_title = get_recipe_metadata().get("title", "Untitled")
2610 title = f"{recipe_title}"
2612 if filename is None:
2613 filename = slugify(recipe_title)
2615 # Add file extension.
2616 plot_filename = f"{filename.rsplit('.', 1)[0]}.png"
2618 # Do the actual plotting.
2619 _plot_and_save_scatter_plot(cube_x, cube_y, plot_filename, title, one_to_one)
2621 # Add list of plots to plot metadata.
2622 plot_index = _append_to_plot_index([plot_filename])
2624 # Make a page to display the plots.
2625 _make_plot_html_page(plot_index)
2627 return iris.cube.CubeList([cube_x, cube_y])
2630def vector_plot(
2631 cube_u: iris.cube.Cube,
2632 cube_v: iris.cube.Cube,
2633 filename: str = None,
2634 sequence_coordinate: str = "time",
2635 **kwargs,
2636) -> iris.cube.CubeList:
2637 """Plot a vector plot based on the input u and v components."""
2638 recipe_title = get_recipe_metadata().get("title", "Untitled")
2640 # Cubes must have a matching sequence coordinate.
2641 try:
2642 # Check that the u and v cubes have the same sequence coordinate.
2643 if cube_u.coord(sequence_coordinate) != cube_v.coord(sequence_coordinate): 2643 ↛ anywhereline 2643 didn't jump anywhere: it always raised an exception.
2644 raise ValueError("Coordinates do not match.")
2645 except (iris.exceptions.CoordinateNotFoundError, ValueError) as err:
2646 raise ValueError(
2647 f"Cubes should have matching {sequence_coordinate} coordinate:\n{cube_u}\n{cube_v}"
2648 ) from err
2650 # Create a plot for each value of the sequence coordinate.
2651 plot_index = []
2652 nplot = np.size(cube_u[0].coord(sequence_coordinate).points)
2653 for cube_u_slice, cube_v_slice in zip(
2654 cube_u.slices_over(sequence_coordinate),
2655 cube_v.slices_over(sequence_coordinate),
2656 strict=True,
2657 ):
2658 # Format the coordinate value in a unit appropriate way.
2659 seq_coord = cube_u_slice.coord(sequence_coordinate)
2660 plot_title, plot_filename = _set_title_and_filename(
2661 seq_coord, nplot, recipe_title, filename
2662 )
2664 # Do the actual plotting.
2665 _plot_and_save_vector_plot(
2666 cube_u_slice,
2667 cube_v_slice,
2668 filename=plot_filename,
2669 title=plot_title,
2670 method="contourf",
2671 )
2672 plot_index.append(plot_filename)
2674 # Add list of plots to plot metadata.
2675 complete_plot_index = _append_to_plot_index(plot_index)
2677 # Make a page to display the plots.
2678 _make_plot_html_page(complete_plot_index)
2680 return iris.cube.CubeList([cube_u, cube_v])
2683def plot_histogram_series(
2684 cubes: iris.cube.Cube | iris.cube.CubeList,
2685 filename: str = None,
2686 sequence_coordinate: str = "time",
2687 stamp_coordinate: str = "realization",
2688 single_plot: bool = False,
2689 **kwargs,
2690) -> iris.cube.Cube | iris.cube.CubeList:
2691 """Plot a histogram plot for each vertical level provided.
2693 A histogram plot can be plotted, but if the sequence_coordinate (i.e. time)
2694 is present then a sequence of plots will be produced using the time slider
2695 functionality to scroll through histograms against time. If a
2696 stamp_coordinate is present then postage stamp plots will be produced. If
2697 stamp_coordinate and single_plot is True, all postage stamp plots will be
2698 plotted in a single plot instead of separate postage stamp plots.
2700 Parameters
2701 ----------
2702 cubes: Cube | iris.cube.CubeList
2703 Iris cube or CubeList of the data to plot. It should have a single dimension other
2704 than the stamp coordinate.
2705 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2706 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2707 filename: str, optional
2708 Name of the plot to write, used as a prefix for plot sequences. Defaults
2709 to the recipe name.
2710 sequence_coordinate: str, optional
2711 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2712 This coordinate must exist in the cube and will be used for the time
2713 slider.
2714 stamp_coordinate: str, optional
2715 Coordinate about which to plot postage stamp plots. Defaults to
2716 ``"realization"``.
2717 single_plot: bool, optional
2718 If True, all postage stamp plots will be plotted in a single plot. If
2719 False, each postage stamp plot will be plotted separately. Is only valid
2720 if stamp_coordinate exists and has more than a single point.
2722 Returns
2723 -------
2724 iris.cube.Cube | iris.cube.CubeList
2725 The original Cube or CubeList (so further operations can be applied).
2726 Plotted data.
2728 Raises
2729 ------
2730 ValueError
2731 If the cube doesn't have the right dimensions.
2732 TypeError
2733 If the cube isn't a Cube or CubeList.
2734 """
2735 recipe_title = get_recipe_metadata().get("title", "Untitled")
2737 cubes = iter_maybe(cubes)
2739 # Internal plotting function.
2740 plotting_func = _plot_and_save_histogram_series
2742 num_models = get_num_models(cubes)
2744 validate_cube_shape(cubes, num_models)
2746 # If several histograms are plotted, check sequence_coordinate
2747 check_sequence_coordinate(cubes, sequence_coordinate)
2749 # Get axis minimum and maximum from levels information.
2750 # If no levels set, derive minima and maxima from data in CubeList.
2751 vmin, vmax = _set_axis_range(cubes)
2753 # Make postage stamp plots if stamp_coordinate exists and has more than a
2754 # single point. If single_plot is True:
2755 # -- all postage stamp plots will be plotted in a single plot instead of
2756 # separate postage stamp plots.
2757 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
2758 # produced per single model only
2759 if num_models == 1:
2760 if ( 2760 ↛ 2764line 2760 didn't jump to line 2764 because the condition on line 2760 was never true
2761 stamp_coordinate in [c.name() for c in cubes[0].coords()]
2762 and cubes[0].coord(stamp_coordinate).shape[0] > 1
2763 ):
2764 if single_plot:
2765 plotting_func = (
2766 _plot_and_save_postage_stamps_in_single_plot_histogram_series
2767 )
2768 else:
2769 plotting_func = _plot_and_save_postage_stamp_histogram_series
2770 cube_iterables = cubes[0].slices_over(sequence_coordinate)
2771 else:
2772 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2774 plot_index = []
2775 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2776 # Create a plot for each value of the sequence coordinate. Allowing for
2777 # multiple cubes in a CubeList to be plotted in the same plot for similar
2778 # sequence values. Passing a CubeList into the internal plotting function
2779 # for similar values of the sequence coordinate. cube_slice can be an
2780 # iris.cube.Cube or an iris.cube.CubeList.
2781 for cube_slice in cube_iterables:
2782 single_cube = cube_slice
2783 if isinstance(cube_slice, iris.cube.CubeList):
2784 single_cube = cube_slice[0]
2786 # Ensure valid stamp coordinate in cube dimensions
2787 if stamp_coordinate == "realization": 2787 ↛ 2790line 2787 didn't jump to line 2790 because the condition on line 2787 was always true
2788 stamp_coordinate = check_stamp_coordinate(single_cube)
2789 # Set plot titles and filename, based on sequence coordinate
2790 seq_coord = single_cube.coord(sequence_coordinate)
2791 # Use time coordinate in title and filename if single histogram output.
2792 if sequence_coordinate == "realization" and nplot == 1: 2792 ↛ 2793line 2792 didn't jump to line 2793 because the condition on line 2792 was never true
2793 seq_coord = single_cube.coord("time")
2794 # Use station name in title and filename if model vs obs comparison
2795 if sequence_coordinate == "station": 2795 ↛ 2796line 2795 didn't jump to line 2796 because the condition on line 2795 was never true
2796 seq_coord = single_cube.coord("Station_Name")
2798 plot_title, plot_filename = _set_title_and_filename(
2799 seq_coord, nplot, recipe_title, filename
2800 )
2802 # Do the actual plotting.
2803 plotting_func(
2804 cube_slice,
2805 filename=plot_filename,
2806 stamp_coordinate=stamp_coordinate,
2807 title=plot_title,
2808 vmin=vmin,
2809 vmax=vmax,
2810 )
2811 plot_index.append(plot_filename)
2813 # Add list of plots to plot metadata.
2814 complete_plot_index = _append_to_plot_index(plot_index)
2816 # Make a page to display the plots.
2817 _make_plot_html_page(complete_plot_index)
2819 return cubes
2822def plot_scatter_series(
2823 cubes: iris.cube.Cube | iris.cube.CubeList,
2824 filename: str = None,
2825 sequence_coordinate: str = "time",
2826 stamp_coordinate: str = "realization",
2827 hexbin: bool = False,
2828 **kwargs,
2829) -> iris.cube.Cube | iris.cube.CubeList:
2830 """Plot a scatter plot for each sequence coordinate provided.
2832 A scatter plot can be plotted, but if the sequence_coordinate (i.e. time)
2833 is present then a sequence of plots will be produced using the time slider
2834 functionality to scroll through scatter against time. If a
2835 stamp_coordinate is present then postage stamp plots will be produced. If
2836 stamp_coordinate and single_plot is True, all postage stamp plots will be
2837 plotted in a single plot instead of separate postage stamp plots.
2839 Parameters
2840 ----------
2841 cubes: Cube | iris.cube.CubeList
2842 Iris cube or CubeList of the data to plot. It should have a single dimension other
2843 than the stamp coordinate.
2844 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2845 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2846 filename: str, optional
2847 Name of the plot to write, used as a prefix for plot sequences. Defaults
2848 to the recipe name.
2849 sequence_coordinate: str, optional
2850 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2851 This coordinate must exist in the cube and will be used for the time
2852 slider.
2853 stamp_coordinate: str, optional
2854 Coordinate about which to plot postage stamp plots. Defaults to
2855 ``"realization"``.
2856 hexbin: bool, optional
2857 If True, generate hexbin comparison plot.
2858 If False, generate point-by-point scatter plot.
2860 Returns
2861 -------
2862 iris.cube.Cube | iris.cube.CubeList
2863 The original Cube or CubeList (so further operations can be applied).
2864 Plotted data.
2866 Raises
2867 ------
2868 ValueError
2869 If the cube doesn't have the right dimensions.
2870 TypeError
2871 If the cube isn't a Cube or CubeList.
2872 """
2873 recipe_title = get_recipe_metadata().get("title", "Untitled")
2875 cubes = iter_maybe(cubes)
2877 # Internal plotting function.
2878 plotting_func = _plot_and_save_scatter_series
2880 num_models = get_num_models(cubes)
2882 validate_cube_shape(cubes, num_models)
2884 check_sequence_coordinate(cubes, sequence_coordinate)
2886 vmin, vmax = _set_axis_range(cubes)
2888 # Require >1 models to compare on scatter plot
2889 if num_models > 1:
2890 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
2891 else:
2892 raise ValueError(
2893 "Scatter plot series requires multiple number of models in input data."
2894 )
2896 plot_index = []
2897 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
2898 # Create a plot for each value of the sequence coordinate. Allowing for
2899 # multiple cubes in a CubeList to be plotted in the same plot for similar
2900 # sequence values. Passing a CubeList into the internal plotting function
2901 # for similar values of the sequence coordinate. cube_slice can be an
2902 # iris.cube.Cube or an iris.cube.CubeList.
2903 for cube_slice in cube_iterables:
2904 single_cube = cube_slice
2905 if isinstance(cube_slice, iris.cube.CubeList): 2905 ↛ 2909line 2905 didn't jump to line 2909 because the condition on line 2905 was always true
2906 single_cube = cube_slice[0]
2908 # Ensure valid stamp coordinate in cube dimensions
2909 if stamp_coordinate == "realization": 2909 ↛ 2912line 2909 didn't jump to line 2912 because the condition on line 2909 was always true
2910 stamp_coordinate = check_stamp_coordinate(single_cube)
2911 # Set plot titles and filename, based on sequence coordinate
2912 seq_coord = single_cube.coord(sequence_coordinate)
2913 # Use time coordinate in title and filename if single histogram output.
2914 if sequence_coordinate == "realization" and nplot == 1:
2915 seq_coord = single_cube.coord("time")
2916 # Use station name in title and filename if model vs obs comparison
2917 if sequence_coordinate == "station":
2918 seq_coord = single_cube.coord("Station_Name")
2920 plot_title, plot_filename = _set_title_and_filename(
2921 seq_coord, nplot, recipe_title, filename
2922 )
2924 # Do the actual plotting.
2925 plotting_func(
2926 cube_slice,
2927 filename=plot_filename,
2928 stamp_coordinate=stamp_coordinate,
2929 title=plot_title,
2930 vmin=vmin,
2931 vmax=vmax,
2932 hexbin=hexbin,
2933 )
2934 plot_index.append(plot_filename)
2936 # Add list of plots to plot metadata.
2937 complete_plot_index = _append_to_plot_index(plot_index)
2939 # Make a page to display the plots.
2940 _make_plot_html_page(complete_plot_index)
2942 return cubes
2945def plot_power_spectrum_series(
2946 cubes: iris.cube.Cube | iris.cube.CubeList,
2947 filename: str = None,
2948 sequence_coordinate: str = "time",
2949 stamp_coordinate: str = "realization",
2950 single_plot: bool = False,
2951 **kwargs,
2952) -> iris.cube.Cube | iris.cube.CubeList:
2953 """Plot a power spectrum plot for each vertical level provided.
2955 A power spectrum plot can be plotted, but if the sequence_coordinate (i.e. time)
2956 is present then a sequence of plots will be produced using the time slider
2957 functionality to scroll through power spectrum against time. If a
2958 stamp_coordinate is present then postage stamp plots will be produced. If
2959 stamp_coordinate and single_plot is True, all postage stamp plots will be
2960 plotted in a single plot instead of separate postage stamp plots.
2962 Parameters
2963 ----------
2964 cubes: Cube | iris.cube.CubeList
2965 Iris cube or CubeList of the data to plot. It should have a single dimension other
2966 than the stamp coordinate.
2967 The cubes should cover the same phenomenon i.e. all cubes contain temperature data.
2968 We do not support different data such as temperature and humidity in the same CubeList for plotting.
2969 filename: str, optional
2970 Name of the plot to write, used as a prefix for plot sequences. Defaults
2971 to the recipe name.
2972 sequence_coordinate: str, optional
2973 Coordinate about which to make a plot sequence. Defaults to ``"time"``.
2974 This coordinate must exist in the cube and will be used for the time
2975 slider.
2976 stamp_coordinate: str, optional
2977 Coordinate about which to plot postage stamp plots. Defaults to
2978 ``"realization"``.
2979 single_plot: bool, optional
2980 If True, all postage stamp plots will be plotted in a single plot. If
2981 False, each postage stamp plot will be plotted separately. Is only valid
2982 if stamp_coordinate exists and has more than a single point.
2984 Returns
2985 -------
2986 iris.cube.Cube | iris.cube.CubeList
2987 The original Cube or CubeList (so further operations can be applied).
2988 Plotted data.
2990 Raises
2991 ------
2992 ValueError
2993 If the cube doesn't have the right dimensions.
2994 TypeError
2995 If the cube isn't a Cube or CubeList.
2996 """
2997 recipe_title = get_recipe_metadata().get("title", "Untitled")
2999 cubes = iter_maybe(cubes)
3001 # Internal plotting function.
3002 plotting_func = _plot_and_save_power_spectrum_series
3004 num_models = get_num_models(cubes)
3006 validate_cube_shape(cubes, num_models)
3008 # If several power spectra are plotted, check sequence_coordinate
3009 check_sequence_coordinate(cubes, sequence_coordinate)
3011 # Make postage stamp plots if stamp_coordinate exists and has more than a
3012 # single point. If single_plot is True:
3013 # -- all postage stamp plots will be plotted in a single plot instead of
3014 # separate postage stamp plots.
3015 # -- model names (hidden in cube attrs) are ignored, that is stamp plots are
3016 # produced per single model only
3017 if num_models == 1: 3017 ↛ 3030line 3017 didn't jump to line 3030 because the condition on line 3017 was always true
3018 if ( 3018 ↛ 3022line 3018 didn't jump to line 3022 because the condition on line 3018 was never true
3019 stamp_coordinate in [c.name() for c in cubes[0].coords()]
3020 and cubes[0].coord(stamp_coordinate).shape[0] > 1
3021 ):
3022 if single_plot:
3023 plotting_func = (
3024 _plot_and_save_postage_stamps_in_single_plot_power_spectrum_series
3025 )
3026 else:
3027 plotting_func = _plot_and_save_postage_stamp_power_spectrum_series
3028 cube_iterables = cubes[0].slices_over(sequence_coordinate)
3029 else:
3030 cube_iterables = _find_matched_slices(cubes, sequence_coordinate)
3032 plot_index = []
3033 nplot = np.size(cubes[0].coord(sequence_coordinate).points)
3034 # Create a plot for each value of the sequence coordinate. Allowing for
3035 # multiple cubes in a CubeList to be plotted in the same plot for similar
3036 # sequence values. Passing a CubeList into the internal plotting function
3037 # for similar values of the sequence coordinate. cube_slice can be an
3038 # iris.cube.Cube or an iris.cube.CubeList.
3039 for cube_slice in cube_iterables:
3040 single_cube = cube_slice
3041 if isinstance(cube_slice, iris.cube.CubeList): 3041 ↛ 3042line 3041 didn't jump to line 3042 because the condition on line 3041 was never true
3042 single_cube = cube_slice[0]
3044 # Set stamp coordinate
3045 if stamp_coordinate == "realization": 3045 ↛ 3048line 3045 didn't jump to line 3048 because the condition on line 3045 was always true
3046 stamp_coordinate = check_stamp_coordinate(single_cube)
3047 # Set plot title and filenames based on sequence values
3048 seq_coord = single_cube.coord(sequence_coordinate)
3049 plot_title, plot_filename = _set_title_and_filename(
3050 seq_coord, nplot, recipe_title, filename
3051 )
3053 # Do the actual plotting.
3054 plotting_func(
3055 cube_slice,
3056 filename=plot_filename,
3057 stamp_coordinate=stamp_coordinate,
3058 title=plot_title,
3059 )
3060 plot_index.append(plot_filename)
3062 # Add list of plots to plot metadata.
3063 complete_plot_index = _append_to_plot_index(plot_index)
3065 # Make a page to display the plots.
3066 _make_plot_html_page(complete_plot_index)
3068 return cubes