clear; clc
warning('off');

par = {'ham37_LeftDG_002','_190401_154740\',3,'Acetylcholine',[0:31];
       'ham37_RightCA1_003','_190401_134302\',1.5,'Acetylcholine',[33:63];
        };
ip = 1;
fps = 100; mvt = 1; otoffset = par{ip,3};
rootstr = par{ip,1};
optitrack = csvread([rootstr,par{ip,2},rootstr,'.csv'],7,1);
idx=find(optitrack(:,2)==0);
optitrack(idx,2:8)=nan(size(idx,1),7);  %%nan means make it empty
optitrack=fillmissing(optitrack,'linear'); %%fillmissing means fill all empty unit of whole matrix with methods
t=optitrack(:,1);
x=optitrack(:,6);
y=optitrack(:,7);
z=optitrack(:,8);
vel=sqrt((x(2:end)-x(1:end-1)).^2+(y(2:end)-y(1:end-1)).^2);
vel=smoothdata(vel,'movmean',1);

session=['.\',rootstr,par{ip,2}];
num_channels = 2;
fileinfo = dir(fullfile(session,'analogin.dat'));
num_samples = fileinfo.bytes/(num_channels * 2); % uint16 = 2 bytes
fid = fopen(fullfile(session,'analogin.dat'), 'r');
v = fread(fid, [num_channels, num_samples], 'uint16'); 
fclose(fid);

sr = 20000;
t_1 = (1:num_samples)/sr;

v = v * 3.3/2^16; % convert to volts
v = downsample(v',20);
tim = downsample(t_1',20);
tmpdata.samplingRate = 20000;
tmpdata.timestamps = tim;
tmpdata.data = v(:,2);
sigfil = bz_Filter(tmpdata,'passband',[0,10],'filter','butter');
idx = sigfil.data(:) > 0;
p = polyfit(sigfil.timestamps(idx), sigfil.data(idx), 3);
signal = sigfil.data(idx) - polyval(p, sigfil.timestamps(idx));
tim = tim(idx);

figure
subplot(12,2,[1:12])
plot(t(2:end)+otoffset,vel,'r'); xlim([min(tim),max(tim)]);
box off; axis off; text(0,max(vel)*1.0,"Speed");
subplot(12,2,[7:12])
plot(tim,signal,'b'); xlim([min(tim),max(tim)]);
box off; axis off; text(0,max(signal)*1.2,[par{ip,4}," (polynomial detrended)"]);
