Compute average sensor array
The following Matlab code will produce the average sensor array from a list of input files, and select the input file which is closest to the average (which I would recommend to use for interpolation). As input files, you can for example specify the on-line averages for your subjects.
% From a list of fiff-files (e.g. on-line average files), compute the one % with average sensor transformation, and find the input file which is % closest to the average (according to both translation and rotation parameters) % file list should be specified in cell array "files{}" % output file should be specified in string "file_out" % both average file "_avg.fif" and file closest to it will be written % if not output file specified, output will be written to current working directory % OH, March 2009
addpath /opt/mne/matlab/toolbox/; % fiff read/write tools
if exist('file_out')~=1, % if not output file specified: write to current working directory
- cwd = pwd; file_out = fullfile( cwd, 'FIFF4TRANS.fif' );
end;
if exist('files')~=1, % you can specify a list of files here (e.g. on-line averages for all subjects in your study)
- files = {'/fullpath/file4subj1.fif', ...
- '/fullpath/file4subj2.fif', ... '/fullpath/file4subj3.fif'};
end;
fid = fopen(file_out, 'a'); % check if output can be written if fid == -1,
- fprintf(1, 'Cannot access output file %s\n', file_out); return;
end; fclose(fid);
nr_files = length(files); epoch_gm = [0]; trans_gm = [0]; for f = 1:nr_files, % read all files and average all that's relevant
- fprintf(1, '%s\n', files{f}); data = fiff_read_evoked(files{f}); epoch_gm = epoch_gm + data.evoked.epochs; % average MEG data epoch trans_gm = trans_gm + data.info.dev_head_t.trans; % average sensor transformation matrix all_trans(:,:,f) = data.info.dev_head_t.trans; % remember individual transformations
end; epoch_gm = epoch_gm/nr_files; trans_gm = trans_gm/nr_files;
data_new = data; data_new.evoked.epoch = epoch_gm; data_new.info.dev_head_t.trans = trans_gm;
[thispath, thisfile,thisext,thisversn] = fileparts(file_out); file_out_avg = fullfile(thispath, [thisfile '_avg' thisext]); fiff_write_evoked(file_out_avg, data_new); % Write average coordinate transformation
% find input file that is closest to average... for f = 1:nr_files,
- diff_trans_t = all_trans(1:3,4,f) - trans_gm(1:3,4); % translation difference diff_trans_r = all_trans(1:3,1:3,f) - trans_gm(1:3,1:3); % rotation difference norm_t(f) = norm( diff_trans_t ); norm_r(f) = norm( diff_trans_r );
end;
[y_t, i_t] = sort( norm_t ); % sort according to translation difference [y_r, i_r] = sort( norm_r ); % sort according to rotation difference
for f = 1:nr_files, % combine the two sorted lists
- rank_t(f) = find( (i_t-f)==0 ); % rank for translation rank_r(f) = find( (i_r-f)==0 ); % rank for rotation avg_rank(f) = (rank_t(f) + rank_r(f))/2; [min_val, min_file] = min( avg_rank ); % best average rank
end;
min_file = files{min_file}; % filename for file with best average rank
fprintf('\nThe winner is... %s\n', min_file);
copyfile(min_file, file_out); % copy file with best average rank to output file
% The end of this