up

src/totalpersistence/totalpersistence.py

@@ -9,14 +9,16 @@
     kercoker_bars,
     matrix_size_from_condensed,
 )
-from IPython import embed
 
 
 def totalpersistence(X, Y, f, maxdim=1, eps=0, tol=1e-11, perturb=1e-7):
     dX = general_position_distance_matrix(X, perturb)
     dY = general_position_distance_matrix(Y, perturb)
 
     data, dgm, dgmX, dgmY = kercoker_via_cone(dX, dY, f, maxdim, eps, tol, perturb)
+
+    distance_bottleneck, matching = persim.bottleneck(dgm_clean, dgm_noisy, matching=True)
+
     return data, dgm, dgmX, dgmY
 
 
@@ -60,7 +62,7 @@ def kercoker_via_cone(dX, dY, f, maxdim=1, cone_eps=0, tol=1e-11):
 
     dgmX = ripser(squareform(dX), distance_matrix=True, maxdim=maxdim)["dgms"]
     dgmY = ripser(squareform(dY), distance_matrix=True, maxdim=maxdim)["dgms"]
-    dgm = ripser(D, maxdim=maxdim, distance_matrix=True)["dgms"]
+    cone_dgm = ripser(D, maxdim=maxdim, distance_matrix=True)["dgms"]
 
-    bars = kercoker_bars(dgm, dgmX, dgmY, cone_eps, tol)
-    return bars, dgm, dgmX, dgmY
+    coker_dgm, ker_dgm = kercoker_bars(cone_dgm, dgmX, dgmY, cone_eps, tol)
+    return coker_dgm, ker_dgm, cone_dgm, dgmX, dgmY

src/totalpersistence/utils.py

@@ -49,15 +49,24 @@ def conematrix(DX, DY, DY_fy, eps):
     return D
 
 
+def format_bars(bars):
+    bars = [np.array(b) for b in bars]
+    lens = list(map(len, bars))
+    for i in range(len(bars)):
+        if all(l == 0 for l in lens[i:]):
+            bars = bars[:i]
+            break
+    return bars
+
+
 def kercoker_bars(dgm, dgmX, dgmY, cone_eps, tol=1e-11):
     """
     Find cokernel and kernel bars in the persistence diagram.
     TODO: optimize
     """
-    data = []
+    coker_dgm = [[] for _ in range(len(dgm))]
+    ker_dgm = [[] for _ in range(len(dgm))]
     for k in range(len(dgm)):
-        coker = []
-        ker = []
         for r in dgm[k]:
             b, d = r
             if d > cone_eps + tol:
@@ -66,29 +75,26 @@ def kercoker_bars(dgm, dgmX, dgmY, cone_eps, tol=1e-11):
                 # d_c = d_y_i
                 m = findclose(b, dgmY[k][:, 0], tol) & findclose(d, dgmY[k][:, 1], tol)
                 if sum(m):
-                    coker.append((b, d))
+                    coker_dgm[k].append((b, d))
 
                 # b_c = b_y_i
                 # d_c = b_x_j
                 if any(findclose(b, dgmY[k][:, 0], tol)) and any(findclose(d, dgmX[k][:, 0], tol)):
-                    coker.append((b, d))
+                    coker_dgm[k].append((b, d))
 
                 # ker
                 if k > 0:
                     # b_c = b_x_i (dim-1)
                     # d_c = d_x_i (dim-1)
                     m = findclose(b, dgmX[k - 1][:, 0], tol) & findclose(d, dgmX[k - 1][:, 1], tol)
                     if sum(m):
-                        ker.append((b, d))
+                        ker_dgm[k - 1].append((b, d))
 
                     # b_c = d_y_i (dim-1)
                     # d_c = d_x_j (dim-1)
                     if any(findclose(b, dgmY[k - 1][:, 1], tol)) and any(findclose(d, dgmX[k - 1][:, 1], tol)):
-                        ker.append((b, d))
-
-        for c in coker:
-            data.append({"dim": k, "set": "coker", "b": c[0], "d": c[1]})
-        for c in ker:
-            data.append({"dim": k - 1, "set": "ker", "b": c[0], "d": c[1]})
+                        ker_dgm[k - 1].append((b, d))
 
-    return data
+    coker_dgm = format_bars(coker_dgm)
+    ker_dgm = format_bars(ker_dgm)
+    return coker_dgm, ker_dgm

tests/test_kercoker.py

@@ -13,26 +13,23 @@ def test_kercoker_via_cone():
     dY = general_position_distance_matrix(Y, perturb=1e-13)
 
     # Run the function
-    result = kercoker_via_cone(dX, dY, f, maxdim=2, cone_eps=0)[0]
-    print(result)
+    coker_dgm, ker_dgm, cone_dgm, dgmX, dgmY = kercoker_via_cone(dX, dY, f, maxdim=2, cone_eps=0)
 
     # Expected output
-    expected = [
-        {"dim": 0, "set": "coker", "b": 0.0, "d": float("inf")},
-        {"dim": 0, "set": "ker", "b": 0.5, "d": 1.0},
-        {"dim": 0, "set": "ker", "b": 0.5, "d": 1.0},
-        {"dim": 0, "set": "ker", "b": 0.0, "d": 1.0},
-        {"dim": 0, "set": "ker", "b": 0.0, "d": 1.0},
-        {"dim": 0, "set": "ker", "b": 0.0, "d": 1.0},
-        {"dim": 1, "set": "ker", "b": 1.0, "d": 1.4142135381698608},
-        {"dim": 1, "set": "ker", "b": 1.0, "d": 1.4142135381698608},
+    # for coker in [array([[ 0., inf]])]
+    # for ker in [array([[ 0.5,  1. ], [ 0.5,  1. ], [ 0. ,  1. ], [ 0. ,  1. ], [ 0. ,  1. ]),
+    # array([[ 1. ,  1.41421356], [ 1. ,  1.41421356]])]
+
+    expected_coker = np.array([[0.0, float("inf")]])
+    expected_ker = [
+        np.array([[0.5, 1.0], [0.5, 1.0], [0.0, 1.0], [0.0, 1.0], [0.0, 1.0]]),
+        np.array([[1.0, 1.41421356], [1.0, 1.41421356]]),
     ]
     tol = 1e-10
 
     # Compare results
-    assert len(result) == len(expected)
-    for r, e in zip(result, expected):
-        assert r["dim"] == e["dim"]
-        assert r["set"] == e["set"]
-        assert abs(r["b"] - e["b"]) < tol
-        assert abs(r["d"] - e["d"]) < tol if e["d"] != float("inf") else r["d"] == float("inf")
+    assert len(coker_dgm) == 1
+    assert np.allclose(coker_dgm[0], expected_coker, atol=tol)
+    assert len(ker_dgm) == 2
+    for r, e in zip(ker_dgm, expected_ker):
+        assert np.allclose(r, e, atol=tol)