Skip to content

Commit

Permalink
improve the IRASA implementation
Browse files Browse the repository at this point in the history
This closes the IRASA issue reported in #1546

Features of the new implementation:
1) This version corrected the computation order of the old one.
2) It combines IRASA alorigthm and ft_specest_mtfft and returns fractal power and original power as users defined in cfg.output.
3) It remained comparable frequency resolution in the fractal and orignal spetrums.

Co-authored-by: Robert Oostenveld <r.oostenveld@gmail.com>
  • Loading branch information
ruiliu-psy and robertoostenveld committed Jan 14, 2021
1 parent 7f71464 commit cb490a3
Show file tree
Hide file tree
Showing 4 changed files with 304 additions and 234 deletions.
19 changes: 8 additions & 11 deletions ft_freqanalysis.m
Expand Up @@ -288,14 +288,11 @@

case 'irasa'
cfg.taper = ft_getopt(cfg, 'taper', 'hanning');
cfg.output = ft_getopt(cfg, 'output', 'pow');
cfg.output = ft_getopt(cfg, 'output', 'fractal');
cfg.pad = ft_getopt(cfg, 'pad', 'nextpow2');
if ~isequal(cfg.taper, 'hanning')
ft_error('the irasa method supports hanning tapers only');
end
if ~isequal(cfg.output, 'pow')
ft_error('the irasa method outputs power only');
end
if ~isequal(cfg.pad, 'nextpow2')
ft_warning('consider using cfg.pad=''nextpow2'' for the irasa method');
end
Expand Down Expand Up @@ -345,7 +342,7 @@
cfg.pad = 'maxperlen';
end
cfg.padtype = ft_getopt(cfg, 'padtype', 'zero');
cfg.output = ft_getopt(cfg, 'output', 'pow');
cfg.output = ft_getopt(cfg, 'output', 'pow'); % the default for irasa is set earlier
cfg.calcdof = ft_getopt(cfg, 'calcdof', 'no');
cfg.channel = ft_getopt(cfg, 'channel', 'all');
cfg.precision = ft_getopt(cfg, 'precision', 'double');
Expand Down Expand Up @@ -384,7 +381,7 @@
end

% Set flags for output
if strcmp(cfg.output, 'pow')
if ismember(cfg.output, {'pow','fractal','original'})
powflg = 1;
csdflg = 0;
fftflg = 0;
Expand Down Expand Up @@ -482,9 +479,9 @@

% options that don't change over trials
if isfield(cfg, 'tapsmofrq')
options = {'pad', cfg.pad, 'padtype', cfg.padtype, 'freqoi', cfg.foi, 'tapsmofrq', cfg.tapsmofrq, 'polyorder', cfg.polyremoval};
options = {'pad', cfg.pad, 'padtype', cfg.padtype, 'freqoi', cfg.foi, 'tapsmofrq', cfg.tapsmofrq, 'polyorder', cfg.polyremoval, 'output', cfg.output};
else
options = {'pad', cfg.pad, 'padtype', cfg.padtype, 'freqoi', cfg.foi, 'polyorder', cfg.polyremoval};
options = {'pad', cfg.pad, 'padtype', cfg.padtype, 'freqoi', cfg.foi, 'polyorder', cfg.polyremoval, 'output', cfg.output};
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Expand Down Expand Up @@ -531,11 +528,11 @@
case 'mtmfft'
[spectrum,ntaper,foi] = ft_specest_mtmfft(dat, time, 'taper', cfg.taper, options{:}, 'feedback', fbopt);
hastime = false;

case 'irasa'
[spectrum,ntaper,foi] = ft_specest_irasa(dat, time, 'taper', cfg.taper, options{:}, 'feedback', fbopt);
[spectrum,ntaper,foi] = ft_specest_irasa(dat, time, options{:}, 'feedback', fbopt);
hastime = false;

case 'wavelet'
[spectrum,foi,toi] = ft_specest_wavelet(dat, time, 'timeoi', cfg.toi, 'width', cfg.width, 'gwidth', cfg.gwidth,options{:}, 'feedback', fbopt);

Expand Down

0 comments on commit cb490a3

Please sign in to comment.