BUG: Raise error for EDF+D/BDF+D files with acquisition gaps

Detect and raise NotImplementedError when loading EDF+D or BDF+D files
that contain actual gaps between data records, instead of silently
treating them as continuous data. This prevents incorrect time alignment.

Fixes #13429

Author: larsoner

doc/changes/dev/13583.bugfix.rst

@@ -0,0 +1 @@
+Raise a :class:`NotImplementedError` when reading EDF+D or BDF+D files with acquisition gaps instead of silently loading them as continuous data, by `Arnav Kumar`_ (:gh:`13583`).

doc/changes/names.inc

@@ -29,6 +29,7 @@
 .. _Antti Rantala: https://github.com/Odingod
 .. _Apoorva Karekal: https://github.com/apoorva6262
 .. _Archit Singhal: https://github.com/architsinghal-mriirs
+.. _Arnav Kumar: https://github.com/Arnav1709
 .. _Arne Pelzer: https://github.com/aplzr
 .. _Ashley Drew: https://github.com/ashdrew
 .. _Asish Panda: https://github.com/kaichogami

mne/io/edf/edf.py

@@ -235,6 +235,30 @@ def __init__(
                 np.ones((len(idx), 1)),
                 None,
             )
+
+            # Check for discontinuous EDF+D files with actual gaps
+            if edf_info.get("discontinuous", False):
+                record_times = _get_tal_record_times(tal_data[0], encoding=encoding)
+                record_length = edf_info["record_length"][0]
+                n_records = edf_info["n_records"]
+                has_gaps, gaps = _check_edf_discontinuity(
+                    record_times, record_length, n_records
+                )
+                if has_gaps:
+                    gap_info = ", ".join(
+                        [f"{onset:.3f}s (duration: {dur:.3f}s)" for onset, dur in gaps]
+                    )
+                    raise NotImplementedError(
+                        "This EDF+D file contains discontinuous data with gaps "
+                        f"between records. Gaps found at: {gap_info}. "
+                        "MNE-Python does not currently support reading EDF+D files "
+                        "with acquisition gaps. The data would be incorrectly "
+                        "treated as continuous, leading to incorrect time alignment. "
+                        "Consider using specialized tools like luna/lunapi for "
+                        "discontinuous EDF+ files, or convert the file to EDF+C "
+                        "format if the gaps are not significant for your analysis."
+                    )
+
             annotations = _read_annotations_edf(
                 tal_data[0],
                 ch_names=info["ch_names"],
@@ -447,6 +471,30 @@ def __init__(
                 np.ones((len(idx), 1)),
                 None,
             )
+
+            # Check for discontinuous BDF+D files with actual gaps
+            if edf_info.get("discontinuous", False):
+                record_times = _get_tal_record_times(tal_data[0], encoding=encoding)
+                record_length = edf_info["record_length"][0]
+                n_records = edf_info["n_records"]
+                has_gaps, gaps = _check_edf_discontinuity(
+                    record_times, record_length, n_records
+                )
+                if has_gaps:
+                    gap_info = ", ".join(
+                        [f"{onset:.3f}s (duration: {dur:.3f}s)" for onset, dur in gaps]
+                    )
+                    raise NotImplementedError(
+                        "This BDF+D file contains discontinuous data with gaps "
+                        f"between records. Gaps found at: {gap_info}. "
+                        "MNE-Python does not currently support reading BDF+D files "
+                        "with acquisition gaps. The data would be incorrectly "
+                        "treated as continuous, leading to incorrect time alignment. "
+                        "Consider using specialized tools like luna/lunapi for "
+                        "discontinuous BDF+ files, or convert the file to BDF+C "
+                        "format if the gaps are not significant for your analysis."
+                    )
+
             annotations = _read_annotations_edf(
                 tal_data[0],
                 ch_names=info["ch_names"],
@@ -1159,9 +1207,16 @@ def _read_edf_header(
         # to determine the subtype (EDF or BDF, which differ in the
         # number of bytes they use for the data records; EDF uses 2 bytes
         # whereas BDF uses 3 bytes).
-        fid.read(44)
+        # However, we still need to check for EDF+D/BDF+D (discontinuous) files.
+        reserved = fid.read(44).decode("latin-1").rstrip()
         subtype = file_type
 
+        # Check for discontinuous EDF+D/BDF+D files
+        if reserved in ("EDF+D", "BDF+D"):
+            edf_info["discontinuous"] = True
+        else:
+            edf_info["discontinuous"] = False
+
         n_records = int(_edf_str(fid.read(8)))
         record_length = float(_edf_str(fid.read(8)))
         record_length = np.array([record_length, 1.0])  # in seconds
@@ -2005,6 +2060,11 @@ def read_raw_edf(
 
     The EDF specification allows storage of subseconds in measurement date.
     However, this reader currently sets subseconds to 0 by default.
+
+    EDF+D (discontinuous) files with actual gaps between data records are not
+    supported and will raise a :class:`NotImplementedError`. EDF+D files that
+    are marked as discontinuous but have no actual gaps (e.g., from some
+    Nihon Kohden systems) will load normally.
     """
     _check_args(input_fname, preload, "edf")
 
@@ -2144,6 +2204,10 @@ def read_raw_bdf(
     If channels named 'status' or 'trigger' are present, they are considered as
     STIM channels by default. Use func:`mne.find_events` to parse events
     encoded in such analog stim channels.
+
+    BDF+D (discontinuous) files with actual gaps between data records are not
+    supported and will raise a :class:`NotImplementedError`. BDF+D files that
+    are marked as discontinuous but have no actual gaps will load normally.
     """
     _check_args(input_fname, preload, "bdf")
 
@@ -2355,3 +2419,96 @@ def _get_annotations_gdf(edf_info, sfreq):
         desc = events[2]
 
     return onset, duration, desc
+
+
+def _get_tal_record_times(annotations, encoding="utf8"):
+    """Extract TAL record onset times from EDF+ annotation data.
+
+    In EDF+ files, each data record contains a Time-stamped Annotation List (TAL)
+    that starts with the onset time of that data record. This function extracts
+    these onset times to detect gaps between records in EDF+D (discontinuous) files.
+
+    Parameters
+    ----------
+    annotations : ndarray (n_chans, n_samples) | str
+        Channel data in EDF+ TAL format or path to annotation file.
+    encoding : str
+        Encoding to use when decoding the TAL data.
+
+    Returns
+    -------
+    record_times : list of float
+        List of onset times for each data record, in seconds.
+    """
+    pat = "([+-]\\d+\\.?\\d*)(\x15(\\d+\\.?\\d*))?(\x14.*?)\x14\x00"
+    if isinstance(annotations, str | Path):
+        with open(annotations, "rb") as annot_file:
+            triggers = re.findall(pat.encode(), annot_file.read())
+            triggers = [tuple(map(lambda x: x.decode(encoding), t)) for t in triggers]
+    else:
+        tals = bytearray()
+        annotations = np.atleast_2d(annotations)
+        for chan in annotations:
+            this_chan = chan.ravel()
+            if this_chan.dtype == INT32:  # BDF
+                this_chan = this_chan.view(dtype=UINT8)
+                this_chan = this_chan.reshape(-1, 4)
+                this_chan = this_chan[:, :3].ravel()
+                tals.extend(this_chan)
+            else:
+                this_chan = chan.astype(np.int64)
+                tals.extend(np.uint8([this_chan % 256, this_chan // 256]).flatten("F"))
+        try:
+            triggers = re.findall(pat, tals.decode(encoding))
+        except UnicodeDecodeError:
+            return []
+
+    # Extract record onset times (first TAL entry of each record has empty description)
+    record_times = []
+    for ev in triggers:
+        onset = float(ev[0])
+        # Check if this is a record timestamp (empty description after \x14)
+        descriptions = ev[3].split("\x14")[1:]
+        # The first TAL in each record has the record onset time
+        # If there's no description, it's the record timestamp
+        if not any(descriptions):
+            record_times.append(onset)
+
+    return record_times
+
+
+def _check_edf_discontinuity(record_times, record_length, n_records, tolerance=1e-6):
+    """Check if an EDF+D file has actual gaps between records.
+
+    Parameters
+    ----------
+    record_times : list of float
+        List of onset times for each data record, extracted from TAL annotations.
+    record_length : float
+        Duration of each data record in seconds.
+    n_records : int
+        Expected number of data records.
+    tolerance : float
+        Tolerance for comparing times (in seconds).
+
+    Returns
+    -------
+    has_gaps : bool
+        True if gaps exist between records.
+    gaps : list of tuple
+        List of (onset, duration) tuples for each gap.
+    """
+    if len(record_times) < 2:
+        return False, []
+
+    gaps = []
+    for i in range(len(record_times) - 1):
+        expected_next = record_times[i] + record_length
+        actual_next = record_times[i + 1]
+        gap = actual_next - expected_next
+
+        if gap > tolerance:
+            # Found a gap
+            gaps.append((expected_next, gap))
+
+    return len(gaps) > 0, gaps