DM-55875: Implement Gomes 2025 with the help of Claude/Fable5.

README.rst

@@ -20,12 +20,14 @@ Overview
 ``treegp`` has some special features compared to other available Gaussian Processes codes:
 
 *   Hyperparameters estimation will scale in O(N log(N)) with the the 2-points correlation function estimation compared to O(N^3) with the classical maximum likelihood.
-    
+
 *   Gaussian process interpolation can be performed around a mean function
-    
+
 *   A tool is provided to compute the mean function (``meanify``)
 
-``treegp`` was originally developed for Point Spread Function interpolation within `Piff <https://github.com/rmjarvis/Piff>`_. There is a specific article that describes the math used in ``treegp`` in the context of modelling astrometric shifts of the Subaru Telescope due to atmospheric turbulences. This article can be found 
+*   The measured anisotropic 2D 2-points correlation function can be used directly as the kernel, without fitting any hyperparameters (``optimizer="empirical-2pcf"``), implementing `Gomes et al. (2025) <https://doi.org/10.3847/1538-3881/ae1a7b>`_. The other optimizers (``"two-pcf"``, ``"anisotropic"``) implement `Léget et al. (2021) <https://arxiv.org/abs/2103.09881>`_.
+
+``treegp`` was originally developed for Point Spread Function interpolation within `Piff <https://github.com/rmjarvis/Piff>`_. There is a specific article that describes the math used in ``treegp`` in the context of modelling astrometric shifts of the Subaru Telescope due to atmospheric turbulences. This article can be found
 `here <https://arxiv.org/abs/2103.09881>`_.
 
 

docs/treegp_gp_interp.rst

@@ -118,6 +118,37 @@ To do the gaussian process interpolation with ``treegp`` it follow this API:
 .. image:: fitted.png
 
 
+Using the measured 2-points correlation function as the kernel
+==============================================================
+
+Instead of fitting the hyperparameters of a parametric kernel on the
+measured 2-points correlation function (`Léget et al 2021
+<https://doi.org/10.1051/0004-6361/202140463>`_), the measured anisotropic 2D
+2-points correlation function can be used directly as the kernel, following
+`Gomes et al 2025 <https://doi.org/10.3847/1538-3881/ae1a7b>`_. The measured
+correlation function is cleaned by apodization and thresholding of its
+Fourier power spectrum, tabulated on a grid, and interpolated at the pair
+separations. No hyperparameters are fitted, so the ``kernel`` argument is
+ignored; the grid is controlled by ``max_sep`` (half width) and
+``pixel_size``:
+
+.. code:: ipython3
+
+    gp = treegp.GPInterpolation(
+                optimizer="empirical-2pcf",
+                normalize=True,
+                max_sep=6.0,
+                pixel_size=0.5,
+            )
+    gp.initialize(x, y, y_err=y_err)
+    # measure, clean, and tabulate the 2-points correlation function
+    gp.solve()
+    y_test, y_test_cov = gp.predict(x_test, return_cov=True)
+
+    # measured and cleaned correlation functions can be inspected with:
+    xi, xi_clean, distance, pixel_size = gp.return_empirical_2pcf()
+
+