Analytic continuation#

Note

Improvements to analytic continuation in AmpForm are currently being developed in Chew-Mandelstam and Analyticity.

Analytic continuation allows one to handle resonances just below threshold (\(m_0 < m_a + m_b\) in Eq. (4)). In practice, this entails using a specific function for \(\rho\) in Eq. (2).

Definitions#

Three usual choices for \(\rho\) are the following:

Hide code cell source
%config InlineBackend.figure_formats = ['svg']

import warnings

import sympy as sp
from IPython.display import Math

warnings.filterwarnings("ignore")

1) Break-up momentum#

The sqrt() or ComplexSqrt of BreakupMomentumSquared:

from ampform.dynamics import BreakupMomentumSquared

s, m_a, m_b = sp.symbols("s, m_a, m_b", nonnegative=True)
q_squared = BreakupMomentumSquared(s, m_a, m_b)
Math(sp.multiline_latex(q_squared, q_squared.doit(), environment="eqnarray"))
\[\displaystyle \begin{eqnarray} q^2\left(s\right) & = & \frac{\left(s - \left(m_{a} - m_{b}\right)^{2}\right) \left(s - \left(m_{a} + m_{b}\right)^{2}\right)}{4 s} \end{eqnarray}\]

2) ‘Normal’ phase space factor#

The ‘normal’ PhaseSpaceFactor (the denominator makes the difference to (2)!):

from ampform.dynamics import PhaseSpaceFactor

rho = PhaseSpaceFactor(s, m_a, m_b)
Math(sp.multiline_latex(rho, rho.evaluate(), environment="eqnarray"))
\[\displaystyle \begin{eqnarray} \rho\left(s\right) & = & \frac{2 \sqrt{q^2\left(s\right)}}{\sqrt{s}} \end{eqnarray}\]

3) ‘Complex’ phase space factor#

A PhaseSpaceFactorComplex that uses ComplexSqrt:

from ampform.dynamics import PhaseSpaceFactorComplex

rho_c = PhaseSpaceFactorComplex(s, m_a, m_b)
Math(sp.multiline_latex(rho_c, rho_c.evaluate(), environment="eqnarray"))
\[\displaystyle \begin{eqnarray} \rho^\mathrm{c}\left(s\right) & = & \frac{2 \sqrt[\mathrm{c}]{q^2\left(s\right)}}{\sqrt{s}} \end{eqnarray}\]

4) ‘Analytic continuation’ of the phase space factor#

The following ‘case-by-case’ analytic continuation for decay products with an equal mass, EqualMassPhaseSpaceFactor:

from ampform.dynamics import EqualMassPhaseSpaceFactor

rho_ac = EqualMassPhaseSpaceFactor(s, m_a, m_b)
Math(sp.multiline_latex(rho_ac, rho_ac.evaluate(), environment="eqnarray"))
\[\begin{split}\displaystyle \begin{eqnarray} \rho^\mathrm{eq}\left(s\right) & = & \begin{cases} \frac{i \log{\left(\left|{\frac{\hat{\rho}\left(s\right) + 1}{\hat{\rho}\left(s\right) - 1}}\right| \right)} \hat{\rho}\left(s\right)}{\pi} + \hat{\rho}\left(s\right) & \text{for}\: s > \left(m_{a} + m_{b}\right)^{2} \\\frac{2 i \operatorname{atan}{\left(\frac{1}{\hat{\rho}\left(s\right)} \right)} \hat{\rho}\left(s\right)}{\pi} & \text{otherwise} \end{cases} \end{eqnarray}\end{split}\]

with

Hide code cell source
from ampform.dynamics import PhaseSpaceFactorAbs

rho_hat = PhaseSpaceFactorAbs(s, m_a, m_b)
Math(sp.multiline_latex(rho_hat, rho_hat.evaluate(), environment="eqnarray"))
\[\displaystyle \begin{eqnarray} \hat{\rho}\left(s\right) & = & \frac{2 \sqrt{\left|{q^2\left(s\right)}\right|}}{\sqrt{s}} \end{eqnarray}\]

(Mind the absolute value.)

5) Chew-Mandelstam for \(S\)-waves#

A PhaseSpaceFactorSWave that uses chew_mandelstam_s_wave():

from ampform.dynamics import PhaseSpaceFactorSWave

rho_cm = PhaseSpaceFactorSWave(s, m_a, m_b)
Math(sp.multiline_latex(rho_cm, rho_cm.evaluate(), environment="eqnarray"))
\[\displaystyle \begin{eqnarray} \rho^\mathrm{CM}\left(s\right) & = &- \frac{i \left(\frac{2 \sqrt[\mathrm{c}]{q^2\left(s\right)}}{\sqrt{s}} \log{\left(\frac{m_{a}^{2} + m_{b}^{2} + 2 \sqrt{s} \sqrt[\mathrm{c}]{q^2\left(s\right)} - s}{2 m_{a} m_{b}} \right)} - \left(m_{a}^{2} - m_{b}^{2}\right) \left(- \frac{1}{\left(m_{a} + m_{b}\right)^{2}} + \frac{1}{s}\right) \log{\left(\frac{m_{a}}{m_{b}} \right)}\right)}{\pi} \end{eqnarray}\]

Visualization#

%matplotlib widget
Hide code cell source
import matplotlib.pyplot as plt
import mpl_interactions.ipyplot as iplt
import numpy as np

import symplot
from ampform.sympy.math import ComplexSqrt

m = sp.Symbol("m", nonnegative=True)
rho_c = PhaseSpaceFactorComplex(m**2, m_a, m_b)
rho_cm = PhaseSpaceFactorSWave(m**2, m_a, m_b)
rho_ac = EqualMassPhaseSpaceFactor(m**2, m_a, m_b)
np_rho_c, sliders = symplot.prepare_sliders(plot_symbol=m, expression=rho_c.doit())
np_rho_ac = sp.lambdify((m, m_a, m_b), rho_ac.doit())
np_rho_cm = sp.lambdify((m, m_a, m_b), rho_cm.doit())
np_breakup_momentum = sp.lambdify(
    (m, m_a, m_b),
    2 * ComplexSqrt(q_squared.subs(s, m**2).doit()),
)
plot_domain = np.linspace(0, 3, 500)
sliders.set_ranges(
    m_a=(0, 2, 200),
    m_b=(0, 2, 200),
)
sliders.set_values(
    m_a=0.3,
    m_b=0.75,
)
Hide code cell source
fig, axes = plt.subplots(
    ncols=2,
    nrows=2,
    figsize=[8, 8],
    sharex=True,
    sharey=True,
)
fig.canvas.toolbar_visible = False

(ax_q, ax_rho), (ax_rho_ac, ax_rho_cm) = axes
for ax in [ax_q, ax_rho, ax_rho_cm, ax_rho_ac]:
    ax.set_xlabel("$m$")
    ax.set_yticks([])
for ax in [ax_rho_cm, ax_rho_ac]:
    ax.set_yticks([])

ylim = (-0.1, 1.4)


def func_imag(func, *args, **kwargs):
    return lambda *args, **kwargs: func(*args, **kwargs).imag


def func_real(func, *args, **kwargs):
    return lambda *args, **kwargs: func(*args, **kwargs).real


q_math = ComplexSqrt(sp.Symbol("q^2")) / (8 * sp.pi)
ax_q.set_title(f"${sp.latex(q_math)}$")
controls = iplt.plot(
    plot_domain,
    func_real(np_breakup_momentum),
    label="real",
    **sliders,
    ylim=ylim,
    ax=ax_q,
    alpha=0.7,
)
iplt.plot(
    plot_domain,
    func_imag(np_breakup_momentum),
    label="imaginary",
    controls=controls,
    ylim=ylim,
    ax=ax_q,
    alpha=0.7,
)

ax_rho.set_title(f"${sp.latex(rho_c)}$")
iplt.plot(
    plot_domain,
    func_real(np_rho_c),
    label="real",
    controls=controls,
    ylim=ylim,
    ax=ax_rho,
    alpha=0.7,
)
iplt.plot(
    plot_domain,
    func_imag(np_rho_c),
    label="imaginary",
    controls=controls,
    ylim=ylim,
    ax=ax_rho,
    alpha=0.7,
)

ax_rho_ac.set_title(R"equal mass $\rho^\mathrm{eq}(m^2)$")
iplt.plot(
    plot_domain,
    func_real(np_rho_ac),
    label="real",
    controls=controls,
    ylim=ylim,
    ax=ax_rho_ac,
    alpha=0.7,
)
iplt.plot(
    plot_domain,
    func_imag(np_rho_ac),
    label="imaginary",
    controls=controls,
    ylim=ylim,
    ax=ax_rho_ac,
    alpha=0.7,
)

ax_rho_cm.set_title(R"Chew-Mandelstam $\rho^\mathrm{CM}(m^2)$")
iplt.plot(
    plot_domain,
    func_real(np_rho_cm),
    label="real",
    controls=controls,
    ylim=ylim,
    ax=ax_rho_cm,
    alpha=0.7,
)
iplt.plot(
    plot_domain,
    func_imag(np_rho_cm),
    label="imaginary",
    controls=controls,
    ylim=ylim,
    ax=ax_rho_cm,
    alpha=0.7,
)

fig.tight_layout()
plt.legend(loc="upper right")
plt.show()
../../../_images/dc1a394c253b6a2aa575fa505077d0297011e42e522c4d19485553e8623090da.svg