From 619fdef014c44a89a7ecef9078905ee44e373b84 Mon Sep 17 00:00:00 2001 From: Mario Mulansky Date: Fri, 18 Dec 2015 14:37:45 +0100 Subject: bugfix for edge correction fixed bug within new edge correction (auxiliary spike was ignored in some cases) added regression test with 10000 random spike train sets --- pyspike/cython/cython_distances.pyx | 10 +++++++--- pyspike/cython/python_backend.py | 21 +++++++++++++++++++-- 2 files changed, 26 insertions(+), 5 deletions(-) (limited to 'pyspike') diff --git a/pyspike/cython/cython_distances.pyx b/pyspike/cython/cython_distances.pyx index 41baa60..6c6e7e5 100644 --- a/pyspike/cython/cython_distances.pyx +++ b/pyspike/cython/cython_distances.pyx @@ -184,7 +184,8 @@ def spike_distance_cython(double[:] t1, double[:] t2, N1 = len(t1) N2 = len(t2) - with nogil: # release the interpreter to allow multithreading + # with nogil: # release the interpreter to allow multithreading + if True: t_last = t_start # t_p1 = t_start # t_p2 = t_start @@ -193,6 +194,7 @@ def spike_distance_cython(double[:] t1, double[:] t2, t_aux1[1] = fmax(t_end, t1[N1-1]+(t1[N1-1]-t1[N1-2])) t_aux2[0] = fmin(t_start, t2[0]-(t2[1]-t2[0])) t_aux2[1] = fmax(t_end, t2[N2-1]+(t2[N2-1]-t2[N2-2])) + # print "aux spikes %.15f, %.15f ; %.15f, %.15f" % (t_aux1[0], t_aux1[1], t_aux2[0], t_aux2[1]) t_p1 = t_start if (t1[0] == t_start) else t_aux1[0] t_p2 = t_start if (t2[0] == t_start) else t_aux2[0] if t1[0] > t_start: @@ -207,7 +209,8 @@ def spike_distance_cython(double[:] t1, double[:] t2, else: t_f1 = t1[1] dt_f1 = get_min_dist_cython(t_f1, t2, N2, 0, t_aux2[0], t_aux2[1]) - dt_p1 = 0.0 + # dt_p1 = t_start-t_p2 # 0.0 + dt_p1 = get_min_dist_cython(t_p1, t2, N2, 0, t_aux2[0], t_aux2[1]) isi1 = t1[1]-t1[0] s1 = dt_p1 index1 = 0 @@ -223,7 +226,8 @@ def spike_distance_cython(double[:] t1, double[:] t2, else: t_f2 = t2[1] dt_f2 = get_min_dist_cython(t_f2, t1, N1, 0, t_aux1[0], t_aux1[1]) - dt_p2 = 0.0 + # dt_p2 = t_start-t_p1 # 0.0 + dt_p2 = get_min_dist_cython(t_p2, t1, N1, 0, t_aux1[0], t_aux1[1]) isi2 = t2[1]-t2[0] s2 = dt_p2 index2 = 0 diff --git a/pyspike/cython/python_backend.py b/pyspike/cython/python_backend.py index a5f7ae4..a007a7f 100644 --- a/pyspike/cython/python_backend.py +++ b/pyspike/cython/python_backend.py @@ -153,6 +153,8 @@ def spike_distance_python(spikes1, spikes2, t_start, t_end): t_p1 = t_start if (t1[0] == t_start) else t_aux1[0] t_p2 = t_start if (t2[0] == t_start) else t_aux2[0] + print "t_aux1", t_aux1, ", t_aux2:", t_aux2 + spike_events[0] = t_start if t1[0] > t_start: t_f1 = t1[0] @@ -163,7 +165,8 @@ def spike_distance_python(spikes1, spikes2, t_start, t_end): s1 = dt_p1 index1 = -1 else: - dt_p1 = 0.0 + # dt_p1 = t_start-t_p2 + dt_p1 = get_min_dist(t_p1, t2, 0, t_aux2[0], t_aux2[1]) t_f1 = t1[1] dt_f1 = get_min_dist(t_f1, t2, 0, t_aux2[0], t_aux2[1]) isi1 = t1[1]-t1[0] @@ -179,13 +182,18 @@ def spike_distance_python(spikes1, spikes2, t_start, t_end): s2 = dt_p2 index2 = -1 else: - dt_p2 = 0.0 + dt_p2 = t_start-t_p1 + dt_p2 = get_min_dist(t_p2, t1, 0, t_aux1[0], t_aux1[1]) t_f2 = t2[1] dt_f2 = get_min_dist(t_f2, t1, 0, t_aux1[0], t_aux1[1]) isi2 = t2[1]-t2[0] s2 = dt_p2 index2 = 0 + print "t_p1:", repr(t_p1), ", t_f1:", repr(t_f1), ", dt_p1:", repr(dt_p1), ", dt_f1:", repr(dt_f1) + print "t_p2:", repr(t_p2), ", t_f2:", repr(t_f2), ", dt_p2:", repr(dt_p2), ", dt_f2:", repr(dt_f2) + print "s1: ", repr(s1), ", s2:", repr(s2) + y_starts[0] = (s1*isi2 + s2*isi1) / (0.5*(isi1+isi2)**2) index = 1 @@ -276,6 +284,11 @@ def spike_distance_python(spikes1, spikes2, t_start, t_end): dt_f2 = dt_p2 isi2 = max(t_end-t2[N2-1], t2[N2-1]-t2[N2-2]) index += 1 + + print "t_p1:", repr(t_p1), ", t_f1:", repr(t_f1), ", dt_p1:", repr(dt_p1), ", dt_f1:", repr(dt_f1) + print "t_p2:", repr(t_p2), ", t_f2:", repr(t_f2), ", dt_p2:", repr(dt_p2), ", dt_f2:", repr(dt_f2) + print "s1: ", repr(s1), ", s2:", repr(s2) + # the last event is the interval end if spike_events[index-1] == t_end: index -= 1 @@ -288,6 +301,10 @@ def spike_distance_python(spikes1, spikes2, t_start, t_end): s2 = dt_f2 # *(t_end-t2[N2-1])/isi2 y_ends[index-1] = (s1*isi2 + s2*isi1) / (0.5*(isi1+isi2)**2) + print "t_p1:", repr(t_p1), ", t_f1:", repr(t_f1), ", dt_p1:", repr(dt_p1), ", dt_f1:", repr(dt_f1) + print "t_p2:", repr(t_p2), ", t_f2:", repr(t_f2), ", dt_p2:", repr(dt_p2), ", dt_f2:", repr(dt_f2) + print "s1: ", repr(s1), ", s2:", repr(s2) + # use only the data added above # could be less than original length due to equal spike times return spike_events[:index+1], y_starts[:index], y_ends[:index] -- cgit v1.2.3