Homec4science
Diffusion iCAPs public version 2 (master) Browse

iCAPs_v2/packages/spm12/spm_bms_partition.m
e5f5005bb2efmaster

Last Change
Raw File

spm_bms_partition.m
View Options

function spm_bms_partition(BMS)
% Compute model partitioning for BMS
% FORMAT spm_bms_partition(BMS)
%
% Input:
% BMS structure (BMS.mat)
%
% Output:
% PPM (images) for each of the subsets defined
% xppm_subsetn.img (RFX) and ppm_subsetn.img (FFX)
%__________________________________________________________________________
% Copyright (C) 2009-2011 Wellcome Trust Centre for Neuroimaging
% Maria Joao Rosa
% $Id: spm_bms_partition.m 4492 2011-09-16 12:11:09Z guillaume $
% Contrast vector
% -------------------------------------------------------------------------
spm_input('Specify contrast vector. Example: [1 1 2 2 3 3]',1,'d');
contrast = spm_input('Contrast vector',2,'e',[]);
% Inference method to plot
% -------------------------------------------------------------------------
method = spm_input('Inference method',3,'b','FFX|RFX',['FFX';'RFX']);
nb_subsets = length(unique(contrast));
max_cont = max(contrast);
nb_models = length(contrast);
switch method
case 'FFX'
str_method = 'ffx';
str_output = 'ppm';
case 'RFX'
str_method = 'rfx';
str_output = 'xppm';
otherwise
error('Unknown inference method.');
end
% Check if ffx exists
% -------------------------------------------------------------------------
if ~isfield(BMS.map,str_method)
msgbox(sprintf('No %s analysis in current BMS.mat.',method));
return
end
% Check number of subsets and nb of models
% -------------------------------------------------------------------------
bms_fields = eval(sprintf('BMS.map.%s.ppm',str_method));
nmodels = size(bms_fields,2);
if nb_models ~= nmodels || nb_subsets == 1 || max_cont ~= nb_subsets
msgbox('Invalid contrast vector!')
return
end
% Get data for each subset
% -------------------------------------------------------------------------
data = cell(1,nb_subsets);
for i = 1:nmodels,
num = contrast(i);
data{num} = [data{num};bms_fields{i}];
end
% Create new images by summing old the ppms
% -------------------------------------------------------------------------
pth = fileparts(BMS.fname);
data_vol = cell(nb_subsets,1);
ftmp = cell(nb_subsets,1);
for j = 1:nb_subsets,
data_vol{j} = spm_vol(char(data{j}));
n_models_sub = size(data{j},1);
ftmp{j} = 'i1';
for jj = 1:n_models_sub-1
ftmp{j} = [ftmp{j},sprintf(' + i%d',jj+1)];
end
fname = fullfile(pth,[sprintf('subset%d_%s',j,str_output) spm_file_ext]);
save_fn{j} = fname;
Vo = calc_im(j,data_vol,fname,ftmp);
end
% Save new BMS structure
% -------------------------------------------------------------------------
bms_struct = eval(sprintf('BMS.map.%s',str_method));
bms_struct.subsets = save_fn;
switch method
case 'FFX'
BMS.map.ffx = bms_struct;
case 'RFX'
BMS.map.rfx = bms_struct;
end
file_name = BMS.fname;
BMS.xSPM = [];
save(file_name,'BMS', spm_get_defaults('mat.format'))
% Return to results
%==========================================================================
spm_input('Done',1,'d');
return;
%==========================================================================
% out = calc_im(j,data_vol,fname,ftmp)
%==========================================================================
% Function to sum the data (taken from spm_imcalc)
%--------------------------------------------------------------------------
function out = calc_im(j,data_vol,fname,ftmp)
Vi_tmp = data_vol{j};
Vi = Vi_tmp(1);
Vo(j) = struct(...
'fname', fname,...
'dim', Vi.dim,...
'dt', [spm_type('float32') spm_platform('bigend')],...
'mat', Vi.mat,...
'descrip', 'spm - algebra');
hold = 1; mask = 0; dmtx = 0;
Vi = data_vol{j};
n = numel(Vi);
Y = zeros(Vo(j).dim(1:3));
f = ftmp{j};
for p = 1:Vo(j).dim(3),
B = spm_matrix([0 0 -p 0 0 0 1 1 1]);
if dmtx, X=zeros(n,prod(Vo(j).dim(1:2))); end
for i = 1:n
M = inv(B*inv(Vo(j).mat)*Vi(i).mat);
d = spm_slice_vol(Vi(i),M,Vo(j).dim(1:2),[hold,NaN]);
if (mask<0), d(isnan(d))=0; end;
if (mask>0) && ~spm_type(Vi(i).dt(1),'nanrep'), d(d==0)=NaN; end
if dmtx, X(i,:) = d(:)'; else eval(['i',num2str(i),'=d;']); end
end
try
eval(['Yp = ' f ';']);
catch
error(['Can''t evaluate "',f,'".']);
end
if prod(Vo(j).dim(1:2)) ~= numel(Yp),
error(['"',f,'" produced incompatible image.']); end
if (mask<0), Yp(isnan(Yp))=0; end
Y(:,:,p) = reshape(Yp,Vo(j).dim(1:2));
end
temp = Vo(j);
temp = spm_write_vol(temp,Y);
out(j) = temp;
c4science · Help