From 453ec1571ceb38dbbc283d3cfcea9de22a1db652 Mon Sep 17 00:00:00 2001 From: Martin Skrodzki Date: Wed, 19 Jun 2019 18:47:41 +0200 Subject: [PATCH] Implement vertex (TM1) to triangle (TM2) distance in traversal processing. --- .../CGAL/Polygon_mesh_processing/distance.h | 2 +- ...traversal_traits_with_Hausdorff_distance.h | 62 ++++++++++++++++--- 2 files changed, 54 insertions(+), 10 deletions(-) diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/distance.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/distance.h index 31ff50ccff7..4c269ce9823 100644 --- a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/distance.h +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/distance.h @@ -922,7 +922,7 @@ double bounded_error_Hausdorff_impl( std::pair hint = tm2_tree.any_reference_point_and_id(); // Build traversal traits for tm1_tree and tm2_tree - Hausdorff_primitive_traits_tm1 traversal_traits_tm1( tm1_tree.traits(), tm2_tree ); + Hausdorff_primitive_traits_tm1 traversal_traits_tm1( tm1_tree.traits(), tm2_tree, tm1, tm2, vpm1, vpm2 ); tm1_tree.traversal( Point_3(0,0,0), traversal_traits_tm1 ); // For each vertex in tm1, store the distance to the closest triangle of tm2 diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/AABB_traversal_traits_with_Hausdorff_distance.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/AABB_traversal_traits_with_Hausdorff_distance.h index 6be21480eb3..86cf640af02 100644 --- a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/AABB_traversal_traits_with_Hausdorff_distance.h +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/AABB_traversal_traits_with_Hausdorff_distance.h @@ -32,15 +32,18 @@ namespace CGAL { /** * @class Hausdorff_primitive_traits_tm2 */ - template + template class Hausdorff_primitive_traits_tm2 { typedef typename AABBTraits::Primitive Primitive; typedef ::CGAL::AABB_node Node; + typedef typename boost::graph_traits::vertex_descriptor vertex_descriptor; + typedef typename boost::graph_traits::halfedge_descriptor halfedge_descriptor; + typename Kernel::Construct_projected_point_3 project_point; public: - Hausdorff_primitive_traits_tm2(const AABBTraits& traits) - : m_traits(traits) { + Hausdorff_primitive_traits_tm2(const AABBTraits& traits, const TriangleMesh& tm1, const TriangleMesh& tm2, const VPM1& vpm1, const VPM2& vpm2) + : m_traits(traits), m_tm1(tm1), m_tm2(tm2), m_vpm1(vpm1), m_vpm2(vpm2) { // Initialize the global bounds with 0., they will only grow. h_local_upper = std::numeric_limits::infinity(); h_local_lower_0 = std::numeric_limits::infinity(); @@ -58,12 +61,39 @@ namespace CGAL { // Have reached a single triangle std::cout << "Reached Triangle in TM2:" << primitive.id() << '\n'; - double distance = std::numeric_limits::infinity(); /* - / TODO Determine the distance accroding to + / Determine the distance accroding to / min_{b \in primitive} ( max_{vertex in query} ( d(vertex, b))) + / + / Here, we only have one triangle in B, i.e. tm2. Thus, it suffices to + / compute the distance of the vertices of the query triangles to the + / primitive triangle and use the maximum of the obtained distances. */ + // The query object is a triangle from TM1, get its vertices + halfedge_descriptor hd = halfedge(query.id(), m_tm1); + vertex_descriptor v0 = source(hd,m_tm1); + vertex_descriptor v1 = target(hd,m_tm1); + vertex_descriptor v2 = target(next(hd, m_tm1), m_tm1); + + // Compute distances of the vertices to the primitive triangle in TM2 + Triangle_from_face_descriptor_map face_to_triangle_map(&m_tm2, m_vpm2); + double v0_dist = squared_distance( + project_point(get(face_to_triangle_map, primitive.id()), get(m_vpm1, v0)), + get(m_vpm1, v0) + ); + double v1_dist = squared_distance( + project_point(get(face_to_triangle_map, primitive.id()), get(m_vpm1, v1)), + get(m_vpm1, v1) + ); + double v2_dist = squared_distance( + project_point(get(face_to_triangle_map, primitive.id()), get(m_vpm1, v2)), + get(m_vpm1, v2) + ); + + // Get the distance as maximizers over all vertices + double distance = approximate_sqrt(std::max(std::max(v0_dist,v1_dist),v2_dist)); + // Update local upper bound if ( distance < h_local_upper ) h_local_upper = distance; @@ -99,6 +129,12 @@ namespace CGAL { private: const AABBTraits& m_traits; + // The two triangle meshes + const TriangleMesh& m_tm1; + const TriangleMesh& m_tm2; + // Their respective vertex-point Maps + const VPM1& m_vpm1; + const VPM2& m_vpm2; // Local Hausdorff bounds for the query triangle double h_local_upper; double h_local_lower; @@ -111,7 +147,7 @@ namespace CGAL { /** * @class Hausdorff_primitive_traits_tm1 */ - template + template class Hausdorff_primitive_traits_tm1 { typedef AABB_face_graph_triangle_primitive TM2_primitive; @@ -122,10 +158,12 @@ namespace CGAL { typedef typename Kernel::Point_3 Point_3; typedef typename Kernel::Triangle_3 Triangle_3; typedef AABB_tree< AABB_traits > TM2_tree; + typedef typename boost::graph_traits::vertex_descriptor vertex_descriptor; + typedef typename boost::graph_traits::face_descriptor face_descriptor; public: - Hausdorff_primitive_traits_tm1(const AABBTraits& traits, const TM2_tree& tree) - : m_traits(traits), tm2_tree(tree) { + Hausdorff_primitive_traits_tm1(const AABBTraits& traits, const TM2_tree& tree, const TriangleMesh& tm1, const TriangleMesh& tm2 , const VPM1& vpm1, const VPM2& vpm2 ) + : m_traits(traits), tm2_tree(tree), m_tm1(tm1), m_tm2(tm2), m_vpm1(vpm1), m_vpm2(vpm2) { // Initialize the global bounds with 0., they will only grow. h_lower = 0.; h_upper = 0.; @@ -142,7 +180,7 @@ namespace CGAL { std::cout << "Reached Triangle in TM1: " << primitive.id() << '\n'; // Call Culling on B with the single triangle found. - Hausdorff_primitive_traits_tm2 traversal_traits_tm2( tm2_tree.traits() ); + Hausdorff_primitive_traits_tm2 traversal_traits_tm2( tm2_tree.traits(), m_tm1, m_tm2, m_vpm1, m_vpm2 ); tm2_tree.traversal(primitive, traversal_traits_tm2); // Update global Hausdorff bounds according to the obtained local bounds @@ -187,6 +225,12 @@ namespace CGAL { private: const AABBTraits& m_traits; + // The two triangle meshes + const TriangleMesh& m_tm1; + const TriangleMesh& m_tm2; + // Their vertex-point-maps + const VPM1 m_vpm1; + const VPM2 m_vpm2; // AABB tree for the second triangle meshes const TM2_tree& tm2_tree; // Global Hausdorff bounds to be taken track of during the traversal