52a09392c1
- modifiche e correzioni varie sugli Zmap.
314 lines
10 KiB
C++
314 lines
10 KiB
C++
//----------------------------------------------------------------------------
|
|
// EgalTech 2016-2016
|
|
//----------------------------------------------------------------------------
|
|
// File : VolZmap.cpp Data : 03.11.16 Versione : 1.6w1
|
|
// Contenuto : Implementazione della classe Volume Zmap : intersezioni.
|
|
//
|
|
//
|
|
//
|
|
// Modifiche : 03.11.16 LM Creazione modulo.
|
|
//
|
|
//
|
|
//----------------------------------------------------------------------------
|
|
|
|
//--------------------------- Include ----------------------------------------
|
|
#include "stdafx.h"
|
|
#include "VolZmap.h"
|
|
#include "GeoConst.h"
|
|
#include "\EgtDev\Include\EgtNumUtils.h"
|
|
|
|
using namespace std ;
|
|
|
|
//----------------------------------------------------------------------------
|
|
bool
|
|
VolZmap::IntersLineBox( const Point3d& ptP, const Vector3d& vtV,
|
|
const Point3d& ptMin, const Point3d& ptMax, double& dU1, double& dU2)
|
|
{
|
|
// Il box è allineato agli assi
|
|
|
|
dU1 = - INFINITO ;
|
|
dU2 = INFINITO ;
|
|
|
|
// confronto con piani YZ (perpendicolari ad asse X)
|
|
if ( vtV.x > EPS_ZERO) {
|
|
dU1 = max( dU1, ( ptMin.x - ptP.x) / vtV.x) ;
|
|
dU2 = min( dU2, ( ptMax.x - ptP.x) / vtV.x) ;
|
|
}
|
|
else if ( vtV.x < - EPS_ZERO) {
|
|
dU1 = max( dU1, ( ptMax.x - ptP.x) / vtV.x) ;
|
|
dU2 = min( dU2, ( ptMin.x - ptP.x) / vtV.x) ;
|
|
}
|
|
else if ( ptP.x < ptMin.x - EPS_SMALL || ptP.x > ptMax.x + EPS_SMALL)
|
|
return false ;
|
|
|
|
// confronto con piani ZX (perpendicolari ad asse Y)
|
|
if ( vtV.y > EPS_ZERO) {
|
|
dU1 = max( dU1, ( ptMin.y - ptP.y) / vtV.y) ;
|
|
dU2 = min( dU2, ( ptMax.y - ptP.y) / vtV.y) ;
|
|
}
|
|
else if ( vtV.y < - EPS_ZERO) {
|
|
dU1 = max( dU1, ( ptMax.y - ptP.y) / vtV.y) ;
|
|
dU2 = min( dU2, ( ptMin.y - ptP.y) / vtV.y) ;
|
|
}
|
|
else if ( ptP.y < ptMin.y - EPS_SMALL || ptP.y > ptMax.y + EPS_SMALL)
|
|
return false ;
|
|
|
|
// confronto con piani XZ (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) ;
|
|
}
|
|
else if ( vtV.z < - EPS_ZERO) {
|
|
dU1 = max( dU1, ( ptMax.z - ptP.z) / vtV.z) ;
|
|
dU2 = min( dU2, ( ptMin.z - ptP.z) / vtV.z) ;
|
|
}
|
|
else if ( ptP.z < ptMin.z - EPS_SMALL || ptP.z > ptMax.z + EPS_SMALL)
|
|
return false ;
|
|
|
|
return ( dU2 >= dU1) ;
|
|
}
|
|
|
|
//----------------------------------------------------------------------------
|
|
bool
|
|
VolZmap::IntersLineZMapBBox( const Point3d& ptP, const Vector3d& vtV, double& dU1, double& dU2)
|
|
{
|
|
// Punti estremi del box dello Zmap
|
|
Point3d ptMin = ORIG ;
|
|
Point3d ptMax = ptMin + m_nNx * m_dStep * X_AX + m_nNy * m_dStep * Y_AX + m_dMaxZ * Z_AX ;
|
|
|
|
if ( IntersLineBox( ptP, vtV, ptMin, ptMax, dU1, dU2) && ( dU1 > 0 || dU2 > 0)) {
|
|
return true ;
|
|
}
|
|
|
|
else {
|
|
return false ;
|
|
}
|
|
}
|
|
|
|
//----------------------------------------------------------------------------
|
|
bool
|
|
VolZmap::IntersLineDexel( const Point3d& ptP, const Vector3d& vtV, unsigned int nI,
|
|
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_nNx + nI ;
|
|
unsigned int nDexelSize = unsigned int( m_ZValues[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 += 2) {
|
|
// estremi del box del singolo intervallo
|
|
Point3d ptE1( dXmin, dYmin, m_ZValues[nDexelPos][nIndex]) ;
|
|
Point3d ptE2( dXmax, dYmax, m_ZValues[nDexelPos][nIndex+1]) ;
|
|
double dt1, dt2 ;
|
|
if ( IntersLineBox( ptP, vtV, ptE1, ptE2, dt1, dt2)) {
|
|
bInters = true ;
|
|
dU1 = min( dU1, dt1) ;
|
|
dU2 = max( dU2, dt2) ;
|
|
}
|
|
}
|
|
|
|
return bInters ;
|
|
}
|
|
|
|
//----------------------------------------------------------------------------
|
|
bool
|
|
VolZmap::GetDepth( const Point3d& ptPGlob, const Vector3d& vtDir, double& dInLength, double& dOutLength)
|
|
{
|
|
// Porto il raggio nel riferimento intrinseco
|
|
Point3d ptP = ptPGlob ;
|
|
ptP.ToLoc( m_LocalFrame) ;
|
|
Vector3d vtV = vtDir ;
|
|
vtV.ToLoc( m_LocalFrame) ;
|
|
vtV.Normalize() ;
|
|
|
|
// Studio dell'intersezione fra semiretta e BBox dello Zmap
|
|
double dU1, dU2 ;
|
|
bool bTest = IntersLineZMapBBox( ptP, vtV, dU1, dU2) ;
|
|
|
|
// Semiretta esterna al box dello Zmap
|
|
if ( ! bTest) {
|
|
dInLength = - 2 ;
|
|
dOutLength = - 2 ;
|
|
return true ;
|
|
}
|
|
|
|
Point3d ptI, ptF ;
|
|
// Una sola intersezione valida ( punto interno, intersezione valida 2)
|
|
if ( dU1 < 0 && dU2 > 0) {
|
|
ptI = ptP ;
|
|
ptF = ptP + dU2 * vtV ;
|
|
}
|
|
// due soluzioni valide ( punto esterno)
|
|
else {
|
|
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 - 1) ;
|
|
int nIj = Clamp( int( floor( ptI.y / m_dStep)), 0, m_nNy - 1) ;
|
|
int nFi = Clamp( int( floor( ptF.x / m_dStep)), 0, m_nNx - 1) ;
|
|
int nFj = Clamp( int( floor( ptF.y / m_dStep)), 0, m_nNy - 1) ;
|
|
|
|
// Inizializzo distanze
|
|
dInLength = INFINITO ;
|
|
dOutLength = - INFINITO ;
|
|
|
|
// Variazioni
|
|
double dDeltaX = ptF.x - ptI.x ;
|
|
double dDeltaY = ptF.y - ptI.y ;
|
|
|
|
// 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) {
|
|
|
|
// Controllo con nuovo i e j corrente (considero il bordo sinistro del dexel)
|
|
double dU1, dU2 ;
|
|
if ( IntersLineDexel( 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 - 1) ;
|
|
|
|
// Analisi del dexel
|
|
if ( j != OldJ) {
|
|
double dU1, dU2 ;
|
|
if ( IntersLineDexel( ptP, vtV, 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( 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 ;
|
|
double dX = ptI.x + dMoveX ;
|
|
int OldI = i ;
|
|
i = Clamp( int( floor( dX / m_dStep)), 0, m_nNx - 1) ;
|
|
|
|
// Analisi del dexel
|
|
if ( i != OldI) {
|
|
double dU1, dU2 ;
|
|
if ( IntersLineDexel( ptP, vtV, i, j, dU1, dU2)) {
|
|
dInLength = min( dInLength, dU1) ;
|
|
dOutLength = max( dOutLength, dU2) ;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Se non abbiamo incontrato materiale
|
|
if ( dInLength > dOutLength - EPS_SMALL) {
|
|
dInLength = - 2 ;
|
|
dOutLength = - 2 ;
|
|
return true ;
|
|
}
|
|
|
|
// Se parto dall'interno
|
|
if ( dInLength < - EPS_SMALL)
|
|
dInLength = - 1 ;
|
|
|
|
return true ;
|
|
}
|
|
|
|
//----------------------------------------------------------------------------
|
|
bool
|
|
VolZmap::AvoidBox( const Frame3d& frBox, const Vector3d& vtDiag)
|
|
{
|
|
// BBox
|
|
BBox3d b3Box( ORIG, ORIG + vtDiag) ;
|
|
|
|
// lo porto nel riferimento intrinseco dello Zmap
|
|
b3Box.LocToLoc( frBox, m_LocalFrame) ;
|
|
|
|
// BBox dello Zmap nel suo riferimento intrinseco
|
|
BBox3d b3Zmap( ORIG, Point3d( m_nNx * m_dStep, m_nNy * m_dStep, m_dMaxZ)) ;
|
|
|
|
// Se non interferiscono, posso uscire
|
|
BBox3d b3Int ;
|
|
if ( ! b3Zmap.FindIntersection( b3Box, b3Int))
|
|
return true ;
|
|
|
|
// Limiti su indici
|
|
int nStI = Clamp( int( b3Int.GetMin().x / m_dStep), 0, m_nNx -1) ;
|
|
int nEnI = Clamp( int( b3Int.GetMax().x / m_dStep), 0, m_nNx -1) ;
|
|
int nStJ = Clamp( int( b3Int.GetMin().y / m_dStep), 0, m_nNy -1) ;
|
|
int nEnJ = Clamp( int( b3Int.GetMax().y / m_dStep), 0, m_nNy -1) ;
|
|
|
|
// Vettore direzione dei dexel nel riferimento del Box
|
|
Vector3d vtK = Z_AX ; vtK.LocToLoc( m_LocalFrame, frBox) ;
|
|
|
|
// Riferimento intrinseco dei dexel nel riferimento del box
|
|
Point3d ptO = ORIG ; ptO.LocToLoc( m_LocalFrame, frBox) ;
|
|
Vector3d vtX = X_AX ; vtX.LocToLoc( m_LocalFrame, frBox) ;
|
|
Vector3d vtY = Y_AX ; vtY.LocToLoc( m_LocalFrame, frBox) ;
|
|
|
|
// Ciclo di intersezione dei dexel con il BBox
|
|
for ( int i = nStI ; i <= nEnI ; ++ i) {
|
|
|
|
for ( int j = nStJ ; j <= nEnJ ; ++ j) {
|
|
|
|
int nPos = j * m_nNx + i ;
|
|
int nSize = int( m_ZValues[nPos].size()) ;
|
|
if ( nSize == 0)
|
|
continue ;
|
|
|
|
Point3d ptC = ptO + ( i + 0.5) * m_dStep * vtX + ( j + 0.5) * m_dStep * vtY ;
|
|
|
|
double dZmin, dZmax ;
|
|
if ( IntersLineBox( ptC, vtK, ORIG, ORIG + vtDiag, dZmin, dZmax)) {
|
|
|
|
for ( int nIndex = 0 ; nIndex < nSize ; nIndex += 2) {
|
|
if ( m_ZValues[nPos].size() == 0)
|
|
continue ;
|
|
if ( ! ( dZmax < m_ZValues[nPos][nIndex] - EPS_SMALL ||
|
|
dZmin > m_ZValues[nPos][nIndex + 1] + EPS_SMALL))
|
|
return false ;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
return true ;
|
|
}
|