%% ACh OXT simultaneous recording: behavior tranist

dataname = ['F:\oxytocin\ACh_OXT004\ACh_OXT004_007_230306_152436', %% corrected
            'F:\oxytocin\ACh_OXT004\ACh_OXT004_011_230325_130235',   %% corrected
            'F:\oxytocin\ACh_OXT007\ACh_OXT007_015_230507_172920',   %% corrected
            'F:\oxytocin\ACh_OXT004\ACh_OXT004_019_230507_140912']   %% corrected no REM
dF_NRtranAll = {}; dF_NRtranAll.ACh = [];dF_NRtranAll.Oxt = [];
dF_RNtranAll = {}; dF_RNtranAll.ACh = [];dF_RNtranAll.Oxt = [];
dF_RWtranAll = {}; dF_RWtranAll.ACh = [];dF_RWtranAll.Oxt = [];
dF_NWtranAll = {}; dF_NWtranAll.ACh = [];dF_NWtranAll.Oxt = [];
for ii = 1: size(dataname,1)
    cd(dataname(ii,:))
    basepath = pwd;
    basename = bz_BasenameFromBasepath(basepath);
    load([basename,'.SleepState.states.mat']);
    
    info = read_Intan_RHD2000_file(basepath);
    digitalin = read_Intan_digitalin_file('digitalin.dat',info.num_board_dig_in_channels)';
    lfp_fbtime = find(digitalin(1,:)>(max(digitalin)+min(digitalin))/2,1);
    delete digitalin
    % sitmulation and signal
    if exist([basename,'_FluoData.mat'])==0
        FluoData = Get_OxtAch_Fluos();
        figure
        plot(FluoData.time,FluoData.Ach_dFF); hold on;
        plot(FluoData.time,FluoData.Oxt_dFF);xlim([0 inf]);
    else load([basename,'_FluoData.mat'])
    end
    %% find NREM-REM
    sleep_ints = bz_INTtoIDX(SleepState.ints,statenames ={'WAKE','','NREM','','REM'});
    idx = find(ischange(sleep_ints.states)==1) ;
    aa = diff(idx);idx_ff =[];
    for k = 1:length(aa)-1
        if aa(k)>20 & aa(k+1)>20
            ff = idx(k+1);
            idx_ff = [idx_ff;ff]
        end
    end
    idx_ff = idx_ff(idx_ff < length(dFoverF)/100 + lfp_fbtime(1)/20000 & idx_ff > lfp_fbtime(1)/20000);

    %% seperate into NREM-REM, NREM-WAKE, REM-WAKE
    NR_tim = []; NW_tim = []; RW_tim = [];
    for j = 1:length(idx_ff)
        if sleep_ints.states(idx_ff(j)-1) == 3  & sleep_ints.states(idx_ff(j)) == 5
            NR_tim = [NR_tim;idx_ff(j)];
        elseif sleep_ints.states(idx_ff(j)-1) == 3  & sleep_ints.states(idx_ff(j)) == 1
            NW_tim = [NW_tim;idx_ff(j)];
        else sleep_ints.states(idx_ff(j)-1) == 5  & sleep_ints.states(idx_ff(j)) == 1
            RW_tim = [RW_tim;idx_ff(j)];
        end
    end
    %% get NREM-REM DATA

    NR_tim_dff = NR_tim -lfp_fbtime(1)/20000;
    dF_NR_transit = [];
    for jj = 1:length(NR_tim_dff)
        idx = round((NR_tim_dff(jj) + [-20 20])*100);
        dF_NR_transit(:,jj) = dFoverF(idx(1):idx(2));
    end
    dF_NR_transit = dF_NR_transit';
    dF_NR_tranAll = [dF_NR_tranAll;dF_NR_transit];
    
   %% get NREM-WAKE
    NW_tim_dff= NW_tim -lfp_fbtime(1)/20000;
    dF_NW_transit = [];
    for jj = 1:length(NW_tim_dff)
        idx = round((NW_tim_dff(jj) + [-20 20])*100);
        dF_NW_transit(:,jj) = dFoverF(idx(1):idx(2));
    end
    dF_NW_transit = dF_NW_transit';
    dF_NW_tranAll = [dF_NW_tranAll;dF_NW_transit];
  %% get REM-WAKE
    RW_tim_dff= RW_tim -lfp_fbtime(1)/20000;
    dF_RW_transit = [];
    for jj = 1:length(RW_tim_dff)
        idx = round((RW_tim_dff(jj) + [-20 20])*100);
        dF_RW_transit(:,jj) = dFoverF(idx(1):idx(2));
    end
    dF_RW_transit = dF_RW_transit';
    dF_RW_tranAll = [dF_RW_tranAll;dF_RW_transit];


end
%% figure NREM-REM
figure;
outsr = 100;
imagesc([1:size(dF_NRtranAll.ACh,2)./outsr]'-20,[1:size(dF_NRtranAll.ACh,1)]',...
    zscore(dF_NRtranAll.ACh,[],2));
zlim = [-8 3];
caxis(zlim)
colormap(NegativeEnhancingColormap(256, zlim, ...
        [0 0 1], [1 0 0],[1 1 1],1.25));
ylabel('Cell Number');colorbar
xlabel('Times from transit (s)'); ylabel('REM rumber #'); title('NREM-REM-ACh')

figure;
outsr = 100;
imagesc([1:size(dF_NRtranAll.Oxt,2)./outsr]'-20,[1:size(dF_NRtranAll.Oxt,1)]',...
    zscore(dF_NRtranAll.Oxt,[],2));
zlim = [-8 3];
caxis(zlim)
colormap(NegativeEnhancingColormap(256, zlim, ...
        [0 0 1], [1 0 0],[1 1 1],1.25));
ylabel('Cell Number');colorbar
xlabel('Times from transit (s)'); ylabel('REM rumber #'); title('NREM-REM-OXT')

%% figure REM-wake
figure;
outsr = 100;
imagesc([1:size(dF_RWtranAll.ACh,2)./outsr]'-20,[1:size(dF_RWtranAll.ACh,1)]',...
    zscore(dF_RWtranAll.ACh,[],2));
zlim = [-8 3];
caxis(zlim)
colormap(NegativeEnhancingColormap(256, zlim, ...
        [0 0 1], [1 0 0],[1 1 1],1.25));
ylabel('Cell Number');colorbar
xlabel('Times from transit (s)'); ylabel('REM rumber #');title('REM-WAKE-ACh')

figure;
outsr = 100;
imagesc([1:size(dF_RWtranAll.Oxt,2)./outsr]'-20,[1:size(dF_RWtranAll.Oxt,1)]',...
    zscore(dF_RWtranAll.Oxt,[],2));
zlim = [-8 3];
caxis(zlim)
colormap(NegativeEnhancingColormap(256, zlim, ...
        [0 0 1], [1 0 0],[1 1 1],1.25));
ylabel('Cell Number');colorbar
xlabel('Times from transit (s)'); ylabel('REM rumber #');title('REM-WAKE-Oxt')

%% NREM-wake
figure;
outsr = 100;
imagesc([1:size(dF_NWtranAll.ACh,2)./outsr]'-50,[1:size(dF_NWtranAll.ACh,1)]',...
    zscore(dF_NWtranAll.ACh,[],2));
zlim = [-8 3];
caxis(zlim)
colormap(NegativeEnhancingColormap(256, zlim, ...
        [0 0 1], [1 0 0],[1 1 1],1.25));
ylabel('Cell Number');colorbar
xlabel('Times from transit (s)'); ylabel('REM rumber #');title('NREM-WAKE')

figure;
outsr = 100;
imagesc([1:size(dF_NWtranAll.Oxt,2)./outsr]'-50,[1:size(dF_NWtranAll.Oxt,1)]',...
    zscore(dF_NWtranAll.Oxt,[],2));
zlim = [-8 3];
caxis(zlim)
colormap(NegativeEnhancingColormap(256, zlim, ...
        [0 0 1], [1 0 0],[1 1 1],1.25));
ylabel('Cell Number');colorbar
xlabel('Times from transit (s)'); ylabel('REM rumber #');title('NREM-WAKE')

save('G:\My Drive\Paper_ACh_OXT\Fig2\Data\Figure2B\ACh_OXT_behaviorTrans.mat','dF_NRtranAll','dF_RWtranAll','dF_NWtranAll')
