%% instanteneous phase values
load('Phase_sample.mat')

e1_phase = phase_vals(1,:);
e2_phase = phase_vals(2,:);
n_phase  = size(phase_vals,2);

figure(2),clf
set(gcf, 'color','w', 'WindowState', 'maximize');
subplot(1,2,1)
% electrode 1 - all trials
polarplot([zeros(1,n_phase); e1_phase] , [zeros(1,n_phase);ones(1,n_phase)],'LineWidth',3,Color='k')
set(gca,'FontSize',20,'LineWidth',3)
title('Electrode 1')

subplot(1,2,2)
% electrode 2 - all trials
polarplot([zeros(1,n_phase); e2_phase] , [zeros(1,n_phase);ones(1,n_phase)],'LineWidth',3,Color='k')
set(gca,'FontSize',20,'LineWidth',3)
title('Electrode 2')

%% Single trial

figure(3),clf
set(gcf, 'color','w', 'WindowState', 'maximize');
subplot(1,2,1)
% electrode 1
polarplot([zeros(1,1); e1_phase(1,100)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k')
set(gca,'FontSize',20,'LineWidth',3)
title('Electrode 1')

subplot(1,2,2)
% electrode 2
polarplot([zeros(1,1); e2_phase(1,100)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k')
set(gca,'FontSize',20,'LineWidth',3)
title('Electrode 2')

%% Difference

diff_phase = e1_phase - e2_phase;

figure(4),clf
set(gcf, 'color','w', 'WindowState', 'maximize');

subplot(1,2,1)
polarplot([zeros(1,1); e2_phase(1,105)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k'), hold on
polarplot([zeros(1,1); e1_phase(1,105)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k')
polarplot([zeros(1,1); diff_phase(1,105)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='m')

set(gca,'FontSize',20,'LineWidth',3)
title('Phase difference: E1(t1) - E2 (t1) ) ')

subplot(1,2,2)
% electrodes - all trials
polarplot([zeros(1,20); diff_phase(1,21:40)] , [zeros(1,20);ones(1,20)],'LineWidth',3,Color='g')
set(gca,'FontSize',20,'LineWidth',3)
title('Phase difference: 20 trials')

%% Average

% Regular average:
avg_phase = angle((e1_phase(1,2) + e2_phase(1,2))/2);
% Average over complex numbers:
cm_avg_phase = angle( exp(1i*e1_phase(1,2)) + exp(1i*e2_phase(1,2)) )/2;

figure(6),clf
set(gcf, 'color','w', 'WindowState', 'maximize');

subplot(1,2,1)
polarplot([zeros(1,1); e2_phase(1,2)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k'), hold on
polarplot([zeros(1,1); e1_phase(1,2)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k')
polarplot([zeros(1,1); avg_phase] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='g')

set(gca,'FontSize',20,'LineWidth',3)
title('Regular average ')

subplot(1,2,2)
polarplot([zeros(1,1); e2_phase(1,2)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k'), hold on
polarplot([zeros(1,1); e1_phase(1,2)] , [zeros(1,1);ones(1,1)],'LineWidth',3,Color='k')
polarplot([zeros(1,1); cm_avg_phase], [zeros(1,1);ones(1,1)],'LineWidth',3,Color='m')

set(gca,'FontSize',20,'LineWidth',3)
title('Complex average ')

%% ISPC


load('Phase_sample.mat')
phase_vals = double(phase_vals);

% Extract the phase values:
e1_phase = phase_vals(1,:);
e2_phase = phase_vals(2,:);
n_phase  = size(phase_vals,2);

% Difference between e1 and e2:
diff_phase = e1_phase - e2_phase;
% Projecting the diffrence in to complex plain 
cmx_phase = exp(1i*diff_phase);
% Average in complex space
mean_cm_phase = mean(cmx_phase,2);
% Vector length (ISPC)
ispc = abs(mean_cm_phase);


figure(1),clf
set(gcf, 'color','w', 'WindowState', 'maximize');

polarplot([zeros(1,n_phase); angle(cmx_phase)] , [zeros(1,n_phase);ones(1,n_phase)],'LineWidth',1,Color='k'), hold on
polarplot([zeros(1,1); angle(mean_cm_phase)] , [zeros(1,1);ispc(1,1)],'LineWidth',5,Color='g')

set(gca,'FontSize',20,'LineWidth',3)
title('ISPC')


%% IPSC in one line of code

ispc = abs(mean(exp(1i*(e1_phase - e2_phase)),2));

figure(2),clf
set(gcf, 'color','w', 'WindowState', 'maximize');

polarplot([zeros(1,n_phase); angle(cmx_phase)] , [zeros(1,n_phase);ones(1,n_phase)],'LineWidth',1,Color='k'), hold on
polarplot([zeros(1,1); angle(mean_cm_phase)] , [zeros(1,1);ispc(1,1)],'LineWidth',5,Color='g')

set(gca,'FontSize',20,'LineWidth',3)
title('ISPC - one line')


%% ISPC over time and frequency

load('TF_phase_sample.mat')
phase_vals = double(phase_vals);

% Extract the phase values:
e1_phase = squeeze(phase_vals(1,:,:,:));
e2_phase = squeeze(phase_vals(2,:,:,:));

% ISPC
tf_ispc = abs(mean(exp(1i*(e1_phase - e2_phase)),3));

figure(3),clf
set(gcf, 'color','w', 'WindowState', 'maximize');
    % indicies and lines
        contourf(timevec,frex,tf_ispc,128,'linecolor','none','fill','on'),hold on    
        plot ([0 0], [frex(1) frex(end)],'LineStyle', '--','color','k','LineWidth',1.5);
        plot ([-1000 -1000], [frex(1) frex(end)],'LineStyle', '--','color','k','LineWidth',1.5);  
    % axis settings
        set(gca,'ylim',ylim,'xlim',xlim,'clim',[.15 .5],'xtick',[-1000,-500,0,500,800],'FontSize',35,...
        'LineWidth',3,'xticklabel',{'C',-500,'T',500,[]})   
        axis square; box off;colormap (hot)


%%

%%EEG_tmp.data = laplacian_perrinX(EEG_tmp.data,[EEG_tmp.chanlocs.X],[EEG_tmp.chanlocs.Y],[EEG_tmp.chanlocs.Z]);