clear all

Subject_ID =  { '01', '02', '05', '06', '08' '10', '11' , '12', '13', '14', '15', '16', '17', '18', '20' , '22', '23', '24', '26', '27', '28', '29', '30'};


path = '/home/marius/Dokumente/Project_2/Matlab_output/';
version = 'version9';


load([path version '/behavioural_data_' version '_' num2str(length(Subject_ID)) '_Subjects.mat']);



%% fMRI GLM regressors:

version = 'MB6';
n_trials = 160;
rating_distance = 10; 

%% Accumulation phase

%Preassign regressor matrices with 1 column per participant

% % Accumulation phase
all_screen = nan(700, length(Subject_ID));
all_signal_update = nan(size(all_screen,1), length(Subject_ID));
all_screen_start = nan(n_trials, length(Subject_ID));
all_screen_end =  nan(n_trials, length(Subject_ID));
all_prior_sumv = nan(size(all_screen));
all_sum_vec = nan(size(all_screen));

all_signal_update_duration = nan(size(all_screen));
all_number_sample_onset = nan(size(all_screen));

all_tip_scale_sign=nan(size(all_screen));
all_number_sample = nan(size(all_screen));
all_sample_length =  nan(n_trials, length(Subject_ID));

% Decision Phase occluded_sample_no_nan = [];


all_shared_ev = nan(n_trials, length(Subject_ID));
all_occluded_ev = nan(n_trials, length(Subject_ID));
all_question_self_onset = nan(n_trials, length(Subject_ID));
all_question_other_onset = nan(n_trials, length(Subject_ID));
all_choice_disc = nan(n_trials, length(Subject_ID));
all_acc_pred_reg = nan(n_trials, length(Subject_ID));

% Onsets

all_sample_onset = nan(700, length(Subject_ID));
all_DecSelf_onset = nan(n_trials, length(Subject_ID));
all_DecOther_onset = nan(n_trials, length(Subject_ID));
all_Outcome_onset = nan(n_trials, length(Subject_ID));



%% Loop through columns (participants)

for g = 1:length(Subject_ID)
    
        
    %% Empty vectors in each iteration (different schedules!)
    
    %% Needed to create regressors
    
    screen = [];
    screen_start = [];
    screen_end = [];
    sum_vec = [];
    prior_sumv  = [];
    signal_update = [];
    number_sample = [];
    start_samples=[];
    tip_scale_sign = [];
    sum_vec_shared = [];
    sum_vec_occluded = [];
    
    % Number of sample in trials
    
    
    sample_length = [];
    number_sample_shared = [];
    number_sample_shared_no_nan = [];
    
    number_sample_occluded = [];
    number_sample_occluded = []; 
    
    
    % For own prior sums
    sum_Sample_signals = [];
    prior_Sample = [];
    sum_priorSample_signals= [];
    sample_update = [];
    
    %For others' prior sums 
      
    prior_Sample_other = [];
    sum_priorSample_signals_other = [];
    
    
    % All prior sums 
    
    prior_sum_other = [];
    prior_sum_self = [];
    prior_sum_diff = [];
    
    % All signal update
    
    signal_update_self = [];
    signal_update_other = [];
    signal_update_diff = [];
    
    %% Phase 1 (Sampling Shared)
    
    % Evidence prior to sample
    
    prior_sum_diff_shared = [];
    prior_sum_other_shared = [];
    prior_sum_self_shared = [];
    
    % Signal update for shared information only (information on current sample)
    
     signal_update_diff_shared = [];
     signal_update_self_shared =[];
     signal_update_other_shared = []; 
    
    
    % Change in outcome at current sample (yes/no)
    
    tip_scale_sign_shared = [];
    tip_scale_sign_shared_no_nan = [];
    
    
    
    %% Phase 2 Decision self
    
    % Evidence shared
    
    shared_ev_self = [];
    shared_ev_self_abs = [];
    shared_ev_self_signed = [];
    
    % Evidence occluded
    
    occluded_self_ev = [];
    occluded_ev_self_abs = [];
    occluded_ev_self_signed = [];
    
    %% Phase 3 Decision other
    
    % Evidence shared
    
    shared_ev_other = [];
    shared_ev_other_abs = [];
    shared_ev_other_signed = [];
    
    % Evidence occluded
    
    occluded_ev_other = [];
    occluded_ev_other_abs = [];
    occluded_ev_other_signed = [];
    
    %% Phase 4 (Outcome)
    
    % Correct answer for self decision
    
    correct_self = [];
    
    % correct answer for other decision
    
    correct_other = [];
    
    % Discrepancy in decision for other and behaviour for other (yes/no)
    
    prediction_disc = [];
    
    % Confidence self (chance to choose red for yourself given the evidence in
    % a trial)
    
    % Confidence self
    
    choices_self_red_single =[];
    evidence_against_blue_self_single = [];
    evidence_for_blue_self_single = [];
    evidence_for_blue_self_single_shift = [];
    bins_self = [];  
    confidence_self = nan(n_trials, 1);
    
    % Confidence other (chance to choose red for the other given the evidence in
    % a trial)
 
    choices_other_red_single =[];
    evidence_against_blue_other_single = [];
    evidence_for_blue_other_single = [];
    evidence_for_blue_other_single_shift = [];
    bins_self = [];  
    confidence_other = nan(n_trials, 1);
    
    % Prediction error other (choice(yes/no) - confidence (chance to choose red)
    
    prediction_error = [];
    
    % Agent ability
    
    agent_ability = [];
    
    % Interaction agent_ability * prediction_error
    
    inter_error_ability = [];
    
    
    
    
    %% Phase 5 (Sampling Occluded)
    
    % Sum of evidence prior to occluded sample
    
    prior_sum_diff_occluded = [];
    prior_sum_other_occluded = [];
    prior_sum_self_occluded = [];
    
    % Signal update of occluded sample (information shown on sample)    

     signal_update_diff_occluded = [];
     signal_update_self_occluded =[];
     signal_update_other_occluded = []; 

    % Change in outcome at current occluded sample (yes/no)
    
    tip_scale_sign_occluded = [];
    tip_scale_sign_occluded_no_nan = [];
    
    
    % Number of occluded sample in trial
    
    occluded_sample = [];
    occluded_sample_no_nan = [];
    
    %Indices of shared/occluded samples in full sample vector: 
    
    sample_index_shared = [];
    sample_index_occluded = [];
    
    
    %% Onset times
    
    Onset = [];
    
    % Phase 1 (onsets of shared samples)
    
    sample_onset = [];
    number_sample_onset_shared = [];
    number_sample_onset_shared_no_nan = [];
    
    % Phase 2 (onset of question self)
    
    question_self_onset = [];
    DecSelf_onset = [];
    
    % Phase 3 (onset of question other)
    question_other_onset = [];
    DecOther_onset = [];
    
    % Phase 4 (onset of feedback other's behaviour)
    
    Outcome_onset = [];
    
    % Phase 5 (onset of occluded sample)
    
    number_sample_onset_occluded = [];
    occluded_sample_onset = [];
    
    
    % Clear Phase Struct every iteration, so no old values are kept in
    % there as well as in the vectors feeding into the Struct
    
    if ~iscell('Phase')
        clear Phase
    end
    
    
    
    %% Create regressors
    
    
    %% Needed to create regressors
    
    % Extract indices of Screens(/Samples)
    
    screen = find(strcmp(behavioural_data.Subject(g).codes, 'Screen'));
    
    % Extract indices of first and last screen in trial
    
    for f = 1:length(behavioural_data.output.indices.welcome(:,g))
        k = min(find(screen > behavioural_data.output.indices.welcome(f,g)));
        o = max(find(screen < behavioural_data.output.indices.question_you(f,g)));
        screen_start(f) = k;
        screen_end(f) = o;
    end
    
    
    Sample_signals_other = behavioural_data.Subject(g).schedule.MainSchedule.SampleSignals;
    
    % Vector of total information until this sample
    % sum_Sample_signals_shared = sum_Sample_signals; 
    
    for o=1:160
        Sample_signals_other(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(o, ~isnan(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(o,:))),o) = nan;
       % prior_Sample_other(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(o, ~isnan(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(o,:))),o) = nan;
    end
    
   prior_Sample_other = [zeros(1,size(behavioural_data.Subject(g).schedule.MainSchedule.SampleSignals,2)); Sample_signals_other(1:end-1,:)];

    % Keep values of prior sums at positions of occluded samples (doesn't
    % change so prior sum stays the same at this point?) 
    
    % sum prior to the current sample for shared information only 
    
    sum_priorSample_signals_other = nancumsum(prior_Sample_other,1);
    
    
    prior_Sample=[zeros(1,size(behavioural_data.Subject(g).schedule.MainSchedule.SampleSignals,2));behavioural_data.Subject(g).schedule.MainSchedule.SampleSignals(1:end-1,:)];
    sum_Sample_signals = cumsum(behavioural_data.Subject(g).schedule.MainSchedule.SampleSignals, 1);

    sum_priorSample_signals=cumsum(prior_Sample, 1);
    sample_update=abs(sum_Sample_signals)-abs(sum_priorSample_signals);
    
    sum_vec(1:length(sum_Sample_signals(~isnan(sum_Sample_signals)))) = sum_Sample_signals(~isnan(sum_Sample_signals));
    
    % Vector of total information prior to the current sample
    
    prior_sumv(1:length(sum_priorSample_signals(~isnan(sum_Sample_signals))))= sum_priorSample_signals(~isnan(sum_Sample_signals));
    
    % Update of signal (information shown on the current sample
    
    
    
    
    % Compute the sample length (3-5) in the current trial in order to extract the number of the sample
    
    for d = 1:length(behavioural_data.output.indices.welcome(:,g))
        sample_length(d) = sum(strcmp({behavioural_data.output.codes{behavioural_data.output.indices.welcome(d, g):behavioural_data.output.indices.question_you(d, g), g}}, 'Screen'));
    end
    
    
    % Number of the sample in the trial (1-5)
    
    for k = 1:length(sample_length)
        for t = 1:sample_length(k)
            if k < 2
                number_sample(t) = t;
            else
                number_sample(find(screen == screen(screen_start(k)))+t-1) = t;
            end
        end
    end
    
    
    % Compute indices of first screen in a trial
    
    start_samples = find(number_sample == 1);
    
    % Compute indices of occluded samples in each trial
    
    for t = 1:length(start_samples)
        for z = 1:size(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples, 2)
            if ~isnan(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(t, z))
                number_sample_occluded = [number_sample_occluded (start_samples(t)+behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(t, z)-1)];
                
            
            end
        end
    end
    
     % Sample number in trial
    
    number_sample_shared = number_sample;
    number_sample_shared(number_sample_occluded) = nan;
    number_sample_shared_no_nan = number_sample_shared(~isnan(number_sample_shared));
    
    %Positions of shared/occluded samples in vector: 
    
    sample_index_shared = find(number_sample_shared(~isnan(number_sample_shared))); 
    sample_index_occluded = find(number_sample_shared(isnan(number_sample_shared)));
    
    % Number of occluded sample in trial
    
    occluded_sample = nan(size(number_sample_shared));
    occluded_sample(number_sample_occluded) = number_sample(number_sample_occluded);
    occluded_sample_no_nan = occluded_sample(~isnan(occluded_sample));
    
    
    

    
    % Correct answer for self decision
    
    correct_self = behavioural_data.decision.self.correct(:,g);
    
    % Correct answer for other decision
    
    correct_other = behavioural_data.decision.other.correct(:,g);
    
    
      for l = 2:length(start_samples)
        if correct_self(l) == 1
            sample_sequence_answer_self(start_samples(l-1):start_samples(l)) = 1; 
        elseif correct_self(l) == 0
            sample_sequence_answer_self(start_samples(l-1):start_samples(l)) = 0; 
        end
        
        if correct_other(l) == 1
             sample_sequence_answer_other(start_samples(l-1):start_samples(l)) = 1; 

        elseif correct_other(l) == 0
             sample_sequence_answer_other(start_samples(l-1):start_samples(l)) = 0; 
        end
    end
    
    % Regressors for Phase 1 and 5 are split after creation to contian only
    % information of shared/occluded samples respectively
    
    % Boolean statement to see whether the total evidence changes sign after a
    % screen (decision change (yes/no))
    
    tip_scale_sign(number_sample(~isnan(number_sample))~=1) =   ((prior_sumv(number_sample(~isnan(number_sample))~=1)>0) & (sum_vec(number_sample(~isnan(number_sample))~=1)<0)) | ((prior_sumv(number_sample(~isnan(number_sample))~=1)<0) & (sum_vec(number_sample(~isnan(number_sample))~=1)>0));
    tip_scale_sign(number_sample(~isnan(number_sample))==1)= 0;
    
    
    %% Phase 1 (Sampling Shared)
    
    % Prior sum of total information

    prior_sum_self = prior_sumv; 
    prior_sum_other(1:length(sum_priorSample_signals_other(~isnan(sum_Sample_signals))))= (sum_priorSample_signals_other(~isnan(sum_Sample_signals)));
    
    % Prior sum of shared information 
    
    prior_sum_other_shared = prior_sum_other(sample_index_shared);
    prior_sum_self_shared = prior_sum_self(sample_index_shared); 
    
    %Difference in prior_sum for self and other
    
    prior_sum_diff = prior_sum_self - prior_sum_other;
    
    prior_sum_diff_shared = prior_sum_diff(sample_index_shared);  
    
    % Signal update for self (total information) and other (shared information) 
    
    for k = 2:length(prior_sumv)
        signal_update_self(k) = (prior_sumv(k) - prior_sumv(k-1)) *(1-2*(prior_sumv(k-1) < 0));
        signal_update_other(k) = (prior_sum_other(k) - prior_sum_other(k-1)) * (1-2*(prior_sum_other(k-1)<0));
    end
    signal_update_self(start_samples) = 0; 
    signal_update_other(start_samples) = 0; 

    
    
    % Diff in signal update between own and other
    
    signal_update_diff = signal_update_self - signal_update_other; 
    
    signal_update_self_shared = signal_update_self(sample_index_shared);
    signal_update_other_shared = signal_update_other(sample_index_shared);
    
    signal_update_diff_shared = signal_update_diff(sample_index_shared);    
    



    % Change in outcome at current sample (yes/no)
    
    tip_scale_sign_shared = tip_scale_sign;
    tip_scale_sign_shared(sample_index_occluded) = nan;
    tip_scale_sign_shared_no_nan = tip_scale_sign_shared(~isnan(tip_scale_sign_shared));
    
    
   
    
    
    
    %% Phase 2 (Decision self)
    
    % Evidence shared
    
    shared_ev_self = behavioural_data.Subject(g).schedule.MainSchedule.OtherEvidence';
    shared_ev_self_abs = abs(shared_ev_self);
    shared_ev_self_signed(correct_self == 1) = abs(shared_ev_self(correct_self == 1)); 
    shared_ev_self_signed(correct_self == 0) = -abs(shared_ev_self(correct_self == 0)); 

        
    % Evidence occluded
    
    occluded_ev_self = abs(behavioural_data.Subject(g).schedule.MainSchedule.OverallEvidence - behavioural_data.Subject(g).schedule.MainSchedule.OtherEvidence);
    %occluded_ev_self_abs = abs(occluded_ev_self); 
    occluded_ev_self_signed(correct_self == 1) = abs(occluded_ev_self(correct_self == 1));
    occluded_ev_self_signed(correct_self == 0) = -abs(occluded_ev_self(correct_self == 0));
    
    
    
    %% Phase 3 (Decision other)
    
    % Evidence shared
    
    shared_ev_other = behavioural_data.Subject(g).schedule.MainSchedule.OtherEvidence';
    shared_ev_other_abs = abs(shared_ev_other);
    shared_ev_other_signed(correct_other == 1) = abs(shared_ev_other(correct_other == 1)); 
    shared_ev_other_signed(correct_other == 0) = -abs(shared_ev_other(correct_other == 0)); 

        
    % Evidence occluded
    
    occluded_ev_other = abs(behavioural_data.Subject(g).schedule.MainSchedule.OverallEvidence - behavioural_data.Subject(g).schedule.MainSchedule.OtherEvidence);
    occluded_ev_other_2 = [];
    shared_ev_other_2 = [];

    for r = 1:size(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples, 1)
        occluded_ev_other_2 = [occluded_ev_other_2 sum(abs(diff(behavioural_data.Subject(g).schedule.MainSchedule.NrDots(:,behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(r,~isnan(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(r,:))),r))))];
       % shared_ev_other_2 = [shared_ev_other_2 sum(abs(diff(behavioural_data.Subject(g).schedule.MainSchedule.NrDots(:,behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(r,1:end, find(~isnan(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(r,:))),r)))))];
    end
    
    if occluded_ev_other_2 ~= occluded_ev_other % || shared_ev_other ~= shared_ev_other_2
        fprintf('Error: Occluded ev Vecotrs are not the same!')
        keyboard;
   
        
    end
    
    occluded_ev_other_signed(correct_other == 1) = abs(occluded_ev_other(correct_other == 1));
    occluded_ev_other_signed(correct_other == 0) = -abs(occluded_ev_other(correct_other == 0));
    
    %% Phase 4 (Outcome)
    
    % Discrepancy in decision for other and behaviour for other (yes/no)
    
    prediction_disc = behavioural_data.decision.other.redchoice(:,g) ~= behavioural_data.Subject(g).schedule.MainSchedule.SimulateOtherResponsesB;
    
    
    % Confidence self (chance to choose red for yourself given the evidence in
    % a trial)
    
    
    evidence_against_blue_self_single = squeeze(sum(behavioural_data.Subject(g).schedule.MainSchedule.NrDots(2,:,:), 2));
    evidence_for_blue_self_single = squeeze(sum(behavioural_data.Subject(g).schedule.MainSchedule.NrDots(1,:,:), 2));
    
    
    % Save vectors in cell array for 'group' case 
    % evidence_for_red_group{g} = evidence_for_red_single;
    % evidence_against_red_group{g} = evidence_against_red_single;
    
    %Save red chocies as vector 
    
    choices_self_red_single = behavioural_data.Subject(g).decision.self.redchoice';
    
    bins_self = linspace(min(evidence_for_blue_self_single-evidence_against_blue_self_single),max(evidence_for_blue_self_single-evidence_against_blue_self_single), 11);
    
    
     %For blue
     
     
   [self_Ts(g),shiftT_self(g),self_equipoints(g)] =  plot_choice_softmax('single',evidence_for_blue_self_single, evidence_against_blue_self_single, choices_self_red_single ,bins_self);
    
   
   
   evidence_for_blue_self_single_shift = evidence_for_blue_self_single-self_equipoints(g);
   
   %P_choice_1 = blue
   
   p_choice_self_1(:, g) = exp((evidence_for_blue_self_single_shift/shiftT_self(g)))./(exp((evidence_for_blue_self_single_shift/shiftT_self(g))) + exp((evidence_against_blue_self_single/shiftT_self(g))));
   
   %P_choice_2 = red
   
   p_choice_self_2(:, g) = exp((evidence_against_blue_self_single/shiftT_self(g)))./(exp((evidence_for_blue_self_single_shift/shiftT_self(g))) + exp((evidence_against_blue_self_single/shiftT_self(g)))); % should be same as 1-p_choice1

   confidence_self(choices_self_red_single == 1) = p_choice_self_2(choices_self_red_single == 1,g);  
   confidence_self(choices_self_red_single == 0) = p_choice_self_1(choices_self_red_single == 0,g);

    
    
    % Own confidence in other chocies 
    
    choices_other_red_single =[];
    evidence_against_blue_other_single = [];
    evidence_for_blue_other_single = [];
    evidence_for_blue_other_single_shift = [];
    bins_other = [];
    
    
    choices_other_red_single = behavioural_data.Subject(g).decision.other.redchoice';
    
    
    obscured_mat = repmat(1:5, 160, 1);
    
    for k = 1:size(obscured_mat, 1)
        obscured_mat(k,behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(k,~isnan(behavioural_data.Subject(g).schedule.MainSchedule.AllOscuredsamples(k,:))))=nan;    
        evidence_against_blue_other_single(k) = squeeze(sum(behavioural_data.Subject(g).schedule.MainSchedule.NrDots(2,obscured_mat(k,~isnan(obscured_mat(k,:))),k), 2));
        evidence_for_blue_other_single(k) = squeeze(sum(behavioural_data.Subject(g).schedule.MainSchedule.NrDots(1,obscured_mat(k,~isnan(obscured_mat(k,:))),k), 2));
    end
    
    bins_other = linspace(min(evidence_for_blue_other_single-evidence_against_blue_other_single),max(evidence_for_blue_other_single-evidence_against_blue_other_single), 8);

    [other_Ts(g), shiftT_other(g), other_equipoints(g)] =  plot_choice_softmax('single',evidence_for_blue_other_single', evidence_against_blue_other_single', choices_other_red_single ,bins_other);
    
   
   
   evidence_for_blue_other_single_shift = evidence_for_blue_other_single'-other_equipoints(g);
   p_choice_other_1(:, g) = exp((evidence_for_blue_other_single_shift/shiftT_other(g)))./(exp((evidence_for_blue_other_single_shift/shiftT_other(g))) + exp((evidence_against_blue_other_single'/shiftT_other(g))));
   p_choice_other_2(:, g) = exp((evidence_against_blue_other_single'/shiftT_other(g)))./(exp((evidence_for_blue_other_single_shift/shiftT_other(g))) + exp((evidence_against_blue_other_single'/shiftT_other(g)))); % should be same as 1-p_choice1
    
 
   confidence_other(choices_other_red_single == 1) = p_choice_other_2(choices_other_red_single == 1,g);  
   confidence_other(choices_other_red_single == 0) = p_choice_other_1(choices_other_red_single == 0,g);
   
    
    % Prediction error other (choice(yes/no) - confidence in choice)
    
    prediction_error = prediction_disc - confidence_other;
    
    % Agent_ability (Level 1 = good, Level 2 = bad)
    
    agent_ability = nan(size(behavioural_data.Subject(g).agent_rating.agent_ID, 1)* rating_distance, 1);
    
    for a = 1:max(behavioural_data.Subject(g).agent_rating.agent_ID)
        if a == 1
        agent_ability(min(find(behavioural_data.Subject(g).agent_rating.agent_ID == a)):max(find(behavioural_data.Subject(g).agent_rating.agent_ID == a))*rating_distance) = ...
            behavioural_data.Subject(g).schedule.MainSchedule.AgentAblLvL(a);
        elseif a > 1
           agent_ability((min(find(behavioural_data.Subject(g).agent_rating.agent_ID == a))*rating_distance-9):max(find(behavioural_data.Subject(g).agent_rating.agent_ID == a))*rating_distance) = ...
            behavioural_data.Subject(g).schedule.MainSchedule.AgentAblLvL(a); 
        end
    end
    
    % Interaction of ability and prediction discrepancy (normalised product)
    
    inter_disc_ability = normalise(prediction_disc) .* normalise(agent_ability); 

    % Interaction of ability and prediction error (normalised product)
    
    inter_error_ability = normalise(prediction_error) .* normalise(agent_ability); 

    
    %% Phase 5 (Sampling Occluded)
    
    
    % Sum of evidence prior to occluded sample for self/other and
    % difference
    
    prior_sum_self_occluded = prior_sum_self(sample_index_occluded); 
    prior_sum_other_occluded = prior_sum_other(sample_index_occluded);
    prior_sum_diff_occluded = prior_sum_diff(sample_index_occluded); 

    % Signal update at occluded samples for self/other and difference
    
    signal_update_other_occluded = signal_update_other(sample_index_occluded); 
    signal_update_self_occluded = signal_update_self(sample_index_occluded); 
    signal_update_diff_occluded = signal_update_diff(sample_index_occluded);

     
    % Change in outcome at current occluded sample (yes/no)
    
    tip_scale_sign_occluded = nan(size(tip_scale_sign));
    tip_scale_sign_occluded(sample_index_occluded) = tip_scale_sign(sample_index_occluded);
    tip_scale_sign_occluded_no_nan = tip_scale_sign_occluded(~isnan(tip_scale_sign_occluded));
    
    
    
  
    
    
    
    %% Onset times
    
    Onset = [];
    
    % Phase 1 (onsets of shared samples)
    
    number_sample_onset_shared = behavioural_data.Subject(g).time(screen) - behavioural_data.Subject(g).time(min(find(strcmp(behavioural_data.Subject(g).event, 'Pulse'))));
    number_sample_onset_shared(number_sample_occluded) = nan;
    number_sample_onset_shared_no_nan = number_sample_onset_shared(~isnan(number_sample_onset_shared));
    
    
    % Phase 2 (onset of question self)
    
    question_self_onset = behavioural_data.output.time(behavioural_data.output.indices.question_you(:,g), g) - behavioural_data.output.time(min(find(strcmp(behavioural_data.output.event_type(:,g), 'Pulse'))),g);
    DecSelf_onset = question_self_onset;
    
    
    % Phase 3 (onset of question other)
    
    question_other_onset = behavioural_data.output.time(behavioural_data.output.indices.question_other(:,g),g) - behavioural_data.output.time(min(find(strcmp(behavioural_data.output.event_type(:,g), 'Pulse'))), g);
    DecOther_onset = question_other_onset;
    
    
    % Phase 4 (onset of feedback other's behaviour)
    
    Outcome_onset = question_other_onset;
    
    % Phase 5 (onset of occluded sample)
    
    number_sample_onset_occluded = behavioural_data.Subject(g).time(screen(number_sample_occluded)) - behavioural_data.Subject(g).time(min(find(strcmp(behavioural_data.Subject(g).event, 'Pulse'))));
    
       
    %% Make Design Matrices
    
    designName='design05';
    
    % Sampling Event
    
    % Exclude occluded samples
    
    
    % Now (22.05.17) absolute, later sign it by answer (error = negative
    % evidence, no error = positive evidence) 
    Phase{1}.reglabels={'PreviousSampleSum';'SignalUpdate'; 'DiffPriorSum'; 'DiffSignalUpdate'; 'choiceChange'; 'sampleNr'};
    Phase{1}.regs=[abs(prior_sum_self_shared)' signal_update_self_shared' prior_sum_diff_shared' signal_update_diff_shared' tip_scale_sign_shared_no_nan' number_sample_shared_no_nan'];
    Phase{1}.Type='SamplingSharedSelfSigned';
    Phase{1}.Duration=0;
    
    % Own Decision Event
    
    Phase{2}.reglabels={'SharedEvidence\_Self';'OccludedEvidence\_Self'};
    Phase{2}.regs=[shared_ev_self_signed' occluded_ev_self_signed'];
    Phase{2}.Type='DecSelf';
    Phase{2}.Duration=0;
    
    % Other Decision Eventc
    Phase{3}.reglabels={'SharedEvidence\_Other';'OccludedEvidence\_Other'};
    Phase{3}.regs=[shared_ev_other_signed' occluded_ev_other_signed'];
    Phase{3}.Type='DecOther';
    Phase{3}.Duration=0;
    
    % Outcome Event  (Temporally very close to the Other decision?)
    
   
    
    Phase{4}.reglabels={'Prediction\_error'; 'ConfidenceOther'; 'ConfidenceSelf'; 'AgentAbility'};
    Phase{4}.regs=[prediction_error confidence_other confidence_self agent_ability];
    Phase{4}.Type='Outcome';
    Phase{4}.Duration=0;
    
    Phase{5}.reglabels={'PreviousSampleSum\_occluded_self\_signed'; 'SignaUpdateOccluded'; 'PriorSumDiffOccluded'; 'SignalUpdateDiffOccluded'; 'Choice\_change\_occluded'; 'SampleNr\_occluded'};
    Phase{5}.regs = [abs(prior_sum_self_occluded)' signal_update_self_occluded' prior_sum_diff_occluded' signal_update_diff_occluded' tip_scale_sign_occluded_no_nan' occluded_sample_no_nan'];
    Phase{5}.Type = 'SamplingOccludedSelfSigned';
    Phase{5}.Duration=0;
   
    

    for p=1:length(Phase)
        for c=1:size(Phase{p}.regs,2)
            Phase{p}.normregs(~isnan(Phase{p}.regs(:,c)),c)=normalise(Phase{p}.regs(~isnan(Phase{p}.regs(:,c)),c));
            Phase{p}.normregs(isnan(Phase{p}.regs(:,c)),c)=0;
        end
    end
    
    Phase{4}.reglabels = [Phase{4}.reglabels; 'InteractionPredAbility'];
    Phase{4}.normregs = [Phase{4}.normregs inter_error_ability];
    
    %% Text Out
    %filepath=['/Users/nkolling/Documents/SocialEvidence/fMRI_Regressors/regs/' designName '/'];
    filepath=['/home/marius/Dokumente/Project_2/fMRI_Regressors/' version '/' designName '/S' Subject_ID{g} '/'];
    
    if ~isdir(filepath)
        mkdir(filepath);
    end
    
    for p=1:length(Phase)
        if strcmp(Phase{p}.Type,'SamplingSharedSelfSigned')
            Onset=number_sample_onset_shared_no_nan;
        elseif strcmp(Phase{p}.Type,'DecSelf')
            Onset=DecSelf_onset;
        elseif strcmp(Phase{p}.Type,'DecOther')
            Onset=DecOther_onset;
        elseif strcmp(Phase{p}.Type,'Outcome')
            Onset=Outcome_onset;
        elseif strcmp(Phase{p}.Type, 'SamplingOccludedSelfSigned')
            Onset=number_sample_onset_occluded;
        elseif strcmp(Phase{p}.Type,'SamplingSharedOtherSigned')
            Onset=number_sample_onset_shared_no_nan;
        elseif strcmp(Phase{p}.Type,'SamplingOccludedOtherSigned')
            Onset=number_sample_onset_occluded;
        end
        
        Duration=Phase{p}.Duration;
        filename=['Phase' num2str(p) '_S' Subject_ID{g} '_' 'main.txt'];
        dlmwrite([filepath filename],[Onset/10000 repmat(Duration,length(Onset),1) ones(length(Onset),1)],'delimiter', '\t');
        
        
        for c=1:size(Phase{p}.normregs,2)
            filename=['Phase' num2str(p) '_reg' num2str(c) '_S' Subject_ID{g} '.txt'];
            dlmwrite([filepath filename],[Onset/10000 repmat(Duration,length(Onset),1) Phase{p}.normregs(:,c)],'delimiter', '\t');
        end
    end
    
    %     %Save vectors in matrices
    %
    %     all_screen(1:length(screen),g) = screen;
    %     all_screen_start(1:length(screen_start),g) = screen_start;
    %     all_screen_end(1:length(screen_end),g) =  screen_end;
    %     all_sum_vec(1:length(sum_vec),g) = sum_vec;
    %
    % %     all_signal_update_duration(1:length(signal_update_duration),g) = signal_update_duration;
    % %     all_number_sample_onset(1:length(number_sample_onset),g) = number_sample_onset;
    %
    %
    %     % Sampling
    %
    %     all_prior_sumv(1:length(prior_sumv),g) = prior_sumv;
    %     all_signal_update(1:length(signal_update),g) = signal_update;
    %     all_tip_scale_sign(1:length(tip_scale_sign),g) = tip_scale_sign;
    %     all_number_sample(1:length(number_sample),g) = number_sample;
    %
    %
    %     % Decision Phase
    %
    %     all_shared_ev(1:length(shared_ev),g) = shared_ev;
    %     all_occluded_ev(1:length(occluded_ev),g) = occluded_ev;
    %
    %
    %     % Outcome
    %
    %     all_prediction_disc(1:length(prediction_disc),g) = prediction_disc;
    %     all_correct_self(1:length(correct_self), g) = correct_self;
    %     all_correct_other(1:length(correct_other), g) = correct_other;
    %
    %
    %     % Onsets
    %
    %     all_sample_onset(1:length(sample_onset),g) = sample_onset;
    %     all_DecSelf_onset(1:length(DecSelf_onset),g) = DecSelf_onset;
    %     all_DecOther_onset(1:length(DecOther_onset),g) = DecOther_onset;
    %     all_Outcome_onset(1:length(Outcome_onset),g) = Outcome_onset;
    %
    %     % all_question_self_onset(1:length(question_self_onset),g) = question_self_onset;
    %     % all_question_other_onset(1:length(question_other_onset),g) = question_other_onset;
    %
    %     % Normalised Vectors as Matrix
    %
    %     % Sampling
    %
    %     all_prior_sumv_norm(1:length(prior_sumv),g) = Phase{1}.normregs(:,1);
    %     all_signal_update_norm(1:length(signal_update),g) = Phase{1}.normregs(:,2);
    %     all_tip_scale_sign_norm(1:length(tip_scale_sign),g) = Phase{1}.normregs(:,3);
    %     all_number_sample_norm(1:length(number_sample),g) = Phase{1}.normregs(:,4);
    %
    %
    %     % Decision Phase
    %
    %     all_shared_ev_norm(1:length(shared_ev),g) = Phase{2}.normregs(:,1);
    %     all_occluded_ev_norm(1:length(occluded_ev),g) = Phase{2}.normregs(:,2);
    %
    %
    %     % Outcome
    %
    %     all_prediction_disc_norm(1:length(prediction_disc),g) = Phase{4}.normregs(:,1);
    %     all_correct_self_norm(1:length(correct_self), g) = Phase{4}.normregs(:,2);
    %     all_correct_other_norm(1:length(correct_other), g) = Phase{4}.normregs(:,3);
    %
    
    
    
    
    % %% Safe regressors in struct one field per Participant
    %
    % % Onsets
    %
    % Regressor_data.Subject{g}.Phase{1}.onset = all_sample_onset;
    % Regressor_data.Subject{g}.Phase{2}.onset = all_DecSelf_onset;
    % Regressor_data.Subject{g}.Phase{3}.onset = all_DecOther_onset;
    % Regressor_data.Subject{g}.Phase{4}.onset = all_Outcome_onset;
    %
    %
    %
    % % Regressors
    %
    % % Phase 1
    %
    % Regressor_data.Subject{g}.Phase{1}.regs.tip_scale_sign = tip_scale_sign;
    % Regressor_data.Subject{g}.Phase{1}.regs.prior_sum = prior_sumv;
    % Regressor_data.Subject{g}.Phase{1}.regs.signal_update = signal_update;
    % Regressor_data.Subject{g}.Phase{1}.regs.number_sample = number_sample;
    %
    %
    % % Phase 2
    %
    % Regressor_data.Subject{g}.Phase{2}.regs.shared_evidence = shared_ev;
    % Regressor_data.Subject{g}.Phase{2}.regs.occluded_evidence = occluded_ev;
    %
    %
    % % Phase 3
    %
    % Regressor_data.Subject{g}.Phase{3}.regs.shared_evidence = shared_ev;
    % Regressor_data.Subject{g}.Phase{3}.regs.occluded_evidence = occluded_ev;
    %
    %
    % % Phase 4
    %
    % Regressor_data.Subject{g}.Phase{4}.regs.self_correct = self_correct;
    % Regressor_data.Subject{g}.Phase{4}.regs.other_correct = other_correct;
    % Regressor_data.Subject{g}.Phase{4}.regs.prediction_disc = prediction_disc;
    %
    % % Regressors normalized
    %
    %
    % % Phase 1
    %
    % Regressor_data.Subject{g}.Phase{1}.normregs.prior_sum = Phase{1}.normregs(:,1);
    % Regressor_data.Subject{g}.Phase{1}.normregs.signal_update = Phase{1}.normregs(:,2);
    % Regressor_data.Subject{g}.Phase{1}.normregs.tip_scale_sign = Phase{1}.normregs(:,3);
    % Regressor_data.Subject{g}.Phase{1}.normregs.number_sample = Phase{1}.normregs(:,4);
    %
    %
    % % Phase 2
    %
    % Regressor_data.Subject{g}.Phase{2}.normregs.shared_evidence = Phase{2}.normregs(:,1);
    % Regressor_data.Subject{g}.Phase{2}.normregs.occluded_evidence = Phase{2}.normregs(:,2);
    %
    %
    % % Phase 3
    %
    % Regressor_data.Subject{g}.Phase{3}.normregs.shared_evidence = Phase{3}.normregs(:,1);
    % Regressor_data.Subject{g}.Phase{3}.normregs.occluded_evidence = Phase{3}.normregs(:,2);
    %
    %
    % % Phase 4
    %
    % Regressor_data.Phase{4}.normregs.self_correct = Phase{4}.normregs(:,1);
    % Regressor_data.Phase{4}.normregs.other_correct = Phase{4}.normregs(:,2);
    % Regressor_data.Phase{4}.normregs.prediction_disc = Phase{4}.normregs(:,3);
    
    
%% Plot regressors to see relationship between them 

% figure('name', ['PriorSum vs update (diff occluded) Subject:', Subject_ID{g}])
% hold on 
% plot(normalise(abs(prior_sum_diff_occluded)))
% plot(normalise(signal_update_diff_occluded))
%     
% 
% figure('name', ['PriorSum vs update (occluded) Subject:', Subject_ID{g}])
% hold on 
% plot(normalise(abs(prior_sum_self_occluded)))
% plot(normalise(signal_update_self_occluded))
% 
% figure('name', ['PriorSum vs update (diff shared) Subject:', Subject_ID{g}])
% hold on 
% plot(prior_sum_diff_shared)
% plot(prior_sum_self_shared)
% 
% figure('name', ['PriorSum vs update shared Subject:', Subject_ID{g}])
% hold on 
% plot(normalise(abs(prior_sum_self_shared)))
% plot(normalise(signal_update_self_shared))
% 
% figure('name', ['PriorSum self vs PriorSumOther shared Subject:', Subject_ID{g}])
% hold on 
% plot(normalise(abs(prior_sum_self_shared)))
% plot(normalise(abs(prior_sum_other_shared)))
% 
% figure('name', ['Update self vs update diff shared Subject:', Subject_ID{g}])
% hold on 
% plot(normalise(signal_update_self_shared))
% plot(normalise(signal_update_diff_shared))
% 
% figure('name', ['Update self vs update other shared Subject:', Subject_ID{g}])
% hold on 
% plot(signal_update_self_occluded)
% plot(signal_update_diff_occluded)


    

%% Compute correlation for regressors and plot correlation as a function of schedule/participant

% Important: Must be the same across schedules! 

corr_prior_diff_occluded = corrcoef(normalise(abs(prior_sum_diff_occluded)), normalise(signal_update_diff_occluded));
corr_prior_diff_shared = corrcoef(normalise(abs(prior_sum_diff_shared)), normalise(signal_update_diff_shared));
corr_update_occluded = corrcoef(normalise(abs(prior_sum_self_occluded)), normalise(signal_update_self_occluded));
corr_update_shared = corrcoef(normalise(abs(prior_sum_self_shared)), normalise(signal_update_self_shared));

corr_update_diff_shared = corrcoef(normalise(signal_update_self_shared), normalise(signal_update_diff_shared));
corr_update_diff_occluded = corrcoef(normalise(signal_update_self_occluded), normalise(signal_update_diff_occluded));

corr_vec_prior_diff_shared(g) = corr_prior_diff_shared(2); 

corr_vec_update_diff_shared(g) = corr_update_diff_shared(2); 
corr_vec_update_diff_occluded(g) = corr_update_diff_occluded(2); 

corr_vec_prior_diff_occluded(g) = corr_prior_diff_occluded(2);
corr_vec_shared(g) = corr_update_shared(2); 
corr_vec_occluded(g) = corr_update_occluded(2);

corr_update_diff_occluded = [];
corr_update_diff_shared = [];
corr_update_diff_occluded = [];
corr_update_diff_shared = [];
    
end


%% Plot correlations over participants/schedules

figure('Name', 'Correlation prior sum diff/update_diff occluded (Phase 5)')
plot(corr_vec_prior_diff_occluded)
hold on 
plot(corr_vec_occluded)
legend('Corr prior sum, update', 'Corr prior sum diff, update diff' )


figure('Name', 'Correlation prior sum diff/update_diff shared (Phase 1)')
plot(corr_vec_update_diff_shared)
hold on 
plot(corr_vec_shared)
legend('Corr prior sum; update', 'Corr prior sum diff, update diff' )

%% Safe in struct/ .mat file as matrix per phase

% Onsets


Regressor_data.Phase{1}.onset = all_sample_onset;
Regressor_data.Phase{2}.onset = all_DecSelf_onset;
Regressor_data.Phase{3}.onset = all_DecOther_onset;
Regressor_data.Phase{4}.onset = all_Outcome_onset;


% Regressors

%Phase 1
Regressor_data.Phase{1}.regs.tip_scale_sign = all_tip_scale_sign;
Regressor_data.Phase{1}.regs.prior_sum = all_prior_sumv;
Regressor_data.Phase{1}.regs.signal_update = all_signal_update;
Regressor_data.Phase{1}.regs.number_sample = all_number_sample;


% Phase 2

Regressor_data.Phase{2}.regs.shared_evidence = all_shared_ev;
Regressor_data.Phase{2}.regs.occluded_evidence = all_occluded_ev;


% Phase 3

Regressor_data.Phase{3}.regs.shared_evidence = all_shared_ev;
Regressor_data.Phase{3}.regs.occluded_evidence = all_occluded_ev;


% Phase 4

% Regressor_data.Phase{4}.regs.self_correct = all_self_correct;
% Regressor_data.Phase{4}.regs.other_correct = all_other_correct;
% Regressor_data.Phase{4}.regs.prediction_disc = all_prediction_disc;


%Regressors normalized


% %Phase 1
% Regressor_data.Phase{1}.normregs.tip_scale_sign = all_tip_scale_sign_norm;
% Regressor_data.Phase{1}.normregs.prior_sum = all_prior_sumv_norm;
% Regressor_data.Phase{1}.normregs.signal_update = all_signal_update_norm;
% Regressor_data.Phase{1}.normregs.number_sample = all_number_sample_norm;
%
%
% % Phase 2
%
% Regressor_data.Phase{2}.normregs.shared_evidence = all_shared_ev_norm;
% Regressor_data.Phase{2}.normregs.occluded_evidence = all_occluded_ev_norm;
%
%
% % Phase 3
%
% Regressor_data.Phase{3}.normregs.shared_evidence = all_shared_ev_norm;
% Regressor_data.Phase{3}.normregs.occluded_evidence = all_occluded_ev_norm;
%
%
% % Phase 4
%
% % Regressor_data.Phase{4}.normregs.self_correct = all_self_correct_norm;
% % Regressor_data.Phase{4}.normregs.other_correct = all_other_correct_norm;
% % Regressor_data.Phase{4}.normregs.prediction_disc = all_prediction_disc_norm;
%







%% Correlations within Phases

for p=1:length(Phase)
    Phase{p}.cormat=corrcoef(Phase{p}.normregs);
    figure('name',['WithinSubject Corr Phase' num2str(p)],'color',[1 1 1]);subplot(1,1,1,'XTickLabel',Phase{p}.reglabels,'XTick',1:length(Phase{p}.reglabels),'YTickLabel',Phase{p}.reglabels,'YTick',1:length(Phase{p}.reglabels));
    hold on;imagesc(Phase{p}.cormat,[-1 1]);subplot(1,1,1,'XTickLabel',Phase{p}.reglabels);
end


% Group regressors of same length in one matrix (Phase 2-4)

normregs = [];
labels = {};
for l=2:4
    normregs = [normregs Phase{l}.normregs];
    labels = [labels; Phase{l}.reglabels];
end

%normregs = [normregs question_self_onset/100000 question_other_onset/10000 Outcome_onset/10000]; 
%labels = [labels; 'Question\_Self\_Onset'; 'Quesiton\_other\_onset'; 'Outcome\_onset'];

%% Plot correlation across phases

cormat_all = corrcoef(normregs)
figure('name',['WithinSubject Corr allPhases'],'color',[1 1 1]);subplot(1,1,1,'XTickLabel',labels,'XTick',1:length(labels),'YTickLabel',labels,'YTick',1:length(labels));
    hold on;imagesc(cormat_all,[-1 1]);subplot(1,1,1,'XTickLabel',labels)
    
    

%keyboard;


