/
ft_determine_units.m
178 lines (148 loc) · 6.28 KB
/
ft_determine_units.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
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
function [obj] = ft_determine_units(obj)
% FT_DETERMINE_UNITS tries to determine the units of a geometrical object by
% looking at its size and by relating this to the approximate size of the
% human head according to the following table:
% from 0.050 to 0.500 -> meter
% from 0.500 to 5.000 -> decimeter
% from 5.000 to 50.000 -> centimeter
% from 50.000 to 500.000 -> millimeter
%
% Use as
% [output] = ft_determine_units(input)
%
% The following input data structures are supported
% electrode or gradiometer array, see FT_DATATYPE_SENS
% volume conduction model, see FT_DATATYPE_HEADMODEL
% source model, see FT_DATATYPE_SOURCE and FT_PREPARE_SOURCEMODEL
% anatomical mri, see FT_DATATYPE_VOLUME
% segmented mri, see FT_DATATYPE_SEGMENTATION
% anatomical or functional atlas, see FT_READ_ATLAS
%
% This function will add the field 'unit' to the output data structure with the
% possible values 'm', 'cm ' or 'mm'.
%
% See also FT_CONVERT_UNITS, FT_DETERMINE_COODSYS, FT_CONVERT_COORDSYS, FT_PLOT_AXES, FT_PLOT_XXX
% Copyright (C) 2005-2020, 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$
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if isstruct(obj) && numel(obj)>1
% deal with a structure array
for i=1:numel(obj)
tmp(i) = ft_determine_units(obj(i));
end
obj = tmp;
return
elseif iscell(obj) && numel(obj)>1
% deal with a cell-array
% this might represent combined EEG, ECoG and/or MEG
for i=1:numel(obj)
obj{i} = ft_determine_units(obj{i});
end
return
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if isfield(obj, 'unit') && ~isempty(obj.unit)
% use the units specified in the object
unit = obj.unit;
elseif isfield(obj, 'bnd') && isfield(obj.bnd, 'unit')
unit = unique({obj.bnd.unit});
if ~all(strcmp(unit, unit{1}))
ft_error('inconsistent units in the individual boundaries');
else
unit = unit{1};
end
% keep one representation of the units rather than keeping it with each boundary
% the units will be reassigned further down
obj.bnd = rmfield(obj.bnd, 'unit');
else
% try to determine the units by looking at the size of the object
if isfield(obj, 'chanpos') && ~isempty(obj.chanpos) && ~all(isnan(obj.chanpos(:)))
siz = norm(idrange(obj.chanpos));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'elecpos') && ~isempty(obj.elecpos) && ~all(isnan(obj.elecpos(:)))
siz = norm(idrange(obj.elecpos));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'coilpos') && ~isempty(obj.coilpos) && ~all(isnan(obj.coilpos(:)))
siz = norm(idrange(obj.coilpos));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'optopos') && ~isempty(obj.optopos) && ~all(isnan(obj.optopos(:)))
siz = norm(idrange(obj.optopos));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'pnt') && ~isempty(cat(1, obj.pnt))
% the obj can be a struct.array, hence the concatenation
siz = norm(idrange(cat(1, obj.pnt)));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'pos') && ~isempty(cat(1, obj.pos))
% the obj can be a struct.array, hence the concatenation
siz = norm(idrange(cat(1, obj.pos)));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'transform') && ~isempty(obj.transform)
% construct the corner points of the volume in voxel and in head coordinates
[pos_voxel, pos_head] = cornerpoints(obj.dim, obj.transform);
siz = norm(idrange(pos_head));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'fid') && isfield(obj.fid, 'pnt') && ~isempty(obj.fid.pnt)
siz = norm(idrange(obj.fid.pnt));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'fid') && isfield(obj.fid, 'pos') && ~isempty(obj.fid.pos)
siz = norm(idrange(obj.fid.pos));
unit = ft_estimate_units(siz);
elseif ft_headmodeltype(obj, 'infinite')
% this is an infinite medium volume conductor, which does not care about units
unit = 'm';
elseif ft_headmodeltype(obj,'singlesphere')
siz = obj.r;
unit = ft_estimate_units(siz);
elseif ft_headmodeltype(obj,'localspheres')
siz = median(obj.r);
unit = ft_estimate_units(siz);
elseif ft_headmodeltype(obj,'concentricspheres')
siz = max(obj.r);
unit = ft_estimate_units(siz);
elseif isfield(obj, 'bnd') && isstruct(obj.bnd) && isfield(obj.bnd(1), 'pnt') && ~isempty(obj.bnd(1).pnt)
siz = norm(idrange(obj.bnd(1).pnt));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'bnd') && isstruct(obj.bnd) && isfield(obj.bnd(1), 'pos') && ~isempty(obj.bnd(1).pos)
siz = norm(idrange(obj.bnd(1).pos));
unit = ft_estimate_units(siz);
elseif isfield(obj, 'nas') && isfield(obj, 'lpa') && isfield(obj, 'rpa')
pnt = [obj.nas; obj.lpa; obj.rpa];
siz = norm(idrange(pnt));
unit = ft_estimate_units(siz);
else
ft_error('cannot determine geometrical units');
end % recognized type of volume conduction model or sensor array
end % determine input units
% add the units to the output object
for i=1:numel(obj)
% obj can be a struct-array, hence the loop
obj(i).unit = unit;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% IDRANGE interdecile range for more robust range estimation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function r = idrange(x)
keeprow=true(size(x,1),1);
for l=1:size(x,2)
keeprow = keeprow & isfinite(x(:,l));
end
sx = sort(x(keeprow,:), 1);
ii = round(interp1([0, 1], [1, size(x(keeprow,:), 1)], [.1, .9])); % indices for 10 & 90 percentile
r = diff(sx(ii, :));