Skip to content

Commit

Permalink
Fixed the parfor issue with non-consecutive channels
Browse files Browse the repository at this point in the history
  • Loading branch information
VisLab committed Jul 31, 2021
1 parent 3ed337e commit 45e553b
Show file tree
Hide file tree
Showing 2 changed files with 130 additions and 128 deletions.
133 changes: 67 additions & 66 deletions PrepPipeline/utilities/blasstLineNoise.m
Original file line number Diff line number Diff line change
@@ -1,66 +1,67 @@
function [signal, lineNoiseOut] = blasstLineNoise(signal, lineNoiseIn)
% Remove sharp spectral peaks from signal using Sleppian filters
%
% Usage:
% signal = cleanLineNoise(signal)
% [signal, lineNoiseOut] = hcleanLineNoise(signal, lineNoiseIn)
%
% Parameters:
% signal Structure with .data and .srate fields
% lineNoiseIn Input structure with fields described below
%
% Structure parameters (lineNoiseIn):
% fPassBand Frequency band used (default [0, Fs/2] = entire band)
% Fs Sampling frequency
% fScanBandWidth +/- bandwidth centered on each f0 to scan for significant
% lines (TM)
% lineFrequencies Line frequencies to be removed (default
% [60, 120, 180, 240, 300])
% lineNoiseChannels Channels to remove line noise from (default
% size(data, 1))
% maximumIterations Maximum times to iterate removal (default = 10)
% p Significance level cutoff (default = 0.01)
% pad FFT padding factor ( -1 corresponds to no padding,
% 0 corresponds to padding to next highest power of 2
% etc.) (default is 0)
% pnts
% tapers Precomputed tapers from dpss
% taperBandWidth Taper bandwidth (default 2 Hz)
% taperWindowSize Taper sliding window length (default 4 sec)
% taperWindowStep Sliding window step size (default 4 sec = no overlap)
% tau Window overlap smoothing factor (default 100)
%
% This function is based on code originally written by Tim Mullen in a
% package called tmullen-cleanline which is based on the chronux_2
% libraries.
%

lineNoiseOut = lineNoiseIn;
%% Remove line frequencies that are greater than Nyquist frequencies
tooLarge = lineNoiseOut.lineFrequencies >= lineNoiseOut.Fs/2;
if any(tooLarge)
warning('cleanLineNoise:LineFrequenciesTooLarge', ...
'Eliminating frequencies greater than half the sampling rate');
lineNoiseOut.lineFrequencies(tooLarge) = [];
lineNoiseOut.lineFrequencies = squeeze(lineNoiseOut.lineFrequencies);
end

%% Set up parameters for blassting the line noise
fRange = lineNoiseOut.fScanBandWidth;
frequencyRanges = repmat(fRange, length(lineNoiseOut.lineFrequencies));
sRate = lineNoiseOut.Fs;
lineFrequencies = lineNoiseOut.lineFrequencies;
maxIterations = lineNoiseOut.maximumIterations;

%% Perform the calculation for each channel separately
data = double(signal.data);
chans = sort(lineNoiseOut.lineNoiseChannels);
parfor ch = chans
data(ch, :) = blasst(squeeze(data(ch, :)), lineFrequencies, ...
frequencyRanges, sRate, ...
'MaximumIterations', maxIterations, ...
'Verbose', 0);
end
signal.data = data;
clear data;

function [signal, lineNoiseOut] = blasstLineNoise(signal, lineNoiseIn)
% Remove sharp spectral peaks from signal using Sleppian filters
%
% Usage:
% signal = cleanLineNoise(signal)
% [signal, lineNoiseOut] = hcleanLineNoise(signal, lineNoiseIn)
%
% Parameters:
% signal Structure with .data and .srate fields
% lineNoiseIn Input structure with fields described below
%
% Structure parameters (lineNoiseIn):
% fPassBand Frequency band used (default [0, Fs/2] = entire band)
% Fs Sampling frequency
% fScanBandWidth +/- bandwidth centered on each f0 to scan for significant
% lines (TM)
% lineFrequencies Line frequencies to be removed (default
% [60, 120, 180, 240, 300])
% lineNoiseChannels Channels to remove line noise from (default
% size(data, 1))
% maximumIterations Maximum times to iterate removal (default = 10)
% p Significance level cutoff (default = 0.01)
% pad FFT padding factor ( -1 corresponds to no padding,
% 0 corresponds to padding to next highest power of 2
% etc.) (default is 0)
% pnts
% tapers Precomputed tapers from dpss
% taperBandWidth Taper bandwidth (default 2 Hz)
% taperWindowSize Taper sliding window length (default 4 sec)
% taperWindowStep Sliding window step size (default 4 sec = no overlap)
% tau Window overlap smoothing factor (default 100)
%
% This function is based on code originally written by Tim Mullen in a
% package called tmullen-cleanline which is based on the chronux_2
% libraries.
%

lineNoiseOut = lineNoiseIn;
%% Remove line frequencies that are greater than Nyquist frequencies
tooLarge = lineNoiseOut.lineFrequencies >= lineNoiseOut.Fs/2;
if any(tooLarge)
warning('cleanLineNoise:LineFrequenciesTooLarge', ...
'Eliminating frequencies greater than half the sampling rate');
lineNoiseOut.lineFrequencies(tooLarge) = [];
lineNoiseOut.lineFrequencies = squeeze(lineNoiseOut.lineFrequencies);
end

%% Set up parameters for blassting the line noise
fRange = lineNoiseOut.fScanBandWidth;
frequencyRanges = repmat(fRange, length(lineNoiseOut.lineFrequencies));
sRate = lineNoiseOut.Fs;
lineFrequencies = lineNoiseOut.lineFrequencies;
maxIterations = lineNoiseOut.maximumIterations;

%% Perform the calculation for each channel separately
signal.data = double(signal.data);
chans = lineNoiseOut.lineNoiseChannels;
data = signal.data(chans, :);
parfor ch = 1:size(data, 1)
data(ch, :) = blasst(squeeze(data(ch, :)), lineFrequencies, ...
frequencyRanges, sRate, ...
'MaximumIterations', maxIterations, ...
'Verbose', 0);
end
signal.data(chans, :) = data;
clear data;

125 changes: 63 additions & 62 deletions PrepPipeline/utilities/cleanLineNoise.m
Original file line number Diff line number Diff line change
@@ -1,62 +1,63 @@
function [signal, lineNoiseOut] = cleanLineNoise(signal, lineNoiseIn)
% Remove sharp spectral peaks from signal using Sleppian filters
%
% Usage:
% signal = cleanLineNoise(signal)
% [signal, lineNoiseOut] = hcleanLineNoise(signal, lineNoiseIn)
%
% Parameters:
% signal Structure with .data and .srate fields
% lineNoiseIn Input structure with fields described below
%
% Structure parameters (lineNoiseIn):
% fPassBand Frequency band used (default [0, Fs/2] = entire band)
% Fs Sampling frequency
% fScanBandWidth +/- bandwidth centered on each f0 to scan for significant
% lines (TM)
% lineFrequencies Line frequencies to be removed (default
% [60, 120, 180, 240, 300])
% lineNoiseChannels Channels to remove line noise from (default
% size(data, 1))
% maximumIterations Maximum times to iterate removal (default = 10)
% p Significance level cutoff (default = 0.01)
% pad FFT padding factor ( -1 corresponds to no padding,
% 0 corresponds to padding to next highest power of 2
% etc.) (default is 0)
% pnts
% tapers Precomputed tapers from dpss
% taperBandWidth Taper bandwidth (default 2 Hz)
% taperWindowSize Taper sliding window length (default 4 sec)
% taperWindowStep Sliding window step size (default 4 sec = no overlap)
% tau Window overlap smoothing factor (default 100)
%
% This function is based on code originally written by Tim Mullen in a
% package called tmullen-cleanline which is based on the chronux_2
% libraries.
%

lineNoiseOut = lineNoiseIn;
%% Remove line frequencies that are greater than Nyquist frequencies
tooLarge = lineNoiseOut.lineFrequencies >= lineNoiseOut.Fs/2;
if any(tooLarge)
warning('cleanLineNoise:LineFrequenciesTooLarge', ...
'Eliminating frequencies greater than half the sampling rate');
lineNoiseOut.lineFrequencies(tooLarge) = [];
lineNoiseOut.lineFrequencies = squeeze(lineNoiseOut.lineFrequencies);
end

%% Set up multi-taper parameters
hbw = lineNoiseOut.taperBandWidth/2; % half-bandwidth
lineNoiseOut.taperTemplate = [hbw, lineNoiseOut.taperWindowSize, 1];
Nwin = round(lineNoiseOut.Fs*lineNoiseOut.taperWindowSize); % number of samples in window
lineNoiseOut.tapers = checkTapers(lineNoiseOut.taperTemplate, Nwin, lineNoiseOut.Fs);

%% Perform the calculation for each channel separately
data = double(signal.data);
chans = sort(lineNoiseOut.lineNoiseChannels);
parfor ch = chans
data(ch, :) = removeLinesMovingWindow(squeeze(data(ch, :)), lineNoiseOut);
end
signal.data = data;
clear data;

function [signal, lineNoiseOut] = cleanLineNoise(signal, lineNoiseIn)
% Remove sharp spectral peaks from signal using Sleppian filters
%
% Usage:
% signal = cleanLineNoise(signal)
% [signal, lineNoiseOut] = hcleanLineNoise(signal, lineNoiseIn)
%
% Parameters:
% signal Structure with .data and .srate fields
% lineNoiseIn Input structure with fields described below
%
% Structure parameters (lineNoiseIn):
% fPassBand Frequency band used (default [0, Fs/2] = entire band)
% Fs Sampling frequency
% fScanBandWidth +/- bandwidth centered on each f0 to scan for significant
% lines (TM)
% lineFrequencies Line frequencies to be removed (default
% [60, 120, 180, 240, 300])
% lineNoiseChannels Channels to remove line noise from (default
% size(data, 1))
% maximumIterations Maximum times to iterate removal (default = 10)
% p Significance level cutoff (default = 0.01)
% pad FFT padding factor ( -1 corresponds to no padding,
% 0 corresponds to padding to next highest power of 2
% etc.) (default is 0)
% pnts
% tapers Precomputed tapers from dpss
% taperBandWidth Taper bandwidth (default 2 Hz)
% taperWindowSize Taper sliding window length (default 4 sec)
% taperWindowStep Sliding window step size (default 4 sec = no overlap)
% tau Window overlap smoothing factor (default 100)
%
% This function is based on code originally written by Tim Mullen in a
% package called tmullen-cleanline which is based on the chronux_2
% libraries.
%

lineNoiseOut = lineNoiseIn;
%% Remove line frequencies that are greater than Nyquist frequencies
tooLarge = lineNoiseOut.lineFrequencies >= lineNoiseOut.Fs/2;
if any(tooLarge)
warning('cleanLineNoise:LineFrequenciesTooLarge', ...
'Eliminating frequencies greater than half the sampling rate');
lineNoiseOut.lineFrequencies(tooLarge) = [];
lineNoiseOut.lineFrequencies = squeeze(lineNoiseOut.lineFrequencies);
end

%% Set up multi-taper parameters
hbw = lineNoiseOut.taperBandWidth/2; % half-bandwidth
lineNoiseOut.taperTemplate = [hbw, lineNoiseOut.taperWindowSize, 1];
Nwin = round(lineNoiseOut.Fs*lineNoiseOut.taperWindowSize); % number of samples in window
lineNoiseOut.tapers = checkTapers(lineNoiseOut.taperTemplate, Nwin, lineNoiseOut.Fs);

%% Perform the calculation for each channel separately
signal.data = double(signal.data);
chans = lineNoiseOut.lineNoiseChannels;
data = signal.data(chans, :);
parfor ch = 1:size(data, 1)
data(ch, :) = removeLinesMovingWindow(squeeze(data(ch, :)), lineNoiseOut);
end
signal.data(chans, :) = data;
clear data;

0 comments on commit 45e553b

Please sign in to comment.