Skip to content

Move plot plasma current into summary.py #4170

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
je-cook merged 1 commit into from
Apr 9, 2026
Merged

Move plot plasma current into summary.py #4170

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
90 changes: 88 additions & 2 deletions process/core/io/plot/summary.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -60,7 +60,6 @@
from process.models.physics.l_h_transition import PlasmaConfinementTransitionModel
from process.models.physics.physics import BetaComponentLimits, BetaNormMaxModel
from process.models.physics.plasma_current import (
PlasmaCurrent,
PlasmaCurrentModel,
PlasmaDiamagneticCurrentModel,
)
Expand Down Expand Up @@ -11724,6 +11723,93 @@ def plot_plasma_pressure_profiles(axis: plt.Axes, mfile: MFile, scan: int):
)


def plot_plasma_current_comparison(axis: plt.Axes, mfile: MFile, scan: int):
"""Function to plot a scatter box plot of different plasma current comparisons.

Parameters
----------
axis :
Axis object to plot to.
mfile :
MFILE data object.
scan :
Scan number to use.
"""
c_plasma_peng_analytic = mfile.get("c_plasma_peng_analytic", scan=scan)
c_plasma_peng_double_null = mfile.get("c_plasma_peng_double_null", scan=scan)
c_plasma_cyclindrical = mfile.get("c_plasma_cyclindrical", scan=scan)
c_plasma_ipdg89 = mfile.get("c_plasma_ipdg89", scan=scan)
c_plasma_todd_empirical_i = mfile.get("c_plasma_todd_empirical_i", scan=scan)
c_plasma_todd_empirical_ii = mfile.get("c_plasma_todd_empirical_ii", scan=scan)
c_plasma_connor_hastie = mfile.get("c_plasma_connor_hastie", scan=scan)
c_plasma_sauter = mfile.get("c_plasma_sauter", scan=scan)
c_plasma_fiesta_st = mfile.get("c_plasma_fiesta_st", scan=scan)

# Data for the box plot
data = {
f"{PlasmaCurrentModel.PENG_ANALYTIC_FIT.full_name}": c_plasma_peng_analytic,
f"{PlasmaCurrentModel.PENG_DIVERTOR_SCALING.full_name}": c_plasma_peng_double_null,
f"{PlasmaCurrentModel.ITER_SCALING.full_name}": c_plasma_cyclindrical,
f"{PlasmaCurrentModel.IPDG89_SCALING.full_name}": c_plasma_ipdg89,
f"{PlasmaCurrentModel.TODD_EMPIRICAL_SCALING_I.full_name}": c_plasma_todd_empirical_i,
f"{PlasmaCurrentModel.TODD_EMPIRICAL_SCALING_II.full_name}": c_plasma_todd_empirical_ii,
f"{PlasmaCurrentModel.CONNOR_HASTIE_MODEL.full_name}": c_plasma_connor_hastie,
f"{PlasmaCurrentModel.SAUTER_SCALING.full_name}": c_plasma_sauter,
f"{PlasmaCurrentModel.FIESTA_ST_SCALING.full_name}": c_plasma_fiesta_st,
}

# Create the violin plot
axis.violinplot(data.values(), showextrema=False)

# Create the box plot
axis.boxplot(
data.values(), showfliers=True, showmeans=True, meanline=True, widths=0.3
)

# Scatter plot for each data point
colors = plt.cm.plasma(np.linspace(0, 1, len(data.values())))
for index, (key, value) in enumerate(data.items()):
axis.scatter(1, value, color=colors[index], label=key, alpha=1.0)
axis.legend(loc="upper left", bbox_to_anchor=(-0.9, 1))

# Calculate average, standard deviation, and median
data_values = list(data.values())
avg_density_limit = np.mean(data_values)
std_density_limit = np.std(data_values)
median_density_limit = np.median(data_values)

# Plot average, standard deviation, and median as text
axis.text(
-0.45,
0.15,
rf"Average: {avg_density_limit * 1e-6:.4f}",
transform=axis.transAxes,
fontsize=9,
)
axis.text(
-0.45,
0.1,
rf"Standard Dev: {std_density_limit * 1e-6:.4f}",
transform=axis.transAxes,
fontsize=9,
)
axis.text(
-0.45,
0.05,
rf"Median: {median_density_limit * 1e-6:.4f}",
transform=axis.transAxes,
fontsize=9,
)

axis.set_title("Plasma Current Comparison")
axis.set_ylabel(r"Plasma Current [MA]")
axis.yaxis.set_major_formatter(plt.FuncFormatter(lambda x, _: f"{x * 1e-6:.1f}"))
axis.set_xlim([0.5, 1.5])
axis.set_xticks([])
axis.set_xticklabels([])
axis.set_facecolor("#f0f0f0")


def plot_plasma_pressure_gradient_profiles(axis: plt.Axes, mfile: MFile, scan: int):
# Get the plasma pressure profiles
n_plasma_profile_elements = int(mfile.get("n_plasma_profile_elements", scan=scan))
Expand Down Expand Up @@ -13718,7 +13804,7 @@ def main_plot(
plot_ebw_ecrh_coupling_graph(figs[13].add_subplot(111), m_file, scan)

plot_bootstrap_comparison(figs[14].add_subplot(221), m_file, scan)
PlasmaCurrent.plot_plasma_current_comparison(figs[14].add_subplot(224), m_file, scan)
plot_plasma_current_comparison(figs[14].add_subplot(224), m_file, scan)
plot_h_threshold_comparison(figs[15].add_subplot(224), m_file, scan)
plot_density_limit_comparison(figs[15].add_subplot(221), m_file, scan)
plot_confinement_time_comparison(figs[16].add_subplot(224), m_file, scan)
Expand Down
89 changes: 0 additions & 89 deletions process/models/physics/plasma_current.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,11 +2,9 @@
from enum import IntEnum
from types import DynamicClassAttribute

import matplotlib.pyplot as plt
import numba as nb
import numpy as np

import process.core.io.mfile as mf
from process.core import constants
from process.core import process_output as po
from process.core.exceptions import ProcessValueError
Expand Down Expand Up @@ -399,93 +397,6 @@ def output_plasma_current_models(self) -> None:
"OP ",
)

@staticmethod
def plot_plasma_current_comparison(axis: plt.Axes, mfile: mf.MFile, scan: int):
"""Function to plot a scatter box plot of different plasma current comparisons.

Parameters
----------
axis :
Axis object to plot to.
mfile :
MFILE data object.
scan :
Scan number to use.
"""
c_plasma_peng_analytic = mfile.get("c_plasma_peng_analytic", scan=scan)
c_plasma_peng_double_null = mfile.get("c_plasma_peng_double_null", scan=scan)
c_plasma_cyclindrical = mfile.get("c_plasma_cyclindrical", scan=scan)
c_plasma_ipdg89 = mfile.get("c_plasma_ipdg89", scan=scan)
c_plasma_todd_empirical_i = mfile.get("c_plasma_todd_empirical_i", scan=scan)
c_plasma_todd_empirical_ii = mfile.get("c_plasma_todd_empirical_ii", scan=scan)
c_plasma_connor_hastie = mfile.get("c_plasma_connor_hastie", scan=scan)
c_plasma_sauter = mfile.get("c_plasma_sauter", scan=scan)
c_plasma_fiesta_st = mfile.get("c_plasma_fiesta_st", scan=scan)

# Data for the box plot
data = {
f"{PlasmaCurrentModel.PENG_ANALYTIC_FIT.full_name}": c_plasma_peng_analytic,
f"{PlasmaCurrentModel.PENG_DIVERTOR_SCALING.full_name}": c_plasma_peng_double_null,
f"{PlasmaCurrentModel.ITER_SCALING.full_name}": c_plasma_cyclindrical,
f"{PlasmaCurrentModel.IPDG89_SCALING.full_name}": c_plasma_ipdg89,
f"{PlasmaCurrentModel.TODD_EMPIRICAL_SCALING_I.full_name}": c_plasma_todd_empirical_i,
f"{PlasmaCurrentModel.TODD_EMPIRICAL_SCALING_II.full_name}": c_plasma_todd_empirical_ii,
f"{PlasmaCurrentModel.CONNOR_HASTIE_MODEL.full_name}": c_plasma_connor_hastie,
f"{PlasmaCurrentModel.SAUTER_SCALING.full_name}": c_plasma_sauter,
f"{PlasmaCurrentModel.FIESTA_ST_SCALING.full_name}": c_plasma_fiesta_st,
}

# Create the violin plot
axis.violinplot(data.values(), showextrema=False)

# Create the box plot
axis.boxplot(
data.values(), showfliers=True, showmeans=True, meanline=True, widths=0.3
)

# Scatter plot for each data point
colors = plt.cm.plasma(np.linspace(0, 1, len(data.values())))
for index, (key, value) in enumerate(data.items()):
axis.scatter(1, value, color=colors[index], label=key, alpha=1.0)
axis.legend(loc="upper left", bbox_to_anchor=(-0.9, 1))

# Calculate average, standard deviation, and median
data_values = list(data.values())
avg_density_limit = np.mean(data_values)
std_density_limit = np.std(data_values)
median_density_limit = np.median(data_values)

# Plot average, standard deviation, and median as text
axis.text(
-0.45,
0.15,
rf"Average: {avg_density_limit * 1e-6:.4f}",
transform=axis.transAxes,
fontsize=9,
)
axis.text(
-0.45,
0.1,
rf"Standard Dev: {std_density_limit * 1e-6:.4f}",
transform=axis.transAxes,
fontsize=9,
)
axis.text(
-0.45,
0.05,
rf"Median: {median_density_limit * 1e-6:.4f}",
transform=axis.transAxes,
fontsize=9,
)

axis.set_title("Plasma Current Comparison")
axis.set_ylabel(r"Plasma Current [MA]")
axis.yaxis.set_major_formatter(plt.FuncFormatter(lambda x, _: f"{x * 1e-6:.1f}"))
axis.set_xlim([0.5, 1.5])
axis.set_xticks([])
axis.set_xticklabels([])
axis.set_facecolor("#f0f0f0")

@staticmethod
def calculate_cyclindrical_plasma_current(
rminor: float, rmajor: float, q95: float, b_plasma_toroidal_on_axis: float
Expand Down
Loading