diff options
Diffstat (limited to 'src')
21 files changed, 467 insertions, 150 deletions
diff --git a/src/Simplex_tree/example/simple_simplex_tree.cpp b/src/Simplex_tree/example/simple_simplex_tree.cpp index 4353939f..47ea7e36 100644 --- a/src/Simplex_tree/example/simple_simplex_tree.cpp +++ b/src/Simplex_tree/example/simple_simplex_tree.cpp @@ -166,10 +166,19 @@ int main(int argc, char* const argv[]) { // ++ GENERAL VARIABLE SET std::cout << "********************************************************************\n"; - // Display the Simplex_tree - Can not be done in the middle of 2 inserts std::cout << "* The complex contains " << simplexTree.num_simplices() << " simplices\n"; std::cout << " - dimension " << simplexTree.dimension() << "\n"; - std::cout << "* Iterator on Simplices in the filtration, with [filtration value]:\n"; + std::cout << "* Iterator on simplices, with [filtration value]:\n"; + for (Simplex_tree::Simplex_handle f_simplex : simplexTree.complex_simplex_range()) { + std::cout << " " + << "[" << simplexTree.filtration(f_simplex) << "] "; + for (auto vertex : simplexTree.simplex_vertex_range(f_simplex)) std::cout << "(" << vertex << ")"; + std::cout << std::endl; + } + + std::cout << "********************************************************************\n"; + // Can not be done in the middle of 2 inserts + std::cout << "* Iterator on simplices sorted by filtration values, with [filtration value]:\n"; for (auto f_simplex : simplexTree.filtration_simplex_range()) { std::cout << " " << "[" << simplexTree.filtration(f_simplex) << "] "; diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree.h b/src/Simplex_tree/include/gudhi/Simplex_tree.h index 76608008..7b39a500 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree.h @@ -1317,6 +1317,8 @@ class Simplex_tree { * \post Some simplex tree functions require the filtration to be valid. `make_filtration_non_decreasing()` * function is not launching `initialize_filtration()` but returns the filtration modification information. If the * complex has changed , please call `initialize_filtration()` to recompute it. + * + * If a simplex has a `NaN` filtration value, it is considered lower than any other defined filtration value. */ bool make_filtration_non_decreasing() { bool modified = false; @@ -1347,7 +1349,9 @@ class Simplex_tree { }); Filtration_value max_filt_border_value = filtration(*max_border); - if (simplex.second.filtration() < max_filt_border_value) { + // Replacing if(f<max) with if(!(f>=max)) would mean that if f is NaN, we replace it with the max of the children. + // That seems more useful than keeping NaN. + if (!(simplex.second.filtration() >= max_filt_border_value)) { // Store the filtration modification information modified = true; simplex.second.assign_filtration(max_filt_border_value); diff --git a/src/Simplex_tree/test/CMakeLists.txt b/src/Simplex_tree/test/CMakeLists.txt index 8b9163f5..cf2b0153 100644 --- a/src/Simplex_tree/test/CMakeLists.txt +++ b/src/Simplex_tree/test/CMakeLists.txt @@ -28,3 +28,9 @@ if (TBB_FOUND) target_link_libraries(Simplex_tree_ctor_and_move_test_unit ${TBB_LIBRARIES}) endif() gudhi_add_boost_test(Simplex_tree_ctor_and_move_test_unit) + +add_executable ( Simplex_tree_make_filtration_non_decreasing_test_unit simplex_tree_make_filtration_non_decreasing_unit_test.cpp ) +if (TBB_FOUND) + target_link_libraries(Simplex_tree_make_filtration_non_decreasing_test_unit ${TBB_LIBRARIES}) +endif() +gudhi_add_boost_test(Simplex_tree_make_filtration_non_decreasing_test_unit) diff --git a/src/Simplex_tree/test/simplex_tree_make_filtration_non_decreasing_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_make_filtration_non_decreasing_unit_test.cpp new file mode 100644 index 00000000..4697ec05 --- /dev/null +++ b/src/Simplex_tree/test/simplex_tree_make_filtration_non_decreasing_unit_test.cpp @@ -0,0 +1,148 @@ +/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT. + * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. + * Author(s): Vincent Rouvreau + * + * Copyright (C) 2020 Inria + * + * Modification(s): + * - YYYY/MM Author: Description of the modification + */ + +#include <iostream> +#include <limits> // for NaN +#include <cmath> // for isNaN + +#define BOOST_TEST_DYN_LINK +#define BOOST_TEST_MODULE "simplex_tree_make_filtration_non_decreasing" +#include <boost/test/unit_test.hpp> +#include <boost/mpl/list.hpp> + +// ^ +// /!\ Nothing else from Simplex_tree shall be included to test includes are well defined. +#include "gudhi/Simplex_tree.h" + +using namespace Gudhi; + +typedef boost::mpl::list<Simplex_tree<>, Simplex_tree<Simplex_tree_options_fast_persistence>> list_of_tested_variants; + +BOOST_AUTO_TEST_CASE_TEMPLATE(make_filtration_non_decreasing, typeST, list_of_tested_variants) { + typeST st; + + st.insert_simplex_and_subfaces({2, 1, 0}, 2.0); + st.insert_simplex_and_subfaces({3, 0}, 2.0); + st.insert_simplex_and_subfaces({3, 4, 5}, 2.0); + + /* Inserted simplex: */ + /* 1 */ + /* o */ + /* /X\ */ + /* o---o---o---o */ + /* 2 0 3\X/4 */ + /* o */ + /* 5 */ + + std::cout << "Check default insertion ensures the filtration values are non decreasing" << std::endl; + BOOST_CHECK(!st.make_filtration_non_decreasing()); + + // Because of non decreasing property of simplex tree, { 0 } , { 1 } and { 0, 1 } are going to be set from value 2.0 + // to 1.0 + st.insert_simplex_and_subfaces({0, 1, 6, 7}, 1.0); + + // Inserted simplex: + // 1 6 + // o---o + // /X\7/ + // o---o---o---o + // 2 0 3\X/4 + // o + // 5 + + std::cout << "Check default second insertion ensures the filtration values are non decreasing" << std::endl; + BOOST_CHECK(!st.make_filtration_non_decreasing()); + + // Copy original simplex tree + typeST st_copy = st; + + // Modify specific values for st to become like st_copy thanks to make_filtration_non_decreasing + st.assign_filtration(st.find({0,1,6,7}), 0.8); + st.assign_filtration(st.find({0,1,6}), 0.9); + st.assign_filtration(st.find({0,6}), 0.6); + st.assign_filtration(st.find({3,4,5}), 1.2); + st.assign_filtration(st.find({3,4}), 1.1); + st.assign_filtration(st.find({4,5}), 1.99); + + std::cout << "Check the simplex_tree is rolled back in case of decreasing filtration values" << std::endl; + BOOST_CHECK(st.make_filtration_non_decreasing()); + BOOST_CHECK(st == st_copy); + + // Other simplex tree + typeST st_other; + st_other.insert_simplex_and_subfaces({2, 1, 0}, 3.0); // This one is different from st + st_other.insert_simplex_and_subfaces({3, 0}, 2.0); + st_other.insert_simplex_and_subfaces({3, 4, 5}, 2.0); + st_other.insert_simplex_and_subfaces({0, 1, 6, 7}, 1.0); + + // Modify specific values for st to become like st_other thanks to make_filtration_non_decreasing + st.assign_filtration(st.find({2}), 3.0); + // By modifying just the simplex {2} + // {0,1,2}, {1,2} and {0,2} will be modified + + std::cout << "Check the simplex_tree is repaired in case of decreasing filtration values" << std::endl; + BOOST_CHECK(st.make_filtration_non_decreasing()); + BOOST_CHECK(st == st_other); + + // Modify specific values for st still to be non-decreasing + st.assign_filtration(st.find({0,1,2}), 10.0); + st.assign_filtration(st.find({0,2}), 9.0); + st.assign_filtration(st.find({0,1,6,7}), 50.0); + st.assign_filtration(st.find({0,1,6}), 49.0); + st.assign_filtration(st.find({0,1,7}), 48.0); + // Other copy simplex tree + typeST st_other_copy = st; + + std::cout << "Check the simplex_tree is not modified in case of non-decreasing filtration values" << std::endl; + BOOST_CHECK(!st.make_filtration_non_decreasing()); + BOOST_CHECK(st == st_other_copy); + +} + +BOOST_AUTO_TEST_CASE_TEMPLATE(make_filtration_non_decreasing_on_nan_values, typeST, list_of_tested_variants) { + typeST st; + + st.insert_simplex_and_subfaces({2, 1, 0}, std::numeric_limits<double>::quiet_NaN()); + st.insert_simplex_and_subfaces({3, 0}, std::numeric_limits<double>::quiet_NaN()); + st.insert_simplex_and_subfaces({3, 4, 5}, std::numeric_limits<double>::quiet_NaN()); + + /* Inserted simplex: */ + /* 1 */ + /* o */ + /* /X\ */ + /* o---o---o---o */ + /* 2 0 3\X/4 */ + /* o */ + /* 5 */ + + std::cout << "SPECIFIC CASE:" << std::endl; + std::cout << "Insertion with NaN values does not ensure the filtration values are non decreasing" << std::endl; + st.make_filtration_non_decreasing(); + + std::cout << "Check all filtration values are NaN" << std::endl; + for (auto f_simplex : st.complex_simplex_range()) { + BOOST_CHECK(std::isnan(st.filtration(f_simplex))); + } + + st.assign_filtration(st.find({0}), 0.); + st.assign_filtration(st.find({1}), 0.); + st.assign_filtration(st.find({2}), 0.); + st.assign_filtration(st.find({3}), 0.); + st.assign_filtration(st.find({4}), 0.); + st.assign_filtration(st.find({5}), 0.); + + std::cout << "Check make_filtration_non_decreasing is modifying the simplicial complex" << std::endl; + BOOST_CHECK(st.make_filtration_non_decreasing()); + + std::cout << "Check all filtration values are now defined" << std::endl; + for (auto f_simplex : st.complex_simplex_range()) { + BOOST_CHECK(!std::isnan(st.filtration(f_simplex))); + } +}
\ No newline at end of file diff --git a/src/Simplex_tree/test/simplex_tree_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_unit_test.cpp index 58bfa8db..e739ad0a 100644 --- a/src/Simplex_tree/test/simplex_tree_unit_test.cpp +++ b/src/Simplex_tree/test/simplex_tree_unit_test.cpp @@ -791,90 +791,6 @@ BOOST_AUTO_TEST_CASE(non_contiguous) { BOOST_CHECK(++i == std::end(b)); } -BOOST_AUTO_TEST_CASE(make_filtration_non_decreasing) { - std::cout << "********************************************************************" << std::endl; - std::cout << "MAKE FILTRATION NON DECREASING" << std::endl; - typedef Simplex_tree<> typeST; - typeST st; - - st.insert_simplex_and_subfaces({2, 1, 0}, 2.0); - st.insert_simplex_and_subfaces({3, 0}, 2.0); - st.insert_simplex_and_subfaces({3, 4, 5}, 2.0); - - /* Inserted simplex: */ - /* 1 */ - /* o */ - /* /X\ */ - /* o---o---o---o */ - /* 2 0 3\X/4 */ - /* o */ - /* 5 */ - - std::cout << "Check default insertion ensures the filtration values are non decreasing" << std::endl; - BOOST_CHECK(!st.make_filtration_non_decreasing()); - - // Because of non decreasing property of simplex tree, { 0 } , { 1 } and { 0, 1 } are going to be set from value 2.0 - // to 1.0 - st.insert_simplex_and_subfaces({0, 1, 6, 7}, 1.0); - - // Inserted simplex: - // 1 6 - // o---o - // /X\7/ - // o---o---o---o - // 2 0 3\X/4 - // o - // 5 - - std::cout << "Check default second insertion ensures the filtration values are non decreasing" << std::endl; - BOOST_CHECK(!st.make_filtration_non_decreasing()); - - // Copy original simplex tree - typeST st_copy = st; - - // Modify specific values for st to become like st_copy thanks to make_filtration_non_decreasing - st.assign_filtration(st.find({0,1,6,7}), 0.8); - st.assign_filtration(st.find({0,1,6}), 0.9); - st.assign_filtration(st.find({0,6}), 0.6); - st.assign_filtration(st.find({3,4,5}), 1.2); - st.assign_filtration(st.find({3,4}), 1.1); - st.assign_filtration(st.find({4,5}), 1.99); - - std::cout << "Check the simplex_tree is rolled back in case of decreasing filtration values" << std::endl; - BOOST_CHECK(st.make_filtration_non_decreasing()); - BOOST_CHECK(st == st_copy); - - // Other simplex tree - typeST st_other; - st_other.insert_simplex_and_subfaces({2, 1, 0}, 3.0); // This one is different from st - st_other.insert_simplex_and_subfaces({3, 0}, 2.0); - st_other.insert_simplex_and_subfaces({3, 4, 5}, 2.0); - st_other.insert_simplex_and_subfaces({0, 1, 6, 7}, 1.0); - - // Modify specific values for st to become like st_other thanks to make_filtration_non_decreasing - st.assign_filtration(st.find({2}), 3.0); - // By modifying just the simplex {2} - // {0,1,2}, {1,2} and {0,2} will be modified - - std::cout << "Check the simplex_tree is repaired in case of decreasing filtration values" << std::endl; - BOOST_CHECK(st.make_filtration_non_decreasing()); - BOOST_CHECK(st == st_other); - - // Modify specific values for st still to be non-decreasing - st.assign_filtration(st.find({0,1,2}), 10.0); - st.assign_filtration(st.find({0,2}), 9.0); - st.assign_filtration(st.find({0,1,6,7}), 50.0); - st.assign_filtration(st.find({0,1,6}), 49.0); - st.assign_filtration(st.find({0,1,7}), 48.0); - // Other copy simplex tree - typeST st_other_copy = st; - - std::cout << "Check the simplex_tree is not modified in case of non-decreasing filtration values" << std::endl; - BOOST_CHECK(!st.make_filtration_non_decreasing()); - BOOST_CHECK(st == st_other_copy); - -} - typedef boost::mpl::list<boost::adjacency_list<boost::setS, boost::vecS, boost::directedS, boost::property<vertex_filtration_t, double>, diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 20e72a5f..22af3ec9 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -56,6 +56,7 @@ if(PYTHONINTERP_FOUND) # Modules that should not be auto-imported in __init__.py set(GUDHI_PYTHON_MODULES_EXTRA "${GUDHI_PYTHON_MODULES_EXTRA}'representations', ") set(GUDHI_PYTHON_MODULES_EXTRA "${GUDHI_PYTHON_MODULES_EXTRA}'wasserstein', ") + set(GUDHI_PYTHON_MODULES_EXTRA "${GUDHI_PYTHON_MODULES_EXTRA}'point_cloud', ") add_gudhi_debug_info("Python version ${PYTHON_VERSION_STRING}") add_gudhi_debug_info("Cython version ${CYTHON_VERSION}") @@ -226,6 +227,7 @@ if(PYTHONINTERP_FOUND) file(COPY "gudhi/persistence_graphical_tools.py" DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/gudhi") file(COPY "gudhi/representations" DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/gudhi/") file(COPY "gudhi/wasserstein.py" DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/gudhi") + file(COPY "gudhi/point_cloud" DESTINATION "${CMAKE_CURRENT_BINARY_DIR}/gudhi") add_custom_command( OUTPUT gudhi.so @@ -404,6 +406,9 @@ if(PYTHONINTERP_FOUND) add_gudhi_py_test(test_representations) endif() + # Time Delay + add_gudhi_py_test(test_time_delay) + # Documentation generation is available through sphinx - requires all modules if(SPHINX_PATH) if(MATPLOTLIB_FOUND) diff --git a/src/python/doc/point_cloud.rst b/src/python/doc/point_cloud.rst index d668428a..c0d4b303 100644 --- a/src/python/doc/point_cloud.rst +++ b/src/python/doc/point_cloud.rst @@ -20,3 +20,11 @@ Subsampling :members: :special-members: :show-inheritance: + +TimeDelayEmbedding +------------------ + +.. autoclass:: gudhi.point_cloud.timedelay.TimeDelayEmbedding + :members: + :special-members: __call__ + diff --git a/src/python/example/alpha_complex_from_points_example.py b/src/python/example/alpha_complex_from_points_example.py index 844d7a82..73faf17c 100755 --- a/src/python/example/alpha_complex_from_points_example.py +++ b/src/python/example/alpha_complex_from_points_example.py @@ -46,8 +46,14 @@ if simplex_tree.find([4]): else: print("[4] Not found...") +# Some insertions, simplex_tree needs to initialize filtrations +simplex_tree.initialize_filtration() + print("dimension=", simplex_tree.dimension()) -print("filtrations=", simplex_tree.get_filtration()) +print("filtrations=") +for simplex_with_filtration in simplex_tree.get_filtration(): + print("(%s, %.2f)" % tuple(simplex_with_filtration)) + print("star([0])=", simplex_tree.get_star([0])) print("coface([0], 1)=", simplex_tree.get_cofaces([0], 1)) diff --git a/src/python/example/rips_complex_from_points_example.py b/src/python/example/rips_complex_from_points_example.py index 59d8a261..c05703c6 100755 --- a/src/python/example/rips_complex_from_points_example.py +++ b/src/python/example/rips_complex_from_points_example.py @@ -22,6 +22,9 @@ rips = gudhi.RipsComplex(points=[[0, 0], [1, 0], [0, 1], [1, 1]], max_edge_lengt simplex_tree = rips.create_simplex_tree(max_dimension=1) -print("filtrations=", simplex_tree.get_filtration()) +print("filtrations=") +for simplex_with_filtration in simplex_tree.get_filtration(): + print("(%s, %.2f)" % tuple(simplex_with_filtration)) + print("star([0])=", simplex_tree.get_star([0])) print("coface([0], 1)=", simplex_tree.get_cofaces([0], 1)) diff --git a/src/python/example/simplex_tree_example.py b/src/python/example/simplex_tree_example.py index 30de00da..34833899 100755 --- a/src/python/example/simplex_tree_example.py +++ b/src/python/example/simplex_tree_example.py @@ -38,8 +38,15 @@ else: print("dimension=", st.dimension()) +print("simplices=") +for simplex_with_filtration in st.get_simplices(): + print("(%s, %.2f)" % tuple(simplex_with_filtration)) + st.initialize_filtration() -print("filtration=", st.get_filtration()) +print("filtration=") +for simplex_with_filtration in st.get_filtration(): + print("(%s, %.2f)" % tuple(simplex_with_filtration)) + print("filtration[1, 2]=", st.filtration([1, 2])) print("filtration[4, 2]=", st.filtration([4, 2])) diff --git a/src/python/gudhi/point_cloud/__init__.py b/src/python/gudhi/point_cloud/__init__.py new file mode 100644 index 00000000..e69de29b --- /dev/null +++ b/src/python/gudhi/point_cloud/__init__.py diff --git a/src/python/gudhi/point_cloud/timedelay.py b/src/python/gudhi/point_cloud/timedelay.py new file mode 100644 index 00000000..f01df442 --- /dev/null +++ b/src/python/gudhi/point_cloud/timedelay.py @@ -0,0 +1,95 @@ +# This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT. +# See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. +# Author(s): Martin Royer, Yuichi Ike, Masatoshi Takenouchi +# +# Copyright (C) 2020 Inria, Copyright (C) 2020 Fujitsu Laboratories Ltd. +# Modification(s): +# - YYYY/MM Author: Description of the modification + +import numpy as np + + +class TimeDelayEmbedding: + """Point cloud transformation class. + Embeds time-series data in the R^d according to [Takens' Embedding Theorem] + (https://en.wikipedia.org/wiki/Takens%27s_theorem) and obtains the + coordinates of each point. + + Parameters + ---------- + dim : int, optional (default=3) + `d` of R^d to be embedded. + delay : int, optional (default=1) + Time-Delay embedding. + skip : int, optional (default=1) + How often to skip embedded points. + + Example + ------- + + Given delay=3 and skip=2, a point cloud which is obtained by embedding + a scalar time-series into R^3 is as follows:: + + time-series = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10] + point cloud = [[1, 4, 7], + [3, 6, 9]] + + Given delay=1 and skip=1, a point cloud which is obtained by embedding + a 2D vector time-series data into R^4 is as follows:: + + time-series = [[0, 1], [2, 3], [4, 5], [6, 7], [8, 9]] + point cloud = [[0, 1, 2, 3], + [2, 3, 4, 5], + [4, 5, 6, 7], + [6, 7, 8, 9]] + """ + + def __init__(self, dim=3, delay=1, skip=1): + self._dim = dim + self._delay = delay + self._skip = skip + + def __call__(self, ts): + """Transform method for single time-series data. + + Parameters + ---------- + ts : Iterable[float] or Iterable[Iterable[float]] + A single time-series data, with scalar or vector values. + + Returns + ------- + point cloud : n x dim numpy arrays + Makes point cloud from a single time-series data. + """ + return self._transform(np.array(ts)) + + def fit(self, ts, y=None): + return self + + def _transform(self, ts): + """Guts of transform method.""" + if ts.ndim == 1: + repeat = self._dim + else: + assert self._dim % ts.shape[1] == 0 + repeat = self._dim // ts.shape[1] + end = len(ts) - self._delay * (repeat - 1) + short = np.arange(0, end, self._skip) + vertical = np.arange(0, repeat * self._delay, self._delay) + return ts[np.add.outer(short, vertical)].reshape(len(short), -1) + + def transform(self, ts): + """Transform method for multiple time-series data. + + Parameters + ---------- + ts : Iterable[Iterable[float]] or Iterable[Iterable[Iterable[float]]] + Multiple time-series data, with scalar or vector values. + + Returns + ------- + point clouds : list of n x dim numpy arrays + Makes point cloud from each time-series data. + """ + return [self._transform(np.array(s)) for s in ts] diff --git a/src/python/gudhi/simplex_tree.pxd b/src/python/gudhi/simplex_tree.pxd index 96d14079..82f155de 100644 --- a/src/python/gudhi/simplex_tree.pxd +++ b/src/python/gudhi/simplex_tree.pxd @@ -21,6 +21,22 @@ cdef extern from "Simplex_tree_interface.h" namespace "Gudhi": cdef cppclass Simplex_tree_options_full_featured: pass + cdef cppclass Simplex_tree_simplex_handle "Gudhi::Simplex_tree_interface<Gudhi::Simplex_tree_options_full_featured>::Simplex_handle": + pass + + cdef cppclass Simplex_tree_simplices_iterator "Gudhi::Simplex_tree_interface<Gudhi::Simplex_tree_options_full_featured>::Complex_simplex_iterator": + Simplex_tree_simplices_iterator() + Simplex_tree_simplex_handle& operator*() + Simplex_tree_simplices_iterator operator++() + bint operator!=(Simplex_tree_simplices_iterator) + + cdef cppclass Simplex_tree_skeleton_iterator "Gudhi::Simplex_tree_interface<Gudhi::Simplex_tree_options_full_featured>::Skeleton_simplex_iterator": + Simplex_tree_skeleton_iterator() + Simplex_tree_simplex_handle& operator*() + Simplex_tree_skeleton_iterator operator++() + bint operator!=(Simplex_tree_skeleton_iterator) + + cdef cppclass Simplex_tree_interface_full_featured "Gudhi::Simplex_tree_interface<Gudhi::Simplex_tree_options_full_featured>": Simplex_tree() double simplex_filtration(vector[int] simplex) @@ -34,8 +50,6 @@ cdef extern from "Simplex_tree_interface.h" namespace "Gudhi": bool find_simplex(vector[int] simplex) bool insert_simplex_and_subfaces(vector[int] simplex, double filtration) - vector[pair[vector[int], double]] get_filtration() - vector[pair[vector[int], double]] get_skeleton(int dimension) vector[pair[vector[int], double]] get_star(vector[int] simplex) vector[pair[vector[int], double]] get_cofaces(vector[int] simplex, int dimension) @@ -43,6 +57,14 @@ cdef extern from "Simplex_tree_interface.h" namespace "Gudhi": void remove_maximal_simplex(vector[int] simplex) bool prune_above_filtration(double filtration) bool make_filtration_non_decreasing() + # Iterators over Simplex tree + pair[vector[int], double] get_simplex_and_filtration(Simplex_tree_simplex_handle f_simplex) + Simplex_tree_simplices_iterator get_simplices_iterator_begin() + Simplex_tree_simplices_iterator get_simplices_iterator_end() + vector[Simplex_tree_simplex_handle].const_iterator get_filtration_iterator_begin() + vector[Simplex_tree_simplex_handle].const_iterator get_filtration_iterator_end() + Simplex_tree_skeleton_iterator get_skeleton_iterator_begin(int dimension) + Simplex_tree_skeleton_iterator get_skeleton_iterator_end(int dimension) cdef extern from "Persistent_cohomology_interface.h" namespace "Gudhi": cdef cppclass Simplex_tree_persistence_interface "Gudhi::Persistent_cohomology_interface<Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_full_featured>>": diff --git a/src/python/gudhi/simplex_tree.pyx b/src/python/gudhi/simplex_tree.pyx index b18627c4..c01cc905 100644 --- a/src/python/gudhi/simplex_tree.pyx +++ b/src/python/gudhi/simplex_tree.pyx @@ -7,6 +7,7 @@ # Modification(s): # - YYYY/MM Author: Description of the modification +from cython.operator import dereference, preincrement from libc.stdint cimport intptr_t from numpy import array as np_array cimport simplex_tree @@ -207,22 +208,34 @@ cdef class SimplexTree: return self.get_ptr().insert_simplex_and_subfaces(csimplex, <double>filtration) - def get_filtration(self): - """This function returns a list of all simplices with their given + def get_simplices(self): + """This function returns a generator with simplices and their given filtration values. + :returns: The simplices. + :rtype: generator with tuples(simplex, filtration) + """ + cdef Simplex_tree_simplices_iterator it = self.get_ptr().get_simplices_iterator_begin() + cdef Simplex_tree_simplices_iterator end = self.get_ptr().get_simplices_iterator_end() + cdef Simplex_tree_simplex_handle sh = dereference(it) + + while it != end: + yield self.get_ptr().get_simplex_and_filtration(dereference(it)) + preincrement(it) + + def get_filtration(self): + """This function returns a generator with simplices and their given + filtration values sorted by increasing filtration values. + :returns: The simplices sorted by increasing filtration values. - :rtype: list of tuples(simplex, filtration) + :rtype: generator with tuples(simplex, filtration) """ - cdef vector[pair[vector[int], double]] filtration \ - = self.get_ptr().get_filtration() - ct = [] - for filtered_complex in filtration: - v = [] - for vertex in filtered_complex.first: - v.append(vertex) - ct.append((v, filtered_complex.second)) - return ct + cdef vector[Simplex_tree_simplex_handle].const_iterator it = self.get_ptr().get_filtration_iterator_begin() + cdef vector[Simplex_tree_simplex_handle].const_iterator end = self.get_ptr().get_filtration_iterator_end() + + while it != end: + yield self.get_ptr().get_simplex_and_filtration(dereference(it)) + preincrement(it) def get_skeleton(self, dimension): """This function returns the (simplices of the) skeleton of a maximum @@ -233,15 +246,12 @@ cdef class SimplexTree: :returns: The (simplices of the) skeleton of a maximum dimension. :rtype: list of tuples(simplex, filtration) """ - cdef vector[pair[vector[int], double]] skeleton \ - = self.get_ptr().get_skeleton(<int>dimension) - ct = [] - for filtered_simplex in skeleton: - v = [] - for vertex in filtered_simplex.first: - v.append(vertex) - ct.append((v, filtered_simplex.second)) - return ct + cdef Simplex_tree_skeleton_iterator it = self.get_ptr().get_skeleton_iterator_begin(dimension) + cdef Simplex_tree_skeleton_iterator end = self.get_ptr().get_skeleton_iterator_end(dimension) + + while it != end: + yield self.get_ptr().get_simplex_and_filtration(dereference(it)) + preincrement(it) def get_star(self, simplex): """This function returns the star of a given N-simplex. diff --git a/src/python/include/Simplex_tree_interface.h b/src/python/include/Simplex_tree_interface.h index 06f31341..4a7062d6 100644 --- a/src/python/include/Simplex_tree_interface.h +++ b/src/python/include/Simplex_tree_interface.h @@ -33,7 +33,10 @@ class Simplex_tree_interface : public Simplex_tree<SimplexTreeOptions> { using Simplex_handle = typename Base::Simplex_handle; using Insertion_result = typename std::pair<Simplex_handle, bool>; using Simplex = std::vector<Vertex_handle>; - using Filtered_simplices = std::vector<std::pair<Simplex, Filtration_value>>; + using Simplex_and_filtration = std::pair<Simplex, Filtration_value>; + using Filtered_simplices = std::vector<Simplex_and_filtration>; + using Skeleton_simplex_iterator = typename Base::Skeleton_simplex_iterator; + using Complex_simplex_iterator = typename Base::Complex_simplex_iterator; public: bool find_simplex(const Simplex& vh) { @@ -82,29 +85,12 @@ class Simplex_tree_interface : public Simplex_tree<SimplexTreeOptions> { Base::initialize_filtration(); } - Filtered_simplices get_filtration() { - Base::initialize_filtration(); - Filtered_simplices filtrations; - for (auto f_simplex : Base::filtration_simplex_range()) { - Simplex simplex; - for (auto vertex : Base::simplex_vertex_range(f_simplex)) { - simplex.insert(simplex.begin(), vertex); - } - filtrations.push_back(std::make_pair(simplex, Base::filtration(f_simplex))); - } - return filtrations; - } - - Filtered_simplices get_skeleton(int dimension) { - Filtered_simplices skeletons; - for (auto f_simplex : Base::skeleton_simplex_range(dimension)) { - Simplex simplex; - for (auto vertex : Base::simplex_vertex_range(f_simplex)) { - simplex.insert(simplex.begin(), vertex); - } - skeletons.push_back(std::make_pair(simplex, Base::filtration(f_simplex))); + Simplex_and_filtration get_simplex_and_filtration(Simplex_handle f_simplex) { + Simplex simplex; + for (auto vertex : Base::simplex_vertex_range(f_simplex)) { + simplex.insert(simplex.begin(), vertex); } - return skeletons; + return std::make_pair(std::move(simplex), Base::filtration(f_simplex)); } Filtered_simplices get_star(const Simplex& simplex) { @@ -135,6 +121,38 @@ class Simplex_tree_interface : public Simplex_tree<SimplexTreeOptions> { Base::initialize_filtration(); pcoh = new Gudhi::Persistent_cohomology_interface<Base>(*this); } + + // Iterator over the simplex tree + Complex_simplex_iterator get_simplices_iterator_begin() { + // this specific case works because the range is just a pair of iterators - won't work if range was a vector + return Base::complex_simplex_range().begin(); + } + + Complex_simplex_iterator get_simplices_iterator_end() { + // this specific case works because the range is just a pair of iterators - won't work if range was a vector + return Base::complex_simplex_range().end(); + } + + typename std::vector<Simplex_handle>::const_iterator get_filtration_iterator_begin() { + // Base::initialize_filtration(); already performed in filtration_simplex_range + // this specific case works because the range is just a pair of iterators - won't work if range was a vector + return Base::filtration_simplex_range().begin(); + } + + typename std::vector<Simplex_handle>::const_iterator get_filtration_iterator_end() { + // this specific case works because the range is just a pair of iterators - won't work if range was a vector + return Base::filtration_simplex_range().end(); + } + + Skeleton_simplex_iterator get_skeleton_iterator_begin(int dimension) { + // this specific case works because the range is just a pair of iterators - won't work if range was a vector + return Base::skeleton_simplex_range(dimension).begin(); + } + + Skeleton_simplex_iterator get_skeleton_iterator_end(int dimension) { + // this specific case works because the range is just a pair of iterators - won't work if range was a vector + return Base::skeleton_simplex_range(dimension).end(); + } }; } // namespace Gudhi diff --git a/src/python/test/test_alpha_complex.py b/src/python/test/test_alpha_complex.py index 3761fe16..77121302 100755 --- a/src/python/test/test_alpha_complex.py +++ b/src/python/test/test_alpha_complex.py @@ -40,7 +40,7 @@ def test_infinite_alpha(): assert simplex_tree.num_simplices() == 11 assert simplex_tree.num_vertices() == 4 - assert simplex_tree.get_filtration() == [ + assert list(simplex_tree.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([2], 0.0), @@ -53,6 +53,7 @@ def test_infinite_alpha(): ([0, 1, 2], 0.5), ([1, 2, 3], 0.5), ] + assert simplex_tree.get_star([0]) == [ ([0], 0.0), ([0, 1], 0.25), @@ -105,7 +106,7 @@ def test_filtered_alpha(): else: assert False - assert simplex_tree.get_filtration() == [ + assert list(simplex_tree.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([2], 0.0), diff --git a/src/python/test/test_euclidean_witness_complex.py b/src/python/test/test_euclidean_witness_complex.py index c18d2484..f3664d39 100755 --- a/src/python/test/test_euclidean_witness_complex.py +++ b/src/python/test/test_euclidean_witness_complex.py @@ -40,7 +40,7 @@ def test_witness_complex(): assert landmarks[1] == euclidean_witness_complex.get_point(1) assert landmarks[2] == euclidean_witness_complex.get_point(2) - assert simplex_tree.get_filtration() == [ + assert list(simplex_tree.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([0, 1], 0.0), @@ -78,13 +78,13 @@ def test_strong_witness_complex(): assert landmarks[1] == euclidean_strong_witness_complex.get_point(1) assert landmarks[2] == euclidean_strong_witness_complex.get_point(2) - assert simplex_tree.get_filtration() == [([0], 0.0), ([1], 0.0), ([2], 0.0)] + assert list(simplex_tree.get_filtration()) == [([0], 0.0), ([1], 0.0), ([2], 0.0)] simplex_tree = euclidean_strong_witness_complex.create_simplex_tree( max_alpha_square=100.0 ) - assert simplex_tree.get_filtration() == [ + assert list(simplex_tree.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([2], 0.0), diff --git a/src/python/test/test_rips_complex.py b/src/python/test/test_rips_complex.py index b02a68e1..b86e7498 100755 --- a/src/python/test/test_rips_complex.py +++ b/src/python/test/test_rips_complex.py @@ -32,7 +32,7 @@ def test_rips_from_points(): assert simplex_tree.num_simplices() == 10 assert simplex_tree.num_vertices() == 4 - assert simplex_tree.get_filtration() == [ + assert list(simplex_tree.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([2], 0.0), @@ -44,6 +44,7 @@ def test_rips_from_points(): ([1, 2], 1.4142135623730951), ([0, 3], 1.4142135623730951), ] + assert simplex_tree.get_star([0]) == [ ([0], 0.0), ([0, 1], 1.0), @@ -95,7 +96,7 @@ def test_rips_from_distance_matrix(): assert simplex_tree.num_simplices() == 10 assert simplex_tree.num_vertices() == 4 - assert simplex_tree.get_filtration() == [ + assert list(simplex_tree.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([2], 0.0), @@ -107,6 +108,7 @@ def test_rips_from_distance_matrix(): ([1, 2], 1.4142135623730951), ([0, 3], 1.4142135623730951), ] + assert simplex_tree.get_star([0]) == [ ([0], 0.0), ([0, 1], 1.0), diff --git a/src/python/test/test_simplex_tree.py b/src/python/test/test_simplex_tree.py index 1822c43b..f7848379 100755 --- a/src/python/test/test_simplex_tree.py +++ b/src/python/test/test_simplex_tree.py @@ -55,7 +55,7 @@ def test_insertion(): assert st.filtration([1]) == 0.0 # skeleton test - assert st.get_skeleton(2) == [ + assert list(st.get_skeleton(2)) == [ ([0, 1, 2], 4.0), ([0, 1], 0.0), ([0, 2], 4.0), @@ -64,7 +64,7 @@ def test_insertion(): ([1], 0.0), ([2], 4.0), ] - assert st.get_skeleton(1) == [ + assert list(st.get_skeleton(1)) == [ ([0, 1], 0.0), ([0, 2], 4.0), ([0], 0.0), @@ -72,12 +72,12 @@ def test_insertion(): ([1], 0.0), ([2], 4.0), ] - assert st.get_skeleton(0) == [([0], 0.0), ([1], 0.0), ([2], 4.0)] + assert list(st.get_skeleton(0)) == [([0], 0.0), ([1], 0.0), ([2], 4.0)] # remove_maximal_simplex test assert st.get_cofaces([0, 1, 2], 1) == [] st.remove_maximal_simplex([0, 1, 2]) - assert st.get_skeleton(2) == [ + assert list(st.get_skeleton(2)) == [ ([0, 1], 0.0), ([0, 2], 4.0), ([0], 0.0), @@ -126,7 +126,8 @@ def test_expansion(): assert st.num_vertices() == 7 assert st.num_simplices() == 17 - assert st.get_filtration() == [ + + assert list(st.get_filtration()) == [ ([2], 0.1), ([3], 0.1), ([2, 3], 0.1), @@ -151,7 +152,7 @@ def test_expansion(): assert st.num_simplices() == 22 st.initialize_filtration() - assert st.get_filtration() == [ + assert list(st.get_filtration()) == [ ([2], 0.1), ([3], 0.1), ([2, 3], 0.1), @@ -248,3 +249,15 @@ def test_make_filtration_non_decreasing(): assert st.filtration([3, 4, 5]) == 2.0 assert st.filtration([3, 4]) == 2.0 assert st.filtration([4, 5]) == 2.0 + +def test_simplices_iterator(): + st = SimplexTree() + + assert st.insert([0, 1, 2], filtration=4.0) == True + assert st.insert([2, 3, 4], filtration=2.0) == True + + for simplex in st.get_simplices(): + print("simplex is: ", simplex[0]) + assert st.find(simplex[0]) == True + print("filtration is: ", simplex[1]) + assert st.filtration(simplex[0]) == simplex[1] diff --git a/src/python/test/test_tangential_complex.py b/src/python/test/test_tangential_complex.py index e650e99c..8668a2e0 100755 --- a/src/python/test/test_tangential_complex.py +++ b/src/python/test/test_tangential_complex.py @@ -37,7 +37,7 @@ def test_tangential(): assert st.num_simplices() == 6 assert st.num_vertices() == 4 - assert st.get_filtration() == [ + assert list(st.get_filtration()) == [ ([0], 0.0), ([1], 0.0), ([2], 0.0), @@ -45,6 +45,7 @@ def test_tangential(): ([3], 0.0), ([1, 3], 0.0), ] + assert st.get_cofaces([0], 1) == [([0, 2], 0.0)] assert point_list[0] == tc.get_point(0) diff --git a/src/python/test/test_time_delay.py b/src/python/test/test_time_delay.py new file mode 100755 index 00000000..1ead9bca --- /dev/null +++ b/src/python/test/test_time_delay.py @@ -0,0 +1,43 @@ +from gudhi.point_cloud.timedelay import TimeDelayEmbedding +import numpy as np + + +def test_normal(): + # Sample array + ts = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10] + # Normal case. + prep = TimeDelayEmbedding() + pointclouds = prep(ts) + assert (pointclouds[0] == np.array([1, 2, 3])).all() + assert (pointclouds[1] == np.array([2, 3, 4])).all() + assert (pointclouds[2] == np.array([3, 4, 5])).all() + assert (pointclouds[3] == np.array([4, 5, 6])).all() + assert (pointclouds[4] == np.array([5, 6, 7])).all() + assert (pointclouds[5] == np.array([6, 7, 8])).all() + assert (pointclouds[6] == np.array([7, 8, 9])).all() + assert (pointclouds[7] == np.array([8, 9, 10])).all() + # Delay = 3 + prep = TimeDelayEmbedding(delay=3) + pointclouds = prep(ts) + assert (pointclouds[0] == np.array([1, 4, 7])).all() + assert (pointclouds[1] == np.array([2, 5, 8])).all() + assert (pointclouds[2] == np.array([3, 6, 9])).all() + assert (pointclouds[3] == np.array([4, 7, 10])).all() + # Skip = 3 + prep = TimeDelayEmbedding(skip=3) + pointclouds = prep(ts) + assert (pointclouds[0] == np.array([1, 2, 3])).all() + assert (pointclouds[1] == np.array([4, 5, 6])).all() + assert (pointclouds[2] == np.array([7, 8, 9])).all() + # Delay = 2 / Skip = 2 + prep = TimeDelayEmbedding(delay=2, skip=2) + pointclouds = prep(ts) + assert (pointclouds[0] == np.array([1, 3, 5])).all() + assert (pointclouds[1] == np.array([3, 5, 7])).all() + assert (pointclouds[2] == np.array([5, 7, 9])).all() + + # Vector series + ts = np.arange(0, 10).reshape(-1, 2) + prep = TimeDelayEmbedding(dim=4) + prep.fit([ts]) + assert (prep.transform([ts])[0] == [[0, 1, 2, 3], [2, 3, 4, 5], [4, 5, 6, 7], [6, 7, 8, 9]]).all() |