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

separate approximation and mesh extraction class

This commit is contained in:
Lingjie Zhu 2017-07-23 17:41:56 +08:00
parent 75d0b80288
commit f708fe81c1
2 changed files with 58 additions and 700 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

736
Surface_mesh_approximation/include/CGAL/internal/Surface_mesh_approximation/VSA.h
View File

@ -30,23 +30,21 @@ namespace internal
{ {
/** /**
* @brief Main class for Variational Shape Approximation algorithm. * @brief Main class for Variational Shape Approximation algorithm.
*
* @tparam TriangleMesh a CGAL TriangleMesh * @tparam TriangleMesh a CGAL TriangleMesh
* @tparam GeomTraits a model of ApproximationGeomTraits * @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @tparam ApproximationTraits a model of ApproximationGeomTraits
*/ */
template <typename TriangleMesh, template <typename TriangleMesh,
typename ApproximationTrait, typename FacetSegmentMap,
typename FacetAreaMap, typename ApproximationTraits>
typename VertexPointPmap>
class VSA class VSA
{ {
// type definitions // type definitions
private: private:
typedef typename ApproximationTrait::GeomTraits GeomTraits; typedef typename ApproximationTraits::GeomTraits GeomTraits;
typedef typename ApproximationTrait::Proxy Proxy; typedef typename ApproximationTraits::Proxy Proxy;
typedef typename ApproximationTrait::ErrorMetric ErrorMetric; typedef typename ApproximationTraits::ErrorMetric ErrorMetric;
typedef typename ApproximationTrait::ProxyFitting ProxyFitting; typedef typename ApproximationTraits::ProxyFitting ProxyFitting;
typedef typename ApproximationTrait::PlaneFitting PlaneFitting;
typedef typename GeomTraits::FT FT; typedef typename GeomTraits::FT FT;
typedef typename GeomTraits::Point_3 Point_3; typedef typename GeomTraits::Point_3 Point_3;
@ -66,39 +64,16 @@ private:
typedef typename boost::graph_traits<TriangleMesh>::edge_descriptor edge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::edge_descriptor edge_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::edge_iterator edge_iterator; typedef typename boost::graph_traits<TriangleMesh>::edge_iterator edge_iterator;
typedef boost::associative_property_map<std::map<vertex_descriptor, int> > VertexStatusPMap;
typedef boost::associative_property_map<std::map<halfedge_descriptor, int> > HalfedgeStatusPMap;
typedef std::vector<halfedge_descriptor> ChordVector;
typedef typename ChordVector::iterator ChordVectorIterator;
public: public:
enum Initialization { enum Initialization {
RandomInit, RandomInit,
IncrementalInit, IncrementalInit,
HierarchicalInit HierarchicalInit
}; };
// The average positioned anchor attached to a vertex.
struct Anchor {
vertex_descriptor vtx; // The associated vertex.
Point_3 pos; // The position of the anchor.
};
// The border cycle of a region.
// One region may have multiple border cycles.
struct Border {
Border(const halfedge_descriptor &h)
: he_head(h), num_anchors(0) {}
halfedge_descriptor he_head; // The heading halfedge of the border cylce.
std::size_t num_anchors; // The number of anchors on the border.
};
// member variables // member variables
private: private:
const TriangleMesh &mesh; const TriangleMesh &mesh;
const VertexPointPmap &vertex_point_pmap;
Construct_vector_3 vector_functor; Construct_vector_3 vector_functor;
Construct_scaled_vector_3 scale_functor; Construct_scaled_vector_3 scale_functor;
Construct_sum_of_vectors_3 sum_functor; Construct_sum_of_vectors_3 sum_functor;
@ -106,20 +81,6 @@ private:
// Proxy. // Proxy.
std::vector<Proxy> proxies; std::vector<Proxy> proxies;
// The attached anchor index of a vertex.
std::map<vertex_descriptor, int> vertex_status_map;
VertexStatusPMap vertex_status_pmap;
// Mesh facet area map, for anchor position average.
const FacetAreaMap area_pmap;
// All anchors.
std::size_t anchor_index;
std::vector<Anchor> anchors;
// All borders cycles.
std::vector<Border> borders;
// The error metric. // The error metric.
ErrorMetric fit_error; ErrorMetric fit_error;
@ -127,50 +88,32 @@ private:
// The proxy fitting functor. // The proxy fitting functor.
ProxyFitting proxy_fitting; ProxyFitting proxy_fitting;
// The proxy plane fitting functor.
PlaneFitting plane_fitting;
//member functions //member functions
public: public:
/** /**
* Initialize and prepare for the approximation. * Initialize and prepare for the approximation.
* @pre @a TriangleMesh.is_pure_triangle() * @pre @a TriangleMesh.is_pure_triangle()
* @param _mesh `CGAL TriangleMesh` on which approximation operate. * @param _mesh `CGAL TriangleMesh` on which approximation operate.
* @param _vertex_point_map vertex point map of the input mesh. * @param _appx_trait approximation traits object.
* @param _traits geometric trait object.
*/ */
VSA(const TriangleMesh &_mesh, VSA(const TriangleMesh &_mesh, const ApproximationTraits &_appx_trait)
const ApproximationTrait &_appx_trait,
const VertexPointPmap &_vertex_point_map,
const FacetAreaMap &_facet_area_map)
: mesh(_mesh), : mesh(_mesh),
vertex_point_pmap(_vertex_point_map),
area_pmap(_facet_area_map),
vertex_status_pmap(vertex_status_map),
fit_error(_appx_trait.construct_fit_error_functor()), fit_error(_appx_trait.construct_fit_error_functor()),
proxy_fitting(_appx_trait.construct_proxy_fitting_functor()), proxy_fitting(_appx_trait.construct_proxy_fitting_functor()) {
plane_fitting(_appx_trait.construct_plane_fitting_functor()) {
GeomTraits traits; GeomTraits traits;
vector_functor = traits.construct_vector_3_object(); vector_functor = traits.construct_vector_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scalar_product_functor = traits.compute_scalar_product_3_object(); scalar_product_functor = traits.compute_scalar_product_3_object();
// initialize all vertex anchor status
enum Vertex_status { NO_ANCHOR = -1 };
BOOST_FOREACH(vertex_descriptor v, vertices(mesh))
vertex_status_map.insert(std::pair<vertex_descriptor, int>(v, static_cast<int>(NO_ANCHOR)));
} }
/** /**
* Partitions the mesh into the designated number of regions, and stores them in @a seg_map. * Partitions the mesh into the designated number of regions, and stores them in @a seg_map.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param number_of_segments number of designated proxies * @param number_of_segments number of designated proxies
* @param number_of_iterations number of iterations when fitting the proxy * @param number_of_iterations number of iterations when fitting the proxy
* @param[out] seg_map facet partition index * @param[out] seg_map facet partition index
*/ */
template<typename FacetSegmentMap>
void partition(const std::size_t number_of_segments, const std::size_t number_of_iterations, FacetSegmentMap &seg_pmap) { void partition(const std::size_t number_of_segments, const std::size_t number_of_iterations, FacetSegmentMap &seg_pmap) {
random_seed(number_of_segments); random_seed(number_of_segments);
//random_seed(number_of_segments / 2); //random_seed(number_of_segments / 2);
@ -182,12 +125,10 @@ public:
/** /**
* Partitions the mesh incrementally into the designated number of regions, and stores them in @a seg_map. * Partitions the mesh incrementally into the designated number of regions, and stores them in @a seg_map.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param number_of_segments number of designated proxies * @param number_of_segments number of designated proxies
* @param number_of_iterations number of iterations when fitting the proxy * @param number_of_iterations number of iterations when fitting the proxy
* @param[out] seg_map facet partition index * @param[out] seg_map facet partition index
*/ */
template<typename FacetSegmentMap>
void partition_incre(const std::size_t number_of_segments, const std::size_t number_of_iterations, FacetSegmentMap &seg_pmap) { void partition_incre(const std::size_t number_of_segments, const std::size_t number_of_iterations, FacetSegmentMap &seg_pmap) {
// random_seed(number_of_segments); // random_seed(number_of_segments);
random_seed(number_of_segments / 2); random_seed(number_of_segments / 2);
@ -208,12 +149,10 @@ public:
/** /**
* Partitions the mesh into the designated number of regions, and stores them in @a seg_map. * Partitions the mesh into the designated number of regions, and stores them in @a seg_map.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param number_of_segments number of designated proxies * @param number_of_segments number of designated proxies
* @param number_of_iterations number of iterations when fitting the proxy * @param number_of_iterations number of iterations when fitting the proxy
* @param[out] seg_map facet partition index * @param[out] seg_map facet partition index
*/ */
template<typename FacetSegmentMap>
void partition_hierarchical(const std::size_t number_of_segments, const std::size_t number_of_iterations, FacetSegmentMap &seg_pmap) { void partition_hierarchical(const std::size_t number_of_segments, const std::size_t number_of_iterations, FacetSegmentMap &seg_pmap) {
hierarchical_seed(number_of_segments, seg_pmap); hierarchical_seed(number_of_segments, seg_pmap);
for (std::size_t i = 0; i < number_of_iterations; ++i) { for (std::size_t i = 0; i < number_of_iterations; ++i) {
@ -222,104 +161,6 @@ public:
} }
} }
/**
* Extracts the approximated triangle mesh from a partition @a seg_pmap, and stores the triangles in @a tris.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
* @param[out] tris indexed triangles
*/
template<typename FacetSegmentMap>
void extract_mesh(const FacetSegmentMap &seg_pmap, std::vector<int> &tris) {
anchor_index = 0;
find_anchors(seg_pmap);
find_edges(seg_pmap);
add_anchors(seg_pmap);
pseudo_CDT(seg_pmap, tris);
std::vector<Plane_3> px_planes(proxies.size());
std::vector<std::list<face_descriptor> > px_facets(proxies.size());
BOOST_FOREACH(face_descriptor f, faces(mesh))
px_facets[seg_pmap[f]].push_back(f);
for (std::size_t i = 0; i < proxies.size(); ++i)
px_planes[i] = plane_fitting(px_facets[i].begin(), px_facets[i].end());
compute_anchor_position(seg_pmap, px_planes);
std::vector<Point_3> vtx;
BOOST_FOREACH(const Anchor &a, anchors)
vtx.push_back(a.pos);
if (is_manifold_surface(tris, vtx))
std::cout << "Manifold surface." << std::endl;
else
std::cout << "Non-manifold surface." << std::endl;
}
/**
* Use a incremental builder to test if the indexed triangle surface is manifold
* @param tris indexed triangles
* @param vtx vertex positions
* @return true if build successfully
*/
bool is_manifold_surface(const std::vector<int> &tris, const std::vector<Point_3> &vtx) {
typedef CGAL::Polyhedron_3<GeomTraits> PolyhedronSurface;
typedef typename PolyhedronSurface::HalfedgeDS HDS;
HDS hds;
CGAL::Polyhedron_incremental_builder_3<HDS> builder(hds, true);
builder.begin_surface(vtx.size(), tris.size() / 3);
BOOST_FOREACH(const Point_3 &v, vtx)
builder.add_vertex(v);
for (std::vector<int>::const_iterator itr = tris.begin(); itr != tris.end(); itr += 3) {
if (builder.test_facet(itr, itr + 3)) {
builder.begin_facet();
builder.add_vertex_to_facet(*itr);
builder.add_vertex_to_facet(*(itr + 1));
builder.add_vertex_to_facet(*(itr + 2));
builder.end_facet();
}
else {
// std::cerr << "test_facet failed" << std::endl;
builder.end_surface();
return false;
}
}
builder.end_surface();
return true;
}
/**
* Collect the anchors.
* @return vector of anchors
*/
std::vector<Anchor> collect_anchors() {
return anchors;
}
/**
* Collect the partition @a seg_pmap approximated borders.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
* @return anchor indexes of each border
*/
template<typename FacetSegmentMap>
std::vector<std::vector<std::size_t> >
collect_borders(const FacetSegmentMap &seg_pmap) {
std::vector<std::vector<std::size_t> > bdrs;
for (typename std::vector<Border>::iterator bitr = borders.begin();
bitr != borders.end(); ++bitr) {
std::vector<std::size_t> bdr;
const halfedge_descriptor he_mark = bitr->he_head;
halfedge_descriptor he = he_mark;
do {
ChordVector chord;
walk_to_next_anchor(he, chord, seg_pmap);
bdr.push_back(vertex_status_pmap[target(he, mesh)]);
} while(he != he_mark);
bdrs.push_back(bdr);
}
return bdrs;
}
private: private:
/** /**
* Random initialize proxies. * Random initialize proxies.
@ -343,11 +184,9 @@ private:
/** /**
* Hierarchical initialize proxies. * Hierarchical initialize proxies.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param initial_px number of proxies * @param initial_px number of proxies
* @param[out] seg_map facet partition index * @param[out] seg_map facet partition index
*/ */
template<typename FacetSegmentMap>
void hierarchical_seed(const std::size_t initial_px, FacetSegmentMap &seg_pmap) { void hierarchical_seed(const std::size_t initial_px, FacetSegmentMap &seg_pmap) {
if (initial_px == 0) if (initial_px == 0)
return; return;
@ -435,10 +274,8 @@ private:
/** /**
* Calculates and updates the best fitting proxies of a given partition @a seg_pmap. * Calculates and updates the best fitting proxies of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index * @param seg_map facet partition index
*/ */
template <typename FacetSegmentMap>
void fitting(const FacetSegmentMap &seg_pmap) { void fitting(const FacetSegmentMap &seg_pmap) {
std::vector<std::list<face_descriptor> > px_facets(proxies.size()); std::vector<std::list<face_descriptor> > px_facets(proxies.size());
BOOST_FOREACH(face_descriptor f, faces(mesh)) BOOST_FOREACH(face_descriptor f, faces(mesh))
@ -450,10 +287,8 @@ private:
/** /**
* Inserts a proxy of a given partition @a seg_pmap at the furthest facet of the region with the maximum fitting error. * Inserts a proxy of a given partition @a seg_pmap at the furthest facet of the region with the maximum fitting error.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index * @param seg_map facet partition index
*/ */
template <typename FacetSegmentMap>
void insert_proxy(const FacetSegmentMap &seg_pmap) { void insert_proxy(const FacetSegmentMap &seg_pmap) {
std::vector<FT> px_error(proxies.size(), FT(0.0)); std::vector<FT> px_error(proxies.size(), FT(0.0));
std::vector<FT> max_facet_error(proxies.size(), FT(0.0)); std::vector<FT> max_facet_error(proxies.size(), FT(0.0));
@ -488,12 +323,10 @@ private:
* Add proxies by diffusing fitting error into current partitions. * Add proxies by diffusing fitting error into current partitions.
* Each partition is added with the number of proxies in proportional to its fitting error. * Each partition is added with the number of proxies in proportional to its fitting error.
* Note that the number of inserted proxies doesn't necessarily equal the requested number. * Note that the number of inserted proxies doesn't necessarily equal the requested number.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param num_proxies_to_be_added added number of proxies * @param num_proxies_to_be_added added number of proxies
* @param seg_map facet partition index * @param seg_map facet partition index
* @return inserted number of proxies * @return inserted number of proxies
*/ */
template<typename FacetSegmentMap>
std::size_t insert_proxy_error_diffusion(const std::size_t num_proxies_to_be_added, const FacetSegmentMap &seg_pmap) { std::size_t insert_proxy_error_diffusion(const std::size_t num_proxies_to_be_added, const FacetSegmentMap &seg_pmap) {
struct ProxyError { struct ProxyError {
ProxyError(const std::size_t &id, const FT &er) ProxyError(const std::size_t &id, const FT &er)
@ -556,280 +389,11 @@ private:
return num_inserted; return num_inserted;
} }
/**
* Finds the anchors of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
*/
template<typename FacetSegmentMap>
void find_anchors(const FacetSegmentMap &seg_pmap) {
anchors.clear();
BOOST_FOREACH(vertex_descriptor vtx, vertices(mesh)) {
std::size_t border_count = 0;
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(vtx, mesh)) {
if (CGAL::is_border_edge(h, mesh))
++border_count;
else if (seg_pmap[face(h, mesh)] != seg_pmap[face(opposite(h, mesh), mesh)])
++border_count;
}
if (border_count >= 3)
attach_anchor(vtx);
}
}
/**
* Finds and approximates the edges connecting the anchors of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
*/
template<typename FacetSegmentMap>
void find_edges(const FacetSegmentMap &seg_pmap) {
// collect candidate halfedges in a set
std::set<halfedge_descriptor> he_candidates;
BOOST_FOREACH(halfedge_descriptor h, halfedges(mesh)) {
if (!CGAL::is_border(h, mesh)
&& (CGAL::is_border(opposite(h, mesh), mesh)
|| seg_pmap[face(h, mesh)] != seg_pmap[face(opposite(h, mesh), mesh)]))
he_candidates.insert(h);
}
// pick up one candidate halfedge each time and traverse the connected border
borders.clear();
while (!he_candidates.empty()) {
halfedge_descriptor he_start = *he_candidates.begin();
walk_to_first_anchor(he_start, seg_pmap);
// no anchor in this connected border, make a new anchor
if (!is_anchor_attached(he_start))
attach_anchor(he_start);
// a new connected border
borders.push_back(Border(he_start));
std::cerr << "#border " << borders.size() << std::endl;
const halfedge_descriptor he_mark = he_start;
do {
ChordVector chord;
walk_to_next_anchor(he_start, chord, seg_pmap);
borders.back().num_anchors += subdivide_chord(chord.begin(), chord.end(), seg_pmap);
std::cerr << "#chord_anchor " << borders.back().num_anchors << std::endl;
for (ChordVectorIterator citr = chord.begin(); citr != chord.end(); ++citr)
he_candidates.erase(*citr);
} while (he_start != he_mark);
}
}
/**
* Adds anchors to the border cycles with only 2 anchors of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
*/
template<typename FacetSegmentMap>
void add_anchors(const FacetSegmentMap &seg_pmap) {
typedef typename std::vector<Border>::iterator BorderIterator;
for (BorderIterator bitr = borders.begin(); bitr != borders.end(); ++bitr) {
if (bitr->num_anchors > 2)
continue;
// 2 initial anchors at least
CGAL_assertion(bitr->num_anchors == 2);
// borders with only 2 initial anchors
const halfedge_descriptor he_mark = bitr->he_head;
Point_3 pt_begin = vertex_point_pmap[target(he_mark, mesh)];
Point_3 pt_end = pt_begin;
halfedge_descriptor he = he_mark;
ChordVector chord;
std::size_t count = 0;
do {
walk_to_next_border_halfedge(he, seg_pmap);
if (!is_anchor_attached(he))
chord.push_back(he);
else {
if (count == 0)
pt_end = vertex_point_pmap[target(he, mesh)];
++count;
}
} while(he != he_mark);
// anchor count may be increased to more than 2 afterwards
// due to the new anchors added by the neighboring border (< 2 anchors)
if (count > 2) {
bitr->num_anchors = count;
continue;
}
FT dist_max(0.0);
halfedge_descriptor he_max;
Vector_3 chord_vec = vector_functor(pt_begin, pt_end);
chord_vec = scale_functor(chord_vec,
FT(1.0 / std::sqrt(CGAL::to_double(chord_vec.squared_length()))));
for (ChordVectorIterator citr = chord.begin(); citr != chord.end(); ++citr) {
Vector_3 vec = vector_functor(pt_begin, vertex_point_pmap[target(*citr, mesh)]);
vec = CGAL::cross_product(chord_vec, vec);
FT dist(std::sqrt(CGAL::to_double(vec.squared_length())));
if (dist > dist_max) {
dist_max = dist;
he_max = *citr;
}
}
attach_anchor(he_max);
// increase border anchors by one
bitr->num_anchors++;
}
}
/**
* Runs the pseudo Constrained Delaunay Triangulation at each region of a given partition @a seg_pmap, and stores the extracted indexed triangles in @a tris.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
* @param tris extracted tirangles, index of anchors
*/
template<typename FacetSegmentMap>
void pseudo_CDT(const FacetSegmentMap &seg_pmap,
std::vector<int> &tris) {
// subgraph attached with vertex anchor status and edge weight
typedef boost::property<boost::vertex_index1_t, int,
boost::property<boost::vertex_index2_t, int> > VertexProperty;
typedef boost::property<boost::edge_weight_t, FT,
boost::property<boost::edge_index_t, int> > EdgeProperty;
typedef boost::subgraph<boost::adjacency_list<
boost::listS, boost::vecS,
boost::undirectedS,
VertexProperty, EdgeProperty> > SubGraph;
typedef typename boost::property_map<SubGraph, boost::vertex_index1_t>::type VertexIndex1Map;
typedef typename boost::property_map<SubGraph, boost::vertex_index2_t>::type VertexIndex2Map;
typedef typename boost::property_map<SubGraph, boost::edge_weight_t>::type EdgeWeightMap;
typedef typename SubGraph::vertex_descriptor sg_vertex_descriptor;
typedef typename SubGraph::edge_descriptor sg_edge_descriptor;
typedef std::vector<sg_vertex_descriptor> VertexVector;
typedef std::map<vertex_descriptor, sg_vertex_descriptor> VertexMap;
typedef boost::associative_property_map<VertexMap> ToSGVertexMap;
VertexMap vmap;
ToSGVertexMap to_sgv_map(vmap);
// mapping the TriangleMesh mesh into a SubGraph
SubGraph gmain;
VertexIndex1Map global_vanchor_map = get(boost::vertex_index1, gmain);
VertexIndex2Map global_vtag_map = get(boost::vertex_index2, gmain);
EdgeWeightMap global_eweight_map = get(boost::edge_weight, gmain);
BOOST_FOREACH(vertex_descriptor v, vertices(mesh)) {
sg_vertex_descriptor sgv = add_vertex(gmain);
global_vanchor_map[sgv] = vertex_status_pmap[v];
global_vtag_map[sgv] = vertex_status_pmap[v];
vmap.insert(std::pair<vertex_descriptor, sg_vertex_descriptor>(v, sgv));
}
BOOST_FOREACH(edge_descriptor e, edges(mesh)) {
vertex_descriptor vs = source(e, mesh);
vertex_descriptor vt = target(e, mesh);
FT len(std::sqrt(CGAL::to_double(
CGAL::squared_distance(vertex_point_pmap[vs], vertex_point_pmap[vt]))));
add_edge(to_sgv_map[vs], to_sgv_map[vt], len, gmain);
}
std::vector<VertexVector> vertex_patches(proxies.size());
BOOST_FOREACH(vertex_descriptor v, vertices(mesh)) {
std::set<std::size_t> px_set;
BOOST_FOREACH(face_descriptor f, faces_around_target(halfedge(v, mesh), mesh)) {
if (f != boost::graph_traits<TriangleMesh>::null_face())
px_set.insert(seg_pmap[f]);
}
BOOST_FOREACH(std::size_t p, px_set)
vertex_patches[p].push_back(to_sgv_map[v]);
}
BOOST_FOREACH(VertexVector &vpatch, vertex_patches) {
// add a super vertex connecting to its boundary anchors into the main graph
const sg_vertex_descriptor superv = add_vertex(gmain);
global_vanchor_map[superv] = 0;
global_vtag_map[superv] = 0;
BOOST_FOREACH(sg_vertex_descriptor v, vpatch) {
if (is_anchor_attached(v, global_vanchor_map))
add_edge(superv, v, FT(0), gmain);
}
vpatch.push_back(superv);
}
// multi-source Dijkstra's shortest path algorithm applied to each proxy patch
BOOST_FOREACH(VertexVector &vpatch, vertex_patches) {
// construct subgraph
SubGraph &glocal = gmain.create_subgraph();
BOOST_FOREACH(sg_vertex_descriptor v, vpatch)
add_vertex(v, glocal);
// most subgraph functions work with local descriptors
VertexIndex1Map local_vanchor_map = get(boost::vertex_index1, glocal);
VertexIndex2Map local_vtag_map = get(boost::vertex_index2, glocal);
EdgeWeightMap local_eweight_map = get(boost::edge_weight, glocal);
const sg_vertex_descriptor source = glocal.global_to_local(vpatch.back());
VertexVector pred(num_vertices(glocal));
boost::dijkstra_shortest_paths(glocal, source,
boost::predecessor_map(&pred[0]).weight_map(local_eweight_map));
// backtrack to the anchor and tag each vertex in the local patch graph
BOOST_FOREACH(sg_vertex_descriptor v, vertices(glocal)) {
sg_vertex_descriptor curr = v;
while (!is_anchor_attached(curr, local_vanchor_map))
curr = pred[curr];
local_vtag_map[v] = local_vanchor_map[curr];
}
}
// tag all boundary chord
BOOST_FOREACH(const Border &bdr, borders) {
const halfedge_descriptor he_mark = bdr.he_head;
halfedge_descriptor he = he_mark;
do {
ChordVector chord;
walk_to_next_anchor(he, chord, seg_pmap);
std::vector<FT> vdist;
vdist.push_back(FT(0));
BOOST_FOREACH(halfedge_descriptor h, chord) {
FT elen = global_eweight_map[edge(
to_sgv_map[source(h, mesh)],
to_sgv_map[target(h, mesh)],
gmain).first];
vdist.push_back(vdist.back() + elen);
}
FT half_chord_len = vdist.back() / FT(2);
const int anchorleft = vertex_status_pmap[source(chord.front(), mesh)];
const int anchorright = vertex_status_pmap[target(chord.back(), mesh)];
typename std::vector<FT>::iterator ditr = vdist.begin() + 1;
for (typename ChordVector::iterator hitr = chord.begin();
hitr != chord.end() - 1; ++hitr, ++ditr) {
if (*ditr < half_chord_len)
global_vtag_map[to_sgv_map[target(*hitr, mesh)]] = anchorleft;
else
global_vtag_map[to_sgv_map[target(*hitr, mesh)]] = anchorright;
}
} while(he != he_mark);
}
// collect triangles
BOOST_FOREACH(face_descriptor f, faces(mesh)) {
halfedge_descriptor he = halfedge(f, mesh);
int i = global_vtag_map[to_sgv_map[source(he, mesh)]];
int j = global_vtag_map[to_sgv_map[target(he, mesh)]];
int k = global_vtag_map[to_sgv_map[target(next(he, mesh), mesh)]];
if (i != j && i != k && j != k) {
tris.push_back(i);
tris.push_back(j);
tris.push_back(k);
}
}
}
/** /**
* Computes fitting error of a given partition @a seg_pmap. * Computes fitting error of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index * @param seg_map facet partition index
* @return total fitting error * @return total fitting error
*/ */
template<typename FacetSegmentMap>
FT fitting_error(const FacetSegmentMap &seg_pmap) { FT fitting_error(const FacetSegmentMap &seg_pmap) {
FT sum_error(0); FT sum_error(0);
BOOST_FOREACH(face_descriptor f, faces(mesh)) BOOST_FOREACH(face_descriptor f, faces(mesh))
@ -839,12 +403,10 @@ private:
/** /**
* Computes fitting error of a given partition @a seg_pmap. * Computes fitting error of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index * @param seg_map facet partition index
* @param px_error vector of error of each proxy * @param px_error vector of error of each proxy
* @return total fitting error * @return total fitting error
*/ */
template<typename FacetSegmentMap>
FT fitting_error(const FacetSegmentMap &seg_pmap, std::vector<FT> &px_error) { FT fitting_error(const FacetSegmentMap &seg_pmap, std::vector<FT> &px_error) {
FT sum_error(0); FT sum_error(0);
BOOST_FOREACH(face_descriptor f, faces(mesh)) { BOOST_FOREACH(face_descriptor f, faces(mesh)) {
@ -855,240 +417,6 @@ private:
} }
return sum_error; return sum_error;
} }
/**
* Computes and updates the proxy areas of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param seg_map facet partition index
* @param proxies_area proxy areas
*/
template<typename FacetSegmentMap>
void compute_proxy_area(const FacetSegmentMap &seg_pmap, std::vector<FT> &proxies_area) {
std::vector<FT> areas(proxies.size(), FT(0));
BOOST_FOREACH(face_descriptor f, faces(mesh)) {
areas[seg_pmap[f]] += area_pmap[f];
}
proxies_area.swap(areas);
}
/**
* Walks along the region border to the first halfedge pointing to a vertex associated with an anchor of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param[in/out] he_start region border halfedge
* @param seg_map facet partition index
*/
template<typename FacetSegmentMap>
void walk_to_first_anchor(halfedge_descriptor &he_start, const FacetSegmentMap &seg_pmap) {
const halfedge_descriptor start_mark = he_start;
while (!is_anchor_attached(he_start)) {
// no anchor attached to the halfedge target
walk_to_next_border_halfedge(he_start, seg_pmap);
if (he_start == start_mark) // back to where started, a circular border
return;
}
}
/**
* Walks along the region border to the next anchor and records the path as @a chord of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param[in/out] he_start starting region border halfedge pointing to a vertex associated with an anchor
* @param[out] chord recorded path chord
* @param seg_map facet partition index
*/
template<typename FacetSegmentMap>
void walk_to_next_anchor(
halfedge_descriptor &he_start,
ChordVector &chord,
const FacetSegmentMap &seg_pmap) {
do {
walk_to_next_border_halfedge(he_start, seg_pmap);
chord.push_back(he_start);
} while (!is_anchor_attached(he_start));
}
/**
* Walks to next border halfedge of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param[in/out] he_start region border halfedge
* @param seg_map facet partition index
*/
template<typename FacetSegmentMap>
void walk_to_next_border_halfedge(halfedge_descriptor &he_start, const FacetSegmentMap &seg_pmap) {
const std::size_t px_idx = seg_pmap[face(he_start, mesh)];
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(he_start, mesh)) {
if (CGAL::is_border(h, mesh) || seg_pmap[face(h, mesh)] != px_idx) {
he_start = opposite(h, mesh);
return;
}
}
}
/**
* Subdivides a chord recursively in range [@a chord_begin, @a chord_end) of a given partition @a seg_pmap.
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
* @param chord_begin begin iterator of the chord
* @param chord_end end iterator of the chord
* @param seg_map facet partition index
* @return the number of anchors of the chord apart from the first one
*/
template<typename FacetSegmentMap>
std::size_t subdivide_chord(
const ChordVectorIterator &chord_begin,
const ChordVectorIterator &chord_end,
const FacetSegmentMap &seg_pmap,
const FT thre = FT(0.2)) {
const std::size_t chord_size = std::distance(chord_begin, chord_end);
// do not subdivide trivial chord
if (chord_size < 4)
return 1;
halfedge_descriptor he_start = *chord_begin;
std::size_t px_left = seg_pmap[face(he_start, mesh)];
std::size_t px_right = px_left;
if (!CGAL::is_border(opposite(he_start, mesh), mesh))
px_right = seg_pmap[face(opposite(he_start, mesh), mesh)];
// suppose the proxy normal angle is acute
FT norm_sin(1.0);
if (!CGAL::is_border(opposite(he_start, mesh), mesh)) {
Vector_3 vec = CGAL::cross_product(proxies[px_left].normal, proxies[px_right].normal);
norm_sin = FT(std::sqrt(CGAL::to_double(scalar_product_functor(vec, vec))));
}
FT criterion = thre + FT(1.0);
ChordVectorIterator he_max;
const ChordVectorIterator chord_last = chord_end - 1;
std::size_t anchor_begin = vertex_status_pmap[source(he_start, mesh)];
std::size_t anchor_end = vertex_status_pmap[target(*chord_last, mesh)];
const Point_3 &pt_begin = vertex_point_pmap[source(he_start, mesh)];
const Point_3 &pt_end = vertex_point_pmap[target(*chord_last, mesh)];
if (anchor_begin == anchor_end) {
// circular chord
CGAL_assertion(chord_size > 2);
// if (chord_size < 3)
// return;
FT dist_max(0.0);
for (ChordVectorIterator citr = chord_begin; citr != chord_last; ++citr) {
FT dist = CGAL::squared_distance(pt_begin, vertex_point_pmap[target(*citr, mesh)]);
dist = FT(std::sqrt(CGAL::to_double(dist)));
if (dist > dist_max) {
he_max = citr;
dist_max = dist;
}
}
}
else {
FT dist_max(0.0);
Vector_3 chord_vec = vector_functor(pt_begin, pt_end);
FT chord_len(std::sqrt(CGAL::to_double(chord_vec.squared_length())));
chord_vec = scale_functor(chord_vec, FT(1.0) / chord_len);
for (ChordVectorIterator citr = chord_begin; citr != chord_last; ++citr) {
Vector_3 vec = vector_functor(pt_begin, vertex_point_pmap[target(*citr, mesh)]);
vec = CGAL::cross_product(chord_vec, vec);
FT dist(std::sqrt(CGAL::to_double(vec.squared_length())));
if (dist > dist_max) {
he_max = citr;
dist_max = dist;
}
}
criterion = dist_max * norm_sin / chord_len;
}
if (criterion > thre) {
// subdivide at the most remote vertex
attach_anchor(*he_max);
std::size_t num0 = subdivide_chord(chord_begin, he_max + 1, seg_pmap);
std::size_t num1 = subdivide_chord(he_max + 1, chord_end, seg_pmap);
return num0 + num1;
}
return 1;
}
/**
* Check if the target vertex of a halfedge is attached with an anchor.
* @param he halfedge
*/
bool is_anchor_attached(const halfedge_descriptor &he) {
return is_anchor_attached(target(he, mesh), vertex_status_pmap);
}
/**
* Check if a vertex is attached with an anchor.
* @tparam VertexAnchorIndexMap `WritablePropertyMap` with `boost::graph_traights<TriangleMesh>::vertex_descriptor` as key and `std::size_t` as value type
* @param v vertex
* @param vertex_anchor_map vertex anchor index map
*/
template<typename VertexAnchorIndexMap>
bool is_anchor_attached(
const typename boost::property_traits<VertexAnchorIndexMap>::key_type &v,
const VertexAnchorIndexMap &vertex_anchor_map) {
return vertex_anchor_map[v] >= 0;
}
/**
* Attachs an anchor to the vertex.
* @param vtx vertex
*/
void attach_anchor(const vertex_descriptor &vtx) {
vertex_status_pmap[vtx] = static_cast<int>(anchor_index++);
}
/**
* Attachs an anchor to the target vertex of the halfedge.
* @param he halfedge
*/
void attach_anchor(const halfedge_descriptor &he) {
vertex_descriptor vtx = target(he, mesh);
attach_anchor(vtx);
}
/**
* Calculate the anchor positions.
* @param
*/
template<typename FacetSegmentMap>
void compute_anchor_position(
const FacetSegmentMap &seg_pmap,
const std::vector<Plane_3> &px_planes) {
std::vector<FT> proxies_area; // The proxy area.
compute_proxy_area(seg_pmap, proxies_area);
anchors = std::vector<Anchor>(anchor_index);
BOOST_FOREACH(vertex_descriptor v, vertices(mesh)) {
if (is_anchor_attached(v, vertex_status_pmap)) {
// construct an anchor from vertex and the incident proxies
std::set<std::size_t> px_set;
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(v, mesh)) {
if (!CGAL::is_border(h, mesh))
px_set.insert(seg_pmap[face(h, mesh)]);
}
// construct an anchor from vertex and the incident proxies
FT avgx(0), avgy(0), avgz(0), sum_area(0);
const Point_3 vtx_pt = vertex_point_pmap[v];
for (std::set<std::size_t>::iterator pxitr = px_set.begin();
pxitr != px_set.end(); ++pxitr) {
std::size_t px_idx = *pxitr;
Point_3 proj = px_planes[px_idx].projection(vtx_pt);
FT area = proxies_area[px_idx];
avgx += proj.x() * area;
avgy += proj.y() * area;
avgz += proj.z() * area;
sum_area += area;
}
Point_3 pos = Point_3(avgx / sum_area, avgy / sum_area, avgz / sum_area);
std::size_t aidx = vertex_status_pmap[v];
anchors[aidx].vtx = v;
anchors[aidx].pos = pos;
}
}
}
}; // end class VSA }; // end class VSA
/** /**
@ -1097,14 +425,15 @@ private:
* @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type * @tparam FacetSegmentMap `WritablePropertyMap` with `boost::graph_traits<TriangleMesh>::face_handle` as key and `std::size_t` as value type
*/ */
template <typename TriangleMesh, template <typename TriangleMesh,
typename ApproximationTrait, typename ApproximationTraits,
typename VertexPointMap, typename VertexPointMap,
typename FacetSegmentMap, typename FacetSegmentMap,
typename FacetAreaMap> typename FacetAreaMap>
class VSA_mesh_extraction class VSA_mesh_extraction
{ {
typedef typename ApproximationTrait::GeomTraits GeomTraits; typedef typename ApproximationTraits::GeomTraits GeomTraits;
typedef typename ApproximationTrait::PlaneFitting PlaneFitting; typedef typename ApproximationTraits::Proxy Proxy;
typedef typename ApproximationTraits::PlaneFitting PlaneFitting;
typedef typename GeomTraits::FT FT; typedef typename GeomTraits::FT FT;
typedef typename GeomTraits::Point_3 Point_3; typedef typename GeomTraits::Point_3 Point_3;
@ -1117,6 +446,8 @@ template <typename TriangleMesh,
typedef typename boost::graph_traits<TriangleMesh>::vertex_descriptor vertex_descriptor; typedef typename boost::graph_traits<TriangleMesh>::vertex_descriptor vertex_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::edge_descriptor edge_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typedef boost::associative_property_map<std::map<vertex_descriptor, int> > VertexAnchorMap; typedef boost::associative_property_map<std::map<vertex_descriptor, int> > VertexAnchorMap;
@ -1143,15 +474,18 @@ public:
/** /**
* Extracts the surface mesh from an approximation partition @a _seg_pmap of mesh @a _mesh. * Extracts the surface mesh from an approximation partition @a _seg_pmap of mesh @a _mesh.
* @param _mesh the approximated triangle mesh. * @param _mesh the approximated triangle mesh.
* @param _appx_trait approximation traits object.
* @param _point_pmap vertex point map of the input mesh.
* @param _seg_pmap approximation partition. * @param _seg_pmap approximation partition.
* @param _area_pmap facet area map of the input mesh.
*/ */
VSA_mesh_extraction(const TriangleMesh &_mesh, VSA_mesh_extraction(const TriangleMesh &_mesh,
const ApproximationTrait &_appx_trait, const ApproximationTraits &_appx_trait,
const VertexPointMap &_point_pmap, const VertexPointMap &_point_pmap,
const FacetSegmentMap &_seg_pmap, const FacetSegmentMap &_seg_pmap,
const FacetAreaMap &_area_pmap) const FacetAreaMap &_area_pmap)
: mesh(_mesh), : mesh(_mesh),
point_pmap(_point_pmap) point_pmap(_point_pmap),
seg_pmap(_seg_pmap), seg_pmap(_seg_pmap),
area_pmap(_area_pmap), area_pmap(_area_pmap),
vanchor_map(vertex_int_map), vanchor_map(vertex_int_map),
@ -1166,13 +500,26 @@ public:
enum Vertex_status { NO_ANCHOR = -1 }; enum Vertex_status { NO_ANCHOR = -1 };
BOOST_FOREACH(vertex_descriptor v, vertices(mesh)) BOOST_FOREACH(vertex_descriptor v, vertices(mesh))
vertex_int_map.insert(std::pair<vertex_descriptor, int>(v, static_cast<int>(NO_ANCHOR))); vertex_int_map.insert(std::pair<vertex_descriptor, int>(v, static_cast<int>(NO_ANCHOR)));
std::size_t num_proxies = 0;
BOOST_FOREACH(face_descriptor f, faces(mesh))
num_proxies = num_proxies < seg_pmap[f] ? seg_pmap[f] : num_proxies;
num_proxies++;
std::vector<std::list<face_descriptor> > px_facets(num_proxies);
BOOST_FOREACH(face_descriptor f, faces(mesh))
px_facets[seg_pmap[f]].push_back(f);
typename ApproximationTraits::ProxyFitting proxy_fitting = _appx_trait.construct_proxy_fitting_functor();
for (std::size_t i = 0; i < num_proxies; ++i)
proxies.push_back(proxy_fitting(px_facets[i].begin(), px_facets[i].end()));
std::cerr << "#proxies " << num_proxies << ", " << proxies.size() << std::endl;
} }
/** /**
* Extracts the approximated triangle mesh in @a tris. * Extracts the approximated triangle mesh in @a tris.
* @param[out] tris indexed triangles * @param[out] tris indexed triangles
*/ */
template<typename FacetSegmentMap>
void extract_mesh(std::vector<int> &tris) { void extract_mesh(std::vector<int> &tris) {
anchor_index = 0; anchor_index = 0;
find_anchors(); find_anchors();
@ -1305,7 +652,7 @@ private:
do { do {
ChordVector chord; ChordVector chord;
walk_to_next_anchor(he_start, chord); walk_to_next_anchor(he_start, chord);
borders.back().num_anchors += subdivide_chord(chord.begin(), chord.end(), seg_pmap); borders.back().num_anchors += subdivide_chord(chord.begin(), chord.end());
std::cerr << "#chord_anchor " << borders.back().num_anchors << std::endl; std::cerr << "#chord_anchor " << borders.back().num_anchors << std::endl;
for (ChordVectorIterator citr = chord.begin(); citr != chord.end(); ++citr) for (ChordVectorIterator citr = chord.begin(); citr != chord.end(); ++citr)
@ -1624,8 +971,8 @@ private:
// subdivide at the most remote vertex // subdivide at the most remote vertex
attach_anchor(*he_max); attach_anchor(*he_max);
std::size_t num0 = subdivide_chord(chord_begin, he_max + 1, seg_pmap); std::size_t num0 = subdivide_chord(chord_begin, he_max + 1);
std::size_t num1 = subdivide_chord(he_max + 1, chord_end, seg_pmap); std::size_t num1 = subdivide_chord(he_max + 1, chord_end);
return num0 + num1; return num0 + num1;
} }
@ -1675,7 +1022,7 @@ private:
* Computes and the proxy fitting planes. * Computes and the proxy fitting planes.
* @param px_planes proxy planes. * @param px_planes proxy planes.
*/ */
void compute_proxy_planes(std::vector<FT> &px_planes) { void compute_proxy_planes(std::vector<Plane_3> &px_planes) {
std::vector<std::list<face_descriptor> > px_facets(proxies.size()); std::vector<std::list<face_descriptor> > px_facets(proxies.size());
BOOST_FOREACH(face_descriptor f, faces(mesh)) BOOST_FOREACH(face_descriptor f, faces(mesh))
px_facets[seg_pmap[f]].push_back(f); px_facets[seg_pmap[f]].push_back(f);
@ -1745,6 +1092,9 @@ private:
Construct_sum_of_vectors_3 sum_functor; Construct_sum_of_vectors_3 sum_functor;
Compute_scalar_product_3 scalar_product_functor; Compute_scalar_product_3 scalar_product_functor;
// proxy
std::vector<Proxy> proxies;
// The attached anchor index of a vertex. // The attached anchor index of a vertex.
std::map<vertex_descriptor, int> vertex_int_map; std::map<vertex_descriptor, int> vertex_int_map;
VertexAnchorMap vanchor_map; VertexAnchorMap vanchor_map;

22
Surface_mesh_approximation/include/CGAL/vsa_mesh_approximation.h
View File

@ -70,11 +70,10 @@ template<typename TriangleMesh,
typedef CGAL::internal::VSA< typedef CGAL::internal::VSA<
TriangleMesh, TriangleMesh,
ApproximationTrait, SegmentPropertyMap,
FacetAreaMap, ApproximationTrait> VSA;
PointPropertyMap> VSA;
VSA algorithm(tm, app_trait, ppmap, area_pmap); VSA algorithm(tm, app_trait);
switch (init) { switch (init) {
case VSA::RandomInit: case VSA::RandomInit:
@ -88,13 +87,22 @@ template<typename TriangleMesh,
break; break;
} }
algorithm.extract_mesh(segment_ids, tris); typedef CGAL::internal::VSA_mesh_extraction<
BOOST_FOREACH(const typename VSA::Anchor &a, algorithm.collect_anchors()) { TriangleMesh,
ApproximationTrait,
PointPropertyMap,
SegmentPropertyMap,
FacetAreaMap> VSA_mesh_extraction;
VSA_mesh_extraction extractor(tm, app_trait, ppmap, segment_ids, area_pmap);
extractor.extract_mesh(tris);
BOOST_FOREACH(const typename VSA_mesh_extraction::Anchor &a, extractor.collect_anchors()) {
vtx.push_back(a.vtx); vtx.push_back(a.vtx);
pos.push_back(a.pos); pos.push_back(a.pos);
} }
bdrs = algorithm.collect_borders(segment_ids); bdrs = extractor.collect_borders();
} }
} }