diff --git a/Color.cpp b/Color.cpp index 2a98971..b972986 100644 --- a/Color.cpp +++ b/Color.cpp @@ -90,8 +90,16 @@ Color GetSurfBackColor( Color cCol) { float fLum = ( cCol.GetRed() + cCol.GetGreen() + cCol.GetBlue()) / 3 ; - return Color( 0.4 * fLum + 0.4 * cCol.GetRed(), - 0.4 * fLum + 0.4 * cCol.GetGreen(), - 0.4 * fLum + 0.4 * cCol.GetBlue(), - cCol.GetAlpha()) ; + // se opache + if ( cCol.GetIntAlpha() > ALPHA_LIM) + return Color( 0.4 * fLum + 0.4 * cCol.GetRed(), + 0.4 * fLum + 0.4 * cCol.GetGreen(), + 0.4 * fLum + 0.4 * cCol.GetBlue(), + cCol.GetAlpha()) ; + // altrimenti semitrasparenti + else + return Color( 0.4 * fLum + 0.6 * cCol.GetRed(), + 0.4 * fLum + 0.6 * cCol.GetGreen(), + 0.4 * fLum + 0.6 * cCol.GetBlue(), + cCol.GetAlpha()) ; } diff --git a/CurveArc.cpp b/CurveArc.cpp index 54ef8d1..7cbc2d7 100644 --- a/CurveArc.cpp +++ b/CurveArc.cpp @@ -16,6 +16,7 @@ #include "CurveArc.h" #include "CurveComposite.h" #include "DistPointArc.h" +#include "PolygonPlane.h" #include "GeoConst.h" #include "GeoObjFactory.h" #include "NgeWriter.h" @@ -782,6 +783,28 @@ CurveArc::GetCenterPoint( Point3d& ptCen) const return true ; } +//---------------------------------------------------------------------------- +bool +CurveArc::GetCentroid( Point3d& ptCen) const +{ + // verifico lo stato + if ( m_nStatus != OK) + return false ; + // approssimo la curva con una polilinea + PolyLine PL ; + if ( ! ApproxWithLines( LIN_TOL_STD, ANG_TOL_STD_DEG, PL)) + return false ; + // calcolo il centro mediante PolygonPlane + Point3d ptP ; + PolygonPlane PolyPlane ; + for ( bool bFound = PL.GetFirstPoint( ptP) ; bFound ; bFound = PL.GetNextPoint( ptP)) + PolyPlane.AddPoint( ptP) ; + if ( PolyPlane.GetCentroid( ptCen)) + return true ; + // se non riuscito, uso il punto medio + return GetMidPoint( ptCen) ; +} + //---------------------------------------------------------------------------- bool CurveArc::GetDir( double dU, Vector3d& vtDir) const diff --git a/CurveArc.h b/CurveArc.h index 8c13ff1..c42602a 100644 --- a/CurveArc.h +++ b/CurveArc.h @@ -64,6 +64,7 @@ class CurveArc : public ICurveArc, public IGeoObjRW virtual bool GetEndPoint( Point3d& ptEnd) const ; virtual bool GetMidPoint( Point3d& ptMid) const ; virtual bool GetCenterPoint( Point3d& ptCen) const ; + virtual bool GetCentroid( Point3d& ptCen) const ; virtual bool GetStartDir( Vector3d& vtDir) const { return GetDir( 0, vtDir) ; } virtual bool GetEndDir( Vector3d& vtDir) const diff --git a/CurveBezier.cpp b/CurveBezier.cpp index 6fdc178..fa1ac2d 100644 --- a/CurveBezier.cpp +++ b/CurveBezier.cpp @@ -18,6 +18,7 @@ #include "DistPointCrvBezier.h" #include "BiArcs.h" #include "GeoConst.h" +#include "PolygonPlane.h" #include "DistPointLine.h" #include "GeoObjFactory.h" #include "NgeWriter.h" @@ -645,6 +646,28 @@ CurveBezier::GetMidPoint( Point3d& ptMid) const return GetPointD1D2( dMid, ptMid) ; } +//---------------------------------------------------------------------------- +bool +CurveBezier::GetCentroid( Point3d& ptCen) const +{ + // verifico lo stato + if ( m_nStatus != OK) + return false ; + // approssimo la curva con una polilinea + PolyLine PL ; + if ( ! ApproxWithLines( LIN_TOL_STD, ANG_TOL_STD_DEG, PL)) + return false ; + // calcolo il centro mediante PolygonPlane + Point3d ptP ; + PolygonPlane PolyPlane ; + for ( bool bFound = PL.GetFirstPoint( ptP) ; bFound ; bFound = PL.GetNextPoint( ptP)) + PolyPlane.AddPoint( ptP) ; + if ( PolyPlane.GetCentroid( ptCen)) + return true ; + // se non riuscito, uso il punto medio + return GetMidPoint( ptCen) ; +} + //---------------------------------------------------------------------------- bool CurveBezier::GetMidDir( Vector3d& vtDir) const diff --git a/CurveBezier.h b/CurveBezier.h index e8d7314..1dd9445 100644 --- a/CurveBezier.h +++ b/CurveBezier.h @@ -66,6 +66,7 @@ class CurveBezier : public ICurveBezier, public IGeoObjRW virtual bool GetMidPoint( Point3d& ptMid) const ; virtual bool GetCenterPoint( Point3d& ptCen) const { return GetMidPoint( ptCen) ; } + virtual bool GetCentroid( Point3d& ptCen) const ; virtual bool GetStartDir( Vector3d& vtDir) const { return ::GetTang( *this, 0, FROM_PLUS, vtDir) ; } virtual bool GetEndDir( Vector3d& vtDir) const diff --git a/CurveComposite.cpp b/CurveComposite.cpp index 8881ee7..a133ee9 100644 --- a/CurveComposite.cpp +++ b/CurveComposite.cpp @@ -18,6 +18,7 @@ #include "CurveLine.h" #include "CurveArc.h" #include "CurveBezier.h" +#include "PolygonPlane.h" #include "GeoConst.h" #include "GeoObjFactory.h" #include "NgeWriter.h" @@ -74,7 +75,7 @@ CurveComposite::AddCurve( const ICurve& cCrv, bool bEndOrStart, double dLinTol) // se curva semplice if ( cCrv.IsSimple()) { // creo una copia della curva - ICurve* pSmplCrv = GetCurve( cCrv.Clone()) ; + ICurve* pSmplCrv = ::GetCurve( cCrv.Clone()) ; if ( pSmplCrv == nullptr) return false ; // inserisco la curva @@ -93,7 +94,7 @@ CurveComposite::AddCurve( const ICurve& cCrv, bool bEndOrStart, double dLinTol) pCrvOri = ( bEndOrStart ? pCrvCompo->GetFirstCurve() : pCrvCompo->GetLastCurve()) ; while ( pCrvOri != nullptr) { // creo una copia della curva - pSmplCrv = GetCurve( pCrvOri->Clone()) ; + pSmplCrv = ::GetCurve( pCrvOri->Clone()) ; if ( pSmplCrv == nullptr) return false ; // inserisco la curva @@ -274,7 +275,7 @@ CurveComposite::FromSplit( const ICurve& cCrv, int nParts) for ( int i = 1 ; i <= nParts ; ++ i) { dStartLen = dEndLen ; dEndLen = dTotLen * i / (double) nParts ; - PtrOwner pSmplCrv( GetCurve( cCrv.Clone())) ; + PtrOwner pSmplCrv( ::GetCurve( cCrv.Clone())) ; if ( IsNull( pSmplCrv)) return false ; if ( ! pSmplCrv->TrimEndAtLen( dEndLen) || @@ -503,7 +504,7 @@ CurveComposite::CopyFrom( const IGeoObj* pGObjSrc) if ( pCC != nullptr) return CopyFrom( *pCC) ; // se sorgente è un'altro tipo di curva - const ICurve* pCrv = GetCurve( pGObjSrc) ; + const ICurve* pCrv = ::GetCurve( pGObjSrc) ; if ( pCrv != nullptr) { Clear() ; pCrv->GetExtrusion( m_VtExtr) ; @@ -944,6 +945,28 @@ CurveComposite::GetMidPoint( Point3d& ptMid) const return GetPointD1D2( dMid, FROM_MINUS, ptMid) ; } +//---------------------------------------------------------------------------- +bool +CurveComposite::GetCentroid( Point3d& ptCen) const +{ + // verifico lo stato + if ( m_nStatus != OK) + return false ; + // approssimo la curva con una polilinea + PolyLine PL ; + if ( ! ApproxWithLines( LIN_TOL_STD, ANG_TOL_STD_DEG, PL)) + return false ; + // calcolo il centro mediante PolygonPlane + Point3d ptP ; + PolygonPlane PolyPlane ; + for ( bool bFound = PL.GetFirstPoint( ptP) ; bFound ; bFound = PL.GetNextPoint( ptP)) + PolyPlane.AddPoint( ptP) ; + if ( PolyPlane.GetCentroid( ptCen)) + return true ; + // se non riuscito, uso il punto medio + return GetMidPoint( ptCen) ; +} + //---------------------------------------------------------------------------- bool CurveComposite::GetMidDir( Vector3d& vtDir) const @@ -1921,6 +1944,20 @@ CurveComposite::LocToLoc( const Frame3d& frOri, const Frame3d& frDest) return true ; } +//---------------------------------------------------------------------------- +const ICurve* +CurveComposite::GetCurve( int nCrv) const +{ + // la curva deve essere validata + if ( m_nStatus != OK) + return nullptr ; + // verifico che l'indice sia nei limiti + if ( nCrv < 0 || nCrv >= int( m_CrvSmplS.size())) + return nullptr ; + // restituisco la curva + return m_CrvSmplS[nCrv] ; +} + //---------------------------------------------------------------------------- const ICurve* CurveComposite::GetFirstCurve( void) const diff --git a/CurveComposite.h b/CurveComposite.h index 0245d76..b66bcf5 100644 --- a/CurveComposite.h +++ b/CurveComposite.h @@ -66,6 +66,7 @@ class CurveComposite : public ICurveComposite, public IGeoObjRW virtual bool GetMidPoint( Point3d& ptMid) const ; virtual bool GetCenterPoint( Point3d& ptCen) const { return GetMidPoint( ptCen) ; } + virtual bool GetCentroid( Point3d& ptCen) const ; virtual bool GetStartDir( Vector3d& vtDir) const { return ::GetTang( *this, 0, FROM_PLUS, vtDir) ; } virtual bool GetEndDir( Vector3d& vtDir) const @@ -125,6 +126,7 @@ class CurveComposite : public ICurveComposite, public IGeoObjRW virtual bool PolygonSide( int nSides, const Point3d& ptStart, const Point3d& ptEnd, const Vector3d& vtN) ; virtual int GetCurveNumber( void) const { return m_nCounter ; } + virtual const ICurve* GetCurve( int nCrv) const ; virtual const ICurve* GetFirstCurve( void) const ; virtual const ICurve* GetNextCurve( void) const ; virtual const ICurve* GetLastCurve( void) const ; diff --git a/CurveLine.h b/CurveLine.h index c3cacea..cda0ebc 100644 --- a/CurveLine.h +++ b/CurveLine.h @@ -65,6 +65,8 @@ class CurveLine : public ICurveLine, public IGeoObjRW virtual bool GetMidPoint( Point3d& ptMid) const ; virtual bool GetCenterPoint( Point3d& ptCen) const { return GetMidPoint( ptCen) ; } + virtual bool GetCentroid( Point3d& ptCen) const + { return GetMidPoint( ptCen) ; } virtual bool GetStartDir( Vector3d& vtDir) const ; virtual bool GetEndDir( Vector3d& vtDir) const { return GetStartDir( vtDir) ; } diff --git a/EgtGeomKernel.rc b/EgtGeomKernel.rc index 64beb71..c27ce80 100644 Binary files a/EgtGeomKernel.rc and b/EgtGeomKernel.rc differ diff --git a/EgtGeomKernel.vcxproj b/EgtGeomKernel.vcxproj index 6e4ec49..058bc95 100644 --- a/EgtGeomKernel.vcxproj +++ b/EgtGeomKernel.vcxproj @@ -289,6 +289,13 @@ copy $(TargetPath) \EgtProg\Dll64 + + false + false + false + false + + @@ -354,6 +361,8 @@ copy $(TargetPath) \EgtProg\Dll64 + + @@ -430,6 +439,7 @@ copy $(TargetPath) \EgtProg\Dll64 + diff --git a/EgtGeomKernel.vcxproj.filters b/EgtGeomKernel.vcxproj.filters index 404c1f4..d7709a1 100644 --- a/EgtGeomKernel.vcxproj.filters +++ b/EgtGeomKernel.vcxproj.filters @@ -252,6 +252,12 @@ File di origine\GeoCreate + + File di origine\GeoInters + + + File di origine\GeoInters + @@ -593,6 +599,15 @@ File di intestazione\Include + + File di intestazione\Include + + + File di intestazione\Include + + + File di intestazione + diff --git a/IdManager.h b/IdManager.h index dbe4f67..6465db5 100644 --- a/IdManager.h +++ b/IdManager.h @@ -48,7 +48,7 @@ class IdManager -- m_nMaxId ; return true ; } int GetNewId( void) - { return ( ++ m_nMaxId) ; } + { return ( m_nMaxId + 1) ; } int GetMaxId( void) const { return m_nMaxId ; } void UpdateMaxId( int nId) diff --git a/IntersLinePlane.cpp b/IntersLinePlane.cpp new file mode 100644 index 0000000..b851794 --- /dev/null +++ b/IntersLinePlane.cpp @@ -0,0 +1,63 @@ +//---------------------------------------------------------------------------- +// EgalTech 2015-2015 +//---------------------------------------------------------------------------- +// File : IntersLinePlane.cpp Data : 18.02.15 Versione : 1.6b7 +// Contenuto : Implementazione della intersezione linea/piano. +// +// +// +// Modifiche : 18.02.15 DS Creazione modulo. +// +// +//---------------------------------------------------------------------------- + +//--------------------------- Include ---------------------------------------- +#include "stdafx.h" +#include "/EgtDev/Include/EGkIntersLinePlane.h" + +//---------------------------------------------------------------------------- +int +IntersLinePlane( const Point3d& ptL1, const Point3d& ptL2, const Plane3d& plPlane, Point3d& ptInt) +{ + Vector3d vtL = ptL2 - ptL1 ; + double dLen = vtL.Len() ; + if ( dLen > EPS_SMALL) + vtL /= dLen ; + else + vtL = V_NULL ; + return IntersLinePlane( ptL1, vtL, dLen, plPlane, ptInt) ; +} + +//---------------------------------------------------------------------------- +int +IntersLinePlane( const Point3d& ptL, const Vector3d& vtL, double dLen, const Plane3d& plPlane, Point3d& ptInt) +{ + Vector3d vtDir = vtL ; + if ( dLen > EPS_SMALL && ! vtDir.Normalize( EPS_ZERO)) + return ILPT_NO ; + // coseno dell'angolo tra direzione linea e normale al piano + double dCosAng = vtDir * plPlane.vtN ; + // posizione del punto rispetto al piano + double dPosL = ( ptL - ORIG) * plPlane.vtN - plPlane.dDist ; + // verifico se linea parallela al piano (ovvero ortogonale alla normale) + if ( fabs( dCosAng) < COS_ORTO_ANG_ZERO) { + // se il punto giace nel piano, vi giace anche la linea + if ( fabs( dPosL) < EPS_SMALL) + return ILPT_INPLANE ; + // alrimenti paralleli e non ci sono intersezioni + else + return ILPT_NO ; + } + // determino il punto di intersezione + double dDistI = - dPosL / dCosAng ; + ptInt = ptL + vtDir * dDistI ; + // intersezione sul primo estremo + if ( fabs( dDistI) < EPS_SMALL) + return ILPT_START ; + else if ( fabs( dDistI - dLen) < EPS_SMALL) + return ILPT_END ; + else if ( dDistI > 0 && dDistI < dLen) + return ILPT_YES ; + else + return ILPT_NO ; +} \ No newline at end of file diff --git a/IntersLineTria.cpp b/IntersLineTria.cpp new file mode 100644 index 0000000..d155694 --- /dev/null +++ b/IntersLineTria.cpp @@ -0,0 +1,65 @@ +//---------------------------------------------------------------------------- +// EgalTech 2015-2015 +//---------------------------------------------------------------------------- +// File : IntersLinePlane.cpp Data : 18.02.15 Versione : 1.6b7 +// Contenuto : Implementazione della intersezione linea/piano. +// +// +// +// Modifiche : 18.02.15 DS Creazione modulo. +// +// +//---------------------------------------------------------------------------- + +//--------------------------- Include ---------------------------------------- +#include "stdafx.h" +#include "ProjPlane.h" +#include "/EgtDev/Include/EGkIntersLineTria.h" +#include "/EgtDev/Include/EGkIntersLinePlane.h" + +//---------------------------------------------------------------------------- +int +IntersLineTria( const Point3d& ptL1, const Point3d& ptL2, const Triangle3d& trTria, Point3d& ptInt) +{ + Vector3d vtL = ptL2 - ptL1 ; + double dLen = vtL.Len() ; + if ( dLen > EPS_SMALL) + vtL /= dLen ; + else + vtL = V_NULL ; + return IntersLineTria( ptL1, vtL, dLen, trTria, ptInt) ; +} + +//---------------------------------------------------------------------------- +int +IntersLineTria( const Point3d& ptL, const Vector3d& vtL, double dLen, const Triangle3d& trTria, Point3d& ptInt) +{ + // determino il piano del triangolo + Plane3d plPlane ; + plPlane.vtN = trTria.GetN() ; + plPlane.dDist = ( trTria.GetP(0) - ORIG) * plPlane.vtN ; + // calcolo l'intersezione tra il segmento di linea e il piano del triangolo + int nRes = IntersLinePlane( ptL, vtL, dLen, plPlane, ptInt) ; + // se non c'è intersezione + if ( nRes == ILPT_NO) + return ILTT_NO ; + // se c'è una intersezione + else if ( nRes == ILPT_START || nRes == ILPT_END || nRes == ILPT_YES) { + int nPTT = PointInTria( ptInt, trTria.GetP(0), trTria.GetP(1), trTria.GetP(2), trTria.GetN()) ; + switch ( nPTT) { + case PTT_OUT : + return ILTT_NO ; + case PTT_VERT : + return ILTT_VERT ; + case PTT_EDGE : + return ILTT_EDGE ; + default : + return ILTT_IN ; + } + } + // se la linea giace nel piano del triangolo + else { + // !!!! DA FARE !!!! + return ILTT_NO ; + } +} diff --git a/PointGrid3d.cpp b/PointGrid3d.cpp index 3912d80..ddf6ad0 100644 --- a/PointGrid3d.cpp +++ b/PointGrid3d.cpp @@ -267,6 +267,7 @@ PointGrid3d::FindNearest( const Point3d& ptTest, INTVECTOR& vnIds) if ( ! m_BBox.GetRadius( dBoxRad)) return false ; int nMaxSteps = int( 2 * ceil( dBoxRad / dDelta)) ; + nMaxSteps = max( nMaxSteps, 1) ; // ciclo di ricerca con raggio crescente IBox iBoxPrev ; IBox ibR[N_SUBIBOX] ; diff --git a/PolygonPlane.cpp b/PolygonPlane.cpp index f52549c..d4e48f3 100644 --- a/PolygonPlane.cpp +++ b/PolygonPlane.cpp @@ -1,13 +1,13 @@ //---------------------------------------------------------------------------- -// EgalTech 2014-2014 +// EgalTech 2014-2015 //---------------------------------------------------------------------------- -// File : PolygonPlane.cpp Data : 12.08.14 Versione : 1.5h3 +// File : PolygonPlane.cpp Data : 23.02.15 Versione : 1.6b7 // Contenuto : Implementazione della classe PolygonPlane. // Calcolo della normale e delle aree con il metodo di Newell. // // // Modifiche : 12.08.14 DS Creazione modulo. -// +// 23.02.15 DS Aggiunta gestione centro geometrico (centroid). // //---------------------------------------------------------------------------- @@ -15,14 +15,22 @@ #include "stdafx.h" #include "PolygonPlane.h" - //---------------------------------------------------------------------------- void PolygonPlane::AddPoint( const Point3d& ptP) { // se è il primo punto (parto da -1 perchè verrà contato alla chiusura) - if ( m_nPntNbr == -1) + if ( m_nPntNbr == -1) { + // inizializzazioni + m_dLenN = 0 ; + m_vtN = V_NULL ; + m_ptMid = ORIG ; + m_dSXy = 0 ; m_dSXz = 0 ; + m_dSYz = 0 ; m_dSYx = 0 ; + m_dSZx = 0 ; m_dSZy = 0 ; + // salvo il punto come primo (per verificare chiusura alla fine) m_ptFirst = ptP ; + } // dal secondo punto posso cominciare ad eseguire le somme parziali else { // Compute normal as being proportional to projected areas of polygon onto the yz, @@ -30,7 +38,18 @@ PolygonPlane::AddPoint( const Point3d& ptP) m_vtN.x += ( m_ptLast.y - ptP.y) * ( m_ptLast.z + ptP.z) ; // projection on yz m_vtN.y += ( m_ptLast.z - ptP.z) * ( m_ptLast.x + ptP.x) ; // projection on xz m_vtN.z += ( m_ptLast.x - ptP.x) * ( m_ptLast.y + ptP.y) ; // projection on xy - m_ptCen += ptP ; + // Very approximate center + m_ptMid += ptP ; + // First moment of Area (momento statico) + double dTmpX = ( m_ptLast.y * ptP.z - ptP.y * m_ptLast.z) ; + m_dSXy += ( m_ptLast.y + ptP.y) * dTmpX ; + m_dSXz += ( m_ptLast.z + ptP.z) * dTmpX ; + double dTmpY = ( m_ptLast.z * ptP.x - ptP.z * m_ptLast.x) ; + m_dSYz += ( m_ptLast.z + ptP.z) * dTmpY ; + m_dSYx += ( m_ptLast.x + ptP.x) * dTmpY ; + double dTmpZ = ( m_ptLast.x * ptP.y - ptP.x * m_ptLast.y) ; + m_dSZx += ( m_ptLast.x + ptP.x) * dTmpZ ; + m_dSZy += ( m_ptLast.y + ptP.y) * dTmpZ ; } // salvo punto e incremento numero di punti m_ptLast = ptP ; @@ -58,11 +77,54 @@ PolygonPlane::Finalize( void) // normalizzo m_vtN /= m_dLenN ; // sistemo baricentro - m_ptCen /= m_nPntNbr ; + m_ptMid /= m_nPntNbr ; } return true ; } +//---------------------------------------------------------------------------- +bool +PolygonPlane::GetCentroid( Point3d& ptCen) +{ + // verifico completamento conti + if ( ! Finalize()) + return false ; + // assegno i dati del centro geometrico + if ( fabs( m_vtN.z) > fabs( m_vtN.x) && + fabs( m_vtN.z) > fabs( m_vtN.y)) { + // calcoli nel piano xy perpendicolare a Z + ptCen.x = m_dSZx / ( 3 * m_vtN.z * m_dLenN) ; + ptCen.y = m_dSZy / ( 3 * m_vtN.z * m_dLenN) ; + ptCen.z = (( m_ptMid.x - ptCen.x) * m_vtN.x + ( m_ptMid.y - ptCen.y) * m_vtN.y + m_ptMid.z * m_vtN.z) / m_vtN.z ; + } + else if ( fabs( m_vtN.x) > fabs( m_vtN.y)) { + // calcoli nel piano yz perpendicolare a X + ptCen.y = m_dSXy / ( 3 * m_vtN.x * m_dLenN) ; + ptCen.z = m_dSXz / ( 3 * m_vtN.x * m_dLenN) ; + ptCen.x = ( m_ptMid.x * m_vtN.x + ( m_ptMid.y - ptCen.y) * m_vtN.y + ( m_ptMid.z - ptCen.z) * m_vtN.z) / m_vtN.x ; + } + else { + // calcoli nel piano zx perpendicolare a Y + ptCen.z = m_dSYz / ( 3 * m_vtN.y * m_dLenN) ; + ptCen.x = m_dSYx / ( 3 * m_vtN.y * m_dLenN) ; + ptCen.y = (( m_ptMid.x - ptCen.x) * m_vtN.x + m_ptMid.y * m_vtN.y + ( m_ptMid.z - ptCen.z) * m_vtN.z) / m_vtN.y ; + } + + return true ; +} + +//---------------------------------------------------------------------------- +bool +PolygonPlane::GetNormal( Vector3d& vtN) +{ + // verifico completamento conti + if ( ! Finalize()) + return false ; + // assegno i dati del centro + vtN = m_vtN ; + return true ; +} + //---------------------------------------------------------------------------- bool PolygonPlane::GetPlane( Plane3d& plPlane) @@ -72,7 +134,7 @@ PolygonPlane::GetPlane( Plane3d& plPlane) return false ; // assegno i dati al piano plPlane.vtN = m_vtN ; - plPlane.dDist = ( ( m_ptCen - ORIG) * plPlane.vtN) ; + plPlane.dDist = ( ( m_ptMid - ORIG) * plPlane.vtN) ; return true ; } diff --git a/PolygonPlane.h b/PolygonPlane.h index 39e06e3..7e7598c 100644 --- a/PolygonPlane.h +++ b/PolygonPlane.h @@ -1,13 +1,13 @@ //---------------------------------------------------------------------------- -// EgalTech 2014-2014 +// EgalTech 2014-2015 //---------------------------------------------------------------------------- -// File : PolygonPlane.h Data : 12.08.14 Versione : 1.5h3 +// File : PolygonPlane.h Data : 23.02.15 Versione : 1.6b7 // Contenuto : Dichiarazione della classe PolygonPlane. // // // // Modifiche : 12.08.14 DS Creazione modulo. -// +// 23.02.15 DS Aggiunta gestione centro geometrico (centroid). // //---------------------------------------------------------------------------- @@ -20,8 +20,13 @@ class PolygonPlane { public : - PolygonPlane( void) : m_nPntNbr( -1), m_dLenN( 0) {} + PolygonPlane( void) + : m_nPntNbr( -1), m_dLenN( 0) {} + void Reset( void) + { m_nPntNbr = -1 ; m_dLenN = 0 ; } void AddPoint( const Point3d& ptP) ; + bool GetCentroid( Point3d& ptCen) ; + bool GetNormal( Vector3d& vtN) ; bool GetPlane( Plane3d& plPlane) ; bool GetArea( double& dArea) ; @@ -34,5 +39,8 @@ class PolygonPlane Point3d m_ptFirst ; // primo punto aggiunto Point3d m_ptLast ; // ultimo punto aggiunto Vector3d m_vtN ; // versore normale - Point3d m_ptCen ; // baricentro + Point3d m_ptMid ; // media dei punti + double m_dSXy, m_dSXz ; // momenti statici dell'area proiettata sul piano yz (perp. a X) + double m_dSYz, m_dSYx ; // momenti statici dell'area proiettata sul piano zx (perp. a Y) + double m_dSZx, m_dSZy ; // momenti statici dell'area proiettata sul piano xy (perp. a Z) } ; \ No newline at end of file diff --git a/ProjPlane.h b/ProjPlane.h new file mode 100644 index 0000000..5320e36 --- /dev/null +++ b/ProjPlane.h @@ -0,0 +1,124 @@ +//---------------------------------------------------------------------------- +// EgalTech 2015-2015 +//---------------------------------------------------------------------------- +// File : ProjPlane.h Data : 20.02.15 Versione : 1.6b7 +// Contenuto : Gestione calcoli su proiezione in piano canonico. +// +// +// +// Modifiche : 20.02.15 DS Creazione modulo. +// +// +//---------------------------------------------------------------------------- + +#pragma once + +#include "/EgtDev/Include/EGkTriangle3d.h" + +//----------------------------------------------------------------------------- +enum PlaneType { PL_NULL = 0, + PL_XY = 1, + PL_YZ = 2, + PL_ZX = 3} ; +enum PntTriaPos { PTT_OUT = 0, // punto esterno al triangolo + PTT_VERT = 1, // punto su un vertice + PTT_EDGE = 2, // punto interno ad un lato + PTT_IN = 3} ; // punto interno + +//---------------------------------------------------------------------------- +// Piano canonico di miglior proiezione (perpendicolare alla componente maggiore della normale) +inline bool +CalcProjPlane( const Vector3d& vtN, int& nPlane, bool& bCCW) +{ + // verifico che la normale non sia nulla + if ( vtN.IsZero()) + return false ; + // proiezione sul piano XY (Nz con valore maggiore) + if ( fabs( vtN.z) > fabs( vtN.x) && + fabs( vtN.z) > fabs( vtN.y)) { + nPlane = PL_XY ; + bCCW = ( vtN.z > 0) ; + } + // proiezione sul piano YZ (Nx con valore maggiore) + else if ( fabs( vtN.x) > fabs( vtN.y)) { + nPlane = PL_YZ ; + bCCW = ( vtN.x > 0) ; + } + // proiezione sul piano ZX (Ny con valore maggiore) + else { + nPlane = PL_ZX ; + bCCW = ( vtN.y > 0) ; + } + return true ; +} + +//---------------------------------------------------------------------------- +// Lunghezza del segmento nel piano di proiezione +inline double +LenSegment( int nPlane, const Point3d& ptA, const Point3d& ptB) +{ + switch ( nPlane) { + default : // PL_XY + return sqrt( ( ptB.x - ptA.x) * ( ptB.x - ptA.x) + ( ptB.y - ptA.y) * ( ptB.y - ptA.y)) ; + case PL_YZ : + return sqrt( ( ptB.y - ptA.y) * ( ptB.y - ptA.y) + ( ptB.z - ptA.z) * ( ptB.z - ptA.z)) ; + case PL_ZX : + return sqrt( ( ptB.z - ptA.z) * ( ptB.z - ptA.z) + ( ptB.x - ptA.x) * ( ptB.x - ptA.x)) ; + } +} + +//---------------------------------------------------------------------------- +// Prodotto vettoriale nel piano di proiezione ( positivo se i 3 punti in ordine CCW) +inline double +TwoArea( int nPlane, const Point3d& ptA, const Point3d& ptB, const Point3d& ptC) +{ + switch ( nPlane) { + default : // PL_XY + return ( ptB.x - ptA.x) * ( ptC.y - ptB.y) - ( ptB.y - ptA.y) * ( ptC.x - ptB.x) ; + case PL_YZ : + return ( ptB.y - ptA.y) * ( ptC.z - ptB.z) - ( ptB.z - ptA.z) * ( ptC.y - ptB.y) ; + case PL_ZX : + return ( ptB.z - ptA.z) * ( ptC.x - ptB.x) - ( ptB.x - ptA.x) * ( ptC.z - ptB.z) ; + } +} + +//---------------------------------------------------------------------------- +inline int +PointInTria( const Point3d& ptP, + const Point3d& ptA, const Point3d& ptB, const Point3d& ptC, const Vector3d& vtN) +{ + // calcolo il piano ottimale di proiezione + int nPlane ; + bool bCCW ; + if ( ! CalcProjPlane( vtN, nPlane, bCCW)) + return PTT_OUT ; + // calcolo la lunghezza dei lati + double dLenAB = LenSegment( nPlane, ptA, ptB) ; + double dLenBC = LenSegment( nPlane, ptB, ptC) ; + double dLenCA = LenSegment( nPlane, ptC, ptA) ; + if ( dLenAB < EPS_SMALL || dLenBC < EPS_SMALL || dLenCA < EPS_SMALL) + return PTT_OUT ; + // determino la posizione del punto rispetto ai lati del triangolo + double dPosAB = TwoArea( nPlane, ptP, ptA, ptB) / dLenAB ; + double dPosBC = TwoArea( nPlane, ptP, ptB, ptC) / dLenBC ; + double dPosCA = TwoArea( nPlane, ptP, ptC, ptA) / dLenCA ; + // se tre valori nulli, triangolo non definito correttamente + if ( fabs( dPosAB) < EPS_SMALL && fabs( dPosBC) < EPS_SMALL && fabs( dPosCA) < EPS_SMALL) + return PTT_OUT ; + // se due valori nulli, sta su un vertice + if ( ( fabs( dPosAB) < EPS_SMALL && fabs( dPosBC) < EPS_SMALL) || + ( fabs( dPosBC) < EPS_SMALL && fabs( dPosCA) < EPS_SMALL) || + ( fabs( dPosCA) < EPS_SMALL && fabs( dPosAB) < EPS_SMALL)) + return PTT_VERT ; + // se un valore nullo e gli altri due con lo stesso segno, sta su un lato + if ( ( fabs( dPosAB) < EPS_SMALL && ( dPosBC * dPosCA > 0)) || + ( fabs( dPosBC) < EPS_SMALL && ( dPosCA * dPosAB > 0)) || + ( fabs( dPosCA) < EPS_SMALL && ( dPosAB * dPosBC > 0))) + return PTT_EDGE ; + // se tre valori con lo stesso segno, sta dentro + if ( ( dPosAB > EPS_SMALL && dPosBC > EPS_SMALL && dPosCA > EPS_SMALL) || + ( dPosAB < - EPS_SMALL && dPosBC < - EPS_SMALL && dPosCA < - EPS_SMALL)) + return PTT_IN ; + // deve essere esterno + return PTT_OUT ; +} diff --git a/SurfTriMesh.cpp b/SurfTriMesh.cpp index c4db10b..4c8c49e 100644 --- a/SurfTriMesh.cpp +++ b/SurfTriMesh.cpp @@ -202,6 +202,18 @@ SurfTriMesh::GetTriangleNum( void) const return ( GetTriangleSize() - nErased) ; } +//---------------------------------------------------------------------------- +bool +SurfTriMesh::GetVertex( int nId, Point3d& ptP) const +{ + // verifico esistenza del vertice + if ( nId < 0 || nId >= GetVertexSize() || m_vVert[nId].nIdTria == SVT_DEL) + return false ; + // recupero i dati + ptP = m_vVert[nId].ptP ; + return true ; +} + //---------------------------------------------------------------------------- int SurfTriMesh::GetFirstVertex( Point3d& ptP) const @@ -226,6 +238,20 @@ SurfTriMesh::GetNextVertex( int nId, Point3d& ptP) const return nId ; } +//---------------------------------------------------------------------------- +bool +SurfTriMesh::GetTriangle( int nId, int nIdVert[3]) const +{ + // verifico esistenza del triangolo + if ( nId < 0 || nId >= GetTriangleSize() || m_vTria[nId].nIdVert[0] == SVT_DEL) + return false ; + // recupero i dati + nIdVert[0] = m_vTria[nId].nIdVert[0] ; + nIdVert[1] = m_vTria[nId].nIdVert[1] ; + nIdVert[2] = m_vTria[nId].nIdVert[2] ; + return true ; +} + //---------------------------------------------------------------------------- int SurfTriMesh::GetFirstTriangle( int nIdVert[3]) const @@ -252,6 +278,21 @@ SurfTriMesh::GetNextTriangle( int nId, int nIdVert[3]) const return nId ; } +//---------------------------------------------------------------------------- +bool +SurfTriMesh::GetTriangle( int nId, Triangle3d& Tria) const +{ + // verifico esistenza del triangolo + if ( nId < 0 || nId >= GetTriangleSize() || m_vTria[nId].nIdVert[0] == SVT_DEL) + return false ; + // recupero i dati + Tria.Set( m_vVert[m_vTria[nId].nIdVert[0]].ptP, + m_vVert[m_vTria[nId].nIdVert[1]].ptP, + m_vVert[m_vTria[nId].nIdVert[2]].ptP, + m_vTria[nId].vtN) ; + return true ; +} + //---------------------------------------------------------------------------- int SurfTriMesh::GetFirstTriangle( Triangle3d& Tria) const @@ -279,6 +320,88 @@ SurfTriMesh::GetNextTriangle( int nId, Triangle3d& Tria) const return nId ; } +//---------------------------------------------------------------------------- +int +SurfTriMesh::GetAllTriaAroundVertex( int nV, INTVECTOR& vT, bool& bCirc) const +{ + const int MAX_VT_SIZE = 512 ; + + // pulisco il vettore risultato + vT.clear() ; + + // verifico esistenza del vertice + if ( nV < 0 || nV >= GetVertexSize() || m_vVert[nV].nIdTria == SVT_DEL) + return false ; + + // recupero il triangolo puntato dal vertice + int nT = m_vVert[nV].nIdTria ; + if ( nT == SVT_NULL) + return 0 ; + vT.push_back( nT) ; + + // cerco i triangoli adiacenti con lo stesso vertice in CCW + int k ; + int nTa = nT ; + do { + if ( FindVertexInTria( nV, nTa, k)) + nTa = m_vTria[nTa].nIdAdjac[Prev(k)] ; + else + nTa = SVT_NULL ; + // per evitare cicli infiniti dovuti a triangoli invertiti + if ( vT.size() >= 2 && nTa == vT[vT.size()-2]) + break ; + // se valido + if ( nTa != nT && nTa != SVT_NULL) + vT.push_back( nTa) ; + } while ( nTa != nT && nTa != SVT_NULL && vT.size() < MAX_VT_SIZE) ; + + // se sono ritornato al triangolo di partenza ho fatto un giro e concluso la ricerca + if ( nTa == nT) { + bCirc = true ; + return int( vT.size()) ; + } + + // altrimenti, devo cercare i triangoli adiacenti con lo stesso vertice in CW + nTa = nT ; + do { + if ( FindVertexInTria( nV, nTa, k)) + nTa = m_vTria[nTa].nIdAdjac[k] ; + else + nTa = SVT_NULL ; + // per evitare cicli infiniti dovuti a triangoli invertiti + if ( vT.size() >= 2 && nTa == vT[vT.size()-2]) + break ; + // se valido + if ( nTa != nT && nTa != SVT_NULL) + vT.insert( vT.begin(), nTa) ; + } while ( nTa != nT && nTa != SVT_NULL && vT.size() < MAX_VT_SIZE) ; + + bCirc = ( nTa == nT) ; + return int( vT.size()) ; +} + +//---------------------------------------------------------------------------- +bool +SurfTriMesh::GetVertexSmoothNormal( int nV, int nT, Vector3d& vtN) const +{ + // verifico esistenza del vertice + if ( nV < 0 || nV >= GetVertexSize() || m_vVert[nV].nIdTria == SVT_DEL) + return false ; + // verifico esistenza del triangolo + if ( nT < 0 || nT >= GetTriangleSize() || m_vTria[nT].nIdVert[0] == SVT_DEL) + return false ; + // recupero la normale del vertice + for ( int i = 0 ; i < 3 ; ++ i) { + if ( nV == m_vTria[nT].nIdVert[i]) { + if ( ! GetTriangleSmoothNormal( nT, i, vtN)) + vtN = m_vTria[nT].vtN ; + return true ; + } + } + // il vertice non appartiene al triangolo + return false ; +} + //---------------------------------------------------------------------------- bool SurfTriMesh::GetTriangleBoundaryEdges( int nId, TriFlags3d& TFlags) const @@ -614,62 +737,6 @@ SurfTriMesh::FindVertexInTria( int nV, int nT, int& nK) const return ( nK != -1) ; } -//---------------------------------------------------------------------------- -int -SurfTriMesh::GetAllTriaAroundVertex( int nV, INTVECTOR& vT, bool& bCirc) const -{ - const int MAX_VT_SIZE = 512 ; - - // pulisco il vettore risultato - vT.clear() ; - - // recupero il triangolo puntato dal vertice - int nT = m_vVert[nV].nIdTria ; - if ( nT == SVT_NULL) - return 0 ; - vT.push_back( nT) ; - - // cerco i triangoli adiacenti con lo stesso vertice in CCW - int k ; - int nTa = nT ; - do { - if ( FindVertexInTria( nV, nTa, k)) - nTa = m_vTria[nTa].nIdAdjac[Prev(k)] ; - else - nTa = SVT_NULL ; - // per evitare cicli infiniti dovuti a triangoli invertiti - if ( vT.size() >= 2 && nTa == vT[vT.size()-2]) - break ; - // se valido - if ( nTa != nT && nTa != SVT_NULL) - vT.push_back( nTa) ; - } while ( nTa != nT && nTa != SVT_NULL && vT.size() < MAX_VT_SIZE) ; - - // se sono ritornato al triangolo di partenza ho fatto un giro e concluso la ricerca - if ( nTa == nT) { - bCirc = true ; - return int( vT.size()) ; - } - - // altrimenti, devo cercare i triangoli adiacenti con lo stesso vertice in CW - nTa = nT ; - do { - if ( FindVertexInTria( nV, nTa, k)) - nTa = m_vTria[nTa].nIdAdjac[k] ; - else - nTa = SVT_NULL ; - // per evitare cicli infiniti dovuti a triangoli invertiti - if ( vT.size() >= 2 && nTa == vT[vT.size()-2]) - break ; - // se valido - if ( nTa != nT && nTa != SVT_NULL) - vT.insert( vT.begin(), nTa) ; - } while ( nTa != nT && nTa != SVT_NULL && vT.size() < MAX_VT_SIZE) ; - - bCirc = ( nTa == nT) ; - return int( vT.size()) ; -} - //---------------------------------------------------------------------------- int SurfTriMesh::NextIndAroundVertex( int nInd, int nSize, bool bCirc) const diff --git a/SurfTriMesh.h b/SurfTriMesh.h index 8b2f1aa..79bef44 100644 --- a/SurfTriMesh.h +++ b/SurfTriMesh.h @@ -132,12 +132,17 @@ class SurfTriMesh : public ISurfTriMesh, public IGeoObjRW { return m_dLinTol ; } virtual double GetSmoothAngle( void) { return m_dSmoothAng ; } + virtual bool GetVertex( int nId, Point3d& ptP) const ; virtual int GetFirstVertex( Point3d& ptP) const ; virtual int GetNextVertex( int nId, Point3d& ptP) const ; + virtual bool GetTriangle( int nId, int nIdVert[3]) const ; virtual int GetFirstTriangle( int nIdVert[3]) const ; virtual int GetNextTriangle( int nId, int nIdVert[3]) const ; + virtual bool GetTriangle( int nId, Triangle3d& Tria) const ; virtual int GetFirstTriangle( Triangle3d& Tria) const ; virtual int GetNextTriangle( int nId, Triangle3d& Tria) const ; + virtual int GetAllTriaAroundVertex( int nV, INTVECTOR& vT, bool& bCirc) const ; + virtual bool GetVertexSmoothNormal( int nV, int nT, Vector3d& vtN) const ; virtual bool GetTriangleBoundaryEdges( int nId, TriFlags3d& TFlags) const ; virtual bool GetTriangleSmoothNormals( int nId, TriNormals3d& TNrms) const ; @@ -171,7 +176,6 @@ class SurfTriMesh : public ISurfTriMesh, public IGeoObjRW int Prev( int i) const { return (( i + 2) % 3) ; } bool FindVertexInTria( int nV, int nT, int& nK) const ; - int GetAllTriaAroundVertex( int nV, INTVECTOR& vT, bool& bCirc) const ; int NextIndAroundVertex( int nInd, int nSize, bool bCirc) const ; int PrevIndAroundVertex( int nInd, int nSize, bool bCirc) const ; bool GetTriangleSmoothNormal( int nT, int nV, Vector3d& vtN) const ; diff --git a/Triangulate.cpp b/Triangulate.cpp index 6173cb0..04f6102 100644 --- a/Triangulate.cpp +++ b/Triangulate.cpp @@ -15,13 +15,13 @@ #include "stdafx.h" #include #include "Triangulate.h" +#include "ProjPlane.h" #include "\EgtDev\Include\EGkPolyLine.h" #include "\EgtDev\Include\EGkPlane3d.h" using namespace std ; //---------------------------------------------------------------------------- -static enum PlaneType { PL_XY = 1, PL_YZ = 2, PL_ZX = 3} ; static enum EarStatus{ EAS_NULL = -1, EAS_NO = 0, EAS_OK = 1} ; //---------------------------------------------------------------------------- @@ -40,20 +40,10 @@ Triangulate::Make( const PolyLine& PL, PNTVECTOR& vPt, INTVECTOR& vTr) Plane3d plPlane ; if ( ! PL.IsClosedAndFlat( plPlane, dArea, 50 * EPS_SMALL)) return false ; + // determino il piano ottimale di proiezione e il relativo senso di rotazione bool bCCW ; - if ( fabs( plPlane.vtN.z) >= fabs( plPlane.vtN.x) && - fabs( plPlane.vtN.z) >= fabs( plPlane.vtN.y)) { - m_nPlane = PL_XY ; - bCCW = ( plPlane.vtN.z > 0) ; - } - else if ( fabs( plPlane.vtN.x) >= fabs( plPlane.vtN.y)) { - m_nPlane = PL_YZ ; - bCCW = ( plPlane.vtN.x > 0) ; - } - else { - m_nPlane = PL_ZX ; - bCCW = ( plPlane.vtN.y > 0) ; - } + if ( ! CalcProjPlane( plPlane.vtN, m_nPlane, bCCW)) + return false ; // riempio il vettore con i vertici del poligono da triangolare vPt.clear() ; vPt.reserve( PL.GetPointNbr() - 1) ; @@ -107,20 +97,10 @@ Triangulate::Make( const POLYLINEVECTOR& vPL, PNTVECTOR& vPt, INTVECTOR& vTr) Plane3d plExtPlane ; if ( ! vPL[0].IsClosedAndFlat( plExtPlane, dArea, 50 * EPS_SMALL)) return false ; + // determino il piano ottimale di proiezione e il relativo senso di rotazione bool bCCW ; - if ( fabs( plExtPlane.vtN.z) >= fabs( plExtPlane.vtN.x) && - fabs( plExtPlane.vtN.z) >= fabs( plExtPlane.vtN.y)) { - m_nPlane = PL_XY ; - bCCW = ( plExtPlane.vtN.z > 0) ; - } - else if ( fabs( plExtPlane.vtN.x) >= fabs( plExtPlane.vtN.y)) { - m_nPlane = PL_YZ ; - bCCW = ( plExtPlane.vtN.x > 0) ; - } - else { - m_nPlane = PL_ZX ; - bCCW = ( plExtPlane.vtN.y > 0) ; - } + if ( ! CalcProjPlane( plExtPlane.vtN, m_nPlane, bCCW)) + return false ; // verifico le altre polilinee for ( int i = 1 ; i < int( vPL.size()) ; ++i) { // deve essere chiusa, giacere nello stesso piano ed essere orientata al contrario della esterna