From e78697e5574c213caf8e9e33cd926b585ee9f81f Mon Sep 17 00:00:00 2001 From: Dario Sassi Date: Thu, 25 Aug 2016 15:59:32 +0000 Subject: [PATCH] =?UTF-8?q?EgtGeomKernel=201.6t3=20:=20-=20migliorata=20ve?= =?UTF-8?q?locit=C3=A0=20di=20esecuzione=20di=20FindNearest=20di=20PointGr?= =?UTF-8?q?id3d=20per=20ChainCurves=20-=20a=20PolyLine=20aggiunte=20GetCon?= =?UTF-8?q?vexHullXY=20e=20GetMinAreaRectangleXY.?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- ChainCurves.cpp | 2 +- EgtGeomKernel.rc | Bin 11718 -> 11718 bytes GeoConst.h | 4 +- IntersLineSurfTm.cpp | 2 +- PointGrid3d.cpp | 151 ++++++++++++++++++++------------- PolyLine.cpp | 167 ++++++++++++++++++++++++++++++++++++- SfrCreate.cpp | 6 +- Simple CDSurfFrMove.cpp | 8 +- StmFromCurves.cpp | 2 +- SurfFlatRegion.cpp | 4 +- SurfFlatRegionBooleans.cpp | 2 +- Triangulate.cpp | 2 +- 12 files changed, 273 insertions(+), 77 deletions(-) diff --git a/ChainCurves.cpp b/ChainCurves.cpp index cd024cb..766a9fe 100644 --- a/ChainCurves.cpp +++ b/ChainCurves.cpp @@ -27,7 +27,7 @@ ChainCurves::Init( bool bAllowInvert, double dToler, int nCrvNbrHint) m_sCrvId.rehash( nCrvNbrHint) ; m_vCrvData.clear() ; m_vCrvData.reserve( nCrvNbrHint) ; - const double DIM_CELL = 5.0 ; + const double DIM_CELL = 10.0 ; m_PointGrid.Init( 2 * nCrvNbrHint, DIM_CELL) ; return true ; } diff --git a/EgtGeomKernel.rc b/EgtGeomKernel.rc index 0bda6a3e44c02b7ee4e182019f2a06291b6400b2..cea901bbf75e042fa2d26f50cbef54be71b2908e 100644 GIT binary patch delta 110 zcmX>WeJpyzA2vqg&G-4vGfl1&(wY2&Q;pGhvZAo=W*@FeEMOU1X1ItO delta 110 zcmX>WeJpyzA2vp#&G-4vGfl1&(wY2&Q;pGRvZAo=W*@FeEMOU1X1ItO= nXmin && nX <= nXmax && nY >= nYmin && nY <= nYmax && @@ -255,72 +255,103 @@ bool PointGrid3d::FindNearest( const Point3d& ptTest, INTVECTOR& vnIds) { const double NEAR_TOL = 100 * EPS_SMALL ; + const int SIZE_LIMIT = 1000 ; // pulisco il risultato vnIds.clear() ; - // raggio di ricerca - double dRad = max( m_BBox.DistFromPoint( ptTest), NEAR_TOL) ; - // delta di incremento per raggio di ricerca - const int N_CELL = 10 ; - double dDelta = N_CELL * m_dCellDim ; - // massimo numero di step di ricerca basato su raggio box di ingombro - double dBoxRad ; - 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] ; - bool bFound = false ; - double dMinDist ; - double dSqMinDist ; - for ( int m = 0 ; m <= nMaxSteps ; ++ m) { - // determino il range di celle sui tre assi - IBox iBox ; - if ( ! Get3dRangeNbr( ptTest, dRad, iBox)) - continue ; - // tolgo il range precedente - SubtractIBox( iBox, iBoxPrev, ibR) ; - // ciclo su tutte le celle dei range rimanenti - for ( int l = 0 ; l < N_SUBIBOX ; ++ l) { - for ( int i = ibR[l].nXmin ; i <= ibR[l].nXmax ; ++ i) { - for ( int j = ibR[l].nYmin ; j <= ibR[l].nYmax ; ++ j) { - for ( int k = ibR[l].nZmin ; k <= ibR[l].nZmax ; ++ k) { - // ciclo sui punti della cella - IPNTI_UMMAP_CRANGE MMrange = m_MMap.equal_range( PointHash( i, j, k)) ; - for ( ; MMrange.first != MMrange.second ; ++ MMrange.first) { - // se distanza inferiore al minimo, aggiorno... - double dSqDist = SqDist( (*MMrange.first).second.first, ptTest) ; - // altro punto con la stessa minima distanza gią trovata - if ( bFound && fabs( dSqDist - dSqMinDist) < 2 * dMinDist * NEAR_TOL + NEAR_TOL * NEAR_TOL) { - // inserisco il punto nel vettore dei risultati - vnIds.push_back( (*MMrange.first).second.second) ; - } - // primo punto o punto con minima distanza pił bassa - else if ( ! bFound || dSqDist < dSqMinDist) { - bFound = true ; - // aggiorno i minimi - dSqMinDist = dSqDist ; - dMinDist = sqrt( dSqMinDist) ; - // pulisco il vettore dei risultati ed inserisco questo - vnIds.clear() ; - vnIds.push_back( (*MMrange.first).second.second) ; + // se ci sono pochi elementi, eseguo direttamente la ricerca sugli stessi + if ( m_MMap.size() < SIZE_LIMIT) { + // ciclo di ricerca sugli elementi + bool bFound = false ; + double dMinDist ; + double dSqMinDist ; + for each ( const auto& PntI in m_MMap) { + // quadrato della distanza + double dSqDist = SqDist( ptTest, PntI.second.first) ; + // altro punto con la stessa minima distanza gią trovata + if ( bFound && fabs( dSqDist - dSqMinDist) < 2 * dMinDist * NEAR_TOL + NEAR_TOL * NEAR_TOL) { + // inserisco il punto nel vettore dei risultati + vnIds.push_back( PntI.second.second) ; + } + // primo punto o punto con minima distanza pił bassa + else if ( ! bFound || dSqDist < dSqMinDist) { + bFound = true ; + // aggiorno i minimi + dSqMinDist = dSqDist ; + dMinDist = sqrt( dSqMinDist) ; + // pulisco il vettore dei risultati ed inserisco questo + vnIds.clear() ; + vnIds.push_back( PntI.second.second) ; + } + } + return bFound ; + } + // altrimenti, verifico cerchi concentrici di celle + else { + // raggio di ricerca + double dRad = max( m_BBox.DistFromPoint( ptTest), NEAR_TOL) ; + // delta di incremento per raggio di ricerca + const int N_CELL = 10 ; + double dDelta = N_CELL * m_dCellDim ; + // massimo numero di step di ricerca basato su raggio box di ingombro + double dBoxRad ; + 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] ; + bool bFound = false ; + double dMinDist ; + double dSqMinDist ; + for ( int m = 0 ; m <= nMaxSteps ; ++ m) { + // determino il range di celle sui tre assi + IBox iBox ; + if ( ! Get3dRangeNbr( ptTest, dRad, iBox)) + continue ; + // tolgo il range precedente + SubtractIBox( iBox, iBoxPrev, ibR) ; + // ciclo su tutte le celle dei range rimanenti + for ( int l = 0 ; l < N_SUBIBOX ; ++ l) { + for ( int i = ibR[l].nXmin ; i <= ibR[l].nXmax ; ++ i) { + for ( int j = ibR[l].nYmin ; j <= ibR[l].nYmax ; ++ j) { + for ( int k = ibR[l].nZmin ; k <= ibR[l].nZmax ; ++ k) { + // ciclo sui punti della cella + IPNTI_UMMAP_CRANGE MMrange = m_MMap.equal_range( PointHash( i, j, k)) ; + for ( ; MMrange.first != MMrange.second ; ++ MMrange.first) { + // se distanza inferiore al minimo, aggiorno... + double dSqDist = SqDist( (*MMrange.first).second.first, ptTest) ; + // altro punto con la stessa minima distanza gią trovata + if ( bFound && fabs( dSqDist - dSqMinDist) < 2 * dMinDist * NEAR_TOL + NEAR_TOL * NEAR_TOL) { + // inserisco il punto nel vettore dei risultati + vnIds.push_back( (*MMrange.first).second.second) ; + } + // primo punto o punto con minima distanza pił bassa + else if ( ! bFound || dSqDist < dSqMinDist) { + bFound = true ; + // aggiorno i minimi + dSqMinDist = dSqDist ; + dMinDist = sqrt( dSqMinDist) ; + // pulisco il vettore dei risultati ed inserisco questo + vnIds.clear() ; + vnIds.push_back( (*MMrange.first).second.second) ; + } } } } } } + // se trovato, inutile continuare perchč ci si allontanerą + if ( bFound) + return true ; + // incremento il raggio di ricerca + dRad += dDelta ; + // salvo box di precedente ricerca + iBoxPrev = iBox ; } - // se trovato, inutile continuare perchč ci si allontanerą - if ( bFound) - return true ; - // incremento il raggio di ricerca - dRad += dDelta ; - // salvo box di precedente ricerca - iBoxPrev = iBox ; - } - return bFound ; + return bFound ; + } } //---------------------------------------------------------------------------- diff --git a/PolyLine.cpp b/PolyLine.cpp index 63af6b2..51d5198 100644 --- a/PolyLine.cpp +++ b/PolyLine.cpp @@ -941,4 +941,169 @@ PolyLine::Flatten( double dZ) } } return true ; -} \ No newline at end of file +} + +//---------------------------------------------------------------------------- +bool +PolyLine::GetConvexHullXY( PNTVECTOR& vConvHull) const +{ + int nSize = int( m_lUPoints.size()) ; + if ( nSize == 0) + return false ; + // inserisco i punti in un array ( considero solo x e y, annullo le z) + PNTVECTOR vPnt( nSize) ; + int k = 0 ; + for each ( const auto& Pnt in m_lUPoints) + vPnt[k++] = Point3d( Pnt.first.x, Pnt.first.y, 0) ; + // ordino secondo le X crescenti + std::sort( vPnt.begin(), vPnt.end(), + []( const Point3d& a, const Point3d& b) { return ( a.x < b.x) ; }) ; + // elimino eventuali punti coincidenti + for ( int i = 0 ; i < nSize ; ++ i) { + for ( int j = i + 1 ; j < nSize ; ++ j) { + if ( ( vPnt[j].x - vPnt[i].x) > EPS_SMALL) + break ; + else if ( AreSamePointXYApprox( vPnt[i], vPnt[j])) { + vPnt.erase( vPnt.begin() + j) ; + -- nSize ; + } + } + } + // applico l'algoritmo di Andrew + int j = 0 ; + vConvHull.resize( 2 * nSize) ; + // costruisco la parte inferiore + for ( int i = 0 ; i < nSize ; ++ i) { + while ( j >= 2 && CrossXY( vConvHull[j-1] - vConvHull[j-2], vPnt[i] - vConvHull[j-2]) <= 0) + -- j ; + vConvHull[j++] = vPnt[i] ; + } + // costruisco la parte superiore + for ( int i = nSize - 2, t = j + 1 ; i >= 0 ; -- i) { + while ( j >= t && CrossXY( vConvHull[j-1] - vConvHull[j-2], vPnt[i] - vConvHull[j-2]) <= 0) + -- j ; + vConvHull[j++] = vPnt[i] ; + } + vConvHull.resize( j - 1) ; + + return true ; +} + +//---------------------------------------------------------------------------- +bool +PolyLine::GetMinAreaRectangleXY( Point3d& ptCen, Vector3d& vtAx, double& dLen, double& dHeight) const +{ + // Convex Hull + PNTVECTOR vConvHull ; + if ( ! GetConvexHullXY( vConvHull)) + return false ; + // all points in ConvexHull are different, minimum a segment + int nCount = int( vConvHull.size()) ; + if ( nCount < 2) + return false ; + + // Starting edge nCount-1 -> 0 + int l = 0, m = 0, n = 0 ; + double dMinArea = INFINITO * INFINITO ; + { + // Edge indexes + int i = 0 ; + int j = nCount - 1 ; + // Edge versor + Vector3d vtE0 = vConvHull[i] - vConvHull[j] ; + vtE0.Normalize() ; + // Edge perpendicular versor + Vector3d vtE1 = Vector3d( -vtE0.y, vtE0.x, 0) ; + // Loop through all points to get maximum extents + double dMin0 = INFINITO, dMax0 = - INFINITO, dMin1 = 0, dMax1 = - INFINITO ; + for ( int k = 0 ; k < nCount ; ++ k) { + // Project points onto axes vtE0 and vtE1 and keep track + // of minimum and maximum values along both axes + Vector3d vtDiff = vConvHull[k] - vConvHull[j] ; + double dSca = ScalarXY( vtDiff, vtE0) ; + if ( dSca < dMin0) { + dMin0 = dSca ; + l = k ; + } + if ( dSca > dMax0) { + dMax0 = dSca ; + m = k ; + } + dSca = ScalarXY( vtDiff, vtE1) ; + if ( dSca > dMax1) { + dMax1 = dSca ; + n = k ; + } + } + // Remember area, center and axes + dMinArea = ( dMax0 - dMin0) * ( dMax1 - dMin1) ; + ptCen = vConvHull[j] + 0.5 * (( dMin0 + dMax0) * vtE0 + ( dMin1 + dMax1) * vtE1) ; + vtAx = vtE0 ; + dLen = dMax0 - dMin0 ; + dHeight = dMax1 - dMin1 ; + } + + // Loop through all other edges (j trails i by 1) + for ( int i = 1, j = 0 ; i < nCount ; j = i, ++ i) { + // Get current edge, normalized + Vector3d vtE0 = vConvHull[i] - vConvHull[j] ; + vtE0.Normalize() ; + // Get an axis vtE1 orthogonal to edge vtE0 + Vector3d vtE1 = Vector3d( -vtE0.y, vtE0.x, 0) ; + // Find new min on vtE0 + double dMin0 = ScalarXY( ( vConvHull[l] - vConvHull[j]), vtE0) ; + for ( int k = 0 ; k < nCount ; ++ k) { + int lnext = ( l + 1) % nCount ; + double dMin0next = ScalarXY( (vConvHull[lnext] - vConvHull[j]), vtE0) ; + if ( dMin0next < dMin0) { + dMin0 = dMin0next ; + l = lnext ; + } + else + break ; + } + // Find new max on vtE0 + double dMax0 = ScalarXY( ( vConvHull[m] - vConvHull[j]), vtE0) ; + for ( int k = 0 ; k < nCount ; ++ k) { + int mnext = ( m + 1) % nCount ; + double dMax0next = ScalarXY( (vConvHull[mnext] - vConvHull[j]), vtE0) ; + if ( dMax0next > dMax0) { + dMax0 = dMax0next ; + m = mnext ; + } + else + break ; + } + // Find new min on vtE1 + double dMin1 = 0 ; + // Find new max on vtE1 + double dMax1 = ScalarXY( ( vConvHull[n] - vConvHull[j]), vtE1) ; + for ( int k = 0 ; k < nCount ; ++ k) { + int nnext = ( n + 1) % nCount ; + double dMax1next = ScalarXY( (vConvHull[nnext] - vConvHull[j]), vtE1) ; + if ( dMax1next > dMax1) { + dMax1 = dMax1next ; + n = nnext ; + } + else + break ; + } + double dArea = ( dMax0 - dMin0) * ( dMax1 - dMin1) ; + // If best so far, remember area, center, and axes + if ( dArea < dMinArea) { + dMinArea = dArea ; + ptCen = vConvHull[j] + 0.5 * (( dMin0 + dMax0) * vtE0 + ( dMin1 + dMax1) * vtE1) ; + vtAx = vtE0 ; + dLen = dMax0 - dMin0 ; + dHeight = dMax1 - dMin1 ; + } + } + + // Axis aligned with max dimension + if ( dHeight > dLen) { + vtAx.Rotate( Z_AX, 0, 1) ; + swap( dLen, dHeight) ; + } + + return true ; +} diff --git a/SfrCreate.cpp b/SfrCreate.cpp index 0e0a423..8ea21b8 100644 --- a/SfrCreate.cpp +++ b/SfrCreate.cpp @@ -134,7 +134,7 @@ GetSurfFlatRegionFromFatCurve( ICurve* pCrv, double dRadius, bool bSquareEnds, b // se curva chiusa if ( pCompo1->IsClosed()) { // fondo le curve allineate - pCompo1->MergeCurves( LIN_TOL_SFR, ANG_TOL_STD_DEG) ; + pCompo1->MergeCurves( LIN_TOL_FINE, ANG_TOL_STD_DEG) ; // ne faccio una copia e la inverto PtrOwner pCompo2( GetCurveComposite( pCompo1->Clone())) ; if ( IsNull( pCompo2) || ! pCompo2->Invert()) @@ -166,7 +166,7 @@ GetSurfFlatRegionFromFatCurve( ICurve* pCrv, double dRadius, bool bSquareEnds, b pCompo1->ExtendEndByLen( dRadius) ; } // fondo le curve allineate - pCompo1->MergeCurves( LIN_TOL_SFR, ANG_TOL_STD_DEG) ; + pCompo1->MergeCurves( LIN_TOL_FINE, ANG_TOL_STD_DEG) ; // ne faccio una copia e la inverto PtrOwner pCompo2( GetCurveComposite( pCompo1->Clone())) ; if ( IsNull( pCompo2) || ! pCompo2->Invert()) @@ -284,7 +284,7 @@ SurfFlatRegionByContours::Prepare( void) } // ordino in senso decrescente sull'area sort( m_vArea.begin(), m_vArea.end(), - []( const INDAREA& a, const INDAREA&b) { return fabs( a.second) > fabs( b.second) ; }) ; + []( const INDAREA& a, const INDAREA& b) { return ( fabs( a.second) > fabs( b.second)) ; }) ; return true ; } diff --git a/Simple CDSurfFrMove.cpp b/Simple CDSurfFrMove.cpp index 2ffa45d..986ee10 100644 --- a/Simple CDSurfFrMove.cpp +++ b/Simple CDSurfFrMove.cpp @@ -267,10 +267,10 @@ MySimpleCDSurfFrMove::TranslateCurveNoCollisionCurve( const ICurve* pCrv1, const // altrimenti confronto le approssimazioni con linee delle curve // determino le due polilinee PolyLine PL1 ; - if ( ! pCrv1->ApproxWithLines( LIN_TOL_SFR, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL1)) + if ( ! pCrv1->ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL1)) return false ; PolyLine PL2 ; - if ( ! pCrv2->ApproxWithLines( LIN_TOL_SFR, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL2)) + if ( ! pCrv2->ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL2)) return false ; // ciclo sulle linee della prima polilinea Point3d ptStart1, ptEnd1 ; @@ -437,10 +437,10 @@ MySimpleCDSurfFrMove::RotateCurveNoCollisionCurve( const ICurve* pCrv1, const IC // altrimenti confronto le approssimazioni con linee delle curve // determino le due polilinee PolyLine PL1 ; - if ( ! pCrv1->ApproxWithLines( LIN_TOL_SFR, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL1)) + if ( ! pCrv1->ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL1)) return false ; PolyLine PL2 ; - if ( ! pCrv2->ApproxWithLines( LIN_TOL_SFR, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL2)) + if ( ! pCrv2->ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, PL2)) return false ; // ciclo sulle linee della prima polilinea Point3d ptStart1, ptEnd1 ; diff --git a/StmFromCurves.cpp b/StmFromCurves.cpp index c56451a..dc24414 100644 --- a/StmFromCurves.cpp +++ b/StmFromCurves.cpp @@ -365,7 +365,7 @@ CalcRegionPolyLines( const CICURVEPVECTOR& vpCurve, double dLinTol, vArea.emplace_back( i, - dArea) ; } sort( vArea.begin(), vArea.end(), - []( const INDAREA& a, const INDAREA&b) { return fabs( a.second) > fabs( b.second) ; }) ; + []( const INDAREA& a, const INDAREA& b) { return ( fabs( a.second) > fabs( b.second)) ; }) ; // sposto le polilinee nel vettore da restituire secondo l'ordine vPL.clear() ; vPL.resize( vPLtmp.size()) ; diff --git a/SurfFlatRegion.cpp b/SurfFlatRegion.cpp index 6c349a5..cbf02af 100644 --- a/SurfFlatRegion.cpp +++ b/SurfFlatRegion.cpp @@ -991,7 +991,7 @@ SurfFlatRegion::GetAuxSurf( void) const ICurve* pLoop = GetMyLoop( i, j) ; while ( pLoop != nullptr) { // approssimo con linee a destra per non avere problemi in punti di contatto tra esterni e interni - if ( ! pLoop->ApproxWithLines( LIN_TOL_SFR, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, vPL[j])) + if ( ! pLoop->ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurve::APL_RIGHT, vPL[j])) return nullptr ; pLoop = GetMyLoop( i, ++j) ; } @@ -1163,7 +1163,7 @@ SurfFlatRegion::MyGetCurveClassification( const ICurve& Crv, CRVCVECTOR& ccClass bFound = inTOnM.GetNext( dMin, dMax) ; } sort( ccClass.begin(), ccClass.end(), - []( const CrvClass& a, const CrvClass& b) { return a.dParS < b.dParS ; }) ; + []( const CrvClass& a, const CrvClass& b) { return ( a.dParS < b.dParS) ; }) ; return true ; } diff --git a/SurfFlatRegionBooleans.cpp b/SurfFlatRegionBooleans.cpp index 404657f..868d208 100644 --- a/SurfFlatRegionBooleans.cpp +++ b/SurfFlatRegionBooleans.cpp @@ -376,7 +376,7 @@ SurfFlatRegion::MyNewSurfFromLoops( PCRV_DEQUE& vpLoop) } // ordino in senso decrescente sull'area sort( vArea.begin(), vArea.end(), - []( const INDAREA& a, const INDAREA&b) { return a.second > b.second ; }) ; + []( const INDAREA& a, const INDAREA& b) { return ( a.second > b.second) ; }) ; // determino quanti sono i nuovi cicli esterni (area positiva) int nLoopCount = int( vArea.size()) ; diff --git a/Triangulate.cpp b/Triangulate.cpp index 04f6102..cc9dbc6 100644 --- a/Triangulate.cpp +++ b/Triangulate.cpp @@ -716,7 +716,7 @@ Triangulate::SortInternalLoops( const POLYLINEVECTOR& vPL, INTVECTOR& vOrd) } // ordino vettore in senso decrescente rispetto al massimo sort( vMax.begin() + 1, vMax.end(), - []( const INDMAX& a, const INDMAX&b) { return a.second > b.second ; }) ; + []( const INDMAX& a, const INDMAX& b) { return ( a.second > b.second) ; }) ; // copio indice nel vettore di ordine vOrd.reserve( vPL.size()) ; for ( int i = 0 ; i < int( vPL.size()) ; ++ i)