{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", " \n", "Disclaimer: \n", " \n", "The main objective of the Jupyter notebooks is to show how to use the models of the QENS library by\n", "\n", " \n", "- building a fitting model: composition of models, convolution with a resolution function \n", "- setting and running the fit \n", "- extracting and displaying information about the results \n", "\n", "\n", "These steps have a minimizer-dependent syntax. That's one of the reasons why different minimizers have been used in the notebooks provided as examples. \n", "But, note that the initial guessed parameters might not be optimal, resulting in a poor fit of the reference data.\n", "\n", "
\n", "\n", "# Chudley-Elliott diffusion with bumps\n", "\n", "## Introduction\n", "\n", "
\n", " \n", "The objective of this notebook is to show how to use the Chudley Elliott diffusion model to perform some \n", "fits using bumps .\n", "
\n", "\n", "### Physical units\n", "For information about unit conversion, please refer to the jupyter notebook called `Convert_units.ipynb` in the `tools` folder.\n", "\n", "The dictionary of units defined in the cell below specify the units of the refined parameters adapted to the convention used in the experimental datafile." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Units of parameters for selected QENS model and experimental data\n", "dict_physical_units = {'omega': \"1/ps\", \n", " 'q': \"1/Angstrom\", \n", " 'D': \"ps.Angstrom^2\", \n", " 'L': \"Angstrom\", \n", " 'scale': \"unit_of_signal/ps\",\n", " 'center': \"1/ps\"}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Import required libraries" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import ipywidgets\n", "import QENSmodels\n", "import bumps.names as bmp\n", "from bumps import fitters\n", "from bumps.formatnum import format_uncertainty_pm\n", "import matplotlib.pyplot as plt\n", "\n", "%matplotlib widget" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setting of fitting\n", "\n", "### Create reference data" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nb_points = 500\n", "xx = np.linspace(-4, 4, nb_points)\n", "q = np.linspace(0.2, 2, 10)\n", "added_noise = np.random.normal(0, 1, nb_points)\n", "\n", "chudley_elliott_noisy = QENSmodels.sqwChudleyElliottDiffusion(\n", " xx, \n", " q, \n", " scale=1., \n", " center=0., \n", " D=4,\n", " L=10\n", ") * (1. + 0.01 * added_noise) + 0.01 * added_noise\n", "\n", "# arbitrary values for errors \n", "err = [0.01 for i in range(nb_points)]" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fig = plt.figure()\n", "\n", "for i in range(len(q)):\n", " plt.plot(xx, chudley_elliott_noisy[i,:], label=f\"q={q[i]:.1f}\")\n", " \n", "plt.legend()\n", "plt.title('Signal')\n", "plt.xlabel(f\"$\\hbar \\omega $\")\n", "plt.grid()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Create fitting model" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "model_all_qs = []\n", "for i in range(len(q)):\n", " # Bumps fitting model\n", " model_q = bmp.Curve(\n", " QENSmodels.sqwChudleyElliottDiffusion, \n", " xx, \n", " chudley_elliott_noisy[i],\n", " err,\n", " name=f'q_{q[i]:.2f}',\n", " q=q[i], \n", " scale=1, \n", " center=0., \n", " D=2., \n", " L=8.\n", " )\n", " \n", " model_q.scale.range(0.1, 1e5)\n", " model_q.center.range(-0.1, 0.1)\n", " model_q.D.range(0.1, 5)\n", " model_q.L.range(0.1, 15)\n", " \n", " # Q-independent parameters\n", " if i == 0:\n", " D_q = model_q.D\n", " L_q = model_q.L\n", " else:\n", " model_q.D = D_q\n", " model_q.L = L_q\n", " \n", " model_all_qs.append(model_q)\n", " \n", "problem = bmp.FitProblem(model_all_qs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Display initial configuration: experimental data, fitting model with initial guesses" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "slider = ipywidgets.IntSlider(value=0, min=0, max=len(q)-1, continuous_update=False)\n", "output = ipywidgets.Output()\n", "\n", "def fig_q(model, ax, q_index=0):\n", " \"\"\"\n", " Plot of experimental data, fitting model and residual for a selected q value\n", " \n", " Parameters\n", " ----------\n", " model: list of bumps.curve.Curves for all q\n", " \n", " ax: matplotlib.axes to be updated when changing the ipywidgets\n", " \n", " q_index: int\n", " index of q to be plotted\n", " \n", " \"\"\"\n", " model = model[q_index]\n", " ax[0].errorbar(model.x,\n", " model.y,\n", " yerr=model.dy,\n", " label='experimental data',\n", " color='C0')\n", " ax[0].plot(model.x,\n", " model.theory(), \n", " label='theory (model)',\n", " color='C1')\n", " ax[0].set_title(f'Model {model.name} - $\\chi^2$={problem.chisq_str()}')\n", " ax[0].legend()\n", " ax[1].plot(model.x, model.residuals(), marker='o', linewidth=0, markersize=3, color='C0')\n", " \n", "\n", "with output:\n", " fig, ax = plt.subplots(nrows=2, ncols=1, sharex=True)\n", " ax[0].grid(); ax[1].grid()\n", " ax[1].set_ylabel('Residual')\n", " ax[1].set_xlabel(f\"$\\hbar \\omega ({dict_physical_units['center']})$\")\n", " fig_q(model_all_qs, ax, 0)\n", " \n", " \n", "def update_profile(change):\n", " \"\"\"\n", " Update plots for a new q-value\n", " \"\"\"\n", " with output:\n", " ax[0].clear(); ax[1].lines.clear()\n", " ax[0].grid()\n", " fig_q(model_all_qs, ax,change['new'])\n", " \n", "slider.observe(update_profile, names=\"value\")\n", "\n", "slider_label = ipywidgets.Label(\"q value to display\")\n", "slider_comp = ipywidgets.HBox([slider_label, slider])\n", "ipywidgets.VBox([slider_comp, output])\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "problem.summarize().splitlines()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Choice of minimizer for bumps" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "options_dict = {} \n", "\n", "for item in fitters.__dict__.keys():\n", " if item.endswith('Fit') and fitters.__dict__[item].id in fitters.FIT_AVAILABLE_IDS:\n", " options_dict[fitters.__dict__[item].name] = fitters.__dict__[item].id\n", "\n", "w_choice_minimizer = ipywidgets.Dropdown(\n", " options=list(options_dict.keys()),\n", " value='Levenberg-Marquardt',\n", " description='Minimizer:',\n", " layout=ipywidgets.Layout(height='40px'))\n", "\n", "w_choice_minimizer" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Setting for running bumps" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "steps_fitting = ipywidgets.IntText(\n", " value=100,\n", " step=100,\n", " description='Number of steps when fitting',\n", " style={'description_width': 'initial'})\n", "\n", "steps_fitting" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Running the fit\n", "\n", "Run the fit using the *minimizer* defined above with a number of *steps* also specified above" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Preview of the settings\n", "print('Initial chisq', problem.chisq_str())" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "problem.summarize().splitlines()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def settings_selected_optimizer(chosen_minimizer):\n", " \"\"\" \n", " List the settings available for the selected optimizer\n", " \n", " This list can be used as arguments for the `fit` function\n", " \"\"\"\n", " \n", " assert type(chosen_minimizer) == ipywidgets.widgets.widget_selection.Dropdown\n", " \n", " for item in fitters.__dict__.keys():\n", " if item.endswith('Fit') and \\\n", " fitters.__dict__[item].id == options_dict[chosen_minimizer.value]:\n", " return [elt[0] for elt in fitters.__dict__[item].settings]" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print((f\"With {w_choice_minimizer.value} optimizer, \"\n", " f\"you can use {settings_selected_optimizer(w_choice_minimizer)} as arguments of `fit`\"))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "result = fitters.fit(\n", " problem,\n", " starts=10,\n", " keep_best=True,\n", " method=options_dict[w_choice_minimizer.value], \n", " steps=int(steps_fitting.value))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Showing the results" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "problem.summarize().splitlines()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Other method to display the results of the fit (chi**2 and parameters' values)\n", "print(\"final chisq\", problem.chisq_str())\n", "for k, v, dv in zip(problem.labels(), result.x, result.dx):\n", " if k in dict_physical_units.keys():\n", " print(k, \":\", format_uncertainty_pm(v, dv), dict_physical_units[k])\n", " else:\n", " print(k, \":\", format_uncertainty_pm(v, dv))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Display final configuration: experimental data, fitting model with output of fitting for the refined parameters" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "slider1 = ipywidgets.IntSlider(value=0, min=0, max=len(q)-1, continuous_update=False)\n", "output1 = ipywidgets.Output()\n", "\n", "with output1:\n", " fig1, ax1 = plt.subplots(nrows=2, ncols=1, sharex=True)\n", " ax1[0].grid(); ax1[1].grid()\n", " ax1[1].set_ylabel('Residual')\n", " ax1[1].set_xlabel(f\"$\\hbar \\omega ({dict_physical_units['center']})$\")\n", " fig_q(model_all_qs, ax1, 0)\n", " \n", "def update_profile1(change):\n", " \"\"\"\n", " Update plots for a new q-value\n", " \"\"\"\n", " with output1:\n", " ax1[0].clear(); ax1[1].lines.clear()\n", " ax1[0].grid()\n", " fig_q(model_all_qs, ax1, change['new'])\n", " \n", "slider1.observe(update_profile1, names=\"value\")\n", "\n", "slider1_label = ipywidgets.Label(\"q value to display\")\n", "slider1_comp = ipywidgets.HBox([slider1_label, slider1])\n", "ipywidgets.VBox([slider1_comp, output1])" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "hide_input": false, "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.8.9" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": true, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 4 }