Skip to content

Commit 30faf58

Browse files
lrineaulrineau
committed
WIP: fix the sub-stars (with a retriangulation)
1 parent 78c510b commit 30faf58

1 file changed

Lines changed: 98 additions & 35 deletions

File tree

‎Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h‎

Lines changed: 98 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -1274,6 +1274,7 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
12741274
using Geom_traits = typename T_3::Geom_traits;
12751275
using Point_3 = typename T_3::Point;
12761276
using Segment_3 = typename Geom_traits::Segment_3;
1277+
using Triangle_3 = typename Geom_traits::Triangle_3;
12771278
using Vector_3 = typename Geom_traits::Vector_3;
12781279
using Locate_type = typename T_3::Locate_type;
12791280
using size_type = typename T_3::size_type;
@@ -1486,21 +1487,24 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
14861487
return std::nullopt;
14871488
}
14881489

1489-
template <typename Offsets, typename Facets>
1490+
template <typename Components, typename Map>
14901491
CGAL::Surface_mesh<Point_3>
14911492
construct_star_component_mesh(std::size_t component_index,
1492-
const Offsets& component_offsets,
1493-
const Facets& facets) const
1493+
const Components& components,
1494+
const Map& map_from_incident_constraint_face_id_to_new_triangles) const
14941495
{
14951496
CGAL::Surface_mesh<Point_3> component_mesh;
14961497
std::vector<Point_3> points;
14971498
std::vector<std::array<std::size_t, 3>> triangles;
14981499
std::size_t vertex_index = 0;
1499-
for(std::size_t facet_id = component_offsets[component_index - 1],
1500-
end = component_offsets[component_index];
1500+
for(std::size_t facet_id = components.component_offsets[component_index - 1],
1501+
end = components.component_offsets[component_index];
15011502
facet_id < end; ++facet_id)
15021503
{
1503-
const auto& [cell, facet_index, _] = facets[facet_id];
1504+
const auto& [cell, facet_index, type] = components.facets[facet_id];
1505+
if(type == Star_components_facets::Type_of_facet::INCIDENT_TO_V) {
1506+
continue;
1507+
}
15041508
if(tr().is_infinite(Facet(cell, facet_index))) {
15051509
continue;
15061510
}
@@ -1512,34 +1516,40 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
15121516
triangles.push_back({vertex_index, vertex_index + 2, vertex_index + 1});
15131517
vertex_index += 3;
15141518
} // end for each facet in the component
1519+
const auto& incident_face_ids = components.component_incident_constraint_face_ids[component_index - 1];
1520+
for(auto face_id : incident_face_ids) {
1521+
const auto& new_triangles = map_from_incident_constraint_face_id_to_new_triangles.at(face_id);
1522+
for(const auto& triangle : new_triangles) {
1523+
points.push_back(triangle[0]);
1524+
points.push_back(triangle[1]);
1525+
points.push_back(triangle[2]);
1526+
triangles.push_back({vertex_index, vertex_index + 1, vertex_index + 2});
1527+
vertex_index += 3;
1528+
}
1529+
}
15151530
CGAL::Polygon_mesh_processing::merge_duplicate_points_in_polygon_soup(points, triangles);
1516-
1517-
std::vector<std::pair<std::size_t, typename CGAL::Surface_mesh<Point_3>::Face_index>>
1518-
triangle_index_to_face_descriptor;
1519-
CGAL::Polygon_mesh_processing::polygon_soup_to_polygon_mesh(
1520-
points,
1521-
triangles,
1522-
component_mesh,
1523-
CGAL::parameters::polygon_to_face_output_iterator(std::back_inserter(triangle_index_to_face_descriptor)));
1531+
CGAL::Polygon_mesh_processing::orient_polygon_soup(points, triangles);
1532+
CGAL::Polygon_mesh_processing::polygon_soup_to_polygon_mesh(points, triangles, component_mesh);
15241533
return component_mesh;
15251534
}
15261535

15271536
std::optional<Point_3> compute_center_point_if_possible(const CGAL::Surface_mesh<Point_3>& component_mesh) const {
15281537
return CGAL::Polygon_mesh_processing::kernel_point(component_mesh, parameters::allow_open_input(true));
15291538
}
15301539

1531-
template <typename Offsets, typename Facets>
1540+
template <typename Components, typename Map>
15321541
CGAL::Surface_mesh<Point_3>
1533-
construct_link_mesh_for_debug(const Offsets& component_offsets,
1534-
const Facets& facets) const
1542+
construct_link_mesh_for_debug(const Components& components,
1543+
const Map& map_from_incident_constraint_face_id_to_new_triangles) const
15351544
{
15361545
using Mesh = CGAL::Surface_mesh<Point_3>;
15371546
using face_descriptor = typename Mesh::Face_index;
15381547
Mesh link_mesh;
15391548
auto [patch_id_map, _] = link_mesh.template add_property_map<face_descriptor, int>("f:patch_id", 0);
15401549

1541-
for(std::size_t component_index = 1; component_index < component_offsets.size(); ++component_index) {
1542-
Mesh component_mesh = construct_star_component_mesh(component_index, component_offsets, facets);
1550+
for(std::size_t component_index = 1; component_index < components.component_offsets.size(); ++component_index) {
1551+
Mesh component_mesh = construct_star_component_mesh(component_index, components,
1552+
map_from_incident_constraint_face_id_to_new_triangles);
15431553

15441554
std::vector<face_descriptor> new_faces;
15451555
auto out_it = boost::make_function_output_iterator(
@@ -1552,9 +1562,12 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
15521562
return link_mesh;
15531563
}
15541564

1555-
template <typename Offsets, typename Facets>
1556-
void dump_link_mesh_of_vertex_to_ply(Vertex_handle v, const Offsets& component_offsets, const Facets& facets) const {
1557-
const auto link_mesh = construct_link_mesh_for_debug(component_offsets, facets);
1565+
template <typename Components, typename Map>
1566+
void dump_link_mesh_of_vertex_to_ply(Vertex_handle v,
1567+
const Components& components,
1568+
const Map& map_from_incident_constraint_face_id_to_new_triangles) const {
1569+
const auto link_mesh =
1570+
construct_link_mesh_for_debug(components, map_from_incident_constraint_face_id_to_new_triangles);
15581571

15591572
std::stringstream filename;
15601573
filename << "dump_link_mesh_vertex" << IO::oformat(v, With_offset_tag{}) << ".ply";
@@ -1759,15 +1772,18 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
17591772
return out;
17601773
} // end collect_star_components_facets(..)
17611774

1762-
template <typename New_CDT_2_face_handles_output_iterator, typename Map_from_hole_edge_to_facet_3d>
1775+
enum Run_mode { NORMAL_RUN, DRY_RUN };
1776+
1777+
template <typename CDT_2, typename New_CDT_2_face_handles_output_iterator, typename Map_from_hole_edge_to_facet_3d>
17631778
std::pair<New_CDT_2_face_handles_output_iterator, Facet>
17641779
fill_cdt_2_hole_with_delaunay(CDT_2& cdt_2,
17651780
Face_index face_index,
17661781
CDT_3_vertex_type vertex_type,
17671782
std::list<CDT_2_edge>& hole_boundary,
17681783
const Map_from_hole_edge_to_facet_3d& map_from_hole_edge_to_facet_3d,
17691784
Face_is_outside_the_face is_outside_the_face,
1770-
New_CDT_2_face_handles_output_iterator new_face_handles_out)
1785+
New_CDT_2_face_handles_output_iterator new_face_handles_out,
1786+
Run_mode run_mode = NORMAL_RUN)
17711787
{
17721788
const bool is_a_half_hole = vertex_type == CDT_3_vertex_type::STEINER_ON_EDGE;
17731789
if(this->debug().move_Steiner_vertices()) {
@@ -1985,7 +2001,16 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
19852001

19862002
*new_face_handles_out++ = fh_2d;
19872003
}; // end of `action_for_each_new_2d_face`
1988-
CGAL::fill_hole_delaunay(cdt_2, hole_boundary, boost::make_function_output_iterator(action_for_each_new_2d_face));
2004+
switch(run_mode) {
2005+
case NORMAL_RUN:
2006+
CGAL::fill_hole_delaunay(cdt_2, hole_boundary, boost::make_function_output_iterator(action_for_each_new_2d_face));
2007+
break;
2008+
case DRY_RUN:
2009+
CGAL::fill_hole_delaunay(cdt_2, hole_boundary, new_face_handles_out);
2010+
break;
2011+
default:
2012+
CGAL_unreachable();
2013+
}
19892014
return {new_face_handles_out, half_hole_front_facet};
19902015
}
19912016

@@ -1994,12 +2019,20 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
19942019
remove_Steiner_vertex_from_cdt_2(Vertex_handle v,
19952020
CDT_3_vertex_type vertex_type,
19962021
Face_index face_id,
1997-
New_CDT_2_face_handles_output_iterator new_face_handles_out)
2022+
New_CDT_2_face_handles_output_iterator new_face_handles_out,
2023+
Run_mode run_mode = NORMAL_RUN)
19982024
{
19992025
std::tuple<New_CDT_2_face_handles_output_iterator, Facet, Facet> result{new_face_handles_out, {}, {}};
20002026
auto& [new_face_handles_out_it, half_hole_front_facet, second_half_hole_front_facet] = result;
2001-
auto& mutable_cdt_2 = non_const_face_cdt_2(face_id);
2002-
const auto& cdt_2 = mutable_cdt_2;
2027+
2028+
using CDT_2_base = typename CDT_2::Triangulation;
2029+
auto& original_mutable_cdt_2 = static_cast<CDT_2_base&>(non_const_face_cdt_2(face_id));
2030+
std::unique_ptr<CDT_2_base> cdt_2_copy_ptr{};
2031+
if(run_mode == DRY_RUN) {
2032+
cdt_2_copy_ptr.reset(new CDT_2_base(original_mutable_cdt_2));
2033+
}
2034+
auto& mutable_cdt_2 = (run_mode == DRY_RUN) ? *cdt_2_copy_ptr : original_mutable_cdt_2;
2035+
const auto & cdt_2 = mutable_cdt_2;
20032036

20042037
typename CDT_2::Locate_type lt;
20052038
int i;
@@ -2089,7 +2122,7 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
20892122
std::for_each(faces_to_delete.begin(), faces_to_delete.end(),
20902123
[&mutable_cdt_2](CDT_2_face_handle f) { mutable_cdt_2.delete_face(f); });
20912124
return fill_cdt_2_hole_with_delaunay(mutable_cdt_2, face_id, vertex_type, hole_boundary, hole_boundary_map,
2092-
is_outside_the_face, new_face_handles_out);
2125+
is_outside_the_face, new_face_handles_out, run_mode);
20932126
};
20942127

20952128
auto set_constrained = [&](CDT_2_edge e, bool is_constrained) {
@@ -2170,7 +2203,8 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
21702203
auto [out_it, facet_side_1] = bore_hole_and_fill_with_delaunay(
21712204
vertex_type, hole_side_1, hole_side_1_faces_to_delete, side_1_hole_is_outside_the_face);
21722205
new_face_handles_out = out_it;
2173-
CGAL_assertion(side_1_hole_is_outside_the_face == OUTSIDE || facet_side_1.first != Cell_handle{});
2206+
CGAL_assertion(run_mode == DRY_RUN ||
2207+
side_1_hole_is_outside_the_face == OUTSIDE || facet_side_1.first != Cell_handle{});
21742208

21752209
// - now the edge opposite to `(side_2_last_face, side_2_last_edge_index)` is on the border of the hole
21762210
// `hole_side_2`, so we can start sewing from there
@@ -2179,7 +2213,8 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
21792213
auto [out_it_2, facet_side_2] = bore_hole_and_fill_with_delaunay(
21802214
vertex_type, hole_side_2, hole_side_2_faces_to_delete, side_2_hole_is_outside_the_face);
21812215
new_face_handles_out = out_it_2;
2182-
CGAL_assertion(side_2_hole_is_outside_the_face == OUTSIDE || facet_side_2.first != Cell_handle{});
2216+
CGAL_assertion(run_mode == DRY_RUN ||
2217+
side_2_hole_is_outside_the_face == OUTSIDE || facet_side_2.first != Cell_handle{});
21832218

21842219
half_hole_front_facet = facet_side_1;
21852220
second_half_hole_front_facet = facet_side_2;
@@ -2264,13 +2299,42 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
22642299
}
22652300
return false;
22662301
}
2302+
const auto v_type = vertex_type(v);
2303+
2304+
std::map<Face_index, std::vector<Triangle_3>> map_from_incident_constraint_face_id_to_new_triangles;
2305+
2306+
auto dry_run_output_iterator = [&](Face_index face_id) {
2307+
auto* triangles_for_this_face_ptr = &map_from_incident_constraint_face_id_to_new_triangles[face_id];
2308+
std::function<void(CDT_2_face_handle)> func = [&,triangles_for_this_face_ptr](CDT_2_face_handle fh_2d) {
2309+
const auto v1 = vertex_3d(fh_2d->vertex(0));
2310+
const auto v2 = vertex_3d(fh_2d->vertex(1));
2311+
const auto v3 = vertex_3d(fh_2d->vertex(2));
2312+
triangles_for_this_face_ptr->emplace_back(v1->point(), v2->point(), v3->point());
2313+
};
2314+
return boost::make_function_output_iterator(func);
2315+
};
2316+
2317+
if(v_type == CDT_3_vertex_type::STEINER_IN_FACE) {
2318+
CGAL_assertion(star_components.component_incident_constraint_face_ids.size() == 2);
2319+
[[maybe_unused]] const auto& [face_id, _] = star_components.face_ids_between_consecutive_components.front();
2320+
remove_Steiner_vertex_from_cdt_2(v, v_type, face_id, dry_run_output_iterator(face_id), DRY_RUN);
2321+
} else {
2322+
CGAL_assertion(v_type == CDT_3_vertex_type::STEINER_ON_EDGE);
2323+
for(std::size_t component_index = 0, end = star_components.face_ids_between_consecutive_components.size();
2324+
component_index < end; ++component_index)
2325+
{
2326+
const auto& [face_id, _] = star_components.face_ids_between_consecutive_components[component_index];
2327+
2328+
remove_Steiner_vertex_from_cdt_2(v, v_type, face_id, dry_run_output_iterator(face_id), DRY_RUN);
2329+
}
2330+
}
22672331

22682332
const auto components_central_points = std::invoke([&]() {
22692333
std::optional<boost::container::small_vector<Point_3, 8>> central_points{std::in_place};
22702334
central_points->reserve(nb_of_star_components);
22712335
for(std::size_t component_index = 1; component_index <= nb_of_star_components; ++component_index) {
2272-
auto component_mesh =
2273-
construct_star_component_mesh(component_index, star_components.component_offsets, star_components.facets);
2336+
auto component_mesh = construct_star_component_mesh(component_index, star_components,
2337+
map_from_incident_constraint_face_id_to_new_triangles);
22742338
auto opt_kernel_center_point = compute_center_point_if_possible(component_mesh);
22752339
if(this->debug().move_Steiner_vertices()) {
22762340
std::cerr << " - incident constraint faces ids of component #" << component_index << ": {"
@@ -2300,7 +2364,7 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
23002364
}); // end compute central points for each component
23012365

23022366
if(this->debug().move_Steiner_vertices()) {
2303-
dump_link_mesh_of_vertex_to_ply(v, star_components.component_offsets, star_components.facets);
2367+
dump_link_mesh_of_vertex_to_ply(v, star_components, map_from_incident_constraint_face_id_to_new_triangles);
23042368
}
23052369

23062370
if(!components_central_points.has_value()) {
@@ -2356,7 +2420,6 @@ class Conforming_constrained_Delaunay_triangulation_3_impl : public Conforming_D
23562420
std::size_t nb_of_new_2d_faces = 0;
23572421
CGAL::Counting_output_iterator output_iterator{&nb_of_new_2d_faces};
23582422

2359-
const auto v_type = vertex_type(v);
23602423
if(v_type == CDT_3_vertex_type::STEINER_IN_FACE) {
23612424
CGAL_assertion(star_components.component_incident_constraint_face_ids.size() == 2);
23622425
const auto& [face_id, orientation] = star_components.face_ids_between_consecutive_components.front();

0 commit comments

Comments
 (0)