songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #1578 from maxGimeno/PMP-add_distance-GF 

Add approximated Hausdorff distance
This commit is contained in:
Laurent Rineau 2017-01-12 11:10:25 +01:00
parent a5e254f0d8 2cac99b2d6
commit 27cb95f022
35 changed files with 2807 additions and 209 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
AABB_tree/doc/AABB_tree/Concepts/AABBGeomTraits.h
View File

@ -65,7 +65,6 @@ A functor object to construct the sphere centered at one point and passing throu
typedef unspecified_type Construct_sphere_3;
/*!
\todo This is not correct! that is not used!
A functor object to compute the point on a geometric primitive which is closest from a query. Provides the operator:
`Point_3 operator()(const Type_2& type_2, const Point_3& p);` where `Type_2` is any type among `Segment_3` and `Triangle_3`. The operator returns the point on `type_2` which is closest to `p`.
*/
@ -76,7 +75,7 @@ A functor object to compare the distance of two points wrt a third one.
Provides the operator:
`CGAL::Comparision_result operator()(const Point_3& p1, const Point_3& p2, const Point_3& p3)`. The operator compare the distance between `p1 and `p2`, and between `p2` and `p3`.
*/
typedef unspecified_type Compare_distance_3
typedef unspecified_type Compare_distance_3;
/*!
A functor object to detect if a point lies inside a sphere or not.

12
Documentation/doc/biblio/geom.bib
View File

@ -151977,6 +151977,7 @@ pages = {179--189}
year={2015}
}
@incollection{bfhhm-epsph-15,
title={Exact Minkowski sums of polygons with holes},
author={Baram, Alon and Fogel, Efi and Halperin, Dan and Hemmer, Michael and Morr, Sebastian},
@ -151985,3 +151986,14 @@ pages = {179--189}
year={2015},
publisher={Springer}
}
@inproceedings{cignoni1998metro,
title={Metro: measuring error on simplified surfaces},
author={Cignoni, Paolo and Rocchini, Claudio and Scopigno, Roberto},
booktitle={Computer Graphics Forum},
volume={17},
number={2},
pages={167--174},
year={1998},
organization={Blackwell Publishers}
}

8
Generator/doc/Generator/CGAL/point_generators_2.h
View File

@ -377,7 +377,9 @@ typedef const Point_2& reference;
\sa `std::random_shuffle`
*/
template< typename Point_2, typename Triangulation >
template< typename Point_2,
typename Triangulation,
typename Creator = Creator_uniform_2<typename Kernel_traits<Point_2>::Kernel::RT,Point_2> >
class Random_points_in_triangle_mesh_2 {
public:
@ -446,7 +448,9 @@ get_default_random() );
\sa `std::random_shuffle`
*/
template< typename Point_2 >
template< typename Point_2,
typename Triangle_2 = typename Kernel_traits<Point_2>::Kernel::Triangle_2,
typename Creator = Creator_uniform_2<typename Kernel_traits<Point_2>::Kernel::RT,Point_2> >
class Random_points_in_triangles_2 {
public:

100
Generator/doc/Generator/CGAL/point_generators_3.h
View File

@ -248,7 +248,71 @@ get_default_random());
/// @}
}; /* end Random_points_in_triangle_3 */
} /* end namespace CGAL */
namespace CGAL{
/*!
The class `Random_points_on_segment_3` is an input iterator creating points uniformly
distributed on a segment. The default `Creator` is
`Creator_uniform_3<Kernel_traits<Point_3>::Kernel::RT,Point_3>`.
\cgalModels `InputIterator`
\cgalModels `PointGenerator`
\sa `CGAL::cpp11::copy_n()`
\sa `CGAL::Counting_iterator`
\sa `std::random_shuffle`
*/
template< typename Point_3, typename Creator >
class Random_points_on_segment_3 {
public:
/// \name Types
/// @{
/*!
*/
typedef std::input_iterator_tag iterator_category;
/*!
*/
typedef Point_3 value_type;
/*!
*/
typedef std::ptrdiff_t difference_type;
/*!
*/
typedef const Point_3* pointer;
/*!
*/
typedef const Point_3& reference;
/*!
creates an input iterator `g` generating points of type `Point_3` uniformly
distributed on the segment from \f$ p\f$ to \f$ q\f$ (excluding \f$ q\f$),
i.e.\ \f$ *g == (1-\lambda)\, p + \lambda q\f$ where \f$ 0 \le\lambda< 1\f$.
A single random number is needed from `rnd` for each point.
The expressions `to_double(p.x())`, `to_double(p.y())`, and `to_double(p.z())` must result
in the respective `double` representation of the coordinates of \f$ p\f$, and similarly for \f$ q\f$.
*/
Random_points_on_segment_3( const Point_3& p, const Point_3& q,
Random& rnd = get_default_random());
/// @}
}; /* end Random_points_on_segment_3 */
} /* end namespace CGAL */
namespace CGAL {
@ -264,7 +328,6 @@ distributed inside a tetrahedron. The default `Creator` is
\sa `CGAL::cpp11::copy_n()`
\sa `CGAL::Counting_iterator`
\sa `CGAL::Random_points_in_disc_2<Point_2, Creator>`
\sa `CGAL::Random_points_in_cube_3<Point_3, Creator>`
\sa `CGAL::Random_points_in_triangle_3<Point_3, Creator>`
\sa `CGAL::Random_points_on_sphere_3<Point_3, Creator>`
@ -344,13 +407,9 @@ The triangle range must be valid and unchanged while the iterator is used.
\sa `CGAL::cpp11::copy_n()`
\sa `CGAL::Counting_iterator`
\sa `CGAL::Points_on_segment_2<Point_2>`
\sa `CGAL::Random_points_in_disc_2<Point_2, Creator>`
\sa `CGAL::Random_points_on_segment_2<Point_2, Creator>`
\sa `CGAL::Random_points_on_square_2<Point_2, Creator>`
\sa `CGAL::Random_points_in_cube_3<Point_3, Creator>`
\sa `CGAL::Random_points_in_triangle_3<Point_2, Creator>`
\sa `CGAL::Random_points_in_tetrahedron_3<Point_2, Creator>`
\sa `CGAL::Random_points_in_triangle_3<Point_3, Creator>`
\sa `CGAL::Random_points_in_tetrahedron_3<Point_3, Creator>`
\sa `CGAL::Random_points_in_triangle_mesh_3<Point_3, TriangleMesh>`
\sa `CGAL::Random_points_in_tetrahedral_mesh_boundary_3<C3T3>`
\sa `CGAL::Random_points_in_tetrahedral_mesh_3<C3T3>`
@ -358,7 +417,10 @@ The triangle range must be valid and unchanged while the iterator is used.
\sa `std::random_shuffle`
*/
template< typename Point_3>
template< typename Point_3,
typename Triangle_3=typename Kernel_traits<Point_3>::Kernel::Triangle_3,
typename Creator = Creator_uniform_3< typename Kernel_traits< Point_3 >::Kernel::RT,
Point_3 > >
class Random_points_in_triangles_3 {
public:
@ -429,8 +491,12 @@ The triangle mesh must be valid and unchanged while the iterator is used.
\sa `std::random_shuffle`
*/
template < class TriangleMesh, class VertexPointMap = typename boost::property_map<TriangleMesh,
CGAL::vertex_point_t>::type> >
template < class TriangleMesh,
class VertexPointMap = typename boost::property_map<TriangleMesh,
CGAL::vertex_point_t>::type>,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename boost::property_traits<VertexPointMap>::value_type >::Kernel::RT,
typename boost::property_traits<VertexPointMap>::value_type > >
class Random_points_in_triangle_mesh_3 {
public:
@ -476,7 +542,6 @@ Similar to the previous constructor using `get(vertex_point, mesh)` as vertex po
*/
Random_points_in_triangle_mesh_3(const TriangleMesh& mesh, Random& rnd = get_default_random() );
/// @}
}; /* end Random_points_in_triangle_mesh_3 */
@ -509,7 +574,11 @@ C3T3 is a model of `Mesh_complex_3_in_triangulation_3`
\sa `std::random_shuffle`
*/
template <class C3T3>
template <class C3T3,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename C3t3::Point >::Kernel::RT,
typename C3t3::Point >
>
class Random_points_in_tetrahedral_mesh_boundary_3 {
public:
@ -582,7 +651,10 @@ C3T3 is a model of `Mesh_complex_3_in_triangulation_3`
\sa `std::random_shuffle`
*/
template <class C3T3>
template <class C3T3,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename C3t3::Point >::Kernel::RT,
typename C3t3::Point > >
class Random_points_in_tetrahedral_mesh_3 {
public:

2
Generator/include/CGAL/internal/Generic_random_point_generator.h
View File

@ -115,7 +115,7 @@ template< class Functor >
struct Apply_approx_sqrt: public Functor
{
template <class T>
double operator()(const T& t) const
typename cpp11::result_of<Functor(T)>::type operator()(const T& t) const
{
return approximate_sqrt( static_cast<const Functor&>(*this)(t) );
}

20
Generator/include/CGAL/point_generators_2.h
View File

@ -505,7 +505,7 @@ class Random_points_in_triangle_2 : public Random_generator_base<P> {
void generate_point();
public:
typedef P result_type;
typedef Random_points_in_triangle_2<P> This;
typedef Random_points_in_triangle_2<P, Creator> This;
typedef typename Kernel_traits<P>::Kernel::Triangle_2 Triangle_2;
Random_points_in_triangle_2() {}
Random_points_in_triangle_2( const This& x,Random& rnd = get_default_random())
@ -561,7 +561,11 @@ public:
}
};
}//end namespace internal
template <class P, class T>
template <class P,
class T,
class Creator =
Creator_uniform_2<typename Kernel_traits<P>::Kernel::RT,P>
>
class Random_points_in_triangle_mesh_2 : public Generic_random_point_generator<
typename T::Face_handle ,
internal::Triangle_from_face_2<T>,
@ -569,11 +573,11 @@ class Random_points_in_triangle_mesh_2 : public Generic_random_point_generator<
public:
typedef Generic_random_point_generator<typename T::Face_handle,
internal::Triangle_from_face_2<T>,
Random_points_in_triangle_2<P>,
Random_points_in_triangle_2<P, Creator>,
P> Base;
typedef typename T::Face_handle Id;
typedef P result_type;
typedef Random_points_in_triangle_mesh_2<P, T> This;
typedef Random_points_in_triangle_mesh_2<P, T, Creator> This;
Random_points_in_triangle_mesh_2( const T& triangulation, Random& rnd = get_default_random())
@ -621,7 +625,9 @@ struct Address_of {
}//namesapce internal
template <class Point_2,
class Triangle_2=typename Kernel_traits<Point_2>::Kernel::Triangle_2>
class Triangle_2=typename Kernel_traits<Point_2>::Kernel::Triangle_2,
class Creator =
Creator_uniform_2<typename Kernel_traits<Point_2>::Kernel::RT,Point_2> >
struct Random_points_in_triangles_2
: public Generic_random_point_generator<const Triangle_2*,
internal::Deref<Triangle_2>,
@ -630,11 +636,11 @@ struct Random_points_in_triangles_2
{
typedef Generic_random_point_generator<const Triangle_2*,
internal::Deref<Triangle_2>,
Random_points_in_triangle_2<Point_2>,
Random_points_in_triangle_2<Point_2, Creator>,
Point_2> Base;
typedef const Triangle_2* Id;
typedef Point_2 result_type;
typedef Random_points_in_triangles_2<Point_2> This;
typedef Random_points_in_triangles_2<Point_2, Triangle_2, Creator> This;
template<typename TriangleRange>
Random_points_in_triangles_2( const TriangleRange& triangles, Random& rnd = get_default_random())

181
Generator/include/CGAL/point_generators_3.h
View File

@ -205,7 +205,7 @@ class Random_points_in_triangle_3 : public Random_generator_base<P> {
void generate_point();
public:
typedef P result_type;
typedef Random_points_in_triangle_3<P> This;
typedef Random_points_in_triangle_3<P, Creator> This;
typedef typename Kernel_traits<P>::Kernel::Triangle_3 Triangle_3;
Random_points_in_triangle_3() {}
Random_points_in_triangle_3( const This& x,Random& rnd = get_default_random())
@ -251,13 +251,78 @@ void Random_points_in_triangle_3<P, Creator>::generate_point() {
}
template < class P, class Creator =
Creator_uniform_3<typename Kernel_traits<P>::Kernel::RT,P> >
class Random_points_on_segment_3 : public Random_generator_base<P> {
P _p;
P _q;
void generate_point();
public:
typedef Random_points_on_segment_3<P,Creator> This;
Random_points_on_segment_3() {}
Random_points_on_segment_3( const P& p,
const P& q,
Random& rnd = CGAL::get_default_random())
// g is an input iterator creating points of type `P' uniformly
// distributed on the segment from p to q except q, i.e. `*g' ==
// \lambda p + (1-\lambda)\, q where 0 <= \lambda < 1 . A single
// random number is needed from `rnd' for each point.
: Random_generator_base<P>( (std::max)( (std::max)( (std::max)(to_double(p.x()), to_double(q.x())),
(std::max)(to_double(p.y()), to_double(q.y()))),
(std::max)(to_double(p.z()), to_double(q.z()))),
rnd) , _p(p), _q(q)
{
generate_point();
}
template <class Segment_3>
Random_points_on_segment_3( const Segment_3& s,
Random& rnd = CGAL::get_default_random())
// g is an input iterator creating points of type `P' uniformly
// distributed on the segment from p to q except q, i.e. `*g' ==
// \lambda p + (1-\lambda)\, q where 0 <= \lambda < 1 . A single
// random number is needed from `rnd' for each point.
: Random_generator_base<P>( (std::max)( (std::max)( (std::max)(to_double(s[0].x()), to_double(s[1].x())),
(std::max)(to_double(s[0].y()), to_double(s[1].y()))),
(std::max)(to_double(s[0].z()), to_double(s[1].z()))),
rnd) , _p(s[0]), _q(s[1])
{
generate_point();
}
const P& source() const { return _p; }
const P& target() const { return _q; }
This& operator++() {
generate_point();
return *this;
}
This operator++(int) {
This tmp = *this;
++(*this);
return tmp;
}
};
template < class P, class Creator >
void
Random_points_on_segment_3<P,Creator>::
generate_point() {
typedef typename Creator::argument_type T;
double la = this->_rnd.get_double();
double mu = 1.0 - la;
Creator creator;
this->d_item = creator(T(mu * to_double(_p.x()) + la * to_double(_q.x())),
T(mu * to_double(_p.y()) + la * to_double(_q.y())),
T(mu * to_double(_p.z()) + la * to_double(_q.z())));
}
template < class P, class Creator =
Creator_uniform_3<typename Kernel_traits<P>::Kernel::RT,P> >
class Random_points_in_tetrahedron_3 : public Random_generator_base<P> {
P _p,_q,_r,_s;
void generate_point();
public:
typedef P result_type;
typedef Random_points_in_tetrahedron_3<P> This;
typedef Random_points_in_tetrahedron_3<P, Creator> This;
typedef typename Kernel_traits<P>::Kernel::Tetrahedron_3 Tetrahedron_3;
Random_points_in_tetrahedron_3() {}
Random_points_in_tetrahedron_3( const This& x,Random& rnd = get_default_random())
@ -306,42 +371,44 @@ void Random_points_in_tetrahedron_3<P, Creator>::generate_point() {
template <class TriangleMesh, class VertexPointMap = typename boost::property_map<TriangleMesh,
CGAL::vertex_point_t>::const_type>
template <class TriangleMesh,
class VertexPointMap = typename boost::property_map<TriangleMesh,
CGAL::vertex_point_t>::const_type,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename boost::property_traits<VertexPointMap>::value_type >::Kernel::RT,
typename boost::property_traits<VertexPointMap>::value_type >
>
struct Random_points_in_triangle_mesh_3
: public Generic_random_point_generator<
typename boost::graph_traits <TriangleMesh>::face_descriptor ,
CGAL::Property_map_to_unary_function<CGAL::Triangle_from_face_descriptor_map<
TriangleMesh, VertexPointMap > >,
Random_points_in_triangle_3<typename boost::property_traits<VertexPointMap>::value_type>,
Random_points_in_triangle_3<typename boost::property_traits<VertexPointMap>::value_type, Creator>,
typename boost::property_traits<VertexPointMap>::value_type>
{
typedef typename boost::property_traits<VertexPointMap>::value_type P;
typedef Generic_random_point_generator<
typename boost::graph_traits <TriangleMesh>::face_descriptor ,
CGAL::Property_map_to_unary_function<typename CGAL::Triangle_from_face_descriptor_map<
CGAL::Property_map_to_unary_function<CGAL::Triangle_from_face_descriptor_map<
TriangleMesh,VertexPointMap> >,
Random_points_in_triangle_3<P> , P> Base;
Random_points_in_triangle_3<P, Creator> , P> Base;
typedef typename CGAL::Triangle_from_face_descriptor_map<
TriangleMesh,VertexPointMap> Pmap;
typedef typename CGAL::Triangle_from_face_descriptor_map<
TriangleMesh,VertexPointMap> Object_from_id_map;
typedef Random_points_in_triangle_3<P> Generator_on_object;
TriangleMesh,VertexPointMap> Object_from_id;
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor Id;
typedef P result_type;
typedef Random_points_in_triangle_mesh_3< TriangleMesh, VertexPointMap> This;
typedef Random_points_in_triangle_mesh_3< TriangleMesh, VertexPointMap, Creator> This;
Random_points_in_triangle_mesh_3( const TriangleMesh& mesh,Random& rnd = get_default_random())
: Base( faces(mesh),
CGAL::Property_map_to_unary_function<Pmap>(Pmap(&mesh, get(vertex_point, mesh))),
CGAL::Property_map_to_unary_function<Object_from_id>(Object_from_id(&mesh, get(vertex_point, mesh))),
internal::Apply_approx_sqrt<typename Kernel_traits<P>::Kernel::Compute_squared_area_3>(),
rnd )
{
}
Random_points_in_triangle_mesh_3( const TriangleMesh& mesh, VertexPointMap vpm, Random& rnd = get_default_random())
: Base( faces(mesh),
CGAL::Property_map_to_unary_function<Pmap>(Pmap(&mesh, vpm)),
CGAL::Property_map_to_unary_function<Object_from_id>(Object_from_id(&mesh, vpm)),
internal::Apply_approx_sqrt<typename Kernel_traits<P>::Kernel::Compute_squared_area_3>(),
rnd )
{
@ -361,6 +428,62 @@ struct Random_points_in_triangle_mesh_3
}
};
template <class EdgeListGraph,
class VertexPointMap = typename boost::property_map<EdgeListGraph,
CGAL::vertex_point_t>::const_type,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename boost::property_traits<VertexPointMap>::value_type >::Kernel::RT,
typename boost::property_traits<VertexPointMap>::value_type >
>
struct Random_points_on_edge_list_graph_3
: public Generic_random_point_generator<
typename boost::graph_traits <EdgeListGraph>::edge_descriptor,
CGAL::Property_map_to_unary_function<CGAL::Segment_from_edge_descriptor_map<
EdgeListGraph, VertexPointMap > >,
Random_points_on_segment_3<typename boost::property_traits<VertexPointMap>::value_type, Creator>,
typename boost::property_traits<VertexPointMap>::value_type>
{
typedef typename boost::property_traits<VertexPointMap>::value_type P;
typedef Generic_random_point_generator<
typename boost::graph_traits <EdgeListGraph>::edge_descriptor,
CGAL::Property_map_to_unary_function<CGAL::Segment_from_edge_descriptor_map<
EdgeListGraph, VertexPointMap > >,
Random_points_on_segment_3<P, Creator> , P> Base;
typedef typename CGAL::Segment_from_edge_descriptor_map<
EdgeListGraph,VertexPointMap> Object_from_id;
typedef typename boost::graph_traits<EdgeListGraph>::edge_descriptor Id;
typedef P result_type;
typedef Random_points_on_edge_list_graph_3< EdgeListGraph, VertexPointMap, Creator> This;
Random_points_on_edge_list_graph_3( const EdgeListGraph& mesh,Random& rnd = get_default_random())
: Base( edges(mesh),
CGAL::Property_map_to_unary_function<Object_from_id>(Object_from_id(&mesh, get(vertex_point, mesh))),
internal::Apply_approx_sqrt<typename Kernel_traits<P>::Kernel::Compute_squared_length_3>(),
rnd )
{
}
Random_points_on_edge_list_graph_3( const EdgeListGraph& mesh, VertexPointMap vpm, Random& rnd = get_default_random())
: Base( edges(mesh),
CGAL::Property_map_to_unary_function<Object_from_id>(Object_from_id(&mesh, vpm)),
internal::Apply_approx_sqrt<typename Kernel_traits<P>::Kernel::Compute_squared_length_3>(),
rnd )
{
}
This& operator++() {
Base::generate_point();
return *this;
}
This operator++(int) {
This tmp = *this;
++(*this);
return tmp;
}
double mesh_length() const
{
return this->sum_of_weights();
}
};
namespace internal
{
@ -408,7 +531,11 @@ public:
};
}//end namespace internal
template <class C3t3>
template <class C3t3,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename C3t3::Point >::Kernel::RT,
typename C3t3::Point >
>
struct Random_points_in_tetrahedral_mesh_boundary_3
: public Generic_random_point_generator<
std::pair<typename C3t3::Triangulation::Cell_handle, int>,
@ -419,11 +546,11 @@ struct Random_points_in_tetrahedral_mesh_boundary_3
typedef Generic_random_point_generator<
std::pair<typename C3t3::Triangulation::Cell_handle, int>,
internal::Triangle_from_face_C3t3<typename C3t3::Triangulation>,
Random_points_in_triangle_3<typename C3t3::Point>,
Random_points_in_triangle_3<typename C3t3::Point, Creator>,
typename C3t3::Point> Base;
typedef std::pair<typename C3t3::Triangulation::Cell_handle, int> Id;
typedef typename C3t3::Point result_type;
typedef Random_points_in_tetrahedral_mesh_boundary_3<C3t3> This;
typedef Random_points_in_tetrahedral_mesh_boundary_3<C3t3, Creator> This;
Random_points_in_tetrahedral_mesh_boundary_3( const C3t3& c3t3,Random& rnd = get_default_random())
@ -445,7 +572,11 @@ struct Random_points_in_tetrahedral_mesh_boundary_3
}
};
template <class C3t3>
template <class C3t3,
class Creator = Creator_uniform_3<
typename Kernel_traits< typename C3t3::Point >::Kernel::RT,
typename C3t3::Point >
>
struct Random_points_in_tetrahedral_mesh_3
: public Generic_random_point_generator<
typename C3t3::Triangulation::Cell_handle,
@ -456,11 +587,11 @@ struct Random_points_in_tetrahedral_mesh_3
typedef Generic_random_point_generator<
typename C3t3::Triangulation::Cell_handle,
internal::Tetrahedron_from_cell_C3t3<typename C3t3::Triangulation>,
Random_points_in_tetrahedron_3<typename C3t3::Point>,
Random_points_in_tetrahedron_3<typename C3t3::Point, Creator>,
typename C3t3::Point> Base;
typedef typename C3t3::Triangulation::Cell_handle Id;
typedef typename C3t3::Point result_type;
typedef Random_points_in_tetrahedral_mesh_3<C3t3> This;
typedef Random_points_in_tetrahedral_mesh_3<C3t3, Creator> This;
Random_points_in_tetrahedral_mesh_3( const C3t3& c3t3,Random& rnd = get_default_random())
@ -484,7 +615,11 @@ struct Random_points_in_tetrahedral_mesh_3
template <class Point_3,
class Triangle_3=typename Kernel_traits<Point_3>::Kernel::Triangle_3>
class Triangle_3=typename Kernel_traits<Point_3>::Kernel::Triangle_3,
class Creator = Creator_uniform_3<
typename Kernel_traits< Point_3 >::Kernel::RT,
Point_3 >
>
struct Random_points_in_triangles_3
: public Generic_random_point_generator<const Triangle_3*,
internal::Deref<Triangle_3>,
@ -493,11 +628,11 @@ struct Random_points_in_triangles_3
{
typedef Generic_random_point_generator<const Triangle_3*,
internal::Deref<Triangle_3>,
Random_points_in_triangle_3<Point_3>,
Random_points_in_triangle_3<Point_3, Creator>,
Point_3> Base;
typedef const Triangle_3* Id;
typedef Point_3 result_type;
typedef Random_points_in_triangles_3<Point_3> This;
typedef Random_points_in_triangles_3<Point_3, Triangle_3, Creator> This;
template<typename TriangleRange>
Random_points_in_triangles_3( const TriangleRange& triangles, Random& rnd = get_default_random())

7
Generator/test/Generator/test_generators.cpp
View File

@ -115,13 +115,16 @@ void test_point_generators_3() {
/* Create test point set. */
std::vector<Point_3> points;
points.reserve(500);
points.reserve(600);
Random_points_in_sphere_3<Point_3,Creator> g1( 100.0);
CGAL::cpp11::copy_n( g1, 100, std::back_inserter(points));
Random_points_on_sphere_3<Point_3,Creator> g2( 100.0);
Random_points_in_cube_3<Point_3,Creator> g3( 100.0);
Random_points_on_segment_3<Point_3,Creator> g4( Point_3(-100,-100, -100),
Point_3( 100, 100, 100));
CGAL::cpp11::copy_n( g2, 100, std::back_inserter(points));
CGAL::cpp11::copy_n( g3, 100, std::back_inserter(points));
CGAL::cpp11::copy_n( g4, 100, std::back_inserter(points));
points_on_cube_grid_3( 50.0, (std::size_t)1,
std::back_inserter(points), Creator());
points_on_cube_grid_3( 50.0, (std::size_t)2,
@ -133,7 +136,7 @@ void test_point_generators_3() {
random_selection( points.begin(), points.end(), 100,
std::back_inserter(points));
assert( points.size() == 500);
assert( points.size() == 600);
for ( std::vector<Point_3>::iterator i = points.begin();
i != points.end(); i++){
assert( i->x() <= 100);

11
Installation/changes.html
View File

@ -278,6 +278,16 @@ and <code>src/</code> directories).
point sets based on the assumption that they sample a curve in 2D
or a surface in 3D.</li>
</ul>
<h3>Polygon Mesh Processing</h3>
<ul>
<li>
Add functions to compute approximated Hausdorff distances between two meshes, a mesh and a point set, or a point set and a mesh:
<code>sample_triangle_mesh()</code>, <code>approximated_Hausdorff_distance()</code>,
<code>approximated_symmetric_Hausdorff_distance()</code>,
<code>max_distance_to_triangle_mesh()</code>,
<code>max_distance_to_point_set()</code>.
</li>
</ul>
<!-- Spatial Searching and Sorting -->
<!-- Geometric Optimization -->
<!-- Interpolation -->
@ -310,6 +320,7 @@ and <code>src/</code> directories).
<li>in a tetrahedral mesh (<code>CGAL::Random_points_in_tetrahedral_mesh_3</code>),</li>
<li>in a 2D triangle mesh (<code>CGAL::Random_points_in_triangle_mesh_2</code>),</li>
<li>in a range of 2D or 3D triangles (<code>CGAL::Random_points_in_triangles_3</code> and <code>CGAL::Random_points_in_triangles_2</code>).</li>
<li>on a 3D segment (<code>CGAL::Random_points_on_segment_3</code>).</li>
</ul>
</li>
</ul>

12
Poisson_surface_reconstruction_3/examples/Poisson_surface_reconstruction_3/poisson_reconstruction_example.cpp
View File

@ -12,6 +12,8 @@
#include <CGAL/IO/read_xyz_points.h>
#include <CGAL/compute_average_spacing.h>
#include <CGAL/Polygon_mesh_processing/distance.h>
#include <vector>
#include <fstream>
@ -104,5 +106,15 @@ int main(void)
CGAL::output_surface_facets_to_polyhedron(c2t3, output_mesh);
out << output_mesh;
/// [PMP_distance_snippet]
// computes the approximation error of the reconstruction
double max_dist =
CGAL::Polygon_mesh_processing::approximate_max_distance_to_point_set(output_mesh,
points,
4000);
std::cout << "Max distance to point_set: " << max_dist << std::endl;
/// [PMP_distance_snippet]
return EXIT_SUCCESS;
}

63
Polygon_mesh_processing/doc/Polygon_mesh_processing/Concepts/PMPDistanceTraits.h Normal file
View File

@ -0,0 +1,63 @@
/// \ingroup PkgPolygonMeshProcessingConcepts
/// \cgalConcept
///
/// The concept `PMPDistanceTraits` is a refinement of the
/// concepts `AABBGeomTraits` and `SpatialSortingTraits_3`. In addition to the types required by
/// those concepts, it also requires types and functors needed by the functions `approximate_max_distance_to_point_set()`,
/// `sample_triangle_mesh()`, `approximate_Hausdorff_distance()` and `max_distance_to_triangle_mesh()`
///
/// \cgalRefines `AABBGeomTraits`
/// \cgalRefines `SpatialSortingTraits_3`
/// \cgalHasModel Any 3D Kernel is a model of this concept.
class PMPDistanceTraits{
public:
/*!
* A number type model of `Field` and `RealEmbeddable`
*/
typedef unspecified_type FT;
///
/*! 3D point type
* It must be default constructible, and can be constructed from 3 objects of type `FT`.
* `bool operator<(Point_3, Point_3)` to lexicographically compare two points must be available.
* Access to Cartesian coordinates must be possible using `Point_3::x()`, `Point_3::y(), Point_3::z()` and
* `FT operator[](int i)` with `0 <= i < 3`.
*
* There must be a specialization of `CGAL::Kernel_traits` such that
* `CGAL::Kernel_traits<Point_3>::%Kernel` is a model implementing this concept.
*/
typedef unspecified_type Point_3;
/// 3D vector type
typedef unspecified_type Vector_3;
/// @name Functors
/// @{
/// Functor for computing squared area of a triangle.
/// It provides `FT operator()(const Point_3&, const Point_3&, const Point_3&) const`
/// and `FT operator()(const Triangle_3&) const` and has `FT` as result_type.
typedef unspecified_type Compute_squared_area_3;
/// Functor for computing squared length of a segment.
/// and `FT operator()(const Segment_3&) const` and has `FT` as result_type.
typedef unspecified_type Compute_squared_length_3;
/// Functor for constructing translated points.
/// It provides `Point_3 operator()(const Point_3 &, const Vector_3 &)`
typedef unspecified_type Construct_translated_point_3;
/// Functor for constructing vectors.
/// It provides `Vector_3 operator()(const Point_3 &, const Point_3 &)`
typedef unspecified_type Construct_vector_3;
/// Functor for constructing scaled vectors.
/// It provides `Vector_3 operator()(const Vector_3 &, const FT &)`
typedef unspecified_type Construct_scaled_vector_3;
/// @}
/// @name Functions
/// @{
Compute_squared_area_3 compute_squared_area_3_object();
Compute_squared_length_3 compute_squared_length_3_object();
Construct_translated_point_3 construct_translated_point_3_object();
Construct_vector_3 construct_vector_3_object();
Construct_scaled_vector_3 construct_scaled_vector_3_object();
/// @}
};

2
Polygon_mesh_processing/doc/Polygon_mesh_processing/Concepts/SelfIntersectionTraits.h → Polygon_mesh_processing/doc/Polygon_mesh_processing/Concepts/PMPSelfIntersectionTraits.h
View File

@ -2,7 +2,7 @@
/// \cgalConcept
///
/// Geometric traits concept for the functions `CGAL::self_intersections()` and `CGAL::does_self_intersect()`.
class SelfIntersectionTraits{
class PMPSelfIntersectionTraits{
public:
/// @name Geometric Types
/// @{

2
Polygon_mesh_processing/doc/Polygon_mesh_processing/Doxyfile.in
View File

@ -18,6 +18,8 @@ ALIASES += "cgalNPTableEnd=</table> </dd> </dl>"
ALIASES += "cgalNPBegin{1}=<tr><td class=\"paramname\">\1 </td><td>"
ALIASES += "cgalNPEnd=</td></tr>"
EXAMPLE_PATH += ${CGAL_Poisson_surface_reconstruction_3_EXAMPLE_DIR}
MACRO_EXPANSION = YES
EXPAND_ONLY_PREDEF = YES
EXPAND_AS_DEFINED = CGAL_PMP_NP_TEMPLATE_PARAMETERS \

106
Polygon_mesh_processing/doc/Polygon_mesh_processing/NamedParameters.txt
View File

@ -210,7 +210,113 @@ set to `true`, this parameter is ignored.
\n
\b Type : `bool` \n
\b Default value is `true`
\cgalNPEnd
\cgalNPBegin{use_random_uniform_sampling} \anchor PMP_use_random_uniform_sampling
Parameter used in `sample_triangle_mesh()` to indicate if points should be picked
in a random uniform way.
\n
\b Type : `bool` \n
\b Default value is `true`
\cgalNPEnd
\cgalNPBegin{use_grid_sampling} \anchor PMP_use_grid_sampling
Parameter used in `sample_triangle_mesh()` to indicate if points should be picked
in on a grid in each face.
\n
\b Type : `bool` \n
\b Default value is `false`
\cgalNPEnd
\cgalNPBegin{use_monte_carlo_sampling} \anchor PMP_use_monte_carlo_sampling
Parameter used in `sample_triangle_mesh()` to indicate if points should be picked
using a Monte-Carlo approach.
\n
\b Type : `bool` \n
\b Default value is `false`
\cgalNPEnd
\cgalNPBegin{sample_edges} \anchor PMP_sample_edges
Parameter used in `sample_triangle_mesh()` to indicate if a dedicated sampling
of edges should be done.
\n
\b Type : `bool` \n
\b Default value is `true`
\cgalNPEnd
\cgalNPBegin{sample_vertices} \anchor PMP_sample_vertices
Parameter used in `sample_triangle_mesh()` to indicate if triangle vertices should
be copied in the output iterator.
\n
\b Type : `bool` \n
\b Default value is `true`
\cgalNPEnd
\cgalNPBegin{sample_faces} \anchor PMP_sample_faces
Parameter used in `sample_triangle_mesh()` to indicate if the interior of faces
should be considered for the sampling.
\n
\b Type : `bool` \n
\b Default value is `true`
\cgalNPEnd
\cgalNPBegin{number_of_points_on_faces} \anchor PMP_number_of_points_on_faces
Parameter used in `sample_triangle_mesh()` to set the number of points picked
using the random uniform method on faces.
\n
\b Type : `std::size_t` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPBegin{number_of_points_on_edges} \anchor PMP_number_of_points_on_edges
Parameter used in `sample_triangle_mesh()` to set the number of points picked
using the random uniform method on edges.
\n
\b Type : `std::size_t` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPBegin{number_of_points_per_face} \anchor PMP_number_of_points_per_face
Parameter used in `sample_triangle_mesh()` to set the number of points picked
per face using the Monte-Carlo method.
\n
\b Type : `std::size_t` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPBegin{number_of_points_per_edge} \anchor PMP_number_of_points_per_edge
Parameter used in `sample_triangle_mesh()` to set the number of points picked
per edge using the Monte-Carlo method.
\n
\b Type : `std::size_t` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPBegin{grid_spacing} \anchor PMP_grid_spacing
Parameter used in `sample_triangle_mesh()` to set the grid spacing when using
the grid sampling method.
\n
\b Type : `double` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPBegin{number_of_points_per_area_unit} \anchor PMP_number_of_points_per_area_unit
Parameter used in `sample_triangle_mesh()` to set the number of points per
area unit to be picked up in faces for the random uniform sampling and
Monte-Carlo methods.
\n
\b Type : `double` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPBegin{number_of_points_per_distance_unit} \anchor PMP_number_of_points_per_distance_unit
Parameter used in `sample_triangle_mesh()` to set the number of points per
distance unit to be picked up on edges for the random uniform sampling and
Monte-Carlo methods.
\n
\b Type : `double` \n
\b Default value is `0`
\cgalNPEnd
\cgalNPTableEnd

13
Polygon_mesh_processing/doc/Polygon_mesh_processing/PackageDescription.txt
View File

@ -40,6 +40,12 @@
/// Functions to orient polygon soups and to stitch geometrically identical boundaries.
/// \ingroup PkgPolygonMeshProcessing
/// \defgroup PMP_distance_grp Distance Functions
/// Functions to compute the distance between meshes, between a mesh and a point set and between a point set and a mesh.
/// \ingroup PkgPolygonMeshProcessing
/*!
\addtogroup PkgPolygonMeshProcessing
@ -126,6 +132,13 @@ and provides a list of the parameters that are used in this package.
- \link measure_grp `CGAL::Polygon_mesh_processing::edge_length()` \endlink
- \link measure_grp `CGAL::Polygon_mesh_processing::face_border_length()` \endlink
## Distance Functions ##
- `CGAL::Polygon_mesh_processing::approximate_Hausdorff_distance()`
- `CGAL::Polygon_mesh_processing::approximate_symmetric_Hausdorff_distance()`
- `CGAL::Polygon_mesh_processing::approximate_max_distance_to_point_set()`
- `CGAL::Polygon_mesh_processing::max_distance_to_triangle_mesh()`
- `CGAL::Polygon_mesh_processing::sample_triangle_mesh()`
## Miscellaneous ##
- `CGAL::Polygon_mesh_slicer`
- `CGAL::Side_of_triangle_mesh`

42
Polygon_mesh_processing/doc/Polygon_mesh_processing/Polygon_mesh_processing.txt
View File

@ -481,11 +481,53 @@ the propagation of a connected component index to cross it.
\cgalExample{Polygon_mesh_processing/connected_components_example.cpp}
\section PMPDistance Approximate Hausdorff Distance
This package provides methods to compute (approximate) distances between meshes and point sets.
The function \link CGAL::Polygon_mesh_processing::approximate_Hausdorff_distance() `approximate_Hausdorff_distance()`\endlink
computes an approximation of the Hausdorff distance from a mesh `tm1` to a mesh `tm2`. Given a
a sampling of `tm1`, it computes the distance to `tm2` of the farthest sample point to `tm2` \cgalCite{cignoni1998metro}.
The symmetric version (\link CGAL::Polygon_mesh_processing::approximate_symmetric_Hausdorff_distance() `approximate_symmetric_Hausdorff_distance()`\endlink)
is the maximum of the two non-symmetric distances. Internally, points are sampled using
\link CGAL::Polygon_mesh_processing::sample_triangle_mesh() `sample_triangle_mesh()`\endlink and the distance
to each sample point is computed using
\link CGAL::Polygon_mesh_processing::max_distance_to_triangle_mesh() `max_distance_to_triangle_mesh()`\endlink.
The quality of the approximation depends on the quality of the sampling and the runtime depends on the number of sample points.
Three sampling methods with different parameters are provided (see \cgalFigureRef{sampling_bunny}).
\cgalFigureBegin{sampling_bunny, pmp_sampling_bunny.jpg}
Sampling of a triangle mesh using different sampling methods. From left to right: (a) Grid sampling, (b) Monte-Carlo sampling with fixed number of points per face and per edge,
(c) Monte-Carlo sampling with a number of points proportional to the area/length, and (d) Uniform random sampling.
The four pictures represent the sampling on the same portion of a mesh, parameters were adjusted so that the total number of points sampled in faces (blue points) and on
edges (red points) are roughly the same.
Note that when using the random uniform sampling some faces/edges may not contain any point, but this method is the only one that allows to exactly match a given number of points.
\cgalFigureEnd
The function \link CGAL::Polygon_mesh_processing::approximate_max_distance_to_point_set() `approximate_max_distance_to_point_set()`\endlink
computes an approximation of the Hausdorff distance from a mesh to a point set.
For each triangle, a lower and upper bound of the Hausdorff distance to the point set are computed.
Triangles are refined until the difference between the bounds is lower than a user-defined precision threshold.
\subsection AHDExample Approximate Hausdorff Distance Example
In the following example, a mesh is isotropically remeshed and the approximate distance between the input and the output is computed.
\cgalExample{Polygon_mesh_processing/hausdorff_distance_remeshing_example.cpp}
\subsection PoissonDistanceExample Max Distance Between Point Set and Surface Example
In \ref Poisson_surface_reconstruction_3/poisson_reconstruction_example.cpp,
a triangulated surface mesh is constructed from a point set using the
\link PkgPoissonSurfaceReconstructionSummary Poisson reconstruction algorithm \endlink,
and the distance between the point set and the reconstructed surface is computed
with the following code:
\snippet Poisson_surface_reconstruction_3/poisson_reconstruction_example.cpp PMP_distance_snippet
\section PMPHistory Implementation History
A first version of this package was started by Ilker %O. Yaz and Sébastien Loriot.
Jane Tournois worked on the finalization of the API, code, and documentation.
*/
} /* namespace CGAL */

1
Polygon_mesh_processing/doc/Polygon_mesh_processing/dependencies
View File

@ -11,3 +11,4 @@ Surface_mesh
Surface_mesh_deformation
AABB_tree
Triangulation_2
Spatial_sorting

2
Polygon_mesh_processing/doc/Polygon_mesh_processing/examples.txt
View File

@ -13,4 +13,6 @@
\example Polygon_mesh_processing/mesh_slicer_example.cpp
\example Polygon_mesh_processing/isotropic_remeshing_example.cpp
\example Polygon_mesh_processing/compute_normals_example_Polyhedron.cpp
\example Polygon_mesh_processing/hausdorff_distance_remeshing_example.cpp
\example Poisson_surface_reconstruction_3/poisson_reconstruction_example.cpp
*/

BIN
Polygon_mesh_processing/doc/Polygon_mesh_processing/fig/pmp_sampling_bunny.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 123 KiB

34
Polygon_mesh_processing/examples/Polygon_mesh_processing/CMakeLists.txt
View File

@ -58,11 +58,7 @@ endif()
# include for local package
include_directories( BEFORE ../../include )
find_package(Eigen3 3.1.0) #(requires 3.1.0 or greater)
if (EIGEN3_FOUND)
# Executables that require Eigen 3.1
include( ${EIGEN3_USE_FILE} )
find_package(Eigen3 3.2.0) #(requires 3.2.0 or greater)
# Creating entries for all .cpp/.C files with "main" routine
@ -70,6 +66,23 @@ find_package(Eigen3 3.1.0) #(requires 3.1.0 or greater)
include( CGAL_CreateSingleSourceCGALProgram )
find_package( TBB )
if( TBB_FOUND )
include( ${TBB_USE_FILE} )
list( APPEND CGAL_3RD_PARTY_LIBRARIES ${TBB_LIBRARIES} )
else()
message( STATUS "NOTICE: Intel TBB was not found. hausdorff_distance_example will use sequential code." )
endif()
create_single_source_cgal_program( "hausdorff_distance_remeshing_example.cpp")
if (EIGEN3_FOUND)
# Executables that require Eigen 3.2
include( ${EIGEN3_USE_FILE} )
create_single_source_cgal_program( "hole_filling_example.cpp" )
create_single_source_cgal_program( "hole_filling_example_SM.cpp" )
create_single_source_cgal_program( "refine_fair_example.cpp")
endif(EIGEN3_FOUND)
create_single_source_cgal_program( "self_intersections_example.cpp" )
create_single_source_cgal_program( "stitch_borders_example.cpp" )
create_single_source_cgal_program( "compute_normals_example_Polyhedron.cpp" CXX_FEATURES cxx_range_for )
@ -83,13 +96,6 @@ create_single_source_cgal_program( "mesh_slicer_example.cpp")
#create_single_source_cgal_program( "remove_degeneracies_example.cpp")
create_single_source_cgal_program( "isotropic_remeshing_example.cpp")
if(NOT (${EIGEN3_VERSION} VERSION_LESS 3.2.0))
create_single_source_cgal_program( "hole_filling_example.cpp" )
create_single_source_cgal_program( "hole_filling_example_SM.cpp" )
create_single_source_cgal_program( "refine_fair_example.cpp")
else()
message(STATUS "NOTICE: Some examples require Eigen 3.2 (or greater) and will not be compiled.")
endif()
if(OpenMesh_FOUND)
create_single_source_cgal_program( "compute_normals_example_OM.cpp" )
@ -111,7 +117,5 @@ create_single_source_cgal_program( "triangulate_faces_example_OM.cpp")
target_link_libraries( triangulate_faces_example_OM ${OPENMESH_LIBRARIES} )
endif(OpenMesh_FOUND)
else(EIGEN3_FOUND)
message(STATUS "NOTICE: Some examples require Eigen 3.1 (or greater) and will not be compiled.")
endif(EIGEN3_FOUND)

34
Polygon_mesh_processing/examples/Polygon_mesh_processing/hausdorff_distance_remeshing_example.cpp Normal file
View File

@ -0,0 +1,34 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/boost/graph/graph_traits_Surface_mesh.h>
#include <CGAL/Polygon_mesh_processing/distance.h>
#include <CGAL/Polygon_mesh_processing/remesh.h>
#if CGAL_LINKED_WITH_TBB
#define TAG CGAL::Parallel_tag
#else
#define TAG CGAL::Sequential_tag
#endif
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef K::Point_3 Point;
typedef CGAL::Surface_mesh<K::Point_3> Mesh;
namespace PMP = CGAL::Polygon_mesh_processing;
int main()
{
Mesh tm1, tm2;
CGAL::make_tetrahedron(Point(.0,.0,.0),
Point(2,.0,.0),
Point(1,1,1),
Point(1,.0,2),
tm1);
tm2=tm1;
CGAL::Polygon_mesh_processing::isotropic_remeshing(tm2.faces(),.05, tm2);
std::cout << "Approximated Hausdorff distance: "
<< CGAL::Polygon_mesh_processing::approximate_Hausdorff_distance
<TAG>(tm1, tm2, PMP::parameters::number_of_points_per_area_unit(4000))
<< std::endl;
}

792
Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/distance.h
View File

@ -16,7 +16,7 @@
// $Id$
//
//
// Author(s) : Sebastien Loriot
// Author(s) : Maxime Gimeno and Sebastien Loriot
#ifndef CGAL_POLYGON_MESH_PROCESSING_DISTANCE_H
#define CGAL_POLYGON_MESH_PROCESSING_DISTANCE_H
@ -26,28 +26,42 @@
#include <CGAL/AABB_tree.h>
#include <CGAL/AABB_traits.h>
#include <CGAL/AABB_triangle_primitive.h>
#include <CGAL/AABB_face_graph_triangle_primitive.h>
#include <CGAL/utility.h>
#include <boost/foreach.hpp>
#include <CGAL/Polygon_mesh_processing/internal/named_function_params.h>
#include <CGAL/Polygon_mesh_processing/internal/named_params_helper.h>
#include <CGAL/point_generators_3.h>
/// \cond SKIP_IN_MANUAL
#include <CGAL/spatial_sort.h>
#include <CGAL/Polygon_mesh_processing/measure.h>
#include <CGAL/Polygon_mesh_processing/internal/mesh_to_point_set_hausdorff_distance.h>
#ifdef CGAL_LINKED_WITH_TBB
#include <tbb/parallel_for.h>
#include <tbb/blocked_range.h>
#include <CGAL/atomic.h>
#endif // CGAL_LINKED_WITH_TBB
#include <boost/foreach.hpp>
#include <boost/unordered_set.hpp>
namespace CGAL{
namespace Polygon_mesh_processing {
namespace internal{
template <class Kernel, class OutputIterator>
OutputIterator
triangle_grid_sampling(const typename Kernel::Point_3& p0,
const typename Kernel::Point_3& p1,
const typename Kernel::Point_3& p2,
double distance,
OutputIterator out)
triangle_grid_sampling( const typename Kernel::Point_3& p0,
const typename Kernel::Point_3& p1,
const typename Kernel::Point_3& p2,
double distance,
OutputIterator out)
{
const double d_p0p1 = CGAL::sqrt( CGAL::squared_distance(p0, p1) );
const double d_p0p2 = CGAL::sqrt( CGAL::squared_distance(p0, p2) );
typename Kernel::Compute_squared_distance_3 squared_distance;
const double d_p0p1 = to_double(approximate_sqrt( squared_distance(p0, p1) ));
const double d_p0p2 = to_double(approximate_sqrt( squared_distance(p0, p2) ));
const double n = (std::max)(std::ceil( d_p0p1 / distance ),
std::ceil( d_p0p2 / distance ));
std::ceil( d_p0p2 / distance ));
for (double i=1; i<n; ++i)
for (double j=1; j<n-i; ++j)
@ -62,121 +76,719 @@ triangle_grid_sampling(const typename Kernel::Point_3& p0,
return out;
}
template <class Kernel, class OutputIterator>
OutputIterator
triangle_grid_sampling(const typename Kernel::Triangle_3& t, double distance, OutputIterator out)
{
return triangle_grid_sampling<Kernel>(t[0], t[1], t[2], distance, out);
}
} //end of namespace internal
template <class Kernel, class TriangleRange, class OutputIterator>
template <class Kernel,
class FaceRange,
class TriangleMesh,
class VertexPointMap,
class OutputIterator>
OutputIterator
sample_triangles(const TriangleRange& triangles, double distance, OutputIterator out)
sample_triangles(const FaceRange& triangles,
const TriangleMesh& tm,
VertexPointMap vpm,
double distance,
OutputIterator out,
bool sample_faces,
bool sample_edges,
bool add_vertices)
{
typedef typename Kernel::Point_3 Point_3;
typedef typename boost::property_traits<VertexPointMap>::reference Point_ref;
typedef typename Kernel::Vector_3 Vector_3;
typedef typename Kernel::Triangle_3 Triangle_3;
typedef boost::graph_traits<TriangleMesh> GT;
typedef typename GT::face_descriptor face_descriptor;
typedef typename GT::halfedge_descriptor halfedge_descriptor;
std::set< std::pair<Point_3, Point_3> > sampled_edges;
boost::unordered_set<typename GT::edge_descriptor> sampled_edges;
boost::unordered_set<typename GT::vertex_descriptor> endpoints;
// sample edges but skip endpoints
BOOST_FOREACH(const Triangle_3& t, triangles)
BOOST_FOREACH(face_descriptor fd, triangles)
{
// sample edges but skip endpoints
halfedge_descriptor hd = halfedge(fd, tm);
for (int i=0;i<3; ++i)
{
const Point_3& p0=t[i];
const Point_3& p1=t[(i+1)%3];
if ( sampled_edges.insert(CGAL::make_sorted_pair(p0, p1)).second )
if (sample_edges && sampled_edges.insert(edge(hd, tm)).second )
{
const double d_p0p1 = CGAL::sqrt( CGAL::squared_distance(p0, p1) );
Point_ref p0 = get(vpm, source(hd, tm));
Point_ref p1 = get(vpm, target(hd, tm));
typename Kernel::Compute_squared_distance_3 squared_distance;
const double d_p0p1 = to_double(approximate_sqrt( squared_distance(p0, p1) ));
const double nb_pts = std::ceil( d_p0p1 / distance );
const Vector_3 step_vec = (p1 - p0) / nb_pts;
const Vector_3 step_vec = typename Kernel::Construct_scaled_vector_3()(
typename Kernel::Construct_vector_3()(p0, p1),
typename Kernel::FT(1)/typename Kernel::FT(nb_pts));
for (double i=1; i<nb_pts; ++i)
{
*out++=p0 + step_vec * i;
*out++=typename Kernel::Construct_translated_point_3()(p0,
typename Kernel::Construct_scaled_vector_3()(step_vec ,
typename Kernel::FT(i)));
}
}
//add endpoints once
if ( add_vertices && endpoints.insert(target(hd, tm)).second )
*out++=get(vpm, target(hd, tm));
hd=next(hd, tm);
}
// sample triangles
if (sample_faces)
{
Point_ref p0 = get(vpm, source(hd, tm));
Point_ref p1 = get(vpm, target(hd, tm));
Point_ref p2 = get(vpm, target(next(hd, tm), tm));
out=internal::triangle_grid_sampling<Kernel>(p0, p1, p2, distance, out);
}
}
return out;
}
#ifdef CGAL_LINKED_WITH_TBB
template <class AABB_tree, class Point_3>
struct Distance_computation{
const AABB_tree& tree;
const std::vector<Point_3>& sample_points;
Point_3 initial_hint;
cpp11::atomic<double>* distance;
Distance_computation(
const AABB_tree& tree,
const Point_3& p,
const std::vector<Point_3>& sample_points,
cpp11::atomic<double>* d)
: tree(tree)
, sample_points(sample_points)
, initial_hint(p)
, distance(d)
{}
void
operator()(const tbb::blocked_range<std::size_t>& range) const
{
Point_3 hint = initial_hint;
double hdist = 0;
for( std::size_t i = range.begin(); i != range.end(); ++i)
{
hint = tree.closest_point(sample_points[i], hint);
typename Kernel_traits<Point_3>::Kernel::Compute_squared_distance_3 squared_distance;
double d = to_double(CGAL::approximate_sqrt( squared_distance(hint,sample_points[i]) ));
if (d>hdist) hdist=d;
}
if (hdist > distance->load(cpp11::memory_order_acquire))
distance->store(hdist, cpp11::memory_order_release);
}
};
#endif
/** \ingroup PMP_distance_grp
* generates points taken on `tm` and outputs them to `out`, the sampling method
* is selected using named parameters.
* @tparam TriangleMesh a model of the concept `FaceListGraph`
* @tparam OutputIterator a model of `OutputIterator`
* holding objects of the same point type as
* the value type of the internal vertex point map of `tm`
*
* @param tm the triangle mesh that will be sampled
* @param out output iterator to be filled with sampled points
* @param np an optional sequence of \ref namedparameters among the ones listed below
*
* \cgalNamedParamsBegin
* \cgalParamBegin{vertex_point_map} the property map with the points
* associated to the vertices of `tm`. If this parameter is omitted,
* an internal property map for `CGAL::vertex_point_t`
* should be available for `TriangleMesh`.
* \cgalParamEnd
* \cgalParamBegin{geom_traits} a model of `PMPDistanceTraits`. \cgalParamEnd
* \cgalParamBegin{use_random_uniform_sampling}
* if `true` is passed (the default), points are generated in a random
* and uniform way on the surface of `tm`, and/or on edges of `tm`.
* For faces, the number of sample points is the value passed to the named
* parameter `number_of_points_on_faces()`. If not set,
* the value passed to the named parameter `number_of_points_per_area_unit()`
* is multiplied by the area of `tm` to get the number of sample points.
* If none of these parameters is set, the number of points sampled is `num_vertices(tm)`.
* For edges, the number of the number of sample points is the value passed to the named
* parameter `number_of_points_on_edges()`. If not set,
* the value passed to the named parameter `number_of_points_per_distance_unit()`
* is multiplied by the sum of the length of edges of `tm` to get the number of sample points.
* If none of these parameters is set, the number of points sampled is `num_vertices(tm)`.
* \cgalParamEnd
* \cgalParamBegin{use_grid_sampling}
* if `true` is passed, points are generated on a grid in each triangle,
* with a minimum of one point per triangle. The distance between
* two consecutive points in the grid is that of the length of the
* smallest non-null edge of `tm` or the value passed to
* the named parameter `grid_spacing()`. Edges are also split using the
* same distance, if requested.
* \cgalParamEnd
* \cgalParamBegin{use_monte_carlo_sampling}
* if `true` is passed, points are generated randomly in each triangle and/or
* on each edge.
* For faces, the number of points per triangle is the value passed to the named
* parameter `number_of_points_per_face()`. If not set, the value passed
* to the named parameter `number_of_points_per_area_unit()` is
* used to pick a number of points per face proportional to the triangle
* area with a minimum of one point per face. If none of these parameters
* is set, 2 divided by the square of the length of the smallest non-null
* edge of `tm` is used as if it was passed to
* `number_of_points_per_area_unit()`.
* For edges, the number of points per edge is the value passed to the named
* parameter `number_of_points_per_edge()`. If not set, the value passed
* to the named parameter `number_of_points_per_distance_unit()` is
* used to pick a number of points per edge proportional to the length of
* the edge with a minimum of one point per face. If none of these parameters
* is set, 1 divided by the length of the smallest non-null edge of `tm`
* is used as if it was passed to `number_of_points_per_distance_unit()`.
* \cgalParamEnd
* \cgalParamBegin{sample_vertices} if `true` is passed (default value),
* vertices of `tm` are put into `out`.\cgalParamEnd
* \cgalParamBegin{sample_edges} if `true` is passed (default value),
* edges of `tm` are sampled.\cgalParamEnd
* \cgalParamBegin{sample_faces} if `true` is passed (default value),
* faces of `tm` are sampled.\cgalParamEnd
* \cgalParamBegin{grid_spacing} a double value used as the grid spacing
* for the grid sampling method.
* \cgalParamEnd
* \cgalParamBegin{number_of_points_on_edges} an unsigned integral value used
* for the random sampling method as the number of points to pick exclusively
* on edges.
* \cgalParamEnd
* \cgalParamBegin{number_of_points_on_faces} an unsigned integral value used
* for the random sampling method as the number of points to pick on the surface.
* \cgalParamEnd
* \cgalParamBegin{number_of_points_per_distance_unit} a double value
* used for the random sampling and the Monte Carlo sampling methods to
* repectively determine the total number of points on edges and the
* number of points per edge.
* \cgalParamEnd
* \cgalParamBegin{number_of_points_per_edge} an unsigned integral value
* used by the Monte-Carlo sampling method as the number of points per edge
* to pick.
* \cgalParamEnd
* \cgalParamBegin{number_of_points_per_face} an unsigned integral value
* used by the Monte-Carlo sampling method as the number of points per face
* to pick.
* \cgalParamEnd
* \cgalParamBegin{number_of_points_per_area_unit} a double value
* used for the random sampling and the Monte Carlo sampling methods to
* repectively determine the total number of points inside faces
* and the number of points per face.
* \cgalParamEnd
* \cgalNamedParamsEnd
*/
template<class OutputIterator, class TriangleMesh, class NamedParameters>
OutputIterator
sample_triangle_mesh(const TriangleMesh& tm,
OutputIterator out,
NamedParameters np)
{
typedef typename GetGeomTraits<TriangleMesh,
NamedParameters>::type Geom_traits;
typedef typename GetVertexPointMap<TriangleMesh,
NamedParameters>::const_type Vpm;
typedef boost::graph_traits<TriangleMesh> GT;
typedef typename GT::face_descriptor face_descriptor;
typedef typename GT::halfedge_descriptor halfedge_descriptor;
typedef typename GT::edge_descriptor edge_descriptor;
using boost::choose_param;
using boost::get_param;
using boost::is_default_param;
Vpm pmap = choose_param(get_param(np, vertex_point),
get_const_property_map(vertex_point, tm));
typedef Creator_uniform_3<typename Geom_traits::FT,
typename Geom_traits::Point_3> Creator;
Geom_traits geomtraits = choose_param(get_param(np, geom_traits), Geom_traits());
bool use_rs = choose_param(get_param(np, random_uniform_sampling), true);
bool use_gs = choose_param(get_param(np, grid_sampling), false);
bool use_ms = choose_param(get_param(np, monte_carlo_sampling), false);
if (use_gs || use_ms)
if (is_default_param(get_param(np, random_uniform_sampling)))
use_rs=false;
bool smpl_vrtcs = choose_param(get_param(np, do_sample_vertices), true);
bool smpl_dgs = choose_param(get_param(np, do_sample_edges), true);
bool smpl_fcs = choose_param(get_param(np, do_sample_faces), true);
double nb_pts_a_u = choose_param(get_param(np, nb_points_per_area_unit), 0.);
double nb_pts_l_u = choose_param(get_param(np, nb_points_per_distance_unit), 0.);
// sample vertices
if (smpl_vrtcs)
{
Property_map_to_unary_function<Vpm> unary(pmap);
out = std::copy(
boost::make_transform_iterator(boost::begin(vertices(tm)), unary),
boost::make_transform_iterator(boost::end(vertices(tm)), unary),
out);
}
// grid sampling
if (use_gs)
{
double grid_spacing_ = choose_param(get_param(np, grid_spacing), 0.);
if (grid_spacing_==0.)
{
// set grid spacing to the shortest edge length
double grid_spacing_ = (std::numeric_limits<double>::max)();
typedef typename boost::graph_traits<TriangleMesh>
::edge_descriptor edge_descriptor;
BOOST_FOREACH(edge_descriptor ed, edges(tm))
{
double el = std::sqrt(
to_double( typename Geom_traits::Compute_squared_distance_3()(
get(pmap, source(ed, tm)), get(pmap, target(ed, tm)) )));
if (el > 0 && el < grid_spacing_)
grid_spacing_ = el;
}
}
out=sample_triangles<Geom_traits>(
faces(tm), tm, pmap, grid_spacing_, out,smpl_fcs, smpl_dgs, false);
}
// monte carlo sampling
if (use_ms)
{
typename Geom_traits::Compute_squared_distance_3 squared_distance;
double min_edge_length = (std::numeric_limits<double>::max)();
std::size_t nb_points_per_face =
choose_param(get_param(np, number_of_points_per_face), 0);
std::size_t nb_points_per_edge =
choose_param(get_param(np, number_of_points_per_edge), 0);
if ((nb_points_per_face == 0 && nb_pts_a_u ==0.) ||
(nb_points_per_edge == 0 && nb_pts_l_u ==0.) )
{
typedef typename boost::graph_traits<TriangleMesh>
::edge_descriptor edge_descriptor;
BOOST_FOREACH(edge_descriptor ed, edges(tm))
{
double el = std::sqrt(
to_double( squared_distance(get(pmap, source(ed, tm)),
get(pmap, target(ed, tm)) )));
if (min_edge_length > 0 && el < min_edge_length)
min_edge_length = el;
}
}
// sample faces
if (smpl_fcs)
{
// set default value
if (nb_points_per_face == 0 && nb_pts_a_u ==0.)
nb_pts_a_u = 2. / CGAL::square(min_edge_length);
BOOST_FOREACH(face_descriptor f, faces(tm))
{
std::size_t nb_points = nb_points_per_face;
if (nb_points == 0)
{
nb_points = (std::max)(
static_cast<std::size_t>(
std::ceil(to_double(
face_area(f,tm,parameters::geom_traits(geomtraits)))*nb_pts_a_u))
,std::size_t(1));
}
// extract triangle face points
typename Geom_traits::Point_3 points[3];
halfedge_descriptor hd(halfedge(f,tm));
for(int i=0; i<3; ++i)
{
points[i] = get(pmap, target(hd, tm));
hd = next(hd, tm);
}
// sample the triangle face
Random_points_in_triangle_3<typename Geom_traits::Point_3, Creator>
g(points[0], points[1], points[2]);
out=CGAL::cpp11::copy_n(g, nb_points, out);
}
}
// sample edges
if (smpl_dgs)
{
if (nb_points_per_edge == 0 && nb_pts_l_u == 0)
nb_pts_l_u = 1. / min_edge_length;
BOOST_FOREACH(edge_descriptor ed, edges(tm))
{
std::size_t nb_points = nb_points_per_edge;
if (nb_points == 0)
{
nb_points = (std::max)(
static_cast<std::size_t>( std::ceil( std::sqrt( to_double(
squared_distance(get(pmap, source(ed, tm)),
get(pmap, target(ed, tm)) )) )*nb_pts_l_u ) ),
std::size_t(1));
}
// now do the sampling of the edge
Random_points_on_segment_3<typename Geom_traits::Point_3, Creator>
g(get(pmap, source(ed,tm)), get(pmap, target(ed,tm)));
out=CGAL::cpp11::copy_n(g, nb_points, out);
}
}
}
// sample triangles
BOOST_FOREACH(const Triangle_3& t, triangles)
out=internal::triangle_grid_sampling<Kernel>(t, distance, out);
//add endpoints
std::set< Point_3 > endpoints;
BOOST_FOREACH(const Triangle_3& t, triangles)
for(int i=0; i<3; ++i)
// random uniform sampling
if (use_rs)
{
// sample faces
if(smpl_fcs)
{
if ( endpoints.insert(t[i]).second ) *out++=t[i];
std::size_t nb_points = choose_param(get_param(np, number_of_points_on_faces), 0);
Random_points_in_triangle_mesh_3<TriangleMesh, Vpm, Creator> g(tm, pmap);
if (nb_points == 0)
{
if (nb_pts_a_u == 0.)
nb_points = num_vertices(tm);
else
nb_points = static_cast<std::size_t>(
std::ceil(g.mesh_area()*nb_pts_a_u) );
}
out = CGAL::cpp11::copy_n(g, nb_points, out);
}
// sample edges
if (smpl_dgs)
{
std::size_t nb_points =
choose_param(get_param(np, number_of_points_on_edges), 0);
Random_points_on_edge_list_graph_3<TriangleMesh, Vpm, Creator> g(tm, pmap);
if (nb_points == 0)
{
if (nb_pts_l_u == 0)
nb_points = num_vertices(tm);
else
nb_points = static_cast<std::size_t>(
std::ceil( g.mesh_length()*nb_pts_a_u) );
}
out = CGAL::cpp11::copy_n(g, nb_points, out);
}
}
return out;
}
template <class Kernel, class TriangleRange1, class TriangleRange2>
double approximated_Hausdorff_distance(
const TriangleRange1& triangles_1,
const TriangleRange2& triangles_2,
double targeted_precision
)
template<class OutputIterator, class TriangleMesh>
OutputIterator
sample_triangle_mesh(const TriangleMesh& tm,
OutputIterator out)
{
std::vector<typename Kernel::Point_3> sample_points;
sample_triangles<Kernel>( triangles_1,
targeted_precision,
std::back_inserter(sample_points) );
// std::ofstream out("/tmp/samples.xyz");
// std::copy(sample_points.begin(), sample_points.end(), std::ostream_iterator<typename Kernel::Point_3>(out,"\n"));
typedef typename TriangleRange2::const_iterator Triangle_iterator;
typedef AABB_triangle_primitive<Kernel, Triangle_iterator> Primitive;
typedef AABB_traits<Kernel, Primitive> Traits;
typedef AABB_tree< Traits > Tree;
Tree tree(triangles_2.begin(), triangles_2.end());
tree.accelerate_distance_queries();
double hdist = 0;
BOOST_FOREACH(const typename Kernel::Point_3& pt, sample_points)
{
double d = CGAL::sqrt( tree.squared_distance(pt) );
// if ( d > 1e-1 ) std::cout << pt << "\n";
if (d>hdist) hdist=d;
}
return hdist;
return sample_triangle_mesh(tm, out, parameters::all_default());
}
template <class Kernel, class TriangleRange1, class TriangleRange2>
double approximated_symmetric_Hausdorff_distance(
const TriangleRange1& triangles_1,
const TriangleRange2& triangles_2,
double targeted_precision
)
template <class Concurrency_tag,
class Kernel,
class PointRange,
class TriangleMesh,
class VertexPointMap>
double approximate_Hausdorff_distance(
const PointRange& original_sample_points,
const TriangleMesh& tm,
VertexPointMap vpm)
{
CGAL_assertion_code( bool is_triangle = is_triangle_mesh(tm) );
CGAL_assertion_msg (is_triangle,
"Mesh is not triangulated. Distance computing impossible.");
#ifdef CGAL_HAUSDORFF_DEBUG
std::cout << "Nb sample points " << sample_points.size() << "\n";
#endif
typedef typename Kernel::Point_3 Point_3;
std::vector<Point_3> sample_points
(boost::begin(original_sample_points), boost::end(original_sample_points) );
spatial_sort(sample_points.begin(), sample_points.end());
typedef AABB_face_graph_triangle_primitive<TriangleMesh> Primitive;
typedef AABB_tree< AABB_traits<Kernel, Primitive> > Tree;
Tree tree( faces(tm).first, faces(tm).second, tm);
tree.accelerate_distance_queries();
tree.build();
Point_3 hint = get(vpm, *vertices(tm).first);
#ifndef CGAL_LINKED_WITH_TBB
CGAL_static_assertion_msg (!(boost::is_convertible<Concurrency_tag, Parallel_tag>::value),
"Parallel_tag is enabled but TBB is unavailable.");
#else
if (boost::is_convertible<Concurrency_tag,Parallel_tag>::value)
{
cpp11::atomic<double> distance(0);
Distance_computation<Tree, Point_3> f(tree, hint, sample_points, &distance);
tbb::parallel_for(tbb::blocked_range<std::size_t>(0, sample_points.size()), f);
return distance;
}
else
#endif
{
double hdist = 0;
BOOST_FOREACH(const Point_3& pt, sample_points)
{
hint = tree.closest_point(pt, hint);
typename Kernel::Compute_squared_distance_3 squared_distance;
typename Kernel::FT dist = squared_distance(hint,pt);
double d = to_double(CGAL::approximate_sqrt(dist));
if(d>hdist)
hdist=d;
}
return hdist;
}
}
template <class Concurrency_tag, class Kernel, class TriangleMesh,
class NamedParameters,
class VertexPointMap >
double approximate_Hausdorff_distance(
const TriangleMesh& tm1,
const TriangleMesh& tm2,
NamedParameters np,
VertexPointMap vpm_2)
{
std::vector<typename Kernel::Point_3> sample_points;
sample_triangle_mesh(
tm1,
std::back_inserter(sample_points),
np);
return approximate_Hausdorff_distance<Concurrency_tag, Kernel>(sample_points, tm2, vpm_2);
}
// documented functions
/**
* \ingroup PMP_distance_grp
* computes the approximate Hausdorff distance from `tm1` to `tm2` by returning
* the distance of the farthest point from `tm2` amongst a sampling of `tm1`
* generated with the function `sample_triangle_mesh()` with
* `tm1` and `np1` as parameter.
*
* A parallel version is provided and requires the executable to be
* linked against the <a href="http://www.threadingbuildingblocks.org">Intel TBB library</a>.
* To control the number of threads used, the user may use the `tbb::task_scheduler_init` class.
* See the <a href="http://www.threadingbuildingblocks.org/documentation">TBB documentation</a>
* for more details.
*
* @tparam Concurrency_tag enables sequential versus parallel algorithm.
* Possible values are `Sequential_tag`
* and `Parallel_tag`.
* @tparam TriangleMesh a model of the concept `FaceListGraph`
* @tparam NamedParameters1 a sequence of \ref namedparameters for `tm1`
* @tparam NamedParameters2 a sequence of \ref namedparameters for `tm2`
*
* @param tm1 the triangle mesh that will be sampled
* @param tm2 the triangle mesh to compute the distance to
* @param np1 optional sequence of \ref namedparameters for `tm1` passed to `sample_triangle_mesh()`.
*
* @param np2 optional sequence of \ref namedparameters for `tm2` among the ones listed below
*
* \cgalNamedParamsBegin
* \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `tm2`
* If this parameter is omitted, an internal property map for CGAL::vertex_point_t should be available in `TriangleMesh`
* and in all places where vertex_point_map is used.
* \cgalParamEnd
* \cgalNamedParamsEnd
* The function `CGAL::Polygon_mesh_processing::params::all_default()` can be used to indicate to use the default values for
* `np1` and specify custom values for `np2`
*/
template< class Concurrency_tag,
class TriangleMesh,
class NamedParameters1,
class NamedParameters2>
double approximate_Hausdorff_distance( const TriangleMesh& tm1,
const TriangleMesh& tm2,
const NamedParameters1& np1,
const NamedParameters2& np2)
{
typedef typename GetGeomTraits<TriangleMesh,
NamedParameters1>::type Geom_traits;
return approximate_Hausdorff_distance<Concurrency_tag, Geom_traits>(
tm1, tm2, np1, choose_param(get_param(np2, vertex_point),
get_const_property_map(vertex_point, tm2)));
}
/**
* \ingroup PMP_distance_grp
* computes the approximate symmetric Hausdorff distance between `tm1` and `tm2`.
* It returns the maximum of `approximate_Hausdorff_distance(tm1, tm2, np1, np2)`
* and `approximate_Hausdorff_distance(tm2, tm1, np2, np1)`.
*/
template< class Concurrency_tag,
class TriangleMesh,
class NamedParameters1,
class NamedParameters2>
double approximate_symmetric_Hausdorff_distance(
const TriangleMesh& tm1,
const TriangleMesh& tm2,
const NamedParameters1& np1,
const NamedParameters2& np2)
{
return (std::max)(
approximated_Hausdorff_distance<Kernel>(triangles_1, triangles_2, targeted_precision),
approximated_Hausdorff_distance<Kernel>(triangles_2, triangles_1, targeted_precision)
approximate_Hausdorff_distance<Concurrency_tag>(tm1,tm2,np1,np2),
approximate_Hausdorff_distance<Concurrency_tag>(tm2,tm1,np2,np1)
);
}
/**
* \ingroup PMP_distance_grp
* returns the distance to `tm` of the point from `points`
* that is the furthest from `tm`.
* @tparam PointRange a range of `Point_3`, model of `Range`.
* @tparam TriangleMesh a model of the concept `FaceListGraph`
* @tparam NamedParameters a sequence of \ref namedparameters
* @param points the range of points of interest
* @param tm the triangle mesh to compute the distance to
* @param np an optional sequence of \ref namedparameters among the ones listed below
*
* \cgalNamedParamsBegin
* \cgalParamBegin{vertex_point_map}
* the property map with the points associated to the vertices of `tm`. If this parameter is omitted,
* an internal property map for `CGAL::vertex_point_t` should be available for the
vertices of `tm` \cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `PMPDistanceTraits`\cgalParamEnd
* \cgalNamedParamsEnd
*/
template< class Concurrency_tag,
class TriangleMesh,
class PointRange,
class NamedParameters>
double max_distance_to_triangle_mesh(const PointRange& points,
const TriangleMesh& tm,
const NamedParameters& np)
{
typedef typename GetGeomTraits<TriangleMesh,
NamedParameters>::type Geom_traits;
///\todo add approximated_Hausdorff_distance(FaceGraph, FaceGraph)
///\todo add approximated_Hausdorff_distance(TriangleRange, FaceGraph)
///\todo add approximated_Hausdorff_distance(FaceGraph, TriangleRange)
///\todo add approximated_Hausdorff_distance(PointRange, FaceGraph)
///\todo add approximated_Hausdorff_distance(PointRange, TriangleRange)
///\todo add barycentric_triangle_sampling(Triangle, PointPerAreaUnit)
///\todo add random_triangle_sampling(Triangle, PointPerAreaUnit)
/// \todo maybe we should use a sampling using epec to avoid being too far from the sampled triangle
return approximate_Hausdorff_distance<Concurrency_tag, Geom_traits>
(points,tm,choose_param(get_param(np, vertex_point),
get_const_property_map(vertex_point, tm)));
}
/*!
*\ingroup PMP_distance_grp
* returns an approximation of the distance between `points` and the point lying on `tm` that is the farthest from `points`
* @tparam PointRange a range of `Point_3`, model of `Range`.
* @tparam TriangleMesh a model of the concept `FaceListGraph`
* @tparam NamedParameters a sequence of \ref namedparameters
* @param tm a triangle mesh
* @param points a range of points
* @param precision for each triangle of `tm`, the distance of its farthest point from `points` is bounded.
* A triangle is subdivided into sub-triangles so that the difference of its distance bounds
* is smaller than `precision`. `precision` must be strictly positive to avoid infinite loops.
* @param np an optional sequence of \ref namedparameters among the ones listed below
*
* \cgalNamedParamsBegin
* \cgalParamBegin{vertex_point_map}
* the property map with the points associated to the vertices of `tm`. If this parameter is omitted,
* an internal property map for `CGAL::vertex_point_t` should be available for the
vertices of `tm` \cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `PMPDistanceTraits`. \cgalParamEnd
* \cgalNamedParamsEnd
*/
template< class TriangleMesh,
class PointRange,
class NamedParameters>
double approximate_max_distance_to_point_set(const TriangleMesh& tm,
const PointRange& points,
const double precision,
const NamedParameters& np)
{
typedef typename GetGeomTraits<TriangleMesh,
NamedParameters>::type Geom_traits;
typedef boost::graph_traits<TriangleMesh> GT;
typedef Orthogonal_k_neighbor_search<Search_traits_3<Geom_traits> > Knn;
typedef typename Knn::Tree Tree;
Tree tree(points.begin(), points.end());
CRefiner<Geom_traits> ref;
BOOST_FOREACH(typename GT::face_descriptor f, faces(tm))
{
typename Geom_traits::Point_3 points[3];
typename GT::halfedge_descriptor hd(halfedge(f,tm));
for(int i=0; i<3; ++i)
{
points[i] = get(choose_param(get_param(np, vertex_point),
get_const_property_map(vertex_point, tm)),
target(hd, tm));
hd = next(hd, tm);
}
ref.add(points[0], points[1], points[2], tree);
}
return to_double(ref.refine(precision, tree));
}
// convenience functions with default parameters
template< class Concurrency_tag,
class TriangleMesh,
class PointRange>
double max_distance_to_triangle_mesh(const PointRange& points,
const TriangleMesh& tm)
{
return max_distance_to_triangle_mesh<Concurrency_tag,
TriangleMesh,
PointRange>
(points, tm, parameters::all_default());
}
template< class TriangleMesh,
class PointRange>
double approximate_max_distance_to_point_set(const TriangleMesh& tm,
const PointRange& points,
const double precision)
{
return approximate_max_distance_to_point_set(tm, points, precision,
parameters::all_default());
}
template< class Concurrency_tag,
class TriangleMesh,
class NamedParameters>
double approximate_Hausdorff_distance(const TriangleMesh& tm1,
const TriangleMesh& tm2,
const NamedParameters& np)
{
return approximate_Hausdorff_distance<Concurrency_tag>(
tm1, tm2, np, parameters::all_default());
}
template< class Concurrency_tag,
class TriangleMesh>
double approximate_Hausdorff_distance(const TriangleMesh& tm1,
const TriangleMesh& tm2)
{
return approximate_Hausdorff_distance<Concurrency_tag>(
tm1, tm2, parameters::all_default(), parameters::all_default());
}
template< class Concurrency_tag,
class TriangleMesh,
class NamedParameters>
double approximate_symmetric_Hausdorff_distance(const TriangleMesh& tm1,
const TriangleMesh& tm2,
const NamedParameters& np)
{
return approximate_symmetric_Hausdorff_distance<Concurrency_tag>(
tm1, tm2, np, parameters::all_default());
}
} } // end of namespace CGAL::Polygon_mesh_processing
template< class Concurrency_tag,
class TriangleMesh>
double approximate_symmetric_Hausdorff_distance(const TriangleMesh& tm1,
const TriangleMesh& tm2)
{
return approximate_symmetric_Hausdorff_distance<Concurrency_tag>(
tm1, tm2, parameters::all_default(), parameters::all_default());
}
}
} // end of namespace CGAL::Polygon_mesh_processing
/// \endcond
#endif //CGAL_POLYGON_MESH_PROCESSING_DISTANCE_H

341
Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/mesh_to_point_set_hausdorff_distance.h Normal file
View File

@ -0,0 +1,341 @@
// Copyright (c) 2016 INRIA Sophia-Antipolis (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
// You can redistribute it and/or modify it under the terms of the GNU
// General Public License as published by the Free Software Foundation,
// either version 3 of the License, or (at your option) any later version.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL$
// $Id$
//
//
// Author(s) : Simon Giraudot and Maxime Gimeno
#ifndef CGAL_MESH_TO_POINT_SET_HAUSDORFF_DISTANCE_H
#define CGAL_MESH_TO_POINT_SET_HAUSDORFF_DISTANCE_H
#include <CGAL/Orthogonal_k_neighbor_search.h>
#include <CGAL/Search_traits_3.h>
#include <CGAL/squared_distance_3.h>
#include <queue>
#include <iostream>
namespace CGAL{
namespace internal{
template <typename Kernel>
class CPointH
{
public:
typedef typename Kernel::Point_3 Point;
typedef typename Kernel::FT FT;
CPointH ()
{
m_point = Point (FT(0.), FT(0.), FT(0.));
m_hausdorff = 0.;
}
CPointH (const Point& p, FT h = 0.)
{
m_point = p;
m_hausdorff = h;
}
CPointH (const CPointH& p)
{
m_point = p ();
m_hausdorff = p.hausdorff ();
}
const Point& operator() () const { return m_point; }
Point& operator() () { return m_point; }
FT hausdorff () const { return m_hausdorff; }
static Point mid_point (const CPointH& a, const CPointH& b)
{
return Point (FT(0.5) * (a().x() + b().x()),
FT(0.5) * (a().y() + b().y()),
FT(0.5) * (a().z() + b().z()));
}
private:
Point m_point;
FT m_hausdorff;
};
template <typename Kernel>
class CRefTriangle
{
private:
typedef typename Kernel::FT FT;
typedef CPointH<Kernel> PointH;
typedef typename PointH::Point Point;
PointH m_point[3];
std::size_t m_edge;
FT m_upper_bound;
FT m_lower_bound;
Point m_bisector;
public:
CRefTriangle (const PointH& a, const PointH& b, const PointH& c) // Triangle
{
m_point[0] = a;
m_point[1] = b;
m_point[2] = c;
m_edge = 0;
typename Kernel::Compute_squared_distance_3 squared_distance;
FT length_max = squared_distance(m_point[1](), m_point[2]());
FT length1 = squared_distance(m_point[2](), m_point[0]());
if (length1 > length_max)
{
m_edge = 1;
length_max = length1;
}
FT length2 = squared_distance(m_point[0](), m_point[1]());
if (length2 > length_max)
m_edge = 2;
m_bisector = PointH::mid_point (m_point[(m_edge+1)%3],
m_point[(m_edge+2)%3]);
m_lower_bound = 0.;
m_upper_bound = 0.;
for (unsigned int i = 0; i < 3; ++ i)
{
if (m_point[i].hausdorff () > m_lower_bound)
m_lower_bound = m_point[i].hausdorff ();
FT up = m_point[i].hausdorff ()
+ CGAL::approximate_sqrt (squared_distance (m_point[i](), m_bisector));
if (up > m_upper_bound)
m_upper_bound = up;
}
}
CRefTriangle (const CRefTriangle& t)
{
m_point[0] = t.points ()[0];
m_point[1] = t.points ()[1];
m_point[2] = t.points ()[2];
m_edge = t.edge ();
m_lower_bound = t.lower_bound ();
m_upper_bound = t.upper_bound ();
m_bisector = t.bisector ();
}
FT lower_bound () const
{
return m_lower_bound;
}
FT upper_bound () const
{
return m_upper_bound;
}
const Point& bisector () const
{
return m_bisector;
}
friend bool operator< (const CRefTriangle& a, const CRefTriangle& b)
{
return a.upper_bound () < b.upper_bound ();
}
const PointH* points () const { return m_point; }
std::size_t edge () const { return m_edge; }
#ifdef CGAL_MTPS_HD_DEBUG
void print () const
{
std::cerr << "[Refinement triangle]" << std::endl
<< " Bounds: " << m_lower_bound << " to "
<< m_upper_bound << std::endl
<< " * " << m_point[0]() << std::endl
<< " * " << m_point[1]() << std::endl
<< " * " << m_point[2]() << std::endl
<< " -> " << m_bisector << std::endl;
}
#endif
};
}//internal
template <typename Kernel>
class CRefiner
{
private:
typedef typename Kernel::Point_3 Point;
typedef typename Kernel::Triangle_3 Triangle;
typedef typename Kernel::FT FT;
typedef internal::CPointH<Kernel> PointH;
typedef internal::CRefTriangle<Kernel> RefTriangle;
typedef CGAL::Orthogonal_k_neighbor_search<CGAL::Search_traits_3<Kernel> > Knn;
typedef typename Knn::Tree Tree;
std::priority_queue<RefTriangle> m_queue;
FT m_lower_bound;
FT m_upper_bound;
Point m_point;
public:
CRefiner ()
{
m_lower_bound = 0.;
m_upper_bound = (std::numeric_limits<FT>::max) ();
}
bool empty ()
{
return m_queue.empty ();
}
void reset (FT lower_bound = 0.)
{
m_queue = std::priority_queue<RefTriangle> ();
m_lower_bound = lower_bound;
m_upper_bound = (std::numeric_limits<FT>::max) ();
}
inline FT uncertainty () const
{
return m_upper_bound - m_lower_bound;
}
inline FT mid () const
{
return ((m_upper_bound + m_lower_bound) / 2.);
}
inline FT lower_bound () const
{
return m_lower_bound;
}
inline FT upper_bound () const
{
return m_upper_bound;
}
inline FT relative_error () const
{
return 0.5 * uncertainty () / mid ();
}
const Point& hausdorff_point () const { return m_point; }
bool add (const Point& a, const Point& b, const Point& c, const Tree& tree)
{
RefTriangle r (PointH (a, CGAL::approximate_sqrt (Knn(tree, a, 1).begin()->second)),
PointH (b, CGAL::approximate_sqrt (Knn(tree, b, 1).begin()->second)),
PointH (c, CGAL::approximate_sqrt (Knn(tree, c, 1).begin()->second)));
if (r.lower_bound () > m_lower_bound)
{
m_lower_bound = r.lower_bound ();
}
if (r.upper_bound () > m_lower_bound)
m_queue.push (r);
return true;
}
bool clean_up_queue ()
{
std::size_t before = m_queue.size ();
std::vector<RefTriangle> to_keep;
while (!(m_queue.empty ()))
{
const RefTriangle& current = m_queue.top ();
if (current.upper_bound () > m_lower_bound)
to_keep.push_back (current);
m_queue.pop ();
}
m_queue = std::priority_queue<RefTriangle> ();
BOOST_FOREACH(RefTriangle& r, to_keep)
m_queue.push (r);
return (m_queue.size () < before);
}
FT refine (FT limit, const Tree& tree, FT upper_bound = 1e30)
{
unsigned int nb_clean = 0;
while (uncertainty () > limit && !(m_queue.empty ()))
{
if (m_queue.size () > 100000)
{
#ifdef CGAL_MTPS_HD_DEBUG
m_queue.top ().print ();
#endif
++ nb_clean;
if (nb_clean > 5)
return m_upper_bound;
if (!clean_up_queue ())
return m_upper_bound;
}
const RefTriangle& current = m_queue.top ();
m_upper_bound = current.upper_bound ();
if (current.lower_bound () > m_lower_bound)
m_lower_bound = current.lower_bound ();
typename Kernel::Compute_squared_area_3 squared_area;
if(squared_area(current.points()[0](),
current.points()[1](),
current.points()[2]()
) < 1e-20)
{
m_queue.pop();
continue;
}
const Point& bisector = current.bisector ();
m_point = bisector;
//squared distance between bisector and its closst point in the mesh
FT hausdorff = CGAL::approximate_sqrt (Knn(tree, bisector, 1).begin()->second);
if (hausdorff > m_lower_bound)
m_lower_bound = hausdorff;
if (m_lower_bound > upper_bound)
return m_upper_bound;
PointH new_point (bisector, hausdorff);
std::size_t i = current.edge ();
PointH p0 (current.points()[i]);
PointH p1 (current.points()[(i+1)%3]);
PointH p2 (current.points()[(i+2)%3]);
m_queue.pop ();
m_queue.push (RefTriangle (new_point, p0, p1));
m_queue.push (RefTriangle (new_point, p0, p2));
}
return m_upper_bound;
}
};
}//CGAL
#endif // CGAL_MESH_TO_POINT_SET_HAUSDORFF_DISTANCE_H

235
Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/internal/named_function_params.h
View File

@ -14,7 +14,7 @@
//
// $URL$
// $Id$
//
//
//
// Author(s) : Jane Tournois
@ -39,6 +39,19 @@ namespace CGAL{
enum relax_constraints_t { relax_constraints };
enum vertex_is_constrained_t { vertex_is_constrained };
enum face_patch_t { face_patch };
enum random_uniform_sampling_t { random_uniform_sampling };
enum grid_sampling_t { grid_sampling };
enum monte_carlo_sampling_t { monte_carlo_sampling };
enum do_sample_edges_t { do_sample_edges };
enum do_sample_vertices_t { do_sample_vertices };
enum do_sample_faces_t { do_sample_faces };
enum number_of_points_on_faces_t { number_of_points_on_faces };
enum number_of_points_per_face_t { number_of_points_per_face };
enum grid_spacing_t { grid_spacing };
enum nb_points_per_area_unit_t { nb_points_per_area_unit };
enum number_of_points_per_edge_t { number_of_points_per_edge };
enum number_of_points_on_edges_t { number_of_points_on_edges };
enum nb_points_per_distance_unit_t{ nb_points_per_distance_unit };
//to be documented
enum face_normal_t { face_normal };
@ -216,6 +229,109 @@ namespace CGAL{
return Params(p, *this);
}
template <typename Boolean>
pmp_bgl_named_params<Boolean, random_uniform_sampling_t, self>
use_random_uniform_sampling(const Boolean& p) const
{
typedef pmp_bgl_named_params<Boolean, random_uniform_sampling_t, self> Params;
return Params(p, *this);
}
template <typename Boolean>
pmp_bgl_named_params<Boolean, grid_sampling_t, self>
use_grid_sampling(const Boolean& p) const
{
typedef pmp_bgl_named_params<Boolean, grid_sampling_t, self> Params;
return Params(p, *this);
}
template <typename Boolean>
pmp_bgl_named_params<Boolean, monte_carlo_sampling_t, self>
use_monte_carlo_sampling(const Boolean& p) const
{
typedef pmp_bgl_named_params<Boolean, monte_carlo_sampling_t, self> Params;
return Params(p, *this);
}
template <typename Boolean>
pmp_bgl_named_params<Boolean, do_sample_edges_t, self>
sample_edges(const Boolean& p) const
{
typedef pmp_bgl_named_params<Boolean, do_sample_edges_t, self> Params;
return Params(p, *this);
}
template <typename Boolean>
pmp_bgl_named_params<Boolean, do_sample_vertices_t, self>
sample_vertices(const Boolean& p) const
{
typedef pmp_bgl_named_params<Boolean, do_sample_vertices_t, self> Params;
return Params(p, *this);
}
template <typename Boolean>
pmp_bgl_named_params<Boolean, do_sample_faces_t, self>
sample_faces(const Boolean& p) const
{
typedef pmp_bgl_named_params<Boolean, do_sample_faces_t, self> Params;
return Params(p, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_on_faces_t, self>
number_of_points_on_faces(const NT& n) const
{
typedef pmp_bgl_named_params<NT, number_of_points_on_faces_t, self> Params;
return Params(n, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_per_face_t, self>
number_of_points_per_face(const NT& n) const
{
typedef pmp_bgl_named_params<NT, number_of_points_per_face_t, self> Params;
return Params(n, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, grid_spacing_t, self>
grid_spacing(const NT& n) const
{
typedef pmp_bgl_named_params<NT, grid_spacing_t, self> Params;
return Params(n, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, nb_points_per_area_unit_t, self>
number_of_points_per_area_unit(const NT& n) const
{
typedef pmp_bgl_named_params<NT, nb_points_per_area_unit_t, self> Params;
return Params(n, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_per_edge_t, self>
number_of_points_per_edge(const NT& n) const
{
typedef pmp_bgl_named_params<NT, number_of_points_per_edge_t, self> Params;
return Params(n, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_on_edges_t, self>
number_of_points_on_edges(const NT& n) const
{
typedef pmp_bgl_named_params<NT, number_of_points_on_edges_t, self> Params;
return Params(n, *this);
}
template<typename NT>
pmp_bgl_named_params<NT, nb_points_per_distance_unit_t, self>
number_of_points_per_distance_unit(const NT& n) const
{
typedef pmp_bgl_named_params<NT, nb_points_per_distance_unit_t, self> Params;
return Params(n, *this);
}
};
namespace Polygon_mesh_processing{
@ -379,6 +495,123 @@ namespace parameters{
return Params(p);
}
//overload
template <typename Boolean>
pmp_bgl_named_params<Boolean, random_uniform_sampling_t>
use_random_uniform_sampling(const Boolean& p)
{
typedef pmp_bgl_named_params<Boolean, random_uniform_sampling_t> Params;
return Params(p);
}
//overload
template <typename Boolean>
pmp_bgl_named_params<Boolean, grid_sampling_t>
use_grid_sampling(const Boolean& p)
{
typedef pmp_bgl_named_params<Boolean, grid_sampling_t> Params;
return Params(p);
}
//overload
template <typename Boolean>
pmp_bgl_named_params<Boolean, monte_carlo_sampling_t>
use_monte_carlo_sampling(const Boolean& p)
{
typedef pmp_bgl_named_params<Boolean, monte_carlo_sampling_t> Params;
return Params(p);
}
//overload
template <typename Boolean>
pmp_bgl_named_params<Boolean, do_sample_edges_t>
sample_edges(const Boolean& p)
{
typedef pmp_bgl_named_params<Boolean, do_sample_edges_t> Params;
return Params(p);
}
//overload
template <typename Boolean>
pmp_bgl_named_params<Boolean, do_sample_vertices_t>
sample_vertices(const Boolean& p)
{
typedef pmp_bgl_named_params<Boolean, do_sample_vertices_t> Params;
return Params(p);
}
//overload
template <typename Boolean>
pmp_bgl_named_params<Boolean, do_sample_faces_t>
sample_faces(const Boolean& p)
{
typedef pmp_bgl_named_params<Boolean, do_sample_faces_t> Params;
return Params(p);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_on_faces_t>
number_of_points_on_faces(const NT& n)
{
typedef pmp_bgl_named_params<NT, number_of_points_on_faces_t> Params;
return Params(n);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_per_face_t>
number_of_points_per_face(const NT& n)
{
typedef pmp_bgl_named_params<NT, number_of_points_per_face_t> Params;
return Params(n);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, grid_spacing_t>
grid_spacing(const NT& n)
{
typedef pmp_bgl_named_params<NT, grid_spacing_t> Params;
return Params(n);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, nb_points_per_area_unit_t>
number_of_points_per_area_unit(const NT& n)
{
typedef pmp_bgl_named_params<NT, nb_points_per_area_unit_t> Params;
return Params(n);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_per_edge_t>
number_of_points_per_edge(const NT& n)
{
typedef pmp_bgl_named_params<NT, number_of_points_per_edge_t> Params;
return Params(n);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, number_of_points_on_edges_t>
number_of_points_on_edges(const NT& n)
{
typedef pmp_bgl_named_params<NT, number_of_points_on_edges_t> Params;
return Params(n);
}
//overload
template<typename NT>
pmp_bgl_named_params<NT, nb_points_per_distance_unit_t>
number_of_points_per_distance_unit(const NT& n)
{
typedef pmp_bgl_named_params<NT, nb_points_per_distance_unit_t> Params;
return Params(n);
}
} //namespace parameters
} //namespace Polygon_mesh_processing

56
Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/measure.h
View File

@ -14,7 +14,7 @@
//
// $URL$
// $Id$
//
//
//
// Author(s) : Andreas Fabri
@ -413,39 +413,43 @@ namespace Polygon_mesh_processing {
#else
typename GetGeomTraits<TriangleMesh, CGAL_PMP_NP_CLASS>::type::FT
#endif
volume(const TriangleMesh& tmesh, const CGAL_PMP_NP_CLASS& np)
volume(const TriangleMesh& tmesh, const CGAL_PMP_NP_CLASS& np)
{
CGAL_assertion(is_triangle_mesh(tmesh));
CGAL_assertion(is_closed(tmesh));
using boost::choose_param;
using boost::get_param;
typename GetVertexPointMap<TriangleMesh, CGAL_PMP_NP_CLASS>::const_type
vpm = choose_param(get_param(np, vertex_point),
get_const_property_map(CGAL::vertex_point, tmesh));
typename GetGeomTraits<TriangleMesh, CGAL_PMP_NP_CLASS>::type::Point_3
origin(0, 0, 0);
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typename GetGeomTraits<TriangleMesh, CGAL_PMP_NP_CLASS>::type::FT volume = 0.;
typename CGAL::Kernel_traits<typename property_map_value<TriangleMesh,
CGAL::vertex_point_t>::type>::Kernel::Compute_volume_3 cv3;
BOOST_FOREACH(face_descriptor f, faces(tmesh))
{
CGAL_assertion(is_triangle_mesh(tmesh));
CGAL_assertion(is_closed(tmesh));
using boost::choose_param;
using boost::get_param;
typename GetVertexPointMap<TriangleMesh, CGAL_PMP_NP_CLASS>::const_type
vpm = choose_param(get_param(np, vertex_point),
get_const_property_map(CGAL::vertex_point, tmesh));
typename GetGeomTraits<TriangleMesh, CGAL_PMP_NP_CLASS>::type::Point_3
origin(0, 0, 0);
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typename GetGeomTraits<TriangleMesh, CGAL_PMP_NP_CLASS>::type::FT volume = 0.;
BOOST_FOREACH(face_descriptor f, faces(tmesh))
{
volume += CGAL::volume(origin,
get(vpm, target(halfedge(f, tmesh), tmesh)),
get(vpm, target(next(halfedge(f, tmesh), tmesh), tmesh)),
get(vpm, target(prev(halfedge(f, tmesh), tmesh), tmesh)));
exact(volume);
}
return volume;
volume += cv3(origin,
get(vpm, target(halfedge(f, tmesh), tmesh)),
get(vpm, target(next(halfedge(f, tmesh), tmesh), tmesh)),
get(vpm, target(prev(halfedge(f, tmesh), tmesh), tmesh)));
exact(volume);
}
return volume;
}
template<typename TriangleMesh>
typename CGAL::Kernel_traits<typename property_map_value<TriangleMesh,
CGAL::vertex_point_t>::type>::Kernel::FT
volume(const TriangleMesh& tmesh)
{
return volume(tmesh,
CGAL::Polygon_mesh_processing::parameters::all_default());
}

6
Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/self_intersections.h
View File

@ -241,7 +241,7 @@ self_intersections( const FaceRange& face_range,
* \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `pmesh`.
* If this parameter is omitted, an internal property map for
* `CGAL::vertex_point_t` should be available in `TriangleMesh`\cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `SelfIntersectionTraits` \cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `PMPSelfIntersectionTraits` \cgalParamEnd
* \cgalNamedParamsEnd
*
* @return `out`
@ -300,7 +300,7 @@ self_intersections(const TriangleMesh& tmesh, OutputIterator out)
* \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `pmesh`.
* If this parameter is omitted, an internal property map for
* `CGAL::vertex_point_t` should be available in `TriangleMesh`\cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `SelfIntersectionTraits` \cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `PMPSelfIntersectionTraits` \cgalParamEnd
* \cgalNamedParamsEnd
*/
@ -386,7 +386,7 @@ OutputIterator self_intersections(const FaceRange& face_range,
* \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `tmesh`.
* If this parameter is omitted, an internal property map for
* `CGAL::vertex_point_t` should be available in `TriangleMesh`\cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `SelfIntersectionTraits` \cgalParamEnd
* \cgalParamBegin{geom_traits} an instance of a geometric traits class, model of `PMPSelfIntersectionTraits` \cgalParamEnd
* \cgalNamedParamsEnd
*
* @return true if `tmesh` self-intersects

31
Polygon_mesh_processing/test/Polygon_mesh_processing/CMakeLists.txt
View File

@ -49,16 +49,26 @@ endif()
# include for local package
include_directories( BEFORE ../../include )
find_package(Eigen3 3.1.0) #(requires 3.1.0 or greater)
find_package(Eigen3 3.2.0) #(requires 3.2.0 or greater)
find_package( TBB )
if( TBB_FOUND )
include( ${TBB_USE_FILE} )
list( APPEND CGAL_3RD_PARTY_LIBRARIES ${TBB_LIBRARIES} )
else()
message( STATUS "NOTICE: Intel TBB was not found. test_pmp_distance will use sequential code." )
endif()
include( CGAL_CreateSingleSourceCGALProgram )
if (EIGEN3_FOUND)
# Executables that require Eigen 3.1
include( ${EIGEN3_USE_FILE} )
# Creating entries for all .cpp/.C files with "main" routine
# ##########################################################
include( CGAL_CreateSingleSourceCGALProgram )
create_single_source_cgal_program("fairing_test.cpp")
create_single_source_cgal_program("triangulate_hole_Polyhedron_3_no_delaunay_test.cpp" )
create_single_source_cgal_program("triangulate_hole_Polyhedron_3_test.cpp")
endif(EIGEN3_FOUND)
create_single_source_cgal_program("connected_component_polyhedron.cpp")
create_single_source_cgal_program("connected_component_surface_mesh.cpp")
@ -77,16 +87,5 @@ if (EIGEN3_FOUND)
create_single_source_cgal_program("measures_test.cpp")
create_single_source_cgal_program("triangulate_faces_test.cpp")
create_single_source_cgal_program("test_pmp_remove_border_edge.cpp")
if(NOT (${EIGEN3_VERSION} VERSION_LESS 3.2.0))
create_single_source_cgal_program("fairing_test.cpp")
create_single_source_cgal_program("triangulate_hole_Polyhedron_3_no_delaunay_test.cpp" )
create_single_source_cgal_program("triangulate_hole_Polyhedron_3_test.cpp")
create_single_source_cgal_program("test_pmp_distance.cpp")
create_single_source_cgal_program("triangulate_hole_polyline_test.cpp")
else()
message(STATUS "NOTICE: Some tests require Eigen 3.2 (or greater) and will not be compiled.")
endif()
else(EIGEN3_FOUND)
message(STATUS "NOTICE: Some examples require Eigen 3.1 (or greater) and will not be compiled.")
endif(EIGEN3_FOUND)

1
Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_distance.cmd Normal file
View File

@ -0,0 +1 @@
data/elephant.off data/blobby_3cc.off

309
Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_distance.cpp Normal file
View File

@ -0,0 +1,309 @@
#include <CGAL/Surface_mesh.h>
#include <CGAL/boost/graph/graph_traits_Surface_mesh.h>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Real_timer.h>
#include <CGAL/boost/graph/property_maps.h>
#include <fstream>
#include <ostream>
struct Custom_traits_Hausdorff
{
// New requirements {
struct FT
{
FT(){}
FT(double){}
FT operator/(FT)const{return FT();}
FT operator-(const FT)const{return FT();}
FT operator+(FT)const{return FT();}
FT operator*(FT)const{return FT();}
bool operator<=(FT)const{return false;}
bool operator>=(FT)const{return false;}
bool operator!=(FT)const{return false;}
bool operator<(FT)const{return false;}
bool operator>(FT)const{return false;}
FT& operator+=(FT){ return *this; }
};
struct Point_3
{
Point_3(){}
Point_3(FT, FT, FT){}
FT operator[](int)const{return FT();}
bool operator<(Point_3)const{return false;}
double x()const{return 0;}
double y()const{return 0;}
double z()const{return 0;}
};
struct Vector_3{};
struct Triangle_3
{
Triangle_3(const Point_3&, const Point_3&, const Point_3&){}
Point_3 operator[](int)const{return Point_3();}
CGAL::Bbox_3 bbox(){return CGAL::Bbox_3();}
};
struct Segment_3
{
Segment_3(const Point_3&, const Point_3&){}
Point_3 operator[](int)const{return Point_3();}
};
struct Compute_squared_area_3
{
typedef FT result_type;
FT operator()(Point_3,Point_3,Point_3)const{return FT();}
FT operator()(Triangle_3)const{return FT();}
};
struct Compute_squared_length_3
{
typedef FT result_type;
FT operator()(Segment_3)const{return FT();}
};
struct Construct_translated_point_3
{
Point_3 operator() (const Point_3 &, const Vector_3 &){return Point_3();}
};
struct Construct_vector_3{
Vector_3 operator() (const Point_3 &, const Point_3 &){return Vector_3();}
};
struct Construct_scaled_vector_3
{
Vector_3 operator() (const Vector_3 &, const FT &)
{return Vector_3();}
};
Compute_squared_area_3 compute_squared_area_3_object(){return Compute_squared_area_3();}
// } end of new requirements
// requirements from AABBGeomTraits {
struct Sphere_3
{};
struct Equal_3
{
bool operator()(Point_3,Point_3) const {return false;}
};
struct Do_intersect_3
{
CGAL::Comparison_result operator()(const Sphere_3& , const CGAL::Bbox_3&){ return CGAL::ZERO;}
};
struct Intersect_3
{
struct result_type{};
};
struct Construct_sphere_3
{
Sphere_3 operator()(const Point_3& , FT ){return Sphere_3();}
};
struct Construct_projected_point_3
{
const Point_3 operator()(Triangle_3, Point_3)const{return Point_3();}
};
struct Compare_distance_3
{
CGAL::Comparison_result operator()(Point_3, Point_3, Point_3)
{return CGAL::ZERO;}
};
struct Has_on_bounded_side_3
{
//documented as Comparision_result
CGAL::Comparison_result operator()(const Point_3&, const Point_3&, const Point_3&)
{return CGAL::ZERO;}
bool operator()(const Sphere_3&, const Point_3&)
{return false;}
};
struct Compute_squared_radius_3
{
FT operator()(const Sphere_3&){return FT();}
};
struct Compute_squared_distance_3
{
FT operator()(const Point_3& , const Point_3& ){return FT();}
};
Compare_distance_3 compare_distance_3_object(){return Compare_distance_3();}
Construct_sphere_3 construct_sphere_3_object(){return Construct_sphere_3();}
Construct_projected_point_3 construct_projected_point_3_object(){return Construct_projected_point_3();}
Compute_squared_distance_3 compute_squared_distance_3_object(){return Compute_squared_distance_3();}
Do_intersect_3 do_intersect_3_object(){return Do_intersect_3();}
Equal_3 equal_3_object(){return Equal_3();}
// } end of requirments from AABBGeomTraits
// requirements from SearchGeomTraits_3 {
struct Iso_cuboid_3{};
struct Construct_iso_cuboid_3{};
struct Construct_min_vertex_3{};
struct Construct_max_vertex_3{};
struct Construct_center_3{};
typedef const FT* Cartesian_const_iterator_3;
struct Construct_cartesian_const_iterator_3{
Construct_cartesian_const_iterator_3(){}
Construct_cartesian_const_iterator_3(const Point_3&){}
const FT* operator()(const Point_3&) const
{ return 0; }
const FT* operator()(const Point_3&, int) const
{ return 0; }
typedef const FT* result_type;
};
// } end of requirements from SearchGeomTraits_3
// requirements from SpatialSortingTraits_3 {
struct Less_x_3{
bool operator()(Point_3, Point_3){return false;}
};
struct Less_y_3{
bool operator()(Point_3, Point_3){return false;}
};
struct Less_z_3{
bool operator()(Point_3, Point_3){return false;}
};
Less_x_3 less_x_3_object()const{return Less_x_3();}
Less_y_3 less_y_3_object()const{return Less_y_3();}
Less_z_3 less_z_3_object()const{return Less_z_3();}
// } end of requirements from SpatialSortingTraits_3
};
namespace CGAL{
template<>struct Kernel_traits<Custom_traits_Hausdorff::Point_3>
{
typedef Custom_traits_Hausdorff Kernel;
};
template<>
struct Algebraic_structure_traits< Custom_traits_Hausdorff::FT >
: public Algebraic_structure_traits_base< Custom_traits_Hausdorff::FT, Field_tag >
{
};
template<>
struct Real_embeddable_traits< Custom_traits_Hausdorff::FT >
: public INTERN_RET::Real_embeddable_traits_base< Custom_traits_Hausdorff::FT , CGAL::Tag_true>
{
class To_double
: public std::unary_function< Custom_traits_Hausdorff::FT, double > {
public:
double operator()( const Custom_traits_Hausdorff::FT& ) const { return 0; }
};
};
}//end CGAL
void exact(Custom_traits_Hausdorff::FT){}
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
#include <CGAL/Polygon_mesh_processing/distance.h>
namespace PMP = CGAL::Polygon_mesh_processing;
typedef CGAL::Surface_mesh<K::Point_3> Mesh;
template <class GeomTraits, class TriangleMesh>
void general_tests(const TriangleMesh& m1,
const TriangleMesh& m2 )
{
if (m1.number_of_vertices()==0 || m2.number_of_vertices()==0) return;
std::vector<typename GeomTraits::Point_3> points;
points.push_back(typename GeomTraits::Point_3(0,0,0));
points.push_back(typename GeomTraits::Point_3(0,1,0));
points.push_back(typename GeomTraits::Point_3(1,0,0));
points.push_back(typename GeomTraits::Point_3(0,0,1));
points.push_back(typename GeomTraits::Point_3(0,2,0));
std::cout << "Symmetric distance between meshes (sequential) "
<< PMP::approximate_symmetric_Hausdorff_distance<CGAL::Sequential_tag>(
m1,m2,
PMP::parameters::number_of_points_per_area_unit(4000),
PMP::parameters::number_of_points_per_area_unit(4000))
<< "\n";
std::cout << "Max distance to point set "
<< PMP::approximate_max_distance_to_point_set(m1,points,4000)
<< "\n";
std::cout << "Max distance to triangle mesh (sequential) "
<< PMP::max_distance_to_triangle_mesh<CGAL::Sequential_tag>(points,m1)
<< "\n";
}
void test_concept()
{
typedef Custom_traits_Hausdorff CK;
CGAL::Surface_mesh<CK::Point_3> m1, m2;
general_tests<CK>(m1,m2);
}
int main(int, char** argv)
{
Mesh m1,m2;
std::ifstream input(argv[1]);
input >> m1;
input.close();
input.open(argv[2]);
input >> m2;
std::cout << "First mesh has " << num_faces(m1) << " faces\n";
std::cout << "Second mesh has " << num_faces(m2) << " faces\n";
CGAL::Real_timer time;
#ifdef CGAL_LINKED_WITH_TBB
time.start();
std::cout << "Distance between meshes (parallel) "
<< PMP::approximate_Hausdorff_distance<CGAL::Parallel_tag>(
m1,m2,PMP::parameters::number_of_points_per_area_unit(4000))
<< "\n";
time.stop();
std::cout << "done in " << time.time() << "s.\n";
#endif
time.reset();
time.start();
std::cout << "Distance between meshes (sequential) "
<< PMP::approximate_Hausdorff_distance<CGAL::Sequential_tag>(
m1,m2,PMP::parameters::number_of_points_per_area_unit(4000))
<< "\n";
time.stop();
std::cout << "done in " << time.time() << "s.\n";
general_tests<K>(m1,m2);
test_concept();
return 0;
}

1
Polyhedron/demo/Polyhedron/Plugins/PCA/Edit_box_plugin.cpp
View File

@ -111,6 +111,7 @@ void Edit_box_plugin::exportToPoly()
break;
}
}
Polyhedron::Point_3 points[8];
for(int i=0; i<8; ++i)
{

3
Polyhedron/demo/Polyhedron/Plugins/PMP/CMakeLists.txt
View File

@ -75,3 +75,6 @@ target_link_libraries(repair_polyhedron_plugin scene_polyhedron_item)
qt5_wrap_ui( isotropicRemeshingUI_FILES Isotropic_remeshing_dialog.ui)
polyhedron_demo_plugin(isotropic_remeshing_plugin Isotropic_remeshing_plugin ${isotropicRemeshingUI_FILES})
target_link_libraries(isotropic_remeshing_plugin scene_polyhedron_item scene_polyhedron_selection_item)
polyhedron_demo_plugin(distance_plugin Distance_plugin)
target_link_libraries(distance_plugin scene_polyhedron_item scene_color_ramp)

467
Polyhedron/demo/Polyhedron/Plugins/PMP/Distance_plugin.cpp Normal file
View File

@ -0,0 +1,467 @@
#include <CGAL/Three/Polyhedron_demo_plugin_interface.h>
#include <CGAL/Three/Scene_interface.h>
#include <QApplication>
#include <QObject>
#include <QAction>
#include <QMainWindow>
#include <QInputDialog>
#include <QMessageBox>
#include <QMap>
#include "Messages_interface.h"
#include "Scene_polyhedron_item.h"
#include "Color_ramp.h"
#include "triangulate_primitive.h"
#include <CGAL/Polygon_mesh_processing/distance.h>
#include <CGAL/Polygon_mesh_processing/compute_normal.h>
#include <boost/iterator/counting_iterator.hpp>
#include <CGAL/Search_traits_3.h>
#include <CGAL/Spatial_sort_traits_adapter_3.h>
#include <CGAL/property_map.h>
#include <boost/container/flat_map.hpp>
using namespace CGAL::Three;
namespace PMP = CGAL::Polygon_mesh_processing;
#ifdef CGAL_LINKED_WITH_TBB
template <class AABB_tree, class Point_3>
struct Distance_computation{
const AABB_tree& tree;
const std::vector<Point_3>& sample_points;
Point_3 initial_hint;
CGAL::cpp11::atomic<double>* distance;
std::vector<double>& output;
Distance_computation(const AABB_tree& tree,
const Point_3 p,
const std::vector<Point_3>& sample_points,
CGAL::cpp11::atomic<double>* d,
std::vector<double>& out )
: tree(tree)
, sample_points(sample_points)
, initial_hint(p)
, distance(d)
, output(out)
{
}
void
operator()(const tbb::blocked_range<std::size_t>& range) const
{
Point_3 hint = initial_hint;
double hdist = 0;
for( std::size_t i = range.begin(); i != range.end(); ++i)
{
hint = tree.closest_point(sample_points[i], hint);
Kernel::FT dist = squared_distance(hint,sample_points[i]);
double d = CGAL::sqrt(dist);
output[i] = d;
if (d>hdist) hdist=d;
}
if (hdist > distance->load(CGAL::cpp11::memory_order_acquire))
distance->store(hdist, CGAL::cpp11::memory_order_release);
}
};
#endif
class Scene_distance_polyhedron_item: public Scene_item
{
Q_OBJECT
public:
Scene_distance_polyhedron_item(Polyhedron* poly, Polyhedron* polyB, QString other_name, int sampling_pts)
:Scene_item(NbOfVbos,NbOfVaos),
poly(poly),
poly_B(polyB),
are_buffers_filled(false),
other_poly(other_name)
{
nb_pts_per_face = sampling_pts;
this->setRenderingMode(FlatPlusEdges);
thermal_ramp.build_thermal();
}
bool supportsRenderingMode(RenderingMode m) const {
return (m == Flat || m == FlatPlusEdges);
}
Scene_item* clone() const {return 0;}
QString toolTip() const {return QString("Item %1 with color indicating distance with %2").arg(this->name()).arg(other_poly);}
void draw(Viewer_interface *viewer) const
{
if(!are_buffers_filled)
{
computeElements();
initializeBuffers(viewer);
compute_bbox();
}
vaos[Facets]->bind();
attribBuffers(viewer, PROGRAM_WITH_LIGHT);
program = getShaderProgram(PROGRAM_WITH_LIGHT);
program->bind();
program->setUniformValue("is_selected", false);
viewer->glDrawArrays(GL_TRIANGLES, 0, static_cast<GLsizei>(nb_pos/3));
program->release();
vaos[Facets]->release();
}
void drawEdges(Viewer_interface* viewer) const
{
vaos[Edges]->bind();
attribBuffers(viewer, PROGRAM_WITHOUT_LIGHT);
program = getShaderProgram(PROGRAM_WITHOUT_LIGHT);
program->bind();
//draw the edges
program->setAttributeValue("colors", QColor(Qt::black));
program->setUniformValue("is_selected", false);
viewer->glDrawArrays(GL_LINES, 0, static_cast<GLsizei>(nb_edge_pos/3));
vaos[Edges]->release();
program->release();
}
void compute_bbox() const {
const Kernel::Point_3& p = *(poly->points_begin());
CGAL::Bbox_3 bbox(p.x(), p.y(), p.z(), p.x(), p.y(), p.z());
for(Polyhedron::Point_iterator it = poly->points_begin();
it != poly->points_end();
++it) {
bbox = bbox + it->bbox();
}
_bbox = Bbox(bbox.xmin(),bbox.ymin(),bbox.zmin(),
bbox.xmax(),bbox.ymax(),bbox.zmax());
}
private:
Polyhedron* poly;
Polyhedron* poly_B;
mutable bool are_buffers_filled;
QString other_poly;
mutable std::vector<float> vertices;
mutable std::vector<float> edge_vertices;
mutable std::vector<float> normals;
mutable std::vector<float> colors;
Color_ramp thermal_ramp;
int nb_pts_per_face;
enum VAOs {
Facets=0,
Edges,
NbOfVaos};
enum VBOs {
Vertices=0,
Edge_vertices,
Normals,
Colors,
NbOfVbos};
mutable std::size_t nb_pos;
mutable std::size_t nb_edge_pos;
mutable QOpenGLShaderProgram *program;
//fills 'out' and returns the hausdorff distance for calibration of the color_ramp.
double compute_distances(const Polyhedron& m, const std::vector<Kernel::Point_3>& sample_points,
std::vector<double>& out)const
{
typedef CGAL::AABB_face_graph_triangle_primitive<Polyhedron> Primitive;
typedef CGAL::AABB_traits<Kernel, Primitive> Traits;
typedef CGAL::AABB_tree< Traits > Tree;
Tree tree( faces(m).first, faces(m).second, m);
tree.accelerate_distance_queries();
tree.build();
typename Traits::Point_3 hint = m.vertices_begin()->point();
#ifndef CGAL_LINKED_WITH_TBB
double hdist = 0;
for(std::size_t i = 0; i<sample_points.size(); ++i)
{
hint = tree.closest_point(sample_points[i], hint);
typename Kernel::FT dist = squared_distance(hint,sample_points[i]);
double d = CGAL::sqrt(dist);
out[i]= d;
if (d>hdist) hdist=d;
}
return hdist;
#else
CGAL::cpp11::atomic<double> distance(0);
Distance_computation<Tree, typename Kernel::Point_3> f(tree, hint, sample_points, &distance, out);
tbb::parallel_for(tbb::blocked_range<std::size_t>(0, sample_points.size()), f);
return distance;
#endif
}
void computeElements()const
{
QApplication::setOverrideCursor(Qt::WaitCursor);
vertices.resize(0);
edge_vertices.resize(0);
normals.resize(0);
colors.resize(0);
typedef Polyhedron::Traits Kernel;
typedef Kernel::Vector_3 Vector;
typedef Polyhedron::Facet_iterator Facet_iterator;
typedef boost::graph_traits<Polyhedron>::face_descriptor face_descriptor;
typedef boost::graph_traits<Polyhedron>::vertex_descriptor vertex_descriptor;
//facets
{
boost::container::flat_map<face_descriptor, Vector> face_normals_map;
boost::associative_property_map< boost::container::flat_map<face_descriptor, Vector> >
nf_pmap(face_normals_map);
boost::container::flat_map<vertex_descriptor, Vector> vertex_normals_map;
boost::associative_property_map< boost::container::flat_map<vertex_descriptor, Vector> >
nv_pmap(vertex_normals_map);
PMP::compute_normals(*poly, nv_pmap, nf_pmap);
std::vector<Kernel::Point_3> total_points(0);
Facet_iterator f = poly->facets_begin();
for(f = poly->facets_begin();
f != poly->facets_end();
f++)
{
Vector nf = get(nf_pmap, f);
f->plane() = Kernel::Plane_3(f->halfedge()->vertex()->point(), nf);
typedef FacetTriangulator<Polyhedron, Polyhedron::Traits, boost::graph_traits<Polyhedron>::vertex_descriptor> FT;
double diagonal;
if(this->diagonalBbox() != std::numeric_limits<double>::infinity())
diagonal = this->diagonalBbox();
else
diagonal = 0.0;
//compute distance with other polyhedron
//sample facet
std::vector<Kernel::Point_3> sampled_points;
std::size_t nb_points = (std::max)((int)std::ceil(nb_pts_per_face * PMP::face_area(f,*poly,PMP::parameters::geom_traits(Kernel()))),
1);
CGAL::Random_points_in_triangle_3<typename Kernel::Point_3> g(f->halfedge()->vertex()->point(), f->halfedge()->next()->vertex()->point(),
f->halfedge()->next()->next()->vertex()->point());
CGAL::cpp11::copy_n(g, nb_points, std::back_inserter(sampled_points));
sampled_points.push_back(f->halfedge()->vertex()->point());
sampled_points.push_back(f->halfedge()->next()->vertex()->point());
sampled_points.push_back(f->halfedge()->next()->next()->vertex()->point());
//triangle facets with sample points for color display
FT triangulation(f,sampled_points,nf,poly,diagonal);
if(triangulation.cdt->dimension() != 2 )
{
qDebug()<<"Error : cdt not right (dimension != 2). Facet not displayed";
return;
}
//iterates on the internal faces to add the vertices to the positions
//and the normals to the appropriate vectors
for(FT::CDT::Finite_faces_iterator
ffit = triangulation.cdt->finite_faces_begin(),
end = triangulation.cdt->finite_faces_end();
ffit != end; ++ffit)
{
if(ffit->info().is_external)
continue;
for (int i = 0; i<3; ++i)
{
total_points.push_back(ffit->vertex(i)->point());
vertices.push_back(ffit->vertex(i)->point().x());
vertices.push_back(ffit->vertex(i)->point().y());
vertices.push_back(ffit->vertex(i)->point().z());
normals.push_back(nf.x());
normals.push_back(nf.y());
normals.push_back(nf.z());
}
}
}
//compute the distances
typedef CGAL::Spatial_sort_traits_adapter_3<Kernel,
CGAL::Pointer_property_map<Kernel::Point_3>::type > Search_traits_3;
std::vector<double> distances(total_points.size());
std::vector<std::size_t> indices;
indices.reserve(total_points.size());
std::copy(boost::counting_iterator<std::size_t>(0),
boost::counting_iterator<std::size_t>(total_points.size()),
std::back_inserter(indices));
spatial_sort(indices.begin(),
indices.end(),
Search_traits_3(CGAL::make_property_map(total_points)));
std::vector<Kernel::Point_3> sorted_points(total_points.size());
for(std::size_t i = 0; i < sorted_points.size(); ++i)
{
sorted_points[i] = total_points[indices[i]];
}
double hausdorff = compute_distances(*poly_B,
sorted_points,
distances);
//compute the colors
colors.resize(sorted_points.size()*3);
for(std::size_t i=0; i<sorted_points.size(); ++i)
{
std::size_t k = indices[i];
double d = distances[i]/hausdorff;
colors[3*k]=thermal_ramp.r(d);
colors[3*k+1]=thermal_ramp.g(d);
colors[3*k+2]=thermal_ramp.b(d);
}
}
//edges
{
//Lines
typedef Kernel::Point_3 Point;
typedef Polyhedron::Edge_iterator Edge_iterator;
Edge_iterator he;
for(he = poly->edges_begin();
he != poly->edges_end();
he++)
{
const Point& a = he->vertex()->point();
const Point& b = he->opposite()->vertex()->point();
{
edge_vertices.push_back(a.x());
edge_vertices.push_back(a.y());
edge_vertices.push_back(a.z());
edge_vertices.push_back(b.x());
edge_vertices.push_back(b.y());
edge_vertices.push_back(b.z());
}
}
}
QApplication::restoreOverrideCursor();
}
void initializeBuffers(Viewer_interface *viewer)const
{
program = getShaderProgram(PROGRAM_WITH_LIGHT, viewer);
program->bind();
vaos[Facets]->bind();
buffers[Vertices].bind();
buffers[Vertices].allocate(vertices.data(),
static_cast<GLsizei>(vertices.size()*sizeof(float)));
program->enableAttributeArray("vertex");
program->setAttributeBuffer("vertex",GL_FLOAT,0,3);
buffers[Vertices].release();
buffers[Normals].bind();
buffers[Normals].allocate(normals.data(),
static_cast<GLsizei>(normals.size()*sizeof(float)));
program->enableAttributeArray("normals");
program->setAttributeBuffer("normals",GL_FLOAT,0,3);
buffers[Normals].release();
buffers[Colors].bind();
buffers[Colors].allocate(colors.data(),
static_cast<GLsizei>(colors.size()*sizeof(float)));
program->enableAttributeArray("colors");
program->setAttributeBuffer("colors",GL_FLOAT,0,3);
buffers[Colors].release();
vaos[Facets]->release();
program->release();
program = getShaderProgram(PROGRAM_WITHOUT_LIGHT, viewer);
program->bind();
vaos[Edges]->bind();
buffers[Edge_vertices].bind();
buffers[Edge_vertices].allocate(edge_vertices.data(),
static_cast<GLsizei>(edge_vertices.size()*sizeof(float)));
program->enableAttributeArray("vertex");
program->setAttributeBuffer("vertex",GL_FLOAT,0,3);
buffers[Edge_vertices].release();
vaos[Facets]->release();
program->release();
nb_pos = vertices.size();
vertices.resize(0);
//"Swap trick" insures that the memory is indeed freed and not kept available
std::vector<float>(vertices).swap(vertices);
nb_edge_pos = edge_vertices.size();
edge_vertices.resize(0);
std::vector<float>(edge_vertices).swap(edge_vertices);
normals.resize(0);
std::vector<float>(normals).swap(normals);
colors.resize(0);
std::vector<float>(colors).swap(colors);
are_buffers_filled = true;
}
};
class DistancePlugin :
public QObject,
public Polyhedron_demo_plugin_interface
{
Q_OBJECT
Q_INTERFACES(CGAL::Three::Polyhedron_demo_plugin_interface)
Q_PLUGIN_METADATA(IID "com.geometryfactory.PolyhedronDemo.PluginInterface/1.0")
typedef Kernel::Point_3 Point_3;
public:
//decides if the plugin's actions will be displayed or not.
bool applicable(QAction*) const
{
return scene->selectionIndices().size() == 2 &&
qobject_cast<Scene_polyhedron_item*>(scene->item(scene->selectionIndices().first())) &&
qobject_cast<Scene_polyhedron_item*>(scene->item(scene->selectionIndices().last()));
}
//the list of the actions of the plugin.
QList<QAction*> actions() const
{
return _actions;
}
//this acts like a constructor for the plugin. It gets the references to the mainwindow and the scene, and connects the action.
void init(QMainWindow* mw, Scene_interface* sc, Messages_interface* mi)
{
//gets the reference to the message interface, to display text in the console widget
this->messageInterface = mi;
//get the references
this->scene = sc;
this->mw = mw;
//creates the action
QAction *actionComputeDistance= new QAction(QString("Compute Distance Between Polyhedra"), mw);
//specifies the subMenu
actionComputeDistance->setProperty("submenuName", "Polygon Mesh Processing");
//links the action
if(actionComputeDistance) {
connect(actionComputeDistance, SIGNAL(triggered()),
this, SLOT(createDistanceItems()));
_actions << actionComputeDistance;
}
}
public Q_SLOTS:
void createDistanceItems()
{
bool ok = false;
nb_pts_per_face = QInputDialog::getInt(mw, tr("Sampling"),
tr("Number of points per face:"),400, 1,2147483647,1, &ok);
if (!ok)
return;
//check the initial conditions
Scene_polyhedron_item* itemA = qobject_cast<Scene_polyhedron_item*>(scene->item(scene->selectionIndices().first()));
Scene_polyhedron_item* itemB = qobject_cast<Scene_polyhedron_item*>(scene->item(scene->selectionIndices().last()));
if(!itemA->polyhedron()->is_pure_triangle() ||
!itemB->polyhedron()->is_pure_triangle() ){
messageInterface->error(QString("Distance not computed. (Both polyhedra must be triangulated)"));
return;
}
QApplication::setOverrideCursor(Qt::WaitCursor);
Scene_distance_polyhedron_item* new_itemA = new Scene_distance_polyhedron_item(itemA->polyhedron(),itemB->polyhedron(), itemB->name(), nb_pts_per_face);
Scene_distance_polyhedron_item* new_itemB = new Scene_distance_polyhedron_item(itemB->polyhedron(),itemA->polyhedron(), itemA->name(), nb_pts_per_face);
itemA->setVisible(false);
itemB->setVisible(false);
new_itemA->setName(QString("%1 to %2").arg(itemA->name()).arg(itemB->name()));
new_itemB->setName(QString("%1 to %2").arg(itemB->name()).arg(itemA->name()));
scene->addItem(new_itemA);
scene->addItem(new_itemB);
QApplication::restoreOverrideCursor();
}
private:
int nb_pts_per_face;
QList<QAction*> _actions;
Messages_interface* messageInterface;
//The reference to the scene
Scene_interface* scene;
//The reference to the main window
QMainWindow* mw;
};
#include "Distance_plugin.moc"

2
Polyhedron/demo/Polyhedron/Plugins/Point_set/Surface_reconstruction_plugin_impl.cpp
View File

@ -237,8 +237,6 @@ Polyhedron* poisson_reconstruct(Point_set& points,
}
}
return output_mesh;
}

113
Polyhedron/demo/Polyhedron/triangulate_primitive.h
View File

@ -7,7 +7,7 @@
#include <CGAL/Three/Scene_item.h>
#include <CGAL/boost/graph/graph_traits_Polyhedron_3.h>
#include <CGAL/boost/graph/graph_traits_Surface_mesh.h>
#include <queue>
#include <QColor>
@ -67,6 +67,23 @@ public:
}
triangulate(idPoints, normal, item_diag);
}
FacetTriangulator(typename boost::graph_traits<Mesh>::face_descriptor fd,
const std::vector<typename Kernel::Point_3>& more_points,
const Vector& normal,
Mesh *poly,
const double item_diag)
{
std::vector<PointAndId> idPoints;
BOOST_FOREACH(halfedge_descriptor he_circ, halfedges_around_face( halfedge(fd, *poly), *poly))
{
PointAndId idPoint;
idPoint.point = get(boost::vertex_point,*poly,source(he_circ, *poly));
idPoint.id = source(he_circ, *poly);
idPoints.push_back(idPoint);
}
triangulate_with_points(idPoints,more_points, normal, item_diag);
}
FacetTriangulator(std::vector<PointAndId > &idPoints,
const Vector& normal,
@ -74,6 +91,13 @@ public:
{
triangulate(idPoints, normal, item_diag);
}
FacetTriangulator(std::vector<PointAndId > &idPoints,
const std::vector<typename Kernel::Point_3>& more_points,
const Vector& normal,
const double item_diag)
{
triangulate_with_points(idPoints, more_points, normal, item_diag);
}
private:
void triangulate( std::vector<PointAndId > &idPoints,
@ -89,11 +113,16 @@ private:
double min_sq_dist = CGAL::square(0.0001*item_diag);
// Iterate the points of the facet and decide if they must be inserted in the CDT
typename Kernel::FT x(0), y(0), z(0);
BOOST_FOREACH(PointAndId idPoint, idPoints)
{
x += idPoint.point.x();
y += idPoint.point.y();
z += idPoint.point.z();
typename CDT::Vertex_handle vh;
//Always insert the first point, then only insert
// if the distance with the precedent is reasonable.
// if the distance with the previous is reasonable.
if(first == 0 || CGAL::squared_distance(idPoint.point, previous->point()) > min_sq_dist)
{
vh = cdt->insert(idPoint.point);
@ -150,6 +179,86 @@ private:
}
}
}
void triangulate_with_points( std::vector<PointAndId > &idPoints,
const std::vector<typename Kernel::Point_3>& more_points,
const Vector& normal,
const double item_diag )
{
P_traits cdt_traits(normal);
cdt = new CDT(cdt_traits);
typename CDT::Vertex_handle previous, first, last_inserted;
//Compute a reasonable precision level used to decide
//if two consecutive points in a facet can be estimated
//equal.
double min_sq_dist = CGAL::square(0.0001*item_diag);
// Iterate the points of the facet and decide if they must be inserted in the CDT
BOOST_FOREACH(PointAndId idPoint, idPoints)
{
typename CDT::Vertex_handle vh;
//Always insert the first point, then only insert
// if the distance with the previous is reasonable.
if(first == 0 || CGAL::squared_distance(idPoint.point, previous->point()) > min_sq_dist)
{
vh = cdt->insert(idPoint.point);
v2v[vh] = idPoint.id;
if(first == 0) {
first = vh;
}
if(previous != 0 && previous != vh) {
double sq_dist = CGAL::squared_distance(previous->point(), vh->point());
if(sq_dist > min_sq_dist)
{
cdt->insert_constraint(previous, vh);
sq_dist = CGAL::squared_distance(previous->point(), first->point());
if(sq_dist > min_sq_dist)
{
last_inserted = previous;
}
}
}
previous = vh;
}
}
double sq_dist = CGAL::squared_distance(previous->point(), first->point());
if(sq_dist > min_sq_dist)
{
cdt->insert_constraint(previous, first);
}
else
{
cdt->insert_constraint(last_inserted, first);
}
BOOST_FOREACH(typename Kernel::Point_3 point, more_points)
{
cdt->insert(point);
}
// sets mark is_external
for(typename CDT::All_faces_iterator
fit2 = cdt->all_faces_begin(),
end = cdt->all_faces_end();
fit2 != end; ++fit2)
{
fit2->info().is_external = false;
}
//check if the facet is external or internal
std::queue<typename CDT::Face_handle> face_queue;
face_queue.push(cdt->infinite_vertex()->face());
while(! face_queue.empty() ) {
typename CDT::Face_handle fh = face_queue.front();
face_queue.pop();
if(fh->info().is_external) continue;
fh->info().is_external = true;
for(int i = 0; i <3; ++i) {
if(!cdt->is_constrained(std::make_pair(fh, i)))
{
face_queue.push(fh->neighbor(i));
}
}
}
}
};
#endif // TRIANGULATE_PRIMITIVE

4
Spatial_searching/include/CGAL/Orthogonal_k_neighbor_search.h
View File

@ -145,7 +145,7 @@ private:
FT diff2 = val - node->upper_low_value();
if ( (diff1 + diff2 >= FT(0.0)) )
{
new_off = 2*val < node->upper_low_value()+node->upper_high_value() ?
new_off = node->upper_low_value()+node->upper_high_value() > val*2?
val - node->upper_high_value():
val - node->upper_low_value();
bestChild = node->lower();
@ -153,7 +153,7 @@ private:
}
else // compute new distance
{
new_off = 2*val < node->lower_low_value()+node->lower_high_value() ?
new_off = node->lower_low_value()+node->lower_high_value() > val*2 ?
val - node->lower_high_value():
val - node->lower_low_value();
bestChild = node->upper();