EgtGeomKernel 1.8f2 :

- migliorie e correzioni a Zmap.
This commit is contained in:
Dario Sassi
2017-06-10 15:43:14 +00:00
parent 17425bc753
commit d9d8b40a02
7 changed files with 287 additions and 367 deletions
+146 -272
View File
@@ -205,18 +205,14 @@ VolZmap::GetDexelLines( int nDir, int nPos1, int nPos2, POLYLINELIST& lstPL) con
bool
VolZmap::GetAllTriangles( TRIA3DLIST& lstTria) const
{
if ( m_nMapNum == 1) {
const int MAX_DIM_CHUNK = 128 ;
for ( int i = 0 ; i < int( m_nNx[0]) ; i += MAX_DIM_CHUNK) {
int nDimChunkX = min( MAX_DIM_CHUNK, int( m_nNx[0]) - i) ;
for ( int j = 0 ; j < int( m_nNy[0]) ; j += MAX_DIM_CHUNK) {
int nDimChunkY = min( MAX_DIM_CHUNK, int( m_nNy[0]) - j) ;
GetChunkPrisms( i, j, nDimChunkX, nDimChunkY, MAX_DIM_CHUNK, lstTria) ;
}
}
INTVECTOR nModifiedBlocks ;
TRIA3DLISTVECTOR vLstTria ;
if ( ! GetTriangles( true, nModifiedBlocks, vLstTria))
return false ;
lstTria.clear() ;
for ( size_t i = 0 ; i < vLstTria.size() ; ++ i) {
lstTria.splice( lstTria.end(), vLstTria[i]) ;
}
else
MarchingCubes( lstTria) ;
return true ;
}
@@ -225,6 +221,17 @@ VolZmap::GetAllTriangles( TRIA3DLIST& lstTria) const
bool
VolZmap::GetTriangles( bool bAllBlocks, INTVECTOR& nModifiedBlocks, TRIA3DLISTVECTOR& vLstTria) const
{
// Se nessun blocco modificato, è richiesta esterna e li considero tutti modificati
bool bSomeModif = false ;
for ( size_t i = 0 ; i < m_nNumBlock ; ++ i) {
if ( m_BlockToUpdate[i]) {
bSomeModif = true ;
break ;
}
}
if ( ! bSomeModif)
bAllBlocks = true ;
// Caso di singola mappa
if ( m_nMapNum == 1) {
@@ -281,82 +288,82 @@ VolZmap::GetTriangles( bool bAllBlocks, INTVECTOR& nModifiedBlocks, TRIA3DLISTVE
TriaMatrix VecTriHold ;
VecTriHold.resize( m_nNumBlock) ;
bool bCalcInterBlock = false ;
// Ciclo sui blocchi
// Calcolo i triangoli sui blocchi
for ( size_t t = 0 ; t < m_nNumBlock ; ++ t) {
// Calcolo i limiti sugli indici dei voxel del blocco
int nIJK[3] ;
// Vettore indici i,j,k del blocco
GetBlockIJKFromN( int( t), nIJK) ;
// Vettore limiti sugli indici dei voxel del blocco
int LimitIndexes[6] ;
GetBlockLimitsIJK( nIJK, LimitIndexes) ;
// Se il blocco deve essere processato, eseguo
// Se il blocco deve essere processato
if ( bAllBlocks || m_BlockToUpdate[t]) {
// processo ...
vLstTria.emplace_back() ;
nModifiedBlocks[t] = int( vLstTria.size()) - 1 ;
m_InterBlockTria[t].clear() ;
ExtMarchingCubes( LimitIndexes, t, vLstTria.back(), VecTriHold[t]) ;
// Flipping fra voxel interni
FlipEdgesII( VecTriHold[t]) ;
#if 1
ExtMarchingCubes( int( t), vLstTria.back(), VecTriHold[t]) ;
// Flipping fra voxel interni
FlipEdgesII( VecTriHold[t]) ;
bCalcInterBlock = true ;
#else
MarchingCubes( int( t), vLstTria.back()) ;
#endif
m_BlockToUpdate[t] = false ;
}
else
nModifiedBlocks[t] = -1 ;
}
// Copio i triangoli di frontiera in una matrice gemella
// di m_InterBlockTria per avere sempre a disposizione
// i triangoli non flippati.
// Calcolo i triangoli di frontiera tra feature di blocchi diversi
// copio i triangoli di frontiera in una matrice gemella
// di m_InterBlockTria per avere sempre a disposizione
// i triangoli non flippati.
TriaMatrix InterBlockTria ;
InterBlockTria = m_InterBlockTria ;
FlipEdgesBB( InterBlockTria) ;
if ( bCalcInterBlock) {
InterBlockTria = m_InterBlockTria ;
FlipEdgesBB( InterBlockTria) ;
}
// ciclo sui blocchi
// Inserisco in lista i triangoli di feature derivanti dai blocchi
for ( size_t t = 0 ; t < m_nNumBlock ; ++ t) {
// ciclo sui voxel del blocco
for ( size_t t1 = 0 ; t1 < VecTriHold[t].size() ; ++ t1) {
// ciclo sulle componenti connesse del voxel
for ( size_t t2 = 0 ; t2 < VecTriHold[t][t1].vCompoTria.size() ; ++ t2) {
// ciclo sui triangoli delle componenti connesse
for ( size_t t3 = 0 ; t3 < VecTriHold[t][t1].vCompoTria[t2].size() ; ++ t3) {
// aggiungo un singolo triangolo alla lista
if ( nModifiedBlocks[t] >= 0)
if ( nModifiedBlocks[t] >= 0) {
// ciclo sui voxel del blocco
for ( size_t t1 = 0 ; t1 < VecTriHold[t].size() ; ++ t1) {
// ciclo sulle componenti connesse del voxel
for ( size_t t2 = 0 ; t2 < VecTriHold[t][t1].vCompoTria.size() ; ++ t2) {
// ciclo sui triangoli delle componenti connesse
for ( size_t t3 = 0 ; t3 < VecTriHold[t][t1].vCompoTria[t2].size() ; ++ t3) {
// aggiungo triangolo alla lista
vLstTria[nModifiedBlocks[t]].emplace_back( VecTriHold[t][t1].vCompoTria[t2][t3]) ;
}
}
}
}
vLstTria.resize( vLstTria.size() + 1) ;
size_t nPos = size_t( vLstTria.size() - 1) ;
for ( size_t t = 0 ; t < m_InterBlockTria.size() ; ++ t) {
for ( size_t t1 = 0 ; t1 < m_InterBlockTria[t].size() ; ++ t1) {
for ( size_t t2 = 0 ; t2 < m_InterBlockTria[t][t1].vCompoTria.size() ; ++ t2) {
for ( size_t t3 = 0 ; t3 < m_InterBlockTria[t][t1].vCompoTria[t2].size() ; ++ t3) {
if ( m_InterBlockTria[t][t1].vCompoTria[t2][t3].GetArea() > EPS_SMALL) {
vLstTria[nPos].emplace_back( InterBlockTria[t][t1].vCompoTria[t2][t3]) ;
}
}
}
}
}
// Inserisco in lista i triangoli di frontiera tra feature di blocchi diversi
if ( bCalcInterBlock) {
vLstTria.resize( vLstTria.size() + 1) ;
size_t nPos = size_t( vLstTria.size() - 1) ;
for ( size_t t = 0 ; t < m_InterBlockTria.size() ; ++ t) {
for ( size_t t1 = 0 ; t1 < m_InterBlockTria[t].size() ; ++ t1) {
for ( size_t t2 = 0 ; t2 < m_InterBlockTria[t][t1].vCompoTria.size() ; ++ t2) {
for ( size_t t3 = 0 ; t3 < m_InterBlockTria[t][t1].vCompoTria[t2].size() ; ++ t3) {
if ( m_InterBlockTria[t][t1].vCompoTria[t2][t3].GetArea() > EPS_SMALL) {
vLstTria[nPos].emplace_back( InterBlockTria[t][t1].vCompoTria[t2][t3]) ;
}
}
}
}
}
nModifiedBlocks.back() = int( nPos) ;
}
nModifiedBlocks[nModifiedBlocks.size()-1] = int( vLstTria.size() - 1) ;
else
nModifiedBlocks.back() = - 1 ;
}
return true ;
@@ -560,106 +567,13 @@ VolZmap::AddDexelSideFace( int nPos, int nPosAdj, const Point3d& ptP, const Poin
return true ;
}
//----------------------------------------------------------------------------
bool
VolZmap::MarchingCubes( TRIA3DLIST& lstTria) const
{
// Ciclo su tutti i voxel dello Zmap
for ( int i = - 1 ; i < int( m_nNx[0]) ; ++ i) {
for ( int j = - 1 ; j < int( m_nNy[0]) ; ++ j) {
for ( int k = - 1 ; k < int( m_nNy[1]) ; ++ k) {
// Indici i,j,k dei vertici
int IndexCorner[8][3] = {
{ i, j, k},
{ i + 1, j, k},
{ i + 1, j + 1, k},
{ i, j + 1, k},
{ i, j, k + 1},
{ i + 1, j, k + 1},
{ i + 1, j + 1, k + 1},
{ i, j + 1, k + 1}
} ;
// Classificazione dei vertici: interni o esterni al materiale
int nIndex = 0 ;
if ( IsThereMat( i, j, k))
nIndex |= ( 1 << 0) ;
if ( IsThereMat( i + 1, j, k))
nIndex |= ( 1 << 1) ;
if ( IsThereMat( i + 1, j + 1, k))
nIndex |= ( 1 << 2) ;
if ( IsThereMat( i, j + 1, k))
nIndex |= ( 1 << 3) ;
if ( IsThereMat( i, j, k + 1))
nIndex |= ( 1 << 4) ;
if ( IsThereMat( i + 1, j, k + 1))
nIndex |= ( 1 << 5) ;
if ( IsThereMat( i + 1, j + 1, k + 1))
nIndex |= ( 1 << 6) ;
if ( IsThereMat( i, j + 1, k + 1))
nIndex |= ( 1 << 7) ;
// Se vi è qualche intersezione fra segmenti e superficie
// continuo altrimenti passo al prossimo voxel
if ( EdgeTable[nIndex] == 0)
continue ;
static int intersections[12][2] = {
{ 0, 1 }, { 1, 2 }, { 2, 3 }, { 3, 0 }, { 4, 5 }, { 5, 6 },
{ 6, 7 }, { 7, 4 }, { 0, 4 }, { 1, 5 }, { 2, 6 }, { 3, 7 }
} ;
Point3d ptIntPoint[12] ;
// Ciclo sui segmenti
for ( int EdgeIndex = 0 ; EdgeIndex < 12 ; ++ EdgeIndex) {
// Se il segmento non attraversa la superficie
// passo al successivo
if ( ! ( EdgeTable[nIndex] & ( 1 << EdgeIndex)))
continue ;
int n1 = intersections[EdgeIndex][0] ;
int n2 = intersections[EdgeIndex][1] ;
// Determino con precisione il punto di intersezione sullo spigolo
IntersPos( IndexCorner[n1], IndexCorner[n2], ptIntPoint[EdgeIndex]) ;
ptIntPoint[EdgeIndex].ToGlob( m_MapFrame[0]) ;
}
// Costruzione dei triangoli
for ( int TriIndex = 0 ; TriangleTableEn[nIndex][0][TriIndex] != - 1 ; TriIndex += 3) {
// Costruzione triangolo
int i0 = TriangleTableEn[nIndex][0][TriIndex + 2] ;
int i1 = TriangleTableEn[nIndex][0][TriIndex + 1] ;
int i2 = TriangleTableEn[nIndex][0][TriIndex] ;
// Il triangolo è pronto
Triangle3d CurrentTriangle ;
CurrentTriangle.Set( ptIntPoint[i0], ptIntPoint[i1], ptIntPoint[i2]) ;
CurrentTriangle.Validate() ;
// Aggiungo triangolo
lstTria.emplace_back( CurrentTriangle) ;
}
}
}
}
return true ;
}
//----------------------------------------------------------------------------
bool
VolZmap::MarchingCubes( int nBlock, TRIA3DLIST& lstTria) const
{
if ( nBlock < 0 || nBlock >= int( m_BlockToUpdate.size()))
if ( nBlock < 0 || nBlock >= int( m_nNumBlock))
return false ;
Point3d ptMapOrig = m_MapFrame[0].Orig() ;
// Calcolo posizione del blocco nel reticolo
int nIBlock = ( nBlock % int( m_nFracLin[0] * m_nFracLin[1])) % int( m_nFracLin[0]) ;
int nJBlock = ( nBlock % int( m_nFracLin[0] * m_nFracLin[1])) / int( m_nFracLin[0]) ;
@@ -667,18 +581,18 @@ VolZmap::MarchingCubes( int nBlock, TRIA3DLIST& lstTria) const
// Calcolo limiti per l'indice i
int nStartI = nIBlock * int( m_nDexNumPBlock) - 1 ;
int nEndI = ( nIBlock + 1 == int( m_nFracLin[0]) ?
int( m_nNx[0]) : ( nIBlock + 1) * int( m_nDexNumPBlock)) ;
int nEndI = ( nIBlock + 1 == int( m_nFracLin[0]) ?
int( m_nNx[0]) : ( nIBlock + 1) * int( m_nDexNumPBlock)) ;
// Calcolo limiti per l'indice j
int nStartJ = nJBlock * int( m_nDexNumPBlock) - 1 ;
int nEndJ = ( nJBlock + 1 == int( m_nFracLin[1]) ?
int( m_nNy[0]) : ( nJBlock + 1) * int( m_nDexNumPBlock)) ;
int nEndJ = ( nJBlock + 1 == int( m_nFracLin[1]) ?
int( m_nNy[0]) : ( nJBlock + 1) * int( m_nDexNumPBlock)) ;
// Calcolo limiti per l'indice k
int nStartK = nKBlock * int( m_nDexNumPBlock) - 1 ;
int nEndK = ( nKBlock + 1 == int( m_nFracLin[2]) ?
int( m_nNy[1]) : ( nKBlock + 1) * int( m_nDexNumPBlock)) ;
int nEndK = ( nKBlock + 1 == int( m_nFracLin[2]) ?
int( m_nNy[1]) : ( nKBlock + 1) * int( m_nDexNumPBlock)) ;
// Ciclo su tutti i voxel dello Zmap
for ( int i = nStartI ; i < nEndI ; ++ i) {
@@ -773,24 +687,27 @@ VolZmap::MarchingCubes( int nBlock, TRIA3DLIST& lstTria) const
//----------------------------------------------------------------------------
bool
VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST& lstTria, TriHolder& triHold) const
VolZmap::ExtMarchingCubes( int nBlock, TRIA3DLIST& lstTria, TriHolder& triHold) const
{
Point3d ptMapOrig = m_MapFrame[0].Orig() ;
if ( nBlock < 0 || nBlock >= int( m_nNumBlock))
return false ;
// Calcolo i limiti sugli indici dei voxel del blocco
int nIJK[3] ;
// Vettore indici i,j,k del blocco
GetBlockIJKFromN( nBlock, nIJK) ;
// Vettore limiti sugli indici dei voxel del blocco
int nLimits[6] ;
GetBlockLimitsIJK( nIJK, nLimits) ;
// Ciclo su tutti i voxel dello Zmap
for ( int i = nLimits[0] ; i < nLimits[1] ; ++ i) {
for ( int j = nLimits[2] ; j < nLimits[3] ; ++ j) {
for ( int k = nLimits[4] ; k < nLimits[5] ; ++ k) {
/*if ( !( i == 2 && j == 14 && k == 24))
continue ;*/
// Riconoscimento dei voxel di frontiera
bool bBoundary = false ;
int nVoxIndexes[3] = { i, j, k} ;
if ( IsAVoxelOnBoundary( nLimits, nVoxIndexes, true))
bBoundary = true ;
bool bBoundary = IsAVoxelOnBoundary( nLimits, nVoxIndexes, true) ;
// Indici i,j,k dei vertici
int IndexCorner[8][3] = {
@@ -804,9 +721,8 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
{ i, j + 1, k + 1}
} ;
int nIndex = 0 ;
// Classificazione dei vertici: interni o esterni al materiale
int nIndex = 0 ;
if ( IsThereMat( i, j, k))
nIndex |= ( 1 << 0) ;
if ( IsThereMat( i + 1, j, k))
@@ -1289,7 +1205,6 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
Point3d ptSol( dUnknownVector( 0) + vtO.x,
dUnknownVector( 1) + vtO.y,
dUnknownVector( 2) + vtO.z) ;
Triangle3d CurrentTriangle ;
TRIA3DVECTOR triContainer ;
@@ -1312,10 +1227,9 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
triContainer.emplace_back( CurrentTriangle) ;
// Controllo sull'inversione delle normali
if ( ( CurrentTriangle.GetN() * CompoVert[nj].vtNorm < - 0.01 &&
CurrentTriangle.GetN() * CompoVert[ni].vtNorm < - 0.01) &&
( ! bInvNormB)) {
if ( ( CurrentTriangle.GetN() * CompoVert[nj].vtNorm < - 0.01 && // - 0.01 passando a - 0.002 si tolgono i 4 restanti
CurrentTriangle.GetN() * CompoVert[ni].vtNorm < - 0.01) && // - 0.01 passando a - 0.002 si tolgono i 4 restanti
! bInvNormB) {
ptSol = ptGravityCenter ;
bInvNormB = true ;
triContainer.resize( 0) ;
@@ -1328,7 +1242,7 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
// interno o esterno al voxel a cui appartiene.
bool bInsideVoxel = IsPointInsideVoxelApprox( i, j, k, ptSol) ;
// Proprietà gemoetriche locali
// Proprietà geometriche locali
// del campo vettoriale: se vero
// procediamo con un ulteriore
// analisi per eliminare triangoli
@@ -1361,23 +1275,22 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
// esito negativo.
if ( nFeatureType == 2 && bLocProp &&
! ( bOutside || bInsideVoxel || bInvNormB)) {
if ( abs( nMin1 - nMin2) == 1 ||
abs( nMin1 - nMin2) == 3) {
abs( nMin1 - nMin2) == 3) {
int nSum = nMin1 + nMin2 ;
int nTriNum = ( nSum == 3 && ( nMin1 == 3 || nMin2 == 3) ?
max( nMin1, nMin2) : min( nMin1, nMin2)) ;
max( nMin1, nMin2) : min( nMin1, nMin2)) ;
double dDot1 = triContainer[nTriNum].GetN() * CompoVert[nMin1].vtNorm ;
double dDot2 = triContainer[nTriNum].GetN() * CompoVert[nMin2].vtNorm ;
if ( ( dDot1 < - 0.2 && dDot2 > - EPS_ZERO) ||
( dDot2 < - 0.2 && dDot1 > - EPS_ZERO)) {
( dDot2 < - 0.2 && dDot1 > - EPS_ZERO)) {
int nNm = dDot1 < - 0.2 ? nMin1 : nMin2 ;
int nNp = dDot1 < - 0.2 ? nMin2 : nMin1 ;
Vector3d vtVV = CompoVert[nNp].ptInt - CompoVert[nNm].ptInt ;
Point3d ptSolZmFrame = ptSol ;
@@ -1409,14 +1322,14 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
for ( int nc = 0 ; nc < nVertComp ; ++ nc) {
int nd = ( nc + 1 < nVertComp) ? nc + 1 : 0 ;
// Il triangolo è pronto
// Il triangolo è pronto
CurrentTriangle.Set( ptSol, CompoVert[nd].ptInt, CompoVert[nc].ptInt) ;
CurrentTriangle.Validate( true) ;
// Aggiungo triangolo al vettore temporaneo
// Aggiungo triangolo al vettore temporaneo
triContainer.emplace_back( CurrentTriangle) ;
}
}
}
}
}
else {
@@ -1431,17 +1344,17 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
if ( CompoVert[nPrev1].vtNorm * CompoVert[nMin1].vtNorm > 0) {
bool bTestOnVert = abs( nPrev1 - nMin2) == 0 || abs( nPrev1 - nMin2) == 3 ;
bool bTestOnVert = abs( nPrev1 - nMin2) == 1 || abs( nPrev1 - nMin2) == 3 ;
nNeighbourIndex = bTestOnVert ? nPrev1 : nNext1 ;
nStartIndex = nMin2 ;
}
else {
bool bTestOnVert = abs( nPrev2 - nMin1) == 0 || abs( nPrev2 - nMin1) == 3 ;
bool bTestOnVert = abs( nPrev2 - nMin1) == 1 || abs( nPrev2 - nMin1) == 3 ;
nNeighbourIndex = bTestOnVert ? nPrev2 : nNext2 ;
nStartIndex = nMin1 ;
}
Vector3d vtVV = CompoVert[nNeighbourIndex].ptInt - CompoVert[nStartIndex].ptInt ;
double dVVLen = vtVV.Len() ;
vtVV.Normalize() ;
@@ -1449,9 +1362,9 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
Vector3d vtVF = ptSol - CompoVert[nStartIndex].ptInt ;
Vector3d vtNewVF = ( vtVF * vtVV) * vtVV ;
double dNewVFLen = vtNewVF.Len() ;
double dNewVFLen = vtNewVF.Len() ;
if ( dNewVFLen > dVVLen) {
if ( dNewVFLen > dVVLen || vtVV * vtNewVF < 0) {
ptSol = CompoVert[nNeighbourIndex].ptInt ;
}
@@ -1459,23 +1372,21 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
ptSol = CompoVert[nStartIndex].ptInt + vtNewVF ;
}
triContainer.resize( 0) ;
for ( int nc = 0 ; nc < nVertComp ; ++ nc) {
int nd = ( nc + 1 < nVertComp) ? nc + 1 : 0 ;
// Il triangolo è pronto
// Il triangolo è pronto
CurrentTriangle.Set( ptSol, CompoVert[nd].ptInt, CompoVert[nc].ptInt) ;
CurrentTriangle.Validate( true) ;
// Aggiungo triangolo al vettore temporaneo
// Aggiungo triangolo al vettore temporaneo
triContainer.emplace_back( CurrentTriangle) ;
}
}
}
}
}
// Valuto normali: questo è ancora un controllo
// sulle normali, se risultano in tutti i punti
@@ -1512,7 +1423,7 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
// Se l'area è non nulla aggiungo il triangolo
// a uno dei due vettori.
if ( triContainer[ni].GetArea() > EPS_SMALL) {
if ( triContainer[ni].GetArea() > EPS_SMALL * EPS_SMALL) {
int nVoxIJK[3] = { i, j, k} ;
@@ -1558,7 +1469,6 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
// nel vector.
int nCurrent = int( triHold.size()) - 1 ;
// Aggiungo vertice della componente
// connessa alla lista dei vertici.
triHold[nCurrent].ptCompoVert.emplace_back( ptSol) ;
@@ -1586,36 +1496,36 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
// la struttura dati del voxel.
if ( nBorderFeatureInVoxel == 1) {
m_InterBlockTria[nBlockNumber].resize( m_InterBlockTria[nBlockNumber].size() + 1) ;
m_InterBlockTria[nBlock].resize( m_InterBlockTria[nBlock].size() + 1) ;
// Questi dati dipendono solo dal voxel,
// quindi sono validi per tutte le
// componenti che vi appartengono.
int nCurrent = int( m_InterBlockTria[nBlockNumber].size()) - 1 ;
int nCurrent = int( m_InterBlockTria[nBlock].size()) - 1 ;
m_InterBlockTria[nBlockNumber][nCurrent].i = i ;
m_InterBlockTria[nBlockNumber][nCurrent].j = j ;
m_InterBlockTria[nBlockNumber][nCurrent].k = k ;
m_InterBlockTria[nBlock][nCurrent].i = i ;
m_InterBlockTria[nBlock][nCurrent].j = j ;
m_InterBlockTria[nBlock][nCurrent].k = k ;
}
// Indice che identifica la posizione del voxel
// nel vector.
int nCurrent = int( m_InterBlockTria[nBlockNumber].size()) - 1 ;
int nCurrent = int( m_InterBlockTria[nBlock].size()) - 1 ;
// Aggiungo vertice della componente
// connessa alla lista dei vertici.
m_InterBlockTria[nBlockNumber][nCurrent].ptCompoVert.emplace_back( ptSol) ;
m_InterBlockTria[nBlock][nCurrent].ptCompoVert.emplace_back( ptSol) ;
int nOldFeatureNum = int( m_InterBlockTria[nBlockNumber][nCurrent].vCompoTria.size()) ;
int nOldFeatureNum = int( m_InterBlockTria[nBlock][nCurrent].vCompoTria.size()) ;
int nNewFeatureNum = nOldFeatureNum + 1 ;
// Aggiungo una componente per il vettore
// dei triangoli della componente connessa.
m_InterBlockTria[nBlockNumber][nCurrent].vCompoTria.resize( nNewFeatureNum) ;
m_InterBlockTria[nBlock][nCurrent].vCompoTria.resize( nNewFeatureNum) ;
for ( int ni = 0 ; ni < int( vBorderTriaTemp.size()) ; ++ ni)
m_InterBlockTria[nBlockNumber][nCurrent].vCompoTria[nOldFeatureNum].emplace_back( vBorderTriaTemp[ni]) ;
m_InterBlockTria[nBlock][nCurrent].vCompoTria[nOldFeatureNum].emplace_back( vBorderTriaTemp[ni]) ;
}
}
// ExtMC non confermato, si passa a MC
@@ -1634,9 +1544,8 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
CurrentTriangle.Set( VecField[i0].ptInt, VecField[i1].ptInt, VecField[i2].ptInt) ;
CurrentTriangle.Validate( true) ;
// Se il triangolo non è degenere lo aggiungo alla lista
if ( CurrentTriangle.GetArea() > EPS_SMALL)
lstTria.emplace_back( CurrentTriangle) ;
// Aggiungo alla lista
lstTria.emplace_back( CurrentTriangle) ;
}
}
}
@@ -1657,9 +1566,8 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
CurrentTriangle.Set( VecField[i0].ptInt, VecField[i1].ptInt, VecField[i2].ptInt) ;
CurrentTriangle.Validate( true) ;
// Se il triangolo non è degenere lo aggiungo alla lista
if ( CurrentTriangle.GetArea() > EPS_SMALL)
lstTria.emplace_back( CurrentTriangle) ;
// Aggiungo alla lista
lstTria.emplace_back( CurrentTriangle) ;
}
}
nTableOffset += nVertComp ;
@@ -1675,28 +1583,18 @@ VolZmap::ExtMarchingCubes( const int nLimits[], size_t nBlockNumber, TRIA3DLIST&
bool
VolZmap::FlipEdgesII( TriHolder& TriHold) const
{
double dSqEps = EPS_SMALL * EPS_SMALL ;
// Numero di voxel in cui si presentano sharp feature
int nVoxelNum = int( TriHold.size()) ;
// Determino estremi del ciclo sui voxel esterno
int nSt1 = 0 ;
int nEn1 = nVoxelNum - 1 ;
// Ciclo su tali voxel
for ( int n1 = nSt1 ; n1 < nEn1 ; ++ n1) {
for ( int n1 = 0 ; n1 < nVoxelNum ; ++ n1) {
// Determino estremi del ciclo sui voxel interno
int nSt2 = nSt1 + 1 ;
int nEn2 = nVoxelNum ;
for ( int n2 = nSt2 ; n2 < nEn2 ; ++ n2) {
for ( int n2 = n1 ; n2 < nVoxelNum ; ++ n2) {
// Se i voxel sono adiacenti proseguo
if ( TriHold[n2].i == TriHold[n1].i + 1 ||
TriHold[n2].j == TriHold[n1].j + 1 ||
TriHold[n2].k == TriHold[n1].k + 1) {
if ( ( TriHold[n2].i >= TriHold[n1].i && TriHold[n2].i <= TriHold[n1].i + 1) ||
( TriHold[n2].j >= TriHold[n1].j && TriHold[n2].j <= TriHold[n1].j + 1) ||
( TriHold[n2].k >= TriHold[n1].k && TriHold[n2].k <= TriHold[n1].k + 1)) {
// Numero delle componenti connesse nei due voxel
int nNumCompo1 = int( TriHold[n1].ptCompoVert.size()) ;
@@ -1707,7 +1605,7 @@ VolZmap::FlipEdgesII( TriHolder& TriHold) const
// Ciclo sulle componenti
for ( ; nCompo1 < nNumCompo1 ; ++ nCompo1) {
int nCompo2 = 0 ;
int nCompo2 = ( n1 == n2 ? nCompo1 + 1 : 0) ;
for ( ; nCompo2 < nNumCompo2 ; ++ nCompo2) {
@@ -1730,13 +1628,13 @@ VolZmap::FlipEdgesII( TriHolder& TriHold) const
Point3d ptP1 = TriHold[n1].vCompoTria[nCompo1][nTri1].GetP( nVert1) ;
Point3d ptP2 = TriHold[n2].vCompoTria[nCompo2][nTri2].GetP( nVert2) ;
if ( SqDist( ptP1, ptP2) < dSqEps) {
if ( AreSamePointEpsilon( ptP1, ptP2, EPS_SMALL)) {
Point3d ptVert1 = TriHold[n1].ptCompoVert[nCompo1] ;
Point3d ptVert2 = TriHold[n2].ptCompoVert[nCompo2] ;
if ( ! ( AreSamePointApprox( ptP1, ptVert1) ||
AreSamePointApprox( ptP2, ptVert2))) {
if ( ! ( AreSamePointEpsilon( ptP1, ptVert1, EPS_SMALL) ||
AreSamePointEpsilon( ptP2, ptVert2, EPS_SMALL))) {
SharedIndex.emplace_back( nVert1) ;
SharedIndex.emplace_back( nVert2) ;
@@ -1753,20 +1651,10 @@ VolZmap::FlipEdgesII( TriHolder& TriHold) const
// Si deve operare la modifica dei triangoli
if ( SharedIndex.size() > 2) {
// verifico che i due lai adiacenti siano controversi
int nProd = ( SharedIndex[2] - SharedIndex[0]) * ( SharedIndex[3] - SharedIndex[1]) ;
if ( nProd != 1 && nProd != - 2 && nProd != 4) {
Triangle3d TriVox1 = TriHold[n1].vCompoTria[nCompo1][nTri1] ;
Triangle3d TriVox2 = TriHold[n2].vCompoTria[nCompo2][nTri2] ;
TriVox1.SetP( SharedIndex[0], TriHold[n2].ptCompoVert[nCompo2]) ;
TriVox2.SetP( SharedIndex[3], TriHold[n1].ptCompoVert[nCompo1]) ;
TriVox1.Validate( true) ;
TriVox2.Validate( true) ;
TriHold[n1].vCompoTria[nCompo1][nTri1].SetP( SharedIndex[0],
TriHold[n2].ptCompoVert[nCompo2]) ;
@@ -1796,28 +1684,26 @@ VolZmap::FlipEdgesII( TriHolder& TriHold) const
//----------------------------------------------------------------------------
bool
VolZmap::FlipEdgesBB( TriaMatrix& InterTria) const {
double dSqEps = EPS_SMALL * EPS_SMALL ;
VolZmap::FlipEdgesBB( TriaMatrix& InterTria) const
{
// Numero di blocchi
size_t nBlocksNum = InterTria.size() ;
// ciclo sui blocchi
for ( size_t tFB = 0 ; tFB < nBlocksNum - 1 ; ++ tFB) {
int nFBijk[3] ;
GetBlockIJKFromN( int( tFB), nFBijk) ;
for ( size_t tLB = tFB + 1 ; tLB < nBlocksNum ; ++ tLB) {
int nFBijk[3] ;
int nLBijk[3] ;
GetBlockIJKFromN( int( tFB), nFBijk) ;
GetBlockIJKFromN( int( tLB), nLBijk) ;
// Se i blocchi non sono adiacenti salto l'iterazione
if ( ! ( abs( nFBijk[0] - nLBijk[0]) <= 1 &&
abs( nFBijk[1] - nLBijk[1]) <= 1 &&
abs( nFBijk[2] - nLBijk[2]) <= 1))
continue ;
// Numero di voxel nei blocchi correnti
@@ -1861,13 +1747,13 @@ VolZmap::FlipEdgesBB( TriaMatrix& InterTria) const {
Point3d ptPF = InterTria[tFB][tVFB].vCompoTria[tCmpF][tTriFB].GetP( nVertF) ;
Point3d ptPL = InterTria[tLB][tVLB].vCompoTria[tCmpL][tTriLB].GetP( nVertL) ;
if ( SqDist( ptPF, ptPL) < dSqEps) {
if ( AreSamePointEpsilon( ptPF, ptPL, EPS_SMALL)) {
Point3d ptVertF = InterTria[tFB][tVFB].ptCompoVert[tCmpF] ;
Point3d ptVertL = InterTria[tLB][tVLB].ptCompoVert[tCmpL] ;
if ( ! ( AreSamePointApprox( ptPF, ptVertF) ||
AreSamePointApprox( ptPL, ptVertL))) {
if ( ! ( AreSamePointEpsilon( ptPF, ptVertF, EPS_SMALL) ||
AreSamePointEpsilon( ptPL, ptVertL, EPS_SMALL))) {
SharedIndex.emplace_back( nVertF) ;
SharedIndex.emplace_back( nVertL) ;
@@ -1885,21 +1771,10 @@ VolZmap::FlipEdgesBB( TriaMatrix& InterTria) const {
// Si deve operare la modifica dei triangoli
if ( SharedIndex.size() > 2) {
// verifico che i due lai adiacenti siano controversi
int nProd = ( SharedIndex[2] - SharedIndex[0]) * ( SharedIndex[3] - SharedIndex[1]) ;
if ( nProd != 1 && nProd != - 2 && nProd != 4) {
Triangle3d TriVoxF = InterTria[tFB][tVFB].vCompoTria[tCmpF][tTriFB] ;
Triangle3d TriVoxL = InterTria[tLB][tVLB].vCompoTria[tCmpL][tTriLB] ;
TriVoxF.SetP( SharedIndex[0], InterTria[tLB][tVLB].ptCompoVert[tCmpL]) ;
TriVoxL.SetP( SharedIndex[3], InterTria[tFB][tVFB].ptCompoVert[tCmpF]) ;
TriVoxF.Validate( true) ;
TriVoxL.Validate( true) ;
InterTria[tFB][tVFB].vCompoTria[tCmpF][tTriFB].SetP( SharedIndex[0],
InterTria[tLB][tVLB].ptCompoVert[tCmpL]) ;
InterTria[tLB][tVLB].vCompoTria[tCmpL][tTriLB].SetP( SharedIndex[3],
@@ -1910,7 +1785,6 @@ VolZmap::FlipEdgesBB( TriaMatrix& InterTria) const {
bModified = true ;
break ;
}
}
}