[BUG] Cross-spectral density missing frequencies (#12633)

mne-tools · May 29, 2024 · 9a222ba · 9a222ba
1 parent 5fe73cb
commit 9a222ba
Show file tree

Hide file tree

Showing 3 changed files with 10 additions and 6 deletions.
diff --git a/doc/changes/devel/12633.bugfix.rst b/doc/changes/devel/12633.bugfix.rst
@@ -0,0 +1 @@
+Fix bug where :func:`mne.time_frequency.csd_multitaper`, :func:`mne.time_frequency.csd_fourier`, :func:`mne.time_frequency.csd_array_multitaper`, and :func:`mne.time_frequency.csd_array_fourier` would return cross-spectral densities with the ``fmin`` and ``fmax`` frequencies missing, by `Thomas Binns`_
diff --git a/mne/time_frequency/csd.py b/mne/time_frequency/csd.py
@@ -810,9 +810,10 @@ def csd_array_fourier(
     n_fft = n_times if n_fft is None else n_fft
 
     # Preparing frequencies of interest
-    # orig_frequencies = fftfreq(n_fft, 1. / sfreq)
     orig_frequencies = rfftfreq(n_fft, 1.0 / sfreq)
-    freq_mask = (orig_frequencies > fmin) & (orig_frequencies < fmax)
+    freq_mask = (
+        (orig_frequencies > 0) & (orig_frequencies >= fmin) & (orig_frequencies <= fmax)
+    )
     frequencies = orig_frequencies[freq_mask]
 
     if len(frequencies) == 0:
@@ -1013,7 +1014,9 @@ def csd_array_multitaper(
 
     # Preparing frequencies of interest
     orig_frequencies = rfftfreq(n_fft, 1.0 / sfreq)
-    freq_mask = (orig_frequencies > fmin) & (orig_frequencies < fmax)
+    freq_mask = (
+        (orig_frequencies > 0) & (orig_frequencies >= fmin) & (orig_frequencies <= fmax)
+    )
     frequencies = orig_frequencies[freq_mask]
 
     if len(frequencies) == 0:

diff --git a/mne/time_frequency/tests/test_csd.py b/mne/time_frequency/tests/test_csd.py
@@ -394,15 +394,15 @@ def _test_fourier_multitaper_parameters(epochs, csd_epochs, csd_array):
         fmin=20,
         fmax=10,
     )
-    raises(ValueError, csd_epochs, epochs, fmin=20, fmax=20.1)
+    raises(ValueError, csd_epochs, epochs, fmin=20.11, fmax=20.19)
     raises(
         ValueError,
         csd_array,
         epochs._data,
         epochs.info["sfreq"],
         epochs.tmin,
-        fmin=20,
-        fmax=20.1,
+        fmin=20.11,
+        fmax=20.19,
     )
     raises(ValueError, csd_epochs, epochs, tmin=0.15, tmax=0.1)
     raises(