diff --git a/gtda/externals/bindings/collapser_bindings.cpp b/gtda/externals/bindings/collapser_bindings.cpp index 90a542fe7..19a7f7a6a 100644 --- a/gtda/externals/bindings/collapser_bindings.cpp +++ b/gtda/externals/bindings/collapser_bindings.cpp @@ -6,6 +6,7 @@ #include #include +#include #include #include @@ -13,7 +14,7 @@ namespace py = pybind11; /* GUDHI Collapser required types */ using Filtration_value = float; -using Vertex_handle = uint32_t; +using Vertex_handle = int32_t; using Filtered_edge = std::tuple; using Filtered_edge_list = std::vector; @@ -29,7 +30,7 @@ using Filtration_values = std::vector; using COO_data = std::tuple; /* Dense matrix input types */ -using Distance_matrix = std::vector>; +using Distance_matrix_np = py::array_t; /* constants */ const Filtration_value filtration_max = @@ -86,13 +87,18 @@ PYBIND11_MODULE(gtda_collapser, m) { "collapses edges while preserving the persistent homology"); m.def("flag_complex_collapse_edges_coo", - [](Row_idx& row, Col_idx& col, Filtration_values& data, + [](py::array_t& row_, py::array_t& col_, + py::array_t& data_, Filtration_value thresh = filtration_max) { Filtered_edge_list graph; + Vertex_handle* row = (Vertex_handle*)row_.request().ptr; + Vertex_handle* col = (Vertex_handle*)col_.request().ptr; + Filtration_value* data = (Filtration_value*)data_.request().ptr; + /* Convert from COO input format to Filtered_edge_list */ /* Applying threshold to the input data */ - int size = data.size(); + int size = data_.size(); for (size_t k = 0; k < size; ++k) if (data[k] <= thresh) graph.push_back(Filtered_edge(row[k], col[k], data[k])); @@ -108,15 +114,15 @@ PYBIND11_MODULE(gtda_collapser, m) { "collapses edges while preserving the persistent homology"); m.def("flag_complex_collapse_edges_dense", - [](Distance_matrix& dm, Filtration_value thresh = filtration_max) { + [](Distance_matrix_np& dm, Filtration_value thresh = filtration_max) { Filtered_edge_list graph; /* Convert from dense format to Filtered edge list */ /* Applying threshold to the input data */ - for (size_t i = 0; i < dm.size(); i++) - for (size_t j = 0; j < dm[i].size(); j++) - if (j > i && (dm[i][j] <= thresh)) - graph.push_back(Filtered_edge(i, j, dm[i][j])); + for (size_t i = 0; i < dm.shape(0); i++) + for (size_t j = 0; j < dm.shape(1); j++) + if (j > i && (*(dm.data(i, j)) <= thresh)) + graph.push_back(Filtered_edge(i, j, *(dm.data(i, j)))); /* Start collapser */ auto vec_triples = diff --git a/gtda/externals/bindings/ripser_bindings.cpp b/gtda/externals/bindings/ripser_bindings.cpp index a6ae892ed..e28b5fe05 100644 --- a/gtda/externals/bindings/ripser_bindings.cpp +++ b/gtda/externals/bindings/ripser_bindings.cpp @@ -7,8 +7,10 @@ #include // PYBIND11 +#include #include #include +#include namespace py = pybind11; @@ -32,21 +34,22 @@ PYBIND11_MODULE(gtda_ripser, m) { .def_readwrite("num_edges", &ripserResults::num_edges); m.def("rips_dm", - [](std::vector D, int N, int modulus, int dim_max, + [](py::array_t& D, int N, int modulus, int dim_max, float threshold, int do_cocycles) { - ripserResults ret = - rips_dm(&D[0], N, modulus, dim_max, threshold, do_cocycles); + ripserResults ret = rips_dm((float*)D.request().ptr, N, modulus, + dim_max, threshold, do_cocycles); return ret; }, "D"_a, "N"_a, "modulus"_a, "dim_max"_a, "threshold"_a, "do_cocycles"_a, "ripser distance matrix"); m.def("rips_dm_sparse", - [](std::vector I, std::vector J, std::vector V, + [](py::array_t& I, py::array_t& J, py::array_t& V, int NEdges, int N, int modulus, int dim_max, float threshold, int do_cocycles) { ripserResults ret = - rips_dm_sparse(&I[0], &J[0], &V[0], NEdges, N, modulus, dim_max, - threshold, do_cocycles); + rips_dm_sparse((int*)I.request().ptr, (int*)J.request().ptr, + (float*)V.request().ptr, NEdges, N, modulus, + dim_max, threshold, do_cocycles); return ret; }, "I"_a, "J"_a, "V"_a, "NEdges"_a, "N"_a, "modulus"_a, "dim_max"_a, diff --git a/gtda/externals/python/ripser_interface.py b/gtda/externals/python/ripser_interface.py index 077b277cb..bc425f659 100644 --- a/gtda/externals/python/ripser_interface.py +++ b/gtda/externals/python/ripser_interface.py @@ -21,10 +21,7 @@ def DRFDM(DParam, maxHomDim, thresh=-1, coeff=2, do_cocycles=0): else: ret = gtda_ripser_coeff.rips_dm(DParam, DParam.shape[0], coeff, maxHomDim, thresh, do_cocycles) - ret_rips = {} - ret_rips.update({"births_and_deaths_by_dim": ret.births_and_deaths_by_dim}) - ret_rips.update({"num_edges": ret.num_edges}) - return ret_rips + return ret def DRFDMSparse(I, J, V, N, maxHomDim, thresh=-1, coeff=2, do_cocycles=0): @@ -34,10 +31,7 @@ def DRFDMSparse(I, J, V, N, maxHomDim, thresh=-1, coeff=2, do_cocycles=0): else: ret = gtda_ripser_coeff.rips_dm_sparse(I, J, V, I.size, N, coeff, maxHomDim, thresh, do_cocycles) - ret_rips = {} - ret_rips.update({"births_and_deaths_by_dim": ret.births_and_deaths_by_dim}) - ret_rips.update({"num_edges": ret.num_edges}) - return ret_rips + return ret def dpoint2pointcloud(X, i, metric): @@ -295,19 +289,19 @@ def ripser(X, maxdim=1, thresh=np.inf, coeff=2, metric="euclidean", ) else: # Only consider strict upper diagonal - idx = np.triu_indices(n_points, 1) - DParam = dm[idx].astype(np.float32) + DParam = dm[np.invert(np.tri(n_points, k=0, dtype=np.bool))].astype( + np.float32).flatten() res = DRFDM(DParam, maxdim, thresh, coeff) # Unwrap persistence diagrams - dgms = res["births_and_deaths_by_dim"] + dgms = res.births_and_deaths_by_dim for dim in range(len(dgms)): N = int(len(dgms[dim]) / 2) dgms[dim] = np.reshape(np.array(dgms[dim]), [N, 2]) ret = { "dgms": dgms, - "num_edges": res["num_edges"], + "num_edges": res.num_edges, "dperm2all": dperm2all, "idx_perm": idx_perm, "r_cover": r_cover,