% script to compute figures of analyzed EMG data
% Open matfile which contains a struct all with timecourse data as well as
% condition info. Likeability 1=positive, -1=negative. Affect 1=positive,
% -1=negative
% Marijn Struiksma December 2015

%% Required directories

% ==> Name your analysis directory where the results are stored:
%sjCode	SjNr	markertrial	incCorr_mean	incZygo_mean	likeConditie_mean	affectConditie_mean	ZLkblCorrPctTrial_mean	ZAfctCorrPctTrial_mean	ZLkblZygoPctTrial_mean	ZAfctZygoPctTrial_mean
wdir = 'M:\jvbproject\Björn\ACE-BRaT\Analyzer na EDFChecker v1.0\analysis\';
% load z-scores file
fileZ = dir([wdir, 'ACE-BRaT_zscores_ordered.txt']);
zscores=load([wdir,fileZ.name]); %loads struct all that contains the segments data: struct segment
% zscores contains the following columns:
% [ppnr trialnr incCor incZygo LikeCond AffectCond zLikeCorr zLikeZygo zAffectCorr zAffectZygo]
% [1 52 1 1 -1 1 0.09 0.72 -0.64 -0.79]
% Trials are included for -1.98 =< z <= 1.98

% load matfile with segment results
%file = dir([wdir, 'timecourse2sjs.mat']);
file = dir([wdir, 'ACE-BRaT_timecourse_60sjs_18-Sep-2015.mat']);
%load([wdir,file.name]); %loads struct all that contains the segments data: struct segment
%segment=segment2sjs; % only necessary for debugging. put in comment for final analysis

% dir where figures should be stored
fdir = 'M:\jvbproject\Björn\ACE-BRaT\Analyzer na EDFChecker v1.0\analysis\figures\';
if (exist(fdir)~=7)
    mkdir(fdir);
end




%% compute mean timecourse per participant per condition
dims=size(segment);
nrsjs=dims(2);
nrtrials=length(segment(1).marker);

AllLikePosCorr=[];
AllLikeNegCorr=[];
AllLikePosEventPosCorr=[];
AllLikePosEventNegCorr=[];
AllLikeNegEventPosCorr=[];
AllLikeNegEventNegCorr=[];

AllLikePosZygo=[];
AllLikeNegZygo=[];
AllLikePosEventPosZygo=[];
AllLikePosEventNegZygo=[];
AllLikeNegEventPosZygo=[];
AllLikeNegEventNegZygo=[];
meanInclusionPr=[];
ppdetails=[];
counter=1;
for iSjs = 1:nrsjs
    %disp('sjnr:');
    %disp(iSjs);
    %disp('trialnr:');
    LikePosCorr=[];
    LikeNegCorr=[];
    LikePosEventPosCorr=[];
    LikePosEventNegCorr=[];
    LikeNegEventPosCorr=[];
    LikeNegEventNegCorr=[];
    
    LikePosZygo=[];
    LikeNegZygo=[];
    LikePosEventPosZygo=[];
    LikePosEventNegZygo=[];
    LikeNegEventPosZygo=[];
    LikeNegEventNegZygo=[];
    inclusionPr=[];
    for iTrials = 1:nrtrials
        %disp(iTrials);
        %make sure the correct zscore info is used
        if ~zscores(counter,1)==iSjs
            disp('sj does not match');
            disp(iSjs);
        end
        if ~zscores(counter,2)==iTrials
            disp('trial does not match');
            disp(iTrials);
        end
        % should trial be included? based on baseline
        inclusion=(segment(iSjs).incCorr(iTrials))*(segment(iSjs).incZygo(iTrials));
        %disp(inclusion)
        % should trial be included? based on zscore
        PzscoresCorrLike=abs(zscores(counter,7));
        inclusionCorrLike=(max(PzscoresCorrLike)<=1.98)*inclusion;
        PzscoresZygoLike=abs(zscores(counter,8));
        inclusionZygoLike=(max(PzscoresZygoLike)<=1.98)*inclusion;
        PzscoresCorrAffect=abs(zscores(counter,9));
        inclusionCorrAffect=(max(PzscoresCorrAffect)<=1.98)*inclusion;
        PzscoresZygoAffect=abs(zscores(counter,10));
        inclusionZygoAffect=(max(PzscoresZygoAffect)<=1.98)*inclusion;
        
        % Likeability condition
        condLike=segment(iSjs).likeConditie(iTrials);
        % Affect event condition
        condAffect=segment(iSjs).affectConditie(iTrials);
        
        % Collect timecourse data
        if inclusionCorrLike==1
            switch condLike
                case 1 %positive
                    LikePosCorr=vertcat(LikePosCorr,segment(iSjs).bin100ms(iTrials).CorrPr);
                case -1 %negative
                    LikeNegCorr=vertcat(LikeNegCorr,segment(iSjs).bin100ms(iTrials).CorrPr);
            end
        end
        
        if inclusionZygoLike==1
            switch condLike
                case 1 %positive
                    LikePosZygo=vertcat(LikePosZygo,segment(iSjs).bin100ms(iTrials).ZygoPr);
                case -1 %negative
                    LikeNegZygo=vertcat(LikeNegZygo,segment(iSjs).bin100ms(iTrials).ZygoPr);
            end
        end
        
        if inclusionCorrAffect==1
            switch condLike
                case 1 %positive
                    switch condAffect
                        case 1 %positive
                            LikePosEventPosCorr=vertcat(LikePosEventPosCorr,segment(iSjs).bin100ms(iTrials).CorrPr);
                        case -1 %negative
                            LikePosEventNegCorr=vertcat(LikePosEventNegCorr,segment(iSjs).bin100ms(iTrials).CorrPr);
                    end
                case -1 %negative
                    switch condAffect
                        case 1 %positive
                            LikeNegEventPosCorr=vertcat(LikeNegEventPosCorr,segment(iSjs).bin100ms(iTrials).CorrPr);
                        case -1 %negative
                            LikeNegEventNegCorr=vertcat(LikeNegEventNegCorr,segment(iSjs).bin100ms(iTrials).CorrPr);
                    end
            end
        end
        
        if inclusionZygoAffect==1
            switch condLike
                case 1 %positive
                    switch condAffect
                        case 1 %positive
                            LikePosEventPosZygo=vertcat(LikePosEventPosZygo,segment(iSjs).bin100ms(iTrials).ZygoPr);
                        case -1 %negative
                            LikePosEventNegZygo=vertcat(LikePosEventNegZygo,segment(iSjs).bin100ms(iTrials).ZygoPr);
                    end
                case -1 %negative
                    switch condAffect
                        case 1 %positive
                            LikeNegEventPosZygo=vertcat(LikeNegEventPosZygo,segment(iSjs).bin100ms(iTrials).ZygoPr);
                        case -1 %negative
                            LikeNegEventNegZygo=vertcat(LikeNegEventNegZygo,segment(iSjs).bin100ms(iTrials).ZygoPr);
                    end
            end
        end
        
        counter=counter+1;
    end %iTrials
    
    % compute means per participant
    AllLikePosCorr=vertcat(AllLikePosCorr,mean(LikePosCorr));
    AllLikeNegCorr=vertcat(AllLikeNegCorr,mean(LikeNegCorr));
    AllLikePosEventPosCorr=vertcat(AllLikePosEventPosCorr,mean(LikePosEventPosCorr));
    AllLikePosEventNegCorr=vertcat(AllLikePosEventNegCorr,mean(LikePosEventNegCorr));
    AllLikeNegEventPosCorr=vertcat(AllLikeNegEventPosCorr,mean(LikeNegEventPosCorr));
    AllLikeNegEventNegCorr=vertcat(AllLikeNegEventNegCorr,mean(LikeNegEventNegCorr));
    
    AllLikePosZygo=vertcat(AllLikePosZygo,mean(LikePosZygo));
    AllLikeNegZygo=vertcat(AllLikeNegZygo,mean(LikeNegZygo));
    AllLikePosEventPosZygo=vertcat(AllLikePosEventPosZygo,mean(LikePosEventPosZygo));
    AllLikePosEventNegZygo=vertcat(AllLikePosEventNegZygo,mean(LikePosEventNegZygo));
    AllLikeNegEventPosZygo=vertcat(AllLikeNegEventPosZygo,mean(LikeNegEventPosZygo));
    AllLikeNegEventNegZygo=vertcat(AllLikeNegEventNegZygo,mean(LikeNegEventNegZygo));
end %iSjs

%% compute mean timecourse and standard deviation per condition
MeanLikePosCorr=[];
MeanLikeNegCorr=[];
MeanLikePosEventPosCorr=[];
MeanLikePosEventNegCorr=[];
MeanLikeNegEventPosCorr=[];
MeanLikeNegEventNegCorr=[];
StdLikePosCorr=[];
StdLikeNegCorr=[];
StdLikePosEventPosCorr=[];
StdLikePosEventNegCorr=[];
StdLikeNegEventPosCorr=[];
StdLikeNegEventNegCorr=[];

MeanLikePosZygo=[];
MeanLikeNegZygo=[];
MeanLikePosEventPosZygo=[];
MeanikePosEventNegZygo=[];
MeanLikeNegEventPosZygo=[];
MeanLikeNegEventNegZygo=[];
StdLikePosZygo=[];
StdLikeNegZygo=[];
StdLikePosEventPosZygo=[];
StdikePosEventNegZygo=[];
StdLikeNegEventPosZygo=[];
StdLikeNegEventNegZygo=[];

MeanLikePosCorr=mean(AllLikePosCorr);
MeanLikeNegCorr=mean(AllLikeNegCorr);
MeanLikePosEventPosCorr=mean(AllLikePosEventPosCorr);
MeanLikePosEventNegCorr=mean(AllLikePosEventNegCorr);
MeanLikeNegEventPosCorr=mean(AllLikeNegEventPosCorr);
MeanLikeNegEventNegCorr=mean(AllLikeNegEventNegCorr);
StdLikePosCorr=std(AllLikePosCorr);
StdLikeNegCorr=std(AllLikeNegCorr);
StdLikePosEventPosCorr=std(AllLikePosEventPosCorr);
StdLikePosEventNegCorr=std(AllLikePosEventNegCorr);
StdLikeNegEventPosCorr=std(AllLikeNegEventPosCorr);
StdLikeNegEventNegCorr=std(AllLikeNegEventNegCorr);

MeanLikePosZygo=mean(AllLikePosZygo);
MeanLikeNegZygo=mean(AllLikeNegZygo);
MeanLikePosEventPosZygo=mean(AllLikePosEventPosZygo);
MeanLikePosEventNegZygo=mean(AllLikePosEventNegZygo);
MeanLikeNegEventPosZygo=mean(AllLikeNegEventPosZygo);
MeanLikeNegEventNegZygo=mean(AllLikeNegEventNegZygo);
StdLikePosZygo=std(AllLikePosZygo);
StdLikeNegZygo=std(AllLikeNegZygo);
StdLikePosEventPosZygo=std(AllLikePosEventPosZygo);
StdikePosEventNegZygo=std(AllLikePosEventNegZygo);
StdLikeNegEventPosZygo=std(AllLikeNegEventPosZygo);
StdLikeNegEventNegZygo=std(AllLikeNegEventNegZygo);


%% create figures
% set the line colors
set(0,'DefaultAxesColorOrder',[1 0 0;0 1 0],...
      'DefaultAxesLineStyleOrder','-|--|:')
x=-2.9:0.1:47;


plot(x,MeanLikeNegCorr(1:1:500),x,MeanLikePosCorr(1:1:500))
title('Corrugator Activation Likeability Exp 1','FontSize',12)
legend('dislikeable','likeable')
v=axis;
axis([-3 47 v(3) v(4)])
xlabel('time*1000 ms')
ylabel('activation as % of baseline')
annotation('line',[0.470 0.470],[0.10 0.92]);
annotation('line',[0.548 0.548],[0.10 0.92]);
annotation('line',[0.829 0.829],[0.10 0.92]);
saveas(gcf,[fdir 'LikeCorr.tiff']);
saveas(gcf,[fdir 'LikeCorr.jpg']);

plot(x,MeanLikePosEventNegCorr(1:1:500),x,MeanLikePosEventPosCorr(1:1:500),x,MeanLikeNegEventNegCorr(1:1:500),x,MeanLikeNegEventPosCorr(1:1:500))
title('Corrugator Activation Likeability Exp 1','FontSize',12)
legend('likeable-neg','likeable-pos','dislikeable-neg','dislikeable-pos')
v=axis;
axis([-3 47 v(3) v(4)])
xlabel('time*1000 ms')
ylabel('activation as % of baseline')
annotation('line',[0.470 0.470],[0.10 0.92]);
annotation('line',[0.548 0.548],[0.10 0.92]);
annotation('line',[0.829 0.829],[0.10 0.92]);
saveas(gcf,[fdir 'LikeEventCorr.tiff']);
saveas(gcf,[fdir 'LikeEventCorr.jpg']);

%% Zygomaticus plot
plot(x,MeanLikeNegZygo(30:1:500),x,MeanLikePosZygo(30:1:500))
title('Zygomaticus Activation Likeability','FontSize',12)
legend('dislikeable','likeable')
%axis tight
axis([0 47 90 220])
xlabel('time*1000 ms')
ylabel('activation as % of baseline')
annotation('line',[0.445 0.445],[0.10 0.92]);
annotation('line',[0.525 0.525],[0.10 0.92]);
annotation('line',[0.822 0.822],[0.10 0.92]);
saveas(gcf,[fdir 'LikeZygo.tiff']);
saveas(gcf,[fdir 'LikeZygo.jpg']);
plot(x,MeanLikePosEventNegZygo(30:1:500),x,MeanLikePosEventPosZygo(30:1:500),x,MeanLikeNegEventNegZygo(30:1:500),x,MeanLikeNegEventPosZygo(30:1:500))
title('Zygomaticus Activation Likeability','FontSize',12)
legend('likeable-neg','likeable-pos','dislikeable-neg','dislikeable-pos')
%axis tight
axis([0 47 90 220])
xlabel('time*1000 ms')
ylabel('activation as % of baseline')
annotation('line',[0.445 0.445],[0.10 0.92]);
annotation('line',[0.525 0.525],[0.10 0.92]);
annotation('line',[0.822 0.822],[0.10 0.92]);
saveas(gcf,[fdir 'LikeEventZygo.tiff']);
saveas(gcf,[fdir 'LikeEventZygo.jpg']);

