"""Lightweight Matplotlib helpers that follow the active Matplotlib style."""
from __future__ import annotations
from types import MappingProxyType
from typing import Any, Mapping, cast
import matplotlib
import numpy as np
from matplotlib.axes import Axes
from matplotlib.collections import PolyCollection
from matplotlib.figure import Figure
from matplotlib.projections.polar import PolarAxes
from shapely import wkt as shapely_wkt
from shapely.geometry import LineString, MultiLineString, MultiPolygon, Polygon
from ._site.base import sector_angles
from ._site.coverage import Turbulence
_SINGLE_FIGSIZE = (6.0, 6.0)
__all__ = ["STYLE_RCPARAMS", "use_style"]
#: Matplotlib readability defaults applied by :func:`use_style`.
STYLE_RCPARAMS: Mapping[str, Any] = MappingProxyType(
{
"figure.dpi": 110.0,
"savefig.dpi": 110.0,
"axes.grid": True,
"grid.alpha": 0.3,
"font.size": 10.0,
"axes.titlesize": 12.0,
"axes.labelsize": 10.0,
"legend.fontsize": 9.0,
"xtick.labelsize": 9.0,
"ytick.labelsize": 9.0,
}
)
[docs]
def use_style() -> None:
"""Apply OOWL's documentation defaults to global Matplotlib settings.
Updates ``matplotlib.rcParams`` with :data:`STYLE_RCPARAMS`, so repeated
calls are idempotent. OOWL's own plots do not need it. It does not import
:mod:`matplotlib.pyplot`, select a backend or enable interactive mode.
"""
matplotlib.rcParams.update(STYLE_RCPARAMS)
# Display backends dispose of pyplot figures without enabling interactive mode.
_DISPLAY_BACKENDS = ("inline", "ipympl", "widget", "nbagg", "webagg")
def _shows_pyplot_figures() -> bool:
"""Whether the active backend draws the figures pyplot is holding."""
if matplotlib.is_interactive():
return True
try:
backend = matplotlib.get_backend().lower()
except Exception: # pragma: no cover - a backend that refuses to resolve
return False
return any(name in backend for name in _DISPLAY_BACKENDS)
def _new_figure(figsize: tuple[float, float], **kwargs) -> Figure:
"""Create a figure, registered with pyplot only if the backend displays it.
Otherwise return a detached :class:`~matplotlib.figure.Figure` so figures
do not leak in pyplot's registry.
"""
if _shows_pyplot_figures():
import matplotlib.pyplot as plt
return plt.figure(figsize=figsize, **kwargs)
return Figure(figsize=figsize, **kwargs)
def _new_axes(figsize: tuple[float, float] = _SINGLE_FIGSIZE, **subplot_kw) -> Axes:
"""Create a single-axes figure and return its axes."""
return _new_figure(figsize).subplots(subplot_kw=subplot_kw)
def _style_axes(ax: Axes, *, equal: bool = False, grid: bool = True) -> None:
"""Tidy one axes without changing global Matplotlib state."""
ax.set_axisbelow(True)
if grid:
ax.grid(True)
if equal:
ax.set_aspect("equal", adjustable="datalim")
for side in ("top", "right"):
ax.spines[side].set_visible(False)
def _cmap(name: str | None):
"""The named colormap, or the style's ``image.cmap`` when ``None``."""
if name is None:
name = matplotlib.rcParams["image.cmap"]
return matplotlib.colormaps[name]
def _domain_geometry(cluster):
"""Scatter domain recorded as WKT in ``metadata["placement"]``, or ``None``."""
placement = cluster.metadata.get("placement")
recorded = placement.get("domain") if isinstance(placement, dict) else None
if not isinstance(recorded, str):
return None
geometry = shapely_wkt.loads(recorded)
return geometry if isinstance(geometry, (Polygon, MultiPolygon)) else None
def _polygons(geometry):
"""The polygons of a Polygon or MultiPolygon, as a tuple."""
if isinstance(geometry, MultiPolygon):
return tuple(geometry.geoms)
if isinstance(geometry, Polygon):
return (geometry,)
raise TypeError(f"Unsupported geometry type: {type(geometry).__name__}")
def _plot_exclusion(ax: Axes, interior, edgecolor, linewidth, zorder: float) -> None:
"""Hatch one ring of excluded ground over whatever lies beneath it."""
ix, iy = interior.xy
# The axes facecolor hides the region fill, so the hatch reads as a hole.
ax.fill(ix, iy, facecolor=ax.get_facecolor(), edgecolor="none", zorder=zorder)
ax.fill(
ix,
iy,
fill=False,
hatch="///",
edgecolor=edgecolor,
linewidth=linewidth,
zorder=zorder + 0.1,
)
def _plot_domain(ax: Axes, geometry) -> None:
"""Outline a scatter domain, with its exclusions hatched."""
edgecolor = matplotlib.rcParams["axes.edgecolor"]
for poly in _polygons(geometry):
x, y = poly.exterior.xy
ax.plot(x, y, color=edgecolor, linewidth=1.0, linestyle="--", zorder=0)
for interior in poly.interiors:
_plot_exclusion(ax, interior, edgecolor, 1.0, zorder=0.1)
def _plot_polygon(ax: Axes, poly: Polygon, regions: bool, **kwargs) -> None:
x, y = poly.exterior.xy
ax.fill(x, y, **kwargs)
z = kwargs["zorder"] if "zorder" in kwargs else 1
edgecolor, linewidth = kwargs["edgecolor"], kwargs["linewidth"]
for interior in poly.interiors:
if regions:
_plot_exclusion(ax, interior, edgecolor, linewidth, zorder=z + 0.1)
continue
# Punch the hole out to the axes background, then edge it like the shape.
ix, iy = interior.xy
ax.fill(ix, iy, facecolor=ax.get_facecolor(), edgecolor="none", zorder=z + 0.1)
ax.plot(
ix,
iy,
color=edgecolor,
linewidth=linewidth,
alpha=kwargs["alpha"],
zorder=z + 0.2,
)
def _plot_boundary(
geom,
ax: Axes | None = None,
title: str | None = None,
facecolor: str | None = None,
edgecolor: str | None = None,
alpha: float | None = None,
linewidth: float | None = None,
regions: bool = True,
) -> Axes:
"""Plot a Polygon, MultiPolygon or ``Boundary`` on equal axes and return them.
``regions`` hatches holes as excluded ground, else they are punched out.
"""
if hasattr(geom, "parts"): # a Boundary instance
geometries = [part.geometry for part in geom.parts]
else:
geometries = [geom]
if ax is None:
ax = _new_axes()
# ``None`` selects the first color of the style's property cycle.
style = dict(
facecolor="C0" if facecolor is None else facecolor,
edgecolor="C0" if edgecolor is None else edgecolor,
alpha=0.25 if alpha is None else alpha,
linewidth=1.5 if linewidth is None else linewidth,
)
for geometry in geometries:
for poly in _polygons(geometry):
_plot_polygon(ax, poly, regions, **style)
_style_axes(ax, equal=True, grid=False)
ax.autoscale_view()
ax.set(xlabel="x [m]", ylabel="y [m]")
if title:
ax.set_title(title)
return ax
def _plot_farm(
layout,
ax: Axes | None = None,
title: str | None = None,
facecolor: str | None = None,
edgecolor: str | None = None,
alpha: float | None = None,
linewidth: float | None = None,
regions: bool = True,
**turbine_kwargs,
) -> Axes:
"""Plot a :class:`~oowl.layouts.base.Farm` boundary and turbines.
Extra keywords go to ``scatter``. Turbines take the first property-cycle
color, or one cycle color per turbine type, unless ``color`` is given.
"""
style: dict[str, Any] = {"s": 20.0, "zorder": 5, **turbine_kwargs}
if style.get("color") is None:
style.pop("color", None)
if layout.boundary is not None:
ax = _plot_boundary(
layout.boundary,
ax=ax,
title=title,
facecolor=facecolor,
edgecolor=edgecolor,
alpha=alpha,
linewidth=linewidth,
regions=regions,
)
if layout.guide is not None:
ax = _plot_string_guide(layout.guide, ax=ax, title=title)
if ax is None:
ax = _new_axes()
pts = np.asarray(layout.positions).reshape(-1, 2)
if pts.size:
turbines = layout.turbines
multiple_types = (
isinstance(turbines, tuple)
and len({id(turbine) for turbine in turbines}) > 1
)
if "color" in style or not multiple_types:
style.setdefault("color", "C0")
ax.scatter(pts[:, 0], pts[:, 1], **style)
else:
groups: dict[int, tuple[object, list[int]]] = {}
for index, turbine in enumerate(turbines):
groups.setdefault(id(turbine), (turbine, []))[1].append(index)
for color_index, (turbine, indices) in enumerate(groups.values()):
selected = pts[indices]
label = getattr(turbine, "name", None) or type(turbine).__name__
typed: dict[str, Any] = {**style, "color": f"C{color_index}"}
ax.scatter(selected[:, 0], selected[:, 1], label=label, **typed)
ax.legend(title="Turbine type", loc="upper left")
_style_axes(ax, equal=True, grid=False)
ax.autoscale_view()
ax.set(xlabel="x [m]", ylabel="y [m]")
if title and layout.boundary is None and layout.guide is None:
ax.set_title(title)
return ax
def _plot_cluster(
cluster,
ax: Axes | None = None,
title: str | None = None,
facecolor: str | None = None,
edgecolor: str | None = None,
alpha: float | None = 0.2,
linewidth: float | None = None,
regions: bool = True,
**turbine_kwargs,
) -> Axes:
"""Plot every member of a :class:`~oowl.Cluster` on shared axes.
Members take successive property-cycle colors unless ``color`` is given.
``regions`` also draws the recorded placement domain and its exclusions.
"""
if ax is None:
ax = _new_axes()
if turbine_kwargs.get("color") is None:
turbine_kwargs.pop("color", None)
if regions:
domain = _domain_geometry(cluster)
if domain is not None:
_plot_domain(ax, domain)
for index, farm in enumerate(cluster.farms):
color = f"C{index}"
_plot_farm(
farm,
ax=ax,
facecolor=color if facecolor is None else facecolor,
edgecolor=color if edgecolor is None else edgecolor,
alpha=alpha,
linewidth=linewidth,
regions=regions,
**{"color": color, **turbine_kwargs},
)
if title:
ax.set_title(title)
return ax
def _plot_wind_rose(
source,
wind_speed=None,
*,
ax: Axes | None = None,
n_sectors: int = 16,
speed_bins=None,
n_speed_bins: int = 6,
cmap: str | None = None,
opening: float = 0.8,
edgecolor="none",
legend: bool = True,
title: str | None = None,
figsize: tuple[float, float] = _SINGLE_FIGSIZE,
) -> Axes:
"""Plot a wind rose from a wind rose, resource, site conditions or arrays.
Bars are percentages, stacked by speed bin when speeds are available. A
resource's bars are percent of all time, so they sum to the share of time
its speed bins cover, and sampled cases' bars are percent of the cases.
The colormap defaults to the style's ``image.cmap``.
"""
theta_deg, table, bins = _wind_rose_table(
source, wind_speed, n_sectors, speed_bins, n_speed_bins
)
cmap_obj = _cmap(cmap)
if bins is None:
# Without speeds, one bar per sector shows its probability.
colors, labels, legend_title = [cmap_obj(0.5)], None, None
else:
colors = [cmap_obj(i) for i in np.linspace(0.0, 1.0, len(bins) - 1)]
labels = [f"{lo:.1f}\u2013{hi:.1f}" for lo, hi in zip(bins[:-1], bins[1:])]
legend_title = "Wind speed [m/s]" if legend else None
return _draw_rose(
ax,
figsize,
theta_deg,
table,
colors,
labels=labels,
opening=opening,
edgecolor=edgecolor,
legend_title=legend_title,
title=title,
)
def _draw_rose(
ax,
figsize,
theta_deg,
table,
colors,
*,
labels,
opening,
edgecolor,
legend_title,
title,
) -> Axes:
"""Draw stacked sector bars on a north-up, clockwise polar axes."""
if ax is None:
ax = _new_axes(figsize, projection="polar")
ax = cast(PolarAxes, ax)
ax.set_theta_zero_location("N")
ax.set_theta_direction(-1)
theta = np.deg2rad(theta_deg)
bar_width = np.deg2rad(360.0 / len(theta_deg)) * opening
_plot_stacked_polar_bars(
ax, theta, table, bar_width, colors, labels=labels, edgecolor=edgecolor
)
ax.set_thetagrids(np.arange(0, 360, 90), labels=["N", "E", "S", "W"])
if legend_title is not None:
ax.legend(
title=legend_title,
loc="upper left",
bbox_to_anchor=(1.02, 1.0),
frameon=False,
)
if title:
ax.set_title(title, pad=20)
return ax
def _plot_ti_rose(
source,
*,
ax: Axes | None = None,
ti_bins=None,
n_ti_bins: int = 8,
u_min: float = 4.0,
u_max: float = 25.0,
n_quadrature: int = 2000,
cmap: str | None = None,
opening: float = 0.8,
edgecolor="none",
legend: bool = True,
title: str | None = None,
figsize: tuple[float, float] = _SINGLE_FIGSIZE,
) -> Axes:
"""Plot an analytical direction x turbulence-intensity rose.
The colormap defaults to the style's ``image.cmap``.
"""
if not hasattr(source, "ti_rose"):
raise TypeError("TI plotting requires an analytical wind resource")
if ti_bins is None:
ti_bins = source._ti_bins(u_min, u_max, n_ti_bins)
bins = np.asarray(ti_bins, dtype=float)
theta_deg = np.asarray(source.theta, dtype=float)
table = (
source.ti_rose(
bins,
u_min=u_min,
u_max=u_max,
n_quadrature=n_quadrature,
).T
* 100.0
)
cmap_obj = _cmap(cmap)
colors = [cmap_obj(i) for i in np.linspace(0.0, 1.0, len(bins) - 1)]
labels = [f"{lo:.3f}\u2013{hi:.3f}" for lo, hi in zip(bins[:-1], bins[1:])]
return _draw_rose(
ax,
figsize,
theta_deg,
table,
colors,
labels=labels,
opening=opening,
edgecolor=edgecolor,
legend_title="Turbulence intensity" if legend else None,
title=title,
)
def _plot_site_conditions(
conditions,
*,
ax=None,
bins: int = 30,
n_sectors: int | None = None,
title: str | None = None,
figsize: tuple[float, float] | None = None,
boundary_color: str | None = None,
scatter_color: str = "C0",
hist_color: str = "C0",
) -> np.ndarray:
"""Plot wind rose, turbulence, stability and speed panels for site conditions.
``ax`` holds the four Axes to draw on, the first polar. The four Axes are
returned. ``boundary_color`` defaults to the style's ``axes.edgecolor``.
"""
if boundary_color is None:
boundary_color = matplotlib.rcParams["axes.edgecolor"]
metadata = conditions.metadata or {}
if n_sectors is None:
n_sectors = int(metadata.get("resource", {}).get("n_sectors", 16))
n_axes = 4
if figsize is None:
figsize = (18.0, 4.5)
if ax is None:
fig = _new_figure(figsize, constrained_layout=True)
axes = np.asarray(
[fig.add_subplot(1, n_axes, 1, projection="polar")]
+ [fig.add_subplot(1, n_axes, index) for index in range(2, n_axes + 1)],
dtype=object,
)
else:
axes = np.asarray(ax).ravel()
if len(axes) != n_axes:
raise ValueError(f"ax must contain exactly {n_axes} Axes")
if getattr(axes[0], "name", None) != "polar":
raise ValueError("the first Axes in ax must use a polar projection")
fig = axes[0].figure
u = np.asarray(conditions.wind_speed, dtype=float)
ti = np.asarray(conditions.turbulence_intensity, dtype=float)
weight = np.asarray(conditions.weight, dtype=float)
u_min = float(metadata.get("u_min", 4.0))
u_max = float(metadata.get("u_max", 25.0))
_plot_wind_rose(
conditions,
ax=axes[0],
n_sectors=n_sectors,
title="Wind rose",
)
legend = axes[0].get_legend()
if legend is not None:
handles, labels = axes[0].get_legend_handles_labels()
legend.remove()
axes[0].legend(
handles,
labels,
title="Wind speed [m/s]",
loc="upper center",
bbox_to_anchor=(0.5, -0.08),
frameon=False,
ncols=2,
)
ax = axes[1]
ax.scatter(u, ti, s=12, alpha=0.35, color=scatter_color, edgecolors="none")
# The Dimitrov band depends on wind speed alone, so it is one envelope.
band = metadata.get("resource", {}).get("turbulence", {})
if band.get("type") == "Dimitrov.Turbulence":
u_curve = np.linspace(u_min, u_max, 300)
sigma_low, sigma_high = Turbulence(
band["minimum_intensity"],
band["reference_intensity"],
band.get("maximum_intensity"),
).bounds(u_curve)
_plot_envelope(
ax, u_curve, sigma_low / u_curve, sigma_high / u_curve, boundary_color
)
ax.set(
xlabel="Wind speed [m/s]",
ylabel="Turbulence intensity",
title="Turbulence samples",
)
_style_axes(ax)
ax = axes[2]
inverse_length = np.asarray(conditions.inverse_obukhov_length, dtype=float)
reference_height = float(np.asarray(conditions.reference_height, dtype=float)[0])
values = reference_height * inverse_length
ax.scatter(
u,
values,
s=12,
alpha=0.35,
color=scatter_color,
edgecolors="none",
)
# Binned percentiles of the drawn samples summarize the stability spread.
order = np.argsort(u)
chunks = np.array_split(order, min(12, max(1, len(order) // 25)))
centers, medians, lower, upper = [], [], [], []
for chunk in chunks:
if not len(chunk):
continue
centers.append(np.mean(u[chunk]))
medians.append(np.median(values[chunk]))
lower.append(np.percentile(values[chunk], 5))
upper.append(np.percentile(values[chunk], 95))
if len(centers) > 1:
ax.plot(centers, medians, color=boundary_color, linewidth=1.4, label="median")
ax.plot(centers, lower, color=boundary_color, linewidth=1.0, ls="--")
ax.plot(
centers,
upper,
color=boundary_color,
linewidth=1.0,
ls="--",
label="5th to 95th percentile",
)
ax.legend(frameon=False)
ylabel = rf"Stability $z/L$ at {reference_height:g} m"
ax.axhline(0.0, color=boundary_color, linewidth=1.0)
ax.set(
xlabel="Wind speed [m/s]",
ylabel=ylabel,
title="Atmospheric stability",
)
_style_axes(ax)
ax = axes[3]
ax.hist(
u,
bins=bins,
weights=weight,
color=hist_color,
alpha=0.85,
edgecolor=ax.get_facecolor(),
linewidth=0.5,
)
ax.set(
xlabel="Wind speed [m/s]",
ylabel="Probability",
title="Wind-speed distribution",
)
_style_axes(ax)
if title:
fig.suptitle(title)
return axes
def _plot_envelope(
ax: Axes,
wind_speed: np.ndarray,
lower: np.ndarray,
upper: np.ndarray,
color: str,
) -> None:
"""Draw the exact sampling envelope over its complete wind-speed range."""
ax.fill_between(wind_speed, lower, upper, color=color, alpha=0.08, linewidth=0)
ax.plot(wind_speed, lower, color=color, linewidth=1.2, linestyle="--")
ax.plot(wind_speed, upper, color=color, linewidth=1.2, linestyle="--")
def _plot_stacked_polar_bars(
ax: Axes,
theta: np.ndarray,
table: np.ndarray,
bar_width: float,
colors,
*,
labels=None,
edgecolor="none",
) -> None:
"""Draw stacked polar bars as one collection per row in *table*."""
table = np.asarray(table, dtype=float)
bottom = np.zeros(len(theta), dtype=float)
theta_offsets = np.linspace(-0.5 * bar_width, 0.5 * bar_width, 5)
labels = labels or [None] * len(table)
for heights, color, label in zip(table, colors, labels):
top = bottom + heights
active = top > bottom
if np.any(active):
polygons = []
for center, inner, outer in zip(theta[active], bottom[active], top[active]):
arc = center + theta_offsets
polygons.append(
np.column_stack(
[
np.concatenate([arc, arc[::-1]]),
np.concatenate(
[
np.full_like(arc, outer),
np.full_like(arc, inner),
]
),
]
)
)
collection = PolyCollection(
polygons,
facecolors=[color],
edgecolors=edgecolor,
linewidths=0,
label=label,
)
ax.add_collection(collection)
bottom = top
y_max = float(bottom.max(initial=0.0))
if y_max > 0.0:
ax.set_ylim(0.0, y_max * 1.05)
def _resolve_speed_bins(speed_bins, n_speed_bins, wind_speed=None) -> np.ndarray:
"""Return speed-bin edges, auto-generating them when *speed_bins* is None."""
if speed_bins is not None:
return np.asarray(speed_bins, dtype=float)
if wind_speed is None:
return np.linspace(0.0, 30.0, n_speed_bins + 1)
u_lo = float(np.min(wind_speed))
u_hi = float(np.max(wind_speed))
if u_hi <= u_lo:
u_hi = u_lo + 1.0
return np.linspace(u_lo, u_hi, n_speed_bins + 1)
def _wind_rose_table(source, wind_speed, n_sectors, speed_bins, n_speed_bins):
"""Return ``(sector_centers_deg, table, speed_bin_edges)`` for a wind rose.
``table`` holds percentages, shape ``(n_rows, n_sectors)``. Edges are ``None``
for a single-row, direction-only rose.
"""
# An analytical resource tabulates its direction-dependent Weibull distributions.
if hasattr(source, "wind_rose"):
theta_deg = np.asarray(source.theta, dtype=float)
bins = (
np.linspace(4.0, 25.0, n_speed_bins + 1)
if speed_bins is None
else np.asarray(speed_bins, dtype=float)
)
return theta_deg, source.wind_rose(bins).T * 100.0, bins
# Otherwise bin samples from SiteConditions or a raw direction array.
weights = None
if hasattr(source, "wind_direction"):
if wind_speed is None:
wind_speed = source.wind_speed
weights = np.asarray(source.weight, dtype=float)
source = source.wind_direction
wind_direction = np.asarray(source, dtype=float)
sector_width = 360.0 / n_sectors
# Sector 0 is centered on north, so edges lie half a sector either side.
sector_idx = (
np.floor((wind_direction % 360.0 + sector_width / 2) / sector_width).astype(int)
% n_sectors
)
if weights is None:
weights = np.ones(wind_direction.shape, dtype=float)
if wind_speed is None:
# Directions alone give one row of sector shares.
bins = None
table = np.zeros((1, n_sectors))
np.add.at(table, (0, sector_idx), weights)
else:
wind_speed = np.asarray(wind_speed, dtype=float)
bins = _resolve_speed_bins(speed_bins, n_speed_bins, wind_speed)
speed_idx = np.clip(np.digitize(wind_speed, bins) - 1, 0, len(bins) - 2)
table = np.zeros((len(bins) - 1, n_sectors))
np.add.at(table, (speed_idx, sector_idx), weights)
total = table.sum()
if total > 0:
table = table * 100.0 / total
return sector_angles(n_sectors), table, bins
def _plot_string_guide(
guide,
ax: Axes | None = None,
title: str | None = None,
color: str = "C0",
linewidth: float = 2.0,
) -> Axes:
"""Plot a string guide while preserving visible breaks between segments."""
if ax is None:
ax = _new_axes()
geometry = guide.geometry
lines = geometry.geoms if isinstance(geometry, MultiLineString) else [geometry]
for line in lines:
if not isinstance(line, LineString):
raise TypeError(f"Unsupported string geometry: {type(line).__name__}")
x, y = line.xy
ax.plot(x, y, color=color, linewidth=linewidth)
_style_axes(ax, equal=True, grid=False)
ax.autoscale_view()
ax.set(xlabel="x [m]", ylabel="y [m]")
if title:
ax.set_title(title)
return ax
