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

Fix Mesh_3

This commit is contained in:
Andreas Fabri 2016-04-07 15:03:13 +02:00 committed by Jane Tournois
parent 468e861a96
commit 59cedfc7dd
16 changed files with 106 additions and 58 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
Cartesian_kernel/include/CGAL/Cartesian/function_objects.h
View File

@ -3008,6 +3008,11 @@ namespace CartesianKernelFunctors {
typedef const Point_3& type;
};
template<typename F>
struct result<F(Point_3)> {
typedef const Point_3& type;
};
Rep // Point_3
operator()(Return_base_tag, Origin o) const
@ -3021,6 +3026,10 @@ namespace CartesianKernelFunctors {
operator()(Return_base_tag, const RT& x, const RT& y, const RT& z, const RT& w) const
{ return Rep(x, y, z, w); }
const Point_3&
operator()(const Point_3 & p) const
{ return p; }
const Point_3&
operator()(const Weighted_point_3 & p) const
{ return p.rep().point(); }

9
Homogeneous_kernel/include/CGAL/Homogeneous/function_objects.h
View File

@ -3148,6 +3148,11 @@ namespace HomogeneousKernelFunctors {
typedef const Point_3& type;
};
template<typename F>
struct result<F(Point_3)> {
typedef const Point_3& type;
};
Rep // Point_3
operator()(Return_base_tag, Origin o) const
{ return Rep(o); }
@ -3165,6 +3170,10 @@ namespace HomogeneousKernelFunctors {
operator()(Return_base_tag, const RT& x, const RT& y, const RT& z, const RT& w) const
{ return Rep(x, y, z, w); }
const Point_3&
operator()(const Point_3 & p) const
{ return p; }
const Point_3&
operator()(const Weighted_point_3 & p) const
{ return p.rep().point(); }

2
Kernel_23/include/CGAL/Weighted_point_2.h
View File

@ -290,7 +290,7 @@ extract(std::istream& is, Weighted_point_2<R>& p, const Cartesian_tag&)
break;
}
if (is)
p = Weighted_point_2<R>(x, y,weight);
p = Weighted_point_2<R>(typename R::Point_2(x, y),weight);
return is;
}

25
Mesh_3/include/CGAL/Mesh_3/C3T3_helpers.h
View File

@ -384,6 +384,7 @@ class C3T3_helpers_base
protected:
typedef typename Tr::Geom_traits Gt;
typedef typename Gt::Point_3 Point_3;
typedef typename Gt::FT FT;
typedef typename Tr::Vertex_handle Vertex_handle;
typedef typename Tr::Cell_handle Cell_handle;
typedef typename Tr::Facet Facet;
@ -2647,6 +2648,7 @@ C3T3_helpers<C3T3,MD>::
rebuild_restricted_delaunay(OutdatedCells& outdated_cells,
Moving_vertices_set& moving_vertices)
{
typename Gt::Equal_3 equal = Gt().equal_3_object();
typename OutdatedCells::iterator first_cell = outdated_cells.begin();
typename OutdatedCells::iterator last_cell = outdated_cells.end();
Update_c3t3 updater(domain_,c3t3_);
@ -2745,7 +2747,7 @@ rebuild_restricted_delaunay(OutdatedCells& outdated_cells,
{
Point_3 new_pos = project_on_surface((*it)->point(),*it);
if ( new_pos != Point_3() )
if ( ! equal(new_pos, Point_3()) )
{
//freezing needs 'erase' to be done before the vertex is actually destroyed
// Update moving vertices (it becomes new_vertex)
@ -2772,6 +2774,7 @@ rebuild_restricted_delaunay(ForwardIterator first_cell,
ForwardIterator last_cell,
Moving_vertices_set& moving_vertices)
{
typename Gt::Equal_3 equal = Gt().equal_3_object();
Update_c3t3 updater(domain_,c3t3_);
// Get facets (returns each canonical facet only once)
@ -2851,7 +2854,7 @@ rebuild_restricted_delaunay(ForwardIterator first_cell,
{
Point_3 new_pos = project_on_surface((it->first)->point(),it->first,it->second);
if ( new_pos != Point_3() )
if ( ! equal(new_pos, Point_3()) )
{
//freezing needs 'erase' to be done before the vertex is actually destroyed
// Update moving vertices (it becomes new_vertex)
@ -3266,6 +3269,9 @@ project_on_surface_aux(const Point_3& p,
typename Gt::Is_degenerate_3 is_degenerate =
Gt().is_degenerate_3_object();
typename Gt::Construct_translated_point_3 translate =
Gt().construct_translated_point_3_object();
typename MD::Construct_intersection construct_intersection =
domain_.construct_intersection_object();
@ -3278,8 +3284,8 @@ project_on_surface_aux(const Point_3& p,
const Vector_3 projection_scaled_vector =
scale(projection_vector, CGAL::sqrt(sq_dist/sq_proj_length));
const Point_3 source = p + projection_scaled_vector;
const Point_3 target = p - projection_scaled_vector;
const Point_3 source = translate(p, projection_scaled_vector);
const Point_3 target = translate(p, - projection_scaled_vector);
const Segment_3 proj_segment(source,target);
@ -3354,10 +3360,14 @@ get_least_square_surface_plane(const Vertex_handle& v,
// Compute least square fitting plane
Plane_3 plane;
Point_3 point;
CGAL::linear_least_squares_fitting_3(surface_point_vector.begin(),
surface_point_vector.end(),
plane,
Dimension_tag<0>());
point,
Dimension_tag<0>(),
tr_.geom_traits(),
Default_diagonalize_traits<FT, 3>());
reference_point = surface_point_vector.front();
@ -3373,16 +3383,17 @@ project_on_surface(const Point_3& p,
const Vertex_handle& v,
Surface_patch_index index) const
{
typename Gt::Equal_3 equal = Gt().equal_3_object();
// return domain_.project_on_surface(p);
// Get plane
Point_3 reference_point(CGAL::ORIGIN);
Plane_3 plane = get_least_square_surface_plane(v,reference_point, index);
if ( reference_point == CGAL::ORIGIN )
if ( equal(reference_point, CGAL::ORIGIN) )
return p;
// Project
if ( p != v->point() )
if ( ! equal(p, v->point()) )
return project_on_surface_aux(p,
v->point(),
plane.orthogonal_vector());

3
Mesh_3/include/CGAL/Mesh_3/Lloyd_move.h
View File

@ -325,7 +325,8 @@ private:
// Fit plane using point-based PCA: compute least square fitting plane
Plane_3 plane;
CGAL::linear_least_squares_fitting_3(first,last,plane,Dimension_tag<0>());
Point_3 point;
CGAL::linear_least_squares_fitting_3(first,last,plane, point, Dimension_tag<0>(), Gt(),Default_diagonalize_traits<FT, 3>());
// Project all points to the plane
std::transform(first, last, first, Project_on_plane(plane));

2
Mesh_3/include/CGAL/Mesh_3/Mesher_3.h
View File

@ -533,7 +533,7 @@ initialize()
float(tbb::task_scheduler_init::default_num_threads())
* Concurrent_mesher_config::get().num_pseudo_infinite_vertices_per_core);
for (int i = 0 ; i < NUM_PSEUDO_INFINITE_VERTICES ; ++i, ++random_point)
r_c3t3_.add_far_point(*random_point + center);
r_c3t3_.add_far_point(r_c3t3_.triangulation().geom_traits().construct_translated_point_3_object()(*random_point, center));
# ifdef CGAL_CONCURRENT_MESH_3_VERBOSE
std::cerr << "done." << std::endl;

4
Mesh_3/include/CGAL/Mesh_3/Robust_weighted_circumcenter_filtered_traits_3.h
View File

@ -172,7 +172,7 @@ public:
// We use power_side_of_power_sphere_3: it is static filtered and
// we know that p,q,r,s are positive oriented
typename Rt::power_side_of_power_sphere_3 power_side_of_power_sphere = Rt().power_side_of_power_sphere_3_object();
typename Rt::Power_side_of_power_sphere_3 power_side_of_power_sphere = Rt().power_side_of_power_sphere_3_object();
// Compute denominator to swith to exact if it is 0
FT num_x, num_y, num_z, den;
@ -202,7 +202,7 @@ public:
return res;
} else {
// Fast output
if ( power_test(p,q,r,s,res) == CGAL::ON_POSITIVE_SIDE )
if ( power_side_of_power_sphere(p,q,r,s,res) == CGAL::ON_POSITIVE_SIDE )
return res;
}
}

4
Mesh_3/include/CGAL/Mesh_3/Sliver_perturber.h
View File

@ -1304,9 +1304,9 @@ perturb_vertex( PVertex pv
{
return;
}
Point_3 p = pv.vertex()->point();
if (!helper_.try_lock_point_no_spin(p) || p != pv.vertex()->point())
if (!helper_.try_lock_point_no_spin(p) || ! Gt().equal_3_object()(p,pv.vertex()->point()))
{
#ifdef CGAL_CONCURRENT_MESH_3_PROFILING
bcounter.increment_branch_2(); // THIS is an early withdrawal!

38
Mesh_3/include/CGAL/Mesh_3/Triangulation_helpers.h
View File

@ -339,7 +339,7 @@ inside_protecting_balls(const Tr& tr,
{
Vertex_handle nv = tr.nearest_power_vertex(p, v->cell());
if(nv->point().weight() > 0)
return CGAL::compare_squared_distance(p, nv->point(),
return Gt().compare_squared_distance_3_object()(p, nv->point(),
nv->point().weight()) != CGAL::LARGER;
return false;
}
@ -353,6 +353,8 @@ Triangulation_helpers<Tr>::
well_oriented(const Tr& tr,
const Cell_vector& cells_tos) const
{
typedef typename Tr::Geom_traits Gt;
typename Gt::Orientation_3 orientation = tr.geom_traits().orientation_3_object();
typename Cell_vector::const_iterator it = cells_tos.begin();
for( ; it != cells_tos.end() ; ++it)
{
@ -362,16 +364,16 @@ well_oriented(const Tr& tr,
int iv = c->index(tr.infinite_vertex());
Cell_handle cj = c->neighbor(iv);
int mj = tr.mirror_index(c, iv);
if(CGAL::orientation(cj->vertex(mj)->point(),
c->vertex((iv+1)&3)->point(),
c->vertex((iv+2)&3)->point(),
c->vertex((iv+3)&3)->point()) != CGAL::NEGATIVE)
if(orientation(cj->vertex(mj)->point(),
c->vertex((iv+1)&3)->point(),
c->vertex((iv+2)&3)->point(),
c->vertex((iv+3)&3)->point()) != CGAL::NEGATIVE)
return false;
}
else if(CGAL::orientation(c->vertex(0)->point(),
c->vertex(1)->point(),
c->vertex(2)->point(),
c->vertex(3)->point()) != CGAL::POSITIVE)
else if(orientation(c->vertex(0)->point(),
c->vertex(1)->point(),
c->vertex(2)->point(),
c->vertex(3)->point()) != CGAL::POSITIVE)
return false;
}
return true;
@ -387,6 +389,8 @@ well_oriented(const Tr& tr,
const Cell_vector& cells_tos,
const Point_getter& pg) const
{
typedef typename Tr::Geom_traits Gt;
typename Gt::Orientation_3 orientation = tr.geom_traits().orientation_3_object();
typename Cell_vector::const_iterator it = cells_tos.begin();
for( ; it != cells_tos.end() ; ++it)
{
@ -396,16 +400,16 @@ well_oriented(const Tr& tr,
int iv = c->index(tr.infinite_vertex());
Cell_handle cj = c->neighbor(iv);
int mj = tr.mirror_index(c, iv);
if(CGAL::orientation(pg(cj->vertex(mj)),
pg(c->vertex((iv+1)&3)),
pg(c->vertex((iv+2)&3)),
pg(c->vertex((iv+3)&3))) != CGAL::NEGATIVE)
if(orientation(pg(cj->vertex(mj)),
pg(c->vertex((iv+1)&3)),
pg(c->vertex((iv+2)&3)),
pg(c->vertex((iv+3)&3))) != CGAL::NEGATIVE)
return false;
}
else if(CGAL::orientation(pg(c->vertex(0)),
pg(c->vertex(1)),
pg(c->vertex(2)),
pg(c->vertex(3))) != CGAL::POSITIVE)
else if(orientation(pg(c->vertex(0)),
pg(c->vertex(1)),
pg(c->vertex(2)),
pg(c->vertex(3))) != CGAL::POSITIVE)
return false;
}
return true;

8
Mesh_3/include/CGAL/Mesh_3/min_dihedral_angle.h
View File

@ -85,10 +85,10 @@ minimum_dihedral_angle(
template <typename K>
typename K::FT
minimum_dihedral_angle(
const typename K::Point_3& p0,
const typename K::Point_3& p1,
const typename K::Point_3& p2,
const typename K::Point_3& p3,
const typename K::Bare_point& p0,
const typename K::Bare_point& p1,
const typename K::Bare_point& p2,
const typename K::Bare_point& p3,
K k = K())
{
typedef typename K::FT FT;

33
Mesh_3/include/CGAL/Mesh_3/vertex_perturbation.h
View File

@ -447,6 +447,9 @@ protected:
typename Gt::Compute_squared_length_3 sq_length =
Gt().compute_squared_length_3_object();
typename Gt::Construct_translated_point_3 translate =
Gt().construct_translated_point_3_object();
// create a helper
typedef C3T3_helpers<C3T3,MeshDomain> C3T3_helpers;
C3T3_helpers helper(c3t3, domain);
@ -456,7 +459,7 @@ protected:
// norm depends on the local size of the mesh
FT sq_norm = this->compute_perturbation_sq_amplitude(v, c3t3, sq_step_size_);
FT step_length = CGAL::sqrt(sq_norm/sq_length(gradient_vector));
Point_3 new_loc = v->point() + step_length * gradient_vector;
Point_3 new_loc = translate(v->point(), step_length * gradient_vector);
Point_3 final_loc = new_loc;
if ( c3t3.in_dimension(v) < 3 )
@ -471,7 +474,7 @@ protected:
v, final_loc)
&& (++i <= max_step_nb_) )
{
new_loc = new_loc + step_length * gradient_vector;
new_loc = translate(new_loc, step_length * gradient_vector);
if ( c3t3.in_dimension(v) == 3 )
final_loc = new_loc;
@ -484,7 +487,7 @@ protected:
while( Th().no_topological_change(c3t3.triangulation(), v, final_loc)
&& (++i <= max_step_nb_) )
{
new_loc = new_loc + step_length * gradient_vector;
new_loc = translate(new_loc, step_length * gradient_vector);
if ( c3t3.in_dimension(v) == 3 )
final_loc = new_loc;
@ -532,7 +535,7 @@ protected:
typedef typename C3T3::Triangulation::Geom_traits Gt;
typedef typename Gt::Vector_3 Vector_3;
typedef typename Gt::Point_3 Point_3;
typedef typename Gt::Construct_point_3 Construct_point_3;
public:
/**
* @brief Constructor
@ -624,12 +627,16 @@ private:
Vector_3 compute_gradient_vector(const Cell_handle& cell,
const Vertex_handle& v) const
{
typedef typename C3T3::Triangulation::Geom_traits Gt;
typename Gt::Construct_translated_point_3 translate =
Gt().construct_translated_point_3_object();
// translate the tet so that cell->vertex((i+3)&3) is 0_{R^3}
unsigned int index = cell->index(v);
Vector_3 translate_to_origin(CGAL::ORIGIN, cell->vertex((index+3)&3)->point()); //p4
const Point_3& p1 = v->point() - translate_to_origin;
const Point_3& p2 = cell->vertex((index+1)&3)->point() - translate_to_origin;
const Point_3& p3 = cell->vertex((index+2)&3)->point() - translate_to_origin;
const Point_3& p1 = translate(v->point(), - translate_to_origin);
const Point_3& p2 = translate(cell->vertex((index+1)&3)->point(), - translate_to_origin);
const Point_3& p3 = translate(cell->vertex((index+2)&3)->point(), - translate_to_origin);
// pre-compute everything
FT sq_p1 = p1.x()*p1.x() + p1.y()*p1.y() + p1.z()*p1.z();
@ -831,7 +838,8 @@ protected:
typedef typename C3T3::Triangulation::Geom_traits Gt;
typedef typename Gt::Vector_3 Vector_3;
typedef typename Gt::Point_3 Point_3;
typedef typename Gt::Construct_point_3 Construct_point_3;
public:
/**
* @brief constructor
@ -1012,7 +1020,8 @@ private:
const Point_3& p3 = ch->vertex(k3)->point();
// compute normal and return it
return normal(p1,p2,p3);
Construct_point_3 cp;
return normal(cp(p1),cp(p2),cp(p3));
}
};
@ -1225,7 +1234,11 @@ private:
bool *could_lock_zone = NULL) const
{
typedef Triangulation_helpers<typename C3T3::Triangulation> Th;
typedef typename C3T3::Triangulation::Geom_traits Gt;
typename Gt::Construct_translated_point_3 translate =
Gt().construct_translated_point_3_object();
modified_vertices.clear();
// Create an helper
@ -1249,7 +1262,7 @@ private:
Vector_3 delta = this->random_vector(sq_norm,on_sphere_);
// always from initial_location!
Point_3 new_location = initial_location + delta;
Point_3 new_location = translate(initial_location,delta);
if ( c3t3.in_dimension(moving_vertex) < 3 )
new_location = helper.project_on_surface(new_location, moving_vertex);

15
Mesh_3/test/Mesh_3/test_mesh_criteria_creation.cpp
View File

@ -99,9 +99,10 @@ int main()
// Test facet criteria
// -----------------------------------
typedef Tr::Geom_traits::FT FT;
FT radius_facet = CGAL::sqrt(CGAL::squared_radius(ch->vertex(k+1)->point(),
ch->vertex(k+2)->point(),
ch->vertex(k+3)->point()));
Tr::Geom_traits::Compute_squared_radius_3 squared_radius;
FT radius_facet = CGAL::sqrt(squared_radius(ch->vertex(k+1)->point(),
ch->vertex(k+2)->point(),
ch->vertex(k+3)->point()));
FT facet_size_ok = radius_facet*FT(10);
FT facet_size_nok = radius_facet/FT(10);
@ -154,10 +155,10 @@ int main()
// -----------------------------------
// Test cell criteria
// -----------------------------------
FT radius_cell = CGAL::sqrt(CGAL::squared_radius(ch->vertex(0)->point(),
ch->vertex(1)->point(),
ch->vertex(2)->point(),
ch->vertex(3)->point()));
FT radius_cell = CGAL::sqrt(squared_radius(ch->vertex(0)->point(),
ch->vertex(1)->point(),
ch->vertex(2)->point(),
ch->vertex(3)->point()));
FT cell_size_ok = radius_cell*FT(10);
FT cell_size_nok = radius_cell/FT(10);

4
Principal_component_analysis/include/CGAL/PCA_util.h
View File

@ -54,7 +54,7 @@ assemble_covariance_matrix_3(InputIterator first,
InputIterator beyond,
typename DiagonalizeTraits::Covariance_matrix& covariance, // covariance matrix
const typename K::Point_3& c, // centroid
const K& , // kernel
const K& k, // kernel
const typename K::Point_3*, // used for indirection
const CGAL::Dimension_tag<0>&,
const DiagonalizeTraits&)
@ -74,7 +74,7 @@ assemble_covariance_matrix_3(InputIterator first,
it++)
{
const Point& p = *it;
Vector d = p - c;
Vector d = k.construct_vector_3_object()(c,p);
covariance[0] += d.x() * d.x();
covariance[1] += d.x() * d.y();
covariance[2] += d.x() * d.z();

4
Principal_component_analysis_LGPL/include/CGAL/centroid.h
View File

@ -405,7 +405,7 @@ template < typename InputIterator,
typename K::Point_3
centroid(InputIterator begin,
InputIterator end,
const K&,
const K& k,
const typename K::Point_3*,
CGAL::Dimension_tag<0>)
{
@ -418,7 +418,7 @@ centroid(InputIterator begin,
unsigned int nb_pts = 0;
while (begin != end)
{
v = v + (*begin++ - ORIGIN);
v = v + k.construct_vector_3_object()(ORIGIN, *begin++);
nb_pts++;
}
return ORIGIN + v / (FT)nb_pts;

2
Voronoi_diagram_2/include/CGAL/Regular_triangulation_adaptation_traits_2.h
View File

@ -35,7 +35,7 @@ struct Regular_triangulation_adaptation_traits_2
: public CGAL_VORONOI_DIAGRAM_2_INS::Adaptation_traits_base_2
<RT2,
CGAL_VORONOI_DIAGRAM_2_INS::Point_accessor
<typename RT2::Geom_traits::Point_2,RT2,Tag_true>,
<typename RT2::Geom_traits::Point_2,RT2,Tag_false>,
CGAL_VORONOI_DIAGRAM_2_INS::Regular_triangulation_Voronoi_point_2<RT2>,
CGAL_VORONOI_DIAGRAM_2_INS::Regular_triangulation_nearest_site_2<RT2> >
{

2
Voronoi_diagram_2/include/CGAL/Voronoi_diagram_2/Regular_triangulation_degeneracy_testers.h
View File

@ -77,7 +77,7 @@ class Regular_triangulation_edge_tester_2
Site_2 s3 = v3->point();
Site_2 s4 = v4->point();
Oriented_side os =
dual.geom_traits().power_test_2_object()(s1,s2,s3,s4);
dual.geom_traits().power_side_of_power_circle_2_object()(s1,s2,s3,s4);
return os == ON_ORIENTED_BOUNDARY;
}