Preparations for the release

.github/workflows/python-ci.yml

@@ -61,7 +61,7 @@ jobs:
         name: unit-tests-job
         flags: unittests
         files: ./coverage-unit.xml
-        fail_ci_if_error: true
+        fail_ci_if_error: false
         verbose: true
         token: ${{ secrets.CODECOV_TOKEN }}
         slug: EasyScience/EasyReflectometryLib

CHANGELOG.md

@@ -1,3 +1,12 @@
+# Version 1.6.0 (1 May 2026)
+
+Add Mighell-based handling of zero-variance points in fitting (issue #256).
+Zero-variance data points are no longer forcibly discarded; instead, a hybrid
+objective applies a Mighell substitution for zero-variance points while using
+standard weighted least squares for the rest. The previous masking behavior is
+available via `objective='legacy_mask'`. New `objective` parameter on
+`MultiFitter`, `fit()`, and `fit_single_data_set_1d()`.
+
 # Version 1.3.3 (17 June 2025)
 
 Added Chi^2 and fit status to fitting results.

docs/src/api/api.rst

@@ -28,6 +28,15 @@ Project provides a higher-level interface for managing models, experiments, and
 
    project
 
+Fitting
+=======
+Fitting helpers and objective functions.
+
+.. toctree::
+   :maxdepth: 1
+
+   fitting
+
 Assemblies
 ==========
 Assemblies are collections of layers that are used to represent a specific physical setup.

docs/src/api/fitting.rst

@@ -0,0 +1,112 @@
+Fitting
+=======
+
+.. currentmodule:: easyreflectometry.fitting
+
+Objective functions and zero-variance handling
+----------------------------------------------
+
+:class:`MultiFitter` supports several objective modes for handling reflectometry
+data during fitting, especially when measured variances are zero or invalid.
+
+The default objective is ``hybrid``. This uses ordinary weighted least squares
+for points with positive variance and applies a Mighell-style substitution only
+to points whose variance is zero or invalid. The older ``legacy_mask`` mode
+drops zero-variance points before fitting. The ``mighell`` mode applies the
+Mighell transform to every point.
+
+Mighell objective
+~~~~~~~~~~~~~~~~~
+
+The full ``mighell`` objective follows the algebraic form of the
+``chi^2_gamma`` statistic described by Mighell for Poisson-distributed count
+data:
+
+.. math::
+
+   \chi^2_\gamma =
+   \sum_i \frac{[n_i + \min(n_i, 1) - m_i]^2}{n_i + 1}
+
+where ``n_i`` are observed counts and ``m_i`` are model values.
+
+In EasyReflectometry this is implemented as a weighted least-squares problem.
+For each observed value ``y_i`` the fitted target is shifted to
+
+.. math::
+
+   y_{\mathrm{eff},i} = y_i + \min(y_i, 1)
+
+and the effective uncertainty is
+
+.. math::
+
+   \sigma_i = \sqrt{y_i + 1}
+
+so the minimized objective is
+
+.. math::
+
+   \sum_i \left(\frac{y_{\mathrm{eff},i} - f_i}{\sigma_i}\right)^2 =
+   \sum_i \frac{[y_i + \min(y_i, 1) - f_i]^2}{y_i + 1}
+
+This is the same algebraic form as Mighell's statistic, with the model value
+``f_i`` replacing ``m_i``.
+
+Scope and interpretation
+~~~~~~~~~~~~~~~~~~~~~~~~
+
+Mighell's statistic was derived for Poisson-distributed count data. In
+reflectometry workflows, the fitted values are usually normalized
+reflectivities or intensities rather than raw counts. They may already have
+been processed, scaled, background-corrected, or otherwise transformed before
+they reach the fitter.
+
+This distinction matters when interpreting the result. The full ``mighell``
+objective is not only a reweighting of residuals; it also changes the fitted
+target from ``y`` to ``y + min(y, 1)``. For values between zero and one, this
+can substantially increase the target value. A fit can therefore have a good
+Mighell objective value while looking poorer against the originally plotted
+reflectivity curve, or while having a worse classical chi-square.
+
+For reflectometry data, ``hybrid`` is generally the recommended compromise:
+it preserves ordinary weighted least-squares behavior where valid variances are
+available, while still allowing zero-variance points to contribute through the
+Mighell-style substitution.
+
+Objective modes
+~~~~~~~~~~~~~~~
+
+``hybrid``
+   Default. Use standard weighted least squares for points with positive
+   variance and apply the Mighell substitution only where variance is zero or
+   invalid.
+
+``mighell``
+   Apply the Mighell transform to all points. The reported objective chi-square
+   is evaluated in transformed objective space and should not be interpreted as
+   a classical chi-square against the original reflectivity values.
+
+``legacy_mask``
+   Remove zero-variance points before fitting and use standard weighted least
+   squares for the remaining points.
+
+``auto``
+   Alias for ``hybrid``.
+
+Fit metrics
+~~~~~~~~~~~
+
+The fitter exposes both objective-space and classical fit metrics after fitting.
+``objective_chi2`` and ``objective_reduced_chi`` describe the minimized
+objective, which may include transformed targets under ``hybrid`` or
+``mighell``. ``classical_chi2`` and ``classical_reduced_chi`` are computed
+against the original observed reflectivity values using only points with
+positive variance.
+
+API reference
+-------------
+
+.. automodule:: easyreflectometry.fitting
+   :members:
+   :undoc-members:
+   :show-inheritance:

docs/src/tutorials/fitting/simple_fitting.ipynb

@@ -30,7 +30,9 @@
    "source": [
     "%matplotlib inline\n",
     "\n",
+    "import matplotlib.pyplot as plt\n",
     "import pooch\n",
+    "import refl1d\n",
     "import refnx\n",
     "\n",
     "import easyreflectometry\n",
@@ -63,7 +65,8 @@
    "outputs": [],
    "source": [
     "print(f'easyreflectometry: {easyreflectometry.__version__}')\n",
-    "print(f'refnx: {refnx.__version__}')"
+    "print(f'refnx: {refnx.__version__}')\n",
+    "print(f'refl1d: {refl1d.__version__}')"
    ]
   },
   {
@@ -395,7 +398,7 @@
     "## Choosing our calculation engine\n",
     "\n",
     "The `easyreflectometry` package enables the calculation of the reflectometry profile using either [*refnx*](https://refnx.readthedocs.io/) or [*Refl1D*](https://refl1d.readthedocs.io/en/latest/).\n",
-    "For this tutorial, we will stick to the current default, which is *refnx*. \n",
+    "We will first run the fit with the current default, *refnx*, and then rebuild the same starting model with *Refl1D* to compare the results.\n",
     "The calculator must be created and associated with the model that we are to fit. "
    ]
   },
@@ -464,7 +467,8 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "analysed = fitter.fit(data)"
+    "analysed_refnx = fitter.fit(data)\n",
+    "analysed = analysed_refnx"
    ]
   },
   {
@@ -513,18 +517,119 @@
     "model"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "id": "2b2fe558",
+   "metadata": {},
+   "source": [
+    "## Repeating the optimisation with Refl1D\n",
+    "\n",
+    "To compare backends fairly, we rebuild the same initial model and then switch the calculator to `refl1d`.\n",
+    "This ensures that the second optimisation starts from the same parameter guesses rather than from the already-optimised `refnx` result."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "341abc6f",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "print(f'refl1d: {refl1d.__version__}')\n",
+    "\n",
+    "si_refl1d = Material(sld=2.07, isld=0, name='Si')\n",
+    "sio2_refl1d = Material(sld=3.47, isld=0, name='SiO2')\n",
+    "film_refl1d = Material(sld=2.0, isld=0, name='Film')\n",
+    "d2o_refl1d = Material(sld=6.36, isld=0, name='D2O')\n",
+    "\n",
+    "si_layer_refl1d = Layer(material=si_refl1d, thickness=0, roughness=0, name='Si layer')\n",
+    "sio2_layer_refl1d = Layer(material=sio2_refl1d, thickness=30, roughness=3, name='SiO2 layer')\n",
+    "film_layer_refl1d = Layer(material=film_refl1d, thickness=250, roughness=3, name='Film Layer')\n",
+    "subphase_refl1d = Layer(material=d2o_refl1d, thickness=0, roughness=3, name='D2O Subphase')\n",
+    "\n",
+    "superphase_refl1d = Multilayer([si_layer_refl1d, sio2_layer_refl1d], name='Si/SiO2 Superphase')\n",
+    "sample_refl1d = Sample(superphase_refl1d, Multilayer(film_layer_refl1d), Multilayer(subphase_refl1d), name='Film Structure')\n",
+    "\n",
+    "resolution_function_refl1d = PercentageFwhm(0.02)\n",
+    "model_refl1d = Model(\n",
+    "    sample=sample_refl1d,\n",
+    "    scale=1,\n",
+    "    background=1e-6,\n",
+    "    resolution_function=resolution_function_refl1d,\n",
+    "    name='Film Model (Refl1D)'\n",
+    ")\n",
+    "\n",
+    "sio2_layer_refl1d.thickness.bounds = (15, 50)\n",
+    "film_layer_refl1d.thickness.bounds = (200, 300)\n",
+    "sio2_layer_refl1d.roughness.bounds = (1, 15)\n",
+    "film_layer_refl1d.roughness.bounds = (1, 15)\n",
+    "subphase_refl1d.roughness.bounds = (1, 15)\n",
+    "film_layer_refl1d.material.sld.bounds = (0.1, 3)\n",
+    "model_refl1d.background.bounds = (1e-8, 1e-5)\n",
+    "model_refl1d.scale.bounds = (0.5, 1.5)\n",
+    "\n",
+    "interface_refl1d = CalculatorFactory()\n",
+    "interface_refl1d.switch('refl1d')\n",
+    "model_refl1d.interface = interface_refl1d\n",
+    "\n",
+    "print(interface_refl1d.current_interface.name)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "2f764d79",
+   "metadata": {},
+   "source": [
+    "We can now fit the same dataset with `refl1d` and compare the fitted curve and reduced chi-squared with the earlier `refnx` result."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "70229cf9",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "fitter_refl1d = MultiFitter(model_refl1d)\n",
+    "analysed_refl1d = fitter_refl1d.fit(data)\n",
+    "\n",
+    "print(f'refnx reduced chi: {analysed_refnx[\"reduced_chi\"]:.4f}')\n",
+    "print(f'refl1d reduced chi: {analysed_refl1d[\"reduced_chi\"]:.4f}')\n",
+    "\n",
+    "qz = data['coords']['Qz_0'].values\n",
+    "reflectivity = data['data']['R_0'].values\n",
+    "uncertainty = data['data']['R_0'].variances**0.5\n",
+    "\n",
+    "plt.figure(figsize=(8, 5))\n",
+    "plt.errorbar(qz, reflectivity, yerr=uncertainty, fmt='o', markersize=3, color='black', alpha=0.6, label='Data')\n",
+    "plt.plot(qz, analysed_refnx['R_0_model'].values, linewidth=2.5, label='refnx fit')\n",
+    "plt.plot(qz, analysed_refl1d['R_0_model'].values, linewidth=2.5, label='refl1d fit')\n",
+    "plt.yscale('log')\n",
+    "plt.xlabel(r'$Q_z$ (1/Å)')\n",
+    "plt.ylabel('Reflectivity')\n",
+    "plt.grid(True, which='both', alpha=0.25)\n",
+    "plt.legend()\n",
+    "plt.show()\n",
+    "\n",
+    "print('refnx model')\n",
+    "print(model)\n",
+    "print()\n",
+    "print('refl1d model')\n",
+    "print(model_refl1d)"
+   ]
+  },
   {
    "cell_type": "markdown",
    "id": "df6db8c3-6515-478b-bd2e-aac967579231",
    "metadata": {},
    "source": [
-    "We note here that the results obtained are very similar to those from the [*refnx* tutorial](https://refnx.readthedocs.io/en/latest/getting_started.html#Fitting-a-neutron-reflectometry-dataset), which is hardly surprising given that we have used the *refnx* engine in this example."
+    "We note here that the results obtained with [*refnx*](https://refnx.readthedocs.io/en/latest/) and [*Refl1D*](https://refl1d.readthedocs.io/en/latest/) are very similar for this simple slab model, with only small backend-dependent differences in the fitted curve and optimised parameters."
    ]
   }
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "easyref",
+   "display_name": ".venv (3.11.9)",
    "language": "python",
    "name": "python3"
   },
@@ -538,7 +643,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.12.9"
+   "version": "3.11.9"
   }
  },
  "nbformat": 4,