diff --git a/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h b/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h index 71119655024..077bba235e2 100644 --- a/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h +++ b/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h @@ -646,25 +646,27 @@ public: // ---------------------------------- using parameters::choose_parameter; using parameters::get_parameter; + namespace p = parameters; - auto mesh_vp_map = choose_parameter(get_parameter(np, internal_np::vertex_point), get(CGAL::vertex_point, mesh)); + auto mesh_vp_map = choose_parameter(get_parameter(np, internal_np::vertex_point), + get_const_property_map(vertex_point, mesh)); using vertex_descriptor = typename boost::graph_traits::vertex_descriptor; std::vector>> patch_edges; - auto v_selected_map = get(CGAL::dynamic_vertex_property_t{}, mesh); + auto tr_vertex_pmap = get(CGAL::dynamic_vertex_property_t(), mesh); - int number_of_patches = 0; - constexpr bool has_plc_face_id = !parameters::is_default_parameter::value; + constexpr bool has_plc_face_id = !p::is_default_parameter::value; if constexpr (has_plc_face_id) { - auto mesh_plc_face_id = parameters::get_parameter(np, internal_np::plc_face_id); + auto v_selected_map = get(CGAL::dynamic_vertex_property_t{}, mesh); + auto mesh_plc_face_id = get_parameter(np, internal_np::plc_face_id); using Patch_id_type = CGAL::cpp20::remove_cvref_t; Patch_id_type max_patch_id{0}; for(auto f : faces(mesh)) { max_patch_id = (std::max)(max_patch_id, get(mesh_plc_face_id, f)); } - number_of_patches = static_cast(max_patch_id + 1); + auto number_of_patches = max_patch_id + 1; patch_edges.resize(number_of_patches); for(auto h : halfedges(mesh)) { if(is_border(h, mesh)) @@ -680,38 +682,20 @@ public: put(v_selected_map, vb, true); } } - } else { - number_of_patches = num_faces(mesh); - } - using Vertex_handle = typename Triangulation::Vertex_handle; - using Cell_handle = typename Triangulation::Cell_handle; - auto tr_vertex_pmap = get(CGAL::dynamic_vertex_property_t(), mesh); - Cell_handle hint_ch{}; - for(auto v : vertices(mesh)) { - if constexpr(has_plc_face_id) { - if(false == get(v_selected_map, v)) continue; - } - auto p = get(mesh_vp_map, v); - auto vh = cdt_impl.insert(p, hint_ch, false); - vh->ccdt_3_data().set_vertex_type(CDT_3_vertex_type::CORNER); - hint_ch = vh->cell(); - put(tr_vertex_pmap, v, vh); - } - - cdt_impl.add_bbox_points_if_not_dimension_3(); - - if constexpr(has_plc_face_id) { + cdt_impl.insert_mesh_vertices(mesh, mesh_vp_map, tr_vertex_pmap, p::vertex_is_constrained_map(v_selected_map)); for(auto&& edges : patch_edges) { cdt_impl.insert_constrained_face_defined_by_its_border_edges(edges, mesh_vp_map, tr_vertex_pmap); } } else { + cdt_impl.insert_mesh_vertices(mesh, mesh_vp_map, tr_vertex_pmap); for(auto f : faces(mesh)) { auto face_vertices = vertices_around_face(halfedge(f, mesh), mesh); auto range_of_vertices = CGAL::make_transform_range_from_property_map(face_vertices, tr_vertex_pmap); cdt_impl.insert_constrained_face(range_of_vertices, false); } } + cdt_impl.restore_Delaunay(); cdt_impl.restore_constrained_Delaunay(); // std::cerr << cdt_3_format("cdt_impl: {} vertices, {} cells\n", cdt_impl.number_of_vertices(), @@ -1353,6 +1337,49 @@ public: Conforming_Dt::restore_Delaunay(insert_in_conflict_visitor); } + /** + * Insert vertices from a polygon mesh into the triangulation. + * + * @tparam PolygonMesh A model of FaceGraph. + * @tparam VertexPointMap A readable property map from vertex descriptor to Point_3. + * @tparam VertexHandleMap A writable property map from vertex descriptor to Vertex_handle. + * @tparam NamedParameters A sequence of named parameters. + * + * @param mesh The polygon mesh. + * @param mesh_vp_map Property map to access vertex coordinates. + * @param tr_vertex_pmap Property map to store the mapping from mesh vertices to triangulation vertices. + * @param np Optional named parameters: + * - `vertex_is_constrained_map`: A readable property map from vertex descriptor to bool. + * If provided, only vertices where the map returns true will be inserted. + */ + template + void insert_mesh_vertices(const PolygonMesh& mesh, + VertexPointMap mesh_vp_map, + VertexHandleMap tr_vertex_pmap, + const NamedParameters& np = parameters::default_values()) + { + using parameters::get_parameter; + using parameters::is_default_parameter; + + constexpr bool has_v_filter = !is_default_parameter::value; + auto v_filter_map = get_parameter(np, internal_np::vertex_is_constrained); + + Cell_handle hint_ch{}; + for(auto v : vertices(mesh)) { + if constexpr(has_v_filter) { + if(false == get(v_filter_map, v)) continue; + } + auto p = get(mesh_vp_map, v); + auto vh = this->insert(p, hint_ch, false); + vh->ccdt_3_data().set_vertex_type(CDT_3_vertex_type::CORNER); + hint_ch = vh->cell(); + put(tr_vertex_pmap, v, vh); + } + + add_bbox_points_if_not_dimension_3(); + } + /** * Insert a constrained face defined by its border edges. * diff --git a/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp b/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp index 754bff240e5..156d1d51aa0 100644 --- a/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp +++ b/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp @@ -9,13 +9,14 @@ // #define CGAL_CDT_2_DEBUG_INTERSECTIONS 1 #define NO_TRY_CATCH 1 // #define CGAL_DEBUG_CDT_3 1 -#include -#include #include #include -#include #include +#include +#include #include +#include +#include #include #include @@ -448,7 +449,7 @@ Min_distance_result compute_minimum_vertex_distance(const CDT& cdt) { } auto [min_sq_distance, min_edge] = min_p; #endif - auto min_distance = CGAL::approximate_sqrt(min_sq_distance); + auto min_distance = CGAL::to_double(CGAL::approximate_sqrt(min_sq_distance)); auto vertices_of_min_edge = cdt.vertices(min_edge); return {min_distance, vertices_of_min_edge}; @@ -707,13 +708,15 @@ Borders_of_patches maybe_merge_facets( Mesh& mesh, const CDT_options& options, MeshPropertyMaps pmaps, - double bbox_max_span) { - + double bbox_max_span) +{ int number_of_patches = 0; Borders_of_patches patch_edges; - if(options.merge_facets) { - CGAL_CDT_3_TASK_BEGIN(merge_facets_task_handle); + if(!options.merge_facets) return patch_edges; + + CGAL_CDT_3_TASK_BEGIN(merge_facets_task_handle); + { auto start_time = std::chrono::high_resolution_clock::now(); if(options.merge_facets_old_method) { @@ -723,19 +726,23 @@ Borders_of_patches maybe_merge_facets( mesh, pmaps, options.coplanar_polygon_max_distance * bbox_max_span, options.coplanar_polygon_max_angle, options.dump_surface_mesh_after_merge_filename); } - if (!options.quiet) { - std::cout << "[timings] detected " << number_of_patches << " patches in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - } patch_edges = extract_patch_edges(mesh, pmaps, number_of_patches); - CGAL_CDT_3_TASK_END(merge_facets_task_handle); - if(!options.dump_patches_after_merge_filename.empty()) { - CGAL_CDT_3_TASK_BEGIN(output_task_handle); - std::ofstream out(options.dump_patches_after_merge_filename); - CGAL::IO::write_PLY(out, mesh, CGAL::parameters::stream_precision(17)); - CGAL_CDT_3_TASK_END(output_task_handle); + if (!options.quiet) { + auto timing = std::chrono::high_resolution_clock::now() - start_time; + std::cout << "[timings] found " << number_of_patches << " patches in " + << std::chrono::duration_cast(timing).count() << " ms\n"; } } + CGAL_CDT_3_TASK_END(merge_facets_task_handle); + + if(options.dump_patches_after_merge_filename.empty()) return patch_edges; + + CGAL_CDT_3_TASK_BEGIN(output_task_handle); + { + std::ofstream out(options.dump_patches_after_merge_filename); + CGAL::IO::write_PLY(out, mesh, CGAL::parameters::stream_precision(17)); + } + CGAL_CDT_3_TASK_END(output_task_handle); return patch_edges; } @@ -792,115 +799,101 @@ int go(Mesh mesh, CDT_options options) { auto bbox_max_span = (std::max)(bbox.x_span(), (std::max)(bbox.y_span(), bbox.z_span())); auto pmaps = setup_mesh_property_maps(mesh); + auto patch_edges = maybe_merge_facets(mesh, options, pmaps, bbox_max_span); + if(!options.dump_patches_borders_prefix.empty()) { CGAL_CDT_3_TASK_BEGIN(output_task_handle); dump_patches_borders(patch_edges, pmaps, options.dump_patches_borders_prefix); CGAL_CDT_3_TASK_END(output_task_handle); } - int exit_code = EXIT_SUCCESS; - - auto [mesh_descriptor_to_vertex_handle_pmap, tr_vertex_pmap_ok] = mesh.add_property_map("tr_vertex"); - assert(tr_vertex_pmap_ok); CGAL_USE(tr_vertex_pmap_ok); + auto mesh_descriptor_to_vertex_handle_pmap = get(CGAL::dynamic_vertex_property_t(), mesh); CGAL_CDT_3_TASK_BEGIN(insert_vertices_task_handle); - auto start_time = std::chrono::high_resolution_clock::now(); - CDT::Cell_handle hint{}; - for(auto v: vertices(mesh)) { - if(options.merge_facets && false == get(pmaps.v_selected_map, v)) continue; - auto vh = cdt.insert(get(pmaps.mesh_vertex_point_map, v), hint, false); - hint = vh->cell(); - put(mesh_descriptor_to_vertex_handle_pmap, v, vh); - } - if(!options.quiet) { - std::cout << "[timings] inserted vertices in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - std::cout << "Number of vertices: " << cdt.number_of_vertices() << "\n\n"; - } - if(cdt.dimension() < 3) { - if(!options.quiet) { - std::cout << "current is 2D... inserting the 8 vertices of an extended bounding box\n"; + { + auto start_time = std::chrono::high_resolution_clock::now(); + if(options.merge_facets) { + cdt.insert_mesh_vertices(mesh, pmaps.mesh_vertex_point_map, mesh_descriptor_to_vertex_handle_pmap, + CGAL::parameters::vertex_is_constrained_map(pmaps.v_selected_map)); + } else { + cdt.insert_mesh_vertices(mesh, pmaps.mesh_vertex_point_map, mesh_descriptor_to_vertex_handle_pmap); } - - auto bbox = CGAL::Polygon_mesh_processing::bbox(mesh); - auto dx = bbox.x_span(); - auto dy = bbox.y_span(); - auto dz = bbox.z_span(); - auto max_span = (std::max)(dx, (std::max)(dy, dz)); - if(dx == 0) dx = max_span; - if(dy == 0) dy = max_span; - if(dz == 0) dz = max_span; - - cdt.insert(Point(bbox.xmin() - dx, bbox.ymin() - dy, bbox.zmin() - dz)); - cdt.insert(Point(bbox.xmin() - dx, bbox.ymax() + dy, bbox.zmin() - dz)); - cdt.insert(Point(bbox.xmin() - dx, bbox.ymin() - dy, bbox.zmax() + dz)); - cdt.insert(Point(bbox.xmin() - dx, bbox.ymax() + dy, bbox.zmax() + dz)); - cdt.insert(Point(bbox.xmax() + dx, bbox.ymin() - dy, bbox.zmin() - dz)); - cdt.insert(Point(bbox.xmax() + dx, bbox.ymax() + dy, bbox.zmin() - dz)); - cdt.insert(Point(bbox.xmax() + dx, bbox.ymin() - dy, bbox.zmax() + dz)); - cdt.insert(Point(bbox.xmax() + dx, bbox.ymax() + dy, bbox.zmax() + dz)); + CDT::Cell_handle hint{}; + for(auto v: vertices(mesh)) { + if(options.merge_facets && false == get(pmaps.v_selected_map, v)) continue; + auto vh = cdt.insert(get(pmaps.mesh_vertex_point_map, v), hint, false); + hint = vh->cell(); + put(mesh_descriptor_to_vertex_handle_pmap, v, vh); + } + if(!options.quiet) { + auto timing = std::chrono::high_resolution_clock::now() - start_time; + std::cout << "[timings] inserted vertices in " + << std::chrono::duration_cast(timing).count() << " ms\n"; + std::cout << "Number of vertices: " << cdt.number_of_vertices() << "\n\n"; + } + cdt.add_bbox_points_if_not_dimension_3(); } CGAL_CDT_3_TASK_END(insert_vertices_task_handle); - start_time = std::chrono::high_resolution_clock::now(); CGAL_CDT_3_TASK_BEGIN(compute_distances_task_handle); + { + auto start_time = std::chrono::high_resolution_clock::now(); - exit_code = validate_minimum_vertex_distances(cdt, options.vertex_vertex_epsilon * bbox_max_span, options); - if(exit_code != EXIT_SUCCESS) { - return exit_code; + auto exit_code = validate_minimum_vertex_distances(cdt, options.vertex_vertex_epsilon * bbox_max_span, options); + if(exit_code != EXIT_SUCCESS) { + return exit_code; + } + + validate_constraint_vertex_distances_or_throw(cdt, mesh, options, patch_edges, mesh_descriptor_to_vertex_handle_pmap); + + if(!options.quiet) { + auto timing = std::chrono::high_resolution_clock::now() - start_time; + std::cout << "[timings] compute distances on " + << std::chrono::duration_cast(timing).count() << " ms\n"; + } } - - validate_constraint_vertex_distances_or_throw(cdt, mesh, options, patch_edges, mesh_descriptor_to_vertex_handle_pmap); CGAL_CDT_3_TASK_END(compute_distances_task_handle); - if(!options.quiet) { - std::cout << "[timings] compute distances on " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - } - CGAL_CDT_3_TASK_BEGIN(conforming_task_handle); - - int poly_id = 0; - auto output_on_exit_scope_guard = CGAL::make_scope_exit(create_output_finalizer(cdt, options)); - start_time = std::chrono::high_resolution_clock::now(); - if(options.merge_facets) { - for(auto& edges: patch_edges) { - cdt.insert_constrained_face_defined_by_its_border_edges(std::move(edges), pmaps.mesh_vertex_point_map, - mesh_descriptor_to_vertex_handle_pmap); - } - } else { - for(auto face_descriptor : faces(mesh)) { - std::vector polygon; - const auto he = halfedge(face_descriptor, mesh); - for(auto vertex_it : CGAL::vertices_around_face(he, mesh)) { - polygon.push_back(get(pmaps.mesh_vertex_point_map, vertex_it)); + { + int poly_id = 0; + auto output_on_exit_scope_guard = CGAL::make_scope_exit(create_output_finalizer(cdt, options)); + auto start_time = std::chrono::high_resolution_clock::now(); + if(options.merge_facets) { + for(auto& edges: patch_edges) { + cdt.insert_constrained_face_defined_by_its_border_edges(std::move(edges), pmaps.mesh_vertex_point_map, + mesh_descriptor_to_vertex_handle_pmap); } - if(cdt.debug_polygon_insertion()) { - std::cerr << "NEW POLYGON #" << poly_id << '\n'; - } - try { + } else { + for(auto face_descriptor : faces(mesh)) { + std::vector polygon; + const auto he = halfedge(face_descriptor, mesh); + for(auto vertex_it : CGAL::vertices_around_face(he, mesh)) { + polygon.push_back(get(pmaps.mesh_vertex_point_map, vertex_it)); + } + if(cdt.debug_polygon_insertion()) { + std::cerr << "NEW POLYGON #" << poly_id << '\n'; + } [[maybe_unused]] auto id = cdt.insert_constrained_polygon(polygon, false); assert(id == poly_id); ++poly_id; - } catch(int error) { - exit_code = error; + // std::ofstream dump("dump.binary.cgal"); + // CGAL::Mesh_3::save_binary_file(dump, cdt); } - // std::ofstream dump("dump.binary.cgal"); - // CGAL::Mesh_3::save_binary_file(dump, cdt); + } // not merge_facets + if(!options.quiet) { + auto timing = std::chrono::high_resolution_clock::now() - start_time; + std::cout << "[timings] registered facets in " + << std::chrono::duration_cast(timing).count() << " ms\n"; + } + start_time = std::chrono::high_resolution_clock::now(); + cdt.restore_Delaunay(); + if(!options.quiet) { + std::cout << "[timings] restored Delaunay (conforming of facets borders) in " << std::chrono::duration_cast( + std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; } - } // not merge_facets - if(!options.quiet) { - std::cout << "[timings] registered facets in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; } - start_time = std::chrono::high_resolution_clock::now(); - cdt.restore_Delaunay(); - if(!options.quiet) { - std::cout << "[timings] restored Delaunay (conforming of facets borders) in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - } - CGAL_CDT_3_TASK_END(conforming_task_handle); if(!options.dump_after_conforming_filename.empty()) { @@ -932,33 +925,35 @@ int go(Mesh mesh, CDT_options options) { if(!options.quiet) { std::cout << "Number of vertices after conforming: " << cdt.number_of_vertices() << "\n\n"; } + CGAL_CDT_3_TASK_BEGIN(validation_task_handle); - CGAL_assertion(!options.call_is_valid || cdt.Base_triantulation::is_valid(true)); - CGAL_assertion(!options.call_is_valid || cdt.is_valid(true)); - CGAL_assertion(!options.call_is_valid || cdt.is_conforming()); + { + CGAL_assertion(!options.call_is_valid || cdt.Base_triantulation::is_valid(true)); + CGAL_assertion(!options.call_is_valid || cdt.is_valid(true)); + CGAL_assertion(!options.call_is_valid || cdt.is_conforming()); + } CGAL_CDT_3_TASK_END(validation_task_handle); - if(exit_code == EXIT_SUCCESS) { - try { - CGAL_CDT_3_TASK_BEGIN(cdt_task_handle); - start_time = std::chrono::high_resolution_clock::now(); - cdt.restore_constrained_Delaunay(); - if(!options.quiet) { - std::cout << "[timings] restored constrained Delaunay in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - std::cout << "Number of vertices after CDT: " << cdt.number_of_vertices() << "\n\n"; - } - CGAL_CDT_3_TASK_END(cdt_task_handle); - } catch(int error) { - exit_code = error; + + CGAL_CDT_3_TASK_BEGIN(cdt_task_handle); + { + auto start_time = std::chrono::high_resolution_clock::now(); + cdt.restore_constrained_Delaunay(); + if(!options.quiet) { + std::cout << "[timings] restored constrained Delaunay in " << std::chrono::duration_cast( + std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; + std::cout << "Number of vertices after CDT: " << cdt.number_of_vertices() << "\n\n"; } } + CGAL_CDT_3_TASK_END(cdt_task_handle); CGAL_CDT_3_TASK_BEGIN(validation_task_handle); - CGAL_assertion(!options.call_is_valid || cdt.is_conforming()); - CGAL_assertion(!options.call_is_valid || cdt.is_valid(true)); + { + CGAL_assertion(!options.call_is_valid || cdt.is_conforming()); + CGAL_assertion(!options.call_is_valid || cdt.is_valid(true)); + } CGAL_CDT_3_TASK_END(validation_task_handle); - return exit_code; + return EXIT_SUCCESS; } int bisect_errors(Mesh mesh, CDT_options options) {