From c41a0e38c2035a35cb911cb356b8e7e966bada39 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ivan=20Pa=C4=91en?= <49401914+ipadjen@users.noreply.github.com> Date: Fri, 6 Oct 2023 19:00:11 -0600 Subject: [PATCH] Merge two tangential relaxation functions into one --- .../tangential_relaxation.h | 326 +----------------- .../remeshing_quality_test.cpp | 2 +- 2 files changed, 11 insertions(+), 317 deletions(-) diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/tangential_relaxation.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/tangential_relaxation.h index e41adb9ca89..83fdd5ff37f 100644 --- a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/tangential_relaxation.h +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/tangential_relaxation.h @@ -116,6 +116,13 @@ struct Allow_all_moves{ * \cgalParamDefault{If not provided, all moves are allowed.} * \cgalParamNEnd * +* \cgalParamNBegin{sizing_function} +* \cgalParamDescription{An object containing sizing field for individual vertices. +* Used to derive smoothing weights.} +* \cgalParamType{A model of `PMPSizingField`.} +* \cgalParamDefault{If not provided, smoothing weights are the same for all vertices.} +* \cgalParamNEnd +* * \cgalNamedParamsEnd * * \todo check if it should really be a triangle mesh or if a polygon mesh is fine @@ -232,12 +239,6 @@ void tangential_relaxation(const VertexRange& vertices, auto gt_barycenter = gt.construct_barycenter_3_object(); auto gt_project = gt.construct_projected_point_3_object(); -// if constexpr (is_default_parameter::value) -// { -// auto gt_area = gt.compute_area_3_object(); -// auto gt_centroid = gt.construct_centroid_3_object(); -// } - // at each vertex, compute vertex normal std::unordered_map vnormals; compute_vertex_normals(tm, boost::make_assoc_property_map(vnormals), np); @@ -287,7 +288,9 @@ void tangential_relaxation(const VertexRange& vertices, const double tri_area = gt_area(get(vpm, v), get(vpm, v1), get(vpm, v2)); const double face_weight = tri_area - / (1. / 3. * (sizing.get_sizing(v) + sizing.get_sizing(v1) + sizing.get_sizing(v2))); + / (1. / 3. * (sizing.get_sizing(v) + + sizing.get_sizing(v1) + + sizing.get_sizing(v2))); weight += face_weight; const Point_3 centroid = gt_centroid(get(vpm, v), get(vpm, v1), get(vpm, v2)); @@ -341,305 +344,6 @@ void tangential_relaxation(const VertexRange& vertices, new_locations.emplace_back(v, qv + (nv * Vector_3(qv, pv)) * nv); } - // perform moves - for(const VP_pair& vp : new_locations) - { - const Point_3 initial_pos = get(vpm, vp.first); // make a copy on purpose - const Vector_3 move(initial_pos, vp.second); - - put(vpm, vp.first, vp.second); - - //check that no inversion happened - double frac = 1.; - while (frac > 0.03 //5 attempts maximum - && ( !check_normals(vp.first) - || !shall_move(vp.first, initial_pos, get(vpm, vp.first)))) //if a face has been inverted - { - frac = 0.5 * frac; - put(vpm, vp.first, initial_pos + frac * move);//shorten the move by 2 - } - if (frac <= 0.02) - put(vpm, vp.first, initial_pos);//cancel move - } - }//end for loop (nit == nb_iterations) - -#ifdef CGAL_PMP_TANGENTIAL_RELAXATION_VERBOSE - std::cout << "\rTangential relaxation : " - << nb_iterations << " iterations done." << std::endl; -#endif -} - -/*! -* \ingroup PMP_meshing_grp -* applies an iterative area-based tangential smoothing to the given range of vertices based on the -* underlying sizing field. -* Each vertex `v` is relocated to its weighted centroid, where weights depend on the area of the -* adjacent triangle and its averaged vertex sizing values. -* The relocation vector is projected back to the tangent plane to the surface at `v`, iteratively. -* The connectivity remains unchanged. -* -* @tparam TriangleMesh model of `FaceGraph` and `VertexListGraph`. -* The descriptor types `boost::graph_traits::%face_descriptor` -* and `boost::graph_traits::%halfedge_descriptor` must be -* models of `Hashable`. -* @tparam VertexRange range of `boost::graph_traits::%vertex_descriptor`, -* model of `Range`. Its iterator type is `ForwardIterator`. -* @tparam SizingFunction model of `PMPSizingField` -* @tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters" -* -* @param vertices the range of vertices which will be relocated by relaxation -* @param tm the triangle mesh to which `vertices` belong -* @param sizing a map containing sizing field for individual vertices. -* Used to derive smoothing weights. -* @param np an optional sequence of \ref bgl_namedparameters "Named Parameters" among the ones listed below -* -* \cgalNamedParamsBegin -* \cgalParamNBegin{vertex_point_map} -* \cgalParamDescription{a property map associating points to the vertices of `tm`} -* \cgalParamType{a class model of `ReadWritePropertyMap` with `boost::graph_traits::%vertex_descriptor` -* as key type and `%Point_3` as value type} -* \cgalParamDefault{`boost::get(CGAL::vertex_point, tm)`} -* \cgalParamExtra{If this parameter is omitted, an internal property map for `CGAL::vertex_point_t` -* must be available in `TriangleMesh`.} -* \cgalParamNEnd -* -* \cgalParamNBegin{geom_traits} -* \cgalParamDescription{an instance of a geometric traits class} -* \cgalParamType{a class model of `Kernel`} -* \cgalParamDefault{a \cgal Kernel deduced from the `Point_3` type, using `CGAL::Kernel_traits`} -* \cgalParamExtra{The geometric traits class must be compatible with the vertex `Point_3` type.} -* \cgalParamExtra{Exact constructions kernels are not supported by this function.} -* \cgalParamNEnd -* -* \cgalParamNBegin{number_of_iterations} -* \cgalParamDescription{the number of smoothing iterations} -* \cgalParamType{unsigned int} -* \cgalParamDefault{`1`} -* \cgalParamNEnd -* -* \cgalParamNBegin{edge_is_constrained_map} -* \cgalParamDescription{a property map containing the constrained-or-not status of each edge of `tm`. -* The endpoints of a constrained edge cannot be moved by relaxation, unless `relax_constraints` is `true`. -* The endpoints of a constrained polyline can never be moved by relaxation. -* \cgalParamType{a class model of `ReadWritePropertyMap` with `boost::graph_traits::%edge_descriptor` -* as key type and `bool` as value type. It must be default constructible.} -* \cgalParamDefault{a default property map where no edges are constrained} -* \cgalParamExtra{Boundary edges are always considered as constrained edges.} -* \cgalParamNEnd -* -* \cgalParamNBegin{vertex_is_constrained_map} -* \cgalParamDescription{a property map containing the constrained-or-not status of each vertex of `tm`. -* A constrained vertex cannot be modified during relaxation.} -* \cgalParamType{a class model of `ReadWritePropertyMap` with `boost::graph_traits::%vertex_descriptor` -* as key type and `bool` as value type. It must be default constructible.} -* \cgalParamDefault{a default property map where no vertices are constrained} -* \cgalParamNEnd -* -* \cgalParamNBegin{relax_constraints} -* \cgalParamDescription{If `true`, the end vertices of the edges set as constrained -* in `edge_is_constrained_map` and boundary edges move along the -* constrained polylines they belong to.} -* \cgalParamType{Boolean} -* \cgalParamDefault{`false`} -* \cgalParamNEnd -* -* \cgalParamNBegin{allow_move_functor} -* \cgalParamDescription{A function object used to determinate if a vertex move should be allowed or not} -* \cgalParamType{Unary functor that provides `bool operator()(vertex_descriptor v, Point_3 src, Point_3 tgt)` returning `true` -* if the vertex `v` can be moved from `src` to `tgt`; `Point_3` being the value type of the vertex point map } -* \cgalParamDefault{If not provided, all moves are allowed.} -* \cgalParamNEnd -* -* \cgalNamedParamsEnd -* -* \todo check if it should really be a triangle mesh or if a polygon mesh is fine -*/ -template -void tangential_relaxation_with_sizing(const VertexRange& vertices, - TriangleMesh& tm, - const SizingFunction& sizing, - const CGAL_NP_CLASS& np = parameters::default_values()) -{ - typedef typename boost::graph_traits::vertex_descriptor vertex_descriptor; - typedef typename boost::graph_traits::halfedge_descriptor halfedge_descriptor; - typedef typename boost::graph_traits::edge_descriptor edge_descriptor; - - using parameters::get_parameter; - using parameters::choose_parameter; - - typedef typename GetGeomTraits::type GT; - GT gt = choose_parameter(get_parameter(np, internal_np::geom_traits), GT()); - - typedef typename GetVertexPointMap::type VPMap; - VPMap vpm = choose_parameter(get_parameter(np, internal_np::vertex_point), - get_property_map(vertex_point, tm)); - - typedef Static_boolean_property_map Default_ECM; - typedef typename internal_np::Lookup_named_param_def < - internal_np::edge_is_constrained_t, - CGAL_NP_CLASS, - Static_boolean_property_map // default (no constraint) - > ::type ECM; - ECM ecm = choose_parameter(get_parameter(np, internal_np::edge_is_constrained), - Default_ECM()); - - typedef typename internal_np::Lookup_named_param_def < - internal_np::vertex_is_constrained_t, - CGAL_NP_CLASS, - Static_boolean_property_map // default (no constraint) - > ::type VCM; - VCM vcm = choose_parameter(get_parameter(np, internal_np::vertex_is_constrained), - Static_boolean_property_map()); - - const bool relax_constraints = choose_parameter(get_parameter(np, internal_np::relax_constraints), false); - const unsigned int nb_iterations = choose_parameter(get_parameter(np, internal_np::number_of_iterations), 1); - - typedef typename GT::Vector_3 Vector_3; - typedef typename GT::Point_3 Point_3; - - auto check_normals = [&](vertex_descriptor v) - { - bool first_run = true; - Vector_3 prev = NULL_VECTOR, first = NULL_VECTOR; - halfedge_descriptor first_h = boost::graph_traits::null_halfedge(); - for (halfedge_descriptor hd : CGAL::halfedges_around_target(v, tm)) - { - if (is_border(hd, tm)) continue; - - Vector_3 n = compute_face_normal(face(hd, tm), tm, np); - if (n == CGAL::NULL_VECTOR) //for degenerate faces - continue; - - if (first_run) - { - first_run = false; - first = n; - first_h = hd; - } - else - { - if (!get(ecm, edge(hd, tm))) - if (to_double(n * prev) <= 0) - return false; - } - prev = n; - } - - if (first_run) - return true; //vertex incident only to degenerate faces - - if (!get(ecm, edge(first_h, tm))) - if (to_double(first * prev) <= 0) - return false; - - return true; - }; - - typedef typename internal_np::Lookup_named_param_def < - internal_np::allow_move_functor_t, - CGAL_NP_CLASS, - internal::Allow_all_moves// default - > ::type Shall_move; - Shall_move shall_move = choose_parameter(get_parameter(np, internal_np::allow_move_functor), - internal::Allow_all_moves()); - - for (unsigned int nit = 0; nit < nb_iterations; ++nit) - { -#ifdef CGAL_PMP_TANGENTIAL_RELAXATION_VERBOSE - std::cout << "\r\t(Tangential relaxation iteration " << (nit + 1) << " / "; - std::cout << nb_iterations << ") "; - std::cout.flush(); -#endif - - typedef std::tuple VNP; - std::vector< VNP > barycenters; - auto gt_barycenter = gt.construct_barycenter_3_object(); - auto gt_centroid = gt.construct_centroid_3_object(); - auto gt_area = gt.compute_area_3_object(); - - // at each vertex, compute vertex normal - std::unordered_map vnormals; - compute_vertex_normals(tm, boost::make_assoc_property_map(vnormals), np); - - // at each vertex, compute centroids of neighbouring faces weighted by - // area and sizing field - for(vertex_descriptor v : vertices) - { - if (get(vcm, v) || CGAL::internal::is_isolated(v, tm)) - continue; - - // collect hedges to detect if we have to handle boundary cases - std::vector interior_hedges, border_halfedges; - for(halfedge_descriptor h : halfedges_around_target(v, tm)) - { - if (is_border_edge(h, tm) || get(ecm, edge(h, tm))) - border_halfedges.push_back(h); - else - interior_hedges.push_back(h); - } - - if (border_halfedges.empty()) - { - const Vector_3& vn = vnormals.at(v); - Vector_3 move = CGAL::NULL_VECTOR; - double weight = 0; - for(halfedge_descriptor h :interior_hedges) - { - // calculate weight - // need v, v1 and v2 - const vertex_descriptor v1 = target(next(h, tm), tm); - const vertex_descriptor v2 = source(h, tm); - - const double tri_area = gt_area(get(vpm, v), get(vpm, v1), get(vpm, v2)); - const double face_weight = tri_area - / (1. / 3. * (sizing.get_sizing(v) + sizing.get_sizing(v1) + sizing.get_sizing(v2))); - weight += face_weight; - - const Point_3 centroid = gt_centroid(get(vpm, v), get(vpm, v1), get(vpm, v2)); - move = move + Vector_3(get(vpm, v), centroid) * face_weight; - } - move = move / weight; //todo ip: what if weight ends up being close to 0? - - barycenters.emplace_back(v, vn, get(vpm, v) + move); - } - else - { - if (!relax_constraints) continue; - Vector_3 vn(NULL_VECTOR); - - if (border_halfedges.size() == 2)// corners are constrained - { - vertex_descriptor ph0 = source(border_halfedges[0], tm); - vertex_descriptor ph1 = source(border_halfedges[1], tm); - double dot = to_double(Vector_3(get(vpm, v), get(vpm, ph0)) - * Vector_3(get(vpm, v), get(vpm, ph1))); - // \todo shouldn't it be an input parameter? - //check squared cosine is < 0.25 (~120 degrees) - if (0.25 < dot*dot / ( squared_distance(get(vpm,ph0), get(vpm, v)) * - squared_distance(get(vpm,ph1), get(vpm, v))) ) - barycenters.emplace_back(v, vn, - gt_barycenter(get(vpm, ph0), 0.25, get(vpm, ph1), 0.25, get(vpm, v), 0.5)); - } - } - } - - // compute moves - typedef std::pair VP_pair; - std::vector< std::pair > new_locations; - new_locations.reserve(barycenters.size()); - for(const VNP& vnp : barycenters) - { - vertex_descriptor v = std::get<0>(vnp); - const Point_3& pv = get(vpm, v); - const Vector_3& nv = std::get<1>(vnp); - const Point_3& qv = std::get<2>(vnp); //barycenter at v - - new_locations.emplace_back(v, qv + (nv * Vector_3(qv, pv)) * nv); - } - // perform moves for(const VP_pair& vp : new_locations) { @@ -680,16 +384,6 @@ void tangential_relaxation(TriangleMesh& tm, tangential_relaxation(vertices(tm), tm, np); } -template -void tangential_relaxation_with_sizing(TriangleMesh& tm, - const SizingFunction& sizing, - const CGAL_NP_CLASS& np = parameters::default_values()) -{ - tangential_relaxation(vertices(tm), tm, sizing, np); -} - } } // CGAL::Polygon_mesh_processing #endif //CGAL_POLYGON_MESH_PROCESSING_TANGENTIAL_RELAXATION_H diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_quality_test.cpp b/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_quality_test.cpp index a6216663596..93b3c4cad05 100644 --- a/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_quality_test.cpp +++ b/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_quality_test.cpp @@ -107,4 +107,4 @@ int main(int argc, char* argv[]) << " deg" << std::endl; return 0; -} \ No newline at end of file +}