% Jon Newman

%% file id
fids = {
    'J:/newmonster/2012_07_02_Xeko-1-minute-clamps/2012_07_02_Xeko-20709_one-min-clamp_07';
    'J:/newmonster/2012_10_11-Herb-I5387-1-minute-clamps/2012_10_11_Herb-I5387_one-min-clamp_02';
    'J:/newmonster/2012_10_22_Herb-I5386-1-minute-clamps/2012_10_22_Herb-I5386_one-min-clamp_10';
    'Z:/jon/2013_02_26_Yam-22565_1-minute-clamp/2013-02-26_yam_22565_one-minute-clamp_06';
    'Z:/jon/2012_11_21_Mandelbrot-21332-1-minute-clamps/2012_11_21_Mandelbrot_21332_one-min-clamp_12';
    };



sym1 = {'o','o','o','o','o','o'};


RT = nan(length(fids),13); % rise-time
TAU = nan(length(fids),13); % time-constant
ST = nan(length(fids),13); % settling time
JIT = nan(length(fids),13); % settling time
% ctl_2 = nan(length(fids),11);

dfr_r = 0.1;
dfr_s = 0.25;
sw = 2.5; %sec
st_time = 1; % second
loop_dt = 0.01;

for f = 1:length(fids)
    
    fid = fids{f};
    load(fid);
    
    % Determine rise-time measures
    for i = 1:length(avg_perf.steps)
        
        targ = avg_perf.steps(i);
        
        b = ts.t{i} > 0;
        t = ts.t{i}(b);
        loop_dt = mean(diff(t));
        
        tru_fr = ts.fr{i}(b);
        sm_fr = smooth(ts.fr{i}(b),round(sw/loop_dt),'lowess');
        
        if size(sm_fr) ~= size(t)
            sm_fr = sm_fr';
        end
        
        
        rt_i = find(abs(sm_fr - targ) < dfr_r & t > 0.1,1,'first');
        [~, mt_i] = max(ts.fr{i});
        
        if ~isempty(rt_i)
            RT(f,i) = loop_dt*rt_i;
        end
        
        ind_in_bound = find(abs(sm_fr - targ) < dfr_s); % & t > 2);
        ind_out_bound = find(abs(sm_fr - targ) > dfr_s);
        
        if  isempty(ind_in_bound)
            continue;
        elseif size(ind_in_bound,2) > 1
            ind_in_bound = ind_in_bound';
            ind_out_bound = ind_out_bound';
        end
        
        ind_enter_bound = [ind_in_bound(1);  ind_in_bound(find(diff(ind_in_bound) > 1)+1)];
        ind_exit_bound = [ind_out_bound(1); ind_out_bound(find(diff(ind_out_bound) > 1)+1)];
        
        
        if ind_exit_bound(1) < ind_enter_bound(1)
            ind_exit_bound(1) = [];
        end
        if length(ind_exit_bound) ~= length(ind_enter_bound)
            ind_exit_bound = [ind_exit_bound; length(sm_fr)];
        end
        
        t_in_bound = ind_exit_bound-ind_enter_bound;
        ind_settled = ind_enter_bound(find(t_in_bound > st_time/loop_dt,1,'first'));
        
        if ~isempty(ind_settled)
            ST(f,i) = loop_dt*ind_settled;
            JIT(f,i) = avg_perf.sfr(i);
            
            
            
        end
        
        close all
        figure
        hold on
        plot(tru_fr,'y')
        plot(sm_fr,'k')
        plot(ind_settled,sm_fr(ind_settled),'r*')
        
        
    end
    
end

JIT(JIT > 1) = nan; % remove the case where the there was a control failure

%%

close all

col = [ [0 100 255];
    [255 205 41]]/255;

col2 = [
    [220 000 000];
    [255 102 000];
    [255 205 041];
    [000 170 000];
    [000 170 212];
    [000 100 255];
    [085 000 212];
    [102 102 102]
    ]/255;

figure('units','centimeter','position',[40 5 15 20]);

ax1 = axes('units','centimeter','position',[5 10 3 3]);
hold on;
set(gca,'TickDir','out','FontName','Arial','FontSize',10,'TickLength',[0.025 0.025],'Xlim',[0 10],'Ylim',[0 20]);
ylabel('Settling time (sec)')
xlabel('Target rate (Hz/unit)')


for i = 1:size(ST,1)
    x = 0:10;
    plot(x,ST(i,1:11),'o','MarkerSize',4,'Color',col2(i,:),'MarkerSize',3,'MarkerFaceColor',[1 1 1]);
    ind_out = find(ST(i,1:11) > 20);
    if ~isempty(ind_out)
       plot(x(ind_out),19*ones(size(ind_out)), 'r^')
        
    end
    
end

plot(0:12,nanmean(ST),'k-','Linewidth',2)


ax2 = axes('units','centimeter','position',[10 10 3 3]);
hold on;
set(gca,'TickDir','out','FontName','Arial','FontSize',10,'TickLength',[0.025 0.025],'Xlim',[0 10],'Ylim',[0 0.2]);
ylabel('SD FR (Hz/unit)')
xlabel('Target rate (Hz/unit)')

for i = 1:size(ST,1)
    
    plot(0:10,JIT(i,1:11),'o','MarkerSize',4,'Color',col2(i,:),'MarkerSize',3,'MarkerFaceColor',[1 1 1])
end

plot(0:12,nanmean(JIT),'k-','Linewidth',2)

%%
saveas(gca,'1-minute-clamp_settling-time_and_jitter_01.fig')
export_fig -pdf -transparent 1-minute-clamp_settling-time_and_jitter_01.pdf

