@@ -1,7 +1,7 @@
 general:
   backend: default                      # The matplotlib backend used for visualization. `default` uses the system default, can specify specific backend (e.g. TKAgg, Qt5Agg, WXAgg).
   dpi: 150                              # Resolution in dots per inch used when saving figures. Lower values reduce file size (e.g. 150 gives ~50% smaller files than 300 with negligible quality loss for diagnostic subplots).
-  imshow_origin: upper                  # The `origin` input of `imshow`, determining if pixel values are ascending or descending on the y-axis.
+  imshow_origin: upper                  # The `origin` input of `imshow`, determining if pixel values are ascending or descending on the y-axis. Presentation only: overlays are reflected with the image so they stay registered, but under `lower` the displayed y axis and reported model y parameters (e.g. centres) disagree in sign. Flipping the data itself at load time is a different operation.
   log10_min_value: 1.0e-4               # If negative values are being plotted on a log10 scale, values below this value are rounded up to it (e.g. to remove negative values).
   log10_max_value: 1.0e99               # If positive values are being plotted on a log10 scale, values above this value are rounded down to it.
   zoom_around_mask: true                # If True, plots of data structures with a mask automatically zoom in the masked region.
 
@@ -21,6 +21,8 @@
     norm_from,
     _apply_contours,
     _conf_imshow_origin,
+    _overlay_yx_for_origin,
+    _vector_yx_for_origin,
     plot_regions,
 )
 
@@ -233,6 +235,45 @@ def plot_array(
         _apply_colorbar(im, ax, cb_unit=cb_unit, is_subplot=not owns_figure)
 
     # --- overlays --------------------------------------------------------------
+    # Every (y, x) overlay is reflected for display when ``origin_imshow`` is
+    # ``"lower"`` so it stays registered to the raster (PyAutoArray#565). When
+    # ``extent`` is ``None`` the raster is drawn in pixel-index coordinates,
+    # where matplotlib already places row ``i`` at ``y = i`` under either
+    # origin, so no reflection is applied.
+    def _for_origin(yx):
+        return _overlay_yx_for_origin(yx, extent, origin_imshow)
+
+    mask = _for_origin(mask)
+    border = _for_origin(border)
+    origin = _for_origin(origin)
+    grid = _for_origin(grid)
+    mesh_grid = _for_origin(mesh_grid)
+    if positions is not None:
+        positions = [_for_origin(np.asarray(pos).reshape(-1, 2)) for pos in positions]
+    if lines is not None:
+        lines = [
+            (
+                _for_origin(np.asarray(line).reshape(-1, 2))
+                if line is not None and len(line) > 0
+                else line
+            )
+            for line in lines
+        ]
+    if regions is not None:
+        regions = [
+            (
+                None
+                if region is None
+                else (
+                    _for_origin(region)
+                    if isinstance(region, np.ndarray) and region.ndim == 2
+                    else [_for_origin(polygon) for polygon in region]
+                )
+            )
+            for region in regions
+        ]
+    vector_yx = _vector_yx_for_origin(vector_yx, extent, origin_imshow)
+
     if array_overlay is not None:
         ax.imshow(
             array_overlay,
@@ -299,6 +340,7 @@ def plot_array(
             vector_yx[:, 2],
         )
 
+    # Patches are drawn as given: they are not origin-aware.
     if patches is not None:
         for patch in patches:
             import copy
 
@@ -16,6 +16,8 @@
     save_figure,
     norm_from,
     _conf_imshow_origin,
+    _overlay_yx_for_origin,
+    numpy_grid,
     plot_regions,
 )
 
@@ -130,10 +132,15 @@ def plot_inversion_reconstruction(
 
     is_subplot = not owns_figure
 
+    # Extent of an ``imshow`` raster drawn on *ax* (uniform rectangular meshes
+    # only); ``None`` for ``pcolormesh`` / ``tripcolor`` panels, which are drawn
+    # in data coordinates and need no origin handling.
+    raster_extent = None
+
     if isinstance(
         mapper.interpolator, (InterpolatorRectangular, InterpolatorRectangularUniform)
     ):
-        _plot_rectangular(
+        raster_extent = _plot_rectangular(
             ax, pixel_values, mapper, norm, colormap, extent, is_subplot=is_subplot
         )
     elif isinstance(
@@ -144,6 +151,39 @@ def plot_inversion_reconstruction(
         )
 
     # --- overlays --------------------------------------------------------------
+    # Over an ``imshow`` raster drawn with ``origin="lower"`` the (y, x)
+    # overlays are reflected for display so they stay registered to the
+    # pixels (PyAutoArray#565); every other case is left unchanged.
+    origin_imshow = _conf_imshow_origin() if raster_extent is not None else None
+    if origin_imshow == "lower":
+
+        def _for_origin(yx):
+            return _overlay_yx_for_origin(yx, raster_extent, origin_imshow)
+
+        if lines is not None:
+            lines = [
+                (
+                    _for_origin(np.asarray(line).reshape(-1, 2))
+                    if line is not None and len(line) > 0
+                    else line
+                )
+                for line in lines
+            ]
+        if regions is not None:
+            regions = [
+                (
+                    None
+                    if region is None
+                    else (
+                        _for_origin(region)
+                        if isinstance(region, np.ndarray) and region.ndim == 2
+                        else [_for_origin(polygon) for polygon in region]
+                    )
+                )
+                for region in regions
+            ]
+        grid = _for_origin(numpy_grid(grid))
+
     if lines is not None:
         for i, line in enumerate(lines):
             if line is not None and len(line) > 0:
@@ -217,6 +257,12 @@ def _plot_rectangular(
     is_subplot
         When ``True`` uses ``labelsize_subplot`` from config for the colorbar
         tick labels (matches the behaviour of :func:`~autoarray.plot.array.plot_array`).
+
+    Returns
+    -------
+    tuple or None
+        The ``imshow`` extent when the uniform path draws a raster (so callers
+        can register overlays to its origin), else ``None``.
     """
     from autoarray.inversion.mesh.interpolator.rectangular_uniform import (
         InterpolatorRectangularUniform,
@@ -258,6 +304,7 @@ def _plot_rectangular(
             origin=_conf_imshow_origin(),
         )
         _apply_colorbar(im, ax, is_subplot=is_subplot)
+        return pix_array.geometry.extent
     else:
         y_edges, x_edges = mapper.mesh_geometry.edges_transformed.T
         Y, X = np.meshgrid(y_edges, x_edges, indexing="ij")
@@ -270,6 +317,7 @@ def _plot_rectangular(
             cmap=colormap,
         )
         _apply_colorbar(im, ax, is_subplot=is_subplot)
+        return None
 
 
 def _plot_delaunay(ax, pixel_values, mapper, norm, colormap, extent, is_subplot=False):
 
@@ -1030,6 +1030,60 @@ def _conf_imshow_origin() -> str:
     return origin
 
 
+def _overlay_yx_for_origin(yx, extent, origin_imshow):
+    """
+    Map ``(y, x)`` overlay coordinates so they stay registered to an ``imshow``
+    raster drawn with ``origin=origin_imshow`` over ``extent``.
+
+    PyAutoArray stores native arrays with row 0 at ``ymax``, which is exactly
+    what ``imshow(origin="upper")`` renders, so under ``"upper"`` the input is
+    returned unchanged. Under ``"lower"`` matplotlib draws row 0 at ``ymin``,
+    reflecting the raster about the extent's y midpoint; overlays given in data
+    ``(y, x)`` coordinates are reflected the same way
+    (``y -> ymin + ymax - y``) so that they sit on the pixels they describe.
+
+    This is a display-only reflection: the internal y convention of the data
+    structures is unchanged, and nothing outside the drawn figure is affected.
+
+    Parameters
+    ----------
+    yx
+        Array of shape ``(..., 2+)`` whose column 0 is ``y`` and column 1 is
+        ``x`` (extra trailing columns are carried through untouched). ``None``
+        passes straight through.
+    extent
+        ``(xmin, xmax, ymin, ymax)`` of the raster the overlay is drawn on.
+        ``None`` means no reflection can be defined and the input is returned
+        unchanged.
+    origin_imshow
+        The ``origin`` passed to ``imshow`` (``"upper"`` or ``"lower"``).
+    """
+    if yx is None or extent is None or origin_imshow != "lower":
+        return yx
+
+    yx = np.array(yx, dtype=float)
+    yx[..., 0] = extent[2] + extent[3] - yx[..., 0]
+    return yx
+
+
+def _vector_yx_for_origin(vector_yx, extent, origin_imshow):
+    """
+    Quiver sibling of :func:`_overlay_yx_for_origin` for ``(y, x, dy, dx)`` rows.
+
+    Under ``"lower"`` the anchor ``y`` is reflected about the extent's y
+    midpoint and the ``dy`` component (column 2) is negated, so each arrow
+    keeps pointing at the same pixels of the reflected raster. Under
+    ``"upper"`` (or with ``vector_yx`` / ``extent`` ``None``) the input is
+    returned unchanged. Display-only, like :func:`_overlay_yx_for_origin`.
+    """
+    if vector_yx is None or extent is None or origin_imshow != "lower":
+        return vector_yx
+
+    vector_yx = _overlay_yx_for_origin(vector_yx, extent, origin_imshow)
+    vector_yx[..., 2] = -vector_yx[..., 2]
+    return vector_yx
+
+
 def _conf_output_format() -> str:
     """Return the default output_format from config (``"show"``, ``"png"``, etc.).
 
 
@@ -0,0 +1,142 @@
+"""
+Overlays drawn by ``plot_array`` must stay registered to the raster whatever the
+``visualize/general.yaml -> general -> imshow_origin`` config is (PyAutoArray#565).
+
+The checks are origin-agnostic: they read the drawn overlay coordinates back off
+the axes and map them through the rendered image's own extent and origin to the
+pixel they land on, which must be a bright pixel of an off-centre block.
+"""
+
+import matplotlib
+
+matplotlib.use("Agg")
+
+import matplotlib.pyplot as plt
+import numpy as np
+import pytest
+
+from autonerves import conf
+
+import autoarray as aa
+import autoarray.plot as aplt
+
+
+def _bright_block_array():
+    values = np.zeros((50, 50))
+    values[10:14, 23:27] = 1.0
+    return aa.Array2D.no_mask(values=values, pixel_scales=0.1)
+
+
+def _native_grid(array):
+    return np.asarray(array.mask.derive_grid.all_false.native)
+
+
+def _pixel_value_under(ax, draw_x, draw_y):
+    image = ax.images[0]
+    data = np.asarray(image.get_array())
+    xmin, xmax, ymin, ymax = image.get_extent()
+    ny, nx = data.shape[:2]
+
+    col = int((draw_x - xmin) / (xmax - xmin) * nx)
+    frac = (draw_y - ymin) / (ymax - ymin)
+    row = int((1.0 - frac) * ny) if image.origin == "upper" else int(frac * ny)
+
+    return data[row, col]
+
+
+@pytest.fixture(name="imshow_origin")
+def make_imshow_origin(request):
+    general = conf.instance["visualize"]["general"]["general"]
+    original = general["imshow_origin"]
+    try:
+        general["imshow_origin"] = request.param
+        yield request.param
+    finally:
+        general["imshow_origin"] = original
+
+
+@pytest.mark.parametrize("imshow_origin", ["upper", "lower"], indirect=True)
+def test__plot_array__positions_marker_lands_on_block(imshow_origin):
+    array = _bright_block_array()
+    grid = _native_grid(array)
+    y, x = float(grid[12, 25, 0]), float(grid[12, 25, 1])
+
+    fig, ax = plt.subplots()
+    try:
+        aplt.plot_array(array=array, positions=aa.Grid2DIrregular([(y, x)]), ax=ax)
+
+        assert ax.images[0].origin == imshow_origin
+
+        # The positions scatter is the only s=20, zorder=5 collection (the auto
+        # mask-edge scatter is s=1).
+        positions = [
+            c
+            for c in ax.collections
+            if c.get_zorder() == 5 and np.allclose(c.get_sizes(), 20)
+        ]
+        assert len(positions) == 1
+        draw_x, draw_y = positions[0].get_offsets()[0]
+
+        assert _pixel_value_under(ax, float(draw_x), float(draw_y)) == 1.0
+    finally:
+        plt.close(fig)
+
+
+@pytest.mark.parametrize("imshow_origin", ["upper", "lower"], indirect=True)
+def test__plot_array__lines_vertices_land_on_block(imshow_origin):
+    array = _bright_block_array()
+    grid = _native_grid(array)
+
+    # An asymmetric polyline whose every vertex sits on a bright pixel.
+    pixels = [(10, 23), (13, 24), (11, 26)]
+    line = np.array([[grid[r, c, 0], grid[r, c, 1]] for r, c in pixels])
+
+    fig, ax = plt.subplots()
+    try:
+        aplt.plot_array(array=array, lines=[line], ax=ax)
+
+        assert ax.images[0].origin == imshow_origin
+
+        xy = ax.lines[-1].get_xydata()
+        assert xy.shape == (3, 2)
+
+        for draw_x, draw_y in xy:
+            assert _pixel_value_under(ax, float(draw_x), float(draw_y)) == 1.0
+    finally:
+        plt.close(fig)
+
+
+@pytest.mark.parametrize("imshow_origin", ["upper", "lower"], indirect=True)
+def test__plot_inversion_reconstruction__uniform_mesh_grid_lands_on_pixel(
+    imshow_origin, rectangular_mapper_7x7_3x3
+):
+    mapper = rectangular_mapper_7x7_3x3
+
+    # One off-centre bright source pixel (row 0, column 1 of the 3x3 mesh).
+    pixel_values = np.zeros(9)
+    pixel_values[1] = 1.0
+
+    mesh = aa.Array2D.no_mask(
+        values=np.zeros(mapper.mesh_geometry.shape),
+        pixel_scales=mapper.mesh_geometry.pixel_scales,
+        origin=mapper.mesh_geometry.origin,
+    )
+    y, x = _native_grid(mesh)[0, 1]
+
+    fig, ax = plt.subplots()
+    try:
+        aplt.plot_inversion_reconstruction(
+            pixel_values=pixel_values,
+            mapper=mapper,
+            grid=np.array([[y, x]]),
+            zoom_to_brightest=False,
+            ax=ax,
+        )
+
+        assert ax.images[0].origin == imshow_origin
+
+        draw_x, draw_y = ax.collections[-1].get_offsets()[0]
+
+        assert _pixel_value_under(ax, float(draw_x), float(draw_y)) == 1.0
+    finally:
+        plt.close(fig)