Egalware/EgtGeomKernel
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

EgtGeomKernel 1.9b4 :

Browse Source 
- modifiche per visualizzazione Zmap.
This commit is contained in:
Dario Sassi
2018-02-16 17:45:59 +00:00
parent 01886cc428
commit 189a942595
2 changed files with 176 additions and 100 deletions
Download Patch File Download Diff File Expand all files Collapse all files
VolZmapGraphics.cpp
+176 -100
View File
@@ -1648,7 +1648,7 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
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) {
// Classificazione dei vertici: interni o esterni al materiale
int nIndex = CalcIndex( i, j, k) ;
@@ -1857,21 +1857,39 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
else if ( nAllConfig[nIndex] == 6) {
// Procedura analoga a quella della configurazione 3
Vector3d vtCmpAvg0, vtCmpAvg1 ;
bool bTest0 = DotTest( CompoVert[0], 4, vtCmpAvg0, 0.7) ;
/*Vector3d vtCmpAvg0, vtCmpAvg1 ;
bool bTest0 = DotTest( CompoVert[0], 4, vtCmpAvg0, 0.7) ;
bool bTest1 = DotTest( CompoVert[1], 3, vtCmpAvg1, 0.7) ;
double dScProd = - 2 ;
if ( bTest0 && bTest1)
dScProd = vtCmpAvg0 * vtCmpAvg1 ;
bool bFirstTest = ! ( bTest0 && bTest1) || ( bTest0 && bTest1 && dScProd > 0.7) ;*/
// Altro metodo
double dDotSum = 0 ;
for ( int nFV = 0 ; nFV < 4 ; ++ nFV) {
for ( int nTV = 0 ; nTV < 3 ; ++ nTV) {
dDotSum += CompoVert[0][nFV].vtNorm * CompoVert[1][nTV].vtNorm ;
}
}
for ( int nFVI = 0 ; nFVI < 3 ; ++ nFVI) {
for ( int nFVJ = nFVI + 1 ; nFVJ < 4 ; ++ nFVJ) {
dDotSum -= CompoVert[0][nFVI].vtNorm * CompoVert[0][nFVJ].vtNorm ;
}
}
for ( int nTVI = 0 ; nTVI < 2 ; ++ nTVI) {
for ( int nTVJ = nTVI + 1 ; nTVJ < 3 ; ++ nTVJ) {
dDotSum -= CompoVert[1][nTVI].vtNorm * CompoVert[1][nTVJ].vtNorm ;
}
}
bool bTestOnSum = dDotSum > - 4 ;
if ( ! ( bTest0 && bTest1) || ( bTest0 && bTest1 && dScProd > 0.7)) {
if ( bTestOnSum) {
int nt = 0 ;
while ( nIndexVsIndex6[nt][0] != nIndex)
++ nt ;
int nRotCase = nIndexVsIndex6[nt][1] ;
// Costruzione dei triangoli
for ( int TriIndex = 0 ; TriIndex < 15 ; TriIndex += 3) {
// Indici vertici
@@ -1885,10 +1903,37 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
CurrentTriangle.SetVertexNorm( 0, VecField[i0].vtNorm) ;
CurrentTriangle.SetVertexNorm( 1, VecField[i1].vtNorm) ;
CurrentTriangle.SetVertexNorm( 2, VecField[i2].vtNorm) ;
// Setto il numero di utensile (conta solo positivo, nullo o negativo)
int nTool0 = Clamp( VecField[i0].nToolFlag, -1, 1) ;
int nTool1 = Clamp( VecField[i1].nToolFlag, -1, 1) ;
int nTool2 = Clamp( VecField[i2].nToolFlag, -1, 1) ;
if ( nTool0 == nTool1 || nTool0 == nTool2)
CurrentTriangle.SetGrade( nTool0) ;
else if ( nTool1 == nTool2)
CurrentTriangle.SetGrade( nTool1) ;
// Valido il triangolo e setto le normali del campo vettoriale ai corrispondenti vertici
if ( CurrentTriangle.Validate( true)) {
Vector3d vtVertNorm = VecField[i0].vtNorm ;
if ( CurrentTriangle.GetN() * vtVertNorm > 0.6)
CurrentTriangle.SetVertexNorm( 0, vtVertNorm) ;
else
CurrentTriangle.SetVertexNorm( 0, CurrentTriangle.GetN()) ;
vtVertNorm = VecField[i1].vtNorm ;
if ( CurrentTriangle.GetN() * vtVertNorm > 0.6)
CurrentTriangle.SetVertexNorm( 1, vtVertNorm) ;
else
CurrentTriangle.SetVertexNorm( 1, CurrentTriangle.GetN()) ;
vtVertNorm = VecField[i2].vtNorm ;
if ( CurrentTriangle.GetN() * vtVertNorm > 0.6)
CurrentTriangle.SetVertexNorm( 2, vtVertNorm) ;
else
CurrentTriangle.SetVertexNorm( 2, CurrentTriangle.GetN()) ;
}
// Riporto il triangolo nel sistema in cui è immerso quello dello Zmap
CurrentTriangle.ToGlob( m_MapFrame) ;
// Aggiungo alla lista
// Aggiungo alla lista
lstTria.emplace_back( CurrentTriangle) ;
}
}
continue ;
}
}
@@ -2079,26 +2124,20 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
// Caso ventaglio con tre vertici di base
if ( nVertComp[nCompCount - 1] == 3) {
// Controllo se esiste almeno un triangolo con normale avente prodotto scalare
// negativo con la normale di almeno uno dei vertici di base del ventaglio.
bool bInversione = false ;
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] ; ++ ni) {
int nj = ( ni + 1 < nVertComp[nCompCount - 1]) ? ni + 1 : 0 ;
double dDI = triContainer[ni].GetN() * CompoVert[nCompCount - 1][ni].vtNorm ;
double dDJ = triContainer[ni].GetN() * CompoVert[nCompCount - 1][nj].vtNorm ;
if ( dDI < - EPS_SMALL || dDJ < - EPS_SMALL) {
bInversione = true ;
break ;
}
}
// Se tale triangolo esiste procedo
if ( bInversione) {
// Conto le coppie di normali con angolo compreso maggiore di 90°
int nNegDot = 0 ;
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] - 1 ; ++ ni) {
@@ -2107,7 +2146,6 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
nNegDot ++ ;
}
}
if ( nNegDot == nVertComp[nCompCount - 1] - 1) {
// Cerco se esiste un punto in cui la normale ha prodotto scalre negativo
// con le normali di entrambi i triangoli che lo contengono
@@ -2122,14 +2160,11 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
break ;
}
}
// Se tale normale esiste
if ( bInversione2) {
// Se la soluzione non cade nel voxel di appartenenza vedo se può
// essere riportata dentro muovendosi lungo la linea di feature.
if ( ! IsPointInsideVoxelApprox( i, j, k, ptSol)) {
Vector3d vtNullSpace( dMatrixV( 0, 2), dMatrixV( 1, 2), dMatrixV( 2, 2)) ;
double dParInt1, dParInt2 ;
Point3d ptVoxMin( ( i * N_DEXVOXRATIO + 0.5) * m_dStep,
@@ -2139,14 +2174,10 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
// Caso in cui può essere riportata dentro: se il voxel in cui cade la feature è pieno
// la riporto muovendola lungo la sua linea
if ( IntersLineBox( ptSol, vtNullSpace, ptVoxMin, ptVoxMax, dParInt1, dParInt2)) {
triContainer.resize( 0) ;
double dPar = abs( dParInt1) < abs( dParInt2) ? dParInt1 + ( dParInt2 - dParInt1) / 100 :
dParInt2 + ( dParInt1 - dParInt2) / 100 ;
ptSol += dPar * vtNullSpace ;
// Costruisco triangoli di prova
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] ; ++ ni) {
int nj = ( ni + 1 < nVertComp[nCompCount - 1]) ? ni + 1 : 0 ;
@@ -2175,18 +2206,14 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
triContainer.emplace_back( CurrentTriangle) ;
}
}
// Caso in cui non può essere riportata dentro:
// si passa alla routine standard se il voxel
// in cui cade è pieno.
else {
int nAdjVoxI, nAdjVoxJ, nAdjVoxK ;
if ( GetPointVoxel( ptSol, nAdjVoxI, nAdjVoxJ, nAdjVoxK)) {
// Classificazione del voxel adiacente
int nAdjIndex = CalcIndex( nAdjVoxI, nAdjVoxJ, nAdjVoxK) ;
// Se il voxel è pieno
if ( EdgeTable[nAdjIndex] != 0)
bDangerInversion = true ;
@@ -2197,84 +2224,133 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
}
}
}
// Ventaglio con base a quattro vertici
// Ventaglio con base a quattro vertici
else if ( nVertComp[nCompCount - 1] == 4) {
// Controllo se esiste almeno un triangolo con normale avente prodotto scalare
// negativo con la normale di almeno uno dei vertici di base del ventaglio.
// Controllo preliminare sulle normali
bool bInversione = false ;
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] ; ++ ni) {
int nj = ( ni + 1 < nVertComp[nCompCount - 1]) ? ni + 1 : 0 ;
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] ; ++ ni) {
int nj = ( ni + 1 < nVertComp[nCompCount - 1]) ? ni + 1 : 0 ;
double dDI = triContainer[ni].GetN() * CompoVert[nCompCount - 1][ni].vtNorm ;
double dDJ = triContainer[ni].GetN() * CompoVert[nCompCount - 1][nj].vtNorm ;
if ( dDI < - EPS_SMALL || dDJ < - EPS_SMALL)
bInversione = true ;
double dDJ = triContainer[ni].GetN() * CompoVert[nCompCount - 1][nj].vtNorm ;
if ( dDI < - EPS_SMALL || dDJ < - EPS_SMALL) {
bInversione = true ;
break ;
}
}
// Se tale triangolo esiste continuo i test
// Si necessita di un ulteriore controllo
if ( bInversione) {
// Se la feature non cade nel suo voxel
if ( ! IsPointInsideVoxelApprox( i, j, k, ptSol)) {
Vector3d vtNullSpace( dMatrixV( 0, 2), dMatrixV( 1, 2), dMatrixV( 2, 2)) ;
double dParInt1, dParInt2 ;
Point3d ptVoxMin( ( i * N_DEXVOXRATIO + 0.5) * m_dStep,
( j * N_DEXVOXRATIO + 0.5) * m_dStep, ( k * N_DEXVOXRATIO + 0.5) * m_dStep) ;
Point3d ptVoxMax( ( ( i + 1) * N_DEXVOXRATIO + 0.5) * m_dStep,
( ( j + 1) * N_DEXVOXRATIO + 0.5) * m_dStep, ( ( k + 1) * N_DEXVOXRATIO + 0.5) * m_dStep) ;
// Se è possibile riportarla dentro e il voxel in cui cade è pieno, la riporto nel suo voxel
// lungo la sua linea
if ( IntersLineBox( ptSol, vtNullSpace, ptVoxMin, ptVoxMax, dParInt1, dParInt2)) {
triContainer.resize( 0) ;
double dPar = abs( dParInt1) < abs( dParInt2) ? dParInt1 + ( dParInt2 - dParInt1) / 100 :
dParInt2 + ( dParInt1 - dParInt2) / 100 ;
ptSol += dPar * vtNullSpace ;
// Costruisco triangoli di prova
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] ; ++ ni) {
int nj = ( ni + 1 < nVertComp[nCompCount - 1]) ? ni + 1 : 0 ;
// Il triangolo è pronto
Triangle3dEx CurrentTriangle ;
CurrentTriangle.Set( ptSol, CompoVert[nCompCount - 1][nj].ptInt, CompoVert[nCompCount - 1][ni].ptInt) ;
CurrentTriangle.SetAttrib( 1, CompoVert[nCompCount - 1][nj].nToolFlag) ;
CurrentTriangle.SetAttrib( 2, CompoVert[nCompCount - 1][ni].nToolFlag) ;
CurrentTriangle.SetVertexNorm( 1, CompoVert[nCompCount - 1][nj].vtNorm) ;
CurrentTriangle.SetVertexNorm( 2, CompoVert[nCompCount - 1][ni].vtNorm) ;
if ( CurrentTriangle.GetVertexNorm( 1) * CurrentTriangle.GetVertexNorm( 2) > 0.5)
CurrentTriangle.SetVertexNorm( 0, 0.5 * ( CurrentTriangle.GetVertexNorm( 1) +
CurrentTriangle.GetVertexNorm( 2))) ;
CurrentTriangle.Validate( true) ;
// Assegno i colori
if ( CurrentTriangle.GetAttrib( 1) < 0 ||
CurrentTriangle.GetAttrib( 2) < 0)
CurrentTriangle.SetGrade( - 1) ;
else if ( CurrentTriangle.GetAttrib( 1) > 0 ||
CurrentTriangle.GetAttrib( 2) > 0)
CurrentTriangle.SetGrade( 1) ;
else
CurrentTriangle.SetGrade( 0) ;
// Aggiungo triangolo al vettore temporaneo
triContainer.emplace_back( CurrentTriangle) ;
}
PolygonPlane Polygon ;
Point3d ptP0 = CompoVert[nCompCount - 1][0].ptInt ;
Point3d ptP1 = CompoVert[nCompCount - 1][1].ptInt ;
Point3d ptP2 = CompoVert[nCompCount - 1][2].ptInt ;
Point3d ptP3 = CompoVert[nCompCount - 1][3].ptInt ;
Polygon.AddPoint( ptP0) ;
Polygon.AddPoint( ptP1) ;
Polygon.AddPoint( ptP2) ;
Polygon.AddPoint( ptP3) ;
Point3d ptCentr ;
Vector3d vtPolyN ;
bool bCentFound = Polygon.GetCentroid( ptCentr) ;
bool bNormFound = Polygon.GetNormal( vtPolyN) ;
if ( bCentFound && bNormFound) {
Point3d ptPoly0 = ptP0 + ( ptP0 - ptCentr) - ( ptP0 - ptCentr) * vtPolyN * vtPolyN ;
Point3d ptPoly1 = ptP1 + ( ptP1 - ptCentr) - ( ptP1 - ptCentr) * vtPolyN * vtPolyN ;
Point3d ptPoly2 = ptP2 + ( ptP2 - ptCentr) - ( ptP2 - ptCentr) * vtPolyN * vtPolyN ;
Point3d ptPoly3 = ptP3 + ( ptP3 - ptCentr) - ( ptP3 - ptCentr) * vtPolyN * vtPolyN ;
Point3d ptSolP = ptSol + ( ptSol - ptCentr) - ( ptSol - ptCentr) * vtPolyN * vtPolyN ;
Vector3d vtSeg0 = ptPoly1 - ptPoly0 ;
Vector3d vtSeg1 = ptPoly2 - ptPoly1 ;
Vector3d vtSeg2 = ptPoly3 - ptPoly2 ;
Vector3d vtSeg3 = ptPoly0 - ptPoly3 ;
Vector3d vtSol0 = ptSolP - ptPoly0 ;
Vector3d vtSol1 = ptSolP - ptPoly1 ;
Vector3d vtSol2 = ptSolP - ptPoly2 ;
Vector3d vtSol3 = ptSolP - ptPoly3 ;
vtSeg0.Normalize() ;
vtSeg1.Normalize() ;
vtSeg2.Normalize() ;
vtSeg3.Normalize() ;
/*vtSol0.Normalize() ;
vtSol1.Normalize() ;
vtSol2.Normalize() ;
vtSol3.Normalize() ;*/
double dIn[4] ;
dIn[0] = ( vtSeg0 ^ vtSol0) * vtPolyN ;
dIn[1] = ( vtSeg1 ^ vtSol1) * vtPolyN ;
dIn[2] = ( vtSeg2 ^ vtSol2) * vtPolyN ;
dIn[3] = ( vtSeg3 ^ vtSol3) * vtPolyN ;
double dThr = 0.1 * N_DEXVOXRATIO * m_dStep ;
bool bSpecial = false ;
for ( int nE = 0 ; nE < 4 ; ++ nE) {
if ( dIn[nE] < dThr) {
int nF = nE + 1 < 4 ? nE + 1 : 0 ;
double dDE = triContainer[nE].GetN() * CompoVert[nCompCount - 1][nE].vtNorm ;
double dDF = triContainer[nE].GetN() * CompoVert[nCompCount - 1][nF].vtNorm ;
double dLim = - 0.5/*EPS_SMALL*/ ;
if ( dDE < dLim || dDF < dLim) {
bSpecial = true ;
break ;
}
}
}
// Se non è possibile riportarla dentro e il voxel in
// cui cade è pieno passo alla routine standard
else {
int nAdjVoxI, nAdjVoxJ, nAdjVoxK ;
if ( GetPointVoxel( ptSol, nAdjVoxI, nAdjVoxJ, nAdjVoxK)) {
// Classificazione del voxel adiacente
int nAdjIndex = CalcIndex( nAdjVoxI, nAdjVoxJ, nAdjVoxK) ;
// Se il voxel è pieno
if ( EdgeTable[nAdjIndex] != 0)
bDangerInversion = true ;
}
if ( bSpecial) {
// Se la feature non cade nel suo voxel
if ( ! IsPointInsideVoxelApprox( i, j, k, ptSol)) {
// Vettore di direzione feature
Vector3d vtNullSpace( dMatrixV( 0, 2), dMatrixV( 1, 2), dMatrixV( 2, 2)) ;
double dParInt1, dParInt2 ;
Point3d ptVoxMin( ( i * N_DEXVOXRATIO + 0.5) * m_dStep,
( j * N_DEXVOXRATIO + 0.5) * m_dStep, ( k * N_DEXVOXRATIO + 0.5) * m_dStep) ;
Point3d ptVoxMax( ( ( i + 1) * N_DEXVOXRATIO + 0.5) * m_dStep,
( ( j + 1) * N_DEXVOXRATIO + 0.5) * m_dStep, ( ( k + 1) * N_DEXVOXRATIO + 0.5) * m_dStep) ;
// Se è possibile riportarla dentro e il voxel in cui cade è pieno, la riporto nel suo voxel
// lungo la sua linea
if ( IntersLineBox( ptSol, vtNullSpace, ptVoxMin, ptVoxMax, dParInt1, dParInt2)) {
triContainer.resize( 0) ;
double dPar = abs( dParInt1) < abs( dParInt2) ? dParInt1 + ( dParInt2 - dParInt1) / 100 :
dParInt2 + ( dParInt1 - dParInt2) / 100 ;
ptSol += dPar * vtNullSpace ;
// Costruisco triangoli di prova
for ( int ni = 0 ; ni < nVertComp[nCompCount - 1] ; ++ ni) {
int nj = ( ni + 1 < nVertComp[nCompCount - 1]) ? ni + 1 : 0 ;
// Il triangolo è pronto
Triangle3dEx CurrentTriangle ;
CurrentTriangle.Set( ptSol, CompoVert[nCompCount - 1][nj].ptInt, CompoVert[nCompCount - 1][ni].ptInt) ;
CurrentTriangle.SetAttrib( 1, CompoVert[nCompCount - 1][nj].nToolFlag) ;
CurrentTriangle.SetAttrib( 2, CompoVert[nCompCount - 1][ni].nToolFlag) ;
CurrentTriangle.SetVertexNorm( 1, CompoVert[nCompCount - 1][nj].vtNorm) ;
CurrentTriangle.SetVertexNorm( 2, CompoVert[nCompCount - 1][ni].vtNorm) ;
if ( CurrentTriangle.GetVertexNorm( 1) * CurrentTriangle.GetVertexNorm( 2) > 0.5)
CurrentTriangle.SetVertexNorm( 0, 0.5 * ( CurrentTriangle.GetVertexNorm( 1) +
CurrentTriangle.GetVertexNorm( 2))) ;
CurrentTriangle.Validate( true) ;
// Assegno i colori
if ( CurrentTriangle.GetAttrib( 1) < 0 ||
CurrentTriangle.GetAttrib( 2) < 0)
CurrentTriangle.SetGrade( - 1) ;
else if ( CurrentTriangle.GetAttrib( 1) > 0 ||
CurrentTriangle.GetAttrib( 2) > 0)
CurrentTriangle.SetGrade( 1) ;
else
CurrentTriangle.SetGrade( 0) ;
// Aggiungo triangolo al vettore temporaneo
triContainer.emplace_back( CurrentTriangle) ;
}
}
// Se non è possibile riportarla dentro e il voxel in
// cui cade è pieno passo alla routine standard
else {
int nAdjVoxI, nAdjVoxJ, nAdjVoxK ;
if ( GetPointVoxel( ptSol, nAdjVoxI, nAdjVoxJ, nAdjVoxK)) {
// Classificazione del voxel adiacente
int nAdjIndex = CalcIndex( nAdjVoxI, nAdjVoxJ, nAdjVoxK) ;
// Se il voxel è pieno
if ( EdgeTable[nAdjIndex] != 0)
bDangerInversion = true ;
}
}
}
}
}
}
@@ -2408,7 +2484,7 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
}
// ExtMC non confermato, si passa a MC
else {
else {
vector<Triangle3dEx> vTria ;
// Costruzione dei triangoli
for ( int TriIndex = 0; TriIndex < ( nVertComp[nCompCount - 1] - 2) * 3 ; TriIndex += 3) {
@@ -2459,7 +2535,7 @@ VolZmap::ExtMarchingCubes( int nBlock, TRIA3DEXLIST& lstTria, TriHolder& triHold
}
// Standard MC
else {
else {
vector<Triangle3dEx> vTria ;
// Costruzione dei triangoli
for ( int TriIndex = 0; TriIndex < ( nVertComp[nCompCount - 1] - 2) * 3 ; TriIndex += 3) {