Split large and eclectic file PMP/repair.h into smaller files

This commit is contained in:
Mael Rouxel-Labbé 2019-07-25 09:25:54 +02:00
parent dac6c04529
commit 3a89e8240e
5 changed files with 3099 additions and 2890 deletions

View File

@ -0,0 +1,430 @@
// Copyright (c) 2015 GeometryFactory (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
// You can redistribute it and/or modify it under the terms of the GNU
// General Public License as published by the Free Software Foundation,
// either version 3 of the License, or (at your option) any later version.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
//
// Author(s) : Sebastien Loriot,
// Mael Rouxel-Labbé
#ifndef CGAL_POLYGON_MESH_PROCESSING_MANIFOLDNESS_H
#define CGAL_POLYGON_MESH_PROCESSING_MANIFOLDNESS_H
#include <CGAL/license/Polygon_mesh_processing/repair.h>
#include <CGAL/algorithm.h>
#include <CGAL/assertions.h>
#include <CGAL/Dynamic_property_map.h>
#include <CGAL/property_map.h>
#include <CGAL/Polygon_mesh_processing/internal/named_function_params.h>
#include <CGAL/Polygon_mesh_processing/internal/named_params_helper.h>
#include <fstream>
#include <iostream>
#include <map>
#include <utility>
#include <vector>
namespace CGAL {
namespace Polygon_mesh_processing {
/// \ingroup PMP_repairing_grp
/// checks whether a vertex of a polygon mesh is non-manifold.
///
/// @tparam PolygonMesh a model of `HalfedgeListGraph`
///
/// @param v a vertex of `pm`
/// @param pm a triangle mesh containing `v`
///
/// \warning This function has linear runtime with respect to the size of the mesh.
///
/// \sa `duplicate_non_manifold_vertices()`
///
/// \return `true` if the vertex is non-manifold, `false` otherwise.
template <typename PolygonMesh>
bool is_non_manifold_vertex(typename boost::graph_traits<PolygonMesh>::vertex_descriptor v,
const PolygonMesh& pm)
{
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
typedef CGAL::dynamic_halfedge_property_t<bool> Halfedge_property_tag;
typedef typename boost::property_map<PolygonMesh, Halfedge_property_tag>::const_type Visited_halfedge_map;
// Dynamic pmaps do not have default initialization values (yet)
Visited_halfedge_map visited_halfedges = get(Halfedge_property_tag(), pm);
for(halfedge_descriptor h : halfedges(pm))
put(visited_halfedges, h, false);
std::size_t incident_null_faces_counter = 0;
for(halfedge_descriptor h : halfedges_around_target(v, pm))
{
put(visited_halfedges, h, true);
if(CGAL::is_border(h, pm))
++incident_null_faces_counter;
}
if(incident_null_faces_counter > 1)
{
// The vertex is the sole connection between two connected components --> non-manifold
return true;
}
for(halfedge_descriptor h : halfedges(pm))
{
if(v == target(h, pm))
{
// Haven't seen that halfedge yet ==> more than one umbrella incident to 'v' ==> non-manifold
if(!get(visited_halfedges, h))
return true;
}
}
return false;
}
namespace internal {
template <typename G>
struct Vertex_collector
{
typedef typename boost::graph_traits<G>::vertex_descriptor vertex_descriptor;
bool has_old_vertex(const vertex_descriptor v) const { return collections.count(v) != 0; }
void collect_vertices(vertex_descriptor v1, vertex_descriptor v2)
{
std::vector<vertex_descriptor>& verts = collections[v1];
if(verts.empty())
verts.push_back(v1);
verts.push_back(v2);
}
template<typename OutputIterator>
void dump(OutputIterator out)
{
typedef std::pair<const vertex_descriptor, std::vector<vertex_descriptor> > Pair_type;
for(const Pair_type& p : collections)
*out++ = p.second;
}
void dump(Emptyset_iterator) { }
std::map<vertex_descriptor, std::vector<vertex_descriptor> > collections;
};
template <typename PolygonMesh, typename VPM, typename ConstraintMap>
typename boost::graph_traits<PolygonMesh>::vertex_descriptor
create_new_vertex_for_sector(typename boost::graph_traits<PolygonMesh>::halfedge_descriptor sector_begin_h,
typename boost::graph_traits<PolygonMesh>::halfedge_descriptor sector_last_h,
PolygonMesh& pm,
const VPM& vpm,
const ConstraintMap& cmap)
{
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
vertex_descriptor old_vd = target(sector_begin_h, pm);
vertex_descriptor new_vd = add_vertex(pm);
put(vpm, new_vd, get(vpm, old_vd));
put(cmap, new_vd, true);
set_halfedge(new_vd, sector_begin_h, pm);
halfedge_descriptor h = sector_begin_h;
do
{
set_target(h, new_vd, pm);
if(h == sector_last_h)
break;
else
h = prev(opposite(h, pm), pm);
}
while(h != sector_begin_h); // for safety
CGAL_assertion(h != sector_begin_h);
return new_vd;
}
template <typename PolygonMesh, typename NamedParameters>
std::size_t make_umbrella_manifold(typename boost::graph_traits<PolygonMesh>::halfedge_descriptor h,
PolygonMesh& pm,
internal::Vertex_collector<PolygonMesh>& dmap,
const NamedParameters& np)
{
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
using boost::get_param;
using boost::choose_param;
typedef typename GetVertexPointMap<PolygonMesh, NamedParameters>::type VertexPointMap;
VertexPointMap vpm = choose_param(get_param(np, internal_np::vertex_point),
get_property_map(vertex_point, pm));
typedef typename boost::lookup_named_param_def<internal_np::vertex_is_constrained_t,
NamedParameters,
Constant_property_map<vertex_descriptor, bool> // default (no constraint pmap)
>::type VerticesMap;
VerticesMap cmap = choose_param(get_param(np, internal_np::vertex_is_constrained),
Constant_property_map<vertex_descriptor, bool>(false));
std::size_t nb_new_vertices = 0;
vertex_descriptor old_v = target(h, pm);
put(cmap, old_v, true); // store the duplicates
// count the number of borders
int border_counter = 0;
halfedge_descriptor ih = h, done = ih, border_h = h;
do
{
if(is_border(ih, pm))
{
border_h = ih;
++border_counter;
}
ih = prev(opposite(ih, pm), pm);
}
while(ih != done);
bool is_non_manifold_within_umbrella = (border_counter > 1);
// if there is a single sector, then simply move the full umbrella to a new vertex, and we're done
if(!is_non_manifold_within_umbrella)
{
// note that since this is marked as a non-manifold vertex, we necessarily need to create
// a new vertex for this umbrella (the main umbrella is not marked as non-manifold)
halfedge_descriptor last_h = opposite(next(h, pm), pm);
vertex_descriptor new_v = create_new_vertex_for_sector(h, last_h, pm, vpm, cmap);
dmap.collect_vertices(old_v, new_v);
nb_new_vertices = 1;
}
// if there is more than one sector, look at each sector and split them away from the main one
else
{
// the first manifold sector, described by two halfedges
halfedge_descriptor sector_start_h = border_h;
CGAL_assertion(is_border(border_h, pm));
bool should_stop = false;
bool is_main_sector = true;
do
{
CGAL_assertion(is_border(sector_start_h, pm));
// collect the sector and split it away if it must be
halfedge_descriptor sector_last_h = sector_start_h;
do
{
halfedge_descriptor next_h = prev(opposite(sector_last_h, pm), pm);
if(is_border(next_h, pm))
break;
sector_last_h = next_h;
}
while(sector_last_h != sector_start_h);
CGAL_assertion(!is_border(sector_last_h, pm));
CGAL_assertion(sector_last_h != sector_start_h);
halfedge_descriptor next_start_h = prev(opposite(sector_last_h, pm), pm);
// there are multiple CCs incident to this particular vertex, and we should create a new vertex
// if it's not the first umbrella around 'old_v' or not the first sector, but not if it's
// the first umbrella and first sector.
bool must_create_new_vertex = (!is_main_sector || dmap.has_old_vertex(old_v));
// In any case, we must set up the next pointer correctly
set_next(sector_start_h, opposite(sector_last_h, pm), pm);
if(must_create_new_vertex)
{
vertex_descriptor new_v = create_new_vertex_for_sector(sector_start_h, sector_last_h, pm, vpm, cmap);
dmap.collect_vertices(old_v, new_v);
++nb_new_vertices;
}
is_main_sector = false;
sector_start_h = next_start_h;
should_stop = (sector_start_h == border_h);
}
while(!should_stop);
}
return nb_new_vertices;
}
} // end namespace internal
/// \ingroup PMP_repairing_grp
/// collects the non-manifold vertices (if any) present in the mesh. A non-manifold vertex `v` is returned
/// via one incident halfedge `h` such that `target(h, pm) = v` for all the umbrellas that `v` apppears in
/// (an <i>umbrella</i> being the set of faces incident to all the halfedges reachable by walking around `v`
/// using `hnext = prev(opposite(h, pm), pm)`, starting from `h`).
///
/// @tparam PolygonMesh a model of `HalfedgeListGraph`
/// @tparam OutputIterator a model of `OutputIterator` holding objects of type
/// `boost::graph_traits<PolygonMesh>::%halfedge_descriptor`
///
/// @param pm a triangle mesh
/// @param out the output iterator that collects halfedges incident to `v`
///
/// \sa `is_non_manifold_vertex()`
/// \sa `duplicate_non_manifold_vertices()`
///
/// \return the output iterator.
template <typename PolygonMesh, typename OutputIterator>
OutputIterator non_manifold_vertices(const PolygonMesh& pm,
OutputIterator out)
{
// Non-manifoldness can appear either:
// - if 'pm' is pinched at a vertex. While traversing the incoming halfedges at this vertex,
// we will meet strictly more than one border halfedge.
// - if there are multiple umbrellas around a vertex. In that case, we will find a non-visited
// halfedge that has for target a vertex that is already visited.
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
typedef CGAL::dynamic_vertex_property_t<bool> Vertex_property_tag;
typedef typename boost::property_map<PolygonMesh, Vertex_property_tag>::const_type Visited_vertex_map;
typedef CGAL::dynamic_halfedge_property_t<bool> Halfedge_property_tag;
typedef typename boost::property_map<PolygonMesh, Halfedge_property_tag>::const_type Visited_halfedge_map;
Visited_vertex_map visited_vertices = get(Vertex_property_tag(), pm);
Visited_halfedge_map visited_halfedges = get(Halfedge_property_tag(), pm);
// Dynamic pmaps do not have default initialization values (yet)
for(vertex_descriptor v : vertices(pm))
put(visited_vertices, v, false);
for(halfedge_descriptor h : halfedges(pm))
put(visited_halfedges, h, false);
for(halfedge_descriptor h : halfedges(pm))
{
// If 'h' is not visited yet, we walk around the target of 'h' and mark these
// halfedges as visited. Thus, if we are here and the target is already marked as visited,
// it means that the vertex is non manifold.
if(!get(visited_halfedges, h))
{
put(visited_halfedges, h, true);
bool is_non_manifold = false;
vertex_descriptor vd = target(h, pm);
if(get(visited_vertices, vd)) // already seen this vertex, but not from this star
is_non_manifold = true;
put(visited_vertices, vd, true);
// While walking the star of this halfedge, if we meet a border halfedge more than once,
// it means the mesh is pinched and we are also in the case of a non-manifold situation
halfedge_descriptor ih = h, done = ih;
int border_counter = 0;
do
{
put(visited_halfedges, ih, true);
if(is_border(ih, pm))
++border_counter;
ih = prev(opposite(ih, pm), pm);
}
while(ih != done);
if(border_counter > 1)
is_non_manifold = true;
if(is_non_manifold)
*out++ = h;
}
}
return out;
}
/// \ingroup PMP_repairing_grp
/// duplicates all the non-manifold vertices of the input mesh.
///
/// @tparam PolygonMesh a model of `HalfedgeListGraph` and `MutableHalfedgeGraph`
/// @tparam NamedParameters a sequence of \ref pmp_namedparameters "Named Parameters"
///
/// @param pm the surface mesh to be repaired
/// @param np optional \ref pmp_namedparameters "Named Parameters" described below
///
/// \cgalNamedParamsBegin
/// \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `pmesh`.
/// The type of this map is model of `ReadWritePropertyMap`.
/// If this parameter is omitted, an internal property map for
/// `CGAL::vertex_point_t` should be available in `PolygonMesh`
/// \cgalParamEnd
/// \cgalParamBegin{vertex_is_constrained_map} a writable property map with `vertex_descriptor`
/// as key and `bool` as `value_type`. `put(pmap, v, true)` will be called for each duplicated
/// vertices, as well as the original non-manifold vertex in the input mesh.
/// \cgalParamEnd
/// \cgalParamBegin{output_iterator} a model of `OutputIterator` with value type
/// `std::vector<vertex_descriptor>`. The first vertex of each vector is a non-manifold vertex
/// of the input mesh, followed by the new vertices that were created to fix this precise
/// non-manifold configuration.
/// \cgalParamEnd
/// \cgalNamedParamsEnd
///
/// \return the number of vertices created.
template <typename PolygonMesh, typename NamedParameters>
std::size_t duplicate_non_manifold_vertices(PolygonMesh& pm,
const NamedParameters& np)
{
using boost::get_param;
using boost::choose_param;
typedef boost::graph_traits<PolygonMesh> GT;
typedef typename GT::halfedge_descriptor halfedge_descriptor;
typedef typename boost::lookup_named_param_def <
internal_np::output_iterator_t,
NamedParameters,
Emptyset_iterator
> ::type Output_iterator;
Output_iterator out = choose_param(get_param(np, internal_np::output_iterator), Emptyset_iterator());
std::vector<halfedge_descriptor> non_manifold_cones;
non_manifold_vertices(pm, std::back_inserter(non_manifold_cones));
internal::Vertex_collector<PolygonMesh> dmap;
std::size_t nb_new_vertices = 0;
if(!non_manifold_cones.empty())
{
for(halfedge_descriptor h : non_manifold_cones)
nb_new_vertices += internal::make_umbrella_manifold(h, pm, dmap, np);
dmap.dump(out);
}
return nb_new_vertices;
}
template <class PolygonMesh>
std::size_t duplicate_non_manifold_vertices(PolygonMesh& pm)
{
return duplicate_non_manifold_vertices(pm, parameters::all_default());
}
} // namespace Polygon_mesh_processing
} // namespace CGAL
#endif // CGAL_POLYGON_MESH_PROCESSING_MANIFOLDNESS_H

View File

@ -0,0 +1,684 @@
// Copyright (c) 2015-2019 GeometryFactory (France).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
// You can redistribute it and/or modify it under the terms of the GNU
// General Public License as published by the Free Software Foundation,
// either version 3 of the License, or (at your option) any later version.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0+
//
// Author(s) : Sebastien Loriot,
// Mael Rouxel-Labbé
#ifndef CGAL_POLYGON_MESH_PROCESSING_REMOVE_SELF_INTERSECTIONS_H
#define CGAL_POLYGON_MESH_PROCESSING_REMOVE_SELF_INTERSECTIONS_H
#include <CGAL/license/Polygon_mesh_processing/repair.h>
#include <CGAL/Polygon_mesh_processing/manifoldness.h>
#include <CGAL/Polygon_mesh_processing/remove_degeneracies.h>
#include <CGAL/Polygon_mesh_processing/self_intersections.h>
#include <CGAL/Polygon_mesh_processing/triangulate_hole.h>
#include <CGAL/assertions.h>
#include <CGAL/boost/graph/named_function_params.h>
#include <CGAL/boost/graph/named_params_helper.h>
#include <CGAL/boost/graph/selection.h>
#include <CGAL/utility.h>
#include <array>
#include <iostream>
#include <iterator>
#include <map>
#include <set>
#include <utility>
#include <vector>
namespace CGAL {
namespace Polygon_mesh_processing {
/// \cond SKIP_IN_MANUAL
template <typename TriangleMesh, typename face_descriptor, typename VertexPointMap>
std::pair<bool, bool>
remove_self_intersections_one_step(TriangleMesh& tmesh,
std::set<face_descriptor>& faces_to_remove,
VertexPointMap& vpmap,
int step,
bool preserve_genus,
bool verbose)
{
std::set<face_descriptor> faces_to_remove_copy = faces_to_remove;
if(verbose)
{
std::cout << "DEBUG: running remove_self_intersections_one_step, step " << step
<< " with " << faces_to_remove.size() << " intersecting faces\n";
}
CGAL_assertion(tmesh.is_valid());
typedef boost::graph_traits<TriangleMesh> graph_traits;
typedef typename graph_traits::vertex_descriptor vertex_descriptor;
typedef typename graph_traits::halfedge_descriptor halfedge_descriptor;
typedef typename graph_traits::edge_descriptor edge_descriptor;
bool something_was_done = false; // indicates if a region was successfully remeshed
bool all_fixed = true; // indicates if all removal went well
// indicates if a removal was not possible because the region handle has
// some boundary cycle of halfedges
bool topology_issue = false;
if(verbose)
{
std::cout << " DEBUG: is_valid in one_step(tmesh)? ";
std::cout.flush();
std::cout << is_valid_polygon_mesh(tmesh) << "\n";
}
if(!faces_to_remove.empty())
{
while(!faces_to_remove.empty())
{
// Process a connected component of faces to remove.
// collect all the faces from the connected component
std::set<face_descriptor> cc_faces;
std::vector<face_descriptor> queue(1, *faces_to_remove.begin()); // temporary queue
cc_faces.insert(queue.back());
while(!queue.empty())
{
face_descriptor top = queue.back();
queue.pop_back();
halfedge_descriptor h = halfedge(top, tmesh);
for(int i=0; i<3; ++i)
{
face_descriptor adjacent_face = face(opposite(h, tmesh), tmesh);
if(adjacent_face!=boost::graph_traits<TriangleMesh>::null_face())
{
if(faces_to_remove.count(adjacent_face) != 0 && cc_faces.insert(adjacent_face).second)
queue.push_back(adjacent_face);
}
h = next(h, tmesh);
}
}
// expand the region to be filled
if(step > 0)
{
expand_face_selection(cc_faces, tmesh, step,
make_boolean_property_map(cc_faces),
Emptyset_iterator());
}
// try to compactify the selection region by also selecting all the faces included
// in the bounding box of the initial selection
std::vector<halfedge_descriptor> stack_for_expension;
Bbox_3 bb;
for(face_descriptor fd : cc_faces)
{
for(halfedge_descriptor h : halfedges_around_face(halfedge(fd, tmesh), tmesh))
{
bb += get(vpmap, target(h, tmesh)).bbox();
face_descriptor nf = face(opposite(h, tmesh), tmesh);
if(nf != boost::graph_traits<TriangleMesh>::null_face() && cc_faces.count(nf) == 0)
{
stack_for_expension.push_back(opposite(h, tmesh));
}
}
}
while(!stack_for_expension.empty())
{
halfedge_descriptor h = stack_for_expension.back();
stack_for_expension.pop_back();
if(cc_faces.count(face(h, tmesh)) == 1)
continue;
if(do_overlap(bb, get(vpmap, target(next(h, tmesh), tmesh)).bbox()))
{
cc_faces.insert(face(h, tmesh));
halfedge_descriptor candidate = opposite(next(h, tmesh), tmesh);
if(face(candidate, tmesh) != boost::graph_traits<TriangleMesh>::null_face())
stack_for_expension.push_back(candidate);
candidate = opposite(prev(h, tmesh), tmesh);
if(face(candidate, tmesh) != boost::graph_traits<TriangleMesh>::null_face())
stack_for_expension.push_back(candidate);
}
}
// remove faces from the set to process
for(face_descriptor f : cc_faces)
faces_to_remove.erase(f);
if(cc_faces.size() == 1)
continue; // it is a triangle nothing better can be done
//Check for non-manifold vertices in the selection and remove them by selecting all incident faces:
// extract the set of halfedges that is on the boundary of the holes to be
// made. In addition, we make sure no hole to be created contains a vertex
// visited more than once along a hole border (pinched surface)
// We save the size of boundary_hedges to make sur halfedges added
// from non_filled_hole are not removed.
bool non_manifold_vertex_remaining_in_selection = false;
do
{
bool non_manifold_vertex_removed = false; //here non-manifold is for the 1D polyline
std::vector<halfedge_descriptor> boundary_hedges;
for(face_descriptor fh : cc_faces)
{
halfedge_descriptor h = halfedge(fh, tmesh);
for(int i=0; i<3; ++i)
{
if(is_border(opposite(h, tmesh), tmesh) || cc_faces.count(face(opposite(h, tmesh), tmesh)) == 0)
boundary_hedges.push_back(h);
h = next(h, tmesh);
}
}
// detect vertices visited more than once along
// a hole border. We then remove all faces incident
// to such a vertex to force the removal of the vertex.
// Actually even if two holes are sharing a vertex, this
// vertex will be removed. It is not needed but since
// we do not yet have one halfedge per hole it is simpler
// and does not harm
std::set<vertex_descriptor> border_vertices;
for(halfedge_descriptor h : boundary_hedges)
{
if(!border_vertices.insert(target(h, tmesh)).second)
{
bool any_face_added = false;
for(halfedge_descriptor hh : halfedges_around_target(h, tmesh))
{
if(!is_border(hh, tmesh))
{
// add the face to the current selection
any_face_added |= cc_faces.insert(face(hh, tmesh)).second;
faces_to_remove.erase(face(hh, tmesh));
}
}
if(any_face_added)
non_manifold_vertex_removed = true;
else
non_manifold_vertex_remaining_in_selection = true;
}
}
if(!non_manifold_vertex_removed)
break;
}
while(true);
if(preserve_genus && non_manifold_vertex_remaining_in_selection)
{
topology_issue = true;
if(verbose)
std::cout << " DEBUG: CC not handled due to the presence at least one non-manifold vertex\n";
continue; // cannot replace a patch containing a nm vertex by a disk
}
// before running this function if preserve_genus=false, we duplicated
// all of them
CGAL_assertion(!non_manifold_vertex_remaining_in_selection);
// Collect halfedges on the boundary of the region to be selected
// (pointing inside the domain to be remeshed)
std::vector<halfedge_descriptor> cc_border_hedges;
for(face_descriptor fd : cc_faces)
{
halfedge_descriptor h = halfedge(fd, tmesh);
for(int i=0; i<3;++i)
{
if(is_border(opposite(h, tmesh), tmesh) || cc_faces.count(face(opposite(h, tmesh), tmesh))== 0)
cc_border_hedges.push_back(h);
h = next(h, tmesh);
}
}
if(!is_selection_a_topological_disk(cc_faces, tmesh))
{
// check if the selection contains cycles of border halfedges
bool only_border_edges = true;
std::set<halfedge_descriptor> mesh_border_hedge;
for(halfedge_descriptor h : cc_border_hedges)
{
if(!is_border(opposite(h, tmesh), tmesh))
only_border_edges = false;
else
mesh_border_hedge.insert(opposite(h, tmesh));
}
int nb_cycles=0;
while(!mesh_border_hedge.empty())
{
// we must count the number of cycle of boundary edges
halfedge_descriptor h_b = *mesh_border_hedge.begin(), h=h_b;
mesh_border_hedge.erase(mesh_border_hedge.begin());
do
{
h = next(h, tmesh);
if(h == h_b)
{
// found a cycle
++nb_cycles;
break;
}
else
{
typename std::set<halfedge_descriptor>::iterator it = mesh_border_hedge.find(h);
if(it == mesh_border_hedge.end())
break; // not a cycle
mesh_border_hedge.erase(it);
}
}
while(true);
}
if(nb_cycles > (only_border_edges ? 1 : 0))
{
if(verbose)
std::cout << " DEBUG: CC not handled due to the presence of "
<< nb_cycles << " of boundary edges\n";
topology_issue = true;
continue;
}
else
{
if(preserve_genus)
{
if(verbose)
std::cout << " DEBUG: CC not handled because it is not a topological disk (preserve_genus=true)\n";
all_fixed = false;
continue;
}
// count the number of cycles of halfedges of the boundary
std::map<vertex_descriptor, vertex_descriptor> bhs;
for(halfedge_descriptor h : cc_border_hedges)
bhs[source(h, tmesh)] = target(h, tmesh);
int nbc=0;
while(!bhs.empty())
{
++nbc;
std::pair<vertex_descriptor, vertex_descriptor > top=*bhs.begin();
bhs.erase(bhs.begin());
do
{
typename std::map<vertex_descriptor, vertex_descriptor>::iterator
it_find = bhs.find(top.second);
if(it_find == bhs.end()) break;
top = *it_find;
bhs.erase(it_find);
}
while(true);
}
if(nbc!=1)
{
if(verbose)
std::cout << " DEBUG: CC not handled because it is not a topological disk("
<< nbc << " boundary cycles)\n";
all_fixed = false;
continue;
}
else
{
if(verbose)
std::cout << " DEBUG: CC that is not a topological disk but has only one boundary cycle(preserve_genus=false)\n";
}
}
}
// sort halfedges so that they describe the sequence
// of halfedges of the hole to be made
CGAL_assertion(cc_border_hedges.size() > 2);
for(std::size_t i=0; i < cc_border_hedges.size()-2; ++i)
{
vertex_descriptor tgt = target(cc_border_hedges[i], tmesh);
for(std::size_t j=i+1; j<cc_border_hedges.size(); ++j)
{
if(tgt == source(cc_border_hedges[j], tmesh))
{
std::swap(cc_border_hedges[i+1], cc_border_hedges[j]);
break;
}
CGAL_assertion(j!=cc_border_hedges.size()-1);
}
}
CGAL_assertion(source(cc_border_hedges.front(), tmesh) == target(cc_border_hedges.back(), tmesh));
// collect vertices and edges inside the current selection cc
std::set<vertex_descriptor> cc_interior_vertices;
std::set<edge_descriptor> cc_interior_edges;
// first collect all vertices and edges incident to the faces to remove
for(face_descriptor fh : cc_faces)
{
for(halfedge_descriptor h : halfedges_around_face(halfedge(fh, tmesh), tmesh))
{
if(halfedge(target(h, tmesh), tmesh) == h) // limit the number of insertions
cc_interior_vertices.insert(target(h, tmesh));
cc_interior_edges.insert(edge(h, tmesh));
}
}
// and then remove those on the boundary
for(halfedge_descriptor h : cc_border_hedges)
{
cc_interior_vertices.erase(target(h, tmesh));
cc_interior_edges.erase(edge(h, tmesh));
}
if(verbose)
{
std::cout << " DEBUG: is_valid(tmesh) in one_step, before mesh changes? ";
std::cout << is_valid_polygon_mesh(tmesh) << std::endl;
}
//try hole_filling.
typedef CGAL::Triple<int, int, int> Face_indices;
typedef typename boost::property_traits<VertexPointMap>::value_type Point;
std::vector<Point> hole_points, third_points;
hole_points.reserve(cc_border_hedges.size());
third_points.reserve(cc_border_hedges.size());
std::vector<vertex_descriptor> border_vertices;
for(halfedge_descriptor h : cc_border_hedges)
{
vertex_descriptor v = source(h, tmesh);
hole_points.push_back(get(vpmap, v));
border_vertices.push_back(v);
third_points.push_back(get(vpmap, target(next(opposite(h, tmesh), tmesh), tmesh))); // TODO fix me for mesh border edges
}
CGAL_assertion(hole_points.size() >= 3);
// try to triangulate the hole using default parameters
//(using Delaunay search space if CGAL_HOLE_FILLING_DO_NOT_USE_DT3 is not defined)
std::vector<Face_indices> patch;
if(hole_points.size()>3)
triangulate_hole_polyline(hole_points, third_points, std::back_inserter(patch));
else
patch.push_back(Face_indices(0,1,2)); // trivial hole filling
if(patch.empty())
{
#ifndef CGAL_HOLE_FILLING_DO_NOT_USE_DT3
if(verbose)
std::cout << " DEBUG: Failed to fill a hole using Delaunay search space.\n";
triangulate_hole_polyline(hole_points, third_points, std::back_inserter(patch),
parameters::use_delaunay_triangulation(false));
#endif
if(patch.empty())
{
if(verbose)
std::cout << " DEBUG: Failed to fill a hole using the whole search space.\n";
all_fixed = false;
continue;
}
}
// make sure that the hole filling is valid, we check that no
// edge already in the mesh is present in patch.
bool non_manifold_edge_found = false;
for(const Face_indices& triangle : patch)
{
std::array<int, 6> edges = make_array(triangle.first, triangle.second,
triangle.second, triangle.third,
triangle.third, triangle.first);
for(int k=0; k<3; ++k)
{
int vi=edges[2*k], vj=edges[2*k+1];
// ignore boundary edges
if(vi+1==vj || (vj==0 && static_cast<std::size_t>(vi)==border_vertices.size()-1))
continue;
halfedge_descriptor h = halfedge(border_vertices[vi], border_vertices[vj], tmesh).first;
if(h != boost::graph_traits<TriangleMesh>::null_halfedge() &&
cc_interior_edges.count(edge(h, tmesh)) == 0)
{
non_manifold_edge_found=true;
break;
}
}
if(non_manifold_edge_found)
break;
}
if(non_manifold_edge_found)
{
if(verbose)
std::cout << " DEBUG: Triangulation produced is non-manifold when plugged into the mesh.\n";
all_fixed = false;
continue;
}
// plug the new triangles in the mesh, reusing previous edges and faces
std::vector<edge_descriptor> edge_stack(cc_interior_edges.begin(), cc_interior_edges.end());
std::vector<face_descriptor> face_stack(cc_faces.begin(), cc_faces.end());
std::map<std::pair<int, int>, halfedge_descriptor> halfedge_map;
int i = 0;
// register border halfedges
for(halfedge_descriptor h : cc_border_hedges)
{
int j = static_cast<int>(std::size_t(i+1)%cc_border_hedges.size());
halfedge_map.insert(std::make_pair(std::make_pair(i, j), h));
set_halfedge(target(h, tmesh), h, tmesh); // update vertex halfedge pointer
CGAL_assertion(border_vertices[i] == source(h, tmesh) && border_vertices[j] == target(h, tmesh));
++i;
}
std::vector<halfedge_descriptor> hedges;
hedges.reserve(4);
face_descriptor f = boost::graph_traits<TriangleMesh>::null_face();
for(const Face_indices& triangle : patch)
{
// get the new face
if(face_stack.empty())
{
f = add_face(tmesh);
}
else
{
f = face_stack.back();
face_stack.pop_back();
}
std::array<int, 4> indices = make_array(triangle.first, triangle.second, triangle.third, triangle.first);
for(int i=0; i<3; ++i)
{
// get the corresponding halfedge (either a new one or an already created)
typename std::map<std::pair<int, int>, halfedge_descriptor >::iterator insert_res =
halfedge_map.insert(std::make_pair(std::make_pair(indices[i], indices[i+1]),
boost::graph_traits<TriangleMesh>::null_halfedge())).first;
if(insert_res->second == boost::graph_traits<TriangleMesh>::null_halfedge())
{
if(edge_stack.empty())
{
insert_res->second = halfedge(add_edge(tmesh), tmesh);
}
else
{
insert_res->second = halfedge(edge_stack.back(), tmesh);
edge_stack.pop_back();
}
halfedge_map[std::make_pair(indices[i+1], indices[i])] = opposite(insert_res->second, tmesh);
}
hedges.push_back(insert_res->second);
}
hedges.push_back(hedges.front());
// update halfedge connections + face pointers
for(int i=0; i<3;++i)
{
set_next(hedges[i], hedges[i+1], tmesh);
set_face(hedges[i], f, tmesh);
set_target(hedges[i], border_vertices[indices[i+1]], tmesh);
}
set_halfedge(f, hedges[0], tmesh);
hedges.clear();
}
// now remove remaining edges,
for(edge_descriptor e : edge_stack)
remove_edge(e, tmesh);
// vertices,
for(vertex_descriptor vh : cc_interior_vertices)
remove_vertex(vh, tmesh);
// and remaning faces
for(face_descriptor f : face_stack)
remove_face(f, tmesh);
if(verbose)
std::cout << " DEBUG: " << cc_faces.size() << " triangles removed, " << patch.size() << " created\n";
CGAL_assertion(is_valid_polygon_mesh(tmesh));
something_was_done = true;
}
}
if(!something_was_done)
{
faces_to_remove.swap(faces_to_remove_copy);
if(verbose)
std::cout << " DEBUG: Nothing was changed during this step, self-intersections won't be recomputed." << std::endl;
}
return std::make_pair(all_fixed, topology_issue);
}
template <typename FaceRange, typename TriangleMesh, typename NamedParameters>
bool remove_self_intersections(const FaceRange& face_range,
TriangleMesh& tmesh,
const NamedParameters& np)
{
typedef boost::graph_traits<TriangleMesh> graph_traits;
typedef typename graph_traits::face_descriptor face_descriptor;
// named parameter extraction
typedef typename GetVertexPointMap<TriangleMesh, NamedParameters>::type VertexPointMap;
VertexPointMap vpm = boost::choose_param(boost::get_param(np, internal_np::vertex_point),
get_property_map(vertex_point, tmesh));
const int max_steps = boost::choose_param(boost::get_param(np, internal_np::number_of_iterations), 7);
bool verbose = boost::choose_param(boost::get_param(np, internal_np::verbosity_level), 0) > 0;
bool preserve_genus = boost::choose_param(boost::get_param(np, internal_np::preserve_genus), true);
if(verbose)
std::cout << "DEBUG: Starting remove_self_intersections, is_valid(tmesh)? " << is_valid_polygon_mesh(tmesh) << "\n";
CGAL_precondition_code(std::set<face_descriptor> degenerate_face_set;)
CGAL_precondition_code(degenerate_faces(face_range, tmesh,
std::inserter(degenerate_face_set, degenerate_face_set.begin()), np);)
CGAL_precondition(degenerate_face_set.empty());
if(!preserve_genus)
duplicate_non_manifold_vertices(tmesh, np);
// Look for self-intersections in the polyhedron and remove them
int step = -1;
bool all_fixed = true; // indicates if the filling of all created holes went fine
bool topology_issue = false; // indicates if some boundary cycles of edges are blocking the fixing
std::set<face_descriptor> faces_to_remove;
while(++step<max_steps)
{
if(faces_to_remove.empty()) // the previous round might have been blocked due to topological constraints
{
typedef std::pair<face_descriptor, face_descriptor> Face_pair;
std::vector<Face_pair> self_inter;
// TODO : possible optimization to reduce the range to check with the bbox
// of the previous patches or something.
self_intersections(face_range, tmesh, std::back_inserter(self_inter));
for(const Face_pair& fp : self_inter)
{
faces_to_remove.insert(fp.first);
faces_to_remove.insert(fp.second);
}
}
if(faces_to_remove.empty() && all_fixed)
{
if(verbose)
std::cout << "DEBUG: There is no more face to remove." << std::endl;
break;
}
std::tie(all_fixed, topology_issue) =
remove_self_intersections_one_step(tmesh, faces_to_remove, vpm, step, preserve_genus, verbose);
if(all_fixed && topology_issue)
{
if(verbose)
std::cout<< "DEBUG: Process stopped because of boundary cycles"
" of boundary edges involved in self-intersections.\n";
return false;
}
}
return step < max_steps;
}
template <typename FaceRange, typename TriangleMesh>
bool remove_self_intersections(const FaceRange& face_range,
TriangleMesh& tmesh)
{
return remove_self_intersections(face_range, tmesh, parameters::all_default());
}
template <typename TriangleMesh, typename CGAL_PMP_NP_TEMPLATE_PARAMETERS>
bool remove_self_intersections(TriangleMesh& tmesh,
const CGAL_PMP_NP_CLASS& np)
{
return remove_self_intersections(faces(tmesh), tmesh, np);
}
template <typename TriangleMesh>
bool remove_self_intersections(TriangleMesh& tmesh)
{
return remove_self_intersections(faces(tmesh), tmesh, parameters::all_default());
}
/// \endcond
} // namespace Polygon_mesh_processing
} // namespace CGAL
#endif // CGAL_POLYGON_MESH_PROCESSING_REMOVE_SELF_INTERSECTIONS_H

View File

@ -32,6 +32,7 @@
#include <CGAL/boost/graph/iterator.h>
#include <CGAL/boost/graph/helpers.h>
#include <boost/range/has_range_iterator.hpp>
#include <boost/graph/graph_traits.hpp>
#include <limits>
@ -94,6 +95,81 @@ bool is_degenerate_edge(typename boost::graph_traits<PolygonMesh>::edge_descript
return is_degenerate_edge(e, pm, parameters::all_default());
}
/// \ingroup PMP_repairing_grp
/// collects the degenerate edges within a given range of edges.
///
/// @tparam EdgeRange a model of `Range` with value type `boost::graph_traits<TriangleMesh>::%edge_descriptor`
/// @tparam TriangleMesh a model of `HalfedgeGraph`
/// @tparam NamedParameters a sequence of \ref pmp_namedparameters "Named Parameters"
///
/// @param edges a subset of edges of `tm`
/// @param tm a triangle mesh
/// @param out an output iterator in which the degenerate edges are written
/// @param np optional \ref pmp_namedparameters "Named Parameters" described below
///
/// \cgalNamedParamsBegin
/// \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `tm`.
/// The type of this map is model of `ReadWritePropertyMap`.
/// If this parameter is omitted, an internal property map for
/// `CGAL::vertex_point_t` should be available in `TriangleMesh`
/// \cgalParamEnd
/// \cgalParamBegin{geom_traits} a geometric traits class instance.
/// The traits class must provide the nested type `Point_3`,
/// and the nested functor `Equal_3` to check whether two points are identical.
/// \cgalParamEnd
/// \cgalNamedParamsEnd
template <typename EdgeRange, typename TriangleMesh, typename OutputIterator, typename NamedParameters>
OutputIterator degenerate_edges(const EdgeRange& edges,
const TriangleMesh& tm,
OutputIterator out,
const NamedParameters& np)
{
typedef typename boost::graph_traits<TriangleMesh>::edge_descriptor edge_descriptor;
for(edge_descriptor ed : edges)
if(is_degenerate_edge(ed, tm, np))
*out++ = ed;
return out;
}
template <typename EdgeRange, typename TriangleMesh, typename OutputIterator>
OutputIterator degenerate_edges(const EdgeRange& edges,
const TriangleMesh& tm,
OutputIterator out,
typename boost::enable_if<
typename boost::has_range_iterator<EdgeRange>
>::type* = 0)
{
return degenerate_edges(edges, tm, out, CGAL::parameters::all_default());
}
/// \ingroup PMP_repairing_grp
/// calls the function `degenerate_edges()` with the range: `edges(tm)`.
///
/// See above for the comprehensive description of the parameters.
///
template <typename TriangleMesh, typename OutputIterator, typename NamedParameters>
OutputIterator degenerate_edges(const TriangleMesh& tm,
OutputIterator out,
const NamedParameters& np
#ifndef DOXYGEN_RUNNING
,
typename boost::disable_if<
boost::has_range_iterator<TriangleMesh>
>::type* = 0
#endif
)
{
return degenerate_edges(edges(tm), tm, out, np);
}
template <typename TriangleMesh, typename OutputIterator>
OutputIterator degenerate_edges(const TriangleMesh& tm, OutputIterator out)
{
return degenerate_edges(edges(tm), tm, out, CGAL::parameters::all_default());
}
/// \ingroup PMP_repairing_grp
/// checks whether a triangle face is degenerate.
/// A triangle face is considered degenerate if the geometric positions of its vertices are collinear.
@ -151,6 +227,83 @@ bool is_degenerate_triangle_face(typename boost::graph_traits<TriangleMesh>::fac
return CGAL::Polygon_mesh_processing::is_degenerate_triangle_face(f, tm, parameters::all_default());
}
/// \ingroup PMP_repairing_grp
/// collects the degenerate faces within a given range of faces.
///
/// @tparam FaceRange a model of `Range` with value type `boost::graph_traits<TriangleMesh>::%face_descriptor`
/// @tparam TriangleMesh a model of `FaceGraph`
/// @tparam NamedParameters a sequence of \ref pmp_namedparameters "Named Parameters"
///
/// @param faces a subset of faces of `tm`
/// @param tm a triangle mesh
/// @param out an output iterator in which the degenerate faces are put
/// @param np optional \ref pmp_namedparameters "Named Parameters" described below
///
/// \cgalNamedParamsBegin
/// \cgalParamBegin{vertex_point_map} the property map with the points associated to the vertices of `tm`.
/// The type of this map is model of `ReadWritePropertyMap`.
/// If this parameter is omitted, an internal property map for
/// `CGAL::vertex_point_t` should be available in `TriangleMesh`
/// \cgalParamEnd
/// \cgalParamBegin{geom_traits} a geometric traits class instance.
/// The traits class must provide the nested functor `Collinear_3`
/// to check whether three points are collinear.
/// \cgalParamEnd
/// \cgalNamedParamsEnd
///
template <typename FaceRange, typename TriangleMesh, typename OutputIterator, typename NamedParameters>
OutputIterator degenerate_faces(const FaceRange& faces,
const TriangleMesh& tm,
OutputIterator out,
const NamedParameters& np)
{
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
for(face_descriptor fd : faces)
{
if(is_degenerate_triangle_face(fd, tm, np))
*out++ = fd;
}
return out;
}
template <typename FaceRange, typename TriangleMesh, typename OutputIterator>
OutputIterator degenerate_faces(const FaceRange& faces,
const TriangleMesh& tm,
OutputIterator out,
typename boost::enable_if<
boost::has_range_iterator<FaceRange>
>::type* = 0)
{
return degenerate_faces(faces, tm, out, CGAL::parameters::all_default());
}
/// \ingroup PMP_repairing_grp
/// calls the function `degenerate_faces()` with the range: `faces(tm)`.
///
/// See above for the comprehensive description of the parameters.
///
template <typename TriangleMesh, typename OutputIterator, typename NamedParameters>
OutputIterator degenerate_faces(const TriangleMesh& tm,
OutputIterator out,
const NamedParameters& np
#ifndef DOXYGEN_RUNNING
,
typename boost::disable_if<
boost::has_range_iterator<TriangleMesh>
>::type* = 0
#endif
)
{
return degenerate_faces(faces(tm), tm, out, np);
}
template <typename TriangleMesh, typename OutputIterator>
OutputIterator degenerate_faces(const TriangleMesh& tm, OutputIterator out)
{
return degenerate_faces(faces(tm), tm, out, CGAL::parameters::all_default());
}
/// \ingroup PMP_repairing_grp
/// checks whether a triangle face is needle.
/// A triangle is said to be a <i>needle</i> if its longest edge is much longer than its shortest edge.