function [XC, SYN, bins] = spk_xcorr(T1,T2,dt, maxlag,rmcenter)
%SPK_XCORR Spike cross-correlation. T1 is the spike times of the reference
%unit. T2 is the spike times of the follower unit. dt is the bin size.
%maxlag is the maximal time lag of the cross correlation function. rmcenter
%is a boolean indicating whether or not to remove the value at the center
%bin.

t1 = T1(T1 > maxlag & T1 < T1(end) - maxlag);

bins = -maxlag:dt:maxlag;
XC = zeros(size(bins'));

for i = 1:length(t1)
    b = t1(i) + bins;
    h = histc(T2,b);
    XC = XC + h(:);
end
XC = XC(1:end-1);

bins = bins(1:end-1)+dt/2;

if rmcenter
    mid = find(bins == 0);
    XC([mid-1 mid]) = 0;
end


end

