getdata

 getdata method - fetches time series data for ROI from images 
 FORMAT [Y multv vXYZ mat] = getdata(roiobj, data_imgs, flags)

 roiobj    - an object of type maroi
 data_imgs - images to fetch data from.  These can be in the form of
             a character array, or an array of type spm_vol (see
             spm_vol.m)

 flags can can be none or more of
              z - use zero masking for images without NaN represenation
              n - nearest neighbour resampling of images
              s - sinc resampling of images (why?)
              (trilinear resampling is the default)
              m - remap images
              l - Leave in columns with missing data

 If the resampling is not set with the flags input, then we use the resampling
 value from the ROI ``spm_hold`` value.

 default flags is empty

 Returns
 Y        - no of images x no of voxels in ROI data matrix
 multv    - weighting values from ROI (which have not been applied)
 vXYZ     - voxel coordinates of ROI from first image in series
 mat      - voxels -> mm mat file, again from first in series

 Matthew Brett 8/11/99, 2/8/01 (JBCP)

 $Id$

0001 function [Y, multv, vXYZ, mat] = getdata(roiobj, data_imgs, flags)
0002 % getdata method - fetches time series data for ROI from images
0003 % FORMAT [Y multv vXYZ mat] = getdata(roiobj, data_imgs, flags)
0004 %
0005 % roiobj    - an object of type maroi
0006 % data_imgs - images to fetch data from.  These can be in the form of
0007 %             a character array, or an array of type spm_vol (see
0008 %             spm_vol.m)
0009 %
0010 % flags can can be none or more of
0011 %              z - use zero masking for images without NaN represenation
0012 %              n - nearest neighbour resampling of images
0013 %              s - sinc resampling of images (why?)
0014 %              (trilinear resampling is the default)
0015 %              m - remap images
0016 %              l - Leave in columns with missing data
0017 %
0018 % If the resampling is not set with the flags input, then we use the resampling
0019 % value from the ROI ``spm_hold`` value.
0020 %
0021 % default flags is empty
0022 %
0023 % Returns
0024 % Y        - no of images x no of voxels in ROI data matrix
0025 % multv    - weighting values from ROI (which have not been applied)
0026 % vXYZ     - voxel coordinates of ROI from first image in series
0027 % mat      - voxels -> mm mat file, again from first in series
0028 %
0029 % Matthew Brett 8/11/99, 2/8/01 (JBCP)
0030 %
0031 % $Id$
0032   
0033 if nargin < 2
0034   error('Need ROI and data images');
0035 end
0036 if nargin < 3
0037   flags = '';
0038 end
0039 if isempty(flags)
0040   flags = ' ';
0041 end
0042 
0043 if ischar(data_imgs)
0044   data_imgs = spm_vol(data_imgs);
0045 elseif ~isstruct(data_imgs)
0046   error('Input data files must be strings or structs')
0047 elseif any(flags == 'm')
0048   for i = 1:length(data_imgs)
0049     data_imgs(i) = spm_vol(data_imgs(i).fname);
0050   end
0051 end
0052 
0053 % resampling = set by ROI hold value by default
0054 if any(flags == 's')
0055   holdval = -11;
0056 elseif any(flags == 'n')
0057   holdval = 0;
0058 else % Not specified, use ROI default resampling value
0059   holdval = spm_hold(roiobj);
0060 end
0061 
0062 % NaN replacement
0063 if any(flags == 'z')
0064   zmask = 1;
0065 else
0066   zmask = 0;
0067 end
0068 
0069 % get real points corresponding to first image in series
0070 [XYZ multv] = realpts(roiobj, mars_space(data_imgs(1)));
0071 dlen = length(multv);
0072 if dlen == 0 % no points in space
0073   Y = [];
0074   return
0075 end
0076 XYZ = [XYZ; ones(size(multv))];
0077 
0078 % check for same dims etc - which could save a bag of time
0079 % Code a bit tricksy here to allow comparison of vector and 4x4 matrices
0080 % without doing loops
0081 %---------------------------------------------------------
0082 nimgs = length(data_imgs);
0083 dims = cat(3,data_imgs(:).dim);
0084 dims = dims(:, 1:3, :); % to allow for SPM2/SPM99 4 element dims
0085 chgflgs = any(diff(dims,1,3)) | any(any(diff(cat(3,data_imgs(:).mat),1,3)));
0086 chgflgs = [1; chgflgs(:)];
0087 
0088 % create return matrix
0089 Y = zeros(nimgs, dlen);
0090 
0091 for i = 1:nimgs
0092   % nan replacement
0093   i_type = mars_vol_utils('type', data_imgs(i));
0094   nanrep =  spm_type(i_type, 'nanrep');
0095 
0096   if chgflgs(i)  % images not the same, (re)get resample points
0097     ixyz = data_imgs(i).mat \ XYZ;
0098   end
0099   if i == 1; % record voxel XYZ for return
0100     vXYZ = ixyz(1:3,:);
0101     mat  = data_imgs(1).mat;
0102   end
0103   % resample data at voxel centres of ROI
0104   data = spm_sample_vol(data_imgs(i), ixyz(1,:),ixyz(2,:),ixyz(3,:),holdval);
0105   % clear out missing values
0106   if ~nanrep & zmask
0107     data(data == 0) = NaN;
0108   end
0109   % return all the values
0110   Y(i, :) = data;
0111 end
0112 
0113 % strip missing data
0114 if ~any(flags == 'l')
0115   % Mask out columns with NaNs
0116   msk = ~any(isnan(Y),1);
0117   if ~all(msk)
0118     Y = Y(:, msk);
0119     multv = multv(msk);
0120   end
0121 end
0122 return