/
ft_removetemplateartifact.m
150 lines (127 loc) · 5.59 KB
/
ft_removetemplateartifact.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
function [data] = ft_removetemplateartifact(cfg, data, template)
% FT_REMOVETEMPLATEARTIFACT removes an artifact from preprocessed data by template
% subtraction. The template can for example be formed by averaging an ECG-triggered
% MEG timecourse.
%
% Use as
% dataclean = ft_removetemplateartifact(cfg, data, template)
% where data is raw data as obtained from FT_PREPROCESSING and template is a averaged
% timelock structure as obtained from FT_TIMELOCKANALYSIS. The configuration should
% be according to
%
% cfg.channel = Nx1 cell-array with selection of channels (default = 'all'), see FT_CHANNELSELECTION for details
% cfg.artifact = Mx2 matrix with sample numbers of the artifact segments, e.g. obtained from FT_ARTIFACT_EOG
%
% To facilitate data-handling and distributed computing you can use
% cfg.inputfile = ...
% cfg.outputfile = ...
% If you specify one of these (or both) the input data will be read from a *.mat
% file on disk and/or the output data will be written to a *.mat file. These mat
% files should contain only a single variable, corresponding with the
% input/output structure.
%
% See also FT_ARTIFACT_ECG, FT_PREPROCESSING, FT_TIMELOCKANALYSIS, FT_REJECTCOMPONENT
% Copyright (C) 2014, Robert Oostenveld
%
% This file is part of FieldTrip, see http://www.fieldtriptoolbox.org
% for the documentation and details.
%
% FieldTrip is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% FieldTrip is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with FieldTrip. If not, see <http://www.gnu.org/licenses/>.
%
% $Id$
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% the initial part deals with parsing the input options and data
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% these are used by the ft_preamble/ft_postamble function and scripts
ft_revision = '$Id$';
ft_nargin = nargin;
ft_nargout = nargout;
% do the general setup of the function
ft_defaults
ft_preamble init
ft_preamble debug
ft_preamble loadvar data template
ft_preamble provenance data template
% the ft_abort variable is set to true or false in ft_preamble_init
if ft_abort
% do not continue function execution in case the outputfile is present and the user indicated to keep it
return
end
% ensure that the input data is valid for this function, this will also do
% backward-compatibility conversions of old data that for example was
% read from an old *.mat file
data = ft_checkdata(data, 'datatype', 'raw', 'feedback', 'yes', 'hassampleinfo', 'yes');
template = ft_checkdata(template, 'datatype', 'timelock');
% check if the input cfg is valid for this function
cfg = ft_checkconfig(cfg, 'forbidden', {'channels'}); % prevent accidental typos, see issue 1729
% set the defaults
cfg.channel = ft_getopt(cfg, 'channel', 'all');
cfg.feedback = ft_getopt(cfg, 'method', 'text');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% the actual computation is done in the middle part
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ensure that the same channels are in the data and template
cfg.channel = ft_channelselection(cfg.channel, data.label);
cfg.channel = ft_channelselection(cfg.channel, template.label);
tmpcfg = keepfields(cfg, {'channel'});
data = ft_selectdata(tmpcfg, data);
template = ft_selectdata(tmpcfg, template);
% restore the provenance information
[cfg, data] = rollback_provenance(cfg, data);
[cfg, template] = rollback_provenance(cfg, template);
ntrial = length(data.trial);
nchan = length(cfg.channel);
nartifact = size(cfg.artifact,1);
ft_progress('init', cfg.feedback, 'removing artifacts');
for i=1:ntrial
datbegsample = data.sampleinfo(i,1);
datendsample = data.sampleinfo(i,2);
dattrllength = datendsample-datbegsample+1;
model = zeros(nchan, dattrllength);
count = 0;
for j=1:nartifact
artbegsample = cfg.artifact(j,1);
artendsample = cfg.artifact(j,2);
if artendsample>=datbegsample && artbegsample<=datendsample
% one of the artifacts overlaps with this trial
count = count + 1;
% express the artifact relative to the trial
artbegsample = artbegsample-datbegsample+1;
artendsample = artendsample-datbegsample+1;
% the artifact might partially fall outside the trial, in that case it needs to be trimmed
artbegtrim = 0;
artendtrim = 0;
if artbegsample<1
artbegtrim = 1 - artbegsample;
artbegsample = 1;
end
if artendsample>dattrllength
artendtrim = artendsample - dattrllength;
artendsample = dattrllength;
end
% insert the trimmed template in the model for this trial
model(:, artbegsample:artendsample) = template.avg(:, (1+artbegtrim):(end-artendtrim));
end
end % for each artifact
ft_progress(i/ntrial, 'removing %d artifacts from trial %d of %d\n', count, i, ntrial);
% remove the artifact model from the actual data
data.trial{i} = data.trial{i} - model;
end % for each trial
ft_progress('close');
% do the general cleanup and bookkeeping at the end of the function
ft_postamble debug
ft_postamble previous data template
ft_postamble provenance data
ft_postamble history data
ft_postamble savevar data