From 5c3f140ba66f4d7717f900af05a88b5863675007 Mon Sep 17 00:00:00 2001
From: Hamza Abdelhedi <hamza.abdelhedii@gmail.com>
Date: Thu, 27 Aug 2026 23:40:20 +0200
Subject: [PATCH 1/8] ENH: Add Forward-based projection reconstruction

---
 doc/changes/dev/newfeature.rst      |  1 +
 mne/_fiff/proj.py                   | 52 +++++++++++++--
 mne/forward/_field_interpolation.py | 50 ++++++++++++---
 mne/tests/test_proj.py              | 99 +++++++++++++++++++++++++++++
 4 files changed, 189 insertions(+), 13 deletions(-)
 create mode 100644 doc/changes/dev/newfeature.rst

diff --git a/doc/changes/dev/newfeature.rst b/doc/changes/dev/newfeature.rst
new file mode 100644
index 00000000000..edeb49875d6
--- /dev/null
+++ b/doc/changes/dev/newfeature.rst
@@ -0,0 +1 @@
+Projection reconstruction can now use an explicit :class:`~mne.Forward` model and integer spatial rank, by `Hamza Abdelhedi`_.
diff --git a/mne/_fiff/proj.py b/mne/_fiff/proj.py
index 16418e2d601..5f9aaf30ba9 100644
--- a/mne/_fiff/proj.py
+++ b/mne/_fiff/proj.py
@@ -13,6 +13,7 @@
 from ..fixes import _safe_svd
 from ..utils import (
     _check_option,
+    _ensure_int,
     _validate_type,
     fill_doc,
     logger,
@@ -561,7 +562,15 @@ def plot_projs_topomap(
         )
         return fig
 
-    def reconstruct_proj(self, *, projs=None, mode="accurate", origin="auto"):
+    def reconstruct_proj(
+        self,
+        *,
+        projs=None,
+        mode="accurate",
+        origin="auto",
+        forward=None,
+        rank=None,
+    ):
         """Apply SSP projectors and reconstruct the resulting signal in sensor space.
 
         Operates in place.
@@ -574,18 +583,46 @@ def reconstruct_proj(self, *, projs=None, mode="accurate", origin="auto"):
             ``None``, all projectors attached to the instance are used.
         mode : str
             Either ``'accurate'`` or ``'fast'``, determines the quality of the
-            Legendre polynomial expansion used for reconstruction.
+            Legendre polynomial expansion used for geometry-based
+            reconstruction. Ignored when ``forward`` is provided.
         origin : array-like, shape (3,) | str
             Origin of the sphere in the head coordinate frame and in meters.
             Can be ``'auto'`` (default), which means a head-digitization-based
-            origin fit.
+            origin fit. Used for geometry-based reconstruction and ignored when
+            ``forward`` is provided.
+        forward : instance of Forward | None
+            Forward model used to construct the reconstruction field mapping.
+            If ``None`` (default), use the geometry-based field mapping model.
+            If provided, ``rank`` must also be specified.
+        rank : int | None
+            Number of spatial modes to retain when reconstructing with
+            ``forward``. This is a sensor-space reconstruction rank, not a
+            number of sources or dipoles. Must be provided together with
+            ``forward``.
+
+        Notes
+        -----
+        When ``forward`` is provided, the reconstruction uses the sensor-space
+        field covariance formed from the Forward gain matrix. ``rank`` specifies
+        the number of sensor-space modes retained in the reconstruction.
 
         Returns
         -------
         self : same type as the input data
             The modified instance.
         """
-        from ..forward import _map_meg_or_eeg_channels
+        from ..forward import Forward, _map_meg_or_eeg_channels
+
+        if forward is None:
+            if rank is not None:
+                raise ValueError("rank can only be used when forward is provided")
+        else:
+            _validate_type(forward, Forward, "forward")
+            if rank is None:
+                raise ValueError("rank must be provided when forward is provided")
+            rank = _ensure_int(rank, "rank")
+            if rank <= 0:
+                raise ValueError(f"rank must be positive, got {rank}")
 
         if projs is None:
             if len(self.info["projs"]) == 0:
@@ -626,7 +663,12 @@ def reconstruct_proj(self, *, projs=None, mode="accurate", origin="auto"):
                         make_eeg_average_ref_proj(info_to, verbose=False)
                     ]
             mapping = _map_meg_or_eeg_channels(
-                info_from, info_to, mode=mode, origin=origin
+                info_from,
+                info_to,
+                mode=mode,
+                origin=origin,
+                forward=forward,
+                rank=rank,
             )
             self.data[..., picks, :] = np.matmul(mapping, self.data[..., picks, :])
         return self
diff --git a/mne/forward/_field_interpolation.py b/mne/forward/_field_interpolation.py
index 684c955e806..2c631be67c3 100644
--- a/mne/forward/_field_interpolation.py
+++ b/mne/forward/_field_interpolation.py
@@ -12,7 +12,7 @@
 
 from .._fiff.constants import FIFF
 from .._fiff.meas_info import _simplify_info
-from .._fiff.pick import pick_info, pick_types
+from .._fiff.pick import pick_channels_forward, pick_info, pick_types
 from .._fiff.proj import _has_eeg_average_ref_proj, make_projector
 from ..bem import _check_origin
 from ..cov import make_ad_hoc_cov
@@ -21,7 +21,14 @@
 from ..fixes import _safe_svd
 from ..surface import get_head_surf, get_meg_helmet_surf
 from ..transforms import _find_trans, transform_surface_to
-from ..utils import _check_fname, _check_option, _pl, _reg_pinv, logger, verbose
+from ..utils import (
+    _check_fname,
+    _check_option,
+    _pl,
+    _reg_pinv,
+    logger,
+    verbose,
+)
 from ._lead_dots import _do_cross_dots, _do_self_dots, _do_surface_dots, _get_legen_fun
 from ._make_forward import _create_eeg_els, _create_meg_coils, _read_coil_defs
 
@@ -29,14 +36,19 @@
 def _setup_dots(mode, info, coils, ch_type):
     """Set up dot products."""
     int_rad = 0.06
-    noise = make_ad_hoc_cov(info, dict(mag=20e-15, grad=5e-13, eeg=1e-6))
+    noise = _make_field_mapping_noise(info)
     # "fast" uses a coarser (n_coeff=50) Legendre series than "accurate" (n_coeff=100)
     n_coeff = 50 if mode == "fast" else 100
     leg_fun, n_fact = _get_legen_fun(ch_type, False, n_coeff)
     return int_rad, noise, leg_fun, n_fact
 
 
-def _compute_mapping_matrix(fmd, info):
+def _make_field_mapping_noise(info):
+    """Create the ad hoc noise covariance used for field mapping."""
+    return make_ad_hoc_cov(info, dict(mag=20e-15, grad=5e-13, eeg=1e-6))
+
+
+def _compute_mapping_matrix(fmd, info, *, rank=None):
     """Do the hairy computations."""
     logger.info("    Preparing the mapping matrix...")
     # assemble a projector and apply it to the data
@@ -54,7 +66,9 @@ def _compute_mapping_matrix(fmd, info):
 
     # SVD is numerically better than the eigenvalue composition even if
     # mat is supposed to be symmetric and positive definite
-    if fmd.get("pinv_method", "tsvd") == "tsvd":
+    if rank is not None:
+        inv, _, _ = _reg_pinv(whitened_dots, reg=0, rank=rank)
+    elif fmd.get("pinv_method", "tsvd") == "tsvd":
         inv, fmd["nest"] = _pinv_trunc(whitened_dots, fmd["miss"])
     else:
         assert fmd["pinv_method"] == "tikhonov", fmd["pinv_method"]
@@ -110,7 +124,9 @@ def _pinv_tikhonov(x, reg):
     return inv, n
 
 
-def _map_meg_or_eeg_channels(info_from, info_to, mode, *, origin, miss=None):
+def _map_meg_or_eeg_channels(
+    info_from, info_to, mode, *, origin, miss=None, forward=None, rank=None
+):
     """Find mapping from one set of channels to another.
 
     Parameters
@@ -133,8 +149,6 @@ def _map_meg_or_eeg_channels(info_from, info_to, mode, *, origin, miss=None):
     mapping : array, shape (n_to, n_from)
         A mapping matrix.
     """
-    assert origin is not None  # should be assured elsewhere
-
     # no need to apply trans because both from and to coils are in device
     # coordinates
     info_kinds = set(ch["kind"] for ch in info_to["chs"])
@@ -150,6 +164,26 @@ def _map_meg_or_eeg_channels(info_from, info_to, mode, *, origin, miss=None):
     )
     kind = "eeg" if info_kinds[0] == FIFF.FIFFV_EEG_CH else "meg"
 
+    if forward is not None:
+        forward = pick_channels_forward(
+            forward, include=info_from["ch_names"], ordered=True
+        )
+        assert forward["sol"]["row_names"] == info_from["ch_names"]
+        lead_field = forward["sol"]["data"]
+        # Form the sensor-space field covariance from the Forward gain matrix.
+        # As with any Gram representation, very weak modes can be numerically unstable.
+        dots = lead_field @ lead_field.T
+        fmd = dict(
+            kind=kind,
+            ch_names=info_from["ch_names"],
+            noise=_make_field_mapping_noise(info_from),
+            self_dots=dots,
+            surface_dots=dots,
+        )
+        return _compute_mapping_matrix(fmd, info_from, rank=rank)
+
+    assert origin is not None  # should be assured elsewhere
+
     #
     # Step 1. Prepare the coil definitions
     #
diff --git a/mne/tests/test_proj.py b/mne/tests/test_proj.py
index 316258c6eb0..a04c36d06b9 100644
--- a/mne/tests/test_proj.py
+++ b/mne/tests/test_proj.py
@@ -19,10 +19,14 @@
     compute_raw_covariance,
     convert_forward_solution,
     create_info,
+    make_forward_solution,
+    make_sphere_model,
+    pick_channels_forward,
     pick_types,
     read_events,
     read_forward_solution,
     read_source_estimate,
+    read_source_spaces,
     sensitivity_map,
 )
 from mne._fiff.proj import (
@@ -209,6 +213,101 @@ def test_reconstruct_proj(raw_orig, events):
         )
 
 
+@pytest.fixture(scope="module")
+def eeg_forward():
+    """Create a small local EEG forward for projection reconstruction tests."""
+    raw = read_raw_fif(raw_fname, preload=False, verbose=False).pick(picks="eeg")
+    raw.pick(raw.ch_names[:8])
+    src = read_source_spaces(base_dir / "small-src.fif.gz", verbose=False)
+    sphere = make_sphere_model(verbose=False)
+    return make_forward_solution(
+        raw.info,
+        trans=None,
+        src=src,
+        bem=sphere,
+        meg=False,
+        eeg=True,
+        mindist=0.0,
+        verbose=False,
+    )
+
+
+def _direct_forward_reconstruction(evoked, forward, rank):
+    """Compute the independent direct lead-field reconstruction."""
+    projector = make_projector(evoked.info["projs"], evoked.ch_names)[0]
+    forward = pick_channels_forward(forward, include=evoked.ch_names, ordered=True)
+    lead_field = forward["sol"]["data"]
+    u, s, vh = np.linalg.svd(projector @ lead_field, full_matrices=False)
+    mapping = lead_field @ (vh[:rank].T / s[:rank]) @ u[:, :rank].T
+    if _has_eeg_average_ref_proj(evoked.info):
+        mapping -= mapping.mean(axis=0)
+    return mapping @ (projector @ evoked.data)
+
+
+def test_reconstruct_proj_forward(raw_orig, eeg_forward):
+    """Test Forward reconstruction and channel handling."""
+    raw = raw_orig.copy().pick(picks="eeg")
+    raw.pick(raw.ch_names[:8])
+    evoked = EvokedArray(raw.get_data()[:, :10], raw.info, tmin=0.0)
+    evoked.add_proj(
+        _make_test_proj(
+            evoked.ch_names,
+            np.arange(1.0, len(evoked.ch_names) + 1.0),
+            "Forward reconstruction",
+        ),
+        verbose=False,
+    )
+    rank = 3
+    for average_ref in (False, True):
+        this_evoked = evoked.copy()
+        if average_ref:
+            this_evoked.set_eeg_reference(projection=True)
+        expected = _direct_forward_reconstruction(this_evoked, eeg_forward, rank)
+        got = this_evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank).data
+        assert_allclose(got, expected, rtol=1e-10, atol=1e-12)
+
+    reordered = pick_channels_forward(
+        eeg_forward, include=eeg_forward.ch_names[::-1], ordered=True
+    )
+    got = evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank).data
+    got_reordered = evoked.copy().reconstruct_proj(forward=reordered, rank=rank).data
+    assert_allclose(got, got_reordered, rtol=1e-10, atol=1e-12)
+
+    evoked_bad = evoked.copy()
+    evoked_bad.info["bads"] = [evoked_bad.ch_names[0]]
+    got_bad = evoked_bad.reconstruct_proj(forward=eeg_forward, rank=rank).data
+    assert_allclose(got_bad[0], evoked_bad.data[0])
+
+
+def test_reconstruct_proj_forward_validation(eeg_forward):
+    """Test validation of the explicit Forward reconstruction arguments."""
+    info = create_info(eeg_forward["info"]["ch_names"], 100.0, "eeg")
+    evoked = EvokedArray(np.zeros((len(info["ch_names"]), 1)), info, tmin=0.0)
+    evoked.add_proj(
+        _make_test_proj(
+            evoked.ch_names,
+            np.ones(len(evoked.ch_names)),
+            "Forward reconstruction",
+        ),
+        verbose=False,
+    )
+    with pytest.raises(ValueError, match="rank can only be used"):
+        evoked.copy().reconstruct_proj(rank=1)
+    with pytest.raises(ValueError, match="rank must be provided"):
+        evoked.copy().reconstruct_proj(forward=eeg_forward)
+    with pytest.raises(TypeError, match="forward must be an instance of Forward"):
+        evoked.copy().reconstruct_proj(forward=[], rank=1)
+    for rank in (0, -1):
+        with pytest.raises(ValueError, match="rank must be positive"):
+            evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank)
+    for rank in (1.5, True):
+        with pytest.raises(TypeError, match="rank must be an int"):
+            evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank)
+    rank = len(evoked.ch_names) + 1
+    with pytest.raises(ValueError, match="Invalid value for the rank parameter"):
+        evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank)
+
+
 @pytest.mark.parametrize("kind", ["raw", "epochs", "evoked"])
 def test_apply_proj_default(kind):
     """Test that ``projs=None`` preserves legacy behavior."""

From 50f59a02a7be7de9874c1dcfb4fa5e50d4a01b33 Mon Sep 17 00:00:00 2001
From: Hamza Abdelhedi <hamza.abdelhedii@gmail.com>
Date: Fri, 28 Aug 2026 14:31:48 +0200
Subject: [PATCH 2/8] DOC: Fix reconstruct_proj docstring section order

---
 mne/_fiff/proj.py | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/mne/_fiff/proj.py b/mne/_fiff/proj.py
index 5f9aaf30ba9..aeda5c694c5 100644
--- a/mne/_fiff/proj.py
+++ b/mne/_fiff/proj.py
@@ -600,16 +600,16 @@ def reconstruct_proj(
             number of sources or dipoles. Must be provided together with
             ``forward``.
 
+        Returns
+        -------
+        self : same type as the input data
+            The modified instance.
+
         Notes
         -----
         When ``forward`` is provided, the reconstruction uses the sensor-space
         field covariance formed from the Forward gain matrix. ``rank`` specifies
         the number of sensor-space modes retained in the reconstruction.
-
-        Returns
-        -------
-        self : same type as the input data
-            The modified instance.
         """
         from ..forward import Forward, _map_meg_or_eeg_channels
 

From 874820f9fd221859f487a45e9526bcd740f194a4 Mon Sep 17 00:00:00 2001
From: "autofix-ci[bot]" <114827586+autofix-ci[bot]@users.noreply.github.com>
Date: Fri, 28 Aug 2026 09:48:41 +0000
Subject: [PATCH 3/8] [autofix.ci] apply automated fixes

---
 doc/changes/dev/{newfeature.rst => 14235.newfeature.rst} | 0
 1 file changed, 0 insertions(+), 0 deletions(-)
 rename doc/changes/dev/{newfeature.rst => 14235.newfeature.rst} (100%)

diff --git a/doc/changes/dev/newfeature.rst b/doc/changes/dev/14235.newfeature.rst
similarity index 100%
rename from doc/changes/dev/newfeature.rst
rename to doc/changes/dev/14235.newfeature.rst

From ddca38fe0100d4da41389cca2f1feed57b2b0204 Mon Sep 17 00:00:00 2001
From: Hamza Abdelhedi <hamza.abdelhedii@gmail.com>
Date: Fri, 28 Aug 2026 18:20:32 +0200
Subject: [PATCH 4/8] ENH: Use standard rank handling for projection
 reconstruction

---
 doc/changes/dev/14235.newfeature.rst |  2 +-
 mne/_fiff/proj.py                    | 38 ++++++-------
 mne/forward/_field_interpolation.py  | 17 +++++-
 mne/tests/test_proj.py               | 80 ++++++++++++++++++++++------
 4 files changed, 99 insertions(+), 38 deletions(-)

diff --git a/doc/changes/dev/14235.newfeature.rst b/doc/changes/dev/14235.newfeature.rst
index edeb49875d6..c5062e3a7f8 100644
--- a/doc/changes/dev/14235.newfeature.rst
+++ b/doc/changes/dev/14235.newfeature.rst
@@ -1 +1 @@
-Projection reconstruction can now use an explicit :class:`~mne.Forward` model and integer spatial rank, by `Hamza Abdelhedi`_.
+Projection reconstruction can now use an explicit :class:`~mne.Forward` model and standard MNE rank handling, by `Hamza Abdelhedi`_.
diff --git a/mne/_fiff/proj.py b/mne/_fiff/proj.py
index aeda5c694c5..aa27a7a8495 100644
--- a/mne/_fiff/proj.py
+++ b/mne/_fiff/proj.py
@@ -13,7 +13,7 @@
 from ..fixes import _safe_svd
 from ..utils import (
     _check_option,
-    _ensure_int,
+    _check_rank,
     _validate_type,
     fill_doc,
     logger,
@@ -593,12 +593,13 @@ def reconstruct_proj(
         forward : instance of Forward | None
             Forward model used to construct the reconstruction field mapping.
             If ``None`` (default), use the geometry-based field mapping model.
-            If provided, ``rank`` must also be specified.
-        rank : int | None
-            Number of spatial modes to retain when reconstructing with
-            ``forward``. This is a sensor-space reconstruction rank, not a
-            number of sources or dipoles. Must be provided together with
-            ``forward``.
+        rank : None | 'info' | dict
+            Rank specification for Forward reconstruction. If ``None``,
+            estimate the rank from the projected Forward field covariance. If
+            ``'info'``, infer the rank from the measurement info and the
+            projectors applied during reconstruction. If a dict, explicitly
+            specify the rank by channel type, e.g. ``{'eeg': M}``. Used only
+            when ``forward`` is provided.
 
         Returns
         -------
@@ -607,22 +608,21 @@ def reconstruct_proj(
 
         Notes
         -----
-        When ``forward`` is provided, the reconstruction uses the sensor-space
-        field covariance formed from the Forward gain matrix. ``rank`` specifies
-        the number of sensor-space modes retained in the reconstruction.
+        When ``forward`` is provided, ``rank`` controls the sensor-space rank
+        of the projected Forward field covariance.
         """
         from ..forward import Forward, _map_meg_or_eeg_channels
 
-        if forward is None:
-            if rank is not None:
-                raise ValueError("rank can only be used when forward is provided")
-        else:
+        if forward is not None:
             _validate_type(forward, Forward, "forward")
-            if rank is None:
-                raise ValueError("rank must be provided when forward is provided")
-            rank = _ensure_int(rank, "rank")
-            if rank <= 0:
-                raise ValueError(f"rank must be positive, got {rank}")
+        rank = _check_rank(rank)
+        if forward is None and rank is not None:
+            raise ValueError("rank can only be used when forward is provided")
+        if forward is not None and rank == "full":
+            raise ValueError(
+                "rank='full' is incompatible with Forward reconstruction "
+                "after projection; use rank='info', None, or an explicit rank dict"
+            )
 
         if projs is None:
             if len(self.info["projs"]) == 0:
diff --git a/mne/forward/_field_interpolation.py b/mne/forward/_field_interpolation.py
index 2c631be67c3..bce98a0339a 100644
--- a/mne/forward/_field_interpolation.py
+++ b/mne/forward/_field_interpolation.py
@@ -15,10 +15,11 @@
 from .._fiff.pick import pick_channels_forward, pick_info, pick_types
 from .._fiff.proj import _has_eeg_average_ref_proj, make_projector
 from ..bem import _check_origin
-from ..cov import make_ad_hoc_cov
+from ..cov import Covariance, make_ad_hoc_cov
 from ..epochs import BaseEpochs, EpochsArray
 from ..evoked import Evoked, EvokedArray
 from ..fixes import _safe_svd
+from ..rank import _compute_rank_int
 from ..surface import get_head_surf, get_meg_helmet_surf
 from ..transforms import _find_trans, transform_surface_to
 from ..utils import (
@@ -143,6 +144,10 @@ def _map_meg_or_eeg_channels(
         Origin of the sphere in the head coordinate frame and in meters.
         Can be ``'auto'``, which means a head-digitization-based origin
         fit.
+    forward : instance of Forward | None
+        Forward model used instead of geometry-based field interpolation.
+    rank : None | 'info' | dict
+        Rank specification for Forward-based reconstruction.
 
     Returns
     -------
@@ -173,6 +178,14 @@ def _map_meg_or_eeg_channels(
         # Form the sensor-space field covariance from the Forward gain matrix.
         # As with any Gram representation, very weak modes can be numerically unstable.
         dots = lead_field @ lead_field.T
+        field_cov = Covariance(
+            dots,
+            info_from["ch_names"],
+            info_from["bads"],
+            info_from["projs"],
+            nfree=1,
+        )
+        rank_int = _compute_rank_int(field_cov, rank=rank, info=info_from)
         fmd = dict(
             kind=kind,
             ch_names=info_from["ch_names"],
@@ -180,7 +193,7 @@ def _map_meg_or_eeg_channels(
             self_dots=dots,
             surface_dots=dots,
         )
-        return _compute_mapping_matrix(fmd, info_from, rank=rank)
+        return _compute_mapping_matrix(fmd, info_from, rank=rank_int)
 
     assert origin is not None  # should be assured elsewhere
 
diff --git a/mne/tests/test_proj.py b/mne/tests/test_proj.py
index a04c36d06b9..4c703ff11d6 100644
--- a/mne/tests/test_proj.py
+++ b/mne/tests/test_proj.py
@@ -263,22 +263,70 @@ def test_reconstruct_proj_forward(raw_orig, eeg_forward):
         if average_ref:
             this_evoked.set_eeg_reference(projection=True)
         expected = _direct_forward_reconstruction(this_evoked, eeg_forward, rank)
-        got = this_evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank).data
+        got = (
+            this_evoked.copy()
+            .reconstruct_proj(forward=eeg_forward, rank={"eeg": 3})
+            .data
+        )
         assert_allclose(got, expected, rtol=1e-10, atol=1e-12)
 
     reordered = pick_channels_forward(
         eeg_forward, include=eeg_forward.ch_names[::-1], ordered=True
     )
-    got = evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank).data
-    got_reordered = evoked.copy().reconstruct_proj(forward=reordered, rank=rank).data
+    got = evoked.copy().reconstruct_proj(forward=eeg_forward, rank={"eeg": 3}).data
+    got_reordered = (
+        evoked.copy().reconstruct_proj(forward=reordered, rank={"eeg": 3}).data
+    )
     assert_allclose(got, got_reordered, rtol=1e-10, atol=1e-12)
 
     evoked_bad = evoked.copy()
     evoked_bad.info["bads"] = [evoked_bad.ch_names[0]]
-    got_bad = evoked_bad.reconstruct_proj(forward=eeg_forward, rank=rank).data
+    got_bad = evoked_bad.reconstruct_proj(forward=eeg_forward, rank={"eeg": 3}).data
     assert_allclose(got_bad[0], evoked_bad.data[0])
 
 
+def test_reconstruct_proj_forward_rank_info(eeg_forward):
+    """Test that rank='info' only uses projectors selected for reconstruction."""
+    ch_names = eeg_forward["info"]["ch_names"]
+    info = create_info(ch_names, 100.0, "eeg")
+    evoked = EvokedArray(np.arange(len(ch_names))[:, np.newaxis], info, tmin=0.0)
+    projs = [
+        _make_test_proj(ch_names, np.eye(len(ch_names))[ii], f"EEG-{ii}")
+        for ii in range(3)
+    ]
+    evoked.add_proj(projs, verbose=False)
+    selected = projs[:2]
+    nproj = make_projector(selected, ch_names)[1]
+    expected_rank = len(ch_names) - nproj
+    got = (
+        evoked.copy()
+        .reconstruct_proj(projs=selected, forward=eeg_forward, rank="info")
+        .data
+    )
+    expected = (
+        evoked.copy()
+        .reconstruct_proj(
+            projs=selected, forward=eeg_forward, rank={"eeg": expected_rank}
+        )
+        .data
+    )
+    assert_allclose(got, expected, rtol=1e-10, atol=1e-12)
+
+
+def test_reconstruct_proj_forward_rank_none(eeg_forward):
+    """Test automatic Forward reconstruction with the default rank."""
+    ch_names = eeg_forward["info"]["ch_names"]
+    info = create_info(ch_names, 100.0, "eeg")
+    evoked = EvokedArray(np.ones((len(ch_names), 1)), info, tmin=0.0)
+    evoked.add_proj(
+        _make_test_proj(ch_names, np.arange(1.0, len(ch_names) + 1.0), "EEG"),
+        verbose=False,
+    )
+    got = evoked.copy().reconstruct_proj(forward=eeg_forward).data
+    assert got.shape == evoked.data.shape
+    assert np.isfinite(got).all()
+
+
 def test_reconstruct_proj_forward_validation(eeg_forward):
     """Test validation of the explicit Forward reconstruction arguments."""
     info = create_info(eeg_forward["info"]["ch_names"], 100.0, "eeg")
@@ -292,20 +340,20 @@ def test_reconstruct_proj_forward_validation(eeg_forward):
         verbose=False,
     )
     with pytest.raises(ValueError, match="rank can only be used"):
-        evoked.copy().reconstruct_proj(rank=1)
-    with pytest.raises(ValueError, match="rank must be provided"):
-        evoked.copy().reconstruct_proj(forward=eeg_forward)
+        evoked.copy().reconstruct_proj(rank={"eeg": 1})
     with pytest.raises(TypeError, match="forward must be an instance of Forward"):
-        evoked.copy().reconstruct_proj(forward=[], rank=1)
-    for rank in (0, -1):
-        with pytest.raises(ValueError, match="rank must be positive"):
-            evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank)
-    for rank in (1.5, True):
-        with pytest.raises(TypeError, match="rank must be an int"):
-            evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank)
+        evoked.copy().reconstruct_proj(forward=[], rank={"eeg": 1})
+    invalid = evoked.copy()
+    with pytest.raises(TypeError, match="rank must be an instance of"):
+        invalid.reconstruct_proj(forward=eeg_forward, rank=1)
+    assert _active_projs(invalid) == [False]
+    with pytest.raises(ValueError, match="rank, if str"):
+        evoked.copy().reconstruct_proj(forward=eeg_forward, rank="bad")
+    with pytest.raises(ValueError, match="rank='full'"):
+        evoked.copy().reconstruct_proj(forward=eeg_forward, rank="full")
     rank = len(evoked.ch_names) + 1
-    with pytest.raises(ValueError, match="Invalid value for the rank parameter"):
-        evoked.copy().reconstruct_proj(forward=eeg_forward, rank=rank)
+    with pytest.raises(ValueError, match=r"rank\['eeg'\]=\d+ exceeds"):
+        evoked.copy().reconstruct_proj(forward=eeg_forward, rank={"eeg": rank})
 
 
 @pytest.mark.parametrize("kind", ["raw", "epochs", "evoked"])

From d504aaa2f0d2d3c0312ab393a4840790bdaf4b05 Mon Sep 17 00:00:00 2001
From: Hamza Abdelhedi <hamza.abdelhedii@gmail.com>
Date: Fri, 28 Aug 2026 18:28:00 +0200
Subject: [PATCH 5/8] TEST: Tighten Forward full-rank regression

---
 mne/tests/test_proj.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/mne/tests/test_proj.py b/mne/tests/test_proj.py
index 4c703ff11d6..0f65b1fff76 100644
--- a/mne/tests/test_proj.py
+++ b/mne/tests/test_proj.py
@@ -349,7 +349,7 @@ def test_reconstruct_proj_forward_validation(eeg_forward):
     assert _active_projs(invalid) == [False]
     with pytest.raises(ValueError, match="rank, if str"):
         evoked.copy().reconstruct_proj(forward=eeg_forward, rank="bad")
-    with pytest.raises(ValueError, match="rank='full'"):
+    with pytest.raises(ValueError, match="rank='full' is incompatible"):
         evoked.copy().reconstruct_proj(forward=eeg_forward, rank="full")
     rank = len(evoked.ch_names) + 1
     with pytest.raises(ValueError, match=r"rank\['eeg'\]=\d+ exceeds"):

From 7b5cfd7b8c0379b013ac9bbaa3d8ea5233415036 Mon Sep 17 00:00:00 2001
From: Hamza Abdelhedi <hamza.abdelhedii@gmail.com>
Date: Sat, 29 Aug 2026 00:23:27 +0200
Subject: [PATCH 6/8] MAINT: Add standard docs and verbosity to
 reconstruct_proj

---
 mne/_fiff/proj.py | 14 +++++++-------
 1 file changed, 7 insertions(+), 7 deletions(-)

diff --git a/mne/_fiff/proj.py b/mne/_fiff/proj.py
index aa27a7a8495..7b8f0e1e8a2 100644
--- a/mne/_fiff/proj.py
+++ b/mne/_fiff/proj.py
@@ -562,6 +562,8 @@ def plot_projs_topomap(
         )
         return fig
 
+    @fill_doc
+    @verbose
     def reconstruct_proj(
         self,
         *,
@@ -570,6 +572,7 @@ def reconstruct_proj(
         origin="auto",
         forward=None,
         rank=None,
+        verbose=None,
     ):
         """Apply SSP projectors and reconstruct the resulting signal in sensor space.
 
@@ -593,13 +596,10 @@ def reconstruct_proj(
         forward : instance of Forward | None
             Forward model used to construct the reconstruction field mapping.
             If ``None`` (default), use the geometry-based field mapping model.
-        rank : None | 'info' | dict
-            Rank specification for Forward reconstruction. If ``None``,
-            estimate the rank from the projected Forward field covariance. If
-            ``'info'``, infer the rank from the measurement info and the
-            projectors applied during reconstruction. If a dict, explicitly
-            specify the rank by channel type, e.g. ``{'eeg': M}``. Used only
-            when ``forward`` is provided.
+        %(rank_none)s
+            For Forward-based reconstruction, ``rank`` is used only when
+            ``forward`` is provided; ``rank='full'`` is not supported.
+        %(verbose)s
 
         Returns
         -------

From 4350e741445259518d13e8f1666354d964cf2fd8 Mon Sep 17 00:00:00 2001
From: Eric Larson <larson.eric.d@gmail.com>
Date: Sat, 29 Aug 2026 11:10:36 +0200
Subject: [PATCH 7/8] MAINT: Collapse _pinv_trunc into _reg_pinv

---
 mne/_fiff/proj.py                   | 10 ++++--
 mne/forward/_field_interpolation.py | 49 +++++++++++++----------------
 2 files changed, 28 insertions(+), 31 deletions(-)

diff --git a/mne/_fiff/proj.py b/mne/_fiff/proj.py
index 7b8f0e1e8a2..fd8328dd490 100644
--- a/mne/_fiff/proj.py
+++ b/mne/_fiff/proj.py
@@ -596,9 +596,13 @@ def reconstruct_proj(
         forward : instance of Forward | None
             Forward model used to construct the reconstruction field mapping.
             If ``None`` (default), use the geometry-based field mapping model.
-        %(rank_none)s
-            For Forward-based reconstruction, ``rank`` is used only when
-            ``forward`` is provided; ``rank='full'`` is not supported.
+        %(rank)s
+            Only used when ``forward`` is provided, where the default ``None``
+            estimates the rank of the projected field covariance. ``'full'`` is
+            not supported, as that covariance is rank-deficient after
+            projection. Without ``forward``, ``rank`` must be ``None``, and the
+            geometry-based field mapping uses its own internal truncation
+            rather than a rank estimated from the data.
         %(verbose)s
 
         Returns
diff --git a/mne/forward/_field_interpolation.py b/mne/forward/_field_interpolation.py
index bce98a0339a..f60f2c0f94d 100644
--- a/mne/forward/_field_interpolation.py
+++ b/mne/forward/_field_interpolation.py
@@ -18,7 +18,6 @@
 from ..cov import Covariance, make_ad_hoc_cov
 from ..epochs import BaseEpochs, EpochsArray
 from ..evoked import Evoked, EvokedArray
-from ..fixes import _safe_svd
 from ..rank import _compute_rank_int
 from ..surface import get_head_surf, get_meg_helmet_surf
 from ..transforms import _find_trans, transform_surface_to
@@ -65,14 +64,26 @@ def _compute_mapping_matrix(fmd, info, *, rank=None):
     whitener = np.diag(1.0 / np.sqrt(noise_cov["data"].ravel()))
     whitened_dots = np.dot(whitener.T, np.dot(proj_dots, whitener))
 
-    # SVD is numerically better than the eigenvalue composition even if
-    # mat is supposed to be symmetric and positive definite
-    if rank is not None:
-        inv, _, _ = _reg_pinv(whitened_dots, reg=0, rank=rank)
-    elif fmd.get("pinv_method", "tsvd") == "tsvd":
-        inv, fmd["nest"] = _pinv_trunc(whitened_dots, fmd["miss"])
+    # whitened_dots is symmetric and positive semi-definite, so _reg_pinv (which
+    # requires square Hermitian input) can do the truncated pseudoinversion
+    if fmd.get("pinv_method", "tsvd") == "tsvd":
+        n = len(whitened_dots)
+        if rank is None:
+            # truncate at most "miss" fraction of the singular value energy
+            s = np.linalg.svd(whitened_dots, compute_uv=False, hermitian=True)
+            varexp = np.cumsum(s)
+            varexp /= varexp[-1]
+            rank = np.where(varexp >= 1.0 - fmd["miss"])[0][0] + 1
+            logger.info(
+                f"    Truncating at {rank}/{n} components to omit less than "
+                f"{fmd['miss']:g} ({1.0 - varexp[rank - 1]:0.2g})"
+            )
+        else:
+            logger.info(f"    Truncating at {rank}/{n} components")
+        inv, _, fmd["nest"] = _reg_pinv(whitened_dots, reg=0, rank=rank)
     else:
         assert fmd["pinv_method"] == "tikhonov", fmd["pinv_method"]
+        assert rank is None, rank  # only the tsvd path supports an explicit rank
         inv, fmd["nest"] = _pinv_tikhonov(whitened_dots, fmd["miss"])
 
     # Sandwich with the whitener
@@ -96,26 +107,6 @@ def _compute_mapping_matrix(fmd, info, *, rank=None):
     return mapping_mat
 
 
-def _pinv_trunc(x, miss):
-    """Compute pseudoinverse, truncating at most "miss" fraction of varexp."""
-    u, s, v = _safe_svd(x, full_matrices=False)
-
-    # Eigenvalue truncation
-    varexp = np.cumsum(s)
-    varexp /= varexp[-1]
-    n = np.where(varexp >= (1.0 - miss))[0][0] + 1
-    logger.info(
-        "    Truncating at %d/%d components to omit less than %g (%0.2g)",
-        n,
-        len(s),
-        miss,
-        1.0 - varexp[n - 1],
-    )
-    s = 1.0 / s[:n]
-    inv = ((u[:, :n] * s) @ v[:n]).T
-    return inv, n
-
-
 def _pinv_tikhonov(x, reg):
     # _reg_pinv requires square Hermitian, which we have here
     inv, _, n = _reg_pinv(x, reg=reg, rank=None)
@@ -147,7 +138,9 @@ def _map_meg_or_eeg_channels(
     forward : instance of Forward | None
         Forward model used instead of geometry-based field interpolation.
     rank : None | 'info' | dict
-        Rank specification for Forward-based reconstruction.
+        Rank specification for Forward-based reconstruction, where ``None``
+        estimates it from the projected field covariance. Must be ``None`` for
+        the geometry-based path, which uses its own truncation instead.
 
     Returns
     -------

From 7927d59cd01ea0896b6d3a16f04110e4ac23e780 Mon Sep 17 00:00:00 2001
From: Eric Larson <larson.eric.d@gmail.com>
Date: Sat, 29 Aug 2026 13:38:46 +0200
Subject: [PATCH 8/8] FIX: Better

---
 mne/_fiff/proj.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/mne/_fiff/proj.py b/mne/_fiff/proj.py
index 0b47f96ab39..0f557cf90b0 100644
--- a/mne/_fiff/proj.py
+++ b/mne/_fiff/proj.py
@@ -562,7 +562,6 @@ def plot_projs_topomap(
         )
         return fig
 
-    @fill_doc
     @verbose
     def reconstruct_proj(
         self,
@@ -597,6 +596,7 @@ def reconstruct_proj(
             Forward model used to construct the reconstruction field mapping.
             If ``None`` (default), use the geometry-based field mapping model.
         %(rank)s
+
             Only used when ``forward`` is provided, where the default ``None``
             estimates the rank of the projected field covariance. ``'full'`` is
             not supported, as that covariance is rank-deficient after