Fix: update raw attributes when appending cHPI channels

doc/changes/dev/14180.bugfix.rst

@@ -0,0 +1 @@
+Fix a bug that caused :meth:`mne.io.Raw.pick` and :meth:`mne.io.Raw.drop_channels` to fail on the result of :func:`mne.preprocessing.maxwell_filter` when movement compensation was enabled with ``head_pos``, by `Christian Brodbeck`_.

mne/preprocessing/maxwell.py

@@ -312,6 +312,11 @@ def maxwell_filter(
     -----
     .. versionadded:: 0.11
 
+    When ``head_pos`` is provided, the returned object contains cHPI result
+    channels that have no counterpart in the source file. It therefore cannot
+    be concatenated using ``raw_sss.append(..., preload=False)``. Use
+    ``preload=True`` or a memory-mapped filename instead.
+
     Some of this code was adapted and relicensed (with BSD form) with
     permission from Jussi Nurminen. These algorithms are based on work
     from :footcite:`TauluKajola2005` and :footcite:`TauluSimola2006`.
@@ -1258,6 +1263,20 @@ def _copy_preload_add_channels(raw, add_channels, copy, info):
         raw.info["chs"].extend(chpi_chs)
         raw.info._update_redundant()
         raw.info._check_consistency()
+        # The remaining per-channel attributes must grow along with info, otherwise
+        # any later channel operation (e.g., raw_sss.drop_channels) indexes them out
+        # of bounds
+        raw._cals = np.concatenate([raw._cals, raw.info._cals[off:]])
+        # The added channels have no counterpart in the source file, but the data are
+        # preloaded, so _read_picks (and _raw_extras) will never be used to read from
+        # disk again -- use indices that are likely to break loudly if they ever are
+        extra_idx = [2147483647] * len(chpi_chs)  # 2 ** 31 - 1
+        raw._read_picks = [np.concatenate([r, extra_idx]) for r in raw._read_picks]
+        assert raw._comp is None  # preloading the data above unsets it
+        if raw._projector is not None:  # identity for the added channels
+            projector = np.eye(raw.info["nchan"])
+            projector[:off, :off] = raw._projector
+            raw._projector = projector
         assert raw._data.shape == (raw.info["nchan"], len(raw.times))
         # Return the pos picks
         pos_picks = np.arange(len(raw.ch_names) - len(chpi_chs), len(raw.ch_names))

mne/preprocessing/tests/test_maxwell.py

@@ -256,6 +256,16 @@ def test_movement_compensation_basic(tmp_path):
     assert_meg_snr(
         raw_sss, read_crop(sss_movecomp_fname, lims), 4.6, 12.4, chpi_med_tol=58
     )
+    # appending the head position channels must leave the instance consistent, i.e.,
+    # picking and dropping channels must still work (and agree with info)
+    want = raw_sss.ch_names[:2] + ["CHPI002"]
+    raw_sss_pick = raw_sss.copy().pick(want)
+    assert raw_sss_pick.ch_names == want
+    assert len(raw_sss_pick._cals) == len(raw_sss_pick._read_picks[0]) == len(want)
+    assert_allclose(raw_sss_pick.get_data(), raw_sss.get_data(want))
+    raw_sss_drop = raw_sss.copy().drop_channels(["CHPI000"])
+    assert raw_sss_drop.ch_names == [c for c in raw_sss.ch_names if c != "CHPI000"]
+    assert_allclose(raw_sss_drop.get_data(), raw_sss.get_data(raw_sss_drop.ch_names))
     # IO
     temp_fname = tmp_path / "test_raw_sss.fif"
     raw_sss.save(temp_fname)