RF: Rework get_matched_empty_room()

doc/whats_new.rst

@@ -24,10 +24,11 @@ Changelog
 
 - :func:`read_raw_bids` now reads in participants tsv data, by `Adam Li`_ (`#392 <https://github.com/mne-tools/mne-bids/pull/392>`_)
 - :func:`mne_bids.utils.get_entity_vals` has gained ``ignore_*`` keyword arguments to exclude specific values from the list of results, e.g. the entities for a particular subject or task, by `Richard Höchenberger`_ (`#409 <https://github.com/mne-tools/mne-bids/pull/409>`_)
-- :func:`mne_bids.write_raw_bids` now uses the 'space' BIDS-entity when writing iEEG electrodes and coordinate frames by `Adam Li`_ (`#390 <https://github.com/mne-tools/mne-bids/pull/390>`_)
-- :code:`convert_ieeg_to_bids` to now use sample ECoG EDF data by `Adam Li`_ (`#390 <https://github.com/mne-tools/mne-bids/pull/390>`_)
-- :func:`mne_bids.write_raw_bids` now writes CoordinateSystemDescription as specified in BIDS Specification if CoordinateSystem is MNE-compatible by `Adam Li`_ (`#416 <https://github.com/mne-tools/mne-bids/pull/416>`_)
-- :func:`mne_bids.write_raw_bids` and :func:`mne_bids.read_raw_bids` now handles scalp EEG if Captrak coordinate system and NAS/LPA/RPA landmarks are present by `Adam Li`_ (`#416 <https://github.com/mne-tools/mne-bids/pull/416>`_)
+- :func:`mne_bids.write_raw_bids` now uses the 'space' BIDS-entity when writing iEEG electrodes and coordinate frames, by `Adam Li`_ (`#390 <https://github.com/mne-tools/mne-bids/pull/390>`_)
+- :code:`convert_ieeg_to_bids` to now use sample ECoG EDF data, by `Adam Li`_ (`#390 <https://github.com/mne-tools/mne-bids/pull/390>`_)
+- :func:`mne_bids.write_raw_bids` now writes CoordinateSystemDescription as specified in BIDS Specification if CoordinateSystem is MNE-compatible, by `Adam Li`_ (`#416 <https://github.com/mne-tools/mne-bids/pull/416>`_)
+- :func:`mne_bids.write_raw_bids` and :func:`mne_bids.read_raw_bids` now handle scalp EEG if Captrak coordinate system and NAS/LPA/RPA landmarks are present, by `Adam Li`_ (`#416 <https://github.com/mne-tools/mne-bids/pull/416>`_)
+- :func:`mne_bids.get_matched_empty_room` now implements an algorithm for discovering empty-room recordings that do not have the recording date set as their session, by `Richard Höchenberger`_ (`#421 <https://github.com/mne-tools/mne-bids/pull/421>`_)
 
 Bug
 ~~~

mne_bids/read.py

@@ -10,6 +10,7 @@
 from datetime import datetime
 import glob
 import json
+import pathlib
 
 import numpy as np
 import mne
@@ -25,8 +26,8 @@
 from mne_bids.utils import (_parse_bids_filename, _extract_landmarks,
                             _find_matching_sidecar, _parse_ext,
                             _get_ch_type_mapping, make_bids_folders,
-                            make_bids_basename, _estimate_line_freq,
-                            _get_kinds_for_sub)
+                            _gen_bids_basename,
+                            _estimate_line_freq, _get_kinds_for_sub)
 
 reader = {'.con': io.read_raw_kit, '.sqd': io.read_raw_kit,
           '.fif': io.read_raw_fif, '.pdf': io.read_raw_bti,
@@ -522,53 +523,128 @@ def get_matched_empty_room(bids_basename, bids_root):
         The basename corresponding to the best-matching empty-room measurement.
         Returns None if none was found.
     """
-    kind = 'meg'
+    kind = 'meg'  # We're only concerned about MEG data here
     bids_fname = _make_bids_fname(bids_basename=bids_basename,
                                   bids_root=bids_root, kind=kind)
     _, ext = _parse_ext(bids_fname)
+    if ext == '.fif':
+        extra_params = dict(allow_maxshield=True)
+    else:
+        extra_params = None
 
     raw = read_raw_bids(bids_basename=bids_basename, bids_root=bids_root,
-                        kind='meg')
+                        kind=kind, extra_params=extra_params)
     if raw.info['meas_date'] is None:
-        raise ValueError('Measurement date not available. Cannot get matching'
-                         ' empty room file')
+        raise ValueError('The provided recording does not have a measurement '
+                         'date set. Cannot get matching empty-room file.')
 
     ref_date = raw.info['meas_date']
     if not isinstance(ref_date, datetime):
         # for MNE < v0.20
         ref_date = datetime.fromtimestamp(raw.info['meas_date'][0])
-    search_path = make_bids_folders(bids_root=bids_root, subject='emptyroom',
-                                    session='**', make_dir=False)
-    search_path = op.join(search_path, '**', '**%s' % ext)
-    er_fnames = glob.glob(search_path)
-
-    best_er_fname = None
-    min_seconds = np.inf
-    for er_fname in er_fnames:
+
+    emptyroom_dir = pathlib.Path(make_bids_folders(bids_root=bids_root,
+                                                   subject='emptyroom',
+                                                   make_dir=False))
+
+    if not emptyroom_dir.exists():
+        return None
+
+    # Find the empty-room recording sessions.
+    emptyroom_session_dirs = [x for x in emptyroom_dir.iterdir()
+                              if x.is_dir() and str(x.name).startswith('ses-')]
+    if not emptyroom_session_dirs:  # No session sub-directories found
+        emptyroom_session_dirs = [emptyroom_dir]
+
+    # Now try to discover all recordings inside the session directories.
+
+    allowed_extensions = list(reader.keys())
+    # `.pdf` is just a "virtual" extension for BTi data (which is stored inside
+    # a dedicated directory that doesn't have an extension)
+    del allowed_extensions[allowed_extensions.index('.pdf')]
+
+    candidate_er_fnames = []
+    for session_dir in emptyroom_session_dirs:
+        dir_contents = glob.glob(op.join(session_dir, kind,
+                                         f'sub-emptyroom_*_{kind}*'))
+        for item in dir_contents:
+            item = pathlib.Path(item)
+            if ((item.suffix in allowed_extensions) or
+                    (not item.suffix and item.is_dir())):  # Hopefully BTi?
+                candidate_er_fnames.append(item.name)
+
+    # Walk through recordings, trying to extract the recording date:
+    # First, from the filename; and if that fails, from `info['meas_date']`.
+    best_er_basename = None
+    min_delta_t = np.inf
+    date_tie = False
+
+    failed_to_get_er_date_count = 0
+    for er_fname in candidate_er_fnames:
         params = _parse_bids_filename(er_fname, verbose=False)
-        dt = datetime.strptime(params['ses'], '%Y%m%d')
-        dt = dt.replace(tzinfo=ref_date.tzinfo)
-        delta_t = dt - ref_date
-        if abs(delta_t.total_seconds()) < min_seconds:
-            min_seconds = abs(delta_t.total_seconds())
-            best_er_fname = er_fname
-
-    if best_er_fname is None:
-        er_basename = None
-    else:
-        params = _parse_bids_filename(best_er_fname, verbose='warning')
-        er_basename = make_bids_basename(
-            subject=params.get('sub', None),
+        er_meas_date = None
+
+        er_basename = _gen_bids_basename(
+            subject='emptyroom',
             session=params.get('ses', None),
             task=params.get('task', None),
             acquisition=params.get('acq', None),
             run=params.get('run', None),
             processing=params.get('proc', None),
             recording=params.get('recording', None),
-            space=params.get('space', None)
+            space=params.get('space', None),
+            # BIDS specification does not enforce use of ses-YYYYMMDD and
+            # task-emptyroom entities.
+            on_invalid_er_session='continue',
+            on_invalid_er_task='continue'
         )
 
-    return er_basename
+        if 'ses' in params:  # Try to extract date from filename.
+            try:
+                er_meas_date = datetime.strptime(params['ses'], '%Y%m%d')
+            except (ValueError, TypeError):
+                # There is a session in the filename, but it doesn't encode a
+                # valid date.
+                pass
+
+        if er_meas_date is None:  # No luck so far! Check info['meas_date']
+            _, ext = _parse_ext(er_fname)
+            if ext == '.fif':
+                extra_params = dict(allow_maxshield=True)
+            else:
+                extra_params = None
+
+            er_raw = read_raw_bids(bids_basename=er_basename,
+                                   bids_root=bids_root,
+                                   kind=kind,
+                                   extra_params=extra_params)
+
+            er_meas_date = er_raw.info['meas_date']
+            if er_meas_date is None:  # There's nothing we can do.
+                failed_to_get_er_date_count += 1
+                continue
+
+        er_meas_date = er_meas_date.replace(tzinfo=ref_date.tzinfo)
+        delta_t = er_meas_date - ref_date
+
+        if abs(delta_t.total_seconds()) == min_delta_t:
+            date_tie = True
+        elif abs(delta_t.total_seconds()) < min_delta_t:
+            min_delta_t = abs(delta_t.total_seconds())
+            best_er_basename = er_basename
+            date_tie = False
+
+    if failed_to_get_er_date_count > 0:
+        msg = (f'Could not retrieve the empty-room measurement date from '
+               f'a total of {failed_to_get_er_date_count} recording(s).')
+        warn(msg)
+
+    if date_tie:
+        msg = ('Found more than one matching empty-room measurement with the '
+               'same recording date. Selecting the first match.')
+        warn(msg)
+
+    return best_er_basename
 
 
 def get_head_mri_trans(bids_basename, bids_root):

mne_bids/tests/test_read.py

@@ -647,7 +647,8 @@ def test_get_matched_empty_room():
         raw.annotations.orig_time = None
     anonymize_info(raw.info)
     write_raw_bids(raw, bids_basename, bids_root, overwrite=True)
-    with pytest.raises(ValueError, match='Measurement date not available'):
+    with pytest.raises(ValueError, match='The provided recording does not '
+                                         'have a measurement date set'):
         get_matched_empty_room(bids_basename=bids_basename,
                                bids_root=bids_root)
 

mne_bids/utils.py

@@ -804,8 +804,8 @@ def _gen_bids_basename(*, subject=None, session=None, task=None,
                        recording=None, space=None, prefix=None, suffix=None,
                        on_invalid_er_session='raise',
                        on_invalid_er_task='raise'):
-    if on_invalid_er_session not in ['raise', 'warn']:
-        msg = (f'on_invalid_er_session must be raise or warn, '
+    if on_invalid_er_session not in ['raise', 'warn', 'continue']:
+        msg = (f'on_invalid_er_session must be raise, warn, or continue, '
                f'but received: {on_invalid_er_session}')
         raise ValueError(msg)
 
@@ -834,25 +834,28 @@ def _gen_bids_basename(*, subject=None, session=None, task=None,
         raise ValueError("At least one parameter must be given.")
 
     if subject == 'emptyroom':
-        if task != 'noise' and on_invalid_er_task != 'continue':
+        if task != 'noise':
             msg = (f'task must be "noise" if subject is "emptyroom", but '
                    f'received: {task}')
             if on_invalid_er_task == 'raise':
                 raise ValueError(msg)
-            else:
+            elif on_invalid_er_task == 'warn':
                 logger.critical(msg)
-
+            else:
+                pass
         try:
             datetime.strptime(session, '%Y%m%d')
         except (ValueError, TypeError):
             msg = (f'empty-room session should be a string of format '
                    f'YYYYMMDD, but received: {session}')
             if on_invalid_er_session == 'raise':
                 raise ValueError(msg)
-            else:
+            elif on_invalid_er_session == 'warn':
                 msg = (f'{msg}. Will proceed anyway, but you should consider '
                        f'fixing your dataset.')
                 logger.critical(msg)
+            else:
+                pass
 
     basename = []
     for key, val in order.items():
@@ -869,6 +872,8 @@ def _gen_bids_basename(*, subject=None, session=None, task=None,
 
     return basename
 
+    return basename
+
 
 def make_bids_basename(subject=None, session=None, task=None,
                        acquisition=None, run=None, processing=None,

mne_bids/write.py

@@ -987,12 +987,12 @@ def write_raw_bids(raw, bids_basename, bids_root, events_data=None,
         subject=subject_id, session=session_id, acquisition=acquisition,
         suffix='electrodes.tsv', prefix=data_path,
         on_invalid_er_task=on_invalid_er_task)
+
     # For the remaining files, we can use make_bids_basename() as usual.
     participants_tsv_fname = make_bids_basename(prefix=bids_root,
                                                 suffix='participants.tsv')
     participants_json_fname = make_bids_basename(prefix=bids_root,
                                                  suffix='participants.json')
-
     sidecar_fname = make_bids_basename(
         subject=subject_id, session=session_id, task=task, run=run,
         acquisition=acquisition, suffix='%s.json' % kind, prefix=data_path)
@@ -1003,6 +1003,7 @@ def write_raw_bids(raw, bids_basename, bids_root, events_data=None,
     channels_fname = make_bids_basename(
         subject=subject_id, session=session_id, task=task, run=run,
         acquisition=acquisition, suffix='channels.tsv', prefix=data_path)
+
     if ext not in ['.fif', '.ds', '.vhdr', '.edf', '.bdf', '.set', '.con',
                    '.sqd']:
         bids_raw_folder = bids_fname.split('.')[0]