Merge branch 'PMP-corefine_non_manifold' into master

This commit is contained in:
Sébastien Loriot 2020-09-30 13:59:30 +02:00
commit 68124448dc
13 changed files with 967 additions and 169 deletions

View File

@ -37,6 +37,20 @@
namespace CGAL { namespace CGAL {
namespace parameters
{
template <class Parameter, class NamedParameters>
struct Is_default
{
typedef typename internal_np::Lookup_named_param_def <
Parameter,
NamedParameters,
internal_np::Param_not_found > ::type NP_type;
static const bool value = boost::is_same<NP_type, internal_np::Param_not_found>::value;
typedef CGAL::Boolean_tag<value> type;
};
} // end of parameters namespace
// forward declarations to avoid dependency to Solver_interface // forward declarations to avoid dependency to Solver_interface
template <typename FT, unsigned int dim> template <typename FT, unsigned int dim>
class Default_diagonalize_traits; class Default_diagonalize_traits;

View File

@ -102,6 +102,7 @@ CGAL_add_named_parameter(halfedges_keeper_t, halfedges_keeper, halfedges_keeper)
CGAL_add_named_parameter(volume_threshold_t, volume_threshold, volume_threshold) CGAL_add_named_parameter(volume_threshold_t, volume_threshold, volume_threshold)
CGAL_add_named_parameter(dry_run_t, dry_run, dry_run) CGAL_add_named_parameter(dry_run_t, dry_run, dry_run)
CGAL_add_named_parameter(do_not_modify_t, do_not_modify, do_not_modify) CGAL_add_named_parameter(do_not_modify_t, do_not_modify, do_not_modify)
CGAL_add_named_parameter(non_manifold_feature_map_t, non_manifold_feature_map, non_manifold_feature_map)
// List of named parameters that we use in the package 'Surface Mesh Simplification' // List of named parameters that we use in the package 'Surface Mesh Simplification'
CGAL_add_named_parameter(get_cost_policy_t, get_cost_policy, get_cost) CGAL_add_named_parameter(get_cost_policy_t, get_cost_policy, get_cost)

View File

@ -99,6 +99,7 @@ void test(const NamedParameters& np)
assert(get_parameter(np, CGAL::internal_np::do_simplify_border).v == 64); assert(get_parameter(np, CGAL::internal_np::do_simplify_border).v == 64);
assert(get_parameter(np, CGAL::internal_np::do_not_modify).v == 65); assert(get_parameter(np, CGAL::internal_np::do_not_modify).v == 65);
assert(get_parameter(np, CGAL::internal_np::maximum_number_of_faces).v == 78910); assert(get_parameter(np, CGAL::internal_np::maximum_number_of_faces).v == 78910);
assert(get_parameter(np, CGAL::internal_np::non_manifold_feature_map).v == 60);
// Named parameters that we use in the package 'Surface Mesh Simplification' // Named parameters that we use in the package 'Surface Mesh Simplification'
assert(get_parameter(np, CGAL::internal_np::get_cost_policy).v == 34); assert(get_parameter(np, CGAL::internal_np::get_cost_policy).v == 34);
@ -207,6 +208,7 @@ void test(const NamedParameters& np)
check_same_type<62>(get_parameter(np, CGAL::internal_np::halfedges_keeper)); check_same_type<62>(get_parameter(np, CGAL::internal_np::halfedges_keeper));
check_same_type<64>(get_parameter(np, CGAL::internal_np::do_simplify_border)); check_same_type<64>(get_parameter(np, CGAL::internal_np::do_simplify_border));
check_same_type<78910>(get_parameter(np, CGAL::internal_np::maximum_number_of_faces)); check_same_type<78910>(get_parameter(np, CGAL::internal_np::maximum_number_of_faces));
check_same_type<60>(get_parameter(np, CGAL::internal_np::non_manifold_feature_map));
// Named parameters that we use in the package 'Surface Mesh Simplification' // Named parameters that we use in the package 'Surface Mesh Simplification'
check_same_type<34>(get_parameter(np, CGAL::internal_np::get_cost_policy)); check_same_type<34>(get_parameter(np, CGAL::internal_np::get_cost_policy));
@ -325,6 +327,7 @@ int main()
.throw_on_self_intersection(A<43>(43)) .throw_on_self_intersection(A<43>(43))
.clip_volume(A<44>(44)) .clip_volume(A<44>(44))
.use_compact_clipper(A<45>(45)) .use_compact_clipper(A<45>(45))
.non_manifold_feature_map(A<60>(60))
.apply_per_connected_component(A<46>(46)) .apply_per_connected_component(A<46>(46))
.output_iterator(A<47>(47)) .output_iterator(A<47>(47))
.erase_all_duplicates(A<48>(48)) .erase_all_duplicates(A<48>(48))

View File

@ -83,6 +83,7 @@ create_single_source_cgal_program( "corefinement_mesh_union.cpp" )
create_single_source_cgal_program( "corefinement_mesh_union_and_intersection.cpp" ) create_single_source_cgal_program( "corefinement_mesh_union_and_intersection.cpp" )
create_single_source_cgal_program( "corefinement_mesh_union_with_attributes.cpp" ) create_single_source_cgal_program( "corefinement_mesh_union_with_attributes.cpp" )
create_single_source_cgal_program( "corefinement_polyhedron_union.cpp" ) create_single_source_cgal_program( "corefinement_polyhedron_union.cpp" )
create_single_source_cgal_program( "corefine_non_manifold.cpp" )
create_single_source_cgal_program( "random_perturbation_SM_example.cpp" ) create_single_source_cgal_program( "random_perturbation_SM_example.cpp" )
create_single_source_cgal_program( "corefinement_LCC.cpp") create_single_source_cgal_program( "corefinement_LCC.cpp")
create_single_source_cgal_program( "detect_features_example.cpp" ) create_single_source_cgal_program( "detect_features_example.cpp" )

View File

@ -0,0 +1,232 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/Polygon_mesh_processing/corefinement.h>
#include <CGAL/Polygon_mesh_processing/stitch_borders.h>
#include <CGAL/Polygon_mesh_processing/triangulate_faces.h>
#include <CGAL/boost/graph/helpers.h>
#include <fstream>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef K::Point_3 P;
typedef CGAL::Surface_mesh<P> Surface_mesh;
namespace PMP = CGAL::Polygon_mesh_processing;
int main()
{
// polyline intersection with a non-manifold edge
{
std::cout << "running polyline test\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,4,0), P(0,4,0), tm1);
CGAL::make_quad(P(0,-4,0), P(4,-4,0), P(4,0,0), P(0,0,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,0,4), P(0,0,4), tm1);
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,0,-4), P(0,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(2,3,-2), P(2,-2,-2), P(2,-2,2), P(2,3,2), tm2);
PMP::triangulate_faces(tm2);
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
std::vector< std::vector<P> > polylines;
PMP::surface_intersection(tm1, tm2, std::back_inserter(polylines), CGAL::parameters::non_manifold_feature_map(nm_map_1));
//dump polylines
std::ofstream output("intersection_polylines.cgal");
for(const std::vector<P>& polyline : polylines)
{
output << polyline.size() << " ";
std::copy(polyline.begin(), polyline.end(),std::ostream_iterator<P>(output," "));
output << "\n";
}
}
// simple case with only one non-manifold edge
{
std::cout << "running corefinement test 1\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,4,0), P(0,4,0), tm1);
CGAL::make_quad(P(0,-4,0), P(4,-4,0), P(4,0,0), P(0,0,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,0,4), P(0,0,4), tm1);
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,0,-4), P(0,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(2,3,-2), P(2,-2,-2), P(2,-2,2), P(2,3,2), tm2);
PMP::triangulate_faces(tm2);
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
PMP::corefine(tm1, tm2, CGAL::parameters::non_manifold_feature_map(nm_map_1));
std::ofstream("t1_tm1_corefined.off") << tm1;
std::ofstream("t1_tm2_corefined.off") << tm2;
}
// edge-edge intersection on a non-manifold edge
{
std::cout << "running corefinement test 2\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,4,0), P(0,4,0), tm1);
CGAL::make_quad(P(0,-4,0), P(4,-4,0), P(4,0,0), P(0,0,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,0,4), P(0,0,4), tm1);
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,0,-4), P(0,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(2,2,-2), P(2,-2,-2), P(2,-2,2), P(2,2,2), tm2); //TODO: test me + also test splitting the diagonal
PMP::triangulate_faces(tm2);
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
PMP::corefine(tm1, tm2, CGAL::parameters::non_manifold_feature_map(nm_map_1));
std::ofstream("t2_tm1_corefined.off") << tm1;
std::ofstream("t2_tm2_corefined.off") << tm2;
}
// coplanar edge and non-manifold edge
{
std::cout << "running corefinement test 3\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(8,0,0), P(8,8,0), P(0,8,0), tm1);
CGAL::make_quad(P(8,0,0), P(12,0,0), P(12,8,0), P(8,8,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(8,0,0), P(8,8,0), P(8,8,4), P(8,0,4), tm1);
CGAL::make_quad(P(8,0,0), P(8,8,0), P(8,8,-4), P(8,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(6,1,0), P(14,3,0), P(13,7,4), P(5,5,0), tm2);
CGAL::make_quad(P(6,1,0), P(5,5,0), P(5,5,2), P(6,1,3), tm2);
PMP::stitch_borders(tm2);
PMP::triangulate_faces(tm2);
std::ofstream("t3_tm1.off") << tm1;
std::ofstream("t3_tm2.off") << tm2;
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
PMP::corefine(tm1, tm2, CGAL::parameters::non_manifold_feature_map(nm_map_1));
std::ofstream("t3_tm1_corefined.off") << tm1;
std::ofstream("t3_tm2_corefined.off") << tm2;
}
//TODO: add more tests nm-edge on vertex, nm-edge vs nm-edge, ...
// coplanar face and non-manifold edge
{
std::cout << "running corefinement test 4\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(8,0,0), P(8,8,0), P(0,8,0), tm1);
CGAL::make_quad(P(8,0,0), P(12,0,0), P(12,8,0), P(8,8,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(8,0,0), P(8,8,0), P(8,8,4), P(8,0,4), tm1);
CGAL::make_quad(P(8,0,0), P(8,8,0), P(8,8,-4), P(8,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(6,1,0), P(14,3,0), P(13,7,0), P(5,5,0), tm2);
CGAL::make_quad(P(6,1,0), P(5,5,0), P(5,5,2), P(6,1,3), tm2);
PMP::stitch_borders(tm2);
PMP::triangulate_faces(tm2);
std::ofstream("t4_tm1.off") << tm1;
std::ofstream("t4_tm2.off") << tm2;
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
PMP::corefine(tm1, tm2, CGAL::parameters::non_manifold_feature_map(nm_map_1));
std::ofstream("t4_tm1_corefined.off") << tm1;
std::ofstream("t4_tm2_corefined.off") << tm2;
}
//
// // coplanar face and non-manifold edge and regular intersection with incident face
{
std::cout << "running corefinement test 5\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(4,0,0), P(4,8,0), P(0,8,0), tm1);
CGAL::make_quad(P(4,0,0), P(12,0,0), P(12,8,0), P(4,8,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(4,0,0), P(4,8,0), P(4,8,4), P(4,0,4), tm1);
CGAL::make_quad(P(4,0,0), P(4,8,0), P(4,8,-4), P(4,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(2,1,0), P(14,3,0), P(13,7,0), P(5,5,0), tm2);
CGAL::make_quad(P(2,1,0), P(5,5,0), P(5,5,2), P(2,1,3), tm2);
PMP::stitch_borders(tm2);
PMP::triangulate_faces(tm2);
std::ofstream("t5_tm1.off") << tm1;
std::ofstream("t5_tm2.off") << tm2;
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
PMP::corefine(tm1, tm2, CGAL::parameters::non_manifold_feature_map(nm_map_1));
std::ofstream("t5_tm1_corefined.off") << tm1;
std::ofstream("t5_tm2_corefined.off") << tm2;
}
// coplanar face and 2 non-manifold edges
{
std::cout << "running corefinement test 6\n";
Surface_mesh tm1;
CGAL::make_quad(P(0,0,0), P(8,0,0), P(8,8,0), P(0,8,0), tm1);
CGAL::make_quad(P(8,0,0), P(12,0,0), P(12,8,0), P(8,8,0), tm1);
PMP::stitch_borders(tm1);
CGAL::make_quad(P(8,0,0), P(8,8,0), P(8,8,4), P(8,0,4), tm1);
CGAL::make_quad(P(8,0,0), P(8,8,0), P(8,8,-4), P(8,0,-4), tm1);
PMP::triangulate_faces(tm1);
Surface_mesh tm2;
CGAL::make_quad(P(6,1,0), P(14,3,0), P(13,7,0), P(5,5,0), tm2);
CGAL::make_quad(P(6,1,0), P(5,5,0), P(5,5,2), P(6,1,3), tm2);
PMP::stitch_borders(tm2);
CGAL::make_quad(P(6,1,0), P(5,5,0), P(5,5,-2), P(6,1,-3), tm2);
CGAL::make_quad(P(6,1,0), P(5,5,0), P(3,5,0), P(3,1,0), tm2);
PMP::triangulate_faces(tm2);
std::ofstream("t6_tm1.off") << tm1;
std::ofstream("t6_tm2.off") << tm2;
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_1(tm1, get(boost::vertex_point, tm1));
PMP::Non_manifold_feature_map<Surface_mesh> nm_map_2(tm2, get(boost::vertex_point, tm2));
#if 0
std::vector< std::vector<P> > polylines;
PMP::surface_intersection(tm1, tm2, std::back_inserter(polylines),
CGAL::parameters::non_manifold_feature_map(nm_map_1),
CGAL::parameters::non_manifold_feature_map(nm_map_2));
//dump polylines
std::ofstream output("intersection_polylines.cgal");
for(const std::vector<P>& polyline : polylines)
{
output << polyline.size() << " ";
std::copy(polyline.begin(), polyline.end(),std::ostream_iterator<P>(output," "));
output << "\n";
}
#else
PMP::corefine(tm1, tm2, CGAL::parameters::non_manifold_feature_map(nm_map_1),
CGAL::parameters::non_manifold_feature_map(nm_map_2));
std::ofstream("t6_tm1_corefined.off") << tm1;
std::ofstream("t6_tm2_corefined.off") << tm2;
#endif
}
}

View File

@ -0,0 +1,114 @@
// Copyright (c) 2019 GeometryFactory (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
// You can redistribute it and/or modify it under the terms of the GNU
// General Public License as published by the Free Software Foundation,
// either version 3 of the License, or (at your option) any later version.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
//
// Author(s) : Sebastien Loriot
#ifndef CGAL_POLYGON_MESH_PROCESSING_NON_MANIFOLD_FEATURE_MAP
#define CGAL_POLYGON_MESH_PROCESSING_NON_MANIFOLD_FEATURE_MAP
namespace CGAL {
namespace Polygon_mesh_processing {
//TODO: right now the base name parameter mechanism will make a deep copy, we probably want to avoid that
template <class PolygonMesh>
struct Non_manifold_feature_map
{
typedef boost::graph_traits<PolygonMesh> Graph_traits;
typedef typename Graph_traits::vertex_descriptor vertex_descriptor;
typedef typename Graph_traits::edge_descriptor edge_descriptor;
typedef typename Graph_traits::halfedge_descriptor halfedge_descriptor;
typedef dynamic_edge_property_t<std::size_t> Edge_to_id_tag;
typedef dynamic_vertex_property_t<std::size_t> Vertex_to_id_tag;
typedef typename boost::property_map<PolygonMesh, Edge_to_id_tag>::type Edge_to_nm_id;
typedef typename boost::property_map<PolygonMesh, Vertex_to_id_tag>::type Vertex_to_nm_id;
Edge_to_nm_id e_nm_id;
Vertex_to_nm_id v_nm_id;
std::vector< std::vector<edge_descriptor> > non_manifold_edges;
std::vector< std::vector<vertex_descriptor> > non_manifold_vertices;
Non_manifold_feature_map()
{}
template <class Vpm>
Non_manifold_feature_map(PolygonMesh& pm, Vpm vpm)
: e_nm_id(get(Edge_to_id_tag(), pm))
, v_nm_id(get(Vertex_to_id_tag(), pm))
{
typedef typename boost::property_traits<Vpm>::value_type Point_3;
// detect non-manifold vertices
std::map<Point_3, std::vector<vertex_descriptor> > vertex_map;
for(vertex_descriptor vd : vertices(pm))
{
put(v_nm_id, vd, std::size_t(-1)); // init map
vertex_map[get(vpm, vd)].push_back(vd);
}
for(std::pair< const Point_3,
std::vector<vertex_descriptor> >& p : vertex_map)
{
if (p.second.size()!=1)
{
for (vertex_descriptor vd : p.second)
put(v_nm_id, vd, non_manifold_vertices.size());
non_manifold_vertices.resize(non_manifold_vertices.size()+1);
// we steal the vertor from the map
p.second.swap(non_manifold_vertices.back());
}
}
// detect non-manifold edges
std::map< std::pair<Point_3, Point_3>, std::vector<edge_descriptor> > edge_map;
for(edge_descriptor ed : edges(pm))
{
put(e_nm_id, ed, std::size_t(-1)); // init map
halfedge_descriptor hd = halfedge(ed, pm);
// an edge can be non-manifold only if both its vertices are non-manifold
if ( get(v_nm_id, source(hd, pm))==std::size_t(-1) ||
get(v_nm_id, target(hd, pm))==std::size_t(-1) ) continue;
const Point_3& src = get(vpm, source(ed, pm));
const Point_3& tgt = get(vpm, target(ed, pm));
// TODO: what to do with null edges?
if (src > tgt)
hd = opposite(hd, pm);
edge_map[ make_sorted_pair(src, tgt) ].push_back(edge(hd,pm));
}
for(std::pair< const std::pair<Point_3, Point_3>,
std::vector<edge_descriptor> >& p : edge_map)
{
if (p.second.size()!=1)
{
for (edge_descriptor ed : p.second)
put(e_nm_id, ed, non_manifold_edges.size());
non_manifold_edges.resize(non_manifold_edges.size()+1);
// we steal the vertor from the map
p.second.swap(non_manifold_edges.back());
}
}
}
};
} } // end of CGAL::Polygon_mesh_processing
#endif

View File

@ -783,14 +783,23 @@ corefine( TriangleMesh& tm1,
> ::type User_visitor; > ::type User_visitor;
User_visitor uv(choose_parameter<User_visitor>(get_parameter(np1, internal_np::graph_visitor))); User_visitor uv(choose_parameter<User_visitor>(get_parameter(np1, internal_np::graph_visitor)));
static const bool handle_non_manifold_features =
!parameters::Is_default<internal_np::non_manifold_feature_map_t, NamedParameters1>::value ||
!parameters::Is_default<internal_np::non_manifold_feature_map_t, NamedParameters2>::value;
// surface intersection algorithm call // surface intersection algorithm call
typedef Corefinement::No_extra_output_from_corefinement<TriangleMesh> Ob; typedef Corefinement::No_extra_output_from_corefinement<TriangleMesh> Ob;
typedef Corefinement::Surface_intersection_visitor_for_corefinement< typedef Corefinement::Surface_intersection_visitor_for_corefinement<
TriangleMesh, Vpm, Ob, Ecm, User_visitor> Algo_visitor; TriangleMesh, Vpm, Ob, Ecm, User_visitor, false, handle_non_manifold_features> Algo_visitor;
Ob ob; Ob ob;
Ecm ecm(tm1,tm2,ecm1,ecm2); Ecm ecm(tm1,tm2,ecm1,ecm2);
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor> Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor>
functor(tm1, tm2, vpm1, vpm2, Algo_visitor(uv,ob,ecm,const_mesh_ptr)); functor(tm1, tm2, vpm1, vpm2, Algo_visitor(uv,ob,ecm,const_mesh_ptr));
// Fill non-manifold feature maps if provided
functor.set_non_manifold_feature_map_1(parameters::get_parameter(np1, internal_np::non_manifold_feature_map));
functor.set_non_manifold_feature_map_2(parameters::get_parameter(np2, internal_np::non_manifold_feature_map));
functor(CGAL::Emptyset_iterator(), throw_on_self_intersection, true); functor(CGAL::Emptyset_iterator(), throw_on_self_intersection, true);
} }

View File

@ -871,7 +871,7 @@ public:
vpm1, vpm2, vpm1, vpm2,
nodes) ) //p1==q2 nodes) ) //p1==q2
{ {
CGAL_assertion( index_p1!=index_p2 || index_p1==Node_id((std::numeric_limits<Node_id>::max)()) ); CGAL_assertion( index_p1!=index_p2 || index_p1==NID );
coplanar_patches_of_tm1.set(patch_id_p1); coplanar_patches_of_tm1.set(patch_id_p1);
coplanar_patches_of_tm2.set(patch_id_q2); coplanar_patches_of_tm2.set(patch_id_q2);
bool q1_is_between_p1p2 = sorted_around_edge( bool q1_is_between_p1p2 = sorted_around_edge(
@ -941,17 +941,17 @@ public:
#endif //CGAL_COREFINEMENT_POLYHEDRA_DEBUG #endif //CGAL_COREFINEMENT_POLYHEDRA_DEBUG
CGAL_assertion( CGAL_assertion(
( index_p1 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm1,p1)): nodes.exact_node(index_p1) ) != ( index_p1 == NID ? nodes.to_exact(get(vpm1,p1)): nodes.exact_node(index_p1) ) !=
( index_q1 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm2,q1)): nodes.exact_node(index_q1) ) ( index_q1 == NID ? nodes.to_exact(get(vpm2,q1)): nodes.exact_node(index_q1) )
&& &&
( index_p2 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm1,p2)): nodes.exact_node(index_p2) ) != ( index_p2 == NID ? nodes.to_exact(get(vpm1,p2)): nodes.exact_node(index_p2) ) !=
( index_q1 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm2,q1)): nodes.exact_node(index_q1) ) ( index_q1 == NID ? nodes.to_exact(get(vpm2,q1)): nodes.exact_node(index_q1) )
&& &&
( index_p1 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm1,p1)): nodes.exact_node(index_p1) ) != ( index_p1 == NID ? nodes.to_exact(get(vpm1,p1)): nodes.exact_node(index_p1) ) !=
( index_q2 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm2,q2)): nodes.exact_node(index_q2) ) ( index_q2 == NID ? nodes.to_exact(get(vpm2,q2)): nodes.exact_node(index_q2) )
&& &&
( index_p2 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm1,p2)): nodes.exact_node(index_p2) ) != ( index_p2 == NID ? nodes.to_exact(get(vpm1,p2)): nodes.exact_node(index_p2) ) !=
( index_q2 == Node_id((std::numeric_limits<Node_id>::max)()) ? nodes.to_exact(get(vpm2,q2)): nodes.exact_node(index_q2) ) ( index_q2 == NID ? nodes.to_exact(get(vpm2,q2)): nodes.exact_node(index_q2) )
); );
bool q1_is_between_p1p2 = sorted_around_edge( bool q1_is_between_p1p2 = sorted_around_edge(

View File

@ -106,7 +106,8 @@ template< class TriangleMesh,
class OutputBuilder_ = Default, class OutputBuilder_ = Default,
class EdgeMarkMapBind_ = Default, class EdgeMarkMapBind_ = Default,
class UserVisitor_ = Default, class UserVisitor_ = Default,
bool doing_autorefinement = false > bool doing_autorefinement = false,
bool handle_non_manifold_features = false >
class Surface_intersection_visitor_for_corefinement{ class Surface_intersection_visitor_for_corefinement{
//default template parameters //default template parameters
typedef typename Default::Get<EdgeMarkMapBind_, typedef typename Default::Get<EdgeMarkMapBind_,
@ -133,13 +134,14 @@ private:
typedef boost::unordered_map<edge_descriptor,Node_ids> On_edge_map; typedef boost::unordered_map<edge_descriptor,Node_ids> On_edge_map;
//to keep the correspondance between node_id and vertex_handle in each mesh //to keep the correspondance between node_id and vertex_handle in each mesh
typedef std::vector<vertex_descriptor> Node_id_to_vertex; typedef std::vector<vertex_descriptor> Node_id_to_vertex;
typedef std::map<TriangleMesh*, Node_id_to_vertex > Mesh_to_map_node; typedef std::map<const TriangleMesh*, Node_id_to_vertex > Mesh_to_map_node;
//to handle coplanar halfedge of polyhedra that are full in the intersection //to handle coplanar halfedge of polyhedra that are full in the intersection
typedef std::multimap<Node_id,halfedge_descriptor> Node_to_target_of_hedge_map; typedef std::multimap<Node_id,halfedge_descriptor> Node_to_target_of_hedge_map;
typedef std::map<TriangleMesh*,Node_to_target_of_hedge_map> typedef std::map<TriangleMesh*,Node_to_target_of_hedge_map>
Mesh_to_vertices_on_intersection_map; Mesh_to_vertices_on_intersection_map;
typedef boost::unordered_map<vertex_descriptor,Node_id> Vertex_to_node_id; typedef boost::unordered_map<vertex_descriptor,Node_id> Vertex_to_node_id;
typedef std::map<TriangleMesh*, Vertex_to_node_id> Mesh_to_vertex_to_node_id; typedef std::map<TriangleMesh*, Vertex_to_node_id> Mesh_to_vertex_to_node_id;
typedef Non_manifold_feature_map<TriangleMesh> NM_features_map;
// typedef for the CDT // typedef for the CDT
typedef typename Intersection_nodes<TriangleMesh, typedef typename Intersection_nodes<TriangleMesh,
VertexPointMap, Predicates_on_constructions_needed>::Exact_kernel EK; VertexPointMap, Predicates_on_constructions_needed>::Exact_kernel EK;
@ -151,7 +153,7 @@ private:
typedef typename CDT::Vertex_handle CDT_Vertex_handle; typedef typename CDT::Vertex_handle CDT_Vertex_handle;
// data members // data members
private: private:
// boost::dynamic_bitset<> non_manifold_nodes; boost::dynamic_bitset<> is_node_on_boundary; // indicate if a vertex is a border vertex in tm1 or tm2
std::vector< std::vector<Node_id> > graph_of_constraints; std::vector< std::vector<Node_id> > graph_of_constraints;
boost::dynamic_bitset<> is_node_of_degree_one; boost::dynamic_bitset<> is_node_of_degree_one;
//nb of intersection points between coplanar faces, see fixes XSL_TAG_CPL_VERT //nb of intersection points between coplanar faces, see fixes XSL_TAG_CPL_VERT
@ -165,6 +167,9 @@ private:
std::map< Node_id,std::set<Node_id> > coplanar_constraints; std::map< Node_id,std::set<Node_id> > coplanar_constraints;
// optional data members to handle non-manifold issues
std::map<const TriangleMesh*, const NM_features_map*> non_manifold_feature_maps;
//data members that require initialization in the constructor //data members that require initialization in the constructor
UserVisitor& user_visitor; UserVisitor& user_visitor;
OutputBuilder& output_builder; OutputBuilder& output_builder;
@ -210,6 +215,70 @@ public:
, const_mesh_ptr(const_mesh_ptr) , const_mesh_ptr(const_mesh_ptr)
{} {}
void
set_non_manifold_feature_map(
const TriangleMesh& tm,
const NM_features_map& nm)
{
non_manifold_feature_maps[&tm] = &nm;
}
void copy_nodes_ids_for_non_manifold_features()
{
static const constexpr std::size_t NM_NID((std::numeric_limits<std::size_t>::max)());
for(const std::pair<const TriangleMesh*, const NM_features_map*>& tm_and_nm :
non_manifold_feature_maps)
{
TriangleMesh* tm_ptr = const_cast<TriangleMesh*>(tm_and_nm.first);
// update nodes on edges
On_edge_map& on_edge_map = on_edge[tm_ptr];
std::vector< std::pair<std::size_t, const Node_ids*> > edges_to_copy;
for (const std::pair<const edge_descriptor, Node_ids>& ed_and_ids : on_edge_map)
{
std::size_t eid = get(tm_and_nm.second->e_nm_id, ed_and_ids.first);
if (eid!=NM_NID)
edges_to_copy.push_back(std::make_pair(eid,&(ed_and_ids.second)));
}
for(const std::pair<const std::size_t, const Node_ids*>& id_and_nodes : edges_to_copy)
{
const std::vector<edge_descriptor>& nm_edges =
tm_and_nm.second->non_manifold_edges[id_and_nodes.first];
CGAL_assertion( on_edge_map.count(nm_edges.front())==1 );
for (std::size_t i=1; i<nm_edges.size(); ++i)
on_edge_map[nm_edges[i]] = *id_and_nodes.second;
}
// update map vertex -> node_id
Vertex_to_node_id& vertex_to_node_id = mesh_to_vertex_to_node_id[tm_ptr];
Node_to_target_of_hedge_map& vertices_on_inter = mesh_to_vertices_on_inter[tm_ptr];
std::vector< std::pair<vertex_descriptor, Node_id> > vertices_to_add;
for (const typename std::pair<const vertex_descriptor, Node_id>& vd_and_id
: vertex_to_node_id)
{
std::size_t vid = get(tm_and_nm.second->v_nm_id, vd_and_id.first);
if (vid!=NM_NID)
vertices_to_add.push_back(std::make_pair(vd_and_id.first,vd_and_id.second));
}
for(const std::pair<vertex_descriptor, Node_id>& vd_and_nid : vertices_to_add)
{
std::size_t vid = get(tm_and_nm.second->v_nm_id, vd_and_nid.first);
for(vertex_descriptor vd : tm_and_nm.second->non_manifold_vertices[vid])
{
if (vd != vd_and_nid.first)
{
vertex_to_node_id.insert(std::make_pair(vd,vd_and_nid.second));
output_builder.set_vertex_id(vd, vd_and_nid.second, *tm_ptr);
vertices_on_inter.insert(std::make_pair(vd_and_nid.second,halfedge(vd,*tm_ptr)));
}
}
}
}
}
template<class Graph_node> template<class Graph_node>
void annotate_graph(std::vector<Graph_node>& graph) void annotate_graph(std::vector<Graph_node>& graph)
{ {
@ -218,14 +287,32 @@ public:
is_node_of_degree_one.resize(nb_nodes); is_node_of_degree_one.resize(nb_nodes);
for(std::size_t node_id=0;node_id<nb_nodes;++node_id) for(std::size_t node_id=0;node_id<nb_nodes;++node_id)
{ {
// if (non_manifold_nodes.test(node_id))
// graph[node_id].make_terminal();
graph_of_constraints[node_id].assign( graph_of_constraints[node_id].assign(
graph[node_id].neighbors.begin(), graph[node_id].neighbors.begin(),
graph[node_id].neighbors.end()); graph[node_id].neighbors.end());
if (graph_of_constraints[node_id].size()==1) if (graph_of_constraints[node_id].size()==1)
is_node_of_degree_one.set(node_id); is_node_of_degree_one.set(node_id);
// mark every vertex on the boundary connected to another vertex by a non-boundary edge
// The logic is somehow equivalent to what was done with the container `non_manifold_nodes`
// that was used to split polylines at certains points. Non-manifold was used in the context
// of the combinatorial map where the import inside the combinatorial map was possible (see broken_bound-[12].off)
if (is_node_on_boundary.test(node_id))
{
// TODO: the commented condition is not restrictive enough. It only tests the vertices, not the edge used
// to go from one vertex to the other. Right now we are marking too many norder vertices as
// terminal than necessary. We cannot use mesh_to_node_id_to_vertex since it will be filled only
// after the call to visitor.finalize() in the main algorithm
// The current version of the code is correct but too many nodes are potentially
// marked as terminal
if (graph_of_constraints[node_id].size()==2 /* && (
!is_node_on_boundary.test(graph_of_constraints[node_id][0]) ||
!is_node_on_boundary.test(graph_of_constraints[node_id][1]) ) */ )
{
graph[node_id].make_terminal();
}
}
} }
} }
@ -261,33 +348,26 @@ public:
CGAL_error_msg("This function should not be called"); CGAL_error_msg("This function should not be called");
} }
// The following code was used to split polylines at certains points. void check_node_on_boundary_edge_case(std::size_t node_id,
// Here manifold was used in the context of the combinatorial map where halfedge_descriptor h,
// the import inside the combinatorial map was possible (see broken_bound-[12].off) const TriangleMesh& tm)
// I keep the code here for now as it could be use to detect non-manifold {
// situation of surfaces if ( is_border_edge(h,tm) )
// void check_node_on_non_manifold_edge( is_node_on_boundary.set(node_id);
// std::size_t node_id, }
// halfedge_descriptor h,
// const TriangleMesh& tm) void check_node_on_boundary_vertex_case(std::size_t node_id,
// { halfedge_descriptor h,
// if ( is_border_edge(h,tm) ) const TriangleMesh& tm)
// non_manifold_nodes.set(node_id); {
// } //we turn around the hedge and check no halfedge is a border halfedge
// for(halfedge_descriptor hc :halfedges_around_target(h,tm))
// void check_node_on_non_manifold_vertex( if ( is_border_edge(hc,tm) )
// std::size_t node_id, {
// halfedge_descriptor h, is_node_on_boundary.set(node_id);
// const TriangleMesh& tm) return;
// { }
// //we turn around the hedge and check no halfedge is a border halfedge }
// for(halfedge_descriptor hc :halfedges_around_target(h,tm))
// if ( is_border_edge(hc,tm) )
// {
// non_manifold_nodes.set(node_id);
// return;
// }
// }
//keep track of the fact that a polyhedron original vertex is a node //keep track of the fact that a polyhedron original vertex is a node
void all_incident_faces_got_a_node_as_vertex( void all_incident_faces_got_a_node_as_vertex(
@ -327,7 +407,7 @@ public:
bool is_target_coplanar, bool is_target_coplanar,
bool is_source_coplanar) bool is_source_coplanar)
{ {
// non_manifold_nodes.resize(node_id+1); is_node_on_boundary.resize(node_id+1, false);
TriangleMesh* tm1_ptr = const_cast<TriangleMesh*>(&tm1); TriangleMesh* tm1_ptr = const_cast<TriangleMesh*>(&tm1);
TriangleMesh* tm2_ptr = const_cast<TriangleMesh*>(&tm2); TriangleMesh* tm2_ptr = const_cast<TriangleMesh*>(&tm2);
@ -344,7 +424,7 @@ public:
case ON_EDGE: //Edge intersected by an edge case ON_EDGE: //Edge intersected by an edge
{ {
on_edge[tm2_ptr][edge(h_2,tm2)].push_back(node_id); on_edge[tm2_ptr][edge(h_2,tm2)].push_back(node_id);
// check_node_on_non_manifold_edge(node_id,h_2,tm2); check_node_on_boundary_edge_case(node_id,h_2,tm2);
} }
break; break;
case ON_VERTEX: case ON_VERTEX:
@ -356,7 +436,7 @@ public:
node_id_to_vertex.resize(node_id+1,Graph_traits::null_vertex()); node_id_to_vertex.resize(node_id+1,Graph_traits::null_vertex());
node_id_to_vertex[node_id]=target(h_2,tm2); node_id_to_vertex[node_id]=target(h_2,tm2);
all_incident_faces_got_a_node_as_vertex(h_2,node_id,*tm2_ptr); all_incident_faces_got_a_node_as_vertex(h_2,node_id,*tm2_ptr);
// check_node_on_non_manifold_vertex(node_id,h_2,tm2); check_node_on_boundary_vertex_case(node_id,h_2,tm2);
output_builder.set_vertex_id(target(h_2, tm2), node_id, tm2); output_builder.set_vertex_id(target(h_2, tm2), node_id, tm2);
} }
break; break;
@ -380,7 +460,7 @@ public:
all_incident_faces_got_a_node_as_vertex(h_1,node_id, *tm1_ptr); all_incident_faces_got_a_node_as_vertex(h_1,node_id, *tm1_ptr);
// register the vertex in the output builder // register the vertex in the output builder
output_builder.set_vertex_id(target(h_1, tm1), node_id, tm1); output_builder.set_vertex_id(target(h_1, tm1), node_id, tm1);
// check_node_on_non_manifold_vertex(node_id,h_1,tm1); check_node_on_boundary_vertex_case(node_id,h_1,tm1);
} }
else{ else{
if ( is_source_coplanar ){ if ( is_source_coplanar ){
@ -394,12 +474,34 @@ public:
all_incident_faces_got_a_node_as_vertex(h_1_opp,node_id, *tm1_ptr); all_incident_faces_got_a_node_as_vertex(h_1_opp,node_id, *tm1_ptr);
// register the vertex in the output builder // register the vertex in the output builder
output_builder.set_vertex_id(source(h_1, tm1), node_id, tm1); output_builder.set_vertex_id(source(h_1, tm1), node_id, tm1);
// check_node_on_non_manifold_vertex(node_id,h_1_opp,tm1); check_node_on_boundary_vertex_case(node_id,h_1_opp,tm1);
} }
else{ else{
//handle intersection on principal edge //handle intersection on principal edge
typename std::map<const TriangleMesh*, const NM_features_map*>::iterator it_find =
non_manifold_feature_maps.find(&tm1);
if ( it_find != non_manifold_feature_maps.end() )
{
// update h_1 if it is not the canonical non-manifold edge
// This is important to make sure intersection points on non-manifold
// edges are all connected for the same edge so that the redistribution
// on other edges does not overwrite some nodes.
// This update might be required in case of EDGE-EDGE intersection or
// COPLANAR intersection.
const NM_features_map& nm_features_map_1 = *it_find->second;
std::size_t eid1 = nm_features_map_1.non_manifold_edges.empty()
? std::size_t(-1)
: get(nm_features_map_1.e_nm_id, edge(h_1, tm1));
if (eid1 != std::size_t(-1))
{
if ( edge(h_1, tm1) != nm_features_map_1.non_manifold_edges[eid1].front() )
h_1 = halfedge(nm_features_map_1.non_manifold_edges[eid1].front(), tm1);
}
}
on_edge[tm1_ptr][edge(h_1,tm1)].push_back(node_id); on_edge[tm1_ptr][edge(h_1,tm1)].push_back(node_id);
// check_node_on_non_manifold_edge(node_id,h_1,tm1); check_node_on_boundary_edge_case(node_id,h_1,tm1);
} }
} }
} }
@ -468,8 +570,9 @@ public:
std::copy(begin,end,std::back_inserter(node_ids_array[it_id->second])); std::copy(begin,end,std::back_inserter(node_ids_array[it_id->second]));
} }
// Used by the autorefinement to re-set the id of nodes on the boundary of a // Used by the autorefinement and non-manifold edge handling to re-set
// face since another vertex (inside a face or on another edge) might have // the id of nodes on the boundary of a face since another vertex
// (inside a face or on another edge) might have
// overwritten the vertex in node_id_to_vertex // overwritten the vertex in node_id_to_vertex
template <class Node_id_to_vertex> template <class Node_id_to_vertex>
void update_node_id_to_vertex_map(Node_id_to_vertex& node_id_to_vertex, void update_node_id_to_vertex_map(Node_id_to_vertex& node_id_to_vertex,
@ -559,7 +662,7 @@ public:
// this condition ensures to consider only graph edges that are in // this condition ensures to consider only graph edges that are in
// the same triangle // the same triangle
if ( !points_on_triangle || it_vh!=id_to_CDT_vh.end() ){ if ( !points_on_triangle || it_vh!=id_to_CDT_vh.end() ){
CGAL_assertion(doing_autorefinement || it_vh!=id_to_CDT_vh.end()); CGAL_assertion(doing_autorefinement || handle_non_manifold_features || it_vh!=id_to_CDT_vh.end());
if (it_vh==id_to_CDT_vh.end()) continue; // needed for autorefinement (interior nodes) if (it_vh==id_to_CDT_vh.end()) continue; // needed for autorefinement (interior nodes)
cdt.insert_constraint(vh,it_vh->second); cdt.insert_constraint(vh,it_vh->second);
constrained_edges.push_back(std::make_pair(id,id_n)); constrained_edges.push_back(std::make_pair(id,id_n));
@ -599,14 +702,16 @@ public:
const VertexPointMap& vpm1, const VertexPointMap& vpm1,
const VertexPointMap& vpm2) const VertexPointMap& vpm2)
{ {
copy_nodes_ids_for_non_manifold_features();
nodes.all_nodes_created(); nodes.all_nodes_created();
TriangleMesh* tm1_ptr = const_cast<TriangleMesh*>(&tm1); TriangleMesh* tm1_ptr = const_cast<TriangleMesh*>(&tm1);
TriangleMesh* tm2_ptr = const_cast<TriangleMesh*>(&tm2); TriangleMesh* tm2_ptr = const_cast<TriangleMesh*>(&tm2);
std::map<TriangleMesh*, VertexPointMap> vpms; std::map<TriangleMesh*, VertexPointMap> vpms;
vpms[tm1_ptr] = vpm1; vpms.insert( std::make_pair(tm1_ptr, vpm1) );
vpms[tm2_ptr] = vpm2; vpms.insert( std::make_pair(tm2_ptr, vpm2) );
vertex_descriptor null_vertex = Graph_traits::null_vertex(); vertex_descriptor null_vertex = Graph_traits::null_vertex();
const Node_id nb_nodes = nodes.size(); const Node_id nb_nodes = nodes.size();
@ -635,15 +740,19 @@ public:
// Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm]; // Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm];
Node_to_target_of_hedge_map& nodes_to_hedge=it->second; Node_to_target_of_hedge_map& nodes_to_hedge=it->second;
// iterate on the vertices that are on the intersection between the input meshes
for(typename Node_to_target_of_hedge_map::iterator for(typename Node_to_target_of_hedge_map::iterator
it_node_2_hedge=nodes_to_hedge.begin(); it_node_2_hedge=nodes_to_hedge.begin();
it_node_2_hedge!=nodes_to_hedge.end(); it_node_2_hedge!=nodes_to_hedge.end();
++it_node_2_hedge) ++it_node_2_hedge)
{ {
Node_id node_id_of_first=it_node_2_hedge->first; Node_id node_id_of_first=it_node_2_hedge->first;
// look for neighbors of the current node in the intersection graph
std::vector<Node_id>& neighbors=graph_of_constraints[node_id_of_first]; std::vector<Node_id>& neighbors=graph_of_constraints[node_id_of_first];
if ( !neighbors.empty() ) if ( !neighbors.empty() )
{ {
// for all neighbors look for input vertices that are also on the intersection
for(Node_id node_id : neighbors) for(Node_id node_id : neighbors)
{ {
//if already done for the opposite //if already done for the opposite
@ -661,7 +770,7 @@ public:
target(it_node_2_hedge_two->second,tm) ) target(it_node_2_hedge_two->second,tm) )
{ {
hedge=opposite(next(hedge,tm),tm); hedge=opposite(next(hedge,tm),tm);
if (tm1_ptr==tm2_ptr && hedge==start) if ((doing_autorefinement || handle_non_manifold_features) && hedge==start)
{ {
++it_node_2_hedge_two; // we are using a multimap and ++it_node_2_hedge_two; // we are using a multimap and
// the halfedge we are looking for // the halfedge we are looking for
@ -825,7 +934,7 @@ public:
f_vertices[1]=it_fb->second.vertices[1]; f_vertices[1]=it_fb->second.vertices[1];
f_vertices[2]=it_fb->second.vertices[2]; f_vertices[2]=it_fb->second.vertices[2];
update_face_indices(f_vertices,f_indices,vertex_to_node_id); update_face_indices(f_vertices,f_indices,vertex_to_node_id);
if (doing_autorefinement) if (doing_autorefinement || handle_non_manifold_features)
it_fb->second.update_node_id_to_vertex_map(node_id_to_vertex, tm); it_fb->second.update_node_id_to_vertex_map(node_id_to_vertex, tm);
} }
else{ else{
@ -871,7 +980,7 @@ public:
{ {
id_to_CDT_vh.insert( id_to_CDT_vh.insert(
std::make_pair(f_indices[ik],triangle_vertices[ik])); std::make_pair(f_indices[ik],triangle_vertices[ik]));
if (doing_autorefinement) if (doing_autorefinement || handle_non_manifold_features)
// update the current vertex in node_id_to_vertex // update the current vertex in node_id_to_vertex
// to match the one of the face // to match the one of the face
node_id_to_vertex[f_indices[ik]]=f_vertices[ik]; node_id_to_vertex[f_indices[ik]]=f_vertices[ik];
@ -1040,7 +1149,7 @@ public:
} }
} }
nodes.finalize(); nodes.finalize(mesh_to_node_id_to_vertex);
// additional operations // additional operations
output_builder(nodes, output_builder(nodes,

View File

@ -521,28 +521,38 @@ struct Patch_container{
typedef typename GT::face_descriptor face_descriptor; typedef typename GT::face_descriptor face_descriptor;
Patch_description<PolygonMesh>& patch=this->operator[](i); Patch_description<PolygonMesh>& patch=this->operator[](i);
out << "OFF\n" << patch.interior_vertices.size() +
patch.shared_edges.size(); std::stringstream ss;
out << " " << patch.faces.size() << " 0\n";
std::map<vertex_descriptor, int> vertexid; std::map<vertex_descriptor, int> vertexid;
int id=0; int id=0;
for(vertex_descriptor vh : patch.interior_vertices) for(vertex_descriptor vh : patch.interior_vertices)
{ {
vertexid[vh]=id++; vertexid[vh]=id++;
out << get(vertex_point, pm, vh) << "\n"; ss << get(vertex_point, pm, vh) << "\n";
} }
for(halfedge_descriptor hh : patch.shared_edges) for(halfedge_descriptor hh : patch.shared_edges)
{ {
vertexid[target(hh, pm)]=id++; if( vertexid.insert( std::make_pair(target(hh,pm), id) ).second )
out << get(vertex_point, pm, target(hh, pm)) << "\n"; {
ss << get(vertex_point, pm, target(hh, pm)) << "\n";
++id;
}
if( vertexid.insert( std::make_pair(source(hh,pm), id) ).second )
{
ss << get(vertex_point, pm, source(hh, pm)) << "\n";
++id;
}
} }
out << "OFF\n" << id << " " << patch.faces.size() << " 0\n";
out << ss.str();
for(face_descriptor f : patch.faces) for(face_descriptor f : patch.faces)
{ {
out << "3 " << vertexid[source(halfedge(f,pm),pm)] << out << "3 " << vertexid.at(source(halfedge(f,pm),pm)) <<
" " << vertexid[target(halfedge(f,pm),pm)] << " " << vertexid.at(target(halfedge(f,pm),pm)) <<
" " << vertexid[target(next(halfedge(f,pm),pm),pm)] << "\n"; " " << vertexid.at(target(next(halfedge(f,pm),pm),pm)) << "\n";
} }
return out; return out;
@ -559,6 +569,16 @@ struct Patch_container{
dump_patch(i, output); dump_patch(i, output);
} }
} }
void dump_patches(std::string prefix)
{
for (std::size_t i=0;i <patches.size() ; ++i)
{
std::stringstream ss;
ss << prefix << "-" << i << ".off";
std::ofstream output(ss.str().c_str());
dump_patch(i, output);
}
}
}; };
template <class PolygonMesh, template <class PolygonMesh,
@ -1637,13 +1657,16 @@ void remove_unused_polylines(
do{ do{
halfedge_descriptor tmp_start = h; halfedge_descriptor starter = h;
while ( !is_border(h, tm) || is_border(opposite(h, tm), tm) ) while ( !is_border(h, tm) || is_border(opposite(h, tm), tm) )
{ {
h = opposite(next(h, tm), tm); h = opposite(next(h, tm), tm);
if (tmp_start==h) break; if (starter==h) break; // it might happen that no update is required because
// it is already closed thanks to another patch from tm'
// (usually incident to coplanar patches)
} }
if( !is_border(h, tm) ) if( !is_border(h, tm) || is_border(opposite(h, tm), tm) )
{ {
// nothing to do: the vertex has already been updated and is now in the middle of a patch kept. // nothing to do: the vertex has already been updated and is now in the middle of a patch kept.
// This function can be called after the stitching of the patches kept, the vertex halfedge // This function can be called after the stitching of the patches kept, the vertex halfedge

View File

@ -23,6 +23,7 @@
#include <CGAL/Polygon_mesh_processing/internal/Corefinement/intersection_of_coplanar_triangles_3.h> #include <CGAL/Polygon_mesh_processing/internal/Corefinement/intersection_of_coplanar_triangles_3.h>
#include <CGAL/Polygon_mesh_processing/internal/Corefinement/intersection_nodes.h> #include <CGAL/Polygon_mesh_processing/internal/Corefinement/intersection_nodes.h>
#include <CGAL/Polygon_mesh_processing/internal/Corefinement/intersect_triangle_and_segment_3.h> #include <CGAL/Polygon_mesh_processing/internal/Corefinement/intersect_triangle_and_segment_3.h>
#include <CGAL/Polygon_mesh_processing/Non_manifold_feature_map.h>
#include <CGAL/utility.h> #include <CGAL/utility.h>
#include <boost/unordered_map.hpp> #include <boost/unordered_map.hpp>
@ -98,9 +99,12 @@ struct Default_surface_intersection_visitor{
// If we implement a predicate only test, we can get rid of it. // If we implement a predicate only test, we can get rid of it.
static const bool Predicates_on_constructions_needed = doing_autorefinement; static const bool Predicates_on_constructions_needed = doing_autorefinement;
static const bool do_need_vertex_graph = false; static const bool do_need_vertex_graph = false;
void set_non_manifold_feature_map(
const TriangleMesh&,
const Non_manifold_feature_map<TriangleMesh>&)
{}
}; };
struct Node_id_set { struct Node_id_set {
typedef std::size_t Node_id; typedef std::size_t Node_id;
@ -187,12 +191,17 @@ class Intersection_of_triangle_meshes
Node_visitor visitor; Node_visitor visitor;
Faces_to_nodes_map f_to_node; //Associate a pair of triangles to their intersection points Faces_to_nodes_map f_to_node; //Associate a pair of triangles to their intersection points
std::vector<Node_id> extra_terminal_nodes; //used only for autorefinement std::vector<Node_id> extra_terminal_nodes; //used only for autorefinement
Non_manifold_feature_map<TriangleMesh> non_manifold_feature_map_1,
non_manifold_feature_map_2;
static const constexpr std::size_t NM_NID = (std::numeric_limits<std::size_t>::max)();
CGAL_assertion_code(bool doing_autorefinement;) CGAL_assertion_code(bool doing_autorefinement;)
// member functions // member functions
void filter_intersections(const TriangleMesh& tm_f, void filter_intersections(const TriangleMesh& tm_f,
const TriangleMesh& tm_e, const TriangleMesh& tm_e,
const VertexPointMap& vpm_f, const VertexPointMap& vpm_f,
const VertexPointMap& vpm_e, const VertexPointMap& vpm_e,
const Non_manifold_feature_map<TriangleMesh>& non_manifold_feature_map,
bool throw_on_self_intersection) bool throw_on_self_intersection)
{ {
std::vector<Box> face_boxes, edge_boxes; std::vector<Box> face_boxes, edge_boxes;
@ -213,15 +222,38 @@ class Intersection_of_triangle_meshes
edge_boxes.reserve(num_edges(tm_e)); edge_boxes.reserve(num_edges(tm_e));
edge_boxes_ptr.reserve(num_edges(tm_e)); edge_boxes_ptr.reserve(num_edges(tm_e));
for(edge_descriptor ed : edges(tm_e)) if (non_manifold_feature_map.non_manifold_edges.empty())
{ // general manifold case
halfedge_descriptor h=halfedge(ed,tm_e); for(edge_descriptor ed : edges(tm_e))
edge_boxes.push_back( Box( {
get(vpm_e,source(h,tm_e)).bbox() + halfedge_descriptor h=halfedge(ed,tm_e);
get(vpm_e,target(h,tm_e)).bbox(), edge_boxes.push_back( Box(
h ) ); get(vpm_e,source(h,tm_e)).bbox() +
edge_boxes_ptr.push_back( &edge_boxes.back() ); get(vpm_e,target(h,tm_e)).bbox(),
} h ) );
edge_boxes_ptr.push_back( &edge_boxes.back() );
}
else
// non-manifold case
for(edge_descriptor ed : edges(tm_e))
{
std::size_t eid=get(non_manifold_feature_map.e_nm_id, ed);
halfedge_descriptor h=halfedge(ed,tm_e);
// insert only one copy of a non-manifold edge
if (eid!=NM_NID)
{
if (non_manifold_feature_map.non_manifold_edges[eid].front()!=ed)
continue;
else
// make sure the halfedge used is consistant with stored one
h = halfedge(non_manifold_feature_map.non_manifold_edges[eid].front(), tm_e);
}
edge_boxes.push_back( Box(
get(vpm_e,source(h,tm_e)).bbox() +
get(vpm_e,target(h,tm_e)).bbox(),
h ) );
edge_boxes_ptr.push_back( &edge_boxes.back() );
}
/// \todo experiments different cutoff values /// \todo experiments different cutoff values
std::ptrdiff_t cutoff = 2 * std::ptrdiff_t( std::ptrdiff_t cutoff = 2 * std::ptrdiff_t(
@ -311,18 +343,34 @@ class Intersection_of_triangle_meshes
get_or_create_node(const Cpl_inter_pt& ipt, get_or_create_node(const Cpl_inter_pt& ipt,
Node_id& current_node, Node_id& current_node,
std::map<Key,Node_id>& coplanar_node_map, std::map<Key,Node_id>& coplanar_node_map,
const TriangleMesh& tm1, const Non_manifold_feature_map<TriangleMesh>& nm_features_map_1,
const TriangleMesh& tm2) const Non_manifold_feature_map<TriangleMesh>& nm_features_map_2,
const TriangleMesh& tm1,
const TriangleMesh& tm2)
{ {
halfedge_descriptor h1=graph_traits::null_halfedge(),h2=h1; halfedge_descriptor h1=graph_traits::null_halfedge(),h2=h1;
switch(ipt.type_1){ switch(ipt.type_1){
case ON_VERTEX: case ON_VERTEX:
h1=halfedge(target(ipt.info_1,tm1),tm1); {
std::size_t vid1 = nm_features_map_1.non_manifold_vertices.empty()
? NM_NID
: get(nm_features_map_1.v_nm_id, target(ipt.info_1,tm1));
if (vid1==NM_NID)
h1=halfedge(target(ipt.info_1,tm1),tm1);
else
h1=halfedge(nm_features_map_1.non_manifold_vertices[vid1][0],tm1);
}
break; break;
case ON_EDGE : case ON_EDGE :
{ {
h1=opposite(ipt.info_1,tm1); std::size_t eid1 = nm_features_map_1.non_manifold_edges.empty()
if (h1>ipt.info_1) h1=ipt.info_1; ? NM_NID
: get(nm_features_map_1.e_nm_id, edge(ipt.info_1,tm1));
if (eid1==NM_NID)
h1=ipt.info_1;
else
h1=halfedge(nm_features_map_1.non_manifold_edges[eid1][0],tm1);
h1=(std::max)(h1, opposite(h1, tm1));
} }
break; break;
case ON_FACE : case ON_FACE :
@ -332,12 +380,26 @@ class Intersection_of_triangle_meshes
} }
switch(ipt.type_2){ switch(ipt.type_2){
case ON_VERTEX: case ON_VERTEX:
h2=halfedge(target(ipt.info_2,tm2),tm2); {
std::size_t vid2 = nm_features_map_2.non_manifold_vertices.empty()
? NM_NID
: get(nm_features_map_2.v_nm_id, target(ipt.info_2,tm2));
if (vid2==NM_NID)
h2=halfedge(target(ipt.info_2,tm2),tm2);
else
h2=halfedge(nm_features_map_2.non_manifold_vertices[vid2][0],tm2);
}
break; break;
case ON_EDGE : case ON_EDGE :
{ {
h2=opposite(ipt.info_2,tm2); std::size_t eid2 = nm_features_map_2.non_manifold_edges.empty()
if (h2>ipt.info_2) h2=ipt.info_2; ? NM_NID
: get(nm_features_map_2.e_nm_id, edge(ipt.info_2,tm2));
if (eid2==NM_NID)
h2=ipt.info_2;
else
h2=halfedge(nm_features_map_2.non_manifold_edges[eid2][0],tm2);
h2=(std::max)(h2, opposite(h2, tm2));
} }
break; break;
case ON_FACE : case ON_FACE :
@ -360,14 +422,14 @@ class Intersection_of_triangle_meshes
} }
void add_intersection_point_to_face_and_all_edge_incident_faces(face_descriptor f_1, void add_intersection_point_to_face_and_all_edge_incident_faces(face_descriptor f_1,
halfedge_descriptor e_2, halfedge_descriptor h_2,
const TriangleMesh& tm1, const TriangleMesh& tm1,
const TriangleMesh& tm2, const TriangleMesh& tm2,
Node_id node_id) Node_id node_id)
{ {
if (!is_border(e_2, tm2)) if (!is_border(h_2, tm2))
{ {
face_descriptor f_2 = face(e_2, tm2); face_descriptor f_2 = face(h_2, tm2);
if(&tm1!=&tm2 || f_1!=f_2) if(&tm1!=&tm2 || f_1!=f_2)
{ {
Face_pair face_pair = &tm1==&tm2 ? make_sorted_pair(f_1,f_2): Face_pair face_pair = &tm1==&tm2 ? make_sorted_pair(f_1,f_2):
@ -378,10 +440,10 @@ class Intersection_of_triangle_meshes
f_to_node[Face_pair_and_int(face_pair,0)].insert(node_id); f_to_node[Face_pair_and_int(face_pair,0)].insert(node_id);
} }
} }
e_2 = opposite(e_2, tm2); h_2 = opposite(h_2, tm2);
if (!is_border(e_2, tm2)) if (!is_border(h_2, tm2))
{ {
face_descriptor f_2 = face(e_2, tm2); face_descriptor f_2 = face(h_2, tm2);
if(&tm1!=&tm2 || f_1!=f_2) if(&tm1!=&tm2 || f_1!=f_2)
{ {
Face_pair face_pair = &tm1==&tm2 ? make_sorted_pair(f_1,f_2): Face_pair face_pair = &tm1==&tm2 ? make_sorted_pair(f_1,f_2):
@ -462,6 +524,7 @@ class Intersection_of_triangle_meshes
halfedge_descriptor f_2, halfedge_descriptor f_2,
const TriangleMesh& tm1, const TriangleMesh& tm1,
const TriangleMesh& tm2, const TriangleMesh& tm2,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_1,
Node_id node_id, Node_id node_id,
bool is_new_node) bool is_new_node)
{ {
@ -472,42 +535,78 @@ class Intersection_of_triangle_meshes
? stm_edge_to_ltm_faces ? stm_edge_to_ltm_faces
: ltm_edge_to_stm_faces; : ltm_edge_to_stm_faces;
for(halfedge_descriptor h_1 : std::vector<vertex_descriptor> tmp_vertices_1(1, target(v_1, tm1));
halfedges_around_target(v_1,tm1))
{ std::size_t vid1 = nm_features_map_1.non_manifold_vertices.empty()
add_intersection_point_to_face_and_all_edge_incident_faces(face(f_2,tm2),h_1,tm2,tm1,node_id); ? NM_NID
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(h_1,tm1)); : get(nm_features_map_1.v_nm_id, target(v_1, tm1));
if (it_ets!=tm1_edge_to_tm2_faces.end()) it_ets->second.erase(face(f_2,tm2));
} const std::vector<vertex_descriptor>& vertices_1 = vid1==NM_NID
? tmp_vertices_1
: nm_features_map_1.non_manifold_vertices[vid1];
for(vertex_descriptor v1 : vertices_1)
for(halfedge_descriptor h_1 :
halfedges_around_target(v1,tm1))
{
add_intersection_point_to_face_and_all_edge_incident_faces(face(f_2,tm2),h_1,tm2,tm1,node_id);
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(h_1,tm1));
if (it_ets!=tm1_edge_to_tm2_faces.end()) it_ets->second.erase(face(f_2,tm2));
}
} }
void handle_coplanar_case_VERTEX_EDGE(halfedge_descriptor v_1, void handle_coplanar_case_VERTEX_EDGE(halfedge_descriptor v_1,
halfedge_descriptor e_2, halfedge_descriptor h_2,
const TriangleMesh& tm1, const TriangleMesh& tm1,
const TriangleMesh& tm2, const TriangleMesh& tm2,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_1,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_2,
Node_id node_id, Node_id node_id,
bool is_new_node) bool is_new_node)
{ {
if(is_new_node) if(is_new_node)
visitor.new_node_added(node_id,ON_VERTEX,e_2,v_1,tm2,tm1,false,false); visitor.new_node_added(node_id,ON_VERTEX,h_2,v_1,tm2,tm1,false,false);
Edge_to_faces& tm1_edge_to_tm2_faces = &tm1 <= &tm2 Edge_to_faces& tm1_edge_to_tm2_faces = &tm1 <= &tm2
? stm_edge_to_ltm_faces ? stm_edge_to_ltm_faces
: ltm_edge_to_stm_faces; : ltm_edge_to_stm_faces;
for(halfedge_descriptor h_1 : std::vector<vertex_descriptor> tmp_vertices_1(1, target(v_1, tm1));
halfedges_around_target(v_1,tm1))
{ std::size_t vid1 = nm_features_map_1.non_manifold_vertices.empty()
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(h_1,tm1)); ? NM_NID
Face_set* fset = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr; : get(nm_features_map_1.v_nm_id, target(v_1, tm1));
cip_handle_case_edge(node_id,fset,h_1,e_2,tm1,tm2);
} const std::vector<vertex_descriptor>& vertices_1 = vid1==NM_NID
? tmp_vertices_1
: nm_features_map_1.non_manifold_vertices[vid1];
std::vector<edge_descriptor> tmp_edges_2(1, edge(h_2, tm2));
std::size_t eid2 = nm_features_map_2.non_manifold_edges.empty()
? NM_NID
: get(nm_features_map_2.e_nm_id, edge(h_2, tm2));
const std::vector<edge_descriptor>& edges_2 = eid2==NM_NID
? tmp_edges_2
: nm_features_map_2.non_manifold_edges[eid2];
for (vertex_descriptor v1 : vertices_1)
for(halfedge_descriptor h_1 :
halfedges_around_target(v1,tm1))
{
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(h_1,tm1));
Face_set* fset = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr;
for (edge_descriptor e2 : edges_2)
cip_handle_case_edge(node_id,fset,h_1,halfedge(e2, tm2),tm1,tm2);
}
} }
void handle_coplanar_case_VERTEX_VERTEX(halfedge_descriptor v_1, void handle_coplanar_case_VERTEX_VERTEX(halfedge_descriptor v_1,
halfedge_descriptor v_2, halfedge_descriptor v_2,
const TriangleMesh& tm1, const TriangleMesh& tm1,
const TriangleMesh& tm2, const TriangleMesh& tm2,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_1,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_2,
Node_id node_id, Node_id node_id,
bool is_new_node) bool is_new_node)
{ {
@ -518,13 +617,33 @@ class Intersection_of_triangle_meshes
? stm_edge_to_ltm_faces ? stm_edge_to_ltm_faces
: ltm_edge_to_stm_faces; : ltm_edge_to_stm_faces;
for(halfedge_descriptor h_1 : std::vector<vertex_descriptor> tmp_vertices_1(1, target(v_1, tm1)),
halfedges_around_target(v_1,tm1)) tmp_vertices_2(1, target(v_2, tm2));
{
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(h_1,tm1)); std::size_t vid1 = nm_features_map_1.non_manifold_vertices.empty()
Face_set* fset = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr; ? NM_NID
cip_handle_case_vertex(node_id,fset,h_1,v_2,tm1,tm2); : get(nm_features_map_1.v_nm_id, target(v_1, tm1));
}
std::size_t vid2 = nm_features_map_2.non_manifold_vertices.empty()
? NM_NID
: get(nm_features_map_2.v_nm_id, target(v_2, tm2));
const std::vector<vertex_descriptor>& vertices_1 = vid1==NM_NID
? tmp_vertices_1
: nm_features_map_1.non_manifold_vertices[vid1];
const std::vector<vertex_descriptor>& vertices_2 = vid2==NM_NID
? tmp_vertices_2
: nm_features_map_2.non_manifold_vertices[vid2];
for (vertex_descriptor v1 : vertices_1)
for(halfedge_descriptor h_1 :
halfedges_around_target(v1,tm1))
{
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(h_1,tm1));
Face_set* fset = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr;
for (vertex_descriptor v2 : vertices_2)
cip_handle_case_vertex(node_id,fset,h_1,halfedge(v2, tm2),tm1,tm2);
}
} }
void compute_intersection_of_coplanar_faces( void compute_intersection_of_coplanar_faces(
@ -532,7 +651,9 @@ class Intersection_of_triangle_meshes
const TriangleMesh& tm1, const TriangleMesh& tm1,
const TriangleMesh& tm2, const TriangleMesh& tm2,
const VertexPointMap& vpm1, const VertexPointMap& vpm1,
const VertexPointMap& vpm2) const VertexPointMap& vpm2,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_1,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_2)
{ {
CGAL_assertion( &tm1 < &tm2 || &tm1==&tm2 ); CGAL_assertion( &tm1 < &tm2 || &tm1==&tm2 );
@ -566,7 +687,7 @@ class Intersection_of_triangle_meshes
Node_id node_id; Node_id node_id;
bool is_new_node; bool is_new_node;
std::tie(node_id, is_new_node) = std::tie(node_id, is_new_node) =
get_or_create_node(ipt,current_node,coplanar_node_map,tm1,tm2); get_or_create_node(ipt,current_node,coplanar_node_map,nm_features_map_1,nm_features_map_2,tm1,tm2);
cpln_nodes.push_back(node_id); cpln_nodes.push_back(node_id);
switch(ipt.type_1){ switch(ipt.type_1){
@ -574,13 +695,13 @@ class Intersection_of_triangle_meshes
{ {
switch(ipt.type_2){ switch(ipt.type_2){
case ON_VERTEX: case ON_VERTEX:
handle_coplanar_case_VERTEX_VERTEX(ipt.info_1,ipt.info_2,tm1,tm2,node_id,is_new_node); handle_coplanar_case_VERTEX_VERTEX(ipt.info_1,ipt.info_2,tm1,tm2,nm_features_map_1,nm_features_map_2,node_id,is_new_node);
break; break;
case ON_EDGE: case ON_EDGE:
handle_coplanar_case_VERTEX_EDGE(ipt.info_1,ipt.info_2,tm1,tm2,node_id,is_new_node); handle_coplanar_case_VERTEX_EDGE(ipt.info_1,ipt.info_2,tm1,tm2,nm_features_map_1,nm_features_map_2,node_id,is_new_node);
break; break;
case ON_FACE: case ON_FACE:
handle_coplanar_case_VERTEX_FACE(ipt.info_1,ipt.info_2,tm1,tm2,node_id,is_new_node); handle_coplanar_case_VERTEX_FACE(ipt.info_1,ipt.info_2,tm1,tm2,nm_features_map_1,node_id,is_new_node);
break; break;
default: CGAL_error_msg("Should not get there!"); default: CGAL_error_msg("Should not get there!");
} }
@ -590,15 +711,33 @@ class Intersection_of_triangle_meshes
{ {
switch(ipt.type_2){ switch(ipt.type_2){
case ON_VERTEX: case ON_VERTEX:
handle_coplanar_case_VERTEX_EDGE(ipt.info_2,ipt.info_1,tm2,tm1,node_id,is_new_node); handle_coplanar_case_VERTEX_EDGE(ipt.info_2,ipt.info_1,tm2,tm1,nm_features_map_2,nm_features_map_1,node_id,is_new_node);
break; break;
case ON_EDGE: case ON_EDGE:
{ {
std::vector<edge_descriptor> tmp_edges_1(1, edge(ipt.info_1,tm1)),
tmp_edges_2(1, edge(ipt.info_2,tm2));
std::size_t eid1 = nm_features_map_1.non_manifold_edges.empty()
? NM_NID
: get(nm_features_map_1.e_nm_id, edge(ipt.info_1, tm1));
std::size_t eid2 = nm_features_map_2.non_manifold_edges.empty()
? NM_NID
: get(nm_features_map_2.e_nm_id, edge(ipt.info_2, tm2));
const std::vector<edge_descriptor>& edges_1 = eid1==NM_NID
? tmp_edges_1
: nm_features_map_1.non_manifold_edges[eid1];
const std::vector<edge_descriptor>& edges_2 = eid2==NM_NID
? tmp_edges_2
: nm_features_map_2.non_manifold_edges[eid2];
if(is_new_node) if(is_new_node)
visitor.new_node_added(node_id,ON_EDGE,ipt.info_1,ipt.info_2,tm1,tm2,false,false); visitor.new_node_added(node_id,ON_EDGE,halfedge(edges_1.front(), tm1),halfedge(edges_2.front(),tm2),tm1,tm2,false,false);
typename Edge_to_faces::iterator it_ets=stm_edge_to_ltm_faces.find(edge(ipt.info_1,tm1)); for(edge_descriptor e1 : edges_1)
Face_set* fset = (it_ets!=stm_edge_to_ltm_faces.end())?&(it_ets->second):nullptr; for(edge_descriptor e2 : edges_2)
cip_handle_case_edge(node_id,fset,ipt.info_1,ipt.info_2,tm1,tm2); {
typename Edge_to_faces::iterator it_ets=stm_edge_to_ltm_faces.find(e1);
Face_set* fset = (it_ets!=stm_edge_to_ltm_faces.end())?&(it_ets->second):nullptr;
cip_handle_case_edge(node_id,fset,halfedge(e1,tm1),halfedge(e2,tm2),tm1,tm2);
}
} }
break; break;
default: CGAL_error_msg("Should not get there!"); default: CGAL_error_msg("Should not get there!");
@ -608,7 +747,7 @@ class Intersection_of_triangle_meshes
case ON_FACE: case ON_FACE:
{ {
CGAL_assertion(ipt.type_2==ON_VERTEX); CGAL_assertion(ipt.type_2==ON_VERTEX);
handle_coplanar_case_VERTEX_FACE(ipt.info_2,ipt.info_1,tm2,tm1,node_id,is_new_node); handle_coplanar_case_VERTEX_FACE(ipt.info_2,ipt.info_1,tm2,tm1,nm_features_map_2,node_id,is_new_node);
} }
break; break;
default: CGAL_error_msg("Should not get there!"); default: CGAL_error_msg("Should not get there!");
@ -662,6 +801,8 @@ class Intersection_of_triangle_meshes
const TriangleMesh& tm2, const TriangleMesh& tm2,
const VertexPointMap& vpm1, const VertexPointMap& vpm1,
const VertexPointMap& vpm2, const VertexPointMap& vpm2,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_1,
const Non_manifold_feature_map<TriangleMesh>& nm_features_map_2,
Node_id& current_node) Node_id& current_node)
{ {
typedef std::tuple<Intersection_type, halfedge_descriptor, bool,bool> Inter_type; typedef std::tuple<Intersection_type, halfedge_descriptor, bool,bool> Inter_type;
@ -671,6 +812,7 @@ class Intersection_of_triangle_meshes
it!=tm1_edge_to_tm2_faces.end();++it) it!=tm1_edge_to_tm2_faces.end();++it)
{ {
edge_descriptor e_1=it->first; edge_descriptor e_1=it->first;
halfedge_descriptor h_1=halfedge(e_1,tm1); halfedge_descriptor h_1=halfedge(e_1,tm1);
Face_set& fset=it->second; Face_set& fset=it->second;
while (!fset.empty()){ while (!fset.empty()){
@ -682,18 +824,93 @@ class Intersection_of_triangle_meshes
//handle degenerate case: one extremity of edge belong to f_2 //handle degenerate case: one extremity of edge belong to f_2
std::vector<halfedge_descriptor> all_edges; std::vector<halfedge_descriptor> all_edges;
if ( std::get<3>(res) ) // is edge target in triangle plane if ( std::get<3>(res) ) // is edge target in triangle plane
std::copy(halfedges_around_target(h_1,tm1).first, {
halfedges_around_target(h_1,tm1).second, if (!nm_features_map_1.non_manifold_edges.empty())
std::back_inserter(all_edges)); {
std::size_t vid1 = get(nm_features_map_1.v_nm_id, target(h_1, tm1));
if (vid1 != NM_NID)
{
for (vertex_descriptor vd : nm_features_map_1.non_manifold_vertices[vid1])
{
std::copy(halfedges_around_target(vd,tm1).first,
halfedges_around_target(vd,tm1).second,
std::back_inserter(all_edges));
}
if (all_edges.front()!=h_1)
{
// restore expected property
typename std::vector<halfedge_descriptor>::iterator pos =
std::find(all_edges.begin(), all_edges.end(), h_1);
CGAL_assertion(pos!=all_edges.end());
std::swap(*pos, all_edges.front());
}
}
else
std::copy(halfedges_around_target(h_1,tm1).first,
halfedges_around_target(h_1,tm1).second,
std::back_inserter(all_edges));
}
else
std::copy(halfedges_around_target(h_1,tm1).first,
halfedges_around_target(h_1,tm1).second,
std::back_inserter(all_edges));
}
else{ else{
if ( std::get<2>(res) ) // is edge source in triangle plane if ( std::get<2>(res) ) // is edge source in triangle plane
{
if (!nm_features_map_1.non_manifold_edges.empty())
{
std::size_t vid1 = get(nm_features_map_1.v_nm_id, source(h_1, tm1));
if (vid1 != NM_NID)
{
for (vertex_descriptor vd : nm_features_map_1.non_manifold_vertices[vid1])
{
std::copy(halfedges_around_source(vd,tm1).first,
halfedges_around_source(vd,tm1).second,
std::back_inserter(all_edges));
}
if (all_edges.front()!=h_1)
{
// restore expected property
typename std::vector<halfedge_descriptor>::iterator pos =
std::find(all_edges.begin(), all_edges.end(), h_1);
CGAL_assertion(pos!=all_edges.end());
std::swap(*pos, all_edges.front());
}
}
else
std::copy(halfedges_around_source(h_1,tm1).first,
halfedges_around_source(h_1,tm1).second,
std::back_inserter(all_edges));
}
else
std::copy(halfedges_around_source(h_1,tm1).first, std::copy(halfedges_around_source(h_1,tm1).first,
halfedges_around_source(h_1,tm1).second, halfedges_around_source(h_1,tm1).second,
std::back_inserter(all_edges)); std::back_inserter(all_edges));
}
else else
{
all_edges.push_back(h_1); all_edges.push_back(h_1);
edge_descriptor e_1 = edge(h_1, tm1);
if (!nm_features_map_1.non_manifold_edges.empty())
{
std::size_t eid1 = get(nm_features_map_1.e_nm_id, e_1);
if (eid1 != NM_NID)
{
CGAL_assertion( nm_features_map_1.non_manifold_edges[eid1][0]==e_1 );
for (std::size_t k=1;
k<nm_features_map_1.non_manifold_edges[eid1].size();
++k)
{
edge_descriptor e_1b = nm_features_map_1.non_manifold_edges[eid1][k];
// note that the orientation of the halfedge pushed back is
// not relevant for how it is used in the following
all_edges.push_back(halfedge(e_1b, tm1));
}
}
}
}
} }
CGAL_precondition(all_edges[0]==h_1 || all_edges[0]==opposite(h_1,tm1)); CGAL_precondition(all_edges[0]==h_1 || all_edges[0]==opposite(h_1,tm1));
// #ifdef USE_DETECTION_MULTIPLE_DEFINED_EDGES // #ifdef USE_DETECTION_MULTIPLE_DEFINED_EDGES
@ -726,7 +943,9 @@ class Intersection_of_triangle_meshes
add_intersection_point_to_face_and_all_edge_incident_faces(f_2,*it_edge,tm2,tm1,node_id); add_intersection_point_to_face_and_all_edge_incident_faces(f_2,*it_edge,tm2,tm1,node_id);
//erase face from the list to test intersection with it_edge //erase face from the list to test intersection with it_edge
if ( it_edge==all_edges.begin() ) if ( it_edge==all_edges.begin() )
{
fset.erase(fset.begin()); fset.erase(fset.begin());
}
else else
{ {
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(*it_edge,tm1)); typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(*it_edge,tm1));
@ -742,15 +961,35 @@ class Intersection_of_triangle_meshes
Node_id node_id=++current_node; Node_id node_id=++current_node;
add_new_node(h_1,f_2,tm1,tm2,vpm1,vpm2,res); add_new_node(h_1,f_2,tm1,tm2,vpm1,vpm2,res);
halfedge_descriptor h_2=std::get<1>(res); halfedge_descriptor h_2=std::get<1>(res);
std::size_t eid2 = nm_features_map_2.non_manifold_edges.empty()
? NM_NID
: get(nm_features_map_2.e_nm_id, edge(h_2, tm2));
if (eid2!=NM_NID)
h_2 = halfedge(nm_features_map_2.non_manifold_edges[eid2].front(), tm2);
visitor.new_node_added(node_id,ON_EDGE,h_1,h_2,tm1,tm2,std::get<3>(res),std::get<2>(res)); visitor.new_node_added(node_id,ON_EDGE,h_1,h_2,tm1,tm2,std::get<3>(res),std::get<2>(res));
for (;it_edge!=all_edges.end();++it_edge){ for (;it_edge!=all_edges.end();++it_edge){
if ( it_edge!=all_edges.begin() ){ if ( it_edge!=all_edges.begin() ){
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(*it_edge,tm1)); typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(*it_edge,tm1));
Face_set* fset_bis = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr; Face_set* fset_bis = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr;
cip_handle_case_edge(node_id,fset_bis,*it_edge,h_2,tm1,tm2); if( eid2 == NM_NID )
cip_handle_case_edge(node_id,fset_bis,*it_edge,h_2,tm1,tm2);
else
{
for (edge_descriptor e2 : nm_features_map_2.non_manifold_edges[eid2])
cip_handle_case_edge(node_id,fset_bis,*it_edge,halfedge(e2, tm2),tm1,tm2);
}
} }
else else
cip_handle_case_edge(node_id,&fset,*it_edge,h_2,tm1,tm2); {
if( eid2 == NM_NID )
cip_handle_case_edge(node_id,&fset,*it_edge,h_2,tm1,tm2);
else
for (edge_descriptor e2 : nm_features_map_2.non_manifold_edges[eid2])
cip_handle_case_edge(node_id,&fset,*it_edge,halfedge(e2, tm2),tm1,tm2);
}
} }
} // end case ON_EDGE } // end case ON_EDGE
break; break;
@ -762,14 +1001,27 @@ class Intersection_of_triangle_meshes
nodes.add_new_node(get(vpm2, target(h_2,tm2))); //we use the original vertex to create the node nodes.add_new_node(get(vpm2, target(h_2,tm2))); //we use the original vertex to create the node
//before it was ON_FACE but do not remember why, probably a bug... //before it was ON_FACE but do not remember why, probably a bug...
visitor.new_node_added(node_id,ON_VERTEX,h_1,h_2,tm1,tm2,std::get<3>(res),std::get<2>(res)); visitor.new_node_added(node_id,ON_VERTEX,h_1,h_2,tm1,tm2,std::get<3>(res),std::get<2>(res));
std::size_t vid2 = nm_features_map_2.non_manifold_vertices.empty()
? NM_NID
: get(nm_features_map_2.v_nm_id, target(h_2, tm2));
for (;it_edge!=all_edges.end();++it_edge){ for (;it_edge!=all_edges.end();++it_edge){
if ( it_edge!=all_edges.begin() ){ if ( it_edge!=all_edges.begin() ){
typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(*it_edge,tm1)); typename Edge_to_faces::iterator it_ets=tm1_edge_to_tm2_faces.find(edge(*it_edge,tm1));
Face_set* fset_bis = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr; Face_set* fset_bis = (it_ets!=tm1_edge_to_tm2_faces.end())?&(it_ets->second):nullptr;
cip_handle_case_vertex(node_id,fset_bis,*it_edge,h_2,tm1,tm2); if( vid2 == NM_NID )
cip_handle_case_vertex(node_id,fset_bis,*it_edge,h_2,tm1,tm2);
else
for (vertex_descriptor vd2 : nm_features_map_2.non_manifold_vertices[vid2])
cip_handle_case_vertex(node_id,fset_bis,*it_edge,halfedge(vd2, tm2),tm1,tm2);
} }
else else
cip_handle_case_vertex(node_id,&fset,*it_edge,h_2,tm1,tm2); if( vid2 == NM_NID )
cip_handle_case_vertex(node_id,&fset,*it_edge,h_2,tm1,tm2);
else
for (vertex_descriptor vd2 : nm_features_map_2.non_manifold_vertices[vid2])
cip_handle_case_vertex(node_id,&fset,*it_edge,halfedge(vd2, tm2),tm1,tm2);
} }
} // end case ON_VERTEX } // end case ON_VERTEX
break; break;
@ -1098,7 +1350,7 @@ class Intersection_of_triangle_meshes
} }
else{ else{
CGAL_assertion(segment.size()==1); CGAL_assertion(segment.size()==1);
isolated_point_seen=true; isolated_point_seen=true; // NOT TRUE CAN BE END POINT OF POLYLINE FALLING ONTO AN INPUT EDGE
} }
} }
@ -1284,6 +1536,20 @@ public:
CGAL_assertion_code( doing_autorefinement=true; ) CGAL_assertion_code( doing_autorefinement=true; )
} }
// setting maps of non manifold features
void set_non_manifold_feature_map_1(internal_np::Param_not_found){}
void set_non_manifold_feature_map_2(internal_np::Param_not_found){}
void set_non_manifold_feature_map_1(const Non_manifold_feature_map<TriangleMesh>& m)
{
non_manifold_feature_map_1=m;
visitor.set_non_manifold_feature_map(nodes.tm1, non_manifold_feature_map_1);
}
void set_non_manifold_feature_map_2(const Non_manifold_feature_map<TriangleMesh>& m)
{
non_manifold_feature_map_2=m;
visitor.set_non_manifold_feature_map(nodes.tm2, non_manifold_feature_map_2);
}
template <class OutputIterator> template <class OutputIterator>
OutputIterator operator()(OutputIterator output, OutputIterator operator()(OutputIterator output,
bool throw_on_self_intersection, bool throw_on_self_intersection,
@ -1296,18 +1562,18 @@ public:
const VertexPointMap& vpm1=nodes.vpm1; const VertexPointMap& vpm1=nodes.vpm1;
const VertexPointMap& vpm2=nodes.vpm2; const VertexPointMap& vpm2=nodes.vpm2;
filter_intersections(tm1, tm2, vpm1, vpm2, throw_on_self_intersection); filter_intersections(tm1, tm2, vpm1, vpm2, non_manifold_feature_map_2, throw_on_self_intersection);
filter_intersections(tm2, tm1, vpm2, vpm1, throw_on_self_intersection); filter_intersections(tm2, tm1, vpm2, vpm1, non_manifold_feature_map_1, throw_on_self_intersection);
Node_id current_node((std::numeric_limits<Node_id>::max)()); Node_id current_node((std::numeric_limits<Node_id>::max)());
CGAL_assertion(current_node+1==0); CGAL_assertion(current_node+1==0);
// TODO: handle non-manifold edges in coplanar
#ifndef DO_NOT_HANDLE_COPLANAR_FACES #ifndef DO_NOT_HANDLE_COPLANAR_FACES
//first handle coplanar triangles //first handle coplanar triangles
if (&tm1<&tm2) if (&tm1<&tm2)
compute_intersection_of_coplanar_faces(current_node, tm1, tm2, vpm1, vpm2); compute_intersection_of_coplanar_faces(current_node, tm1, tm2, vpm1, vpm2, non_manifold_feature_map_1, non_manifold_feature_map_2);
else else
compute_intersection_of_coplanar_faces(current_node, tm2, tm1, vpm2, vpm1); compute_intersection_of_coplanar_faces(current_node, tm2, tm1, vpm2, vpm1, non_manifold_feature_map_2, non_manifold_feature_map_1);
visitor.set_number_of_intersection_points_from_coplanar_faces(current_node+1); visitor.set_number_of_intersection_points_from_coplanar_faces(current_node+1);
if (!coplanar_faces.empty()) if (!coplanar_faces.empty())
visitor.input_have_coplanar_faces(); visitor.input_have_coplanar_faces();
@ -1322,8 +1588,11 @@ public:
? stm_edge_to_ltm_faces ? stm_edge_to_ltm_faces
: ltm_edge_to_stm_faces; : ltm_edge_to_stm_faces;
compute_intersection_points(tm1_edge_to_tm2_faces, tm1, tm2, vpm1, vpm2, current_node); compute_intersection_points(tm1_edge_to_tm2_faces, tm1, tm2, vpm1, vpm2, non_manifold_feature_map_1, non_manifold_feature_map_2, current_node);
compute_intersection_points(tm2_edge_to_tm1_faces, tm2, tm1, vpm2, vpm1, current_node); compute_intersection_points(tm2_edge_to_tm1_faces, tm2, tm1, vpm2, vpm1, non_manifold_feature_map_2, non_manifold_feature_map_1, current_node);
nodes.check_no_duplicates();
if (!build_polylines){ if (!build_polylines){
visitor.finalize(nodes,tm1,tm2,vpm1,vpm2); visitor.finalize(nodes,tm1,tm2,vpm1,vpm2);
return output; return output;
@ -1369,7 +1638,7 @@ public:
CGAL_assertion(current_node+1==0); CGAL_assertion(current_node+1==0);
//first handle coplanar triangles //first handle coplanar triangles
compute_intersection_of_coplanar_faces(current_node, tm, tm, vpm, vpm); compute_intersection_of_coplanar_faces(current_node, tm, tm, vpm, vpm, non_manifold_feature_map_1, non_manifold_feature_map_1);
if (!coplanar_faces.empty()) if (!coplanar_faces.empty())
visitor.input_have_coplanar_faces(); visitor.input_have_coplanar_faces();
@ -1377,7 +1646,7 @@ public:
//compute intersection points of segments and triangles. //compute intersection points of segments and triangles.
//build the nodes of the graph and connectivity infos //build the nodes of the graph and connectivity infos
compute_intersection_points(stm_edge_to_ltm_faces, tm, tm, vpm, vpm, current_node); compute_intersection_points(stm_edge_to_ltm_faces, tm, tm, vpm, vpm, non_manifold_feature_map_1, non_manifold_feature_map_1, current_node);
if (!build_polylines){ if (!build_polylines){
visitor.finalize(nodes,tm,tm,vpm,vpm); visitor.finalize(nodes,tm,tm,vpm,vpm);

View File

@ -125,7 +125,13 @@ public:
} }
void all_nodes_created(){} void all_nodes_created(){}
void finalize() {} template <class Mesh_to_map_node>
void finalize(const Mesh_to_map_node&) {}
void check_no_duplicates()
{
CGAL_assertion(nodes.size() == std::set<Point_3>(nodes.begin(), nodes.end()).size());
}
}; // end specialization }; // end specialization
// Intersection_nodes<Polyhedron,Kernel,No_predicates_on_constructions,false> // Intersection_nodes<Polyhedron,Kernel,No_predicates_on_constructions,false>
@ -158,7 +164,6 @@ private:
Exact_to_double exact_to_double; Exact_to_double exact_to_double;
Exact_kernel ek; Exact_kernel ek;
Exact_kernel::Intersect_3 exact_intersection; Exact_kernel::Intersect_3 exact_intersection;
std::vector<vertex_descriptor> tm1_vertices, tm2_vertices;
const bool doing_autorefinement; const bool doing_autorefinement;
public: public:
const TriangleMesh &tm1, &tm2; const TriangleMesh &tm1, &tm2;
@ -268,39 +273,45 @@ public:
} }
void all_nodes_created() void all_nodes_created()
{}
void call_put(const VertexPointMap& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh&)
{ {
tm1_vertices.resize(enodes.size(), GT::null_vertex()); put(vpm, vd, exact_to_double(enodes[i])); // Note this call is useless and only useful to see something in debug for intermediate results
tm2_vertices.resize(enodes.size(), GT::null_vertex());
} }
void call_put(const VertexPointMap& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh& tm) template <class Node_id_to_vertex>
void finalize(const std::map<const TriangleMesh*, Node_id_to_vertex>& mesh_to_node_id_to_vertex)
{ {
put(vpm, vd, exact_to_double(enodes[i])); if (!doing_autorefinement)
if (&tm1==&tm)
{ {
if ( tm1_vertices[i] == GT::null_vertex() ) const Node_id_to_vertex& tm1_vertices = mesh_to_node_id_to_vertex.find(&tm1)->second;
const Node_id_to_vertex& tm2_vertices = mesh_to_node_id_to_vertex.find(&tm2)->second;
for (std::size_t i=0, e=enodes.size(); i!=e; ++i)
{ {
tm1_vertices[i] = vd; Point_3 pt = exact_to_double(enodes[i]);
return; if ( tm1_vertices[i] != GT::null_vertex() )
put(vpm1, tm1_vertices[i], pt);
if ( tm2_vertices[i] != GT::null_vertex() )
put(vpm2, tm2_vertices[i], pt);
}
}
else{
const Node_id_to_vertex& tm1_vertices = mesh_to_node_id_to_vertex.find(&tm1)->second;
for (std::size_t i=0, e=enodes.size(); i!=e; ++i)
{
Point_3 pt = exact_to_double(enodes[i]);
if ( tm1_vertices[i] != GT::null_vertex() )
put(vpm1, tm1_vertices[i], pt);
} }
if (doing_autorefinement)
tm2_vertices[i] = vd;
} }
else
tm2_vertices[i] = vd;
} }
void finalize() void check_no_duplicates()
{ {
for (std::size_t i=0, e=enodes.size(); i!=e; ++i) CGAL_assertion(enodes.size() == std::set<typename Exact_kernel::Point_3>(enodes.begin(), enodes.end()).size());
{
Point_3 pt = exact_to_double(enodes[i]);
if ( tm1_vertices[i] != GT::null_vertex() )
put(vpm1, tm1_vertices[i], pt);
if ( tm2_vertices[i] != GT::null_vertex() )
put(vpm2, tm2_vertices[i], pt);
}
} }
}; // end specialization }; // end specialization
// Intersection_nodes<Polyhedron,Kernel,Predicates_on_constructions,false> // Intersection_nodes<Polyhedron,Kernel,Predicates_on_constructions,false>
@ -413,9 +424,16 @@ public:
} }
void all_nodes_created(){} void all_nodes_created(){}
void finalize() {}
template <class Node_id_to_vertex>
void finalize(const std::map<const TriangleMesh*, Node_id_to_vertex>&)
{}
void check_no_duplicates()
{
CGAL_assertion(nodes.size() == std::set<Point_3>(nodes.begin(), nodes.end()).size());
}
}; // end specialization }; // end specialization

View File

@ -1774,6 +1774,11 @@ surface_intersection(const TriangleMesh& tm1,
Corefinement::Intersection_of_triangle_meshes<TriangleMesh,Vpm> Corefinement::Intersection_of_triangle_meshes<TriangleMesh,Vpm>
functor(tm1, tm2, vpm1, vpm2); functor(tm1, tm2, vpm1, vpm2);
// Fill non-manifold feature maps if provided
functor.set_non_manifold_feature_map_1(parameters::get_parameter(np1, internal_np::non_manifold_feature_map));
functor.set_non_manifold_feature_map_2(parameters::get_parameter(np2, internal_np::non_manifold_feature_map));
return functor(polyline_output, throw_on_self_intersection, true); return functor(polyline_output, throw_on_self_intersection, true);
} }