EgtGeomKernel :

- correzioni a Zmap.
This commit is contained in:
Dario Sassi
2018-02-10 10:23:28 +00:00
parent 1aa106305a
commit 6cf3040131
3 changed files with 45 additions and 166 deletions
+33 -148
View File
@@ -808,111 +808,6 @@ VolZmap::IntersLineCylinder( const Point3d& ptLineSt, const Vector3d& vtLineDir,
return true ;
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersZLineCylinder( const Point3d& ptLine,
const Point3d& ptBase, const Point3d& ptTop, double dCylR,
double& dInfZ, double& dSupZ)
{
// NB: Le coordinate sono espresse nel sistema griglia
// La funzione restituisce true in caso di intersezione, false altrimenti.
double dSqRad = dCylR * dCylR ;
// Cilindro verticale
if ( AreSamePointXYApprox( ptBase, ptTop)) {
// Intersezione
if ( SqDistXY( ptLine, ptBase) < dSqRad) {
dInfZ = min( ptBase.z, ptTop.z) ;
dSupZ = max( ptBase.z, ptTop.z) ;
return true ;
}
// Non vi è intersezione
else
return false ;
}
// Cilindro non verticale
else {
// Studio delle simmetrie
Point3d ptS = ( ptBase.z < ptTop.z ? ptBase : ptTop) ;
Point3d ptE = ( ptBase.z < ptTop.z ? ptTop : ptBase) ;
Vector3d vtAx = ptE - ptS ;
Vector3d vtV1( vtAx.x, vtAx.y, 0) ;
double dLenXY = vtV1.LenXY() ;
double dSZ = ptS.z ;
double dEZ = ptE.z ;
double dDeltaZ = dEZ - dSZ ;
Vector3d vtL( ptLine.x - ptS.x, ptLine.y - ptS.y, 0) ;
// vtV1 e vtV2 formano un sistema ortonormale
// sul piano e insieme a ptSxy formano un sistema
// di riferimento bidimensionale
vtV1.Normalize() ;
Vector3d vtV2 = vtV1 ;
vtV2.Rotate( Z_AX, 90) ;
double dLen = vtAx.Len() ;
// Sono seno e coseno dell'angolo complementare
// rispetto a quello formato dal vettore movimento
// con il piano, per questo motivo si ha dCos con
// dDeltaZ e dSin con dLenXY
double dCos = dDeltaZ / dLen ;
double dSin = dLenXY / dLen ;
// Nuove coordinate piane del punto
double dLocX1 = vtL * vtV1 ;
double dLocX2 = vtL * vtV2 ;
double dSqRoot = sqrt( dSqRad - dLocX2 * dLocX2) ;
double dX1_0 = dCos * dSqRoot ;
if ( dLocX1 >= - dX1_0 && dLocX1 <= dLenXY + dX1_0 &&
abs( dLocX2) < dCylR) {
// Minimi
if ( dLocX1 < dX1_0) {
double dDotS = vtAx * ( ptS - ORIG) ;
// Qui usiamo ptLine perché servono coordinate griglia
dInfZ = ( dDotS - vtAx.x * ptLine.x - vtAx.y * ptLine.y) / vtAx.z ;
}
else {
double dZ0 = - dSin * dSqRoot ;
dInfZ = dSZ + dZ0 + ( dLocX1 - dX1_0) * dDeltaZ / dLenXY ;
}
// Massimi
if ( dLocX1 <= dLenXY - dX1_0) {
double dZ0 = dSin * dSqRoot ;
dSupZ = dSZ + dZ0 + ( dLocX1 + dX1_0) * dDeltaZ / dLenXY ;
}
else {
double dDotE = vtAx * ( ptE - ORIG) ;
// Qui usiamo ptLine perché servono coordinate griglia
dSupZ = ( dDotE - vtAx.x * ptLine.x - vtAx.y * ptLine.y) / vtAx.z ;
}
return true ;
}
return false ;
}
}
//----------------------------------------------------------------------------
bool
VolZmap::IntersLineConus( const Point3d& ptLineSt, const Vector3d& vtLineDir,
@@ -1115,24 +1010,23 @@ VolZmap::IntersLineConus( const Point3d& ptLineSt, const Vector3d& vtLineDir,
}
//----------------------------------------------------------------------------
// NB: L'origine del sistema di riferimento deve essere
// nel centro della circonferenza di base, la cui traslazione obliqua
// genera il cilindro ellittico, e l'asse z deve essere l'asse
// di simmetria di tale circonferenza.
// La funzione restituisce true in caso di intersezione, false altrimenti.
// 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 z e x del sistema di riferimento CircFrame.
bool
VolZmap::IntersLineEllipticalCylinder( const Vector3d& vtLineDir, const Point3d& ptLineSt,
const Frame3d& CircFrame, double dSqRad, double dLongMvLen, double dOrtMvLen,
const Frame3d& CircFrame, double dRad, double dLongMvLen, double dOrtMvLen,
bool bTapLow, bool bTapUp,
Point3d& ptInt1, Point3d& ptInt2, Vector3d& vtN1, Vector3d& vtN2)
{
// NB: L'origine del sistema di riferimento deve essere
// nel centro della circonferenza di base, la cui traslazione obliqua
// genera il cilindro ellittico, e l'asse z deve essere l'asse
// di simmetria di tale circonferenza.
// La funzione restituisce true in caso di intersezione, false altrimenti.
// 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 z e x del sistema di riferimento CircFrame.
double dObCoef = dOrtMvLen / dLongMvLen ;
double dSqCoef = dObCoef * dObCoef ;
// Quadrato del raggio
double dSqRad = dRad * dRad ;
// Punto e vettore individuanti la retta
Point3d ptP = ptLineSt ;
Vector3d vtV = vtLineDir ;
@@ -1152,7 +1046,7 @@ VolZmap::IntersLineEllipticalCylinder( const Vector3d& vtLineDir, const Point3d&
Point3d ptOLsCirc( dOrtMvLen, 0, dLongMvLen) ;
ptInt1 = ptP - ( ptP.z / vtV.z) * vtV ;
ptInt2 = ptP - ( ( ptP.z - dLongMvLen) / vtV.z) * vtV ;
double dSafeSqRad = dSqRad - 2 * sqrt( dSqRad) * EPS_SMALL ;
double dSafeSqRad = dSqRad - 2 * dRad * EPS_SMALL ;
if ( ( ptInt1 - ORIG).SqLenXY() < dSafeSqRad &&
( ptInt2 - ptOLsCirc).SqLenXY() < dSafeSqRad) {
vtN1 = Z_AX ;
@@ -1169,10 +1063,18 @@ VolZmap::IntersLineEllipticalCylinder( const Vector3d& vtLineDir, const Point3d&
vector <double> vdCoef(3) ;
vector <double> vdRoots ;
// Coefficiente angolare della retta di movimento nel
// piano ZX del sistema di riferimento del movimento
// e suo quadrato
double dObCoef = dOrtMvLen / dLongMvLen ;
double dSqCoef = dObCoef * dObCoef ;
// Setto i coeficienti dell'equazione
vdCoef[0] = dSqCoef * ptP.z * ptP.z + ptP.x * ptP.x + ptP.y * ptP.y - 2 * dObCoef * ptP.z * ptP.x - dSqRad ;
vdCoef[1] = 2 * ( dSqCoef * vtV.z * ptP.z + vtV.x * ptP.x + vtV.y * ptP.y - dObCoef * ( vtV.z * ptP.x + vtV.x * ptP.z)) ;
vdCoef[2] = dSqCoef * vtV.z * vtV.z + vtV.x * vtV.x + vtV.y * vtV.y - 2 * dObCoef * vtV.z * vtV.x ;
vdCoef[0] = dSqCoef * ptP.z * ptP.z + ptP.x * ptP.x + ptP.y * ptP.y
- 2 * dObCoef * ptP.z * ptP.x - dSqRad ;
vdCoef[1] = 2 * ( dSqCoef * vtV.z * ptP.z + vtV.x * ptP.x + vtV.y * ptP.y
- dObCoef * ( vtV.z * ptP.x + vtV.x * ptP.z)) ;
vdCoef[2] = dSqCoef * vtV.z * vtV.z + vtV.x * vtV.x + vtV.y * vtV.y
- 2 * dObCoef * vtV.z * vtV.x ;
// Numero di soluzioni
int nRoot = PolynomialRoots( 2, vdCoef, vdRoots) ;
@@ -1180,31 +1082,13 @@ VolZmap::IntersLineEllipticalCylinder( const Vector3d& vtLineDir, const Point3d&
// coincidenti) oppure nessuna o infinite se la la retta
// appartiene alla superficie
// Sperimentale
if ( nRoot == 0 || nRoot == 1) {
if ( abs( vtV.z) > EPS_ZERO) {
ptInt1 = ptP - ( ptP.z / vtV.z) * vtV ;
ptInt2 = ptP + ( ( dLongMvLen - ptP.z) / vtV.z) * vtV ;
if ( ptInt1.x * ptInt1.x + ptInt1.y * ptInt1.y < dSqRad &&
( ptInt2.x - dOrtMvLen) * ( ptInt2.x - dOrtMvLen) + ptInt2.y * ptInt2.y < dSqRad) {
ptInt1.ToGlob( CircFrame) ;
ptInt2.ToGlob( CircFrame) ;
vtN1 = Z_AX ;
vtN2 = - Z_AX ;
vtN1.ToGlob( CircFrame) ;
vtN2.ToGlob( CircFrame) ;
return true ;
}
// Nessuna intersezione
else
return false ;
}
// Nessuna intersezione
else
return false ;
}// Fine sperimentale
if ( nRoot == 2) {
// Se ci sono intersezioni con i tappi o l'equazione
// degenera in una di primo grado, le eventuali
// soluzioni sono già state trovate.
if ( nRoot == 0 || nRoot == 1)
return false ;
// Due soluzioni trovate
else if ( nRoot == 2) {
// Flag per i tappi
double dEpsLow = ( bTapLow ? - EPS_SMALL : EPS_SMALL) ;
double dEpsUp = ( bTapUp ? EPS_SMALL : - EPS_SMALL) ;
@@ -1285,8 +1169,9 @@ VolZmap::IntersLineEllipticalCylinder( const Vector3d& vtLineDir, const Point3d&
ptInt2.ToGlob( CircFrame) ;
vtN1.ToGlob( CircFrame) ;
vtN2.ToGlob( CircFrame) ;
return true ;
}
return true ;
return false ;
}
//----------------------------------------------------------------------------