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

570 lines
20 KiB
C++
Raw Blame History

// ======================================================================
//
// Copyright (c) 2005-2006 by Le-Jeng Shiue. All Rights Reserved.
//
// This file is part of CGAL (www.cgal.org); you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public License as
// published by the Free Software Foundation; version 2.1 of the License.
// See the file LICENSE.LGPL distributed with CGAL.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL$
// $Id$
//
//
// Author(s): Le-Jeng Shiue <Andy.Shiue@gmail.com>
// <1998-2005 SurfLab, CISE, University of Florida>
// <2005- RapidMind Inc.>
//
// ======================================================================
#ifndef CGAL_POLYHEDRON_SUBDIVISION_IMPL_H_02102006
#define CGAL_POLYHEDRON_SUBDIVISION_IMPL_H_02102006
#include <CGAL/basic.h>
#include <vector>
#include <CGAL/circulator.h>
#include <CGAL/Polyhedron_decorator_3.h>
CGAL_BEGIN_NAMESPACE
// ======================================================================
namespace Subdivision_method_3 {
namespace Private {
// ======================================================================
template <class Poly, template <typename> class Mask>
void PQQ_1step(Poly& p, Mask<Poly> mask) {
typedef Polyhedron_decorator_3<Poly> PD;
typedef typename Poly::Halfedge_data_structure HDS;
typedef typename Poly::Vertex Vertex;
typedef typename Poly::Halfedge Halfedge;
typedef typename Poly::Facet Facet;
typedef typename Poly::Vertex_handle Vertex_handle;
typedef typename Poly::Halfedge_handle Halfedge_handle;
typedef typename Poly::Facet_handle Facet_handle;
typedef typename Poly::Vertex_iterator Vertex_iterator;
typedef typename Poly::Halfedge_iterator Halfedge_iterator;
typedef typename Poly::Edge_iterator Edge_iterator;
typedef typename Poly::Facet_iterator Facet_iterator;
typedef typename Poly::Halfedge_around_facet_circulator
Halfedge_around_facet_circulator;
typedef typename Poly::Halfedge_around_vertex_circulator
Halfedge_around_vertex_circulator;
typedef typename Poly::Traits Traits;
typedef typename Traits::Kernel Kernel;
typedef typename Kernel::Point_3 Point;
p.normalize_border();
// Build a new vertices buffer has the following structure
//
// 0 1 ... e_begin ... f_begin ... (end_of_buffer)
// 0 ... e_begin-1 : store the positions of the vertex-vertices
// e_begin ... f_begin-1 : store the positions of the edge-vertices
// f_begin ... (end) : store the positions of the face-vertices
// The index of the vertices buffer should 1-1 map to the distance
// of the corresponding iterator to the begin of the iterator.
size_t num_vertex = p.size_of_vertices();
size_t num_edge = p.size_of_halfedges()/2;
size_t num_facet = p.size_of_facets();
// If Polyhedron is using vector, we need to reserve the memory to prevent
// the CGAL_assertion.
// This function for polyhedron using list is VOID.
p.reserve(num_vertex+num_edge+num_facet, 4*2*num_edge, 4*num_edge/2);
Point* vertex_point_buffer = new Point[num_vertex + num_edge + num_facet];
Point* edge_point_buffer = vertex_point_buffer + num_vertex;
Point* face_point_buffer = edge_point_buffer + num_edge;
std::vector<bool> v_onborder(num_vertex);
Facet_iterator fitr = p.facets_begin();
for (size_t i = 0; i < num_facet; i++, ++fitr)
mask.facet_node(fitr, face_point_buffer[i]);
size_t sb = p.size_of_border_edges();
Edge_iterator eitr = p.edges_begin();
for (size_t i = 0; i < num_edge-sb; i++, ++eitr)
mask.edge_node(eitr, edge_point_buffer[i]);
for (size_t i = num_edge-sb; i < num_edge; i++, ++eitr) {
int v = (int) std::distance(p.vertices_begin(), eitr->vertex());
v_onborder[v] = true;
mask.border_node(eitr, edge_point_buffer[i], vertex_point_buffer[v]);
}
Vertex_iterator vitr = p.vertices_begin();
for (size_t i = 0; i < num_vertex; i++, ++vitr)
if (!v_onborder[i]) mask.vertex_node(vitr, vertex_point_buffer[i]);
// Build the connectivity using insert_vertex() and insert_edge()
// 1. insert_vertex() to all edges and set them to new positions
// 2. insert_edge() between 2 randomly selected neighboring new inserted
// vertices
// 3. insert_vertex() to the new inserted edge and set them to new positions
// 4. insert_edge() between all other new inserted vertices of step 1 and
// the new inserted vertex of step 3
// Step 1.
eitr = p.edges_begin();
for (size_t i = 0; i < num_edge; i++, ++eitr) {
Vertex_handle vh = PD::insert_vertex(p, eitr);
vh->point() = edge_point_buffer[i];
}
fitr = p.facets_begin();
// TODO: the topoloy modification can be done by a template function
// and that gives the user a chance to create new topological masks.
for (size_t i = 0; i < num_facet; i++, ++fitr) {
// Step 2.
Halfedge_around_facet_circulator hcir_begin = fitr->facet_begin();
Halfedge_around_facet_circulator hcir = hcir_begin;
Halfedge_handle e1 = ++hcir; // e1 points to the newly inserted vertex
++hcir; // Skips one original vertex
Halfedge_handle e2 = ++hcir; // points to the next newly inserted vertex
++hcir; // Must move the cir before inserts the new edge !!
Halfedge_handle newe = PD::insert_edge(p, e1, e2);
// Step 3.
Halfedge_handle newv = PD::insert_vertex_return_edge(p, newe);
newv = newv->opposite()->prev(); // change newv to the larger face and
// still points to the newly inserted
// vertex
// Update the geometry data of the newly inserted face-vertices
newv->vertex()->point() = face_point_buffer[i];
// Step 4.
while (hcir != hcir_begin) {
e1 = ++hcir;
++hcir; // Must move the cir before inserts the new edge !!
PD::insert_edge(p, e1, newv);
}
}
// Update the geometry data of the newly inserted vertices by the
// vertices buffer
vitr = p.vertices_begin();
for (size_t i = 0; i < num_vertex; i++, ++vitr)
vitr->point() = vertex_point_buffer[i];
delete []vertex_point_buffer;
}
// ======================================================================
template <class Poly, template <typename> class Mask>
void PTQ_1step(Poly& p, Mask<Poly> mask) {
typedef Polyhedron_decorator_3<Poly> PD;
typedef typename Poly::Halfedge_data_structure HDS;
typedef typename Poly::Vertex Vertex;
typedef typename Poly::Halfedge Halfedge;
typedef typename Poly::Facet Facet;
typedef typename Poly::Vertex_handle Vertex_handle;
typedef typename Poly::Halfedge_handle Halfedge_handle;
typedef typename Poly::Facet_handle Facet_handle;
typedef typename Poly::Vertex_iterator Vertex_iterator;
typedef typename Poly::Halfedge_iterator Halfedge_iterator;
typedef typename Poly::Edge_iterator Edge_iterator;
typedef typename Poly::Facet_iterator Facet_iterator;
typedef typename Poly::Halfedge_around_facet_circulator
Halfedge_around_facet_circulator;
typedef typename Poly::Halfedge_around_vertex_circulator
Halfedge_around_vertex_circulator;
typedef typename Poly::Traits Traits;
typedef typename Traits::Kernel Kernel;
typedef typename Kernel::Point_3 Point;
p.normalize_border();
// Build a new vertices buffer has the following structure
//
// 0 1 ... e_begin ... f_begin ... (end_of_buffer)
// 0 ... e_begin-1 : store the positions of the vertex-vertices
// e_begin ... (end) : store the positions of the edge-vertices
// The index of the vertices buffer should 1-1 map to the distance
// of the corresponding iterator to the begin of the iterator.
size_t num_vertex = p.size_of_vertices();
size_t num_edge = p.size_of_halfedges()/2;
size_t num_facet = p.size_of_facets();
// If Polyhedron is using vector, we need to reserve the memory to prevent
// the CGAL_assertion.
// This function for polyhedron using list is VOID.
p.reserve(num_vertex+num_edge, 2*2*num_edge, 4*num_edge/2);
Point* vertex_point_buffer = new Point[num_vertex + num_edge];
Point* edge_point_buffer = vertex_point_buffer + num_vertex;
std::vector<bool> v_onborder(num_vertex);
size_t sb = p.size_of_border_edges();
Edge_iterator eitr = p.edges_begin();
for (size_t i = 0; i < num_edge-sb; i++, ++eitr)
mask.edge_node(eitr, edge_point_buffer[i]);
for (size_t i = num_edge-sb; i < num_edge; i++, ++eitr) {
int v = (int) std::distance(p.vertices_begin(), eitr->vertex());
v_onborder[v] = true;
mask.border_node(eitr, edge_point_buffer[i], vertex_point_buffer[v]);
}
Vertex_iterator vitr = p.vertices_begin();
for (size_t i = 0; i < num_vertex; i++, ++vitr)
if (!v_onborder[i]) mask.vertex_node(vitr, vertex_point_buffer[i]);
// Build the connectivity using insert_vertex() and insert_edge()
// 1. insert_vertex() to all edges and set them to new positions
// 2. insert_edge() between 2 randomly selected neighboring new inserted
// vertices
// 3. insert_vertex() to the new inserted edge and set them to new positions
// 4. insert_edge() between all other new inserted vertices of step 1 and
// the new inserted vertex of step 3
// Step 1.
eitr = p.edges_begin();
for (size_t i = 0; i < num_edge; i++, ++eitr) {
Vertex_handle vh = PD::insert_vertex(p, eitr);
vh->point() = edge_point_buffer[i];
}
Facet_iterator fitr = p.facets_begin();
for (size_t i = 0; i < num_facet; i++, ++fitr) {
// Step 2.
Halfedge_around_facet_circulator hcir_begin = fitr->facet_begin();
Halfedge_around_facet_circulator hcir = hcir_begin;
// After linsub, the facet valence = 6
CGAL_assertion(circulator_size(hcir)==6);
Halfedge_handle e1 = ++hcir;
++hcir;
Halfedge_handle e2 = ++hcir;
++hcir;
Halfedge_handle e3 = ++hcir;
e2 = PD::insert_edge(p, e1, e2);
e3 = PD::insert_edge(p, e2, e3);
PD::insert_edge(p, e3, e1);
}
// Update the geometry data of the newly inserted vertices by the
// vertices buffer
vitr = p.vertices_begin();
for (size_t i = 0; i < num_vertex; i++, ++vitr)
vitr->point() = vertex_point_buffer[i];
delete []vertex_point_buffer;
}
// ======================================================================
//#define _EULER_DQQ_SPLITTING // should be CGAL_ prefixed
//#define _EULER_DQQ_TILTING // Tilting is faster // commented out since unused and not CGAL_ prefixed
template <class Poly, template <typename> class Mask>
void DQQ_1step(Poly& p, Mask<Poly> mask) {
typedef Polyhedron_decorator_3<Poly> PD;
typedef typename Poly::Halfedge_data_structure HDS;
typedef typename Poly::Vertex Vertex;
typedef typename Poly::Halfedge Halfedge;
typedef typename Poly::Facet Facet;
typedef typename Poly::Vertex_handle Vertex_handle;
typedef typename Poly::Halfedge_handle Halfedge_handle;
typedef typename Poly::Facet_handle Facet_handle;
typedef typename Poly::Vertex_iterator Vertex_iterator;
typedef typename Poly::Halfedge_iterator Halfedge_iterator;
typedef typename Poly::Edge_iterator Edge_iterator;
typedef typename Poly::Facet_iterator Facet_iterator;
typedef typename Poly::Halfedge_around_facet_circulator
Halfedge_around_facet_circulator;
typedef typename Poly::Halfedge_around_vertex_circulator
Halfedge_around_vertex_circulator;
typedef typename Poly::Traits Traits;
typedef typename Traits::Kernel Kernel;
typedef typename Kernel::Point_3 Point;
p.normalize_border();
size_t num_v = p.size_of_vertices();
size_t num_e = p.size_of_halfedges()/2;
size_t num_f = p.size_of_facets();
size_t num_be = p.size_of_border_edges();
Point* point_buffer = new Point[num_e*2];
//
#ifdef _EULER_DQQ_SPLITTING
//
// Splitting
//! Splitting is not implemented to support border
// build the point_buffer
Facet_iterator fitr, fitr_end = p.facets_end();
int pi = 0;
for (fitr = p.facets_begin(); fitr != fitr_end; ++fitr) {
Halfedge_around_facet_circulator cir = fitr->facet_begin();
do {
mask.corner_node(cir, point_buffer[pi++]);
} while (--cir != fitr->facet_begin());
}
// If Polyhedron is using vector, we need to reserve the memory to prevent
// the CGAL_assertion. This function for polyhedron using list is VOID.
p.reserve(num_v+num_e+num_f, 2*num_e, (2+4+2)*num_e);
// Build the connectivity using insert_vertex() and insert_edge()
// 1. create barycentric centers of each facet
fitr = p.facets_begin();
pi = 0;
for (size_t i = 0; i < num_f; i++) {
Facet_handle fh = fitr;
++fitr;
Vertex_handle vh = (p.create_center_vertex(fh->facet_begin()))->vertex();
// 1.1 add vertex on each new edges
Halfedge_around_vertex_circulator vcir = vh->vertex_begin();
int vn = circulator_size(vcir);
for (int j = 0; j < vn; ++j) {
Halfedge_handle e = vcir;
++vcir;
Vertex_handle v = PD::insert_vertex(p, e);
v->point() = point_buffer[pi++];
}
// 1.2 connect new vertices surround each barycentric center
for (int j = 0; j < vn; ++j) {
Halfedge_handle e1 = vcir->prev();
++vcir;
Halfedge_handle e2 = vcir->opposite();
PD::insert_edge(p, e1, e2);
}
// 1.3 remove the barycentric centers
p.erase_center_vertex(vcir);
}
// 2. remove old edges
Edge_iterator eitr = p.edges_begin();
for (size_t i = 0; i < num_e; ++i) {
Halfedge_handle eh = eitr;
++eitr;
p.join_facet(eh);
}
// 3. connect new vertices surround old vertices and then remove
// old vertices.
Vertex_iterator vitr = p.vertices_begin();
for (size_t i = 0; i < num_v; ++i) {
Halfedge_around_vertex_circulator vcir = vitr->vertex_begin();
int vn = circulator_size(vcir);
for (int j = 0; j < vn; ++j) {
Halfedge_handle e1 = vcir->prev();
++vcir;
Halfedge_handle e2 = vcir->opposite();
PD::insert_edge(p, e1, e2);
}
++vitr;
p.erase_center_vertex(vcir);
}
//
#else
//
// Tilting
// build the point_buffer
Vertex_iterator vitr, vitr_end = p.vertices_end();
int pi = 0;
for (vitr = p.vertices_begin(); vitr != vitr_end; ++vitr) {
Halfedge_around_vertex_circulator cir = vitr->vertex_begin();
do {
if (!cir->is_border()) mask.corner_node(cir, point_buffer[pi++]);
} while (++cir != vitr->vertex_begin());
}
// If Polyhedron is using vector, we need to reserve the memory to prevent
// the CGAL_assertion. This function for polyhedron using list is VOID.
p.reserve(num_v+num_e+num_f, 2*num_e, (2+4+2)*num_e);
// Build the connectivity using insert_vertex() and insert_edge()
vitr = p.vertices_begin();
pi = 0;
for (size_t i = 0; i < num_v; ++i) {
Vertex_handle vh = vitr;
++vitr;
Halfedge_around_vertex_circulator vcir = vh->vertex_begin();
size_t vn = circulator_size(vcir);
for (size_t j = 0; j < vn; ++j) {
Halfedge_handle e = vcir;
++vcir;
if (!e->is_border()) {
Vertex_handle v = PD::insert_vertex(p, e);
v->point() = point_buffer[pi++];
}
}
vcir = vh->vertex_begin();
for (size_t j = 0; j < vn; ++j) {
if (!vcir->is_border()) {
Halfedge_handle e1 = vcir->prev();
++vcir;
if (!vcir->is_border()) {
Halfedge_handle e2 = vcir->opposite();
PD::insert_edge(p, e1, e2);
}
} else ++vcir;
}
//p.erase_center_vertex(vh->vertex_begin());
}
Edge_iterator eitr = p.edges_begin();
for (size_t i = 0; i < num_e; ++i) {
Halfedge_handle eh = eitr;
++eitr;
if (!eh->is_border_edge()) {
PD::insert_edge(p, eh->prev()->prev(), eh);
eh = eh->opposite();
PD::insert_edge(p, eh->prev()->prev(), eh);
p.join_facet(eh);
} else {
if (eh->is_border()) {
eh = eh->opposite();
PD::insert_edge(p, eh, eh->prev()->prev());
} else
PD::insert_edge(p, eh->prev()->prev(), eh);
}
}
// after this point, the original border edges are in front!
eitr = p.edges_begin();
for (size_t i = 0; i < num_be; ++i) {
Halfedge_handle eh = eitr;
++eitr;
if (eh->is_border()) eh = eh->opposite();
Halfedge_handle ehe = eh;
eh = eh->prev()->opposite();
while (!eh->is_border()) {
p.erase_facet(ehe);
ehe = eh;
eh = eh->prev()->opposite();
}
p.erase_facet(ehe);
}
vitr = p.vertices_begin();
for (size_t i = 0; i < num_v-num_be; ++i) {
Vertex_handle vh = vitr;
++vitr;
p.erase_center_vertex(vh->vertex_begin());
}
#endif //_EULER_DQQ_SPLITTING
delete []point_buffer;
}
// ======================================================================
template <class Poly, template <typename> class Mask>
void Sqrt3_1step(Poly& p, Mask<Poly> mask) {
typedef Polyhedron_decorator_3<Poly> PD;
typedef typename Poly::Halfedge_data_structure HDS;
typedef typename Poly::Vertex Vertex;
typedef typename Poly::Halfedge Halfedge;
typedef typename Poly::Facet Facet;
typedef typename Poly::Vertex_handle Vertex_handle;
typedef typename Poly::Halfedge_handle Halfedge_handle;
typedef typename Poly::Facet_handle Facet_handle;
typedef typename Poly::Vertex_iterator Vertex_iterator;
typedef typename Poly::Halfedge_iterator Halfedge_iterator;
typedef typename Poly::Edge_iterator Edge_iterator;
typedef typename Poly::Facet_iterator Facet_iterator;
typedef typename Poly::Halfedge_around_facet_circulator
Halfedge_around_facet_circulator;
typedef typename Poly::Halfedge_around_vertex_circulator
Halfedge_around_vertex_circulator;
typedef typename Poly::Traits Traits;
typedef typename Traits::Kernel Kernel;
typedef typename Kernel::Point_3 Point;
//
p.normalize_border();
//
size_t num_v = p.size_of_vertices();
size_t num_e = p.size_of_halfedges()/2;
size_t num_f = p.size_of_facets();
p.reserve(num_v+num_f, (num_e+3*num_f)*2, 3*num_f);
// prepare the smoothed center points
Point* cpt = new Point[num_f];
Facet_iterator fitr = p.facets_begin();
for (size_t i = 0; i < num_f; ++i, ++fitr) {
//ASSERTION_MSG(circulator_size(fitr->facet_begin())==3, "(ERROR) Non-triangle facet!");
mask.facet_node(fitr, cpt[i]);
}
// smooth the vertex points
Vertex_iterator vitr = p.vertices_begin();
for (size_t i = 0; i < num_v; ++i, ++vitr)
mask.vertex_node(vitr, vitr->point());
// insert the facet points
fitr = p.facets_begin();
for (size_t i = 0; i < num_f; ++i, ++fitr) {
Halfedge_handle center = p.create_center_vertex(fitr->halfedge());
center->vertex()->point() = cpt[i];
}
delete []cpt;
// flip the old edges except the border edges
Edge_iterator eitr = p.edges_begin();
for (size_t i = 0; i < num_e; ++i) {
Halfedge_handle e = eitr;
++eitr; // move to next edge before flip since flip destroys current edge
if (!e->is_border_edge()) {
Halfedge_handle h = p.join_facet(e);
p.split_facet(h->prev(), h->next());
}
}
// TODO: border ...
CGAL_postcondition(p.is_valid());
}
}
}
CGAL_END_NAMESPACE
#endif //CGAL_POLYHEDRON_SUBDIVISION_H_01292002
Reference in New Issue View Git Blame Copy Permalink