function figure_neuropriming(fld)

fprintf('\n======================================================\n');
fprintf('-- Creating Figure neuropriming --\n');
fprintf('======================================================\n');

%% settings
yellow = [247 148 29]./255;
blue = [27 117 188]./255;
cols = [blue; yellow];
    
snrthres = 2.5;
smoothfact = 15;
mindays = 3;
falpha = 0.4; % set to 1 for converting to illustrator

%% plots
% load traces
monkeys = {'M1','M2'};
wm = [];
f1 = figure; set(gcf,'Position',[200   500   1200   600]);
f2 = figure; set(gcf,'Position',[200   500   1200   600]);
f1b = figure; set(gcf,'Position',[100 100  600 1200]);
f2b = figure; set(gcf,'Position',[100 100  600 1200]);
for mi = 1:2
    monkey = monkeys{mi};
    savedir = fullfile(fld.basedir, 'results','figure_neuropriming');
    load(fullfile(savedir, [monkey '_averages_snr' num2str(snrthres) '_mindays' num2str(mindays) '.mat']),...
        'targmod','distmod','targmag','distmag','tracesLUT','env_t');

    % smooth first
    q = []; q2 = [];
    for i = 1:size(targmod)
        q(i,:) = smooth(targmod(i,:),smoothfact);
        q2(i,:) = smooth(distmod(i,:),smoothfact);
    end
    targmod = q;
    distmod = q2;

    %% plot the effect of shape change on target modulation
    figure(f1);  subplot(2,3,mi+1);
    [D, q] = make_plot(env_t,targmod,tracesLUT.shapeswitch,'Shape');
    DShape{mi,1}=D;
    DShape_q{mi,1}=q;
    ylabel('target modulation');
    title(monkey);
    disp(['---- ' monkey ' | shape ----']);
    [~,p,~,st] = ttest(D(:,1),D(:,2));
    disp([monkey ', target vs non-target, p=' num2str(p) ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);
    
    figure(f1b); subplot(2,3,mi+1); hold on;
    bar([1 2],[mean(D(:,1)) mean(D(:,2))]);
    errorbar([1 2],[mean(D(:,1)) mean(D(:,2))],...
        [std(D(:,1))./sqrt(size(D,1)) std(D(:,2))./sqrt(size(D,1))],...
        'LineStyle','none','Color','k')
    x = ones(size(D,1),1); x = x + (randi(10,length(x),1)-5)/30;
    cx = [];
    for l=1:2; cx(:,l) = x + 1*(l-1); end
    incr = D(:,2)>D(:,1);
    plot(cx(incr,:)',D(incr,:)','Color',[.4 .4 .4]); hold all;
    plot(cx(~incr,:)',D(~incr,:)','Color',[.8 .8 .8]);
    for l=1:size(D,2)
        scatter(cx(incr,l),D(incr,l),20,[.4 .4 .4],'f');
        scatter(cx(~incr,l),D(~incr,l),20,[.8 .8 .8],'f');
    end
    title(monkey);
    set(gca,'ylim',[-0.1 0.2]);

    figure(f1); subplot(2,3,mi+4); 
    [D, q] = make_plot(env_t,distmod,tracesLUT.shapeswitch,'Shape');
    DShape{mi,2}=D; 
    DShape_q{mi,2}=q;
    ylabel('distractor modulation');
    title(monkey);
    [~,p,~,st] = ttest(D(:,1),D(:,2));
    disp([monkey ', distractor vs non-target, p=' num2str(p) ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);
    
    figure(f1b); subplot(2,3,mi+4); hold on;
    bar([1 2],[mean(D(:,1)) mean(D(:,2))]);
    errorbar([1 2],[mean(D(:,1)) mean(D(:,2))],...
        [std(D(:,1))./sqrt(size(D,1)) std(D(:,2))./sqrt(size(D,1))],...
        'LineStyle','none','Color','k')
    x = ones(size(D,1),1); x = x + (randi(10,length(x),1)-5)/30;
    cx = [];
    for l=1:2; cx(:,l) = x + 1*(l-1); end
    incr = D(:,2)>D(:,1);
    plot(cx(incr,:)',D(incr,:)','Color',[.8 .8 .8]); hold all;
    plot(cx(~incr,:)',D(~incr,:)','Color',[.4 .4 .4]);
    for l=1:size(D,2)
        scatter(cx(incr,l),D(incr,l),20,[.8 .8 .8],'f');
        scatter(cx(~incr,l),D(~incr,l),20,[.4 .4 .4],'f');
    end
    title(monkey);
    set(gca,'ylim',[-0.1 0.05]);


    %% do the same on color
    figure(f2); subplot(2,3,mi+1);
    [D, q] = make_plot(env_t,targmod,tracesLUT.colorswitch,'Color');
    DColor{mi,1}=D;
    DColor_q{mi,1}=q;
    ylabel('target modulation');
    title(monkey);
    disp(['---- ' monkey ' | color ----']);
    [~,p,~,st] = ttest(D(:,1),D(:,2));
    disp([monkey ', target vs non-target, p=' num2str(p) ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);
    
    figure(f2b); subplot(2,3,mi+1); hold on;
    bar([1 2],[mean(D(:,1)) mean(D(:,2))]);
    errorbar([1 2],[mean(D(:,1)) mean(D(:,2))],...
        [std(D(:,1))./sqrt(size(D,1)) std(D(:,2))./sqrt(size(D,1))],...
        'LineStyle','none','Color','k')
    x = ones(size(D,1),1); x = x + (randi(10,length(x),1)-5)/30;
    cx = [];
    for l=1:2; cx(:,l) = x + 1*(l-1); end
    incr = D(:,2)>D(:,1);
    plot(cx(incr,:)',D(incr,:)','Color',[.4 .4 .4]); hold all;
    plot(cx(~incr,:)',D(~incr,:)','Color',[.8 .8 .8]);
    for l=1:size(D,2)
        scatter(cx(incr,l),D(incr,l),20,[.4 .4 .4],'f');
        scatter(cx(~incr,l),D(~incr,l),20,[.8 .8 .8],'f');
    end
    title(monkey);
    set(gca,'ylim',[-0.1 0.2]);
    
    figure(f2); subplot(2,3,mi+4);
    [D, q] = make_plot(env_t,distmod,tracesLUT.colorswitch,'Color');
    DColor{mi,2}=D;
    DColor_q{mi,2}=q;
    ylabel('distractor modulation');
    title(monkey);
    [~,p,~,st] = ttest(D(:,1),D(:,2));
    disp([monkey ', distractor vs non-target, p=' num2str(p) ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);
    
    figure(f2b); subplot(2,3,mi+4); hold on;
    bar([1 2],[mean(D(:,1)) mean(D(:,2))]);
    errorbar([1 2],[mean(D(:,1)) mean(D(:,2))],...
        [std(D(:,1))./sqrt(size(D,1)) std(D(:,2))./sqrt(size(D,1))],...
        'LineStyle','none','Color','k')
    x = ones(size(D,1),1); x = x + (randi(10,length(x),1)-5)/30;
    cx = [];
    for l=1:2; cx(:,l) = x + 1*(l-1); end
    incr = D(:,2)>D(:,1);
    plot(cx(incr,:)',D(incr,:)','Color',[.8 .8 .8]); hold all;
    plot(cx(~incr,:)',D(~incr,:)','Color',[.4 .4 .4]);
    for l=1:size(D,2)
        scatter(cx(incr,l),D(incr,l),20,[.8 .8 .8],'f');
        scatter(cx(~incr,l),D(~incr,l),20,[.4 .4 .4],'f');
    end
    title(monkey);
    set(gca,'ylim',[-0.1 0.05]);
end

%% Now do the same thing on the pooled data
%% plot the effect of shape change on target modulation
%f0 = figure; set(gcf,'Position',[216   492   411   606]);
figure(f1);  subplot(2,3,1);
DD = [DShape{1,1};DShape{2,1}];
q = [DShape_q{1,1};DShape_q{2,1}];
make_plot_pooled(q, env_t, 'Shape')
ylabel('target modulation');
title('BOTH monkeys');
disp('---- BOTH | shape ----');
[~,p,~,st] = ttest(DD(:,1),DD(:,2));
disp(['BOTH, target vs non-target, p=' num2str(p) ...
    ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);
suptitle('Effect of Shape-change on neural activity');

%f0b = figure; set(gcf,'Position',[100 100 250 500]);
figure(f1b); subplot(2,3,1); hold on;
bar([1 2],[mean(DD(:,1)) mean(DD(:,2))]);
errorbar([1 2],[mean(DD(:,1)) mean(DD(:,2))],...
    [std(DD(:,1))./sqrt(size(DD,1)) std(DD(:,2))./sqrt(size(DD,1))],...
    'LineStyle','none','Color','k')
x = ones(size(DD,1),1); x = x + (randi(10,length(x),1)-5)/30;
cx = [];
for l=1:2; cx(:,l) = x + 1*(l-1); end
incr = DD(:,2)>DD(:,1);
plot(cx(incr,:)',DD(incr,:)','Color',[.4 .4 .4]); hold all;
plot(cx(~incr,:)',DD(~incr,:)','Color',[.8 .8 .8]);
for l=1:size(DD,2)
    scatter(cx(incr,l),DD(incr,l),20,[.4 .4 .4],'f');
    scatter(cx(~incr,l),DD(~incr,l),20,[.8 .8 .8],'f');
end
title('BOTH');
set(gca,'ylim',[-0.1 0.2]);
suptitle('Effect of Shape-change on neural activity');

figure(f1); subplot(2,3,4);
DD = [DShape{1,2};DShape{2,2}];
q = [DShape_q{1,2};DShape_q{2,2}];
make_plot_pooled(q, env_t, 'Shape')
ylabel('distractor modulation');
title('BOTH monkeys');
[~,p,~,st] = ttest(DD(:,1),DD(:,2));
disp(['BOTH, distractor vs non-target, p=' num2str(p) ...
    ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);

figure(f1b); subplot(2,3,4); hold on;
bar([1 2],[mean(DD(:,1)) mean(DD(:,2))]);
errorbar([1 2],[mean(DD(:,1)) mean(DD(:,2))],...
    [std(DD(:,1))./sqrt(size(DD,1)) std(DD(:,2))./sqrt(size(DD,1))],...
    'LineStyle','none','Color','k')
x = ones(size(DD,1),1); x = x + (randi(10,length(x),1)-5)/30;
cx = [];
for l=1:2; cx(:,l) = x + 1*(l-1); end
incr = DD(:,2)>DD(:,1);
plot(cx(incr,:)',DD(incr,:)','Color',[.8 .8 .8]); hold all;
plot(cx(~incr,:)',DD(~incr,:)','Color',[.4 .4 .4]);
for l=1:size(DD,2)
    scatter(cx(incr,l),DD(incr,l),20,[.8 .8 .8],'f');
    scatter(cx(~incr,l),DD(~incr,l),20,[.4 .4 .4],'f');
end
title('BOTH');
set(gca,'ylim',[-0.1 0.05]);


%% do the same on color
%f00 = figure; set(gcf,'Position',[216   492   411   606]);
figure(f2);  subplot(2,3,1);
DD = [DColor{1,1};DColor{2,1}];
q = [DColor_q{1,1};DColor_q{2,1}];
make_plot_pooled(q, env_t, 'Color')
ylabel('target modulation');
title('BOTH');
disp('---- BOTH | color ----');
[~,p,~,st] = ttest(DD(:,1),DD(:,2));
disp(['BOTH, target vs non-target, p=' num2str(p) ...
    ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);
suptitle('Effect of Color-change on neural activity');

%f00b = figure; set(gcf,'Position',[100 100 250 500]);
figure(f2b); subplot(2,3,1); hold on;
bar([1 2],[mean(DD(:,1)) mean(DD(:,2))]);
errorbar([1 2],[mean(DD(:,1)) mean(DD(:,2))],...
    [std(DD(:,1))./sqrt(size(DD,1)) std(DD(:,2))./sqrt(size(DD,1))],...
    'LineStyle','none','Color','k')
x = ones(size(DD,1),1); x = x + (randi(10,length(x),1)-5)/30;
cx = [];
for l=1:2; cx(:,l) = x + 1*(l-1); end
incr = DD(:,2)>DD(:,1);
plot(cx(incr,:)',DD(incr,:)','Color',[.4 .4 .4]); hold all;
plot(cx(~incr,:)',DD(~incr,:)','Color',[.8 .8 .8]);
for l=1:size(DD,2)
    scatter(cx(incr,l),DD(incr,l),20,[.4 .4 .4],'f');
    scatter(cx(~incr,l),DD(~incr,l),20,[.8 .8 .8],'f');
end
title('BOTH');
set(gca,'ylim',[-0.1 0.2]);
suptitle('Effect of Color-change on neural activity');

figure(f2); subplot(2,3,4);
DD = [DColor{1,2};DColor{2,2}];
q = [DColor_q{1,2};DColor_q{2,2}];
make_plot_pooled(q, env_t, 'Color')
ylabel('distractor modulation');
title('BOTH');
[~,p,~,st] = ttest(DD(:,1),DD(:,2));
disp(['BOTH, distractor vs non-target, p=' num2str(p) ...
    ', t=' num2str(st.tstat) ', df=' num2str(st.df)]);

figure(f2b); subplot(2,3,4); hold on;
bar([1 2],[mean(DD(:,1)) mean(DD(:,2))]);
errorbar([1 2],[mean(DD(:,1)) mean(DD(:,2))],...
    [std(DD(:,1))./sqrt(size(DD,1)) std(DD(:,2))./sqrt(size(DD,1))],...
    'LineStyle','none','Color','k')
x = ones(size(DD,1),1); x = x + (randi(10,length(x),1)-5)/30;
cx = [];
for l=1:2; cx(:,l) = x + 1*(l-1); end
incr = DD(:,2)>DD(:,1);
plot(cx(incr,:)',DD(incr,:)','Color',[.8 .8 .8]); hold all;
plot(cx(~incr,:)',DD(~incr,:)','Color',[.4 .4 .4]);
for l=1:size(DD,2)
    scatter(cx(incr,l),DD(incr,l),20,[.8 .8 .8],'f');
    scatter(cx(~incr,l),DD(~incr,l),20,[.4 .4 .4],'f');
end
title('BOTH');
set(gca,'ylim',[-0.1 0.05]);

end


function [D , q] = make_plot(env_t,data,selector,featurename)
    yellow = [247 148 29]./255;
    blue = [27 117 188]./255;
    cols = [blue; yellow];
    
    falpha = 0.4;

    hold off
    p = [];
    win_idx = find(env_t>0.15 & env_t<0.20);
    D = [];
    for featureswitch = [0 1]
        incl = selector==featureswitch;
        q{featureswitch+1} = data(incl,:);
        mn = mean(q{featureswitch+1});
        sem = std(q{featureswitch+1})./sqrt(sum(incl));
        
        % also prepare for some stats
        D(:,featureswitch+1) = mean(data(incl,win_idx),2);

        t2 = [env_t, fliplr(env_t)];
        inBetween = [mn+sem, fliplr(mn-sem)];
        fill(t2*1e3, inBetween, 'g','LineStyle','none','FaceColor',...
            cols(featureswitch+1,:),'FaceAlpha',falpha); hold all
        p(end+1) = plot(env_t*1e3,mn,'Color',cols(featureswitch+1,:));
    end
    legend(p,{[featurename '-same'],[featurename '-swap']},'Location','southwest');
    xlabel('time(ms)');
    plot([env_t(1)*1e3 env_t(end)*1e3],[0 0],'k-')
    set(gca,'xlim',[-50 250]);
end

function make_plot_pooled(q, env_t, featurename)
    yellow = [247 148 29]./255;
    blue = [27 117 188]./255;
    cols = [blue; yellow];
    falpha = 0.4;

    hold off
    p = [];
    for featureswitch = [0 1]
        mn = mean(q{featureswitch+1});
        sem = std(q{featureswitch+1})./sqrt(sum(size(q{featureswitch+1},1)));
        t2 = [env_t, fliplr(env_t)];
        inBetween = [mn+sem, fliplr(mn-sem)];
        fill(t2*1e3, inBetween, 'g','LineStyle','none','FaceColor',...
            cols(featureswitch+1,:),'FaceAlpha',falpha); hold all
        p(end+1) = plot(env_t*1e3,mn,'Color',cols(featureswitch+1,:));
    end
    legend(p,{[featurename '-same'],[featurename '-swap']},'Location','southwest');
    xlabel('time(s)');
    plot([env_t(1)*1e3 env_t(end)*1e3],[0 0],'k-')
    set(gca,'xlim',[-50 250]);
end