Skip to content

Commit

Permalink
Revert #465 due to different handling of COO matrices by ripser follo…
Browse files Browse the repository at this point in the history
…wing #530 (#532)

* Revert #465 following #530

* Add copyright notice inherited from ripser.py following @MonkeyBreaker's comment
  • Loading branch information
ulupo authored Nov 16, 2020
2 parents 5675d2a + 21c33e0 commit 8c29934
Showing 1 changed file with 45 additions and 56 deletions.
101 changes: 45 additions & 56 deletions gtda/externals/python/ripser_interface.py
Original file line number Diff line number Diff line change
@@ -1,3 +1,23 @@
"""
MIT License
Copyright (c) 2018 Christopher Tralie and Nathaniel Saul
Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
"""

import gc
from warnings import warn

Expand All @@ -9,13 +29,6 @@
from ..modules import gtda_ripser, gtda_ripser_coeff, gtda_collapser


def _lexsort_coo_data(row, col, data):
lex_sort_idx = np.lexsort((col, row))
row, col, data = \
row[lex_sort_idx], col[lex_sort_idx], data[lex_sort_idx]
return row, col, data


def DRFDM(DParam, maxHomDim, thresh=-1, coeff=2, do_cocycles=0):
if coeff == 2:
ret = gtda_ripser.rips_dm(DParam, DParam.shape[0], coeff, maxHomDim,
Expand Down Expand Up @@ -204,26 +217,20 @@ def ripser(X, maxdim=1, thresh=np.inf, coeff=2, metric="euclidean",
"""
if n_perm and sparse.issparse(X):
raise Exception(
"Greedy permutation is not supported for sparse distance matrices"
)
raise Exception("Greedy permutation is not supported for sparse "
"distance matrices")
if n_perm and n_perm > X.shape[0]:
raise Exception(
"Number of points in greedy permutation is greater"
" than number of points in the point cloud"
)
raise Exception("Number of points in greedy permutation is greater "
"than number of points in the point cloud")
if n_perm and n_perm < 0:
raise Exception(
"Should be a strictly positive number of points in the greedy "
"permutation"
)
raise Exception("There should be a strictly positive number of points "
"in the greedy permutation")

idx_perm = np.arange(X.shape[0])
r_cover = 0.0
if n_perm:
idx_perm, lambdas, dperm2all = get_greedy_perm(
X, n_perm=n_perm, metric=metric
)
idx_perm, lambdas, dperm2all = \
get_greedy_perm(X, n_perm=n_perm, metric=metric)
r_cover = lambdas[-1]
dm = dperm2all[:, idx_perm]
else:
Expand All @@ -234,7 +241,6 @@ def ripser(X, maxdim=1, thresh=np.inf, coeff=2, metric="euclidean",
dperm2all = None

n_points = max(dm.shape)
sort_coo = True
if (dm.diagonal() != 0).any():
if collapse_edges:
warn("Edge collapses are not supported when any of the diagonal "
Expand All @@ -246,11 +252,9 @@ def ripser(X, maxdim=1, thresh=np.inf, coeff=2, metric="euclidean",
# format, because currently that's the only format that handles
# nonzero births
dm = sparse.coo_matrix(dm)
sort_coo = False

if sparse.issparse(dm) or collapse_edges:
if collapse_edges:
sort_coo = True
if not sparse.issparse(dm):
row, col, data = \
gtda_collapser.flag_complex_collapse_edges_dense(dm,
Expand All @@ -262,31 +266,17 @@ def ripser(X, maxdim=1, thresh=np.inf, coeff=2, metric="euclidean",
coo.col,
coo.data,
thresh)
else:
if sparse.isspmatrix_coo(dm):
# If the matrix is already COO, we need to order the row and
# column indices lexicographically to avoid errors. See
# https://github.com/scikit-tda/ripser.py/issues/103
row, col, data = dm.row, dm.col, dm.data
else:
coo = dm.tocoo()
row, col, data = coo.row, coo.col, coo.data
sort_coo = False

if sort_coo:
row, col, data = _lexsort_coo_data(np.asarray(row),
np.asarray(col),
np.asarray(data))

res = DRFDMSparse(
row.astype(dtype=np.int32, order="C"),
col.astype(dtype=np.int32, order="C"),
np.array(data, dtype=np.float32, order="C"),
n_points,
maxdim,
thresh,
coeff
)
elif sparse.issparse(dm):
coo = dm.tocoo()
row, col, data = coo.row, coo.col, coo.data

res = DRFDMSparse(np.asarray(row, dtype=np.int32, order="C"),
np.asarray(col, dtype=np.int32, order="C"),
np.asarray(data, dtype=np.float32, order="C"),
n_points,
maxdim,
thresh,
coeff)
else:
# Only consider strict upper diagonal
DParam = squareform(dm, checks=False).astype(np.float32)
Expand All @@ -301,11 +291,10 @@ def ripser(X, maxdim=1, thresh=np.inf, coeff=2, metric="euclidean",
N = int(len(dgms[dim]) / 2)
dgms[dim] = np.reshape(np.array(dgms[dim]), [N, 2])

ret = {
"dgms": dgms,
"num_edges": res.num_edges,
"dperm2all": dperm2all,
"idx_perm": idx_perm,
"r_cover": r_cover,
}
ret = {"dgms": dgms,
"num_edges": res.num_edges,
"dperm2all": dperm2all,
"idx_perm": idx_perm,
"r_cover": r_cover}

return ret

0 comments on commit 8c29934

Please sign in to comment.