
clear all;

datafile = {'F:\oxytocin\OXYT010\OXYT010_004_220126_112539',
            'F:\oxytocin\OXYT012\OXYT012_001_220323_130702',
            'F:\oxytocin\OXYT011\OXYT011_005_220224_160927',
            'F:\oxytocin\OXYT013\OXYT013_001_220411_204913',
            'F:\oxytocin\OXYT014\OXYT014_001_220429_155806'}


for ii = 1:length(datafile)
    cd(datafile{ii,:})
    
    basepath = pwd;
    basename = bz_BasenameFromBasepath(basepath);
    load([basename,'.SleepScoreLFP.LFP.mat']);
    load([basename,'.SleepState.states.mat'])

    % get alighment
    info = read_Intan_RHD2000_file(basepath);
    analogin = read_Intan_analogin_file('analogin.dat',info.num_board_adc_channels);
    lfp_fbtime = find(analogin(1,:)>0.6,1);lfp_fbtime(1,2) = find(analogin(1,:)>0.6,1,'last');
    delete analogin
    % get name 
    if exist('dFoverF.mat')==0
    
        [filepath,name,ext] = fileparts(basename);
        [k y] = strtok(name,'_');
        name = [k,y(1:4)];
        % get syn
        load([name,ext],[name,'_AChSti']);
        ACh_sti = eval([name,'_AChSti']);
        fbtime = ACh_sti.times(find(ACh_sti.values(:,1)>0.5,1));fbtime(1,2) = ACh_sti.times(find(ACh_sti.values(:,1)>0.5,1,'last'));
        fbIDX = round(fbtime*20000);
        % signal
        load([name,ext],[name,'_AChSig1']);
        ACh_signal = eval([name,'_AChSig1']);
        ACh_signal.data = downsample(ACh_signal.values((fbIDX(1,1)):(fbIDX(1,2))),200);
        ACh_signal.time = downsample(ACh_signal.times((fbIDX(1,1)):(fbIDX(1,2))),200)- fbtime(1,1);
        
        p = [];
        p = polyfit(ACh_signal.time,ACh_signal.data,5);
        Ybaseline = polyval(p,ACh_signal.time);
        
        figure
        plot(ACh_signal.time, ACh_signal.data); hold on;
        plot(ACh_signal.time,Ybaseline);
        dFoverF = smoothdata((ACh_signal.data - Ybaseline)./Ybaseline,'movmean',10);
        hold on;
        plot(ACh_signal.time, dFoverF);xlim([0 inf])
    else load('dFoverF.mat'); figure;plot(dFoverF);
    end
    %get LFPs
    thCh = SleepScoreLFP.THchanID;
    THLFPs = bz_GetLFP(thCh,'basepath',basepath,'noPrompts',true);


    datpath = 'G:\My Drive\Paper_ACh_OXT\Fig1\Data\Figure1F\OXT_theta'
    save(fullfile(datpath,[basename,'_OXT_theta.mat']),'THLFPs','lfp_fbtime','dFoverF','SleepState');
end
