diff --git a/.cspell.json b/.cspell.json index 54b417d7..d40af59a 100644 --- a/.cspell.json +++ b/.cspell.json @@ -62,6 +62,7 @@ "ampform", "arange", "asdot", + "aslatex", "atfi", "autoupdate", "axhline", @@ -138,6 +139,7 @@ "nbformat", "nbins", "nbody", + "nbsp", "nbsphinx", "ncalls", "ncols", diff --git a/docs/amplitude-analysis.ipynb b/docs/amplitude-analysis.ipynb index f49a4cbf..52ff4bf4 100644 --- a/docs/amplitude-analysis.ipynb +++ b/docs/amplitude-analysis.ipynb @@ -61,7 +61,9 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "tags": [] + }, "source": [ "# Amplitude analysis" ] @@ -70,17 +72,17 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "While TensorWaves can handle {doc}`general mathematical expressions `, it was originally created to perform **amplitude analysis** / **Partial Wave Analysis** (PWA), that is, to fit amplitude models to four-momenta data samples.\n", + "While TensorWaves can handle [general mathematical expressions](./usage.ipynb), it was originally created to perform **amplitude analysis** / **Partial Wave Analysis** (PWA), that is, to fit amplitude models to four-momenta data samples.\n", "\n", "This notebook shows how to do an amplitude analysis with the [ComPWA](https://github.com/ComPWA) packages [QRules](https://qrules.rtfd.io), [AmpForm](https://AmpForm.rtfd.io), and [TensorWaves](https://tensorwaves.rtfd.io). The ComPWA workflow generally consists of three stages:\n", "\n", - "1. {ref}`Create an amplitude model ` with {mod}`qrules` and {mod}`ampform`.\n", - "2. {ref}`Generate hit-and-miss data samples ` with this amplitude model.\n", - "3. {ref}`Fit model to the data samples `.\n", + "1. [Create an amplitude model](compwa-step-1) with {mod}`qrules` and {mod}`ampform`.\n", + "2. [Generate hit-and-miss data samples](compwa-step-2) with this amplitude model.\n", + "3. [Fit model to the data samples](compwa-step-3).\n", "\n", ":::{note}\n", "\n", - "This notebook show several tricks that _can_ be helpful when doing an amplitude analysis. Most are optional though―they only serve to illustrate some tips that can be adopted and worked out further for specific analyses.\n", + "This notebook shows several tricks that _can_ be helpful when doing an amplitude analysis. **Most steps are optional** though—they only serve to illustrate some tips that can be adopted and worked out further for specific analyses.\n", "\n", ":::" ] @@ -88,7 +90,8 @@ { "cell_type": "markdown", "metadata": { - "jp-MarkdownHeadingCollapsed": true + "jp-MarkdownHeadingCollapsed": true, + "tags": [] }, "source": [ "(compwa-step-1)=\n", @@ -101,7 +104,7 @@ "jp-MarkdownHeadingCollapsed": true }, "source": [ - "Whether {ref}`generating data ` or {ref}`fitting a model `, TensorWaves takes mathematical expressions as input. When that expression is an amplitude model, it is most convenient to formulate it with {mod}`qrules` and {mod}`ampform`." + "Whether [generating data](compwa-step-2) or [fitting a model](compwa-step-3), TensorWaves takes mathematical expressions as input. When that expression is an amplitude model, it is most convenient to formulate it with {mod}`qrules` and {mod}`ampform`." ] }, { @@ -117,7 +120,7 @@ "class: dropdown\n", "---\n", "\n", - "As {ref}`compwa-step-3` serves to illustrate usage only, we make the amplitude model here a bit simpler by not allowing $\\omega$ resonances (which are narrow and therefore hard to fit). For this reason, we can also limit the {class}`~qrules.settings.InteractionType` to {attr}`~qrules.settings.InteractionType.STRONG`.\n", + "As [](compwa-step-3) serves to illustrate usage only, we make the amplitude model here a bit simpler by not allowing $\\omega$ resonances (which are narrow and therefore hard to fit). For this reason, we can also limit the {class}`~qrules.settings.InteractionType` to {attr}`~qrules.settings.InteractionType.STRONG`.\n", "\n", ":::" ] @@ -187,17 +190,51 @@ "cell_type": "code", "execution_count": null, "metadata": { - "tags": [ - "full-width" - ] + "tags": [] }, "outputs": [], "source": [ "import ampform\n", "\n", "model_builder = ampform.get_builder(reaction)\n", - "model = model_builder.formulate()\n", - "dict(model.parameter_defaults)" + "model_no_dynamics = model_builder.formulate()\n", + "model_no_dynamics.intensity" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "full-width", + "hide-input" + ] + }, + "outputs": [], + "source": [ + "from ampform.io import aslatex\n", + "from IPython.display import Math\n", + "\n", + "Math(aslatex(model_no_dynamics.amplitudes))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "jupyter": { + "source_hidden": true + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "Math(aslatex(model_no_dynamics.parameter_defaults))" ] }, { @@ -214,7 +251,9 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "tags": [] + }, "outputs": [], "source": [ "from ampform.dynamics.builder import (\n", @@ -239,13 +278,21 @@ "cell_type": "code", "execution_count": null, "metadata": { + "jupyter": { + "source_hidden": true + }, "tags": [ - "full-width" + "hide-input" ] }, "outputs": [], "source": [ - "sorted(model.parameter_defaults, key=str)" + "sorted_parameter_defaults = {\n", + " symbol: model.parameter_defaults[symbol]\n", + " for symbol in sorted(model.parameter_defaults, key=str)\n", + "}\n", + "src = aslatex(sorted_parameter_defaults)\n", + "Math(src)" ] }, { @@ -272,13 +319,8 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "In the next steps, we use this {class}`~ampform.helicity.HelicityModel` as a template for a computational function to {ref}`generate data ` and to {ref}`perform a fit `." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ + "In the next steps, we use this {class}`~ampform.helicity.HelicityModel` as a template for a computational function to [generate data](compwa-step-2) and to [perform a fit](compwa-step-3).\n", + "\n", ":::{tip}\n", "\n", "See more advanced examples on {doc}`AmpForm's usage page `.\n", @@ -307,9 +349,9 @@ ":::\n", "::::\n", "\n", - "In this section, we use the {class}`~ampform.helicity.HelicityModel` that we created with {mod}`ampform` in {ref}`the previous step ` to generate a data sample via hit & miss Monte Carlo. We do this with the {mod}`.data` module.\n", + "In this section, we use the {class}`~ampform.helicity.HelicityModel` that we created with {mod}`ampform` in [the previous step](compwa-step-1) to generate a data sample via hit & miss Monte Carlo. We do this with the {mod}`.data` module.\n", "\n", - "First, we {func}`~pickle.load` the {class}`~ampform.helicity.HelicityModel` that was created in the previous step. This does not have to be done if the model has been generated in the same script or notebook, but can be useful if the model was generated elsewhere." + "**Optionally**, we can {func}`~pickle.load` the {class}`~ampform.helicity.HelicityModel` that was created in the previous step. This does not have to be done if the model has been generated in the same script or notebook (like in this notebook), but can be useful if the model was generated elsewhere." ] }, { @@ -323,7 +365,7 @@ "from ampform.helicity import HelicityModel\n", "\n", "with open(\"helicity_model.pickle\", \"rb\") as model_file:\n", - " model: HelicityModel = pickle.load(model_file)" + " imported_model: HelicityModel = pickle.load(model_file)" ] }, { @@ -339,19 +381,20 @@ }, "outputs": [], "source": [ - "reaction_info = model.reaction_info\n", - "initial_state = next(iter(reaction_info.initial_state.values()))\n", + "initial_state, *_ = imported_model.reaction_info.initial_state.values()\n", "print(\"Initial state:\")\n", "print(\" \", initial_state.name)\n", "print(\"Final state:\")\n", - "for i, p in reaction_info.final_state.items():\n", - " print(f\" {i}: {p.name}\")" + "for i, p in imported_model.reaction_info.final_state.items():\n", + " print(f\" {i}: {p.name}\")\n", + "del initial_state" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ + "(compwa-step-2.1)=\n", "### 2.1 Generate phase space sample" ] }, @@ -361,7 +404,7 @@ "source": [ "The {class}`~qrules.transition.ReactionInfo` class defines the constraints of the phase space. As such, we have enough information to generate a **phase-space sample** for this particle reaction. We do this with a {class}`.TFPhaseSpaceGenerator` class, which is a {class}`.DataGenerator` for a {obj}`.DataSample` of **four-momenta** arrays (using {obj}`tensorflow ` and the [`phasespace`](https://phasespace.readthedocs.io) package as a back-end). We also need to construct a {class}`.RealNumberGenerator` that can generate random numbers. {class}`.TFUniformRealNumberGenerator` is the natural choice here.\n", "\n", - "As opposed to the main {ref}`amplitude-analysis:Step 2: Generate data` of the main usage example page, we will generate a **deterministic** data sample. This can be done by feeding a {class}`.RealNumberGenerator` with a specific {attr}`~.RealNumberGenerator.seed` and giving that generator to the {meth}`.TFPhaseSpaceGenerator.generate` method:" + "As opposed to the main [](compwa-step-2) of the main usage example page, we will generate a **deterministic** data sample. This can be done by feeding a {class}`.RealNumberGenerator` with a specific {attr}`~.RealNumberGenerator.seed` and giving that generator to the {meth}`.TFPhaseSpaceGenerator.generate` method:" ] }, { @@ -392,8 +435,8 @@ "\n", "rng = TFUniformRealNumberGenerator(seed=0)\n", "phsp_generator = TFPhaseSpaceGenerator(\n", - " initial_state_mass=reaction_info.initial_state[-1].mass,\n", - " final_state_masses={i: p.mass for i, p in reaction_info.final_state.items()},\n", + " initial_state_mass=reaction.initial_state[-1].mass,\n", + " final_state_masses={i: p.mass for i, p in reaction.final_state.items()},\n", ")\n", "phsp_momenta = phsp_generator.generate(100_000, rng)" ] @@ -430,13 +473,47 @@ "raw_mimetype": "text/restructuredtext" }, "source": [ - "The resulting phase space sample is a {obj}`dict` of final state IDs to an {obj}`~numpy.array` of four-momenta. In the last step, we converted this sample in such a way that it is rendered as an understandable {class}`pandas.DataFrame`." + "The resulting phase space sample is a {obj}`dict` of final state IDs to an {obj}`~numpy.array` of four-momenta. In the last step, we converted this sample in such a way that it is rendered as an understandable {class}`pandas.DataFrame`.\n", + "\n", + ":::{admonition} Four-momentum arrays\n", + "{doc}`Kinematic expressions ` from AmpForm that involve four-momenta should be formatted as $\\left(E, p_x, p_y, p_z\\right)$. In addition, the shape of input arrays should be `(n, 4)` with `n` the number of events.\n", + ":::\n", + "\n", + ":::{warning}\n", + "When using the helicity formalism, the sum of the four-momenta should be in the rest frame, that is $\\sum_i p_i = \\left(m_A, 0, 0, 0\\right)$ with $m_A$ the mass of the decaying particle $A$. This is because the helicity formalisms boosts through the decay chain starting from particle $A$. **Take care to boost your experimental data into the rest frame**, optionally following the {doc}`kinematics classes provided by AmpForm `.\n", + ":::" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "tags": [] + }, + "outputs": [], + "source": [ + "import numpy as np\n", + "\n", + "p = np.array(list(phsp_momenta.values()))\n", + "p.shape" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "tags": [] + }, + "outputs": [], + "source": [ + "p.sum(axis=0).round(decimals=14)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ + "(compwa-step-2.2)=\n", "### 2.2 Generate intensity-based sample" ] }, @@ -448,7 +525,7 @@ "\n", "An intensity-based sample is generated over a phase space sample using the {class}`.IntensityDistributionGenerator`. Its usage is similar to {class}`.TFPhaseSpaceGenerator`, but now you have to provide a {obj}`.Function` as well as a {obj}`.DataTransformer` that is used to transform the four-momentum phase space sample to a data sample that can be understood by the {obj}`.Function`.\n", "\n", - "Now, recall that in {ref}`compwa-step-1`, we used the helicity formalism to mathematically express the reaction in terms of an amplitude model. TensorWaves needs to convert this {obj}`~ampform.helicity.HelicityModel` to a {obj}`.Function` object that can perform fast computations. This can be done with {func}`.create_parametrized_function`:" + "Now, recall that in [](compwa-step-1), we used the helicity formalism to mathematically express the reaction in terms of an amplitude model. TensorWaves needs to convert this {obj}`~ampform.helicity.HelicityModel` to a {obj}`.Function` object that can perform fast computations. This can be done with {func}`.create_parametrized_function`:" ] }, { @@ -464,10 +541,9 @@ "from tensorwaves.function.sympy import create_parametrized_function\n", "\n", "unfolded_expression = model.expression.doit()\n", - "intensity = create_parametrized_function(\n", + "intensity_func = create_parametrized_function(\n", " expression=unfolded_expression,\n", " parameters=model.parameter_defaults,\n", - " max_complexity=200,\n", " backend=\"numpy\",\n", ")" ] @@ -478,26 +554,16 @@ "source": [ ":::{tip}\n", "\n", - "We made use of {func}`.fast_lambdify` here by specifying `max_complexity`. See {ref}`usage/faster-lambdify:Specifying complexity`.\n", + "If {func}`.create_parametrized_function` takes a long time, have a look at {doc}`usage/faster-lambdify`.\n", + "\n", + ":::\n", "\n", - ":::" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ ":::{seealso}\n", "\n", "{ref}`usage/basics:Hit & miss`\n", "\n", - ":::" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ + ":::\n", + "\n", "A problem is that {class}`.ParametrizedBackendFunction` takes a {obj}`.DataSample` with kinematic variables for the helicity formalism as input, not a set of four-momenta. We therefore need to construct a {class}`.DataTransformer` to transform these four-momenta to function variables. In this case, we work with the helicity formalism, so we construct a {class}`.SympyDataTransformer`:" ] }, @@ -522,14 +588,9 @@ "\n", "{ref}`usage:Generate and transform data`\n", "\n", - ":::" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "That's it, now we have enough info to create an intensity-based data sample. Notice how the structure of the output data is the same as the {ref}`phase-space sample we generated previously `:" + ":::\n", + "\n", + "That's it, now we have enough info to create an intensity-based data sample. Notice how the structure of the output data is the same as the [phase-space sample we generated previously](compwa-step-2.1):" ] }, { @@ -548,12 +609,12 @@ ")\n", "\n", "weighted_phsp_generator = TFWeightedPhaseSpaceGenerator(\n", - " initial_state_mass=reaction_info.initial_state[-1].mass,\n", - " final_state_masses={i: p.mass for i, p in reaction_info.final_state.items()},\n", + " initial_state_mass=reaction.initial_state[-1].mass,\n", + " final_state_masses={i: p.mass for i, p in reaction.final_state.items()},\n", ")\n", "data_generator = IntensityDistributionGenerator(\n", " domain_generator=weighted_phsp_generator,\n", - " function=intensity,\n", + " function=intensity_func,\n", " domain_transformer=helicity_transformer,\n", ")\n", "data_momenta = data_generator.generate(10_000, rng)\n", @@ -583,6 +644,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ + "(compwa-step-2.3)=\n", "### 2.3 Visualize kinematic variables" ] }, @@ -608,6 +670,10 @@ "cell_type": "markdown", "metadata": {}, "source": [ + "::::{note}\n", + "Check the remark about four-momentum format and rest frame of the decaying particle [here](compwa-step-2.1).\n", + "::::\n", + "\n", "The {obj}`.DataSample` is a mapping of kinematic variables names to a 1-dimensional array of values. The numbers you see here are final state IDs as defined in the {class}`~ampform.helicity.HelicityModel` member of the {class}`~ampform.helicity.HelicityModel`:" ] }, @@ -621,7 +687,7 @@ }, "outputs": [], "source": [ - "for state_id, particle in reaction_info.final_state.items():\n", + "for state_id, particle in reaction.final_state.items():\n", " print(f\"ID {state_id}:\", particle.name)" ] }, @@ -634,13 +700,8 @@ "class: dropdown\n", "---\n", "By default, {mod}`tensorwaves` only generates invariant masses of the {class}`Topologies ` that are of relevance to the decay problem. In this case, we only have resonances $f_0 \\to \\pi^0\\pi^0$. If you are interested in more invariant mass combinations, you can do so with the method {meth}`~ampform.kinematics.HelicityAdapter.register_topology`.\n", - "````" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ + "````\n", + "\n", "The {obj}`.DataSample` can easily be converted to a {class}`pandas.DataFrame`:" ] }, @@ -686,7 +747,7 @@ "from matplotlib import cm\n", "\n", "resonances = sorted(\n", - " reaction_info.get_intermediate_particles(),\n", + " reaction.get_intermediate_particles(),\n", " key=lambda p: p.mass,\n", ")\n", "evenly_spaced_interval = np.linspace(0, 1, len(resonances))\n", @@ -711,7 +772,7 @@ "source": [ ":::{seealso}\n", "\n", - "{ref}`amplitude-analysis:Extract intensity components`\n", + "[](#extract-intensity-components)\n", "\n", ":::" ] @@ -748,7 +809,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "In the {ref}`next step `, we illustrate how to {meth}`~.Minuit2.optimize` the intensity model to these data samples." + "In the [next step](compwa-step-3), we illustrate how to {meth}`~.Minuit2.optimize` the intensity model to these data samples." ] }, { @@ -766,7 +827,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "As explained in the {ref}`previous step `, a {class}`.ParametrizedFunction` can compute a list of intensities (real numbers) for an input {obj}`.DataSample`. At this stage, we want to optimize the parameters of this {class}`.ParametrizedFunction`, so that it matches the distribution of our data sample. This is what we call 'fitting'." + "As explained in the [previous step](compwa-step-2), a {class}`.ParametrizedFunction` can compute a list of intensities (real numbers) for an input {obj}`.DataSample`. At this stage, we want to optimize the parameters of this {class}`.ParametrizedFunction`, so that it matches the distribution of our data sample. This is what we call 'fitting'." ] }, { @@ -796,6 +857,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ + "(compwa-step-3.1)=\n", "### 3.1 Prepare parametrized function" ] }, @@ -803,11 +865,11 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "In principle, we can use the same {class}`.ParametrizedFunction` as the one that we created in {ref}`amplitude-analysis:2.2 Generate intensity-based sample`. However, when fitting such a function to a data distribution, an {class}`.Optimizer` will evaluate this function numerous times, so it is smart to apply some optimizations to the underlying expression tree before-hand.\n", + "In principle, we can use the same {class}`.ParametrizedFunction` as the one that we created in [](compwa-step-2.2). However, when fitting such a function to a data distribution, an {class}`.Optimizer` will evaluate this function numerous times, so it is smart to apply some optimizations to the underlying expression tree before-hand.\n", "\n", ":::{tip}\n", "\n", - "The sections below illustrate some tricks for how to simplify the expression tree underneath a {class}`.ParametrizedFunction`. **Most of this can also be achieved with {func}`.create_cached_function`**, which is illustrated in {doc}`/usage/caching` and under {ref}`compwa-create_cached_function`. But note that it is still smart to {ref}`cast complex-valued data `.\n", + "The sections below illustrate some tricks for how to simplify the expression tree underneath a {class}`.ParametrizedFunction`. **Most of this can also be achieved with {func}`.create_cached_function`**, which is illustrated in {doc}`/usage/caching` and under [](compwa-create_cached_function). But note that it is still smart to [cast complex-valued data](#cast-complex-valued-data).\n", "\n", ":::" ] @@ -887,7 +949,7 @@ "source": [ ":::{admonition} Complex-valued parameters\n", "\n", - "If initial parameter values are {obj}`complex`, the parameter is split into a real and an imaginary part during the fit. See also {ref}`amplitude-analysis:Covariance matrix`.\n", + "If initial parameter values are {obj}`complex`, the parameter is split into a real and an imaginary part during the fit. See also [](#covariance-matrix).\n", "\n", ":::" ] @@ -946,8 +1008,8 @@ }, "outputs": [], "source": [ - "old = __\n", - "new = _\n", + "old = sp.count_ops(unfolded_expression)\n", + "new = sp.count_ops(optimized_expression)\n", "assert old > new" ] }, @@ -1001,8 +1063,8 @@ }, "outputs": [], "source": [ - "old = ___\n", - "new = _\n", + "old = sp.count_ops(unfolded_expression)\n", + "new = sp.count_ops(optimized_expression)\n", "assert old > new" ] }, @@ -1047,8 +1109,8 @@ }, "outputs": [], "source": [ - "old = __\n", - "new = _\n", + "old = sp.count_ops(unfolded_expression)\n", + "new = sp.count_ops(optimized_expression)\n", "assert old > new" ] }, @@ -1077,10 +1139,7 @@ "\n", "def safe_downcast_to_real(data: DataSample) -> DataSample:\n", " # using isrealobj instead of real_if_close to keep same array backend\n", - " return {\n", - " key: array.real if np.isrealobj(array) else array\n", - " for key, array in data.items()\n", - " }" + " return {k: v.real if np.isrealobj(v) else v for k, v in data.items()}" ] }, { @@ -1122,7 +1181,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Note that the intensities computed by the optimized function are indeed the same as the original intensity function that was created in {ref}`amplitude-analysis:2.2 Generate intensity-based sample` and that it is much faster!" + "Note that the intensities computed by the optimized function are indeed the same as the original intensity function that was created in [](compwa-step-2.2) and that it is much faster!" ] }, { @@ -1134,7 +1193,7 @@ "# JIT-compile functions and test equality\n", "np.testing.assert_array_almost_equal(\n", " optimized_function(data_real),\n", - " intensity(data_real),\n", + " intensity_func(data_real),\n", " decimal=13,\n", ")" ] @@ -1153,7 +1212,7 @@ "metadata": {}, "outputs": [], "source": [ - "%timeit -n1 intensity(data)\n", + "%timeit -n1 intensity_func(data)\n", "%timeit -n1 optimized_function(data_real)" ] }, @@ -1221,9 +1280,9 @@ "(compwa-create_cached_function)=\n", "#### Simplified procedure: {func}`.create_cached_function`\n", "\n", - "As noted under {ref}`amplitude-analysis:3.1 Prepare parametrized function`, most of what is described in this section can be achieved with the function {func}`.create_cached_function`. The idea is described on {doc}`/usage/caching`, but in this section, we show how this translates to amplitude analysis.\n", + "As noted under [](compwa-step-3.1), most of what is described in this section can be achieved with the function {func}`.create_cached_function`. The idea is described on {doc}`/usage/caching`, but in this section, we show how this translates to amplitude analysis.\n", "\n", - "First, note that {func}`.create_cached_function` works with mappings and iterable of {class}`sympy.Symbol ` and not with the {obj}`str` mappings that we defined in {ref}`amplitude-analysis:Determine free parameters`. We can convert the parameter names back to {class}`~sympy.core.symbol.Symbol`s as follows:" + "First, note that {func}`.create_cached_function` works with mappings and iterable of {class}`sympy.Symbol ` and not with the {obj}`str` mappings that we defined in [](#determine-free-parameters). We can convert the parameter names back to {class}`~sympy.core.symbol.Symbol`s as follows:" ] }, { @@ -1259,7 +1318,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "This gives us all the information we need to create a cached function, convert the data samples to cached data samples and construct an {class}`.Estimator` with these transformed items (compare {ref}`amplitude-analysis:3.2 Define estimator`):" + "This gives us all the information we need to create a cached function, convert the data samples to cached data samples and construct an {class}`.Estimator` with these transformed items (compare [](compwa-step-3.2)):" ] }, { @@ -1290,7 +1349,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Note that, just like in {ref}`amplitude-analysis:Create computational function`, the computed intensities of both the original intensity function and the cached function are indeed the same:" + "Note that, just like in [](#create-computational-function), the computed intensities of both the original intensity function and the cached function are indeed the same:" ] }, { @@ -1301,7 +1360,7 @@ "source": [ "np.testing.assert_array_almost_equal(\n", " cached_function(cached_data),\n", - " intensity(data_real),\n", + " intensity_func(data_real),\n", " decimal=13,\n", ")" ] @@ -1320,7 +1379,7 @@ "metadata": {}, "outputs": [], "source": [ - "%timeit -n1 intensity(data_real)\n", + "%timeit -n1 intensity_func(data_real)\n", "%timeit -n1 cached_function(cached_data)" ] }, @@ -1328,6 +1387,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ + "(compwa-step-3.2)=\n", "### 3.2 Define estimator" ] }, @@ -1337,7 +1397,7 @@ "source": [ "To perform a fit, you need to define an {class}`.Estimator`. This is a measure for the discrepancy between the {class}`.ParametrizedFunction` and the data distribution to which you fit it. In PWA, we usually use an **unbinned negative log likelihood estimator** ({class}`.UnbinnedNLL`).\n", "\n", - "Generally, the {class}`.ParametrizedFunction` is not normalized with regards to the data sample, while a log likelihood estimator requires a normalized function. This is where the {ref}`phase space data ` comes into play again: the {class}`.ParametrizedFunction` is evaluated over the phase space data, so that its output can be used as a normalization factor.\n", + "Generally, the {class}`.ParametrizedFunction` is not normalized with regards to the data sample, while a log likelihood estimator requires a normalized function. This is where the [phase space data](compwa-step-2.1) comes into play again: the {class}`.ParametrizedFunction` is evaluated over the phase space data, so that its output can be used as a normalization factor.\n", "\n", "```{margin}\n", "If you want to correct for the efficiency of the detector, you should use a *detector-reconstructed* phase space sample.\n", @@ -1371,6 +1431,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ + "(compwa-step-3.3)=\n", "### 3.3 Optimize fit parameters" ] }, @@ -1380,7 +1441,7 @@ "source": [ "Starting the fit itself is quite simple: just create an {mod}`.optimizer` instance of your choice and call its {meth}`~.Optimizer.optimize` method to start the fitting process. The {meth}`~.Optimizer.optimize` method requires a mapping of parameter names to their initial values. **Only the parameters listed in the mapping are optimized.**\n", "\n", - "Let's have a look at our {ref}`first guess for the parameter values `. Recall that a {class}`.ParametrizedFunction` object computes the intensity for a certain {obj}`.DataSample`. This can be seen nicely when we use these intensities as weights on the phase space sample and plot it together with the original data sample. Here, we look at the invariant mass distribution projection of the final states `1` and `2`, which, {ref}`as we saw before `, is the final state particle pair $\\pi^0\\pi^0$.\n", + "Let's have a look at our [first guess for the parameter values](#determine-free-parameters). Recall that a {class}`.ParametrizedFunction` object computes the intensity for a certain {obj}`.DataSample`. This can be seen nicely when we use these intensities as weights on the phase space sample and plot it together with the original data sample. Here, we look at the invariant mass distribution projection of the final states `1` and `2`, which, [as we saw before](compwa-step-2.3), is the final state particle pair $\\pi^0\\pi^0$.\n", "\n", "Don't forget to use {meth}`~.ParametrizedFunction.update_parameters` first!" ] @@ -1402,9 +1463,8 @@ "import numpy as np\n", "from matplotlib import cm\n", "\n", - "reaction_info = model.reaction_info\n", "resonances = sorted(\n", - " reaction_info.get_intermediate_particles(),\n", + " reaction.get_intermediate_particles(),\n", " key=lambda p: p.mass,\n", ")\n", "\n", @@ -1562,7 +1622,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "In {ref}`amplitude-analysis:3.3 Optimize fit parameters`, we initialized {obj}`.Minuit2` with a {class}`.Loadable` callback. Such callback classes offer the possibility to {meth}`~.Loadable.load_latest_parameters`, so you can pick up the optimize process in case it crashes or if you pause it. Loading the latest parameters goes as follows:" + "In [](compwa-step-3.3), we initialized {obj}`.Minuit2` with a {class}`.Loadable` callback. Such callback classes offer the possibility to {meth}`~.Loadable.load_latest_parameters`, so you can pick up the optimize process in case it crashes or if you pause it. Loading the latest parameters goes as follows:" ] }, { @@ -1585,7 +1645,7 @@ "source": [ ":::{seealso} \n", "\n", - "{ref}`usage/basics:Callbacks` and {ref}`this example ` of a (custom) plotting callback.\n", + "{ref}`usage/basics:Callbacks` and [this example](./usage.ipynb#optimize-parameters) of a (custom) plotting callback.\n", "\n", ":::" ] @@ -1777,7 +1837,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Just like in {ref}`amplitude-analysis:2.2 Generate intensity-based sample`, these _intensity components_ can each be expressed in a computational backend. We do not have to parametrize this function now that we have already optimized the parameters, so we can just substitute the {class}`~sympy.core.symbol.Symbol`s in all expression beforehand and use {func}`.create_function` instead:" + "Just like in [](compwa-step-2.2), these _intensity components_ can each be expressed in a computational backend. We do not have to parametrize this function now that we have already optimized the parameters, so we can just substitute the {class}`~sympy.core.symbol.Symbol`s in all expression beforehand and use {func}`.create_function` instead:" ] }, { @@ -1837,7 +1897,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "The result is a {class}`.PositionalArgumentFunction` that can be plotted just like in {ref}`amplitude-analysis:Plot optimized model`:" + "The result is a {class}`.PositionalArgumentFunction` that can be plotted just like in [](#plot-optimized-model):" ] }, { @@ -1911,8 +1971,8 @@ " for name, sub_expression in model.components.items()\n", " if name.startswith(\"I\")\n", "}\n", - "initial_state_mass = reaction_info.initial_state[-1].mass\n", - "final_state_masses = {i: p.mass for i, p in reaction_info.final_state.items()}\n", + "initial_state_mass = reaction.initial_state[-1].mass\n", + "final_state_masses = {i: p.mass for i, p in reaction.final_state.items()}\n", "\n", "masses = []\n", "for sub_function in sub_intensity_functions.values():\n", diff --git a/docs/conf.py b/docs/conf.py index 89f4fff3..6cbefabc 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -357,6 +357,7 @@ def get_tensorflow_url() -> str: "smartquotes", "substitution", ] +myst_heading_anchors = 4 BINDER_LINK = ( f"https://mybinder.org/v2/gh/ComPWA/{REPO_NAME}/{BRANCH}?filepath=docs/usage" )