From a0753e0e09341009ad3eb40f67982412029591d0 Mon Sep 17 00:00:00 2001 From: Remco de Boer Date: Wed, 2 Mar 2022 22:38:26 +0100 Subject: [PATCH] docs: create TR-015 from AmpForm's spin-alignment (#119) * build: install ComPWA/ampform@98de70f from main branch * ci: ignore vscode.dev * ci: run linkcheck on Python 3.8 --- .cspell.json | 3 + .github/workflows/linkcheck.yml | 6 +- .pre-commit-config.yaml | 2 +- docs/bibliography.bib | 43 ++ docs/conf.py | 3 +- docs/report/015.ipynb | 1013 +++++++++++++++++++++++++++++++ 6 files changed, 1065 insertions(+), 5 deletions(-) create mode 100644 docs/report/015.ipynb diff --git a/.cspell.json b/.cspell.json index 6c58c5f5..eee00604 100644 --- a/.cspell.json +++ b/.cspell.json @@ -116,6 +116,7 @@ "axvline", "azim", "bdist", + "bgcolor", "boldsymbol", "celltoolbar", "clim", @@ -146,6 +147,7 @@ "facecolors", "figsize", "filterwarnings", + "fontcolor", "fontsize", "getitem", "getsource", @@ -175,6 +177,7 @@ "linkcheck", "linspace", "livereveal", + "marangotto", "markdownlint", "mathbb", "mathcal", diff --git a/.github/workflows/linkcheck.yml b/.github/workflows/linkcheck.yml index 21dc2fa9..77271724 100644 --- a/.github/workflows/linkcheck.yml +++ b/.github/workflows/linkcheck.yml @@ -17,14 +17,14 @@ jobs: runs-on: ubuntu-20.04 steps: - uses: actions/checkout@v2 - - name: Set up Python 3.7 + - name: Set up Python 3.8 uses: actions/setup-python@v2 with: - python-version: "3.7" + python-version: "3.8" - name: Install dependencies run: | python -m pip install --upgrade pip - pip install -c .constraints/py3.7.txt -e .[doc] + pip install -c .constraints/py3.8.txt -e .[doc] sudo apt-get -y install graphviz pandoc - name: Check external links working-directory: docs diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 1c1a57e0..fffd607f 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -41,7 +41,7 @@ repos: )$ - repo: https://github.com/ComPWA/repo-maintenance - rev: 0.0.111 + rev: 0.0.112 hooks: - id: check-dev-files args: diff --git a/docs/bibliography.bib b/docs/bibliography.bib index d1c59e0b..ccaacb47 100644 --- a/docs/bibliography.bib +++ b/docs/bibliography.bib @@ -96,6 +96,20 @@ @misc{kutschkeAngularDistributionCookbook1996 url = {https://home.fnal.gov/~kutschke/Angdist/angdist.ps} } +@article{marangottoHelicityAmplitudesGeneric2020, + title = {Helicity {{Amplitudes}} for {{Generic Multibody Particle Decays Featuring Multiple Decay Chains}}}, + author = {Marangotto, Daniele}, + editor = {Vagnozzi, Sunny}, + year = {2020}, + month = dec, + journal = {Advances in High Energy Physics}, + volume = {2020}, + pages = {1--15}, + issn = {1687-7365, 1687-7357}, + doi = {10.1155/2020/6674595}, + url = {https://www.hindawi.com/journals/ahep/2020/6674595/} +} + @book{martinCleanCodeHandbook2009, title = {Clean {{Code}}: {{A Handbook}} of {{Agile Software Craftsmanship}}}, shorttitle = {Clean Code}, @@ -116,6 +130,23 @@ @misc{meyerMatrixTutorial2008 url = {http://www.curtismeyer.com/talks/PWA_Munich_KMatrix.pdf} } +@article{mikhasenkoDalitzplotDecompositionThreebody2020, + title = {Dalitz-Plot Decomposition for Three-Body Decays}, + author = {Mikhasenko, M. and Albaladejo, M. and Bibrzycki, Ł. and {Fernandez-Ramirez}, C. and Mathieu, V. and Mitchell, S. and Pappagallo, M. and Pilloni, A. and Winney, D. and Skwarnicki, T. and Szczepaniak, A. P.}, + year = {2020}, + month = feb, + journal = {Physical Review D}, + volume = {101}, + number = {3}, + eprint = {1910.04566}, + eprinttype = {arxiv}, + pages = {034033}, + issn = {2470-0010, 2470-0029}, + doi = {10.1103/PhysRevD.101.034033}, + url = {https://journals.aps.org/prd/abstract/10.1103/PhysRevD.101.034033}, + archiveprefix = {arXiv} +} + @book{percivalTestDrivenDevelopmentPython2017, title = {Test-{{Driven Development}} with {{Python}}: {{Obey}} the {{Testing Goat}}: {{Using Django}}, {{Selenium}}, and {{JavaScript}}}, shorttitle = {Test-Driven Development with {{Python}}}, @@ -168,4 +199,16 @@ @book{slatkinEffectivePython902019 annotation = {OCLC: 1127093006} } +@article{wangNovelMethodTest2020, + title = {A Novel Method to Test Particle Ordering and Final State Alignment in Helicity Formalism}, + author = {Wang, Mengzhen and Jiang, Yi and Liu, Yinrui and Qian, Wenbin and Lyu, Xiaorui and Zhang, Liming}, + year = {2020}, + month = dec, + journal = {arXiv}, + eprint = {2012.03699}, + eprinttype = {arxiv}, + url = {http://arxiv.org/abs/2012.03699}, + archiveprefix = {arXiv} +} + diff --git a/docs/conf.py b/docs/conf.py index 65cd0e49..79b9b4c7 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -208,8 +208,8 @@ def get_minor_version(package_name: str) -> str: f"https://ipython.readthedocs.io/en/{get_version('IPython')}", None, ), + "ampform": ("https://ampform--245.org.readthedocs.build/en/245", None), "attrs": (f"https://www.attrs.org/en/{get_version('attrs')}", None), - "ampform": ("https://ampform.readthedocs.io/en/stable", None), "expertsystem": ("https://expertsystem.readthedocs.io/en/stable", None), "graphviz": ("https://graphviz.readthedocs.io/en/stable", None), "ipywidgets": ( @@ -251,6 +251,7 @@ def get_minor_version(package_name: str) -> str: linkcheck_anchors = False linkcheck_ignore = [ "http://127.0.0.1:8000", + "https://open.vscode.dev", ] # Settings for myst_nb diff --git a/docs/report/015.ipynb b/docs/report/015.ipynb new file mode 100644 index 00000000..c97570aa --- /dev/null +++ b/docs/report/015.ipynb @@ -0,0 +1,1013 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "hideCode": true, + "hideOutput": true, + "hidePrompt": true, + "jupyter": { + "source_hidden": true + }, + "slideshow": { + "slide_type": "skip" + }, + "tags": [ + "remove-cell" + ] + }, + "outputs": [], + "source": [ + "%config InlineBackend.figure_formats = ['svg']\n", + "import os\n", + "\n", + "STATIC_WEB_PAGE = {\"EXECUTE_NB\", \"READTHEDOCS\"}.intersection(os.environ)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "```{autolink-concat}\n", + "```" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "tags": [] + }, + "source": [ + "# [TR-015] Spin alignment" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + ":::{note}\n", + "\n", + "This report has been implemented through [ampform#245](https://github.com/ComPWA/ampform/pull/245). For details on how to use it, see [this notebook](https://ampform--245.org.readthedocs.build/en/245/usage/helicity/spin-alignment.html).\n", + "\n", + ":::" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "```{autolink-skip}\n", + "```" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input", + "remove-output" + ] + }, + "outputs": [], + "source": [ + "%pip install -q git+https://github.com/ComPWA/ampform@98de70f qrules[viz]==0.9.7 sympy==1.9" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-cell" + ] + }, + "outputs": [], + "source": [ + "import logging\n", + "import warnings\n", + "\n", + "import ampform\n", + "import graphviz\n", + "import qrules\n", + "import sympy as sp\n", + "from ampform.helicity import formulate_wigner_d\n", + "from IPython.display import Math\n", + "\n", + "LOGGER = logging.getLogger()\n", + "LOGGER.setLevel(logging.ERROR)\n", + "warnings.filterwarnings(\"ignore\")\n", + "\n", + "\n", + "def show_transition(transition, **kwargs):\n", + " if \"size\" not in kwargs:\n", + " kwargs[\"size\"] = 5\n", + " dot = qrules.io.asdot(transition, **kwargs)\n", + " display(graphviz.Source(dot))" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "tags": [] + }, + "source": [ + "## Helicity formalism" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Imagine we want to formulate the amplitude for the following **single** {class}`~qrules.transition.StateTransition`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "full_reaction = qrules.generate_transitions(\n", + " initial_state=\"J/psi(1S)\",\n", + " final_state=[\"K0\", (\"Sigma+\", [+0.5]), (\"p~\", [+0.5])],\n", + " allowed_intermediate_particles=[\"Sigma(1660)~-\", \"N(1650)+\"],\n", + " allowed_interaction_types=\"strong\",\n", + " formalism=\"helicity\",\n", + ")\n", + "graphs = full_reaction.to_graphs()\n", + "single_transition_reaction = full_reaction.from_graphs(\n", + " [graphs[0]], formalism=full_reaction.formalism\n", + ")\n", + "transition = single_transition_reaction.transitions[0]\n", + "show_transition(transition)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The specific {attr}`~qrules.transition.State.spin_projection`s for each {attr}`~qrules.transition.State.particle` only make sense _given a specific reference frame_. AmpForm's {class}`~ampform.helicity.HelicityAmplitudeBuilder` interprets these projections as the **helicity** $\\lambda=\\vec{S}\\cdot\\vec{p}$ of each particle _in the rest frame of the parent particle_. For example, the helicity $\\lambda_2=+\\tfrac{1}{2}$ of $\\bar p$ is the helicity as measured in the rest frame of resonance $\\bar\\Sigma(1660)^-$. The reason is that these helicities are needed when formulating the two-particle state for the decay node $\\bar\\Sigma(1660)^- \\to K^0\\bar p$ (see {doc}`ampform:usage/helicity/formalism`).\n", + "\n", + "Ignoring dynamics and coefficients, the {class}`~ampform.helicity.HelicityModel` for this single transition is rather simple:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "builder = ampform.get_builder(single_transition_reaction)\n", + "model = builder.formulate()\n", + "model.expression.subs(model.parameter_defaults).subs(1.0, 1)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The two Wigner-$D$ functions come from the two **two-body decay nodes** that appear in the {class}`~qrules.transition.StateTransition` above. They were formulated as follows:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "sp.Mul(\n", + " formulate_wigner_d(transition, node_id=0),\n", + " formulate_wigner_d(transition, node_id=1),\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now, as {func}`~ampform.helicity.formulate_wigner_d` explains, the numbers that appear in the Wigner-$D$ functions here are computed from the helicities of the decay products. But there's a subtle problem: these helicities are _assumed to be in the rest frame of each parent particle_. For the first node, this is fine, because the parent particle rest frame matches that of the initial state in the {class}`~qrules.transition.StateTransition` above. In the second node, however, we are in a different rest frame. This can result in phase differences for the different amplitudes.\n", + "\n", + "If there is a single decay {class}`~qrules.topology.Topology` in the {class}`~qrules.transition.ReactionInfo` object for which we are formulating an amplitude model, the problem we identified here can be ignored. The reason is that the phase difference for each {class}`~qrules.transition.StateTransition` (with each an identical decay {class}`~qrules.topology.Topology`) is the same and does not introduce interference effects within the coherent sum. It again becomes a problem, however, when we are formulating an amplitude model _with different topologies_. An example would be the following reaction:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "show_transition(full_reaction, collapse_graphs=True)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n", + "When formulating the amplitude model for this reaction, the {class}`~ampform.helicity.HelicityAmplitudeBuilder` implements the 'standard' helicity formalism as described in {cite}`richmanExperimenterGuideHelicity1984, kutschkeAngularDistributionCookbook1996, chungSpinFormalismsUpdated2014` and simply sums over the different amplitudes to get the full amplitude:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input", + "full-width" + ] + }, + "outputs": [], + "source": [ + "builder = ampform.get_builder(full_reaction)\n", + "model = builder.formulate()\n", + "latex = sp.multiline_latex(\n", + " sp.Symbol(\"I\"),\n", + " model.expression.subs(model.parameter_defaults).subs(1.0, 1),\n", + " environment=\"eqnarray\",\n", + ")\n", + "Math(latex)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "As pointed out in {cite}`marangottoHelicityAmplitudesGeneric2020, mikhasenkoDalitzplotDecompositionThreebody2020, wangNovelMethodTest2020`, this is wrong because of the mismatch in reference frames for the helicities." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Aligning reference frames" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In the rest of this document, we follow {cite}`marangottoHelicityAmplitudesGeneric2020` to align all amplitudes in the different topologies back to the initial state reference frame $A$, so that they can be correctly summed up. Specifically, we want to formulate a new, correctly aligned amplitude $\\mathcal{A}^{A\\to 0,1,\\dots}_{m_A,m_0,m_1,\\dots}$ from the original amplitudes $\\mathcal{A}^{A\\to R,S,i,...\\to 0,1,\\dots}_{\\lambda_A,\\lambda_0,\\lambda_1,\\dots}$ by applying Eq.(45) and Eq.(47) for generic, multi-body decays. Here, the $\\lambda$ values are the helicities in the parent rest frame of each two-body decay and the $m$ are the canonical[^canonical] spin projections in the rest frame of the mother particle that is the same no matter the {class}`~qrules.topology.Topology`.\n", + "\n", + "[^canonical]: The canonical rest frame differs from the 'helicity' rest frame in that it is reached by a direct Lorentz boost without rotating the coordinate system in such a way that the $z$-axis aligns with the momentum of one of the decay products.\n", + "\n", + "Just as in {cite}`marangottoHelicityAmplitudesGeneric2020`, we test the implementation with 1-to-3 body decays. We use the notation from {func}`~ampform.helicity.naming.get_boost_chain_suffix` to indicate resonances $R,S,U$. This results in the following figure for the two alignments sums of Equations (45) and (46) in {cite}`marangottoHelicityAmplitudesGeneric2020`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "dot1 = \"\"\"\n", + "digraph {\n", + " bgcolor=none\n", + " rankdir=LR\n", + " edge [arrowhead=none]\n", + " node [shape=none, width=0]\n", + " A\n", + " 0 [fontcolor=red]\n", + " 1 [fontcolor=green, label=<1>]\n", + " 2 [fontcolor=blue, label=<2>]\n", + " { rank=same A }\n", + " { rank=same 0, 1, 2 }\n", + " N0 [label=\"\"]\n", + " N1 [label=\"\"]\n", + " A -> N0 [style=dotted]\n", + " N0 -> N1 [label=\"R = 01\", fontcolor=orange]\n", + " N1 -> 0\n", + " N0 -> 2 [style=dashed]\n", + " N1 -> 1 [style=dashed]\n", + "}\n", + "\"\"\"\n", + "dot2 = \"\"\"\n", + "digraph {\n", + " bgcolor=none\n", + " rankdir=LR\n", + " edge [arrowhead=none]\n", + " node [shape=none, width=0]\n", + " A\n", + " 0 [label=0, fontcolor=red]\n", + " 1 [label=1, fontcolor=green, label=<1>]\n", + " 2 [label=2, fontcolor=blue, label=<2>]\n", + " { rank=same A }\n", + " { rank=same 0, 1, 2 }\n", + " N0 [label=\"\"]\n", + " N1 [label=\"\"]\n", + " A -> N0 [style=dotted]\n", + " N0 -> N1 [label=\"S = 02\", fontcolor=violet]\n", + " N1 -> 0\n", + " N0 -> 1 [style=dashed]\n", + " N1 -> 2 [style=dashed]\n", + "}\n", + "\"\"\"\n", + "display(*map(graphviz.Source, [dot1, dot2]))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The dashed edges and bars above the state IDs indicate \"opposite helicity\" states. The helicity of an **opposite helicity state** gets a minus sign in the Wigner-$D$ function for a two-body state as formulated by {func}`~ampform.helicity.formulate_wigner_d` (see {ref}`report/015:Helicity formalism`) and therefore needs to be defined consistently. AmpForm does this with {func}`~ampform.helicity.decay.is_opposite_helicity_state`.\n", + "\n", + "Opposite helicity states are also of importance in the spin alignment procedure sketched by {cite}`marangottoHelicityAmplitudesGeneric2020`. The Wigner-$D$ functions that appear in Equations (45) and (46) from {cite}`marangottoHelicityAmplitudesGeneric2020`, operate on the spin of the final state, but the angles in the Wigner-$D$ function are taken from the sibling state:" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "$$\n", + "\\begin{eqnarray}\n", + "\\mathcal{A}^{A \\to {\\color{orange}R},2 \\to 0,1,2}_{m_A,m_0,m_1,m_2}\n", + "&=&\n", + " \\sum_{\\lambda_0^{01},\\mu_0^{01},\\nu_0^{01}}\n", + " {\\color{red}{D^{s_0}_{m_0,\\nu_0^{01}}}}\\!\\left({\\color{red}{\\alpha_0^{01}, \\beta_0^{01}, \\gamma_0^{01}}}\\right)\n", + " {\\color{red}{D^{s_0}_{\\nu_0^{01},\\mu_0^{01}}}}\\!\\left({\\color{orange}{\\phi_{_{01}}, \\theta_{_{01}}}}, 0\\right)\n", + " {\\color{red}{D^{s_0}_{\\mu_0^{01},\\lambda_0^{01}}}}\\!\\left({\\color{red}{\\phi_0^{01}, \\theta_0^{01}}}\\right) \\\\\n", + "&\\times&\n", + " \\sum_{\\lambda_1^{01},\\mu_1^{01},\\nu_1^{01}}\n", + " {\\color{green}{D^{s_1}_{m_1,\\nu_1^{01}}}}\\!\\left({\\color{green}{\\alpha_1^{01}, \\beta_1^{01}, \\gamma_1^{01}}}\\right)\n", + " {\\color{green}{D^{s_1}_{\\nu_1^{01},\\mu_1^{01}}}}\\!\\left({\\color{orange}{\\phi_{_{01}}, \\theta_{_{01}}}}, 0\\right)\n", + " {\\color{green}{D^{s_1}_{\\mu_1^{01},\\lambda_1^{01}}}}\\!\\left({\\color{red}{\\phi_0^{01}, \\theta_0^{01}}}\\right) \\\\\n", + "&\\times&\n", + " \\sum_{\\lambda_2^{01}}\n", + " {\\color{blue}{D^{s_2}_{m_2,\\lambda_2^{01}}}}\\!\\left({\\color{orange}{\\phi_{_{01}}, \\theta_{_{01}}}}, 0\\right) \\\\\n", + "&\\times&\n", + " \\mathcal{A}^{A \\to {\\color{orange}R},2 \\to 0,1,2}_{m_A,\\lambda_0^{01},\\bar\\lambda_1^{01},\\bar\\lambda_2^{01}}\n", + "\\end{eqnarray}\n", + "$$ (alignment-R)\n", + "\n", + "$$\n", + "\\begin{eqnarray}\n", + "\\mathcal{A}^{A \\to {\\color{violet}S},1 \\to 0,1,2}_{m_A,m_0,m_1,m_2}\n", + "&=&\n", + " \\sum_{\\lambda_0^{02},\\mu_0^{02},\\nu_0^{02}}\n", + " {\\color{red}{D^{s_0}_{m_0,\\nu_0^{02}}}}\\!\\left({\\color{red}{\\alpha_0^{02}, \\beta_0^{02}, \\gamma_0^{02}}}\\right)\n", + " {\\color{red}{D^{s_0}_{\\nu_0^{02},\\mu_0^{02}}}}\\!\\left({\\color{violet}{\\phi_{_{02}}, \\theta_{_{02}}}}, 0\\right)\n", + " {\\color{red}{D^{s_0}_{\\mu_0^{02},\\lambda_0^{02}}}}\\!\\left({\\color{red}{\\phi_0^{02}, \\theta_0^{02}}}\\right) \\\\\n", + "&\\times&\n", + " \\sum_{\\lambda_1^{02}}\n", + " {\\color{green}{D^{s_1}_{m_1,\\lambda_1^{02}}}}\\!\\left({\\color{violet}{\\phi_{_{02}}, \\theta_{_{02}}}}, 0\\right) \\\\\n", + "&\\times&\n", + " \\sum_{\\lambda_2^{02},\\mu_2^{02},\\nu_2^{02}}\n", + " {\\color{blue}{D^{s_2}_{m_2,\\nu_2^{02}}}}\\!\\left({\\color{blue}{\\alpha_2^{02}, \\beta_2^{02}, \\gamma_2^{02}}}\\right)\n", + " {\\color{blue}{D^{s_2}_{\\nu_2^{02},\\mu_2^{02}}}}\\!\\left({\\color{violet}{\\phi_{_{02}}, \\theta_{_{02}}}}, 0\\right)\n", + " {\\color{blue}{D^{s_2}_{\\mu_2^{02},\\lambda_2^{02}}}}\\!\\left({\\color{red}{\\phi_0^{02}, \\theta_0^{02}}}\\right) \\\\\n", + "&\\times&\n", + " \\mathcal{A}^{A \\to {\\color{violet}S},2 \\to 0,1,2}_{m_A,\\lambda_0^{02},\\bar\\lambda_1^{02},\\bar\\lambda_2^{02}}\n", + "\\end{eqnarray}\n", + "$$ (alignment-S)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This procedure also allows us to formulate the alignment summation for $\\mathcal{A}^{A \\to {\\color{turquoise}U},0 \\to 0,1,2}_{m_A,m_0,m_1,m_2}$:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "dot3 = \"\"\"\n", + "digraph {\n", + " bgcolor=none\n", + " rankdir=LR\n", + " edge [arrowhead=none]\n", + " node [shape=none, width=0]\n", + " 0 [shape=none, label=0, fontcolor=red]\n", + " 1 [shape=none, label=1, fontcolor=green]\n", + " 2 [shape=none, label=2, fontcolor=blue, label=<2>]\n", + " A [shape=none, label=A]\n", + " { rank=same A }\n", + " { rank=same 0, 1, 2 }\n", + " N0 [label=\"\"]\n", + " N1 [label=\"\"]\n", + " A -> N0 [style=dotted]\n", + " N0 -> N1 [label=12>, fontcolor=turquoise, style=dashed]\n", + " N0 -> 0\n", + " N1 -> 1\n", + " N1 -> 2 [style=dashed]\n", + "}\n", + "\"\"\"\n", + "graphviz.Source(dot3)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "$$\n", + "\\begin{eqnarray}\n", + "\\mathcal{A}^{A \\to {\\color{turquoise}U},0 \\to 0,1,2}_{m_A,m_0,m_1,m_2}\n", + "&=&\n", + " \\sum_{\\lambda_0^{12}}\n", + " {\\color{red}{D^{s_0}_{m_0,\\lambda_0^{12}}}}\\!\\left({\\color{red}{\\phi_0, \\theta_0}}, 0\\right) \\\\\n", + "&\\times&\n", + " \\sum_{\\lambda_1^{12},\\mu_1^{12},\\nu_1^{12}}\n", + " {\\color{green}{D^{s_0}_{m_0,\\nu_1^{12}}}}\\!\\left({\\color{green}{\\alpha_1^{12}, \\beta_1^{12}, \\gamma_1^{12}}}\\right)\n", + " {\\color{green}{D^{s_0}_{\\nu_1^{12},\\mu_1^{12}}}}\\!\\left({\\color{red}{\\phi_0, \\theta_0}}, 0\\right)\n", + " {\\color{green}{D^{s_0}_{\\mu_1^{12},\\lambda_1^{12}}}}\\!\\left({\\color{green}{\\phi_1^{12}, \\theta_1^{12}}}\\right) \\\\\n", + "&\\times&\n", + " \\sum_{\\lambda_2^{12},\\mu_2^{12},\\nu_2^{12}}\n", + " {\\color{blue}{D^{s_2}_{m_2,\\nu_2^{12}}}}\\!\\left({\\color{blue}{\\alpha_2^{12}, \\beta_2^{12}, \\gamma_2^{12}}}\\right)\n", + " {\\color{blue}{D^{s_2}_{\\nu_2^{12},\\mu_2^{12}}}}\\!\\left({\\color{red}{\\phi_0, \\theta_0}}, 0\\right)\n", + " {\\color{blue}{D^{s_2}_{\\mu_2^{12},\\lambda_2^{12}}}}\\!\\left({\\color{green}{\\phi_1^{12}, \\theta_1^{12}}}\\right) \\\\\n", + "&\\times&\n", + " \\mathcal{A}^{A \\to {\\color{turquoise}U},2 \\to 0,1,2}_{m_A,\\lambda_1^{12},\\bar\\lambda_1^{12},\\bar\\lambda_2^{12}}\n", + "\\end{eqnarray}\n", + "$$ (alignment-U)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Finally, the total intensity can be computed from these amplitudes by incoherently summing over the initial and final state canonical spin projections (see [Equation (47)](https://downloads.hindawi.com/journals/ahep/2020/6674595.pdf#page=7) in {cite}`marangottoHelicityAmplitudesGeneric2020`):\n", + "\n", + "$$\n", + "I = \\sum_{m_A,m_0,m_1,m_2}\\left|\n", + " \\mathcal{A}^{A \\to {\\color{orange}R},2 \\to 0,1,2}_{m_A,m_0,m_1,m_2} +\n", + " \\mathcal{A}^{A \\to {\\color{violet}S},1 \\to 0,1,2}_{m_A,m_0,m_1,m_2} +\n", + " \\mathcal{A}^{A \\to {\\color{turquoise}U},0 \\to 0,1,2}_{m_A,m_0,m_1,m_2}\n", + "\\right|^2\n", + "$$ (total-intensity)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### $J/\\psi \\to K^0 \\Sigma^+ \\bar{p}$" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-cell" + ] + }, + "outputs": [], + "source": [ + "from ampform.helicity import (\n", + " formulate_helicity_rotation_chain,\n", + " formulate_rotation_chain,\n", + " formulate_spin_alignment,\n", + ")\n", + "\n", + "\n", + "def show_all_spin_matrices(transition, functor, cleanup: bool) -> None:\n", + " for i in transition.final_states:\n", + " state = transition.states[i]\n", + " particle_name = state.particle.latex\n", + " s = sp.Rational(state.particle.spin)\n", + " m = sp.Rational(state.spin_projection)\n", + " display(\n", + " Math(\n", + " Rf\"|s_{i},m_{i}\\rangle=|{s},{m}\\rangle \\quad ({particle_name})\"\n", + " )\n", + " )\n", + " if functor is formulate_rotation_chain:\n", + " args = (transition, i)\n", + " else:\n", + " args = (transition, i, state.spin_projection)\n", + " summation = functor(*args)\n", + " if cleanup:\n", + " summation = summation.cleanup()\n", + " display(summation)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In this section, we test some of the functions from the {mod}`~ampform.helicity` and {mod}`~ampform.kinematics` modules to see if they reproduce Equations {eq}`alignment-R`, {eq}`alignment-S`, and {eq}`alignment-U`. We perform this test on the channel $J/\\psi \\to K^0 \\Sigma^+ \\bar{p}$ with resonances generated for each of the three allowed three-body topologies. The transition that corresponds to Equation {eq}`alignment-R` is shown below.\n", + "\n", + "The first step is to use {func}`~ampform.helicity.formulate_helicity_rotation_chain` to generate the Wigner-$D$ functions for all **helicity rotations** for each final state. These helicity rotations \"undo\" all rotations that came from each Lorentz boosts when boosting from initial state $J/\\psi$ to each final state:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "transition_r = full_reaction.transitions[-1]\n", + "show_transition(transition_r)\n", + "show_all_spin_matrices(\n", + " transition_r, formulate_helicity_rotation_chain, cleanup=True\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The function {func}`~ampform.helicity.formulate_rotation_chain` goes one step further. It adds a **Wigner rotation** to the generated list of helicity rotation Wigner-$D$ functions in case there are resonances in between the initial state and rotated final state. If there are no resonances in between (here, state `2`, the $\\bar p$), there is only one helicity rotation and there is no need for a Wigner rotation." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "show_all_spin_matrices(transition_r, formulate_rotation_chain, cleanup=False)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**These are indeed all the terms that we see in Equation {eq}`alignment-R`!**\n", + "\n", + "To create all sum combinations for all final states, we can use {func}`~ampform.helicity.formulate_spin_alignment`. This should give the sum of Eq.(45):" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "tags": [ + "full-width" + ] + }, + "outputs": [], + "source": [ + "alignment_summation = formulate_spin_alignment(transition_r)\n", + "alignment_summation.cleanup()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Finally, here are the generated spin alignment terms for the other two decay chains. Notice that the first is indeed the same as {eq}`alignment-S`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "reaction_s = qrules.generate_transitions(\n", + " initial_state=\"J/psi(1S)\",\n", + " final_state=[\"K0\", (\"Sigma+\", [+0.5]), (\"p~\", [+0.5])],\n", + " allowed_intermediate_particles=[\"N(1650)+\"],\n", + " allowed_interaction_types=\"strong\",\n", + " formalism=\"helicity\",\n", + ")\n", + "transition_s = reaction_s.transitions[0]\n", + "show_transition(transition_s)\n", + "show_all_spin_matrices(transition_s, formulate_rotation_chain, cleanup=False)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "...and that the second matches Equation {eq}`alignment-U`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "reaction_u = qrules.generate_transitions(\n", + " initial_state=\"J/psi(1S)\",\n", + " final_state=[\"K0\", (\"Sigma+\", [+0.5]), (\"p~\", [+0.5])],\n", + " allowed_intermediate_particles=[\"K*(1680)~0\"],\n", + " allowed_interaction_types=\"strong\",\n", + " formalism=\"helicity\",\n", + ")\n", + "transition_u = reaction_u.transitions[0]\n", + "show_transition(transition_u)\n", + "show_all_spin_matrices(transition_u, formulate_rotation_chain, cleanup=False)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Compute Wigner rotation angles" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now it's still a matter of computing the values for the angles $\\alpha,\\beta,\\gamma$ in the Wigner rotation matrices. These angles represents the difference between the canonical spin frame as attained by a direct boost from the initial state versus a chain of boosts through each resonance. See Equation (36) in {cite}`marangottoHelicityAmplitudesGeneric2020`.\n", + "\n", + "The {mod}`~ampform.kinematics` module can generate an expression for the chain of Lorentz boosts from the initial state to the final state with {func}`~ampform.kinematics.compute_boost_chain`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from ampform.kinematics import (\n", + " compute_boost_chain,\n", + " create_four_momentum_symbols,\n", + ")\n", + "\n", + "dot = qrules.io.asdot(transition_u)\n", + "topology = transition_u.topology\n", + "display(graphviz.Source(dot))\n", + "momenta = create_four_momentum_symbols(topology)\n", + "for state_id in topology.outgoing_edge_ids:\n", + " boosts = compute_boost_chain(topology, momenta, state_id)\n", + " display(sp.Array(boosts))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This chain of Lorentz boosts needs to be 'undo' with a direct Lorentz boost back to the initial state. A contraction of inverse Lorentz boost with the chain of Lorentz boosts can be generated with {func}`~ampform.kinematics.compute_wigner_rotation_matrix`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from ampform.kinematics import compute_wigner_rotation_matrix\n", + "\n", + "for state_id in topology.outgoing_edge_ids:\n", + " expr = compute_wigner_rotation_matrix(topology, momenta, state_id)\n", + " display(expr)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The result of this matrix product is the rotation matrix for the Wigner rotation. The function {func}`~ampform.kinematics.compute_wigner_angles` computes the required Euler angles from this rotation matrix by implementing Equations (B.2-3) from {cite}`marangottoHelicityAmplitudesGeneric2020`:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from ampform.kinematics import compute_wigner_angles\n", + "\n", + "angles = {}\n", + "for state_id in topology.outgoing_edge_ids:\n", + " angle_definitions = compute_wigner_angles(topology, momenta, state_id)\n", + " for name, expr in angle_definitions.items():\n", + " angle_symbol = sp.Symbol(name, real=True)\n", + " angles[angle_symbol] = expr" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "latex_lines = [R\"\\begin{eqnarray}\"]\n", + "for symbol, expr in angles.items():\n", + " latex_lines.append(Rf\"{sp.latex(symbol)}&=&{sp.latex(expr)}\\\\\")\n", + "latex_lines.append(R\"\\end{eqnarray}\")\n", + "Math(\"\\n\".join(latex_lines))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + ":::{note}\n", + "\n", + "In the topology underlying {eq}`alignment-U`, the Wigner rotation matrix with angles $\\alpha_0^{12}, \\beta_0^{12}, \\gamma_0^{12}$ is simply the identity matrix. This is the reason why it can be omitted in {func}`~ampform.helicity.formulate_rotation_chain` and we only have one helicity rotation.\n", + "\n", + ":::" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Computational code" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Classes like {class}`~ampform.kinematics.BoostMatrix` have been split up into smaller {class}`~ampform.sympy.UnevaluatedExpression` classes so that lambdification to NumPy code results in relatively small and fast code, when using `cse=True` in {func}`~sympy.utilities.lambdify.lambdify` (see {class}`~ampform.sympy.NumPyPrintable`)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "import inspect\n", + "\n", + "beta = sp.Symbol(\"beta_1^12\", real=True)\n", + "beta_expr = angles[beta]\n", + "\n", + "func = sp.lambdify(momenta.values(), beta_expr.doit(), cse=True)\n", + "src = inspect.getsource(func)\n", + "n_characters = len(src)\n", + "latex = sp.latex(beta)\n", + "latex += Rf\":\\quad\\text{{{n_characters:,} characters in generated code}}\"\n", + "Math(latex)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Test on data sample" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + ":::{tip}\n", + "\n", + "A test with a larger data distribution is being developed in [TR-013](https://compwa-org--111.org.readthedocs.build/report/013.html).\n", + "\n", + ":::" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The following phase space mini-sample of four-momenta has been generated for the decay $J/\\psi \\to K^0 \\Sigma^+ \\bar{p}$ with the [`tensorwaves.data`](https://tensorwaves.readthedocs.io/en/stable/api/tensorwaves.data.html) module." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "\n", + "phsp = {\n", + " \"p0\": np.array(\n", + " [\n", + " [0.63140486, 0.13166435, -0.35734744, 0.07760603],\n", + " [0.65169531, 0.37242432, 0.12027178, 0.15467675],\n", + " [0.60647425, -0.22286205, -0.175258, 0.19952806],\n", + " [0.72744323, 0.05529811, 0.30502402, -0.43064999],\n", + " [0.76778868, -0.43557036, 0.35491651, -0.16185017],\n", + " ]\n", + " ),\n", + " \"p1\": np.array(\n", + " [\n", + " [1.37017117, 0.173769668, 0.355893315, -0.553093198],\n", + " [1.34556663, -5.272033e-04, -0.3074542, -0.54901747],\n", + " [1.41660182, 0.634007973, -0.0457976658, -0.433700564],\n", + " [1.38592340, 0.138369075, -0.258624859, 0.648189682],\n", + " [1.37858847, 0.551405385, -0.338705615, 0.259105737],\n", + " ]\n", + " ),\n", + " \"p2\": np.array(\n", + " [\n", + " [1.09532397, -0.30543402, 0.00145413, 0.47548716],\n", + " [1.09963805, -0.37189712, 0.18718247, 0.39434072],\n", + " [1.07382393, -0.41114592, 0.22105567, 0.2341725],\n", + " [0.98353336, -0.19366719, -0.04639917, -0.21753969],\n", + " [0.95052285, -0.11583502, -0.01621089, -0.09725557],\n", + " ]\n", + " ),\n", + "}" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "matrix_expr = compute_wigner_rotation_matrix(topology, momenta, state_id=2)\n", + "matrix_expr" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "matrix_func = sp.lambdify(momenta.values(), matrix_expr.doit(), cse=True)\n", + "matrix_array = matrix_func(*phsp.values())\n", + "np.round(matrix_array, decimals=2).real" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "latex_lines = [R\"\\begin{eqnarray}\"]\n", + "for angle_symbol, angle_expr in angles.items():\n", + " angle_func = sp.lambdify(momenta.values(), angle_expr.doit(), cse=True)\n", + " angle_array = angle_func(*phsp.values())\n", + " rounded_values = np.round(angle_array, decimals=2).real\n", + " latex_lines.append(\n", + " Rf\"{sp.latex(angle_symbol)}&=&{sp.latex(sp.Array(rounded_values))}\\\\\"\n", + " )\n", + "latex_lines.append(R\"\\end{eqnarray}\")\n", + "Math(\"\\n\".join(latex_lines))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "dot = qrules.io.asdot(transition_u, collapse_graphs=True)\n", + "graphviz.Source(dot)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + ":::{note}\n", + "\n", + "The {obj}`~numpy.NAN` values above come from the fact that the inverse boost on a boost results in negative values under the square root of $\\gamma=\\sqrt{1-\\beta^2}$. This can be ignored, because the Wigner rotation is simply omitted when formulating the chain of rotation matrices, as noted in {ref}`report/015:Compute Wigner rotation angles`.\n", + "\n", + ":::" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Four-body decay" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The algorithm for computing Euler angles for the Wigner rotation works an arbitrary number of final states. Here, we illustrate this by formulating an expression for the Wigner rotation matrix in a four-body decay." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from qrules.topology import create_isobar_topologies\n", + "\n", + "topology_4body = create_isobar_topologies(4)[1]\n", + "dot = qrules.io.asdot(topology_4body)\n", + "graphviz.Source(dot)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "momenta_4body = create_four_momentum_symbols(topology_4body)\n", + "compute_wigner_rotation_matrix(topology_4body, momenta_4body, state_id=3)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "name": "python", + "version": "3.8.12" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +}