eyeHeadVelocity.m

filename = 'ExampleData.csv';
data = csvread(filename,1);
t = data(:,1);
velocityThresh = 10; % eye-balled this
minpeakdist = 100; % in samples
extractRange = -20:100; % padding not handled
t_compiled = t(1:numel(extractRange));

eyeHorizontal = rad2deg(data(:,2));
eyeVelocity = filterVelocity(diff(eyeHorizontal));
[eyeLocs_pos,eyePks_pos] = peakseek(eyeVelocity,minpeakdist,velocityThresh);
eyeLocs_pos = eyeLocs_pos + 1; % fix diff
eyeCompiled_pos = extractSignal(eyeHorizontal,eyeLocs_pos,extractRange);

[eyeLocs_neg,eyePks_neg] = peakseek(-eyeVelocity,minpeakdist,velocityThresh);
eyeLocs_neg = eyeLocs_neg + 1; % fix diff
eyePks_neg = -eyePks_neg; % convert sign
eyeCompiled_neg = extractSignal(eyeHorizontal,eyeLocs_neg,extractRange);

headHorizontal = rad2deg(data(:,3));
headVelocity = filterVelocity(diff(headHorizontal));
[headLocs_pos,headPks_pos] = peakseek(headVelocity,minpeakdist,velocityThresh);
headLocs_pos = headLocs_pos + 1; % fix diff
headCompiled_pos = extractSignal(headHorizontal,headLocs_pos,extractRange);

[headLocs_neg,headPks_neg] = peakseek(-headVelocity,minpeakdist,velocityThresh);
headLocs_neg = headLocs_neg + 1; % fix diff
headPks_neg = -headPks_neg; % convert sign
headCompiled_neg = extractSignal(headHorizontal,headLocs_neg,extractRange);

close all;
move_ylim = [-30 30];
vel_ylim = [-50 50];
h = figure('position',[0 0 1400,800]);

subplot(221);
yyaxis left;
plot(t,eyeHorizontal);
ylabel('movement');
ylim(move_ylim);
yticks(sort([0,ylim]));
yyaxis right;
plot(t(2:end),eyeVelocity);
hold on;
plot(t,repmat(velocityThresh,size(t)),'k:');
plot(t,-repmat(velocityThresh,size(t)),'k:');
text(0.5,velocityThresh,'thresh','verticalAlignment','bottom');
plot(t(eyeLocs_pos),eyePks_pos,'k*');
plot(t(eyeLocs_neg),eyePks_neg,'r*');
xlim([min(t) max(t)]);
ylim(vel_ylim);
ylabel('velocity');
yticks(sort([0,ylim]));
xlabel('time (s)');
title('Eye');

subplot(223);
yyaxis left;
plot(t,headHorizontal);
ylabel('movement');
ylim(move_ylim);
yticks(sort([0,ylim]));
yyaxis right;
plot(t(2:end),headVelocity);
hold on;
plot(t,repmat(velocityThresh,size(t)),'k:');
plot(t,-repmat(velocityThresh,size(t)),'k:');
text(0.5,velocityThresh,'thresh','verticalAlignment','bottom');
plot(t(headLocs_pos),headPks_pos,'k*');
plot(t(headLocs_neg),headPks_neg,'r*');
xlim([min(t) max(t)]);
ylabel('velocity');
ylim(vel_ylim);
yticks(sort([0,ylim]));
xlabel('time (s)');
title('Head');

subplot(222);
plot(t_compiled,eyeCompiled_pos,'k');
hold on;
plot(t_compiled,eyeCompiled_neg,'r');
xlim([min(t_compiled) max(t_compiled)]);
ylabel('movement');
ylim(move_ylim);
yticks(sort([0,ylim]));
xlabel('time (s)');
title('Eye Compiled');

subplot(224);
plot(t_compiled,headCompiled_pos,'k');
hold on;
plot(t_compiled,headCompiled_neg,'r');
xlim([min(t_compiled) max(t_compiled)]);
ylabel('movement');
ylim(move_ylim);
yticks(sort([0,ylim]));
xlabel('time (s)');
title('Head Compiled');

set(gcf,'color','w');

saveas(h,[mfilename,'.fig']);
saveas(h,[mfilename,'.png']);

function compiled = extractSignal(data,dataLocs,dataRange)
    % if dataRange is out of array bounds this will throw an error
    compiled = zeros(numel(dataLocs),numel(dataRange));
    for iLoc = 1:numel(dataLocs)
        compiled(iLoc,:) = data(dataLocs(iLoc) + dataRange - 1);
    end
end

function data = filterVelocity(data)
    % just cleans the data a bit, makes it 'peakier'
    data = smooth(sign(data) .* data.^2);
end