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

Use dynamic_property_map

This commit is contained in:
Andreas Fabri 2017-11-27 12:55:20 +00:00
parent 22661170f8
commit 0cb55dfd8d
4 changed files with 173 additions and 122 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Surface_mesh_approximation/examples/Surface_mesh_approximation/CMakeLists.txt
View File

@ -47,6 +47,8 @@ create_single_source_cgal_program( "vsa_approximation_2_example.cpp" )
create_single_source_cgal_program( "vsa_approximation_example.cpp" ) create_single_source_cgal_program( "vsa_approximation_example.cpp" )
create_single_source_cgal_program( "vsa_approximation_SM_example.cpp" )
create_single_source_cgal_program( "vsa_class_interface_example.cpp" ) create_single_source_cgal_program( "vsa_class_interface_example.cpp" )
create_single_source_cgal_program( "vsa_isotropic_metric_example.cpp" ) create_single_source_cgal_program( "vsa_isotropic_metric_example.cpp" )

2
Surface_mesh_approximation/examples/Surface_mesh_approximation/vsa_segmentation_example.cpp
View File

@ -7,7 +7,7 @@
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef CGAL::Polyhedron_3<Kernel> Polyhedron; typedef CGAL::Polyhedron_3<Kernel> Polyhedron;
typedef typename boost::graph_traits<Polyhedron>::face_descriptor face_descriptor; typedef boost::graph_traits<Polyhedron>::face_descriptor face_descriptor;
typedef boost::unordered_map<face_descriptor, std::size_t> Facet_index_map; typedef boost::unordered_map<face_descriptor, std::size_t> Facet_index_map;
typedef boost::associative_property_map<Facet_index_map> Facet_proxy_pmap; typedef boost::associative_property_map<Facet_index_map> Facet_proxy_pmap;

143
Surface_mesh_approximation/include/CGAL/vsa_approximation.h
View File

@ -14,14 +14,12 @@
#include <CGAL/vsa_metrics.h> #include <CGAL/vsa_metrics.h>
#include <CGAL/Default.h> #include <CGAL/Default.h>
#include <boost/unordered_map.hpp>
#include <boost/graph/graph_traits.hpp> #include <boost/graph/graph_traits.hpp>
#include <boost/graph/adjacency_list.hpp> #include <boost/graph/adjacency_list.hpp>
#include <boost/graph/dijkstra_shortest_paths.hpp> #include <boost/graph/dijkstra_shortest_paths.hpp>
#include <boost/graph/subgraph.hpp> #include <boost/graph/subgraph.hpp>
#include <boost/foreach.hpp> #include <boost/foreach.hpp>
#include <boost/optional.hpp> #include <boost/optional.hpp>
#include <vector> #include <vector>
#include <queue> #include <queue>
#include <iterator> #include <iterator>
@ -102,7 +100,12 @@ private:
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor; typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
// internal typedefs // internal typedefs
typedef boost::associative_property_map<boost::unordered_map<vertex_descriptor, std::size_t> > Vertex_anchor_map;
typedef CGAL::internal::vertex_property_t<std::size_t> Vertex_anchor_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Vertex_anchor_tag >::type Vertex_anchor_map;
typedef CGAL::internal::face_property_t<std::size_t> Face_proxy_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_proxy_tag >::type Fproxy_map;
typedef std::vector<halfedge_descriptor> Chord_vector; typedef std::vector<halfedge_descriptor> Chord_vector;
typedef typename Chord_vector::iterator Chord_vector_iterator; typedef typename Chord_vector::iterator Chord_vector_iterator;
@ -231,10 +234,9 @@ private:
Compute_scalar_product_3 scalar_product_functor; Compute_scalar_product_3 scalar_product_functor;
// The facet proxy index map. // The facet proxy index map.
boost::unordered_map<face_descriptor, std::size_t> internal_fidx_map; Fproxy_map fproxy_map;
boost::associative_property_map<boost::unordered_map<face_descriptor, std::size_t> > fproxy_map;
// The attached anchor index of a vertex. // The attached anchor index of a vertex.
boost::unordered_map<vertex_descriptor, std::size_t> internal_vidx_map;
Vertex_anchor_map vanchor_map; Vertex_anchor_map vanchor_map;
// Proxies. // Proxies.
@ -261,8 +263,8 @@ public:
m_pmesh(NULL), m_pmesh(NULL),
fit_error(NULL), fit_error(NULL),
proxy_fitting(NULL), proxy_fitting(NULL),
fproxy_map(internal_fidx_map), fproxy_map(),
vanchor_map(internal_vidx_map), vanchor_map(),
average_edge_length(0.0) { average_edge_length(0.0) {
Geom_traits traits; Geom_traits traits;
vector_functor = traits.construct_vector_3_object(); vector_functor = traits.construct_vector_3_object();
@ -282,14 +284,30 @@ public:
point_pmap(_point_pmap), point_pmap(_point_pmap),
fit_error(NULL), fit_error(NULL),
proxy_fitting(NULL), proxy_fitting(NULL),
fproxy_map(internal_fidx_map), fproxy_map(),
vanchor_map(internal_vidx_map), vanchor_map(),
average_edge_length(0.0) { average_edge_length(0.0) {
Geom_traits traits; Geom_traits 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();
vanchor_map = CGAL::internal::add_property(Vertex_anchor_tag("VSA-vertex_anchor"),
*(const_cast<TriangleMesh*>(m_pmesh)));
fproxy_map = CGAL::internal::add_property(Face_proxy_tag("VSA-fproxy"),
*(const_cast<TriangleMesh*>(m_pmesh)));
}
~Mesh_approximation()
{
if(m_pmesh){
CGAL::internal::remove_property(vanchor_map, *(const_cast<TriangleMesh*>(m_pmesh)));
CGAL::internal::remove_property(fproxy_map, *(const_cast<TriangleMesh*>(m_pmesh)));
}
} }
/*! /*!
@ -334,10 +352,11 @@ public:
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
internal_fidx_map[f] = CGAL_VSA_INVALID_TAG; internal_fidx_map[f] = CGAL_VSA_INVALID_TAG;
internal_vidx_map.clear(); CGAL::internal::remove_property(vanchor_map, *m_pmesh);
vanchor_map = CGAL::internal::add_property(Vertex_anchor_tag("VSA-vertex_anchor"),
*(const_cast<TriangleMesh*>(p_mesh)));
BOOST_FOREACH(vertex_descriptor v, vertices(*m_pmesh)) BOOST_FOREACH(vertex_descriptor v, vertices(*m_pmesh))
internal_vidx_map.insert( put(vanchor_map, v, CGAL_VSA_INVALID_TAG);
std::pair<vertex_descriptor, std::size_t>(v, CGAL_VSA_INVALID_TAG));
// compute average edge length // compute average edge length
FT sum(0.0); FT sum(0.0);
@ -426,7 +445,7 @@ public:
for (std::size_t i = 0; i < nb_iterations; ++i) { for (std::size_t i = 0; i < nb_iterations; ++i) {
// tag the whole surface // tag the whole surface
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
fproxy_map[f] = CGAL_VSA_INVALID_TAG; put(fproxy_map,f, CGAL_VSA_INVALID_TAG);
partition(proxies.begin(), proxies.end()); partition(proxies.begin(), proxies.end());
fit(proxies.begin(), proxies.end()); fit(proxies.begin(), proxies.end());
@ -485,7 +504,7 @@ public:
BOOST_FOREACH(Proxy_wrapper &pxw, proxies) BOOST_FOREACH(Proxy_wrapper &pxw, proxies)
pxw.err = FT(0.0); pxw.err = FT(0.0);
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
std::size_t pxidx = fproxy_map[f]; std::size_t pxidx = get(fproxy_map,f);
proxies[pxidx].err += (*fit_error)(f, proxies[pxidx].px); proxies[pxidx].err += (*fit_error)(f, proxies[pxidx].px);
} }
@ -549,7 +568,7 @@ public:
static_cast<double>(num_faces(*m_pmesh)) / static_cast<double>(num_proxies)); static_cast<double>(num_faces(*m_pmesh)) / static_cast<double>(num_proxies));
std::vector<FT> px_size(proxies.size(), FT(0.0)); std::vector<FT> px_size(proxies.size(), FT(0.0));
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
px_size[fproxy_map[f]] += FT(1.0); px_size[get(fproxy_map,f)] += FT(1.0);
FT residual(0.0); FT residual(0.0);
for (std::size_t i = 0; i < proxies.size(); ++i) { for (std::size_t i = 0; i < proxies.size(); ++i) {
FT to_add = (residual + px_size[i]) / avg_facet; FT to_add = (residual + px_size[i]) / avg_facet;
@ -585,7 +604,7 @@ public:
std::size_t num_added = 0; std::size_t num_added = 0;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
const std::size_t px_id = fproxy_map[f]; const std::size_t px_id = get(fproxy_map, f);
if (proxies[px_id].seed == f) if (proxies[px_id].seed == f)
continue; continue;
@ -631,7 +650,7 @@ public:
bool found = false; bool found = false;
face_descriptor tele_to; face_descriptor tele_to;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
if (fproxy_map[f] == px_worst && f != proxies[px_worst].seed) { if (get(fproxy_map, f) == px_worst && f != proxies[px_worst].seed) {
// teleport to anywhere but the seed // teleport to anywhere but the seed
tele_to = f; tele_to = f;
found = true; found = true;
@ -652,16 +671,16 @@ public:
// update merged proxies // update merged proxies
std::list<face_descriptor> merged_patch; std::list<face_descriptor> merged_patch;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
std::size_t &px_idx = fproxy_map[f]; std::size_t px_idx = get(fproxy_map,f);
if (px_idx == px_enlarged || px_idx == px_merged) { if (px_idx == px_enlarged || px_idx == px_merged) {
px_idx = px_enlarged; put(fproxy_map, f, px_enlarged);
merged_patch.push_back(f); merged_patch.push_back(f);
} }
} }
proxies[px_enlarged] = fit_new_proxy(merged_patch.begin(), merged_patch.end(), px_enlarged); proxies[px_enlarged] = fit_new_proxy(merged_patch.begin(), merged_patch.end(), px_enlarged);
// replace the merged proxy position to the newly teleported proxy // replace the merged proxy position to the newly teleported proxy
proxies[px_merged] = fit_new_proxy(tele_to, px_merged); proxies[px_merged] = fit_new_proxy(tele_to, px_merged);
fproxy_map[tele_to] = px_merged; put(fproxy_map, tele_to, px_merged);
num_teleported++; num_teleported++;
// coarse re-fitting // coarse re-fitting
@ -692,10 +711,10 @@ public:
FT err_sum(0.0); FT err_sum(0.0);
std::list<face_descriptor> merged_patch; std::list<face_descriptor> merged_patch;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
std::size_t px_idx = fproxy_map[f]; std::size_t px_idx = get(fproxy_map, f);
if (px_idx == px1 || px_idx == px0) { if (px_idx == px1 || px_idx == px0) {
err_sum += (*fit_error)(f, proxies[px_idx].px); err_sum += (*fit_error)(f, proxies[px_idx].px);
fproxy_map[f] = px0; set(fproxy_map, f, px0);
merged_patch.push_back(f); merged_patch.push_back(f);
} }
} }
@ -707,8 +726,8 @@ public:
proxies[i].idx = i; proxies[i].idx = i;
// keep facet proxy map valid // keep facet proxy map valid
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
if (fproxy_map[f] > px1) if (get(fproxy_map,f) > px1)
--fproxy_map[f]; put(fproxy_map, f, get(fproxy_map,f)-1);
} }
FT err_merged(0.0); FT err_merged(0.0);
@ -734,7 +753,7 @@ public:
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(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
px_facets[fproxy_map[f]].push_back(f); px_facets[get(fproxy_map,f)].push_back(f);
// find best merge // find best merge
MergedPair merged_set; MergedPair merged_set;
@ -743,8 +762,8 @@ public:
BOOST_FOREACH(edge_descriptor e, edges(*m_pmesh)) { BOOST_FOREACH(edge_descriptor e, edges(*m_pmesh)) {
if (CGAL::is_border(e, *m_pmesh)) if (CGAL::is_border(e, *m_pmesh))
continue; continue;
std::size_t pxi = fproxy_map[face(halfedge(e, *m_pmesh), *m_pmesh)]; std::size_t pxi = get(fproxy_map, face(halfedge(e, *m_pmesh), *m_pmesh));
std::size_t pxj = fproxy_map[face(opposite(halfedge(e, *m_pmesh), *m_pmesh), *m_pmesh)]; std::size_t pxj = get(fproxy_map, face(opposite(halfedge(e, *m_pmesh), *m_pmesh), *m_pmesh));
if (pxi == pxj) if (pxi == pxj)
continue; continue;
if (pxi > pxj) if (pxi > pxj)
@ -804,7 +823,7 @@ public:
// collect confined proxy area // collect confined proxy area
std::vector<face_descriptor> confined_area; std::vector<face_descriptor> confined_area;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
if (fproxy_map[f] == px_idx) if (get(fproxy_map, f) == px_idx)
confined_area.push_back(f); confined_area.push_back(f);
// not enough facets to split // not enough facets to split
if (n > confined_area.size()) if (n > confined_area.size())
@ -820,8 +839,8 @@ public:
if (count >= n) if (count >= n)
break; break;
if (fproxy_map[f] == px_idx && f != proxies[px_idx].seed) { if (get(fproxy_map,f) == px_idx && f != proxies[px_idx].seed) {
fproxy_map[f] = proxies.size(); put(fproxy_map, f, proxies.size());
proxies.push_back(fit_new_proxy(f, proxies.size())); proxies.push_back(fit_new_proxy(f, proxies.size()));
++count; ++count;
// copy // copy
@ -832,7 +851,7 @@ public:
// relaxation on confined area and proxies // relaxation on confined area and proxies
for (std::size_t i = 0; i < nb_relaxations; ++i) { for (std::size_t i = 0; i < nb_relaxations; ++i) {
BOOST_FOREACH(face_descriptor f, confined_area) BOOST_FOREACH(face_descriptor f, confined_area)
fproxy_map[f] = CGAL_VSA_INVALID_TAG; put(fproxy_map, f, CGAL_VSA_INVALID_TAG);
partition(confined_proxies.begin(), confined_proxies.end()); partition(confined_proxies.begin(), confined_proxies.end());
fit(confined_proxies.begin(), confined_proxies.end()); fit(confined_proxies.begin(), confined_proxies.end());
@ -866,7 +885,7 @@ public:
const bool pca_plane = false) { const bool pca_plane = false) {
// initialize all vertex anchor status // initialize all vertex anchor status
BOOST_FOREACH(vertex_descriptor v, vertices(*m_pmesh)) BOOST_FOREACH(vertex_descriptor v, vertices(*m_pmesh))
internal_vidx_map[v] = CGAL_VSA_INVALID_TAG; put(vanchor_map, v, CGAL_VSA_INVALID_TAG);
anchors.clear(); anchors.clear();
borders.clear(); borders.clear();
tris.clear(); tris.clear();
@ -891,7 +910,7 @@ public:
template <typename FacetProxyMap> template <typename FacetProxyMap>
void get_proxy_map(FacetProxyMap &facet_proxy_map) const { void get_proxy_map(FacetProxyMap &facet_proxy_map) const {
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
facet_proxy_map[f] = fproxy_map[f]; facet_proxy_map[f] = get(fproxy_map, f);
} }
/*! /*!
@ -906,7 +925,7 @@ public:
return; return;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
if (fproxy_map[f] == px_idx) if (get(fproxy_map,f) == px_idx)
*out_itr++ = f; *out_itr++ = f;
} }
@ -988,7 +1007,7 @@ public:
do { do {
Chord_vector chord; Chord_vector chord;
walk_to_next_anchor(he, chord); walk_to_next_anchor(he, chord);
bdr.push_back(vanchor_map[target(he, *m_pmesh)]); bdr.push_back(get(vanchor_map, target(he, *m_pmesh)));
} while (he != he_mark); } while (he != he_mark);
*out_itr++ = bdr; *out_itr++ = bdr;
} }
@ -1173,11 +1192,11 @@ private:
std::priority_queue<Facet_to_integrate> facet_pqueue; std::priority_queue<Facet_to_integrate> facet_pqueue;
for (ProxyWrapperIterator pxw_itr = beg; pxw_itr != end; ++pxw_itr) { for (ProxyWrapperIterator pxw_itr = beg; pxw_itr != end; ++pxw_itr) {
face_descriptor f = pxw_itr->seed; face_descriptor f = pxw_itr->seed;
fproxy_map[f] = pxw_itr->idx; put(fproxy_map, f, pxw_itr->idx);
BOOST_FOREACH(face_descriptor fadj, faces_around_face(halfedge(f, *m_pmesh), *m_pmesh)) { BOOST_FOREACH(face_descriptor fadj, faces_around_face(halfedge(f, *m_pmesh), *m_pmesh)) {
if (fadj != boost::graph_traits<TriangleMesh>::null_face() if (fadj != boost::graph_traits<TriangleMesh>::null_face()
&& fproxy_map[fadj] == CGAL_VSA_INVALID_TAG) { && get(fproxy_map, fadj) == CGAL_VSA_INVALID_TAG) {
facet_pqueue.push(Facet_to_integrate( facet_pqueue.push(Facet_to_integrate(
fadj, pxw_itr->idx, (*fit_error)(fadj, pxw_itr->px))); fadj, pxw_itr->idx, (*fit_error)(fadj, pxw_itr->px)));
} }
@ -1187,11 +1206,11 @@ private:
while (!facet_pqueue.empty()) { while (!facet_pqueue.empty()) {
const Facet_to_integrate c = facet_pqueue.top(); const Facet_to_integrate c = facet_pqueue.top();
facet_pqueue.pop(); facet_pqueue.pop();
if (fproxy_map[c.f] == CGAL_VSA_INVALID_TAG) { if (get(fproxy_map,c.f) == CGAL_VSA_INVALID_TAG) {
fproxy_map[c.f] = c.px; put(fproxy_map, c.f, c.px);
BOOST_FOREACH(face_descriptor fadj, faces_around_face(halfedge(c.f, *m_pmesh), *m_pmesh)) { BOOST_FOREACH(face_descriptor fadj, faces_around_face(halfedge(c.f, *m_pmesh), *m_pmesh)) {
if (fadj != boost::graph_traits<TriangleMesh>::null_face() if (fadj != boost::graph_traits<TriangleMesh>::null_face()
&& fproxy_map[fadj] == CGAL_VSA_INVALID_TAG) { && get(fproxy_map,fadj) == CGAL_VSA_INVALID_TAG) {
facet_pqueue.push(Facet_to_integrate( facet_pqueue.push(Facet_to_integrate(
fadj, c.px, (*fit_error)(fadj, proxies[c.px].px))); fadj, c.px, (*fit_error)(fadj, proxies[c.px].px)));
} }
@ -1210,7 +1229,7 @@ private:
void fit(const ProxyWrapperIterator beg, const ProxyWrapperIterator end) { void fit(const ProxyWrapperIterator beg, const ProxyWrapperIterator end) {
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(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
px_facets[fproxy_map[f]].push_back(f); px_facets[get(fproxy_map,f)].push_back(f);
// update proxy parameters and seed // update proxy parameters and seed
for (ProxyWrapperIterator pxw_itr = beg; pxw_itr != end; ++pxw_itr) { for (ProxyWrapperIterator pxw_itr = beg; pxw_itr != end; ++pxw_itr) {
@ -1242,7 +1261,7 @@ private:
face_descriptor fworst; face_descriptor fworst;
bool first = true; bool first = true;
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) {
std::size_t px_idx = fproxy_map[f]; std::size_t px_idx = get(fproxy_map,f);
if (px_idx != px_worst || f == proxies[px_idx].seed) if (px_idx != px_worst || f == proxies[px_idx].seed)
continue; continue;
@ -1257,7 +1276,7 @@ private:
if (first) if (first)
return false; return false;
fproxy_map[fworst] = proxies.size(); put(fproxy_map,fworst, proxies.size());
proxies.push_back(fit_new_proxy(fworst, proxies.size())); proxies.push_back(fit_new_proxy(fworst, proxies.size()));
return true; return true;
@ -1345,7 +1364,7 @@ private:
proxies.clear(); proxies.clear();
proxies.push_back(fit_new_proxy(*(faces(*m_pmesh).first), proxies.size())); proxies.push_back(fit_new_proxy(*(faces(*m_pmesh).first), proxies.size()));
BOOST_FOREACH(face_descriptor f, faces(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
fproxy_map[f] = 0; put(fproxy_map,f, 0);
} }
/*! /*!
@ -1356,7 +1375,7 @@ private:
// fit proxy planes, areas, normals // fit proxy planes, areas, normals
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(*m_pmesh)) BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
px_facets[fproxy_map[f]].push_back(f); px_facets[get(fproxy_map,f)].push_back(f);
BOOST_FOREACH(const std::list<face_descriptor> &px_patch, px_facets) { BOOST_FOREACH(const std::list<face_descriptor> &px_patch, px_facets) {
Plane_3 fit_plane = if_pca_plane ? Plane_3 fit_plane = if_pca_plane ?
@ -1392,7 +1411,7 @@ private:
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(vtx, *m_pmesh)) { BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(vtx, *m_pmesh)) {
if (CGAL::is_border_edge(h, *m_pmesh)) if (CGAL::is_border_edge(h, *m_pmesh))
++border_count; ++border_count;
else if (fproxy_map[face(h, *m_pmesh)] != fproxy_map[face(opposite(h, *m_pmesh), *m_pmesh)]) else if (get(fproxy_map,face(h, *m_pmesh)) != get(fproxy_map, face(opposite(h, *m_pmesh), *m_pmesh)))
++border_count; ++border_count;
} }
if (border_count >= 3) if (border_count >= 3)
@ -1415,7 +1434,7 @@ private:
BOOST_FOREACH(halfedge_descriptor h, halfedges(*m_pmesh)) { BOOST_FOREACH(halfedge_descriptor h, halfedges(*m_pmesh)) {
if (!CGAL::is_border(h, *m_pmesh) if (!CGAL::is_border(h, *m_pmesh)
&& (CGAL::is_border(opposite(h, *m_pmesh), *m_pmesh) && (CGAL::is_border(opposite(h, *m_pmesh), *m_pmesh)
|| fproxy_map[face(h, *m_pmesh)] != fproxy_map[face(opposite(h, *m_pmesh), *m_pmesh)])) || get(fproxy_map, face(h, *m_pmesh)) != get(fproxy_map, face(opposite(h, *m_pmesh), *m_pmesh))))
he_candidates.insert(h); he_candidates.insert(h);
} }
@ -1542,8 +1561,8 @@ private:
EdgeWeightMap global_eweight_map = get(boost::edge_weight, gmain); EdgeWeightMap global_eweight_map = get(boost::edge_weight, gmain);
BOOST_FOREACH(vertex_descriptor v, vertices(*m_pmesh)) { BOOST_FOREACH(vertex_descriptor v, vertices(*m_pmesh)) {
sg_vertex_descriptor sgv = add_vertex(gmain); sg_vertex_descriptor sgv = add_vertex(gmain);
global_vanchor_map[sgv] = vanchor_map[v]; global_vanchor_map[sgv] = get(vanchor_map, v);
global_vtag_map[sgv] = vanchor_map[v]; global_vtag_map[sgv] = get(vanchor_map,v);
vmap.insert(std::pair<vertex_descriptor, sg_vertex_descriptor>(v, sgv)); vmap.insert(std::pair<vertex_descriptor, sg_vertex_descriptor>(v, sgv));
} }
BOOST_FOREACH(edge_descriptor e, edges(*m_pmesh)) { BOOST_FOREACH(edge_descriptor e, edges(*m_pmesh)) {
@ -1559,7 +1578,7 @@ private:
std::set<std::size_t> px_set; std::set<std::size_t> px_set;
BOOST_FOREACH(face_descriptor f, faces_around_target(halfedge(v, *m_pmesh), *m_pmesh)) { BOOST_FOREACH(face_descriptor f, faces_around_target(halfedge(v, *m_pmesh), *m_pmesh)) {
if (f != boost::graph_traits<TriangleMesh>::null_face()) if (f != boost::graph_traits<TriangleMesh>::null_face())
px_set.insert(fproxy_map[f]); px_set.insert(get(fproxy_map, f));
} }
BOOST_FOREACH(std::size_t p, px_set) BOOST_FOREACH(std::size_t p, px_set)
vertex_patches[p].push_back(to_sgv_map[v]); vertex_patches[p].push_back(to_sgv_map[v]);
@ -1624,8 +1643,8 @@ private:
} }
FT half_chord_len = vdist.back() / FT(2.0); FT half_chord_len = vdist.back() / FT(2.0);
const std::size_t anchorleft = vanchor_map[source(chord.front(), *m_pmesh)]; const std::size_t anchorleft = get(vanchor_map, source(chord.front(), *m_pmesh));
const std::size_t anchorright = vanchor_map[target(chord.back(), *m_pmesh)]; const std::size_t anchorright = get(vanchor_map, target(chord.back(), *m_pmesh));
typename std::vector<FT>::iterator ditr = vdist.begin() + 1; typename std::vector<FT>::iterator ditr = vdist.begin() + 1;
for (typename Chord_vector::iterator hitr = chord.begin(); for (typename Chord_vector::iterator hitr = chord.begin();
hitr != chord.end() - 1; ++hitr, ++ditr) { hitr != chord.end() - 1; ++hitr, ++ditr) {
@ -1684,9 +1703,9 @@ private:
* @param[in/out] he_start region border halfedge * @param[in/out] he_start region border halfedge
*/ */
void walk_to_next_border_halfedge(halfedge_descriptor &he_start) const { void walk_to_next_border_halfedge(halfedge_descriptor &he_start) const {
const std::size_t px_idx = fproxy_map[face(he_start, *m_pmesh)]; const std::size_t px_idx = get(fproxy_map, face(he_start, *m_pmesh));
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(he_start, *m_pmesh)) { BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(he_start, *m_pmesh)) {
if (CGAL::is_border(h, *m_pmesh) || fproxy_map[face(h, *m_pmesh)] != px_idx) { if (CGAL::is_border(h, *m_pmesh) || get(fproxy_map, face(h, *m_pmesh)) != px_idx) {
he_start = opposite(h, *m_pmesh); he_start = opposite(h, *m_pmesh);
return; return;
} }
@ -1712,8 +1731,8 @@ private:
const std::size_t chord_size = std::distance(chord_begin, chord_end); const std::size_t chord_size = std::distance(chord_begin, chord_end);
const halfedge_descriptor he_first = *chord_begin; const halfedge_descriptor he_first = *chord_begin;
const halfedge_descriptor he_last = *(chord_end - 1); const halfedge_descriptor he_last = *(chord_end - 1);
const std::size_t anchor_first = vanchor_map[source(he_first, *m_pmesh)]; const std::size_t anchor_first = get(vanchor_map, source(he_first, *m_pmesh));
const std::size_t anchor_last = vanchor_map[target(he_last, *m_pmesh)]; const std::size_t anchor_last = get(vanchor_map, target(he_last, *m_pmesh));
// do not subdivide trivial non-circular chord // do not subdivide trivial non-circular chord
if ((anchor_first != anchor_last) && (chord_size < 4)) if ((anchor_first != anchor_last) && (chord_size < 4))
@ -1763,10 +1782,10 @@ private:
if (with_dihedral_angle) { if (with_dihedral_angle) {
// suppose the proxy normal angle is acute // suppose the proxy normal angle is acute
std::size_t px_left = fproxy_map[face(he_first, *m_pmesh)]; std::size_t px_left = get(fproxy_map, face(he_first, *m_pmesh));
std::size_t px_right = px_left; std::size_t px_right = px_left;
if (!CGAL::is_border(opposite(he_first, *m_pmesh), *m_pmesh)) if (!CGAL::is_border(opposite(he_first, *m_pmesh), *m_pmesh))
px_right = fproxy_map[face(opposite(he_first, *m_pmesh), *m_pmesh)]; px_right = get(fproxy_map, face(opposite(he_first, *m_pmesh), *m_pmesh));
FT norm_sin(1.0); FT norm_sin(1.0);
if (!CGAL::is_border(opposite(he_first, *m_pmesh), *m_pmesh)) { if (!CGAL::is_border(opposite(he_first, *m_pmesh), *m_pmesh)) {
Vector_3 vec = CGAL::cross_product( Vector_3 vec = CGAL::cross_product(
@ -1813,7 +1832,7 @@ private:
bool is_anchor_attached( bool is_anchor_attached(
const typename boost::property_traits<VertexAnchorIndexMap>::key_type &v, const typename boost::property_traits<VertexAnchorIndexMap>::key_type &v,
const VertexAnchorIndexMap &vertex_anchor_map) const { const VertexAnchorIndexMap &vertex_anchor_map) const {
return vertex_anchor_map[v] != CGAL_VSA_INVALID_TAG; return get(vertex_anchor_map, v) != CGAL_VSA_INVALID_TAG;
} }
/*! /*!
@ -1821,7 +1840,7 @@ private:
* @param vtx vertex * @param vtx vertex
*/ */
void attach_anchor(const vertex_descriptor &vtx) { void attach_anchor(const vertex_descriptor &vtx) {
vanchor_map[vtx] = anchors.size(); put(vanchor_map, vtx, anchors.size());
anchors.push_back(Anchor(vtx, compute_anchor_position(vtx))); anchors.push_back(Anchor(vtx, compute_anchor_position(vtx)));
} }
@ -1844,7 +1863,7 @@ private:
std::set<std::size_t> px_set; std::set<std::size_t> px_set;
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(v, *m_pmesh)) { BOOST_FOREACH(halfedge_descriptor h, halfedges_around_target(v, *m_pmesh)) {
if (!CGAL::is_border(h, *m_pmesh)) if (!CGAL::is_border(h, *m_pmesh))
px_set.insert(fproxy_map[face(h, *m_pmesh)]); px_set.insert(get(fproxy_map, face(h, *m_pmesh)));
} }
// construct an anchor from vertex and the incident proxies // construct an anchor from vertex and the incident proxies

146
Surface_mesh_approximation/include/CGAL/vsa_metrics.h
View File

@ -72,8 +72,10 @@ class L21_metric
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor; typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, Vector_3> > FacetNormalMap; typedef CGAL::internal::face_property_t<Vector_3> Face_normal_tag;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, FT> > FacetAreaMap; typedef CGAL::internal::face_property_t<FT> Face_area_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_normal_tag >::type FacetNormalMap;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_area_tag >::type FacetAreaMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -81,12 +83,16 @@ public:
// constructor // constructor
L21_metric(const TriangleMesh &tm, const VertexPointMap &point_pmap) L21_metric(const TriangleMesh &tm, const VertexPointMap &point_pmap)
: normal_pmap(facet_normals), area_pmap(facet_areas) { : normal_pmap(), area_pmap() {
GeomTraits traits; GeomTraits traits;
scalar_product_functor = traits.compute_scalar_product_3_object(); scalar_product_functor = traits.compute_scalar_product_3_object();
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
normal_pmap = CGAL::internal::add_property(Face_normal_tag("VSA-face_normal"),
const_cast<TriangleMesh&>(tm));
area_pmap = CGAL::internal::add_property(Face_area_tag("VSA-face_area"),
const_cast<TriangleMesh&>(tm));
// construct internal facet normal & area map // construct internal facet normal & area map
BOOST_FOREACH(face_descriptor f, faces(tm)) { BOOST_FOREACH(face_descriptor f, faces(tm)) {
const halfedge_descriptor he = halfedge(f, tm); const halfedge_descriptor he = halfedge(f, tm);
@ -94,23 +100,21 @@ public:
const Point_3 &p1 = point_pmap[target(he, tm)]; const Point_3 &p1 = point_pmap[target(he, tm)];
const Point_3 &p2 = point_pmap[target(next(he, tm), tm)]; const Point_3 &p2 = point_pmap[target(next(he, tm), tm)];
Vector_3 normal = CGAL::unit_normal(p0, p1, p2); Vector_3 normal = CGAL::unit_normal(p0, p1, p2);
facet_normals.insert(std::pair<face_descriptor, Vector_3>(f, normal)); put(normal_pmap, f, normal);
FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2)))); FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2))));
facet_areas.insert(std::pair<face_descriptor, FT>(f, area)); put(area_pmap, f, area);
} }
} }
// returns L21 error of a facet f to a proxy px. // returns L21 error of a facet f to a proxy px.
FT operator()(const face_descriptor &f, const Proxy &px) const { FT operator()(const face_descriptor &f, const Proxy &px) const {
Vector_3 v = sum_functor(normal_pmap[f], scale_functor(px.normal, FT(-1))); Vector_3 v = sum_functor(get(normal_pmap, f), scale_functor(px.normal, FT(-1)));
return area_pmap[f] * scalar_product_functor(v, v); return get(area_pmap, f) * scalar_product_functor(v, v);
} }
private: private:
boost::unordered_map<face_descriptor, Vector_3> facet_normals; FacetNormalMap normal_pmap;
boost::unordered_map<face_descriptor, FT> facet_areas; FacetAreaMap area_pmap;
const FacetNormalMap normal_pmap;
const FacetAreaMap area_pmap;
Construct_scaled_vector_3 scale_functor; Construct_scaled_vector_3 scale_functor;
Compute_scalar_product_3 scalar_product_functor; Compute_scalar_product_3 scalar_product_functor;
Construct_sum_of_vectors_3 sum_functor; Construct_sum_of_vectors_3 sum_functor;
@ -138,8 +142,11 @@ class L21_metric<TriangleMesh,
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap; typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, Vector_3> > FacetNormalMap;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, FT> > FacetAreaMap; typedef CGAL::internal::face_property_t<Vector_3> Face_normal_tag;
typedef CGAL::internal::face_property_t<FT> Face_area_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_normal_tag >::type FacetNormalMap;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_area_tag >::type FacetAreaMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -147,12 +154,17 @@ public:
// constructor // constructor
L21_metric(const TriangleMesh &tm) L21_metric(const TriangleMesh &tm)
: normal_pmap(facet_normals), area_pmap(facet_areas) { : normal_pmap(), area_pmap() {
GeomTraits traits; GeomTraits traits;
scalar_product_functor = traits.compute_scalar_product_3_object(); scalar_product_functor = traits.compute_scalar_product_3_object();
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
normal_pmap = CGAL::internal::add_property(Face_normal_tag("VSA-face_normal"),
const_cast<TriangleMesh&>(tm));
area_pmap = CGAL::internal::add_property(Face_area_tag("VSA-face_area"),
const_cast<TriangleMesh&>(tm));
// construct internal facet normal & area map // construct internal facet normal & area map
VertexPointMap point_pmap = get(boost::vertex_point, const_cast<TriangleMesh &>(tm)); VertexPointMap point_pmap = get(boost::vertex_point, const_cast<TriangleMesh &>(tm));
BOOST_FOREACH(face_descriptor f, faces(tm)) { BOOST_FOREACH(face_descriptor f, faces(tm)) {
@ -161,23 +173,21 @@ public:
const Point_3 &p1 = point_pmap[target(he, tm)]; const Point_3 &p1 = point_pmap[target(he, tm)];
const Point_3 &p2 = point_pmap[target(next(he, tm), tm)]; const Point_3 &p2 = point_pmap[target(next(he, tm), tm)];
Vector_3 normal = CGAL::unit_normal(p0, p1, p2); Vector_3 normal = CGAL::unit_normal(p0, p1, p2);
facet_normals.insert(std::pair<face_descriptor, Vector_3>(f, normal)); put(normal_pmap, f, normal);
FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2)))); FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2))));
facet_areas.insert(std::pair<face_descriptor, FT>(f, area)); put(area_pmap, f, area);
} }
} }
// returns L21 error of a facet f to a proxy px. // returns L21 error of a facet f to a proxy px.
FT operator()(const face_descriptor &f, const Proxy &px) const { FT operator()(const face_descriptor &f, const Proxy &px) const {
Vector_3 v = sum_functor(normal_pmap[f], scale_functor(px.normal, FT(-1))); Vector_3 v = sum_functor(get(normal_pmap, f), scale_functor(px.normal, FT(-1)));
return area_pmap[f] * scalar_product_functor(v, v); return get(area_pmap, f) * scalar_product_functor(v, v);
} }
private: private:
boost::unordered_map<face_descriptor, Vector_3> facet_normals; FacetNormalMap normal_pmap;
boost::unordered_map<face_descriptor, FT> facet_areas; FacetAreaMap area_pmap;
const FacetNormalMap normal_pmap;
const FacetAreaMap area_pmap;
Construct_scaled_vector_3 scale_functor; Construct_scaled_vector_3 scale_functor;
Compute_scalar_product_3 scalar_product_functor; Compute_scalar_product_3 scalar_product_functor;
Construct_sum_of_vectors_3 sum_functor; Construct_sum_of_vectors_3 sum_functor;
@ -204,7 +214,8 @@ class L21_metric<TriangleMesh,
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor; typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, Vector_3> > FacetNormalMap; typedef CGAL::internal::face_property_t<Vector_3> Face_normal_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_normal_tag >::type FacetNormalMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -212,12 +223,15 @@ public:
// constructor // constructor
L21_metric(const TriangleMesh &tm, const VertexPointMap &point_pmap) L21_metric(const TriangleMesh &tm, const VertexPointMap &point_pmap)
: normal_pmap(facet_normals) { : normal_pmap() {
GeomTraits traits; GeomTraits traits;
scalar_product_functor = traits.compute_scalar_product_3_object(); scalar_product_functor = traits.compute_scalar_product_3_object();
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
normal_pmap = CGAL::internal::add_property(Face_normal_tag("VSA-face_normal"),
const_cast<TriangleMesh&>(tm));
// construct internal facet normal map // construct internal facet normal map
BOOST_FOREACH(face_descriptor f, faces(tm)) { BOOST_FOREACH(face_descriptor f, faces(tm)) {
const halfedge_descriptor he = halfedge(f, tm); const halfedge_descriptor he = halfedge(f, tm);
@ -225,20 +239,19 @@ public:
const Point_3 &p1 = point_pmap[target(he, tm)]; const Point_3 &p1 = point_pmap[target(he, tm)];
const Point_3 &p2 = point_pmap[target(next(he, tm), tm)]; const Point_3 &p2 = point_pmap[target(next(he, tm), tm)];
Vector_3 normal = CGAL::unit_normal(p0, p1, p2); Vector_3 normal = CGAL::unit_normal(p0, p1, p2);
facet_normals.insert(std::pair<face_descriptor, Vector_3>(f, normal)); put(normal_pmap, f, normal);
FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2)))); FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2))));
} }
} }
// returns L21 error of a facet f to a proxy px. // returns L21 error of a facet f to a proxy px.
FT operator()(const face_descriptor &f, const Proxy &px) const { FT operator()(const face_descriptor &f, const Proxy &px) const {
Vector_3 v = sum_functor(normal_pmap[f], scale_functor(px.normal, FT(-1))); Vector_3 v = sum_functor(get(normal_pmap,f), scale_functor(px.normal, FT(-1)));
return scalar_product_functor(v, v); return scalar_product_functor(v, v);
} }
private: private:
boost::unordered_map<face_descriptor, Vector_3> facet_normals; FacetNormalMap normal_pmap;
const FacetNormalMap normal_pmap;
Construct_scaled_vector_3 scale_functor; Construct_scaled_vector_3 scale_functor;
Compute_scalar_product_3 scalar_product_functor; Compute_scalar_product_3 scalar_product_functor;
Construct_sum_of_vectors_3 sum_functor; Construct_sum_of_vectors_3 sum_functor;
@ -265,7 +278,9 @@ class L21_metric<TriangleMesh,
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap; typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, Vector_3> > FacetNormalMap;
typedef CGAL::internal::face_property_t<Vector_3> Face_normal_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_normal_tag >::type FacetNormalMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -273,12 +288,15 @@ public:
// constructor // constructor
L21_metric(const TriangleMesh &tm) L21_metric(const TriangleMesh &tm)
: normal_pmap(facet_normals) { : normal_pmap() {
GeomTraits traits; GeomTraits traits;
scalar_product_functor = traits.compute_scalar_product_3_object(); scalar_product_functor = traits.compute_scalar_product_3_object();
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
normal_pmap = CGAL::internal::add_property(Face_normal_tag("VSA-face_normal"),
const_cast<TriangleMesh&>(tm));
// construct internal facet normal map // construct internal facet normal map
VertexPointMap point_pmap = get(boost::vertex_point, const_cast<TriangleMesh &>(tm)); VertexPointMap point_pmap = get(boost::vertex_point, const_cast<TriangleMesh &>(tm));
BOOST_FOREACH(face_descriptor f, faces(tm)) { BOOST_FOREACH(face_descriptor f, faces(tm)) {
@ -287,19 +305,18 @@ public:
const Point_3 &p1 = point_pmap[target(he, tm)]; const Point_3 &p1 = point_pmap[target(he, tm)];
const Point_3 &p2 = point_pmap[target(next(he, tm), tm)]; const Point_3 &p2 = point_pmap[target(next(he, tm), tm)];
Vector_3 normal = CGAL::unit_normal(p0, p1, p2); Vector_3 normal = CGAL::unit_normal(p0, p1, p2);
facet_normals.insert(std::pair<face_descriptor, Vector_3>(f, normal)); put(normal_pmap, f, normal);
} }
} }
// returns L21 error of a facet f to a proxy px. // returns L21 error of a facet f to a proxy px.
FT operator()(const face_descriptor &f, const Proxy &px) const { FT operator()(const face_descriptor &f, const Proxy &px) const {
Vector_3 v = sum_functor(normal_pmap[f], scale_functor(px.normal, FT(-1))); Vector_3 v = sum_functor(get(normal_pmap,f), scale_functor(px.normal, FT(-1)));
return scalar_product_functor(v, v); return scalar_product_functor(v, v);
} }
private: private:
boost::unordered_map<face_descriptor, Vector_3> facet_normals; FacetNormalMap normal_pmap;
const FacetNormalMap normal_pmap;
Construct_scaled_vector_3 scale_functor; Construct_scaled_vector_3 scale_functor;
Compute_scalar_product_3 scalar_product_functor; Compute_scalar_product_3 scalar_product_functor;
Construct_sum_of_vectors_3 sum_functor; Construct_sum_of_vectors_3 sum_functor;
@ -335,8 +352,10 @@ class L21_proxy_fitting
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor; typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, Vector_3> > FacetNormalMap; typedef CGAL::internal::face_property_t<Vector_3> Face_normal_tag;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, FT> > FacetAreaMap; typedef CGAL::internal::face_property_t<FT> Face_area_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_normal_tag >::type FacetNormalMap;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_area_tag >::type FacetAreaMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -344,11 +363,16 @@ public:
// constructor. // constructor.
L21_proxy_fitting(const TriangleMesh &tm, const VertexPointMap &point_pmap) L21_proxy_fitting(const TriangleMesh &tm, const VertexPointMap &point_pmap)
: normal_pmap(facet_normals), area_pmap(facet_areas) { : normal_pmap(), area_pmap() {
GeomTraits traits; GeomTraits traits;
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
normal_pmap = CGAL::internal::add_property(Face_normal_tag("VSA-face_normal"),
const_cast<TriangleMesh&>(tm));
area_pmap = CGAL::internal::add_property(Face_area_tag("VSA-face_area"),
const_cast<TriangleMesh&>(tm));
// construct internal facet normal & area map // construct internal facet normal & area map
BOOST_FOREACH(face_descriptor f, faces(tm)) { BOOST_FOREACH(face_descriptor f, faces(tm)) {
const halfedge_descriptor he = halfedge(f, tm); const halfedge_descriptor he = halfedge(f, tm);
@ -356,9 +380,9 @@ public:
const Point_3 &p1 = point_pmap[target(he, tm)]; const Point_3 &p1 = point_pmap[target(he, tm)];
const Point_3 &p2 = point_pmap[target(next(he, tm), tm)]; const Point_3 &p2 = point_pmap[target(next(he, tm), tm)];
Vector_3 normal = CGAL::unit_normal(p0, p1, p2); Vector_3 normal = CGAL::unit_normal(p0, p1, p2);
facet_normals.insert(std::pair<face_descriptor, Vector_3>(f, normal)); put(normal_pmap, f, normal);
FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2)))); FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2))));
facet_areas.insert(std::pair<face_descriptor, FT>(f, area)); put(area_pmap, f, area);
} }
} }
@ -371,7 +395,7 @@ public:
Vector_3 norm = CGAL::NULL_VECTOR; Vector_3 norm = CGAL::NULL_VECTOR;
for (FacetIterator fitr = beg; fitr != end; ++fitr) { for (FacetIterator fitr = beg; fitr != end; ++fitr) {
norm = sum_functor(norm, norm = sum_functor(norm,
scale_functor(normal_pmap[*fitr], area_pmap[*fitr])); scale_functor(get(normal_pmap, *fitr), get(area_pmap, *fitr)));
} }
norm = scale_functor(norm, norm = scale_functor(norm,
FT(1.0 / std::sqrt(CGAL::to_double(norm.squared_length())))); FT(1.0 / std::sqrt(CGAL::to_double(norm.squared_length()))));
@ -384,10 +408,8 @@ public:
} }
private: private:
boost::unordered_map<face_descriptor, Vector_3> facet_normals; FacetNormalMap normal_pmap;
boost::unordered_map<face_descriptor, FT> facet_areas; FacetAreaMap area_pmap;
const FacetNormalMap normal_pmap;
const FacetAreaMap area_pmap;
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;
}; };
@ -411,8 +433,11 @@ class L21_proxy_fitting<TriangleMesh,
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap; typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, Vector_3> > FacetNormalMap;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, FT> > FacetAreaMap; typedef CGAL::internal::face_property_t<Vector_3> Face_normal_tag;
typedef CGAL::internal::face_property_t<FT> Face_area_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_normal_tag >::type FacetNormalMap;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_area_tag >::type FacetAreaMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -420,11 +445,15 @@ public:
// constructor. // constructor.
L21_proxy_fitting(const TriangleMesh &tm) L21_proxy_fitting(const TriangleMesh &tm)
: normal_pmap(facet_normals), area_pmap(facet_areas) { : normal_pmap(), area_pmap() {
GeomTraits traits; GeomTraits traits;
sum_functor = traits.construct_sum_of_vectors_3_object(); sum_functor = traits.construct_sum_of_vectors_3_object();
scale_functor = traits.construct_scaled_vector_3_object(); scale_functor = traits.construct_scaled_vector_3_object();
normal_pmap = CGAL::internal::add_property(Face_normal_tag("VSA-face_normal"),
const_cast<TriangleMesh&>(tm));
area_pmap = CGAL::internal::add_property(Face_area_tag("VSA-face_area"),
const_cast<TriangleMesh&>(tm));
// construct internal facet normal & area map // construct internal facet normal & area map
VertexPointMap point_pmap = get(boost::vertex_point, const_cast<TriangleMesh &>(tm)); VertexPointMap point_pmap = get(boost::vertex_point, const_cast<TriangleMesh &>(tm));
BOOST_FOREACH(face_descriptor f, faces(tm)) { BOOST_FOREACH(face_descriptor f, faces(tm)) {
@ -433,9 +462,9 @@ public:
const Point_3 &p1 = point_pmap[target(he, tm)]; const Point_3 &p1 = point_pmap[target(he, tm)];
const Point_3 &p2 = point_pmap[target(next(he, tm), tm)]; const Point_3 &p2 = point_pmap[target(next(he, tm), tm)];
Vector_3 normal = CGAL::unit_normal(p0, p1, p2); Vector_3 normal = CGAL::unit_normal(p0, p1, p2);
facet_normals.insert(std::pair<face_descriptor, Vector_3>(f, normal)); put(normal_pmap, f, normal);
FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2)))); FT area(std::sqrt(CGAL::to_double(CGAL::squared_area(p0, p1, p2))));
facet_areas.insert(std::pair<face_descriptor, FT>(f, area)); put(area_pmap, f, area);
} }
} }
@ -448,7 +477,7 @@ public:
Vector_3 norm = CGAL::NULL_VECTOR; Vector_3 norm = CGAL::NULL_VECTOR;
for (FacetIterator fitr = beg; fitr != end; ++fitr) { for (FacetIterator fitr = beg; fitr != end; ++fitr) {
norm = sum_functor(norm, norm = sum_functor(norm,
scale_functor(normal_pmap[*fitr], area_pmap[*fitr])); scale_functor(get(normal_pmap, *fitr), get(area_pmap, *fitr)));
} }
norm = scale_functor(norm, norm = scale_functor(norm,
FT(1.0 / std::sqrt(CGAL::to_double(norm.squared_length())))); FT(1.0 / std::sqrt(CGAL::to_double(norm.squared_length()))));
@ -461,10 +490,8 @@ public:
} }
private: private:
boost::unordered_map<face_descriptor, Vector_3> facet_normals; FacetNormalMap normal_pmap;
boost::unordered_map<face_descriptor, FT> facet_areas; FacetAreaMap area_pmap;
const FacetNormalMap normal_pmap;
const FacetAreaMap area_pmap;
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;
}; };
@ -496,7 +523,8 @@ class L2_metric
typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor; typedef typename boost::graph_traits<TriangleMesh>::face_descriptor face_descriptor;
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, FT> > FacetAreaMap; typedef CGAL::internal::face_property_t<FT> Face_area_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_area_tag >::type FacetAreaMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -528,7 +556,7 @@ public:
const FT d1(std::sqrt(CGAL::to_double(sq_d1))); const FT d1(std::sqrt(CGAL::to_double(sq_d1)));
const FT d2(std::sqrt(CGAL::to_double(sq_d2))); const FT d2(std::sqrt(CGAL::to_double(sq_d2)));
return (sq_d0 + sq_d1 + sq_d2 + d0 * d1 + d1 * d2 + d2 * d0) * area_pmap[f] / FT(6.0); return (sq_d0 + sq_d1 + sq_d2 + d0 * d1 + d1 * d2 + d2 * d0) * get(area_pmap, f) / FT(6.0);
} }
private: private:
@ -554,7 +582,9 @@ class L2_metric<TriangleMesh,
typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor; typedef typename boost::graph_traits<TriangleMesh>::halfedge_descriptor halfedge_descriptor;
typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap; typedef typename boost::property_map<TriangleMesh, boost::vertex_point_t>::type VertexPointMap;
typedef boost::associative_property_map<boost::unordered_map<face_descriptor, FT> > FacetAreaMap;
typedef CGAL::internal::face_property_t<FT> Face_area_tag;
typedef typename CGAL::internal::dynamic_property_map<TriangleMesh, Face_area_tag >::type FacetAreaMap;
public: public:
// type required by `ErrorMetric` concept // type required by `ErrorMetric` concept
@ -587,7 +617,7 @@ public:
const FT d1(std::sqrt(CGAL::to_double(sq_d1))); const FT d1(std::sqrt(CGAL::to_double(sq_d1)));
const FT d2(std::sqrt(CGAL::to_double(sq_d2))); const FT d2(std::sqrt(CGAL::to_double(sq_d2)));
return (sq_d0 + sq_d1 + sq_d2 + d0 * d1 + d1 * d2 + d2 * d0) * area_pmap[f] / FT(6); return (sq_d0 + sq_d1 + sq_d2 + d0 * d1 + d1 * d2 + d2 * d0) * get(area_pmap, f) / FT(6);
} }
private: private: