clear; clc;
warning('off');

sesslist = {'LAMP3_004_190612_151313';
            'LAMP3_006_190613_194519';
            'LAMP3_011_190617_142002';
            'LAMP3_015_190618_172714'};

iplist = [1,2,3,4]; outsr = 100;
datapath = 'C:\Users\Liang\Dropbox\AchProject\Fig2 data\Ripple_Ach\';

%% load scatter data
plotdata = zeros(0);
for ip = iplist
    session = sesslist{ip,1};
    load([datapath,session,filesep,'scatdata_ripple.mat'],'scatdata_ripple');
    plotdata = [plotdata; scatdata_ripple];
%     load([datapath,session,filesep,'scatdata.mat'],'scatdata');
%     scatdata(scatdata(:,2)<pow2db(0.4e4) | scatdata(:,2)>pow2db(4e4),:) = [];
%     plotdata = [plotdata; (scatdata-nanmean(scatdata))./nanstd(scatdata)];
%     plotdata = [plotdata; scatdata];    
    
end

figure;
subplot(2,1,1)
scatter(plotdata(:,1)*100,plotdata(:,2),16,'filled');
xlabel('dFoverF (%)'); ylabel('Ripple rate (Hz)')
title('dFoverF vs. Ripple rate');
xlim([-15,15]); ylim([0.01,1]);
[r,p] = corrcoef(plotdata);
legend(sprintf('r = %.3f, p = %.2e',r(1,2),p(1,2)))


subplot(2,1,2)
ca2 = [128 193 219]./255;    % Blue theme
cl2 = [ 52 148 186]./255;
achedge = linspace(-15,15,7);
achc = achedge(1:end-1) + 0.5*diff(achedge(1:2));
[N,~,binX] = histcounts(plotdata(:,1)*100,achedge);
histmean = zeros(0);
for ib = 1:length(achedge)-1
    histmean{ib,1} = plotdata(binX == ib,2);
    iosr.statistics.boxPlot(ib,histmean{ib,1},'theme','colorall',...
        'showScatter',true,'scatterMarker','.','symbolMarker','none',...
        'scatterColor',cl2 ,'scatterSize', 64, 'scatterAlpha', 0.5,...
        'boxWidth',0.4,'boxColor',ca2,'showViolin',false); hold on;
end
xlim([1,length(achedge)-1]+0.25*0.5*diff(achedge(1:2))*[-1,1]);
xticks(1:length(achedge)-1); xticklabels(achc);
xlabel('dFoverF (%)'); ylabel('Ripple rate (Hz)')
title('dFoverF vs. Ripple rate');
% set(gca,'YScale','log'); ylim([0.01,1])


%% load scatter data
plotdata = zeros(0);
for ip = iplist
    session = sesslist{ip,1};
    load([datapath,session,filesep,'ripple_peri_event.mat'],'ripple_peri_event');
    plotdata = [plotdata; ripple_peri_event];
    clear scatdata;
end

figure;
imagesc([1:size(plotdata,2)./outsr]'-30,[1:size(plotdata,1)]',...
        zscore(plotdata,[],2));
clim([-3,3]); cb = colorbar();
cb.Label.String = 'z-scored dFoverF';
xlabel('Times from ripple peaks (s)'); ylabel('Ripple events #');

figure
options.color_area = [128 193 219]./255;    % Blue theme
options.color_line = [ 52 148 186]./255;
options.alpha      = 0.5;
options.line_width = 2;
options.error      = 'sem';
options.x_axis     = [1:size(plotdata,2)]./outsr-30;
h2 = plot_areaerrorbar(plotdata*100,options);
xlabel('Times from ripple peaks (s)'); ylabel('dFoverF (%)');
