%
% cat all maps
% longer rip window
% larger spat bins
clearvars;
load('day5.PCs5cm.mat');
load('day5.WmazeTrials.mat');
load('day5.spkSWRTask.mat')
load('day5.cellMetrics.mat')
load('day5.spkCLpulses.mat')

%%
for unit = 1:size(ripSpkTask.timesRelat,1)
    for rip = 1:size(ripSpkTask.timesRelat,2)
        %edges = 0:0.01:ripSpkTask.ripDuration(rip);
        edges = -0.05:0.01:ripSpkTask.ripDuration(rip)+0.05;
        %ripSpkTask.timesBin{unit,rip} = histc(ripSpkTask.timesRelat{unit,rip},edges);
        ripSpkTask.timesBinLonger{unit,rip} = histc(ripSpkTask.timesLonger{unit,rip},edges);
    end
end
for rip = 1:size(ripSpkTask.timesRelat,2)
        edges = -0.05:0.01:ripSpkTask.ripDuration(rip)+0.05;
    for unit = 1:size(ripSpkTask.timesRelat,1)
        if ~isempty(ripSpkTask.timesBinLonger{unit,rip})
           ripSpkTask.timesBinLongerCat{rip}(unit,:) = ripSpkTask.timesBinLonger{unit,rip};
        else
           ripSpkTask.timesBinLongerCat{rip}(unit,:) = zeros(1,length(edges));
        end
    end
end
for rip = 1:size(ripSpkTask.timesRelat,2)
    spikesperbin{rip,1} = ripSpkTask.timesBinLongerCat{rip};
end 

clear spikesperbinConv
for rip = 1:size(ripSpkTask.timesRelat,2)
    for unit = 1:size(ripSpkTask.timesRelat,1)
        spikesperbinConv{rip,1}(unit,:) = conv(spikesperbin{rip}(unit,:),ones(2,1),'same');
    end
end
% cat all rate maps
for i = 1%:2
    for unit = 1:size(ripSpkTask.timesRelat,1)
        placefieldmaps{1}(unit,:) = cat(2,curve{unit}{1}.rate,curve{unit}{2}.rate,curve{unit}{3}.rate,curve{unit}{4}.rate);
%         placefieldmaps{1}(unit,:) = cat(2,curve{unit}{1}.rate,curve{unit}{2}.rate);
%         placefieldmaps{2}(unit,:) = cat(2,curve{unit}{3}.rate,curve{unit}{4}.rate);        
       %[val, maxT(unit)] = max(placefieldmaps(unit,:));
    end
end
%sort place fields
[a,b]=sort(maxT);
placefieldmaps_sorted=placefieldmaps(b,:);
figure;
imagesc(placefieldmaps_sorted);

tbin= 0.02; Xsyze=1; Ysyze=size(placefieldmaps{1},2);clims=[0 0.05];
for rip = 1:size(ripSpkTask.timesRelat,2)
    if ripSpkTask.nCellsEvent(rip) > 0.2
     [e r] = size(spikesperbinConv{rip});
    for i = 1:length(placefieldmaps)
    logpfmaps{i} = log(placefieldmaps{i}*tbin);
    logpfmaps{i}(find(logpfmaps{i} == -Inf)) = -450;

    %logpfmaps(find(logpfmaps == -Inf)) = 0;
    sumpfmaps{i} = tbin*sum(placefieldmaps{i});
    firstlook{i} = spikesperbinConv{rip}' * logpfmaps{i}-repmat(sumpfmaps{i},r,1);
    ppost{i} = exp(firstlook{i});
    ppost_norm{i} = bsxfun(@times, ppost{i}, 1./sum(ppost{i},2));
    firstlooktranslate{i}{rip} = reshape(ppost_norm{i},r,Xsyze,Ysyze);
    end
    titles = {'allMaps'};%{'inbound';'outbound'};

    figure;
    for i = 1:length(placefieldmaps)
    %subplot(2,1,i)
    imagesc(squeeze(firstlooktranslate{i}{rip}),clims);hold on;title(titles{i});
    ylabel('rip time');xlabel('position');
     plot([25.6 25.6],ylim,'k');hold on;plot([38.4 38.4],ylim,'k');hold on;
     plot([89.6 89.6],ylim,'k');hold on;plot([102.4 102.4],ylim,'k');hold on;
    end
    plot([ size(placefieldmaps{1},2)/4  size(placefieldmaps{1},2)/4],ylim,'r');hold on;
    plot([ size(placefieldmaps{1},2)/2  size(placefieldmaps{1},2)/2],ylim,'r');hold on;
    plot([ size(placefieldmaps{1},2)*(3/4)  size(placefieldmaps{1},2)*(3/4)],ylim,'r');hold on;
     plot([153.6 153.6],ylim,'k');hold on;plot([166.4 166.4],ylim,'k');hold on;
     plot([217.6 217.6],ylim,'k');hold on;plot([230.4 230.4],ylim,'k');hold on;    
    %set(gcf, 'Units', 'Normalized', 'OuterPosition', [0.1, 0.1, 0.7, 0.5]);
    saveas(gcf,['decoding/decodingSWR_rip' num2str(rip) '_4maps'],'png')%pause;
    close all;
    end
end

%% pulses CL
%pc=find(cellMetrics.isPC==1); 
for unit = 1:size(pulSpkCL.timesRelat,1)
    for rip = 1:size(pulSpkCL.timesRelat,2)
        %edges = 0:0.01:ripSpkTask.ripDuration(rip);
        edges = -0.11:0.01:0.1;
        %ripSpkTask.timesBin{unit,rip} = histc(ripSpkTask.timesRelat{unit,rip},edges);
        pulSpkCL.timesBinLonger{unit,rip} = histc(pulSpkCL.timesLonger{unit,rip},edges);
    end
end
for rip = 1:size(pulSpkCL.timesRelat,2)
        edges = -0.11:0.01:0.1;
        for unit = 1:size(pulSpkCL.timesRelat,1)
            if ~isempty(pulSpkCL.timesBinLonger{unit,rip})
               pulSpkCL.timesBinLongerCat{rip}(unit,:) = pulSpkCL.timesBinLonger{unit,rip};
            else
               pulSpkCL.timesBinLongerCat{rip}(unit,:) = zeros(1,length(edges));
            end
        end
end
for rip = 1:size(pulSpkCL.timesRelat,2)
    spikesperbinCL{rip,1} = pulSpkCL.timesBinLongerCat{rip};
end 

clear spikesperbinCLConv placefieldmaps
for rip = 1:size(pulSpkCL.timesRelat,2)
    for unit = 1:size(pulSpkCL.timesRelat,1)
        spikesperbinCLConv{rip,1}(unit,:) = conv(spikesperbinCL{rip}(unit,:),ones(2,1),'same');
    end
end
% cat all rate maps
for i = 1%:2
    for unit = 1:size(pulSpkCL.timesRelat,1)
        placefieldmaps{i}(unit,:) = cat(2,curve{unit}{1}.rate,curve{unit}{2}.rate,curve{unit}{3}.rate,curve{unit}{4}.rate);
        %placefieldmaps{1}(unit,:) = cat(2,curve{unit}{1}.rate,curve{unit}{2}.rate);
        %placefieldmaps{2}(unit,:) = cat(2,curve{unit}{3}.rate,curve{unit}{4}.rate);
       %[val, maxT(unit)] = max(placefieldmaps(unit,:));
    end
end
%sort place fields
[a,b]=sort(maxT);
placefieldmaps_sorted=placefieldmaps(b,:);
figure;
imagesc(placefieldmaps_sorted);

tbin= 0.02; Xsyze=1; Ysyze=size(placefieldmaps{1},2);clims=[0 0.05];
for rip = 1:size(pulSpkCL.timesRelat,2)
    if pul2SpkCL.nCellsEvent1(rip) > 0.05 && pul2SpkCL.nCellsEvent2(rip) > 0.12
    [e r] = size(spikesperbinCLConv{rip});
    for i = 1:length(placefieldmaps)
    logpfmaps{i} = log(placefieldmaps{i}*tbin);
    logpfmaps{i}(find(logpfmaps{i} == -Inf)) = -450;

    %logpfmaps(find(logpfmaps == -Inf)) = 0;
    sumpfmaps{i} = tbin*sum(placefieldmaps{i});
    firstlook{i} = spikesperbinCLConv{rip}' * logpfmaps{i}-repmat(sumpfmaps{i},r,1);
    ppost{i} = exp(firstlook{i});
    ppost_norm{i} = bsxfun(@times, ppost{i}, 1./sum(ppost{i},2));
    firstlooktranslate{i}{rip} = reshape(ppost_norm{i},r,Xsyze,Ysyze);
    end
        titles = {'allMaps'};%{'inbound';'outbound'};

    figure;
    for i = 1:length(placefieldmaps)
    %subplot(2,1,i)
    imagesc(squeeze(firstlooktranslate{i}{rip}),clims);hold on;title(titles{i});
    ylabel('rip time');xlabel('position');
     plot([25.6 25.6],ylim,'k');hold on;plot([38.4 38.4],ylim,'k');hold on;
     plot([89.6 89.6],ylim,'k');hold on;plot([102.4 102.4],ylim,'k');hold on;
     plot(xlim,[11 11],'--r');
    end
    plot([ size(placefieldmaps{1},2)/4  size(placefieldmaps{1},2)/4],ylim,'r');hold on;
    plot([ size(placefieldmaps{1},2)/2  size(placefieldmaps{1},2)/2],ylim,'r');hold on;
    plot([ size(placefieldmaps{1},2)*(3/4)  size(placefieldmaps{1},2)*(3/4)],ylim,'r');hold on;
     plot([153.6 153.6],ylim,'k');hold on;plot([166.4 166.4],ylim,'k');hold on;
     plot([217.6 217.6],ylim,'k');hold on;plot([230.4 230.4],ylim,'k');hold on;    
    %set(gcf, 'Units', 'Normalized', 'OuterPosition', [0.1, 0.1, 0.7, 0.5]);
    saveas(gcf,['decoding/decodingCL_pul' num2str(rip) '_v3'],'png')%pause;
    close all;
    end
end
clear logpfmaps sumpfmaps firstlook ppost ppost_norm firstlooktranslate

% good = 5
%% Test for spk in trial
for unit = 1:size(ripSpkTask.timesRelat,1)
    for trial = 1:size(trials.int,1)
        edges = trials.int(trial,1):0.2:trials.int(trial,2);
        temp = Restrict(spikes.times{unit},trials.int(trial,:));
        timesBin{unit,trial} = histc(temp,edges);
    end
end
for trial = 1:size(trials.int,1)
        edges = trials.int(trial,1):0.2:trials.int(trial,2);
    for unit = 1:size(ripSpkTask.timesRelat,1)
        if ~isempty(timesBin{unit,trial})
           timesRelatCat{trial}(unit,:) = timesBin{unit,trial};
        else
           timesRelatCat{trial}(unit,:) = zeros(1,length(edges));
        end
    end
end

spikesperbinT = timesRelatCat{6};

% for unit = 1:68
% spikesperbinTConv(unit,:) = conv(spikesperbinT(unit,:),ones(2,1),'same');
% end


tbin= 0.2; Xsyze=1; Ysyze=1280;%320;
[e r] = size(spikesperbinT);
logpfmaps = log(placefieldmaps*tbin);
logpfmaps(find(logpfmaps == -Inf)) = -450;

%logpfmaps(find(logpfmaps == -Inf)) = 0;
sumpfmaps = tbin*sum(placefieldmaps);
firstlook = spikesperbinT' * logpfmaps -repmat(sumpfmaps,r,1);
ppost = exp(firstlook);
ppost_norm = bsxfun(@times, ppost, 1./sum(ppost,2));
firstlooktranslate = reshape(ppost_norm,r,Xsyze,Ysyze);

figure;
imagesc(squeeze(firstlooktranslate));
