phynexis-deps/cgal
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Fix project_on_surface's projection direction

Browse Source 
See also https://github.com/CGAL/cgal/pull/5209
This commit is contained in:
Mael Rouxel-Labbé
2022-11-29 12:31:44 +01:00
parent c255b51b4b
commit cdc3bd22cf
3 changed files with 28 additions and 52 deletions
Download Patch File Download Diff File Expand all files Collapse all files
Mesh_3/include/CGAL/Mesh_3/C3T3_helpers.h
+9 -15
View File
@@ -3452,17 +3452,10 @@ get_least_square_surface_plane(const Vertex_handle& v,
(patch_index == Surface_patch_index() ||
c3t3_.surface_patch_index(f) == patch_index) )
{
ref_facet = f;
if(ref_facet.first == Cell_handle())
ref_facet = f;
// In the case of a periodic triangulation, the incident facets of a point
// do not necessarily have the same offsets. Worse, the surface centers
// might not have the same offset as their facet. Thus, no solution except
// calling a function 'get_closest_triangle(p, t)' that simply returns t
// for a non-periodic triangulation, and checks all possible offsets for
// periodic triangulations
Triangle t = c3t3_.triangulation().triangle(f);
Triangle ct = tr_.get_closest_triangle(cp(position), t);
const Triangle& ct = tr_.get_incident_triangle(f, v);
triangles.push_back(ct);
}
}
@@ -3474,7 +3467,6 @@ get_least_square_surface_plane(const Vertex_handle& v,
// Compute least square fitting plane
Plane_3 plane;
Bare_point point;
CGAL::linear_least_squares_fitting_3(triangles.begin(),
triangles.end(),
plane,
@@ -3483,12 +3475,13 @@ get_least_square_surface_plane(const Vertex_handle& v,
tr_.geom_traits(),
Default_diagonalize_traits<FT, 3>());
return std::make_pair(plane,
ref_facet.first->get_facet_surface_center(ref_facet.second));
// The surface center of a facet might have an offset in periodic triangulations
const Bare_point& ref_facet_scp = ref_facet.first->get_facet_surface_center(ref_facet.second);
const Bare_point& ref_point = tr_.get_closest_point(cp(position), ref_facet_scp);
return std::make_pair(plane, ref_point);
}
template <typename C3T3, typename MD>
typename C3T3_helpers<C3T3,MD>::Bare_point
C3T3_helpers<C3T3,MD>::
@@ -3510,6 +3503,7 @@ project_on_surface_if_possible(const Vertex_handle& v,
const Bare_point& p,
Surface_patch_index index) const
{
// @todo should call below if it's available...
// return domain_.project_on_surface(p);
typename Gt::Construct_point_3 cp = tr_.geom_traits().construct_point_3_object();
Mesh_3/include/CGAL/Mesh_triangulation_3.h
+4 -2
View File
@@ -62,12 +62,14 @@ public:
typedef typename Base::Triangle Triangle;
typedef typename Base::Vertex_handle Vertex_handle;
typedef typename Base::Facet Facet;
typedef typename Base::Cell_handle Cell_handle;
typedef typename Geom_traits::Vector_3 Vector;
using Base::geom_traits;
using Base::point;
using Base::triangle;
static std::string io_signature() { return Get_io_signature<Base>()(); }
@@ -83,9 +85,9 @@ public:
return q;
}
const Triangle& get_closest_triangle(const Bare_point& /*p*/, const Triangle& t) const
Triangle get_incident_triangle(const Facet& f, const Vertex_handle) const
{
return t;
return triangle(f);
}
void set_point(const Vertex_handle v,
Periodic_3_mesh_3/include/CGAL/Periodic_3_mesh_triangulation_3.h
+15 -35
View File
@@ -468,45 +468,25 @@ public:
return cwp(get_closest_point(cp(wp), cp(wq)), cw(wq));
}
Triangle get_closest_triangle(const Bare_point& p, const Triangle& t) const
// returns the triangle corresponding to f, with a geometric shift
// so that it is incident to ref_v's canonical position
Triangle get_incident_triangle(const Facet& f, const Vertex_handle ref_v) const
{
typename Geom_traits::Construct_vector_3 cv = geom_traits().construct_vector_3_object();
typename Geom_traits::Construct_point_3 cp = geom_traits().construct_point_3_object();
typename Geom_traits::Construct_translated_point_3 tr = geom_traits().construct_translated_point_3_object();
typename Geom_traits::Compute_squared_distance_3 csd = geom_traits().compute_squared_distance_3_object();
typename Geom_traits::Construct_vector_3 cv = geom_traits().construct_vector_3_object();
typename Geom_traits::Construct_triangle_3 ct = geom_traits().construct_triangle_3_object();
// It doesn't matter which point we use to canonicalize the triangle since we have to look
// at all the neighboring copies anyway since we do not have control of 'p'.
Bare_point canon_p0 = canonicalize_point(t[0]);
Vector_3 move_to_canonical = cv(t[0], canon_p0);
const std::array<Bare_point, 3> ct = { canon_p0,
tr(t[1], move_to_canonical),
tr(t[2], move_to_canonical) };
CGAL_precondition(f.first->has_vertex(ref_v));
const int ref_v_pos = f.first->index(ref_v);
const Bare_point& ref_p = cp(point(ref_v));
const Bare_point& ref_p_in_f = cp(point(f.first, ref_v_pos));
Vector_3 move_to_canonical = cv(ref_p_in_f, ref_p);
FT min_sq_dist = std::numeric_limits<FT>::infinity();
Triangle rt;
for(int i = -1; i < 2; ++i) {
for(int j = -1; j < 2; ++j) {
for(int k = -1; k < 2; ++k) {
const Triangle tt(
construct_point(std::make_pair(ct[0], Offset(i, j, k))),
construct_point(std::make_pair(ct[1], Offset(i, j, k))),
construct_point(std::make_pair(ct[2], Offset(i, j, k))));
const FT sq_dist = csd(p, tt);
if(sq_dist == FT(0))
return rt;
if(sq_dist < min_sq_dist) {
rt = tt;
min_sq_dist = sq_dist;
}
}
}
}
return rt;
const int s = f.second;
return ct(tr(cp(point(f.first, (s+1)%4)), move_to_canonical),
tr(cp(point(f.first, (s+2)%4)), move_to_canonical),
tr(cp(point(f.first, (s+3)%4)), move_to_canonical));
}
// Warning: This is a periodic version that computes the smallest possible