Commit b6f9a130 authored by Julian Kosciessa's avatar Julian Kosciessa
Browse files

implement Chebychev filter definition for bandpass version

parent e241ced1
...@@ -200,12 +200,21 @@ for s = 1:numel(timescales) % loop through timescales ...@@ -200,12 +200,21 @@ for s = 1:numel(timescales) % loop through timescales
case 'bp' case 'bp'
fs = data.fsample; fs = data.fsample;
nyquist = fs/2; nyquist = fs/2;
fcLowPass = (1./sc).*nyquist + ((1./sc).*nyquist)./4; fcLowPass = (1./sc).*nyquist + .05*((1./sc).*nyquist);
fcHighPass = (1./sc).*nyquist - ((1./sc).*nyquist)./4; fcHighPass = (1./(sc+1)).*nyquist - .05*((1./(sc+1)).*nyquist);
if sc > 1 % don't define low-pass for first scale, as its upper frequency limit is specified anyway if sc > 1 % don't define low-pass for first scale, as its upper frequency limit is specified anyway
[B,A] = butter(6,fcLowPass/nyquist, 'low'); % define low-pass filter: https://de.mathworks.com/help/signal/ref/butter.html if fcLowPass-fcHighPass > .05*nyquist
end [B,A] = cheby1(4,1,fcLowPass/nyquist, 'low'); % define low-pass filter: https://de.mathworks.com/help/signal/ref/butter.html
[D,C] = butter(6,fcHighPass/nyquist, 'high'); % define high-pass filter [D,C] = cheby1(4,1,fcHighPass/nyquist,'high'); % define high-pass filter
else
[B,A]=butter(10,fcLowPass/nyquist, 'low'); % Lowpass
[D,C]=butter(10,fcHighPass/nyquist,'high'); % Highpass
end
else
[D,C]= butter(10,fcHighPass/nyquist,'high'); % use Butterworth highpass
fcLowPass = nyquist;
end
cfg.freq(1,s) = fcLowPass; cfg.freq(1,s) = fcLowPass;
cfg.freq(2,s) = fcHighPass; cfg.freq(2,s) = fcHighPass;
...@@ -314,7 +323,7 @@ for s = 1:numel(timescales) % loop through timescales ...@@ -314,7 +323,7 @@ for s = 1:numel(timescales) % loop through timescales
cg_data = {}; cg_data = {};
switch coarsegrainmethod switch coarsegrainmethod
case 'filtskip' case 'filtskip'
if strcmp(filtmethod, 'hp') || strcmp(filtmethod, 'bp') if strcmp(filtmethod, 'hp')
nloops = 1; % keep original sampling rate nloops = 1; % keep original sampling rate
stepSize = 1; stepSize = 1;
else else
......
...@@ -165,12 +165,21 @@ for s = 1:numel(timescales) % loop through timescales ...@@ -165,12 +165,21 @@ for s = 1:numel(timescales) % loop through timescales
case 'bp' case 'bp'
fs = data.fsample; fs = data.fsample;
nyquist = fs/2; nyquist = fs/2;
fcLowPass = (1./sc).*nyquist + ((1./sc).*nyquist)./4; fcLowPass = (1./sc).*nyquist + .05*((1./sc).*nyquist);
fcHighPass = (1./sc).*nyquist - ((1./sc).*nyquist)./4; fcHighPass = (1./(sc+1)).*nyquist - .05*((1./(sc+1)).*nyquist);
if sc > 1 % don't define low-pass for first scale, as its upper frequency limit is specified anyway if sc > 1 % don't define low-pass for first scale, as its upper frequency limit is specified anyway
[B,A] = butter(6,fcLowPass/nyquist, 'low'); % define low-pass filter: https://de.mathworks.com/help/signal/ref/butter.html if fcLowPass-fcHighPass > .05*nyquist
end [B,A] = cheby1(4,1,fcLowPass/nyquist, 'low'); % define low-pass filter: https://de.mathworks.com/help/signal/ref/butter.html
[D,C] = butter(6,fcHighPass/nyquist, 'high'); % define high-pass filter [D,C] = cheby1(4,1,fcHighPass/nyquist,'high'); % define high-pass filter
else
[B,A]=butter(10,fcLowPass/nyquist, 'low'); % Lowpass
[D,C]=butter(10,fcHighPass/nyquist,'high'); % Highpass
end
else
[D,C]= butter(10,fcHighPass/nyquist,'high'); % use Butterworth highpass
fcLowPass = nyquist;
end
cfg.freq(1,s) = fcLowPass; cfg.freq(1,s) = fcLowPass;
cfg.freq(2,s) = fcHighPass; cfg.freq(2,s) = fcHighPass;
...@@ -239,7 +248,7 @@ for s = 1:numel(timescales) % loop through timescales ...@@ -239,7 +248,7 @@ for s = 1:numel(timescales) % loop through timescales
cg_data = {}; cg_data = {};
switch coarsegrainmethod switch coarsegrainmethod
case 'filtskip' case 'filtskip'
if strcmp(filtmethod, 'hp') || strcmp(filtmethod, 'bp') if strcmp(filtmethod, 'hp')
nloops = 1; % keep original sampling rate for hp option nloops = 1; % keep original sampling rate for hp option
stepSize = 1; stepSize = 1;
else else
...@@ -253,7 +262,7 @@ for s = 1:numel(timescales) % loop through timescales ...@@ -253,7 +262,7 @@ for s = 1:numel(timescales) % loop through timescales
end end
clear resamp_x; clear resamp_x;
case 'pointavg' % original point averaging coarse graining, no loop over starting points case 'pointavg' % original point averaging coarse graining, no loop over starting points
if sc == 1 || strcmp(filtmethod, 'hp') % no coarse graining for native sampling rate or high-pass entropy if sc == 1 % no coarse graining for native sampling rate
cg_data{1} = data_sel.trial; %only keep trial data cg_data{1} = data_sel.trial; %only keep trial data
nloops = 1; % no loop across starting points nloops = 1; % no loop across starting points
else % coarse-grain time series at this time scale else % coarse-grain time series at this time scale
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment