%%%%%%%%%%%%%%
% Function for selecting a specific interval of the CARISMA B PSD in
% the B time series and splitting it into sections for analysis later
% with FFT.
%%%%%%%%%%%%%%

function [sampwin set1 set2 subset1 subset2 Ntitlestr] = samplewindow(time, B_output, spikeloc, smoothing, smoothspan, logging, units, npt, plottoggle)

global fs overlap

global fntsize res xsize ysize lmarg bmarg widthv heightv xstr titlestr fstr pstr colpl1 colpl2

% Define axes
if strcmp(plottoggle,'on')
ax1 = subplot(4,1,1);
ax2 = subplot(4,2,3);
ax3 = subplot(4,2,4);
ax4 = subplot(4,1,3);
ax5 = subplot(4,2,7);
ax6 = subplot(4,2,8);
end

% Define spike location
N = spikeloc;

sampstart = (N-1)*(npt-overlap); % start of vector overlaps with previous window
sampend = sampstart + npt - 1; % end of vector is npt points later

sampwintime = time(sampstart:sampend); % extract interval vector for time axis

sampwindata = B_output(sampstart:sampend); % extract interval vector for data axis

sampwin = {sampwintime sampwindata}; % create cell with datetime time vector and data vector

formatout = 'HH:MM:SS.FFF'; % format for addition of time interval to title string

Ntitlestr = strcat(titlestr, {' '}, datestr(sampwintime(1),formatout), '-', ...
    datestr(sampwintime(end),formatout)); % title string that specifies the time interval

set1 = {sampwintime(1:length(sampwintime)/2) ...
    sampwindata(1:length(sampwintime)/2)}; % create cell of first set
set2 = {sampwintime(length(sampwintime)/2+1:end) ...
    sampwindata(length(sampwintime)/2+1:end)}; % create cell of second set

% Calculate FFT of sets 1 & 2
Lset1 = length(set1{1}); % length of subset
set1{3} = fs*(0:(Lset1/2))/Lset1; % frequency vector
sampfftset1 = fft(set1{2},length(set1{2})); % fft of time series y values
P2sampset1 = abs(sampfftset1/Lset1); % two sided spectrum of the fft
P1sampset1 = P2sampset1(1:Lset1/2+1); % single sided spectrum based on P2samp1
P1sampset1(2:end-1) = 2*P1sampset1(2:end-1); % single sided spectrum based on P2samp1
if strcmp(logging,'log') % toggle log of output
    set1{4} = log(P1sampset1); % take log of single sided spectrum and save in the cell
        pstr = 'log|P1(f)|';
else
    set1{4} = P1sampset1; % save in the cell
        pstr = '|P1(f)|';
end
set1{5} = smooth(set1{4},smoothspan,'moving'); % smoothing fft of time series y values with a moving average
Lset2 = length(set2{1}); % length of subset
set2{3} = fs*(0:(Lset2/2))/Lset2; % frequency vector
sampfftset2 = fft(set2{2},length(set2{2})); % fft of time series y values
P2sampset2 = abs(sampfftset2/Lset2); % two sided spectrum of the fft
P1sampset2 = P2sampset2(1:Lset2/2+1); % single sided spectrum based on P2samp1
P1sampset2(2:end-1) = 2*P1sampset2(2:end-1); % single sided spectrum based on P2samp1
if strcmp(logging,'log') % toggle log of output
    set2{4} = log(P1sampset2); % take log of single sided spectrum and save in the cell
else
    set2{4} = P1sampset2; % save in the cell
end
set2{5} = smooth(set2{4},smoothspan,'moving'); % smoothing fft of time series y values with a moving average
if strcmp(plottoggle,'on')
% Plot sets 1 & 2 time series
axes(ax1);
plot(set1{1},set1{2},'b',set2{1},set2{2},'r');
title(ax1,Ntitlestr, 'fontsize', fntsize, 'fontname', 'helvetica');
axis tight;
xlim(ax1,[sampwintime(1) sampwintime(end)]);
yL1 = get(gca,'YLim');
line([set1{1}(end) set1{1}(end)],yL1,'Color','k','LineStyle','--');
xlabel(xstr, 'fontsize', fntsize, 'fontname', 'helvetica'); % label axes and title
ylabel(strcat('B (',units,'/s)'));
legend(ax1,'Set 1','Set 2');

% Plot sets 1 & 2 FFT's
axes(ax2);
if strcmp(smoothing,'smooth')
    plot(set1{3},set1{5},'b');
    title('Smoothed single sided amplitude spectrum of B, set 1', 'fontsize', fntsize, 'fontname', 'helvetica');
else
    plot(set1{3},set1{4},'b');
    title('Single sided amplitude spectrum of B, set 1', 'fontsize', fntsize, 'fontname', 'helvetica');
end
axis tight;
xlabel(fstr, 'fontsize', fntsize, 'fontname', 'helvetica');
ylabel(pstr, 'fontsize', fntsize, 'fontname', 'helvetica');

axes(ax3);
if strcmp(smoothing,'smooth')
    plot(set2{3},set2{5},'r');
    title('Smoothed single sided amplitude spectrum of B, set 2', 'fontsize', fntsize, 'fontname', 'helvetica');
else
    plot(set2{3},set2{4},'r');
    title('Single sided amplitude spectrum of B, set 2', 'fontsize', fntsize, 'fontname', 'helvetica');
end
axis tight
xlabel(fstr, 'fontsize', fntsize, 'fontname', 'helvetica');
ylabel(pstr, 'fontsize', fntsize, 'fontname', 'helvetica');
end
% Prep plot for subsets
colpl1 = {[1 0 1] [0.5 0 0] [1 0 0] [1 0.9 0] [0.4 0.8 0] [1 0.5 0] [0 0 1] [0.5 0.5 0.5]}; % Predefined colour pattern for set 1
colpl2 = {[1 0 1] [0.5 0 0] [1 0 0] [1 0.9 0] [0.4 0.8 0] [1 0.5 0] [0 0 1] [0.5 0.5 0.5]}; % Predefined colour pattern for set 2
for i = 1:8
    subset1{i,1} = set1{1}((i-1)*length(set1{1})/8+1:i*length(set1{1})/8); % time vector
    subset1{i,2} = set1{2}((i-1)*length(set1{2})/8+1:i*length(set1{2})/8); % time series y vector
    L1 = length(subset1{i,1}); % length of subset
    subset1{i,3} = fs*(0:(L1/2))/L1; % frequency vector
    sampfft1 = fft(subset1{i,2},length(subset1{i,2})); % fft of time series y values
    P2samp1 = abs(sampfft1/L1); % two sided spectrum of the fft
    P1samp1 = P2samp1(1:L1/2+1); % single sided spectrum based on P2samp1
    P1samp1(2:end-1) = 2*P1samp1(2:end-1); % single sided spectrum based on P2samp1
    if strcmp(logging,'log')
        subset1{i,4} = log(P1samp1); % take log of single sided spectrum and save in the cell
    else
        subset1{i,4} = P1samp1; % save in the cell
    end
    subset1{i,5} = smooth(subset1{i,4},smoothspan,'moving'); % smoothing fft of time series y values with a moving average
    subset2{i,1} = set2{1}((i-1)*length(set2{1})/8+1:i*length(set2{1})/8); % time vector
    subset2{i,2} = set2{2}((i-1)*length(set2{2})/8+1:i*length(set2{2})/8); % time series y vector
    L2 = length(subset2{i,1}); % length of subset
    subset2{i,3} = fs*(0:(L2/2))/L2; % frequency vector
    sampfft2 = fft(subset2{i,2},length(subset2{i,2})); % fft of time series y values
    P2samp2 = abs(sampfft2/L2); % two sided spectrum of the fft
    P1samp2 = P2samp2(1:L2/2+1); % single sided spectrum based on P2samp1
    P1samp2(2:end-1) = 2*P1samp2(2:end-1); % single sided spectrum based on P2samp1
    if strcmp(logging,'log')
        subset2{i,4} = log(P1samp2); % take log of single sided spectrum and save in the cell
    else
        subset2{i,4} = P1samp2; % save in the cell
    end
    subset2{i,5} = smooth(subset2{i,4},smoothspan,'moving'); % smoothing fft of time series y values with a moving average
    if strcmp(plottoggle,'on')
    axes(ax4);
    hold on;
    p1(i) = plot(subset1{i,1},subset1{i,2},'color',colpl1{i});
    p2(i) = plot(subset2{i,1},subset2{i,2},'color',colpl2{i});
    axes(ax5);
    hold on;
    if strcmp(smoothing,'smooth')
        plot(subset1{i,3},subset1{i,5},'color',colpl1{i});
        title('Smoothed single sided amplitude spectrum of B, subset 1-8', 'fontsize', fntsize, 'fontname', 'helvetica');
    else
        plot(subset1{i,3},subset1{i,4},'color',colpl1{i});
        title('Single sided amplitude spectrum of B, subset 1-8', 'fontsize', fntsize, 'fontname', 'helvetica');
    end
    axes(ax6);
    hold on;
    if strcmp(smoothing,'smooth')
        plot(subset2{i,3},subset2{i,5},'color',colpl2{i});
        title('Smoothed single sided amplitude spectrum of B, subset 9-16', 'fontsize', fntsize, 'fontname', 'helvetica');
    else
        plot(subset2{i,3},subset2{i,4},'color',colpl2{i});
        title('Single sided amplitude spectrum of B, subset 9-16', 'fontsize', fntsize, 'fontname', 'helvetica');
    end
    end
end

if strcmp(plottoggle,'on')
axes(ax4);
title(Ntitlestr, 'fontsize', fntsize, 'fontname', 'helvetica');
axis tight;
xlim([sampwintime(1) sampwintime(end)]);
yL2 = get(gca,'YLim');
line([set1{1}(end) set1{1}(end)],yL2,'Color','k','LineStyle','--');
xlabel(xstr, 'fontsize', fntsize, 'fontname', 'helvetica'); % label axes and title
ylabel(strcat('B (',units,'/s)'), 'fontsize', fntsize, 'fontname', 'helvetica');
legend(ax4,p1,'Subset 1/9','Subset 2/10','Subset 3/11','Subset 4/12','Subset 5/13','Subset 6/14','Subset 7/15','Subset 8/16','Orientation','horizontal');

axes(ax5);
axis tight;
xlabel(fstr, 'fontsize', fntsize, 'fontname', 'helvetica');
ylabel(pstr, 'fontsize', fntsize, 'fontname', 'helvetica');

axes(ax6);
axis tight;
xlabel(fstr, 'fontsize', fntsize, 'fontname', 'helvetica');
ylabel(pstr, 'fontsize', fntsize, 'fontname', 'helvetica');
end
end