% Jon Newman
% Ga Tech

% 17 May 2010

% Find bursts in extracellular data, creates a plot showing location and
% amplitude on ASDR plot and saves a new data structure with this info.

tankID = 'C:\TDTDATA\20100420\20100420-25';
[directory file] = queryfid(tankID);
fid = strcat([directory '\' file]);
load(fid);

spktime = data.extrarec.spiketime;

%% Find bursts
asdrbin = 0.1; % seconds
[ASDR bins] = asdr(spktime,asdrbin);
ASDR = smooth(ASDR,4);
FRthresh = 0; hardthresh = 15; windsize = 30;
bursts = findbursts(ASDR,bins,windsize,FRthresh,hardthresh); % find bursts compared to 100 samples of background

%global firing rate
frsmooth = 1000;
firerate = smooth(ASDR,frsmooth); % smoothed firing rate
firerate = firerate(200:end-200);
fireratetime = bins(200:end-200)/3600;

if ~isempty(bursts)
    burstrate = [1./smooth(diff(bursts(:,1)),10);nan]; % smoothed burst rate
    burstrate = burstrate(5:end-5);
    burstratetime = bursts(5:end-5,1)/3600;
else
    burstrate = [];
end


h = figure();
subplot(211)
hold on
plot(bins/3600,ASDR,'k');
if ~isempty(bursts)
    plot(bursts(:,1)/3600,bursts(:,2),'*','color',[0 0.5 0.1]);
end
xlabel('Time (hrs)')
ylabel('ASDR (s^-1)')
title('FR and Burst Analysis')

if isempty(burstrate)
    subplot(212)
    plot(fireratetime,firerate,'b','LineWidth',2);
    ylabel('Smoothed AWFR (Hz)')
    xlabel('Time (hrs)')
else
    subplot(212)
    [AX,H1,H2] = plotyy(fireratetime,firerate,burstratetime,burstrate);
    
    set(get(AX(1),'Ylabel'),'String','Smoothed AWFR (Hz)')
    set(get(AX(2),'Ylabel'),'String','Smoothed AWBR (Hz)')
    set(get(AX(1),'Xlabel'),'String','Time (hrs)')
    
    set(H1,'LineWidth',2)
    set(H2,'LineWidth',2)
end

%% Modify data file

data.extrarec.bursts.bursts = bursts;
data.extrarec.bursts.FRthresh = FRthresh;
data.extrarec.bursts.hardthresh = hardthresh;
data.extrarec.bursts.windsize = windsize;
data.extrarec.burstrate.AWBR = burstrate;
data.extrarec.burstrate.smoothwind = 10;
data.extrarec.frrate.AWFR = firerate;
data.extrarec.frrate.smoothwind= frsmooth;
data.extrarec.ASDR.binsize = asdrbin;
data.extrarec.ASDR.ASDR = [bins' ASDR];

save(strcat([directory '\' file]),'data')

figdir = strcat([fid(1:end-4) '_FIG']);
figfid = strcat([file(1:end-4) '_BurstAnalysis']);
if ~exist(figdir,'dir')
    mkdir(figdir);
end
saveas(h, strcat([figdir '\' figfid '.fig']))
saveas(h, strcat([figdir '\' figfid '.jpg']))
