From 9b784512d2086ba2f1d80b72215d44cec1c49b4c Mon Sep 17 00:00:00 2001 From: Jacques-Olivier Lachaud Date: Tue, 8 Oct 2024 18:04:34 +0200 Subject: [PATCH] Add covering for triangles --- src/DGtal/geometry/volumes/DigitalConvexity.h | 3 +- .../geometry/volumes/testDigitalConvexity.cpp | 70 +++++++++++++++++++ 2 files changed, 71 insertions(+), 2 deletions(-) diff --git a/src/DGtal/geometry/volumes/DigitalConvexity.h b/src/DGtal/geometry/volumes/DigitalConvexity.h index 1484f6547c..5a57b9517a 100644 --- a/src/DGtal/geometry/volumes/DigitalConvexity.h +++ b/src/DGtal/geometry/volumes/DigitalConvexity.h @@ -468,8 +468,7 @@ namespace DGtal /// i.e. all the cells intersected by the triangle belong to \a /// cells. /// - /// @note This method is supposed to be faster than the others, - /// but is limited to 3D triangles. + /// @note This method is limited to 3D triangles. bool isFullyCovered( const Point& a, const Point& b, const Point& c, const LatticeSet& cells ) const; diff --git a/tests/geometry/volumes/testDigitalConvexity.cpp b/tests/geometry/volumes/testDigitalConvexity.cpp index 084e695d63..43297325fd 100644 --- a/tests/geometry/volumes/testDigitalConvexity.cpp +++ b/tests/geometry/volumes/testDigitalConvexity.cpp @@ -1064,3 +1064,73 @@ SCENARIO( "DigitalConvexity< Z3 > full covering of triangles", "[full_cover][3d] REQUIRE( P.size() == 7 ); } } + + +SCENARIO( "DigitalConvexity< Z3 > full subconvexity and full covering of triangles", "[subconvexity][3d][full_cover]" ) +{ + typedef KhalimskySpaceND<3,int> KSpace; + typedef KSpace::Point Point; + typedef KSpace::Vector Vector; + typedef KSpace::Space Space; + typedef HyperRectDomain< Space > Domain; + typedef DigitalConvexity< KSpace > DConvexity; + + Domain domain( Point( -20, -20, -20 ), Point( 20, 20, 20 ) ); + DConvexity dconv ( Point( -21, -21, -21 ), Point( 21, 21, 21 ) ); + + WHEN( "Computing many tetrahedra" ) { + const unsigned int nb = 50; + unsigned int nb_total = 0; + unsigned int nb_ok_tri = 0; + unsigned int nb_subconvex = 0; + unsigned int nb_covered = 0; + for ( unsigned int l = 0; l < nb; ++l ) + { + const Point a { (rand() % 10 - 10), (rand() % 10 - 10), (rand() % 10 - 10) }; + const Point b { (rand() % 20 ), (rand() % 20 - 10), (rand() % 20 - 10) }; + const Point c { (rand() % 20 - 10), (rand() % 20 ), (rand() % 20 - 10) }; + const Point d { (rand() % 20 - 10), (rand() % 20 - 10), (rand() % 20 ) }; + const std::vector pts = { a, b, c, d }; + const bool fulldim = dconv.isSimplexFullDimensional(pts.cbegin(), pts.cend()); + if ( ! fulldim ) continue; + auto simplex = dconv.makeSimplex( pts.cbegin(), pts.cend() ); + auto cover = dconv.makeCellCover( simplex, 0, 3 ); + auto ls = dconv.StarCvxH( pts ); + { + for ( unsigned int i = 0; i < 100; i++ ) + { + const Point p { (rand() % 10 - 5), (rand() % 10 - 5), (rand() % 10 - 5) }; + const Point q { (rand() % 10 - 5), (rand() % 10 - 5), (rand() % 10 - 5) }; + const Point r { (rand() % 10 - 5), (rand() % 10 - 5), (rand() % 10 - 5) }; + const Vector n = ( q - p ).crossProduct( r - p ); + if ( n == Vector::zero ) continue; + auto tri1 = dconv.makeSimplex( { p, q, r } ); + bool ok1 = dconv.isFullySubconvex( p, q, r, ls ); + bool ok2 = dconv.isFullyCovered( p, q, r, ls ); + nb_subconvex += ok1 ? 1 : 0; + nb_covered += ok2 ? 1 : 0; + bool ok = ok1 == ok2; + if ( ! ok ) + { + std::cout << "***** FULL SUBCONVEXITY ERROR ON TRIANGLE ****" << std::endl; + std::cout << "splx v =" << a << b << c << d << std::endl; + std::cout << "simplex=" << simplex << std::endl; + std::cout << "tri v =" << p << q << r << std::endl; + std::cout << "tri1=" << tri1 << std::endl; + std::cout << "3 points are fully subconvex: " << ( ok1 ? "YES" : "NO" ) << std::endl; + std::cout << "3 points are fully covered by star: " << ( ok2 ? "YES" : "NO" ) << std::endl; + } + nb_ok_tri += ( ok ) ? 1 : 0; + nb_total += 1; + } + } + } + THEN( "The number of subconvex and covered triangles should be equal." ) { + REQUIRE( nb_subconvex == nb_covered ); + } + THEN( "Full subconvexity and covering should agree on every subset." ) { + REQUIRE( nb_ok_tri == nb_total ); + } + } + +}