

bz_LFPfromDat(pwd);
[ripples] = bz_FindRipples(pwd,0 ,'thresholds',[1 4])



%%
%mouse1_180412, UID = 20 = axo axo
spikes = bz_GetSpikes('getWaveforms',false,'noprompts',true);
basepath = pwd;
basename = bz_BasenameFromBasepath(basepath);





ripFil = [basepath '/' basename '.evt.rip'];
LFPFil = [basepath '/' basename '.lfp'];
rip_evs = LoadEvents(ripFil);

rip_ts = rip_evs.time(cellfun(@any,regexp(rip_evs.description,'peak')));

rip_st = rip_evs.time(cellfun(@any,regexp(rip_evs.description,'start')));

rip_end = rip_evs.time(cellfun(@any,regexp(rip_evs.description,'stop')));




%%
ph = [];
ts = [];

for i = 1:length(rip_ts)
    
    
    st = rip_ts(i);
    pr = round(1250*(st - rip_st(i)));
     po = round(1250*(rip_end(i) - st));
    data = LoadBinary(LFPFil,'frequency',1250,'channels',14,'start',st-1,'duration',2,'nchannels',32);
    
    dataf = BandpassFilter(double(data),1250,[80 220]);
    pht = InstPhase(dataf);
 [wavespec] = bz_WaveSpec(double(data),'samplingRate',1250,'frange',[100 120], 'nfreqs' ,1);
%get phase at max power

InstPhase = atan2(imag(wavespec.data), real(wavespec.data))';
%[a,b] = max(abs(wavespec.data)');

%idx = sub2ind(size(InstPhase),b,1:size(InstPhase,2));
%InstPhase = InstPhase(1,:);
 
ph =  [ph pht((-pr:po)+1251)'];
ts = [ts st+((-pr:po)+1251)/1250];
end

%%
close all
[a,b] = histc(spikes.times{20},ts);
b = b(a(b>0)==1);

bin_spk = histc(ph(b),-pi:pi/10:pi);
bin_ph = histc(ph,-pi:pi/10:pi);
rip_rate = (bin_spk./bin_ph)*1250;

bar(-pi:pi/10:pi,rip_rate)

%%

loner = rip_st(diff([0;rip_st])>.5);
rip_peth = [];
for i = 1:length(spikes.times)
rip_peth(i,:) = CrossCorr(loner,spikes.times{i},.01,100);
end

%%
bar(-500:10:500,ok)


%%
preax =CrossCorr(spikes.times{35},spikes.times{20},.001,100);
prenonax =CrossCorr(spikes.times{35},spikes.times{6},.001,100);
%%

stims = LoadEvents([basename '.evt.ait']);
on = stims.time(cellfun(@any,regexp(stims.description,'pulse_on')));
off = stims.time(cellfun(@any,regexp(stims.description,'pulse_off')));
%%
figure
res = .02;
dur = .5;
ok = CrossCorr(on,rip_st,res,2*(dur/res));
ok = ok/length(ok) / dur;
bar(-dur:res:dur,ok)
hold on
plot([0 0;0 0],[0 0;2 2])

 
 %%
 
 %16, 20
stim_peth = [];
for i = 1:length(spikes.times)
stim_peth(i,:) = CrossCorr(on,spikes.times{i},.01,100);
end


 stim_peth = (stim_peth/length(on)*100);
%%
 
 [status, interval] = InIntervals(spikes.times{20},[on off]);

 spks = spikes.times{20}(~status);
 
 ok = CrossCorr(spks,spks,.001,200);
 ok1 = ok([1:100 102:end]);
 
 %%
  [status, interval] = InIntervals(spikes.times{16},[on off]);

 spks = spikes.times{16}(~status);
 
 ok = CrossCorr(spks,spks,.001,200);
 ok = ok([1:100 102:end]);
 