Merge branch 'master' into Polygon_mesh_processing-triangulate_hole_with_cdt2-danston

This commit is contained in:
Dmitry Anisimov
2020-12-04 18:32:25 +01:00
committed by GitHub
970 changed files with 53100 additions and 26360 deletions
@@ -33,7 +33,14 @@
#include <boost/property_map/property_map.hpp>
#include <array>
#include <algorithm>
#include <map>
#include <set>
#include <type_traits>
#include <unordered_map>
#include <utility>
#include <vector>
namespace CGAL{
namespace Polygon_mesh_processing {
@@ -420,10 +427,9 @@ generic_clip_impl(
NamedParameters1>::type Vpm;
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters2>::type Vpm2;
CGAL_USE_TYPE(Vpm2);
CGAL_assertion_code(
static const bool same_vpm = (boost::is_same<Vpm,Vpm2>::value); )
CGAL_static_assertion(same_vpm);
CGAL_static_assertion((std::is_same<typename boost::property_traits<Vpm>::value_type,
typename boost::property_traits<Vpm>::value_type>::value));
Vpm vpm1 = choose_parameter(get_parameter(np1, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm1));
@@ -477,17 +483,17 @@ generic_clip_impl(
// surface intersection algorithm call
typedef Corefinement::Generic_clip_output_builder<TriangleMesh,
Vpm,
Vpm, Vpm2,
Algo_ecm1,
FaceIndexMap1,
Default> Ob;
typedef Corefinement::Surface_intersection_visitor_for_corefinement<
TriangleMesh, Vpm, Ob, Ecm_in, User_visitor> Algo_visitor;
TriangleMesh, Vpm, Vpm2, Ob, Ecm_in, User_visitor> Algo_visitor;
Ecm_in ecm_in(tm1,tm2,ecm1,ecm2);
Ob ob(tm1, tm2, vpm1, vpm2, algo_ecm1, fid_map1, use_compact_clipper);
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor >
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Vpm2, Algo_visitor >
functor(tm1, tm2, vpm1, vpm2, Algo_visitor(uv,ob,ecm_in,&tm2));
functor(CGAL::Emptyset_iterator(), false, true);
}
@@ -663,6 +669,14 @@ clip(TriangleMesh& tm,
* \cgalParamDefault{`true`}
* \cgalParamNEnd
*
* \cgalParamNBegin{allow_self_intersections}
* \cgalParamDescription{If `true`, self-intersections are accepted for `tm`.}
* \cgalParamType{Boolean}
* \cgalParamDefault{`false`}
* \cgalParamExtra{If this option is set to `true`, `tm` is no longer required to be without self-intersection.
* Setting this option to `true` will automatically set `throw_on_self_intersection` to `false`
* and `clip_volume` to `false`.}
* \cgalParamNEnd
* \cgalNamedParamsEnd
*
* @return `true` if the output surface mesh is manifold.
@@ -678,6 +692,10 @@ bool clip(TriangleMesh& tm,
#endif
const NamedParameters& np)
{
using parameters::get_parameter;
using parameters::choose_parameter;
namespace PMP = CGAL::Polygon_mesh_processing;
namespace params = PMP::parameters;
if(boost::begin(faces(tm))==boost::end(faces(tm))) return true;
CGAL::Bbox_3 bbox = ::CGAL::Polygon_mesh_processing::bbox(tm);
@@ -701,7 +719,8 @@ bool clip(TriangleMesh& tm,
break;
}
return clip(tm, clipper, np, parameters::all_default());
const bool do_not_modify = choose_parameter(get_parameter(np, internal_np::allow_self_intersections), false);
return clip(tm, clipper, np, params::do_not_modify(do_not_modify));
}
/**
@@ -752,6 +771,15 @@ bool clip(TriangleMesh& tm,
* \cgalParamType{Boolean}
* \cgalParamDefault{`true`}
* \cgalParamNEnd
*
* \cgalParamNBegin{allow_self_intersections}
* \cgalParamDescription{If `true`, self-intersections are accepted for `tm`.}
* \cgalParamType{Boolean}
* \cgalParamDefault{`false`}
* \cgalParamExtra{If this option is set to `true`, `tm` is no longer required to be without self-intersection.
* Setting this option to `true` will automatically set `throw_on_self_intersection` to `false`
* and `clip_volume` to `false`.}
* \cgalParamNEnd
* \cgalNamedParamsEnd
*
* @return `true` if the output surface mesh is manifold.
@@ -767,6 +795,11 @@ bool clip(TriangleMesh& tm,
#endif
const NamedParameters& np)
{
using parameters::get_parameter;
using parameters::choose_parameter;
namespace PMP = CGAL::Polygon_mesh_processing;
namespace params = PMP::parameters;
if(boost::begin(faces(tm))==boost::end(faces(tm))) return true;
TriangleMesh clipper;
@@ -775,7 +808,8 @@ bool clip(TriangleMesh& tm,
clipper);
triangulate_faces(clipper);
return clip(tm, clipper, np, parameters::all_default());
const bool do_not_modify = choose_parameter(get_parameter(np, internal_np::allow_self_intersections), false);
return clip(tm, clipper, np, params::do_not_modify(do_not_modify));
}
/*!
@@ -818,6 +852,8 @@ bool clip(TriangleMesh& tm,
* will be thrown if at least one self-intersection is found.}
* \cgalParamType{Boolean}
* \cgalParamDefault{`false`}
* \cgalParamNEnd
*
* \cgalParamNBegin{do_not_modify}
* \cgalParamDescription{(`np_s` only) if `true`, `splitter` will not be modified.}
* \cgalParamType{Boolean}
@@ -908,6 +944,13 @@ void split(TriangleMesh& tm,
* \cgalParamDefault{`false`}
* \cgalParamNEnd
*
* \cgalParamNBegin{allow_self_intersections}
* \cgalParamDescription{If `true`, self-intersections are accepted for `tm`.}
* \cgalParamType{Boolean}
* \cgalParamDefault{`false`}
* \cgalParamExtra{If this option is set to `true`, `tm` is no longer required to be without self-intersection.
* Setting this option to `true` will automatically set `throw_on_self_intersection` to `false`.}
* \cgalParamNEnd
* \cgalNamedParamsEnd
*/
template <class TriangleMesh,
@@ -937,7 +980,11 @@ void split(TriangleMesh& tm,
CGAL::Oriented_side os = PMP::internal::clip_to_bbox(plane, bbox, splitter, PMP::parameters::all_default());
if(os == CGAL::ON_ORIENTED_BOUNDARY)
return split(tm, splitter, np, params::all_default());
{
const bool do_not_modify = choose_parameter(get_parameter(np, internal_np::allow_self_intersections), false);
return split(tm, splitter, np, params::do_not_modify(do_not_modify));
}
//else nothing to do, no intersection.
}
@@ -998,6 +1045,15 @@ void split(TriangleMesh& tm,
* \cgalParamDefault{`true`}
* \cgalParamNEnd
*
* * \cgalParamNBegin{allow_self_intersections}
* \cgalParamDescription{If `true`, self-intersections are accepted for `tm`.}
* \cgalParamType{Boolean}
* \cgalParamDefault{`false`}
* \cgalParamExtra{If this option is set to `true`, `tm` is no longer required to be without self-intersection.
* Setting this option to `true` will automatically set `throw_on_self_intersection` to `false`
* and `clip_volume` to `false`.}
* \cgalParamNEnd
*
* \cgalNamedParamsEnd
*/
template <class TriangleMesh,
@@ -1010,14 +1066,18 @@ void split(TriangleMesh& tm,
#endif
const NamedParameters& np)
{
using parameters::get_parameter;
using parameters::choose_parameter;
namespace PMP = CGAL::Polygon_mesh_processing;
namespace params = PMP::parameters;
TriangleMesh splitter;
make_hexahedron(iso_cuboid[0], iso_cuboid[1], iso_cuboid[2], iso_cuboid[3],
iso_cuboid[4], iso_cuboid[5], iso_cuboid[6], iso_cuboid[7],
splitter);
triangulate_faces(splitter);
split(tm, splitter, np, parameters::all_default());
const bool do_not_modify = choose_parameter(get_parameter(np, internal_np::allow_self_intersections), false);
return split(tm, splitter, np, params::do_not_modify(do_not_modify));
}
/// \cond SKIP_IN_MANUAL
@@ -201,22 +201,19 @@ corefine_and_compute_boolean_operations(
// Vertex point maps
//for input meshes
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters1>::type Vpm;
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters2>::type Vpm2;
CGAL_USE_TYPE(Vpm2);
CGAL_assertion_code(
static const bool same_vpm = (boost::is_same<Vpm,Vpm2>::value); )
CGAL_static_assertion(same_vpm);
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters1>::type VPM1;
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters2>::type VPM2;
Vpm vpm1 = choose_parameter(get_parameter(np1, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm1));
CGAL_static_assertion((std::is_same<typename boost::property_traits<VPM1>::value_type,
typename boost::property_traits<VPM2>::value_type>::value));
Vpm vpm2 = choose_parameter(get_parameter(np2, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm2));
VPM1 vpm1 = choose_parameter(get_parameter(np1, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm1));
typedef typename boost::property_traits<Vpm>::value_type Point_3;
VPM2 vpm2 = choose_parameter(get_parameter(np2, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm2));
typedef typename boost::property_traits<VPM1>::value_type Point_3;
// for output meshes: here we have to use a trick so that if for a specific output
// that is not requested, the default vpm does not have the same value type as the
@@ -227,24 +224,24 @@ corefine_and_compute_boolean_operations(
Corefinement::TweakedGetVertexPointMap<Point_3, NamedParametersOut1, TriangleMesh>,
Corefinement::TweakedGetVertexPointMap<Point_3, NamedParametersOut2, TriangleMesh>,
Corefinement::TweakedGetVertexPointMap<Point_3, NamedParametersOut3, TriangleMesh>
> Vpm_out_tuple_helper;
> VPM_out_tuple_helper;
typedef std::tuple<
boost::optional< typename std::tuple_element<0, Vpm_out_tuple_helper>::type::type >,
boost::optional< typename std::tuple_element<1, Vpm_out_tuple_helper>::type::type >,
boost::optional< typename std::tuple_element<2, Vpm_out_tuple_helper>::type::type >,
boost::optional< typename std::tuple_element<3, Vpm_out_tuple_helper>::type::type >
> Vpm_out_tuple;
boost::optional< typename std::tuple_element<0, VPM_out_tuple_helper>::type::type >,
boost::optional< typename std::tuple_element<1, VPM_out_tuple_helper>::type::type >,
boost::optional< typename std::tuple_element<2, VPM_out_tuple_helper>::type::type >,
boost::optional< typename std::tuple_element<3, VPM_out_tuple_helper>::type::type >
> VPM_out_tuple;
Vpm_out_tuple vpm_out_tuple(
VPM_out_tuple vpm_out_tuple(
Corefinement::get_vpm<Point_3>(std::get<0>(nps_out), output[0],
typename std::tuple_element<0, Vpm_out_tuple_helper>::type::Use_default_tag()),
typename std::tuple_element<0, VPM_out_tuple_helper>::type::Use_default_tag()),
Corefinement::get_vpm<Point_3>(std::get<1>(nps_out), output[1],
typename std::tuple_element<1, Vpm_out_tuple_helper>::type::Use_default_tag()),
typename std::tuple_element<1, VPM_out_tuple_helper>::type::Use_default_tag()),
Corefinement::get_vpm<Point_3>(std::get<2>(nps_out), output[2],
typename std::tuple_element<2, Vpm_out_tuple_helper>::type::Use_default_tag()),
typename std::tuple_element<2, VPM_out_tuple_helper>::type::Use_default_tag()),
Corefinement::get_vpm<Point_3>(std::get<3>(nps_out), output[3],
typename std::tuple_element<3, Vpm_out_tuple_helper>::type::Use_default_tag())
typename std::tuple_element<3, VPM_out_tuple_helper>::type::Use_default_tag())
);
if (&tm1==&tm2)
@@ -375,8 +372,9 @@ corefine_and_compute_boolean_operations(
// surface intersection algorithm call
typedef Corefinement::Face_graph_output_builder<TriangleMesh,
Vpm,
Vpm_out_tuple,
VPM1,
VPM2,
VPM_out_tuple,
FaceIndexMap1,
FaceIndexMap2,
Default,
@@ -385,7 +383,8 @@ corefine_and_compute_boolean_operations(
User_visitor> Ob;
typedef Corefinement::Surface_intersection_visitor_for_corefinement<
TriangleMesh, Vpm, Ob, Ecm_in, User_visitor> Algo_visitor;
TriangleMesh, VPM1, VPM2, Ob, Ecm_in, User_visitor> Algo_visitor;
Ecm_in ecm_in(tm1,tm2,ecm1,ecm2);
Edge_mark_map_tuple ecms_out(ecm_out_0, ecm_out_1, ecm_out_2, ecm_out_3);
Ob ob(tm1, tm2, vpm1, vpm2, fid_map1, fid_map2, ecm_in, vpm_out_tuple, ecms_out, uv, output);
@@ -402,7 +401,7 @@ corefine_and_compute_boolean_operations(
ob.setup_for_clipping_a_surface(use_compact_clipper);
}
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor >
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, VPM1, VPM2, Algo_visitor >
functor(tm1, tm2, vpm1, vpm2, Algo_visitor(uv,ob,ecm_in));
functor(CGAL::Emptyset_iterator(), throw_on_self_intersection, true);
@@ -735,20 +734,17 @@ corefine( TriangleMesh& tm1,
choose_parameter(get_parameter(np1, internal_np::throw_on_self_intersection), false);
// Vertex point maps
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters1>::type Vpm;
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters2>::type Vpm2;
CGAL_USE_TYPE(Vpm2);
CGAL_assertion_code(
static const bool same_vpm = (boost::is_same<Vpm,Vpm2>::value);)
CGAL_static_assertion(same_vpm);
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters1>::type VPM1;
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters2>::type VPM2;
Vpm vpm1 = choose_parameter(get_parameter(np1, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm1));
CGAL_static_assertion((std::is_same<typename boost::property_traits<VPM1>::value_type,
typename boost::property_traits<VPM2>::value_type>::value));
Vpm vpm2 = choose_parameter(get_parameter(np2, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm2));
VPM1 vpm1 = choose_parameter(get_parameter(np1, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm1));
VPM2 vpm2 = choose_parameter(get_parameter(np2, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm2));
// Edge is-constrained maps
typedef typename internal_np::Lookup_named_param_def <
@@ -786,10 +782,11 @@ corefine( TriangleMesh& tm1,
// surface intersection algorithm call
typedef Corefinement::No_extra_output_from_corefinement<TriangleMesh> Ob;
typedef Corefinement::Surface_intersection_visitor_for_corefinement<
TriangleMesh, Vpm, Ob, Ecm, User_visitor> Algo_visitor;
TriangleMesh, VPM1, VPM2, Ob, Ecm, User_visitor> Algo_visitor;
Ob ob;
Ecm ecm(tm1,tm2,ecm1,ecm2);
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor>
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, VPM1, VPM2, Algo_visitor>
functor(tm1, tm2, vpm1, vpm2, Algo_visitor(uv,ob,ecm,const_mesh_ptr));
functor(CGAL::Emptyset_iterator(), throw_on_self_intersection, true);
}
@@ -852,9 +849,9 @@ autorefine( TriangleMesh& tm,
using parameters::get_parameter;
// Vertex point maps
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::type Vpm;
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::type VPM;
Vpm vpm = choose_parameter(get_parameter(np, internal_np::vertex_point),
VPM vpm = choose_parameter(get_parameter(np, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm));
// Edge is-constrained maps
@@ -877,10 +874,10 @@ autorefine( TriangleMesh& tm,
// surface intersection algorithm call
typedef Corefinement::No_extra_output_from_corefinement<TriangleMesh> Ob;
typedef Corefinement::Surface_intersection_visitor_for_corefinement<
TriangleMesh, Vpm, Ob, Ecm, User_visitor,true> Algo_visitor;
TriangleMesh, VPM, VPM, Ob, Ecm, User_visitor,true> Algo_visitor;
Ob ob;
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor>
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, VPM, VPM, Algo_visitor>
functor(tm, vpm, Algo_visitor(uv,ob,ecm) );
functor(CGAL::Emptyset_iterator(), true);
@@ -945,8 +942,8 @@ autorefine_and_remove_self_intersections( TriangleMesh& tm,
using parameters::get_parameter;
// Vertex point maps
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::type Vpm;
Vpm vpm = choose_parameter(get_parameter(np, internal_np::vertex_point),
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::type VPM;
VPM vpm = choose_parameter(get_parameter(np, internal_np::vertex_point),
get_property_map(boost::vertex_point, tm));
// Face index map
@@ -972,16 +969,16 @@ autorefine_and_remove_self_intersections( TriangleMesh& tm,
// surface intersection algorithm call
typedef Corefinement::Output_builder_for_autorefinement<TriangleMesh,
Vpm,
VPM,
Fid_map,
Ecm,
Default > Ob;
typedef Corefinement::Surface_intersection_visitor_for_corefinement<
TriangleMesh, Vpm, Ob, Ecm, User_visitor,true> Algo_visitor;
TriangleMesh, VPM, VPM, Ob, Ecm, User_visitor,true> Algo_visitor;
Ob ob(tm, vpm, fid_map, ecm);
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, Vpm, Algo_visitor>
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, VPM, VPM, Algo_visitor>
functor(tm, vpm, Algo_visitor(uv,ob,ecm) );
functor(CGAL::Emptyset_iterator(), true);
@@ -58,7 +58,8 @@ namespace PMP=Polygon_mesh_processing;
namespace params=PMP::parameters;
template <class TriangleMesh,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class VpmOutTuple,
class FaceIdMap1,
class FaceIdMap2,
@@ -72,8 +73,9 @@ class Face_graph_output_builder
typedef typename Default::Get<
Kernel_,
typename Kernel_traits<
typename boost::property_traits<VertexPointMap>::value_type
>::Kernel >::type Kernel;
typename boost::property_traits<VertexPointMap1>::value_type
>::Kernel >::type Kernel;
typedef typename Default::Get<EdgeMarkMapBind_,
Ecm_bind<TriangleMesh, No_mark<TriangleMesh> >
>::type EdgeMarkMapBind;
@@ -111,8 +113,8 @@ class Face_graph_output_builder
//Data members
TriangleMesh &tm1, &tm2;
// property maps of input meshes
const VertexPointMap vpm1;
const VertexPointMap vpm2;
const VertexPointMap1& vpm1;
const VertexPointMap2& vpm2;
FaceIdMap1 fids1;
FaceIdMap2 fids2;
EdgeMarkMapBind& marks_on_input_edges;
@@ -346,8 +348,8 @@ public:
Face_graph_output_builder(TriangleMesh& tm1,
TriangleMesh& tm2,
const VertexPointMap vpm1,
const VertexPointMap vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
FaceIdMap1 fids1,
FaceIdMap2 fids2,
EdgeMarkMapBind& marks_on_input_edges,
@@ -1098,10 +1100,6 @@ public:
if (!is_tm2_closed)
patch_status_not_set_tm1.reset();
typedef Side_of_triangle_mesh<TriangleMesh,
Kernel,
VertexPointMap> Inside_poly_test;
#ifdef CGAL_COREFINEMENT_POLYHEDRA_DEBUG
#warning stop using next_marked_halfedge_around_target_vertex and create lists of halfedges instead?
#endif
@@ -1111,7 +1109,7 @@ public:
CGAL::Bounded_side in_tm2 = is_tm2_inside_out
? ON_UNBOUNDED_SIDE : ON_BOUNDED_SIDE;
Inside_poly_test inside_tm2(tm2, vpm2);
Side_of_triangle_mesh<TriangleMesh, Kernel, VertexPointMap2> inside_tm2(tm2, vpm2);
for(face_descriptor f : faces(tm1))
{
@@ -1173,7 +1171,7 @@ public:
CGAL::Bounded_side in_tm1 = is_tm1_inside_out
? ON_UNBOUNDED_SIDE : ON_BOUNDED_SIDE;
Inside_poly_test inside_tm1(tm1, vpm1);
Side_of_triangle_mesh<TriangleMesh, Kernel, VertexPointMap1> inside_tm1(tm1, vpm1);
for(face_descriptor f : faces(tm2))
{
const std::size_t f_id = get(fids2, f);
@@ -38,7 +38,8 @@ namespace PMP=Polygon_mesh_processing;
namespace params=PMP::parameters;
template <class TriangleMesh,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class Ecm1,
class FaceIdMap1,
class Kernel_=Default>
@@ -48,7 +49,7 @@ class Generic_clip_output_builder
typedef typename Default::Get<
Kernel_,
typename Kernel_traits<
typename boost::property_traits<VertexPointMap>::value_type
typename boost::property_traits<VertexPointMap2>::value_type
>::Kernel >::type Kernel;
// graph_traits typedefs
@@ -76,8 +77,8 @@ class Generic_clip_output_builder
//Data members
TriangleMesh &tm1, &tm2;
// property maps of input meshes
const VertexPointMap vpm1;
const VertexPointMap vpm2;
const VertexPointMap1 vpm1;
const VertexPointMap2 vpm2;
Ecm1 ecm1;
FaceIdMap1 fids1;
bool use_compact_clipper;
@@ -105,8 +106,8 @@ public:
Generic_clip_output_builder(TriangleMesh& tm1,
TriangleMesh& tm2,
const VertexPointMap vpm1,
const VertexPointMap vpm2,
const VertexPointMap1 vpm1,
const VertexPointMap2 vpm2,
const Ecm1& ecm1,
FaceIdMap1 fids1,
bool use_compact_clipper)
@@ -167,7 +168,7 @@ public:
typedef Side_of_triangle_mesh<TriangleMesh,
Kernel,
VertexPointMap> Inside_poly_test;
VertexPointMap2> Inside_poly_test;
CGAL::Bounded_side in_tm2 = is_tm2_inside_out
? ON_UNBOUNDED_SIDE : ON_BOUNDED_SIDE;
@@ -102,7 +102,8 @@ struct No_extra_output_from_corefinement
// A visitor for Intersection_of_triangle_meshes that can be used to corefine
// two meshes
template< class TriangleMesh,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class OutputBuilder_ = Default,
class EdgeMarkMapBind_ = Default,
class UserVisitor_ = Default,
@@ -141,8 +142,9 @@ private:
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 for the CDT
typedef typename Intersection_nodes<TriangleMesh,
VertexPointMap, Predicates_on_constructions_needed>::Exact_kernel EK;
typedef Intersection_nodes<TriangleMesh, VertexPointMap1, VertexPointMap2,
Predicates_on_constructions_needed> INodes;
typedef typename INodes::Exact_kernel EK;
typedef Triangulation_2_projection_traits_3<EK> CDT_traits;
typedef Triangulation_vertex_base_with_info_2<Node_id,CDT_traits> Vb;
typedef Constrained_triangulation_face_base_2<CDT_traits> Fb;
@@ -365,7 +367,8 @@ public:
}
}
if (tm1_ptr==const_mesh_ptr) return;
if (tm1_ptr==const_mesh_ptr)
return;
CGAL_assertion(!is_target_coplanar || !is_source_coplanar); //coplanar edge are not forwarded
@@ -406,16 +409,16 @@ public:
//sort node ids so that we can split the hedge
//consecutively
template <class Node_vector>
template <class VPM, class Node_vector>
void sort_vertices_along_hedge(std::vector<std::size_t>& node_ids,
halfedge_descriptor hedge,
const TriangleMesh& tm,
const VertexPointMap& vpm,
const VPM& vpm,
const Node_vector& nodes)
{
std::sort(node_ids.begin(),
node_ids.end(),
Less_along_a_halfedge<TriangleMesh,VertexPointMap,Node_vector>
Less_along_a_halfedge<TriangleMesh, VPM, Node_vector>
(hedge, tm, vpm, nodes)
);
}
@@ -490,6 +493,8 @@ public:
};
typedef boost::unordered_map<face_descriptor,Face_boundary> Face_boundaries;
//update the id of input mesh vertex that are also a node
void update_face_indices(
std::array<vertex_descriptor,3>& f_vertices,
@@ -592,22 +597,366 @@ public:
return vh;
}
void finalize(Intersection_nodes<TriangleMesh,
VertexPointMap, Predicates_on_constructions_needed>& nodes,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2)
template <class OnEdgeMapIterator, class VPM>
void split_halfedges(OnEdgeMapIterator it,
const VPM& vpm,
INodes& nodes,
std::map<TriangleMesh*, Face_boundaries>& mesh_to_face_boundaries)
{
TriangleMesh& tm=*it->first;
CGAL_assertion(&tm!=const_mesh_ptr);
On_edge_map& on_edge_map=it->second;
On_face_map& on_face_map=on_face[&tm];
Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm];
for(typename On_edge_map::iterator it2=on_edge_map.begin();
it2!=on_edge_map.end();
++it2)
{
//the edge to be split
halfedge_descriptor hedge=halfedge(it2->first,tm);
//indices of the nodes to be inserted
Node_ids& node_ids=it2->second;
CGAL_assertion( std::set<Node_id>(node_ids.begin(), node_ids.end())
.size()==node_ids.size() );
//sort nodes along the egde to allow consecutive splits
sort_vertices_along_hedge(node_ids,hedge,tm,vpm,nodes);
//save original face and nodes for face of hedge (1)
if ( !is_border(hedge,tm) ){
face_descriptor f=face(hedge,tm);
typename Face_boundaries::iterator it_face = face_boundaries.find(f);
if (it_face==face_boundaries.end())
it_face=face_boundaries.insert(std::make_pair(f,Face_boundary(hedge,tm))).first;
it_face->second.copy_node_ids(hedge,node_ids.begin(),node_ids.end());
}
//save original face and nodes for face of hedge->opposite (2)
typename Face_boundaries::iterator opposite_original_info=face_boundaries.end();
halfedge_descriptor hedge_opp = opposite(hedge,tm);
if ( !is_border(hedge_opp,tm) ){
face_descriptor f=face(hedge_opp,tm);
opposite_original_info=face_boundaries.find(f);
if (opposite_original_info==face_boundaries.end())
opposite_original_info=face_boundaries.insert(std::make_pair(f,Face_boundary(hedge_opp,tm))).first;
opposite_original_info->second.copy_node_ids(hedge_opp,node_ids.rbegin(),node_ids.rend());
}
typename Mesh_to_map_node::iterator it_map=mesh_to_node_id_to_vertex.find(&tm);
CGAL_assertion(it_map!=mesh_to_node_id_to_vertex.end());
//a map to identify the vertex in the polyhedron corresponding to an intersection point
Node_id_to_vertex& node_id_to_vertex=it_map->second;
CGAL_assertion_code(vertex_descriptor original_vertex=source(hedge,tm);)
//We need an edge incident to the source vertex of hedge. This is the first opposite edge created.
bool first=true;
halfedge_descriptor hedge_incident_to_src=Graph_traits::null_halfedge();
bool hedge_is_marked = call_get(marks_on_edges,tm,edge(hedge,tm));
//do split the edges
CGAL_assertion_code(vertex_descriptor expected_src=source(hedge,tm));
for(std::size_t node_id : node_ids)
{
halfedge_descriptor hnew = Euler::split_edge(hedge, tm);
CGAL_assertion(expected_src==source(hnew,tm));
vertex_descriptor vnew=target(hnew,tm);
// user_visitor.new_vertex_added(node_id, vnew, tm); // NODE_VISITOR_TAG
nodes.call_put(vpm, vnew, node_id, tm);
// register the new vertex in the output builder
output_builder.set_vertex_id(vnew, node_id, tm);
node_id_to_vertex[node_id]=vnew;
if (first){
first=false;
hedge_incident_to_src=next(opposite(hedge,tm),tm);
}
//update marker tags. If the edge was marked, then the resulting edges in the split must be marked
if ( hedge_is_marked )
call_put(marks_on_edges,tm,edge(hnew,tm),true);
CGAL_assertion_code(expected_src=vnew);
}
CGAL_assertion(target(hedge_incident_to_src,tm)==original_vertex);
CGAL_assertion(face(hedge_incident_to_src,tm)==face(hedge_opp,tm));
//save original face and nodes for face of hedge->opposite (2)
if ( !is_border(hedge_opp,tm) ){
CGAL_assertion(opposite_original_info!=face_boundaries.end());
opposite_original_info->second.update_original_halfedge(
hedge_opp,hedge_incident_to_src,tm);
}
//insert the two incident faces in on_face map so that they will be triangulated.
if (!is_border(hedge,tm)) on_face_map[face(hedge,tm)];
if (!is_border(hedge_opp,tm)) on_face_map[face(hedge_opp,tm)];
}
}
template <class OnFaceMapIterator, class VPM>
void triangulate_intersected_faces(OnFaceMapIterator it,
const VPM& vpm,
INodes& nodes,
std::map<TriangleMesh*, Face_boundaries>& mesh_to_face_boundaries)
{
TriangleMesh& tm=*it->first;
CGAL_assertion(&tm!=const_mesh_ptr);
On_face_map& on_face_map=it->second;
Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm];
Node_id_to_vertex& node_id_to_vertex=mesh_to_node_id_to_vertex[&tm];
Vertex_to_node_id& vertex_to_node_id=mesh_to_vertex_to_node_id[&tm];
const Node_id nb_nodes = nodes.size();
for (typename On_face_map::iterator it=on_face_map.begin();
it!=on_face_map.end();++it)
{
face_descriptor f = it->first; //the face to be triangulated
Node_ids& node_ids = it->second; // ids of nodes in the interior of f
typename Face_boundaries::iterator it_fb=face_boundaries.find(f);
std::map<Node_id,typename CDT::Vertex_handle> id_to_CDT_vh;
//associate an edge of the triangulation to a halfedge in a given polyhedron
std::map<std::pair<Node_id,Node_id>,halfedge_descriptor> edge_to_hedge;
// the vertices of f
std::array<vertex_descriptor,3> f_vertices;
// the node_id of an input vertex or a fake id (>=nb_nodes)
std::array<Node_id,3> f_indices = {{nb_nodes,nb_nodes+1,nb_nodes+2}};
if (it_fb!=face_boundaries.end()){ //the boundary of the triangle face was refined
f_vertices[0]=it_fb->second.vertices[0];
f_vertices[1]=it_fb->second.vertices[1];
f_vertices[2]=it_fb->second.vertices[2];
update_face_indices(f_vertices,f_indices,vertex_to_node_id);
if (doing_autorefinement)
it_fb->second.update_node_id_to_vertex_map(node_id_to_vertex, tm);
}
else{
CGAL_assertion( is_triangle(halfedge(f,tm),tm) );
halfedge_descriptor h0=halfedge(f,tm), h1=next(h0,tm), h2=next(h1,tm);
f_vertices[0]=target(h0,tm); //nb_nodes
f_vertices[1]=target(h1,tm); //nb_nodes+1
f_vertices[2]=target(h2,tm); //nb_nodes+2
update_face_indices(f_vertices,f_indices,vertex_to_node_id);
edge_to_hedge[std::make_pair( f_indices[2],f_indices[0] )] = h0;
edge_to_hedge[std::make_pair( f_indices[0],f_indices[1] )] = h1;
edge_to_hedge[std::make_pair( f_indices[1],f_indices[2] )] = h2;
}
typename EK::Point_3 p = nodes.to_exact(get(vpm,f_vertices[0])),
q = nodes.to_exact(get(vpm,f_vertices[1])),
r = nodes.to_exact(get(vpm,f_vertices[2]));
///TODO use a positive normal and remove all work around to guarantee that triangulation of coplanar patches are compatible
CDT_traits traits(typename EK::Construct_normal_3()(p,q,r));
CDT cdt(traits);
// insert triangle points
std::array<CDT_Vertex_handle,3> triangle_vertices;
//we can do this to_exact because these are supposed to be input points.
triangle_vertices[0]=cdt.insert_outside_affine_hull(p);
triangle_vertices[1]=cdt.insert_outside_affine_hull(q);
triangle_vertices[2]=cdt.tds().insert_dim_up(cdt.infinite_vertex(), false);
triangle_vertices[2]->set_point(r);
triangle_vertices[0]->info()=f_indices[0];
triangle_vertices[1]->info()=f_indices[1];
triangle_vertices[2]->info()=f_indices[2];
node_id_to_vertex[nb_nodes ]=f_vertices[0];
node_id_to_vertex[nb_nodes+1]=f_vertices[1];
node_id_to_vertex[nb_nodes+2]=f_vertices[2];
//if one of the triangle input vertex is also a node
for (int ik=0;ik<3;++ik){
if ( f_indices[ik]<nb_nodes )
{
id_to_CDT_vh.insert(
std::make_pair(f_indices[ik],triangle_vertices[ik]));
if (doing_autorefinement)
// update the current vertex in node_id_to_vertex
// to match the one of the face
node_id_to_vertex[f_indices[ik]]=f_vertices[ik];
}
}
//insert points on edges
if (it_fb!=face_boundaries.end()) //if f not a triangle?
{
// collect infinite faces incident to the initial triangle
typename CDT::Face_handle infinite_faces[3];
for (int i=0;i<3;++i)
{
int oi=-1;
CGAL_assertion_code(bool is_edge = )
cdt.is_edge(triangle_vertices[i], triangle_vertices[(i+1)%3], infinite_faces[i], oi);
CGAL_assertion(is_edge);
CGAL_assertion( cdt.is_infinite( infinite_faces[i]->vertex(oi) ) );
}
// In this loop, for each original edge of the triangle, we insert
// the constrained edges and we recover the halfedge_descriptor
// corresponding to these constrained (they are already in tm)
Face_boundary& f_boundary=it_fb->second;
for (int i=0;i<3;++i){
//handle case of halfedge starting at triangle_vertices[i]
// and ending at triangle_vertices[(i+1)%3]
const Node_ids& ids_on_edge=f_boundary.node_ids_array[i];
CDT_Vertex_handle previous=triangle_vertices[i];
Node_id prev_index=f_indices[i];// node-id of the mesh vertex
halfedge_descriptor hedge = next(f_boundary.halfedges[(i+2)%3],tm);
CGAL_assertion( source(hedge,tm)==f_boundary.vertices[i] );
if (!ids_on_edge.empty()){ //is there at least one node on this edge?
// fh must be an infinite face
// The points must be ordered from fh->vertex(cw(infinite_vertex)) to fh->vertex(ccw(infinite_vertex))
for(Node_id id : ids_on_edge)
{
CDT_Vertex_handle vh=insert_point_on_ch_edge(cdt,infinite_faces[i],nodes.exact_node(id));
vh->info()=id;
id_to_CDT_vh.insert(std::make_pair(id,vh));
edge_to_hedge[std::make_pair(prev_index,id)]=hedge;
previous=vh;
hedge=next(hedge,tm);
prev_index=id;
}
}
else{
CGAL_assertion_code(halfedge_descriptor hd=f_boundary.halfedges[i]);
CGAL_assertion( target(hd,tm) == f_boundary.vertices[(i+1)%3] );
CGAL_assertion( source(hd,tm) == f_boundary.vertices[ i ] );
}
CGAL_assertion(hedge==f_boundary.halfedges[i]);
edge_to_hedge[std::make_pair(prev_index,f_indices[(i+1)%3])] =
it_fb->second.halfedges[i];
}
}
//insert point inside face
for(Node_id node_id : node_ids)
{
CDT_Vertex_handle vh=cdt.insert(nodes.exact_node(node_id));
vh->info()=node_id;
id_to_CDT_vh.insert(std::make_pair(node_id,vh));
}
std::vector<std::pair<Node_id,Node_id> > constrained_edges;
// insert constraints that are interior to the triangle (in the case
// no edges are collinear in the meshes)
insert_constrained_edges(node_ids,cdt,id_to_CDT_vh,constrained_edges);
// insert constraints between points that are on the boundary
// (not a contrained on the triangle boundary)
if (it_fb!=face_boundaries.end()) //is f not a triangle ?
{
for (int i=0;i<3;++i)
{
Node_ids& ids=it_fb->second.node_ids_array[i];
insert_constrained_edges(ids,cdt,id_to_CDT_vh,constrained_edges,1);
}
}
//insert coplanar edges for endpoints of triangles
for (int i=0;i<3;++i){
Node_id nindex=triangle_vertices[i]->info();
if ( nindex < nb_nodes )
insert_constrained_edges_coplanar_case(nindex,cdt,id_to_CDT_vh);
}
//XSL_TAG_CPL_VERT
//collect edges incident to a point that is the intersection of two
// coplanar faces. This ensure that triangulations are compatible.
if (it_fb!=face_boundaries.end()) //is f not a triangle ?
{
for (typename CDT::Finite_vertices_iterator
vit=cdt.finite_vertices_begin(),
vit_end=cdt.finite_vertices_end();vit_end!=vit;++vit)
{
//skip original vertices (that are not nodes) and non-coplanar face
// issued vertices (this is working because intersection points
// between coplanar facets are the first inserted)
if (vit->info() >= nb_nodes ||
vit->info() >= number_coplanar_vertices) continue;
// \todo no need to insert constrained edges (they also are constrained
// in the other mesh)!!
typename std::map< Node_id,std::set<Node_id> >::iterator res =
coplanar_constraints.insert(
std::make_pair(vit->info(),std::set<Node_id>())).first;
//turn around the vertex and get incident edge
typename CDT::Edge_circulator start=cdt.incident_edges(vit);
typename CDT::Edge_circulator curr=start;
do{
if (cdt.is_infinite(*curr) ) continue;
typename CDT::Edge mirror=cdt.mirror_edge(*curr);
if ( cdt.is_infinite( curr->first->vertex(curr->second) ) ||
cdt.is_infinite( mirror.first->vertex(mirror.second) ) )
continue; // skip edges that are on the boundary of the triangle
// (these are already constrained)
//insert edges in the set of constraints
CDT_Vertex_handle vh=vit;
int nindex = curr->first->vertex((curr->second+1)%3)==vh
? (curr->second+2)%3
: (curr->second+1)%3;
CDT_Vertex_handle vn=curr->first->vertex(nindex);
if ( vit->info() > vn->info() || vn->info()>=nb_nodes)
continue; //take only one out of the two edges + skip input
CGAL_assertion(vn->info()<nb_nodes);
res->second.insert( vn->info() );
}while(start!=++curr);
}
}
// import the triangle in `cdt` in the face `f` of `tm`
triangulate_a_face(f, tm, nodes, node_ids, node_id_to_vertex,
edge_to_hedge, cdt, vpm, output_builder, user_visitor);
// TODO Here we do the update only for internal edges.
// Update for border halfedges could be done during the split
//3) mark halfedges that are common to two polyhedral surfaces
//recover halfedges inserted that are on the intersection
typedef std::pair<Node_id,Node_id> Node_id_pair;
for(const Node_id_pair& node_id_pair : constrained_edges)
{
typename std::map<Node_id_pair,halfedge_descriptor>
::iterator it_poly_hedge=edge_to_hedge.find(node_id_pair);
//we cannot have an assertion here in case an edge or part of an edge is a constraints.
//Indeed, the graph_of_constraints report an edge 0,1 and 1,0 for example while only one of the two
//is defined as one of them defines an adjacent face
//CGAL_assertion(it_poly_hedge!=edge_to_hedge.end());
if( it_poly_hedge!=edge_to_hedge.end() ){
call_put(marks_on_edges,tm,edge(it_poly_hedge->second,tm),true);
output_builder.set_edge_per_polyline(tm,node_id_pair,it_poly_hedge->second);
}
else{
//WARNING: in few case this is needed if the marked edge is on the border
//to optimize it might be better to only use sorted pair. TAG_SLXX1
Node_id_pair opposite_pair(node_id_pair.second,node_id_pair.first);
it_poly_hedge=edge_to_hedge.find(opposite_pair);
CGAL_assertion( it_poly_hedge!=edge_to_hedge.end() );
call_put(marks_on_edges,tm,edge(it_poly_hedge->second,tm),true);
output_builder.set_edge_per_polyline(tm,opposite_pair,it_poly_hedge->second);
}
}
}
}
void finalize(INodes& nodes,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2)
{
nodes.all_nodes_created();
TriangleMesh* tm1_ptr = const_cast<TriangleMesh*>(&tm1);
TriangleMesh* tm2_ptr = const_cast<TriangleMesh*>(&tm2);
std::map<TriangleMesh*, VertexPointMap> vpms;
vpms[tm1_ptr] = vpm1;
vpms[tm2_ptr] = vpm2;
vertex_descriptor null_vertex = Graph_traits::null_vertex();
const Node_id nb_nodes = nodes.size();
// we reserve nb_nodes+3 because we use the last three entries for the
@@ -617,7 +966,6 @@ public:
//store for each triangle face which boundary is intersected by the other surface,
//original vertices (and halfedges in the refined mesh pointing on these vertices)
typedef boost::unordered_map<face_descriptor,Face_boundary> Face_boundaries;
std::map<TriangleMesh*,Face_boundaries> mesh_to_face_boundaries;
//0) For each triangle mesh, collect original vertices that belongs to the intersection.
@@ -632,6 +980,7 @@ public:
{
TriangleMesh& tm=*it->first;
CGAL_assertion(&tm!=const_mesh_ptr);
// Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm];
Node_to_target_of_hedge_map& nodes_to_hedge=it->second;
@@ -700,95 +1049,10 @@ public:
for (typename std::map<TriangleMesh*,On_edge_map>::iterator
it=on_edge.begin(); it!=on_edge.end(); ++it)
{
TriangleMesh& tm=*it->first;
CGAL_assertion(&tm!=const_mesh_ptr);
const VertexPointMap& vpm=vpms[&tm];
On_edge_map& on_edge_map=it->second;
On_face_map& on_face_map=on_face[&tm];
Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm];
for(typename On_edge_map::iterator it2=on_edge_map.begin();
it2!=on_edge_map.end();
++it2)
{
//the edge to be split
halfedge_descriptor hedge=halfedge(it2->first,tm);
//indices of the nodes to be inserted
Node_ids& node_ids=it2->second;
CGAL_assertion( std::set<Node_id>(node_ids.begin(), node_ids.end())
.size()==node_ids.size() );
//sort nodes along the egde to allow consecutive splits
sort_vertices_along_hedge(node_ids,hedge,tm,vpm,nodes);
//save original face and nodes for face of hedge (1)
if ( !is_border(hedge,tm) ){
face_descriptor f=face(hedge,tm);
typename Face_boundaries::iterator it_face = face_boundaries.find(f);
if (it_face==face_boundaries.end())
it_face=face_boundaries.insert(std::make_pair(f,Face_boundary(hedge,tm))).first;
it_face->second.copy_node_ids(hedge,node_ids.begin(),node_ids.end());
}
//save original face and nodes for face of hedge->opposite (2)
typename Face_boundaries::iterator opposite_original_info=face_boundaries.end();
halfedge_descriptor hedge_opp = opposite(hedge,tm);
if ( !is_border(hedge_opp,tm) ){
face_descriptor f=face(hedge_opp,tm);
opposite_original_info=face_boundaries.find(f);
if (opposite_original_info==face_boundaries.end())
opposite_original_info=face_boundaries.insert(std::make_pair(f,Face_boundary(hedge_opp,tm))).first;
opposite_original_info->second.copy_node_ids(hedge_opp,node_ids.rbegin(),node_ids.rend());
}
typename Mesh_to_map_node::iterator it_map=mesh_to_node_id_to_vertex.find(&tm);
CGAL_assertion(it_map!=mesh_to_node_id_to_vertex.end());
//a map to identify the vertex in the polyhedron corresponding to an intersection point
Node_id_to_vertex& node_id_to_vertex=it_map->second;
CGAL_assertion_code(vertex_descriptor original_vertex=source(hedge,tm);)
//We need an edge incident to the source vertex of hedge. This is the first opposite edge created.
bool first=true;
halfedge_descriptor hedge_incident_to_src=Graph_traits::null_halfedge();
bool hedge_is_marked = call_get(marks_on_edges,tm,edge(hedge,tm));
//do split the edges
CGAL_assertion_code(vertex_descriptor expected_src=source(hedge,tm));
for(std::size_t node_id : node_ids)
{
halfedge_descriptor hnew = Euler::split_edge(hedge, tm);
CGAL_assertion(expected_src==source(hnew,tm));
vertex_descriptor vnew=target(hnew,tm);
// user_visitor.new_vertex_added(node_id, vnew, tm); // NODE_VISITOR_TAG
nodes.call_put(vpm, vnew, node_id, tm);
// register the new vertex in the output builder
output_builder.set_vertex_id(vnew, node_id, tm);
node_id_to_vertex[node_id]=vnew;
if (first){
first=false;
hedge_incident_to_src=next(opposite(hedge,tm),tm);
}
//update marker tags. If the edge was marked, then the resulting edges in the split must be marked
if ( hedge_is_marked )
call_put(marks_on_edges,tm,edge(hnew,tm),true);
CGAL_assertion_code(expected_src=vnew);
}
CGAL_assertion(target(hedge_incident_to_src,tm)==original_vertex);
CGAL_assertion(face(hedge_incident_to_src,tm)==face(hedge_opp,tm));
//save original face and nodes for face of hedge->opposite (2)
if ( !is_border(hedge_opp,tm) ){
CGAL_assertion(opposite_original_info!=face_boundaries.end());
opposite_original_info->second.update_original_halfedge(
hedge_opp,hedge_incident_to_src,tm);
}
//insert the two incident faces in on_face map so that they will be triangulated.
if (!is_border(hedge,tm)) on_face_map[face(hedge,tm)];
if (!is_border(hedge_opp,tm)) on_face_map[face(hedge_opp,tm)];
}
if(it->first == tm1_ptr)
split_halfedges(it, vpm1, nodes, mesh_to_face_boundaries);
else
split_halfedges(it, vpm2, nodes, mesh_to_face_boundaries);
}
//2)triangulation of the triangle faces containing intersection point in their interior
@@ -796,248 +1060,10 @@ public:
for (typename std::map<TriangleMesh*,On_face_map>::iterator
it=on_face.begin(); it!=on_face.end(); ++it)
{
TriangleMesh& tm=*it->first;
CGAL_assertion(&tm!=const_mesh_ptr);
const VertexPointMap& vpm=vpms[&tm];
On_face_map& on_face_map=it->second;
Face_boundaries& face_boundaries=mesh_to_face_boundaries[&tm];
Node_id_to_vertex& node_id_to_vertex=mesh_to_node_id_to_vertex[&tm];
Vertex_to_node_id& vertex_to_node_id=mesh_to_vertex_to_node_id[&tm];
for (typename On_face_map::iterator it=on_face_map.begin();
it!=on_face_map.end();++it)
{
face_descriptor f = it->first; //the face to be triangulated
Node_ids& node_ids = it->second; // ids of nodes in the interior of f
typename Face_boundaries::iterator it_fb=face_boundaries.find(f);
std::map<Node_id,typename CDT::Vertex_handle> id_to_CDT_vh;
//associate an edge of the triangulation to a halfedge in a given polyhedron
std::map<std::pair<Node_id,Node_id>,halfedge_descriptor> edge_to_hedge;
// the vertices of f
std::array<vertex_descriptor,3> f_vertices;
// the node_id of an input vertex or a fake id (>=nb_nodes)
std::array<Node_id,3> f_indices = {{nb_nodes,nb_nodes+1,nb_nodes+2}};
if (it_fb!=face_boundaries.end()){ //the boundary of the triangle face was refined
f_vertices[0]=it_fb->second.vertices[0];
f_vertices[1]=it_fb->second.vertices[1];
f_vertices[2]=it_fb->second.vertices[2];
update_face_indices(f_vertices,f_indices,vertex_to_node_id);
if (doing_autorefinement)
it_fb->second.update_node_id_to_vertex_map(node_id_to_vertex, tm);
}
else{
CGAL_assertion( is_triangle(halfedge(f,tm),tm) );
halfedge_descriptor h0=halfedge(f,tm), h1=next(h0,tm), h2=next(h1,tm);
f_vertices[0]=target(h0,tm); //nb_nodes
f_vertices[1]=target(h1,tm); //nb_nodes+1
f_vertices[2]=target(h2,tm); //nb_nodes+2
update_face_indices(f_vertices,f_indices,vertex_to_node_id);
edge_to_hedge[std::make_pair( f_indices[2],f_indices[0] )] = h0;
edge_to_hedge[std::make_pair( f_indices[0],f_indices[1] )] = h1;
edge_to_hedge[std::make_pair( f_indices[1],f_indices[2] )] = h2;
}
typename EK::Point_3 p = nodes.to_exact(get(vpm,f_vertices[0])),
q = nodes.to_exact(get(vpm,f_vertices[1])),
r = nodes.to_exact(get(vpm,f_vertices[2]));
///TODO use a positive normal and remove all work around to guarantee that triangulation of coplanar patches are compatible
CDT_traits traits(typename EK::Construct_normal_3()(p,q,r));
CDT cdt(traits);
// insert triangle points
std::array<CDT_Vertex_handle,3> triangle_vertices;
//we can do this to_exact because these are supposed to be input points.
triangle_vertices[0]=cdt.insert_outside_affine_hull(p);
triangle_vertices[1]=cdt.insert_outside_affine_hull(q);
triangle_vertices[2]=cdt.tds().insert_dim_up(cdt.infinite_vertex(), false);
triangle_vertices[2]->set_point(r);
triangle_vertices[0]->info()=f_indices[0];
triangle_vertices[1]->info()=f_indices[1];
triangle_vertices[2]->info()=f_indices[2];
node_id_to_vertex[nb_nodes ]=f_vertices[0];
node_id_to_vertex[nb_nodes+1]=f_vertices[1];
node_id_to_vertex[nb_nodes+2]=f_vertices[2];
//if one of the triangle input vertex is also a node
for (int ik=0;ik<3;++ik){
if ( f_indices[ik]<nb_nodes )
{
id_to_CDT_vh.insert(
std::make_pair(f_indices[ik],triangle_vertices[ik]));
if (doing_autorefinement)
// update the current vertex in node_id_to_vertex
// to match the one of the face
node_id_to_vertex[f_indices[ik]]=f_vertices[ik];
}
}
//insert points on edges
if (it_fb!=face_boundaries.end()) //if f not a triangle?
{
// collect infinite faces incident to the initial triangle
typename CDT::Face_handle infinite_faces[3];
for (int i=0;i<3;++i)
{
int oi=-1;
CGAL_assertion_code(bool is_edge = )
cdt.is_edge(triangle_vertices[i], triangle_vertices[(i+1)%3], infinite_faces[i], oi);
CGAL_assertion(is_edge);
CGAL_assertion( cdt.is_infinite( infinite_faces[i]->vertex(oi) ) );
}
// In this loop, for each original edge of the triangle, we insert
// the constrained edges and we recover the halfedge_descriptor
// corresponding to these constrained (they are already in tm)
Face_boundary& f_boundary=it_fb->second;
for (int i=0;i<3;++i){
//handle case of halfedge starting at triangle_vertices[i]
// and ending at triangle_vertices[(i+1)%3]
const Node_ids& ids_on_edge=f_boundary.node_ids_array[i];
CDT_Vertex_handle previous=triangle_vertices[i];
Node_id prev_index=f_indices[i];// node-id of the mesh vertex
halfedge_descriptor hedge = next(f_boundary.halfedges[(i+2)%3],tm);
CGAL_assertion( source(hedge,tm)==f_boundary.vertices[i] );
if (!ids_on_edge.empty()){ //is there at least one node on this edge?
// fh must be an infinite face
// The points must be ordered from fh->vertex(cw(infinite_vertex)) to fh->vertex(ccw(infinite_vertex))
for(Node_id id : ids_on_edge)
{
CDT_Vertex_handle vh=insert_point_on_ch_edge(cdt,infinite_faces[i],nodes.exact_node(id));
vh->info()=id;
id_to_CDT_vh.insert(std::make_pair(id,vh));
edge_to_hedge[std::make_pair(prev_index,id)]=hedge;
previous=vh;
hedge=next(hedge,tm);
prev_index=id;
}
}
else{
CGAL_assertion_code(halfedge_descriptor hd=f_boundary.halfedges[i]);
CGAL_assertion( target(hd,tm) == f_boundary.vertices[(i+1)%3] );
CGAL_assertion( source(hd,tm) == f_boundary.vertices[ i ] );
}
CGAL_assertion(hedge==f_boundary.halfedges[i]);
edge_to_hedge[std::make_pair(prev_index,f_indices[(i+1)%3])] =
it_fb->second.halfedges[i];
}
}
//insert point inside face
for(Node_id node_id : node_ids)
{
CDT_Vertex_handle vh=cdt.insert(nodes.exact_node(node_id));
vh->info()=node_id;
id_to_CDT_vh.insert(std::make_pair(node_id,vh));
}
std::vector<std::pair<Node_id,Node_id> > constrained_edges;
// insert constraints that are interior to the triangle (in the case
// no edges are collinear in the meshes)
insert_constrained_edges(node_ids,cdt,id_to_CDT_vh,constrained_edges);
// insert constraints between points that are on the boundary
// (not a contrained on the triangle boundary)
if (it_fb!=face_boundaries.end()) //is f not a triangle ?
{
for (int i=0;i<3;++i)
{
Node_ids& ids=it_fb->second.node_ids_array[i];
insert_constrained_edges(ids,cdt,id_to_CDT_vh,constrained_edges,1);
}
}
//insert coplanar edges for endpoints of triangles
for (int i=0;i<3;++i){
Node_id nindex=triangle_vertices[i]->info();
if ( nindex < nb_nodes )
insert_constrained_edges_coplanar_case(nindex,cdt,id_to_CDT_vh);
}
//XSL_TAG_CPL_VERT
//collect edges incident to a point that is the intersection of two
// coplanar faces. This ensure that triangulations are compatible.
if (it_fb!=face_boundaries.end()) //is f not a triangle ?
{
for (typename CDT::Finite_vertices_iterator
vit=cdt.finite_vertices_begin(),
vit_end=cdt.finite_vertices_end();vit_end!=vit;++vit)
{
//skip original vertices (that are not nodes) and non-coplanar face
// issued vertices (this is working because intersection points
// between coplanar facets are the first inserted)
if (vit->info() >= nb_nodes ||
vit->info() >= number_coplanar_vertices) continue;
// \todo no need to insert constrained edges (they also are constrained
// in the other mesh)!!
typename std::map< Node_id,std::set<Node_id> >::iterator res =
coplanar_constraints.insert(
std::make_pair(vit->info(),std::set<Node_id>())).first;
//turn around the vertex and get incident edge
typename CDT::Edge_circulator start=cdt.incident_edges(vit);
typename CDT::Edge_circulator curr=start;
do{
if (cdt.is_infinite(*curr) ) continue;
typename CDT::Edge mirror=cdt.mirror_edge(*curr);
if ( cdt.is_infinite( curr->first->vertex(curr->second) ) ||
cdt.is_infinite( mirror.first->vertex(mirror.second) ) )
continue; // skip edges that are on the boundary of the triangle
// (these are already constrained)
//insert edges in the set of constraints
CDT_Vertex_handle vh=vit;
int nindex = curr->first->vertex((curr->second+1)%3)==vh
? (curr->second+2)%3
: (curr->second+1)%3;
CDT_Vertex_handle vn=curr->first->vertex(nindex);
if ( vit->info() > vn->info() || vn->info()>=nb_nodes)
continue; //take only one out of the two edges + skip input
CGAL_assertion(vn->info()<nb_nodes);
res->second.insert( vn->info() );
}while(start!=++curr);
}
}
// import the triangle in `cdt` in the face `f` of `tm`
triangulate_a_face(f, tm, nodes, node_ids, node_id_to_vertex,
edge_to_hedge, cdt, vpm, output_builder, user_visitor);
// TODO Here we do the update only for internal edges.
// Update for border halfedges could be done during the split
//3) mark halfedges that are common to two polyhedral surfaces
//recover halfedges inserted that are on the intersection
typedef std::pair<Node_id,Node_id> Node_id_pair;
for(const Node_id_pair& node_id_pair : constrained_edges)
{
typename std::map<Node_id_pair,halfedge_descriptor>
::iterator it_poly_hedge=edge_to_hedge.find(node_id_pair);
//we cannot have an assertion here in case an edge or part of an edge is a constraints.
//Indeed, the graph_of_constraints report an edge 0,1 and 1,0 for example while only one of the two
//is defined as one of them defines an adjacent face
//CGAL_assertion(it_poly_hedge!=edge_to_hedge.end());
if( it_poly_hedge!=edge_to_hedge.end() ){
call_put(marks_on_edges,tm,edge(it_poly_hedge->second,tm),true);
output_builder.set_edge_per_polyline(tm,node_id_pair,it_poly_hedge->second);
}
else{
//WARNING: in few case this is needed if the marked edge is on the border
//to optimize it might be better to only use sorted pair. TAG_SLXX1
Node_id_pair opposite_pair(node_id_pair.second,node_id_pair.first);
it_poly_hedge=edge_to_hedge.find(opposite_pair);
CGAL_assertion( it_poly_hedge!=edge_to_hedge.end() );
call_put(marks_on_edges,tm,edge(it_poly_hedge->second,tm),true);
output_builder.set_edge_per_polyline(tm,opposite_pair,it_poly_hedge->second);
}
}
}
if(it->first == tm1_ptr)
triangulate_intersected_faces(it, vpm1, nodes, mesh_to_face_boundaries);
else
triangulate_intersected_faces(it, vpm2, nodes, mesh_to_face_boundaries);
}
nodes.finalize();
@@ -583,7 +583,8 @@ next_marked_halfedge_around_target_vertex(
template <class PolygonMesh,
class EdgeMap,
class VertexMap,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class VertexPointMapOut,
class IntersectionEdgeMap>
void import_polyline(
@@ -598,8 +599,8 @@ void import_polyline(
VertexMap& pm1_to_output_vertices,
const IntersectionEdgeMap& intersection_edges1,
const IntersectionEdgeMap& intersection_edges2,
const VertexPointMap& vpm1,
const VertexPointMap& /*vpm2*/,
const VertexPointMap1& vpm1,
const VertexPointMap2& /*vpm2*/,
const VertexPointMapOut& vpm_out,
std::vector<typename boost::graph_traits<PolygonMesh>
::edge_descriptor>& output_shared_edges)
@@ -1004,7 +1005,8 @@ void append_patches_to_triangle_mesh(
template < class TriangleMesh,
class IntersectionEdgeMap,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class VertexPointMapOut,
class EdgeMarkMap1,
class EdgeMarkMap2,
@@ -1024,8 +1026,8 @@ void fill_new_triangle_mesh(
const IntersectionPolylines& polylines,
const IntersectionEdgeMap& intersection_edges1,
const IntersectionEdgeMap& intersection_edges2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
const VertexPointMapOut& vpm_out,
const EdgeMarkMap1& edge_mark_map1,
const EdgeMarkMap2& edge_mark_map2,
@@ -1261,7 +1263,8 @@ template <class TriangleMesh,
class PatchContainer2,
class IntersectionPolylines,
class EdgeMap,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class EdgeMarkMapIn1,
class EdgeMarkMapIn2,
class EdgeMarkMapOut,
@@ -1275,8 +1278,8 @@ void compute_inplace_operation_delay_removal_and_insideout(
PatchContainer2& patches_of_tm2,
bool reverse_patch_orientation_tm2,
const IntersectionPolylines& polylines,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
EdgeMarkMapIn1&,
const EdgeMarkMapIn2& edge_mark_map2,
const EdgeMarkMapOut& edge_mark_map_out1,
@@ -1420,7 +1423,8 @@ remove_patches(TriangleMesh& tm,
template <class TriangleMesh,
class PatchContainer1,
class PatchContainer2,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class EdgeMarkMapIn1,
class EdgeMarkMapIn2,
class EdgeMarkMapOut1,
@@ -1434,8 +1438,8 @@ void compute_inplace_operation(
PatchContainer2& patches_of_tm2,
bool reverse_patch_orientation_tm1,
bool reverse_patch_orientation_tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
EdgeMarkMapIn1& edge_mark_map_in1,
const EdgeMarkMapIn2& edge_mark_map_in2,
EdgeMarkMapOut1& edge_mark_map_out1,
@@ -1527,7 +1531,8 @@ template <class TriangleMesh,
class PatchContainer1,
class PatchContainer2,
class IntersectionPolylines,
class VertexPointMap,
class VertexPointMap1,
class VertexPointMap2,
class EdgeMarkMapIn1,
class EdgeMarkMapIn2,
class EdgeMarkMapOut1,
@@ -1541,8 +1546,8 @@ void compute_inplace_operation(
PatchContainer2& patches_of_tm2,
bool reverse_patch_orientation_tm1,
bool reverse_patch_orientation_tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
const EdgeMarkMapIn1& edge_mark_map_in1,
const EdgeMarkMapIn2& edge_mark_map_in2,
const EdgeMarkMapOut1& edge_mark_map_out1,
@@ -78,7 +78,7 @@ find_intersection(const Point_3& p, const Point_3& q, //segment
}
template<class TriangleMesh, class VertexPointMap>
template<class TriangleMesh, class VertexPointMap1, class VertexPointMap2>
std::tuple<Intersection_type,
typename boost::graph_traits<TriangleMesh>::halfedge_descriptor,
bool,bool>
@@ -87,23 +87,26 @@ intersection_type(
typename boost::graph_traits<TriangleMesh>::face_descriptor f_2,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2)
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2)
{
typedef boost::graph_traits<TriangleMesh> GT;
typedef typename GT::halfedge_descriptor halfedge_descriptor;
typedef std::tuple<Intersection_type,halfedge_descriptor,bool,bool> result_type;
typedef typename boost::property_traits<VertexPointMap>::reference Point_ref;
typedef typename boost::property_traits<VertexPointMap>::value_type Point_3;
typedef typename boost::property_traits<VertexPointMap1>::reference Point_ref1;
typedef typename boost::property_traits<VertexPointMap2>::reference Point_ref2;
typedef typename boost::property_traits<VertexPointMap1>::value_type Point_3;
typedef typename Kernel_traits<Point_3>::Kernel Kernel;
CGAL_static_assertion((std::is_same<Point_3, typename boost::property_traits<VertexPointMap2>::value_type>::value));
halfedge_descriptor h_2=halfedge(f_2,tm2);
Point_ref a = get(vpm2, target(h_2,tm2) );
Point_ref b = get(vpm2, target(next(h_2,tm2),tm2) );
Point_ref c = get(vpm2, source(h_2,tm2) );
Point_ref p = get(vpm1, source(h_1,tm1) );
Point_ref q = get(vpm1, target(h_1,tm1) );
Point_ref2 a = get(vpm2, target(h_2,tm2) );
Point_ref2 b = get(vpm2, target(next(h_2,tm2),tm2) );
Point_ref2 c = get(vpm2, source(h_2,tm2) );
Point_ref1 p = get(vpm1, source(h_1,tm1) );
Point_ref1 q = get(vpm1, target(h_1,tm1) );
const Orientation abcp = orientation(a,b,c,p);
const Orientation abcq = orientation(a,b,c,q);
@@ -57,7 +57,7 @@ public:
halfedge_descriptor fh = face_box.info();
halfedge_descriptor eh = edge_box.info();
edge_to_faces[eh].insert(face(fh, tm_faces));
edge_to_faces[edge(eh,tm_edges)].insert(face(fh, tm_faces));
}
void operator()( const Box* face_box_ptr, const Box* edge_box_ptr) const
@@ -67,15 +67,15 @@ public:
};
template<class TriangleMesh,
class VertexPointMap,
class VertexPointMapF, class VertexPointMapE,
class EdgeToFaces,
class CoplanarFaceSet>
class Collect_face_bbox_per_edge_bbox_with_coplanar_handling {
protected:
const TriangleMesh& tm_faces;
const TriangleMesh& tm_edges;
const VertexPointMap& vpmap_tmf;
const VertexPointMap& vpmap_tme;
const VertexPointMapF& vpmap_tmf;
const VertexPointMapE& vpmap_tme;
EdgeToFaces& edge_to_faces;
CoplanarFaceSet& coplanar_faces;
@@ -83,7 +83,7 @@ protected:
typedef typename Graph_traits::face_descriptor face_descriptor;
typedef typename Graph_traits::halfedge_descriptor halfedge_descriptor;
typedef typename boost::property_traits<VertexPointMap>::reference Point;
typedef typename boost::property_traits<VertexPointMapF>::reference Point;
typedef CGAL::Box_intersection_d::ID_FROM_BOX_ADDRESS Box_policy;
typedef CGAL::Box_intersection_d::Box_with_info_d<double, 3, halfedge_descriptor, Box_policy> Box;
@@ -92,8 +92,8 @@ public:
Collect_face_bbox_per_edge_bbox_with_coplanar_handling(
const TriangleMesh& tm_faces,
const TriangleMesh& tm_edges,
const VertexPointMap& vpmap_tmf,
const VertexPointMap& vpmap_tme,
const VertexPointMapF& vpmap_tmf,
const VertexPointMapE& vpmap_tme,
EdgeToFaces& edge_to_faces,
CoplanarFaceSet& coplanar_faces)
: tm_faces(tm_faces)
@@ -83,10 +83,10 @@ struct Default_surface_intersection_visitor{
void start_new_polyline(std::size_t,std::size_t){}
void add_node_to_polyline(std::size_t){}
void input_have_coplanar_faces(){}
template<class T,class VertexPointMap>
template<class T,class VPM1,class VPM2>
void finalize(T&,
const TriangleMesh&, const TriangleMesh&,
const VertexPointMap&, const VertexPointMap&)
const VPM1, const VPM2)
{}
void new_node_added_triple_face(std::size_t /* node_id */,
face_descriptor /* f1 */,
@@ -144,7 +144,7 @@ struct Node_id_set {
};
template< class TriangleMesh,
class VertexPointMap,
class VertexPointMap1, class VertexPointMap2,
class Node_visitor=Default_surface_intersection_visitor<TriangleMesh>
>
class Intersection_of_triangle_meshes
@@ -175,7 +175,7 @@ class Intersection_of_triangle_meshes
// may contain several segments.
typedef std::map< Face_pair_and_int, Node_id_set > Faces_to_nodes_map;
typedef Intersection_nodes<TriangleMesh,
VertexPointMap,
VertexPointMap1, VertexPointMap2,
Predicates_on_constructions_needed> Node_vector;
// data members
@@ -188,11 +188,13 @@ class Intersection_of_triangle_meshes
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
CGAL_assertion_code(bool doing_autorefinement;)
// member functions
template <class VPMF, class VPME>
void filter_intersections(const TriangleMesh& tm_f,
const TriangleMesh& tm_e,
const VertexPointMap& vpm_f,
const VertexPointMap& vpm_e,
const VPMF& vpm_f,
const VPME& vpm_e,
bool throw_on_self_intersection)
{
std::vector<Box> face_boxes, edge_boxes;
@@ -237,7 +239,7 @@ class Intersection_of_triangle_meshes
Callback callback(tm_f, tm_e, edge_to_faces);
#else
typedef Collect_face_bbox_per_edge_bbox_with_coplanar_handling<
TriangleMesh, VertexPointMap, Edge_to_faces, Coplanar_face_set>
TriangleMesh, VPMF, VPME, Edge_to_faces, Coplanar_face_set>
Callback;
Callback callback(tm_f, tm_e, vpm_f, vpm_e, edge_to_faces, coplanar_faces);
#endif
@@ -258,8 +260,9 @@ class Intersection_of_triangle_meshes
}
// for autorefinement
template <class VPM>
void filter_intersections(const TriangleMesh& tm,
const VertexPointMap& vpm)
const VPM& vpm)
{
std::vector<Box> face_boxes, edge_boxes;
std::vector<Box*> face_boxes_ptr, edge_boxes_ptr;
@@ -296,7 +299,7 @@ class Intersection_of_triangle_meshes
Edge_to_faces& edge_to_faces = stm_edge_to_ltm_faces;
typedef Collect_face_bbox_per_edge_bbox_with_coplanar_handling_one_mesh<
TriangleMesh, VertexPointMap, Edge_to_faces, Coplanar_face_set>
TriangleMesh, VPM, Edge_to_faces, Coplanar_face_set>
Callback;
Callback callback(tm, vpm, edge_to_faces, coplanar_faces);
@@ -527,12 +530,12 @@ class Intersection_of_triangle_meshes
}
}
void compute_intersection_of_coplanar_faces(
Node_id& current_node,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2)
template <typename VPM1, typename VPM2>
void compute_intersection_of_coplanar_faces(Node_id& current_node,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VPM1& vpm1,
const VPM2& vpm2)
{
CGAL_assertion( &tm1 < &tm2 || &tm1==&tm2 );
@@ -637,12 +640,13 @@ class Intersection_of_triangle_meshes
//add a new node in the final graph.
//it is the intersection of the triangle with the segment
template <typename VPM1, typename VPM2>
void add_new_node(halfedge_descriptor h_1,
face_descriptor f_2,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VPM1& vpm1,
const VPM2& vpm2,
std::tuple<Intersection_type,
halfedge_descriptor,
bool,bool> inter_res)
@@ -657,11 +661,12 @@ class Intersection_of_triangle_meshes
}
}
template <typename VPM1, typename VPM2>
void compute_intersection_points(Edge_to_faces& tm1_edge_to_tm2_faces,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VPM1& vpm1,
const VPM2& vpm2,
Node_id& current_node)
{
typedef std::tuple<Intersection_type, halfedge_descriptor, bool,bool> Inter_type;
@@ -821,8 +826,9 @@ class Intersection_of_triangle_meshes
}
};
template <class VPM>
void detect_intersections_in_the_graph(const TriangleMesh& tm,
const VertexPointMap& vpm,
const VPM& vpm,
Node_id& current_node)
{
boost::unordered_map<face_descriptor,
@@ -1082,7 +1088,7 @@ class Intersection_of_triangle_meshes
template <class Output_iterator>
void construct_polylines(Output_iterator out){
typedef typename boost::property_traits<VertexPointMap>::value_type Point_3;
typedef typename boost::property_traits<VertexPointMap1>::value_type Point_3;
std::size_t nb_nodes=nodes.size();
std::vector<Graph_node> graph(nb_nodes);
//counts the number of time each node has been seen
@@ -1262,8 +1268,8 @@ class Intersection_of_triangle_meshes
public:
Intersection_of_triangle_meshes(const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
const Node_visitor& v=Node_visitor())
: nodes(tm1, tm2, vpm1, vpm2)
, visitor(v)
@@ -1275,7 +1281,7 @@ public:
// for autorefinement
Intersection_of_triangle_meshes(const TriangleMesh& tm,
const VertexPointMap& vpm,
const VertexPointMap1& vpm,
const Node_visitor& v=Node_visitor())
: nodes(tm, tm, vpm, vpm)
, visitor(v)
@@ -1293,8 +1299,8 @@ public:
const TriangleMesh& tm1=nodes.tm1;
const TriangleMesh& tm2=nodes.tm2;
const VertexPointMap& vpm1=nodes.vpm1;
const VertexPointMap& vpm2=nodes.vpm2;
const VertexPointMap1& vpm1=nodes.vpm1;
const VertexPointMap2& vpm2=nodes.vpm2;
filter_intersections(tm1, tm2, vpm1, vpm2, throw_on_self_intersection);
filter_intersections(tm2, tm1, vpm2, vpm1, throw_on_self_intersection);
@@ -1308,6 +1314,7 @@ public:
compute_intersection_of_coplanar_faces(current_node, tm1, tm2, vpm1, vpm2);
else
compute_intersection_of_coplanar_faces(current_node, tm2, tm1, vpm2, vpm1);
visitor.set_number_of_intersection_points_from_coplanar_faces(current_node+1);
if (!coplanar_faces.empty())
visitor.input_have_coplanar_faces();
@@ -1324,6 +1331,7 @@ public:
compute_intersection_points(tm1_edge_to_tm2_faces, tm1, tm2, vpm1, vpm2, current_node);
compute_intersection_points(tm2_edge_to_tm1_faces, tm2, tm1, vpm2, vpm1, current_node);
if (!build_polylines){
visitor.finalize(nodes,tm1,tm2,vpm1,vpm2);
return output;
@@ -1361,7 +1369,7 @@ public:
CGAL_assertion(doing_autorefinement);
const TriangleMesh& tm=nodes.tm1;
const VertexPointMap& vpm=nodes.vpm1;
const VertexPointMap1& vpm=nodes.vpm1;
filter_intersections(tm, vpm);
@@ -28,23 +28,26 @@ namespace Corefinement {
// polylines. Different specializations are available depending whether
// predicates on constructions are needed.
template <class TriangleMesh,
class VertexPointMap,
class VertexPointMap1, class VertexPointMap2,
bool Predicates_on_constructions_needed,
bool Has_exact_constructions=
bool Has_exact_constructions =
!boost::is_floating_point<
typename Kernel_traits<
typename boost::property_traits<VertexPointMap>::value_type
typename boost::property_traits<VertexPointMap1>::value_type
>::Kernel::FT
>::value >
class Intersection_nodes;
//Store only the double version of the intersection points.
template <class TriangleMesh,
class VertexPointMap>
class Intersection_nodes<TriangleMesh,VertexPointMap,false,false>
class VertexPointMap1, class VertexPointMap2>
class Intersection_nodes<TriangleMesh, VertexPointMap1, VertexPointMap2, false, false>
{
//typedefs
typedef typename boost::property_traits<VertexPointMap>::value_type Point_3;
typedef typename boost::property_traits<VertexPointMap1>::value_type Point_3;
CGAL_static_assertion((std::is_same<typename boost::property_traits<VertexPointMap1>::value_type,
typename boost::property_traits<VertexPointMap2>::value_type>::value));
typedef typename Kernel_traits<Point_3>::Kernel Input_kernel;
typedef std::vector <Point_3> Nodes_vector;
typedef CGAL::Exact_predicates_exact_constructions_kernel Exact_kernel;
@@ -67,12 +70,13 @@ class Intersection_nodes<TriangleMesh,VertexPointMap,false,false>
public:
const TriangleMesh &tm1, &tm2;
VertexPointMap vpm1, vpm2;
const VertexPointMap1& vpm1;
const VertexPointMap2& vpm2;
Intersection_nodes(const TriangleMesh& tm1_,
const TriangleMesh& tm2_,
const VertexPointMap& vpm1_,
const VertexPointMap& vpm2_)
const VertexPointMap1& vpm1_,
const VertexPointMap2& vpm2_)
: tm1(tm1_)
, tm2(tm2_)
, vpm1(vpm1_)
@@ -97,12 +101,13 @@ public:
//add a new node in the final graph.
//it is the intersection of the triangle with the segment
template <class VPM_A, class VPM_B> // VertexPointMap1 or VertexPointMap2
void add_new_node(halfedge_descriptor h_a,
face_descriptor f_b,
const TriangleMesh& tm_a,
const TriangleMesh& tm_b,
const VertexPointMap vpm_a,
const VertexPointMap& vpm_b)
const VPM_A& vpm_a,
const VPM_B& vpm_b)
{
halfedge_descriptor h_b = halfedge(f_b, tm_b);
add_new_node(
@@ -114,12 +119,8 @@ public:
to_exact( get(vpm_a, target(h_a,tm_a)) ) ) );
}
void add_new_node(halfedge_descriptor edge_1, face_descriptor face_2)
{
add_new_node(edge_1, face_2, tm1, tm2, vpm1, vpm2);
}
void call_put(const VertexPointMap& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh&)
template <class VPM> // VertexPointMap1 or VertexPointMap2
void call_put(const VPM& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh&)
{
put(vpm, vd, nodes[i]);
}
@@ -132,14 +133,18 @@ public:
// second specializations: store an exact copy of the points so
// that we can answer exactly predicates
template <class TriangleMesh, class VertexPointMap>
class Intersection_nodes<TriangleMesh,VertexPointMap,true,false>
template <class TriangleMesh, class VertexPointMap1, class VertexPointMap2>
class Intersection_nodes<TriangleMesh, VertexPointMap1, VertexPointMap2, true, false>
{
//typedefs
public:
typedef CGAL::Exact_predicates_exact_constructions_kernel Exact_kernel;
private:
typedef typename boost::property_traits<VertexPointMap>::value_type Point_3;
typedef typename boost::property_traits<VertexPointMap1>::value_type Point_3;
CGAL_static_assertion((std::is_same<typename boost::property_traits<VertexPointMap1>::value_type,
typename boost::property_traits<VertexPointMap2>::value_type>::value));
typedef typename Kernel_traits<Point_3>::Kernel Input_kernel;
typedef Cartesian_converter<Input_kernel,Exact_kernel> Double_to_exact;
@@ -160,14 +165,16 @@ private:
Exact_kernel::Intersect_3 exact_intersection;
std::vector<vertex_descriptor> tm1_vertices, tm2_vertices;
const bool doing_autorefinement;
public:
const TriangleMesh &tm1, &tm2;
VertexPointMap vpm1, vpm2;
const VertexPointMap1& vpm1;
const VertexPointMap2& vpm2;
Intersection_nodes(const TriangleMesh& tm1_,
const TriangleMesh& tm2_,
const VertexPointMap& vpm1_,
const VertexPointMap& vpm2_)
const VertexPointMap1& vpm1_,
const VertexPointMap2& vpm2_)
: doing_autorefinement(&tm1_ == &tm2_)
, tm1(tm1_)
, tm2(tm2_)
@@ -211,12 +218,13 @@ public:
//add a new node in the final graph.
//it is the intersection of the triangle with the segment
template <class VPM_A, class VPM_B> // VertexPointMap1 or VertexPointMap2
void add_new_node(halfedge_descriptor h_a,
face_descriptor f_b,
const TriangleMesh& tm_a,
const TriangleMesh& tm_b,
const VertexPointMap vpm_a,
const VertexPointMap& vpm_b)
const VPM_A vpm_a,
const VPM_B vpm_b)
{
halfedge_descriptor h_b = halfedge(f_b, tm_b);
add_new_node(
@@ -229,11 +237,12 @@ public:
}
// use to resolve intersection of 3 faces in autorefinement only
template <class VPM>
void add_new_node(halfedge_descriptor h1,
halfedge_descriptor h2,
halfedge_descriptor h3,
const TriangleMesh& tm,
const VertexPointMap& vpm)
const VPM& vpm)
{
// TODO Far from optimal!
typedef Exact_kernel::Plane_3 Plane_3;
@@ -257,11 +266,6 @@ public:
add_new_node(*pt);
}
void add_new_node(halfedge_descriptor edge_1, face_descriptor face_2)
{
add_new_node(edge_1, face_2, tm1, tm2, vpm1, vpm2);
}
//the point is an input
void add_new_node(const Point_3& p){
enodes.push_back(to_exact(p));
@@ -273,7 +277,8 @@ public:
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 VPM> // VertexPointMap1 or VertexPointMap2
void call_put(const VPM& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh& tm)
{
put(vpm, vd, exact_to_double(enodes[i]));
if (&tm1==&tm)
@@ -306,11 +311,16 @@ public:
//Third specialization: The kernel already has exact constructions.
template <class TriangleMesh,class VertexPointMap,bool Predicates_on_constructions_needed>
class Intersection_nodes<TriangleMesh,VertexPointMap,Predicates_on_constructions_needed,true>
template <class TriangleMesh, class VertexPointMap1, class VertexPointMap2,
bool Predicates_on_constructions_needed>
class Intersection_nodes<TriangleMesh, VertexPointMap1, VertexPointMap2,
Predicates_on_constructions_needed, true>
{
//typedefs
typedef typename boost::property_traits<VertexPointMap>::value_type Point_3;
typedef typename boost::property_traits<VertexPointMap1>::value_type Point_3;
CGAL_static_assertion((std::is_same<typename boost::property_traits<VertexPointMap1>::value_type,
typename boost::property_traits<VertexPointMap2>::value_type>::value));
typedef typename Kernel_traits<Point_3>::Kernel Input_kernel;
typedef std::vector <Point_3> Nodes_vector;
@@ -326,12 +336,13 @@ public:
typedef Input_kernel Exact_kernel;
const TriangleMesh &tm1, &tm2;
VertexPointMap vpm1, vpm2;
const VertexPointMap1& vpm1;
const VertexPointMap2& vpm2;
Intersection_nodes(const TriangleMesh& tm1_,
const TriangleMesh& tm2_,
const VertexPointMap& vpm1_,
const VertexPointMap& vpm2_)
const VertexPointMap1& vpm1_,
const VertexPointMap2& vpm2_)
: tm1(tm1_)
, tm2(tm2_)
, vpm1(vpm1_)
@@ -346,11 +357,12 @@ public:
size_t size() const {return nodes.size();}
const Point_3& exact_node(std::size_t i) const {return nodes[i];}
template <class VPM>
void add_new_node(halfedge_descriptor h1,
halfedge_descriptor h2,
halfedge_descriptor h3,
const TriangleMesh& tm,
const VertexPointMap& vpm)
const VPM& vpm)
{
// TODO Far from optimal!
typedef typename Exact_kernel::Plane_3 Plane_3;
@@ -376,12 +388,13 @@ public:
//add a new node in the final graph.
//it is the intersection of the triangle with the segment
template <class VPM_A, class VPM_B> // VertexPointMap1 or VertexPointMap2
void add_new_node(halfedge_descriptor h_a,
face_descriptor f_b,
const TriangleMesh& tm_a,
const TriangleMesh& tm_b,
const VertexPointMap vpm_a,
const VertexPointMap& vpm_b)
const VPM_A& vpm_a,
const VPM_B& vpm_b)
{
halfedge_descriptor h_b=halfedge(f_b,tm_b);
@@ -394,12 +407,6 @@ public:
get(vpm_a, target(h_a,tm_a)) ) );
}
void add_new_node(halfedge_descriptor edge_1, face_descriptor face_2)
{
add_new_node(edge_1, face_2, tm1, tm2, vpm1, vpm2);
}
void add_new_node(const Point_3& p)
{
nodes.push_back(p);
@@ -407,7 +414,8 @@ public:
const Point_3& to_exact(const Point_3& p) const { return p; }
void call_put(const VertexPointMap& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh&)
template <class VPM> // VertexPointMap1 or VertexPointMap2
void call_put(const VPM& vpm, vertex_descriptor vd, std::size_t i, TriangleMesh&)
{
put(vpm, vd, nodes[i]);
}
@@ -27,10 +27,15 @@ namespace CGAL{
namespace Polygon_mesh_processing {
namespace Corefinement{
template <class TriangleMesh, class VertexPointMap>
struct Intersect_coplanar_faces_3{
template <class TriangleMesh, class VertexPointMap1, class VertexPointMap2>
struct Intersect_coplanar_faces_3
{
// typedefs
typedef typename boost::property_traits<VertexPointMap>::value_type Point;
typedef typename boost::property_traits<VertexPointMap1>::value_type Point;
CGAL_static_assertion((std::is_same<typename boost::property_traits<VertexPointMap1>::value_type,
typename boost::property_traits<VertexPointMap1>::value_type>::value));
typedef typename CGAL::Kernel_traits<Point>::Kernel Input_kernel;
typedef CGAL::Exact_predicates_exact_constructions_kernel Exact_kernel;
@@ -40,12 +45,14 @@ struct Intersect_coplanar_faces_3{
typedef Coplanar_intersection<TriangleMesh, Exact_kernel> Inter_pt_info;
// data members
const TriangleMesh &tm1, &tm2;
const VertexPointMap &vpm1, &vpm2;
const VertexPointMap1& vpm1;
const VertexPointMap2& vpm2;
// constructor
Intersect_coplanar_faces_3(const TriangleMesh& tm1_,
const TriangleMesh& tm2_,
const VertexPointMap& vpm1_,
const VertexPointMap& vpm2_)
const VertexPointMap1& vpm1_,
const VertexPointMap2& vpm2_)
: tm1(tm1_), tm2(tm2_), vpm1(vpm1_), vpm2(vpm2_)
{}
@@ -282,14 +289,14 @@ struct Intersect_coplanar_faces_3{
}
};
template <class TriangleMesh, class VertexPointMap, class Exact_kernel>
template <class TriangleMesh, class VertexPointMap1, class VertexPointMap2, class Exact_kernel>
void intersection_coplanar_faces(
typename boost::graph_traits<TriangleMesh>::face_descriptor f1,
typename boost::graph_traits<TriangleMesh>::face_descriptor f2,
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const VertexPointMap& vpm1,
const VertexPointMap& vpm2,
const VertexPointMap1& vpm1,
const VertexPointMap2& vpm2,
std::list< Coplanar_intersection<TriangleMesh, Exact_kernel> >& inter_pts)
{
typedef boost::graph_traits<TriangleMesh> GT;
@@ -297,7 +304,7 @@ void intersection_coplanar_faces(
halfedge_descriptor h1=halfedge(f1,tm1), h2=halfedge(f2,tm2);
Intersect_coplanar_faces_3<TriangleMesh, VertexPointMap>
Intersect_coplanar_faces_3<TriangleMesh, VertexPointMap1, VertexPointMap2>
intersect_cpln(tm1, tm2, vpm1, vpm2);
// We will add in `inter_pts` the initial triangle of h1
@@ -122,15 +122,15 @@ bool are_triangles_coplanar_same_side(
}
template <class Node_id, class Node_vector, class vertex_descriptor, class Vpm>
template <class Node_id, class Node_vector, class vertex_descriptor, class VPMP, class VPMQ>
bool are_triangles_coplanar_same_side(Node_id o_prime_index,
Node_id o_index,
Node_id p_index,
Node_id q_index,
vertex_descriptor p,
vertex_descriptor q,
const Vpm& vpm_p,
const Vpm& vpm_q,
const VPMP& vpm_p,
const VPMQ& vpm_q,
const Node_vector& nodes)
{
const Node_id NID((std::numeric_limits<Node_id>::max)());
@@ -142,7 +142,7 @@ bool are_triangles_coplanar_same_side(Node_id o_prime_index,
);
}
template <class Node_id, class Node_vector, class vertex_descriptor, class Vpm>
template <class Node_id, class Node_vector, class vertex_descriptor, class VPMP, class VPMQ>
bool sorted_around_edge( Node_id o_prime_index,
Node_id o_index,
Node_id p1_index,
@@ -151,8 +151,8 @@ bool sorted_around_edge( Node_id o_prime_index,
vertex_descriptor p1,
vertex_descriptor p2,
vertex_descriptor q,
const Vpm& vpm_p,
const Vpm& vpm_q,
const VPMP& vpm_p,
const VPMQ& vpm_q,
const Node_vector& nodes)
{
const Node_id NID((std::numeric_limits<Node_id>::max)());
@@ -323,7 +323,6 @@ namespace internal {
{
halfedge_status_pmap_ = get(CGAL::dynamic_halfedge_property_t<Halfedge_status>(),
pmesh);
CGAL_assertion(CGAL::is_triangle_mesh(mesh_));
CGAL_assertion_code(input_mesh_is_valid_ = CGAL::is_valid_polygon_mesh(pmesh));
CGAL_warning_msg(input_mesh_is_valid_,
"The input mesh is not a valid polygon mesh. "
@@ -594,7 +593,6 @@ namespace internal {
#endif
#ifdef CGAL_PMP_REMESHING_DEBUG
CGAL_expensive_assertion(is_triangle_mesh(mesh_));
debug_status_map();
debug_self_intersections();
#endif
@@ -808,7 +806,6 @@ namespace internal {
#endif
#ifdef CGAL_PMP_REMESHING_DEBUG
CGAL_expensive_assertion(is_triangle_mesh(mesh_));
debug_status_map();
debug_self_intersections();
CGAL_assertion(PMP::remove_degenerate_faces(mesh_,
@@ -1038,7 +1035,6 @@ namespace internal {
}
CGAL_assertion(!input_mesh_is_valid_ || is_valid_polygon_mesh(mesh_));
CGAL_assertion(is_triangle_mesh(mesh_));
}//end for loop (nit == nb_iterations)
#ifdef CGAL_PMP_REMESHING_DEBUG
@@ -1072,7 +1068,6 @@ namespace internal {
put(vpmap_, v, proj);
}
CGAL_assertion(!input_mesh_is_valid_ || is_valid_polygon_mesh(mesh_));
CGAL_assertion(is_triangle_mesh(mesh_));
#ifdef CGAL_PMP_REMESHING_DEBUG
debug_self_intersections();
#endif
@@ -1101,7 +1096,6 @@ namespace internal {
put(vpmap_, v, proj(v));
}
CGAL_assertion(is_valid(mesh_));
CGAL_assertion(is_triangle_mesh(mesh_));
#ifdef CGAL_PMP_REMESHING_DEBUG
debug_self_intersections();
#endif
@@ -1315,7 +1309,7 @@ private:
{
if (is_on_patch_border(next(hopp, mesh_)) && is_on_patch_border(prev(hopp, mesh_)))
return false;
else if (next_on_patch_border(h) == hopp)
else if (next_on_patch_border(h) == hopp && prev_on_patch_border(h) == hopp)
return false; //isolated patch border
else
return true;
@@ -0,0 +1,207 @@
// Copyright (c) 2020 GeometryFactory (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
//
// Author(s) : Sébastien Loriot
#ifndef CGAL_POLYGON_MESH_PROCESSING_SIMPLIFY_POLYLINE_H
#define CGAL_POLYGON_MESH_PROCESSING_SIMPLIFY_POLYLINE_H
#include <CGAL/license/Polygon_mesh_processing/repair.h>
#include <CGAL/boost/graph/Named_function_parameters.h>
#include <CGAL/boost/graph/named_params_helper.h>
#include <type_traits>
namespace CGAL {
namespace Polygon_mesh_processing {
namespace experimental {
enum Polyline_simplification_algorithms { DOUGLAS_PEUCKER, ITERATIVE };
template <typename PointRangeIn, typename PointRangeOut,
typename NamedParametersIn, typename NamedParametersOut>
void simplify_polyline(const PointRangeIn& input,
PointRangeOut& output,
const double max_squared_frechet_distance,
const NamedParametersIn& np_in,
const NamedParametersOut& np_out)
{
using parameters::choose_parameter;
using parameters::get_parameter;
static_assert(std::is_same<
typename std::iterator_traits<typename PointRangeIn::const_iterator>::value_type,
typename std::iterator_traits<typename PointRangeOut::const_iterator>::value_type >::value, "");
typedef typename GetPointMap<PointRangeIn, NamedParametersIn>::type Point_map_in;
typedef typename GetPointMap<PointRangeOut, NamedParametersOut>::type Point_map_out;
typedef typename Point_set_processing_3::GetK<PointRangeIn, NamedParametersIn>::Kernel Kernel;
Point_map_in in_pm = choose_parameter<Point_map_in>(get_parameter(np_in, internal_np::point_map));
Point_map_out out_pm = choose_parameter<Point_map_out>(get_parameter(np_out, internal_np::point_map));
const Polyline_simplification_algorithms algorithm =
choose_parameter(get_parameter(np_in, internal_np::algorithm), DOUGLAS_PEUCKER);
switch(algorithm)
{
case ITERATIVE:
{
const bool is_closed = input.front()==input.back();
std::size_t nb_points = is_closed ? input.size()-1 : input.size();
// skip points in the input range that do not contains any information
if (nb_points<=2)
{
output.reserve(input.size());
for (const auto& p : input)
{
output.push_back(p);
put(out_pm, output.back(), get(in_pm, p));
}
return;
}
auto is_valid_approx = [&input, &in_pm, max_squared_frechet_distance](
std::size_t b, std::size_t e,
const typename Kernel::Line_3& line)
{
typename Kernel::Compare_squared_distance_3 compare_squared_distance;
for (std::size_t i=b+1; i<e; ++i)
{
if (compare_squared_distance(get(in_pm, input[i]), line, max_squared_frechet_distance) == LARGER)
return false;
}
return true;
};
std::size_t bi=0;
while(bi!=nb_points)
{
std::size_t ei=bi+2;
output.push_back(input[bi]);
put(out_pm, output.back(), get(in_pm, input[bi]));
while(ei<nb_points)
{
typename Kernel::Line_3 sl(get(in_pm,input[bi]), get(in_pm,input[ei]));
if (is_valid_approx(bi,ei,sl))
++ei; // we skip ei-1
else
{
bi=ei-1; // ei-1 shall not be skipt
break;
}
}
if(ei>=nb_points) break;
}
output.push_back(input[nb_points-1]);
put(out_pm, output.back(), get(in_pm, input[nb_points-1]));
if (is_closed)
{
output.push_back(input.back());
put(out_pm, output.back(), get(in_pm, input.back()));
}
return;
}
case DOUGLAS_PEUCKER:
{
const bool is_closed = input.front()==input.back();
std::size_t nb_points = is_closed ? input.size()-1 : input.size();
if (nb_points<=2)
{
output.reserve(input.size());
for (const auto& p : input)
{
output.push_back(p);
put(out_pm, output.back(), get(in_pm, p));
}
return;
}
std::vector< std::pair<std::size_t, std::size_t> > ranges;
ranges.push_back(std::make_pair(0, nb_points-1));
std::vector<bool> kept(input.size(), false);
if (is_closed) kept[nb_points]=true;
while( !ranges.empty() )
{
std::size_t rb, re;
std::tie(rb, re) = ranges.back();
ranges.pop_back();
kept[rb]=true;
kept[re]=true;
if (rb+1==re) continue;
typename Kernel::Line_3 line(get(in_pm, input[rb]), get(in_pm, input[re]));
double max_d = max_squared_frechet_distance;
std::size_t max_i = 0;
for (std::size_t i=rb; i<re; ++i)
{
double d = squared_distance(line, get(in_pm, input[i]));
if (d > max_d)
{
max_d = d;
max_i = i;
}
}
if (max_i != 0)
{
ranges.push_back( std::make_pair(max_i, re) );
ranges.push_back( std::make_pair(rb, max_i) );
}
}
std::size_t nb_kept=0;
for (std::size_t i=0; i<input.size(); ++i)
if (kept[i]) ++nb_kept;
output.reserve(nb_kept);
for (std::size_t i=0; i<input.size(); ++i)
if (kept[i])
{
output.push_back(input[i]);
put(out_pm, output.back(), get(in_pm, input[i]));
}
//TODO if is_closed-==true, shall we add en extra step to see if we can remove output.front() and output[output.size()-2] (inital endpoints)
}
}
}
template <typename PointRangeIn, typename PointRangeOut>
void simplify_polyline(const PointRangeIn& input,
PointRangeOut& output,
const double max_squared_frechet_distance)
{
simplify_polyline(input, output, max_squared_frechet_distance,
parameters::all_default(), parameters::all_default());
}
template <typename PointRangeIn, typename PointRangeOut, typename NamedParametersIn>
void simplify_polyline(const PointRangeIn& input,
PointRangeOut& output,
const double max_squared_frechet_distance,
const NamedParametersIn& np_in)
{
simplify_polyline(input, output, max_squared_frechet_distance, np_in, parameters::all_default());
}
} } } // end of CGAL::Polygon_mesh_processing::experimental namespace
#endif // CGAL_POLYGON_MESH_PROCESSING_SIMPLIFY_POLYLINE_H
@@ -1758,21 +1758,18 @@ surface_intersection(const TriangleMesh& tm1,
const bool throw_on_self_intersection =
parameters::choose_parameter(parameters::get_parameter(np1, internal_np::throw_on_self_intersection), false);
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters1>::const_type Vpm;
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters2>::const_type Vpm2;
CGAL_USE_TYPE(Vpm2);
CGAL_assertion_code(
static const bool same_vpm = (boost::is_same<Vpm,Vpm2>::value);)
CGAL_static_assertion(same_vpm);
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters1>::const_type VPM1;
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters2>::const_type VPM2;
Vpm vpm1 = parameters::choose_parameter(parameters::get_parameter(np1, internal_np::vertex_point),
get_const_property_map(CGAL::vertex_point, tm1));
Vpm vpm2 = parameters::choose_parameter(parameters::get_parameter(np2, internal_np::vertex_point),
get_const_property_map(CGAL::vertex_point, tm2));
CGAL_static_assertion((std::is_same<typename boost::property_traits<VPM1>::value_type,
typename boost::property_traits<VPM2>::value_type>::value));
Corefinement::Intersection_of_triangle_meshes<TriangleMesh,Vpm>
VPM1 vpm1 = parameters::choose_parameter(parameters::get_parameter(np1, internal_np::vertex_point),
get_const_property_map(CGAL::vertex_point, tm1));
VPM2 vpm2 = parameters::choose_parameter(parameters::get_parameter(np2, internal_np::vertex_point),
get_const_property_map(CGAL::vertex_point, tm2));
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, VPM1, VPM2>
functor(tm1, tm2, vpm1, vpm2);
return functor(polyline_output, throw_on_self_intersection, true);
}
@@ -1814,17 +1811,14 @@ surface_self_intersection(const TriangleMesh& tm,
const NamedParameters& np)
{
// Vertex point maps
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters>::const_type Vpm;
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::const_type VPM;
Vpm vpm = parameters::choose_parameter(parameters::get_parameter(np, internal_np::vertex_point),
get_const_property_map(CGAL::vertex_point, tm));
VPM vpm = parameters::choose_parameter(parameters::get_parameter(np, internal_np::vertex_point),
get_const_property_map(CGAL::vertex_point, tm));
// surface intersection algorithm call
typedef Corefinement::Default_surface_intersection_visitor<TriangleMesh,
true> Visitor;
Corefinement::Intersection_of_triangle_meshes<TriangleMesh,Vpm, Visitor>
functor(tm, vpm);
typedef Corefinement::Default_surface_intersection_visitor<TriangleMesh, true> Visitor;
Corefinement::Intersection_of_triangle_meshes<TriangleMesh, VPM, VPM, Visitor> functor(tm, vpm);
polyline_output=functor(polyline_output, true);
return polyline_output;
@@ -132,11 +132,6 @@ get_descriptor_from_location(const std::pair<typename boost::graph_traits<Triang
std::array<FT, 3> >& loc,
const TriangleMesh& tm);
template <typename FT, typename TriangleMesh>
Face_location<TriangleMesh, FT>
locate_in_face(typename boost::graph_traits<TriangleMesh>::halfedge_descriptor hd,
const FT t,
const TriangleMesh& tm);
// end of forward declarations
namespace internal {
@@ -1694,13 +1689,14 @@ locate_with_AABB_tree(const typename internal::Location_traits<TriangleMesh, Nam
{
typedef typename internal::Location_traits<TriangleMesh, NamedParameters>::Point Point;
typedef internal::Point_to_Point_3<TriangleMesh, Point> P_to_P3;
typedef typename boost::property_traits<Point3VPM>::value_type Point_3;
CGAL_static_assertion((std::is_same<Point_3, typename P_to_P3::Point_3>::value));
typedef typename GetGeomTraits<TriangleMesh, NamedParameters>::type Geom_traits;
typedef typename CGAL::AABB_face_graph_triangle_primitive<TriangleMesh, Point3VPM> Primitive;
typedef typename CGAL::AABB_traits<Geom_traits, Primitive> AABB_traits;
typedef typename Primitive::Point Point_3;
CGAL_static_assertion((std::is_same<Point_3, typename P_to_P3::Point_3>::value));
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::const_type VertexPointMap;
typedef internal::Point_to_Point_3_VPM<TriangleMesh, VertexPointMap> WrappedVPM;
@@ -17,6 +17,7 @@
#include <CGAL/algorithm.h>
#include <CGAL/boost/graph/iterator.h>
#include <CGAL/boost/graph/Named_function_parameters.h>
#include <CGAL/boost/graph/named_params_helper.h>
#include <CGAL/Container_helper.h>
#include <CGAL/Dynamic_property_map.h>
@@ -290,6 +290,17 @@ void polygon_soup_to_polygon_mesh(const PointRange& points,
converter(out, vpm);
}
template<typename PolygonMesh,
typename PointRange, typename PolygonRange,
typename NamedParameters_PS>
void polygon_soup_to_polygon_mesh(const PointRange& points,
const PolygonRange& polygons,
PolygonMesh& out,
const NamedParameters_PS& np_ps)
{
return polygon_soup_to_polygon_mesh(points, polygons, out, np_ps, parameters::all_default());
}
template<typename PolygonMesh, typename PointRange, typename PolygonRange>
void polygon_soup_to_polygon_mesh(const PointRange& points,
const PolygonRange& polygons,
@@ -256,7 +256,7 @@ std::size_t remove_connected_components_of_negligible_size(TriangleMesh& tmesh,
// Volumes make no sense for CCs that are not closed
std::vector<bool> cc_closeness(num, true);
std::vector<FT> component_volumes(num);
std::vector<FT> component_volumes(num, FT(0));
if(use_volumes)
{
@@ -183,13 +183,15 @@ bool is_collapse_geometrically_valid(typename boost::graph_traits<TriangleMesh>:
*/
template <class TriangleMesh, typename VPM, typename Traits>
boost::optional<double> get_collapse_volume(typename boost::graph_traits<TriangleMesh>::halfedge_descriptor h,
const TriangleMesh& tmesh,
const VPM& vpm,
const Traits& gt)
boost::optional<typename Traits::FT>
get_collapse_volume(typename boost::graph_traits<TriangleMesh>::halfedge_descriptor h,
const TriangleMesh& tmesh,
const VPM& vpm,
const Traits& gt)
{
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename Traits::FT FT;
typedef typename boost::property_traits<VPM>::reference Point_ref;
typedef typename Traits::Vector_3 Vector_3;
@@ -204,8 +206,8 @@ boost::optional<double> get_collapse_volume(typename boost::graph_traits<Triangl
Point_ref removed= get(vpm, target(h, tmesh));
// init volume with incident triangles (reversed orientation
double delta_vol = volume(removed, kept, get(vpm, target(next(h, tmesh), tmesh)), origin) +
volume(kept, removed, get(vpm, target(next(opposite(h, tmesh), tmesh), tmesh)), origin);
FT delta_vol = volume(removed, kept, get(vpm, target(next(h, tmesh), tmesh)), origin) +
volume(kept, removed, get(vpm, target(next(opposite(h, tmesh), tmesh), tmesh)), origin);
// consider triangles incident to the vertex removed
halfedge_descriptor stop = prev(opposite(h, tmesh), tmesh);
@@ -248,13 +250,14 @@ get_best_edge_orientation(typename boost::graph_traits<TriangleMesh>::edge_descr
const Traits& gt)
{
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename Traits::FT FT;
halfedge_descriptor h = halfedge(e, tmesh), ho = opposite(h, tmesh);
CGAL_assertion(!get(vcm, source(h, tmesh)) || !get(vcm, target(h, tmesh)));
boost::optional<double> dv1 = get_collapse_volume(h, tmesh, vpm, gt);
boost::optional<double> dv2 = get_collapse_volume(ho, tmesh, vpm, gt);
boost::optional<FT> dv1 = get_collapse_volume(h, tmesh, vpm, gt);
boost::optional<FT> dv2 = get_collapse_volume(ho, tmesh, vpm, gt);
// the resulting point of the collapse of a halfedge is the target of the halfedge before collapse
if(get(vcm, source(h, tmesh)))
@@ -288,6 +291,7 @@ bool should_flip(typename boost::graph_traits<TriangleMesh>::edge_descriptor e,
{
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename Traits::FT FT;
typedef typename boost::property_traits<VPM>::reference Point_ref;
typedef typename Traits::Vector_3 Vector_3;
@@ -318,16 +322,16 @@ bool should_flip(typename boost::graph_traits<TriangleMesh>::edge_descriptor e,
const Vector_3 v23 = gt.construct_vector_3_object()(p2, p3);
const Vector_3 v30 = gt.construct_vector_3_object()(p3, p0);
const double p1p3 = gt.compute_scalar_product_3_object()(
gt.construct_cross_product_vector_3_object()(v12, v23),
gt.construct_cross_product_vector_3_object()(v30, v01));
const FT p1p3 = gt.compute_scalar_product_3_object()(
gt.construct_cross_product_vector_3_object()(v12, v23),
gt.construct_cross_product_vector_3_object()(v30, v01));
const Vector_3 v21 = gt.construct_opposite_vector_3_object()(v12);
const Vector_3 v03 = gt.construct_opposite_vector_3_object()(v30);
const double p0p2 = gt.compute_scalar_product_3_object()(
gt.construct_cross_product_vector_3_object()(v01, v21),
gt.construct_cross_product_vector_3_object()(v23, v03));
const FT p0p2 = gt.compute_scalar_product_3_object()(
gt.construct_cross_product_vector_3_object()(v01, v21),
gt.construct_cross_product_vector_3_object()(v23, v03));
return p0p2 <= p1p3;
}
@@ -1694,37 +1698,23 @@ bool remove_degenerate_faces(const FaceRange& face_range,
// Ignore faces with null edges
if(!all_removed)
{
std::map<edge_descriptor, bool> are_degenerate_edges;
for(face_descriptor fd : degenerate_face_set)
typename std::set<face_descriptor>::iterator it = degenerate_face_set.begin();
while(it != degenerate_face_set.end())
{
for(halfedge_descriptor hd : halfedges_around_face(halfedge(fd, tmesh), tmesh))
bool has_degenerate_edge = false;
for(halfedge_descriptor hd : halfedges_around_face(halfedge(*it, tmesh), tmesh))
{
edge_descriptor ed = edge(hd, tmesh);
std::pair<typename std::map<edge_descriptor, bool>::iterator, bool> is_insert_successful =
are_degenerate_edges.insert(std::make_pair(ed, false));
bool is_degenerate = false;
if(is_insert_successful.second)
const edge_descriptor ed = edge(hd, tmesh);
if(is_degenerate_edge(ed, tmesh, np))
{
// did not previously exist in the map, so actually have to check if it is degenerate
if(traits.equal_3_object()(get(vpmap, target(ed, tmesh)), get(vpmap, source(ed, tmesh))))
is_degenerate = true;
}
is_insert_successful.first->second = is_degenerate;
if(is_degenerate)
{
halfedge_descriptor h = halfedge(ed, tmesh);
if(!is_border(h, tmesh))
degenerate_face_set.erase(face(h, tmesh));
h = opposite(h, tmesh);
if(!is_border(h, tmesh))
degenerate_face_set.erase(face(h, tmesh));
has_degenerate_edge = true;
it = degenerate_face_set.erase(it);
break;
}
}
if(!has_degenerate_edge)
++it;
}
}
@@ -348,13 +348,6 @@ std::size_t remove_invalid_polygons_in_polygon_soup(PointRange& /*points*/,
} // end namespace internal
template <typename PointRange, typename PolygonRange>
std::size_t remove_degenerate_polygons_in_polygon_soup(PointRange& points,
PolygonRange& polygons)
{
return remove_degenerate_polygons_in_polygon_soup(points, polygons, CGAL::parameters::all_default());
}
/// \ingroup PMP_repairing_grp
///
/// removes the isolated points from a polygon soup.
@@ -15,8 +15,7 @@
#include <CGAL/license/Polygon_mesh_processing/repair.h>
#include <CGAL/Polygon_mesh_processing/detect_features.h>
#include <CGAL/Polygon_mesh_processing/distance.h>
#include <CGAL/Polygon_mesh_processing/connected_components.h>
#include <CGAL/Polygon_mesh_processing/manifoldness.h>
#include <CGAL/Polygon_mesh_processing/orient_polygon_soup.h>
#include <CGAL/Polygon_mesh_processing/polygon_soup_to_polygon_mesh.h>
@@ -403,7 +403,7 @@ is_needle_triangle_face(typename boost::graph_traits<TriangleMesh>::face_descrip
if(min_sq_length == 0)
return min_h;
const FT sq_threshold = threshold * threshold;
const FT sq_threshold = square(FT(threshold));
if(max_sq_length / min_sq_length >= sq_threshold)
{
CGAL_assertion(min_h != boost::graph_traits<TriangleMesh>::null_halfedge());
@@ -482,7 +482,7 @@ is_cap_triangle_face(typename boost::graph_traits<TriangleMesh>::face_descriptor
typedef typename Traits::FT FT;
typedef typename Traits::Vector_3 Vector_3;
const FT sq_threshold = threshold * threshold;
const FT sq_threshold = square(FT(threshold));
const halfedge_descriptor h0 = halfedge(f, tm);
std::array<FT, 3> sq_lengths;
File diff suppressed because it is too large Load Diff