diff options
Diffstat (limited to 'src/Tangential_complex')
3 files changed, 54 insertions, 9 deletions
diff --git a/src/Tangential_complex/doc/Intro_tangential_complex.h b/src/Tangential_complex/doc/Intro_tangential_complex.h index f4fc8ac7..501f4a8b 100644 --- a/src/Tangential_complex/doc/Intro_tangential_complex.h +++ b/src/Tangential_complex/doc/Intro_tangential_complex.h @@ -35,9 +35,11 @@ namespace tangential_complex { \section tangentialdefinition Definition -A Tangential Delaunay complex is a <a target="_blank" href="https://en.wikipedia.org/wiki/Simplicial_complex">simplicial complex</a> +A Tangential Delaunay complex is a +<a target="_blank" href="https://en.wikipedia.org/wiki/Simplicial_complex">simplicial complex</a> designed to reconstruct a \f$k\f$-dimensional smooth manifold embedded in \f$d\f$-dimensional Euclidean space. -The input is a point sample coming from an unknown manifold, which means that the points lie close to a structure of "small" intrinsic dimension. +The input is a point sample coming from an unknown manifold, which means that the points lie close to a structure of +"small" intrinsic dimension. The running time depends only linearly on the extrinsic dimension \f$ d \f$ and exponentially on the intrinsic dimension \f$ k \f$. @@ -46,17 +48,19 @@ An extensive description of the Tangential complex can be found in \cite tangent \subsection whatisthetc What is a Tangential Complex? Let us start with the description of the Tangential complex of a simple example, with \f$ k=1 \f$ and \f$ d=2 \f$. -The input data is 4 points \f$ P \f$ located on a curve embedded in 2D. +The point set \f$ \mathscr P \f$ is located on a closed curve embedded in 2D. +Only 4 points will be displayed (more are required for PCA) to simplify the figures. \image html "tc_example_01.png" "The input" -For each point \f$ p \f$, estimate its tangent subspace \f$ T_p \f$ (e.g. using PCA). +For each point \f$ P \f$, estimate its tangent subspace \f$ T_P \f$ using PCA. \image html "tc_example_02.png" "The estimated normals" -Let us add the Voronoi diagram of the points in orange. For each point \f$ p \f$, construct its star in the Delaunay triangulation of \f$ P \f$ restricted to \f$ T_p \f$. +Let us add the Voronoi diagram of the points in orange. For each point \f$ P \f$, construct its star in the Delaunay +triangulation of \f$ \mathscr P \f$ restricted to \f$ T_P \f$. \image html "tc_example_03.png" "The Voronoi diagram" The Tangential Delaunay complex is the union of those stars. In practice, neither the ambient Voronoi diagram nor the ambient Delaunay triangulation is computed. -Instead, local \f$ k \f$-dimensional regular triangulations are computed with a limited number of points as we only need the star of each point. -More details can be found in \cite tangentialcomplex2014. +Instead, local \f$ k \f$-dimensional regular triangulations are computed with a limited number of points as we only +need the star of each point. More details can be found in \cite tangentialcomplex2014. \subsection inconsistencies Inconsistencies @@ -65,7 +69,7 @@ An inconsistency occurs when a simplex is not in the star of all its vertices. Let us take the same example. \image html "tc_example_07_before.png" "Before" -Let us slightly move the tangent subspace \f$ T_q \f$ +Let us slightly move the tangent subspace \f$ T_Q \f$ \image html "tc_example_07_after.png" "After" Now, the star of \f$ Q \f$ contains \f$ QP \f$, but the star of \f$ P \f$ does not contain \f$ QP \f$. We have an inconsistency. \image html "tc_example_08.png" "After" diff --git a/src/Tangential_complex/include/gudhi/Tangential_complex.h b/src/Tangential_complex/include/gudhi/Tangential_complex.h index 37cdf1b4..4a78127c 100644 --- a/src/Tangential_complex/include/gudhi/Tangential_complex.h +++ b/src/Tangential_complex/include/gudhi/Tangential_complex.h @@ -322,7 +322,11 @@ class Tangential_complex { for (std::size_t i = 0; i < m_points.size(); ++i) m_are_tangent_spaces_computed[i] = true; } - /// Computes the tangential complex. + /** \brief Computes the tangential complex. + * \exception std::invalid_argument In debug mode, if the computed star dimension is too low. Try to set a bigger + * maximal edge length value with `Tangential_complex::set_max_squared_edge_length` if + * this happens. + */ void compute_tangential_complex() { #ifdef GUDHI_TC_PERFORM_EXTRA_CHECKS std::cerr << red << "WARNING: GUDHI_TC_PERFORM_EXTRA_CHECKS is defined. " @@ -1984,6 +1988,13 @@ class Tangential_complex { return os; } + /** \brief Sets the maximal possible squared edge length for the edges in the triangulations. + * + * @param[in] max_squared_edge_length Maximal possible squared edge length. + * + * If the maximal edge length value is too low `Tangential_complex::compute_tangential_complex` will throw an + * exception in debug mode. + */ void set_max_squared_edge_length(FT max_squared_edge_length) { m_max_squared_edge_length = max_squared_edge_length; } private: diff --git a/src/Tangential_complex/test/test_tangential_complex.cpp b/src/Tangential_complex/test/test_tangential_complex.cpp index 4e2d4f65..103b8b30 100644 --- a/src/Tangential_complex/test/test_tangential_complex.cpp +++ b/src/Tangential_complex/test/test_tangential_complex.cpp @@ -126,3 +126,33 @@ BOOST_AUTO_TEST_CASE(test_mini_tangential) { BOOST_CHECK(stree.num_vertices() == 4); BOOST_CHECK(stree.num_simplices() == 6); } + +#ifdef GUDHI_DEBUG +BOOST_AUTO_TEST_CASE(test_basic_example_throw) { + typedef CGAL::Epick_d<CGAL::Dynamic_dimension_tag> Kernel; + typedef Kernel::FT FT; + typedef Kernel::Point_d Point; + typedef Kernel::Vector_d Vector; + typedef tc::Tangential_complex<Kernel, CGAL::Dynamic_dimension_tag,CGAL::Parallel_tag> TC; + + const int INTRINSIC_DIM = 2; + const int AMBIENT_DIM = 3; + const int NUM_POINTS = 1000; + + Kernel k; + + // Generate points on a 2-sphere + CGAL::Random_points_on_sphere_d<Point> generator(AMBIENT_DIM, 3.); + std::vector<Point> points; + points.reserve(NUM_POINTS); + for (int i = 0; i < NUM_POINTS; ++i) + points.push_back(*generator++); + + // Compute the TC + TC tc(points, INTRINSIC_DIM, k); + tc.set_max_squared_edge_length(0.01); + std::cout << "test_basic_example_throw - set_max_squared_edge_length(0.01) to make GUDHI_CHECK fail" << std::endl; + BOOST_CHECK_THROW(tc.compute_tangential_complex(), std::invalid_argument); + +} +#endif |