EgtGeomKernel 1.6t3 :

- migliorata velocità di esecuzione di FindNearest di PointGrid3d per ChainCurves - a PolyLine aggiunte GetConvexHullXY e GetMinAreaRectangleXY.
2016-08-25 15:59:32 +00:00
parent 8e7b54e3d5
commit e78697e557
12 changed files with 273 additions and 77 deletions
@@ -941,4 +941,169 @@ PolyLine::Flatten( double dZ)
      }
   }
   return true ;
-}
+}
+
+//----------------------------------------------------------------------------
+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 ;
+}