Phase space factor widget

Phase space factor widget#

Below, we formulate different versions of the phase space factor and investigate their behavior at the real, physical axis.

%matplotlib widget

Hide code cell source

import warnings

import ipywidgets as w
import matplotlib.pyplot as plt
import numpy as np
import sympy as sp
from ipympl.backend_nbagg import Canvas
from IPython.display import display
from matplotlib_inline.backend_inline import set_matplotlib_formats

from ampform.dynamics.phasespace import (
    EqualMassPhaseSpaceFactor,
    PhaseSpaceFactorComplex,
    PhaseSpaceFactorSWave,
)
from ampform.kinematics.phasespace import BreakupMomentumComplex
from ampform.sympy.math import ComplexSqrt
from ampform.sympy.slider import create_slider

set_matplotlib_formats("svg")
warnings.filterwarnings("ignore")


def hide_toolbars(canvas: Canvas) -> None:
    canvas.header_visible = False
    canvas.footer_visible = False
    canvas.toolbar_visible = False

Hide code cell source

m, m1, m2 = sp.symbols("m m1 m2", nonnegative=True)
rho_c = PhaseSpaceFactorComplex(m**2, m1, m2)
rho_cm = PhaseSpaceFactorSWave(m**2, m1, m2)
rho_ac = EqualMassPhaseSpaceFactor(m**2, m1, m2)
q = BreakupMomentumComplex(m**2, m1, m2)
rho_c_func = sp.lambdify((m, m1, m2), rho_c.doit())
rho_ac_func = sp.lambdify((m, m1, m2), rho_ac.doit())
rho_cm_func = sp.lambdify((m, m1, m2), rho_cm.doit())
q_func = sp.lambdify((m, m1, m2), 2 * q.doit())

Hide code cell source

sliders = {str(s): create_slider(s, min=0, max=2, step=0.01) for s in [m1, m2]}
sliders["m1"].value = 0.3
sliders["m2"].value = 0.75
UI = w.VBox(list(sliders.values()))

Hide code cell source

x = np.linspace(0, 3, num=500)
fig, axes = plt.subplots(
    figsize=(8, 5),
    ncols=2,
    nrows=2,
    sharex=True,
    sharey=True,
)
fig.subplots_adjust(bottom=0.08, hspace=0.1, left=0.01, right=0.99, top=1, wspace=0.05)
hide_toolbars(fig.canvas)
fig.patch.set_facecolor("none")
for ax in axes.flatten():
    ax.patch.set_facecolor("none")
    ax.spines["bottom"].set_position("zero")
    ax.spines["left"].set_position("zero")
    ax.spines["right"].set_visible(False)
    ax.spines["top"].set_visible(False)

ax1, ax2, ax3, ax4 = axes.ravel()
for ax in axes[-1]:
    ax.set_xlabel("$m$")
for ax in axes.flatten():
    ax.set_ylim(-0.1, 1.4)
    ax.set_yticks([])
    ax.set_yticks([])

ylim = (-0.1, 1.4)
q_math = ComplexSqrt(sp.Symbol("q^2")) / (8 * sp.pi)
ax1.set_title(f"${sp.latex(q_math)}$", y=0.85)
ax2.set_title(f"${sp.latex(rho_c)}$", y=0.85)
ax3.set_title(R"equal mass $\rho^\mathrm{eq}(m^2)$", y=0.85)
ax4.set_title(R"Chew–Mandelstam $\rho^\mathrm{CM}(m^2)$", y=0.85)

LINES = None


def plot(m1, m2):
    global LINES
    if LINES is None:
        style = dict(alpha=0.7)
        LINES = [
            ax1.plot(x, q_func(x, m1, m2).real, **style, label="real")[0],
            ax1.plot(x, q_func(x, m1, m2).imag, **style, label="imag")[0],
            ax2.plot(x, rho_c_func(x, m1, m2).real, **style, label="real")[0],
            ax2.plot(x, rho_c_func(x, m1, m2).imag, **style, label="imag")[0],
            ax3.plot(x, rho_ac_func(x, m1, m2).real, **style, label="real")[0],
            ax3.plot(x, rho_ac_func(x, m1, m2).imag, **style, label="imag")[0],
            ax4.plot(x, rho_cm_func(x, m1, m2).real, **style, label="real")[0],
            ax4.plot(x, rho_cm_func(x, m1, m2).imag, **style, label="imag")[0],
        ]
    else:
        LINES[0].set_ydata(q_func(x, m1, m2).real)
        LINES[1].set_ydata(q_func(x, m1, m2).imag)
        LINES[2].set_ydata(rho_c_func(x, m1, m2).real)
        LINES[3].set_ydata(rho_c_func(x, m1, m2).imag)
        LINES[4].set_ydata(rho_ac_func(x, m1, m2).real)
        LINES[5].set_ydata(rho_ac_func(x, m1, m2).imag)
        LINES[6].set_ydata(rho_cm_func(x, m1, m2).real)
        LINES[7].set_ydata(rho_cm_func(x, m1, m2).imag)


output = w.interactive_output(plot, controls=sliders)
ax2.legend(loc="upper right")
plt.show()
display(UI, output)