diff --git a/Lab/demo/Lab/Plugins/Mesh_3/Mesh_3_plugin.cpp b/Lab/demo/Lab/Plugins/Mesh_3/Mesh_3_plugin.cpp index d84bd9fc6f6..2aadc761957 100644 --- a/Lab/demo/Lab/Plugins/Mesh_3/Mesh_3_plugin.cpp +++ b/Lab/demo/Lab/Plugins/Mesh_3/Mesh_3_plugin.cpp @@ -750,7 +750,9 @@ void Mesh_3_plugin::mesh_3(const Mesh_type mesh_type, } } - if(mesh_type != Mesh_type::SURFACE_ONLY && !material_ids_valid) + if(mesh_type != Mesh_type::SURFACE_ONLY + && !material_ids_valid + && bounding_sm_item != nullptr) { sm_items.removeAll(make_not_null(bounding_sm_item)); } diff --git a/Lab/demo/Lab/Plugins/Mesh_3/Mesh_function.h b/Lab/demo/Lab/Plugins/Mesh_3/Mesh_function.h index 68c8300f9ee..ec38a369ea0 100644 --- a/Lab/demo/Lab/Plugins/Mesh_3/Mesh_function.h +++ b/Lab/demo/Lab/Plugins/Mesh_3/Mesh_function.h @@ -171,8 +171,10 @@ log() const .arg(detect_connected_components); res << QString("use weights: %1").arg(weights_ptr != nullptr); } - res << QString("use aabb tree: %1").arg(use_sizing_field_with_aabb_tree); - res << QString("manifold: %1").arg(manifold); + if(use_sizing_field_with_aabb_tree) + res << QString("use sizing field with aabb tree: %1").arg(use_sizing_field_with_aabb_tree); + if(manifold) + res << QString("manifold: %1").arg(manifold); return res; } diff --git a/Mesh_3/include/CGAL/Mesh_3/Protect_edges_sizing_field.h b/Mesh_3/include/CGAL/Mesh_3/Protect_edges_sizing_field.h index 4a0e3f131de..e84e7dd933f 100644 --- a/Mesh_3/include/CGAL/Mesh_3/Protect_edges_sizing_field.h +++ b/Mesh_3/include/CGAL/Mesh_3/Protect_edges_sizing_field.h @@ -152,7 +152,7 @@ public: Protect_edges_sizing_field(C3T3& c3t3, const MeshDomain& domain, SizingFunction size=SizingFunction(), - const FT minimal_size = FT(-1), + const FT minimal_size = FT(0), const Distance_Function edge_distance = Distance_Function(), const std::size_t maximal_number_of_vertices = 0, Mesh_error_code* error_code = 0 @@ -273,8 +273,7 @@ private: /// Returns `true` if the edge `(va,vb)` is a not good enough approximation /// of its feature. - bool approx_is_too_large(const Vertex_handle& va, - const Vertex_handle& vb, + bool approx_is_too_large(const Edge& e, const bool is_edge_in_complex) const; /// Returns `true` if the balls of `va` and `vb` intersect. @@ -458,10 +457,12 @@ private: dim = -1 - dim; const FT s = field(p, dim, index); - if(s <= FT(0)) { + if(s < minimal_size_) + { std::stringstream msg; msg.precision(17); - msg << "Error: the field is null at "; + msg << "Error: the field is smaller than minimal size (" + << minimal_size_ << ") at "; if(dim == 0) msg << "corner ("; else msg << "point ("; msg << p << ")"; @@ -493,7 +494,7 @@ private: bool use_minimal_size() const { - return minimal_size_ != FT(-1); + return minimal_size_ != FT(0); } Weight minimal_weight() const { @@ -645,7 +646,10 @@ insert_corners() #endif // Get weight (the ball radius is given by the 'query_size' function) - FT w = CGAL::square(query_size(p, 0, p_index)); + const FT query_weight = CGAL::square(query_size(p, 0, p_index)); + FT w = use_minimal_size() + ? (std::min)(minimal_weight_, query_weight) + : query_weight; #if CGAL_MESH_3_PROTECTION_DEBUG & 1 std::cerr << "Weight from sizing field: " << w << std::endl; @@ -716,7 +720,10 @@ insert_point(const Bare_point& p, const Weight& w, int dim, const Index& index, typename GT::Construct_weighted_point_3 cwp = c3t3_.triangulation().geom_traits().construct_weighted_point_3_object(); - const Weighted_point wp = cwp(p,w*weight_modifier); + const FT wwm = use_minimal_size() + ? (std::max)(w * weight_modifier, minimal_weight()) + : w * weight_modifier; + const Weighted_point wp = cwp(p, wwm); typename Tr::Locate_type lt; int li, lj; @@ -734,7 +741,7 @@ insert_point(const Bare_point& p, const Weight& w, int dim, const Index& index, std::cerr << "SPECIAL "; std::cerr << "protecting ball "; if(v == Vertex_handle()) - std::cerr << cwp(p,w*weight_modifier); + std::cerr << cwp(p, wwm); else std::cerr << disp_vert(v); @@ -860,13 +867,11 @@ smart_insert_point(const Bare_point& p, Weight w, int dim, const Index& index, std::back_inserter(cells_in_conflicts), CGAL::Emptyset_iterator()); - for(typename std::vector::const_iterator - it = cells_in_conflicts.begin(), - end = cells_in_conflicts.end(); it != end; ++it) + for(Cell_handle cit : cells_in_conflicts) { for(int i=0, d=tr.dimension(); i<=d; ++i) { - const Vertex_handle v = (*it)->vertex(i); + const Vertex_handle v = cit->vertex(i); if(c3t3_.triangulation().is_infinite(v)) continue; if(!vertices_in_conflict_zone_set.insert(v).second) @@ -1023,21 +1028,20 @@ insert_balls_on_edges() domain_.get_curves(std::back_inserter(input_features)); // Iterate on edges - for ( typename Input_features::iterator fit = input_features.begin(), - end = input_features.end() ; fit != end ; ++fit ) + for (const Feature_tuple& ft : input_features) { if(forced_stop()) break; - const Curve_index& curve_index = std::get<0>(*fit); + const Curve_index& curve_index = std::get<0>(ft); if ( ! is_treated(curve_index) ) { #if CGAL_MESH_3_PROTECTION_DEBUG & 1 std::cerr << "\n** treat curve #" << curve_index << std::endl; #endif - const Bare_point& p = std::get<1>(*fit).first; - const Bare_point& q = std::get<2>(*fit).first; + const Bare_point& p = std::get<1>(ft).first; + const Bare_point& q = std::get<2>(ft).first; - const Index& p_index = std::get<1>(*fit).second; - const Index& q_index = std::get<2>(*fit).second; + const Index& p_index = std::get<1>(ft).second; + const Index& q_index = std::get<2>(ft).second; Vertex_handle vp,vq; if ( ! domain_.is_loop(curve_index) ) @@ -1178,12 +1182,12 @@ insert_balls(const Vertex_handle& vp, const Weighted_point& vp_wp = c3t3_.triangulation().point(vp); #if ! defined(CGAL_NO_PRECONDITIONS) - if(sp <= 0) { - std::stringstream msg;; + if(sp < minimal_size_) { + std::stringstream msg; msg.precision(17); - msg << "Error: the mesh sizing field is null at point ("; - msg << cp(vp_wp) << ")!"; - CGAL_precondition_msg(sp > 0, msg.str().c_str()); + msg << "Error: the mesh sizing field is smaller than minimal size "; + msg << " at point (" << cp(vp_wp) << ")!"; + CGAL_precondition_msg(sp > minimal_size_, msg.str().c_str()); } #endif // ! CGAL_NO_PRECONDITIONS @@ -1428,7 +1432,7 @@ refine_balls() if( // topology condition non_adjacent_but_intersect(va, vb, is_edge_in_complex) // approximation condition - || (use_distance_field() && approx_is_too_large(va, vb, is_edge_in_complex))) + || (use_distance_field() && approx_is_too_large(*eit, is_edge_in_complex))) { using CGAL::Mesh_3::internal::distance_divisor; @@ -1480,14 +1484,11 @@ refine_balls() new_sizes.clear(); // Update size of balls - for ( typename std::vector >::iterator - it = new_sizes_copy.begin(), - end = new_sizes_copy.end(); - it != end ; ++it ) + for (const std::pair& it : new_sizes_copy) { if(forced_stop()) break; - const Vertex_handle v = it->first; - const FT new_size = it->second; + const Vertex_handle v = it.first; + const FT new_size = it.second; // Set size of the ball to new value if(use_minimal_size() && new_size < minimal_size_) { if(!is_special(v)) { @@ -1557,34 +1558,27 @@ do_balls_intersect(const Vertex_handle& va, const Vertex_handle& vb) const template bool Protect_edges_sizing_field:: -approx_is_too_large(const Vertex_handle& va, const Vertex_handle& vb, const bool is_edge_in_complex) const +approx_is_too_large(const Edge& e, const bool is_edge_in_complex) const { if ( ! is_edge_in_complex ) - { return false; - } - typedef typename Kernel::Point_3 Point_3; - Curve_index curve_index = domain_.curve_index((va->in_dimension() < vb->in_dimension()) ? vb->index() : va->index()); + const Bare_point& pa = e.first->vertex(e.second)->point().point(); + const Bare_point& pb = e.first->vertex(e.third)->point().point(); - const Point_3& pa = va->point().point(); - const Point_3& pb = vb->point().point(); - const Point_3& segment_middle = CGAL::midpoint(pa, pb); - // Obtain the geodesic middle point - FT signed_geodesic_distance = domain_.signed_geodesic_distance(pa, pb, curve_index); - Point_3 geodesic_middle; - if (signed_geodesic_distance >= FT(0)) - { - geodesic_middle = domain_.construct_point_on_curve(pa, curve_index, signed_geodesic_distance / 2); - } - else - { - geodesic_middle = domain_.construct_point_on_curve(pb, curve_index, -signed_geodesic_distance / 2); - } - // Compare distance to the parameter's distance - FT squared_evaluated_distance = CGAL::squared_distance(segment_middle, geodesic_middle); - FT segment_middle_edge_distance = query_distance(segment_middle, 1, curve_index); - return squared_evaluated_distance > CGAL::square(segment_middle_edge_distance); + // Construct the geodesic middle point + const Curve_index curve_index = c3t3_.curve_index(e); + const FT signed_geodesic_distance = domain_.signed_geodesic_distance(pa, pb, curve_index); + const Bare_point geodesic_middle = (signed_geodesic_distance >= FT(0)) + ? domain_.construct_point_on_curve(pa, curve_index, signed_geodesic_distance / 2) + : domain_.construct_point_on_curve(pb, curve_index, -signed_geodesic_distance / 2); + + const Bare_point edge_middle = CGAL::midpoint(pa, pb); + const FT squared_evaluated_distance = CGAL::squared_distance(edge_middle, geodesic_middle); + + // Compare distance to the distance field from criteria + const FT max_distance_to_curve = query_distance(edge_middle, 1, curve_index); + return squared_evaluated_distance > CGAL::square(max_distance_to_curve); } template diff --git a/Mesh_3/include/CGAL/Mesh_domain_with_polyline_features_3.h b/Mesh_3/include/CGAL/Mesh_domain_with_polyline_features_3.h index 61eb0682d83..ddb1ba1109f 100644 --- a/Mesh_3/include/CGAL/Mesh_domain_with_polyline_features_3.h +++ b/Mesh_3/include/CGAL/Mesh_domain_with_polyline_features_3.h @@ -252,8 +252,8 @@ public: : (- negative_distance); } else { return (pit <= qit) - ? curve_segment_length(p, q, CGAL::POSITIVE) - : ( - curve_segment_length(p, q, CGAL::NEGATIVE) ); + ? curve_segment_length(p, q, CGAL::POSITIVE, pit, qit) + : ( - curve_segment_length(p, q, CGAL::NEGATIVE, pit, qit) ); } } diff --git a/Mesh_3/include/CGAL/Mesh_edge_criteria_3.h b/Mesh_3/include/CGAL/Mesh_edge_criteria_3.h index 8d6215cd512..f743eb3f0fa 100644 --- a/Mesh_3/include/CGAL/Mesh_edge_criteria_3.h +++ b/Mesh_3/include/CGAL/Mesh_edge_criteria_3.h @@ -177,7 +177,13 @@ public: /// Returns size of tuple (p,dim,index) FT sizing_field(const Point_3& p, const int dim, const Index& index) const - { return (*p_size_)(p,dim,index); } + { + const FT s = (*p_size_)(p, dim, index); + if (min_length_bound_ == FT(0)) + return s; + else + return (std::max)(s, min_length_bound_); + } FT distance_field(const Point_3& p, const int dim, const Index& index) const {