function batch_discrim_gen_all
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Information about this function here
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% User modified parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
display_all = 0;
display_final = 1;

param.ABS_or_Tvalue = 'ABS';
param.n_voxels = 50;
param.save_jpg = 0;
param.filestring = ['lin_discrim' '_' param.ABS_or_Tvalue '_' int2str(param.n_voxels)];
param.filestring_results_total = ['prediction_results'];
param.basedir= ['D:\V1load_motion\'];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

param.subjects = {'Williamson_Kirk', 'JALKANEN_Lauri', 'Bays_Paul', 'Schulvinik_Marieke'};
param.ROIs = {'V1', 'V2', 'V3', 'V123'};
param.localizers = {'in', 'out'};
param.conditions = {'low', 'high', 'left', 'right'};

ss_values = 1:length(param.subjects);
g_values = 1:length(param.ROIs);     % ROIS
loc_values = 1:length(param.localizers);
cond_values = 1:length(param.conditions);

result_total = cell(length(param.subjects),length(param.localizers));

for subjtodo= ss_values
    for loc = loc_values
        res = cell(length(param.ROIs),length(param.conditions));
        for cond = cond_values
            %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
            if strcmp(param.conditions(cond),'low')
                cond1 = 1; cond2 = 2;
            elseif strcmp(param.conditions(cond),'high')
                cond1 = 3; cond2 = 4;
            elseif strcmp(param.conditions(cond),'left')
                cond1 = 1; cond2 = 3;
            elseif strcmp(param.conditions(cond),'right')
                cond1 = 2; cond2 = 4;
            end
            if display_all == 1      figure,   end

            for g = g_values

                region_todo = param.ROIs{g};
                subjectname = param.subjects{subjtodo}
                basedir2= [param.basedir subjectname '\'];
                [filsubjectnumber filehdr dir_data localizer_100plus] = filsubjectnumbersMOTION(subjectname,param.basedir,region_todo);
                n_runs = size(filehdr{1},1); % in the experiment. Currently the program is set up for different classification of different runs.

                dir_data{1} = {[basedir2 'prediction\' region_todo '\vectors_' param.localizers{loc} '_']};

                m = zeros(length(dir_data),param.n_voxels);
                s = zeros(length(dir_data),param.n_voxels);
                for grp=1:n_runs
                    fn=strcat(dir_data{1},int2str(grp));
                    load(fn{1});
                    fprintf('Subject %d loop %d\n',1, grp);
                    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                    % get indexes for ordering
                    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

                    %------------------------------------------------------
                    %  this section to correct for diff sizes of train
                    %  vectors
                    for i = 1:4
                        ntrain(i)=size(vectors_train{i},2);
                    end
                    compltrain= min(ntrain);
                    for i = 1:4
                        vectors_train{i}=[vectors_train{i}(:,1:compltrain)];
                    end
                    %  this section to correct for diff sizes of test
                    %  vectors
                    for i = 1:4
                        ntest(i)=size(vectors_test{i},2);
                    end
                    compltest= min(ntest);
                    for i = 1:4
                        vectors_test{i}=[vectors_test{i}(:,1:compltest)];
                    end
                    %--------------------------------------------------

                    %  the voxel-selection should be equally based on both conditions therefore
                    %  we create a temporary vectors_all for the selection
                    % however never use test data to select voxels!
                    vectors_all = {[vectors_train{1}' vectors_train{3}']'; [vectors_train{2}' vectors_train{4}']'}';

                    ntrain1=size(vectors_all{1},2);
                    ntrain2=size(vectors_all{2},2);
                    sz1=ntrain1;  % gives number of voxels
                    sz2=ntrain2;  % gives number of voxels
                    compl= min(ntrain1,ntrain2); % in case one has less voxels
                    dat1=[(vectors_all{1}(:,1:compl))];  % load training data for condition 1
                    dat2=[(vectors_all{2}(:,1:compl))];  % load training data for conditon 2
                    m1=mean(dat1); %mean of condition 1 for each voxel
                    m2=mean(dat2);%mean on condition 2 fo each voxel

                    % get indexes for sorting by combined mean activation
                    ampl=(abs(m1)+abs(m2));
                    [Y1,I1]=sortrows(ampl');
                    indexesM=flipud(I1); % gets the order of voxels in order of mean
                    indexes=indexesM(1:param.n_voxels); % ordering for mean
                    clear resultsvector vectors_hilo
                    testindexes=[];

                    for i=1:param.n_voxels;

                        testindexes=[testindexes indexes(i)];

                        ntrain=size(vectors_train{cond1},1);
                        nsample=size(vectors_test{cond1},1);
                        groups=[ones(ntrain,1);ones(ntrain,1)+1];  %column vector of groups 1 and 2
                        trainingtot=[(vectors_train{cond1});(vectors_train{cond2})];  %combine both train data for con1 and con2
                        sampletot=[(vectors_test{cond1});(vectors_test{cond2})];    %combine both test data for con1 and con2


                        training=trainingtot(:,testindexes);
                        sample=sampletot(:,testindexes);
                        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                        %  HERE IS THE ACTUAL CLASSIFICATION
                        [c,err]=classify(sample,training,groups,'linear');
                        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                        classres=c';
                        classerr=err;
                        pcorrect=sum(([ones(nsample,1);(ones(nsample,1)+1)]')==classres)/(2*nsample);
                        resultsvector(i)=pcorrect;

                    end
                    total_error(grp,:)=resultsvector;

                    clear vectors_test vectors_train vectors_all
                end %group
                m(1, :)  = mean(total_error);
                s(1, :)  = std(total_error)/sqrt(n_runs-1);
                if display_all == 1
                    subplot(2,2,g);errorbar(1:param.n_voxels,m(1,:),s(1,:),s(1,:),'b');whitebg
                    title([num2str(subjtodo) ' ' region_todo ' ' param.localizers{loc} ' ' param.conditions(cond)],'FontSize',10);
                    set(gca, 'Ylim', [0.45,0.75]);
                end

                res{g, cond} = {m s};

            end  %ROI
        end % cond
        save ([ basedir2 'prediction\' param.filestring '_' param.localizers{loc}], 'res');

        result_total{subjtodo,loc} = res;

    end %loc

end  % subj

% save all results into one big file
save ([param.basedir '_prediction\' param.filestring_results_total], 'result_total');


if display_final ==1
    make_stats_lin_discrimMOTION(param, ss_values)
end