diff options
Diffstat (limited to 'src/python/gudhi')
| -rw-r--r-- | src/python/gudhi/__init__.py | 1 | ||||
| -rw-r--r-- | src/python/gudhi/__init__.py.in | 40 | ||||
| -rw-r--r-- | src/python/gudhi/alpha_complex.pyx | 116 | ||||
| -rw-r--r-- | src/python/gudhi/bottleneck.pyx | 49 | ||||
| -rw-r--r-- | src/python/gudhi/cubical_complex.pyx | 188 | ||||
| -rw-r--r-- | src/python/gudhi/euclidean_strong_witness_complex.pyx | 92 | ||||
| -rw-r--r-- | src/python/gudhi/euclidean_witness_complex.pyx | 92 | ||||
| -rw-r--r-- | src/python/gudhi/nerve_gic.pyx | 412 | ||||
| -rw-r--r-- | src/python/gudhi/off_reader.pyx | 38 | ||||
| -rw-r--r-- | src/python/gudhi/periodic_cubical_complex.pyx | 190 | ||||
| -rw-r--r-- | src/python/gudhi/persistence_graphical_tools.py | 423 | ||||
| -rw-r--r-- | src/python/gudhi/reader_utils.pyx | 87 | ||||
| -rw-r--r-- | src/python/gudhi/rips_complex.pyx | 103 | ||||
| -rw-r--r-- | src/python/gudhi/simplex_tree.pxd | 56 | ||||
| -rw-r--r-- | src/python/gudhi/simplex_tree.pyx | 508 | ||||
| -rw-r--r-- | src/python/gudhi/strong_witness_complex.pyx | 78 | ||||
| -rw-r--r-- | src/python/gudhi/subsampling.pyx | 130 | ||||
| -rw-r--r-- | src/python/gudhi/tangential_complex.pyx | 173 | ||||
| -rw-r--r-- | src/python/gudhi/witness_complex.pyx | 78 |
19 files changed, 2854 insertions, 0 deletions
diff --git a/src/python/gudhi/__init__.py b/src/python/gudhi/__init__.py new file mode 100644 index 00000000..fde749eb --- /dev/null +++ b/src/python/gudhi/__init__.py @@ -0,0 +1 @@ +# Fake empty __init__.py for cython to accept this directory as a Python package diff --git a/src/python/gudhi/__init__.py.in b/src/python/gudhi/__init__.py.in new file mode 100644 index 00000000..947aa3c9 --- /dev/null +++ b/src/python/gudhi/__init__.py.in @@ -0,0 +1,40 @@ +from importlib import import_module + +"""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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "GUDHI Editorial Board" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "https://gudhi.inria.fr/licensing/" +__version__ = "@GUDHI_VERSION@" +# This variable is used by doctest to find files +__root_source_dir__ = "@CMAKE_SOURCE_DIR@" +__debug_info__ = @GUDHI_PYTHON_DEBUG_INFO@ + +from sys import exc_info +from importlib import import_module + +__all__ = [@GUDHI_PYTHON_MODULES@] + +__available_modules__ = '' +__missing_modules__ = '' + +# try to import * from gudhi.__module_name__ +for __module_name__ in __all__: + try: + __module__ = import_module('gudhi.' + __module_name__) + try: + __to_import__ = __module__.__all__ + except AttributeError: + __to_import__ = [name for name in __module__.__dict__ if not name.startswith('_')] + globals().update({name: __module__.__dict__[name] for name in __to_import__}) + __available_modules__ += __module_name__ + ";" + except: + __missing_modules__ += __module_name__ + ";" diff --git a/src/python/gudhi/alpha_complex.pyx b/src/python/gudhi/alpha_complex.pyx new file mode 100644 index 00000000..6d6309db --- /dev/null +++ b/src/python/gudhi/alpha_complex.pyx @@ -0,0 +1,116 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp.string cimport string +from libcpp cimport bool +from libc.stdint cimport intptr_t +import os + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "GPL v3" + +cdef extern from "Alpha_complex_interface.h" namespace "Gudhi": + cdef cppclass Alpha_complex_interface "Gudhi::alpha_complex::Alpha_complex_interface": + Alpha_complex_interface(vector[vector[double]] points) + # bool from_file is a workaround for cython to find the correct signature + Alpha_complex_interface(string off_file, bool from_file) + vector[double] get_point(int vertex) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square) + +# AlphaComplex python interface +cdef class AlphaComplex: + """AlphaComplex is a simplicial complex constructed from the finite cells + of a Delaunay Triangulation. + + The filtration value of each simplex is computed as the square of the + circumradius of the simplex if the circumsphere is empty (the simplex is + then said to be Gabriel), and as the minimum of the filtration values of + the codimension 1 cofaces that make it not Gabriel otherwise. + + All simplices that have a filtration value strictly greater than a given + alpha squared value are not inserted into the complex. + + .. note:: + + When Alpha_complex is constructed with an infinite value of alpha, the + complex is a Delaunay complex. + + """ + + cdef Alpha_complex_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, points=None, off_file=''): + """AlphaComplex constructor. + + :param points: A list of points in d-Dimension. + :type points: list of list of double + + Or + + :param off_file: An OFF file style name. + :type off_file: string + """ + + # The real cython constructor + def __cinit__(self, points=None, off_file=''): + if off_file is not '': + if os.path.isfile(off_file): + self.thisptr = new Alpha_complex_interface(str.encode(off_file), True) + else: + print("file " + off_file + " not found.") + else: + if points is None: + # Empty Alpha construction + points=[] + self.thisptr = new Alpha_complex_interface(points) + + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if AlphaComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def get_point(self, vertex): + """This function returns the point corresponding to a given vertex. + + :param vertex: The vertex. + :type vertex: int + :rtype: list of float + :returns: the point. + """ + cdef vector[double] point = self.thisptr.get_point(vertex) + return point + + def create_simplex_tree(self, max_alpha_square=float('inf')): + """ + :param max_alpha_square: The maximum alpha square threshold the + simplices shall not exceed. Default is set to infinity, and + there is very little point using anything else since it does + not save time. + :type max_alpha_square: float + :returns: A simplex tree created from the Delaunay Triangulation. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, max_alpha_square) + return stree diff --git a/src/python/gudhi/bottleneck.pyx b/src/python/gudhi/bottleneck.pyx new file mode 100644 index 00000000..4b378cbc --- /dev/null +++ b/src/python/gudhi/bottleneck.pyx @@ -0,0 +1,49 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +import os + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "GPL v3" + +cdef extern from "Bottleneck_distance_interface.h" namespace "Gudhi::persistence_diagram": + double bottleneck(vector[pair[double, double]], vector[pair[double, double]], double) + double bottleneck(vector[pair[double, double]], vector[pair[double, double]]) + +def bottleneck_distance(diagram_1, diagram_2, e=None): + """This function returns the point corresponding to a given vertex. + + :param diagram_1: The first diagram. + :type diagram_1: vector[pair[double, double]] + :param diagram_2: The second diagram. + :type diagram_2: vector[pair[double, double]] + :param e: If `e` is 0, this uses an expensive algorithm to compute the + exact distance. + If `e` is not 0, it asks for an additive `e`-approximation, and + currently also allows a small multiplicative error (the last 2 or 3 + bits of the mantissa may be wrong). This version of the algorithm takes + advantage of the limited precision of `double` and is usually a lot + faster to compute, whatever the value of `e`. + + Thus, by default, `e` is the smallest positive double. + :type e: float + :rtype: float + :returns: the bottleneck distance. + """ + if e is None: + # Default value is the smallest double value (not 0, 0 is for exact version) + return bottleneck(diagram_1, diagram_2) + else: + # Can be 0 for exact version + return bottleneck(diagram_1, diagram_2, e) diff --git a/src/python/gudhi/cubical_complex.pyx b/src/python/gudhi/cubical_complex.pyx new file mode 100644 index 00000000..0dc133d1 --- /dev/null +++ b/src/python/gudhi/cubical_complex.pyx @@ -0,0 +1,188 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp.string cimport string +from libcpp cimport bool +import os + +from numpy import array as np_array + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Cubical_complex_interface.h" namespace "Gudhi": + cdef cppclass Bitmap_cubical_complex_base_interface "Gudhi::Cubical_complex::Cubical_complex_interface<>": + Bitmap_cubical_complex_base_interface(vector[unsigned] dimensions, vector[double] top_dimensional_cells) + Bitmap_cubical_complex_base_interface(string perseus_file) + int num_simplices() + int dimension() + +cdef extern from "Persistent_cohomology_interface.h" namespace "Gudhi": + cdef cppclass Cubical_complex_persistence_interface "Gudhi::Persistent_cohomology_interface<Gudhi::Cubical_complex::Cubical_complex_interface<>>": + Cubical_complex_persistence_interface(Bitmap_cubical_complex_base_interface * st, bool persistence_dim_max) + vector[pair[int, pair[double, double]]] get_persistence(int homology_coeff_field, double min_persistence) + vector[int] betti_numbers() + vector[int] persistent_betti_numbers(double from_value, double to_value) + vector[pair[double,double]] intervals_in_dimension(int dimension) + +# CubicalComplex python interface +cdef class CubicalComplex: + """The CubicalComplex is an example of a structured complex useful in + computational mathematics (specially rigorous numerics) and image + analysis. + """ + cdef Bitmap_cubical_complex_base_interface * thisptr + + cdef Cubical_complex_persistence_interface * pcohptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, dimensions=None, top_dimensional_cells=None, + perseus_file=''): + """CubicalComplex constructor from dimensions and + top_dimensional_cells or from a Perseus-style file name. + + :param dimensions: A list of number of top dimensional cells. + :type dimensions: list of int + :param top_dimensional_cells: A list of cells filtration values. + :type top_dimensional_cells: list of double + + Or + + :param perseus_file: A Perseus-style file name. + :type perseus_file: string + """ + + # The real cython constructor + def __cinit__(self, dimensions=None, top_dimensional_cells=None, + perseus_file=''): + if (dimensions is not None) and (top_dimensional_cells is not None) and (perseus_file is ''): + self.thisptr = new Bitmap_cubical_complex_base_interface(dimensions, top_dimensional_cells) + elif (dimensions is None) and (top_dimensional_cells is None) and (perseus_file is not ''): + if os.path.isfile(perseus_file): + self.thisptr = new Bitmap_cubical_complex_base_interface(str.encode(perseus_file)) + else: + print("file " + perseus_file + " not found.") + else: + print("CubicalComplex can be constructed from dimensions and " + "top_dimensional_cells or from a Perseus-style file name.") + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + if self.pcohptr != NULL: + del self.pcohptr + + def __is_defined(self): + """Returns true if CubicalComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def __is_persistence_defined(self): + """Returns true if Persistence pointer is not NULL. + """ + return self.pcohptr != NULL + + def num_simplices(self): + """This function returns the number of all cubes in the complex. + + :returns: int -- the number of all cubes in the complex. + """ + return self.thisptr.num_simplices() + + def dimension(self): + """This function returns the dimension of the complex. + + :returns: int -- the complex dimension. + """ + return self.thisptr.dimension() + + def persistence(self, homology_coeff_field=11, min_persistence=0): + """This function returns the persistence of the complex. + + :param homology_coeff_field: The homology coefficient field. Must be a + prime number + :type homology_coeff_field: int. + :param min_persistence: The minimum persistence value to take into + account (strictly greater than min_persistence). Default value is + 0.0. + Sets min_persistence to -1.0 to see all values. + :type min_persistence: float. + :returns: list of pairs(dimension, pair(birth, death)) -- the + persistence of the complex. + """ + if self.pcohptr != NULL: + del self.pcohptr + if self.thisptr != NULL: + self.pcohptr = new Cubical_complex_persistence_interface(self.thisptr, True) + cdef vector[pair[int, pair[double, double]]] persistence_result + if self.pcohptr != NULL: + persistence_result = self.pcohptr.get_persistence(homology_coeff_field, min_persistence) + return persistence_result + + def betti_numbers(self): + """This function returns the Betti numbers of the complex. + + :returns: list of int -- The Betti numbers ([B0, B1, ..., Bn]). + + :note: betti_numbers function requires persistence function to be + launched first. + + :note: betti_numbers function always returns [1, 0, 0, ...] as infinity + filtration cubes are not removed from the complex. + """ + cdef vector[int] bn_result + if self.pcohptr != NULL: + bn_result = self.pcohptr.betti_numbers() + return bn_result + + def persistent_betti_numbers(self, from_value, to_value): + """This function returns the persistent Betti numbers of the complex. + + :param from_value: The persistence birth limit to be added in the + numbers (persistent birth <= from_value). + :type from_value: float. + :param to_value: The persistence death limit to be added in the + numbers (persistent death > to_value). + :type to_value: float. + + :returns: list of int -- The persistent Betti numbers ([B0, B1, ..., + Bn]). + + :note: persistent_betti_numbers function requires persistence + function to be launched first. + """ + cdef vector[int] pbn_result + if self.pcohptr != NULL: + pbn_result = self.pcohptr.persistent_betti_numbers(<double>from_value, <double>to_value) + return pbn_result + + def persistence_intervals_in_dimension(self, dimension): + """This function returns the persistence intervals of the complex in a + specific dimension. + + :param dimension: The specific dimension. + :type dimension: int. + :returns: The persistence intervals. + :rtype: numpy array of dimension 2 + + :note: intervals_in_dim function requires persistence function to be + launched first. + """ + cdef vector[pair[double,double]] intervals_result + if self.pcohptr != NULL: + intervals_result = self.pcohptr.intervals_in_dimension(dimension) + else: + print("intervals_in_dim function requires persistence function" + " to be launched first.") + return np_array(intervals_result) diff --git a/src/python/gudhi/euclidean_strong_witness_complex.pyx b/src/python/gudhi/euclidean_strong_witness_complex.pyx new file mode 100644 index 00000000..5d6e4fb9 --- /dev/null +++ b/src/python/gudhi/euclidean_strong_witness_complex.pyx @@ -0,0 +1,92 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libc.stdint cimport intptr_t + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "GPL v3" + +cdef extern from "Euclidean_strong_witness_complex_interface.h" namespace "Gudhi": + cdef cppclass Euclidean_strong_witness_complex_interface "Gudhi::witness_complex::Euclidean_strong_witness_complex_interface": + Euclidean_strong_witness_complex_interface(vector[vector[double]] landmarks, vector[vector[double]] witnesses) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square, + unsigned limit_dimension) + vector[double] get_point(unsigned vertex) + +# EuclideanStrongWitnessComplex python interface +cdef class EuclideanStrongWitnessComplex: + """Constructs strong witness complex for given sets of witnesses and + landmarks in Euclidean space. + """ + + cdef Euclidean_strong_witness_complex_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, landmarks=None, witnesses=None): + """WitnessComplex constructor. + + :param landmarks: A list of landmarks (in the point cloud). + :type landmarks: list of list of double + + :param witnesses: The point cloud. + :type witnesses: list of list of double + """ + + # The real cython constructor + def __cinit__(self, landmarks=None, witnesses=None): + if landmarks is not None and witnesses is not None: + self.thisptr = new Euclidean_strong_witness_complex_interface(landmarks, witnesses) + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if WitnessComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def create_simplex_tree(self, max_alpha_square, limit_dimension = -1): + """ + :param max_alpha_square: The maximum alpha square threshold the + simplices shall not exceed. Default is set to infinity. + :type max_alpha_square: float + :returns: A simplex tree created from the Delaunay Triangulation. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + if limit_dimension is not -1: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square, limit_dimension) + else: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square) + return stree + + def get_point(self, vertex): + """This function returns the point corresponding to a given vertex. + + :param vertex: The vertex. + :type vertex: int. + :returns: The point. + :rtype: list of float + """ + cdef vector[double] point = self.thisptr.get_point(vertex) + return point + diff --git a/src/python/gudhi/euclidean_witness_complex.pyx b/src/python/gudhi/euclidean_witness_complex.pyx new file mode 100644 index 00000000..2531919b --- /dev/null +++ b/src/python/gudhi/euclidean_witness_complex.pyx @@ -0,0 +1,92 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libc.stdint cimport intptr_t + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "GPL v3" + +cdef extern from "Euclidean_witness_complex_interface.h" namespace "Gudhi": + cdef cppclass Euclidean_witness_complex_interface "Gudhi::witness_complex::Euclidean_witness_complex_interface": + Euclidean_witness_complex_interface(vector[vector[double]] landmarks, vector[vector[double]] witnesses) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square, + unsigned limit_dimension) + vector[double] get_point(unsigned vertex) + +# EuclideanWitnessComplex python interface +cdef class EuclideanWitnessComplex: + """Constructs (weak) witness complex for given sets of witnesses and + landmarks in Euclidean space. + """ + + cdef Euclidean_witness_complex_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, landmarks=None, witnesses=None): + """WitnessComplex constructor. + + :param landmarks: A list of landmarks (in the point cloud). + :type landmarks: list of list of double + + :param witnesses: The point cloud. + :type witnesses: list of list of double + """ + + # The real cython constructor + def __cinit__(self, landmarks=None, witnesses=None): + if landmarks is not None and witnesses is not None: + self.thisptr = new Euclidean_witness_complex_interface(landmarks, witnesses) + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if WitnessComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def create_simplex_tree(self, max_alpha_square, limit_dimension = -1): + """ + :param max_alpha_square: The maximum alpha square threshold the + simplices shall not exceed. Default is set to infinity. + :type max_alpha_square: float + :returns: A simplex tree created from the Delaunay Triangulation. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + if limit_dimension is not -1: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square, limit_dimension) + else: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square) + return stree + + def get_point(self, vertex): + """This function returns the point corresponding to a given vertex. + + :param vertex: The vertex. + :type vertex: int. + :returns: The point. + :rtype: list of float + """ + cdef vector[double] point = self.thisptr.get_point(vertex) + return point + diff --git a/src/python/gudhi/nerve_gic.pyx b/src/python/gudhi/nerve_gic.pyx new file mode 100644 index 00000000..2b230b8c --- /dev/null +++ b/src/python/gudhi/nerve_gic.pyx @@ -0,0 +1,412 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp.string cimport string +from libcpp cimport bool +import os +from libc.stdint cimport intptr_t + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2018 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2018 Inria" +__license__ = "GPL v3" + +cdef extern from "Nerve_gic_interface.h" namespace "Gudhi": + cdef cppclass Nerve_gic_interface "Gudhi::cover_complex::Nerve_gic_interface": + Nerve_gic_interface() + double compute_confidence_level_from_distance(double distance) + double compute_distance_from_confidence_level(double alpha) + void compute_distribution(int N) + double compute_p_value() + vector[pair[double, double]] compute_PD() + void find_simplices() + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree) + bool read_point_cloud(string off_file_name) + double set_automatic_resolution() + void set_color_from_coordinate(int k) + void set_color_from_file(string color_file_name) + void set_color_from_range(vector[double] color) + void set_cover_from_file(string cover_file_name) + void set_cover_from_function() + void set_cover_from_Euclidean_Voronoi(int m) + void set_function_from_coordinate(int k) + void set_function_from_file(string func_file_name) + void set_function_from_range(vector[double] function) + void set_gain(double g) + double set_graph_from_automatic_euclidean_rips(int N) + void set_graph_from_file(string graph_file_name) + void set_graph_from_OFF() + void set_graph_from_euclidean_rips(double threshold) + void set_mask(int nodemask) + void set_resolution_with_interval_length(double resolution) + void set_resolution_with_interval_number(int resolution) + void set_subsampling(double constant, double power) + void set_type(string type) + void set_verbose(bool verbose) + vector[int] subpopulation(int c) + void write_info() + void plot_DOT() + void plot_OFF() + void set_point_cloud_from_range(vector[vector[double]] cloud) + void set_distances_from_range(vector[vector[double]] distance_matrix) + +# CoverComplex python interface +cdef class CoverComplex: + """Cover complex data structure. + + The data structure is a simplicial complex, representing a Graph Induced + simplicial Complex (GIC) or a Nerve, and whose simplices are computed with + a cover C of a point cloud P, which often comes from the preimages of + intervals covering the image of a function f defined on P. These intervals + are parameterized by their resolution (either their length or their number) + and their gain (percentage of overlap). To compute a GIC, one also needs a + graph G built on top of P, whose cliques with vertices belonging to + different elements of C correspond to the simplices of the GIC. + """ + + cdef Nerve_gic_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self): + """CoverComplex constructor. + """ + + # The real cython constructor + def __cinit__(self): + self.thisptr = new Nerve_gic_interface() + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if CoverComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def set_point_cloud_from_range(self, cloud): + """ Reads and stores the input point cloud from a vector stored in memory. + + :param cloud: Input vector containing the point cloud. + :type cloud: vector[vector[double]] + """ + return self.thisptr.set_point_cloud_from_range(cloud) + + def set_distances_from_range(self, distance_matrix): + """ Reads and stores the input distance matrix from a vector stored in memory. + + :param distance_matrix: Input vector containing the distance matrix. + :type distance_matrix: vector[vector[double]] + """ + return self.thisptr.set_distances_from_range(distance_matrix) + + def compute_confidence_level_from_distance(self, distance): + """Computes the confidence level of a specific bottleneck distance + threshold. + + :param distance: Bottleneck distance. + :type distance: double + :rtype: double + :returns: Confidence level. + """ + return self.thisptr.compute_confidence_level_from_distance(distance) + + def compute_distance_from_confidence_level(self, alpha): + """Computes the bottleneck distance threshold corresponding to a + specific confidence level. + + :param alpha: Confidence level. + :type alpha: double + :rtype: double + :returns: Bottleneck distance. + """ + return self.thisptr.compute_distance_from_confidence_level(alpha) + + def compute_distribution(self, N=100): + """Computes bootstrapped distances distribution. + + :param N: Loop number (default value is 100). + :type alpha: int + """ + self.thisptr.compute_distribution(N) + + def compute_p_value(self): + """Computes the p-value, i.e. the opposite of the confidence level of + the largest bottleneck distance preserving the points in the + persistence diagram of the output simplicial complex. + + :rtype: double + :returns: p-value. + """ + return self.thisptr.compute_p_value() + + def compute_PD(self): + """Computes the extended persistence diagram of the complex. + """ + return self.thisptr.compute_PD() + + def create_simplex_tree(self): + """ + :returns: A simplex tree created from the Cover complex. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr) + return stree + + def find_simplices(self): + """Computes the simplices of the simplicial complex. + """ + self.thisptr.find_simplices() + + def read_point_cloud(self, off_file): + """Reads and stores the input point cloud from .(n)OFF file. + + :param off_file: Name of the input .OFF or .nOFF file. + :type off_file: string + :rtype: bool + :returns: Read file status. + """ + if os.path.isfile(off_file): + return self.thisptr.read_point_cloud(str.encode(off_file)) + else: + print("file " + off_file + " not found.") + return False + + def set_automatic_resolution(self): + """Computes the optimal length of intervals (i.e. the smallest interval + length avoiding discretization artifacts—see :cite:`Carriere17c`) for a + functional cover. + + :rtype: double + :returns: reso interval length used to compute the cover. + """ + return self.thisptr.set_automatic_resolution() + + def set_color_from_coordinate(self, k=0): + """Computes the function used to color the nodes of the simplicial + complex from the k-th coordinate. + + :param k: Coordinate to use (start at 0). Default value is 0. + :type k: int + """ + return self.thisptr.set_color_from_coordinate(k) + + def set_color_from_file(self, color_file_name): + """Computes the function used to color the nodes of the simplicial + complex from a file containing the function values. + + :param color_file_name: Name of the input color file. + :type color_file_name: string + """ + if os.path.isfile(color_file_name): + self.thisptr.set_color_from_file(str.encode(color_file_name)) + else: + print("file " + color_file_name + " not found.") + + def set_color_from_range(self, color): + """Computes the function used to color the nodes of the simplicial + complex from a vector stored in memory. + + :param color: Input vector of values. + :type color: vector[double] + """ + self.thisptr.set_color_from_range(color) + + def set_cover_from_file(self, cover_file_name): + """Creates the cover C from a file containing the cover elements of + each point (the order has to be the same as in the input file!). + + :param cover_file_name: Name of the input cover file. + :type cover_file_name: string + """ + if os.path.isfile(cover_file_name): + self.thisptr.set_cover_from_file(str.encode(cover_file_name)) + else: + print("file " + cover_file_name + " not found.") + + def set_cover_from_function(self): + """Creates a cover C from the preimages of the function f. + """ + self.thisptr.set_cover_from_function() + + def set_cover_from_Voronoi(self, m=100): + """Creates the cover C from the Voronoï cells of a subsampling of the + point cloud. + + :param m: Number of points in the subsample. Default value is 100. + :type m: int + """ + self.thisptr.set_cover_from_Euclidean_Voronoi(m) + + def set_function_from_coordinate(self, k): + """Creates the function f from the k-th coordinate of the point cloud. + + :param k: Coordinate to use (start at 0). + :type k: int + """ + self.thisptr.set_function_from_coordinate(k) + + def set_function_from_file(self, func_file_name): + """Creates the function f from a file containing the function values. + + :param func_file_name: Name of the input function file. + :type func_file_name: string + """ + if os.path.isfile(func_file_name): + self.thisptr.set_function_from_file(str.encode(func_file_name)) + else: + print("file " + func_file_name + " not found.") + + def set_function_from_range(self, function): + """Creates the function f from a vector stored in memory. + + :param function: Input vector of values. + :type function: vector[double] + """ + self.thisptr.set_function_from_range(function) + + def set_gain(self, g = 0.3): + """Sets a gain from a value stored in memory. + + :param g: Gain (default value is 0.3). + :type g: double + """ + self.thisptr.set_gain(g) + + def set_graph_from_automatic_rips(self, N=100): + """Creates a graph G from a Rips complex whose threshold value is + automatically tuned with subsampling—see. + + :param N: Number of subsampling iteration (the default reasonable value + is 100, but there is no guarantee on how to choose it). + :type N: int + :rtype: double + :returns: Delta threshold used for computing the Rips complex. + """ + return self.thisptr.set_graph_from_automatic_euclidean_rips(N) + + def set_graph_from_file(self, graph_file_name): + """Creates a graph G from a file containing the edges. + + :param graph_file_name: Name of the input graph file. The graph file + contains one edge per line, each edge being represented by the IDs of + its two nodes. + :type graph_file_name: string + """ + if os.path.isfile(graph_file_name): + self.thisptr.set_graph_from_file(str.encode(graph_file_name)) + else: + print("file " + graph_file_name + " not found.") + + def set_graph_from_OFF(self): + """Creates a graph G from the triangulation given by the input OFF + file. + """ + self.thisptr.set_graph_from_OFF() + + def set_graph_from_rips(self, threshold): + """Creates a graph G from a Rips complex. + + :param threshold: Threshold value for the Rips complex. + :type threshold: double + """ + self.thisptr.set_graph_from_euclidean_rips(threshold) + + def set_mask(self, nodemask): + """Sets the mask, which is a threshold integer such that nodes in the + complex that contain a number of data points which is less than or + equal to this threshold are not displayed. + + :param nodemask: Threshold. + :type nodemask: int + """ + self.thisptr.set_mask(nodemask) + + def set_resolution_with_interval_length(self, resolution): + """Sets a length of intervals from a value stored in memory. + + :param resolution: Length of intervals. + :type resolution: double + """ + self.thisptr.set_resolution_with_interval_length(resolution) + + def set_resolution_with_interval_number(self, resolution): + """Sets a number of intervals from a value stored in memory. + + :param resolution: Number of intervals. + :type resolution: int + """ + self.thisptr.set_resolution_with_interval_number(resolution) + + def set_subsampling(self, constant, power): + """Sets the constants used to subsample the data set. These constants + are explained in :cite:`Carriere17c`. + + :param constant: Constant. + :type constant: double + :param power: Power. + :type resolution: double + """ + self.thisptr.set_subsampling(constant, power) + + def set_type(self, type): + """Specifies whether the type of the output simplicial complex. + + :param type: either "GIC" or "Nerve". + :type type: string + """ + self.thisptr.set_type(str.encode(type)) + + def set_verbose(self, verbose): + """Specifies whether the program should display information or not. + + :param verbose: true = display info, false = do not display info. + :type verbose: boolean + """ + self.thisptr.set_verbose(verbose) + + def subpopulation(self, c): + """Returns the data subset corresponding to a specific node of the + created complex. + + :param c: ID of the node. + :type c: int + :rtype: vector[int] + :returns: Vector of IDs of data points. + """ + return self.thisptr.subpopulation(c) + + def write_info(self): + """Creates a .txt file called SC.txt describing the 1-skeleton, which can + then be plotted with e.g. KeplerMapper. + """ + return self.thisptr.write_info() + + def plot_dot(self): + """Creates a .dot file called SC.dot for neato (part of the graphviz + package) once the simplicial complex is computed to get a visualization of + its 1-skeleton in a .pdf file. + """ + return self.thisptr.plot_DOT() + + def plot_off(self): + """Creates a .off file called SC.off for 3D visualization, which contains + the 2-skeleton of the GIC. This function assumes that the cover has been + computed with Voronoi. If data points are in 1D or 2D, the remaining + coordinates of the points embedded in 3D are set to 0. + """ + return self.thisptr.plot_OFF() diff --git a/src/python/gudhi/off_reader.pyx b/src/python/gudhi/off_reader.pyx new file mode 100644 index 00000000..9efd97ff --- /dev/null +++ b/src/python/gudhi/off_reader.pyx @@ -0,0 +1,38 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.string cimport string +import os + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Off_reader_interface.h" namespace "Gudhi": + vector[vector[double]] read_points_from_OFF_file(string off_file) + +def read_off(off_file=''): + """Read points from OFF file. + + :param off_file: An OFF file style name. + :type off_file: string + + :returns: The point set. + :rtype: vector[vector[double]] + """ + if off_file is not '': + if os.path.isfile(off_file): + return read_points_from_OFF_file(str.encode(off_file)) + else: + print("file " + off_file + " not found.") + return [] + diff --git a/src/python/gudhi/periodic_cubical_complex.pyx b/src/python/gudhi/periodic_cubical_complex.pyx new file mode 100644 index 00000000..724fadd4 --- /dev/null +++ b/src/python/gudhi/periodic_cubical_complex.pyx @@ -0,0 +1,190 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp.string cimport string +from libcpp cimport bool +import os + +from numpy import array as np_array + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Cubical_complex_interface.h" namespace "Gudhi": + cdef cppclass Periodic_cubical_complex_base_interface "Gudhi::Cubical_complex::Cubical_complex_interface<Gudhi::cubical_complex::Bitmap_cubical_complex_periodic_boundary_conditions_base<double>>": + Periodic_cubical_complex_base_interface(vector[unsigned] dimensions, vector[double] top_dimensional_cells, vector[bool] periodic_dimensions) + Periodic_cubical_complex_base_interface(string perseus_file) + int num_simplices() + int dimension() + +cdef extern from "Persistent_cohomology_interface.h" namespace "Gudhi": + cdef cppclass Periodic_cubical_complex_persistence_interface "Gudhi::Persistent_cohomology_interface<Gudhi::Cubical_complex::Cubical_complex_interface<Gudhi::cubical_complex::Bitmap_cubical_complex_periodic_boundary_conditions_base<double>>>": + Periodic_cubical_complex_persistence_interface(Periodic_cubical_complex_base_interface * st, bool persistence_dim_max) + vector[pair[int, pair[double, double]]] get_persistence(int homology_coeff_field, double min_persistence) + vector[int] betti_numbers() + vector[int] persistent_betti_numbers(double from_value, double to_value) + vector[pair[double,double]] intervals_in_dimension(int dimension) + +# PeriodicCubicalComplex python interface +cdef class PeriodicCubicalComplex: + """The PeriodicCubicalComplex is an example of a structured complex useful + in computational mathematics (specially rigorous numerics) and image + analysis. + """ + cdef Periodic_cubical_complex_base_interface * thisptr + + cdef Periodic_cubical_complex_persistence_interface * pcohptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, dimensions=None, top_dimensional_cells=None, + periodic_dimensions=None, perseus_file=''): + """PeriodicCubicalComplex constructor from dimensions and + top_dimensional_cells or from a Perseus-style file name. + + :param dimensions: A list of number of top dimensional cells. + :type dimensions: list of int + :param top_dimensional_cells: A list of cells filtration values. + :type top_dimensional_cells: list of double + :param periodic_dimensions: A list of top dimensional cells periodicity value. + :type periodic_dimensions: list of boolean + + Or + + :param perseus_file: A Perseus-style file name. + :type perseus_file: string + """ + + # The real cython constructor + def __cinit__(self, dimensions=None, top_dimensional_cells=None, + periodic_dimensions=None, perseus_file=''): + if (dimensions is not None) and (top_dimensional_cells is not None) and (periodic_dimensions is not None) and (perseus_file is ''): + self.thisptr = new Periodic_cubical_complex_base_interface(dimensions, top_dimensional_cells, periodic_dimensions) + elif (dimensions is None) and (top_dimensional_cells is None) and (periodic_dimensions is None) and (perseus_file is not ''): + if os.path.isfile(perseus_file): + self.thisptr = new Periodic_cubical_complex_base_interface(str.encode(perseus_file)) + else: + print("file " + perseus_file + " not found.") + else: + print("CubicalComplex can be constructed from dimensions and " + "top_dimensional_cells or from a Perseus-style file name.") + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + if self.pcohptr != NULL: + del self.pcohptr + + def __is_defined(self): + """Returns true if PeriodicCubicalComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def __is_persistence_defined(self): + """Returns true if Persistence pointer is not NULL. + """ + return self.pcohptr != NULL + + def num_simplices(self): + """This function returns the number of all cubes in the complex. + + :returns: int -- the number of all cubes in the complex. + """ + return self.thisptr.num_simplices() + + def dimension(self): + """This function returns the dimension of the complex. + + :returns: int -- the complex dimension. + """ + return self.thisptr.dimension() + + def persistence(self, homology_coeff_field=11, min_persistence=0): + """This function returns the persistence of the complex. + + :param homology_coeff_field: The homology coefficient field. Must be a + prime number + :type homology_coeff_field: int. + :param min_persistence: The minimum persistence value to take into + account (strictly greater than min_persistence). Default value is + 0.0. + Sets min_persistence to -1.0 to see all values. + :type min_persistence: float. + :returns: list of pairs(dimension, pair(birth, death)) -- the + persistence of the complex. + """ + if self.pcohptr != NULL: + del self.pcohptr + if self.thisptr != NULL: + self.pcohptr = new Periodic_cubical_complex_persistence_interface(self.thisptr, True) + cdef vector[pair[int, pair[double, double]]] persistence_result + if self.pcohptr != NULL: + persistence_result = self.pcohptr.get_persistence(homology_coeff_field, min_persistence) + return persistence_result + + def betti_numbers(self): + """This function returns the Betti numbers of the complex. + + :returns: list of int -- The Betti numbers ([B0, B1, ..., Bn]). + + :note: betti_numbers function requires persistence function to be + launched first. + + :note: betti_numbers function always returns [1, 0, 0, ...] as infinity + filtration cubes are not removed from the complex. + """ + cdef vector[int] bn_result + if self.pcohptr != NULL: + bn_result = self.pcohptr.betti_numbers() + return bn_result + + def persistent_betti_numbers(self, from_value, to_value): + """This function returns the persistent Betti numbers of the complex. + + :param from_value: The persistence birth limit to be added in the + numbers (persistent birth <= from_value). + :type from_value: float. + :param to_value: The persistence death limit to be added in the + numbers (persistent death > to_value). + :type to_value: float. + + :returns: list of int -- The persistent Betti numbers ([B0, B1, ..., + Bn]). + + :note: persistent_betti_numbers function requires persistence + function to be launched first. + """ + cdef vector[int] pbn_result + if self.pcohptr != NULL: + pbn_result = self.pcohptr.persistent_betti_numbers(<double>from_value, <double>to_value) + return pbn_result + + def persistence_intervals_in_dimension(self, dimension): + """This function returns the persistence intervals of the complex in a + specific dimension. + + :param dimension: The specific dimension. + :type dimension: int. + :returns: The persistence intervals. + :rtype: numpy array of dimension 2 + + :note: intervals_in_dim function requires persistence function to be + launched first. + """ + cdef vector[pair[double,double]] intervals_result + if self.pcohptr != NULL: + intervals_result = self.pcohptr.intervals_in_dimension(dimension) + else: + print("intervals_in_dim function requires persistence function" + " to be launched first.") + return np_array(intervals_result) diff --git a/src/python/gudhi/persistence_graphical_tools.py b/src/python/gudhi/persistence_graphical_tools.py new file mode 100644 index 00000000..181bc8ea --- /dev/null +++ b/src/python/gudhi/persistence_graphical_tools.py @@ -0,0 +1,423 @@ +from os import path +from math import isfinite +import numpy as np + +from gudhi.reader_utils import read_persistence_intervals_in_dimension +from gudhi.reader_utils import read_persistence_intervals_grouped_by_dimension + +""" 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, Bertrand Michel + + Copyright (C) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau, Bertrand Michel" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + + +def __min_birth_max_death(persistence, band=0.0): + """This function returns (min_birth, max_death) from the persistence. + + :param persistence: The persistence to plot. + :type persistence: list of tuples(dimension, tuple(birth, death)). + :param band: band + :type band: float. + :returns: (float, float) -- (min_birth, max_death). + """ + # Look for minimum birth date and maximum death date for plot optimisation + max_death = 0 + min_birth = persistence[0][1][0] + for interval in reversed(persistence): + if float(interval[1][1]) != float("inf"): + if float(interval[1][1]) > max_death: + max_death = float(interval[1][1]) + if float(interval[1][0]) > max_death: + max_death = float(interval[1][0]) + if float(interval[1][0]) < min_birth: + min_birth = float(interval[1][0]) + if band > 0.0: + max_death += band + return (min_birth, max_death) + + +""" +Only 13 colors for the palette +""" +palette = [ + "#ff0000", + "#00ff00", + "#0000ff", + "#00ffff", + "#ff00ff", + "#ffff00", + "#000000", + "#880000", + "#008800", + "#000088", + "#888800", + "#880088", + "#008888", +] + + +def plot_persistence_barcode( + persistence=[], + persistence_file="", + alpha=0.6, + max_intervals=1000, + max_barcodes=1000, + inf_delta=0.1, + legend=False, +): + """This function plots the persistence bar code from persistence values list + or from a :doc:`persistence file <fileformats>`. + + :param persistence: Persistence intervals values list grouped by dimension. + :type persistence: list of tuples(dimension, tuple(birth, death)). + :param persistence_file: A :doc:`persistence file <fileformats>` style name + (reset persistence if both are set). + :type persistence_file: string + :param alpha: barcode transparency value (0.0 transparent through 1.0 + opaque - default is 0.6). + :type alpha: float. + :param max_intervals: maximal number of intervals to display. + Selected intervals are those with the longest life time. Set it + to 0 to see all. Default value is 1000. + :type max_intervals: int. + :param inf_delta: Infinity is placed at :code:`((max_death - min_birth) x + inf_delta)` above :code:`max_death` value. A reasonable value is + between 0.05 and 0.5 - default is 0.1. + :type inf_delta: float. + :param legend: Display the dimension color legend (default is False). + :type legend: boolean. + :returns: A matplotlib object containing horizontal bar plot of persistence + (launch `show()` method on it to display it). + """ + try: + import matplotlib.pyplot as plt + import matplotlib.patches as mpatches + + if persistence_file is not "": + if path.isfile(persistence_file): + # Reset persistence + persistence = [] + diag = read_persistence_intervals_grouped_by_dimension( + persistence_file=persistence_file + ) + for key in diag.keys(): + for persistence_interval in diag[key]: + persistence.append((key, persistence_interval)) + else: + print("file " + persistence_file + " not found.") + return None + + if max_barcodes is not 1000: + print("Deprecated parameter. It has been replaced by max_intervals") + max_intervals = max_barcodes + + if max_intervals > 0 and max_intervals < len(persistence): + # Sort by life time, then takes only the max_intervals elements + persistence = sorted( + persistence, + key=lambda life_time: life_time[1][1] - life_time[1][0], + reverse=True, + )[:max_intervals] + + persistence = sorted(persistence, key=lambda birth: birth[1][0]) + + (min_birth, max_death) = __min_birth_max_death(persistence) + ind = 0 + delta = (max_death - min_birth) * inf_delta + # Replace infinity values with max_death + delta for bar code to be more + # readable + infinity = max_death + delta + axis_start = min_birth - delta + # Draw horizontal bars in loop + for interval in reversed(persistence): + if float(interval[1][1]) != float("inf"): + # Finite death case + plt.barh( + ind, + (interval[1][1] - interval[1][0]), + height=0.8, + left=interval[1][0], + alpha=alpha, + color=palette[interval[0]], + linewidth=0, + ) + else: + # Infinite death case for diagram to be nicer + plt.barh( + ind, + (infinity - interval[1][0]), + height=0.8, + left=interval[1][0], + alpha=alpha, + color=palette[interval[0]], + linewidth=0, + ) + ind = ind + 1 + + if legend: + dimensions = list(set(item[0] for item in persistence)) + plt.legend( + handles=[ + mpatches.Patch(color=palette[dim], label=str(dim)) + for dim in dimensions + ], + loc="lower right", + ) + plt.title("Persistence barcode") + # Ends plot on infinity value and starts a little bit before min_birth + plt.axis([axis_start, infinity, 0, ind]) + return plt + + except ImportError: + print("This function is not available, you may be missing matplotlib.") + + +def plot_persistence_diagram( + persistence=[], + persistence_file="", + alpha=0.6, + band=0.0, + max_intervals=1000, + max_plots=1000, + inf_delta=0.1, + legend=False, +): + """This function plots the persistence diagram from persistence values + list or from a :doc:`persistence file <fileformats>`. + + :param persistence: Persistence intervals values list grouped by dimension. + :type persistence: list of tuples(dimension, tuple(birth, death)). + :param persistence_file: A :doc:`persistence file <fileformats>` style name + (reset persistence if both are set). + :type persistence_file: string + :param alpha: plot transparency value (0.0 transparent through 1.0 + opaque - default is 0.6). + :type alpha: float. + :param band: band (not displayed if :math:`\leq` 0. - default is 0.) + :type band: float. + :param max_intervals: maximal number of intervals to display. + Selected intervals are those with the longest life time. Set it + to 0 to see all. Default value is 1000. + :type max_intervals: int. + :param inf_delta: Infinity is placed at :code:`((max_death - min_birth) x + inf_delta)` above :code:`max_death` value. A reasonable value is + between 0.05 and 0.5 - default is 0.1. + :type inf_delta: float. + :param legend: Display the dimension color legend (default is False). + :type legend: boolean. + :returns: A matplotlib object containing diagram plot of persistence + (launch `show()` method on it to display it). + """ + try: + import matplotlib.pyplot as plt + import matplotlib.patches as mpatches + + if persistence_file is not "": + if path.isfile(persistence_file): + # Reset persistence + persistence = [] + diag = read_persistence_intervals_grouped_by_dimension( + persistence_file=persistence_file + ) + for key in diag.keys(): + for persistence_interval in diag[key]: + persistence.append((key, persistence_interval)) + else: + print("file " + persistence_file + " not found.") + return None + + if max_plots is not 1000: + print("Deprecated parameter. It has been replaced by max_intervals") + max_intervals = max_plots + + if max_intervals > 0 and max_intervals < len(persistence): + # Sort by life time, then takes only the max_intervals elements + persistence = sorted( + persistence, + key=lambda life_time: life_time[1][1] - life_time[1][0], + reverse=True, + )[:max_intervals] + + (min_birth, max_death) = __min_birth_max_death(persistence, band) + delta = (max_death - min_birth) * inf_delta + # Replace infinity values with max_death + delta for diagram to be more + # readable + infinity = max_death + delta + axis_start = min_birth - delta + + # line display of equation : birth = death + x = np.linspace(axis_start, infinity, 1000) + # infinity line and text + plt.plot(x, x, color="k", linewidth=1.0) + plt.plot(x, [infinity] * len(x), linewidth=1.0, color="k", alpha=alpha) + plt.text(axis_start, infinity, r"$\infty$", color="k", alpha=alpha) + # bootstrap band + if band > 0.0: + plt.fill_between(x, x, x + band, alpha=alpha, facecolor="red") + + # Draw points in loop + for interval in reversed(persistence): + if float(interval[1][1]) != float("inf"): + # Finite death case + plt.scatter( + interval[1][0], + interval[1][1], + alpha=alpha, + color=palette[interval[0]], + ) + else: + # Infinite death case for diagram to be nicer + plt.scatter( + interval[1][0], infinity, alpha=alpha, color=palette[interval[0]] + ) + + if legend: + dimensions = list(set(item[0] for item in persistence)) + plt.legend( + handles=[ + mpatches.Patch(color=palette[dim], label=str(dim)) + for dim in dimensions + ] + ) + + plt.title("Persistence diagram") + plt.xlabel("Birth") + plt.ylabel("Death") + # Ends plot on infinity value and starts a little bit before min_birth + plt.axis([axis_start, infinity, axis_start, infinity + delta]) + return plt + + except ImportError: + print("This function is not available, you may be missing matplotlib.") + + +def plot_persistence_density( + persistence=[], + persistence_file="", + nbins=300, + bw_method=None, + max_intervals=1000, + dimension=None, + cmap=None, + legend=False, +): + """This function plots the persistence density from persistence + values list or from a :doc:`persistence file <fileformats>`. Be + aware that this function does not distinguish the dimension, it is + up to you to select the required one. This function also does not handle + degenerate data set (scipy correlation matrix inversion can fail). + + :param persistence: Persistence intervals values list grouped by dimension. + :type persistence: list of tuples(dimension, tuple(birth, death)). + :param persistence_file: A :doc:`persistence file <fileformats>` + style name (reset persistence if both are set). + :type persistence_file: string + :param nbins: Evaluate a gaussian kde on a regular grid of nbins x + nbins over data extents (default is 300) + :type nbins: int. + :param bw_method: The method used to calculate the estimator + bandwidth. This can be 'scott', 'silverman', a scalar constant + or a callable. If a scalar, this will be used directly as + kde.factor. If a callable, it should take a gaussian_kde + instance as only parameter and return a scalar. If None + (default), 'scott' is used. See + `scipy.stats.gaussian_kde documentation + <http://scipy.github.io/devdocs/generated/scipy.stats.gaussian_kde.html>`_ + for more details. + :type bw_method: str, scalar or callable, optional. + :param max_intervals: maximal number of points used in the density + estimation. + Selected intervals are those with the longest life time. Set it + to 0 to see all. Default value is 1000. + :type max_intervals: int. + :param dimension: the dimension to be selected in the intervals + (default is None to mix all dimensions). + :type dimension: int. + :param cmap: A matplotlib colormap (default is + matplotlib.pyplot.cm.hot_r). + :type cmap: cf. matplotlib colormap. + :param legend: Display the color bar values (default is False). + :type legend: boolean. + :returns: A matplotlib object containing diagram plot of persistence + (launch `show()` method on it to display it). + """ + try: + import matplotlib.pyplot as plt + from scipy.stats import kde + + if persistence_file is not "": + if dimension is None: + # All dimension case + dimension = -1 + if path.isfile(persistence_file): + persistence_dim = read_persistence_intervals_in_dimension( + persistence_file=persistence_file, only_this_dim=dimension + ) + print(persistence_dim) + else: + print("file " + persistence_file + " not found.") + return None + + if len(persistence) > 0: + persistence_dim = np.array( + [ + (dim_interval[1][0], dim_interval[1][1]) + for dim_interval in persistence + if (dim_interval[0] == dimension) or (dimension is None) + ] + ) + + persistence_dim = persistence_dim[np.isfinite(persistence_dim[:, 1])] + if max_intervals > 0 and max_intervals < len(persistence_dim): + # Sort by life time, then takes only the max_intervals elements + persistence_dim = np.array( + sorted( + persistence_dim, + key=lambda life_time: life_time[1] - life_time[0], + reverse=True, + )[:max_intervals] + ) + + # Set as numpy array birth and death (remove undefined values - inf and NaN) + birth = persistence_dim[:, 0] + death = persistence_dim[:, 1] + + # line display of equation : birth = death + x = np.linspace(death.min(), birth.max(), 1000) + plt.plot(x, x, color="k", linewidth=1.0) + + # Evaluate a gaussian kde on a regular grid of nbins x nbins over data extents + k = kde.gaussian_kde([birth, death], bw_method=bw_method) + xi, yi = np.mgrid[ + birth.min() : birth.max() : nbins * 1j, + death.min() : death.max() : nbins * 1j, + ] + zi = k(np.vstack([xi.flatten(), yi.flatten()])) + + # default cmap value cannot be done at argument definition level as matplotlib is not yet defined. + if cmap is None: + cmap = plt.cm.hot_r + # Make the plot + plt.pcolormesh(xi, yi, zi.reshape(xi.shape), cmap=cmap) + + if legend: + plt.colorbar() + + plt.title("Persistence density") + plt.xlabel("Birth") + plt.ylabel("Death") + return plt + + except ImportError: + print( + "This function is not available, you may be missing matplotlib and/or scipy." + ) diff --git a/src/python/gudhi/reader_utils.pyx b/src/python/gudhi/reader_utils.pyx new file mode 100644 index 00000000..147fae71 --- /dev/null +++ b/src/python/gudhi/reader_utils.pyx @@ -0,0 +1,87 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.string cimport string +from libcpp.map cimport map +from libcpp.pair cimport pair + +from os import path +from numpy import array as np_array + +""" 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) 2017 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2017 Inria" +__license__ = "MIT" + +cdef extern from "Reader_utils_interface.h" namespace "Gudhi": + vector[vector[double]] read_matrix_from_csv_file(string off_file, char separator) + map[int, vector[pair[double, double]]] read_pers_intervals_grouped_by_dimension(string filename) + vector[pair[double, double]] read_pers_intervals_in_dimension(string filename, int only_this_dim) + +def read_lower_triangular_matrix_from_csv_file(csv_file='', separator=';'): + """Read lower triangular matrix from a CSV style file. + + :param csv_file: A CSV file style name. + :type csv_file: string + :param separator: The value separator in the CSV file. Default value is ';' + :type separator: char + + :returns: The lower triangular matrix. + :rtype: vector[vector[double]] + """ + if csv_file is not '': + if path.isfile(csv_file): + return read_matrix_from_csv_file(str.encode(csv_file), ord(separator[0])) + print("file " + csv_file + " not set or not found.") + return [] + +def read_persistence_intervals_grouped_by_dimension(persistence_file=''): + """Reads a file containing persistence intervals. + Each line might contain 2, 3 or 4 values: [[field] dimension] birth death + The return value is an `map[dim, vector[pair[birth, death]]]` + where `dim` is an `int`, `birth` a `double`, and `death` a `double`. + Note: the function does not check that birth <= death. + + :param persistence_file: A persistence file style name. + :type persistence_file: string + + :returns: The persistence pairs grouped by dimension. + :rtype: map[int, vector[pair[double, double]]] + """ + if persistence_file is not '': + if path.isfile(persistence_file): + return read_pers_intervals_grouped_by_dimension(str.encode(persistence_file)) + print("file " + persistence_file + " not set or not found.") + return [] + +def read_persistence_intervals_in_dimension(persistence_file='', only_this_dim=-1): + """Reads a file containing persistence intervals. + Each line of persistence_file might contain 2, 3 or 4 values: + [[field] dimension] birth death + Note: the function does not check that birth <= death. + + :param persistence_file: A persistence file style name. + :type persistence_file: string + :param only_this_dim: The specific dimension. Default value is -1. + If `only_this_dim` = -1, dimension is ignored and all lines are returned. + If `only_this_dim` is >= 0, only the lines where dimension = + `only_this_dim` (or where dimension is not specified) are returned. + :type only_this_dim: int. + + :returns: The persistence intervals. + :rtype: numpy array of dimension 2 + """ + if persistence_file is not '': + if path.isfile(persistence_file): + return np_array(read_pers_intervals_in_dimension(str.encode( + persistence_file), only_this_dim)) + print("file " + persistence_file + " not set or not found.") + return [] diff --git a/src/python/gudhi/rips_complex.pyx b/src/python/gudhi/rips_complex.pyx new file mode 100644 index 00000000..f2cd6a8d --- /dev/null +++ b/src/python/gudhi/rips_complex.pyx @@ -0,0 +1,103 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp.string cimport string +from libcpp cimport bool +from libc.stdint cimport intptr_t + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Rips_complex_interface.h" namespace "Gudhi": + cdef cppclass Rips_complex_interface "Gudhi::rips_complex::Rips_complex_interface": + Rips_complex_interface() + void init_points(vector[vector[double]] values, double threshold) + void init_matrix(vector[vector[double]] values, double threshold) + void init_points_sparse(vector[vector[double]] values, double threshold, double sparse) + void init_matrix_sparse(vector[vector[double]] values, double threshold, double sparse) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, int dim_max) + +# RipsComplex python interface +cdef class RipsComplex: + """The data structure is a one skeleton graph, or Rips graph, containing + edges when the edge length is less or equal to a given threshold. Edge + length is computed from a user given point cloud with a given distance + function, or a distance matrix. + """ + + cdef Rips_complex_interface thisref + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, points=None, distance_matrix=None, + max_edge_length=float('inf'), sparse=None): + """RipsComplex constructor. + + :param max_edge_length: Rips value. + :type max_edge_length: float + + :param points: A list of points in d-Dimension. + :type points: list of list of double + + Or + + :param distance_matrix: A distance matrix (full square or lower + triangular). + :type points: list of list of double + + And in both cases + + :param sparse: If this is not None, it switches to building a sparse + Rips and represents the approximation parameter epsilon. + :type sparse: float + """ + + # The real cython constructor + def __cinit__(self, points=None, distance_matrix=None, + max_edge_length=float('inf'), sparse=None): + if sparse is not None: + if distance_matrix is not None: + self.thisref.init_matrix_sparse(distance_matrix, + max_edge_length, + sparse) + else: + if points is None: + # Empty Rips construction + points=[] + self.thisref.init_points_sparse(points, max_edge_length, sparse) + else: + if distance_matrix is not None: + self.thisref.init_matrix(distance_matrix, max_edge_length) + else: + if points is None: + # Empty Rips construction + points=[] + self.thisref.init_points(points, max_edge_length) + + + def create_simplex_tree(self, max_dimension=1): + """ + :param max_dimension: graph expansion for rips until this given maximal + dimension. + :type max_dimension: int + :returns: A simplex tree created from the Delaunay Triangulation. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + self.thisref.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_dimension) + return stree diff --git a/src/python/gudhi/simplex_tree.pxd b/src/python/gudhi/simplex_tree.pxd new file mode 100644 index 00000000..5f86cfe2 --- /dev/null +++ b/src/python/gudhi/simplex_tree.pxd @@ -0,0 +1,56 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp cimport bool +from libcpp.string cimport string + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Simplex_tree_interface.h" namespace "Gudhi": + cdef cppclass Simplex_tree_options_full_featured: + pass + + 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) + void assign_simplex_filtration(vector[int] simplex, double filtration) + void initialize_filtration() + int num_vertices() + int num_simplices() + void set_dimension(int dimension) + int dimension() + int upper_bound_dimension() + 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) + void expansion(int max_dim) + void remove_maximal_simplex(vector[int] simplex) + bool prune_above_filtration(double filtration) + bool make_filtration_non_decreasing() + +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>>": + Simplex_tree_persistence_interface(Simplex_tree_interface_full_featured * st, bool persistence_dim_max) + vector[pair[int, pair[double, double]]] get_persistence(int homology_coeff_field, double min_persistence) + vector[int] betti_numbers() + vector[int] persistent_betti_numbers(double from_value, double to_value) + vector[pair[double,double]] intervals_in_dimension(int dimension) + void write_output_diagram(string diagram_file_name) + vector[pair[vector[int], vector[int]]] persistence_pairs() diff --git a/src/python/gudhi/simplex_tree.pyx b/src/python/gudhi/simplex_tree.pyx new file mode 100644 index 00000000..9f490271 --- /dev/null +++ b/src/python/gudhi/simplex_tree.pyx @@ -0,0 +1,508 @@ +from libc.stdint cimport intptr_t +from numpy import array as np_array +cimport simplex_tree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +# SimplexTree python interface +cdef class SimplexTree: + """The simplex tree is an efficient and flexible data structure for + representing general (filtered) simplicial complexes. The data structure + is described in Jean-Daniel Boissonnat and Clément Maria. The Simplex + Tree: An Efficient Data Structure for General Simplicial Complexes. + Algorithmica, pages 1–22, 2014. + + This class is a filtered, with keys, and non contiguous vertices version + of the simplex tree. + """ + # unfortunately 'cdef public Simplex_tree_interface_full_featured* thisptr' is not possible + # Use intptr_t instead to cast the pointer + cdef public intptr_t thisptr + + # Get the pointer casted as it should be + cdef Simplex_tree_interface_full_featured* get_ptr(self): + return <Simplex_tree_interface_full_featured*>(self.thisptr) + + cdef Simplex_tree_persistence_interface * pcohptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self): + """SimplexTree constructor. + """ + + # The real cython constructor + def __cinit__(self): + self.thisptr = <intptr_t>(new Simplex_tree_interface_full_featured()) + + def __dealloc__(self): + cdef Simplex_tree_interface_full_featured* ptr = self.get_ptr() + if ptr != NULL: + del ptr + if self.pcohptr != NULL: + del self.pcohptr + + def __is_defined(self): + """Returns true if SimplexTree pointer is not NULL. + """ + return self.get_ptr() != NULL + + def __is_persistence_defined(self): + """Returns true if Persistence pointer is not NULL. + """ + return self.pcohptr != NULL + + def filtration(self, simplex): + """This function returns the filtration value for a given N-simplex in + this simplicial complex, or +infinity if it is not in the complex. + + :param simplex: The N-simplex, represented by a list of vertex. + :type simplex: list of int. + :returns: The simplicial complex filtration value. + :rtype: float + """ + return self.get_ptr().simplex_filtration(simplex) + + def assign_filtration(self, simplex, filtration): + """This function assigns the simplicial complex filtration value for a + given N-simplex. + + :param simplex: The N-simplex, represented by a list of vertex. + :type simplex: list of int. + :param filtration: The simplicial complex filtration value. + :type filtration: float + """ + self.get_ptr().assign_simplex_filtration(simplex, filtration) + + def initialize_filtration(self): + """This function initializes and sorts the simplicial complex + filtration vector. + + .. note:: + + This function must be launched before + :func:`persistence()<gudhi.SimplexTree.persistence>`, + :func:`betti_numbers()<gudhi.SimplexTree.betti_numbers>`, + :func:`persistent_betti_numbers()<gudhi.SimplexTree.persistent_betti_numbers>`, + or :func:`get_filtration()<gudhi.SimplexTree.get_filtration>` + after :func:`inserting<gudhi.SimplexTree.insert>` or + :func:`removing<gudhi.SimplexTree.remove_maximal_simplex>` + simplices. + """ + self.get_ptr().initialize_filtration() + + def num_vertices(self): + """This function returns the number of vertices of the simplicial + complex. + + :returns: The simplicial complex number of vertices. + :rtype: int + """ + return self.get_ptr().num_vertices() + + def num_simplices(self): + """This function returns the number of simplices of the simplicial + complex. + + :returns: the simplicial complex number of simplices. + :rtype: int + """ + return self.get_ptr().num_simplices() + + def dimension(self): + """This function returns the dimension of the simplicial complex. + + :returns: the simplicial complex dimension. + :rtype: int + + .. note:: + + This function is not constant time because it can recompute + dimension if required (can be triggered by + :func:`remove_maximal_simplex()<gudhi.SimplexTree.remove_maximal_simplex>` + or + :func:`prune_above_filtration()<gudhi.SimplexTree.prune_above_filtration>` + methods). + """ + return self.get_ptr().dimension() + + def upper_bound_dimension(self): + """This function returns a valid dimension upper bound of the + simplicial complex. + + :returns: an upper bound on the dimension of the simplicial complex. + :rtype: int + """ + return self.get_ptr().upper_bound_dimension() + + def set_dimension(self, dimension): + """This function sets the dimension of the simplicial complex. + + :param dimension: The new dimension value. + :type dimension: int. + + .. note:: + + This function must be used with caution because it disables + dimension recomputation when required + (this recomputation can be triggered by + :func:`remove_maximal_simplex()<gudhi.SimplexTree.remove_maximal_simplex>` + or + :func:`prune_above_filtration()<gudhi.SimplexTree.prune_above_filtration>` + ). + """ + self.get_ptr().set_dimension(<int>dimension) + + def find(self, simplex): + """This function returns if the N-simplex was found in the simplicial + complex or not. + + :param simplex: The N-simplex to find, represented by a list of vertex. + :type simplex: list of int. + :returns: true if the simplex was found, false otherwise. + :rtype: bool + """ + cdef vector[int] csimplex + for i in simplex: + csimplex.push_back(i) + return self.get_ptr().find_simplex(csimplex) + + def insert(self, simplex, filtration=0.0): + """This function inserts the given N-simplex and its subfaces with the + given filtration value (default value is '0.0'). If some of those + simplices are already present with a higher filtration value, their + filtration value is lowered. + + :param simplex: The N-simplex to insert, represented by a list of + vertex. + :type simplex: list of int. + :param filtration: The filtration value of the simplex. + :type filtration: float. + :returns: true if the simplex was not yet in the complex, false + otherwise (whatever its original filtration value). + :rtype: bool + """ + cdef vector[int] csimplex + for i in simplex: + csimplex.push_back(i) + 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 + filtration values. + + :returns: The simplices sorted by increasing filtration values. + :rtype: list of 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 + + def get_skeleton(self, dimension): + """This function returns the (simplices of the) skeleton of a maximum + given dimension. + + :param dimension: The skeleton dimension value. + :type dimension: int. + :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 + + def get_star(self, simplex): + """This function returns the star of a given N-simplex. + + :param simplex: The N-simplex, represented by a list of vertex. + :type simplex: list of int. + :returns: The (simplices of the) star of a simplex. + :rtype: list of tuples(simplex, filtration) + """ + cdef vector[int] csimplex + for i in simplex: + csimplex.push_back(i) + cdef vector[pair[vector[int], double]] star \ + = self.get_ptr().get_star(csimplex) + ct = [] + for filtered_simplex in star: + v = [] + for vertex in filtered_simplex.first: + v.append(vertex) + ct.append((v, filtered_simplex.second)) + return ct + + def get_cofaces(self, simplex, codimension): + """This function returns the cofaces of a given N-simplex with a + given codimension. + + :param simplex: The N-simplex, represented by a list of vertex. + :type simplex: list of int. + :param codimension: The codimension. If codimension = 0, all cofaces + are returned (equivalent of get_star function) + :type codimension: int. + :returns: The (simplices of the) cofaces of a simplex + :rtype: list of tuples(simplex, filtration) + """ + cdef vector[int] csimplex + for i in simplex: + csimplex.push_back(i) + cdef vector[pair[vector[int], double]] cofaces \ + = self.get_ptr().get_cofaces(csimplex, <int>codimension) + ct = [] + for filtered_simplex in cofaces: + v = [] + for vertex in filtered_simplex.first: + v.append(vertex) + ct.append((v, filtered_simplex.second)) + return ct + + def remove_maximal_simplex(self, simplex): + """This function removes a given maximal N-simplex from the simplicial + complex. + + :param simplex: The N-simplex, represented by a list of vertex. + :type simplex: list of int. + + .. note:: + + Be aware that removing is shifting data in a flat_map + (:func:`initialize_filtration()<gudhi.SimplexTree.initialize_filtration>` to be done). + + .. note:: + + The dimension of the simplicial complex may be lower after calling + remove_maximal_simplex than it was before. However, + :func:`upper_bound_dimension()<gudhi.SimplexTree.upper_bound_dimension>` + method will return the old value, which + remains a valid upper bound. If you care, you can call + :func:`dimension()<gudhi.SimplexTree.dimension>` + to recompute the exact dimension. + """ + self.get_ptr().remove_maximal_simplex(simplex) + + def prune_above_filtration(self, filtration): + """Prune above filtration value given as parameter. + + :param filtration: Maximum threshold value. + :type filtration: float. + :returns: The filtration modification information. + :rtype: bool + + + .. note:: + + Some simplex tree functions require the filtration to be valid. + prune_above_filtration function is not launching + :func:`initialize_filtration()<gudhi.SimplexTree.initialize_filtration>` + but returns the filtration modification + information. If the complex has changed , please call + :func:`initialize_filtration()<gudhi.SimplexTree.initialize_filtration>` + to recompute it. + + .. note:: + + Note that the dimension of the simplicial complex may be lower + after calling + :func:`prune_above_filtration()<gudhi.SimplexTree.prune_above_filtration>` + than it was before. However, + :func:`upper_bound_dimension()<gudhi.SimplexTree.upper_bound_dimension>` + will return the old value, which remains a + valid upper bound. If you care, you can call + :func:`dimension()<gudhi.SimplexTree.dimension>` + method to recompute the exact dimension. + """ + return self.get_ptr().prune_above_filtration(filtration) + + def expansion(self, max_dim): + """Expands the Simplex_tree containing only its one skeleton + until dimension max_dim. + + The expanded simplicial complex until dimension :math:`d` + attached to a graph :math:`G` is the maximal simplicial complex of + dimension at most :math:`d` admitting the graph :math:`G` as + :math:`1`-skeleton. + The filtration value assigned to a simplex is the maximal filtration + value of one of its edges. + + The Simplex_tree must contain no simplex of dimension bigger than + 1 when calling the method. + + :param max_dim: The maximal dimension. + :type max_dim: int. + """ + self.get_ptr().expansion(max_dim) + + def make_filtration_non_decreasing(self): + """This function ensures that each simplex has a higher filtration + value than its faces by increasing the filtration values. + + :returns: True if any filtration value was modified, + False if the filtration was already non-decreasing. + :rtype: bool + + + .. note:: + + Some simplex tree functions require the filtration to be valid. + make_filtration_non_decreasing function is not launching + :func:`initialize_filtration()<gudhi.SimplexTree.initialize_filtration>` + but returns the filtration modification + information. If the complex has changed , please call + :func:`initialize_filtration()<gudhi.SimplexTree.initialize_filtration>` + to recompute it. + """ + return self.get_ptr().make_filtration_non_decreasing() + + def persistence(self, homology_coeff_field=11, min_persistence=0, persistence_dim_max = False): + """This function returns the persistence of the simplicial complex. + + :param homology_coeff_field: The homology coefficient field. Must be a + prime number. Default value is 11. + :type homology_coeff_field: int. + :param min_persistence: The minimum persistence value to take into + account (strictly greater than min_persistence). Default value is + 0.0. + Sets min_persistence to -1.0 to see all values. + :type min_persistence: float. + :param persistence_dim_max: If true, the persistent homology for the + maximal dimension in the complex is computed. If false, it is + ignored. Default is false. + :type persistence_dim_max: bool + :returns: The persistence of the simplicial complex. + :rtype: list of pairs(dimension, pair(birth, death)) + """ + if self.pcohptr != NULL: + del self.pcohptr + self.pcohptr = new Simplex_tree_persistence_interface(self.get_ptr(), persistence_dim_max) + cdef vector[pair[int, pair[double, double]]] persistence_result + if self.pcohptr != NULL: + persistence_result = self.pcohptr.get_persistence(homology_coeff_field, min_persistence) + return persistence_result + + def betti_numbers(self): + """This function returns the Betti numbers of the simplicial complex. + + :returns: The Betti numbers ([B0, B1, ..., Bn]). + :rtype: list of int + + :note: betti_numbers function requires + :func:`persistence()<gudhi.SimplexTree.persistence>` + function to be launched first. + """ + cdef vector[int] bn_result + if self.pcohptr != NULL: + bn_result = self.pcohptr.betti_numbers() + else: + print("betti_numbers function requires persistence function" + " to be launched first.") + return bn_result + + def persistent_betti_numbers(self, from_value, to_value): + """This function returns the persistent Betti numbers of the + simplicial complex. + + :param from_value: The persistence birth limit to be added in the + numbers (persistent birth <= from_value). + :type from_value: float. + :param to_value: The persistence death limit to be added in the + numbers (persistent death > to_value). + :type to_value: float. + + :returns: The persistent Betti numbers ([B0, B1, ..., Bn]). + :rtype: list of int + + :note: persistent_betti_numbers function requires + :func:`persistence()<gudhi.SimplexTree.persistence>` + function to be launched first. + """ + cdef vector[int] pbn_result + if self.pcohptr != NULL: + pbn_result = self.pcohptr.persistent_betti_numbers(<double>from_value, <double>to_value) + else: + print("persistent_betti_numbers function requires persistence function" + " to be launched first.") + return pbn_result + + def persistence_intervals_in_dimension(self, dimension): + """This function returns the persistence intervals of the simplicial + complex in a specific dimension. + + :param dimension: The specific dimension. + :type dimension: int. + :returns: The persistence intervals. + :rtype: numpy array of dimension 2 + + :note: intervals_in_dim function requires + :func:`persistence()<gudhi.SimplexTree.persistence>` + function to be launched first. + """ + cdef vector[pair[double,double]] intervals_result + if self.pcohptr != NULL: + intervals_result = self.pcohptr.intervals_in_dimension(dimension) + else: + print("intervals_in_dim function requires persistence function" + " to be launched first.") + return np_array(intervals_result) + + def persistence_pairs(self): + """This function returns a list of persistence birth and death simplices pairs. + + :returns: A list of persistence simplices intervals. + :rtype: list of pair of list of int + + :note: persistence_pairs function requires + :func:`persistence()<gudhi.SimplexTree.persistence>` + function to be launched first. + """ + cdef vector[pair[vector[int],vector[int]]] persistence_pairs_result + if self.pcohptr != NULL: + persistence_pairs_result = self.pcohptr.persistence_pairs() + else: + print("persistence_pairs function requires persistence function" + " to be launched first.") + return persistence_pairs_result + + def write_persistence_diagram(self, persistence_file=''): + """This function writes the persistence intervals of the simplicial + complex in a user given file name. + + :param persistence_file: The specific dimension. + :type persistence_file: string. + + :note: intervals_in_dim function requires + :func:`persistence()<gudhi.SimplexTree.persistence>` + function to be launched first. + """ + if self.pcohptr != NULL: + if persistence_file != '': + self.pcohptr.write_output_diagram(str.encode(persistence_file)) + else: + print("persistence_file must be specified") + else: + print("intervals_in_dim function requires persistence function" + " to be launched first.") diff --git a/src/python/gudhi/strong_witness_complex.pyx b/src/python/gudhi/strong_witness_complex.pyx new file mode 100644 index 00000000..e757abea --- /dev/null +++ b/src/python/gudhi/strong_witness_complex.pyx @@ -0,0 +1,78 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libc.stdint cimport intptr_t + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Strong_witness_complex_interface.h" namespace "Gudhi": + cdef cppclass Strong_witness_complex_interface "Gudhi::witness_complex::Strong_witness_complex_interface": + Strong_witness_complex_interface(vector[vector[pair[size_t, double]]] nearest_landmark_table) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square, + unsigned limit_dimension) + +# StrongWitnessComplex python interface +cdef class StrongWitnessComplex: + """Constructs (strong) witness complex for a given table of nearest + landmarks with respect to witnesses. + """ + + cdef Strong_witness_complex_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, nearest_landmark_table=None): + """StrongWitnessComplex constructor. + + :param nearest_landmark_table: A list of lists of nearest landmarks and their distances. + `nearest_landmark_table[w][k]==(l,d)` means that l is the k-th nearest landmark to + witness w, and d is the (squared) distance between l and w. + :type nearest_landmark_table: list of list of pair of int and float + """ + + # The real cython constructor + def __cinit__(self, nearest_landmark_table=None): + if nearest_landmark_table is not None: + self.thisptr = new Strong_witness_complex_interface(nearest_landmark_table) + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if StrongWitnessComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def create_simplex_tree(self, max_alpha_square = float('inf'), limit_dimension = -1): + """ + :param max_alpha_square: The maximum relaxation parameter. + Default is set to infinity. + :type max_alpha_square: float + :returns: A simplex tree created from the Delaunay Triangulation. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + if limit_dimension is not -1: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square, limit_dimension) + else: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square) + return stree diff --git a/src/python/gudhi/subsampling.pyx b/src/python/gudhi/subsampling.pyx new file mode 100644 index 00000000..1135c1fb --- /dev/null +++ b/src/python/gudhi/subsampling.pyx @@ -0,0 +1,130 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.string cimport string +from libcpp cimport bool +import os + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "GPL v3" + +cdef extern from "Subsampling_interface.h" namespace "Gudhi::subsampling": + vector[vector[double]] subsampling_n_farthest_points(vector[vector[double]] points, unsigned nb_points) + vector[vector[double]] subsampling_n_farthest_points(vector[vector[double]] points, unsigned nb_points, unsigned starting_point) + vector[vector[double]] subsampling_n_farthest_points_from_file(string off_file, unsigned nb_points) + vector[vector[double]] subsampling_n_farthest_points_from_file(string off_file, unsigned nb_points, unsigned starting_point) + vector[vector[double]] subsampling_n_random_points(vector[vector[double]] points, unsigned nb_points) + vector[vector[double]] subsampling_n_random_points_from_file(string off_file, unsigned nb_points) + vector[vector[double]] subsampling_sparsify_points(vector[vector[double]] points, double min_squared_dist) + vector[vector[double]] subsampling_sparsify_points_from_file(string off_file, double min_squared_dist) + +def choose_n_farthest_points(points=None, off_file='', nb_points=0, starting_point = ''): + """Subsample by a greedy strategy of iteratively adding the farthest point + from the current chosen point set to the subsampling. + The iteration starts with the landmark `starting point`. + + :param points: The input point set. + :type points: vector[vector[double]]. + + Or + + :param off_file: An OFF file style name. + :type off_file: string + + :param nb_points: Number of points of the subsample. + :type nb_points: unsigned. + :param starting_point: The iteration starts with the landmark `starting \ + point`,which is the index of the poit to start with. If not set, this \ + index is choosen randomly. + :type starting_point: unsigned. + :returns: The subsample point set. + :rtype: vector[vector[double]] + """ + if off_file is not '': + if os.path.isfile(off_file): + if starting_point is '': + return subsampling_n_farthest_points_from_file(str.encode(off_file), + nb_points) + else: + return subsampling_n_farthest_points_from_file(str.encode(off_file), + nb_points, + starting_point) + else: + print("file " + off_file + " not found.") + else: + if points is None: + # Empty points + points=[] + if starting_point is '': + return subsampling_n_farthest_points(points, nb_points) + else: + return subsampling_n_farthest_points(points, nb_points, + starting_point) + +def pick_n_random_points(points=None, off_file='', nb_points=0): + """Subsample a point set by picking random vertices. + + :param points: The input point set. + :type points: vector[vector[double]]. + + Or + + :param off_file: An OFF file style name. + :type off_file: string + + :param nb_points: Number of points of the subsample. + :type nb_points: unsigned. + :returns: The subsample point set. + :rtype: vector[vector[double]] + """ + if off_file is not '': + if os.path.isfile(off_file): + return subsampling_n_random_points_from_file(str.encode(off_file), + nb_points) + else: + print("file " + off_file + " not found.") + else: + if points is None: + # Empty points + points=[] + return subsampling_n_random_points(points, nb_points) + +def sparsify_point_set(points=None, off_file='', min_squared_dist=0.0): + """Outputs a subset of the input points so that the squared distance + between any two points is greater than or equal to min_squared_dist. + + :param points: The input point set. + :type points: vector[vector[double]]. + + Or + + :param off_file: An OFF file style name. + :type off_file: string + + :param min_squared_dist: Minimum squared distance separating the output \ + points. + :type min_squared_dist: float. + :returns: The subsample point set. + :rtype: vector[vector[double]] + """ + if off_file is not '': + if os.path.isfile(off_file): + return subsampling_sparsify_points_from_file(str.encode(off_file), + min_squared_dist) + else: + print("file " + off_file + " not found.") + else: + if points is None: + # Empty points + points=[] + return subsampling_sparsify_points(points, min_squared_dist) diff --git a/src/python/gudhi/tangential_complex.pyx b/src/python/gudhi/tangential_complex.pyx new file mode 100644 index 00000000..3a945fe2 --- /dev/null +++ b/src/python/gudhi/tangential_complex.pyx @@ -0,0 +1,173 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libcpp.string cimport string +from libcpp cimport bool +from libc.stdint cimport intptr_t +import os + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "GPL v3" + +cdef extern from "Tangential_complex_interface.h" namespace "Gudhi": + cdef cppclass Tangential_complex_interface "Gudhi::tangential_complex::Tangential_complex_interface": + Tangential_complex_interface(int intrisic_dim, vector[vector[double]] points) + # bool from_file is a workaround for cython to find the correct signature + Tangential_complex_interface(int intrisic_dim, string off_file, bool from_file) + void compute_tangential_complex() except + + vector[double] get_point(unsigned vertex) + unsigned number_of_vertices() + unsigned number_of_simplices() + unsigned number_of_inconsistent_simplices() + unsigned number_of_inconsistent_stars() + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree) + void fix_inconsistencies_using_perturbation(double max_perturb, double time_limit) + void set_max_squared_edge_length(double max_squared_edge_length) + +# TangentialComplex python interface +cdef class TangentialComplex: + """The class Tangential_complex represents a tangential complex. After the + computation of the complex, an optional post-processing called perturbation + can be run to attempt to remove inconsistencies. + """ + + cdef Tangential_complex_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, intrisic_dim, points=None, off_file=''): + """TangentialComplex constructor. + + :param intrisic_dim: Intrinsic dimension of the manifold. + :type intrisic_dim: integer + + :param points: A list of points in d-Dimension. + :type points: list of list of double + + Or + + :param off_file: An OFF file style name. + :type off_file: string + """ + + # The real cython constructor + def __cinit__(self, intrisic_dim, points=None, off_file=''): + if off_file is not '': + if os.path.isfile(off_file): + self.thisptr = new Tangential_complex_interface(intrisic_dim, str.encode(off_file), True) + else: + print("file " + off_file + " not found.") + else: + if points is None: + # Empty tangential construction + points=[] + self.thisptr = new Tangential_complex_interface(intrisic_dim, points) + + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if TangentialComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def compute_tangential_complex(self): + """This function computes the tangential complex. + + Raises: + ValueError: In debug mode, if the computed star dimension is too + low. Try to set a bigger maximal edge length value with + :func:`~gudhi.Tangential_complex.set_max_squared_edge_length` + if this happens. + """ + self.thisptr.compute_tangential_complex() + + def get_point(self, vertex): + """This function returns the point corresponding to a given vertex. + + :param vertex: The vertex. + :type vertex: int. + :returns: The point. + :rtype: list of float + """ + cdef vector[double] point = self.thisptr.get_point(vertex) + return point + + def num_vertices(self): + """ + :returns: The number of vertices. + :rtype: unsigned + """ + return self.thisptr.number_of_vertices() + + def num_simplices(self): + """ + :returns: Total number of simplices in stars (including duplicates that appear in several stars). + :rtype: unsigned + """ + return self.thisptr.number_of_simplices() + + def num_inconsistent_simplices(self): + """ + :returns: The number of inconsistent simplices. + :rtype: unsigned + """ + return self.thisptr.number_of_inconsistent_simplices() + + def num_inconsistent_stars(self): + """ + :returns: The number of stars containing at least one inconsistent simplex. + :rtype: unsigned + """ + return self.thisptr.number_of_inconsistent_stars() + + def create_simplex_tree(self): + """Exports the complex into a simplex tree. + + :returns: A simplex tree created from the complex. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr) + return stree + + def fix_inconsistencies_using_perturbation(self, max_perturb, time_limit=-1.0): + """Attempts to fix inconsistencies by perturbing the point positions. + + :param max_perturb: Maximum length of the translations used by the + perturbation. + :type max_perturb: double + :param time_limit: Time limit in seconds. If -1, no time limit is set. + :type time_limit: double + """ + self.thisptr.fix_inconsistencies_using_perturbation(max_perturb, + time_limit) + + def set_max_squared_edge_length(self, max_squared_edge_length): + """Sets the maximal possible squared edge length for the edges in the + triangulations. + + :param max_squared_edge_length: Maximal possible squared edge length. + :type max_squared_edge_length: double + + If the maximal edge length value is too low + :func:`~gudhi.Tangential_complex.compute_tangential_complex` + will throw an exception in debug mode. + """ + self.thisptr.set_max_squared_edge_length(max_squared_edge_length) diff --git a/src/python/gudhi/witness_complex.pyx b/src/python/gudhi/witness_complex.pyx new file mode 100644 index 00000000..baa70b7a --- /dev/null +++ b/src/python/gudhi/witness_complex.pyx @@ -0,0 +1,78 @@ +from cython cimport numeric +from libcpp.vector cimport vector +from libcpp.utility cimport pair +from libc.stdint cimport intptr_t + +from gudhi.simplex_tree cimport * +from gudhi.simplex_tree import SimplexTree + +""" 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) 2016 Inria + + Modification(s): + - YYYY/MM Author: Description of the modification +""" + +__author__ = "Vincent Rouvreau" +__copyright__ = "Copyright (C) 2016 Inria" +__license__ = "MIT" + +cdef extern from "Witness_complex_interface.h" namespace "Gudhi": + cdef cppclass Witness_complex_interface "Gudhi::witness_complex::Witness_complex_interface": + Witness_complex_interface(vector[vector[pair[size_t, double]]] nearest_landmark_table) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square) + void create_simplex_tree(Simplex_tree_interface_full_featured* simplex_tree, double max_alpha_square, + unsigned limit_dimension) + +# WitnessComplex python interface +cdef class WitnessComplex: + """Constructs (weak) witness complex for a given table of nearest landmarks + with respect to witnesses. + """ + + cdef Witness_complex_interface * thisptr + + # Fake constructor that does nothing but documenting the constructor + def __init__(self, nearest_landmark_table=None): + """WitnessComplex constructor. + + :param nearest_landmark_table: A list of lists of nearest landmarks and their distances. + `nearest_landmark_table[w][k]==(l,d)` means that l is the k-th nearest landmark to + witness w, and d is the (squared) distance between l and w. + :type nearest_landmark_table: list of list of pair of int and float + """ + + # The real cython constructor + def __cinit__(self, nearest_landmark_table=None): + if nearest_landmark_table is not None: + self.thisptr = new Witness_complex_interface(nearest_landmark_table) + + def __dealloc__(self): + if self.thisptr != NULL: + del self.thisptr + + def __is_defined(self): + """Returns true if WitnessComplex pointer is not NULL. + """ + return self.thisptr != NULL + + def create_simplex_tree(self, max_alpha_square = float('inf'), limit_dimension = -1): + """ + :param max_alpha_square: The maximum relaxation parameter. + Default is set to infinity. + :type max_alpha_square: float + :returns: A simplex tree created from the Delaunay Triangulation. + :rtype: SimplexTree + """ + stree = SimplexTree() + cdef intptr_t stree_int_ptr=stree.thisptr + if limit_dimension is not -1: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square, limit_dimension) + else: + self.thisptr.create_simplex_tree(<Simplex_tree_interface_full_featured*>stree_int_ptr, + max_alpha_square) + return stree |