% Modified version of Ferath's script to threshhold and print SPM5 second
% level results.
% Prints rendered activations as well as a list of coordinates.
% The first 11 lines are the ones you need to (or are most likely to
% want to) change.
% djm 02/07

dataroot = ''; % path to data, and where the results are printed 
outputdir=''; % name of folder, in above directory, where analysis is
cond_dir={''}; % cell array of folder names in outputdir, containing the SPMs for each contrast
outputfilesuffix='_RFX_Rendered.ps'  
%{set the following variable to 'yes' to do a 2-tailed test and show
% rendered activations and deactivations (assumes that positive and negative 
% second level contrasts have already been calculated; Glass brains will 
% only show positive activations). Set it to 'no' to to a 1-tailed test 
% and only calculate positive activations. %}
ShowDeactivations = 'no'; 
Threshold=0.05; % p threshold
ExtentThreshold=0; % minimum number of contiguous voxels 
addpath /imaging/dm01/MoreTools/spm2Batch %{the program needs files cls_getRes.m and cls_getSPM2.m. For CBU people this line will find them, otherwise you need to download them and point to their path.%}
MCC='FDR'; % method for correcting for multiple comparisons : 'FWE','FDR', or 'none'
Brain='/imaging/local/spm/spm5/rend/render_smooth_average.mat'%render_single_subj.mat' %render_smooth_average.mat'; % the brain on which you want to render the results
Style=1; % Rendering Style: NaN=old, 1=normal, <1=brighter

%%%%%%% Do it:
if exist(fullfile(dataroot,strcat(outputdir,outputfilesuffix)))==2
    delete(fullfile(dataroot,strcat(outputdir,outputfilesuffix)));
end;
spm_defaults;
global defaults
defaults.modality='FMRI';
try q=defaults.units{3}; catch; defaults.units={'mm','mm','mm'}; end; % for some reason I had to add this line 15/01/07
    
%--------------------------------------------------------------------
%- xSPM is a structure for the parameters

xSPM = struct( ...
    'swd', '', ... % full path to SPM.mat file
    'Ic', [], ... % no of contrast (or contrasts for conjunction)
    'Im', [],... % no of contrast to mask with. Empty for no masking
    'pm', [],... % masking contrast uncorrected p
    'Ex', [],... % whether masking is inclusive or exclusive
    'title', '',... % if empty results in default contrast title
    'Mcp', MCC,... % Mutiple comp method: FWE|FDR|none
    'u', Threshold,... % threshold (corrected or uncorrected, as above)
    'k', ExtentThreshold); % extent threshold

clear SPM

for i = 1:size(cond_dir,2)
    condir = fullfile(dataroot,outputdir,cond_dir{i});
    cd(condir);
    xSPM.spmmat = fullfile(condir,'SPM.mat'); % Get SPM path
    if lower(ShowDeactivations(1))=='y'
        xSPM.u=Threshold/2; % 2-tailed
        TailText=strcat('2-tailed, p<',num2str(Threshold));
    else
        xSPM.u=Threshold; % 1-tailed
        TailText=strcat('1-tailed, p<',num2str(Threshold));
    end
    xSPM.k=ExtentThreshold;
    bck_xSPM=xSPM;
    delete(gcf) % this is a bit annoying, but for some reason I had problems without it 
    nc=0;
    for i=2:-1:1 
        xSPM=bck_xSPM; % Reset xSPM to default
        xSPM.Ic=i; %Set the current contrast Index
        try
            [hReg,xSPM,SPM] = csl_getRes(xSPM);
            nc=nc+1;
            dat(nc) = struct( 'XYZ', xSPM.XYZ,...
                't', xSPM.Z',...
                'mat', xSPM.M,...
                'dim', xSPM.DIM);
        catch
            ShowDeactivations = 'no'; % only found one contrast 
            xSPM.u=Threshold; % 1-tailed
            TailText=strcat('1-tailed, p<',num2str(Threshold));
        end
    end
    clear flipdat
    flipdat(1)=dat(2);
    if lower(ShowDeactivations(1))=='y'
        flipdat(2)=dat(1);
    end 
    try
        spm_render(flipdat,Style,Brain); 
        ann1=annotation('textbox',[0 .94 1 .06],'Color','r','String',{strcat('RFX results (',TailText,' correction=',MCC,', red=positive)'),strcat('Data:...',condir),strcat('Date:...',datestr(now,0))},'edge','none');
    catch
        try
            spm_render(flipdat(1),Style,Brain);
            ann1=annotation('textbox',[0 .94 1 .06],'Color','r','String',{strcat('RFX results (',TailText,' correction=',MCC,', red=positive)'),strcat('Data:...',condir),strcat('Date:...',datestr(now,0))},'edge','none');
        catch
            try
                spm_render(flipdat(2),Style,Brain);
                ann1=annotation('textbox',[0 .94 1 .06],'Color','r','String',{strcat('RFX results (',TailText,' correction=',MCC,', RED=NEGATIVE)'),strcat('Data:...',condir),strcat('Date:...',datestr(now,0))},'edge','none');
            catch
                ann1=annotation('textbox',[0 .94 1 .06],'Color','r','String',{strcat('RFX results (',TailText,' correction=',MCC,') No significant effects'),strcat('Data:...',condir),strcat('Date:...',datestr(now,0))},'edge','none');
            end
        end
    end
    spm_print(fullfile(dataroot,strcat(outputdir,outputfilesuffix)));
    TabDat=spm_list('List',xSPM,hReg);
    tot=0;
    for r=1:size(TabDat.dat,1)
        if ~isempty(TabDat.dat{r,4}); tot=tot+TabDat.dat{r,4};end
    end
    if tot>0
        delete(ann1);
        ann1=annotation('textbox',[0 .94 1 .06],'Color','r','String',{strcat('RFX results (',TailText,' correction=',MCC,')'),strcat('Data:...',condir),strcat('Date:...',datestr(now,0))},'edge','none');
        ann2=annotation('textbox',[0 0 1 .02],'Color','r','String',strcat('First page of table only. Total number of significant voxels=',num2str(tot)),'edge','none');
        spm_print(fullfile(dataroot,strcat(outputdir,outputfilesuffix)));
        delete(ann2);
    end
    xSPM= bck_xSPM; % Reset xSPM to default
    delete(ann1);
end
fprintf('\nFinished.\nResults saved to: %s.\r',fullfile(dataroot,strcat(outputdir,outputfilesuffix)));