clear; clc;

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

mvlist = [1,2,5,10,30];

for imv = 1:size(mvlist,2)
    mvt = mvlist(imv);
    fps = 100; start = 100; stop = 6000;
    optictrack = csvread('ham37_02_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*mvt);

    % Open a file
    sr1 = 20000; sr2 = 80000;
    fhand1 = CEDS64Open('ACh_400Hz_1.0V_1.5ms.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*mvt);
    signal_cc = interp1(1:size(signal_d,1),signal_d,1:size(vel,1),'linear')';
    signal_cc = fillmissing(signal_cc,'linear');

    clf;
    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', mvt,cor(1,2)));

    savefig(sprintf('ACh002_%dsec_smooth',mvt));  
end

unloadlibrary ceds64int;