diff --git a/cpp_package/examples/signal_processing/src/downsampling.cpp b/cpp_package/examples/signal_processing/src/downsampling.cpp index 5473344cf..3abbd6497 100644 --- a/cpp_package/examples/signal_processing/src/downsampling.cpp +++ b/cpp_package/examples/signal_processing/src/downsampling.cpp @@ -45,6 +45,14 @@ int main (int argc, char *argv[]) int filtered_size = 0; std::vector eeg_channels = BoardShim::get_eeg_channels (board_id); + // Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. + // The operation changes data in-place and may be applied before other signal processing. + if (eeg_channels.size () >= 2) + { + std::vector reference_channels = {eeg_channels[0], eeg_channels[1]}; + DataFilter::reference (data, eeg_channels, reference_channels); + } + for (int i = 0; i < eeg_channels.size (); i++) { std::cout << "Data from :" << eeg_channels[i] << " before downsampling " << std::endl; diff --git a/cpp_package/src/data_filter.cpp b/cpp_package/src/data_filter.cpp index 561fcc124..84e179eec 100644 --- a/cpp_package/src/data_filter.cpp +++ b/cpp_package/src/data_filter.cpp @@ -7,6 +7,65 @@ #include "data_handler.h" +void DataFilter::reference (BrainFlowArray &data, + const std::vector &channels_to_reference, const std::vector &reference_channels) +{ + const int rows = data.get_size (0); + const int cols = data.get_size (1); + if ((rows <= 0) || (cols <= 0) || channels_to_reference.empty () || + reference_channels.empty ()) + { + throw BrainFlowException ( + "data and channel lists must be non-empty", + (int)BrainFlowExitCodes::INVALID_ARGUMENTS_ERROR); + } + + for (const int channel : channels_to_reference) + { + if ((channel < 0) || (channel >= rows)) + { + throw BrainFlowException ( + "channel index is out of range", + (int)BrainFlowExitCodes::INVALID_ARGUMENTS_ERROR); + } + } + for (const int channel : reference_channels) + { + if ((channel < 0) || (channel >= rows)) + { + throw BrainFlowException ( + "reference channel index is out of range", + (int)BrainFlowExitCodes::INVALID_ARGUMENTS_ERROR); + } + } + + // Compute the complete reference signal before mutating any channel. This keeps the result + // correct when a reference channel is also present in channels_to_reference. + std::vector reference_signal (cols, 0.0); + for (const int channel : reference_channels) + { + const double *channel_data = data.get_address (channel); + for (int sample = 0; sample < cols; sample++) + { + reference_signal[sample] += channel_data[sample]; + } + } + const double scale = 1.0 / static_cast (reference_channels.size ()); + for (double &sample : reference_signal) + { + sample *= scale; + } + + for (const int channel : channels_to_reference) + { + double *channel_data = data.get_address (channel); + for (int sample = 0; sample < cols; sample++) + { + channel_data[sample] -= reference_signal[sample]; + } + } +} + double DataFilter::get_oxygen_level (double *ppg_ir, double *ppg_red, int data_len, int sampling_rate, double coef1, double coef2, double coef3) { diff --git a/cpp_package/src/inc/data_filter.h b/cpp_package/src/inc/data_filter.h index de66122de..89d7cb7c0 100644 --- a/cpp_package/src/inc/data_filter.h +++ b/cpp_package/src/inc/data_filter.h @@ -28,6 +28,16 @@ class DataFilter /// write user defined string to BrainFlow logger static void log_message (int log_level, const char *format, ...); + /** + * re-reference selected channels in-place + * @param data input 2d array, rows are channels and columns are samples + * @param channels_to_reference channel rows from which to subtract the reference signal + * @param reference_channels channel rows whose sample-wise mean defines the reference signal + */ + static void reference (BrainFlowArray &data, + const std::vector &channels_to_reference, + const std::vector &reference_channels); + /// perform low pass filter in-place static void perform_lowpass (double *data, int data_len, int sampling_rate, double cutoff, int order, int filter_type, double ripple); diff --git a/csharp_package/brainflow/brainflow/data_filter.cs b/csharp_package/brainflow/brainflow/data_filter.cs index b8f14f429..0bcffc16d 100644 --- a/csharp_package/brainflow/brainflow/data_filter.cs +++ b/csharp_package/brainflow/brainflow/data_filter.cs @@ -76,6 +76,61 @@ public static void set_log_file (string log_file) throw new BrainFlowError (res); } } + + /// + /// Re-reference selected channels in-place using the sample-wise mean of reference channels. + /// + /// Rows are channels and columns are samples. + /// Rows from which to subtract the reference signal. + /// Rows whose mean defines the reference signal. + public static void reference (double[,] data, int[] channels_to_reference, int[] reference_channels) + { + if ((data == null) || (channels_to_reference == null) || (reference_channels == null) || + (data.GetLength (0) == 0) || (data.GetLength (1) == 0) || + (channels_to_reference.Length == 0) || (reference_channels.Length == 0)) + { + throw new BrainFlowError ((int)BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR); + } + + int rows = data.GetLength (0); + int cols = data.GetLength (1); + foreach (int channel in channels_to_reference) + { + if ((channel < 0) || (channel >= rows)) + { + throw new BrainFlowError ((int)BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR); + } + } + foreach (int channel in reference_channels) + { + if ((channel < 0) || (channel >= rows)) + { + throw new BrainFlowError ((int)BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR); + } + } + + // Snapshot the reference before changing any row so overlapping channel lists are safe. + double[] reference_signal = new double[cols]; + foreach (int channel in reference_channels) + { + for (int sample = 0; sample < cols; sample++) + { + reference_signal[sample] += data[channel, sample]; + } + } + for (int sample = 0; sample < cols; sample++) + { + reference_signal[sample] /= reference_channels.Length; + } + + foreach (int channel in channels_to_reference) + { + for (int sample = 0; sample < cols; sample++) + { + data[channel, sample] -= reference_signal[sample]; + } + } + } // accord GetRow returns a copy instead pointer, so we can not easily update data in place like in other bindings /// diff --git a/csharp_package/brainflow/examples/downsampling/downsampling.cs b/csharp_package/brainflow/examples/downsampling/downsampling.cs index e903a9525..813dd34af 100644 --- a/csharp_package/brainflow/examples/downsampling/downsampling.cs +++ b/csharp_package/brainflow/examples/downsampling/downsampling.cs @@ -22,6 +22,15 @@ static void Main (string[] args) board_shim.stop_stream (); double[,] unprocessed_data = board_shim.get_board_data (); int[] eeg_channels = BoardShim.get_eeg_channels (board_id); + + // Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. + // The operation changes unprocessed_data in-place. + if (eeg_channels.Length >= 2) + { + int[] reference_channels = {eeg_channels[0], eeg_channels[1]}; + DataFilter.reference (unprocessed_data, eeg_channels, reference_channels); + } + board_shim.release_session (); for (int i = 0; i < eeg_channels.Length; i++) diff --git a/docs/SwiftAPIParity.rst b/docs/SwiftAPIParity.rst index b07f2922e..aec28f469 100644 --- a/docs/SwiftAPIParity.rst +++ b/docs/SwiftAPIParity.rst @@ -29,8 +29,8 @@ Implemented: - Filters, noise removal, and detrending: lowpass, highpass, bandpass, bandstop, environmental noise removal, rolling filter, and detrend. -- Transforms and features: downsampling, wavelet transform/inverse/denoising, CSP, windowing, - FFT/IFFT, PSD/Welch, band powers, and ICA. +- Transforms and features: channel re-referencing, downsampling, wavelet + transform/inverse/denoising, CSP, windowing, FFT/IFFT, PSD/Welch, band powers, and ICA. - Helpers: standard deviation, railed percentage, oxygen level, heart rate, peak detection, nearest power of two, file I/O, reshape helpers, logging, and version. diff --git a/java_package/brainflow/src/main/java/brainflow/DataFilter.java b/java_package/brainflow/src/main/java/brainflow/DataFilter.java index cff1fc7ba..9a9708fa2 100644 --- a/java_package/brainflow/src/main/java/brainflow/DataFilter.java +++ b/java_package/brainflow/src/main/java/brainflow/DataFilter.java @@ -192,6 +192,74 @@ public static void set_log_file (String log_file) throws BrainFlowError } } + /** + * Re-reference selected channels in-place using the sample-wise mean of reference channels. + * + * @param data rows are channels and columns are samples + * @param channels_to_reference rows from which to subtract the reference signal + * @param reference_channels rows whose mean defines the reference signal + */ + public static void reference (double[][] data, int[] channels_to_reference, int[] reference_channels) + throws BrainFlowError + { + if ((data == null) || (channels_to_reference == null) || (reference_channels == null) + || (data.length == 0) || (data[0] == null) || (data[0].length == 0) + || (channels_to_reference.length == 0) || (reference_channels.length == 0)) + { + throw new BrainFlowError ("Invalid reference arguments", + BrainFlowExitCode.INVALID_ARGUMENTS_ERROR.get_code ()); + } + + final int rows = data.length; + final int cols = data[0].length; + for (double[] row : data) + { + if ((row == null) || (row.length != cols)) + { + throw new BrainFlowError ("Data must be a rectangular matrix", + BrainFlowExitCode.INVALID_ARGUMENTS_ERROR.get_code ()); + } + } + for (int channel : channels_to_reference) + { + if ((channel < 0) || (channel >= rows)) + { + throw new BrainFlowError ("Channel index is out of range", + BrainFlowExitCode.INVALID_ARGUMENTS_ERROR.get_code ()); + } + } + for (int channel : reference_channels) + { + if ((channel < 0) || (channel >= rows)) + { + throw new BrainFlowError ("Reference channel index is out of range", + BrainFlowExitCode.INVALID_ARGUMENTS_ERROR.get_code ()); + } + } + + // Snapshot the reference before changing any row so overlapping channel lists are safe. + double[] reference_signal = new double[cols]; + for (int channel : reference_channels) + { + for (int sample = 0; sample < cols; sample++) + { + reference_signal[sample] += data[channel][sample]; + } + } + for (int sample = 0; sample < cols; sample++) + { + reference_signal[sample] /= reference_channels.length; + } + + for (int channel : channels_to_reference) + { + for (int sample = 0; sample < cols; sample++) + { + data[channel][sample] -= reference_signal[sample]; + } + } + } + /** * calc stddev */ diff --git a/java_package/brainflow/src/main/java/brainflow/examples/Downsampling.java b/java_package/brainflow/src/main/java/brainflow/examples/Downsampling.java index 1e574cdde..643263335 100644 --- a/java_package/brainflow/src/main/java/brainflow/examples/Downsampling.java +++ b/java_package/brainflow/src/main/java/brainflow/examples/Downsampling.java @@ -30,6 +30,14 @@ public static void main (String[] args) throws Exception board_shim.release_session (); int[] eeg_channels = BoardShim.get_eeg_channels (board_id); + // Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. + // The operation changes data in-place. + if (eeg_channels.length >= 2) + { + int[] reference_channels = {eeg_channels[0], eeg_channels[1]}; + DataFilter.reference (data, eeg_channels, reference_channels); + } + for (int i = 0; i < eeg_channels.length; i++) { System.out.println ("Original data:"); diff --git a/julia_package/brainflow/src/data_filter.jl b/julia_package/brainflow/src/data_filter.jl index 3c5e6ad42..3b4a19253 100644 --- a/julia_package/brainflow/src/data_filter.jl +++ b/julia_package/brainflow/src/data_filter.jl @@ -141,6 +141,43 @@ end WaveletType = Union{WaveletTypes, Integer} +""" + reference(data, channels_to_reference, reference_channels) + +Re-reference selected channel rows in-place using the sample-wise mean of +`reference_channels`. Channel indices follow Julia's one-based convention. +""" +function reference(data::AbstractMatrix{<:AbstractFloat}, channels_to_reference, + reference_channels) + rows, cols = size(data) + if rows == 0 || cols == 0 || isempty(channels_to_reference) || isempty(reference_channels) + throw(BrainFlowError("Data and channel lists must be non-empty", Integer(INVALID_ARGUMENTS_ERROR))) + end + if !all(channel -> channel isa Integer && 1 <= channel <= rows, channels_to_reference) + throw(BrainFlowError("Channel index is out of range", Integer(INVALID_ARGUMENTS_ERROR))) + end + if !all(channel -> channel isa Integer && 1 <= channel <= rows, reference_channels) + throw(BrainFlowError("Reference channel index is out of range", Integer(INVALID_ARGUMENTS_ERROR))) + end + + # Snapshot the reference before changing any row so overlapping channel lists are safe. + reference_signal = zeros(eltype(data), cols) + for channel in reference_channels + for sample in 1:cols + reference_signal[sample] += data[channel, sample] + end + end + reference_signal ./= length(reference_channels) + + for channel in channels_to_reference + for sample in 1:cols + data[channel, sample] -= reference_signal[sample] + end + end + return nothing +end + + @brainflow_rethrow function perform_lowpass(data, sampling_rate::Integer, cutoff::Float64, order::Integer, filter_type::FiltType, ripple::Float64) ccall((:perform_lowpass, DATA_HANDLER_INTERFACE), Cint, (Ptr{Float64}, Cint, Cint, Float64, Cint, Cint, Float64), diff --git a/julia_package/brainflow/test/downsampling.jl b/julia_package/brainflow/test/downsampling.jl index 0f7e508b5..a148a9847 100644 --- a/julia_package/brainflow/test/downsampling.jl +++ b/julia_package/brainflow/test/downsampling.jl @@ -16,9 +16,15 @@ BrainFlow.release_session(board_shim) eeg_channels = BrainFlow.get_eeg_channels(BrainFlow.SYNTHETIC_BOARD) sampling_rate = BrainFlow.get_sampling_rate(BrainFlow.SYNTHETIC_BOARD) +# Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. +# The operation changes data in-place. +if length(eeg_channels) >= 2 + BrainFlow.reference(data, eeg_channels, eeg_channels[1:2]) +end + data_first_channel = data[eeg_channels[1], :] println("Original Data First Channel") println(data_first_channel) downsampled_data = BrainFlow.perform_downsampling(data_first_channel, 3, BrainFlow.EACH) println("After Downsampling") -println(downsampled_data) \ No newline at end of file +println(downsampled_data) diff --git a/julia_package/brainflow/test/julia_tests.jl b/julia_package/brainflow/test/julia_tests.jl index 75526647f..399644ea2 100644 --- a/julia_package/brainflow/test/julia_tests.jl +++ b/julia_package/brainflow/test/julia_tests.jl @@ -43,6 +43,17 @@ end @test Int32(BrainFlow.BLACKMAN_HARRIS) == 3 end +@testset "channel reference" begin + data = [1.0 2.0; 3.0 4.0; 5.0 6.0] + BrainFlow.reference(data, [1, 3], [1, 2]) + @test data == [-1.0 -1.0; 3.0 4.0; 3.0 3.0] + + @test_throws BrainFlow.BrainFlowError BrainFlow.reference(copy(data), Int[], [1]) + @test_throws BrainFlow.BrainFlowError BrainFlow.reference(copy(data), [1], Int[]) + @test_throws BrainFlow.BrainFlowError BrainFlow.reference(copy(data), [4], [1]) + @test_throws BrainFlow.BrainFlowError BrainFlow.reference(copy(data), [1], [4]) +end + @testset "presets" begin presets = BrainFlow.get_board_presets(BrainFlow.CYTON_BOARD) @test presets[1] == Int32(BrainFlow.DEFAULT_PRESET) @@ -52,4 +63,4 @@ end params = BrainFlowModelParams(BrainFlow.RESTFULNESS, BrainFlow.DEFAULT_CLASSIFIER) @test params.metric == BrainFlow.RESTFULNESS @test params.classifier == BrainFlow.DEFAULT_CLASSIFIER -end \ No newline at end of file +end diff --git a/matlab_package/brainflow/DataFilter.m b/matlab_package/brainflow/DataFilter.m index e5764ed1b..8a64e31f0 100644 --- a/matlab_package/brainflow/DataFilter.m +++ b/matlab_package/brainflow/DataFilter.m @@ -78,6 +78,33 @@ function disable_data_logger() DataFilter.set_log_level(int32(6)) end + function referenced_data = reference(data, channels_to_reference, reference_channels) + % Re-reference selected channel rows using the sample-wise mean of reference channels. + % Channel indices follow MATLAB's one-based convention. + if ~isnumeric(data) || ~ismatrix(data) || isempty(data) || ... + isempty(channels_to_reference) || isempty(reference_channels) + error('BrainFlow:InvalidArguments', 'Data and channel lists must be non-empty'); + end + rows = size(data, 1); + target_indices = channels_to_reference(:); + reference_indices = reference_channels(:); + valid_targets = all(isfinite(target_indices)) && ... + all(target_indices == fix(target_indices)) && ... + all(target_indices >= 1) && all(target_indices <= rows); + valid_references = all(isfinite(reference_indices)) && ... + all(reference_indices == fix(reference_indices)) && ... + all(reference_indices >= 1) && all(reference_indices <= rows); + if ~valid_targets || ~valid_references + error('BrainFlow:InvalidArguments', 'Channel index is out of range'); + end + + % Snapshot the reference before changing any row so overlapping lists are safe. + reference_signal = mean(data(reference_indices, :), 1); + referenced_data = data; + referenced_data(target_indices, :) = ... + referenced_data(target_indices, :) - reference_signal; + end + function filtered_data = perform_lowpass(data, sampling_rate, cutoff, order, filter_type, ripple) % perform lowpass filtering task_name = 'perform_lowpass'; @@ -462,4 +489,4 @@ function write_file(data, file_name, file_mode) end -end \ No newline at end of file +end diff --git a/matlab_package/brainflow/examples/Downsampling.m b/matlab_package/brainflow/examples/Downsampling.m index 2d7e9200b..e814aa491 100644 --- a/matlab_package/brainflow/examples/Downsampling.m +++ b/matlab_package/brainflow/examples/Downsampling.m @@ -12,7 +12,13 @@ board_shim.release_session(); eeg_channels = BoardShim.get_eeg_channels(int32(BoardIds.SYNTHETIC_BOARD), preset); +% Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. +% DataFilter.reference returns a referenced copy of the matrix. +if length(eeg_channels) >= 2 + data = DataFilter.reference(data, eeg_channels, eeg_channels(1:2)); +end + % apply downsampling to the first eeg channel % first_eeg_channel = eeg_channels(1); original_data = data(first_eeg_channel, :); -downsampled_data = DataFilter.perform_downsampling(original_data, 3, int32(AggOperations.MEAN)); \ No newline at end of file +downsampled_data = DataFilter.perform_downsampling(original_data, 3, int32(AggOperations.MEAN)); diff --git a/nodejs_package/brainflow/data_filter.ts b/nodejs_package/brainflow/data_filter.ts index 66ea6a94d..bb9d911a5 100644 --- a/nodejs_package/brainflow/data_filter.ts +++ b/nodejs_package/brainflow/data_filter.ts @@ -160,6 +160,65 @@ export class DataFilter return out[0].substring(0, len[0]); } + /** + * Re-reference selected channels in-place using the sample-wise mean of reference channels. + * Rows are channels and columns are samples. + */ + public static reference(data: number[][], channelsToReference: number[], + referenceChannels: number[]): void + { + if (!Array.isArray(data) || data.length === 0 || !Array.isArray(data[0]) || + data[0].length === 0 || !Array.isArray(channelsToReference) || + !Array.isArray(referenceChannels) || channelsToReference.length === 0 || + referenceChannels.length === 0) + { + throw new BrainFlowError ( + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR, 'Invalid reference arguments'); + } + + const rows = data.length; + const cols = data[0].length; + if (!data.every((row) => Array.isArray(row) && row.length === cols)) + { + throw new BrainFlowError ( + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR, 'Data must be a rectangular matrix'); + } + if (!channelsToReference.every( + (channel) => Number.isInteger(channel) && channel >= 0 && channel < rows)) + { + throw new BrainFlowError ( + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR, 'Channel index is out of range'); + } + if (!referenceChannels.every( + (channel) => Number.isInteger(channel) && channel >= 0 && channel < rows)) + { + throw new BrainFlowError (BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR, + 'Reference channel index is out of range'); + } + + // Snapshot the reference before changing any row so overlapping channel lists are safe. + const referenceSignal = new Array(cols).fill(0); + for (const channel of referenceChannels) + { + for (let sample = 0; sample < cols; sample++) + { + referenceSignal[sample] += data[channel][sample]; + } + } + for (let sample = 0; sample < cols; sample++) + { + referenceSignal[sample] /= referenceChannels.length; + } + + for (const channel of channelsToReference) + { + for (let sample = 0; sample < cols; sample++) + { + data[channel][sample] -= referenceSignal[sample]; + } + } + } + // signal processing methods public static getRailedPercentage(data: number[], gain: number): number { diff --git a/nodejs_package/tests/downsampling.ts b/nodejs_package/tests/downsampling.ts index 5a941af70..0233c203e 100644 --- a/nodejs_package/tests/downsampling.ts +++ b/nodejs_package/tests/downsampling.ts @@ -16,6 +16,14 @@ async function runExample (): Promise const data = board.getCurrentBoardData(10); board.releaseSession(); const eegChannels = BoardShim.getEegChannels(boardId); + + // Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. + // The operation changes data in-place. + if (eegChannels.length >= 2) + { + DataFilter.reference(data, eegChannels, eegChannels.slice(0, 2)); + } + const oldData = data[eegChannels[0]]; console.info(oldData); const newData = DataFilter.performDownsampling(oldData, 3, AggOperations.MEAN); diff --git a/python_package/brainflow/data_filter.py b/python_package/brainflow/data_filter.py index d98a3dda2..cb68fcf6a 100644 --- a/python_package/brainflow/data_filter.py +++ b/python_package/brainflow/data_filter.py @@ -590,6 +590,41 @@ def set_log_file(cls, log_file: str) -> None: if res != BrainFlowExitCodes.STATUS_OK.value: raise BrainFlowError('unable to redirect logs to a file', res) + @classmethod + def reference(cls, data, channels_to_reference: List[int], reference_channels: List[int]) -> None: + """Re-reference selected channels in-place. + + The sample-wise mean of ``reference_channels`` is computed from the original data and + subtracted from every row in ``channels_to_reference``. Computing the complete reference + signal first makes overlapping target and reference channel lists safe. + + :param data: 2-D C-contiguous float64 array; rows are channels and columns are samples + :param channels_to_reference: rows from which to subtract the reference signal + :param reference_channels: rows whose mean defines the reference signal + """ + check_memory_layout_row_major(data, 2) + if data.shape[0] == 0 or data.shape[1] == 0 or data.dtype != numpy.float64: + raise BrainFlowError('data must be a non-empty float64 matrix', + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR.value) + if (channels_to_reference is None or reference_channels is None or + len(channels_to_reference) == 0 or len(reference_channels) == 0): + raise BrainFlowError('channel lists must be non-empty', + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR.value) + + rows = data.shape[0] + for channel in channels_to_reference: + if not isinstance(channel, (int, numpy.integer)) or channel < 0 or channel >= rows: + raise BrainFlowError('channel index is out of range', + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR.value) + for channel in reference_channels: + if not isinstance(channel, (int, numpy.integer)) or channel < 0 or channel >= rows: + raise BrainFlowError('reference channel index is out of range', + BrainFlowExitCodes.INVALID_ARGUMENTS_ERROR.value) + + reference_signal = numpy.mean(data[reference_channels, :], axis=0) + for channel in channels_to_reference: + data[channel, :] -= reference_signal + @classmethod def perform_lowpass(cls, data, sampling_rate: int, cutoff: float, order: int, filter_type: int, ripple: float) -> None: diff --git a/python_package/examples/tests/downsampling.py b/python_package/examples/tests/downsampling.py index 16c861e84..7c7699edc 100644 --- a/python_package/examples/tests/downsampling.py +++ b/python_package/examples/tests/downsampling.py @@ -19,6 +19,11 @@ def main(): board.release_session() eeg_channels = BoardShim.get_eeg_channels(BoardIds.SYNTHETIC_BOARD.value) + # Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. + # The operation changes data in-place. + if len(eeg_channels) >= 2: + DataFilter.reference(data, eeg_channels, eeg_channels[:2]) + # demo for downsampling, it just aggregates data for count, channel in enumerate(eeg_channels): print('Original data for channel %d:' % channel) diff --git a/r_package/examples/downsampling.R b/r_package/examples/downsampling.R index 89c886d7a..4c5ad9117 100644 --- a/r_package/examples/downsampling.R +++ b/r_package/examples/downsampling.R @@ -9,8 +9,18 @@ board_shim$stop_stream() data <- board_shim$get_current_board_data(as.integer(250)) board_shim$release_session() -# need to convert to numpy array manually -numpy_data <- np$array(data[2,]) +# Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. +# R delegates DataFilter operations to the Python binding, so keep the matrix as a NumPy array. +numpy_matrix <- np$array(data) +eeg_channels <- brainflow_python$BoardShim$get_eeg_channels( + brainflow_python$BoardIds$SYNTHETIC_BOARD$value) +if (length(eeg_channels) >= 2) +{ + brainflow_python$DataFilter$reference( + numpy_matrix, eeg_channels, eeg_channels[1:2]) +} + +numpy_data <- np$array(numpy_matrix[2,]) print(numpy_data) brainflow_python$DataFilter$perform_downsampling(numpy_data, as.integer(3), brainflow_python$AggOperations$EACH$value) -print(numpy_data) \ No newline at end of file +print(numpy_data) diff --git a/rust_package/brainflow/src/data_filter.rs b/rust_package/brainflow/src/data_filter.rs index 8bdb06ee6..0ebb61524 100644 --- a/rust_package/brainflow/src/data_filter.rs +++ b/rust_package/brainflow/src/data_filter.rs @@ -57,6 +57,44 @@ pub fn set_log_file>(log_file: S) -> Result<()> { Ok(check_brainflow_exit_code(res)?) } +/// Re-reference selected channels in-place using the sample-wise mean of reference channels. +/// Rows are channels and columns are samples. +pub fn reference( + data: &mut Array2, + channels_to_reference: &[usize], + reference_channels: &[usize], +) -> Result<()> { + let (rows, cols) = data.dim(); + if rows == 0 + || cols == 0 + || channels_to_reference.is_empty() + || reference_channels.is_empty() + || channels_to_reference.iter().any(|&channel| channel >= rows) + || reference_channels.iter().any(|&channel| channel >= rows) + { + return Err(Error::BrainFlowError(BrainFlowError::InvalidArgumentsError)); + } + + // Snapshot the reference before changing any row so overlapping channel lists are safe. + let mut reference_signal = vec![0.0; cols]; + for &channel in reference_channels { + for sample in 0..cols { + reference_signal[sample] += data[[channel, sample]]; + } + } + let scale = 1.0 / reference_channels.len() as f64; + for value in &mut reference_signal { + *value *= scale; + } + + for &channel in channels_to_reference { + for sample in 0..cols { + data[[channel, sample]] -= reference_signal[sample]; + } + } + Ok(()) +} + /// Apply low pass filter to provided data. pub fn perform_lowpass( data: &mut [f64], @@ -926,4 +964,23 @@ mod tests { let read_data = read_file(filename).unwrap(); assert_eq!(data, read_data); } + + #[test] + fn reference_uses_original_reference_signal_for_overlapping_channels() { + let mut data = array![[1.0, 2.0], [3.0, 4.0], [5.0, 6.0]]; + + reference(&mut data, &[0, 2], &[0, 1]).unwrap(); + + assert_eq!(data, array![[-1.0, -1.0], [3.0, 4.0], [3.0, 3.0]]); + } + + #[test] + fn reference_rejects_empty_or_out_of_range_channel_lists() { + let mut data = array![[1.0, 2.0], [3.0, 4.0]]; + + assert!(reference(&mut data, &[], &[0]).is_err()); + assert!(reference(&mut data, &[0], &[]).is_err()); + assert!(reference(&mut data, &[2], &[0]).is_err()); + assert!(reference(&mut data, &[0], &[2]).is_err()); + } } diff --git a/swift_package/Sources/BrainFlow/DataFilter.swift b/swift_package/Sources/BrainFlow/DataFilter.swift index ec21ef5a9..79fcc1138 100644 --- a/swift_package/Sources/BrainFlow/DataFilter.swift +++ b/swift_package/Sources/BrainFlow/DataFilter.swift @@ -44,6 +44,43 @@ public enum DataFilter { try getVersion(function: \.get_version_data_handler) } + /// Re-reference selected channels in-place using the sample-wise mean of reference channels. + /// Rows are channels and columns are samples. + public static func reference( + data: inout [[Double]], + channels_to_reference: [Int], + reference_channels: [Int] + ) throws { + let (rows, cols) = try BrainFlowArray.validateRectangular(data) + guard !channels_to_reference.isEmpty, !reference_channels.isEmpty else { + throw invalidArguments("Channel lists must be non-empty") + } + guard channels_to_reference.allSatisfy({ $0 >= 0 && $0 < rows }) else { + throw invalidArguments("Channel index is out of range") + } + guard reference_channels.allSatisfy({ $0 >= 0 && $0 < rows }) else { + throw invalidArguments("Reference channel index is out of range") + } + + // Snapshot the reference before changing any row so overlapping channel lists are safe. + var referenceSignal = [Double](repeating: 0.0, count: cols) + for channel in reference_channels { + for sample in 0..= 2 { + // Re-reference every EEG channel using the sample-wise mean of the first two EEG channels. + // The operation changes referencedData in-place. + try DataFilter.reference( + data: &referencedData, + channels_to_reference: sample.eegChannels, + reference_channels: Array(sample.eegChannels.prefix(2)) + ) + } + + guard let firstEEGChannel = sample.eegChannels.first else { + throw BrainFlowError("No EEG channel found", BrainFlowExitCodes.GENERAL_ERROR.rawValue) + } + let downsampled = try DataFilter.perform_downsampling( + data: referencedData[firstEEGChannel], + period: 4, + operation: AggOperations.MEAN + ) print(downsampled) } }