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

EgtGeomKernel :

Browse Source 
- migliorie alla costruzione degli Zmap
- migliorie a GetDepth di Zmap.
- sistemato calcolo BBox esatto di curva Bezier.
This commit is contained in:
Dario Sassi
2017-08-16 16:33:59 +00:00
parent c17e44a3e1
commit 9da43cfd66
5 changed files with 378 additions and 766 deletions
Download Patch File Download Diff File Expand all files Collapse all files
VolTriZmapCalculus.cpp
+111 -266
View File
@@ -57,7 +57,7 @@ VolZmap::IntersLineBox( const Point3d& ptP, const Vector3d& vtV,
else if ( ptP.y < ptMin.y - EPS_SMALL || ptP.y > ptMax.y + EPS_SMALL)
return false ;
// confronto con piani XZ (perpendicolari ad asse Z)
// confronto con piani XY (perpendicolari ad asse Z)
if ( vtV.z > EPS_ZERO) {
dU1 = max( dU1, ( ptMin.z - ptP.z) / vtV.z) ;
dU2 = min( dU2, ( ptMax.z - ptP.z) / vtV.z) ;
@@ -72,205 +72,23 @@ VolZmap::IntersLineBox( const Point3d& ptP, const Vector3d& vtV,
return ( dU2 >= dU1) ;
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineVoxel( const Point3d& ptP, const Vector3d& vtV, int nIndI, int nIndJ, int nIndK,
int& nFaceF, int& nFaceL, double& dUF, double& dUL) const
{
// Controllo sull'ammissibilità del voxel
if ( nIndI < - 1 || nIndI >= int( m_nNx[0]) ||
nIndJ < - 1 || nIndJ >= int( m_nNy[0]) ||
nIndK < - 1 || nIndK >= int( m_nNy[1]))
return false ;
Point3d ptInt, ptIntF, ptIntL ;
double dSqEps = EPS_SMALL * EPS_SMALL ;
int nIntNum = 0 ;
double dU1 = ( ( nIndJ + 0.5) * m_dStep - ptP.y) / vtV.y ;
double dU2 = ( ( nIndI + 0.5) * m_dStep - ptP.x) / vtV.x ;
double dU3 = ( ( nIndJ + 1.5) * m_dStep - ptP.y) / vtV.y ;
double dU4 = ( ( nIndI + 1.5) * m_dStep - ptP.x) / vtV.x ;
double dU5 = ( ( nIndK + 0.5) * m_dStep - ptP.z) / vtV.z ;
double dU6 = ( ( nIndK + 1.5) * m_dStep - ptP.z) / vtV.z ;
// Intersezione con le facce 1 e 3
if ( abs( vtV.y) > EPS_ZERO) {
// Intersezione con la prima faccia
ptInt = ptP + dU1 * vtV ;
if ( ptInt.x >= ( nIndI + 0.5) * m_dStep &&
ptInt.x <= ( nIndI + 1.5) * m_dStep &&
ptInt.z >= ( nIndK + 0.5) * m_dStep &&
ptInt.z <= ( nIndK + 1.5) * m_dStep) {
dUF = dU1 ;
nFaceF = 1 ;
ptIntF = ptInt ;
++ nIntNum ;
}
// Intersezione con la terza faccia
ptInt = ptP + dU3 * vtV ;
if ( ptInt.x >= ( nIndI + 0.5) * m_dStep &&
ptInt.x <= ( nIndI + 1.5) * m_dStep &&
ptInt.z >= ( nIndK + 0.5) * m_dStep &&
ptInt.z <= ( nIndK + 1.5) * m_dStep) {
if ( nIntNum == 0) {
dUF = dU3 ;
nFaceF = 3 ;
ptIntF = ptInt ;
++ nIntNum ;
}
else if ( nIntNum == 1 &&
( ptIntF - ptInt).SqLen() > dSqEps) {
dUL = dU3 ;
nFaceL = 3 ;
ptIntL = ptInt ;
++ nIntNum ;
}
}
}
// Intersezione con le facce 2 e 4
if ( abs( vtV.x) > EPS_ZERO) {
// Intersezione con la seconda faccia
ptInt = ptP + dU2 * vtV ;
if ( ptInt.y >= ( nIndJ + 0.5) * m_dStep &&
ptInt.y <= ( nIndJ + 1.5) * m_dStep &&
ptInt.z >= ( nIndK + 0.5) * m_dStep &&
ptInt.z <= ( nIndK + 1.5) * m_dStep) {
if ( nIntNum == 0) {
dUF = dU2 ;
nFaceF = 2 ;
ptIntF = ptInt ;
++ nIntNum ;
}
else if ( nIntNum == 1 &&
( ptIntF - ptInt).SqLen() > dSqEps) {
dUL = dU2 ;
nFaceL = 2 ;
ptIntL = ptInt ;
++ nIntNum ;
}
}
// Intersezione con la quarta faccia
ptInt = ptP + dU4 * vtV ;
if ( ptInt.y >= ( nIndJ + 0.5) * m_dStep &&
ptInt.y <= ( nIndJ + 1.5) * m_dStep &&
ptInt.z >= ( nIndK + 0.5) * m_dStep &&
ptInt.z <= ( nIndK + 1.5) * m_dStep) {
if ( nIntNum == 0) {
dUF = dU4 ;
nFaceF = 4 ;
ptIntF = ptInt ;
++ nIntNum ;
}
else if ( nIntNum == 1 &&
( ptIntF - ptInt).SqLen() > dSqEps) {
dUL = dU4 ;
nFaceL = 4 ;
ptIntL = ptInt ;
++ nIntNum ;
}
}
}
// Intersezione con le facce 5 e 6
if ( abs( vtV.z) > EPS_ZERO) {
// Intersezione con la quinta faccia
ptInt = ptP + dU5 * vtV ;
if ( ptInt.x >= ( nIndI + 0.5) * m_dStep &&
ptInt.x <= ( nIndI + 1.5) * m_dStep &&
ptInt.y >= ( nIndJ + 0.5) * m_dStep &&
ptInt.y <= ( nIndJ + 1.5) * m_dStep) {
if ( nIntNum == 0) {
dUF = dU5 ;
nFaceF = 5 ;
ptIntF = ptInt ;
++ nIntNum ;
}
else if ( nIntNum == 1 &&
( ptIntF - ptInt).SqLen() > dSqEps) {
dUL = dU5 ;
nFaceL = 5 ;
ptIntL = ptInt ;
++ nIntNum ;
}
}
// Intersezione con la sesta faccia
ptInt = ptP + dU6 * vtV ;
if ( ptInt.x >= ( nIndI + 0.5) * m_dStep &&
ptInt.x <= ( nIndI + 1.5) * m_dStep &&
ptInt.y >= ( nIndJ + 0.5) * m_dStep &&
ptInt.y <= ( nIndJ + 1.5) * m_dStep) {
if ( nIntNum == 0) {
dUF = dU6 ;
nFaceF = 6 ;
ptIntF = ptInt ;
++ nIntNum ;
}
else if ( nIntNum == 1 &&
( ptIntF - ptInt).SqLen() > dSqEps) {
dUL = dU6 ;
nFaceL = 6 ;
ptIntL = ptInt ;
++ nIntNum ;
}
}
}
if ( dUF > dUL) {
swap( dUF, dUL) ;
swap( nFaceF, nFaceL) ;
}
return ( nIntNum == 2) ;
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineZMapBBox( unsigned int nGrid, const Point3d& ptP, const Vector3d& vtV, double& dU1, double& dU2)
VolZmap::IntersLineZMapBBox( const Point3d& ptP, const Vector3d& vtV, unsigned int nGrid, double& dU1, double& dU2)
{
// Punti estremi del box dello Zmap
Point3d ptMin = ORIG ;
Point3d ptMax = ptMin + Vector3d( m_nNx[nGrid] * m_dStep, m_nNy[nGrid] * m_dStep, m_dMaxZ[nGrid]) ;
return ( IntersLineBox( ptP, vtV, ptMin, ptMax, dU1, dU2) && ( dU1 > 0 || dU2 > 0)) ;
return IntersLineBox( ptP, vtV, ptMin, ptMax, dU1, dU2) ;
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineDexel( unsigned int nGrid, const Point3d& ptP, const Vector3d& vtV, unsigned int nI,
unsigned int nJ, double& dU1, double& dU2)
VolZmap::IntersLineDexel( const Point3d& ptP, const Vector3d& vtV, unsigned int nGrid, unsigned int nI, unsigned int nJ,
double& dU1, double& dU2)
{
// Determino l'indice del dexel e il doppio del numero di suo intervalli
// Determino l'indice del dexel e il numero di suoi intervalli
unsigned int nDexelPos = nJ * m_nNx[nGrid] + nI ;
unsigned int nDexelSize = unsigned int( m_Values[nGrid][nDexelPos].size()) ;
@@ -292,11 +110,11 @@ VolZmap::IntersLineDexel( unsigned int nGrid, const Point3d& ptP, const Vector3d
// estremi del box del singolo intervallo
Point3d ptE1( dXmin, dYmin, m_Values[nGrid][nDexelPos][nIndex].dMin) ;
Point3d ptE2( dXmax, dYmax, m_Values[nGrid][nDexelPos][nIndex].dMax) ;
double dt1, dt2 ;
if ( IntersLineBox( ptP, vtV, ptE1, ptE2, dt1, dt2)) {
double dP1, dP2 ;
if ( IntersLineBox( ptP, vtV, ptE1, ptE2, dP1, dP2)) {
bInters = true ;
dU1 = min( dU1, dt1) ;
dU2 = max( dU2, dt2) ;
dU1 = min( dU1, dP1) ;
dU2 = max( dU2, dP2) ;
}
}
@@ -304,21 +122,68 @@ VolZmap::IntersLineDexel( unsigned int nGrid, const Point3d& ptP, const Vector3d
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersRayDexel( const Point3d& ptP, const Vector3d& vtV, unsigned int nGrid, unsigned int nI, unsigned int nJ,
double& dU1, double& dU2)
{
// Determino l'indice del dexel e il numero di suoi intervalli
unsigned int nDexelPos = nJ * m_nNx[nGrid] + nI ;
unsigned int nDexelSize = unsigned int( m_Values[nGrid][nDexelPos].size()) ;
// Se non c'è materiale non devo fare alcunché
if ( nDexelSize == 0)
return false ;
// Determino estremi nel piano XY intrinseco del dexel
double dXmin = nI * m_dStep ;
double dYmin = nJ * m_dStep ;
double dXmax = ( nI + 1) * m_dStep ;
double dYmax = ( nJ + 1) * m_dStep ;
// ciclo sugli intervalli
dU1 = INFINITO ;
dU2 = - INFINITO ;
bool bInters = false ;
for ( unsigned int nIndex = 0 ; nIndex < nDexelSize ; nIndex += 1) {
// estremi del box del singolo intervallo
Point3d ptE1( dXmin, dYmin, m_Values[nGrid][nDexelPos][nIndex].dMin) ;
Point3d ptE2( dXmax, dYmax, m_Values[nGrid][nDexelPos][nIndex].dMax) ;
double dP1, dP2 ;
if ( IntersLineBox( ptP, vtV, ptE1, ptE2, dP1, dP2)) {
if ( dP2 >= 0) {
bInters = true ;
if ( dP1 >= 0)
dU1 = min( dU1, dP1) ;
else
dU1 = - 1 ;
dU2 = max( dU2, dP2) ;
}
}
}
return bInters ;
}
//----------------------------------------------------------------------------
// InLength = distanza di ingresso (se -1 il punto è interno, se -2 il punto è esterno e il raggio non interseca lo Zmap)
// OutLength = distanza di uscita
bool
VolZmap::GetDepth( const Point3d& ptPGlob, const Vector3d& vtDir, double& dInLength, double& dOutLength)
{
int nGrid = 0 ;
// Porto il raggio nel riferimento intrinseco
Point3d ptP = ptPGlob ;
ptP.ToLoc( m_MapFrame[0]) ;
ptP.ToLoc( m_MapFrame[nGrid]) ;
Vector3d vtV = vtDir ;
vtV.ToLoc( m_MapFrame[0]) ;
vtV.ToLoc( m_MapFrame[nGrid]) ;
vtV.Normalize() ;
// Studio dell'intersezione fra semiretta e BBox dello Zmap
// Intersezione fra semiretta e BBox dello Zmap
double dU1, dU2 ;
bool bTest = IntersLineZMapBBox( 0, ptP, vtV, dU1, dU2) ;
bool bTest = IntersLineZMapBBox( ptP, vtV, nGrid, dU1, dU2) && ( dU1 > 0 || dU2 > 0) ;
// Semiretta esterna al box dello Zmap
// Semiretta esterna al box dello Zmap quindi esterna anche allo Zmap
if ( ! bTest) {
dInLength = - 2 ;
dOutLength = - 2 ;
@@ -336,85 +201,65 @@ VolZmap::GetDepth( const Point3d& ptPGlob, const Vector3d& vtDir, double& dInLen
ptI = ptP + dU1 * vtV ;
ptF = ptP + dU2 * vtV ;
}
// Determinazione degli indici i j dei punti ptI e ptF
int nIi = Clamp( int( floor( ptI.x / m_dStep)), 0, m_nNx[0] - 1) ;
int nIj = Clamp( int( floor( ptI.y / m_dStep)), 0, m_nNy[0] - 1) ;
int nFi = Clamp( int( floor( ptF.x / m_dStep)), 0, m_nNx[0] - 1) ;
int nFj = Clamp( int( floor( ptF.y / m_dStep)), 0, m_nNy[0] - 1) ;
// Inizializzo distanze
dInLength = INFINITO ;
dOutLength = - INFINITO ;
// Variazioni
double dDeltaX = ptF.x - ptI.x ;
double dDeltaY = ptF.y - ptI.y ;
// Determino asse di spostamento maggiore
bool bOnX = ( abs( ptF.y - ptI.y) <= abs( ptF.x - ptI.x)) ;
// se inclinazione da asse X minore di 45 gradi (in assoluto)
if ( abs( dDeltaY) <= abs( dDeltaX)) {
// mi muovo lungo X (i)
int nIncrI = ( nFi >= nIi ? 1 : - 1) ;
for ( int i = nIi, j = nIj ;
i != nFi + nIncrI ;
i += nIncrI) {
// Movimento crescente su asse di movimento principale
if ( ( bOnX && ptF.x < ptI.x) || ( ! bOnX && ptF.y < ptI.y) )
swap( ptI, ptF) ;
// Controllo con nuovo i e j corrente (considero il bordo sinistro del dexel)
double dU1, dU2 ;
if ( IntersLineDexel( 0, ptP, vtV, i, j, dU1, dU2)) {
dInLength = min( dInLength, dU1) ;
dOutLength = max( dOutLength, dU2) ;
}
// Mi sposto sul bordo destro del dexel
double dMoveX = ( ( i + max( nIncrI, 0)) * m_dStep - ptI.x) ;
double dMoveY = dMoveX * dDeltaY / dDeltaX ;
double dY = ptI.y + dMoveY ;
int OldJ = j ;
j = Clamp( int( floor( dY / m_dStep)), 0, m_nNy[0] - 1) ;
// Pendenza
double dDeltaX = ( ptF.x - ptI.x) ;
double dDeltaY = ( ptF.y - ptI.y) ;
double dM = 0 ;
if ( bOnX && abs( dDeltaX) > EPS_SMALL)
dM = dDeltaY / dDeltaX ;
else if ( ! bOnX && abs( dDeltaY) > EPS_SMALL)
dM = dDeltaX / dDeltaY ;
// Analisi del dexel
if ( j != OldJ) {
double dU1, dU2 ;
if ( IntersLineDexel( 0, ptP, vtV, i, j, dU1, dU2)) {
dInLength = min( dInLength, dU1) ;
dOutLength = max( dOutLength, dU2) ;
}
}
// Determinazione degli indici i j dei punti ptI e ptF
int nIi = Clamp( int( floor( ptI.x / m_dStep)), 0, m_nNx[nGrid] - 1) ;
int nIj = Clamp( int( floor( ptI.y / m_dStep)), 0, m_nNy[nGrid] - 1) ;
int nFi = Clamp( int( floor( ptF.x / m_dStep)), 0, m_nNx[nGrid] - 1) ;
int nFj = Clamp( int( floor( ptF.y / m_dStep)), 0, m_nNy[nGrid] - 1) ;
// mi muovo
int i = nIi ; int j = nIj ;
while ( ( bOnX && i <= nFi) || ( ! bOnX && j <= nFj)) {
// Eseguo controllo
double dU1, dU2 ;
if ( IntersRayDexel( ptP, vtV, nGrid, i, j, dU1, dU2)) {
dInLength = min( dInLength, dU1) ;
dOutLength = max( dOutLength, dU2) ;
}
}
// altrimenti
else {
// mi muovo lungo Y (j)
int nIncrJ = ( nFj >= nIj ? 1 : - 1) ;
for ( int i = nIi, j = nIj ;
j != nFj + nIncrJ ;
j += nIncrJ) {
// Controllo con nuovo j e i corrente (considero il bordo sotto del dexel)
double dU1, dU2 ;
if ( IntersLineDexel( 0, ptP, vtV, i, j, dU1, dU2)) {
dInLength = min( dInLength, dU1) ;
dOutLength = max( dOutLength, dU2) ;
}
// Mi sposto sul bordo sopra del dexel
double dMoveY = ( ( j + max( nIncrJ, 0)) * m_dStep - ptI.y) ;
double dMoveX = dMoveY * dDeltaX / dDeltaY ;
// Calcolo spostamento (a destra o sopra)
if ( bOnX) {
double dMoveX = ( ( i + 1) * m_dStep - ptI.x) ;
double dMoveY = dMoveX * dM ;
double dY = ptI.y + dMoveY ;
int nNewJ = Clamp( int( floor( dY / m_dStep)), 0, m_nNy[nGrid] - 1) ;
if ( nNewJ != j)
j = nNewJ ;
else
++ i ;
}
else {
double dMoveY = ( ( j + 1) * m_dStep - ptI.y) ;
double dMoveX = dMoveY * dM ;
double dX = ptI.x + dMoveX ;
int OldI = i ;
i = Clamp( int( floor( dX / m_dStep)), 0, m_nNx[0] - 1) ;
// Analisi del dexel
if ( i != OldI) {
double dU1, dU2 ;
if ( IntersLineDexel( 0, ptP, vtV, i, j, dU1, dU2)) {
dInLength = min( dInLength, dU1) ;
dOutLength = max( dOutLength, dU2) ;
}
}
}
int nNewI = Clamp( int( floor( dX / m_dStep)), 0, m_nNx[nGrid] - 1) ;
if ( nNewI != i)
i = nNewI ;
else
++ j ;
}
}
// Se non abbiamo incontrato materiale