0001 function [fid,PRT,tocomp] = prt_init_fs(PRT, in, mids,mask,precmask,headers)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058
0059 fs_exists = false;
0060 if ~(prt_checkAlphaNumUnder(in.fs_name))
0061 beep
0062 disp('Feature set name should be entered in alphanumeric format only')
0063 disp('Please correct')
0064 return
0065 end
0066 if isfield(PRT,'fs')
0067 if any(strcmpi(in.fs_name,{PRT.fs(:).fs_name}))
0068 fid = find(strcmpi(in.fs_name,{PRT.fs(:).fs_name}));
0069 fs_exists = true;
0070 else
0071 fid = length(PRT.fs)+1;
0072 end
0073 else
0074 fid = 1;
0075 end
0076
0077
0078 if nargout == 1
0079 if fs_exists
0080
0081 disp(['Feature set ''',in.fs_name,''' found in PRT.mat.']);
0082 else
0083 error('prt_init_fs:fsNotFoundinPRT',...
0084 ['Feature set ''',in.fs_name,''' not found in PRT.mat.']);
0085 end
0086 else
0087
0088
0089 [pathName fileName]=fileparts(in.fname);
0090 if fs_exists
0091 warning('prt_init_fs:overwriteFsInPRT',...
0092 ['Feature set ''',in.fs_name,''' found in PRT.mat. Overwriting ...']);
0093 else
0094
0095 disp(['Feature set ''',in.fs_name,''' not found in PRT.mat. Creating...'])
0096 end
0097
0098 PRT.fs(fid).fs_name = in.fs_name;
0099 PRT.fs(fid).k_file = in.fs_name;
0100 PRT.fs(fid).id_col_names = {'group','subject','modality','condition','block','scan'};
0101 PRT.fs(fid).fas=struct('im',[],'ifa',[]);
0102 n_vox=0;
0103 n_mods=length(mids);
0104 for m = 1:n_mods
0105 PRT.fs(fid).modality(m).mod_name = in.mod(mids(m)).mod_name;
0106 PRT.fs(fid).modality(m).detrend = in.mod(mids(m)).detrend;
0107 PRT.fs(fid).modality(m).param_dt = in.mod(mids(m)).param_dt;
0108 PRT.fs(fid).modality(m).mode = in.mod(mids(m)).mode;
0109
0110 vm = spm_vol(mask{m});
0111 vm = spm_read_vols(vm);
0112 PRT.fs(fid).modality(m).feat_idx_img = find(vm>0);
0113 if isempty(find(vm>0))
0114 error('prt_prepare_data:NoVoxelinMask',...
0115 ['Mask of modality ',num2str(m),' does not contain any voxel >0'])
0116 end
0117 mid=mids(m);
0118 if m==1
0119 n_vox=length(find(vm>0));
0120 end
0121 if n_vox~=length(find(vm>0))
0122 error('prt_prepare_data:MasksNotConsistent',...
0123 'Masks access areas of different sizes across modalities')
0124 end
0125
0126 if ~isempty(precmask{m})
0127 vm = spm_vol(precmask{m});
0128 vm = spm_read_vols(vm);
0129 if isempty(find(vm>0))
0130 error('prt_prepare_data:NoVoxelinMask',...
0131 ['2nd level mask of modality ',num2str(m),' does not contain any voxel >0'])
0132 end
0133 [d,PRT.fs(fid).modality(m).idfeat_fas] = intersect(PRT.fs(fid).modality(m).feat_idx_img, find(vm>0));
0134 else
0135 PRT.fs(fid).modality(m).idfeat_fas=[];
0136 end
0137 PRT.fs(fid).modality(m).normalise=struct('type',[],'scaling',[]);
0138 end
0139
0140 indm = zeros(n_mods,1);
0141 szm = zeros(n_mods,1);
0142
0143
0144
0145 n = 0;
0146 for gid = 1:length(PRT.group)
0147 for sid = 1:length(PRT.group(gid).subject);
0148 for m = 1:n_mods
0149 mid = mids(m);
0150 if strcmpi(in.mod(mid).mode,'all_scans');
0151 n = n + size(PRT.group(gid).subject(sid).modality(mid).scans,1);
0152 elseif strcmpi(in.mod(mid).mode,'all_cond')
0153 if ~isfield(PRT.group(gid).subject(sid).modality(mid).design,'conds')
0154 error('prt_model:fsIsAllCondModelisAllScans',...
0155 ['''All conditions'' selected for modality ', num2str(m)...
0156 ' but no design was specified. This syntax is invalid, '...
0157 'Please use ''All Scans'' instead.']);
0158 end
0159 for cid = 1:length(PRT.group(gid).subject(sid).modality(mid).design.conds)
0160 n = n + length(PRT.group(gid).subject(sid).modality(mid).design.conds(cid).scans);
0161 end
0162 end
0163 end
0164 end
0165 end
0166 PRT.fs(fid).id_mat = zeros(n,length(PRT.fs(fid).id_col_names));
0167 PRT.fs(fid).fas.im = zeros(n,1);
0168 PRT.fs(fid).fas.ifa= zeros(n,1);
0169
0170
0171
0172 sample_range = 0;
0173 for gid = 1:length(PRT.group)
0174 for sid = 1:length(PRT.group(gid).subject);
0175 for m = 1:n_mods
0176 mid = mids(m);
0177
0178 if strcmpi(in.mod(mid).mode,'all_scans')
0179 n_vols_s = size(PRT.group(gid).subject(sid).modality(mid).scans,1);
0180 all_scans = 1:n_vols_s;
0181
0182
0183 sample_range = (1:n_vols_s)+max(sample_range);
0184 PRT.fs(fid).id_mat(sample_range,1) = gid;
0185 PRT.fs(fid).id_mat(sample_range,2) = sid;
0186 PRT.fs(fid).id_mat(sample_range,3) = mid;
0187
0188 if isfield(PRT.group(gid).subject(sid).modality(mid).design,'conds')
0189 conds = PRT.group(gid).subject(sid).modality(mid).design.conds;
0190 for cid = 1:length(conds)
0191 scans = PRT.group(gid).subject(sid).modality(mid).design.conds(cid).scans;
0192 blocks = PRT.group(gid).subject(sid).modality(mid).design.conds(cid).blocks;
0193
0194 PRT.fs(fid).id_mat(sample_range(scans),4) = cid;
0195 PRT.fs(fid).id_mat(sample_range(scans),5) = blocks;
0196
0197 end
0198
0199 PRT.fs(fid).id_mat(sample_range,6) = 1:length(all_scans);
0200 else
0201 scans = 1:size(PRT.group(gid).subject(sid).modality(mid).scans,1);
0202 PRT.fs(fid).id_mat(sample_range,6) = scans;
0203 end
0204
0205 sctoadd=(1:n_vols_s)+indm(m);
0206 PRT.fs(fid).fas.ifa(sample_range)=sctoadd';
0207 PRT.fs(fid).fas.im(sample_range)=mid*ones(n_vols_s,1);
0208
0209 indm(m)=n_vols_s+max(indm(m));
0210 elseif strcmpi(in.mod(mid).mode,'all_cond')
0211 conds = PRT.group(gid).subject(sid).modality(mid).design.conds;
0212 n_vols_s = size(PRT.group(gid).subject(sid).modality(mid).scans,1);
0213
0214
0215 for cid = 2:4
0216 scans = PRT.group(gid).subject(sid).modality(mid).design.conds(cid).scans;
0217 blocks = PRT.group(gid).subject(sid).modality(mid).design.conds(cid).blocks;
0218 n_vol_s_c = length(scans);
0219 if n_vol_s_c==0
0220 sample_range = 1+max(sample_range);
0221 PRT.fs(fid).id_mat(sample_range,5) = 0;
0222 PRT.fs(fid).id_mat(sample_range,6) = 0;
0223 else
0224 sample_range = (1:n_vol_s_c)+max(sample_range);
0225 PRT.fs(fid).id_mat(sample_range,5) = blocks;
0226 PRT.fs(fid).id_mat(sample_range,6) = scans;
0227
0228 sctoadd=scans+indm(m);
0229 PRT.fs(fid).fas.ifa(sample_range)=sctoadd';
0230 PRT.fs(fid).fas.im(sample_range)=mid*ones(n_vol_s_c,1);
0231 end
0232
0233
0234
0235 PRT.fs(fid).id_mat(sample_range,1) = gid;
0236 PRT.fs(fid).id_mat(sample_range,2) = sid;
0237 PRT.fs(fid).id_mat(sample_range,3) = mid;
0238 PRT.fs(fid).id_mat(sample_range,4) = cid;
0239
0240
0241
0242
0243
0244 end
0245
0246 indm(m)=n_vols_s+max(indm(m));
0247 end
0248 szm(m)=szm(m)+size(PRT.group(gid).subject(sid).modality(mid).scans,1);
0249 end
0250 end
0251 end
0252
0253
0254
0255 if ~isfield(PRT,'fas');
0256
0257 for m = 1:length(PRT.masks)
0258
0259 PRT.fas(m)=struct('mod_name',[],'dat',[],'detrend',[],'param_dt',[],'hdr',[]);
0260 PRT.fas(m).mod_name = PRT.masks(m).mod_name;
0261 end
0262 end
0263 tocomp=zeros(1,length(in.mod));
0264 prt_dir=fileparts(in.fname);
0265 for i=1:n_mods
0266
0267 if mids(i)>length(PRT.fas) ||...
0268 isempty(PRT.fas(mids(i)).dat) || exist(PRT.fas(mids(i)).dat.fname)==0 ||...
0269 PRT.fas(mids(i)).detrend ~= in.mod(mids(i)).detrend || ...
0270 (isempty(PRT.fas(mids(i)).param_dt) && ~isempty(in.mod(mids(i)).param_dt)) || ...
0271 (~isempty(PRT.fas(mids(i)).param_dt) && isempty(in.mod(mids(i)).param_dt)) || ...
0272 ((~isempty(PRT.fas(mids(i)).param_dt) && ~isempty(in.mod(mids(i)).param_dt)) && ...
0273 PRT.fas(mids(i)).param_dt~=in.mod(mids(i)).param_dt)
0274
0275 if mids(i)>length(PRT.fas) || isempty(PRT.fas(mids(i)).dat)
0276 disp(['File array does not exist for modality ''',...
0277 char(in.mod(mids(i)).mod_name),'''. Creating...'])
0278 elseif PRT.fas(mids(i)).detrend ~= in.mod(mids(i)).detrend ...
0279 && any(strcmpi(fieldnames(PRT.fas(mids(i)).dat),'fname')) ...
0280 && exist(PRT.fas(mids(i)).dat.fname,'file')
0281
0282 warning('prt_init_fs:overwriteFileArray',...
0283 ['File array already exists for modality ''',...
0284 char(in.mod(mids(i)).mod_name),''', but parameters ',...
0285 'have changed. Re-creating ...']);
0286
0287 delete(PRT.fas(mids(i)).dat.fname);
0288 end
0289
0290 tocomp(mids(i))=1;
0291
0292 PRT.fas(mids(i)).detrend = in.mod(mids(i)).detrend;
0293 PRT.fas(mids(i)).param_dt = in.mod(mids(i)).param_dt;
0294 PRT.fas(mids(i)).hdr = headers{i};
0295 PRT.fas(mids(i)).idfeat_img = PRT.fs(fid).modality(i).feat_idx_img;
0296 datname=[prt_dir,filesep,'Feature_set_',char(in.mod(mids(i)).mod_name),'.dat'];
0297 PRT.fas(mids(i)).dat = file_array(...
0298 datname, ...
0299 [szm(i),n_vox],...
0300 spm_type('float64'), ...
0301 0, ...
0302 1);
0303 else
0304 disp(['Using existing file array for modality ''', ...
0305 char(in.mod(mids(i)).mod_name),'''.'])
0306 end
0307
0308
0309 if in.mod(mids(i)).normalise==2
0310 try
0311 load(in.mod(mids(i)).matnorm);
0312 catch
0313 error('prt_prepare_data:ScalingMatUnloadable',...
0314 'Could not load the .mat file containing the scaling')
0315 end
0316 try
0317 szin=max(size(scaling));
0318 catch
0319 error('prt_prepare_data:ScalingNotinFile',...
0320 'This file does not contain the "scaling" field required')
0321 end
0322 if szin~=szm(i)
0323 error('prt_prepare_data:Scalingdimensionwrong',...
0324 'The dimension of the .mat file does not correspond to the number of scans in that modality')
0325 end
0326 PRT.fs(fid).modality(i).normalise.type=2;
0327 PRT.fs(fid).modality(i).normalise.scaling=reshape(scaling,1,szm(i));
0328 elseif in.mod(mids(i)).normalise==1
0329 PRT.fs(fid).modality(i).normalise.type=1;
0330 else
0331 PRT.fs(fid).modality(i).normalise.type=0;
0332 end
0333 end
0334
0335 PRT.fs(fid).modality=rmfield(PRT.fs(fid).modality,'feat_idx_img');
0336 end