3 """A simple implementation of persistent homology on 2D images."""
5 __author__
= "Stefan Huber <shuber@sthu.org>"
15 def iter_neighbors(p
, w
, h
):
19 neigh
= [(y
+j
, x
+i
) for i
in [-1, 0, 1] for j
in [-1, 0, 1]]
21 # neigh = [(y-1, x), (y+1, x), (y, x-1), (y, x+1)]
36 # Get indices orderd by value from high to low. As a tie-breaker, we use
37 # (value, index) as key.
38 indices
= [(i
, j
) for i
in range(h
) for j
in range(w
)]
39 # We add p as a secondary key to have an unambiguous total order below when
40 # we enumerate neighboring cells of cells and consistency regarding
42 indices
.sort(key
=lambda p
: (get(im
, p
), p
), reverse
=True)
44 # Maintains the growing sets
45 uf
= union_find
.UnionFind()
49 def get_comp_birth(p
):
52 # Process pixels from high to low
53 for i
, p
in enumerate(indices
):
55 ni
= [uf
[q
] for q
in iter_neighbors(p
, w
, h
) if q
in uf
]
56 # Sort by (value, index) as key. Note that this is the same sorting
57 # order as for indices. Otherwise we have an inconsistence notion of
58 # the "older" component!
59 nc
= sorted([(get_comp_birth(q
), q
) for q
in set(ni
)], reverse
=True)
62 groups0
[p
] = (v
, v
, None)
70 # Merge all others with oldp
72 if uf
[q
] not in groups0
:
73 # print(i, ": Merge", uf[q], "with", oldp, "via", p)
74 groups0
[uf
[q
]] = (bl
, bl
-v
, p
)
77 groups0
= [(k
, groups0
[k
][0], groups0
[k
][1], groups0
[k
][2]) for k
in groups0
]
78 groups0
.sort(key
=lambda g
: g
[2], reverse
=True)