clear; clc;

cedpath = 'C:\CEDMATLAB\CEDS64ML';
addpath(cedpath);
CEDS64LoadLib(cedpath);

%%calculate all
mvlist = 1;
fps = 100; start = 100; stop = 3000;
optictrack = csvread('ham37_03_003.csv',7,1);
idx = find(optictrack(:,2)==0);
optictrack(idx,2:8) = nan(size(idx,1),7);
optictrack = fillmissing(optictrack,'linear');
t = optictrack(:,1);
x = optictrack(:,6);
y = optictrack(:,8);
vel = sqrt((x(2:end) - x(1:end-1)).^2 + (y(2:end) - y(1:end-1)).^2)*fps;
vel = smoothdata(vel,'movmean',fps*mvlist);

 % Open a file
sr1 = 20000; sr2 = 80000;
fhand1 = CEDS64Open('ham37_03-003_2019-03-12.smrx');
maxTimeTicks = CEDS64MaxTime(fhand1) + 1;

ic = 2;
[fRead, fVals, fTime] = CEDS64ReadWaveF(fhand1, ic, maxTimeTicks, 0, maxTimeTicks);
driver = fVals;

ic = 3;
[fRead, fVals, fTime] = CEDS64ReadWaveF(fhand1, ic, maxTimeTicks, 0, maxTimeTicks);
signal = fVals;

ic = 4;
[iRead, vi64T] = CEDS64ReadEvents(fhand1, ic, maxTimeTicks, 0, maxTimeTicks);
dSecs = CEDS64TicksToSecs(fhand1, vi64T );



tim = CEDS64TicksToSecs(fhand1, [1:maxTimeTicks]');
tim = downsample(tim(start*sr2+1:stop*sr2),800);
tim_cc = interp1(1:size(tim,1),tim,1:size(vel,1),'linear')';
signal = downsample(signal(start*sr1+1:stop*sr1),200); sr = 100;
signal_d = smoothdata(detrend(signal),'movmean',sr*mvlist);
signal_cc = interp1(1:size(signal_d,1),signal_d,1:size(vel,1),'linear')';
signal_cc = fillmissing(signal_cc,'linear');



figure;
subplot(4,2,[1,2,3,4]);
plot(t(2:end)+dSecs(1),vel); hold on;
plot(tim_cc,signal_cc); xlim([start,stop]);
cor = corrcoef(signal_cc,vel);
title(sprintf('Average windows: %d sec.; Corrcoef: %.3f', mvlist,cor(1,2)));

subplot(4,2,5);
plot(t(2:end)+dSecs(1),vel); hold on;
plot(tim_cc,signal_cc); xlim([250,550]);ylim([-0.1,0.2]);
cor = corrcoef(signal_cc(tim_cc>300&tim_cc<500),vel(tim_cc>300&tim_cc<500));
title(sprintf('Average windows: %d sec.; Corrcoef: %.3f', mvlist,cor(1,2)));

subplot(4,2,6);
plot(t(2:end)+dSecs(1),vel); hold on;
plot(tim_cc,signal_cc); xlim([1700,2000]);ylim([-0.1,0.2]);
index=tim_cc>1700&tim_cc<2000;
cor = corrcoef(signal_cc(index),vel(index));
title(sprintf('Average windows: %d sec.; Corrcoef: %.3f', mvlist,cor(1,2)));

subplot(4,2,7);
plot(t(2:end)+dSecs(1),vel); hold on;
plot(tim_cc,signal_cc); xlim([4200,4500]);ylim([-0.1,0.2]);
cor = corrcoef(signal_cc(tim_cc>4500&tim_cc<4800),vel(tim_cc>4500&tim_cc<4800));
title(sprintf('Average windows: %d sec.; Corrcoef: %.3f', mvlist,cor(1,2)));

subplot(4,2,8);
plot(t(2:end)+dSecs(1),vel); hold on;
plot(tim_cc,signal_cc); xlim([4700,5000]);ylim([-0.1,0.2]);
index=tim_cc>4700&tim_cc<5000;
cor = corrcoef(signal_cc(index),vel(index));
title(sprintf('Average windows: %d sec.; Corrcoef: %.3f', mvlist,cor(1,2)));