from mhkit.river.resource import exceedance_probability
from mhkit.river.graphics import _xy_plot
from mhkit.utils import convert_to_dataset
import matplotlib.patheffects as pe
import matplotlib.pyplot as plt
from matplotlib import gridspec
import pandas as pd
import xarray as xr
import numpy as np
import matplotlib
import calendar
[docs]
def plot_spectrum(S, ax=None):
"""
Plots wave amplitude spectrum versus omega
Parameters
------------
S: pandas Series, pandas DataFrame, xarray DataArray, or xarray Dataset
Spectral density [m^2/Hz] indexed frequency [Hz]
ax : matplotlib axes object
Axes for plotting. If None, then a new figure is created.
Returns
---------
ax : matplotlib pyplot axes
"""
S = convert_to_dataset(S)
frequency_dimension = list(S.dims)[0]
f = S[frequency_dimension]
for var in S.data_vars:
ax = _xy_plot(
f * 2 * np.pi,
S[var] / (2 * np.pi),
fmt="-",
xlabel="omega [rad/s]",
ylabel="Spectral density [m$^2$s/rad]",
ax=ax,
)
return ax
[docs]
def plot_elevation_timeseries(eta, ax=None):
"""
Plot wave surface elevation time-series
Parameters
----------
eta: pandas Series, pandas DataFrame, xarray DataArray, or xarray Dataset
Wave surface elevation [m] indexed by time [datetime or s]
ax : matplotlib axes object
Axes for plotting. If None, then a new figure is created.
Returns
-------
ax : matplotlib pyplot axes
"""
eta = convert_to_dataset(eta)
time_dimension = list(eta.dims)[0]
t = eta[time_dimension]
for var in eta.data_vars:
ax = _xy_plot(t, eta[var], fmt="-", xlabel="Time", ylabel="$\\eta$ [m]", ax=ax)
return ax
[docs]
def plot_matrix(M, xlabel="Te", ylabel="Hm0", zlabel=None, show_values=True, ax=None):
"""
Plots values in the matrix as a scatter diagram
Parameters
------------
M: pandas Series, pandas DataFrame, xarray DataArray
Matrix with numeric labels for x and y axis, and numeric entries.
An example would be the average capture length matrix generated by
mhkit.device.wave, or something similar.
xlabel: string (optional)
Title of the x-axis
ylabel: string (optional)
Title of the y-axis
zlabel: string (optional)
Colorbar label
show_values : bool (optional)
Show values on the scatter diagram
ax : matplotlib axes object
Axes for plotting. If None, then a new figure is created.
Returns
---------
ax : matplotlib pyplot axes
"""
try:
M = pd.DataFrame(M)
except:
pass
if not isinstance(M, pd.DataFrame):
raise TypeError(f"M must be of type pd.DataFrame. Got: {type(M)}")
if ax is None:
plt.figure()
ax = plt.gca()
im = ax.imshow(M, origin="lower", aspect="auto")
# Add colorbar
cbar = plt.colorbar(im)
if zlabel:
cbar.set_label(zlabel, rotation=270, labelpad=15)
# Set x and y label
ax.set_xlabel(xlabel)
ax.set_ylabel(ylabel)
# Show values in the plot
if show_values:
for i, col in enumerate(M.columns):
for j, index in enumerate(M.index):
if not np.isnan(M.loc[index, col]):
ax.text(
i, j, format(M.loc[index, col], ".2f"), ha="center", va="center"
)
# Reset x and y ticks
ax.set_xticks(np.arange(len(M.columns)))
ax.set_yticks(np.arange(len(M.index)))
ax.set_xticklabels(M.columns)
ax.set_yticklabels(M.index)
return ax
[docs]
def plot_chakrabarti(H, lambda_w, D, ax=None):
"""
Plots, in the style of Chakrabarti (2005), relative importance of viscous,
inertia, and diffraction phemonena
Chakrabarti, Subrata. Handbook of Offshore Engineering (2-volume set).
Elsevier, 2005.
Examples:
**Using floats**
>>> plt.figure()
>>> D = 5
>>> H = 8
>>> lambda_w = 200
>>> wave.graphics.plot_chakrabarti(H, lambda_w, D)
**Using numpy array**
>>> plt.figure()
>>> D = np.linspace(5,15,5)
>>> H = 8*np.ones_like(D)
>>> lambda_w = 200*np.ones_like(D)
>>> wave.graphics.plot_chakrabarti(H, lambda_w, D)
**Using pandas DataFrame**
>>> plt.figure()
>>> D = np.linspace(5,15,5)
>>> H = 8*np.ones_like(D)
>>> lambda_w = 200*np.ones_like(D)
>>> df = pd.DataFrame([H.flatten(),lambda_w.flatten(),D.flatten()], index=['H','lambda_w','D']).transpose()
>>> wave.graphics.plot_chakrabarti(df.H, df.lambda_w, df.D)
Parameters
----------
H: int, float, numpy array, pandas Series, or xarray DataArray
Wave height [m]
lambda_w: int, float, numpy array, pandas Series, or xarray DataArray
Wave length [m]
D: int, float, numpy array, pandas Series, or xarray DataArray
Characteristic length [m]
ax : matplotlib axes object (optional)
Axes for plotting. If None, then a new figure is created.
Returns
-------
ax : matplotlib pyplot axes
"""
if not isinstance(H, (np.ndarray, float, int, np.int64, pd.Series, xr.DataArray)):
raise TypeError(
f"H must be of type float, int, np.int64, np.ndarray, pd.Series, or xr.DataArray. Got: {type(H)}"
)
if not isinstance(
lambda_w, (np.ndarray, float, int, np.int64, pd.Series, xr.DataArray)
):
raise TypeError(
f"lambda_w must be of type float, int, np.int64, np.ndarray, pd.Series, or xr.DataArray. Got: {type(lambda_w)}"
)
if not isinstance(D, (np.ndarray, float, int, np.int64, pd.Series, xr.DataArray)):
raise TypeError(
f"D must be of type float, int, np.int64, np.ndarray, pd.Series, or xr.DataArray. Got: {type(D)}"
)
if any(
[
isinstance(H, (np.ndarray, pd.Series, xr.DataArray)),
isinstance(lambda_w, (np.ndarray, pd.Series, xr.DataArray)),
isinstance(D, (np.ndarray, pd.Series, xr.DataArray)),
]
):
n_H = H.squeeze().shape
n_lambda_w = lambda_w.squeeze().shape
n_D = D.squeeze().shape
if not (n_H == n_lambda_w and n_H == n_D):
raise ValueError("D, H, and lambda_w must be same shape")
if isinstance(H, np.ndarray):
mvals = pd.DataFrame(H.reshape(len(H), 1), columns=["H"])
mvals["lambda_w"] = lambda_w
mvals["D"] = D
elif isinstance(H, (pd.Series, xr.DataArray)):
mvals = pd.DataFrame(H)
mvals["lambda_w"] = lambda_w
mvals["D"] = D
else:
H = np.array([H])
lambda_w = np.array([lambda_w])
D = np.array([D])
mvals = pd.DataFrame(H.reshape(len(H), 1), columns=["H"])
mvals["lambda_w"] = lambda_w
mvals["D"] = D
if ax is None:
plt.figure()
ax = plt.gca()
ax.set_xscale("log")
ax.set_yscale("log")
for index, row in mvals.iterrows():
H = row.H
D = row.D
lambda_w = row.lambda_w
KC = H / D
Diffraction = np.pi * D / lambda_w
label = f"$H$ = {H:g}, $\\lambda_w$ = {lambda_w:g}, $D$ = {D:g}"
ax.plot(Diffraction, KC, "o", label=label)
if (
np.any(KC >= 10 or KC <= 0.02)
or np.any(Diffraction >= 50)
or np.any(lambda_w >= 1000)
):
ax.autoscale(enable=True, axis="both", tight=True)
else:
ax.set_xlim((0.01, 10))
ax.set_ylim((0.01, 50))
graphScale = list(ax.get_xlim())
if graphScale[0] >= 0.01:
graphScale[0] = 0.01
# deep water breaking limit (H/lambda_w = 0.14)
x = np.logspace(1, np.log10(graphScale[0]), 2)
y_breaking = 0.14 * np.pi / x
ax.plot(x, y_breaking, "k-")
graphScale = list(ax.get_xlim())
ax.text(
1,
7,
"wave\nbreaking\n$H/\\lambda_w > 0.14$",
ha="center",
va="center",
fontstyle="italic",
fontsize="small",
clip_on="True",
)
# upper bound of low drag region
ldv = 20
y_small_drag = 20 * np.ones_like(graphScale)
graphScale[1] = 0.14 * np.pi / ldv
ax.plot(graphScale, y_small_drag, "k--")
ax.text(
0.0125,
30,
"drag",
ha="center",
va="top",
fontstyle="italic",
fontsize="small",
clip_on="True",
)
# upper bound of small drag region
sdv = 1.5
y_small_drag = sdv * np.ones_like(graphScale)
graphScale[1] = 0.14 * np.pi / sdv
ax.plot(graphScale, y_small_drag, "k--")
ax.text(
0.02,
7,
"inertia \n& drag",
ha="center",
va="center",
fontstyle="italic",
fontsize="small",
clip_on="True",
)
# upper bound of negligible drag region
ndv = 0.25
graphScale[1] = 0.14 * np.pi / ndv
y_small_drag = ndv * np.ones_like(graphScale)
ax.plot(graphScale, y_small_drag, "k--")
ax.text(
8e-2,
0.7,
"large\ninertia",
ha="center",
va="center",
fontstyle="italic",
fontsize="small",
clip_on="True",
)
ax.text(
8e-2,
6e-2,
"all\ninertia",
ha="center",
va="center",
fontstyle="italic",
fontsize="small",
clip_on="True",
)
# left bound of diffraction region
drv = 0.5
graphScale = list(ax.get_ylim())
graphScale[1] = 0.14 * np.pi / drv
x_diff_reg = drv * np.ones_like(graphScale)
ax.plot(x_diff_reg, graphScale, "k--")
ax.text(
2,
6e-2,
"diffraction",
ha="center",
va="center",
fontstyle="italic",
fontsize="small",
clip_on="True",
)
if index > 0:
ax.legend(fontsize="xx-small", ncol=2)
ax.set_xlabel("Diffraction parameter, $\\frac{\\pi D}{\\lambda_w}$")
ax.set_ylabel("KC parameter, $\\frac{H}{D}$")
plt.tight_layout()
[docs]
def plot_environmental_contour(x1, x2, x1_contour, x2_contour, **kwargs):
"""
Plots an overlay of the x1 and x2 variables to the calculate
environmental contours.
Parameters
----------
x1: list, np.ndarray, pd.Series, xr.DataArray
x-axis data
x2: list, np.ndarray, pd.Series, xr.DataArray
x-axis data
x1_contour: list, np.ndarray, pd.Series, xr.DataArray
Calculated x1 contour values
x2_contour: list, np.ndarray, pd.Series, xr.DataArray
Calculated x2 contour values
**kwargs : optional
x_label: string (optional)
x-axis label. Default None.
y_label: string (optional)
y-axis label. Default None.
data_label: string (optional)
Legend label for x1, x2 data (e.g. 'Buoy 46022').
Default None.
contour_label: string or list of strings (optional)
Legend label for x1_contour, x2_contour countor data
(e.g. '100-year contour'). Default None.
ax : matplotlib axes object (optional)
Axes for plotting. If None, then a new figure is created.
Default None.
markers: string
string or list of strings to use as marker types
Returns
-------
ax : matplotlib pyplot axes
"""
try:
x1 = x1.values
except:
pass
try:
x2 = x2.values
except:
pass
if not isinstance(x1, np.ndarray):
raise TypeError(f"x1 must be of type np.ndarray. Got: {type(x1)}")
if not isinstance(x2, np.ndarray):
raise TypeError(f"x2 must be of type np.ndarray. Got: {type(x2)}")
try:
x1_contour = x1_contour.values
except:
pass
try:
x2_contour = x2_contour.values
except:
pass
if not isinstance(x1_contour, (np.ndarray, list)):
raise TypeError(
f"x1_contour must be of type np.ndarray or list. Got: {type(x1_contour)}"
)
if not isinstance(x2_contour, (np.ndarray, list)):
raise TypeError(
f"x2_contour must be of type np.ndarray or list. Got: {type(x2_contour)}"
)
x_label = kwargs.get("x_label", None)
y_label = kwargs.get("y_label", None)
data_label = kwargs.get("data_label", None)
contour_label = kwargs.get("contour_label", None)
ax = kwargs.get("ax", None)
markers = kwargs.get("markers", "-")
if not isinstance(data_label, (str, type(None))):
raise TypeError(
f"If specified, data_label must be of type str. Got: {type(data_label)}"
)
if not isinstance(contour_label, (str, list, type(None))):
raise TypeError(
f"If specified, contour_label be of type str. Got: {type(contour_label)}"
)
if isinstance(markers, str):
markers = [markers]
if not isinstance(markers, list) or not all(
[isinstance(marker, (str)) for marker in markers]
):
raise TypeError(
f"markers must be of type str or list of strings. Got: {markers}"
)
if not len(x2_contour) == len(x1_contour):
raise ValueError(
f"contour must be of equal dimension got {len(x2_contour)} and {len(x1_contour)}"
)
if isinstance(x1_contour, np.ndarray):
N_contours = 1
x2_contour = [x2_contour]
x1_contour = [x1_contour]
elif isinstance(x1_contour, list):
N_contours = len(x1_contour)
if contour_label != None:
if isinstance(contour_label, str):
contour_label = [contour_label]
N_c_labels = len(contour_label)
if not N_c_labels == N_contours:
raise ValueError(
"If specified, the number of contour labels must"
" be equal to number the number of contour years."
f" Got: {N_c_labels} and {N_contours}"
)
else:
contour_label = [None] * N_contours
if len(markers) == 1:
markers = markers * N_contours
if not len(markers) == N_contours:
raise ValueError(
"Markers must be same length as N contours specified."
f"Got: {len(markers)} and {len(x1_contour)}"
)
for i in range(N_contours):
contour1 = np.array(x1_contour[i]).T
contour2 = np.array(x2_contour[i]).T
ax = _xy_plot(contour1, contour2, markers[i], label=contour_label[i], ax=ax)
plt.plot(x1, x2, "bo", alpha=0.1, label=data_label)
plt.legend(loc="lower right")
plt.xlabel(x_label)
plt.ylabel(y_label)
plt.tight_layout()
return ax
[docs]
def plot_avg_annual_energy_matrix(
Hm0,
Te,
J,
time_index=None,
Hm0_bin_size=None,
Te_bin_size=None,
Hm0_edges=None,
Te_edges=None,
):
"""
Creates an average annual energy matrix with frequency of occurance.
Parameters
----------
Hm0: array-like
Significant wave height
Te: array-like
Energy period
J: array-like
Energy flux
time_index: DateTime Index
time to index by. Optional default None. If None Passed parameters must be series indexed by Datetime.
Hm0_bin_size: float, int
Creates edges of bin using this discrtization. Optional default None. If not passed must pass Hm0_edges.
Te_bin_size: float, int
Creates edges of bin using this discrtization. Optional default None. If not passed must pass Te_edges.
Hm0_edges: array-like
Defines the Hm0 bin edges to use. Optional default None.
Te_edges: array-like
Defines the Te bin edges to use. Optional default None.
Returns
-------
fig: Figure
Average annual energy table plot
"""
fig = plt.figure()
if isinstance(time_index, type(None)):
data = pd.DataFrame(dict(Hm0=Hm0, Te=Te, J=J))
else:
data = pd.DataFrame(dict(Hm0=Hm0, Te=Te, J=J), index=time_index)
years = data.index.year.unique()
if isinstance(Hm0_edges, type(None)):
Hm0_max = data.Hm0.max()
Hm0_edges = np.arange(0, Hm0_max + Hm0_bin_size, Hm0_bin_size)
if isinstance(Te_edges, type(None)):
Te_max = data.Te.max()
Te_edges = np.arange(0, Te_max + Te_bin_size, Te_bin_size)
# Dict for number of hours each sea state occurs
hist_counts = {}
hist_J = {}
# Create hist of counts, and weghted by J for each year
for year in years:
year_data = data.loc[str(year)].copy(deep=True)
# Get the counts of each bin
counts, xedges, yedges = np.histogram2d(
year_data.Te,
year_data.Hm0,
bins=(Te_edges, Hm0_edges),
)
# Get centers for number of counts plot location
xcenters = xedges[:-1] + np.diff(xedges)
ycenters = yedges[:-1] + np.diff(yedges)
year_data["xbins"] = np.digitize(year_data.Te, xcenters)
year_data["ybins"] = np.digitize(year_data.Hm0, ycenters)
total_year_J = year_data.J.sum()
H = counts.copy()
for i in range(len(xcenters)):
for j in range(len(ycenters)):
bin_J = year_data[
(year_data.xbins == i) & (year_data.ybins == j)
].J.sum()
H[i][j] = bin_J / total_year_J
# Save in results dict
hist_counts[year] = counts
hist_J[year] = H
# Calculate avg annual
avg_annual_counts_hist = sum(hist_counts.values()) / len(years)
avg_annual_J_hist = sum(hist_J.values()) / len(years)
# Create a mask of non-zero weights to hide from imshow
Hmasked = np.ma.masked_where(~(avg_annual_J_hist > 0), avg_annual_J_hist)
plt.imshow(
Hmasked.T,
interpolation="none",
vmin=0.005,
origin="lower",
aspect="auto",
extent=[xedges[0], xedges[-1], yedges[0], yedges[-1]],
)
# Plot number of counts as text on the hist of annual avg J
for xi in range(len(xcenters)):
for yi in range(len(ycenters)):
if avg_annual_counts_hist[xi][yi] != 0:
plt.text(
xedges[xi],
yedges[yi],
int(np.ceil(avg_annual_counts_hist[xi][yi])),
fontsize=10,
color="white",
path_effects=[pe.withStroke(linewidth=1, foreground="k")],
)
plt.xlabel("Wave Energy Period (s)")
plt.ylabel("Significant Wave Height (m)")
cbar = plt.colorbar()
cbar.set_label("Mean Normalized Annual Energy")
plt.tight_layout()
return fig
[docs]
def monthly_cumulative_distribution(J):
"""
Creates a cumulative distribution of energy flux as described in
IEC TS 62600-101.
Parameters
----------
J: pd.Series, xr.DataArray
Energy Flux with DateTime index
Returns
-------
ax: axes
Figure of monthly cumulative distribution
"""
J = pd.Series(J)
cumSum = {}
months = J.index.month.unique()
for month in months:
F = exceedance_probability(J[J.index.month == month])
cumSum[month] = 1 - F / 100
cumSum[month].sort_values("F", inplace=True)
plt.figure(figsize=(12, 8))
for month in months:
plt.semilogx(
J.loc[cumSum[month].index],
cumSum[month].F,
"--",
label=calendar.month_abbr[month],
)
F = exceedance_probability(J)
F.sort_values("F", inplace=True)
ax = plt.semilogx(
J.loc[F.index], 1 - F["F"] / 100, "k-", fillstyle="none", label="All"
)
plt.grid()
plt.xlabel("Energy Flux")
plt.ylabel("Cumulative Distribution")
plt.legend()
return ax
[docs]
def plot_compendium(Hs, Tp, Dp, buoy_title=None, ax=None):
"""
Create subplots showing: Significant Wave Height (Hs), Peak Period (Tp),
and Direction (Dp) using OPeNDAP service from CDIP THREDDS Server.
See http://cdip.ucsd.edu/themes/cdip?pb=1&bl=cdip?pb=1&d2=p70&u3=s:100:st:1:v:compendium:dt:201204 for example Compendium plot.
Developed based on: http://cdip.ucsd.edu/themes/media/docs/documents/html_pages/compendium.html
Parameters
----------
Hs: pandas Series or xarray DataArray
significant wave height
Tp: pandas Series or xarray DataArray
significant wave height
Dp: pandas Series or xarray DataArray
significant wave height
buoy_title: string (optional)
Buoy title from the CDIP THREDDS Server
ax : matplotlib axes object (optional)
Axes for plotting. If None, then a new figure is created.
Returns
-------
ax : matplotlib pyplot axes
"""
Hs = pd.Series(Hs)
Tp = pd.Series(Tp)
Dp = pd.Series(Dp)
if not isinstance(Hs, pd.Series):
raise TypeError(f"Hs must be of type pd.Series. Got: {type(Hs)}")
if not isinstance(Tp, pd.Series):
raise TypeError(f"Tp must be of type pd.Series. Got: {type(Tp)}")
if not isinstance(Dp, pd.Series):
raise TypeError(f"Dp must be of type pd.Series. Got: {type(Dp)}")
if not isinstance(buoy_title, (str, type(None))):
raise TypeError(
f"If specified, buoy_title must be of type string. Got: {type(buoy_title)}"
)
f, (pHs, pTp, pDp) = plt.subplots(3, 1, sharex=True, figsize=(15, 10))
pHs.plot(Hs.index, Hs, "b")
pTp.plot(Tp.index, Tp, "b")
pDp.scatter(Dp.index, Dp, color="blue", s=5)
pHs.tick_params(axis="x", which="major", labelsize=12, top="off")
pHs.set_ylim(0, 8)
pHs.tick_params(axis="y", which="major", labelsize=12, right="off")
pHs.set_ylabel("Hs [m]", fontsize=18)
pHs.grid(color="b", linestyle="--")
pHs2 = pHs.twinx()
pHs2.set_ylim(0, 25)
pHs2.set_ylabel("Hs [ft]", fontsize=18)
# Peak Period, Tp
pTp.set_ylim(0, 28)
pTp.set_ylabel("Tp [s]", fontsize=18)
pTp.grid(color="b", linestyle="--")
# Direction, Dp
pDp.set_ylim(0, 360)
pDp.set_ylabel("Dp [deg]", fontsize=18)
pDp.grid(color="b", linestyle="--")
pDp.set_xlabel("Day", fontsize=18)
# Set x-axis tick interval to every 5 days
degrees = 70
days = matplotlib.dates.DayLocator(interval=5)
daysFmt = matplotlib.dates.DateFormatter("%Y-%m-%d")
plt.gca().xaxis.set_major_locator(days)
plt.gca().xaxis.set_major_formatter(daysFmt)
plt.setp(pDp.xaxis.get_majorticklabels(), rotation=degrees)
# Set Titles
month_name_start = Hs.index.month_name()[0][:3]
year_start = Hs.index.year[0]
month_name_end = Hs.index.month_name()[-1][:3]
year_end = Hs.index.year[-1]
plt.suptitle(buoy_title, fontsize=30)
plt.title(f"{Hs.index[0].date()} to {Hs.index[-1].date()}", fontsize=20)
ax = f
return ax
[docs]
def plot_boxplot(Hs, buoy_title=None):
"""
Create plot of monthly-averaged boxes of Significant Wave Height (Hs)
data.
Developed based on:
http://cdip.ucsd.edu/themes/media/docs/documents/html_pages/annualHs_plot.html
Parameters
------------
Hs: pandas Series or xarray DataArray
Spectral density [m^2/Hz] indexed frequency [Hz]
buoy_title: string (optional)
Buoy title from the CDIP THREDDS Server
ax : matplotlib axes object (optional)
Axes for plotting. If None, then a new figure is created.
Returns
---------
ax : matplotlib pyplot axes
"""
Hs = pd.Series(Hs)
if not isinstance(Hs, pd.Series):
raise TypeError(f"Hs must be of type pd.Series. Got: {type(Hs)}")
if not isinstance(buoy_title, (str, type(None))):
raise TypeError(
f"If specified, buoy_title must be of type string. Got: {type(buoy_title)}"
)
months = Hs.index.month
means = Hs.groupby(months).mean()
monthlengths = Hs.groupby(months).count()
fig = plt.figure(figsize=(10, 12))
gs = gridspec.GridSpec(2, 1, height_ratios=[4, 1])
boxprops = dict(color="k")
whiskerprops = dict(linestyle="--", color="k")
flierprops = dict(marker="+", color="r", markeredgecolor="r", markerfacecolor="r")
medianprops = dict(linewidth=2.5, color="firebrick")
meanprops = dict(linewidth=2.5, marker="_", markersize=25)
bp = plt.subplot(gs[0, :])
Hs_months = Hs.to_frame().groupby(months)
bp = Hs_months.boxplot(
subplots=False,
boxprops=boxprops,
whiskerprops=whiskerprops,
flierprops=flierprops,
medianprops=medianprops,
showmeans=True,
meanprops=meanprops,
)
# Add values of monthly means as text
for i, mean in enumerate(means):
bp.annotate(
np.round(mean, 2),
(means.index[i], mean),
fontsize=12,
horizontalalignment="center",
verticalalignment="bottom",
color="g",
)
# Create a second row of x-axis labels for top subplot
newax = bp.twiny()
newax.tick_params(which="major", direction="in", pad=-18)
newax.set_xlim(bp.get_xlim())
newax.xaxis.set_ticks_position("top")
newax.xaxis.set_label_position("top")
newax.set_xticks(np.arange(1, 13, 1))
newax.set_xticklabels(monthlengths, fontsize=10)
# Sample 'legend' boxplot, to go underneath actual boxplot
bp_sample2 = np.random.normal(2.5, 0.5, 500)
bp2 = plt.subplot(gs[1, :])
meanprops = dict(linewidth=2.5, marker="|", markersize=25)
bp2_example = bp2.boxplot(
bp_sample2, vert=False, flierprops=flierprops, medianprops=medianprops
)
sample_mean = 2.3
bp2.scatter(sample_mean, 1, marker="|", color="g", linewidths=1.0, s=200)
for line in bp2_example["medians"]:
xm, ym = line.get_xydata()[0]
for line in bp2_example["boxes"]:
xb, yb = line.get_xydata()[0]
for line in bp2_example["whiskers"]:
xw, yw = line.get_xydata()[0]
bp2.annotate("Median", [xm - 0.1, ym - 0.3 * ym], fontsize=10, color="firebrick")
bp2.annotate("Mean", [sample_mean - 0.1, 0.65], fontsize=10, color="g")
bp2.annotate("25%ile", [xb - 0.05 * xb, yb - 0.15 * yb], fontsize=10)
bp2.annotate("75%ile", [xb + 0.26 * xb, yb - 0.15 * yb], fontsize=10)
bp2.annotate("Outliers", [xw + 0.3 * xw, yw - 0.3 * yw], fontsize=10, color="r")
if buoy_title:
plt.suptitle(buoy_title, fontsize=30, y=0.97)
bp.set_title("Significant Wave Height by Month", fontsize=20, y=1.01)
bp2.set_title("Sample Boxplot", fontsize=10, y=1.02)
# Set axes labels and ticks
months_text = [m[:3] for m in Hs.index.month_name().unique()]
bp.set_xticklabels(months_text, fontsize=12)
bp.set_ylabel("Significant Wave Height, Hs (m)", fontsize=14)
bp.tick_params(axis="y", which="major", labelsize=12, right="off")
bp.tick_params(axis="x", which="major", labelsize=12, top="off")
# Plot horizontal gridlines onto top subplot
bp.grid(axis="x", color="b", linestyle="-", alpha=0.25)
# Remove tickmarks from bottom subplot
bp2.axes.get_xaxis().set_visible(False)
bp2.axes.get_yaxis().set_visible(False)
ax = fig
return ax
[docs]
def plot_directional_spectrum(
spectrum,
color_level_min=None,
fill=True,
nlevels=11,
name="Elevation Variance",
units="m^2",
):
"""
Create a contour polar plot of a directional spectrum.
Parameters
------------
spectrum: xarray.DataArray
Spectral data indexed frequency [Hz] and wave direction [deg].
color_level_min: float (optional)
Minimum color bar level.
fill: bool
Whether to use `contourf` (filled) instead of `contour` (lines).
nlevels: int
Number of contour levels to plot.
name: str
Name of the (integral) spectrum variable.
units: str
Units of the (integral) spectrum variable.
Returns
---------
ax : matplotlib pyplot axes
"""
if not isinstance(spectrum, xr.DataArray):
raise TypeError(f"spectrum must be of type xr.DataArray. Got: {type(spectrum)}")
if not isinstance(color_level_min, (type(None), float)):
raise TypeError(
f"If specified, color_level_min must be of type float. Got: {type(color_level_min)}"
)
if not isinstance(fill, bool):
raise TypeError(f"If specified, fill must be of type bool. Got: {type(fill)}")
if not isinstance(nlevels, int):
raise TypeError(
f"If specified, nlevels must be of type int. Got: {type(nlevels)}"
)
if not isinstance(name, str):
raise TypeError(f"If specified, name must be of type string. Got: {type(name)}")
if not isinstance(units, str):
raise TypeError(
f"If specified, units must be of type string. Got: {type(units)}"
)
a, f = np.meshgrid(np.deg2rad(spectrum.direction), spectrum.frequency)
_, ax = plt.subplots(subplot_kw=dict(projection="polar"))
tmp = np.floor(np.min(spectrum.data) * 10) / 10
color_level_min = tmp if (color_level_min is None) else color_level_min
color_level_max = np.ceil(np.max(spectrum.data) * 10) / 10
levels = np.linspace(color_level_min, color_level_max, nlevels)
if fill:
c = ax.contourf(a, f, spectrum, levels=levels)
else:
c = ax.contour(a, f, spectrum, levels=levels)
cbar = plt.colorbar(c)
cbar.set_label(f"Spectrum [{units}/Hz/deg]", rotation=270, labelpad=20)
ax.set_title(f"{name} Spectrum")
ylabels = ax.get_yticklabels()
ylabels = [ilabel.get_text() for ilabel in ax.get_yticklabels()]
ylabels = [ilabel + "Hz" for ilabel in ylabels]
ticks_loc = ax.get_yticks()
ax.set_yticks(ticks_loc)
ax.set_yticklabels(ylabels)
return ax