%%capture
%config Completer.use_jedi = False
%config InlineBackend.figure_formats = ['svg']
import os
STATIC_WEB_PAGE = {"EXECUTE_NB", "READTHEDOCS"}.intersection(os.environ)
# Install on Google Colab
import subprocess
import sys
from IPython import get_ipython
install_packages = "google.colab" in str(get_ipython())
if install_packages:
for package in ["ampform[doc]", "graphviz"]:
subprocess.check_call(
[sys.executable, "-m", "pip", "install", package]
)
{{ run_interactive }}
By default, the dynamic terms in an amplitude model are set to $1$ by the {class}.HelicityAmplitudeBuilder
. The method {meth}.set_dynamics
can then be used to set dynamics lineshapes for specific resonances. The {mod}.dynamics.builder
module provides some tools to set standard lineshapes (see below), but it is also possible to set {doc}custom dynamics </usage/dynamics/custom>
.
The standard lineshapes provided by AmpForm are illustrated below. For more info, have a look at the following pages:
{toctree}
:maxdepth: 2
dynamics/custom
dynamics/analytic-continuation
dynamics/k-matrix
%matplotlib widget
import logging
import warnings
import matplotlib.pyplot as plt
import mpl_interactions.ipyplot as iplt
import numpy as np
import sympy as sp
from IPython.display import Math
import symplot
logging.basicConfig()
logging.getLogger().setLevel(logging.ERROR)
warnings.filterwarnings("ignore")
AmpForm uses Blatt-Weisskopf functions $B_L$ as barrier factors (also called form factors, see {class}.BlattWeisskopfSquared
):
from ampform.dynamics import BlattWeisskopfSquared
L = sp.Symbol("L", integer=True)
z = sp.Symbol("z", real=True)
ff2 = BlattWeisskopfSquared(L, z)
Math(sp.multiline_latex(ff2, ff2.doit(), environment="eqnarray"))
The Blatt-Weisskopf form factor is used to 'dampen' the breakup-momentum at threshold and when going to infinity. A usual choice for $z$ is therefore $z=q^2d^2$ with $q^2$ the {class}.BreakupMomentumSquared
and $d$ the impact parameter (also called meson radius):
from ampform.dynamics import BreakupMomentumSquared
m, m_a, m_b, d = sp.symbols("m, m_a, m_b, d")
s = m ** 2
q_squared = BreakupMomentumSquared(s, m_a, m_b)
ff2 = BlattWeisskopfSquared(L, z=q_squared * d ** 2)
np_blatt_weisskopf, sliders = symplot.prepare_sliders(
plot_symbol=m,
expression=ff2.doit(),
)
np_breakup_momentum = sp.lambdify((m, L, d, m_a, m_b), q_squared.doit())
plot_domain = np.linspace(0.01, 4, 500)
sliders.set_ranges(
L=(0, 8),
m_a=(0, 2, 200),
m_b=(0, 2, 200),
d=(0, 5),
)
sliders.set_values(
L=1,
m_a=0.3,
m_b=0.2,
d=3,
)
fig, ax = plt.subplots(figsize=(8, 5), tight_layout=True)
ax.set_xlabel("$m$")
ax.axhline(0, c="black", linewidth=0.5)
controls = iplt.plot(
plot_domain,
np_blatt_weisskopf,
**sliders,
ylim="auto",
label=R"$B_L^2\left(q(m)\right)$",
linestyle="dotted",
ax=ax,
)
iplt.plot(
plot_domain,
np_breakup_momentum,
controls=controls,
ylim="auto",
label="$q^2(m)$",
linestyle="dashed",
ax=ax,
)
iplt.title(
"Effect of Blatt-Weisskopf factor\n$L={L}, m_a={m_a:.1f}, m_b={m_b:.1f}$",
controls=controls,
)
plt.legend()
plt.show()
if STATIC_WEB_PAGE:
from IPython.display import SVG
output_file = "blatt-weisskopf.svg"
plt.savefig(output_file)
display(SVG(output_file))
AmpForm has two types of relativistic Breit-Wigner functions. Both are compared below ― for more info, see the links to the API.
The 'normal' {func}.relativistic_breit_wigner
looks as follows:
from ampform.dynamics import relativistic_breit_wigner
m, m0, w0 = sp.symbols("m, m0, Gamma0")
rel_bw = relativistic_breit_wigner(s=m ** 2, mass0=m0, gamma0=w0)
rel_bw
The relativistic Breit-Wigner can be adapted slightly, so that its amplitude goes to zero at threshold ($m_0 = m_a + m_b$) and that it becomes normalizable. This is done with {ref}form factors <usage/dynamics:Form factor>
and can be obtained with the function {func}.relativistic_breit_wigner_with_ff
:
from ampform.dynamics import PhaseSpaceFactor # noqa: F401
from ampform.dynamics import BreakupMomentumSquared, ComplexSqrt
def breakup_momentum(s: sp.Symbol, m_a: sp.Symbol, m_b: sp.Symbol) -> sp.Expr:
q_squared = BreakupMomentumSquared(s, m_a, m_b)
return ComplexSqrt(q_squared)
from ampform.dynamics import (
PhaseSpaceFactorAnalytic,
relativistic_breit_wigner_with_ff,
)
rel_bw_with_ff = relativistic_breit_wigner_with_ff(
s=s,
mass0=m0,
gamma0=w0,
m_a=m_a,
m_b=m_b,
angular_momentum=L,
meson_radius=1,
phsp_factor=PhaseSpaceFactorAnalytic,
)
rel_bw_with_ff
Here, $\Gamma(m)$ is the {class}.EnergyDependentWidth
(also called running width or mass-dependent width), defined as:
from ampform.dynamics import EnergyDependentWidth
L = sp.Symbol("L", integer=True)
width = EnergyDependentWidth(
s=s,
mass0=m0,
gamma0=w0,
m_a=m_a,
m_b=m_b,
angular_momentum=L,
meson_radius=1,
phsp_factor=PhaseSpaceFactorAnalytic,
)
Math(sp.multiline_latex(width, width.evaluate(), environment="eqnarray"))
It is possible to choose different formulations for the phase space factor $\rho$, see {doc}/usage/dynamics/analytic-continuation
.
The following shows the effect of {doc}/usage/dynamics/analytic-continuation
a on relativistic Breit-Wigner:
from ampform.dynamics import PhaseSpaceFactorComplex
# Two types of relativistic Breit-Wigners
rel_bw_with_ff = relativistic_breit_wigner_with_ff(
s=s,
mass0=m0,
gamma0=w0,
m_a=m_a,
m_b=m_b,
angular_momentum=L,
meson_radius=d,
phsp_factor=PhaseSpaceFactorComplex,
)
rel_bw_with_ff_ac = relativistic_breit_wigner_with_ff(
s=s,
mass0=m0,
gamma0=w0,
m_a=m_a,
m_b=m_b,
angular_momentum=L,
meson_radius=d,
phsp_factor=PhaseSpaceFactorAnalytic,
)
# Lambdify
np_rel_bw_with_ff, sliders = symplot.prepare_sliders(
plot_symbol=m,
expression=rel_bw_with_ff.doit(),
)
np_rel_bw_with_ff_ac = sp.lambdify(
args=(m, w0, L, d, m0, m_a, m_b),
expr=rel_bw_with_ff_ac.doit(),
)
np_rel_bw = sp.lambdify(
args=(m, w0, L, d, m0, m_a, m_b),
expr=rel_bw.doit(),
)
# Set sliders
plot_domain = np.linspace(start=0, stop=4, num=500)
sliders.set_ranges(
m0=(0, 4, 200),
Gamma0=(0, 1, 100),
L=(0, 8),
m_a=(0, 2, 200),
m_b=(0, 2, 200),
d=(0, 5),
)
sliders.set_values(
m0=1.8,
Gamma0=0.6,
L=1,
m_a=0.6,
m_b=0.6,
d=1,
)
fig, axes = plt.subplots(
nrows=2,
figsize=(8, 6),
sharex=True,
# gridspec_kw={"height_ratios": [1, 1.8]},
)
ax_ff, ax_ac = axes
ax_ac.set_xlabel("$m$")
for ax in axes:
ax.axhline(0, c="gray", linewidth=0.5)
rho_c = PhaseSpaceFactorComplex(m, m_a, m_b)
ax_ff.set_title(
f"'Complex' phase space factor: ${sp.latex(rho_c.evaluate())}$"
)
ax_ac.set_title("Analytic continuation of the phase space factor")
ylim = "auto" # (-0.6, 1.2)
controls = iplt.plot(
plot_domain,
lambda *args, **kwargs: np_rel_bw_with_ff(*args, **kwargs).real,
label="real",
**sliders,
ylim=ylim,
ax=ax_ff,
)
iplt.plot(
plot_domain.real,
lambda *args, **kwargs: np_rel_bw_with_ff(*args, **kwargs).imag,
label="imaginary",
controls=controls,
ylim=ylim,
ax=ax_ff,
)
iplt.plot(
plot_domain.real,
lambda *args, **kwargs: np.abs(np_rel_bw_with_ff(*args, **kwargs)) ** 2,
label="absolute",
controls=controls,
ylim=ylim,
ax=ax_ff,
c="black",
linestyle="dotted",
)
iplt.plot(
plot_domain.real,
lambda *args, **kwargs: np_rel_bw_with_ff_ac(*args, **kwargs).real,
label="real",
controls=controls,
ylim=ylim,
ax=ax_ac,
)
iplt.plot(
plot_domain.real,
lambda *args, **kwargs: np_rel_bw_with_ff_ac(*args, **kwargs).imag,
label="imaginary",
controls=controls,
ylim=ylim,
ax=ax_ac,
)
iplt.plot(
plot_domain.real,
lambda *args, **kwargs: np.abs(np_rel_bw_with_ff_ac(*args, **kwargs)) ** 2,
label="absolute",
controls=controls,
ylim=ylim,
ax=ax_ac,
c="black",
linestyle="dotted",
)
for ax in axes:
iplt.axvline(
controls["m0"],
ax=ax,
c="red",
label=f"${sp.latex(m0)}$",
alpha=0.3,
)
iplt.axvline(
lambda m_a, m_b, **kwargs: m_a + m_b,
controls=controls,
ax=ax,
c="black",
alpha=0.3,
label=f"${sp.latex(m_a)} + {sp.latex(m_b)}$",
)
ax_ac.legend(loc="upper right")
fig.tight_layout()
plt.show()
if STATIC_WEB_PAGE:
from IPython.display import SVG
output_file = "relativistic-breit-wigner-with-form-factor.svg"
plt.savefig(output_file)
display(SVG(output_file))