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

EgtGeomKernel 1.6x5 :

Browse Source 
- aggiunto marching cube a Zmap.
This commit is contained in:
Dario Sassi
2017-02-02 08:31:15 +00:00
parent acb2cea3b8
commit 0ebaa582c4
9 changed files with 3150 additions and 3756 deletions
Download Patch File Download Diff File Expand all files Collapse all files
VolTriZmapCalculus.cpp
+463 -145
View File
@@ -78,7 +78,7 @@ VolZmap::IntersLineZMapBBox( unsigned int nGrid, const Point3d& ptP, const Vecto
{
// Punti estremi del box dello Zmap
Point3d ptMin = ORIG ;
Point3d ptMax = ptMin + Vector3d( m_nVNx[nGrid] * m_dStep, m_nVNy[nGrid] * m_dStep, m_dVMaxZ[nGrid]) ;
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)) ;
}
@@ -89,8 +89,8 @@ VolZmap::IntersLineDexel( unsigned int nGrid, const Point3d& ptP, const Vector3d
unsigned int nJ, double& dU1, double& dU2)
{
// Determino l'indice del dexel e il doppio del numero di suo intervalli
unsigned int nDexelPos = nJ * m_nVNx[nGrid] + nI ;
unsigned int nDexelSize = unsigned int( m_TriZValues[nGrid][nDexelPos].size()) ;
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)
@@ -108,8 +108,8 @@ VolZmap::IntersLineDexel( unsigned int nGrid, const Point3d& ptP, const Vector3d
bool bInters = false ;
for ( unsigned int nIndex = 0 ; nIndex < nDexelSize ; nIndex += 2) {
// estremi del box del singolo intervallo
Point3d ptE1( dXmin, dYmin, m_TriZValues[nGrid][nDexelPos][nIndex]) ;
Point3d ptE2( dXmax, dYmax, m_TriZValues[nGrid][nDexelPos][nIndex+1]) ;
Point3d ptE1( dXmin, dYmin, m_Values[nGrid][nDexelPos][nIndex].dZVal) ;
Point3d ptE2( dXmax, dYmax, m_Values[nGrid][nDexelPos][nIndex+1].dZVal) ;
double dt1, dt2 ;
if ( IntersLineBox( ptP, vtV, ptE1, ptE2, dt1, dt2)) {
bInters = true ;
@@ -156,10 +156,10 @@ VolZmap::GetDepth( const Point3d& ptPGlob, const Vector3d& vtDir, double& dInLen
}
// Determinazione degli indici i j dei punti ptI e ptF
int nIi = Clamp( int( floor( ptI.x / m_dStep)), 0, m_nVNx[0] - 1) ;
int nIj = Clamp( int( floor( ptI.y / m_dStep)), 0, m_nVNy[0] - 1) ;
int nFi = Clamp( int( floor( ptF.x / m_dStep)), 0, m_nVNx[0] - 1) ;
int nFj = Clamp( int( floor( ptF.y / m_dStep)), 0, m_nVNy[0] - 1) ;
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 ;
@@ -189,7 +189,7 @@ VolZmap::GetDepth( const Point3d& ptPGlob, const Vector3d& vtDir, double& dInLen
double dMoveY = dMoveX * dDeltaY / dDeltaX ;
double dY = ptI.y + dMoveY ;
int OldJ = j ;
j = Clamp( int( floor( dY / m_dStep)), 0, m_nVNy[0] - 1) ;
j = Clamp( int( floor( dY / m_dStep)), 0, m_nNy[0] - 1) ;
// Analisi del dexel
if ( j != OldJ) {
@@ -222,7 +222,7 @@ VolZmap::GetDepth( const Point3d& ptPGlob, const Vector3d& vtDir, double& dInLen
double dMoveX = dMoveY * dDeltaX / dDeltaY ;
double dX = ptI.x + dMoveX ;
int OldI = i ;
i = Clamp( int( floor( dX / m_dStep)), 0, m_nVNx[0] - 1) ;
i = Clamp( int( floor( dX / m_dStep)), 0, m_nNx[0] - 1) ;
// Analisi del dexel
if ( i != OldI) {
@@ -260,7 +260,7 @@ VolZmap::AvoidBox( const Frame3d& frBox, const Vector3d& vtDiag)
b3Box.LocToLoc( frBox, m_MapFrame[0]) ;
// BBox dello Zmap nel suo riferimento intrinseco
BBox3d b3Zmap( ORIG, Point3d( m_nVNx[0] * m_dStep, m_nVNy[0] * m_dStep, m_dVMaxZ[0])) ;
BBox3d b3Zmap( ORIG, Point3d( m_nNx[0] * m_dStep, m_nNy[0] * m_dStep, m_dMaxZ[0])) ;
// Se non interferiscono, posso uscire
BBox3d b3Int ;
@@ -268,10 +268,10 @@ VolZmap::AvoidBox( const Frame3d& frBox, const Vector3d& vtDiag)
return true ;
// Limiti su indici
int nStI = Clamp( int( b3Int.GetMin().x / m_dStep), 0, m_nVNx[0] -1) ;
int nEnI = Clamp( int( b3Int.GetMax().x / m_dStep), 0, m_nVNx[0] -1) ;
int nStJ = Clamp( int( b3Int.GetMin().y / m_dStep), 0, m_nVNy[0] -1) ;
int nEnJ = Clamp( int( b3Int.GetMax().y / m_dStep), 0, m_nVNy[0] -1) ;
int nStI = Clamp( int( b3Int.GetMin().x / m_dStep), 0, m_nNx[0] -1) ;
int nEnI = Clamp( int( b3Int.GetMax().x / m_dStep), 0, m_nNx[0] -1) ;
int nStJ = Clamp( int( b3Int.GetMin().y / m_dStep), 0, m_nNy[0] -1) ;
int nEnJ = Clamp( int( b3Int.GetMax().y / m_dStep), 0, m_nNy[0] -1) ;
// Vettore direzione dei dexel nel riferimento del Box
Vector3d vtK = Z_AX ; vtK.LocToLoc( m_MapFrame[0], frBox) ;
@@ -286,8 +286,8 @@ VolZmap::AvoidBox( const Frame3d& frBox, const Vector3d& vtDiag)
for ( int j = nStJ ; j <= nEnJ ; ++ j) {
int nPos = j * m_nVNx[0] + i ;
int nSize = int( m_TriZValues[0][nPos].size()) ;
int nPos = j * m_nNx[0] + i ;
int nSize = int( m_Values[0][nPos].size()) ;
if ( nSize == 0)
continue ;
@@ -297,8 +297,8 @@ VolZmap::AvoidBox( const Frame3d& frBox, const Vector3d& vtDiag)
if ( IntersLineBox( ptC, vtK, ORIG, ORIG + vtDiag, dZmin, dZmax)) {
for ( int nIndex = 0 ; nIndex < nSize ; nIndex += 2) {
if ( ! ( dZmax < m_TriZValues[0][nPos][nIndex] - EPS_SMALL ||
dZmin > m_TriZValues[0][nPos][nIndex + 1] + EPS_SMALL))
if ( ! ( dZmax < m_Values[0][nPos][nIndex].dZVal - EPS_SMALL ||
dZmin > m_Values[0][nPos][nIndex + 1].dZVal + EPS_SMALL))
return false ;
}
}
@@ -340,7 +340,7 @@ VolZmap::IntersLineCylinder( const Point3d& ptLineSt, const Vector3d& vtLineDir,
int nRoot = PolynomialRoots( 2, vdCoef, vdRoots) ;
// Nessuna soluzione
if ( nRoot == 0) {
if ( nRoot == 0 || nRoot == 1) {
if ( abs( vtV.x) > EPS_ZERO) {
@@ -375,26 +375,30 @@ VolZmap::IntersLineCylinder( const Point3d& ptLineSt, const Vector3d& vtLineDir,
if ( ptInt1.x > ptInt2.x)
swap( ptInt1, ptInt2) ;
if ( ptInt1.x < 0 && ptInt2.x >= 0 && ptInt2.x <= dL) {
if ( ptInt1.x < dL + EPS_SMALL) {
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
if ( ptInt1.x > - EPS_SMALL) {
if ( ptInt2.x > dL)
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else {
if ( ptInt2.x > dL) {
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt2.x > 0)
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
else
return false ;
}
}
else if ( ptInt1.x < 0 && ptInt2.x > dL) {
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x >= 0 && ptInt2.x <= dL) {
;
}
else if ( ptInt1.x >= 0 && ptInt1.x < dL && ptInt2.x >= dL) {
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
// Intersezioni esterne alla regione di interesse
// (quella compresa fra 0 e dL)
else
return false ;
@@ -409,7 +413,7 @@ VolZmap::IntersLineCylinder( const Point3d& ptLineSt, const Vector3d& vtLineDir,
//----------------------------------------------------------------------------
bool
VolZmap::IntersZLineCylinder( const Point3d& ptLine,
const Point3d& ptBase, const Point3d& ptTop, const Vector3d& vtDir, double dCylR,
const Point3d& ptBase, const Point3d& ptTop, double dCylR,
double& dInfZ, double& dSupZ)
{
// NB: Le coordinate sono espresse nel sistema griglia
@@ -422,7 +426,7 @@ VolZmap::IntersZLineCylinder( const Point3d& ptLine,
if ( AreSamePointXYApprox( ptBase, ptTop)) {
// Intersezione
if ( SqDistXY( ptLine, ptBase) <= dSqRad) {
if ( SqDistXY( ptLine, ptBase) < dSqRad) {
dInfZ = min( ptBase.z, ptTop.z) ;
dSupZ = max( ptBase.z, ptTop.z) ;
return true ;
@@ -440,7 +444,9 @@ VolZmap::IntersZLineCylinder( const Point3d& ptLine,
Point3d ptS = ( ptBase.z < ptTop.z ? ptBase : ptTop) ;
Point3d ptE = ( ptBase.z < ptTop.z ? ptTop : ptBase) ;
Vector3d vtV1 = ptE - ptS ; vtV1.z = 0 ;
Vector3d vtAx = ptE - ptS ;
Vector3d vtV1( vtAx.x, vtAx.y, 0) ;
double dLenXY = vtV1.LenXY() ;
double dSZ = ptS.z ;
@@ -456,7 +462,7 @@ VolZmap::IntersZLineCylinder( const Point3d& ptLine,
Vector3d vtV2 = vtV1 ;
vtV2.Rotate( Z_AX, 90) ;
double dLen = ( ptE - ptS).Len() ;
double dLen = vtAx.Len() ;
// Sono seno e coseno dell'angolo complementare
// rispetto a quello formato dal vettore movimento
@@ -478,9 +484,9 @@ VolZmap::IntersZLineCylinder( const Point3d& ptLine,
// Minimi
if ( dLocX1 < dX1_0) {
double dDotS = vtDir * ( ptS - ORIG) ;
double dDotS = vtAx * ( ptS - ORIG) ;
// Qui usiamo ptLine perché servono coordinate griglia
dInfZ = ( dDotS - vtDir.x * ptLine.x - vtDir.y * ptLine.y) / vtDir.z ;
dInfZ = ( dDotS - vtAx.x * ptLine.x - vtAx.y * ptLine.y) / vtAx.z ;
}
else {
@@ -490,17 +496,17 @@ VolZmap::IntersZLineCylinder( const Point3d& ptLine,
}
// Massimi
if ( dLocX1 < dLenXY - dX1_0) {
if ( dLocX1 <= dLenXY - dX1_0) {
double dZ0 = dSin * dSqRoot ;
dSupZ = dSZ + dZ0 + ( dLocX1 - dX1_0) * dDeltaZ / dLenXY ;
dSupZ = dSZ + dZ0 + ( dLocX1 + dX1_0) * dDeltaZ / dLenXY ;
}
else {
double dDotE = vtDir * ( ptE - ORIG) ;
double dDotE = vtAx * ( ptE - ORIG) ;
// Qui usiamo ptLine perché servono coordinate griglia
dSupZ = ( dDotE - vtDir.x * ptLine.x - vtDir.y * ptLine.y) / vtDir.z ;
dSupZ = ( dDotE - vtAx.x * ptLine.x - vtAx.y * ptLine.y) / vtAx.z ;
}
return true ;
@@ -533,6 +539,8 @@ VolZmap::IntersLineConus( const Point3d& ptLineSt, const Vector3d& vtLineDir,
DBLVECTOR vdRoots ;
double dSqTan = dTan * dTan ;
double dMinRad = dTan * dl ;
double dMaxRad = dTan * dL ;
vdCoef[0] = dSqTan * ptP.x * ptP.x - ptP.y * ptP.y - ptP.z * ptP.z ;
vdCoef[1] = 2 * ( dSqTan * ptP.x * vtV.x - ptP.y * vtV.y - ptP.z * vtV.z) ;
@@ -550,21 +558,33 @@ VolZmap::IntersLineConus( const Point3d& ptLineSt, const Vector3d& vtLineDir,
if ( nRoot == 1) {
ptInt1 = ptP + vdRoots[0] * vtV ;
if ( ptInt1.x >= dl && ptInt1.x < dL) {
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x >= 0 && ptInt1.x < dl) {
if ( ptInt1.x < dL + EPS_SMALL) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
if ( ptInt1.x > dl)
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
// Riporto le coordinate nel sistema di riferimento
// griglia
ptInt1.ToGlob( ConusFrame) ;
ptInt2.ToGlob( ConusFrame) ;
else if ( ptInt1.x > - EPS_SMALL) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
if ( ptInt2.y * ptInt2.y + ptInt2.z * ptInt2.z > dMaxRad * dMaxRad)
return false ;
}
else
return false ;
ptInt1.ToGlob( ConusFrame) ;
ptInt2.ToGlob( ConusFrame) ;
return true ;
}
else
return false ;
}
// Due soluzioni: la retta interseca due volte la
// superficie laterale
@@ -573,54 +593,67 @@ VolZmap::IntersLineConus( const Point3d& ptLineSt, const Vector3d& vtLineDir,
ptInt1 = ptP + vdRoots[0] * vtV ;
ptInt2 = ptP + vdRoots[1] * vtV ;
if ( ptInt1.x > ptInt2.x) {
if ( ptInt1.x > ptInt2.x)
swap( ptInt1, ptInt2) ;
}
if ( ptInt1.x < 0 && ptInt2.x > 0 && ptInt2.x < dl) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x < 0 && ptInt2.x >= dl && ptInt2.x < dL) {
ptInt1 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x > 0 && ptInt1.x < dl && ptInt2.x >= dl && ptInt2.x < dL) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x > 0 && ptInt1.x < dl && ptInt2.x >= dL) {
if ( ptInt1.x < dL + EPS_SMALL) {
if ( ptInt1.x > dl - EPS_SMALL) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x >= dl && ptInt1.x < dL && ptInt2.x < dL) {
if ( ptInt2.x > dL)
;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt1.x > - EPS_SMALL) {
if ( ptInt2.x > dL) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt2.x > dl - EPS_SMALL)
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
else
return false ;
}
else {
if ( ptInt2.x < 0)
return false ;
else if ( ptInt2.x < dl) {
ptInt1 = ptP + ( ( dl - ptP.x) / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptInt2.x < dL)
ptInt1 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
else
return false ;
}
ptInt1.ToGlob( ConusFrame) ;
ptInt2.ToGlob( ConusFrame) ;
return true ;
}
else if ( ptInt1.x >= dl && ptInt1.x < dL && ptInt2.x >= dL) {
ptInt2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
// Intersezioni esterne alla regione di interesse
// (quella compresa fra dl e dL)
else
return false ;
// Riporto le coordinate nel sistema di riferimento
// griglia
ptInt1.ToGlob( ConusFrame) ;
ptInt2.ToGlob( ConusFrame) ;
}
return true ;
return false ;
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineEllipticalCylinder( const Frame3d & CircFrame, const Vector3d & vtLineDir, const Point3d & ptLineSt,
std::vector <Point3d> & ptInters, double dObCoef, double dSqRad, double dL)
VolZmap::IntersLineEllipticalCylinder( const Vector3d& vtLineDir, const Point3d& ptLineSt,
const Frame3d& CircFrame, double dSqRad, double dLongMvLen, double dOrtMvLen,
Point3d& ptInt1, Point3d& ptInt2)
{
// NB: L'origine del sistema di riferimento deve essere
// nel centro della circonferenza di base, la cui tralsazione obliqua
@@ -628,21 +661,20 @@ VolZmap::IntersLineEllipticalCylinder( const Frame3d & CircFrame, const Vector3d
// di simmetria di tale circonferenza.
// La funzione restituisce true in caso di intersezione,
// false altrimenti.
// NB: Il Parametro dObCoef è il coeffociente angolare della retta movimento
// rispetto all'asse x, e dSqRad è il quadrato del raggio della circonferenza.
// NB: dSqRad è il quadrato del raggio della circonferenza la cui
// traslazione obliqua genera il cilindro ellittico, dLongMvLen e
// dOrtMvLen sono rispettivamente le lunghezze delle proiezioni del
// movimento su x e y del sistema di riferimento CircFrame.
double dObCoef = dOrtMvLen / dLongMvLen ;
double dSqCoef = dObCoef * dObCoef ;
Point3d ptP = ptLineSt ;
Vector3d vtV = vtLineDir ;
// Sistema di riferimanto grigia
Frame3d GridFrame ;
GridFrame.Set( ORIG, X_AX, Y_AX, Z_AX) ;
// Trasformazione delle coordinate
ptP.LocToLoc( GridFrame, CircFrame) ;
vtV.LocToLoc( GridFrame, CircFrame) ;
ptP.ToLoc( CircFrame) ;
vtV.ToLoc( CircFrame) ;
std::vector <double> vdCoef(3) ;
std::vector <double> vdRoots ;
@@ -653,26 +685,20 @@ VolZmap::IntersLineEllipticalCylinder( const Frame3d & CircFrame, const Vector3d
int nRoot = PolynomialRoots( 2, vdCoef, vdRoots) ;
Point3d ptR1, ptR2 ;
// Nessuna soluzione
if ( nRoot == 0) {
if ( nRoot == 0 || nRoot == 1) {
if ( abs( vtV.x) > EPS_ZERO) {
ptR1 = ptP - ( ptP.x / vtV.x) * vtV ;
ptR2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dLongMvLen - ptP.x) / vtV.x) * vtV ;
if ( ptR1.y * ptR1.y + ptR1.z * ptR1.z < dSqRad &&
ptR1.y * ptR1.y + ptR1.z * ptR1.z < dSqRad) {
if ( ptInt1.y * ptInt1.y + ptInt1.z * ptInt1.z < dSqRad &&
( ptInt2.y - dOrtMvLen) * ( ptInt2.y - dOrtMvLen) + ptInt2.z * ptInt2.z < dSqRad) {
ptR1.LocToLoc( CircFrame, GridFrame) ;
ptR2.LocToLoc( CircFrame, GridFrame) ;
ptInters.resize(2) ;
ptInters[0] = ptR1 ;
ptInters[0] = ptR2 ;
ptInt1.ToGlob( CircFrame) ;
ptInt2.ToGlob( CircFrame) ;
return true ;
}
@@ -689,48 +715,340 @@ VolZmap::IntersLineEllipticalCylinder( const Frame3d & CircFrame, const Vector3d
// coincidenti) oppure nessuna o infinite se la la retta
// appartiene alla superficie
ptInters.resize(2) ;
if ( nRoot == 2) {
ptR1 = ptP + vdRoots[0] * vtV ;
ptR2 = ptP + vdRoots[1] * vtV ;
ptInt1 = ptP + vdRoots[0] * vtV ;
ptInt2 = ptP + vdRoots[1] * vtV ;
if ( ptR1.x > ptR2.x) {
if ( ptInt1.x > ptInt2.x)
Point3d ptTemp = ptR1 ;
ptR1 = ptR2 ;
ptR2 = ptTemp ;
}
swap( ptInt1, ptInt2) ;
if ( ptR1.x >= 0 && ptR1.x < dL &&
ptR2.x > dL) {
ptR1 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptR1.x >= 0 && ptR2.x <= dL) {
if ( ptInt1.x < dLongMvLen + EPS_SMALL) {
if ( ptInt1.x > - EPS_SMALL) {
;
}
else if ( ptR1.x < 0 && ptR2.x > dL) {
if ( ptInt2.x > dLongMvLen)
ptR1 = ptP - ( ptP.x / vtV.x) * vtV ;
ptR2 = ptP + ( ( dL - ptP.x) / vtV.x) * vtV ;
}
else if ( ptR1.x < 0 && ptR2.x >= 0 && ptR2.x <= dL) {
ptInt2 = ptP + ( ( dLongMvLen - ptP.x) / vtV.x) * vtV ;
}
else {
ptR1 = ptP - ( ptP.x / vtV.x) * vtV ;
if ( ptInt2.x > dLongMvLen) {
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
ptInt2 = ptP + ( ( dLongMvLen - ptP.x) / vtV.x) * vtV ;
if ( ptInt1.y * ptInt1.y + ptInt1.z * ptInt1.z > dSqRad &&
ptInt2.y * ptInt2.y + ptInt2.z * ptInt2.z > dSqRad)
return false ;
}
else if ( ptInt2.x > 0)
ptInt1 = ptP - ( ptP.x / vtV.x) * vtV ;
else
return false ;
}
}
else
return false ;
// Riporto le coordinate nel sistema di riferimento
// griglia
ptR1.LocToLoc( CircFrame, GridFrame) ;
ptR2.LocToLoc( CircFrame, GridFrame) ;
ptInters[0] = ptR1 ;
ptInters[1] = ptR2 ;
ptInt1.ToGlob( CircFrame) ;
ptInt2.ToGlob( CircFrame) ;
}
return true ;
}
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineMyPolyhedron( const Point3d& ptLineSt, const Vector3d& vtLineDir,
const Frame3d& PolyFrame, double dLenX, double dLenY, double dLenZ, double dDeltaX,
Point3d& ptInt1, Point3d& ptInt2)
{
double SqIndet = EPS_SMALL * EPS_SMALL ;
Point3d ptP = ptLineSt ;
Vector3d vtV = vtLineDir ;
// Trasformazione delle coordinate
ptP.ToLoc( PolyFrame) ;
vtV.ToLoc( PolyFrame) ;
// Facce 1 e 2 parallele a XY
// Facce 3 e 4 parallele a XZ
// Facce 5 e 6 oblique
Point3d ptFacet135( 0, 0, dLenZ /2) ;
Point3d ptFacet246( dLenX + dDeltaX, dLenY, - dLenZ / 2) ;
// Servono per descrivere i piani obliqui
Vector3d vtFacet5 = ptFacet135 - ptP ;
Vector3d vtFacet6 = ptFacet246 - ptP ;
Vector3d vtOb( dLenY, - dDeltaX, 0) ;
Point3d ptI1 = ptP + ( ( ptFacet135.z - ptP.z) / vtV.z) * vtV ;
Point3d ptI2 = ptP + ( ( ptFacet246.z - ptP.z) / vtV.z) * vtV ;
Point3d ptI3 = ptP + ( ( ptFacet135.y - ptP.y) / vtV.y) * vtV ;
Point3d ptI4 = ptP + ( ( ptFacet246.y - ptP.y) / vtV.y) * vtV ;
Point3d ptI5 = ptP + ( ( vtFacet5 * vtOb) / ( vtV * vtOb)) * vtV ;
Point3d ptI6 = ptP + ( ( vtFacet6 * vtOb) / ( vtV * vtOb)) * vtV ;
int nIntNum = 0 ;
// Intersezione con la prima faccia
if ( ptI1.y >= 0 && ptI1.y <= dLenY &&
ptI1.x * dLenY >= dDeltaX * ptI1.y && ( ptI1.x - dLenX) * dLenY <= dDeltaX * ptI1.y) {
ptInt1 = ptI1 ;
++ nIntNum ;
}
// Intersezione con la seconda faccia
if ( ptI2.y >= 0 && ptI2.y <= dLenY &&
ptI2.x * dLenY >= dDeltaX * ptI2.y && ( ptI2.x - dLenX) * dLenY <= dDeltaX * ptI2.y) {
if ( nIntNum == 0) {
ptInt1 = ptI2 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI2).SqLen() > SqIndet) {
ptInt2 = ptI2 ;
++ nIntNum ;
}
}
// Intersezione con la terza faccia
if ( nIntNum < 2 &&
ptI3.x >= 0 && ptI3.x <= dLenX &&
ptI3.z >= - ptFacet135.z && ptI3.z <= ptFacet135.z) {
if ( nIntNum == 0) {
ptInt1 = ptI3 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI3).SqLen() > SqIndet) {
ptInt2 = ptI3 ;
++ nIntNum ;
}
}
// Intersezione con la quarta faccia
if ( nIntNum < 2 &&
ptI4.x >= dDeltaX && ptI4.x <= dLenX + dDeltaX &&
ptI4.z >= - ptFacet135.z && ptI4.z <= ptFacet135.z) {
if ( nIntNum == 0) {
ptInt1 = ptI4 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI4).SqLen() > SqIndet) {
ptInt2 = ptI4 ;
++ nIntNum ;
}
}
// Intersezione con la quinta faccia
if ( nIntNum < 2 &&
ptI5.y >= 0 && ptI5.y <= dLenY &&
ptI5.z >= - ptFacet135.z && ptI5.z <= ptFacet135.z) {
if ( nIntNum == 0) {
ptInt1 = ptI5 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI5).SqLen() > SqIndet) {
ptInt2 = ptI5 ;
++ nIntNum ;
}
}
// Intersezione con la sesta faccia
if ( nIntNum < 2 &&
ptI6.y >= 0 && ptI6.y <= dLenY &&
ptI6.z >= - ptFacet135.z && ptI6.z <= ptFacet135.z) {
if ( nIntNum == 0) {
ptInt1 = ptI6;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI6).SqLen() > SqIndet) {
ptInt2 = ptI6;
++ nIntNum ;
}
}
if ( nIntNum == 2) {
ptInt1.ToGlob( PolyFrame) ;
ptInt2.ToGlob( PolyFrame) ;
return true ;
}
else
return false ;
}
/*
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineMyPolyhedron( const Point3d& ptLineSt, const Vector3d& vtLineDir,
const Frame3d& PolyFrame, double dLenX, double dLenY, double dLenZ, double dDeltaX,
Point3d& ptInt1, Point3d& ptInt2)
{
double SqIndet = EPS_SMALL * EPS_SMALL ;
Point3d ptP = ptLineSt ;
Vector3d vtV = vtLineDir ;
// Trasformazione delle coordinate
ptP.ToLoc( PolyFrame) ;
vtV.ToLoc( PolyFrame) ;
// Facce 1 e 2 parallele a XY
// Facce 3 e 4 parallele a XZ
// Facce 5 e 6 oblique
Point3d ptFacet135( 0, 0, dLenZ /2) ;
Point3d ptFacet246( dLenX + dDeltaX, dLenY, - dLenZ / 2) ;
// Servono per descrivere i piani obliqui
Vector3d vtFacet5 = ptFacet135 - ptP ;
Vector3d vtFacet6 = ptFacet246 - ptP ;
Vector3d vtOb( dLenY, - dDeltaX, 0) ;
Point3d ptI1 = ptP + ( ( ptFacet135.z - ptP.z) / vtV.z) * vtV ;
Point3d ptI2 = ptP + ( ( ptFacet246.z - ptP.z) / vtV.z) * vtV ;
Point3d ptI3 = ptP + ( ( ptFacet135.y - ptP.y) / vtV.y) * vtV ;
Point3d ptI4 = ptP + ( ( ptFacet246.y - ptP.y) / vtV.y) * vtV ;
Point3d ptI5 = ptP + ( ( vtFacet5 * vtOb) / ( vtV * vtOb)) * vtV ;
Point3d ptI6 = ptP + ( ( vtFacet6 * vtOb) / ( vtV * vtOb)) * vtV ;
int nIntNum = 0 ;
// Intersezione con la prima faccia
if ( ptI1.y > - EPS_SMALL && ptI1.y < dLenY + EPS_SMALL &&
ptI1.x * dLenY > dDeltaX * ptI1.y - EPS_SMALL &&
( ptI1.x - dLenX) * dLenY < dDeltaX * ptI1.y + EPS_SMALL) {
ptInt1 = ptI1 ;
++ nIntNum ;
}
// Intersezione con la seconda faccia
if ( ptI2.y > - EPS_SMALL && ptI2.y <dLenY + EPS_SMALL &&
ptI2.x * dLenY > dDeltaX * ptI2.y - EPS_SMALL &&
( ptI2.x - dLenX) * dLenY < dDeltaX * ptI2.y + EPS_SMALL) {
if ( nIntNum == 0) {
ptInt1 = ptI2 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI2).SqLen() > SqIndet) {
ptInt2 = ptI2 ;
++ nIntNum ;
}
}
// Intersezione con la terza faccia
if ( nIntNum < 2 &&
ptI3.x > - EPS_SMALL && ptI3.x < dLenX + EPS_SMALL &&
ptI3.z > - ptFacet135.z - EPS_SMALL &&
ptI3.z < ptFacet135.z + EPS_SMALL) {
if ( nIntNum == 0) {
ptInt1 = ptI3 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI3).SqLen() > SqIndet) {
ptInt2 = ptI3 ;
++ nIntNum ;
}
}
// Intersezione con la quarta faccia
if ( nIntNum < 2 &&
ptI4.x > dDeltaX - EPS_SMALL && ptI4.x < dLenX + dDeltaX + EPS_SMALL &&
ptI4.z > - ptFacet135.z - EPS_SMALL &&
ptI4.z < ptFacet135.z + EPS_SMALL) {
if ( nIntNum == 0) {
ptInt1 = ptI4 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI4).SqLen() > SqIndet) {
ptInt2 = ptI4 ;
++ nIntNum ;
}
}
// Intersezione con la quinta faccia
if ( nIntNum < 2 &&
ptI5.y > - EPS_SMALL && ptI5.y < dLenY + EPS_SMALL &&
ptI5.z > - ptFacet135.z - EPS_SMALL &&
ptI5.z < ptFacet135.z + EPS_SMALL) {
if ( nIntNum == 0) {
ptInt1 = ptI5 ;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI5).SqLen() > SqIndet) {
ptInt2 = ptI5 ;
++ nIntNum ;
}
}
// Intersezione con la sesta faccia
if ( nIntNum < 2 &&
ptI6.y > - EPS_SMALL && ptI6.y < dLenY + EPS_SMALL &&
ptI6.z > - ptFacet135.z - EPS_SMALL &&
ptI6.z < ptFacet135.z + EPS_SMALL) {
if ( nIntNum == 0) {
ptInt1 = ptI6;
++ nIntNum ;
}
else if ( ( ptInt1 - ptI6).SqLen() > SqIndet) {
ptInt2 = ptI6;
++ nIntNum ;
}
}
if ( nIntNum == 2) {
ptInt1.ToGlob( PolyFrame) ;
ptInt2.ToGlob( PolyFrame) ;
return true ;
}
else
return false ;
} */