MRG, ENH: Add forward sensitivity maps to reports

doc/changes/dev/14107.newfeature.rst

@@ -0,0 +1 @@
+Added the ``sensitivity`` parameter to :meth:`mne.Report.add_forward` to include forward sensitivity maps in reports, by :newcontrib:`Mariam Husain`.

doc/changes/names.inc

@@ -213,6 +213,7 @@
 .. _Mainak Jas: https://jasmainak.github.io
 .. _Maksym Balatsko: https://github.com/mbalatsko
 .. _Marcin Koculak: https://github.com/mkoculak
+.. _Mariam Husain: https://github.com/mariam-hedgie
 .. _Marian Dovgialo: https://github.com/mdovgialo
 .. _Marijn van Vliet: https://github.com/wmvanvliet
 .. _Mark Alexander Henney: https://github.com/henneysq

mne/report/report.py

@@ -41,8 +41,8 @@
 from ..minimum_norm import InverseOperator, read_inverse_operator
 from ..parallel import parallel_func
 from ..preprocessing.ica import read_ica
-from ..proj import read_proj
-from ..source_estimate import SourceEstimate, read_source_estimate
+from ..proj import read_proj, sensitivity_map
+from ..source_estimate import _BaseSourceEstimate, read_source_estimate
 from ..source_space._source_space import _ensure_src
 from ..surface import dig_mri_distances
 from ..transforms import _find_trans
@@ -1583,6 +1583,7 @@ def add_forward(
         subject=None,
         subjects_dir=None,
         plot=False,
+        sensitivity=False,
         tags=("forward-solution",),
         section=None,
         replace=False,
@@ -1604,6 +1605,10 @@ def add_forward(
             If True, plot the source space of the forward solution.
 
             .. versionadded:: 1.10
+        sensitivity : bool
+            If True, render sensitivity maps for all available sensor types.
+
+            .. versionadded:: 1.13
         %(tags_report)s
         %(section_report)s
 
@@ -1626,6 +1631,7 @@ def add_forward(
             tags=tags,
             replace=replace,
             plot=plot,
+            sensitivity=sensitivity,
         )
 
     @fill_doc
@@ -3709,6 +3715,7 @@ def _add_forward(
         title,
         image_format,
         plot,
+        sensitivity,
         section,
         tags,
         replace,
@@ -3720,9 +3727,28 @@ def _add_forward(
         subject = self.subject if subject is None else subject
         subject = forward["src"][0]["subject_his_id"] if subject is None else subject
 
-        # XXX Todo
-        # Render sensitivity maps
         sensitivity_maps_html = ""
+        if sensitivity:
+            ch_types = forward["info"].get_channel_types(unique=True)
+            for ch_type in ("grad", "mag", "eeg"):
+                if ch_type not in ch_types:
+                    continue
+                stc = sensitivity_map(forward, ch_type=ch_type)
+                html_partial = self._render_stc(
+                    stc=stc,
+                    title=f"{ch_type.upper()} sensitivity",
+                    subject=subject,
+                    subjects_dir=subjects_dir,
+                    n_time_points=1,
+                    image_format=image_format,
+                    tags=tags,
+                    stc_plot_kwargs=None,
+                )
+                sensitivity_maps_html += html_partial(
+                    id_=self._get_dom_id(
+                        section=section, title=f"{title}-{ch_type}-sensitivity"
+                    )
+                )
         source_space_html = ""
         if plot:
             source_space_html = self._src_html(
@@ -4508,7 +4534,38 @@ def _add_stc(
         replace,
     ):
         """Render STC."""
-        if isinstance(stc, SourceEstimate):
+        html_partial = self._render_stc(
+            stc=stc,
+            title=title,
+            subject=subject,
+            subjects_dir=subjects_dir,
+            n_time_points=n_time_points,
+            image_format=image_format,
+            tags=tags,
+            stc_plot_kwargs=stc_plot_kwargs,
+        )
+        self._add_or_replace(
+            title=title,
+            section=section,
+            tags=tags,
+            html_partial=html_partial,
+            replace=replace,
+        )
+
+    def _render_stc(
+        self,
+        *,
+        stc,
+        title,
+        subject,
+        subjects_dir,
+        n_time_points,
+        image_format,
+        tags,
+        stc_plot_kwargs,
+    ):
+        """Render an STC as embeddable report HTML."""
+        if isinstance(stc, _BaseSourceEstimate):
             if subject is None:
                 subject = self.subject  # supplied during Report init
                 if not subject:
@@ -4543,96 +4600,40 @@ def _add_stc(
             )
         t_zero_idx = np.abs(times).argmin()  # index of time closest to zero
 
-        # Plot using 3d backend if available, and use Matplotlib
-        # otherwise.
-        import matplotlib.pyplot as plt
+        if get_3d_backend() is None:
+            raise RuntimeError(
+                "A 3D backend is required to render source estimates in a report."
+            )
 
         stc_plot_kwargs = _handle_default("report_stc_plot_kwargs", stc_plot_kwargs)
         stc_plot_kwargs.update(subject=subject, subjects_dir=subjects_dir)
-        # we need to set the size based on the min (img_max_width can be None)
         if self.img_max_width is not None:
             stc_plot_kwargs["size"] = (
                 stc_plot_kwargs["size"][0],
                 min(stc_plot_kwargs["size"][1], self.img_max_width),
             )
-        if get_3d_backend() is not None:
-            brain = stc.plot(**stc_plot_kwargs)
-            brain._renderer.plotter.subplot(0, 0)
-            backend_is_3d = True
-        else:
-            backend_is_3d = False
+
+        brain = stc.plot(**stc_plot_kwargs)
+        brain._renderer.plotter.subplot(0, 0)
 
         figs = []
         for t in times:
-            with warnings.catch_warnings():
-                warnings.filterwarnings(
-                    action="ignore",
-                    message="More than 20 figures have been opened",
-                    category=RuntimeWarning,
-                )
-
-                if backend_is_3d:
-                    brain.set_time(t)
-                    figs.append(brain.screenshot(time_viewer=True, mode="rgb"))
-                else:
-                    fig_lh = plt.figure(layout="constrained")
-                    fig_rh = plt.figure(layout="constrained")
-
-                    brain_lh = stc.plot(
-                        views="lat",
-                        hemi="lh",
-                        initial_time=t,
-                        backend="matplotlib",
-                        subject=subject,
-                        subjects_dir=subjects_dir,
-                        figure=fig_lh,
-                    )
-                    brain_rh = stc.plot(
-                        views="lat",
-                        hemi="rh",
-                        initial_time=t,
-                        subject=subject,
-                        subjects_dir=subjects_dir,
-                        backend="matplotlib",
-                        figure=fig_rh,
-                    )
-                    _constrain_fig_resolution(
-                        fig_lh,
-                        max_width=stc_plot_kwargs["size"][0],
-                        max_res=self.img_max_res,
-                    )
-                    _constrain_fig_resolution(
-                        fig_rh,
-                        max_width=stc_plot_kwargs["size"][0],
-                        max_res=self.img_max_res,
-                    )
-                    figs.append(brain_lh)
-                    figs.append(brain_rh)
-                    plt.close(fig_lh)
-                    plt.close(fig_rh)
-
-        if backend_is_3d:
-            brain.close()
-        else:
-            brain_lh.close()
-            brain_rh.close()
+            brain.set_time(t)
+            figs.append(brain.screenshot(time_viewer=True, mode="rgb"))
+        brain.close()
 
         captions = [f"Time point: {round(t, 3):0.3f} s" for t in times]
-        self._add_slider(
+        return self._render_slider(
             figs=figs,
             imgs=None,
             captions=captions,
             title=title,
             image_format=image_format,
             start_idx=t_zero_idx,
-            section=section,
             tags=tags,
-            replace=replace,
-            own_figure=False,  # prevent rescaling
+            klass="stc",
+            own_figure=False,
         )
-        for fig in figs:
-            if not isinstance(fig, np.ndarray):
-                plt.close(fig)
 
     @_use_agg
     def _add_bem(

mne/report/tests/test_report.py

@@ -4,6 +4,7 @@
 
 import base64
 import glob
+import inspect
 import os
 import pickle
 import re
@@ -35,6 +36,7 @@
     _ALLOWED_IMAGE_FORMATS,
     CONTENT_ORDER,
 )
+from mne.source_estimate import SourceEstimate, VolSourceEstimate
 from mne.utils import Bunch, _record_warnings
 from mne.utils._testing import assert_object_equal
 from mne.viz import plot_alignment
@@ -540,6 +542,100 @@ def test_add_forward(renderer_interactive_pyvistaqt):
     assert len(report.html) == 1
     assert report.html[0].count("<img") == 0
 
+    report = Report(subjects_dir=subjects_dir, image_format="png")
+    report.add_forward(
+        forward=fwd_fname,
+        subjects_dir=subjects_dir,
+        title="Forward solution",
+        plot=False,
+        sensitivity=True,
+    )
+    assert len(report.html) == 1
+    assert report.html[0].count("<img") == 3
+
+
+def test_add_forward_sensitivity_parameter():
+    """Test sensitivity maps can be requested when adding a forward solution."""
+    assert "sensitivity" in inspect.signature(Report.add_forward).parameters
+
+
+class _FakeBrain:
+    """Minimal stand-in for a 3D source-estimate plot."""
+
+    def __init__(self):
+        self._renderer = Bunch(plotter=Bunch(subplot=lambda *_: None))
+
+    def close(self):
+        pass
+
+    def screenshot(self, **kwargs):
+        return np.zeros((2, 2, 3), np.uint8)
+
+    def set_time(self, time):
+        pass
+
+
+def _fake_stc_plot(*args, **kwargs):
+    """Return a lightweight source-estimate plot for report rendering tests."""
+    return _FakeBrain()
+
+
+def test_render_volume_stc(monkeypatch):
+    """Test rendering an in-memory volume source estimate."""
+    stc = VolSourceEstimate(
+        data=np.ones((1, 1)),
+        vertices=[np.array([0])],
+        tmin=0,
+        tstep=1,
+        subject="sample",
+    )
+    report = Report()
+    monkeypatch.setattr(report_mod, "get_3d_backend", lambda: "pyvista")
+    monkeypatch.setattr(VolSourceEstimate, "plot", _fake_stc_plot)
+    monkeypatch.setattr(
+        report_mod,
+        "read_source_estimate",
+        lambda **kwargs: pytest.fail("An in-memory STC must not be read from disk"),
+    )
+    html_partial = report._render_stc(
+        stc=stc,
+        title="Volume STC",
+        subject="sample",
+        subjects_dir=None,
+        n_time_points=1,
+        image_format="png",
+        tags=(),
+        stc_plot_kwargs=None,
+    )
+    assert callable(html_partial)
+
+
+def test_render_stc_requires_3d_backend(monkeypatch):
+    """Test rendering source estimates requires a 3D backend."""
+    stc = SourceEstimate(
+        data=np.ones((2, 3)),
+        vertices=[np.array([0]), np.array([0])],
+        tmin=-1,
+        tstep=1,
+        subject="sample",
+    )
+    report = Report()
+    monkeypatch.setattr(report_mod, "get_3d_backend", lambda: None)
+    monkeypatch.setattr(SourceEstimate, "plot", _fake_stc_plot)
+    monkeypatch.setattr(report, "_render_slider", lambda **kwargs: None)
+
+    with pytest.raises(RuntimeError, match="3D backend"):
+        report._render_stc(
+            stc=stc,
+            title="Surface STC",
+            subject="sample",
+            subjects_dir=None,
+            n_time_points=3,
+            image_format="png",
+            tags=(),
+            stc_plot_kwargs=None,
+        )
+
 
 @testing.requires_testing_data
 def test_render_mri_without_bem(tmp_path):