Rev 150 | Rev 183 | Go to most recent revision | Blame | Compare with Previous | Last modification | View Log | RSS feed
#include <queue>
#include <vector>
#include "triangulate.h"
#include "HMesh/VertexCirculator.h"
#include "HMesh/FaceCirculator.h"
using namespace std;
using namespace CGLA;
using namespace HMesh;
namespace HMesh
{
void get_candidates(const VertexIter& v,
vector<HalfEdgeIter>& candidates)
{
vector<VertexIter> vertices;
vector<HalfEdgeIter> hedges;
VertexCirculator vc(v);
for(;!vc.end();++vc)
{
vertices.push_back(vc.get_vertex());
hedges.push_back(vc.get_halfedge());
}
size_t N = vc.no_steps();
vector<VertexIter> vertices_check(vertices);
assert(N==vertices.size());
for(int i=N-1;i>=0;--i)
{
HalfEdgeIter h = hedges[i]->next;
while(h->vert != vertices[(i+N-1)%N])
{
if(find(vertices_check.begin(),vertices_check.end(),
h->vert) == vertices_check.end())
candidates.push_back(h);
h = h->next;
// TODO : does not compile in debug mode in visual studio
//assert(h!=v);
}
}
}
float curv(const Vec3f& p, std::vector<Vec3f>& vec)
{
int N = vec.size();
std::vector<Vec3f> normals;
for(int i=0;i<N;++i)
{
Vec3f norm = normalize(cross(vec[i]-p, vec[(i+1)%N]-p));
normals.push_back(norm);
}
float alpha = 0;
for(int i=0;i<N;++i)
{
alpha += (vec[i]-p).length()*acos(dot(normals[i],normals[(i+1)%N]));
}
return alpha;
}
float get_badness(const VertexIter& v, const VertexIter& n)
{
vector<HalfEdgeIter> hedges;
VertexCirculator vc(v);
for(;!vc.end();++vc)
hedges.push_back(vc.get_halfedge());
vector<Vec3f> one_ring;
vector<Vec3f> one_ring_n;
int N = vc.no_steps();
for(int i=N-1;i>=0;--i)
{
HalfEdgeIter h = hedges[i];
h = h->next;
while(h->vert != v)
{
one_ring.push_back(h->vert->pos);
if(h->vert != n)
one_ring_n.push_back(h->vert->pos);
h = h->next;
}
}
return curv(v->pos,one_ring)-curv(v->pos,one_ring_n);
}
const CGLA::Vec3f get_normal(const VertexIter& v)
{
vector<HalfEdgeIter> hedges;
VertexCirculator vc(v);
for(;!vc.end();++vc)
hedges.push_back(vc.get_halfedge());
vector<Vec3f> one_ring;
int N = vc.no_steps();
for(int i=N-1;i>=0;--i)
{
HalfEdgeIter h = hedges[i];
h = h->next;
while(h->vert != v)
{
one_ring.push_back(h->vert->pos);
h = h->next;
}
}
Vec3f norm(0);
N= one_ring.size();
Vec3f p0 = v->pos;
for(int i=0;i<N;++i)
{
Vec3f p1 = one_ring[i];
Vec3f p2 = one_ring[(i+1)%N];
Vec3f e0 = normalize(p1-p0);
Vec3f e1 = normalize(p2-p0);
norm += cross(e0,e1)*acos(dot(e0,e1));
}
return normalize(norm);
}
void safe_triangulate(Manifold& m)
{
vector<FaceIter> fv;
for(FaceIter fi=m.faces_begin(); fi != m.faces_end(); ++fi)
fv.push_back(fi);
for(size_t i=0;i<fv.size();++i)
m.safe_triangulate(fv[i]);
}
void triangulate(Manifold& m)
{
vector<FaceIter> fv;
for(FaceIter fi=m.faces_begin(); fi != m.faces_end(); ++fi)
fv.push_back(fi);
for(size_t i=0;i<fv.size();++i)
m.triangulate(fv[i]);
}
struct PotentialEdge
{
int time_tag;
float badness;
FaceIter f;
VertexIter v0, v1;
};
bool operator>(const PotentialEdge& e0, const PotentialEdge& e1)
{
return e0.badness>e1.badness;
}
typedef std::priority_queue<PotentialEdge,
std::vector<PotentialEdge>,
std::greater<PotentialEdge> >
PotentialEdgeQueue;
void insert_potential_edges(const VertexIter& v,
PotentialEdgeQueue& pot_edges)
{
vector<Vec3f> one_ring;
vector<HalfEdgeIter> candidates;
get_candidates(v, candidates);
Vec3f p0 = v->pos;
Vec3f norm = normal(v);
int n = candidates.size();
for(int i=0;i<n;++i)
{
VertexIter v_n = candidates[i]->vert;
Vec3f edir = normalize(v_n->pos-p0);
Vec3f norm_n = normal(v_n);
float bad = sqr(dot(edir, norm));
float bad_n = sqr(dot(edir, norm_n));
PotentialEdge pot;
/* So if the edge between two vertices is not orthogonal to
their normals, the badness is increased. Badness is also
increased if the normals are very different. */
pot.badness = bad+bad_n - dot(norm_n,norm);
pot.time_tag = v->touched;
pot.v0 = v;
pot.v1 = candidates[i]->vert;
pot.f = candidates[i]->face;
pot_edges.push(pot);
}
}
void create_candidates(Manifold& m, PotentialEdgeQueue& pot_edges)
{
for(VertexIter v= m.vertices_begin(); v!= m.vertices_end();++v)
{
v->touched = 0;
insert_potential_edges(v, pot_edges);
}
}
void curvature_triangulate(Manifold& m)
{
PotentialEdgeQueue pot_edges;
// Create candidates
create_candidates(m,pot_edges);
while(!pot_edges.empty())
{
const PotentialEdge& pot_edge = pot_edges.top();
// Record all the vertices of the face. We need to
// recompute the candidates.
std::vector<VertexIter> reeval_vec;
FaceCirculator fc(pot_edge.f);
while(!fc.end())
{
reeval_vec.push_back(fc.get_vertex());
++fc;
}
VertexIter v0 = pot_edge.v0;
// Make sure that the vertex has not been touched since
// we created the potential edge. If the vertex has been
// touched the candidate edge has either (a) been processed,
// (b) received a lower priority or (c) become invalid.
if(pot_edge.time_tag == v0->touched)
{
VertexIter v1 = pot_edge.v1;
vector<Vec3f> one_ring;
vector<HalfEdgeIter> candidates;
m.split_face(pot_edge.f, v0, v1);
// Recompute priorities.
for(size_t i=0;i<reeval_vec.size();++i)
{
VertexIter& v = reeval_vec[i];
v->touched++;
insert_potential_edges(v, pot_edges);
}
}
pot_edges.pop();
}
}
}