% Spatial Preprocessing for fMRI data
% Christian Kaul, July 2007
% presented @ Matlab for Cogntive Neuroscience 12/07/07
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Specify each subject and the number of runs acquired for them.
% this script assumes a folder structure like this:
%
% exp_name/
% exp_name/s1
% exp_name/s1/structural           % the structural folder contains
                                   % structural, fieldmaps, segmentation,...
% exp_name/s1/structural/ROI       %the ROI folder contains all the area defs, ergo put your amygdala def here
%
% exp_name/s1/localizer_data
% exp_name/s1/localizer_data/run1  % contains functional files
% exp_name/s1/localizer_data/run2
% exp_name/s1/localizer_stats      % contains the SPM, betas, con-images
%
% exp_name/s1/experiment
% exp_name/s1/experiment/run1      % contains functional files
% exp_name/s1/experiment/run2
% exp_name/s1/experiment_stats     % contains the SPM, betas, con-images


% this batch script with set perform the following chain of pre-processing
% step:
% 1) Realignment of all functional runs  (realign: estimate & reslice)
% 2) Coregister all epi images to the space of the structural.
% other sections can be added as desired
%


% Specify root directory
origDir  = 'D:\V1load_motion2';
% which group of data: experiment, meridian or motion localizer
dataFolder = 'experiment\run';
structuralFolder = 'structural';
% subjects and the number of volumes collected for each run and the name of
% the structural file
subjects = { ...
    %  '1',  ones(1,48)*48-3, 'ED-MA04903.img'; ...
    %  '2',  ones(1,48)*48-3, 'GV-MA04904.img'; ...
     '3',  ones(1,48)*48-3, 'HR-sMA04830.img'; ...
    %  '4',  ones(1,48)*48-3, 'DB-meansMA04771.img'; ...
    %  '5',  ones(1,48)*48-3, 'AD.img'; ...
%              '7',  [205 205]-5,         'ADs1.img'; ...
%              '8',  [205 205 205]-5,     'AMs1.img'; ...
%              '9',  [205 205]-5,         'meanGRs1.img'; ...
%              '10', [205 205]-5,         'VWs1.img'; ...
             };
% =========================================================================
cd(origDir);
spm fmri
% =========================================================================
% =========================================================================
% =========================================================================
% =========================== Realignment =================================
% =========================================================================
clc;
disp(['Realign: Estimate and Reslice for MEAN and ALL images ']);
for s0 = 1 : size(subjects,1)
    disp(['Realign job specification for Subject : ', subjects{s0,1}]);
    runs = size(subjects{s0,2},2);
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % first specify the f* scans for each run
    for run = 1 : runs
        %% define epi files in the run
        epiDir = [origDir '\s' subjects{s0,1} '\' dataFolder int2str(run)];
        %%% select scans and assign to job
        f   = spm_select('List', epiDir, '^fM.*\.img$');
        fs  = cellstr([repmat([epiDir '\'],size(f,1),1) f]);
        jobs{1}.spatial{1}.realign{1}.estwrite.data{run} = fs;
        %%%% clear temporary variables for next run
        f = []; fs = [];
    end
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % specify Estimation Options
    % all default values taken from SPM5 batch
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.quality = 0.9;
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.sep     = 4.0;
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.fwhm    = 5.0;
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.rtm     = 1.0;
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.interp  = 2.0;
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.wrap    = [0 0 0];
    jobs{1}.spatial{1}.realign{1}.estwrite.eoptions.weight  = {};
    % specify Reslice Options
    jobs{1}.spatial{1}.realign{1}.estwrite.roptions.which   = [2 1]; % mean and all
    jobs{1}.spatial{1}.realign{1}.estwrite.roptions.interp  = 4.0;
    jobs{1}.spatial{1}.realign{1}.estwrite.roptions.wrap    = [0 0 0];
    jobs{1}.spatial{1}.realign{1}.estwrite.roptions.mask    = 1.0;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % save and run job
    save realign.mat jobs
    disp(['RUNNING Realign for subject ' subjects{s0,1}]);
    spm_jobman('run','realign.mat');
    clear jobs
end
% =========================================================================
% ========================== Coregistration ===============================
% =========================================================================
clc;
disp(['Coregister: Estimate (without Reslice)']);
for s0 = 1 : size(subjects,1)
    disp(['Coregister job specification for Subject : ', subjects{s0,1}]);
    runs = size(subjects{s0,2},2);
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % first specify the refernce image: structural scan
    stDir = [origDir '\s' subjects{s0,1} '\' structuralFolder '\'];
    stFilePath = [stDir subjects{s0,3} ',1'];
    jobs{1}.spatial{1}.coreg{1}.estimate.ref{1} = stFilePath;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % then specify the source image: mean scan
    sourceDir = [origDir '\s' subjects{s0,1} '\' dataFolder '1'];
    f  = spm_select('List', sourceDir, '^mean.*\.img$');
    fs = [sourceDir '\', f ',1'];
    jobs{1}.spatial{1}.coreg{1}.estimate.source{1} = fs;
    f = []; fs = [];
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % then specify the OTHER FILES: rf* scans for
    conCat_fs = [];
    for run = 1 : runs
        %% define epi files in the run
        epiDir = [origDir '\s' subjects{s0,1} '\' dataFolder int2str(run)];
        %%% select scans and assign to job
        f   = spm_select('List', epiDir, '^rfM.*\.img$');
        fs  = cellstr([repmat([epiDir '\'],size(f,1),1) f repmat([',1'],size(f,1),1)]);
        conCat_fs = [conCat_fs; fs];
        %%%% clear temporary variables for next run
        f = []; fs = [];
    end
    jobs{1}.spatial{1}.coreg{1}.estimate.other = conCat_fs;
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % then add estimation options
    jobs{1}.spatial{1}.coreg{1}.estimate.eoptions.cost_fun = 'nmi';
    jobs{1}.spatial{1}.coreg{1}.estimate.eoptions.sep      = [4 2];
    jobs{1}.spatial{1}.coreg{1}.estimate.eoptions.tol      = [ 0.0200 0.0200 0.0200 0.0010 0.0010 0.0010 ...
                                                                0.0100 0.0100 0.0100 0.0010 0.0010 0.0010];
    jobs{1}.spatial{1}.coreg{1}.estimate.eoptions.fwhm     = [7 7];
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % save and run job
    save coreg.mat jobs
    disp(['RUNNING Coreg for subject ' subjects{s0,1}]);
    spm_jobman('run','coreg.mat');
    clear jobs
end