7.3. Phase space sample#
Import Python libraries
import matplotlib.pyplot as plt
import sympy as sp
from ampform.kinematics.phasespace import (
Kallen,
Kibble,
compute_third_mandelstam,
is_within_phasespace,
)
from matplotlib.colors import LogNorm
from polarimetry.data import (
create_mass_symbol_mapping,
create_phase_space_filter,
generate_meshgrid_sample,
generate_phasespace_sample,
)
from polarimetry.io import display_doit, display_latex, mute_jax_warnings
from polarimetry.lhcb import load_model_builder
from polarimetry.lhcb.particle import load_particles
from polarimetry.plot import use_mpl_latex_fonts
mute_jax_warnings()
7.3.1. Definition#
See also
AmpForm’s Kinematics page.
Show code cell source
m0, mi, mj, mk = sp.symbols("m0 m_(i:k)", nonnegative=True)
σi, σj, σk = sp.symbols("sigma_(i:k)", nonnegative=True)
is_within_phasespace(σi, σj, m0, mi, mj, mk)
\[\begin{split}\displaystyle \begin{cases} 1 & \text{for}\: \phi\left(\sigma_{i}, \sigma_{j}\right) \leq 0 \\\text{NaN} & \text{otherwise} \end{cases}\end{split}\]
Show code cell source
display_doit(Kibble(σi, σj, σk, m0, mi, mj, mk))
\[\begin{split}\displaystyle \begin{array}{rcl}
\phi\left(\sigma_{i}, \sigma_{j}\right) &=& \lambda\left(\lambda\left(\sigma_{j}, m_{j}^{2}, m_{0}^{2}\right), \lambda\left(\sigma_{k}, m_{k}^{2}, m_{0}^{2}\right), \lambda\left(\sigma_{i}, m_{i}^{2}, m_{0}^{2}\right)\right) \\
\end{array}\end{split}\]
Show code cell source
display_doit(Kallen(*sp.symbols("x:z")))
\[\begin{split}\displaystyle \begin{array}{rcl}
\lambda\left(x, y, z\right) &=& x^{2} - 2 x y - 2 x z + y^{2} - 2 y z + z^{2} \\
\end{array}\end{split}\]
m1, m2, m3 = sp.symbols("m1:4")
display_latex({σk: compute_third_mandelstam(σi, σj, m0, m1, m2, m3)})
\[\begin{split}\displaystyle \begin{array}{rcl}
\sigma_{k} &=& m_{0}^{2} + m_{1}^{2} + m_{2}^{2} + m_{3}^{2} - \sigma_{i} - \sigma_{j} \\
\end{array}\end{split}\]
7.3.2. Visualization#
Show code cell source
model_choice = 0
model_file = "../../data/model-definitions.yaml"
particles = load_particles("../../data/particle-definitions.yaml")
amplitude_builder = load_model_builder(model_file, particles, model_choice)
decay = amplitude_builder.decay
display_latex(create_mass_symbol_mapping(decay))
\[\begin{split}\displaystyle \begin{array}{rcl}
m_{0} &=& 2.28646 \\
m_{1} &=& 0.938272046 \\
m_{2} &=& 0.13957018 \\
m_{3} &=& 0.49367700000000003 \\
\end{array}\end{split}\]
Show code cell source
def plot_phsp_boundary(ax, x_mandelstam: int, y_mandelstam: int):
ax.set_xticks([])
ax.set_xlabel(Rf"$\sigma_{x_mandelstam}$ [GeV$^2$]")
ax.set_ylabel(Rf"$\sigma_{y_mandelstam}$ [GeV$^2$]")
phsp = generate_meshgrid_sample(decay, resolution, x_mandelstam, y_mandelstam)
phsp_filter = create_phase_space_filter(
decay, x_mandelstam, y_mandelstam, outside_value=0
)
xyz = (
phsp[f"sigma{x_mandelstam}"],
phsp[f"sigma{y_mandelstam}"],
phsp_filter(phsp),
)
ax.contourf(*xyz, colors="lightgray", levels=[0.5, 1])
ax.contour(*xyz, colors="black", levels=[0.5, 1])
%config InlineBackend.figure_formats = ['svg']
resolution = 500
use_mpl_latex_fonts()
plt.rc("font", size=18)
fig, ax = plt.subplots(figsize=(6, 5), tight_layout=True)
plot_phsp_boundary(ax, 1, 2)
fig.savefig("../_images/phase-space-boundary.svg")
plt.show()
Show code cell source
def plot_phsp_distribution(ax, x_mandelstam: int, y_mandelstam: int):
ax.set_xticks([])
ax.set_xlabel(Rf"$\sigma_{x_mandelstam}$")
ax.set_ylabel(Rf"$\sigma_{y_mandelstam}$")
ax.hist2d(
phsp[f"sigma{x_mandelstam}"],
phsp[f"sigma{y_mandelstam}"],
bins=500,
norm=LogNorm(),
)
%config InlineBackend.figure_formats = ['png']
phsp = generate_phasespace_sample(decay, n_events=10_000_000, seed=0)
fig, axes = plt.subplots(ncols=3, dpi=200, figsize=(15, 4.9), tight_layout=True)
plot_phsp_distribution(axes[0], 1, 2)
plot_phsp_distribution(axes[1], 2, 3)
plot_phsp_distribution(axes[2], 3, 1)
fig.tight_layout()
plt.show()
