search

mne-tools / mne-python

MNE: Magnetoencephalography (MEG) and Electroencephalography (EEG) in Python
https://mne.tools
BSD 3-Clause "New" or "Revised" License
2.66k stars 1.31k forks source link

_extract_sampling_rate function for io.read_raw_snirf incorrectly estimates sample rate due to floating point error. #12508

Open kdarti opened 5 months ago

kdarti commented 5 months ago

Description of the problem

I have a snirf file with a sample rate of 50hz, the /nirs/data1/time array is correctly incremented by 0.02 per sample, as expected.

However, when mne.io_read_raw_snirf() is used to read the file, the sample rate is estimated to be not 50Hz:

image

This is due to floating errors introduced in this section:

https://github.com/mne-tools/mne-python/blob/415e7f68ed71135baff0ea857ca4fab5a3690bf8/mne/io/snirf/_snirf.py#L542-L558

More precisely, by .np.diff() on line 541.

Easiest fix would be to use Decimal to avoid these floating point errors. using float64 should also work, at least in this case.

Steps to reproduce

import mne

mne.io.read_raw_snirf(r"WR11 Baseline SNIRF.snirf")

Link to data

WR11 Baseline SNIRF.zip

Expected results

Sample rate of 50Hz

Actual results

49.999996Hz

Additional information

Platform Windows-10-10.0.19041-SP0 Python 3.8.3 (default, Jul 2 2020, 17:30:36) [MSC v.1916 64 bit (AMD64)] Executable C:\Users\kdahlslatt\Anaconda3\python.exe CPU Intel64 Family 6 Model 158 Stepping 10, GenuineIntel (12 cores) Memory 15.7 GB

Core ├☑ mne 1.6.1 (latest release) ├☑ numpy 1.24.3 (MKL 2023.1-Product with 6 threads) ├☑ scipy 1.10.1 ├☑ matplotlib 3.7.1 (backend=Qt5Agg) C:\Users\kdahlslatt\Anaconda3\lib\site-packages\paramiko\transport.py:219: CryptographyDeprecationWarning: Blowfish has been deprecated "class": algorithms.Blowfish, ├☑ pooch 1.8.0 └☑ jinja2 3.1.2

Numerical (optional) ├☑ sklearn 1.3.0 ├☑ numba 0.57.1 ├☑ nibabel 3.2.1 ├☑ nilearn 0.9.1 ├☑ pandas 1.5.3 └☐ unavailable dipy, openmeeg, cupy

Visualization (optional) ├☑ pyvista 0.37.0 (OpenGL 4.5.0 - Build 30.0.101.1404 via Intel(R) UHD Graphics 630) ├☑ vtk 9.2.5 ├☑ qtpy 2.2.0 (PyQt5=5.15.2) ├☑ ipympl 0.9.3 ├☑ pyqtgraph 0.11.0 ├☑ ipywidgets 8.0.4 └☐ unavailable pyvistaqt, mne-qt-browser, trame_client, trame_server, trame_vtk, trame_vuetify

Ecosystem (optional) └☐ unavailable mne-bids, mne-nirs, mne-features, mne-connectivity, mne-icalabel, mne-bids-pipeline

larsoner commented 5 months ago

@kdarti agreed ideally we'd get 50 here, are you up for trying some fix?

drammock commented 5 months ago

IIRC, @rob-luke and I wrote this code; if you're up for making a PR yourself to fix this, feel free to ping me for reviews!

kdarti commented 5 months ago

I tried simply using np.float64, but there are still floating point errors. The option of using Decimal package isn't great because it is very very slow (in a quick test, 750x slower with default precision). The remaining option would be to upscale to int and then later downscale, but I don't really have time to do this unfortunately.