load('PRh axons.mat', 'Mouse16D')
%load('200116D_200129_trace.mat')
%%
fs = 30; %30.3017
preStim = 59; %frame %2 sec
postStim = 90; %frame %3 sec

% Before learning
pre_stimframe = [Mouse16D.pre10uA_0127.stimframe];
%pre_stimframe = CleanStimframe(pre_stimframe);
pre_trace = session1(1).normtrace;

prePSTCa = PSTCaTrace(pre_stimframe, pre_trace, preStim, postStim);

% After learning
post_stimframe = [Mouse16D.pre10uA_0128.stimframe];
%post_stimframe = CleanStimframe(post_stimframe);
post_trace = session2(1).normtrace;

postPSTCa = PSTCaTrace(post_stimframe, post_trace, preStim, postStim);


% Reward only
reward_stimframe = [Mouse16D.reward_0128_2.stimframe];
%reward_stimframe = CleanStimframe(reward_stimframe);

reward_trace = session2(4).normtrace;

rewardPSTCa = PSTCaTrace(reward_stimframe, reward_trace, preStim, postStim);


%%

edges = -2:1/fs:3;
figure
boundedline(edges(1:end-1), mean(prePSTCa,2), std(prePSTCa,[],2)/sqrt(size(prePSTCa,2)), 'k-', ...,
    edges(1:end-1), mean(postPSTCa,2), std(postPSTCa,[],2)/sqrt(size(postPSTCa,2)), 'r-')
xlim([-1 2])
legend('prelearning', 'postlearning')
ylabel('dF/F0 (%)')
xlabel('Time (s)')
title('Mouse 16D')


% Reward only
figure
boundedline(edges(1:end-1), mean(rewardPSTCa,2), std(rewardPSTCa,[],2)/sqrt(size(rewardPSTCa,2)), 'g-')
legend('Reward only')
ylabel('dF/F0 (%)')
xlabel('Time (s)')
title('Mouse 16D')
xlim([-1 2])



%%

Mouse16D.pre_trace = pre_trace;
Mouse16D.pre_stimframe = pre_stimframe;
Mouse16D.prePSTCa = prePSTCa;
Mouse16D.post_trace = post_trace;
Mouse16D.post_stimframe = post_stimframe;
Mouse16D.postPSTCa = postPSTCa;
Mouse16D.reward_trace = reward_trace;
Mouse16D.reward_stimframe = reward_stimframe;
Mouse16D.rewardPSTCa = rewardPSTCa;



%% Statistics
% Mean calcium amplitudes (between t = 0 and 1.5 s after stimulus) before and after learning were compared by Mann-Whitney U test.

[p h] = ranksum(mean(Mouse15C.prePSTCa(61:105,:)), mean(Mouse15C.postPSTCa(61:105,:)))

[p h] = ranksum(mean(Mouse16D.prePSTCa(61:105,:)), mean(Mouse16D.postPSTCa(61:105,:)))


[p h] = ranksum(mean(Mouse31F.prePSTCa(61:105,:)), mean(Mouse31F.postPSTCa(61:105,:)))



%% Correction fit
prePSTCa = Mouse15C.prePSTCa;
[fx, Corr_prePSTCa] = CorrectionFit(prePSTCa,edges);

postPSTCa = Mouse15C.postPSTCa;
[fx, Corr_postPSTCa] = CorrectionFit(postPSTCa,edges);


figure
boundedline(edges(1:end-1), mean(Corr_prePSTCa,2), std(Corr_prePSTCa,[],2)/sqrt(size(Corr_prePSTCa,2)), 'k-', ...,
    edges(1:end-1), mean(Corr_postPSTCa,2), std(Corr_postPSTCa,[],2)/sqrt(size(Corr_postPSTCa,2)), 'r-')
xlim([-1 2])



[p h] = ranksum(mean(Corr_prePSTCa(61:105,:)), mean(corr_postPSTCa(61:105,:)))

Mouse15C.Corr_prePSTCa = Corr_prePSTCa;
Mouse15C.Corr_postPSTCa = Corr_postPSTCa;

%% Pooled
prePSTCa = [Mouse15C.prePSTCa Mouse16D.prePSTCa Mouse31F.prePSTCa];
postPSTCa = [Mouse15C.postPSTCa Mouse16D.postPSTCa Mouse31F.postPSTCa];
[p h] = ranksum(mean(prePSTCa(61:105,:)), mean(postPSTCa(61:105,:)))

rewardPSTCa = [Mouse15C.rewardPSTCa Mouse16D.rewardPSTCa Mouse31F.rewardPSTCa];
[p h] = signrank(mean(rewardPSTCa(31:60,:)), mean(rewardPSTCa(61:105,:)))