clear
clc
close all
subject_options;
do_localizer=0;
subid=[1:26];


for sb=subid
    if exist([matpath,'gmf_',subjname{sb},'_ICA_lp30.mat'],'file')
        if do_localizer
            load([matpath,'Other/gmf_',subjname{sb},'.mat']);
        else
            load([matpath,'gmf_',subjname{sb},'_ICA_lp30.mat']);
        end
        
        if strcmp(subjname{sb},'sub14')||strcmp(subjname{sb},'sub24')
            id(sb)=0;
        else
            id(sb)=1;
        end
        
        exp(sb,:)=gmf_exp.avg;
        unexp(sb,:)=gmf_unexp.avg;
        
    end
end
%%
id = logical(id);
%%
for sb = subid
    clear Xhat X
    load(['/home/predatt/biahan/Project/MEG/mat/',subjname{sb},'_decode_exp_20ms_lp.mat'])
    m0_exp(:,:,sb) = mean((Xhat(:,:,X==0)),3);
    m1_exp(:,:,sb) = mean((Xhat(:,:,X==1)),3);
    
    [h,p,subid,t_exp]=ttest2(Xhat(:,:,X==0),Xhat(:,:,X==1),'dim',3);
    diff_exp(:,:,sb) = squeeze(t_exp.tstat);
    
    clear Xhat X
    load(['/home/predatt/biahan/Project/MEG/mat/',subjname{sb},'_decode_unexp_20ms_lp.mat'])
    m0_unexp(:,:,sb) = mean((Xhat(:,:,X==0)),3);
    m1_unexp(:,:,sb) = mean((Xhat(:,:,X==1)),3);
    [h,p,subid,t_unexp]=ttest2(Xhat(:,:,X==0),Xhat(:,:,X==1),'dim',3);
    diff_unexp(:,:,sb) = squeeze(t_unexp.tstat);
end
%%


exp_diff = dig_func(m1_exp)-dig_func(m0_exp);
unexp_diff = dig_func(m1_unexp)-dig_func(m0_unexp);


mean_exp = mean(exp_diff(71:111,:),1)*100;
mean_unexp = mean(unexp_diff(71:111,:),1)*100;


b=largest_cluster(diff_exp(:,:,id),diff_unexp(:,:,id),time(T));


b=logical(b);
dd=repmat(b,1,1,26);

exp_mean = m0_exp(dd)-m1_exp(dd);
unexep_mean = m0_unexp(dd)-m1_unexp(dd);

exp_mean = reshape(exp_mean,[],26);
unexep_mean = reshape(unexep_mean,[],26);

exp_mean_mean = mean(exp_mean,1);
unexp_mean_mean = mean(unexep_mean,1);

erf_unexp = mean(unexp(:,1261:1500),2);
erf_exp = mean(exp(:,1261:1500),2);


decodingdiff_unexp_exp= unexp_mean_mean(id)-exp_mean_mean(id);
erfdiff_unexp_exp = erf_unexp(id)-erf_exp(id);


figure
plot(-decodingdiff_unexp_exp*100,erfdiff_unexp_exp,'o')
lsline
ylabel('Field Strength Difference (T)')
xlabel('Decoding Signal Difference (a.u.)')


[r,p] = corr(-decodingdiff_unexp_exp'.*100,erfdiff_unexp_exp);


set(findall(gcf,'-property','FontSize'),'FontSize',16)
set(gcf, 'Units', 'centimeters', 'Position', [0, 0, 18, 10], 'PaperUnits', 'centimeters', 'PaperSize', [20, 15])
box off
pbaspect([1.8 1.2 1]);