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