From c1773919bc8ac8444ee6008e6eb19f5f8b8a3642 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Fri, 10 Jul 2026 10:25:32 +0100 Subject: [PATCH] fix: rectangular-mesh plot edges follow the mapper's node convention MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The bilinear interpolation mapper places reconstruction value (row r, col c) at node U_y=(n_y-r-1)/(n_y-3), U_x=(c-1)/(n_x-3) in unit CDF space (guard ring at the border, rows flipped by the mapper's row = n - i), but edges_transformed drew a uniform [0,1] partition — pcolormesh plots put flux up to ~1.5 mesh pixels from where the mapper scatters it (the "small offsets" of issue #372, confirmed by a delta-function reproduction). The edges are now the node midpoints pushed through the same CDF machinery (guard-safe denominator for degenerate n<=3 meshes; border guard cells clamp to the data span). Both consumers inherit the fix: the pcolormesh plot path and rectangular_rotated's warped-centre reconstruction; spline and linear CDFs share the upstream index convention. Plot geometry only — the mapper is untouched, so inversions and likelihoods are byte-identical. PyAutoLabs/PyAutoArray#372 Co-Authored-By: Claude Fable 5 --- .../mesh/mesh_geometry/rectangular.py | 29 ++++++++--- .../mesh_geometry/test_rectangular.py | 50 +++++++++++++++++++ 2 files changed, 71 insertions(+), 8 deletions(-) diff --git a/autoarray/inversion/mesh/mesh_geometry/rectangular.py b/autoarray/inversion/mesh/mesh_geometry/rectangular.py index f5c6fb255..5fd310946 100644 --- a/autoarray/inversion/mesh/mesh_geometry/rectangular.py +++ b/autoarray/inversion/mesh/mesh_geometry/rectangular.py @@ -510,16 +510,29 @@ def areas_for_magnification(self): @property def edges_transformed(self): """ - A class packing the ndarrays describing the neighbors of every pixel in the rectangular pixelization (see - `Neighbors` for a complete description of the neighboring scheme). - - The neighbors of a rectangular pixelization are computed by exploiting the uniform and symmetric nature of the - rectangular grid, as described in the method `rectangular_neighbors_from`. + The source-plane cell edges of every mesh pixel, transformed through + the adaptive CDF — one more edge than pixels per axis, shape + ``(n + 1, 2)`` packing (y_edges, x_edges). + + The interpolation mapper is node-based: reconstruction value + ``(row r, col c)`` lives at node ``U_y = (n_y - r - 1)/(n_y - 3)``, + ``U_x = (c - 1)/(n_x - 3)`` in unit CDF space (a guard ring occupies + the border; rows are flipped by the mapper's ``row = n - i``). The + cell edges are therefore the node midpoints — a uniform ``[0, 1]`` + partition drew every cell up to ~1.5 mesh pixels away from where the + mapper scatters its flux (issue #372). Guard-node edges fall outside + the CDF domain and clamp to the data span in the inverse interp, + squashing the border guard cells to the region they actually cover. """ - # edges defined in 0 -> 1 space, there is one more edge than pixel centers on each side - edges_y = self._xp.linspace(1, 0, self.shape_native[0] + 1) - edges_x = self._xp.linspace(0, 1, self.shape_native[1] + 1) + # Node-midpoint edges in unit CDF space (see docstring); the guard + # denominator keeps the degenerate n <= 3 meshes finite (their + # interpolator maps every point to the single interior node anyway). + n_y, n_x = self.shape_native + rows = self._xp.arange(n_y + 1) + cols = self._xp.arange(n_x + 1) + edges_y = (n_y - rows - 0.5) / max(n_y - 3, 1) + edges_x = (cols - 1.5) / max(n_x - 3, 1) edges_reshaped = self._xp.stack([edges_y, edges_x]).T diff --git a/test_autoarray/inversion/pixelization/mesh_geometry/test_rectangular.py b/test_autoarray/inversion/pixelization/mesh_geometry/test_rectangular.py index eadeb26eb..d04fb74ba 100644 --- a/test_autoarray/inversion/pixelization/mesh_geometry/test_rectangular.py +++ b/test_autoarray/inversion/pixelization/mesh_geometry/test_rectangular.py @@ -200,3 +200,53 @@ def test__edges_transformed(mask_2d_7x7): ), abs=1e-8, ) + + +def test__edges_transformed__aligned_with_interpolation_node_convention(): + """ + Regression test for issue #372: the pcolormesh plotting path draws value + (row r, col c) inside the cell bounded by edges_transformed — that cell + must be centred on the interpolation mapper's node for that value, not on + a uniform [0, 1] partition (which shifted plots by ~1.5 mesh pixels). + + A delta function scattered through the mapper must land in plotted cells + whose weighted centroid matches the input point to sub-cell precision. + """ + from autoarray.inversion.mesh.interpolator.rectangular import ( + adaptive_rectangular_mappings_weights_via_interpolation_from, + adaptive_rectangular_transformed_grid_from, + ) + + n = 10 + rng = np.random.default_rng(0) + data_grid = rng.uniform(-1.0, 1.0, (5000, 2)) + test_point = np.array([[0.3, -0.2]]) + + flat_indices, weights = ( + adaptive_rectangular_mappings_weights_via_interpolation_from( + source_grid_size=n, + data_grid=data_grid, + data_grid_over_sampled=test_point, + ) + ) + + # The node-midpoint unit-space edges (what edges_transformed now builds), + # pushed through the same CDF transform. + rows = np.arange(n + 1) + edges_y = (n - rows - 0.5) / (n - 3) + edges_x = (rows - 1.5) / (n - 3) + edges = np.stack([edges_y, edges_x]).T + edges_t = adaptive_rectangular_transformed_grid_from(data_grid, edges) + y_edges, x_edges = edges_t.T + + centroid_y = 0.0 + centroid_x = 0.0 + for flat, weight in zip(flat_indices[0], weights[0]): + r, c = flat // n, flat % n + centroid_y += weight * 0.5 * (y_edges[r] + y_edges[r + 1]) + centroid_x += weight * 0.5 * (x_edges[c] + x_edges[c + 1]) + + # Half a mesh cell in these units is ~0.15; the pre-fix uniform edges + # missed by ~0.4 in y. + assert centroid_y == pytest.approx(test_point[0, 0], abs=0.1) + assert centroid_x == pytest.approx(test_point[0, 1], abs=0.1)