# The ImpedanceFitter is a package to fit impedance spectra to
# equivalent-circuit models using open-source software.
#
# Copyright (C) 2018, 2019 Leonard Thiele, leonard.thiele[AT]uni-rostock.de
# Copyright (C) 2018, 2019, 2020 Julius Zimmermann,
# julius.zimmermann[AT]uni-rostock.de
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <https://www.gnu.org/licenses/>.
import logging
import corner
import matplotlib.pyplot as plt
import numpy as np
from scipy.constants import epsilon_0 as e0
from .utils import (
_return_resistance_capacitance,
get_labels,
return_diel_properties,
return_dielectric_modulus,
)
logger = logging.getLogger(__name__)
[docs]
def plot_complex_permittivity(
omega,
Z,
c0,
Z_comp=None,
title="",
show=True,
save=False,
logscale="permittivity",
labels=None,
append=False,
):
"""
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array
c0: double
unit capacitance of device
Z_comp: :class:`numpy.ndarray`, complex, optional
complex-valued impedance array.
Might be used to compare the properties of two data sets.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True)
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_dielectric_properties.pdf`.
logscale: str, optional
Decide what you want to plot using log scale.
Possible are `permittivity`, `loss` and `both`
labels: list, optional
Give custom labels. Needs to be a list of length 2.
append: bool, optional
Decide if you want to show plot or add line to existing plot.
"""
eps_r, cond_fit = return_diel_properties(omega, Z, c0)
if labels is None:
labels = [r"$Z_1$", r"$Z_2$"]
if not len(labels) == 2:
raise ValueError("You need to provide lables as a list containing 2 strings!")
if Z_comp is not None:
eps_r2, cond_fit2 = return_diel_properties(omega, Z_comp, c0)
plt.figure()
plt.suptitle("complex permittivity", y=1.05)
plt.subplot(211)
plt.title("Relative permittivity")
plt.ylabel("Relative permittivity")
plt.xlabel("Frequency / Hz")
if logscale == "permittivity" or logscale == "both":
plt.yscale("log")
plt.xscale("log")
plt.plot(omega / (2.0 * np.pi), eps_r, label=labels[0])
if Z_comp is not None:
plt.plot(omega / (2.0 * np.pi), eps_r2, label=labels[1])
plt.legend()
plt.subplot(212)
plt.title("Dielectric loss")
plt.ylabel("Dielectric loss")
plt.xlabel("Frequency / Hz")
if logscale == "loss" or logscale == "both":
plt.yscale("log")
plt.xscale("log")
plt.plot(omega / (2.0 * np.pi), cond_fit / (e0 * omega), label=labels[0])
if Z_comp is not None:
plt.plot(omega / (2.0 * np.pi), cond_fit2 / (e0 * omega), label=labels[1])
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_complex_permittivity.pdf")
if show:
plt.show()
else:
plt.close()
[docs]
def plot_dielectric_modulus(
omega,
Z,
c0,
Z_comp=None,
title="",
show=True,
save=False,
logscale=None,
labels=None,
append=False,
):
"""
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array
c0: double
unit capacitance of device
Z_comp: :class:`numpy.ndarray`, complex, optional
complex-valued impedance array.
Might be used to compare the properties of two data sets.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True)
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_dielectric_properties.pdf`.
logscale: str, optional
Decide what you want to plot using log scale.
Possible are `ReM`, `ImM` and `both`
labels: list, optional
Give custom labels. Needs to be a list of length 2.
append: bool, optional
Whether to leave figure active for later plotting
"""
ReM, ImM = return_dielectric_modulus(omega, Z, c0)
if labels is None:
labels = [r"$Z_1$", r"$Z_2$"]
if not len(labels) == 2:
raise ValueError("You need to provide lables as a list containing 2 strings!")
if Z_comp is not None:
ReM2, ImM2 = return_dielectric_modulus(omega, Z_comp, c0)
plt.figure()
plt.suptitle("Real part", y=1.05)
plt.subplot(211)
plt.ylabel("Re M")
plt.xlabel("Frequency / Hz")
if logscale == "ReM" or logscale == "both":
plt.yscale("log")
plt.xscale("log")
plt.plot(omega / (2.0 * np.pi), ReM, label=labels[0])
if Z_comp is not None:
plt.plot(omega / (2.0 * np.pi), ReM2, label=labels[1])
plt.legend()
plt.subplot(212)
plt.title("Imaginary part")
plt.ylabel("Im M")
plt.xlabel("Frequency / Hz")
if logscale == "ImM" or logscale == "both":
plt.yscale("log")
plt.xscale("log")
plt.plot(omega / (2.0 * np.pi), ImM, label=labels[0])
if Z_comp is not None:
plt.plot(omega / (2.0 * np.pi), ImM2, label=labels[1])
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_dielectric_modulus.pdf")
if show:
plt.show()
else:
plt.close()
# ruff: noqa: C901
[docs]
def plot_dielectric_properties(
omega,
Z,
c0,
Z_comp=None,
title="",
show=True,
save=False,
logscale="permittivity",
labels=None,
append=False,
markers=[None, None],
legend=True,
limits=None,
**plotkwargs,
):
"""
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array
c0: double
unit capacitance of device
Z_comp: :class:`numpy.ndarray`, complex, optional
complex-valued impedance array.
Might be used to compare the properties of two data sets.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True)
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_dielectric_properties.pdf`.
logscale: str, optional
Decide what you want to plot using log scale.
Possible are `permittivity`, `conductivity` and `both`
labels: list, optional
Give custom labels. Needs to be a list of length 2.
append: bool, optional
Decide if you want to show plot or add line to existing plot.
legend: bool, optional
Switch legend on/off
markers: list
Two entries to choose custom markers for each property
plotkwargs: kwargs
further keyword arguments passed to matplotlib
limits: list, optional
pass lower and upper limit for y-axis
"""
eps_r, cond_fit = return_diel_properties(omega, Z, c0)
axes = []
plt.figure("dielectricproperties")
axes = plt.gcf().axes
if len(axes) < 2:
plt.suptitle(title, y=1.05)
plt.subplot(211)
else:
plt.sca(axes[0])
if labels is None:
labels = [r"$Z_1$", r"$Z_2$"]
if not len(labels) == 2:
raise ValueError("You need to provide lables as a list containing 2 strings!")
if Z_comp is not None:
eps_r2, cond_fit2 = return_diel_properties(omega, Z_comp, c0)
plt.title("Relative permittivity")
plt.ylabel("Relative permittivity")
plt.xlabel("Frequency / Hz")
if limits:
plt.ylim(limits[0])
if logscale == "permittivity" or logscale == "both":
plt.yscale("log")
plt.xscale("log")
plt.plot(
omega / (2.0 * np.pi), eps_r, label=labels[0], marker=markers[0], **plotkwargs
)
if Z_comp is not None:
plt.plot(
omega / (2.0 * np.pi),
eps_r2,
linestyle="--",
label=labels[1],
marker=markers[1],
**plotkwargs,
)
if legend:
plt.legend()
if len(axes) < 2:
plt.subplot(212)
else:
plt.sca(axes[1])
plt.title("Conductivity")
plt.ylabel(r"Conductivity / S$\cdot$m$^{-1}$")
if limits:
plt.ylim(limits[1])
plt.xlabel("Frequency / Hz")
if logscale == "conductivity" or logscale == "both":
plt.yscale("log")
plt.xscale("log")
plt.plot(
omega / (2.0 * np.pi),
cond_fit,
label=labels[0],
marker=markers[0],
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
omega / (2.0 * np.pi),
cond_fit2,
linestyle="--",
label=labels[1],
marker=markers[1],
**plotkwargs,
)
if legend:
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_dielectric_properties.pdf")
if show:
plt.show()
else:
plt.close()
# ruff: noqa: C901
[docs]
def plot_comparison_dielectric_properties(
omega,
Z,
c0,
Z_comp,
title="",
show=True,
save=False,
residual="relative",
label=None,
append=False,
marker=None,
legend=True,
limits=None,
**plotkwargs,
):
"""
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array
c0: double
unit capacitance of device
Z_comp: :class:`numpy.ndarray`, complex
complex-valued impedance array to compare the properties of two data sets.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True)
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_dielectric_properties.pdf`.
logscale: str, optional
Decide what you want to plot using log scale.
Possible are `permittivity`, `conductivity` and `both`
label: str, optional
Give custom label. Needs to be a str.
append: bool, optional
Decide if you want to show plot or add line to existing plot.
legend: bool, optional
Switch legend on/off
limits: list, optional
pass lower and upper limit for y-axis
plotkwargs: kwargs
further keyword arguments passed to matplotlib
residual: str
Plot relative difference w.r.t. conductivity and
permittivity if `relative`.
Plot relative difference w.r.t. absolute value if `absolute`.
marker: list, optional
Two entries for marker for conductivity and permittivity curve
"""
eps_r, cond_fit = return_diel_properties(omega, Z, c0)
eps_r2, cond_fit2 = return_diel_properties(omega, Z_comp, c0)
if residual == "relative":
diff_eps = 100.0 * np.abs((eps_r - eps_r2) / eps_r2)
diff_sigma = 100.0 * np.abs((cond_fit - cond_fit2) / cond_fit2)
ylabel = "Relative difference / %"
elif residual == "absolute":
diff_eps = eps_r - eps_r2
diff_sigma = cond_fit - cond_fit2
ylabel = r"Absolute Difference"
else:
raise RuntimeError(
f"residual must be either `relative` or `absolute` but not {residual}"
)
axes = []
plt.figure("comparisondielectricproperties")
axes = plt.gcf().axes
if len(axes) < 2:
plt.suptitle(title, y=1.05)
plt.subplot(211)
else:
plt.sca(axes[0])
plt.title("Relative permittivity")
plt.ylabel(ylabel)
plt.xlabel("Frequency / Hz")
if limits:
plt.ylim(limits[0])
plt.xscale("log")
plt.plot(omega / (2.0 * np.pi), diff_eps, label=label, marker=marker, **plotkwargs)
if legend:
plt.legend()
if len(axes) < 2:
plt.subplot(212)
else:
plt.sca(axes[1])
plt.title("Conductivity")
if residual == "absolute":
ylabel += r" / S$\cdot$m$^{-1}$"
plt.ylabel(ylabel)
if limits:
plt.ylim(limits[1])
plt.xlabel("Frequency / Hz")
plt.xscale("log")
plt.plot(
omega / (2.0 * np.pi), diff_sigma, label=label, marker=marker, **plotkwargs
)
if legend:
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(
str(title).replace(" ", "_")
+ f"_comparison_{residual}_dielectric_properties.pdf"
)
if show:
plt.show()
else:
plt.close()
[docs]
def plot_cole_cole(
omega,
Z,
c0,
Z_comp=None,
append=False,
legend=True,
markers=[None, None],
title="",
show=True,
save=False,
labels=None,
limits=None,
):
"""
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array
c0: double
unit capacitance of device
Z_comp: :class:`numpy.ndarray`, complex, optional
complex-valued impedance array.
Might be used to compare the properties of two data sets.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True)
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_dielectric_properties.pdf`.
logscale: str, optional
Decide what you want to plot using log scale.
Possible are `permittivity`, `conductivity` and `both`
labels: list, optional
Give custom labels. Needs to be a list of length 2.
limits: list, optional
pass lower and upper limit for y-axis
append: bool, optional
Decide if you want to show plot or add line to existing plot.
markers: list
Three entries to choose custom markers for each impedance curve
legend: str, optional
Choose if a legend should be shown. Recommended to switch to False
when using large datasets.
"""
axes = []
plt.figure("colecoleplot")
axes = plt.gcf().axes
if len(axes) == 1:
plt.sca(axes[0])
eps_r, cond_fit = return_diel_properties(omega, Z, c0)
epsc_fit = eps_r - 1j * cond_fit / (e0 * omega)
if labels is None:
labels = [r"$Z_1$", r"$Z_2$"]
if not len(labels) == 2:
raise ValueError("You need to provide lables as a list containing 2 strings!")
if Z_comp is not None:
eps_r2, cond_fit2 = return_diel_properties(omega, Z_comp, c0)
epsc_fit2 = eps_r2 - 1j * cond_fit2 / (e0 * omega)
plt.title("Cole-Cole plot")
plt.xlabel(r"Re $\varepsilon$")
plt.ylabel(r"-Im $\varepsilon$")
plt.plot(epsc_fit.real, -epsc_fit.imag, label=labels[0], marker=markers[0])
if limits is not None:
plt.xlim(limits[0])
plt.ylim(limits[1])
if Z_comp is not None:
plt.plot(epsc_fit2.real, -epsc_fit2.imag, label=labels[1], marker=markers[1])
if legend:
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_cole_cole_plot.pdf")
if show:
plt.show()
else:
plt.close()
return
# ruff: noqa: C901
[docs]
def plot_bode(
omega,
Z,
title="",
Z_fit=None,
show=True,
save=False,
Z_comp=None,
labels=["Data", "Best fit", "Init fit"],
append=False,
legend=True,
markers=[None, "^", "v"],
):
"""Bode plot of impedance.
Plots phase and log of magnitude over log of frequency.
Parameters
----------
omega: :class:`numpy.ndarray`, double
Frequency array
Z: :class:`numpy.ndarray`, complex
Impedance array, experimental data or data to compare to.
Z_fit: :class:`numpy.ndarray`, complex
Impedance array, fit result. If provided, the difference
between data and fit will be shown.
title: str
Title of plot.
show: bool, optional
Show figure (default is True).
save: bool, optional
Save figure to pdf (default is False).
Name of figure starts with `title`.
Z_comp: :class:`numpy.ndarray`, complex, optional
Complex-valued impedance array.
Might be used to compare the properties of two data sets.
labels: list
List of labels for three plots. Must have length 3 always.
Is ordered like: `[Z, Z_fit, Z_comp]`
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_bode_plot.pdf`.
append: bool, optional
Decide if you want to show plot or add line to existing plot.
legend: str, optional
Choose if a legend should be shown. Recommended to switch to False
when using large datasets.
markers: list
Three entries to choose custom markers for each impedance curve
"""
axes = []
plt.figure("bodeimpedance")
axes = plt.gcf().axes
if len(axes) < 2:
plt.suptitle(title, y=1.05)
plt.subplot(211)
else:
plt.sca(axes[0])
# plot real part of impedance
plt.xscale("log")
plt.yscale("log")
plt.ylabel(r"|Z| / $\Omega$")
plt.xlabel("Frequency / Hz")
plt.plot(omega / (2.0 * np.pi), np.abs(Z), label=labels[0], marker=markers[0])
if Z_fit is not None:
plt.plot(
omega / (2.0 * np.pi), np.abs(Z_fit), label=labels[1], marker=markers[1]
)
if Z_comp is not None:
plt.plot(
omega / (2.0 * np.pi), np.abs(Z_comp), label=labels[2], marker=markers[2]
)
if legend:
plt.legend()
if len(axes) < 2:
plt.subplot(212)
else:
plt.sca(axes[1])
plt.xscale("log")
plt.ylabel("Phase / °")
plt.xlabel("Frequency / Hz")
plt.plot(
omega / (2.0 * np.pi), np.angle(Z, deg=True), label=labels[0], marker=markers[0]
)
if Z_fit is not None:
plt.plot(
omega / (2.0 * np.pi),
np.angle(Z_fit, deg=True),
label=labels[1],
marker=markers[1],
)
if Z_comp is not None:
plt.plot(
omega / (2.0 * np.pi),
np.angle(Z_comp, deg=True),
label=labels[2],
marker=markers[2],
)
if legend:
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_bode_plot.pdf")
if show:
plt.show()
else:
plt.close()
[docs]
def plot_resistance_capacitance(
omega,
Z,
title="",
Z_fit=None,
show=True,
save=False,
Z_comp=None,
labels=["Data", "Best fit", "Init fit"],
append=False,
legend=True,
):
"""R-C plot of impedance.
Plots phase and log of magnitude over log of frequency.
Parameters
----------
omega: :class:`numpy.ndarray`, double
Frequency array
Z: :class:`numpy.ndarray`, complex
Impedance array, experimental data or data to compare to.
Z_fit: :class:`numpy.ndarray`, complex
Impedance array, fit result. If provided, the difference
between data and fit will be shown.
title: str
Title of plot.
show: bool, optional
Show figure (default is True).
save: bool, optional
Save figure to pdf (default is False).
Name of figure starts with `title`.
Z_comp: :class:`numpy.ndarray`, complex, optional
Complex-valued impedance array.
Might be used to compare the properties of two data sets.
labels: list
List of labels for three plots. Must have length 3 always.
Is ordered like: `[Z, Z_fit, Z_comp]`
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_bode_plot.pdf`.
append: bool, optional
Decide if you want to show plot or add line to existing plot.
legend: str, optional
Choose if a legend should be shown. Recommended to switch to False
when using large datasets.
"""
axes = []
plt.figure("rcimpedance")
axes = plt.gcf().axes
if len(axes) < 2:
plt.suptitle(title, y=1.05)
plt.subplot(211)
else:
plt.sca(axes[0])
# plot real part of impedance
plt.xscale("log")
plt.yscale("log")
plt.ylabel(r"R / $\Omega$")
R, C = _return_resistance_capacitance(omega, Z)
plt.xlabel("Frequency / Hz")
plt.plot(omega / (2.0 * np.pi), R, label=labels[0])
if Z_fit is not None:
R_fit, C_fit = _return_resistance_capacitance(omega, Z_fit)
plt.plot(omega / (2.0 * np.pi), R_fit, "^", label=labels[1])
if Z_comp is not None:
R_comp, C_comp = _return_resistance_capacitance(omega, Z_comp)
plt.plot(omega / (2.0 * np.pi), C_comp, "v", label=labels[2])
if legend:
plt.legend()
if len(axes) < 2:
plt.subplot(212)
else:
plt.sca(axes[1])
plt.xscale("log")
plt.yscale("log")
plt.ylabel("C / F")
plt.xlabel("Frequency / Hz")
plt.plot(omega / (2.0 * np.pi), C, label=labels[0])
if Z_fit is not None:
plt.plot(omega / (2.0 * np.pi), C_fit, "^", label=labels[1])
if Z_comp is not None:
plt.plot(omega / (2.0 * np.pi), C_comp, "v", label=labels[2])
if legend:
plt.legend()
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_rc_plot.pdf")
if show:
plt.show()
else:
plt.close()
# ruff: noqa: C901
[docs]
def plot_impedance(
omega,
Z,
title="",
Z_fit=None,
show=True,
save=False,
Z_comp=None,
labels=["Data", "Best fit", "Init fit"],
residual="parts",
sign=False,
Zlog=False,
append=False,
limits_residual=None,
Nyquist_conventional=False,
omega_fit=None,
omega_comp=None,
legend=True,
compare=True,
**plotkwargs,
):
"""Plot the `result` and compare it to data `Z`.
Generates 4 subplots showing the real and imaginary parts over
the frequency; a Nyquist plot of real and negative imaginary part
and the relative differences of real and imaginary part as well
as absolute value of impedance.
Parameters
----------
omega: :class:`numpy.ndarray`, double
Frequency array
Z: :class:`numpy.ndarray`, complex
Impedance array, experimental data or data to compare to.
Z_fit: :class:`numpy.ndarray`, complex
Impedance array, fit result. If provided, the difference
between data and fit will be shown.
title: str
Title of plot.
show: bool, optional
Show figure (default is True).
save: bool, optional
Save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_impedance_overview.pdf`.
Z_comp: :class:`numpy.ndarray`, complex, optional
Complex-valued impedance array.
Might be used to compare the properties of two data sets.
labels: list
List of labels for three plots. Must have length 3 always.
Is ordered like: `[Z, Z_fit, Z_comp]`
residual: str
Plot relative difference w.r.t. real and imaginary part if `parts`.
Plot relative difference w.r.t. absolute value if `absolute`.
Plot difference (residual) if `diff`.
sign: bool, optional
Use sign of residual. Default is False, i.e. absolute value is plotted.
Zlog: bool, optional
Log-scale of impedance
append: bool, optional
Decide if you want to show plot or add line to existing plot.
omega_fit: :class:`numpy.ndarray`, double, optional
Frequency array, provide only if fitted impedance was evaluated at
different frequencies than the experimental data
omega_comp: :class:`numpy.ndarray`, double, optional
Frequency array, provide only if fitted impedance was evaluated at
different frequencies than the experimental data
legend: str, optional
Choose if a legend should be shown. Recommended to switch to False
when using large datasets.
compare: bool, optional
Choose if the difference between fit and data should be computed.
Nyquist_conventional: bool, optional
Choose if the Nyquist plot should use the same limits for real
and imaginary part.
limits_residual: list, optional
List with entries `[bottom, top]` for y-axis of residual plot.
plotkwargs: kwargs
further keyword arguments passed to matplotlib
"""
if omega_fit is None:
omega_fit = omega
if omega_comp is None:
omega_comp = omega
axes = []
plt.figure("impedance")
axes = plt.gcf().axes
if len(axes) < 3:
plt.suptitle(title, y=1.05)
plt.subplot(221)
else:
plt.sca(axes[0])
# use logscale
if Zlog:
plt.yscale("log")
# plot real part of impedance
plt.xscale("log")
plt.title("Impedance real part")
plt.ylabel(r"Re Z / $\Omega$")
plt.xlabel("Frequency / Hz")
line_Z = plt.plot(omega / (2.0 * np.pi), Z.real, label=labels[0], **plotkwargs)
if Z_fit is not None:
line_Zfit = plt.plot(
omega_fit / (2.0 * np.pi),
Z_fit.real,
linestyle="--",
label=labels[1],
**plotkwargs,
)
if Z_comp is not None:
line_Zcomp = plt.plot(
omega_comp / (2.0 * np.pi),
Z_comp.real,
linestyle="-.",
label=labels[2],
**plotkwargs,
)
if legend:
plt.legend()
# plot imaginary part of impedance
if len(axes) < 3:
plt.subplot(222)
else:
plt.sca(axes[1])
plt.title("Impedance imaginary part")
plt.xscale("log")
plt.ylabel(r"Im Z / $\Omega$")
plt.xlabel("Frequency / Hz")
if Zlog:
plt.yscale("log")
if np.all(np.less_equal(Z.imag, 0)):
plt.ylabel(r"-Im Z / $\Omega$")
plt.plot(
omega / (2.0 * np.pi),
-Z.imag,
label=labels[0],
color=line_Z[0].get_color(),
**plotkwargs,
)
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
-Z_fit.imag,
"--",
label=labels[1],
color=line_Zfit[0].get_color(),
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
-Z_comp.imag,
"-.",
label=labels[2],
color=line_Zcomp[0].get_color(),
**plotkwargs,
)
elif np.all(np.greater_equal(Z.imag, 0)):
plt.plot(omega / (2.0 * np.pi), Z.imag, label=labels[0], **plotkwargs)
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
Z_fit.imag,
"--",
label=labels[1],
color=line_Zfit[0].get_color(),
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
Z_comp.imag,
"-.",
label=labels[2],
color=line_Zcomp[0].get_color(),
**plotkwargs,
)
elif np.where(Z.imag < 0)[0].size > np.where(Z.imag > 0)[0].size:
plt.ylabel(r"-Im Z / $\Omega$")
plt.plot(
omega / (2.0 * np.pi),
-Z.imag,
label=labels[0],
color=line_Z[0].get_color(),
**plotkwargs,
)
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
-Z_fit.imag,
"--",
label=labels[1],
color=line_Zfit[0].get_color(),
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
-Z_comp.imag,
"-.",
label=labels[2],
color=line_Zcomp[0].get_color(),
**plotkwargs,
)
else:
plt.plot(
omega / (2.0 * np.pi),
Z.imag,
label=labels[0],
color=line_Z[0].get_color(),
**plotkwargs,
)
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
Z_fit.imag,
"--",
label=labels[1],
color=line_Zfit[0].get_color(),
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
Z_comp.imag,
"-.",
label=labels[2],
color=line_Zcomp[0].get_color(),
**plotkwargs,
)
else:
plt.plot(omega / (2.0 * np.pi), Z.imag, label=labels[0], **plotkwargs)
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
Z_fit.imag,
"--",
label=labels[1],
color=line_Zfit[0].get_color(),
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
Z_comp.imag,
"-.",
label=labels[2],
color=line_Zcomp[0].get_color(),
**plotkwargs,
)
if legend:
plt.legend()
# plot real vs negative imaginary part
if len(axes) < 3:
plt.subplot(223)
else:
plt.sca(axes[2])
plt.title("Nyquist plot")
plt.ylabel(r"-Im Z / $\Omega$")
plt.xlabel(r"Re Z / $\Omega$")
plt.plot(Z.real, -Z.imag, "o", label=labels[0], **plotkwargs)
if Zlog:
plt.xscale("log")
plt.yscale("log")
# set limits of Nyquist plot
if Nyquist_conventional:
Zmin_real = np.min(Z.real)
Zmin_imag = np.min(-Z.imag)
Zmax_real = np.max(Z.real)
Zmax_imag = np.max(-Z.imag)
if Z_fit is not None:
Zmin_real = np.min([Zmin_real, np.min(Z_fit.real)])
Zmin_imag = np.min([Zmin_imag, np.min(-Z_fit.imag)])
Zmax_real = np.max([Zmax_real, np.max(Z_fit.real)])
Zmax_imag = np.max([Zmax_imag, np.max(-Z_fit.imag)])
if Z_comp is not None:
Zmin_real = np.min([Zmin_real, np.min(Z_comp.real)])
Zmin_imag = np.min([Zmin_imag, np.min(-Z_comp.imag)])
Zmax_real = np.max([Zmax_real, np.max(Z_comp.real)])
Zmax_imag = np.max([Zmax_imag, np.max(-Z_comp.imag)])
Zmin = np.min([Zmin_real, Zmin_imag])
Zmax = np.max([Zmax_real, Zmax_imag])
plt.xlim(left=0.8 * Zmin, right=1.2 * Zmax)
plt.ylim(bottom=0.8 * Zmin, top=1.2 * Zmax)
if Z_fit is not None:
plt.plot(
Z_fit.real,
-Z_fit.imag,
"^",
label=labels[1],
color=line_Zfit[0].get_color(),
**plotkwargs,
)
if Z_comp is not None:
plt.plot(
Z_comp.real,
-Z_comp.imag,
"v",
label=labels[2],
color=line_Zcomp[0].get_color(),
**plotkwargs,
)
if legend:
plt.legend()
if Z_fit is not None and np.all(omega == omega_fit) and compare:
plot_compare_to_data(
omega,
Z,
Z_fit,
subplot=224,
residual=residual,
sign=sign,
limits=limits_residual,
legend=legend,
**plotkwargs,
)
plt.tight_layout()
if append:
return
if save:
cleantitle = str(title).replace(" ", "_").replace("/", "")
plt.savefig(cleantitle + "_impedance_overview.pdf")
if show:
plt.show()
else:
plt.close()
# ruff: noqa: C901
[docs]
def plot_compare_to_data(
omega,
Z,
Z_fit,
subplot=None,
title="",
show=True,
save=False,
residual="parts",
sign=False,
limits=None,
impedance_threshold=1.0,
legend=True,
**plotkwargs,
):
"""
Plots the difference of the fitted function to the data.
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array, experimental data or data to compare to
Z_fit: :class:`numpy.ndarray`, complex
impedance array, fit result
subplot: optional
decide whether it is a new figure or a subplot.
Default is None (yields new figure).
Otherwise it can be an integer to denote the subfigure.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True). Only has an effect when `subplot` is None.
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_relative_difference_to_data.pdf` or `_difference_to_data.pdf`.
relative: str, optional
Plot relative difference if True, else plot residual (i.e. just difference).
sign: bool, optional
Use sign of residual. Default is False, i.e. absolute value is plotted.
residual: str
Plot relative difference w.r.t. real and imaginary part if `parts`.
Plot relative difference w.r.t. absolute value if `absolute`.
Plot difference (residual) if `diff`.
limits: list, optional
List with entries `[bottom, top]` for y-axis of residual plot.
impedance_threshold: double, optional
Threshold for impedance around 0, which is disregarded in the relative
differences plot. Default is that impedances, with an absolute value less than
0 are not considered.
legend: str, optional
Choose if a legend should be shown. Recommended to switch to False
when using large datasets.
plotkwargs: kwargs
further keyword arguments passed to matplotlib
Notes
-----
When computing the relative difference, impedances between -1 and 1 Ohm are not
considered since they might lead to a blow up of the relative difference
(close to division by 0).
Instead of this quantitative measure, qualitative checks should be done.
"""
if subplot is None:
plt.figure()
else:
show = False
plt.subplot(subplot)
if residual == "parts":
close_to_zero_real = np.where(np.isclose(Z.real, 0.0, atol=impedance_threshold))
close_to_zero_imag = np.where(np.isclose(Z.imag, 0.0, atol=impedance_threshold))
diff_real = 100.0 * (Z.real - Z_fit.real) / Z.real
diff_imag = 100.0 * (Z.imag - Z_fit.imag) / Z.imag
diff_real[close_to_zero_real] = np.nan
diff_imag[close_to_zero_imag] = np.nan
diff_abs = 100.0 * np.abs((Z - Z_fit) / Z)
label = "Relative difference / %"
elif residual == "absolute":
diff_real = 100.0 * (Z.real - Z_fit.real) / np.abs(Z)
diff_imag = 100.0 * (Z.imag - Z_fit.imag) / np.abs(Z)
label = "Relative difference / %"
elif residual == "diff":
diff_real = Z.real - Z_fit.real
diff_imag = Z.imag - Z_fit.imag
diff_abs = np.abs(Z - Z_fit)
label = r"Difference / $\Omega$"
else:
raise RuntimeError(
f"""residual must be either `parts`, `absolute` or `diff`,
but not {residual}"""
)
if not sign:
diff_real = np.abs(diff_real)
diff_imag = np.abs(diff_imag)
plt.xscale("log")
if title is not None:
plt.title((str(title) + "relative difference to data").capitalize())
plt.xlabel("Frequency / Hz")
plt.ylabel(label)
plt.plot(omega / (2.0 * np.pi), diff_real, label="Real part", **plotkwargs)
plt.plot(
omega / (2.0 * np.pi),
diff_imag,
label="Imag part",
linestyle="--",
**plotkwargs,
)
if residual != "absolute":
plt.plot(
omega / (2.0 * np.pi),
diff_abs,
label="Abs value",
linestyle="-.",
**plotkwargs,
)
if limits is not None:
if not len(limits) == 2:
raise ValueError("You need to provide upper and lower limit!")
plt.ylim(limits)
if legend:
plt.legend()
if subplot is None:
plt.tight_layout()
if save:
if residual != "diff":
plt.savefig(
str(title).replace(" ", "_") + "_relative_difference_to_data.pdf"
)
else:
plt.savefig(str(title).replace(" ", "_") + "_difference_to_data.pdf")
if show:
plt.show()
elif not show and subplot is None:
plt.close()
[docs]
def emcee_plot(res, clustered=False, **corner_kwargs):
"""
Create corner plot.
Parameters
----------
res: dict
Dictionary containing values and flatchain
clustered: bool
decide which flatchain to use
corner_kwargs: dict, optional
Dictionary with further corner plot options
"""
params = [p for p in res.params if res.params[p].vary]
truths = [res.params[p].value for p in params]
ifLabels = get_labels(params)
labels = [ifLabels[p] for p in res.var_names]
if not clustered:
flatchain = res.flatchain
else:
flatchain = res.new_flatchain
plot = corner.corner(flatchain, labels=labels, truths=truths, **corner_kwargs)
return plot
[docs]
def plot_uncertainty(omega, Zdata, Z, Z1, Z2, sigma, show=True, model=None):
"""Plot best fit with uncertainty interval.
Parameters
----------
omega: :class:`numpy.ndarray`, float
frequency array
Zdata: :class:`numpy.ndarray`, complex
impedance array of experimental data
Z: :class:`numpy.ndarray`, complex
impedance array of best fit
Z1: :class:`numpy.ndarray`, complex
impedance array of upper uncertainty limit
Z2: :class:`numpy.ndarray`, complex
impedance array of lower uncertainty limit
sigma: double
confidence level
show: bool, optional
show figure (default is True)
model: int, optional
numbering of model for sequential plotting
"""
plt.figure()
if model is not None:
plt.suptitle(rf"${sigma} \sigma$ results", y=1.05)
else:
plt.suptitle(rf"${sigma} \sigma$ results, model {model}", y=1.05)
# plot real part of impedance
plt.subplot(211)
plt.xscale("log")
plt.title("Impedance real part")
plt.ylabel(r"Re Z / $\Omega$")
plt.xlabel("Frequency / Hz")
plt.plot(omega / (2.0 * np.pi), Z.real, "^", label="Best fit")
plt.plot(omega / (2.0 * np.pi), Zdata.real, "r", label="Data")
plt.fill_between(
omega / (2.0 * np.pi),
Z1.real,
Z2.real,
color="#888888",
label=rf"${sigma} \sigma$ interval",
)
plt.legend()
# plot imaginary part of impedance
plt.subplot(212)
plt.title("Impedance imaginary part")
plt.xscale("log")
plt.ylabel(r"Im Z / $\Omega$")
plt.xlabel("Frequency / Hz")
plt.plot(omega / (2.0 * np.pi), Z.imag, "^", label="Best fit")
plt.plot(omega / (2.0 * np.pi), Zdata.imag, "r", label="Data")
plt.fill_between(
omega / (2.0 * np.pi),
Z1.imag,
Z2.imag,
color="#888888",
label=rf"${sigma} \sigma$ interval",
)
plt.legend()
plt.tight_layout()
if show:
plt.show()
else:
plt.close()
# ruff: noqa: C901
[docs]
def plot_admittance(
omega,
Z,
title="",
Z_fit=None,
show=True,
save=False,
Z_comp=None,
labels=["Data", "Best fit", "Init fit"],
residual="parts",
sign=False,
Zlog=False,
append=False,
limits_residual=None,
omega_fit=None,
omega_comp=None,
legend=True,
compare=True,
):
"""Plot the admittance and compare it to data 1/`Z`.
Generates 4 subplots showing the real and imaginary parts over
the frequency; a Nyquist plot of real and negative imaginary part
and the relative differences of real and imaginary part as well as
absolute value of admittance.
Parameters
----------
omega: :class:`numpy.ndarray`, double
Frequency array
Z: :class:`numpy.ndarray`, complex
Impedance array, experimental data or data to compare to.
Z_fit: :class:`numpy.ndarray`, complex
Impedance array, fit result. If provided, the difference
between data and fit will be shown.
title: str
Title of plot.
show: bool, optional
Show figure (default is True).
save: bool, optional
Save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_admittance_overview.pdf`.
Z_comp: :class:`numpy.ndarray`, complex, optional
Complex-valued impedance array.
Might be used to compare the properties of two data sets.
labels: list
List of labels for three plots. Must have length 3 always.
Is ordered like: `[Z, Z_fit, Z_comp]`
residual: str
Plot relative difference w.r.t. real and imaginary part if `parts`.
Plot relative difference w.r.t. absolute value if `absolute`.
Plot difference (residual) if `diff`.
sign: bool, optional
Use sign of residual. Default is False, i.e. absolute value is plotted.
Zlog: bool, optional
Log-scale of impedance
append: bool, optional
Decide if you want to show plot or add line to existing plot.
omega_fit: :class:`numpy.ndarray`, double, optional
Frequency array, provide only if fitted impedance was evaluated at
different frequencies than the experimental data
omega_comp: :class:`numpy.ndarray`, double, optional
Frequency array, provide only if fitted impedance was evaluated at
different frequencies than the experimental data
legend: str, optional
Choose if a legend should be shown. Recommended to switch to False
when using large datasets.
compare: bool, optional
Choose if the difference between fit and data should be computed.
limits_residual: list, optional
List with entries `[bottom, top]` for y-axis of residual plot.
"""
if omega_fit is None:
omega_fit = omega
if omega_comp is None:
omega_comp = omega
axes = []
plt.figure("admittance")
axes = plt.gcf().axes
if len(axes) < 3:
plt.suptitle(title, y=1.05)
plt.subplot(221)
else:
plt.sca(axes[0])
# use logscale
if Zlog:
plt.yscale("log")
# plot real part of admittance
plt.xscale("log")
plt.title("Admittance real part")
plt.ylabel(r"Re Y / S")
plt.xlabel("Frequency / Hz")
plt.plot(omega / (2.0 * np.pi), (1.0 / Z).real, label=labels[0])
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
(1.0 / Z_fit).real,
linestyle="--",
label=labels[1],
lw=3,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
(1.0 / Z_comp).real,
linestyle="-.",
label=labels[2],
lw=3,
)
if legend:
plt.legend()
# plot imaginary part of impedance
if len(axes) < 3:
plt.subplot(222)
else:
plt.sca(axes[1])
plt.title("Admittance imaginary part")
plt.xscale("log")
plt.ylabel(r"Im Y / S")
plt.xlabel("Frequency / Hz")
if Zlog:
plt.yscale("log")
if np.all(np.less_equal(Z.imag, 0)):
plt.ylabel(r"Im Y / S")
plt.plot(omega / (2.0 * np.pi), (1.0 / Z).imag, label=labels[0])
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi), (1.0 / Z_fit).imag, "--", label=labels[1]
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
(1.0 / Z_comp).imag,
"-.",
label=labels[2],
)
elif np.all(np.greater_equal(Z.imag, 0)):
plt.plot(omega / (2.0 * np.pi), (1.0 / Z).imag, label=labels[0])
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi), (1.0 / Z_fit).imag, "--", label=labels[1]
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
(1.0 / Z_comp).imag,
"-.",
label=labels[2],
)
elif np.where(Z.imag < 0).size > np.where(Z.imag > 0).size:
plt.ylabel(r"Im Y / S")
plt.plot(omega / (2.0 * np.pi), (1.0 / Z).imag, label=labels[0])
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi), (1.0 / Z_fit).imag, "--", label=labels[1]
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
(1.0 / Z_comp).imag,
"-.",
label=labels[2],
)
else:
plt.plot(omega / (2.0 * np.pi), (1.0 / Z).imag, label=labels[0])
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi), (1.0 / Z_fit).imag, "--", label=labels[1]
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
(1.0 / Z_comp).imag,
"-.",
label=labels[2],
)
else:
plt.plot(omega / (2.0 * np.pi), (1.0 / Z).imag, label=labels[0])
if Z_fit is not None:
plt.plot(
omega_fit / (2.0 * np.pi),
(1.0 / Z_fit).imag,
"--",
label=labels[1],
lw=3,
)
if Z_comp is not None:
plt.plot(
omega_comp / (2.0 * np.pi),
(1.0 / Z_comp).imag,
"-.",
label=labels[2],
lw=3,
)
if legend:
plt.legend()
# plot real vs negative imaginary part
if len(axes) < 3:
plt.subplot(223)
else:
plt.sca(axes[2])
plt.title("Nyquist plot")
plt.ylabel(r"Im Y / S")
plt.xlabel(r"Re Y / S")
if Zlog:
plt.xscale("log")
plt.yscale("log")
plt.plot((1.0 / Z).real, (1.0 / Z).imag, "o", label=labels[0])
if Z_fit is not None:
plt.plot((1.0 / Z_fit).real, (1.0 / Z_fit).imag, "^", label=labels[1])
if Z_comp is not None:
plt.plot((1.0 / Z_comp).real, (1.0 / Z_comp).imag, "v", label=labels[2])
if legend:
plt.legend()
if Z_fit is not None and np.all(omega == omega_fit) and compare:
plot_compare_to_data(
omega,
Z,
Z_fit,
subplot=224,
residual=residual,
sign=sign,
limits=limits_residual,
legend=legend,
)
plt.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_admittance_overview.pdf")
if show:
plt.show()
else:
plt.close()
[docs]
def plot_dielectric_dispersion(
omega,
Z,
c0,
Z_comp=None,
title="",
show=True,
save=False,
logscale="permittivity",
labels=None,
append=False,
**plotkwargs,
):
"""
Parameters
----------
omega: :class:`numpy.ndarray`, double
frequency array
Z: :class:`numpy.ndarray`, complex
impedance array
c0: double
unit capacitance of device
Z_comp: :class:`numpy.ndarray`, complex, optional
complex-valued impedance array.
Might be used to compare the properties of two data sets.
title: str, optional
title of plot. Default is an empty string.
show: bool, optional
show figure (default is True)
save: bool, optional
save figure to pdf (default is False). Name of figure starts with `title`
and ends with `_dielectric_properties.pdf`.
logscale: str, optional
Decide what you want to plot using log scale.
Possible are `permittivity`, `conductivity` and `both`
labels: list, optional
Give custom labels. Needs to be a list of length 2.
append: bool, optional
Decide if you want to show plot or add line to existing plot.
plotkwargs: kwargs
further keyword arguments passed to matplotlib
"""
eps_r, cond_fit = return_diel_properties(omega, Z, c0)
fig, ax1 = plt.subplots()
plt.title(title, y=1.05)
if labels is None:
labels = [r"$Z_1$", r"$Z_2$"]
if not len(labels) == 2:
raise ValueError("You need to provide lables as a list containing 2 strings!")
if Z_comp is not None:
eps_r2, cond_fit2 = return_diel_properties(omega, Z_comp, c0)
ax1.set_ylabel(r"Relative permittivity")
ax1.set_xlabel("Frequency / Hz")
if logscale == "permittivity" or logscale == "both":
ax1.set_yscale("log")
ax1.set_xscale("log")
ax1.plot(omega / (2.0 * np.pi), eps_r, label=labels[0], **plotkwargs)
if Z_comp is not None:
ax1.plot(omega / (2.0 * np.pi), eps_r2, label=labels[1], **plotkwargs)
ax2 = ax1.twinx()
ax2.set_ylabel(r"Dielectric loss")
if logscale == "conductivity" or logscale == "both":
ax2.set_yscale("log")
ax2.plot(omega / (2.0 * np.pi), cond_fit / (e0 * omega), ls="-.", **plotkwargs)
if Z_comp is not None:
plt.plot(omega / (2.0 * np.pi), cond_fit2 / (e0 * omega), ls="-.", **plotkwargs)
ax1.legend()
fig.tight_layout()
if append:
return
if save:
plt.savefig(str(title).replace(" ", "_") + "_dielectric_dispersion.pdf")
if show:
plt.show()
else:
plt.close()
[docs]
def plot_time_domain_signals_with_impedance(
t,
frequencies,
voltage,
current,
impedance,
impedance_expected=None,
save_file="impedance_time_domain.pdf",
t_zoom_range=(0.05, 0.25),
current_scale=50,
):
"""Plot impedance with original time domain signals."""
# TODO more documentation
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(15, 10))
# Plot high-frequency input signals (zoomed to pulse region)
t_zoom = t[(t >= t_zoom_range[0]) & (t <= t_zoom_range[1])]
v_zoom = voltage[(t >= t_zoom_range[0]) & (t <= t_zoom_range[1])]
i_zoom = current[(t >= t_zoom_range[0]) & (t <= t_zoom_range[1])]
ax1.plot(t_zoom * 1000, v_zoom, "b-", label="Voltage", linewidth=2)
ax1.plot(
t_zoom * 1000,
i_zoom * current_scale,
"r-",
label=f"Current x{current_scale}",
linewidth=2,
)
ax1.set_xlabel("Time / ms")
ax1.set_ylabel("Amplitude")
ax1.legend()
ax1.grid(True)
# Plot full signal overview
ax2.plot(t, voltage, "b-", label="Voltage", alpha=0.7)
ax2.plot(
t, current * current_scale, "r-", label=f"Current x{current_scale}", alpha=0.7
)
ax2.set_xlabel("Time / s")
ax2.set_ylabel("Amplitude")
ax2.legend()
ax2.grid(True)
# Plot impedance magnitude
if impedance_expected is not None:
ax3.loglog(
frequencies, np.abs(impedance_expected), "blue", lw=3, label="Reconstructed"
)
ax3.loglog(frequencies, np.abs(impedance), "ro-", markersize=6, label="Fitted")
ax3.set_xlabel("Frequency / Hz")
ax3.set_ylabel(r"|Z| / $\Omega$")
ax3.grid(True)
# Plot impedance phase
if impedance_expected is not None:
ax4.semilogx(
frequencies,
np.angle(impedance_expected, deg=True),
"blue",
lw=3,
label="Expected",
)
ax4.semilogx(
frequencies, np.angle(impedance, deg=True), "ro-", markersize=6, label="Fitted"
)
ax4.set_xlabel("Frequency / Hz")
ax4.set_ylabel(r"Phase / °")
ax4.grid(True)
plt.tight_layout()
plt.savefig(save_file)
plt.close()