Files
EgtGeomKernel/Triangulate.cpp
T
Dario Sassi 77e74ccf4e EgtGeomKernel 1.5h3 :
- aggiunta IsFlat a tutte le Curve 
- aggiunta ApproxWithArcs a tutte le Curve 
- aggiunto oggetto PolyArc (raccolta ordinata di linee e archi con bulge)
- aggiunto oggetto PointsPCA per stima componenti principali di un insieme di punti
- FromSpheriical e FromPolar di Vector3d sono diventati funzioni e aggiunto FromUprightOrtho
- aggiunte Invert e a Vector3d.
2014-08-15 17:36:08 +00:00

336 lines
11 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
//----------------------------------------------------------------------------
// EgalTech 2014-2014
//----------------------------------------------------------------------------
// File : Triangulate.cpp Data : 23.06.14 Versione : 1.5f6
// Contenuto : Implementazione della classe Triangulate.
//
//
//
// Modifiche : 08.04.14 DS Creazione modulo.
// 23.06.14 DS Corretto errore con contorni CW.
//
//----------------------------------------------------------------------------
//--------------------------- Include ----------------------------------------
#include "stdafx.h"
#include "Triangulate.h"
#include "\EgtDev\Include\EGkPolyLine.h"
#include "\EgtDev\Include\EGkPlane3d.h"
using namespace std ;
//----------------------------------------------------------------------------
// In : PolyLine
// Out : PNTVECTOR (Point3d Vector) : points of the polyline
// INTVECTOR (int Vector) : 3*T indices of above points for T triangles
//----------------------------------------------------------------------------
bool
Triangulate::Make( const PolyLine& PL, PNTVECTOR& vPt, INTVECTOR& vTr)
{
// verific che la polilinea sia chiusa e calcolo il piano medio del poligono
double dArea ;
Plane3d plPlane ;
if ( ! PL.IsClosedAndFlat( plPlane, dArea, 50 * EPS_SMALL))
return false ;
bool bCCW ;
if ( fabs( plPlane.vtN.z) >= fabs( plPlane.vtN.x) &&
fabs( plPlane.vtN.z) >= fabs( plPlane.vtN.y)) {
m_nPlane = PL_XY ;
bCCW = ( plPlane.vtN.z > 0) ;
}
else if ( fabs( plPlane.vtN.x) >= fabs( plPlane.vtN.y)) {
m_nPlane = PL_YZ ;
bCCW = ( plPlane.vtN.x > 0) ;
}
else {
m_nPlane = PL_ZX ;
bCCW = ( plPlane.vtN.y > 0) ;
}
// riempio il vettore con i vertici del poligono da triangolare
vPt.clear() ;
vPt.reserve( PL.GetPointNbr() - 1) ;
// salto il primo punto (coincide con l'ultimo)
Point3d ptP ;
if ( ! PL.GetFirstPoint( ptP))
return false ;
// inserisco i punti
while ( PL.GetNextPoint( ptP))
vPt.push_back( ptP) ;
// creo il vettore degli indici del Poligono
INTVECTOR vPol ;
int n = int( vPt.size()) ;
vPol.reserve( n) ;
// se orientato correttamente (componente di N > 0)
if ( bCCW) {
for ( int i = 0 ; i < n ; ++ i)
vPol.push_back( i) ;
}
// altrimenti lo prendo al contrario
else {
for ( int i = n - 1 ; i >= 0 ; -- i)
vPol.push_back( i) ;
}
// eseguo la triangolazione
if ( ! MakeByEC2( vPt, vPol, vTr))
return false ;
// se era CW, devo invertire il senso dei triangoli
if ( ! bCCW)
reverse( vTr.begin(), vTr.end()) ;
return true ;
}
//----------------------------------------------------------------------------
// Triangulate the CCW n-gon specified by the vertices vPt (Pt[n] != Pt[0])
// Ear Clipping algorithm
//----------------------------------------------------------------------------
bool
Triangulate::MakeByEC( const PNTVECTOR& vPt, const INTVECTOR& vPol, INTVECTOR& vTr)
{
// Clear triangle vector
vTr.clear() ;
// At least 3 points
int n = int( vPol.size()) ;
if ( n < 3)
return false ;
// Preallocate triangle vector ( #triangles = n - 2)
vTr.reserve( 3 * ( n - 2)) ;
// Set up previous and next links to effectively form a double-linked vertex list
INTVECTOR vPrev( n) ;
INTVECTOR vNext( n) ;
for ( int j = 0 ; j < n ; ++ j) {
vPrev[j] = j - 1 ;
vNext[j] = j + 1 ;
}
vPrev[0] = n - 1 ;
vNext[n-1] = 0 ;
// Start at vertex 0
int i = 0 ;
int nCount = n ;
// Keep removing vertices until just a triangle left
while ( n > 3) {
// To avoid infinite loop
if ( nCount <= 0)
return false ;
-- nCount ;
// Test if current vertex, v[i], is an ear
bool bIsEar = TestTriangle( vPt, vPol, vPrev, vNext, i) ;
// If current vertex v[i] is an ear, delete it and visit the previous vertex
if ( bIsEar) {
// Triangle (v[prev[i]], v[i], v[next[i]]) is an ear
vTr.push_back( vPol[vPrev[i]]) ;
vTr.push_back( vPol[i]) ;
vTr.push_back( vPol[vNext[i]]) ;
// Delete vertex v[i] by redirecting next and previous links
// of neighboring verts past it. Decrement vertex count
vNext[vPrev[i]] = vNext[i] ;
vPrev[vNext[i]] = vPrev[i] ;
n--;
// Visit the previous vertex next
i = vPrev[i] ;
// Reset Count
nCount = n ;
}
else {
// Current vertex is not an ear; visit the next vertex
i = vNext[i] ;
}
}
// Last triangle is an ear
vTr.push_back( vPol[vPrev[i]]) ;
vTr.push_back( vPol[i]) ;
vTr.push_back( vPol[vNext[i]]) ;
return true ;
}
//----------------------------------------------------------------------------
// Triangulate the CCW n-gon specified by the vertices vPt (Pt[n] != Pt[0])
// Ear Clipping algorithm enhanced
//----------------------------------------------------------------------------
bool
Triangulate::MakeByEC2( const PNTVECTOR& vPt, const INTVECTOR& vPol, INTVECTOR& vTr)
{
// Clear triangle vector
vTr.clear() ;
// At least 3 points
int n = int( vPol.size()) ;
if ( n < 3)
return false ;
// Preallocate triangle vector ( #triangles = n - 2)
vTr.reserve( 3 * ( n - 2)) ;
// Set up previous and next links to effectively form a double-linked vertex list
INTVECTOR vPrev( n) ;
INTVECTOR vNext( n) ;
for ( int j = 0 ; j < n ; ++ j) {
vPrev[j] = j - 1 ;
vNext[j] = j + 1 ;
}
vPrev[0] = n - 1 ;
vNext[n-1] = 0 ;
// Start at vertex 0
int i = 0 ;
int nCount = n ;
// Keep removing vertices until just a triangle left
while ( n > 3) {
// To avoid infinite loop
if ( nCount <= 0)
return false ;
-- nCount ;
// Test if current vertex, v[i], is an ear
bool bIsEar = TestTriangle( vPt, vPol, vPrev, vNext, i) ;
if ( bIsEar) {
// Save square distance of diagonal
double dSqDist = SqDist(vPt[vPol[vPrev[i]]], vPt[vPol[vNext[i]]]) ;
// Try with next or next^2
int j = vNext[i] ;
double dSqDist1 = INFINITO ;
if ( TestTriangle( vPt, vPol, vPrev, vNext, j))
dSqDist1 = SqDist( vPt[vPol[vPrev[j]]], vPt[vPol[vNext[j]]]) ;
else {
j = vNext[j] ;
if ( TestTriangle( vPt, vPol, vPrev, vNext, j))
dSqDist1 = SqDist( vPt[vPol[vPrev[j]]], vPt[vPol[vNext[j]]]) ;
}
// Try with prev or prev^2
int k = vPrev[i] ;
double dSqDist2 = INFINITO ;
if ( TestTriangle( vPt, vPol, vPrev, vNext, k))
dSqDist2 = SqDist( vPt[vPol[vPrev[k]]], vPt[vPol[vNext[k]]]) ;
else {
k = vPrev[k] ;
if ( TestTriangle( vPt, vPol, vPrev, vNext, k))
dSqDist2 = SqDist( vPt[vPol[vPrev[k]]], vPt[vPol[vNext[k]]]) ;
}
// Choose the better (EPS_ZERO to have a difference)
if ( dSqDist < dSqDist1 + EPS_ZERO && dSqDist < dSqDist2 + EPS_ZERO)
; // i is the better
else if ( dSqDist1 < dSqDist2 + EPS_ZERO)
i = j ;
else
i = k ;
}
// If current vertex v[i] is an ear, delete it and visit the previous vertex
if ( bIsEar) {
// Triangle (v[prev[i]], v[i], v[next[i]]) is an ear
vTr.push_back( vPol[vPrev[i]]) ;
vTr.push_back( vPol[i]) ;
vTr.push_back( vPol[vNext[i]]) ;
// Delete vertex v[i] by redirecting next and previous links
// of neighboring verts past it. Decrement vertex count
vNext[vPrev[i]] = vNext[i] ;
vPrev[vNext[i]] = vPrev[i] ;
n--;
// Visit the previous vertex next
i = vPrev[i] ;
// Reset Count
nCount = n ;
}
else {
// Current vertex is not an ear; visit the next vertex
i = vNext[i] ;
}
}
// Last triangle is an ear
vTr.push_back( vPol[vPrev[i]]) ;
vTr.push_back( vPol[i]) ;
vTr.push_back( vPol[vNext[i]]) ;
return true ;
}
//----------------------------------------------------------------------------
bool
Triangulate::TestTriangle( const PNTVECTOR& vPt, const INTVECTOR& vPol,
const INTVECTOR& vPrev, INTVECTOR& vNext, int i)
{
// Test if current vertex, v[i], is an ear
bool bIsEar = true ;
// An ear must be convex (here counterclockwise)
if ( TriangleIsCCW( vPt[vPol[vPrev[i]]], vPt[vPol[i]], vPt[vPol[vNext[i]]])) {
// Loop over all vertices not part of the tentative ear
int k = vNext[vNext[i]] ;
do {
// If vertex k is inside the ear triangle, then this is not an ear
if ( TestPointInTriangle( vPt[vPol[k]], vPt[vPol[vPrev[i]]], vPt[vPol[i]], vPt[vPol[vNext[i]]])) {
bIsEar = false ;
break ;
}
k = vNext[k] ;
} while (k != vPrev[i]) ;
}
else {
// The ear triangle is clockwise so v[i] is not an ear
bIsEar = false ;
}
return bIsEar ;
}
//----------------------------------------------------------------------------
bool
Triangulate::TriangleIsCCW( const Point3d& ptA, const Point3d& ptB, const Point3d& ptC, double dToler)
{
double dV11 ;
double dV12 ;
double dV21 ;
double dV22 ;
switch ( m_nPlane) {
default : // PL_XY
dV11 = ptB.x - ptA.x ;
dV12 = ptB.y - ptA.y ;
dV21 = ptC.x - ptB.x ;
dV22 = ptC.y - ptB.y ;
break ;
case PL_YZ :
dV11 = ptB.y - ptA.y ;
dV12 = ptB.z - ptA.z ;
dV21 = ptC.y - ptB.y ;
dV22 = ptC.z - ptB.z ;
break ;
case PL_ZX :
dV11 = ptB.z - ptA.z ;
dV12 = ptB.x - ptA.x ;
dV21 = ptC.z - ptB.z ;
dV22 = ptC.x - ptB.x ;
break ;
}
return ( ( dV11 * dV22 - dV12 * dV21) > dToler * dToler) ;
}
//----------------------------------------------------------------------------
// Test if point p is contained in triangle (a, b, c)
bool
Triangulate::TestPointInTriangle( const Point3d& ptP, const Point3d& ptA, const Point3d& ptB, const Point3d& ptC)
{
// if P is on a vertex is considered outside
if ( AreSamePointApprox( ptP, ptA))
return false ;
if ( AreSamePointApprox( ptP, ptB))
return false ;
if ( AreSamePointApprox( ptP, ptC))
return false ;
// if P is on the right of at least one edge is outside
if ( TriangleIsCCW( ptA, ptP, ptB))
return false ;
if ( TriangleIsCCW( ptB, ptP, ptC))
return false ;
if ( TriangleIsCCW( ptC, ptP, ptA))
return false ;
return true ;
}