19311527649

Committed 12 Nov 2025 08:51PM UTC coverage: 86.671% (+1.9%) from 84.795%

Build # 19311527649

Build Type

Pull #458

github

Committed by

web-flow

Commit Message

Merge pull request #455 from EIT-ALIVE/tests/330-open-source-test-data

Migrate to open source test data and fix revealed issues

Pull Request Pull Request #458: Release 1.8.5

Run Details

781 of 976 branches covered (80.02%)

Branch coverage included in aggregate %.

73 of 82 new or added lines in 4 files covered. (89.02%)

9 existing lines in 2 files now uncovered.

2828 of 3188 relevant lines covered (88.71%)

0.89 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

91.76

/eitprocessing/datahandling/loading/draeger.py

from __future__ import annotations

import math
import mmap
import sys
import warnings
from functools import partial
from typing import TYPE_CHECKING, NamedTuple
from warnings import catch_warnings

import numpy as np
import scipy as sp

from eitprocessing.datahandling.continuousdata import ContinuousData
from eitprocessing.datahandling.datacollection import DataCollection
from eitprocessing.datahandling.eitdata import EITData, Vendor
from eitprocessing.datahandling.event import Event
from eitprocessing.datahandling.intervaldata import IntervalData
from eitprocessing.datahandling.loading import load_eit_data
from eitprocessing.datahandling.loading.binreader import BinReader
from eitprocessing.datahandling.sparsedata import SparseData

if TYPE_CHECKING:
    from pathlib import Path

    from numpy.typing import NDArray

load_draeger_data = partial(load_eit_data, vendor=Vendor.DRAEGER)
NAN_VALUE_INDICATOR = -1e30
SAMPLE_FREQUENCY_ESTIMATION_PRECISION = 4


def load_from_single_path(  # noqa: PLR0915
    path: Path,
    sample_frequency: float | None = None,
    first_frame: int = 0,
    max_frames: int | None = None,
) -> dict[str, DataCollection]:
    """Load Dräger EIT data from path."""
    file_size = path.stat().st_size

    frame_size: int
    medibus_fields: list

    # iterate over the supported file formats to find the frame size that matches the file size
    for _file_format_data in _bin_file_formats.values():
        frame_size = _file_format_data["frame_size"]
        if file_size % frame_size == 0:
            # if the file size is an integer multiple of the frame size, assume this is the correct format
            medibus_fields = _file_format_data["medibus_fields"]
            break
    else:
        msg = (
            f"File size {file_size} of file {path!s} does not match the supported *.bin file formats.\n"
            "Currently this package does not support loading files containing "
            "esophageal pressure or other non-standard data. "
            "Make sure this is a valid and uncorrupted Dräger data file."
        )
        raise OSError(msg)
    total_frames = file_size // frame_size

    if (f0 := first_frame) > (fn := total_frames):
        msg = f"Invalid input: `first_frame` ({f0}) is larger than the total number of frames in the file ({fn})."
        raise ValueError(msg)

    n_frames = min(total_frames - first_frame, max_frames if max_frames is not None else sys.maxsize)
    if n_frames < 1:
        msg = f"No frames to load with `{first_frame=}` and `{max_frames=}`."
        raise ValueError(msg)

    if max_frames and max_frames != n_frames:
        msg = (
            f"The number of frames requested ({max_frames}) is larger "
            f"than the available number ({n_frames}) of frames after "
            f"the first frame selected ({first_frame}, total frames: "
            f"{total_frames}).\n {n_frames} frames will be loaded."
        )
        warnings.warn(msg, RuntimeWarning, stacklevel=2)

    # We need to load 1 frame before first actual frame to check if there is an event marker. Data for the pre-first
    # (dummy) frame will be removed from self at the end of this function.
    load_dummy_frame = first_frame > 0
    first_frame_to_load = first_frame - 1 if load_dummy_frame else 0

    pixel_impedance = np.zeros((n_frames, 32, 32))
    time = np.zeros((n_frames,))
    events: list[tuple[float, Event]] = []
    phases: list[tuple[float, int]] = []
    medibus_data = np.zeros((len(medibus_fields), n_frames), dtype=np.float32)

    with path.open("br") as fo, mmap.mmap(fo.fileno(), length=0, access=mmap.ACCESS_READ) as fh:
        fh.seek(first_frame_to_load * frame_size)
        reader = BinReader(fh)
        previous_marker = None

        first_index = -1 if load_dummy_frame else 0
        for index in range(first_index, n_frames):
            previous_marker = _read_frame(
                reader,
                index,
                time,
                pixel_impedance,
                medibus_data,
                len(medibus_fields),
                events,
                phases,
                previous_marker,
            )

    # time wraps around the number of seconds in a day
    time = np.unwrap(time, period=24 * 60 * 60)

    if not np.all(np.diff(time) > 0):
        msg = "The time axis is not strictly monotonically increasing."
        raise ValueError(msg)

    sample_frequency = _estimate_sample_frequency(time, sample_frequency)

    eit_data = EITData(
        vendor=Vendor.DRAEGER,
        path=path,
        sample_frequency=sample_frequency,
        nframes=n_frames,
        time=time,
        label="raw",
        pixel_impedance=pixel_impedance,
    )
    eitdata_collection = DataCollection(EITData, raw=eit_data)

    (
        continuousdata_collection,
        sparsedata_collection,
    ) = _convert_medibus_data(medibus_data, medibus_fields, time, sample_frequency)
    intervaldata_collection = DataCollection(IntervalData)
    # TODO: move some medibus data to sparse / interval
    # TODO: move phases and events to sparse / interval

    continuousdata_collection.add(
        ContinuousData(
            label="global_impedance_(raw)",
            name="Global impedance (raw)",
            unit="a.u.",
            category="impedance",
            derived_from=[eit_data],
            time=eit_data.time,
            values=eit_data.calculate_global_impedance(),
            sample_frequency=sample_frequency,
        ),
    )
    sparsedata_collection.add(
        SparseData(
            label="minvalues_(draeger)",
            name="Minimum values detected by Draeger device.",
            unit=None,
            category="minvalue",
            derived_from=[eit_data],
            time=np.array([t for t, d in phases if d == -1]),
        ),
    )
    sparsedata_collection.add(
        SparseData(
            label="maxvalues_(draeger)",
            name="Maximum values detected by Draeger device.",
            unit=None,
            category="maxvalue",
            derived_from=[eit_data],
            time=np.array([t for t, d in phases if d == 1]),
        ),
    )
    if events:
        time_, events_ = zip(*events, strict=True)
        time = np.array(time_)
        events = list(events_)
    else:
        time, events = np.array([]), []
    sparsedata_collection.add(
        SparseData(
            label="events_(draeger)",
            name="Events loaded from Draeger data",
            unit=None,
            category="event",
            derived_from=[eit_data],
            time=time,
            values=events,
        ),
    )

    return {
        "eitdata_collection": eitdata_collection,
        "continuousdata_collection": continuousdata_collection,
        "sparsedata_collection": sparsedata_collection,
        "intervaldata_collection": intervaldata_collection,
    }


def _estimate_sample_frequency(time: np.ndarray, sample_frequency: float | None) -> float:
    """Estimate the sample frequency from the time axis, and check with provided sample frequency."""
    with catch_warnings():
        warnings.filterwarnings("ignore", category=RuntimeWarning)
        unrounded_estimated_sample_frequency = 1 / sp.stats.linregress(np.arange(len(time)), time).slope

    if np.isnan(unrounded_estimated_sample_frequency):
        msg = (
            "Could not estimate sample frequency from time axis, "
            f"which could be due to too few data points ({len(time)})."
        )
        if sample_frequency is not None:
            warnings.warn(msg, RuntimeWarning, stacklevel=2)
            return float(sample_frequency)

        raise ValueError(msg)

    # Rounds to the number of digits, rather than the number of decimals
    estimated_sample_frequency = round(
        unrounded_estimated_sample_frequency,
        -math.ceil(np.log10(abs(unrounded_estimated_sample_frequency))) + SAMPLE_FREQUENCY_ESTIMATION_PRECISION,
    )

    if sample_frequency is None:
        return float(estimated_sample_frequency)
    if not isinstance(sample_frequency, (int, float)):
        msg = f"Provided sample frequency has invalid type {type(sample_frequency)}; should be int or float."
        raise TypeError(msg)

    if not np.isclose(
        sample_frequency, unrounded_estimated_sample_frequency, rtol=10**-SAMPLE_FREQUENCY_ESTIMATION_PRECISION, atol=0
    ):
        msg = (
            "Provided sample frequency "
            f"({sample_frequency:.{SAMPLE_FREQUENCY_ESTIMATION_PRECISION + 2}f} Hz) "
            "does not match the estimated sample frequency "
            f"({unrounded_estimated_sample_frequency:.{SAMPLE_FREQUENCY_ESTIMATION_PRECISION + 2}f} Hz) "
            f"within {SAMPLE_FREQUENCY_ESTIMATION_PRECISION} digits. "
            "Note that the estimate might not be as accurate for very short signals."
        )
        warnings.warn(msg, RuntimeWarning, stacklevel=2)

    return float(sample_frequency)


def _convert_medibus_data(
    medibus_data: NDArray,
    medibus_fields: list,
    time: NDArray,
    sample_frequency: float,
) -> tuple[DataCollection, DataCollection]:
    continuousdata_collection = DataCollection(ContinuousData)
    sparsedata_collection = DataCollection(SparseData)

    for field_info, data in zip(medibus_fields, medibus_data, strict=True):
        data[data < NAN_VALUE_INDICATOR] = np.nan
        if field_info.continuous:
            continuous_data = ContinuousData(
                label=field_info.signal_name,
                name=field_info.signal_name,
                description=f"Continuous {field_info.signal_name} data loaded from file",
                unit=field_info.unit,
                time=time,
                values=data,
                category=field_info.signal_name,
                sample_frequency=sample_frequency,
            )
            continuous_data.lock()
            continuousdata_collection.add(continuous_data)

        else:
            # TODO parse sparse data
            ...

    return continuousdata_collection, sparsedata_collection


def _read_frame(
    reader: BinReader,
    index: int,
    time: NDArray,
    pixel_impedance: NDArray,
    medibus_data: NDArray,
    n_medibus_fields: int,
    events: list,
    phases: list,
    previous_marker: int | None,
) -> int:
    """Read frame by frame data from DRAEGER files.

    This method adds the loaded data to the provided arrays `time` and
    `pixel_impedance` and the provided lists `events` and `phases` when the
    index is non-negative. When the index is negative, no data is saved. In
    any case, the event marker is returned.
    """
    frame_time = reader.float64() * 24 * 60 * 60
    _ = reader.float32()
    frame_pixel_impedance = reader.npfloat32(length=1024)
    frame_pixel_impedance = np.reshape(frame_pixel_impedance, (32, 32), "C")
    min_max_flag = reader.int32()
    event_marker = reader.int32()
    event_text = reader.string(length=30)
    timing_error = reader.int32()

    frame_medibus_data = reader.npfloat32(length=n_medibus_fields)

    if index < 0:
        # do not keep any loaded data, just return the event marker
        return event_marker

    time[index] = frame_time
    pixel_impedance[index, :, :] = frame_pixel_impedance
    medibus_data[:, index] = frame_medibus_data

    # The event marker stays the same until the next event occurs.
    # Therefore, check whether the event marker has changed with
    # respect to the most recent event. If so, create a new event.
    if ((previous_marker is not None) and (event_marker > previous_marker)) or (index == 0 and event_text):
        events.append((frame_time, Event(event_marker, event_text)))
    if timing_error:
        warnings.warn("A timing error was encountered during loading.", RuntimeWarning, stacklevel=2)
        # TODO: expand on what timing errors are in some documentation.
    if min_max_flag in (1, -1):
        phases.append((frame_time, min_max_flag))

    return event_marker


class _MedibusField(NamedTuple):
    signal_name: str
    unit: str
    continuous: bool


_bin_file_formats = {
    "original": {
        "frame_size": 4358,
        "medibus_fields": [
            _MedibusField("airway pressure", "mbar", True),
            _MedibusField("flow", "L/min", True),
            _MedibusField("volume", "mL", True),
            _MedibusField("CO2 (%)", "%", True),
            _MedibusField("CO2 (kPa)", "kPa", True),
            _MedibusField("CO2 (mmHg)", "mmHg", True),
            _MedibusField("dynamic compliance", "mL/mbar", False),
            _MedibusField("resistance", "mbar/L/s", False),
            _MedibusField("r^2", "", False),
            _MedibusField("spontaneous inspiratory time", "s", False),
            _MedibusField("minimal pressure", "mbar", False),
            _MedibusField("P0.1", "mbar", False),
            _MedibusField("mean pressure", "mbar", False),
            _MedibusField("plateau pressure", "mbar", False),
            _MedibusField("PEEP", "mbar", False),
            _MedibusField("intrinsic PEEP", "mbar", False),
            _MedibusField("mandatory respiratory rate", "/min", False),
            _MedibusField("mandatory minute volume", "L/min", False),
            _MedibusField("peak inspiratory pressure", "mbar", False),
            _MedibusField("mandatory tidal volume", "L", False),
            _MedibusField("spontaneous tidal volume", "L", False),
            _MedibusField("trapped volume", "mL", False),
            _MedibusField("mandatory expiratory tidal volume", "mL", False),
            _MedibusField("spontaneous expiratory tidal volume", "mL", False),
            _MedibusField("mandatory inspiratory tidal volume", "mL", False),
            _MedibusField("tidal volume", "mL", False),
            _MedibusField("spontaneous inspiratory tidal volume", "mL", False),
            _MedibusField("negative inspiratory force", "mbar", False),
            _MedibusField("leak minute volume", "L/min", False),
            _MedibusField("leak percentage", "%", False),
            _MedibusField("spontaneous respiratory rate", "/min", False),
            _MedibusField("percentage of spontaneous minute volume", "%", False),
            _MedibusField("spontaneous minute volume", "L/min", False),
            _MedibusField("minute volume", "L/min", False),
            _MedibusField("airway temperature", "degrees C", False),
            _MedibusField("rapid shallow breating index", "1/min/L", False),
            _MedibusField("respiratory rate", "/min", False),
            _MedibusField("inspiratory:expiratory ratio", "", False),
            _MedibusField("CO2 flow", "mL/min", False),
            _MedibusField("dead space volume", "mL", False),
            _MedibusField("percentage dead space of expiratory tidal volume", "%", False),
            _MedibusField("end-tidal CO2", "%", False),
            _MedibusField("end-tidal CO2", "kPa", False),
            _MedibusField("end-tidal CO2", "mmHg", False),
            _MedibusField("fraction inspired O2", "%", False),
            _MedibusField("spontaneous inspiratory:expiratory ratio", "", False),
            _MedibusField("elastance", "mbar/L", False),
            _MedibusField("time constant", "s", False),
            _MedibusField(
                "ratio between upper 20% pressure range and total dynamic compliance",
                "",
                False,
            ),
            _MedibusField("end-inspiratory pressure", "mbar", False),
            _MedibusField("expiratory tidal volume", "mL", False),
            _MedibusField("time at low pressure", "s", False),
        ],
    },
    "pressure_pod": {
        "frame_size": 4382,
        "medibus_fields": [
            _MedibusField("airway pressure", "mbar", True),
            _MedibusField("flow", "L/min", True),
            _MedibusField("volume", "mL", True),
            _MedibusField("CO2 (%)", "%", True),
            _MedibusField("CO2 (kPa)", "kPa", True),
            _MedibusField("CO2 (mmHg)", "mmHg", True),
            _MedibusField("dynamic compliance", "mL/mbar", False),
            _MedibusField("resistance", "mbar/L/s", False),
            _MedibusField("r^2", "", False),
            _MedibusField("spontaneous inspiratory time", "s", False),
            _MedibusField("minimal pressure", "mbar", False),
            _MedibusField("P0.1", "mbar", False),
            _MedibusField("mean pressure", "mbar", False),
            _MedibusField("plateau pressure", "mbar", False),
            _MedibusField("PEEP", "mbar", False),
            _MedibusField("intrinsic PEEP", "mbar", False),
            _MedibusField("mandatory respiratory rate", "/min", False),
            _MedibusField("mandatory minute volume", "L/min", False),
            _MedibusField("peak inspiratory pressure", "mbar", False),
            _MedibusField("mandatory tidal volume", "L", False),
            _MedibusField("spontaneous tidal volume", "L", False),
            _MedibusField("trapped volume", "mL", False),
            _MedibusField("mandatory expiratory tidal volume", "mL", False),
            _MedibusField("spontaneous expiratory tidal volume", "mL", False),
            _MedibusField("mandatory inspiratory tidal volume", "mL", False),
            _MedibusField("tidal volume", "mL", False),
            _MedibusField("spontaneous inspiratory tidal volume", "mL", False),
            _MedibusField("negative inspiratory force", "mbar", False),
            _MedibusField("leak minute volume", "L/min", False),
            _MedibusField("leak percentage", "%", False),
            _MedibusField("spontaneous respiratory rate", "/min", False),
            _MedibusField("percentage of spontaneous minute volume", "%", False),
            _MedibusField("spontaneous minute volume", "L/min", False),
            _MedibusField("minute volume", "L/min", False),
            _MedibusField("airway temperature", "degrees C", False),
            _MedibusField("rapid shallow breating index", "1/min/L", False),
            _MedibusField("respiratory rate", "/min", False),
            _MedibusField("inspiratory:expiratory ratio", "", False),
            _MedibusField("CO2 flow", "mL/min", False),
            _MedibusField("dead space volume", "mL", False),
            _MedibusField("percentage dead space of expiratory tidal volume", "%", False),
            _MedibusField("end-tidal CO2", "%", False),
            _MedibusField("end-tidal CO2", "kPa", False),
            _MedibusField("end-tidal CO2", "mmHg", False),
            _MedibusField("fraction inspired O2", "%", False),
            _MedibusField("spontaneous inspiratory:expiratory ratio", "", False),
            _MedibusField("elastance", "mbar/L", False),
            _MedibusField("time constant", "s", False),
            _MedibusField(
                "ratio between upper 20% pressure range and total dynamic compliance",
                "",
                False,
            ),
            _MedibusField("end-inspiratory pressure", "mbar", False),
            _MedibusField("expiratory tidal volume", "mL", False),
            _MedibusField("high pressure", "mbar", False),
            _MedibusField("low pressure", "mbar", False),
            _MedibusField("time at low pressure", "s", False),
            _MedibusField("airway pressure (pod)", "mbar", True),
            _MedibusField("esophageal pressure (pod)", "mbar", True),
            _MedibusField("transpulmonary pressure (pod)", "mbar", True),
            _MedibusField("gastric pressure/auxiliary pressure (pod)", "mbar", True),
        ],
    },
}

1	from __future__ import annotations	1✔
2
3	import math	1✔
4	import mmap	1✔
5	import sys	1✔
6	import warnings	1✔
7	from functools import partial	1✔
8	from typing import TYPE_CHECKING, NamedTuple	1✔
9	from warnings import catch_warnings	1✔
10
11	import numpy as np	1✔
12	import scipy as sp	1✔
13
14	from eitprocessing.datahandling.continuousdata import ContinuousData	1✔
15	from eitprocessing.datahandling.datacollection import DataCollection	1✔
16	from eitprocessing.datahandling.eitdata import EITData, Vendor	1✔
17	from eitprocessing.datahandling.event import Event	1✔
18	from eitprocessing.datahandling.intervaldata import IntervalData	1✔
19	from eitprocessing.datahandling.loading import load_eit_data	1✔
20	from eitprocessing.datahandling.loading.binreader import BinReader	1✔
21	from eitprocessing.datahandling.sparsedata import SparseData	1✔
22
23	if TYPE_CHECKING:
24	from pathlib import Path
25
26	from numpy.typing import NDArray
27
28	load_draeger_data = partial(load_eit_data, vendor=Vendor.DRAEGER)	1✔
29	NAN_VALUE_INDICATOR = -1e30	1✔
30	SAMPLE_FREQUENCY_ESTIMATION_PRECISION = 4	1✔
31
32
33	def load_from_single_path( # noqa: PLR0915	1✔
34	path: Path,
35	sample_frequency: float \| None = None,
36	first_frame: int = 0,
37	max_frames: int \| None = None,
38	) -> dict[str, DataCollection]:
39	"""Load Dräger EIT data from path."""
40	file_size = path.stat().st_size	1✔
41
42	frame_size: int
43	medibus_fields: list
44
45	# iterate over the supported file formats to find the frame size that matches the file size
46	for _file_format_data in _bin_file_formats.values():	1!
47	frame_size = _file_format_data["frame_size"]	1✔
48	if file_size % frame_size == 0:	1✔
49	# if the file size is an integer multiple of the frame size, assume this is the correct format
50	medibus_fields = _file_format_data["medibus_fields"]	1✔
51	break	1✔
52	else:
UNCOV 53	msg = (	×
54	f"File size {file_size} of file {path!s} does not match the supported *.bin file formats.\n"
55	"Currently this package does not support loading files containing "
56	"esophageal pressure or other non-standard data. "
57	"Make sure this is a valid and uncorrupted Dräger data file."
58	)
UNCOV 59	raise OSError(msg)	×
60	total_frames = file_size // frame_size	1✔
61
62	if (f0 := first_frame) > (fn := total_frames):	1!
UNCOV 63	msg = f"Invalid input: `first_frame` ({f0}) is larger than the total number of frames in the file ({fn})."	×
UNCOV 64	raise ValueError(msg)	×
65
66	n_frames = min(total_frames - first_frame, max_frames if max_frames is not None else sys.maxsize)	1✔
67	if n_frames < 1:	1✔
68	msg = f"No frames to load with `{first_frame=}` and `{max_frames=}`."	1✔
69	raise ValueError(msg)	1✔
70
71	if max_frames and max_frames != n_frames:	1✔
72	msg = (	1✔
73	f"The number of frames requested ({max_frames}) is larger "
74	f"than the available number ({n_frames}) of frames after "
75	f"the first frame selected ({first_frame}, total frames: "
76	f"{total_frames}).\n {n_frames} frames will be loaded."
77	)
78	warnings.warn(msg, RuntimeWarning, stacklevel=2)	1✔
79
80	# We need to load 1 frame before first actual frame to check if there is an event marker. Data for the pre-first
81	# (dummy) frame will be removed from self at the end of this function.
82	load_dummy_frame = first_frame > 0	1✔
83	first_frame_to_load = first_frame - 1 if load_dummy_frame else 0	1✔
84
85	pixel_impedance = np.zeros((n_frames, 32, 32))	1✔
86	time = np.zeros((n_frames,))	1✔
87	events: list[tuple[float, Event]] = []	1✔
88	phases: list[tuple[float, int]] = []	1✔
89	medibus_data = np.zeros((len(medibus_fields), n_frames), dtype=np.float32)	1✔
90
91	with path.open("br") as fo, mmap.mmap(fo.fileno(), length=0, access=mmap.ACCESS_READ) as fh:	1✔
92	fh.seek(first_frame_to_load * frame_size)	1✔
93	reader = BinReader(fh)	1✔
94	previous_marker = None	1✔
95
96	first_index = -1 if load_dummy_frame else 0	1✔
97	for index in range(first_index, n_frames):	1✔
98	previous_marker = _read_frame(	1✔
99	reader,
100	index,
101	time,
102	pixel_impedance,
103	medibus_data,
104	len(medibus_fields),
105	events,
106	phases,
107	previous_marker,
108	)
109
110	# time wraps around the number of seconds in a day
111	time = np.unwrap(time, period=24 * 60 * 60)	1✔
112
113	if not np.all(np.diff(time) > 0):	1!
NEW 114	msg = "The time axis is not strictly monotonically increasing."	×
NEW 115	raise ValueError(msg)	×
116
117	sample_frequency = _estimate_sample_frequency(time, sample_frequency)	1✔
118
119	eit_data = EITData(	1✔
120	vendor=Vendor.DRAEGER,
121	path=path,
122	sample_frequency=sample_frequency,
123	nframes=n_frames,
124	time=time,
125	label="raw",
126	pixel_impedance=pixel_impedance,
127	)
128	eitdata_collection = DataCollection(EITData, raw=eit_data)	1✔
129
130	(	1✔
131	continuousdata_collection,
132	sparsedata_collection,
133	) = _convert_medibus_data(medibus_data, medibus_fields, time, sample_frequency)
134	intervaldata_collection = DataCollection(IntervalData)	1✔
135	# TODO: move some medibus data to sparse / interval
136	# TODO: move phases and events to sparse / interval
137
138	continuousdata_collection.add(	1✔
139	ContinuousData(
140	label="global_impedance_(raw)",
141	name="Global impedance (raw)",
142	unit="a.u.",
143	category="impedance",
144	derived_from=[eit_data],
145	time=eit_data.time,
146	values=eit_data.calculate_global_impedance(),
147	sample_frequency=sample_frequency,
148	),
149	)
150	sparsedata_collection.add(	1✔
151	SparseData(
152	label="minvalues_(draeger)",
153	name="Minimum values detected by Draeger device.",
154	unit=None,
155	category="minvalue",
156	derived_from=[eit_data],
157	time=np.array([t for t, d in phases if d == -1]),
158	),
159	)
160	sparsedata_collection.add(	1✔
161	SparseData(
162	label="maxvalues_(draeger)",
163	name="Maximum values detected by Draeger device.",
164	unit=None,
165	category="maxvalue",
166	derived_from=[eit_data],
167	time=np.array([t for t, d in phases if d == 1]),
168	),
169	)
170	if events:	1✔
171	time_, events_ = zip(*events, strict=True)	1✔
172	time = np.array(time_)	1✔
173	events = list(events_)	1✔
174	else:
175	time, events = np.array([]), []	1✔
176	sparsedata_collection.add(	1✔
177	SparseData(
178	label="events_(draeger)",
179	name="Events loaded from Draeger data",
180	unit=None,
181	category="event",
182	derived_from=[eit_data],
183	time=time,
184	values=events,
185	),
186	)
187
188	return {	1✔
189	"eitdata_collection": eitdata_collection,
190	"continuousdata_collection": continuousdata_collection,
191	"sparsedata_collection": sparsedata_collection,
192	"intervaldata_collection": intervaldata_collection,
193	}
194
195
196	def _estimate_sample_frequency(time: np.ndarray, sample_frequency: float \| None) -> float:	1✔
197	"""Estimate the sample frequency from the time axis, and check with provided sample frequency."""
198	with catch_warnings():	1✔
199	warnings.filterwarnings("ignore", category=RuntimeWarning)	1✔
200	unrounded_estimated_sample_frequency = 1 / sp.stats.linregress(np.arange(len(time)), time).slope	1✔
201
202	if np.isnan(unrounded_estimated_sample_frequency):	1✔
203	msg = (	1✔
204	"Could not estimate sample frequency from time axis, "
205	f"which could be due to too few data points ({len(time)})."
206	)
207	if sample_frequency is not None:	1✔
208	warnings.warn(msg, RuntimeWarning, stacklevel=2)	1✔
209	return float(sample_frequency)	1✔
210
211	raise ValueError(msg)	1✔
212
213	# Rounds to the number of digits, rather than the number of decimals
214	estimated_sample_frequency = round(	1✔
215	unrounded_estimated_sample_frequency,
216	-math.ceil(np.log10(abs(unrounded_estimated_sample_frequency))) + SAMPLE_FREQUENCY_ESTIMATION_PRECISION,
217	)
218
219	if sample_frequency is None:	1✔
220	return float(estimated_sample_frequency)	1✔
221	if not isinstance(sample_frequency, (int, float)):	1!
NEW 222	msg = f"Provided sample frequency has invalid type {type(sample_frequency)}; should be int or float."	×
NEW 223	raise TypeError(msg)	×
224
225	if not np.isclose(	1✔
226	sample_frequency, unrounded_estimated_sample_frequency, rtol=10**-SAMPLE_FREQUENCY_ESTIMATION_PRECISION, atol=0
227	):
228	msg = (	1✔
229	"Provided sample frequency "
230	f"({sample_frequency:.{SAMPLE_FREQUENCY_ESTIMATION_PRECISION + 2}f} Hz) "
231	"does not match the estimated sample frequency "
232	f"({unrounded_estimated_sample_frequency:.{SAMPLE_FREQUENCY_ESTIMATION_PRECISION + 2}f} Hz) "
233	f"within {SAMPLE_FREQUENCY_ESTIMATION_PRECISION} digits. "
234	"Note that the estimate might not be as accurate for very short signals."
235	)
236	warnings.warn(msg, RuntimeWarning, stacklevel=2)	1✔
237
238	return float(sample_frequency)	1✔
239
240
241	def _convert_medibus_data(	1✔
242	medibus_data: NDArray,
243	medibus_fields: list,
244	time: NDArray,
245	sample_frequency: float,
246	) -> tuple[DataCollection, DataCollection]:
247	continuousdata_collection = DataCollection(ContinuousData)	1✔
248	sparsedata_collection = DataCollection(SparseData)	1✔
249
250	for field_info, data in zip(medibus_fields, medibus_data, strict=True):	1✔
251	data[data < NAN_VALUE_INDICATOR] = np.nan	1✔
252	if field_info.continuous:	1✔
253	continuous_data = ContinuousData(	1✔
254	label=field_info.signal_name,
255	name=field_info.signal_name,
256	description=f"Continuous {field_info.signal_name} data loaded from file",
257	unit=field_info.unit,
258	time=time,
259	values=data,
260	category=field_info.signal_name,
261	sample_frequency=sample_frequency,
262	)
263	continuous_data.lock()	1✔
264	continuousdata_collection.add(continuous_data)	1✔
265
266	else:
267	# TODO parse sparse data
268	...	1✔
269
270	return continuousdata_collection, sparsedata_collection	1✔
271
272
273	def _read_frame(	1✔
274	reader: BinReader,
275	index: int,
276	time: NDArray,
277	pixel_impedance: NDArray,
278	medibus_data: NDArray,
279	n_medibus_fields: int,
280	events: list,
281	phases: list,
282	previous_marker: int \| None,
283	) -> int:
284	"""Read frame by frame data from DRAEGER files.
285
286	This method adds the loaded data to the provided arrays `time` and
287	`pixel_impedance` and the provided lists `events` and `phases` when the
288	index is non-negative. When the index is negative, no data is saved. In
289	any case, the event marker is returned.
290	"""
291	frame_time = reader.float64() * 24 * 60 * 60	1✔
292	_ = reader.float32()	1✔
293	frame_pixel_impedance = reader.npfloat32(length=1024)	1✔
294	frame_pixel_impedance = np.reshape(frame_pixel_impedance, (32, 32), "C")	1✔
295	min_max_flag = reader.int32()	1✔
296	event_marker = reader.int32()	1✔
297	event_text = reader.string(length=30)	1✔
298	timing_error = reader.int32()	1✔
299
300	frame_medibus_data = reader.npfloat32(length=n_medibus_fields)	1✔
301
302	if index < 0:	1✔
303	# do not keep any loaded data, just return the event marker
304	return event_marker	1✔
305
306	time[index] = frame_time	1✔
307	pixel_impedance[index, :, :] = frame_pixel_impedance	1✔
308	medibus_data[:, index] = frame_medibus_data	1✔
309
310	# The event marker stays the same until the next event occurs.
311	# Therefore, check whether the event marker has changed with
312	# respect to the most recent event. If so, create a new event.
313	if ((previous_marker is not None) and (event_marker > previous_marker)) or (index == 0 and event_text):	1✔
314	events.append((frame_time, Event(event_marker, event_text)))	1✔
315	if timing_error:	1!
316	warnings.warn("A timing error was encountered during loading.", RuntimeWarning, stacklevel=2)	×
317	# TODO: expand on what timing errors are in some documentation.
318	if min_max_flag in (1, -1):	1✔
319	phases.append((frame_time, min_max_flag))	1✔
320
321	return event_marker	1✔
322
323
324	class _MedibusField(NamedTuple):	1✔
325	signal_name: str	1✔
326	unit: str	1✔
327	continuous: bool	1✔
328
329
330	_bin_file_formats = {	1✔
331	"original": {
332	"frame_size": 4358,
333	"medibus_fields": [
334	_MedibusField("airway pressure", "mbar", True),
335	_MedibusField("flow", "L/min", True),
336	_MedibusField("volume", "mL", True),
337	_MedibusField("CO2 (%)", "%", True),
338	_MedibusField("CO2 (kPa)", "kPa", True),
339	_MedibusField("CO2 (mmHg)", "mmHg", True),
340	_MedibusField("dynamic compliance", "mL/mbar", False),
341	_MedibusField("resistance", "mbar/L/s", False),
342	_MedibusField("r^2", "", False),
343	_MedibusField("spontaneous inspiratory time", "s", False),
344	_MedibusField("minimal pressure", "mbar", False),
345	_MedibusField("P0.1", "mbar", False),
346	_MedibusField("mean pressure", "mbar", False),
347	_MedibusField("plateau pressure", "mbar", False),
348	_MedibusField("PEEP", "mbar", False),
349	_MedibusField("intrinsic PEEP", "mbar", False),
350	_MedibusField("mandatory respiratory rate", "/min", False),
351	_MedibusField("mandatory minute volume", "L/min", False),
352	_MedibusField("peak inspiratory pressure", "mbar", False),
353	_MedibusField("mandatory tidal volume", "L", False),
354	_MedibusField("spontaneous tidal volume", "L", False),
355	_MedibusField("trapped volume", "mL", False),
356	_MedibusField("mandatory expiratory tidal volume", "mL", False),
357	_MedibusField("spontaneous expiratory tidal volume", "mL", False),
358	_MedibusField("mandatory inspiratory tidal volume", "mL", False),
359	_MedibusField("tidal volume", "mL", False),
360	_MedibusField("spontaneous inspiratory tidal volume", "mL", False),
361	_MedibusField("negative inspiratory force", "mbar", False),
362	_MedibusField("leak minute volume", "L/min", False),
363	_MedibusField("leak percentage", "%", False),
364	_MedibusField("spontaneous respiratory rate", "/min", False),
365	_MedibusField("percentage of spontaneous minute volume", "%", False),
366	_MedibusField("spontaneous minute volume", "L/min", False),
367	_MedibusField("minute volume", "L/min", False),
368	_MedibusField("airway temperature", "degrees C", False),
369	_MedibusField("rapid shallow breating index", "1/min/L", False),
370	_MedibusField("respiratory rate", "/min", False),
371	_MedibusField("inspiratory:expiratory ratio", "", False),
372	_MedibusField("CO2 flow", "mL/min", False),
373	_MedibusField("dead space volume", "mL", False),
374	_MedibusField("percentage dead space of expiratory tidal volume", "%", False),
375	_MedibusField("end-tidal CO2", "%", False),
376	_MedibusField("end-tidal CO2", "kPa", False),
377	_MedibusField("end-tidal CO2", "mmHg", False),
378	_MedibusField("fraction inspired O2", "%", False),
379	_MedibusField("spontaneous inspiratory:expiratory ratio", "", False),
380	_MedibusField("elastance", "mbar/L", False),
381	_MedibusField("time constant", "s", False),
382	_MedibusField(
383	"ratio between upper 20% pressure range and total dynamic compliance",
384	"",
385	False,
386	),
387	_MedibusField("end-inspiratory pressure", "mbar", False),
388	_MedibusField("expiratory tidal volume", "mL", False),
389	_MedibusField("time at low pressure", "s", False),
390	],
391	},
392	"pressure_pod": {
393	"frame_size": 4382,
394	"medibus_fields": [
395	_MedibusField("airway pressure", "mbar", True),
396	_MedibusField("flow", "L/min", True),
397	_MedibusField("volume", "mL", True),
398	_MedibusField("CO2 (%)", "%", True),
399	_MedibusField("CO2 (kPa)", "kPa", True),
400	_MedibusField("CO2 (mmHg)", "mmHg", True),
401	_MedibusField("dynamic compliance", "mL/mbar", False),
402	_MedibusField("resistance", "mbar/L/s", False),
403	_MedibusField("r^2", "", False),
404	_MedibusField("spontaneous inspiratory time", "s", False),
405	_MedibusField("minimal pressure", "mbar", False),
406	_MedibusField("P0.1", "mbar", False),
407	_MedibusField("mean pressure", "mbar", False),
408	_MedibusField("plateau pressure", "mbar", False),
409	_MedibusField("PEEP", "mbar", False),
410	_MedibusField("intrinsic PEEP", "mbar", False),
411	_MedibusField("mandatory respiratory rate", "/min", False),
412	_MedibusField("mandatory minute volume", "L/min", False),
413	_MedibusField("peak inspiratory pressure", "mbar", False),
414	_MedibusField("mandatory tidal volume", "L", False),
415	_MedibusField("spontaneous tidal volume", "L", False),
416	_MedibusField("trapped volume", "mL", False),
417	_MedibusField("mandatory expiratory tidal volume", "mL", False),
418	_MedibusField("spontaneous expiratory tidal volume", "mL", False),
419	_MedibusField("mandatory inspiratory tidal volume", "mL", False),
420	_MedibusField("tidal volume", "mL", False),
421	_MedibusField("spontaneous inspiratory tidal volume", "mL", False),
422	_MedibusField("negative inspiratory force", "mbar", False),
423	_MedibusField("leak minute volume", "L/min", False),
424	_MedibusField("leak percentage", "%", False),
425	_MedibusField("spontaneous respiratory rate", "/min", False),
426	_MedibusField("percentage of spontaneous minute volume", "%", False),
427	_MedibusField("spontaneous minute volume", "L/min", False),
428	_MedibusField("minute volume", "L/min", False),
429	_MedibusField("airway temperature", "degrees C", False),
430	_MedibusField("rapid shallow breating index", "1/min/L", False),
431	_MedibusField("respiratory rate", "/min", False),
432	_MedibusField("inspiratory:expiratory ratio", "", False),
433	_MedibusField("CO2 flow", "mL/min", False),
434	_MedibusField("dead space volume", "mL", False),
435	_MedibusField("percentage dead space of expiratory tidal volume", "%", False),
436	_MedibusField("end-tidal CO2", "%", False),
437	_MedibusField("end-tidal CO2", "kPa", False),
438	_MedibusField("end-tidal CO2", "mmHg", False),
439	_MedibusField("fraction inspired O2", "%", False),
440	_MedibusField("spontaneous inspiratory:expiratory ratio", "", False),
441	_MedibusField("elastance", "mbar/L", False),
442	_MedibusField("time constant", "s", False),
443	_MedibusField(
444	"ratio between upper 20% pressure range and total dynamic compliance",
445	"",
446	False,
447	),
448	_MedibusField("end-inspiratory pressure", "mbar", False),
449	_MedibusField("expiratory tidal volume", "mL", False),
450	_MedibusField("high pressure", "mbar", False),
451	_MedibusField("low pressure", "mbar", False),
452	_MedibusField("time at low pressure", "s", False),
453	_MedibusField("airway pressure (pod)", "mbar", True),
454	_MedibusField("esophageal pressure (pod)", "mbar", True),
455	_MedibusField("transpulmonary pressure (pod)", "mbar", True),
456	_MedibusField("gastric pressure/auxiliary pressure (pod)", "mbar", True),
457	],
458	},
459	}

EIT-ALIVE / eitprocessing / 19311527649

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous