diff --git a/IntersLineSurfBez.cpp b/IntersLineSurfBez.cpp index 34e45fe..7ee2c8e 100644 --- a/IntersLineSurfBez.cpp +++ b/IntersLineSurfBez.cpp @@ -17,46 +17,11 @@ #include "/EgtDev/Include/EGkIntersLineSurfTm.h" #include "/EgtDev/Include/EGkIntersLineSurfBez.h" #include "/EgtDev/Include/EGkSurfBezier.h" -#include "/EgtDev/Extern/Eigen/Dense" #include "DistPointLine.h" #include "CurveLine.h" using namespace std ; -//---------------------------------------------------------------------------- -// raffino i punti di intersezione e recupero le rispettive coordinate nello spazio parametrico -bool -FindParametricCoord( const ISurfTriMesh* pSurfTm ,int nIL, int nT, const Point3d& ptI, Point3d& ptSP) { - // recupero i dati dei vertici del triangolo che fa intersezione - int nVert[3] ; - pSurfTm->GetTriangle( nT, nVert) ; - double dU0, dV0, dU1, dV1,dU2, dV2 ; - pSurfTm->GetVertexParam( nVert[0], dU0, dV0) ; - pSurfTm->GetVertexParam( nVert[1], dU1, dV1) ; - pSurfTm->GetVertexParam( nVert[2], dU2, dV2) ; - Point3d pt0, pt1, pt2 ; - pSurfTm->GetVertex( nVert[0], pt0) ; - pSurfTm->GetVertex( nVert[1], pt1) ; - pSurfTm->GetVertex( nVert[2], pt2) ; - // calcolo approssimativamente le coordinate nello spazio parametrico del punto di intersezione - // quindi prima calcolo la composizione lineare tra i vertici del triangolo per ottenere il punto di intersezione - Eigen::Matrix3d mA ; - mA.col(0) << pt0.x, pt0.y , pt0.z ; - mA.col(1) << pt1.x, pt1.y , pt1.z ; - mA.col(2) << pt2.x, pt2.y , pt2.z ; - Eigen::Vector3d b ( ptI.x, ptI.y, ptI.z) ; - Eigen::Vector3d x = mA.fullPivLu().solve(b) ; - // applico questa composizione alle loro coordinate parametriche - Eigen::Matrix3d mB ; - mB.col(0) << dU0, dV0, 0 ; - mB.col(1) << dU1, dV1, 0 ; - mB.col(2) << dU2, dV2, 0 ; - Eigen::Vector3d ptParam = mB * x ; - ptSP.x = ptParam.x() ; - ptSP.y = ptParam.y() ; - return true ; -} - //---------------------- bool RefineIntersNewton( const Point3d& ptL, const Vector3d& vtL, double dLen, bool bFinite, @@ -162,7 +127,10 @@ IntersLineSurfBz( const Point3d& ptL, const Vector3d& vtL, double dLen, const IS ILSBIVECTOR& vInfo, bool bFinite) { PtrOwner pCL( CreateCurveLine()) ; - pCL->SetPVL(ptL, vtL, dLen) ; + if ( bFinite) + pCL->SetPVL(ptL, vtL, dLen) ; + else + pCL->SetPVL(ptL, vtL, 1e6) ; // verifico linea Vector3d vtDir = vtL ; if ( ! vtDir.Normalize( EPS_ZERO)) @@ -186,7 +154,7 @@ IntersLineSurfBz( const Point3d& ptL, const Vector3d& vtL, double dLen, const IS Point3d ptI, ptI2 ; // devo trovare le intersezioni Point3d ptSP, ptSP2 ; // coordinate parametriche delle soluzioni - FindParametricCoord( pSurfTm, InfoTm.nILTT, InfoTm.nT, InfoTm.ptI, ptSP) ; + pSurfBz->UnprojectPointFromStm( InfoTm.nT, InfoTm.ptI, ptSP, InfoTm.nILTT) ; Point3d ptIBz, ptIBz2 ; if ( ! RefineIntersNewton( ptL, vtL, dLen, bFinite, pSurfBz, ptSP, ptIBz)) { /////// posso provare anche a rilanciare newton con un punto di partenza diverso oppure con una direzione di avvicinamento diversa/////////////////////////////////// @@ -205,8 +173,8 @@ IntersLineSurfBz( const Point3d& ptL, const Vector3d& vtL, double dLen, const IS double dCos2 = 0 ; // eventualmente ripeto tutto per ptI2 ( se ho un'intersezione con sovrapposizione) if ( InfoTm.nILTT == ILTT_SEGM || InfoTm.nILTT == ILTT_SEGM_ON_EDGE ) { - FindParametricCoord( pSurfTm, InfoTm.nILTT, InfoTm.nT, InfoTm.ptI2, ptSP2) ; - if ( !RefineIntersNewton(ptL, vtL, dLen, bFinite, pSurfBz, ptSP2, ptIBz2) ) { + pSurfBz->UnprojectPointFromStm( InfoTm.nT, InfoTm.ptI2, ptSP2, InfoTm.nILTT) ; + if ( ! RefineIntersNewton(ptL, vtL, dLen, bFinite, pSurfBz, ptSP2, ptIBz2) ) { int nVert[3] ; pSurfTm->GetTriangle( InfoTm.nT, nVert) ; double dU0, dV0 ; diff --git a/SurfBezier.cpp b/SurfBezier.cpp index 265236c..e8326a5 100644 --- a/SurfBezier.cpp +++ b/SurfBezier.cpp @@ -29,6 +29,13 @@ #include "/EgtDev/Include/EGkUiUnits.h" #include "/EgtDev/Include/EgtNumUtils.h" #include "/EgtDev/Include/EgtPointerOwner.h" +#include "/EgtDev/Include/EGkIntersPlaneSurfTm.h" +#include "/EgtDev/Include/EGkChainCurves.h" +#include "/EgtDev/Include/EGkIntersLineSurfBez.h" +#include "/EgtDev/Include/EGkDistPointSurfTm.h" +#include "/EgtDev/Extern/Eigen/Dense" +#include "/EgtDev/Include/EGkCurveComposite.h" +#include "/EgtDev/Include/EGkGeoObjSave.h" using namespace std ; @@ -129,15 +136,28 @@ SurfBezier::SetControlPoint( int nInd, const Point3d& ptCtrl, double dW) //---------------------------------------------------------------------------- bool -SurfBezier::SetTrimRegion( const ISurfFlatRegion& sfrTrimReg) +SurfBezier::SetTrimRegion( ISurfFlatRegion& sfrTrimReg, bool bIntersectOrSubtract) { // verifico la regione passata if ( &sfrTrimReg == nullptr || ! sfrTrimReg.IsValid()) return false ; + // se la normale ha z negativa ribalto la superficie, sennò le operazioni di intersect e subtract non funzionano + if ( sfrTrimReg.GetNormVersor().z < 0) + sfrTrimReg.Invert() ; // limito la regione allo spazio parametrico della superficie PtrOwner< ISurfFlatRegion> pSfrTrim( GetSurfFlatRegionRectangle( SBZ_TREG_COEFF * m_nSpanU, SBZ_TREG_COEFF * m_nSpanV)) ; - if ( IsNull( pSfrTrim) || ! pSfrTrim->Intersect( sfrTrimReg) || ! pSfrTrim->IsValid()) - return false ; + // bIntersectOrSubtract == true per ottenere lo spazio parametrico trimmato devo fare l'INTERSEZIONE tra il rettangolo totale e l'area passata + if ( bIntersectOrSubtract) { + if ( IsNull( pSfrTrim) || ! pSfrTrim->Intersect( sfrTrimReg) || ! pSfrTrim->IsValid()) + return false ; + } + // bIntersectOrSubtract == false per ottenere lo spazio parametrico trimmato devo fare la SOTTRAZIONE tra il rettangolo totale e l'area passata + else { + if ( IsNull( pSfrTrim) || ! pSfrTrim->Subtract( sfrTrimReg) || ! pSfrTrim->IsValid()) + return false ; + } + ResetAuxSurf() ; + delete m_pSTM ; // assegno la regione di trim m_bTrimmed = true ; delete m_pTrimReg ; @@ -1411,82 +1431,6 @@ SurfBezier::GetCurveOnVApproxLen( double dU) const return 0 ; return dLen ; } -// -////---------------------------------------------------------------------------- -//const SurfTriMesh* -//SurfBezier::GetAuxSurf( void) const -//{ -// // la superficie deve essere validata -// if ( m_nStatus != OK) { -// ResetAuxSurf() ; -// return nullptr ; -// } -// // se gi� calcolata, la restituisco -// if ( m_pSTM != nullptr) -// return m_pSTM ; -// // costruttore della superficie -// StmFromTriangleSoup stmSoup ; -// if ( ! stmSoup.Start()) -// return nullptr ; -// // definisco il numero degli step in U e in V -// double dMaxLenU = 0 ; -// for ( int j = 0 ; j <= m_nDegV * m_nSpanV ; ++ j) -// dMaxLenU = max( dMaxLenU, GetCurveOnUApproxLen( double( j) / m_nDegV)) ; -// int nStepU = GetSteps( m_nDegU, m_nSpanU, dMaxLenU, 2) ; -// double dMaxLenV = 0 ; -// for ( int i = 0 ; i <= m_nDegU * m_nSpanU ; ++ i) -// dMaxLenV = max( dMaxLenV, GetCurveOnVApproxLen( double( i) / m_nDegU)) ; -// int nStepV = GetSteps( m_nDegV, m_nSpanV, dMaxLenV, 2) ; -// // prima curva isoparametrica (potrebbe essere un solo punto) -// PolyLine PL1 ; -// GetCurveOnU( 0, nStepU, PL1) ; -// bool bSingle1 = ( PL1.GetPointNbr() == 1) ; -// // ciclo sulle isoparametriche -// for ( int i = 1 ; i <= nStepV ; ++ i) { -// // seconda curva isoparametrica (con tanti punti quanti la prima, oppure uno solo) -// double dV = double( i) * m_nSpanV / nStepV ; -// PolyLine PL2 ; -// GetCurveOnU( dV, nStepU, PL2) ; -// bool bSingle2 = ( PL2.GetPointNbr() == 1) ; -// // inserisco i triangoli della striscia nel costruttore della TriMesh -// Point3d ptP1c, ptP2c ; -// Point3d ptP1n, ptP2n ; -// bool bNext = PL1.GetFirstPoint( ptP1c) && PL2.GetFirstPoint( ptP2c) ; -// if ( bNext) { -// if ( bSingle1 && bSingle2) -// bNext = false ; -// if ( bSingle1) -// ptP1n = ptP1c ; -// else -// bNext = bNext && PL1.GetNextPoint( ptP1n) ; -// if ( bSingle2) -// ptP2n = ptP2c ; -// else -// bNext = bNext && PL2.GetNextPoint( ptP2n) ; -// } -// while ( bNext) { -// // eventuale primo triangolo (con base sui correnti e vertice su P2 successivo) -// if ( ! AreSamePointApprox( ptP1c, ptP2c)) -// stmSoup.AddTriangle( ptP2c, ptP1c, ptP2n) ; -// // eventuale secondo triangolo (con vertice su P1 corrente e base sui successivi) -// if ( ! AreSamePointApprox( ptP1n, ptP2n)) -// stmSoup.AddTriangle( ptP1c, ptP1n, ptP2n) ; -// // passo alla successiva coppia -// ptP1c = ptP1n ; -// ptP2c = ptP2n ; -// bNext = ( bSingle1 || PL1.GetNextPoint( ptP1n)) && ( bSingle2 || PL2.GetNextPoint( ptP2n)) ; -// } -// // salvo isoparametrica PL2 in PL1 -// PL1.GetUPointList().swap( PL2.GetUPointList()) ; -// bSingle1 = bSingle2 ; -// } -// // la completo -// if ( ! stmSoup.End()) -// return nullptr ; -// // la salvo -// m_pSTM = GetBasicSurfTriMesh( stmSoup.GetSurf()) ; -// return m_pSTM ; -//} //---------------------------------------------------------------------------- const SurfTriMesh* @@ -1531,6 +1475,9 @@ SurfBezier::GetAuxSurf( void) const //} //// per usare i polygon basic/////////////////// + // salvo i bordi in 3d, che servono in caso si voglia trimmare la superficie DOPO aver costruito la trimesh ausiliaria + Tree.GetEdge3D( m_vPLEdges) ; + // qui non sarebbe male stampare un messaggio di errore nel log se avevo un'area da disegnare ma non sono usciti dei poligoni if ( int(vvPL.size()) == 0) LOG_DBG_ERR( GetEGkLogger(), "ERROR : Bezier Surface couldn't be triangulated, hence wasn't drawn") ; @@ -1618,3 +1565,347 @@ SurfBezier::ResetTrimRegion( void) delete( m_pTrimReg) ; m_pTrimReg = nullptr ; } + +//---------------------------------------------------------------------------- +ICurveComposite* +SurfBezier::UnprojectCurveFromStm( const ICurveComposite* pCC) const +{ + // do per scontato che la compo sia una spezzata, visto che arriva dall'intersezione tra un piano e una trimesh + + const ICurve* pCrv0 = pCC->GetCurve( 0) ; + PolyLine pl ; + Point3d pt3D, pt2D ; pCrv0->GetStartPoint( pt3D) ; + DistPointSurfTm distPtStm0( pt3D, *GetAuxSurf()) ; + // aggiungo il primo punto + int nTriaIndex ; distPtStm0.GetMinDistTriaIndex( nTriaIndex) ; + if ( ! UnprojectPointFromStm( nTriaIndex, pt3D, pt2D)) + return nullptr ; + pl.AddUPoint( 0, pt2D) ; + // aggiungo tutti i successivi + for ( int i = 1 ; i < int( pCC->GetCurveCount()) ; ++i) { + const ICurve* pCrv = pCC->GetCurve( i) ; + pCrv->GetEndPoint( pt3D) ; + DistPointSurfTm distPtStm( pt3D, *GetAuxSurf()) ; + distPtStm.GetMinDistTriaIndex( nTriaIndex) ; + if ( ! UnprojectPointFromStm( nTriaIndex, pt3D, pt2D)) + return nullptr ; + pl.AddUPoint( i, pt2D) ; + } + PtrOwner pCC2D ( CreateCurveComposite()) ; + pCC2D->FromPolyLine( pl) ; + return Release( pCC2D) ; +} + + +//---------------------------------------------------------------------------- +bool +SurfBezier::Cut( const Plane3d& plPlane, bool bSaveOnEq) +{ + // faccio l'intersezione della trimesh ausiliaria con il piano posso ottenere: punti, curve 3d e triangoli( coplanari al piano di taglio) + // i punti li escludo + // le curve 3d le trasformo in curve 2d e le aggiungo alle curve di trim + // accorpo eventuali triangoli adiacenti ed estraggo i loop delle regioni ottenute; questi vengono poi portati in 2d e aggiunti alle curve di trim + + PNTVECTOR vPnt ; + BIPNTVECTOR vBPnt ; + TRIA3DVECTOR vTria ; + IntersPlaneSurfTm( plPlane, *GetAuxSurf(), vPnt, vBPnt, vTria) ; + + // recupero la superficie di trim a cui devo aggiungere tutt i nuovi loop di trim + PtrOwner sfrTrimReg ( GetTrimRegion()) ; + if ( IsNull( sfrTrimReg)) + sfrTrimReg.Set( CreateBasicSurfFlatRegion()) ; + + // concateno le curve 3d + ChainCurves chainC ; + double dToler = EPS_SMALL ; + chainC.Init( false, dToler, int( vBPnt.size())) ; + for ( int i = 0 ; i < int( vBPnt.size()) ; ++ i) { + Vector3d vtDir = vBPnt[i].second - vBPnt[i].first ; + vtDir.Normalize() ; + if ( ! chainC.AddCurve( i + 1, vBPnt[i].first, vtDir, vBPnt[i].second, vtDir)) + return false ; + } + // GESTIONE DELLE CURVE OTTENUTE DALL'INTERSEZIONE + + // recupero i percorsi concatenati + Point3d ptNear = ( vBPnt.empty() ? ORIG : vBPnt[0].first) ; + INTVECTOR vId ; + // Devo controllare se i trim intersecano il bordo dello spazio parametrico ed eventualmente spezzare le curve di trim in 2D tenendo conto della periodicità + // costruisco le curveCompo che rappresetano gli edge + ICRVCOMPOPVECTOR vCCEdge ; + for ( int i= 0 ; i < int( m_vPLEdges.size()); ++i) { + vCCEdge.emplace_back() ; + vCCEdge.back()->FromPolyLine( m_vPLEdges[i]) ; + } + // separo tra loop chiusi, interni allo spazio parametrico e loop passanti che tagliano lo spazio intersecando i bordi + ICRVCOMPOPOVECTOR vpCCOpen ; + ICRVCOMPOPOVECTOR vpCCClosed ; + while ( chainC.GetChainFromNear( ptNear, false, vId)) { + // creo una curva composita + PtrOwner pCrvCompo( CreateCurveComposite()) ; + if ( IsNull( pCrvCompo)) + return false ; + // recupero gli estremi dei segmenti, creo le linee e le inserisco nella composita + bool bAdded = true ; + for ( int i = 0 ; i < int( vId.size()) ; ++ i) { + // creo un segmento di retta + PtrOwner pLine( CreateCurveLine()) ; + if ( IsNull( pLine)) + return false ; + // recupero gli estremi (non vanno mai invertiti per opzione di concatenamento) + int nInd = abs( vId[i]) - 1 ; + Point3d ptStart = ( bAdded ? vBPnt[nInd].first : ptNear) ; + Point3d ptEnd = vBPnt[nInd].second ; + // provo ad accodarlo alla composita + bAdded = ( Dist( ptStart, ptEnd) > dToler / 2 && + pLine->Set( ptStart, ptEnd) && + pCrvCompo->AddCurve( Release( pLine), true, dToler)) ; + ptNear = ( bAdded ? ptEnd : ptStart) ; + } + // se lunghezza curva inferiore a 5 volte la tolleranza, la ignoro e sposto il punto finale nel punto di fine della curva che sto ignorando + double dCrvLen ; + if ( ! pCrvCompo->GetLength( dCrvLen) || dCrvLen < 5. * dToler) + continue ; + // se curva chiusa entro 5 volte la tolleranza ma considerata aperta, la chiudo bene + Point3d ptStart, ptEnd ; + if ( pCrvCompo->GetStartPoint( ptStart) && + pCrvCompo->GetEndPoint( ptEnd) && + AreSamePointEpsilon( ptStart, ptEnd, 5. * dToler) && + ! AreSamePointApprox( ptStart, ptEnd)) { + // porto il punto finale a coincidere esattamente con l'inizio + pCrvCompo->ModifyEnd( ptStart) ; + } + // assegno estrusione come direzione normale al piano + pCrvCompo->SetExtrusion( plPlane.GetVersN()) ; + // unisco segmenti allineati + pCrvCompo->MergeCurves( 0.5 * dToler, ANG_TOL_STD_DEG) ; + + // porto la curva nello spazio parametrico + PtrOwner pCC( UnprojectCurveFromStm( pCrvCompo)) ; + + // le curve aperte le tengo da parte per giuntarle alla fine col bordo + if ( ! pCC->IsClosed()) + vpCCOpen.emplace_back( Release( pCC)) ; + // le curve chiuse le metto tutte insieme subito + else + vpCCClosed.emplace_back( Release( pCC)) ; // da controllare che funzioni anche se passo i loop con orientamenti sbagliati////////// + } + + // creo il nuovo bordo usando i loop aperti + POLYLINEMATRIX vPolygons ; + if ( int(vpCCOpen.size()) != 0 ) { + Tree Tree(Point3d( 0, 0, 0), Point3d( m_nSpanU * SBZ_TREG_COEFF, m_nSpanV * SBZ_TREG_COEFF)) ; + POLYLINEVECTOR vPL ; + for ( int i = 0 ; i < int(vpCCOpen.size()); ++i ) { + PolyLine pl ; + vpCCOpen[i]->ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurveComposite::APL_STD, pl) ; + vPL.emplace_back( pl) ; + } + Tree.AddCutsToRoot( vPL) ; + Tree.CreateCellContour( vPolygons) ; + } + + //GESTIONE DEI TRIANGOLI RISULTANTI DALL'INTERSEZIONE + StmFromTriangleSoup StmFts ; + if ( ! StmFts.Start()) + return GDB_ID_NULL ; + for ( int i = 0 ; i < int( vTria.size()) ; ++ i) + // inserisco il triangolo nella nuova superficie + StmFts.AddTriangle( vTria[i]) ; + // valido la superficie e calcolo le adiacenze + if ( ! StmFts.End()) + return GDB_ID_NULL ; + // se superficie con triangoli + PtrOwner pNewStm( StmFts.GetSurf()) ; + POLYLINEVECTOR vPLTria ; + if ( ! IsNull( pNewStm) && pNewStm->GetTriangleCount() > 0) { + pNewStm->GetLoops( vPLTria) ; + } + + SurfFlatRegionByContours sfrContour ; + // aggiungo i tagli sul bordo + for ( int i = 0 ; i < int( vPolygons.size()); ++i ) { + for ( int j = 0 ; j < int( vPolygons[i].size()); ++j) { + PtrOwner pCC( CreateCurveComposite()) ; + pCC->FromPolyLine( vPolygons[i][j]) ; + sfrContour.AddCurve( Release( pCC)) ; + } + } + // aggiungo i loop chiusi + for ( int i = 0 ; i < int( vpCCClosed.size()); ++i ) + sfrContour.AddCurve( Release( vpCCClosed[i])) ; + // aggiungo loop derivati dai triangoli + for ( int i = 0 ; i < int( vPLTria.size()); ++i ) { + PtrOwner pCC( CreateCurveComposite()) ; + pCC->FromPolyLine( vPLTria[i]) ; + sfrContour.AddCurve( Release( pCC)) ; + } + PtrOwner pSFR( sfrContour.GetSurf()) ; + if ( IsNull( pSFR) || ! pSFR->IsValid()) + return false ; + + // se la superficie ha normale con z negativa la inverto + if ( pSFR->GetNormVersor().z < 0) + pSFR->Invert() ; + + // verifico se la superficie che ho ottenuto è corretta o devo prendere il complementare ( rispetto allo spazio parametrico totale) + // per verificarlo prendo un punto su questa superficie e verifico dove sta il suo corrispettivo 3D rispetto al piano di taglio + int nChunkMin = 0 , nLoopMin = 0 ; + // sono nello spazio parametrico, quindi le aree delle curve possono essere molto grandi + double dAreaMin = 1e30 ; + bool bPos = false ; + for ( int c = 0 ; c < int( pSFR->GetChunkCount()); ++c) { + PtrOwner pSurf( pSFR->CloneChunk( c)) ; + for ( int l = 0 ; l < pSurf->GetLoopCount( 0); ++l) { + PtrOwner pCrv( pSurf->GetLoop( 0, l)) ; + double dArea ; pCrv->GetAreaXY( dArea) ; + if ( abs( dArea) < dAreaMin) { + nChunkMin = c ; + nLoopMin = l ; + dAreaMin = abs( dArea) ; + bPos = dArea > 0 ; + } + } + } + + // aggiorno la superficie di trim + Point3d ptStart ; + Vector3d vtDir ; + PtrOwner pCrv( pSFR->GetLoop( nChunkMin, nLoopMin)) ; + pCrv->GetStartPoint( ptStart) ; + pCrv->GetStartDir( vtDir) ; + vtDir.Rotate( Z_AX, bPos? 90 : -90) ; + PtrOwner pCL( CreateCurveLine()) ; + pCL->SetPVL( ptStart, vtDir, 1e6) ; + IntersCurveCurve icc( *Release( pCL), *Release(pCrv)) ; + // verifico di guardare verso l'interno ( il numero di intersezioni deve essere pari visto che partivo da un punto sulla curva) + if ( icc.GetIntersCount()%2 != 0) + return false ; // o sennò un altro tentativo prima di rinunciare + IntCrvCrvInfo iccInfo ; + icc.GetIntCrvCrvInfo( 1, iccInfo) ; + Point3d ptI = iccInfo.IciA[0].ptI ; + Point3d ptToCheck = ( ptStart + ptI) / 2 ; + Point3d pt3D ; GetPointD1D2( ptToCheck.x / SBZ_TREG_COEFF, ptToCheck.y / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, pt3D) ; + double dDist = DistPointPlane( pt3D, plPlane) ; + + // ho due casi in cui devo invertire(prendere il complementare rispetto allo spazio parametrico) la superficie: + // 1. se il punto è sopra il piano ed era dentro una curva CCW + // 2. se il punto è sotto il piano ed era interno ad una curva CW + + PtrOwner pSrfFR_Copy( pSFR->Clone()) ; + SaveGeoObj( Release(pSrfFR_Copy), "D:\\Temp\\inters\\failed_trim.nge", GDB_SV_BIN) ; + if ( ( dDist > 0 && bPos) || ( dDist < 0 && ! bPos)) + SetTrimRegion( *Release( pSFR), false) ; + else + SetTrimRegion( *Release( pSFR)) ; + + return true ; +} + +//---------------------------------------------------------------------------- +bool +SurfBezier::UnprojectPointFromStm( int nT, const Point3d& ptI, Point3d& ptSP, int nIL) const +{ + // dato un punto sulla trimesh ausiliaria, ne ricavo le coordinate parametriche + const ISurfTriMesh* pSurfTm = GetAuxSurf() ; + // recupero i dati dei vertici del triangolo che fa intersezione + int nVert[3] ; + pSurfTm->GetTriangle( nT, nVert) ; + double dU0, dV0, dU1, dV1,dU2, dV2 ; + pSurfTm->GetVertexParam( nVert[0], dU0, dV0) ; + pSurfTm->GetVertexParam( nVert[1], dU1, dV1) ; + pSurfTm->GetVertexParam( nVert[2], dU2, dV2) ; + Point3d pt0, pt1, pt2 ; + pSurfTm->GetVertex( nVert[0], pt0) ; + pSurfTm->GetVertex( nVert[1], pt1) ; + pSurfTm->GetVertex( nVert[2], pt2) ; + // se l'intersezione era su un vertice restituisco le coordinate parametriche del vertice + if ( nIL == 3 ) { + if ( AreSamePointApprox(ptI, pt0)) + ptSP.Set( dU0, dV0, 0) ; + else if ( AreSamePointApprox(ptI, pt1)) + ptSP.Set( dU1, dV1, 0) ; + else if ( AreSamePointApprox(ptI, pt2)) + ptSP.Set( dU2, dV2, 0) ; + return true ; + } + // calcolo approssimativamente le coordinate nello spazio parametrico del punto di intersezione + // quindi prima calcolo la composizione lineare tra i vertici del triangolo per ottenere il punto di intersezione + Eigen::Matrix3d mA ; + mA.col(0) << pt0.x, pt0.y , pt0.z ; + mA.col(1) << pt1.x, pt1.y , pt1.z ; + mA.col(2) << pt2.x, pt2.y , pt2.z ; + Eigen::Vector3d b ( ptI.x, ptI.y, ptI.z) ; + Eigen::Vector3d x = mA.fullPivLu().solve(b) ; + // applico questa composizione alle loro coordinate parametriche + Eigen::Matrix3d mB ; + mB.col(0) << dU0, dV0, 0 ; + mB.col(1) << dU1, dV1, 0 ; + mB.col(2) << dU2, dV2, 0 ; + Eigen::Vector3d ptParam = mB * x ; + ptSP.x = ptParam.x() ; + ptSP.y = ptParam.y() ; + return true ; +} + +//---------------------------------------------------------------------------- +bool +SurfBezier::UnprojectPoint( const Point3d& pt3D, Point3d& ptParam) const +{ + // dato il punto pt3D sulla superficie di Bezier si cercano le coordinate parametriche ( ptParam) , iterativamente con Newton + // trovato un primo candidato ptParam, ne calcolo l'immagine sulla superficie ( ptBez) e ne calcolo la distanza con il punto pt3D + // ripeto cercando di avvicinarmi il più possibile + // per trovare il primo punto trovo il triangolo della trimesh ausiliaria più vicino e il punto più vicino + DistPointSurfTm dptSurfTm( pt3D, *GetAuxSurf()) ; + int nTriaIndex ; dptSurfTm.GetMinDistTriaIndex( nTriaIndex) ; + Point3d ptI ; dptSurfTm.GetMinDistPoint( ptI) ; + UnprojectPointFromStm( nTriaIndex, ptI, ptParam) ; + Point3d ptBez ; + GetPointD1D2( ptParam.x / SBZ_TREG_COEFF, ptParam.y / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptBez) ; + // usando un algoritmo di newton cerco di avvicinarmi il più possibile al punto + double dDistNew = Dist( pt3D, ptBez) ; + double dDistPre = dDistNew ; + + int nCount = 0 ; + double dh = EPS_SMALL ; + // metodo di newton in più dimensioni + // vario sia il parametro U che il parametro V e verifico se la distanza dalla retta diminuisce per scostamenti positivi o negativi. + while ( dDistNew > EPS_SMALL && nCount < 100) { + dDistPre = dDistNew ; + Point3d ptIBzNew1 ; + GetPointD1D2( ( ptParam.x + dh) / SBZ_TREG_COEFF, ptParam.y / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptIBzNew1) ; + dDistNew = Dist( pt3D, ptIBzNew1) ; + double dfdU = ( dDistNew - dDistPre) / dh ; + Point3d ptIBzNew2 ; + GetPointD1D2( ptParam.x / SBZ_TREG_COEFF, ( ptParam.y + dh) / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptIBzNew2) ; + dDistNew = Dist( pt3D, ptIBzNew2) ; + double dfdV = ( dDistNew - dDistPre) / dh ; + //// opzione 0 + ////scelgo h1 e h2 separatamente e in modo da annullare f(x) + //// opzione 1 + //// valore fisso + //double dr = EPS_SMALL ; + //if ( dDistPre > 1) + // dr = 1 ; + //else if ( dDistPre > 0.1) + // dr = 0.1 ; + //else if ( dDistPre > 0.01) + // dr = 0.01 ; + //// opzione 2 + //// valore direttamente vincolato + //double dr = dDistPre ; + //// opzione 3 + //// valuto la deformazione locale in base allo spostamento del punto sulla bezier // non serve + //double dh1 = Dist( ptIBz, ptIBzNew1) ; + //double dh2 = Dist( ptIBz, ptIBzNew2) ; + // potrei valutare il nuovo spostamento in base all'ultima variazione di dDist + // potrei anche vedere se sto uscendo dal triangolo ( definito nello spazio parametrico) + // mi avvicino cercando di annullare la distanza in un colpo solo + double dr = - dDistPre / ( dfdU + dfdV) ; + GetPointD1D2(( ptParam.x + dr * dfdU) / SBZ_TREG_COEFF, ( ptParam.y + dr * dfdV) / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptBez) ; + dDistNew = Dist( pt3D, ptBez) ; + ++nCount ; + } + return nCount != 99 ; +} diff --git a/SurfBezier.h b/SurfBezier.h index 9f357d1..7e150f4 100644 --- a/SurfBezier.h +++ b/SurfBezier.h @@ -87,7 +87,7 @@ class SurfBezier : public ISurfBezier, public IGeoObjRW bool SetControlPoint( int nIndU, int nIndV, const Point3d& ptCtrl, double dW) override { return SetControlPoint( GetInd( nIndU, nIndV), ptCtrl, dW) ; } bool SetControlPoint( int nInd, const Point3d& ptCtrl, double dW) override ; - bool SetTrimRegion( const ISurfFlatRegion& sfrTrimReg) override ; + bool SetTrimRegion( ISurfFlatRegion& sfrTrimReg, bool bIntersectOrSubtract = true) override ; SurfFlatRegion* GetTrimRegion( void) const override ; bool GetInfo( int& nDegU, int& nDegV, int& nSpanU, int& nSpanV, bool& bIsRat, bool& bTrimmed) const override ; const Point3d& GetControlPoint( int nIndU, int nIndV, bool* pbOk) const override @@ -112,6 +112,10 @@ class SurfBezier : public ISurfBezier, public IGeoObjRW bool GetControlCurveOnV( int nIndU, PolyLine& plCtrlV) const override ; const SurfTriMesh* GetAuxSurf( void) const override ; bool GetLeaves( std::vector>& vLeaves) const override ; + bool UnprojectPointFromStm( int nT, const Point3d& ptI, Point3d& ptSP, int nIL = 5) const override ; + bool UnprojectPoint( const Point3d& pt3D, Point3d& ptParam) const override ; + ICurveComposite* UnprojectCurveFromStm( const ICurveComposite* pCC) const override ; + bool Cut( const Plane3d& plPlane, bool bSaveOnEq = false) override ; public : // IGeoObjRW int GetNgeId( void) const override ; @@ -175,6 +179,7 @@ class SurfBezier : public ISurfBezier, public IGeoObjRW SurfFlatRegion* m_pTrimReg ; // eventuale regione di trim int m_nTempProp[2] ; // vettore proprietà temporanee double m_dTempParam[2] ; // vettore parametri temporanei + mutable POLYLINEVECTOR m_vPLEdges ; // vettore delle polyline degli edge nello spazio 3D } ; //----------------------------------------------------------------------------- diff --git a/Tree.cpp b/Tree.cpp index 9bdfaca..2c63934 100644 --- a/Tree.cpp +++ b/Tree.cpp @@ -49,6 +49,14 @@ Tree::~Tree( void) { } +//---------------------------------------------------------------------------- +Tree::Tree( const Point3d ptBl, const Point3d ptTr) + : m_pSrfBz( nullptr), m_bTrimmed( false), m_bBilinear( false), m_bMulti( false), m_bClosedU( false), m_bClosedV( false), m_bSplitPatches( true), m_bTestMode( false) +{ + Cell cRoot( ptBl, ptTr) ; + m_mTree.insert( pair< int, Cell>( -1, cRoot)) ; +} + //---------------------------------------------------------------------------- void Tree::SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches, const Point3d& ptMin, const Point3d& ptMax) @@ -931,6 +939,7 @@ Tree::Balance() void Tree::GetTopNeigh( int nId, INTVECTOR& vTopNeighs) const { + // le celle restituite sono ordinate per x crescente if ( vTopNeighs.empty()) { if ( m_mTree.at( nId).m_nTop == -2) return ; @@ -1000,6 +1009,7 @@ Tree::GetTopNeigh( int nId, INTVECTOR& vTopNeighs) const vector vCells ; for ( int k : vTopNeighs) vCells.push_back( m_mTree.at( k)) ; + // le celle restituite sono ordinate per x crescente sort( vCells.begin(), vCells.end(), Cell::minorX) ; vTopNeighs.clear() ; for ( Cell c : vCells) @@ -1011,6 +1021,7 @@ Tree::GetTopNeigh( int nId, INTVECTOR& vTopNeighs) const void Tree::GetBottomNeigh( int nId, INTVECTOR& vBottomNeighs) const { + // le celle restituite sono ordinate per x crescente if ( vBottomNeighs.empty()) { if ( m_mTree.at( nId).m_nBottom == -2) return ; @@ -1080,6 +1091,7 @@ Tree::GetBottomNeigh( int nId, INTVECTOR& vBottomNeighs) const vector vCells ; for ( int k : vBottomNeighs) vCells.push_back( m_mTree.at( k)) ; + // le celle restituite sono ordinate per x crescente sort( vCells.begin(), vCells.end(), Cell::minorX) ; vBottomNeighs.clear() ; for ( Cell c : vCells) @@ -1090,6 +1102,7 @@ Tree::GetBottomNeigh( int nId, INTVECTOR& vBottomNeighs) const void Tree::GetLeftNeigh( int nId, INTVECTOR& vLeftNeighs) const { + // le celle restituite sono ordinate per y crescente if ( vLeftNeighs.empty()) { if ( m_mTree.at( nId).m_nLeft == -2) return ; @@ -1159,6 +1172,7 @@ Tree::GetLeftNeigh( int nId, INTVECTOR& vLeftNeighs) const vector vCells ; for ( int k : vLeftNeighs) vCells.push_back( m_mTree.at( k)) ; + // le celle restituite sono ordinate per y crescente sort( vCells.begin(), vCells.end(), Cell::minorY) ; vLeftNeighs.clear() ; for ( Cell c : vCells) @@ -1169,6 +1183,7 @@ Tree::GetLeftNeigh( int nId, INTVECTOR& vLeftNeighs) const void Tree::GetRightNeigh( int nId, INTVECTOR& vRightNeighs) const { + // le celle restituite sono ordinate per y crescente if ( vRightNeighs.empty()) { if ( m_mTree.at( nId).m_nRight == -2) return ; @@ -1238,6 +1253,7 @@ Tree::GetRightNeigh( int nId, INTVECTOR& vRightNeighs) const vector vCells ; for ( int k : vRightNeighs) vCells.push_back( m_mTree.at( k)) ; + // le celle restituite sono ordinate per y crescente sort( vCells.begin(), vCells.end(), Cell::minorY) ; vRightNeighs.clear() ; for ( Cell c : vCells) @@ -1418,9 +1434,16 @@ Tree::GetPolygons( POLYLINEMATRIX& vPolygons) //---------------------------------------------------------------------------- bool -Tree::GetPolygonsBasic( POLYLINEVECTOR& vPolygons) +Tree::GetPolygonsBasic( POLYLINEVECTOR& vPolygons, INTVECTOR& vCells) { - if ( m_vPolygons.empty()) { + // condizioni per il calcolo dei poligoni di base + if ( m_vPolygons.empty() || // se non li ho mai calcolati + ( ! m_vPolygons.empty() && ! vCells.empty()) || // se ho già calcolato dei poligoni ma ne sto chiedendo di un altro gruppo di celle + ( vCells.empty() && m_vPolygons.size() != m_vnLeaves.size())) { // se sto chiedendo i poligoni di tutte le celle, ma i poligoni già calcolati sono meno delle celle + // se non ho dato un elenco di celle in input do per scontato che la chiamata sia per calcolare i poligoni di tutte le foglie + if ( vCells.empty()) + vCells = m_vnLeaves ; + PNTVECTOR vVertices ; INTVECTOR vNeigh ; // setto le celle che sono sul LeftEdge e sul TopEdge @@ -1440,7 +1463,7 @@ Tree::GetPolygonsBasic( POLYLINEVECTOR& vPolygons) bool bBottomRight , bTopLeft ; // scorro lungo tutte le celle leaves e oltre agli angoli della cella aggiungo alla polyline della cella anche i vertici, delle celle adiacenti, // che sono sui lati della cella corrente - for ( int nId : m_vnLeaves) { + for ( int nId : vCells) { vVertices.clear() ; vNeigh.clear() ; vVertices.push_back( m_mTree.at( nId).GetBottomLeft()) ; @@ -3521,3 +3544,153 @@ Tree::OnWhichEdge( int nId, const Point3d& ptToAssign, int& nEdge) const return false ; return true ; } + +//---------------------------------------------------------------------------- +bool +Tree::GetEdge3D( POLYLINEVECTOR& vPLEdges) +{ + INTMATRIX vEdges ; // le righe sono gli edge a partire dallo 0, le colonne sono le celle che compongono quell'edge + // recupero le celle sui quattro bordi + vEdges.emplace_back() ; + GetRootNeigh( 0, vEdges[0]) ; + // le celle sui bordi orizzontali sono ordinate per x o y crescente, ma i io voglio costruire gli edge in senso antiorario a partire dal ptTR, + // quindi devo invertire gli Edge 0 e 1 + reverse( vEdges[0].begin(), vEdges[0].end()) ; + vEdges.emplace_back() ; + GetRootNeigh( 1, vEdges[1]) ; + reverse( vEdges[1].begin(), vEdges[1].end()) ; + vEdges.emplace_back() ; + GetRootNeigh( 2, vEdges[2]) ; + vEdges.emplace_back() ; + GetRootNeigh( 3, vEdges[3]) ; + + // recupero i poligoni base delle celle sui bordi + POLYLINEMATRIX mPL ; + mPL.emplace_back() ; + GetPolygonsBasic( mPL[0], vEdges[0]) ; + mPL.emplace_back() ; + GetPolygonsBasic( mPL[1], vEdges[1]) ; + mPL.emplace_back() ; + GetPolygonsBasic( mPL[2], vEdges[2]) ; + mPL.emplace_back() ; + GetPolygonsBasic( mPL[3], vEdges[3]) ; + + // scorro sui gruppi di polyline che rappresentano i poligoni delle celle lungo un lato + for ( int i = 0 ; i < int( mPL.size()) ; ++i) { + vPLEdges.emplace_back() ; + int nPtCount = 0 ; + // scorro sui poligoni delle celle di un lato + for ( int c = 0 ; c < int( mPL[i].size()) ; ++c) { + Point3d pt ; mPL[i][c].GetFirstPoint( pt) ; + // a seconda del lato controllo di stare scorrendo il poligono prendendo solo i punti su quel lato + if ( i == 0 && AreSamePointApprox(pt, m_mTree.at(vEdges[0][c]).GetTopRight()) ) { + // scorro fino alla fine di quel lato + while ( ! AreSamePointApprox(pt, m_mTree.at(vEdges[0][c]).GetTopLeft()) && mPL[i][c].GetNextPoint( pt)) { + vPLEdges.back().AddUPoint( nPtCount, pt) ; + ++ nPtCount ; + } + } + else if ( i == 1 && AreSamePointApprox(pt, m_mTree.at(vEdges[1][c]).GetTopLeft()) ) { + // scorro fino alla fine di quel lato + while ( ! AreSamePointApprox(pt, m_mTree.at(vEdges[1][c]).GetBottomLeft()) && mPL[i][c].GetNextPoint( pt)) { + vPLEdges.back().AddUPoint( nPtCount, pt) ; + ++ nPtCount ; + } + } + else if ( i == 2 && AreSamePointApprox(pt, m_mTree.at(vEdges[2][c]).GetBottomLeft()) ) { + // scorro fino alla fine di quel lato + while ( ! AreSamePointApprox(pt, m_mTree.at(vEdges[2][c]).GetBottomRight()) && mPL[i][c].GetNextPoint( pt)) { + vPLEdges.back().AddUPoint( nPtCount, pt) ; + ++ nPtCount ; + } + } + else if ( i == 3 && AreSamePointApprox(pt, m_mTree.at(vEdges[3][c]).GetBottomRight()) ) { + // scorro fino alla fine di quel lato + while ( ! AreSamePointApprox(pt, m_mTree.at(vEdges[3][c]).GetTopRight()) && mPL[i][c].GetNextPoint( pt)) { + vPLEdges.back().AddUPoint( nPtCount, pt) ; + ++ nPtCount ; + } + } + } + } + + return true ; +} + +//---------------------------------------------------------------------------- +bool +Tree::AddCutsToRoot( POLYLINEVECTOR& vCuts) +{ + // questa funzione dà per scontato di stare lavorando sulla root di uno spazio parametrico di cui si sta solamente calcolando + // il Contour. La funzione da chiamare dopo questa è la CreateCellContour. + // i tagli sono sempre delle linee che tagliano da parte a parte la cella, quindi sono curve aperte + if ( m_mTree.size() > 1) + return false ; + int nRoot = -1 ; + for ( int i = 0 ; i < int( vCuts.size()); ++i) { + PolyLine pl = vCuts.at( i) ; + m_mTree.at( nRoot).m_vInters.emplace_back() ; + Point3d pt ; pl.GetFirstPoint( pt) ; + int nEdgeIn ; OnWhichEdge( nRoot, pt, nEdgeIn) ; + m_mTree.at( nRoot).m_vInters.back().nIn = nEdgeIn ; + PNTVECTOR vInters ; + vInters.emplace_back( pt) ; + while ( pl.GetNextPoint( pt)) + vInters.emplace_back( pt) ; + pl.GetLastPoint( pt) ; + int nEdgeOut ; OnWhichEdge( nRoot, pt, nEdgeOut) ; + m_mTree.at( nRoot).m_vInters.back().nOut = nEdgeOut ; + } + // chiamo una funzione per renderli coerenti + AdjustCuts() ; + return true ; +} + +//---------------------------------------------------------------------------- +bool +Tree::AdjustCuts( void) +{ + // li riordino per ordine di quali taglio incontrerei percorrendo il bordo della cella a partire da ptTR + sort( m_mTree.at( -1).m_vInters.begin(), m_mTree.at( -1).m_vInters.end(), [](Inters& a, Inters& b){ return Inters::FirstEncounter(a,b) ;}) ; + // ora controllo che le intersezioni che trovo siano ingressi alternati ad uscite, sennò inverto l'intersezione + bool bPreviousWasStart = m_mTree.at( -1).m_vInters.at(0).bSortedbyStart ; + for ( int i = 0 ; i < int( m_mTree.at( -1).m_vInters.size()); ++i) { + if ( m_mTree.at( -1).m_vInters.at(i).bSortedbyStart == bPreviousWasStart) { + reverse( m_mTree.at( -1).m_vInters.at(i).vpt.begin(), m_mTree.at( -1).m_vInters.at(i).vpt.end()) ; + bPreviousWasStart = m_mTree.at( -1).m_vInters.at(i).bSortedbyStart ; + } + else + bPreviousWasStart = ! m_mTree.at( -1).m_vInters.at(i).bSortedbyStart ; + } + return true ; +} + +//---------------------------------------------------------------------------- +bool +Tree::CreateCellContour( POLYLINEMATRIX& vPolygons) +{ + // questa funzione è pensata per essere chiamata dopo la AddCutsToRoot, per creare il poligono di un'unica cella a cui sono stati aggiunti dei tagli + if ( m_mTree.size() > 1) + return false ; + int nRoot = -1 ; + // preparo tutto per poter chiamare la createCellPolygon + m_vnLeaves.push_back( nRoot) ; + INTVECTOR vToCheck( (int) m_mTree.at(nRoot).m_vInters.size()) ; + generate_n( vToCheck.begin(), (int) m_mTree.at(nRoot).m_vInters.size(), generator()) ; + int nPoly = 0 ; + INTVECTOR vnParentChunk ; + PolyLine pl ; + pl.AddUPoint(0, m_mTree.at(nRoot).GetTopRight()) ; + pl.AddUPoint(1, m_mTree.at(nRoot).GetTopLeft()) ; + pl.AddUPoint(2, m_mTree.at(nRoot).GetBottomLeft()) ; + pl.AddUPoint(3, m_mTree.at(nRoot).GetBottomRight()) ; + pl.Close() ; + // ora posso creare il poligono della cella con i tagli + while( (int)vToCheck.size() != 0) { + int nPolyBefore = nPoly ; + CreateCellPolygons( 0, vPolygons, vToCheck, nPoly, vnParentChunk, pl) ; + if ( nPolyBefore == nPoly) + break ; + } + return true ; +} \ No newline at end of file diff --git a/Tree.h b/Tree.h index ea60bb6..0ae3288 100644 --- a/Tree.h +++ b/Tree.h @@ -27,11 +27,12 @@ struct Inters { int nOut ; bool bCCW ; int nChunk ; + bool bSortedbyStart ; // riordino le intersezioni per lato in senso antiorario dal top - // se ho pi� intersezioni che entrano in un lato le riordino considerando che percorro i lati in senso antiorario a partire da ptTR + // se ho più intersezioni che entrano in un lato le riordino considerando che percorro i lati in senso antiorario a partire da ptTR bool operator < ( Inters& b) { - // trovo in che ordine stanno i due strat, tenendo conto anche della possibilit� che siano vertici + // trovo in che ordine stanno i due start, tenendo conto anche della possibilità che siano vertici INTVECTOR vEdges = { 7, 0, 4, 1, 5, 2, 6, 3} ; const auto iter1 = find( vEdges.begin(), vEdges.end(), nIn) ; int nPos1 = std::distance( vEdges.begin(), iter1) ; @@ -52,7 +53,7 @@ struct Inters { pl.Close() ; pl.GetAreaXY( dAreaB) ; } - // se nIn � un vertice sistemo il valore + // se nIn è un vertice sistemo il valore int nEdgeIn = nIn ; if ( nIn > 3) nEdgeIn = nIn - 4 ; @@ -63,6 +64,67 @@ struct Inters { ( bEqIn && nEdgeIn == 2 && vpt[0].x < b.vpt[0].x) || ( bEqIn && nEdgeIn == 3 && vpt[0].y < b.vpt[0].y)) ; } + + static bool FirstEncounter ( Inters& a, Inters& b) + { + // riordino in base al lato toccato, o dall'uscita o dall'ingresso, che viene prima. + // ottengo l'ordine che avrei percorrendo il bordo da ptTR e considerando i loop che incontro, indipendentemente se li incontro nel punto di uscita o ingresso + // nell'intersezione salvo se il taglio è stato ordinato guardando l'ingresso o l'uscita + INTVECTOR vEdges = { 7, 0, 4, 1, 5, 2, 6, 3} ; + // trovo i lati di ingresso e uscita + const auto iter1 = find( vEdges.begin(), vEdges.end(), a.nIn) ; + int nPos1 = std::distance( vEdges.begin(), iter1) ; + const auto iter2 = find( vEdges.begin(), vEdges.end(), a.nOut) ; + int nPos2 = std::distance( vEdges.begin(), iter2) ; + const auto iter3 = find( vEdges.begin(), vEdges.end(), b.nIn) ; + int nPos3 = std::distance( vEdges.begin(), iter3) ; + const auto iter4 = find( vEdges.begin(), vEdges.end(), b.nOut) ; + int nPos4 = std::distance( vEdges.begin(), iter4) ; + int nFirstA = 0 ; + int nFirstB = 0 ; + // salvo l'indice del primo punto dell'intersezione che ho incontrato scorrendo il bordo da ptTR + // salvo il lato che viene prima confrontando ingresso e uscita + if ( nPos2 < nPos1) { + nPos1 = nPos2 ; + nFirstA = int( a.vpt.size()) ; + } + // se ingresso e uscita sono sullo stesso lato allora confronto le coordinate per capire se viene prima l'ingresso o l'uscita + else if ( nPos2 == nPos1 ) { + if ( nPos1 == 0 ) + nFirstA = a.vpt[0].x > a.vpt.back().x ? 0 : int( a.vpt.size()) ; + else if ( nPos1 == 1 ) + nFirstA = a.vpt[0].y > a.vpt.back().y ? 0 : int( a.vpt.size()) ; + else if ( nPos1 == 2 ) + nFirstA = a.vpt[0].x < a.vpt.back().x ? 0 : int( a.vpt.size()) ; + else if ( nPos1 == 3 ) + nFirstA = a.vpt[0].y < a.vpt.back().y ? 0 : int( a.vpt.size()) ; + } + if ( nPos4 < nPos3) { + nPos3 = nPos4 ; + nFirstB = int( b.vpt.size()) ; + } + else if ( nPos4 == nPos3 ) { + if ( nPos3 == 0 ) + nFirstB = b.vpt[0].x > b.vpt.back().x ? 0 : int( b.vpt.size()) ; + else if ( nPos3 == 1 ) + nFirstB = b.vpt[0].y > b.vpt.back().y ? 0 : int( b.vpt.size()) ; + else if ( nPos3 == 2 ) + nFirstB = b.vpt[0].x < b.vpt.back().x ? 0 : int( b.vpt.size()) ; + else if ( nPos3 == 3 ) + nFirstB = b.vpt[0].y < b.vpt.back().y ? 0 : int( b.vpt.size()) ; + } + a.bSortedbyStart = nFirstA == 0 ; + b.bSortedbyStart = nFirstB == 0 ; + // se sono diversi ritorno il confronto + if ( nPos1 != nPos3) + return nPos1 < nPos3 ; + // se sono uguali devo valutare il punto di intersezione + return ( nPos1 == 0 && a.vpt[nFirstA].x > b.vpt[nFirstB].x) || + ( nPos1 == 1 && a.vpt[nFirstA].y > b.vpt[nFirstB].y) || + ( nPos1 == 2 && a.vpt[nFirstA].x < b.vpt[nFirstB].x) || + ( nPos1 == 3 && a.vpt[nFirstA].y < b.vpt[nFirstB].y) ; + } + bool operator == ( Inters& b) { return AreSamePointExact( vpt[0], b.vpt[0]) ; @@ -73,8 +135,8 @@ struct Inters { } } ; // nIn e nOut sono flag che indicano da quale lato ho l'ingresso e l'uscita a partire dal lato top in senso antiorario -// oltre il 3 sono le celle adiacenti in diagonale al vertice-> 4 corrisponde al ptTl e da l� in senso antiorario -// -1 se la curva � sempre dentro la cella +// oltre il 3 sono le celle adiacenti in diagonale al vertice-> 4 corrisponde al ptTl e da lì in senso antiorario +// -1 se la curva è sempre dentro la cella //---------------------------------------------------------------------------- class Cell @@ -116,6 +178,10 @@ class Cell { return m_ptPbl ; } Point3d GetTopRight( void) const { return m_ptPtr ; } + Point3d GetTopLeft( void) const + { return Point3d( m_ptPbl.x, m_ptPtr.y) ; } + Point3d GetBottomRight( void) const + { return Point3d( m_ptPtr.x, m_ptPbl.y); } double GetSplitValue( void) const { return m_dSplit ; } bool IsSplitVert( void) const // se true la cella verrebbe splittata verticalmente, senn� orizzontalmente @@ -138,24 +204,24 @@ class Cell int m_nLeft ; // cella adiacente al lato left int m_nRight ; // cella adiacente al lato right int m_nParent ; // cella genitore - int m_nDepth ; // profondit� della cella rispetto a root - double m_dSplit ; // parametro a cui � stata splittata la cella + int m_nDepth ; // profondità della cella rispetto a root + double m_dSplit ; // parametro a cui è stata splittata la cella int m_nChild1 ; // prima cella figlio int m_nChild2 ; // seconda cella figlio int m_nFlag ; // falg che indica la caratterizzazione della cella rispetto ai loop di trim // 0 esterna, 1 intersecata, 2 contiene un loop, 3 intersecata e contenente un loop, 4 contenuta in un loop - int m_nFlag2 ; // falg che indica se la cella � stata attraversata durante l'ultima fase del labelling - int m_nRightEdgeIn ; // 0 right edge fuori, 1 right edge dentro, 2 met� e met� - bool m_bOnLeftEdge ; // flag che indica se la cella � sul lato sinistro ( per superfici chiuse sul parametro U) - bool m_bOnTopEdge ; // flag che indica se la cella � sul lato top ( per superfici chiuse sul parametro V) + int m_nFlag2 ; // falg che indica se la cella è stata attraversata durante l'ultima fase del labelling + int m_nRightEdgeIn ; // 0 right edge fuori, 1 right edge dentro, 2 metà e metà + bool m_bOnLeftEdge ; // flag che indica se la cella è sul lato sinistro ( per superfici chiuse sul parametro U) + bool m_bOnTopEdge ; // flag che indica se la cella è sul lato top ( per superfici chiuse sul parametro V) std::vector m_vInters ; // vettore delle intersezioni della cella con i loop di trim - // ogni elemento del vettore � l'insieme dei punti che caratterizza un atrtaversamento della cella + // ogni elemento del vettore è l'insieme dei punti che caratterizza un attraversamento della cella private : Point3d m_ptPbl ; // punto bottom left Point3d m_ptPtr ; // punto top right - bool m_bProcessed ; // flag che indica se la cella � stata processata - bool m_bSplitVert ; // flag che indica in quale direzione � stata divisa la cella + bool m_bProcessed ; // flag che indica se la cella è stata processata + bool m_bSplitVert ; // flag che indica in quale direzione è stata divisa la cella } ; //---------------------------------------------------------------------------- @@ -165,16 +231,21 @@ class Tree ~Tree( void) ; Tree( void) ; Tree ( const SurfBezier* pSrfBz, bool bSplitPatches = true, const Point3d& ptMin = ORIG, const Point3d& ptMax = ORIG) ; + Tree( const Point3d ptBl, const Point3d ptTr) ; // creatore da usare solo nel caso in cui si voglia aggiungere tagli ad un'unica cella e del risultato ottenere il contorno void SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches = true, const Point3d& ptMin = ORIG, const Point3d& ptMax = ORIG) ; bool GetIndependentTrees( BIPNTVECTOR& vTrees) ; // calcolo la suddivisione della superficie solo sulle singole bbox dei loop di trim ( unendo quelli vicini) bool BuildTree( double dLinTol = LIN_TOL_STD, double dSideMin = 1, double dSideMax = INFINITO) ; // dSideMax � il massimo per la dimensione maggiore di un triangolo della trimesh // dSideMin � lunghezza minima del lato di una cella nello spazio reale bool BuildTree_test( double dLinTol = LIN_TOL_STD, double dSideMin = 1, double dSideMax = INFINITO) ; bool GetPolygons( POLYLINEMATRIX& vPolygons) ; - bool GetPolygonsBasic( POLYLINEVECTOR& vPolygons) ; // restituisce il poligono corrispondente ad ogni cella foglia dell'albero + bool GetPolygonsBasic( POLYLINEVECTOR& vPolygons, INTVECTOR& vCells = {}) ; // restituisce il poligono corrispondente ad ogni cella foglia dell'albero // ad ogni poligono sono stati aggiunti tutti i vertici dei vicini posizionati sui suoi lati bool GetLeaves ( std::vector& vLeaves) const ; // restituisce gli indici delle foglie nell'albero + bool GetEdge3D ( POLYLINEVECTOR& vPLEdges) ; // restituisce gli edge 3D come curve composite void SetTestMode( void) { m_bTestMode = true ;} ; // attivando la test mode, per la costruzione dell'albero viene usata la funzione BuiltTree_test e viene corretta di conseguenza la FindCell + // funzioni da usare per ricostruire tagli che vanno aggiunti allo spazio parametrico + bool AddCutsToRoot( POLYLINEVECTOR& vCuts) ; // aggiunge i tagli al tree + bool CreateCellContour( POLYLINEMATRIX& vPolygons) ; // crea il nuovo contorno esterno, tenendo conto dei tagli private : bool Split( int nId, double dSplitValue) ; // funzione di split di una cella al parametro indicato nella direzione data da bVert @@ -205,6 +276,8 @@ class Tree bool CategorizeCell( int nId) ; // categorizza la cella in base al flag m_nFlag (dentro, fuori, intersecata) bool CheckIfBetween( const Inters& inA, const Inters& inB) const ; // / controllo se inB è compreso tra l'end e lo start di inA (in senso CCW) bool OnWhichEdge( int nId, const Point3d& ptToAssign, int& nEdge) const ; // indica a quale edge o vertice il punto è vicino entro EPS_SMALL + bool AdjustCuts( void) ; + private : const SurfBezier* m_pSrfBz ; // superficie di bezier