From 7a6fcc81e8dd65fa16825ee7fbbcce57b8f4d927 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 23 Apr 2025 15:52:59 +0200 Subject: [PATCH 01/13] plotting routine Af(A) plane for pdf_set_nuclear added --- pdfplotter/pdf_set_nuclear.py | 110 +++++++++++++++++++++++++++++++++- 1 file changed, 109 insertions(+), 1 deletion(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index b042d26..44ce72f 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -6,9 +6,14 @@ import numpy as np import numpy.typing as npt import pandas as pd +import matplotlib.pyplot as plt +import sympy as sp +import math +import seaborn as sns +from pdfplotter import util from pdfplotter.pdf_set import PDFSet - +from pdfplotter import flavors class NuclearPDFSet(PDFSet): @@ -153,3 +158,106 @@ def get( raise ValueError(f"Multiple PDFSets found for Z = {Z}") else: return pdf_set.iloc[0]["pdf_set"] + + + def plot_Af_plane( + self, x_vals: npt.ArrayLike, observables:Sequence[sp.Basic], Q: float = 2, colors:Sequence=[], logx:bool=True, name:str| None = None, + plot_unc:bool=False, plotPRatio:bool=False, plot_args:dict | None=None, legend_args:dict |None =None, label_args:dict={"fontsize":15}, title_args:dict|None={"fontsize":15}, + notation_args:dict={"fontsize":12, + "xy":(0.97,0.96)} + )->Sequence[plt.axes]: + + + + my_sets=self.pdf_sets + my_data={} + + for x in x_vals: + if x not in self.get(A=my_sets["A"][0]).x_values: + raise ValueError(f"Chosen x value {x} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialising") + + if colors==[]: + colors=sns.color_palette("viridis", n_colors=len(x_vals)) + + if len(colors)!=len(x_vals): + raise ValueError("No. of colors must match no. of x-values") + + for obs in observables: + data_obs={} + list_x=[] + list_central=[] + list_unc1=[] + list_unc2=[] + list_A=[] + for A in my_sets["A"]: + for x in self.get(A=A).x_values: + list_x.append(x) + list_central.append(self.get(A=A).get_central(x=x,Q=Q,observable=obs)) + unc1=self.get(A=A).get_uncertainties(x=x,Q=Q,observable=obs)[0] + unc2=self.get(A=A).get_uncertainties(x=x,Q=Q,observable=obs)[1] + if math.isnan(unc1): + list_unc1.append(self.get(A=A).get_central(x=x,Q=Q,observable=obs)) + else: + list_unc1.append(unc1) + if math.isnan(unc2): + list_unc2.append(self.get(A=A).get_central(x=x,Q=Q,observable=obs)) + else: + list_unc2.append(unc2) + i=0 + while i < len(self.get(A=A).x_values): + list_A.append(A) + i+=1 + + data_obs["A"]=list_A + data_obs["x"]=list_x + data_obs["central"]=list_central + data_obs["unc1"]=list_unc1 + data_obs["unc2"]=list_unc2 + + dataframe_obs=pd.DataFrame(data_obs) + my_data[obs]=dataframe_obs + + fig, axs =plt.subplots(1,len(observables), figsize=(9*len(observables),5)) + + for i,(obs, ax_i) in enumerate(zip(observables, axs)): + + if not plotPRatio: + for x,col in zip(x_vals,colors): + ax_i.plot(my_data[obs].query(f"x=={x}")["A"], my_data[obs].query(f"x=={x}")["central"],color=col,label=f"x={x}",**plot_args) + + if plot_unc: + ax_i.fill_between(my_data[obs].query(f"x=={x}")["A"],my_data[obs].query(f"x=={x}")["unc1"], my_data[obs].query(f"x=={x}")["unc2"],color=col,alpha=0.3) + + else: + for x,col in zip(x_vals,colors): + ax_i.plot(my_data[obs].query(f"x=={x}")["A"], np.array(my_data[obs].query(f"x=={x}")["central"])/np.array(my_data[obs].query(f"A=={1} & x=={x}")["central"]),color=col,label=f"x={x}",**plot_args) + + if plot_unc: + ax_i.fill_between(my_data[obs].query(f"x=={x}")["A"],my_data[obs].query(f"x=={x}")["unc1"], my_data[obs].query(f"x=={x}")["unc2"],color=col,alpha=0.3) + ax_i.set_xlabel("$A$", **label_args) + #ax_i.set_title(f"${util.to_str(obs)}$") + ax_i.annotate( + f"${util.to_str(obs,Q=Q)}$", + xycoords="axes fraction", + va="top", + ha="right", + bbox=dict( + facecolor=(1,1,1), + edgecolor=(0.8,0.8,0.8), + lw=0.9, + boxstyle="round, pad=0.2" + ), **notation_args + ) + if logx: + ax_i.set_xscale("log") + axs[-1].legend(**legend_args) + if plotPRatio: + axs[0].set_ylabel("$R_x(A)$", **label_args) + else: + axs[0].set_ylabel("$f_x(A)$", **label_args) + if name: + fig.suptitle(f"{name}, Q= {Q} GeV", **title_args) + else: + fig.suptitle(f"Q= {Q} GeV", **title_args) + + return axs \ No newline at end of file From 0a1609a04abbd2a837743d7ee22dc9b5811db5ea Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Thu, 24 Apr 2025 14:46:49 +0200 Subject: [PATCH 02/13] Add possibility to change keyword arguments for labels, title, central and uncertainties. Format code with black formatter --- pdfplotter/pdf_set_nuclear.py | 327 ++++++++++++++++++++++++++-------- 1 file changed, 252 insertions(+), 75 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 44ce72f..d81e95c 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -1,7 +1,7 @@ from __future__ import annotations from itertools import zip_longest -from typing import Sequence +from typing import Sequence, Any, Literal import numpy as np import numpy.typing as npt @@ -11,10 +11,12 @@ import math import seaborn as sns from pdfplotter import util +from pdfplotter.util import update_kwargs from pdfplotter.pdf_set import PDFSet from pdfplotter import flavors + class NuclearPDFSet(PDFSet): _pdf_sets: pd.DataFrame @@ -159,105 +161,280 @@ def get( else: return pdf_set.iloc[0]["pdf_set"] - def plot_Af_plane( - self, x_vals: npt.ArrayLike, observables:Sequence[sp.Basic], Q: float = 2, colors:Sequence=[], logx:bool=True, name:str| None = None, - plot_unc:bool=False, plotPRatio:bool=False, plot_args:dict | None=None, legend_args:dict |None =None, label_args:dict={"fontsize":15}, title_args:dict|None={"fontsize":15}, - notation_args:dict={"fontsize":12, - "xy":(0.97,0.96)} - )->Sequence[plt.axes]: - + self, + ax: plt.Axes | npt.NDArray[plt.Axes], # pyright: ignore[reportInvalidTypeForm] + x_vals: float | list[float], + observables: ( + sp.Basic + | npt.NDArray[sp.Basic] # pyright: ignore[reportInvalidTypeForm] + | list[sp.Basic] + ), + Q: float = 2, + colors: list = [], + logx: bool = True, + title: str | list[str] | None = None, + plot_unc: bool = False, + plot_ratio: bool = False, + pdf_label: Literal["ylabel", "annotate"] = "annotate", + kwargs_theory: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_legend: dict[str, Any] = {}, + kwargs_xlabel: dict[str, Any] = {}, + kwargs_ylabel: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_title: dict[str, Any] = {}, + kwargs_notation: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_uncertainty: dict[str, Any] | list[dict[str, Any] | None] = {}, + ) -> None: + """Plot nuclear PDFs in the A-f plane for different values of x. + + Parameters + ---------- + ax : matplotlib.axes.Axes | numpy.ndarray[matplotlib.axes.Axes] + The axes to plot on. + x_vals : float | list[float] + The x values to plot for. + observables : sympy.Basic | numpy.ndarray[sympy.Basic] | list[sympy.Basic] + The observables to plot. + Q : float, optional + The scale at which to plot the PDFs, by default 2. + colors : list, optional + The colors to use for the different x values, by default [] and seaborn color palette viridis is used if == []. + logx : bool, optional + If True, use a logarithmic scale for the x axis, by default True. + title : str | list[str] | None, optional + The title of the plot, by default None. If a list is passed, the titles are set for each subplot. If a single string is passed, it is set for the first subplot. + plot_unc : bool, optional + If True, plot the uncertainties, by default False. + plot_ratio : bool, optional + If True, plot the ratio of the PDFs to the Proon PDF, by default False. + pdf_label : str, optional + The label for the PDF, by default "annotate". If "ylabel", the label is set as the y-axis label. If "annotate", the label is set as an annotation in the top right corner of the plot. + kwargs_theory : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the plot function for the central PDF, by default {}. + kwargs_legend : dict[str, Any], optional + The keyword arguments to pass to the legend function, by default {}. + kwargs_xlabel : dict[str, Any], optional + The keyword arguments to pass to the xlabel function, by default {}. + kwargs_ylabel : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the ylabel function, by default {}. + kwargs_title : dict[str, Any], optional + The keyword arguments to pass to the title function, by default {}. + kwargs_notation : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the annotation function, by default {}. + kwargs_uncertainty : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the fill_between function for the uncertainties, by default {}. + + """ + my_sets = self.pdf_sets + my_data = {} - my_sets=self.pdf_sets - my_data={} + if not isinstance(x_vals, list): + x_vals = [x_vals] for x in x_vals: if x not in self.get(A=my_sets["A"][0]).x_values: - raise ValueError(f"Chosen x value {x} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialising") + raise ValueError( + f"Chosen x value {x} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialising" + ) - if colors==[]: - colors=sns.color_palette("viridis", n_colors=len(x_vals)) - - if len(colors)!=len(x_vals): + if isinstance(observables, np.ndarray): + observables = list(observables.flatten()) + + if not isinstance(observables, list): + observables = [observables] + + if colors == []: + colors = sns.color_palette("viridis", n_colors=len(x_vals)) + + if len(colors) != len(x_vals): raise ValueError("No. of colors must match no. of x-values") for obs in observables: - data_obs={} - list_x=[] - list_central=[] - list_unc1=[] - list_unc2=[] - list_A=[] + data_obs = {} + list_x = [] + list_central = [] + list_unc1 = [] + list_unc2 = [] + list_A = [] for A in my_sets["A"]: for x in self.get(A=A).x_values: list_x.append(x) - list_central.append(self.get(A=A).get_central(x=x,Q=Q,observable=obs)) - unc1=self.get(A=A).get_uncertainties(x=x,Q=Q,observable=obs)[0] - unc2=self.get(A=A).get_uncertainties(x=x,Q=Q,observable=obs)[1] + list_central.append( + self.get(A=A).get_central(x=x, Q=Q, observable=obs) + ) + unc1 = self.get(A=A).get_uncertainties(x=x, Q=Q, observable=obs)[0] + unc2 = self.get(A=A).get_uncertainties(x=x, Q=Q, observable=obs)[1] if math.isnan(unc1): - list_unc1.append(self.get(A=A).get_central(x=x,Q=Q,observable=obs)) + list_unc1.append( + self.get(A=A).get_central(x=x, Q=Q, observable=obs) + ) else: list_unc1.append(unc1) if math.isnan(unc2): - list_unc2.append(self.get(A=A).get_central(x=x,Q=Q,observable=obs)) + list_unc2.append( + self.get(A=A).get_central(x=x, Q=Q, observable=obs) + ) else: list_unc2.append(unc2) - i=0 + i = 0 while i < len(self.get(A=A).x_values): list_A.append(A) - i+=1 - - data_obs["A"]=list_A - data_obs["x"]=list_x - data_obs["central"]=list_central - data_obs["unc1"]=list_unc1 - data_obs["unc2"]=list_unc2 - - dataframe_obs=pd.DataFrame(data_obs) - my_data[obs]=dataframe_obs - - fig, axs =plt.subplots(1,len(observables), figsize=(9*len(observables),5)) - - for i,(obs, ax_i) in enumerate(zip(observables, axs)): - - if not plotPRatio: - for x,col in zip(x_vals,colors): - ax_i.plot(my_data[obs].query(f"x=={x}")["A"], my_data[obs].query(f"x=={x}")["central"],color=col,label=f"x={x}",**plot_args) + i += 1 + + data_obs["A"] = list_A + data_obs["x"] = list_x + data_obs["central"] = list_central + data_obs["unc1"] = list_unc1 + data_obs["unc2"] = list_unc2 + + dataframe_obs = pd.DataFrame(data_obs) + my_data[obs] = dataframe_obs + + # fig, ax = plt.subplots(1, len(observables), figsize=(9 * len(observables), 5)) + + if not isinstance(ax, np.ndarray): + ax = np.array([ax]) + + for m, (obs_m, ax_m) in enumerate(zip(observables, ax.flat)): + ax_m: plt.Axes + + if not plot_ratio: + for j, (x, col) in enumerate(zip(x_vals, colors)): + kwargs_default = { + "color": col, + "label": f"x={x}", + "linewidth": 1.5, + } + kwargs = update_kwargs( + kwargs_default, + kwargs_theory, + i=j, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x}")["A"], + my_data[obs_m].query(f"x=={x}")["central"], + **kwargs, + ) if plot_unc: - ax_i.fill_between(my_data[obs].query(f"x=={x}")["A"],my_data[obs].query(f"x=={x}")["unc1"], my_data[obs].query(f"x=={x}")["unc2"],color=col,alpha=0.3) - + kwargs_uncertainty_default = { + "color": col, + "alpha": 0.3, + } + kwargs_u = update_kwargs( + kwargs_uncertainty_default, + kwargs_uncertainty, + i=j, + ) + ax_m.fill_between( + my_data[obs_m].query(f"x=={x}")["A"], + my_data[obs_m].query(f"x=={x}")["unc1"], + my_data[obs_m].query(f"x=={x}")["unc2"], + **kwargs_u, + ) + else: - for x,col in zip(x_vals,colors): - ax_i.plot(my_data[obs].query(f"x=={x}")["A"], np.array(my_data[obs].query(f"x=={x}")["central"])/np.array(my_data[obs].query(f"A=={1} & x=={x}")["central"]),color=col,label=f"x={x}",**plot_args) - + for j, (x, col) in enumerate(zip(x_vals, colors)): + kwargs_default = { + "color": col, + "label": f"x={x}", + "linewidth": 1.5, + } + kwargs = update_kwargs( + kwargs_default, + kwargs_theory, + i=j, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x}")["A"], + np.array(my_data[obs_m].query(f"x=={x}")["central"]) + / np.array(my_data[obs_m].query(f"A=={1} & x=={x}")["central"]), + **kwargs, + ) + if plot_unc: - ax_i.fill_between(my_data[obs].query(f"x=={x}")["A"],my_data[obs].query(f"x=={x}")["unc1"], my_data[obs].query(f"x=={x}")["unc2"],color=col,alpha=0.3) - ax_i.set_xlabel("$A$", **label_args) - #ax_i.set_title(f"${util.to_str(obs)}$") - ax_i.annotate( - f"${util.to_str(obs,Q=Q)}$", - xycoords="axes fraction", - va="top", - ha="right", - bbox=dict( - facecolor=(1,1,1), - edgecolor=(0.8,0.8,0.8), - lw=0.9, - boxstyle="round, pad=0.2" - ), **notation_args + kwargs_uncertainty_default = { + "color": col, + "alpha": 0.3, + } + kwargs_uncertainty = update_kwargs( + kwargs_uncertainty_default, + kwargs_uncertainty, + i=j, + ) + ax_m.fill_between( + my_data[obs_m].query(f"x=={x}")["A"], + my_data[obs_m].query(f"x=={x}")["unc1"], + my_data[obs_m].query(f"x=={x}")["unc2"], + **kwargs_uncertainty, + ) + kwargs_xlabel_default = { + "fontsize": 14, + "xlabel": "$A$", + } + kwargs_x = update_kwargs( + kwargs_xlabel_default, + kwargs_xlabel, + ) + + ax_m.set_xlabel(**kwargs_x) + + kwargs_legend_default = { + "fontsize": 12, + "loc": "upper left", + "bbox_to_anchor": (1, 1), + "frameon": False, + } + kwargs_legend = update_kwargs( + kwargs_legend_default, + kwargs_legend, ) + + if m == len(ax.flat) - 1: + ax_m.legend(**kwargs_legend) + + if pdf_label == "annotate": + kwargs_notation_default = { + "fontsize": 12, + "xy": (0.97, 0.96), + "xycoords": "axes fraction", + "va": "top", + "ha": "right", + "bbox": dict( + facecolor=(1, 1, 1), + edgecolor=(0.8, 0.8, 0.8), + lw=0.9, + boxstyle="round, pad=0.2", + ), + } + kwargs_n = update_kwargs(kwargs_notation_default, kwargs_notation, i=m) + ax_m.annotate(f"${util.to_str(obs_m, Q=Q)}$", **kwargs_n) + + if pdf_label == "ylabel": + kwargs_ylabel_default = { + "fontsize": 14, + "ylabel": f"${util.to_str(obs_m, Q=Q)}$", + } + kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel, i=m) + ax_m.set_ylabel(**kwargs_y) + if logx: - ax_i.set_xscale("log") - axs[-1].legend(**legend_args) - if plotPRatio: - axs[0].set_ylabel("$R_x(A)$", **label_args) - else: - axs[0].set_ylabel("$f_x(A)$", **label_args) - if name: - fig.suptitle(f"{name}, Q= {Q} GeV", **title_args) - else: - fig.suptitle(f"Q= {Q} GeV", **title_args) + ax_m.set_xscale("log") + + if title: + kwargs_title_default = { + "fontsize": 20, + "y": 1.05, + "loc": "center" + } + kwargs_title = update_kwargs( + kwargs_title_default, + kwargs_title, + ) - return axs \ No newline at end of file + if isinstance(title, list): + for k, title_k in enumerate(title): + ax.flatten()[k].set_title(f"{title_k}", **kwargs_title) + else: + ax.flatten()[0].set_title(f"{title}", **kwargs_title) From 62a07a67f6779e041ab2a103825c2e118ba57bdc Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Fri, 25 Apr 2025 11:06:59 +0200 Subject: [PATCH 03/13] Fix the plotting of uncertainities in Af plotting function --- pdfplotter/pdf_set_nuclear.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index d81e95c..1027f7c 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -358,17 +358,18 @@ def plot_Af_plane( "color": col, "alpha": 0.3, } - kwargs_uncertainty = update_kwargs( + kwargs_u = update_kwargs( kwargs_uncertainty_default, kwargs_uncertainty, i=j, ) ax_m.fill_between( my_data[obs_m].query(f"x=={x}")["A"], - my_data[obs_m].query(f"x=={x}")["unc1"], - my_data[obs_m].query(f"x=={x}")["unc2"], - **kwargs_uncertainty, + my_data[obs_m].query(f"x=={x}")["unc1"]/ np.array(my_data[obs_m].query(f"A=={1} & x=={x}")["central"]), + my_data[obs_m].query(f"x=={x}")["unc2"]/ np.array(my_data[obs_m].query(f"A=={1} & x=={x}")["central"]), + **kwargs_u, ) + kwargs_xlabel_default = { "fontsize": 14, "xlabel": "$A$", From 2e415feb815aabf16c16912301540af31141236b Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Tue, 29 Apr 2025 12:54:10 +0200 Subject: [PATCH 04/13] Change to standardised naming of the input parameters. Add Q2 option. Change default colour palette from seaborn to matplotlib tab. Add the option to plot uncertainty edges and pass uncertaintiy_edges_kwargs. Change the name of the plot function to plot_A_dep_vs_A_xfixed --- pdfplotter/pdf_set_nuclear.py | 172 ++++++++++++++++++++++++---------- 1 file changed, 123 insertions(+), 49 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 1027f7c..7e1a5cf 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -161,16 +161,17 @@ def get( else: return pdf_set.iloc[0]["pdf_set"] - def plot_Af_plane( + def plot_A_dep_vs_A_xfixed( self, ax: plt.Axes | npt.NDArray[plt.Axes], # pyright: ignore[reportInvalidTypeForm] - x_vals: float | list[float], + x: float | list[float], observables: ( sp.Basic | npt.NDArray[sp.Basic] # pyright: ignore[reportInvalidTypeForm] - | list[sp.Basic] + | list[sp.Basic] ), - Q: float = 2, + Q: float | None = None, + Q2: float | None = None, colors: list = [], logx: bool = True, title: str | list[str] | None = None, @@ -184,22 +185,25 @@ def plot_Af_plane( kwargs_title: dict[str, Any] = {}, kwargs_notation: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_uncertainty: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_uncertainty_edges: dict[str, Any] | list[dict[str, Any] | None] = {}, ) -> None: """Plot nuclear PDFs in the A-f plane for different values of x. - + Parameters ---------- ax : matplotlib.axes.Axes | numpy.ndarray[matplotlib.axes.Axes] The axes to plot on. - x_vals : float | list[float] + x : float | list[float] The x values to plot for. observables : sympy.Basic | numpy.ndarray[sympy.Basic] | list[sympy.Basic] The observables to plot. Q : float, optional - The scale at which to plot the PDFs, by default 2. + The scale at which to plot the PDFs + Q2 : float, optional + The Q^2 scale at which to plot the PDFs. Either Q or Q2 has to be passed. colors : list, optional - The colors to use for the different x values, by default [] and seaborn color palette viridis is used if == []. - logx : bool, optional + The colors to use for the different x values, by default [], tab color palette is used if == []. + logx : bool, optional If True, use a logarithmic scale for the x axis, by default True. title : str | list[str] | None, optional The title of the plot, by default None. If a list is passed, the titles are set for each subplot. If a single string is passed, it is set for the first subplot. @@ -209,7 +213,7 @@ def plot_Af_plane( If True, plot the ratio of the PDFs to the Proon PDF, by default False. pdf_label : str, optional The label for the PDF, by default "annotate". If "ylabel", the label is set as the y-axis label. If "annotate", the label is set as an annotation in the top right corner of the plot. - kwargs_theory : dict[str, Any] | list[dict[str, Any] | None], optional + kwargs_theory : dict[str, Any] | list[dict[str, Any] | None], optional The keyword arguments to pass to the plot function for the central PDF, by default {}. kwargs_legend : dict[str, Any], optional The keyword arguments to pass to the legend function, by default {}. @@ -223,21 +227,41 @@ def plot_Af_plane( The keyword arguments to pass to the annotation function, by default {}. kwargs_uncertainty : dict[str, Any] | list[dict[str, Any] | None], optional The keyword arguments to pass to the fill_between function for the uncertainties, by default {}. - + kwargs_uncertainty_edges : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the plot function for the uncertainty edges, by default {}. """ my_sets = self.pdf_sets my_data = {} - if not isinstance(x_vals, list): - x_vals = [x_vals] + if not isinstance(x, list): + x = [x] - for x in x_vals: - if x not in self.get(A=my_sets["A"][0]).x_values: + for x_i in x: + if x_i not in self.get(A=my_sets["A"][0]).x_values: raise ValueError( - f"Chosen x value {x} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialising" + f"Chosen x value {x_i} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialization" ) + if Q is None and Q2 is None: + raise ValueError("Please pass either `Q` or `Q2`") + + elif Q is not None and Q2 is not None: + raise ValueError("Please pass either `Q` or `Q2`, not both") + + elif Q is not None: + if Q not in self.get(A=my_sets["A"][0]).Q_values and Q not in np.sqrt(np.array(self.get(A=my_sets["A"][0]).Q2_values)): + raise ValueError( + f"Chosen Q value {Q} was not used for defining nuclear pdf set. \n Please choose Q that was used in initialization" + ) + else: + if ( + Q2 not in self.get(A=my_sets["A"][0]).Q2_values + and Q2 not in np.array(self.get(A=my_sets["A"][0]).Q_values) ** 2 + ): + raise ValueError( + f"Chosen Q2 value {Q2} was not used for defining nuclear pdf set. \n Please choose Q2 that was used in initialization" + ) if isinstance(observables, np.ndarray): observables = list(observables.flatten()) @@ -245,9 +269,9 @@ def plot_Af_plane( observables = [observables] if colors == []: - colors = sns.color_palette("viridis", n_colors=len(x_vals)) + colors = plt.get_cmap("tab10", lut=len(x)).colors - if len(colors) != len(x_vals): + if len(colors) != len(x): raise ValueError("No. of colors must match no. of x-values") for obs in observables: @@ -258,22 +282,26 @@ def plot_Af_plane( list_unc2 = [] list_A = [] for A in my_sets["A"]: - for x in self.get(A=A).x_values: - list_x.append(x) + for x_i in self.get(A=A).x_values: + list_x.append(x_i) list_central.append( - self.get(A=A).get_central(x=x, Q=Q, observable=obs) + self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs) ) - unc1 = self.get(A=A).get_uncertainties(x=x, Q=Q, observable=obs)[0] - unc2 = self.get(A=A).get_uncertainties(x=x, Q=Q, observable=obs)[1] + unc1 = self.get(A=A).get_uncertainties( + x=x_i, Q=Q, Q2=Q2, observable=obs + )[0] + unc2 = self.get(A=A).get_uncertainties( + x=x_i, Q=Q, Q2=Q2, observable=obs + )[1] if math.isnan(unc1): list_unc1.append( - self.get(A=A).get_central(x=x, Q=Q, observable=obs) + self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs) ) else: list_unc1.append(unc1) if math.isnan(unc2): list_unc2.append( - self.get(A=A).get_central(x=x, Q=Q, observable=obs) + self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs) ) else: list_unc2.append(unc2) @@ -300,10 +328,10 @@ def plot_Af_plane( ax_m: plt.Axes if not plot_ratio: - for j, (x, col) in enumerate(zip(x_vals, colors)): + for j, (x_j, col) in enumerate(zip(x, colors)): kwargs_default = { "color": col, - "label": f"x={x}", + "label": f"x={x_j}", "linewidth": 1.5, } kwargs = update_kwargs( @@ -312,8 +340,8 @@ def plot_Af_plane( i=j, ) ax_m.plot( - my_data[obs_m].query(f"x=={x}")["A"], - my_data[obs_m].query(f"x=={x}")["central"], + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["central"], **kwargs, ) @@ -328,17 +356,35 @@ def plot_Af_plane( i=j, ) ax_m.fill_between( - my_data[obs_m].query(f"x=={x}")["A"], - my_data[obs_m].query(f"x=={x}")["unc1"], - my_data[obs_m].query(f"x=={x}")["unc2"], + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc1"], + my_data[obs_m].query(f"x=={x_j}")["unc2"], **kwargs_u, ) - + kwargs_uncertainty_edges_default = { + "color": col, + "linestyle": "solid", + } + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + i=j, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc1"], + **kwargs_u_e, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc2"], + **kwargs_u_e, + ) else: - for j, (x, col) in enumerate(zip(x_vals, colors)): + for j, (x_j, col) in enumerate(zip(x, colors)): kwargs_default = { "color": col, - "label": f"x={x}", + "label": f"x={x_j}", "linewidth": 1.5, } kwargs = update_kwargs( @@ -347,9 +393,11 @@ def plot_Af_plane( i=j, ) ax_m.plot( - my_data[obs_m].query(f"x=={x}")["A"], - np.array(my_data[obs_m].query(f"x=={x}")["central"]) - / np.array(my_data[obs_m].query(f"A=={1} & x=={x}")["central"]), + my_data[obs_m].query(f"x=={x_j}")["A"], + np.array(my_data[obs_m].query(f"x=={x_j}")["central"]) + / np.array( + my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] + ), **kwargs, ) @@ -364,12 +412,42 @@ def plot_Af_plane( i=j, ) ax_m.fill_between( - my_data[obs_m].query(f"x=={x}")["A"], - my_data[obs_m].query(f"x=={x}")["unc1"]/ np.array(my_data[obs_m].query(f"A=={1} & x=={x}")["central"]), - my_data[obs_m].query(f"x=={x}")["unc2"]/ np.array(my_data[obs_m].query(f"A=={1} & x=={x}")["central"]), + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc1"] + / np.array( + my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] + ), + my_data[obs_m].query(f"x=={x_j}")["unc2"] + / np.array( + my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] + ), **kwargs_u, ) - + kwargs_uncertainty_edges_default = { + "color": col, + "linestyle": "solid", + "alpha": 1, + } + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + i=j, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc1"]/ np.array( + my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] + ), + **kwargs_u_e, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc2"]/ np.array( + my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] + ), + **kwargs_u_e, + ) + kwargs_xlabel_default = { "fontsize": 14, "xlabel": "$A$", @@ -410,12 +488,12 @@ def plot_Af_plane( ), } kwargs_n = update_kwargs(kwargs_notation_default, kwargs_notation, i=m) - ax_m.annotate(f"${util.to_str(obs_m, Q=Q)}$", **kwargs_n) + ax_m.annotate(f"${util.to_str(obs_m, Q=Q,Q2=Q2)}$", **kwargs_n) if pdf_label == "ylabel": kwargs_ylabel_default = { "fontsize": 14, - "ylabel": f"${util.to_str(obs_m, Q=Q)}$", + "ylabel": f"${util.to_str(obs_m,Q=Q,Q2=Q2)}$", } kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel, i=m) ax_m.set_ylabel(**kwargs_y) @@ -424,11 +502,7 @@ def plot_Af_plane( ax_m.set_xscale("log") if title: - kwargs_title_default = { - "fontsize": 20, - "y": 1.05, - "loc": "center" - } + kwargs_title_default = {"fontsize": 20, "y": 1.05, "loc": "center"} kwargs_title = update_kwargs( kwargs_title_default, kwargs_title, From c52676595ee1b8272ea83e52d2f7790a328484b9 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 30 Apr 2025 17:06:40 +0200 Subject: [PATCH 05/13] Add horizontal A lines. Add possibility to have a legend or a colorbar for x. --- pdfplotter/pdf_set_nuclear.py | 212 ++++++++++++++++++++++++---------- pdfplotter/util.py | 6 + 2 files changed, 158 insertions(+), 60 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 7e1a5cf..54b3b25 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -1,20 +1,25 @@ from __future__ import annotations -from itertools import zip_longest +from itertools import zip_longest, cycle from typing import Sequence, Any, Literal import numpy as np import numpy.typing as npt import pandas as pd import matplotlib.pyplot as plt +import matplotlib.cm as cm +import matplotlib.colors as cls +from matplotlib import ticker as mticker import sympy as sp import math import seaborn as sns +from sympy import rotations from pdfplotter import util -from pdfplotter.util import update_kwargs +from pdfplotter.util import tick_formatter_exp_to_int, update_kwargs +from pdfplotter.util import log_tick_formatter from pdfplotter.pdf_set import PDFSet -from pdfplotter import flavors +from pdfplotter import elements class NuclearPDFSet(PDFSet): @@ -159,7 +164,9 @@ def get( elif pdf_set.shape[0] > 1: raise ValueError(f"Multiple PDFSets found for Z = {Z}") else: - return pdf_set.iloc[0]["pdf_set"] + return pdf_set.iloc[0][ + "pdf_set" + ] # pyright: ignore[reportInvalidTypeForm] def plot_A_dep_vs_A_xfixed( self, @@ -172,18 +179,23 @@ def plot_A_dep_vs_A_xfixed( ), Q: float | None = None, Q2: float | None = None, - colors: list = [], + A_lines: float | list[float] | None = None, + colors: list[str] | str | cycle = [], + cmap: str = "viridis", + labels_Bjx: Literal["colorbar", "legend","none"] = "legend", logx: bool = True, title: str | list[str] | None = None, plot_unc: bool = False, plot_ratio: bool = False, - pdf_label: Literal["ylabel", "annotate"] = "annotate", + pdf_label: Literal["ylabel", "annotate"] | None = "annotate", + plot_legend: bool = True, + legend_labels: Literal["PDFSet", "x", "Both"] = "x", kwargs_theory: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_legend: dict[str, Any] = {}, kwargs_xlabel: dict[str, Any] = {}, kwargs_ylabel: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_title: dict[str, Any] = {}, - kwargs_notation: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_annotate: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_uncertainty: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_uncertainty_edges: dict[str, Any] | list[dict[str, Any] | None] = {}, ) -> None: @@ -212,7 +224,7 @@ def plot_A_dep_vs_A_xfixed( plot_ratio : bool, optional If True, plot the ratio of the PDFs to the Proon PDF, by default False. pdf_label : str, optional - The label for the PDF, by default "annotate". If "ylabel", the label is set as the y-axis label. If "annotate", the label is set as an annotation in the top right corner of the plot. + The label for the PDF, by default "annotate". If "ylabel", the label is set as the y-axis label. If "annotate", the label is set as an anannotate in the top right corner of the plot. kwargs_theory : dict[str, Any] | list[dict[str, Any] | None], optional The keyword arguments to pass to the plot function for the central PDF, by default {}. kwargs_legend : dict[str, Any], optional @@ -223,8 +235,8 @@ def plot_A_dep_vs_A_xfixed( The keyword arguments to pass to the ylabel function, by default {}. kwargs_title : dict[str, Any], optional The keyword arguments to pass to the title function, by default {}. - kwargs_notation : dict[str, Any] | list[dict[str, Any] | None], optional - The keyword arguments to pass to the annotation function, by default {}. + kwargs_annotate : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the anannotate function, by default {}. kwargs_uncertainty : dict[str, Any] | list[dict[str, Any] | None], optional The keyword arguments to pass to the fill_between function for the uncertainties, by default {}. kwargs_uncertainty_edges : dict[str, Any] | list[dict[str, Any] | None], optional @@ -250,7 +262,9 @@ def plot_A_dep_vs_A_xfixed( raise ValueError("Please pass either `Q` or `Q2`, not both") elif Q is not None: - if Q not in self.get(A=my_sets["A"][0]).Q_values and Q not in np.sqrt(np.array(self.get(A=my_sets["A"][0]).Q2_values)): + if Q not in self.get(A=my_sets["A"][0]).Q_values and Q not in np.sqrt( + np.array(self.get(A=my_sets["A"][0]).Q2_values) + ): raise ValueError( f"Chosen Q value {Q} was not used for defining nuclear pdf set. \n Please choose Q that was used in initialization" ) @@ -268,11 +282,12 @@ def plot_A_dep_vs_A_xfixed( if not isinstance(observables, list): observables = [observables] - if colors == []: - colors = plt.get_cmap("tab10", lut=len(x)).colors + if isinstance(colors, str): + colors = len(x) * [colors] - if len(colors) != len(x): - raise ValueError("No. of colors must match no. of x-values") + elif isinstance(colors, list) and colors != []: + if len(colors) != len(x): + raise ValueError("No. of colors must match no. of x-values") for obs in observables: data_obs = {} @@ -327,12 +342,26 @@ def plot_A_dep_vs_A_xfixed( for m, (obs_m, ax_m) in enumerate(zip(observables, ax.flat)): ax_m: plt.Axes + if labels_Bjx == "legend": + if colors == []: + colors = cycle(plt.rcParams["axes.prop_cycle"].by_key()["color"]) + else: + if colors == []: + colors = cycle( + [cm.get_cmap(cmap, lut=len(x))(i) for i in range(len(x))] + ) + if not plot_ratio: - for j, (x_j, col) in enumerate(zip(x, colors)): + for j, x_j in enumerate(x): + if isinstance(colors, str): + col = colors + elif isinstance(colors, list): + col = colors[j] + else: + col = next(colors) kwargs_default = { "color": col, "label": f"x={x_j}", - "linewidth": 1.5, } kwargs = update_kwargs( kwargs_default, @@ -348,7 +377,7 @@ def plot_A_dep_vs_A_xfixed( if plot_unc: kwargs_uncertainty_default = { "color": col, - "alpha": 0.3, + "alpha": 0.2, } kwargs_u = update_kwargs( kwargs_uncertainty_default, @@ -363,13 +392,19 @@ def plot_A_dep_vs_A_xfixed( ) kwargs_uncertainty_edges_default = { "color": col, - "linestyle": "solid", + "linewidth": 0.5, } - kwargs_u_e = update_kwargs( - kwargs_uncertainty_edges_default, - kwargs_uncertainty_edges, - i=j, - ) + if isinstance(kwargs_uncertainty_edges, list): + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + i=j, + ) + else: + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + ) ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], my_data[obs_m].query(f"x=={x_j}")["unc1"], @@ -382,10 +417,15 @@ def plot_A_dep_vs_A_xfixed( ) else: for j, (x_j, col) in enumerate(zip(x, colors)): + if isinstance(colors, str): + col = colors + elif isinstance(colors, list): + col = colors[j] + else: + col = next(colors) kwargs_default = { "color": col, "label": f"x={x_j}", - "linewidth": 1.5, } kwargs = update_kwargs( kwargs_default, @@ -404,7 +444,7 @@ def plot_A_dep_vs_A_xfixed( if plot_unc: kwargs_uncertainty_default = { "color": col, - "alpha": 0.3, + "alpha": 0.2, } kwargs_u = update_kwargs( kwargs_uncertainty_default, @@ -425,31 +465,63 @@ def plot_A_dep_vs_A_xfixed( ) kwargs_uncertainty_edges_default = { "color": col, - "linestyle": "solid", - "alpha": 1, + "linewidth": 0.5, } - kwargs_u_e = update_kwargs( - kwargs_uncertainty_edges_default, - kwargs_uncertainty_edges, - i=j, - ) + if isinstance(kwargs_uncertainty_edges, list): + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + i=j, + ) + else: + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + ) ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], - my_data[obs_m].query(f"x=={x_j}")["unc1"]/ np.array( + my_data[obs_m].query(f"x=={x_j}")["unc1"] + / np.array( my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] ), **kwargs_u_e, ) ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], - my_data[obs_m].query(f"x=={x_j}")["unc2"]/ np.array( + my_data[obs_m].query(f"x=={x_j}")["unc2"] + / np.array( my_data[obs_m].query(f"A=={1} & x=={x_j}")["central"] ), **kwargs_u_e, ) + if logx: + ax_m.set_xscale("log") + if A_lines is not None: + if not isinstance(A_lines, list): + A_lines = [A_lines] + + for A_line in A_lines: + ax_m.axvline( + x=A_line, + color="black", + linestyle="--", + linewidth=0.8, + ) + #transform = ax_m.get_xaxis_transform() + #ax_m.annotate( + # f"{elements.element_to_str(A=A_line,long=True)}", + # xy=(A_line, 0.03), + # xycoords=transform, + # rotation=90 + #) + ax_m.set_xticks(A_lines,labels=[f"{A_line} {elements.element_to_str(A=A_line,long=True)}" for A_line in A_lines],ha="left",rotation=-30) + ax_m.xaxis.set_tick_params(which="minor", size=0) + else: + ax_m.xaxis.set_major_formatter( + mticker.FuncFormatter(tick_formatter_exp_to_int) + ) kwargs_xlabel_default = { - "fontsize": 14, "xlabel": "$A$", } kwargs_x = update_kwargs( @@ -459,22 +531,36 @@ def plot_A_dep_vs_A_xfixed( ax_m.set_xlabel(**kwargs_x) - kwargs_legend_default = { - "fontsize": 12, - "loc": "upper left", - "bbox_to_anchor": (1, 1), - "frameon": False, - } - kwargs_legend = update_kwargs( - kwargs_legend_default, - kwargs_legend, - ) + if labels_Bjx == "colorbar": + if m == len(ax.flat) - 1: + norm = cls.LogNorm(vmin=min(x), vmax=(max(x))) + sm = cm.ScalarMappable(cmap=cmap, norm=norm) + sm.set_array([]) + if len(x) <= 7: + ticks = list(x) + else: + ticks = np.linspace(min(x), max(x), 7) + cbar = plt.colorbar(sm, ax=ax_m, ticks=ticks) + #cbar.ax.yaxis.set_major_locator(mticker.FixedLocator(ticks)) + cbar.ax.set_yticklabels([f"{(x_i)}" for x_i in ticks]) + cbar.set_label("$x$", labelpad=15) + elif labels_Bjx == "legend": + if m == len(ax.flat) - 1: + if plot_legend: + kwargs_legend_default = { + "loc": "upper left", + "bbox_to_anchor": (1, 1), + "frameon": False, + } + kwargs_legend = update_kwargs( + kwargs_legend_default, + kwargs_legend, + ) - if m == len(ax.flat) - 1: - ax_m.legend(**kwargs_legend) + ax_m.legend(**kwargs_legend) if pdf_label == "annotate": - kwargs_notation_default = { + kwargs_annotate_default = { "fontsize": 12, "xy": (0.97, 0.96), "xycoords": "axes fraction", @@ -487,29 +573,35 @@ def plot_A_dep_vs_A_xfixed( boxstyle="round, pad=0.2", ), } - kwargs_n = update_kwargs(kwargs_notation_default, kwargs_notation, i=m) + kwargs_n = update_kwargs(kwargs_annotate_default, kwargs_annotate, i=m) ax_m.annotate(f"${util.to_str(obs_m, Q=Q,Q2=Q2)}$", **kwargs_n) if pdf_label == "ylabel": kwargs_ylabel_default = { - "fontsize": 14, "ylabel": f"${util.to_str(obs_m,Q=Q,Q2=Q2)}$", } kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel, i=m) ax_m.set_ylabel(**kwargs_y) - if logx: - ax_m.set_xscale("log") + + if title: - kwargs_title_default = {"fontsize": 20, "y": 1.05, "loc": "center"} - kwargs_title = update_kwargs( - kwargs_title_default, - kwargs_title, - ) if isinstance(title, list): + for k, title_k in enumerate(title): - ax.flatten()[k].set_title(f"{title_k}", **kwargs_title) + kwargs_title_default = {"y": 1.05, "loc": "center","label":f"{title_k}"} + kwargs_title = update_kwargs( + kwargs_title_default, + kwargs_title, + i=k + ) + ax.flatten()[k].set_title(**kwargs_title) else: - ax.flatten()[0].set_title(f"{title}", **kwargs_title) + kwargs_title_default = {"y": 1.05, "loc": "center","label":f"{title}"} + kwargs_title = update_kwargs( + kwargs_title_default, + kwargs_title, + ) + ax.flatten()[0].set_title(**kwargs_title) diff --git a/pdfplotter/util.py b/pdfplotter/util.py index 16903fc..756b8cc 100644 --- a/pdfplotter/util.py +++ b/pdfplotter/util.py @@ -201,3 +201,9 @@ def update_kwargs( return kwargs else: raise ValueError("kwargs_user must be a dict or a list") + +def log_tick_formatter(val,pos=None): + return r"$10^{{{:.0f}}}$".format(val) + +def tick_formatter_exp_to_int(val,pos=None): + return r"${{{:.0f}}}$".format(val) \ No newline at end of file From e4c8806f13f4982ba47d9a829d71edd2dbeedc0f Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Fri, 2 May 2025 15:38:04 +0200 Subject: [PATCH 06/13] Add a formatting for colorbar y axis. --- pdfplotter/pdf_set_nuclear.py | 105 ++++++++++++++++++++++++---------- pdfplotter/util.py | 8 ++- 2 files changed, 81 insertions(+), 32 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 54b3b25..86e633c 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -10,12 +10,19 @@ import matplotlib.cm as cm import matplotlib.colors as cls from matplotlib import ticker as mticker + +# from mpl_toolkits.axes_grid1.inset_locator import inset_axes import sympy as sp import math import seaborn as sns from sympy import rotations from pdfplotter import util -from pdfplotter.util import tick_formatter_exp_to_int, update_kwargs +from pdfplotter.util import ( + log_tick_formatter_sci, + log_tick_formatter_sci, + update_kwargs, + tick_formatter_exp_to_int, +) from pdfplotter.util import log_tick_formatter from pdfplotter.pdf_set import PDFSet @@ -182,17 +189,17 @@ def plot_A_dep_vs_A_xfixed( A_lines: float | list[float] | None = None, colors: list[str] | str | cycle = [], cmap: str = "viridis", - labels_Bjx: Literal["colorbar", "legend","none"] = "legend", + labels_Bjx: Literal["colorbar", "legend", "none"] = "legend", logx: bool = True, title: str | list[str] | None = None, plot_unc: bool = False, plot_ratio: bool = False, - pdf_label: Literal["ylabel", "annotate"] | None = "annotate", + pdf_label: Literal["ylabel", "annotate", "none"] | None = "annotate", plot_legend: bool = True, legend_labels: Literal["PDFSet", "x", "Both"] = "x", kwargs_theory: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_legend: dict[str, Any] = {}, - kwargs_xlabel: dict[str, Any] = {}, + kwargs_xlabel: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_ylabel: dict[str, Any] | list[dict[str, Any] | None] = {}, kwargs_title: dict[str, Any] = {}, kwargs_annotate: dict[str, Any] | list[dict[str, Any] | None] = {}, @@ -508,14 +515,22 @@ def plot_A_dep_vs_A_xfixed( linestyle="--", linewidth=0.8, ) - #transform = ax_m.get_xaxis_transform() - #ax_m.annotate( + # transform = ax_m.get_xaxis_transform() + # ax_m.annotate( # f"{elements.element_to_str(A=A_line,long=True)}", # xy=(A_line, 0.03), # xycoords=transform, # rotation=90 - #) - ax_m.set_xticks(A_lines,labels=[f"{A_line} {elements.element_to_str(A=A_line,long=True)}" for A_line in A_lines],ha="left",rotation=-30) + # ) + ax_m.set_xticks( + A_lines, + labels=[ + f"{A_line} {elements.element_to_str(A=A_line,long=True)}" + for A_line in A_lines + ], + ha="left", + rotation=-30, + ) ax_m.xaxis.set_tick_params(which="minor", size=0) else: ax_m.xaxis.set_major_formatter( @@ -524,26 +539,52 @@ def plot_A_dep_vs_A_xfixed( kwargs_xlabel_default = { "xlabel": "$A$", } - kwargs_x = update_kwargs( - kwargs_xlabel_default, - kwargs_xlabel, - ) + if isinstance(kwargs_xlabel, list): + kwargs_x = update_kwargs(kwargs_xlabel_default, kwargs_xlabel, i=m) + else: + kwargs_x = update_kwargs( + kwargs_xlabel_default, + kwargs_xlabel, + ) ax_m.set_xlabel(**kwargs_x) if labels_Bjx == "colorbar": if m == len(ax.flat) - 1: - norm = cls.LogNorm(vmin=min(x), vmax=(max(x))) + norm = cls.LogNorm(vmin=(min(x)), vmax=(max(x))) sm = cm.ScalarMappable(cmap=cmap, norm=norm) sm.set_array([]) + cax = ax_m.inset_axes([1.05, 0, 0.07, 1]) + cax.set_yscale("log") + cax.yaxis.set_tick_params(which="minor", size=0) + cbar = plt.colorbar(sm, cax=cax) + cbar.set_label("$x$", labelpad=15) if len(x) <= 7: - ticks = list(x) + ticks = x + cbar.ax.set_yticks( + ticks, + labels=[ + f"{x_i:.0e}" if x_i != round(x_i, 1) else f"{x_i}" + for x_i in ticks + ], + ) else: - ticks = np.linspace(min(x), max(x), 7) - cbar = plt.colorbar(sm, ax=ax_m, ticks=ticks) - #cbar.ax.yaxis.set_major_locator(mticker.FixedLocator(ticks)) - cbar.ax.set_yticklabels([f"{(x_i)}" for x_i in ticks]) - cbar.set_label("$x$", labelpad=15) + ticks = np.arange( + int(np.min(np.log10(x))), int(np.max(np.log10(x))) + 1 + ) + + cbar.ax.set_yticks( + 10.0 ** (ticks), + labels=[ + ( + f"{10.0**(x_i):.1e}" + if x_i not in [-1, 0] + else (f"{0.1}" if x_i == -1 else f"{1}") + ) + for x_i in ticks + ], + ) + elif labels_Bjx == "legend": if m == len(ax.flat) - 1: if plot_legend: @@ -580,28 +621,30 @@ def plot_A_dep_vs_A_xfixed( kwargs_ylabel_default = { "ylabel": f"${util.to_str(obs_m,Q=Q,Q2=Q2)}$", } - kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel, i=m) + if isinstance(kwargs_ylabel, list): + kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel, i=m) + else: + kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel) ax_m.set_ylabel(**kwargs_y) - - - if title: if isinstance(title, list): for k, title_k in enumerate(title): - kwargs_title_default = {"y": 1.05, "loc": "center","label":f"{title_k}"} + kwargs_title_default = { + "y": 1.05, + "loc": "center", + "label": f"{title_k}", + } kwargs_title = update_kwargs( - kwargs_title_default, - kwargs_title, - i=k + kwargs_title_default, kwargs_title, i=k ) ax.flatten()[k].set_title(**kwargs_title) else: - kwargs_title_default = {"y": 1.05, "loc": "center","label":f"{title}"} + kwargs_title_default = {"y": 1.05, "loc": "center", "label": f"{title}"} kwargs_title = update_kwargs( - kwargs_title_default, - kwargs_title, - ) + kwargs_title_default, + kwargs_title, + ) ax.flatten()[0].set_title(**kwargs_title) diff --git a/pdfplotter/util.py b/pdfplotter/util.py index 756b8cc..9e89e65 100644 --- a/pdfplotter/util.py +++ b/pdfplotter/util.py @@ -206,4 +206,10 @@ def log_tick_formatter(val,pos=None): return r"$10^{{{:.0f}}}$".format(val) def tick_formatter_exp_to_int(val,pos=None): - return r"${{{:.0f}}}$".format(val) \ No newline at end of file + return r"${{{:.0f}}}$".format(val) + +def log_tick_formatter_sci(val,pos=None): + if val!=round(val,1): + return f"${val:.1e}$".format(val) + else: + return f"${val:.1f}$".format(val) From 401238ecd22a02960320c60b59702243fadd51c3 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 21 May 2025 11:49:50 +0200 Subject: [PATCH 07/13] Add plotting routine for A-dep ratio in xslices using an offset. Comment out the A must be greater than one error for construct nuclear PDFs --- pdfplotter/pdf_set_nuclear.py | 410 +++++++++++++++++++++++++++++++++- 1 file changed, 409 insertions(+), 1 deletion(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 86e633c..4e16546 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -602,6 +602,7 @@ def plot_A_dep_vs_A_xfixed( if pdf_label == "annotate": kwargs_annotate_default = { + "text":f"${util.to_str(obs_m, Q=Q,Q2=Q2)}$", "fontsize": 12, "xy": (0.97, 0.96), "xycoords": "axes fraction", @@ -615,7 +616,7 @@ def plot_A_dep_vs_A_xfixed( ), } kwargs_n = update_kwargs(kwargs_annotate_default, kwargs_annotate, i=m) - ax_m.annotate(f"${util.to_str(obs_m, Q=Q,Q2=Q2)}$", **kwargs_n) + ax_m.annotate(**kwargs_n) if pdf_label == "ylabel": kwargs_ylabel_default = { @@ -648,3 +649,410 @@ def plot_A_dep_vs_A_xfixed( kwargs_title, ) ax.flatten()[0].set_title(**kwargs_title) + + + def plot_A_dep_vs_A_xslices( + self, + ax: plt.Axes | npt.NDArray[plt.Axes], # pyright: ignore[reportInvalidTypeForm] + x: float | list[float], + observables: ( + sp.Basic + | npt.NDArray[sp.Basic] # pyright: ignore[reportInvalidTypeForm] + | list[sp.Basic] + ), + Q: float | None = None, + Q2: float | None = None, + A_lines: float | list[float] | None = None, + offset:float=1, + colors: list[str] | str | cycle = [], + unc_conv: Literal["sym", "asym"]="asym", + labels_Bjx: Literal["colorbar", "legend", "none"] = "legend", + logx: bool = True, + title: str | list[str] | None = None, + plot_unc: bool = False, + ratio_to: PDFSet | None =None, + pdf_label: Literal["ylabel", "annotate", "none"] | None = "annotate", + plot_legend: bool = True, + kwargs_theory: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_legend: dict[str, Any] = {}, + kwargs_xlabel: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_ylabel: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_title: dict[str, Any] = {}, + kwargs_annotate: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_uncertainty: dict[str, Any] | list[dict[str, Any] | None] = {}, + kwargs_uncertainty_edges: dict[str, Any] | list[dict[str, Any] | None] = {}, + ) -> None: + """Plot nuclear PDFs in the A-f plane for different values of x. + + Parameters + ---------- + ax : matplotlib.axes.Axes | numpy.ndarray[matplotlib.axes.Axes] + The axes to plot on. + x : float | list[float] + The x values to plot for. + observables : sympy.Basic | numpy.ndarray[sympy.Basic] | list[sympy.Basic] + The observables to plot. + Q : float, optional + The scale at which to plot the PDFs + Q2 : float, optional + The Q^2 scale at which to plot the PDFs. Either Q or Q2 has to be passed. + colors : list, optional + The colors to use for the different x values, by default [], tab color palette is used if == []. + logx : bool, optional + If True, use a logarithmic scale for the x axis, by default True. + title : str | list[str] | None, optional + The title of the plot, by default None. If a list is passed, the titles are set for each subplot. If a single string is passed, it is set for the first subplot. + plot_unc : bool, optional + If True, plot the uncertainties, by default False. + plot_ratio : bool, optional + If True, plot the ratio of the PDFs to the Proon PDF, by default False. + pdf_label : str, optional + The label for the PDF, by default "annotate". If "ylabel", the label is set as the y-axis label. If "annotate", the label is set as an anannotate in the top right corner of the plot. + kwargs_theory : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the plot function for the central PDF, by default {}. + kwargs_legend : dict[str, Any], optional + The keyword arguments to pass to the legend function, by default {}. + kwargs_xlabel : dict[str, Any], optional + The keyword arguments to pass to the xlabel function, by default {}. + kwargs_ylabel : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the ylabel function, by default {}. + kwargs_title : dict[str, Any], optional + The keyword arguments to pass to the title function, by default {}. + kwargs_annotate : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the anannotate function, by default {}. + kwargs_uncertainty : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the fill_between function for the uncertainties, by default {}. + kwargs_uncertainty_edges : dict[str, Any] | list[dict[str, Any] | None], optional + The keyword arguments to pass to the plot function for the uncertainty edges, by default {}. + """ + + my_sets = self.pdf_sets + my_data = {} + + if not isinstance(x, list): + x = [x] + + for x_i in x: + if x_i not in self.get(A=my_sets["A"][0]).x_values: + raise ValueError( + f"Chosen x value {x_i} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialization" + ) + + if Q is None and Q2 is None: + raise ValueError("Please pass either `Q` or `Q2`") + + elif Q is not None and Q2 is not None: + raise ValueError("Please pass either `Q` or `Q2`, not both") + + elif Q is not None: + if Q not in self.get(A=my_sets["A"][0]).Q_values and Q not in np.sqrt( + np.array(self.get(A=my_sets["A"][0]).Q2_values) + ): + raise ValueError( + f"Chosen Q value {Q} was not used for defining nuclear pdf set. \n Please choose Q that was used in initialization" + ) + else: + if ( + Q2 not in self.get(A=my_sets["A"][0]).Q2_values + and Q2 not in np.array(self.get(A=my_sets["A"][0]).Q_values) ** 2 + ): + raise ValueError( + f"Chosen Q2 value {Q2} was not used for defining nuclear pdf set. \n Please choose Q2 that was used in initialization" + ) + if isinstance(observables, np.ndarray): + observables = list(observables.flatten()) + + if not isinstance(observables, list): + observables = [observables] + + if isinstance(colors, str): + colors = len(x) * [colors] + + elif isinstance(colors, list) and colors != []: + if len(colors) != len(x): + raise ValueError("No. of colors must match no. of x-values") + #if ratio_to==None and 1 in self.pdf_sets["A"]: + # ratio_to=self.get(A=1) + for obs in observables: + data_obs = {} + list_x = [] + list_central = [] + list_unc1 = [] + list_unc2 = [] + list_A = [] + for A in my_sets["A"]: + for x_i in self.get(A=A).x_values: + list_x.append(x_i) + list_central.append( + self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to) + ) + unc1 = self.get(A=A).get_uncertainties( + x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to, convention=unc_conv + )[0] + unc2 = self.get(A=A).get_uncertainties( + x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to, convention=unc_conv + )[1] + if math.isnan(unc1): + list_unc1.append( + self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to) + ) + else: + list_unc1.append(unc1) + if math.isnan(unc2): + list_unc2.append( + self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to) + ) + else: + list_unc2.append(unc2) + i = 0 + while i < len(self.get(A=A).x_values): + list_A.append(A) + i += 1 + + data_obs["A"] = list_A + data_obs["x"] = list_x + data_obs["central"] = list_central + data_obs["unc1"] = list_unc1 + data_obs["unc2"] = list_unc2 + + dataframe_obs = pd.DataFrame(data_obs) + my_data[obs] = dataframe_obs + + # fig, ax = plt.subplots(1, len(observables), figsize=(9 * len(observables), 5)) + + if not isinstance(ax, np.ndarray): + ax = np.array([ax]) + + for m, (obs_m, ax_m) in enumerate(zip(observables, ax.flat)): + ax_m: plt.Axes + + if colors == []: + colors = cycle(plt.rcParams["axes.prop_cycle"].by_key()["color"]) + + for j, (x_j, col) in enumerate(zip(x, colors)): + if isinstance(colors, str): + col = colors + elif isinstance(colors, list): + col = colors[j] + else: + col = next(colors) + kwargs_default = { + "color": col, + "label": f"x={x_j}", + } + kwargs = update_kwargs( + kwargs_default, + kwargs_theory, + i=j, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + np.array(my_data[obs_m].query(f"x=={x_j}")["central"]) + #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + +offset*j, + **kwargs, + ) + ax_m.plot( + [min(my_data[obs_m].query(f"x=={x_j}")["A"]),max(my_data[obs_m].query(f"x=={x_j}")["A"])], + [j*offset+1, j*offset+1], + linestyle="--", + color="grey" + ) + ax_m.plot( + [min(my_data[obs_m].query(f"x=={x_j}")["A"]),max(my_data[obs_m].query(f"x=={x_j}")["A"])], + [j*offset+1-0.3, j*offset+1-0.3], + linestyle="--", + color="grey", + linewidth=0.7 + ) + ax_m.plot( + [min(my_data[obs_m].query(f"x=={x_j}")["A"]),max(my_data[obs_m].query(f"x=={x_j}")["A"])], + [j*offset+1+0.3, j*offset+1+0.3], + linestyle="--", + color="grey", + linewidth=0.7 + ) + if plot_unc: + kwargs_uncertainty_default = { + "color": col, + "alpha": 0.2, + } + kwargs_u = update_kwargs( + kwargs_uncertainty_default, + kwargs_uncertainty, + i=j, + ) + ax_m.fill_between( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc1"] + #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + +offset*j, + my_data[obs_m].query(f"x=={x_j}")["unc2"] + #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + +offset*j, + **kwargs_u, + ) + kwargs_uncertainty_edges_default = { + "color": col, + "linewidth": 0.5, + } + if isinstance(kwargs_uncertainty_edges, list): + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + i=j, + ) + else: + kwargs_u_e = update_kwargs( + kwargs_uncertainty_edges_default, + kwargs_uncertainty_edges, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc1"] + #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + +offset*j, + **kwargs_u_e, + ) + ax_m.plot( + my_data[obs_m].query(f"x=={x_j}")["A"], + my_data[obs_m].query(f"x=={x_j}")["unc2"] + #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + +offset*j, + **kwargs_u_e, + ) + if logx: + ax_m.set_xscale("log") + + if A_lines is not None: + if not isinstance(A_lines, list): + A_lines = [A_lines] + + for A_line in A_lines: + ax_m.axvline( + x=A_line, + color="black", + linestyle="--", + linewidth=0.8, + ) + + ax_m.set_xticks( + A_lines, + labels=[ + f"{A_line} {elements.element_to_str(A=A_line,long=True)}" + for A_line in A_lines + ], + ha="left", + rotation=-30, + ) + ax_m.xaxis.set_tick_params(which="minor", size=0) + else: + ax_m.xaxis.set_major_formatter( + mticker.FuncFormatter(tick_formatter_exp_to_int) + ) + kwargs_xlabel_default = { + "xlabel": "$A$", + } + if isinstance(kwargs_xlabel, list): + kwargs_x = update_kwargs(kwargs_xlabel_default, kwargs_xlabel, i=m) + else: + kwargs_x = update_kwargs( + kwargs_xlabel_default, + kwargs_xlabel, + ) + + ax_m.set_xlabel(**kwargs_x) + + if ratio_to !=None: + locs=[] + vals=[] + for j in range(len(x)): + locs+=[j*offset+1-0.3,j*offset+1,j*offset+1+0.3] + vals+=[0.7,1,1.3] + ax_m.set_yticks(locs,vals) + + if labels_Bjx == "legend": + if m == len(ax.flat) - 1: + if plot_legend: + kwargs_legend_default = { + "loc": "upper left", + "bbox_to_anchor": (1, 1), + "frameon": False, + } + kwargs_legend = update_kwargs( + kwargs_legend_default, + kwargs_legend, + ) + + ax_m.legend(**kwargs_legend) + + if pdf_label == "annotate": + if ratio_to !=None: + kwargs_annotate_default = { + "text":f"${util.to_str(obs_m, Q=Q, Q2=Q2, R=True)}$", + "fontsize": 12, + "xy": (0.97, 0.96), + "xycoords": "axes fraction", + "va": "top", + "ha": "right", + "bbox": dict( + facecolor=(1, 1, 1), + edgecolor=(0.8, 0.8, 0.8), + lw=0.9, + boxstyle="round, pad=0.2", + ), + } + else: + kwargs_annotate_default = { + "text":f"${util.to_str(obs_m, Q=Q, Q2=Q2, R=False)}$", + "fontsize": 12, + "xy": (0.97, 0.96), + "xycoords": "axes fraction", + "va": "top", + "ha": "right", + "bbox": dict( + facecolor=(1, 1, 1), + edgecolor=(0.8, 0.8, 0.8), + lw=0.9, + boxstyle="round, pad=0.2", + ), + } + kwargs_n = update_kwargs(kwargs_annotate_default, kwargs_annotate, i=m) + ax_m.annotate(**kwargs_n) + + if pdf_label == "ylabel": + if ratio_to !=None: + kwargs_ylabel_default = { + "ylabel": f"${util.to_str(obs_m,Q=Q,Q2=Q2, R=True)}$", + } + else: + kwargs_ylabel_default = { + "ylabel": f"${util.to_str(obs_m,Q=Q,Q2=Q2, R=False)}$", + } + if isinstance(kwargs_ylabel, list): + kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel, i=m) + else: + kwargs_y = update_kwargs(kwargs_ylabel_default, kwargs_ylabel) + ax_m.set_ylabel(**kwargs_y) + + if title: + + if isinstance(title, list): + + for k, title_k in enumerate(title): + kwargs_title_default = { + "y": 1.05, + "loc": "center", + "label": f"{title_k}", + } + kwargs_title = update_kwargs( + kwargs_title_default, kwargs_title, i=k + ) + ax.flatten()[k].set_title(**kwargs_title) + else: + kwargs_title_default = {"y": 1.05, "loc": "center", "label": f"{title}"} + kwargs_title = update_kwargs( + kwargs_title_default, + kwargs_title, + ) + ax.flatten()[0].set_title(**kwargs_title) \ No newline at end of file From 02b0687ed4f394082d1f1610be90bc1f6a356521 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Thu, 22 May 2025 12:47:16 +0200 Subject: [PATCH 08/13] Fix the offset. Add option to plot ratio to whatever PDF you want --- pdfplotter/pdf_set_nuclear.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 4e16546..fcbc171 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -825,7 +825,7 @@ def plot_A_dep_vs_A_xslices( for m, (obs_m, ax_m) in enumerate(zip(observables, ax.flat)): ax_m: plt.Axes - + ax_m.set_ylim(0,float(len(x))*offset+0.3) if colors == []: colors = cycle(plt.rcParams["axes.prop_cycle"].by_key()["color"]) From 1c40bb86c9db8b8b0ba9c0d667b7862e475d613e Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Thu, 12 Jun 2025 14:33:37 +0200 Subject: [PATCH 09/13] Add formatting and change x_values to x Q values to Q etc, similarly to changes in PDFSet --- pdfplotter/pdf_set_nuclear.py | 142 ++++++++++++++++++++-------------- 1 file changed, 83 insertions(+), 59 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index fcbc171..dea5f8a 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -257,7 +257,7 @@ def plot_A_dep_vs_A_xfixed( x = [x] for x_i in x: - if x_i not in self.get(A=my_sets["A"][0]).x_values: + if x_i not in self.get(A=my_sets["A"][0]).x: raise ValueError( f"Chosen x value {x_i} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialization" ) @@ -269,16 +269,16 @@ def plot_A_dep_vs_A_xfixed( raise ValueError("Please pass either `Q` or `Q2`, not both") elif Q is not None: - if Q not in self.get(A=my_sets["A"][0]).Q_values and Q not in np.sqrt( - np.array(self.get(A=my_sets["A"][0]).Q2_values) + if Q not in self.get(A=my_sets["A"][0]).Q and Q not in np.sqrt( + np.array(self.get(A=my_sets["A"][0]).Q2) ): raise ValueError( f"Chosen Q value {Q} was not used for defining nuclear pdf set. \n Please choose Q that was used in initialization" ) else: if ( - Q2 not in self.get(A=my_sets["A"][0]).Q2_values - and Q2 not in np.array(self.get(A=my_sets["A"][0]).Q_values) ** 2 + Q2 not in self.get(A=my_sets["A"][0]).Q2 + and Q2 not in np.array(self.get(A=my_sets["A"][0]).Q) ** 2 ): raise ValueError( f"Chosen Q2 value {Q2} was not used for defining nuclear pdf set. \n Please choose Q2 that was used in initialization" @@ -304,7 +304,7 @@ def plot_A_dep_vs_A_xfixed( list_unc2 = [] list_A = [] for A in my_sets["A"]: - for x_i in self.get(A=A).x_values: + for x_i in self.get(A=A).x: list_x.append(x_i) list_central.append( self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs) @@ -328,7 +328,7 @@ def plot_A_dep_vs_A_xfixed( else: list_unc2.append(unc2) i = 0 - while i < len(self.get(A=A).x_values): + while i < len(self.get(A=A).x): list_A.append(A) i += 1 @@ -602,7 +602,7 @@ def plot_A_dep_vs_A_xfixed( if pdf_label == "annotate": kwargs_annotate_default = { - "text":f"${util.to_str(obs_m, Q=Q,Q2=Q2)}$", + "text": f"${util.to_str(obs_m, Q=Q,Q2=Q2)}$", "fontsize": 12, "xy": (0.97, 0.96), "xycoords": "axes fraction", @@ -650,7 +650,6 @@ def plot_A_dep_vs_A_xfixed( ) ax.flatten()[0].set_title(**kwargs_title) - def plot_A_dep_vs_A_xslices( self, ax: plt.Axes | npt.NDArray[plt.Axes], # pyright: ignore[reportInvalidTypeForm] @@ -663,14 +662,14 @@ def plot_A_dep_vs_A_xslices( Q: float | None = None, Q2: float | None = None, A_lines: float | list[float] | None = None, - offset:float=1, + offset: float = 1, colors: list[str] | str | cycle = [], - unc_conv: Literal["sym", "asym"]="asym", + unc_conv: Literal["sym", "asym"] = "asym", labels_Bjx: Literal["colorbar", "legend", "none"] = "legend", logx: bool = True, title: str | list[str] | None = None, plot_unc: bool = False, - ratio_to: PDFSet | None =None, + ratio_to: PDFSet | None = None, pdf_label: Literal["ylabel", "annotate", "none"] | None = "annotate", plot_legend: bool = True, kwargs_theory: dict[str, Any] | list[dict[str, Any] | None] = {}, @@ -733,7 +732,7 @@ def plot_A_dep_vs_A_xslices( x = [x] for x_i in x: - if x_i not in self.get(A=my_sets["A"][0]).x_values: + if x_i not in self.get(A=my_sets["A"][0]).x: raise ValueError( f"Chosen x value {x_i} was not used for defining nuclear pdf set. \n Pleas choose x that was used in initialization" ) @@ -745,16 +744,16 @@ def plot_A_dep_vs_A_xslices( raise ValueError("Please pass either `Q` or `Q2`, not both") elif Q is not None: - if Q not in self.get(A=my_sets["A"][0]).Q_values and Q not in np.sqrt( - np.array(self.get(A=my_sets["A"][0]).Q2_values) + if Q not in self.get(A=my_sets["A"][0]).Q and Q not in np.sqrt( + np.array(self.get(A=my_sets["A"][0]).Q2) ): raise ValueError( f"Chosen Q value {Q} was not used for defining nuclear pdf set. \n Please choose Q that was used in initialization" ) else: if ( - Q2 not in self.get(A=my_sets["A"][0]).Q2_values - and Q2 not in np.array(self.get(A=my_sets["A"][0]).Q_values) ** 2 + Q2 not in self.get(A=my_sets["A"][0]).Q2 + and Q2 not in np.array(self.get(A=my_sets["A"][0]).Q) ** 2 ): raise ValueError( f"Chosen Q2 value {Q2} was not used for defining nuclear pdf set. \n Please choose Q2 that was used in initialization" @@ -771,7 +770,7 @@ def plot_A_dep_vs_A_xslices( elif isinstance(colors, list) and colors != []: if len(colors) != len(x): raise ValueError("No. of colors must match no. of x-values") - #if ratio_to==None and 1 in self.pdf_sets["A"]: + # if ratio_to==None and 1 in self.pdf_sets["A"]: # ratio_to=self.get(A=1) for obs in observables: data_obs = {} @@ -781,31 +780,47 @@ def plot_A_dep_vs_A_xslices( list_unc2 = [] list_A = [] for A in my_sets["A"]: - for x_i in self.get(A=A).x_values: + for x_i in self.get(A=A).x: list_x.append(x_i) list_central.append( - self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to) + self.get(A=A).get_central( + x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to + ) ) unc1 = self.get(A=A).get_uncertainties( - x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to, convention=unc_conv + x=x_i, + Q=Q, + Q2=Q2, + observable=obs, + ratio_to=ratio_to, + convention=unc_conv, )[0] unc2 = self.get(A=A).get_uncertainties( - x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to, convention=unc_conv + x=x_i, + Q=Q, + Q2=Q2, + observable=obs, + ratio_to=ratio_to, + convention=unc_conv, )[1] if math.isnan(unc1): list_unc1.append( - self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to) + self.get(A=A).get_central( + x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to + ) ) else: list_unc1.append(unc1) if math.isnan(unc2): list_unc2.append( - self.get(A=A).get_central(x=x_i, Q=Q, Q2=Q2, observable=obs,ratio_to=ratio_to) + self.get(A=A).get_central( + x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to + ) ) else: list_unc2.append(unc2) i = 0 - while i < len(self.get(A=A).x_values): + while i < len(self.get(A=A).x): list_A.append(A) i += 1 @@ -825,9 +840,9 @@ def plot_A_dep_vs_A_xslices( for m, (obs_m, ax_m) in enumerate(zip(observables, ax.flat)): ax_m: plt.Axes - ax_m.set_ylim(0,float(len(x))*offset+0.3) + ax_m.set_ylim(0, float(len(x)) * offset + 0.3) if colors == []: - colors = cycle(plt.rcParams["axes.prop_cycle"].by_key()["color"]) + colors = cycle(plt.rcParams["axes.prop_cycle"].by_key()["color"]) for j, (x_j, col) in enumerate(zip(x, colors)): if isinstance(colors, str): @@ -848,29 +863,38 @@ def plot_A_dep_vs_A_xslices( ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], np.array(my_data[obs_m].query(f"x=={x_j}")["central"]) - #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) - +offset*j, + # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + + offset * j, **kwargs, ) ax_m.plot( - [min(my_data[obs_m].query(f"x=={x_j}")["A"]),max(my_data[obs_m].query(f"x=={x_j}")["A"])], - [j*offset+1, j*offset+1], + [ + min(my_data[obs_m].query(f"x=={x_j}")["A"]), + max(my_data[obs_m].query(f"x=={x_j}")["A"]), + ], + [j * offset + 1, j * offset + 1], linestyle="--", - color="grey" + color="grey", ) ax_m.plot( - [min(my_data[obs_m].query(f"x=={x_j}")["A"]),max(my_data[obs_m].query(f"x=={x_j}")["A"])], - [j*offset+1-0.3, j*offset+1-0.3], + [ + min(my_data[obs_m].query(f"x=={x_j}")["A"]), + max(my_data[obs_m].query(f"x=={x_j}")["A"]), + ], + [j * offset + 1 - 0.3, j * offset + 1 - 0.3], linestyle="--", color="grey", - linewidth=0.7 + linewidth=0.7, ) ax_m.plot( - [min(my_data[obs_m].query(f"x=={x_j}")["A"]),max(my_data[obs_m].query(f"x=={x_j}")["A"])], - [j*offset+1+0.3, j*offset+1+0.3], + [ + min(my_data[obs_m].query(f"x=={x_j}")["A"]), + max(my_data[obs_m].query(f"x=={x_j}")["A"]), + ], + [j * offset + 1 + 0.3, j * offset + 1 + 0.3], linestyle="--", color="grey", - linewidth=0.7 + linewidth=0.7, ) if plot_unc: kwargs_uncertainty_default = { @@ -885,11 +909,11 @@ def plot_A_dep_vs_A_xslices( ax_m.fill_between( my_data[obs_m].query(f"x=={x_j}")["A"], my_data[obs_m].query(f"x=={x_j}")["unc1"] - #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) - +offset*j, + # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + + offset * j, my_data[obs_m].query(f"x=={x_j}")["unc2"] - #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) - +offset*j, + # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + + offset * j, **kwargs_u, ) kwargs_uncertainty_edges_default = { @@ -910,15 +934,15 @@ def plot_A_dep_vs_A_xslices( ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], my_data[obs_m].query(f"x=={x_j}")["unc1"] - #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) - +offset*j, + # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + + offset * j, **kwargs_u_e, ) ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], my_data[obs_m].query(f"x=={x_j}")["unc2"] - #/self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) - +offset*j, + # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + + offset * j, **kwargs_u_e, ) if logx: @@ -962,14 +986,14 @@ def plot_A_dep_vs_A_xslices( ) ax_m.set_xlabel(**kwargs_x) - - if ratio_to !=None: - locs=[] - vals=[] + + if ratio_to != None: + locs = [] + vals = [] for j in range(len(x)): - locs+=[j*offset+1-0.3,j*offset+1,j*offset+1+0.3] - vals+=[0.7,1,1.3] - ax_m.set_yticks(locs,vals) + locs += [j * offset + 1 - 0.3, j * offset + 1, j * offset + 1 + 0.3] + vals += [0.7, 1, 1.3] + ax_m.set_yticks(locs, vals) if labels_Bjx == "legend": if m == len(ax.flat) - 1: @@ -987,9 +1011,9 @@ def plot_A_dep_vs_A_xslices( ax_m.legend(**kwargs_legend) if pdf_label == "annotate": - if ratio_to !=None: + if ratio_to != None: kwargs_annotate_default = { - "text":f"${util.to_str(obs_m, Q=Q, Q2=Q2, R=True)}$", + "text": f"${util.to_str(obs_m, Q=Q, Q2=Q2, R=True)}$", "fontsize": 12, "xy": (0.97, 0.96), "xycoords": "axes fraction", @@ -1004,7 +1028,7 @@ def plot_A_dep_vs_A_xslices( } else: kwargs_annotate_default = { - "text":f"${util.to_str(obs_m, Q=Q, Q2=Q2, R=False)}$", + "text": f"${util.to_str(obs_m, Q=Q, Q2=Q2, R=False)}$", "fontsize": 12, "xy": (0.97, 0.96), "xycoords": "axes fraction", @@ -1016,12 +1040,12 @@ def plot_A_dep_vs_A_xslices( lw=0.9, boxstyle="round, pad=0.2", ), - } + } kwargs_n = update_kwargs(kwargs_annotate_default, kwargs_annotate, i=m) ax_m.annotate(**kwargs_n) if pdf_label == "ylabel": - if ratio_to !=None: + if ratio_to != None: kwargs_ylabel_default = { "ylabel": f"${util.to_str(obs_m,Q=Q,Q2=Q2, R=True)}$", } @@ -1055,4 +1079,4 @@ def plot_A_dep_vs_A_xslices( kwargs_title_default, kwargs_title, ) - ax.flatten()[0].set_title(**kwargs_title) \ No newline at end of file + ax.flatten()[0].set_title(**kwargs_title) From 11f1fe81d7b9d3c1d4528c137642626b65b4c236 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 2 Jul 2025 12:04:15 +0200 Subject: [PATCH 10/13] Fix data type issue during plotting the uncertainties --- pdfplotter/pdf_set_nuclear.py | 7 ++----- 1 file changed, 2 insertions(+), 5 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index dea5f8a..48d69c0 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -863,7 +863,6 @@ def plot_A_dep_vs_A_xslices( ax_m.plot( my_data[obs_m].query(f"x=={x_j}")["A"], np.array(my_data[obs_m].query(f"x=={x_j}")["central"]) - # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + offset * j, **kwargs, ) @@ -908,11 +907,9 @@ def plot_A_dep_vs_A_xslices( ) ax_m.fill_between( my_data[obs_m].query(f"x=={x_j}")["A"], - my_data[obs_m].query(f"x=={x_j}")["unc1"] - # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + np.array(my_data[obs_m].query(f"x=={x_j}")["unc1"], dtype=float) + offset * j, - my_data[obs_m].query(f"x=={x_j}")["unc2"] - # /self.get(A=1).get_central(x=x_j, Q=Q, Q2=Q2, observable=obs_m) + np.array(my_data[obs_m].query(f"x=={x_j}")["unc2"], dtype=float) + offset * j, **kwargs_u, ) From 1b8641326a3ba3d87c1e5bba856b3c46ba638654 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 6 Aug 2025 12:23:08 +0200 Subject: [PATCH 11/13] Add possibiliy to disable plotign the uncertainties for specfic observables/x-values --- pdfplotter/pdf_set_nuclear.py | 14 +++++++++++--- 1 file changed, 11 insertions(+), 3 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 48d69c0..8f80591 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -668,7 +668,7 @@ def plot_A_dep_vs_A_xslices( labels_Bjx: Literal["colorbar", "legend", "none"] = "legend", logx: bool = True, title: str | list[str] | None = None, - plot_unc: bool = False, + plot_unc: bool | list[list[bool]] = False, ratio_to: PDFSet | None = None, pdf_label: Literal["ylabel", "annotate", "none"] | None = "annotate", plot_legend: bool = True, @@ -764,6 +764,9 @@ def plot_A_dep_vs_A_xslices( if not isinstance(observables, list): observables = [observables] + if not isinstance(plot_unc, list): + plot_unc=len(observables)*[plot_unc] + if isinstance(colors, str): colors = len(x) * [colors] @@ -843,8 +846,11 @@ def plot_A_dep_vs_A_xslices( ax_m.set_ylim(0, float(len(x)) * offset + 0.3) if colors == []: colors = cycle(plt.rcParams["axes.prop_cycle"].by_key()["color"]) - for j, (x_j, col) in enumerate(zip(x, colors)): + + if not isinstance(plot_unc[m], list): + plot_unc[m]=len(x)*[plot_unc[m]] + if isinstance(colors, str): col = colors elif isinstance(colors, list): @@ -895,7 +901,9 @@ def plot_A_dep_vs_A_xslices( color="grey", linewidth=0.7, ) - if plot_unc: + + + if plot_unc[m][j]: kwargs_uncertainty_default = { "color": col, "alpha": 0.2, From 1b1a0a044ec0ee87f2e59eb9fe7d56d7c0357bf4 Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 6 Aug 2025 15:52:16 +0200 Subject: [PATCH 12/13] Possibility to change ticks above and below 1 --- pdfplotter/pdf_set_nuclear.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 8f80591..660c7d9 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -663,6 +663,7 @@ def plot_A_dep_vs_A_xslices( Q2: float | None = None, A_lines: float | list[float] | None = None, offset: float = 1, + sub_tick_offset: float = 0.3, colors: list[str] | str | cycle = [], unc_conv: Literal["sym", "asym"] = "asym", labels_Bjx: Literal["colorbar", "legend", "none"] = "legend", @@ -886,7 +887,7 @@ def plot_A_dep_vs_A_xslices( min(my_data[obs_m].query(f"x=={x_j}")["A"]), max(my_data[obs_m].query(f"x=={x_j}")["A"]), ], - [j * offset + 1 - 0.3, j * offset + 1 - 0.3], + [j * offset + 1 - sub_tick_offset, j * offset + 1 - sub_tick_offset], linestyle="--", color="grey", linewidth=0.7, @@ -896,7 +897,7 @@ def plot_A_dep_vs_A_xslices( min(my_data[obs_m].query(f"x=={x_j}")["A"]), max(my_data[obs_m].query(f"x=={x_j}")["A"]), ], - [j * offset + 1 + 0.3, j * offset + 1 + 0.3], + [j * offset + 1 + sub_tick_offset, j * offset + 1 + sub_tick_offset], linestyle="--", color="grey", linewidth=0.7, @@ -996,8 +997,8 @@ def plot_A_dep_vs_A_xslices( locs = [] vals = [] for j in range(len(x)): - locs += [j * offset + 1 - 0.3, j * offset + 1, j * offset + 1 + 0.3] - vals += [0.7, 1, 1.3] + locs += [j * offset + 1 - sub_tick_offset, j * offset + 1, j * offset + 1 + sub_tick_offset] + vals += [1-sub_tick_offset, 1, 1+sub_tick_offset] ax_m.set_yticks(locs, vals) if labels_Bjx == "legend": From da790639d43448e97316ea3d19eb885ddce7791d Mon Sep 17 00:00:00 2001 From: Katrin-Greve Date: Wed, 22 Oct 2025 14:25:31 +0200 Subject: [PATCH 13/13] Fix bug in nuclear modification factor --- pdfplotter/pdf_set_nuclear.py | 21 ++++++++++++++------- 1 file changed, 14 insertions(+), 7 deletions(-) diff --git a/pdfplotter/pdf_set_nuclear.py b/pdfplotter/pdf_set_nuclear.py index 660c7d9..9f65b25 100644 --- a/pdfplotter/pdf_set_nuclear.py +++ b/pdfplotter/pdf_set_nuclear.py @@ -130,6 +130,8 @@ def __init__( pdf_sets_dict, ) + self._confidence_level=confidence_level + @property def pdf_sets(self) -> pd.DataFrame: """The `PDFSet` objects""" @@ -670,7 +672,8 @@ def plot_A_dep_vs_A_xslices( logx: bool = True, title: str | list[str] | None = None, plot_unc: bool | list[list[bool]] = False, - ratio_to: PDFSet | None = None, + #ratio_to: PDFSet | None = None, + ratio_to: str | None=None, pdf_label: Literal["ylabel", "annotate", "none"] | None = "annotate", plot_legend: bool = True, kwargs_theory: dict[str, Any] | list[dict[str, Any] | None] = {}, @@ -783,12 +786,16 @@ def plot_A_dep_vs_A_xslices( list_unc1 = [] list_unc2 = [] list_A = [] - for A in my_sets["A"]: + for A, Z in zip(my_sets["A"], my_sets["Z"]): + if ratio_to: + ratio_to_pdf=PDFSet(x=np.array(x), Q=Q, Q2=Q2, name=ratio_to, A=A, Z=Z, construct_full_nuclear_pdfs=True, confidence_level=self.confidence_level) + else: + ratio_to_pdf=None for x_i in self.get(A=A).x: list_x.append(x_i) list_central.append( self.get(A=A).get_central( - x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to + x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to_pdf ) ) unc1 = self.get(A=A).get_uncertainties( @@ -796,7 +803,7 @@ def plot_A_dep_vs_A_xslices( Q=Q, Q2=Q2, observable=obs, - ratio_to=ratio_to, + ratio_to=ratio_to_pdf, convention=unc_conv, )[0] unc2 = self.get(A=A).get_uncertainties( @@ -804,13 +811,13 @@ def plot_A_dep_vs_A_xslices( Q=Q, Q2=Q2, observable=obs, - ratio_to=ratio_to, + ratio_to=ratio_to_pdf, convention=unc_conv, )[1] if math.isnan(unc1): list_unc1.append( self.get(A=A).get_central( - x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to + x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to_pdf ) ) else: @@ -818,7 +825,7 @@ def plot_A_dep_vs_A_xslices( if math.isnan(unc2): list_unc2.append( self.get(A=A).get_central( - x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to + x=x_i, Q=Q, Q2=Q2, observable=obs, ratio_to=ratio_to_pdf ) ) else: