Get values from a fiber group structure
val = fefgGet(fg,param,varargin)
Copyright (2013-2014), Franco Pestilli, Stanford University, pestillifranco@gmail.com.
------------
Returns the length of each fiber in the fiber group
flen = fefgGet(fg,'length')
------------
Compute the Westin shape parameters (linearity and planarity)
for all the fibers in the fiber group.
Requires a dt structure.
westin = fefgGet(fg,'westinshape',dt)
westin = fefgGet(fg,'westinshape',eigenvals)
westin = fefgGet(fg,'westinshape',dtFileName)
-------------
Compute eigenvalues for the fibers in a fiber group.
It requires a dt structure.
eigenvals = fefgGet(fg,'eigenvals',dt)
eigenvals = fefgGet(fg,'eigenvals',dtFileName)
-------------
Compute the 6 parameters for the tensors at each node in each fiber.
It requires a dt structure.
dt6 = fefgGet(fg,'eigenvals',dt)
dt6 = fefgGet(fg,'eigenvals',dtFileName)
-------------
Compute the radial and axial ADC for all the fibers in the fiber group.
Requires a dt structure.
adc = fefgGet(fg,'adc',dt)
adc = fefgGet(fg,'adc',eigenvals)
adc = fefgGet(fg,'adc',dtFileName)
-------------
Compute the FA for all the fibers in the fiber group.
Requires a dt structure.
fa = fefgGet(fg,'fa',dt)
fa = fefgGet(fg,'fa',eigenvals)
fa = fefgGet(fg,'fa',dtFileName)
-------------
Compute the Mean Diffusivity for all the fibers in the fiber group.
Requires a dt structure.
md = fefgGet(fg,'md',dt)
md = fefgGet(fg,'md',eigenvals)
md = fefgGet(fg,'md',dtFileName)
-------------
Compute the Radial Diffusivity for all the fibers in the fiber group.
Requires a dt structure.
rd = fefgGet(fg,'rd',dt)
rd = fefgGet(fg,'rd',eigenvals)
rd = fefgGet(fg,'rd',dtFileName)
-------------
Compute the Axial Diffusivity for all the fibers in the fiber group.
Requires a dt structure.
ad = fefgGet(fg,'ad',dt)
ad = fefgGet(fg,'ad',eigenvals)
ad = fefgGet(fg,'ad',dtFileName)
-------------
Compute the radial and axial ADC for all the fibers in the fiber group.
Requires a dt structure.
adc = fefgGet(fg,'adc',dt)
adc = fefgGet(fg,'adc',eigenvals)
adc = fefgGet(fg,'adc',dtFileName)
-------------
Parameters
General
'name'
'type'
'colorrgb'
'thickness'
'visible'
Fiber related
'nfibers'- Number of fibers in this group
'nodes per fiber' - Number of nodes per fiber.
'fibers' - Fiber coordinates
'fibernames'
'fiberindex'
Compute the Mean Diffusivity for all the fibers in the fiber group.
Requires a dt structure.
fa = fefgGet(fg,'fa',dt)
fa = fefgGet(fg,'fa',eigenvals)
ROI and image coord related
'unique image coords'
'nodes to imagecoords' -
'voxel2fiber node pairs' - For each roi coord, an Nx2 matrix of
(fiber number,node number)
'nodes in voxels' - Nodes inside the voxels of roi coords
'voxels in fg' - Cell array of the roiCoords touched by each fiber
'voxels2fibermatrix' - Binary matrix (voxels by fibers). 1s when a
fiber is in a voxel of the roiCoords (which are, sadly, implicit).
Tensor and tractography related
'tensors' - Tensors for each node
See also: fgCreate; fgSet
Copyright (2013-2014), Franco Pestilli, Stanford University, pestillifranco@gmail.com.0001 function val = fefgGet(fg,param,varargin) 0002 % Get values from a fiber group structure 0003 % 0004 % val = fefgGet(fg,param,varargin) 0005 % 0006 % Copyright (2013-2014), Franco Pestilli, Stanford University, pestillifranco@gmail.com. 0007 % 0008 %------------ 0009 % Returns the length of each fiber in the fiber group 0010 % flen = fefgGet(fg,'length') 0011 %------------ 0012 % Compute the Westin shape parameters (linearity and planarity) 0013 % for all the fibers in the fiber group. 0014 % Requires a dt structure. 0015 % westin = fefgGet(fg,'westinshape',dt) 0016 % westin = fefgGet(fg,'westinshape',eigenvals) 0017 % westin = fefgGet(fg,'westinshape',dtFileName) 0018 %------------- 0019 % Compute eigenvalues for the fibers in a fiber group. 0020 % It requires a dt structure. 0021 % eigenvals = fefgGet(fg,'eigenvals',dt) 0022 % eigenvals = fefgGet(fg,'eigenvals',dtFileName) 0023 %------------- 0024 % Compute the 6 parameters for the tensors at each node in each fiber. 0025 % It requires a dt structure. 0026 % dt6 = fefgGet(fg,'eigenvals',dt) 0027 % dt6 = fefgGet(fg,'eigenvals',dtFileName) 0028 %------------- 0029 % Compute the radial and axial ADC for all the fibers in the fiber group. 0030 % Requires a dt structure. 0031 % adc = fefgGet(fg,'adc',dt) 0032 % adc = fefgGet(fg,'adc',eigenvals) 0033 % adc = fefgGet(fg,'adc',dtFileName) 0034 %------------- 0035 % Compute the FA for all the fibers in the fiber group. 0036 % Requires a dt structure. 0037 % fa = fefgGet(fg,'fa',dt) 0038 % fa = fefgGet(fg,'fa',eigenvals) 0039 % fa = fefgGet(fg,'fa',dtFileName) 0040 %------------- 0041 % Compute the Mean Diffusivity for all the fibers in the fiber group. 0042 % Requires a dt structure. 0043 % md = fefgGet(fg,'md',dt) 0044 % md = fefgGet(fg,'md',eigenvals) 0045 % md = fefgGet(fg,'md',dtFileName) 0046 %------------- 0047 % Compute the Radial Diffusivity for all the fibers in the fiber group. 0048 % Requires a dt structure. 0049 % rd = fefgGet(fg,'rd',dt) 0050 % rd = fefgGet(fg,'rd',eigenvals) 0051 % rd = fefgGet(fg,'rd',dtFileName) 0052 %------------- 0053 % Compute the Axial Diffusivity for all the fibers in the fiber group. 0054 % Requires a dt structure. 0055 % ad = fefgGet(fg,'ad',dt) 0056 % ad = fefgGet(fg,'ad',eigenvals) 0057 % ad = fefgGet(fg,'ad',dtFileName) 0058 %------------- 0059 % Compute the radial and axial ADC for all the fibers in the fiber group. 0060 % Requires a dt structure. 0061 % adc = fefgGet(fg,'adc',dt) 0062 % adc = fefgGet(fg,'adc',eigenvals) 0063 % adc = fefgGet(fg,'adc',dtFileName) 0064 %------------- 0065 % 0066 % Parameters 0067 % General 0068 % 'name' 0069 % 'type' 0070 % 'colorrgb' 0071 % 'thickness' 0072 % 'visible' 0073 % 0074 % Fiber related 0075 % 'nfibers'- Number of fibers in this group 0076 % 'nodes per fiber' - Number of nodes per fiber. 0077 % 'fibers' - Fiber coordinates 0078 % 'fibernames' 0079 % 'fiberindex' 0080 % Compute the Mean Diffusivity for all the fibers in the fiber group. 0081 % Requires a dt structure. 0082 % fa = fefgGet(fg,'fa',dt) 0083 % fa = fefgGet(fg,'fa',eigenvals) 0084 % 0085 % ROI and image coord related 0086 % 'unique image coords' 0087 % 'nodes to imagecoords' - 0088 % 'voxel2fiber node pairs' - For each roi coord, an Nx2 matrix of 0089 % (fiber number,node number) 0090 % 'nodes in voxels' - Nodes inside the voxels of roi coords 0091 % 'voxels in fg' - Cell array of the roiCoords touched by each fiber 0092 % 'voxels2fibermatrix' - Binary matrix (voxels by fibers). 1s when a 0093 % fiber is in a voxel of the roiCoords (which are, sadly, implicit). 0094 % 0095 % Tensor and tractography related 0096 % 'tensors' - Tensors for each node 0097 % 0098 % 0099 % See also: fgCreate; fgSet 0100 % 0101 % Copyright (2013-2014), Franco Pestilli, Stanford University, pestillifranco@gmail.com. 0102 val = []; 0103 0104 switch strrep(lower(param),' ','') 0105 0106 % Basic fiber parameters 0107 case 'name' 0108 val = fg.name; 0109 case 'type' % Should always be fibergroup 0110 val = fg.type; 0111 0112 % Fiber visualization settings. 0113 case 'colorrgb' 0114 val = fg.colorRgb; 0115 case 'thickness' 0116 val = fg.thickness; 0117 case 'visible' 0118 val = fg.visible; 0119 0120 % Simple fiber properties -- 0121 case {'fibers'} 0122 % val = fefgGet(fg,'fibers',fList); 0123 % 0124 % Returns a 3xN matrix of fiber coordinates corresponding to the 0125 % fibers specified in the integer vector, fList. This differs from 0126 % the dtiH (mrDiffusion) representation, where fiber coordinates 0127 % are stored as a set of cell arrays for each fiber. 0128 if ~isempty(varargin) 0129 list = varargin{1}; 0130 val = cell(length(list),1); 0131 for ii=1:length(list) 0132 val{ii} = fg.fibers{ii}; 0133 end 0134 else 0135 val = fg.fibers; 0136 end 0137 case 'fibernames' 0138 val = fg.fiberNames; 0139 case 'fiberindex' 0140 val = fg.fiberIndex; 0141 case 'nfibers' 0142 val = length(fg.fibers); 0143 case {'nodesperfiber','nsamplesperfiber','nfibersamples'} 0144 % fefgGet(fg,'n samples per fiber ') 0145 % How many samples per fiber. This is about equal to 0146 % their length in mm, though we need to write the fiber lengths 0147 % routine to actually calculate this. 0148 nFibers = fefgGet(fg,'n fibers'); 0149 val = zeros(1,nFibers); 0150 for ii=1:nFibers 0151 val(ii) = length(fg.fibers{ii}); 0152 end 0153 0154 % Fiber group (subgroup) properties. 0155 % These are used when we classify fibers into subgroups. We should 0156 % probably clean up this organization which is currently 0157 % 0158 % subgroup - length of fibers, an index of group identity 0159 % subgroupNames() 0160 % .subgroupIndex - Probably should go away and the index should 0161 % just be 0162 % .subgroupName - Probably should be moved up. 0163 % 0164 case {'ngroups','nsubgroups'} 0165 val = length(fg.subgroupNames); 0166 case {'groupnames'} 0167 val = cell(1,fefgGet(fg,'n groups')); 0168 for ii=1:nGroups 0169 val{ii} = fg.subgroupNames(ii).subgroupName; 0170 end 0171 0172 % DTI properties 0173 case 'tensors' 0174 val = fg.tensors; 0175 0176 % Fiber to coord calculations 0177 case {'imagecoords'} 0178 % c = fefgGet(fgAcpc,'image coords',fgList,xForm); 0179 % c = fefgGet(fgAcpc,'image coords',fgList,xForm); 0180 % 0181 % Return the image coordinates of a specified list of fibers 0182 % Returns a matrix that is fgList by 3 of the image coordinates for 0183 % each node of each fiber. 0184 % 0185 % Fiber coords are represented at fine resolution in ACPC space. 0186 % These coordinates are rounded and in image space 0187 if ~isempty(varargin) 0188 fList = varargin{1}; 0189 if length(varargin) > 1 0190 xForm = varargin{2}; 0191 % Put the fiber coordinates into image space 0192 fg = dtiXformFiberCoords(fg,xForm); 0193 end 0194 else 0195 % In this case, the fiber coords should already be in image 0196 % space. 0197 nFibers = fefgGet(fg,'n fibers'); 0198 fList = 1:nFibers; 0199 end 0200 0201 % Pull out the coordinates and ceil them. These are in image 0202 % space. 0203 nFibers = length(fList); 0204 val = cell(1,nFibers); 0205 if nFibers == 1 0206 val = ceil(fg.fibers{fList(1)}')+1; 0207 0208 else 0209 feOpenLocalCluster 0210 parfor ii=1:nFibers 0211 val{ii} = ceil(fg.fibers{fList(ii)}')+1; 0212 end 0213 end 0214 0215 case {'uniqueimagecoords'} 0216 % coords = fefgGet(fgIMG,'unique image coords'); 0217 % 0218 % The fg input must be in IMG space. 0219 % 0220 % Returns the unique image coordinates of all the fibers as an Nx3 0221 % matrix of integers. 0222 % val = round(horzcat(fg.fibers{:})'); 0223 val = ceil(horzcat(fg.fibers{:})')+1; 0224 val = unique(val,'rows'); 0225 0226 case {'allimagecoords'} 0227 % coords = fefgGet(fgIMG,'all image coords'); 0228 % 0229 % The fg input must be in IMG space. 0230 % 0231 % Returns all image coordinates of all the fibers as an Nx3 0232 % matrix of integers. 0233 % val = round(horzcat(fg.fibers{:})'); 0234 val = ceil(horzcat(fg.fibers{:})')+1; 0235 0236 case {'uniqueacpccoords'} 0237 % coords = fefgGet(fg,'unique acpc coords'); 0238 % 0239 % The fg input must be in ACPC space. 0240 % 0241 % Returns the unique ACPC coordinates of all the fibers as an Nx3 0242 % matrix of integers. 0243 val = fefgGet(fg, 'uniqueimagecoords') ; 0244 0245 case {'nodes2voxels'} 0246 % nodes2voxels = fefgGet(fgImg,'nodes2voxels',roiCoords) 0247 % 0248 % The roiCoords are a matrix of Nx3 coordinates. They describe a 0249 % region of interest, typically in image space or possibly in acpc 0250 % space. 0251 % 0252 % We return a cell array that is a mapping of fiber nodes to voxels in 0253 % the roi. The roi is specified as an Nx3 matrix of coordinates. 0254 % The returned cell array, nodes2voxels, has the same number of 0255 % cells as there are fibers. 0256 % 0257 % Unlike the fiber group cells, which have a 3D coordinate of each 0258 % node, this cell array has an integer that indexes the row of 0259 % roiCoords that contains the node. If a node is not in any of the 0260 % roiCoords, the entry in node2voxels{ii} for that node is zero. 0261 % This means that node is outside the 'roiCoords'. 0262 % 0263 % Once again: The cell nodes2voxels{ii} specifies whether each 0264 % node in the iith fiber is inside a voxel in the roiCoords. The 0265 % value specifies the row in roiCoords that contains the node. 0266 % 0267 if isempty(varargin), error('roiCoords required'); 0268 else 0269 roiCoords = varargin{1}; 0270 end 0271 fprintf('[%s] Computing nodes-to-voxels..',mfilename) 0272 0273 % Find the roiCoord for each node in each fiber. 0274 nFiber = fefgGet(fg,'n fibers'); 0275 val = cell(nFiber,1); 0276 0277 feOpenLocalCluster 0278 parfor ii=1:nFiber 0279 % if ~mod(ii,200), fprintf('%d ',ii); end 0280 % Node coordinates in image space 0281 nodeCoords = fefgGet(fg,'image coords',ii); 0282 0283 % The values in loc are the row of the coords matrix that contains 0284 % that sample point in a fiber. For example, if the number 100 is 0285 % in the 10th position of loc, then the 10th sample point in the 0286 % fiber passes through the voxel in row 100 of coords. 0287 %keyboard 0288 [~, val{ii}] = ismember(nodeCoords, roiCoords, 'rows'); % This operation is slow. 0289 end 0290 fprintf(' took: %2.3fminutes.\n',toc/60) 0291 0292 case {'voxel2fibernodepairs','v2fn'} 0293 % voxel2FNpairs = fefgGet(fgImg,'voxel 2 fibernode pairs',roiCoords); 0294 % voxel2FNpairs = fefgGet(fgImg,'voxel 2 fibernode pairs',roiCoords,nodes2voxels); 0295 % 0296 % The return is a cell array whose size is the number of voxels. 0297 % The cell is a Nx2 matrix of the (fiber, node) pairs that pass 0298 % through it. 0299 % 0300 % The value N is the number of nodes in the voxel. The first 0301 % column is the fiber number. The second column reports the indexes 0302 % of the nodes for each fiber in each voxel. 0303 fprintf('[%s] Computing fibers/nodes pairing in each voxel..\n',mfilename) 0304 0305 if (length(varargin) < 1), error('Requires the roiCoords.'); 0306 else 0307 roiCoords = varargin{1}; 0308 nCoords = size(roiCoords,1); 0309 end 0310 if length(varargin) < 2 0311 % We assume the fg and the ROI coordinates are in the same 0312 % coordinate frame. 0313 tic,nodes2voxels = fefgGet(fg,'nodes 2 voxels',roiCoords); 0314 else nodes2voxels = varargin{2}; 0315 end 0316 0317 nFibers = fefgGet(fg,'nFibers'); 0318 voxelsInFG = fefgGet(fg,'voxels in fg',nodes2voxels); 0319 0320 tic,roiNodesInFG = fefgGet(fg,'nodes in voxels',nodes2voxels); 0321 tic, val = cell(1,nCoords); 0322 for thisFiber=1:nFibers 0323 voxelsInFiber = voxelsInFG{thisFiber}; % A few voxels, in a list 0324 nodesInFiber = roiNodesInFG{thisFiber}; % The corresponding nodes 0325 0326 % Then add a row for each (fiber,node) pairs that pass through 0327 % the voxels for this fiber. 0328 for jj=1:length(voxelsInFiber) 0329 thisVoxel = voxelsInFiber(jj); 0330 % Print out roi coord and fiber coord to verify match 0331 % roiCoords(thisVoxel,:) 0332 % fg.fibers{thisFiber}(:,nodesInFiber(jj)) 0333 % Would horzcat be faster? 0334 val{thisVoxel} = cat(1,val{thisVoxel},[thisFiber,nodesInFiber(jj)]); 0335 end 0336 end 0337 fprintf('[%s] fiber/node pairing completed in: %2.3fs.\n',mfilename, toc) 0338 0339 case {'voxel2fibers','fiberdensity'} 0340 % voxel2FNpairs = fefgGet(fgImg,'fiber density',roiCoords); 0341 % 0342 % The return is a cell array whose size is the number of voxels. 0343 % The cell is a Nx1 matrix of fiber counts. How many fibers in each 0344 % voxel. 0345 % 0346 fprintf('[%s] Computing fiber density in each voxel...\n',mfilename) 0347 0348 if (length(varargin) < 1), 0349 roiCoords = fefgGet(fg,'uniqueimagecoords'); 0350 fprintf('[%s] Computing white-matter volume ROI from fibers coordinates.\n',mfilename) 0351 fprintf(' Assuming fiber coordinates in IMG space.\n') 0352 else 0353 roiCoords = varargin{1}; 0354 end 0355 0356 if length(varargin) < 2 0357 % We assume the fg and the ROI coordinates are in the same 0358 % coordinate frame. 0359 tic,nodes2voxels= fefgGet(fg,'nodes 2 voxels',roiCoords); 0360 else nodes2voxels = varargin{2}; 0361 end 0362 0363 nCoords = size(roiCoords,1); 0364 nFibers = fefgGet(fg,'nFibers'); 0365 voxelsInFG = fefgGet(fg,'voxels in fg',nodes2voxels); 0366 0367 tic, fibersInVox = cell(1,nCoords); 0368 for thisFiber=1:nFibers 0369 voxelsInFiber = voxelsInFG{thisFiber}; % A few voxels, in a list 0370 0371 % Then add a row for each (fiber,node) pairs that pass through 0372 % the voxels for this fiber. 0373 for jj=1:length(voxelsInFiber) 0374 thisVoxel = voxelsInFiber(jj); 0375 % Print out roi coord and fiber coord to verify match 0376 % roiCoords(thisVoxel,:) 0377 % fg.fibers{thisFiber}(:,nodesInFiber(jj)) 0378 % Would horzcat be faster? 0379 fibersInVox{thisVoxel} = cat(1,fibersInVox{thisVoxel},thisFiber); 0380 end 0381 end 0382 0383 % Now create the fiber density, the unique fibers in each voxel 0384 val = zeros(length(fibersInVox),1); 0385 for ivx = 1:length(fibersInVox) 0386 val(ivx) = length(unique(fibersInVox{ivx})); 0387 end 0388 0389 fprintf('[%s] fiber density completed in: %2.3fs.\n',mfilename, toc) 0390 0391 case {'uniquefibersinvox'} 0392 % uniquefvx = fefgGet(fgImg,'uniquefibrsinvox',roiCoords); 0393 % 0394 % The return is a vector whose size is the number of voxels containing 0395 % the unique fibers in each voxel 0396 % 0397 fprintf('[%s] Computing the unique fibers in each voxel...\n',mfilename) 0398 tic 0399 if (length(varargin) < 1), error('Requires the roiCoords.'); 0400 else 0401 roiCoords = varargin{1}; 0402 nCoords = size(roiCoords,1); 0403 end 0404 0405 % We assume the fg and the ROI coordinates are in the same 0406 % coordinate frame. 0407 nodes2voxels = fefgGet(fg,'nodes 2 voxels',roiCoords); 0408 0409 nFibers = fefgGet(fg,'nFibers'); 0410 voxelsInFG = fefgGet(fg,'voxels in fg',nodes2voxels); 0411 0412 fibersInVox = cell(1,nCoords); 0413 for thisFiber=1:nFibers 0414 voxelsInFiber = voxelsInFG{thisFiber}; % A few voxels, in a list 0415 0416 % Then add a row for each (fiber,node) pairs that pass through 0417 % the voxels for this fiber. 0418 if ~isempty(voxelsInFiber) 0419 for jj=1:length(voxelsInFiber) 0420 thisVoxel = voxelsInFiber(jj); 0421 fibersInVox{thisVoxel} = cat(1,fibersInVox{thisVoxel},thisFiber); 0422 end 0423 end 0424 end 0425 0426 % Now create find the unique fibers in each voxel 0427 val = cell(length(fibersInVox),1); 0428 for ivx = 1:length(fibersInVox) 0429 val{ivx} = (unique(fibersInVox{ivx})); 0430 end 0431 fprintf('[%s] done computing unique fibers in each voxel: %2.3fs.\n',mfilename, toc) 0432 0433 case {'nodesinvoxels'} 0434 % nodesInVoxels = fefgGet(fg,'nodes in voxels',nodes2voxels); 0435 % 0436 % This cell array is a modified form of nodes2voxels (see above). 0437 % In that cell array every node in every fiber has a number 0438 % referring to its row in roiCoords, or a 0 when the node is not in 0439 % any roiCoord voxel. 0440 % 0441 % This cell array differs only in that the 0s removed. This 0442 % is used to simplify certain calculations. 0443 % 0444 if (length(varargin) <1), error('Requires nodes2voxels cell array.'); end 0445 tic,fprintf('[%s] Computing nodes-in-voxels..',mfilename) 0446 nodes2voxels = varargin{1}; 0447 nFibers = fefgGet(fg,'nFibers'); 0448 val = cell(1,nFibers); 0449 0450 feOpenLocalCluster 0451 % For each fiber, this is a list of the nodes that pass through 0452 % a voxel in the roiCoords 0453 parfor ii = 1:nFibers 0454 % For each fiber, this is a list of the nodes that pass through 0455 % a voxel in the roiCoords 0456 lst = (nodes2voxels{ii} ~= 0); 0457 val{ii} = find(lst); 0458 end 0459 fprintf(' took: %2.3fs.\n',toc) 0460 0461 case 'voxelsinfg' 0462 % voxelsInFG = fefgGet(fgImg,'voxels in fg',nodes2voxels); 0463 % 0464 % A cell array length n-fibers. Each cell has a list of the voxels 0465 % (rows of roiCoords) for a fiber. 0466 % 0467 % This routine eliminates the 0's in the nodes2voxels lists. 0468 % 0469 if length(varargin) < 1, error('Requires nodes2voxels cell array.'); end 0470 tic,fprintf('[%s] Computing voxels-in-fg..',mfilename) 0471 0472 nodes2voxels = varargin{1}; 0473 nFibers = fefgGet(fg,'nFibers'); 0474 val = cell(1,nFibers); 0475 0476 feOpenLocalCluster 0477 parfor ii = 1:nFibers 0478 % These are the nodes that pass through a voxel in the 0479 % roiCoords 0480 lst = (nodes2voxels{ii} ~= 0); 0481 val{ii} = nodes2voxels{ii}(lst); 0482 end 0483 fprintf(' took: %2.3fs.\n',toc) 0484 0485 case {'voxels2fibermatrix','v2fm'} 0486 % v2fm = fefgGet(fgImg,'voxels 2 fiber matrix',roiCoords); 0487 % Or, 0488 % v2fnPairs = fefgGet(fgImg,'v2fn',roiCoords); 0489 % v2fm = fefgGet(fgImg,'voxels 2 fiber matrix',roiCoords, v2fnPairs); 0490 % 0491 % mrvNewGraphWin; imagesc(v2fm) 0492 % 0493 % Returns a binary matrix of size Voxels by Fibers. 0494 % When voxel ii has at least one node from fiber jj, there is a one 0495 % in v2fm(ii,jj). Otherwise, the entry is zero. 0496 % 0497 0498 % Check that the fg is in the image coordspace: 0499 if isfield(fg, 'coordspace') && ~strcmp(fg.coordspace, 'img') 0500 error('Fiber group is not in the image coordspace, please xform'); 0501 end 0502 0503 if isempty(varargin), error('roiCoords required'); 0504 else 0505 roiCoords = varargin{1}; 0506 nCoords = size(roiCoords,1); 0507 if length(varargin) < 2 0508 v2fnPairs = fefgGet(fg,'v2fn',roiCoords); 0509 else 0510 v2fnPairs = varargin{2}; 0511 end 0512 end 0513 0514 % Allocate matrix of voxels by fibers 0515 val = zeros(nCoords,fefgGet(fg,'n fibers')); 0516 0517 % For each coordinate, find the fibers. Set those entries to 1. 0518 for ii=1:nCoords 0519 if ~isempty(v2fnPairs{ii}) 0520 f = unique(v2fnPairs{ii}(:,1)); 0521 end 0522 val(ii,f) = 1; 0523 end 0524 0525 case {'fibersinroi','fginvoxels','fibersinvoxels'} 0526 % fList = fefgGet(fgImg,'fibersinroi',roiCoords); 0527 % 0528 % v2fn = fefgGet(fgImg,'v2fn',roiCoords); 0529 % fList = fefgGet(fgImg,'fibersinroi',roiCoords,v2fn); 0530 % 0531 % Returns an integer vector of the fibers with at least 0532 % one node in a region of interest. 0533 % 0534 % The fg and roiCoords should be in the same coordinate frame. 0535 % 0536 if isempty(varargin), error('roiCoords required'); 0537 elseif length(varargin) == 1 0538 roiCoords = varargin{1}; 0539 v2fnPairs = fefgGet(fg,'v2fn',roiCoords); 0540 elseif length(varargin) > 1 0541 roiCoords = varargin{1}; 0542 v2fnPairs = varargin{2}; 0543 end 0544 0545 val = []; nCoords = size(roiCoords,1); 0546 for ii=1:nCoords 0547 if ~isempty(v2fnPairs{ii}) 0548 val = cat(1,val,v2fnPairs{ii}(:,1)); 0549 end 0550 end 0551 val = sort(unique(val),'ascend'); 0552 0553 case {'coordspace','fibercoordinatespace','fcspace'} 0554 % In some cases, the fg might contain information telling us in which 0555 % coordinate space its coordinates are set. This information is set 0556 % as a struct. Each entry in the struct can be either a 4x4 xform 0557 % matrix from the fiber coordinates to that space (with eye(4) for 0558 % the space in which the coordinates are defined), or (if the xform 0559 % is not know) an empty matrix. 0560 0561 cspace_fields = fields(fg.coordspace); 0562 val = []; 0563 for f=1:length(cspace_fields) 0564 this_field = cspace_fields{f}; 0565 if isequal(getfield(fg.coordspace, this_field), eye(4)) 0566 val = this_field; 0567 end 0568 end 0569 0570 case {'length'} 0571 % Returns the length of each fiber in the fiber group 0572 % flen = fefgGet(fg,'length') 0573 0574 % Measure the step size of the first fiber. They *should* all be the same! 0575 stepSize = mean(sqrt(sum(diff(fg.fibers{1},1,2).^2))); 0576 0577 % Check that they are 0578 %stepSize2 = mean(sqrt(sum(diff(fg.fibers{2},1,2).^2))); 0579 %assertElementsAlmostEqual(stepSize,stepSize2,'relative',.001,.001); 0580 0581 % Estimate the length for the fibers, as well mean, min and max 0582 val = cellfun('length',fg.fibers)*stepSize; 0583 0584 case {'dt6'} 0585 % Compute the 6 parameters for the tensors at each node in each fiber. 0586 % It requires a dt structure. 0587 % dt6 = fefgGet(fg,'dt6',dt) 0588 % dt6 = fefgGet(fg,'dt6',dtFileName) 0589 0590 if ~ischar(varargin{1}) 0591 if isstruct(varargin{1}) 0592 % A dti structure was passed. 0593 dt = varargin{1}; 0594 else 0595 error('[%s] A ''dt'' structure is necessary.', mfilename) 0596 end 0597 else % The string should be a path to a dt.mat file. 0598 dt = dtiLoadDt6(varargin{1}); 0599 end 0600 0601 % Get the dt6 tensors for each node in each fiber. 0602 nFibers = fefgGet(fg,'nFibers'); 0603 val = cell(1,nFibers); 0604 xform = inv(dt.xformToAcpc); 0605 dt6 = dt.dt6; clear dt 0606 feOpenLocalCluster 0607 parfor ii = 1:nFibers 0608 % Get the fiber coordinates. 0609 coords = fg.fibers{ii}'; % Assumed in ACPC 0610 0611 [val1,val2,val3,val4,val5,val6] = ... 0612 dtiGetValFromTensors(dt6, coords, xform,'dt6','nearest'); 0613 0614 % Build a vector of tesnsors' eigenvalues 0615 val{ii} = [val1,val2,val3,val4,val5,val6]; 0616 end 0617 0618 0619 case {'eigenvals'} 0620 % Compute eigenvalues for the fibers in a fiber group. 0621 % It requires a dt structure. 0622 % eigenvals = fefgGet(fg,'eigenvals',dt) 0623 % eigenvals = fefgGet(fg,'eigenvals',dt6FileName) 0624 if ~ischar(varargin{1}) 0625 if ( ~isstruct(varargin{1}) ) 0626 dt6 = fefgGet(fg,'dt6',varargin{1}); 0627 elseif ( size(varargin{1},1)==6 ) 0628 dt6 = varargin{1}; 0629 else 0630 error('[%s] Either a ''dt'' structure or a ''dt6'' vector of tensor coefficients is necessary.\n',mfilename) 0631 end 0632 else % The string should be a path to a dt6.mat file. 0633 dt6 = fefgGet(fg,'dt6',dtiLoadDt6(varargin{1})); 0634 end 0635 0636 nFibers = fefgGet(fg,'nFibers'); 0637 val = cell(1,nFibers); 0638 0639 feOpenLocalCluster 0640 parfor ii = 1:nFibers 0641 0642 % We now have the dt6 data from all of the fibers. We extract the 0643 % directions into vec and the eigenvalues into val. The units of val are 0644 % um^2/sec or um^2/msec .. somebody answer this here, please. 0645 [~,val{ii}] = dtiEig(dt6{ii}); 0646 end 0647 0648 for ii = 1:nFibers 0649 % Some of the ellipsoid fits are wrong and we get negative eigenvalues. 0650 % These are annoying. If they are just a little less than 0, then clipping 0651 % to 0 is not an entirely unreasonable thing. Maybe we should check for the 0652 % magnitude of the error? 0653 nonPD = find(any(val{ii}<0,2), 1); 0654 if(~isempty(nonPD)) 0655 %fprintf('\n NOTE: %d fiber points had negative eigenvalues. These will be clipped to 0..\n',length(nonPD)); 0656 val{ii}(val{ii}<0) = 0; 0657 end 0658 0659 threeZeroVals = (sum(val{ii}, 2) == 0); 0660 %if ~isempty (threeZeroVals) 0661 % fprintf('\n NOTE: %d of these fiber points had all three negative eigenvalues. These will be excluded from analyses\n', length(threeZeroVals)); 0662 %end 0663 val{ii}(threeZeroVals, :)=[]; 0664 end 0665 0666 case {'fa'} 0667 % Compute the FA for all the fibers in the fiber group. 0668 % Requires a dt structure. 0669 % fa = fefgGet(fg,'fa',dt) 0670 % fa = fefgGet(fg,'fa',eigenvals) 0671 % fa = fefgGet(fg,'fa',dtFileName) 0672 if (~isstruct(varargin{1}) || ischar(varargin{1})) && ~iscell(varargin{1}) 0673 % A dti structue was passed. 0674 eigenvals = fefgGet(fg,'eigenvals',varargin{1}); 0675 else 0676 % We assume that eigenvals were passed already 0677 eigenvals = varargin{1}; 0678 end 0679 0680 nFibers = fefgGet(fg,'nFibers'); 0681 val = cell(1,nFibers); 0682 0683 feOpenLocalCluster 0684 parfor ii = 1:nFibers 0685 [val{ii},~,~,~] = dtiComputeFA(eigenvals{ii}); 0686 end 0687 0688 case {'md'} 0689 % Compute the Mean Diffusivity for all the fibers in the fiber group. 0690 % Requires a dt structure. 0691 % md = fefgGet(fg,'md',dt) 0692 % md = fefgGet(fg,'md',eigenvals) 0693 % md = fefgGet(fg,'md',dtFileName) 0694 0695 if (~isstruct(varargin{1}) || ischar(varargin{1})) && ~iscell(varargin{1}) 0696 % A dti structue was passed. 0697 eigenvals = fefgGet(fg,'eigenvals',varargin{1}); 0698 else 0699 % We assume that eigenvals were passed already 0700 eigenvals = varargin{1}; 0701 end 0702 nFibers = fefgGet(fg,'nFibers'); 0703 val = cell(1,nFibers); 0704 0705 feOpenLocalCluster 0706 parfor ii = 1:nFibers 0707 [~,val{ii},~,~] = dtiComputeFA(eigenvals{ii}); 0708 end 0709 0710 case {'ad'} 0711 % Compute the Axial Diffusivity for all the fibers in the fiber group. 0712 % Requires a dt structure. 0713 % ad = fefgGet(fg,'ad',dt) 0714 % ad = fefgGet(fg,'ad',eigenvals) 0715 % ad = fefgGet(fg,'ad',dtFileName) 0716 0717 if (~isstruct(varargin{1}) || ischar(varargin{1})) && ~iscell(varargin{1}) 0718 % A dti structue was passed. 0719 eigenvals = fefgGet(fg,'eigenvals',varargin{1}); 0720 else 0721 % We assume that eigenvals were passed already 0722 eigenvals = varargin{1}; 0723 end 0724 nFibers = fefgGet(fg,'nFibers'); 0725 val = cell(1,nFibers); 0726 0727 feOpenLocalCluster 0728 parfor ii = 1:nFibers 0729 [~,~,~,val{ii}] = dtiComputeFA(eigenvals{ii}); 0730 end 0731 0732 case {'rd'} 0733 % Compute the Radial Diffusivity for all the fibers in the fiber group. 0734 % Requires a dt structure. 0735 % rd = fefgGet(fg,'rd',dt) 0736 % rd = fefgGet(fg,'rd',eigenvals) 0737 % rd = fefgGet(fg,'rd',dtFileName) 0738 0739 if (~isstruct(varargin{1}) || ischar(varargin{1})) && ~iscell(varargin{1}) 0740 % A dti structue was passed. 0741 eigenvals = fefgGet(fg,'eigenvals',varargin{1}); 0742 else 0743 % We assume that eigenvals were passed already 0744 eigenvals = varargin{1}; 0745 end 0746 nFibers = fefgGet(fg,'nFibers'); 0747 val = cell(1,nFibers); 0748 0749 feOpenLocalCluster 0750 parfor ii = 1:nFibers 0751 [~,~,val{ii},~] = dtiComputeFA(eigenvals{ii}); 0752 end 0753 0754 case {'adc'} 0755 % Compute the radial and axial ADC for all the fibers in the fiber group. 0756 % Requires a dt structure. 0757 % adc = fefgGet(fg,'adc',dt) 0758 % adc = fefgGet(fg,'adc',eigenvals) 0759 % adc = fefgGet(fg,'adc',dtFileName) 0760 0761 if (~isstruct(varargin{1}) || ischar(varargin{1})) && ~iscell(varargin{1}) 0762 % A dti structue was passed. 0763 eigenvals = fefgGet(fg,'eigenvals',varargin{1}); 0764 else 0765 % We assume that eigenvals were passed already 0766 eigenvals = varargin{1}; 0767 end 0768 0769 nFibers = fefgGet(fg,'nFibers'); 0770 val = cell(1,nFibers); 0771 feOpenLocalCluster 0772 parfor ii = 1:nFibers 0773 [~,~,val{ii}.radial, val{ii}.axial] = dtiComputeFA(eigenvals{ii}); 0774 end 0775 0776 case {'westinshape'} 0777 % Compute the Westin shape parameters (linearity and planarity) 0778 % for all the fibers in the fiber group. 0779 % Requires a dt structure. 0780 % westin = fefgGet(fg,'westinshape',dt) 0781 % westin = fefgGet(fg,'westinshape',eigenvals) 0782 % westin = fefgGet(fg,'westinshape',dtFileName) 0783 0784 if (~isstruct(varargin{1}) || ischar(varargin{1})) && ~iscell(varargin{1}) 0785 % A dt structue was passed. 0786 eigenvals = fefgGet(fg,'eigenvals',varargin{1}); 0787 else 0788 % We assume that eigenvals were passed already 0789 eigenvals = varargin{1}; 0790 end 0791 0792 % This was the originial formulation described in: 0793 % C-F. Westin, S. Peled, H. Gubbjartsson, R. Kikinis, and F.A. Jolesz. 0794 % Geometrical diffusion measures for MRI from tensor basis analysis. 0795 % In Proceedings 5th Annual ISMRM, 1997. 0796 nFibers = fefgGet(fg,'nFibers'); 0797 val = cell(1,nFibers); 0798 0799 feOpenLocalCluster 0800 parfor ii = 1:nFibers 0801 [val{ii}.linearity, val{ii}.planarity] = dtiComputeWestinShapes(eigenvals{ii}); 0802 end 0803 0804 otherwise 0805 error('Unknown fg parameter: "%s"\n',param); 0806 end 0807 0808 return