fixed spike-sync bug

fixed ugly bugs in code for computing multi-variate spike sync profile and
multi-variate spike sync value.

5 files changed, 47 insertions, 43 deletions

diff --git a/pyspike/DiscreteFunc.py b/pyspike/DiscreteFunc.py
index fe97bc2..caad290 100644
--- a/pyspike/DiscreteFunc.py
+++ b/pyspike/DiscreteFunc.py
@@ -170,10 +170,6 @@ expected."
                     start_ind, end_ind = get_indices(ival)
                     value += np.sum(self.y[start_ind:end_ind])
                     multiplicity += np.sum(self.mp[start_ind:end_ind])
-        if multiplicity == 0.0:
-            # empty profile, return spike sync of 1
-            value = 1.0
-            multiplicity = 1.0
         return (value, multiplicity)
 
     def avrg(self, interval=None, normalize=True):
@@ -190,7 +186,10 @@ expected."
         """
         val, mp = self.integral(interval)
         if normalize:
-            return val/mp
+            if mp > 0:
+                return val/mp
+            else:
+                return 1.0
         else:
             return val
 
diff --git a/pyspike/cython/cython_add.pyx b/pyspike/cython/cython_add.pyx
index 8da1e53..25f1181 100644
--- a/pyspike/cython/cython_add.pyx
+++ b/pyspike/cython/cython_add.pyx
@@ -182,8 +182,8 @@ def add_discrete_function_cython(double[:] x1, double[:] y1, double[:] mp1,
     cdef int index1 = 0
     cdef int index2 = 0
     cdef int index = 0
-    cdef int N1 = len(y1)
-    cdef int N2 = len(y2)
+    cdef int N1 = len(y1)-1
+    cdef int N2 = len(y2)-1
     x_new[0] = x1[0]
     while (index1+1 < N1) and (index2+1 < N2):
         if x1[index1+1] < x2[index2+1]:
@@ -206,20 +206,21 @@ def add_discrete_function_cython(double[:] x1, double[:] y1, double[:] mp1,
             y_new[index] = y1[index1] + y2[index2]
             mp_new[index] = mp1[index1] + mp2[index2]
     # one array reached the end -> copy the contents of the other to the end
-    if index1+1 < len(y1):
-        x_new[index+1:index+1+len(x1)-index1-1] = x1[index1+1:]
-        y_new[index+1:index+1+len(x1)-index1-1] = y1[index1+1:]
-        mp_new[index+1:index+1+len(x1)-index1-1] = mp1[index1+1:]
-        index += len(x1)-index1-1
-    elif index2+1 < len(y2):
-        x_new[index+1:index+1+len(x2)-index2-1] = x2[index2+1:]
-        y_new[index+1:index+1+len(x2)-index2-1] = y2[index2+1:]
-        mp_new[index+1:index+1+len(x2)-index2-1] = mp2[index2+1:]
-        index += len(x2)-index2-1
-    # else:  # both arrays reached the end simultaneously
-    #     x_new[index+1] = x1[-1]
-    #     y_new[index+1] = y1[-1] + y2[-1]
-    #     mp_new[index+1] = mp1[-1] + mp2[-1]
+    if index1+1 < N1:
+        x_new[index+1:index+1+N1-index1] = x1[index1+1:]
+        y_new[index+1:index+1+N1-index1] = y1[index1+1:]
+        mp_new[index+1:index+1+N1-index1] = mp1[index1+1:]
+        index += N1-index1
+    elif index2+1 < N2:
+        x_new[index+1:index+1+N2-index2] = x2[index2+1:]
+        y_new[index+1:index+1+N2-index2] = y2[index2+1:]
+        mp_new[index+1:index+1+N2-index2] = mp2[index2+1:]
+        index += N2-index2
+    else:  # both arrays reached the end simultaneously
+        x_new[index+1] = x1[index1+1]
+        y_new[index+1] = y1[index1+1] + y2[index2+1]
+        mp_new[index+1] = mp1[index1+1] + mp2[index2+1]
+        index += 1
 
     y_new[0] = y_new[1]
     mp_new[0] = mp_new[1]
diff --git a/pyspike/cython/cython_distances.pyx b/pyspike/cython/cython_distances.pyx
index 7dc9cb9..ac5f226 100644
--- a/pyspike/cython/cython_distances.pyx
+++ b/pyspike/cython/cython_distances.pyx
@@ -428,9 +428,4 @@ def coincidence_value_cython(double[:] spikes1, double[:] spikes2,
             mp += 2
             coinc += 2
 
-    if coinc == 0 and mp == 0:
-        # empty spike trains -> spike sync = 1 by definition
-        coinc = 1
-        mp = 1
-
     return coinc, mp
diff --git a/pyspike/cython/python_backend.py b/pyspike/cython/python_backend.py
index a1e3a65..6b7209a 100644
--- a/pyspike/cython/python_backend.py
+++ b/pyspike/cython/python_backend.py
@@ -560,7 +560,9 @@ def add_discrete_function_python(x1, y1, mp1, x2, y2, mp2):
     index1 = 0
     index2 = 0
     index = 0
-    while (index1+1 < len(y1)) and (index2+1 < len(y2)):
+    N1 = len(x1)-1
+    N2 = len(x2)-1
+    while (index1+1 < N1) and (index2+1 < N2):
         if x1[index1+1] < x2[index2+1]:
             index1 += 1
             index += 1
@@ -581,20 +583,21 @@ def add_discrete_function_python(x1, y1, mp1, x2, y2, mp2):
             y_new[index] = y1[index1] + y2[index2]
             mp_new[index] = mp1[index1] + mp2[index2]
     # one array reached the end -> copy the contents of the other to the end
-    if index1+1 < len(y1):
-        x_new[index+1:index+1+len(x1)-index1-1] = x1[index1+1:]
-        y_new[index+1:index+1+len(x1)-index1-1] = y1[index1+1:]
-        mp_new[index+1:index+1+len(x1)-index1-1] = mp1[index1+1:]
-        index += len(x1)-index1-1
-    elif index2+1 < len(y2):
-        x_new[index+1:index+1+len(x2)-index2-1] = x2[index2+1:]
-        y_new[index+1:index+1+len(x2)-index2-1] = y2[index2+1:]
-        mp_new[index+1:index+1+len(x2)-index2-1] = mp2[index2+1:]
-        index += len(x2)-index2-1
-    # else:  # both arrays reached the end simultaneously
-    #     x_new[index+1] = x1[-1]
-    #     y_new[index+1] = y1[-1] + y2[-1]
-    #     mp_new[index+1] = mp1[-1] + mp2[-1]
+    if index1+1 < N1:
+        x_new[index+1:index+1+N1-index1] = x1[index1+1:]
+        y_new[index+1:index+1+N1-index1] = y1[index1+1:]
+        mp_new[index+1:index+1+N1-index1] = mp1[index1+1:]
+        index += N1-index1
+    elif index2+1 < N2:
+        x_new[index+1:index+1+N2-index2] = x2[index2+1:]
+        y_new[index+1:index+1+N2-index2] = y2[index2+1:]
+        mp_new[index+1:index+1+N2-index2] = mp2[index2+1:]
+        index += N2-index2
+    else:  # both arrays reached the end simultaneously
+        x_new[index+1] = x1[-1]
+        y_new[index+1] = y1[-1] + y2[-1]
+        mp_new[index+1] = mp1[-1] + mp2[-1]
+        index += 1
 
     y_new[0] = y_new[1]
     mp_new[0] = mp_new[1]
diff --git a/pyspike/spike_sync.py b/pyspike/spike_sync.py
index 3dc29ff..802b98c 100644
--- a/pyspike/spike_sync.py
+++ b/pyspike/spike_sync.py
@@ -119,7 +119,10 @@ def spike_sync(spike_train1, spike_train2, interval=None, max_tau=None):
 
     """
     c, mp = _spike_sync_values(spike_train1, spike_train2, interval, max_tau)
-    return 1.0*c/mp
+    if mp == 0:
+        return 1.0
+    else:
+        return 1.0*c/mp
 
 
 ############################################################
@@ -187,7 +190,10 @@ def spike_sync_multi(spike_trains, indices=None, interval=None, max_tau=None):
         coincidence += c
         mp += m
 
-    return coincidence/mp
+    if mp == 0.0:
+        return 1.0
+    else:
+        return coincidence/mp
 
 
 ############################################################