EgtGeomKernel :

- codice di test per Zmap.

CDeBoxClosedSurfTm.cpp

@@ -1,13 +1,14 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2016-2020
+//  EgalTech 2016-2024
 //----------------------------------------------------------------------------
-// File : CDBoxSurfTm.cpp      Data : 09.01.20        Versione : 2.2a2
+// File : CDeBoxClosedSurfTm.cpp      Data : 24.03.24        Versione : 2.6c2
 // Contenuto : Implementazione della verifica di collisione tra
 //             BoundingBox e Closed SurftriMesh.
 //
 //
 // Modifiche : 05.10.16 DS Creazione modulo.
 //             09.01.20 DS Cambio nome alla funzione.
+//             24.03.24 DS Aggiunta TestCylSurfTm.
 //
 //----------------------------------------------------------------------------
 
@@ -66,3 +67,39 @@ CDeBoxClosedSurfTm( const Frame3d& frBox, const Vector3d& vtDiag, const ISurfTri
   // Se il box è interno c'è collisione
    return DistBoxCenSurfCalc.IsPointInside() ;
 }
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra il box e la superficie : restituisce true in caso di interferenza.
+//----------------------------------------------------------------------------
+bool
+TestBoxSurfTm( const Frame3d& frBox, const Vector3d& vtDiag, const ISurfTriMesh& Stm, double dSafeDist)
+{
+  // Se il box non è ben definito non ha senso proseguire
+   if ( vtDiag.IsSmall())
+      return true ;
+  // Se superficie non valida, non ha senso proseguire
+   if ( ! Stm.IsValid())
+      return true ;
+  // Recupero BBox del poliedro
+   BBox3d b3Poly = Stm.GetAllTriaBox() ;
+  // Calcolo il BBox del parallelepipedo
+   BBox3d b3BoxL( ORIG, ORIG + vtDiag) ;
+   if ( dSafeDist > EPS_SMALL)
+      b3BoxL.Expand( dSafeDist) ;
+   BBox3d b3Box = GetToGlob( b3BoxL, frBox) ;
+  // Se i BBox non interferiscono, non c'è collisione
+   if ( ! b3Poly.Overlaps( b3Box)  ||  ! b3Poly.Overlaps( frBox, b3BoxL))
+      return false ;
+  // Verifico se il parallelepipedo interferisce con i triangoli del poliedro presenti nel suo BBox
+   INTVECTOR vT ;
+   Stm.GetAllTriaOverlapBox( b3Box, vT) ;
+   for ( int nT : vT) {
+      Triangle3d Tria ;
+      if ( Stm.GetTriangle( nT, Tria)) {
+         if ( CDeBoxTria( frBox, vtDiag, Tria, dSafeDist))
+            return true ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CDeClosedSurfTmClosedSurfTm.cpp

@@ -1,12 +1,12 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2020-2020
+//  EgalTech 2020-2024
 //----------------------------------------------------------------------------
-// File :      CDeSurfTmSurfTm.h      Data : 14.06.23    Versione : 2.5f3 
+// File : CDeClosedSurfTmClosedSurfTm.h    Data : 24.03.24  Versione : 2.6c2 
 // Contenuto : Implementazione funzione verifica collisione tra
 //             SurfTm e SurfTm.
 //
 // Modifiche : 13.11.20 LM Creazione modulo.
-//
+//             24.03.24 DS Aggiunta TestSurfTmSurfTm.
 //
 //----------------------------------------------------------------------------
 
@@ -122,3 +122,89 @@ CDeClosedSurfTmClosedSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& Surf
    DistPointSurfTm DistPoinBSrfA( ptPointB, *pSrfA) ;
    return ( DistPoinASrfB.IsPointInside()  ||  DistPoinBSrfA.IsPointInside()) ;
 }
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra le due superfici : restituisce true in caso di interferenza.
+//----------------------------------------------------------------------------
+bool
+TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double dSafeDist)
+{
+  // Recupero le superfici base
+   const SurfTriMesh* pSrfA = GetBasicSurfTriMesh( &SurfA) ;
+   const SurfTriMesh* pSrfB = GetBasicSurfTriMesh( &SurfB) ;
+  // Se le superfici non sono valide, non ha senso proseguire
+   if ( pSrfA == nullptr ||  ! pSrfA->IsValid()  ||
+        pSrfB == nullptr ||  ! pSrfB->IsValid())
+      return true ;
+  // Se i box delle superfici non si intersecano, ho finito.
+   BBox3d b3BoxA, b3BoxB ;
+   pSrfA->GetLocalBBox( b3BoxA) ;
+   pSrfB->GetLocalBBox( b3BoxB) ;
+  // Se è necessario, espando il box di B di una costante additiva pari alla distanza di sicurezza.
+   if ( dSafeDist > EPS_SMALL)
+      b3BoxB.Expand( dSafeDist) ;
+  // Se i box non si sovrappongono, non c'è collisione. Ho finito.
+   if ( ! b3BoxA.Overlaps( b3BoxB))
+      return false ;
+  // Recupero i triangoli di A che interferiscono col box di B
+   INTVECTOR vTriaIndex ;
+   pSrfA->GetAllTriaOverlapBox( b3BoxB, vTriaIndex) ;
+  // Ciclo sui triangoli della superficie A che interferiscono col box della superficie B.
+   for ( int nTA : vTriaIndex) {
+      Triangle3d trTriaA ;
+      if ( ! ( pSrfA->GetTriangle( nTA, trTriaA)  &&  trTriaA.Validate()))
+         continue ;
+      BBox3d b3BoxTriaA ;
+      trTriaA.GetLocalBBox( b3BoxTriaA) ;
+     // Se è necessario, espando il box di una costante additiva pari alla distanza di sicurezza.
+      if ( dSafeDist > EPS_SMALL)
+         b3BoxTriaA.Expand( dSafeDist) ;
+     // Recupero i triangoli di B che interferiscono col box del triangolo di A
+      INTVECTOR vNearTria ;
+      pSrfB->GetAllTriaOverlapBox( b3BoxTriaA, vNearTria) ;
+     // Settare tutti i triangoli come già processati.
+     // Al termine della chiamata i TempInt dei triangoli valgono 0.
+      pSrfB->ResetTempInts() ;
+     // Ciclo sui triangoli della superficie B che cadono nel box del triangolo corrente della Superficie A.
+      for ( int nTB : vNearTria) {
+        // Recupero il triangolo corrente della superficie B.
+        // Se triangolo non valido salto al successivo.
+         Triangle3d trTriaB ;
+         if ( ! ( pSrfB->GetTriangle( nTB, trTriaB)  &&  trTriaB.Validate()))
+            continue ;   
+        // Se necessario considero l'offset
+         if ( dSafeDist > EPS_SMALL) {
+            int nAdjTriaId[3] ;
+            pSrfB->GetTriangleAdjacencies( nTB, nAdjTriaId) ;
+           // Ciclo sui vertici del triangolo.
+            for ( int nVB = 0 ; nVB < 3 ; ++ nVB) {
+              // Se il triangolo adiacente al triangolo corrente su questo edge
+              // non è stato processato, processo il vertice e l'edge.
+               int nAdjTriaTempFlag ;
+               if ( ! ( pSrfB->GetTempInt( nAdjTriaId[nVB], nAdjTriaTempFlag)  ||  nAdjTriaTempFlag == 0))
+                  continue ;
+              // Processo il vertice: se c'è collisione fra triangolo A e sfera ho finito.
+               if ( CDeSimpleSpheTria( trTriaB.GetP( nVB), dSafeDist, trTriaA))
+                  return true ;
+              // Processo l'edge: se c'è collisione fra triangolo A e cilindro ho finito.
+               Vector3d vtEdgeV = trTriaB.GetP( ( nVB + 1) % 3) - trTriaB.GetP( nVB) ;
+               double dEdgeLen = vtEdgeV.Len() ;
+               vtEdgeV /= dEdgeLen ;
+               Frame3d frCyl ;
+               frCyl.Set( trTriaB.GetP( nVB), vtEdgeV) ;
+               if ( CDeSimpleCylTria( frCyl, dSafeDist, dEdgeLen, trTriaA))
+                  return true ;  
+            }
+           // Traslo il triangolo
+            trTriaB.Translate( dSafeDist * trTriaB.GetN()) ;
+         }
+        // Processo il triangolo: se i due triangoli collidono ho finito.
+         if ( CDeTriaTria( trTriaA, trTriaB))
+            return true ;
+        // Segno il triangolo come processato: nTemp = 1
+         pSrfB->SetTempInt( nTB, 1) ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CDeConeFrustumClosedSurfTm.cpp

@@ -1,12 +1,13 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2020-2020
+//  EgalTech 2020-2024
 //----------------------------------------------------------------------------
-// File : CDCylSurfTm.cpp      Data : 09.11.20        Versione : 
+// File : CDeConeFrustumClosedSurfTm.cpp  Data : 24.031.24  Versione : 2.6c2 
 // Contenuto : Implementazione della verifica di collisione tra
 //             Cone e Closed SurftriMesh.
 //
 //
 // Modifiche : 09.11.20 LM Creazione modulo.
+//             24.03.24 DS Aggiunta TestConeFrustumSurfTm.
 //
 //
 //----------------------------------------------------------------------------
@@ -69,4 +70,44 @@ CDeConeFrustumClosedSurfTm( const Frame3d& frCone, double dBaseRad, double dTopR
    DistPointSurfTm DistConeCenSurfCalc( ptConeCen, Stm) ;
   // Se il tronco di cono è interno c'è collisione
    return DistConeCenSurfCalc.IsPointInside() ;
-}
+}
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra il tronco di cono e la superficie : restituisce true in caso di interferenza.
+//----------------------------------------------------------------------------
+bool
+TestConeFrustumSurfTm( const Frame3d& frCone, double dBaseRad, double dTopRad, double dHeight,
+                       const ISurfTriMesh& Stm, double dSafeDist)
+{
+  // Se il tronco di cono non è ben definito non ha senso proseguire
+   if ( max( dBaseRad, dTopRad) < EPS_SMALL  ||  dHeight < EPS_SMALL)
+      return true ;
+  // Se superficie non valida, non ha senso proseguire
+   if ( ! Stm.IsValid())
+      return true ;
+  // Recupero BBox della trimesh
+   BBox3d b3Surf = Stm.GetAllTriaBox() ;
+  // Calcolo il BBox del tronco di cono
+   double dMaxRad = max( dBaseRad, dTopRad) ;
+   BBox3d b3ConeL( Point3d( -dMaxRad, -dMaxRad, 0),
+                   Point3d( dMaxRad, dMaxRad, dHeight)) ;
+   if ( dSafeDist > EPS_SMALL)
+      b3ConeL.Expand( dSafeDist) ;
+   BBox3d b3Cone = GetToGlob( b3ConeL, frCone) ;
+  // Se i BBox non interferiscono, non c'è collisione
+   if ( ! b3Surf.Overlaps( b3Cone)  ||  ! b3Surf.Overlaps( frCone, b3ConeL))
+      return false ;
+  // Recupero i triangoli che interferiscono con il box del cono
+   INTVECTOR vT ;
+   Stm.GetAllTriaOverlapBox( b3Cone, vT) ;
+  // Verifico se il tronco di cono interferisce con i triangoli del poliedro presenti nel suo BBox
+   for ( int nT : vT) {
+      Triangle3d trTria ;
+      if ( Stm.GetTriangle( nT, trTria)) {
+         if ( CDeConeFrustumTria( frCone, dBaseRad, dTopRad, dHeight, trTria, dSafeDist))
+            return true ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CDeConeFrustumTria.h

@@ -1,7 +1,7 @@
 //----------------------------------------------------------------------------
 //  EgalTech 2020-2020
 //----------------------------------------------------------------------------
-// File : CDeConTria.cpp   Data : 28.10.20     Versione : 2.2k1
+// File : CDeConeFrustumTria.h   Data : 28.10.20     Versione : 2.2k1
 // Contenuto : Dichiarazione della verifica di collisione tra 
 //             Cone e Triangle3d.
 //

CDeConvexTorusClosedSurfTm.cpp

@@ -1,12 +1,12 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2020-2020
+//  EgalTech 2020-2024
 //----------------------------------------------------------------------------
-// File : CDeConvexTorusTria.cpp      Data : 18.11.20    Versione : 
+// File : CDeConvexTorusTria.cpp      Data : 24.03.24    Versione : 2.6c2
 // Contenuto : Implementazione funzione verifica collisione tra
 //             toro convesso e SurfTriMesh.
 //
 // Modifiche : 18.11.20 LM Creazione modulo.
-//
+//             24.03.24 DS Aggiunta TestConvexTorusSurfTm.
 //
 //----------------------------------------------------------------------------
 
@@ -66,3 +66,42 @@ CDeConvexTorusClosedSurfTm( const Frame3d& frTorus, double dRad1, double dRad2,
   // Se il toro è interno c'è collisione
    return ( DistConeOrigSurfCalc.IsPointInside()) ;
 }
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra il toro convesso e la superficie : restituisce true in caso di interferenza
+//----------------------------------------------------------------------------
+bool
+TestConvexTorusSurfTm( const Frame3d& frTorus, double dRad1, double dRad2,
+                       const ISurfTriMesh& Stm, double dSafeDist)
+{
+  // I raggi devono essere non nulli
+   if ( dRad1 < EPS_SMALL  ||  dRad2 < EPS_SMALL)
+      return true ;
+  // Se superficie non valida, non ha senso proseguire
+   if ( ! Stm.IsValid())
+      return true ;
+  // Box della superficie
+   BBox3d b3Surf = Stm.GetAllTriaBox() ;
+  // Box del toro (sempre completo)
+   BBox3d b3TorusL( Point3d( -dRad1 - dRad2, -dRad1 - dRad2, -dRad2),
+                    Point3d(  dRad1 + dRad2,  dRad1 + dRad2,  dRad2)) ;
+   if ( dSafeDist > EPS_SMALL)
+      b3TorusL.Expand( dSafeDist) ;
+   BBox3d b3Torus = GetToGlob( b3TorusL, frTorus) ;
+  // Se i BBox non interferiscono, non c'è collisione
+   if ( ! b3Surf.Overlaps( b3Torus)  ||  ! b3Surf.Overlaps( frTorus, b3TorusL))
+      return false ;
+  // Recupero i triangoli che interferiscono con il box del toro
+   INTVECTOR vT ;
+   Stm.GetAllTriaOverlapBox( b3Torus, vT) ;
+  // Verifico se il toro interferisce con i triangoli del poliedro presenti nel suo BBox
+   for ( int nT : vT) {
+      Triangle3d trTria ;
+      if ( Stm.GetTriangle( nT, trTria)) {
+         if ( CDeConvexTorusTria( frTorus, dRad1, dRad2, CT_TOT, trTria, dSafeDist))
+            return true ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CDeCylClosedSurfTm.cpp

@@ -1,13 +1,13 @@
 //----------------------------------------------------------------------------
 //  EgalTech 2020-2024
 //----------------------------------------------------------------------------
-// File : CDCylSurfTm.cpp      Data : 15.02.24        Versione : 2.6b2
+// File : CDeCylClosedSurfTm.cpp      Data : 24.03.24        Versione : 2.6c2
 // Contenuto : Implementazione della verifica di collisione tra
 //             Cylinder e Closed SurftriMesh.
 //
 //
 // Modifiche : 09.01.20 DS Creazione modulo.
-//
+//             24.03.24 DS Aggiunta TestCylSurfTm.
 //
 //----------------------------------------------------------------------------
 
@@ -73,3 +73,46 @@ CDeCylClosedSurfTm( const Frame3d& frCyl, double dR, double dH, const ISurfTriMe
   // Se il cilindro è interno c'è collisione
    return ( DistCylCenSurfCalc.IsPointInside()) ;
 }
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra il cilindro e la superficie : restituisce true in caso di interferenza
+//----------------------------------------------------------------------------
+bool
+TestCylSurfTm( const Frame3d& frCyl, double dR, double dH, const ISurfTriMesh& Stm, double dSafeDist)
+{
+  // Il cilindro deve essere ben definito
+   if ( dR < EPS_SMALL  ||  dH < EPS_SMALL)
+      return true ;
+  // Se superficie non valida, non ha senso proseguire
+   if ( ! Stm.IsValid())
+      return true ;
+  // Recupero BBox della superficie poligonale
+   BBox3d b3Poly = Stm.GetAllTriaBox() ;
+  // Sistemazioni cilindro
+   Frame3d frMyCyl = frCyl ;
+   if ( dH < 0) {
+      frMyCyl.Translate( dH * frMyCyl.VersZ()) ;
+      dH = -dH ;
+   }
+  // Calcolo il BBox del cilindro
+   BBox3d b3CylL( Point3d( -dR, -dR, 0),
+                  Point3d(  dR,  dR, dH)) ;
+   if ( dSafeDist > EPS_SMALL)
+      b3CylL.Expand( dSafeDist) ;
+   BBox3d b3Cyl = GetToGlob( b3CylL, frMyCyl) ;
+  // Se i BBox non interferiscono, non c'è interferenza
+   if ( ! b3Poly.Overlaps( b3Cyl)  ||  ! b3Poly.Overlaps( frMyCyl, b3CylL))
+      return false ;
+  // Verifico se il cilindro interferisce con i triangoli del poliedro presenti nel suo BBox
+   INTVECTOR vT ;
+   Stm.GetAllTriaOverlapBox( b3Cyl, vT) ;
+   for ( int nT : vT) {
+      Triangle3d Tria ;
+      if ( Stm.GetTriangle( nT, Tria)) {
+         if ( CDeCylTria( frMyCyl, dR, dH, Tria, dSafeDist))
+            return true ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CDeRectPrismoidClosedSurfTm.cpp

@@ -1,12 +1,12 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2020-2020
+//  EgalTech 2020-2024
 //----------------------------------------------------------------------------
-// File : CDePyramidClosedSurfTm.h      Data : 09.11.20    Versione : 
+// File : CDeRectPrismoidClosedSurfTm.h    Data : 24.03.24  Versione : 2.6c2
 // Contenuto : Implementazione funzione verifica collisione tra
 //             Prismoide a basi rettangolari e Closed SurfTriMesh.
 //
 // Modifiche : 09.11.20 LM Creazione modulo.
-//
+//             24.03.24 DS Aggiunta TestRectPrismoidSurfTm.
 //
 //----------------------------------------------------------------------------
 
@@ -74,3 +74,48 @@ CDeRectPrismoidClosedSurfTm( const Frame3d& frPrismoid, double dLenghtBaseX, dou
   // C'è collisione se il tronco di piramide è interno
    return ( DistPyrCenSurfCalc.IsPointInside()) ;
 }
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra il Prismoide a basi rettangolari e la superficie : restituisce true in caso di interferenza.
+//----------------------------------------------------------------------------
+bool
+TestRectPrismoidSurfTm( const Frame3d& frPrismoid, double dLenghtBaseX, double dLenghtBaseY,
+                        double dLenghtTopX, double dLenghtTopY, double dHeight,
+                        const ISurfTriMesh& Stm, double dSafeDist)
+{
+  // Se il tronco di piramide non è definito non ha senso proseguire.
+   if ( max( dLenghtBaseX, dLenghtTopX) < EPS_SMALL  ||
+        max( dLenghtBaseY, dLenghtTopY) < EPS_SMALL  ||
+        dHeight < EPS_SMALL)
+      return true ;
+  // Se superficie non valida, non ha senso proseguire
+   if ( ! Stm.IsValid())
+      return true ;
+  // Recupero BBox della trimesh
+   BBox3d b3Surf = Stm.GetAllTriaBox() ;
+  // Calcolo il BBox del tronco di piramide
+   double dMaxLenX = max( dLenghtBaseX, dLenghtTopX) ;
+   double dMaxLenY = max( dLenghtBaseY, dLenghtTopY) ;
+   BBox3d b3PyrL( Point3d( -dMaxLenX / 2, -dMaxLenY / 2, 0.),
+                  Point3d(  dMaxLenX / 2,  dMaxLenY / 2, dHeight)) ;
+   if ( dSafeDist > EPS_SMALL)
+      b3PyrL.Expand( dSafeDist) ;
+   BBox3d b3Pyr = GetToGlob( b3PyrL, frPrismoid) ;
+  // Se i BBox non interferiscono, non c'è collisione
+   if ( ! b3Surf.Overlaps( b3Pyr)  ||  ! b3Surf.Overlaps( frPrismoid, b3PyrL))
+      return false ;
+  // Recupero i triangoli che interferiscono con il box del tronco di piramide.
+   INTVECTOR vT ;
+   Stm.GetAllTriaOverlapBox( b3Pyr, vT) ;
+  // Verifico se il tronco di piramide interferisce con i triangoli del poliedro presenti nel suo BBox
+   for ( int nT : vT) {
+      Triangle3d trTria ; 
+      if ( Stm.GetTriangle( nT, trTria)) {
+         if ( CDeRectPrismoidTria( frPrismoid, dLenghtBaseX, dLenghtBaseY, dLenghtTopX, dLenghtTopY, dHeight,
+                                    trTria, dSafeDist))
+            return true ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CDeSpheClosedSurfTm.cpp

@@ -1,13 +1,13 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2020-2020
+//  EgalTech 2020-2024
 //----------------------------------------------------------------------------
-// File : CDCylSurfTm.cpp      Data : 09.01.20        Versione : 2.2a2
+// File : CDCylSurfTm.cpp      Data : 24.03.24        Versione : 2.6c2
 // Contenuto : Implementazione della verifica di collisione tra
 //             Sphere e Closed SurftriMesh.
 //
 //
 // Modifiche : 09.01.20 DS Creazione modulo.
-//
+//             24.03.24 DS Aggiunta TestSpheSurfTm.
 //
 //----------------------------------------------------------------------------
 
@@ -60,3 +60,38 @@ CDeSpheClosedSurfTm( const Point3d& ptCen, double dR, const ISurfTriMesh& Stm, d
   // C'è collisione se la sfera è interna.
    return ( DistCenSurfCalc.IsPointInside()) ;
 }
+
+//----------------------------------------------------------------------------
+// Verifica l'interferenza tra la sfera e la superficie : restituisce true in caso di interferenza
+//----------------------------------------------------------------------------
+bool
+TestSpheSurfTm( const Point3d& ptCen, double dR, const ISurfTriMesh& Stm, double dSafeDist)
+{
+  // Il raggio deve essere non nullo
+   if ( dR < EPS_SMALL)
+      return true ;
+  // Se superficie non valida, non ha senso proseguire
+   if ( ! Stm.IsValid())
+      return true ;
+  // Recupero BBox del poliedro
+   BBox3d b3Poly = Stm.GetAllTriaBox() ;
+  // Calcolo il BBox della sfera
+   BBox3d b3Sphe( ptCen, dR) ;
+   if ( dSafeDist > EPS_SMALL)
+      b3Sphe.Expand( dSafeDist) ;
+  // Se i BBox non interferiscono, non c'è collisione
+   if ( ! b3Sphe.Overlaps( b3Poly))
+      return false ;
+  // Verifico se la sfera interferisce con i triangoli del poliedro presenti nel suo BBox
+   INTVECTOR vT ;
+   Stm.GetAllTriaOverlapBox( b3Sphe, vT) ;
+   for ( int nT : vT) {
+      Triangle3d Tria ;
+      if ( Stm.GetTriangle( nT, Tria)) {
+         if ( CDeSpheTria( ptCen, dR, Tria, dSafeDist))
+            return true ;
+      }
+   }
+  // Non c'è interferenza
+   return false ;
+}

CurveAux.cpp

@@ -524,6 +524,39 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
 
   // se periodica
    if ( cnData.bPeriodic  ||  ! cnData.bClamped) {
+      // se la curva è peridica verifco che effettivamente ci sia un numero di punti ripetituti uguale al grado della curva
+      if ( cnData.bPeriodic ) {
+         bool bRepetead = true ;
+         for ( int i = 0 ; i < cnData.nDeg ; ++i) {
+            if ( ! AreSamePointApprox( cnData.vCP[i], cnData.vCP.end()[-cnData.nDeg + i]) ) {
+               bRepetead = false ;
+               break ;
+            }
+         }
+         if ( ! bRepetead){
+            // se il primo e l'ultimo punto non coincidono allora aggiungo il primo punto in fondo al vettore dei punti di controllo
+            if ( ! AreSamePointApprox( cnData.vCP[0], cnData.vCP.back()) ) {
+               cnData.vCP.push_back( cnData.vCP[0]) ;
+               if ( cnData.bRat)
+                  cnData.vW.push_back( cnData.vW[0]) ;
+            }
+            // se effettivamente ho dei nodi in più da togliere allora li tolgo
+            if ( int( cnData.vU.size()) !=  int(cnData.vCP.size()) + cnData.nDeg - 1) {
+               // devo poi anche togliere i nodi di troppo // presuppongo che la convenzione sia che i nodi di troppo sono alla fine del vettore dei nodi
+               cnData.vU = DBLVECTOR( cnData.vU.begin(), cnData.vU.end() - cnData.nDeg) ;
+               // controlloeventualmente anche i nodi extra
+               // se ne ho due in più ne tolgo uno in cima e uno in fondo
+               if ( cnData.vU.size() == int( cnData.vCP.size()) + cnData.nDeg + 1 )   // significa che ci sono due nodi extra, uno all'inizio e uno alla fine, da togliere
+                  cnData.vU = vector<double>( cnData.vU.begin() + 1, cnData.vU.end() - 1) ;
+               // se ne ho solo uno in più lo tolgo in cima
+               else if ( cnData.vU.size() == int( cnData.vCP.size()) + cnData.nDeg)
+                  cnData.vU = vector<double>( cnData.vU.begin() + 1, cnData.vU.end()) ;
+            }
+            cnData.bPeriodic = false ;
+            return true ;
+         }
+      }
+
      // va trasformata in non-periodica (clamped)
      // bisogna aumentare la molteplicità dei nodi u_p-1 e u_(m-p+1) fino ad arrivare al grado della nurbs
      // e poi scartare nodi e punti fuori dalla regione clamped ( al di fuori della regione u_p-1 -> u_(m-p+1))

EgtGeomKernel.vcxproj

@@ -319,6 +319,9 @@ copy $(TargetPath) \EgtProg\Dll64</Command>
     <ClCompile Include="IntersPlaneVolZmap.cpp" />
     <ClCompile Include="IntersLineSurfBez.cpp" />
     <ClCompile Include="PolygonElevation.cpp" />
+    <ClCompile Include="Quaternion.cpp" />
+    <ClCompile Include="RotationMinimizingFrame.cpp" />
+    <ClCompile Include="RotationXplaneFrame.cpp" />
     <ClCompile Include="SbzStandard.cpp" />
     <ClCompile Include="Voronoi.cpp" />
     <ClInclude Include="..\Include\EGkCDeClosedSurfTmClosedSurfTm.h" />
@@ -331,6 +334,9 @@ copy $(TargetPath) \EgtProg\Dll64</Command>
     <ClInclude Include="..\Include\EGkIntersPlaneBox.h" />
     <ClInclude Include="..\Include\EGkIntersPlaneVolZmap.h" />
     <ClInclude Include="..\Include\EGkPolygonElevation.h" />
+    <ClInclude Include="..\Include\EGkQuaternion.h" />
+    <ClInclude Include="..\Include\EGkRotationMinimizingFrame.h" />
+    <ClInclude Include="..\Include\EGkRotationXplaneFrame.h" />
     <ClInclude Include="..\Include\EGkSubtractProjectedFacesOnStmFace.h" />
     <ClInclude Include="CDeBoxTria.h" />
     <ClInclude Include="CDeCapsTria.h" />
@@ -462,7 +468,6 @@ copy $(TargetPath) \EgtProg\Dll64</Command>
     <ClCompile Include="VolZmapGraphics.cpp" />
     <ClCompile Include="VolZmapVolume.cpp" />
     <ClCompile Include="VolZmap.cpp" />
-    <ClInclude Include="IntersLineSurfBez.h" />
     <ClInclude Include="Voronoi.h" />
   </ItemGroup>
   <ItemGroup>

EgtGeomKernel.vcxproj.filters

@@ -528,6 +528,15 @@
     <ClCompile Include="SbzStandard.cpp">
       <Filter>File di origine\GeoCreate</Filter>
     </ClCompile>
+    <ClCompile Include="RotationMinimizingFrame.cpp">
+      <Filter>File di origine\Base</Filter>
+    </ClCompile>
+    <ClCompile Include="RotationXplaneFrame.cpp">
+      <Filter>File di origine\Base</Filter>
+    </ClCompile>
+    <ClCompile Include="Quaternion.cpp">
+      <Filter>File di origine\Base</Filter>
+    </ClCompile>
   </ItemGroup>
   <ItemGroup>
     <ClInclude Include="stdafx.h">
@@ -1196,6 +1205,15 @@
     <ClInclude Include="..\Include\EGkIntersLineVolZmap.h">
       <Filter>File di intestazione\Include</Filter>
     </ClInclude>
+    <ClInclude Include="..\Include\EGkRotationMinimizingFrame.h">
+      <Filter>File di intestazione\Include</Filter>
+    </ClInclude>
+    <ClInclude Include="..\Include\EGkRotationXplaneFrame.h">
+      <Filter>File di intestazione\Include</Filter>
+    </ClInclude>
+    <ClInclude Include="..\Include\EGkQuaternion.h">
+      <Filter>File di intestazione\Include</Filter>
+    </ClInclude>
   </ItemGroup>
   <ItemGroup>
     <ResourceCompile Include="EgtGeomKernel.rc">

ExtDimension.cpp

@@ -28,6 +28,7 @@
 #include "/EgtDev/Include/EgtNumUtils.h"
 #include "/EgtDev/Include/EGkCurve.h"
 #include <new>
+#include <algorithm>
 
 using namespace std ;
 
@@ -319,8 +320,23 @@ ExtDimension::SetAngularEx( const Point3d& ptV1, const Point3d& ptP1,
       dLenDir = m_vtDir.Len() ;
    }
   // assegno gli altri punti notevoli della quotatura
-   m_ptP1 = ptV1 - ( ptV1 - m_ptP6) * m_vtN * m_vtN ;
-   m_ptP2 = ptV2 - ( ptV2 - m_ptP6) * m_vtN * m_vtN ;
+  // se i ptV concidono con l'intersezione tra le linee le sposto leggermente verso l'altro punto di quel lato
+   Point3d ptV1New = ptV1 ;
+   Point3d ptV2New = ptV2 ;
+   if ( AreSamePointApprox( ptV1, ptV)) {
+      Vector3d vtDir1 ; Line1.GetStartDir( vtDir1) ;
+      double dLen ; Line1.GetLength( dLen) ;
+      dLen = max( min( 1., dLen/ 10), 0.01) ;
+      ptV1New = ptV1 + vtDir1 * dLen ;
+   }
+   if ( AreSamePointApprox( ptV2, ptV)) {
+      Vector3d vtDir2 ; Line2.GetStartDir( vtDir2) ;
+      double dLen ; Line2.GetLength( dLen) ;
+      dLen = max( min( 1., dLen/ 10), 0.01) ;
+      ptV2New = ptV2 + vtDir2 * dLen ;
+   }
+   m_ptP1 = ptV1New - ( ptV1New - m_ptP6) * m_vtN * m_vtN ;
+   m_ptP2 = ptV2New - ( ptV2New - m_ptP6) * m_vtN * m_vtN ;
    m_ptP5 = m_ptP6 + vtLine1 * dLenDir ;
    double dLen1 = Dist( m_ptP6, m_ptP1) ;
    m_ptP3 = m_ptP6 + vtLine1 * ( dLenDir + m_dExtLineLen * ( dLen1 < dLenDir ? 1 : -1)) ;

ExtDimension.h

@@ -108,6 +108,7 @@ class ExtDimension : public IExtDimension, public IGeoObjRW
                         { return m_dTextHeight ; }
       bool GetMidPoint( Point3d& ptMid) const override ;
       bool GetCenterPoint( Point3d& ptCen) const override ;
+      const std::string& GetSubType( void) const override ;
       bool ApproxWithLines( double dLinTol, double dAngTolDeg, POLYLINELIST& lstPL) const override ;
 
    public :                      // IGeoObjRW
@@ -130,7 +131,6 @@ class ExtDimension : public IExtDimension, public IGeoObjRW
 
    private :
       bool CopyFrom( const ExtDimension& gpSrc) ;
-      const std::string& GetSubType( void) const ;
       bool Update( void) const ;
       double GetTextHalfDist( const Point3d& ptText, bool bUseRot = true) const ;
       bool GetArrowHead( const Point3d& ptTip, const Vector3d& vtDir, PolyLine& PL) const ;

FilletChamfer.cpp

@@ -18,6 +18,7 @@
 #include "/EgtDev/Include/EGkFilletChamfer.h"
 #include "/EgtDev/Include/EGkDistPointCurve.h"
 #include "/EgtDev/Include/EGkIntersCurves.h"
+#include "/EgtDev/Include/EgkOffsetCurve.h"
 #include "/EgtDev/Include/EgtPointerOwner.h"
 
 using namespace std ;
@@ -40,13 +41,21 @@ CalcForFillet( const ICurve& cCrv1, const Point3d& ptNear1,
    if ( ! dPC1.GetSideAtMinDistPoint( 0, vtNorm, nSide1))
       return false ;
    double dOffs1 = ( nSide1 == MDS_RIGHT ? dRadius : - dRadius) ;
-  // calcolo l'offset nel piano locale e dal lato opportuno di una copia della curva 1
+  // calcolo gli offset nel piano locale e dal lato opportuno di una copia della curva 1
    PtrOwner<ICurve> pCopy1( cCrv1.Clone()) ;
    if ( IsNull( pCopy1))
       return false ;
    pCopy1->ToLoc( frIntr) ;
    pCopy1->SetExtrusion( Z_AX) ;
-   if ( ! pCopy1->SimpleOffset( dOffs1, ICurve::OFF_FILLET))
+   OffsetCurve OffsCrv1 ;
+   OffsCrv1.Make( pCopy1, dOffs1, ICurve::OFF_FILLET) ;
+   ICURVEPOVECTOR vOffs1 ;  
+   ICurve* pCrv = OffsCrv1.GetLongerCurve() ;
+   while ( pCrv != nullptr) {
+      vOffs1.emplace_back( pCrv) ;
+      pCrv = OffsCrv1.GetLongerCurve() ;
+   }
+   if ( vOffs1.empty())
       return false ;
 
   // determino il lato di offset della curva 2
@@ -54,23 +63,50 @@ CalcForFillet( const ICurve& cCrv1, const Point3d& ptNear1,
    if ( ! dPC2.GetSideAtMinDistPoint( 0, vtNorm, nSide2))
       return false ;
    double dOffs2 = ( nSide2 == MDS_RIGHT ? dRadius : - dRadius) ;
-  // calcolo l'offset nel piano locale e dal lato opportuno di una copia della curva 2
+  // calcolo gli offset nel piano locale e dal lato opportuno di una copia della curva 2
    PtrOwner<ICurve> pCopy2( cCrv2.Clone()) ;
    if ( IsNull( pCopy2))
       return false ;
    pCopy2->ToLoc( frIntr) ;
    pCopy2->SetExtrusion( Z_AX) ;
-   if ( ! pCopy2->SimpleOffset( dOffs2, ICurve::OFF_FILLET))
+   OffsetCurve OffsCrv2 ;
+   OffsCrv2.Make( pCopy2, dOffs2, ICurve::OFF_FILLET) ;
+   ICURVEPOVECTOR vOffs2 ;
+   pCrv = OffsCrv2.GetLongerCurve() ;
+   while ( pCrv != nullptr) {
+      vOffs2.emplace_back( pCrv) ;
+      pCrv = OffsCrv2.GetLongerCurve() ;
+   }
+   if ( vOffs2.empty())
       return false ;
 
-  // calcolo l'intersezione tra le due curve
-   Point3d ptInt1, ptInt2 ;
+  // calcolo le intersezioni tra tutte le curve di offset e seleziono quella più vicina ai punti passati
+   Point3d ptInt1 = P_INVALID, ptInt2 = P_INVALID ;
    Point3d ptNearI = Media( ptNear1, ptNear2) ;
    ptNearI.ToLoc( frIntr) ;
-   IntersCurveCurve intCC( *pCopy1, *pCopy2) ;
-   if ( ! intCC.GetIntersPointNearTo( 0, ptNearI, ptInt1)  ||
-        ! intCC.GetIntersPointNearTo( 1, ptNearI, ptInt2))
-      return false ;
+   double dMinDist = INFINITO ;
+   for ( int i = 0 ; i < int( vOffs1.size()) ; i++) {
+      for ( int j = 0 ; j < int( vOffs2.size()) ; j ++) {
+         IntersCurveCurve intCC( *vOffs1[i], *vOffs2[j]) ;
+         if ( intCC.GetIntersCount() == 0)
+            continue ;
+         Point3d ptInt1Curr, ptInt2Curr ;
+         if ( ! intCC.GetIntersPointNearTo( 0, ptNearI, ptInt1Curr)  ||
+              ! intCC.GetIntersPointNearTo( 1, ptNearI, ptInt2Curr))
+            return false ;
+         Point3d ptCenCurr = Media( ptInt1Curr, ptInt2Curr) ;
+         double dDist = Dist( ptNearI, ptCenCurr) ;
+         if ( dDist < dMinDist - EPS_SMALL) {
+            dMinDist = dDist ;
+            ptInt1 = ptInt1Curr ;
+            ptInt2 = ptInt2Curr ;
+         }
+      }  
+   }
+  // se non sono state trovate intersezioni esco
+   if ( ! ptInt1.IsValid()  ||  ! ptInt2.IsValid())
+      return false ;   
+   
    ptInt1.ToGlob( frIntr) ;
    ptInt2.ToGlob( frIntr) ;
    ptCen = Media( ptInt1, ptInt2) ;
@@ -105,7 +141,7 @@ CreateFillet( const ICurve& cCrv1, const Point3d& ptNear1,
               const ICurve& cCrv2, const Point3d& ptNear2,
               const Vector3d& vtNorm, double dRadius, double& dPar1, double& dPar2)
 {
-  // verifico validità parametri ricevuti
+  // verifico validità parametri ricevuti
    if ( &cCrv1 == nullptr  ||  &ptNear1 == nullptr  ||
         &cCrv2 == nullptr  ||  &ptNear2 == nullptr  ||
         &vtNorm == nullptr  ||  &dPar1 == nullptr  ||  &dPar2 == nullptr)
@@ -119,7 +155,7 @@ CreateFillet( const ICurve& cCrv1, const Point3d& ptNear1,
                          ptCen, ptTg1, ptTg2, nSide1, nSide2, dSinA, dTgPar1, dTgPar2))
       return nullptr ;
 
-  // se tangenti parallele al contatto con fillet, ricalcolo con raggio più piccolo
+  // se tangenti parallele al contatto con fillet, ricalcolo con raggio più piccolo
    bool bParallel = ( abs( dSinA) < EPS_SMALL) ;
    if ( bParallel) {
       Point3d ptQQQ, ptQQ1, ptQQ2 ;
@@ -165,7 +201,7 @@ CreateChamfer( const ICurve& cCrv1, const Point3d& ptNear1,
                const ICurve& cCrv2, const Point3d& ptNear2,
                const Vector3d& vtNorm, double dDist, double& dPar1, double& dPar2)
 {
-  // verifico validità parametri ricevuti
+  // verifico validità parametri ricevuti
    if ( &cCrv1 == nullptr  ||  &ptNear1 == nullptr  ||
         &cCrv2 == nullptr  ||  &ptNear2 == nullptr  ||
         &vtNorm == nullptr  ||  &dPar1 == nullptr  ||  &dPar2 == nullptr)

IntersLineSurfBez.cpp

@@ -22,13 +22,14 @@
 
 using namespace std ;
 
-//----------------------
-bool
+//----------------------------------------------------------------------------
+static bool
 RefineIntersNewton( const Point3d& ptL, const Vector3d& vtL, double dLen, bool bFinite,
-   const ISurfBezier* pSurfBz, Point3d& ptSP, Point3d& ptIBz) {
-   // la funzione raffina la posisione del punto ptSP, minimizzando la distanza dalla retta e restituisce il punto di intersezione ptIBz
+                    const ISurfBezier* pSurfBz, Point3d& ptSP, Point3d& ptIBz)
+{
+  // la funzione raffina la posisione del punto ptSP, minimizzando la distanza dalla retta e restituisce il punto di intersezione ptIBz
    pSurfBz->GetPointD1D2( ptSP.x / SBZ_TREG_COEFF, ptSP.y / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptIBz) ;
-   // usando un algoritmo di newton cerco di avvicinarmi il più possibile alla retta
+  // usando un algoritmo di newton cerco di avvicinarmi il più possibile alla retta
    DistPointLine dpl( ptIBz, ptL, vtL, dLen, bFinite) ;
    double dDistNew = 0, dDistPre = 0 ;
    dpl.GetDist(dDistNew) ;
@@ -36,8 +37,8 @@ RefineIntersNewton( const Point3d& ptL, const Vector3d& vtL, double dLen, bool b
    int nCount = 0 ;
    double dh = EPS_SMALL ;
    pSurfBz->GetPointD1D2( ptSP.x, ptSP.y, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptIBz) ;
-   // 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.
+  // 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 ;
@@ -50,48 +51,21 @@ RefineIntersNewton( const Point3d& ptL, const Vector3d& vtL, double dLen, bool b
       DistPointLine dplNewV( ptIBzNew2, ptL, vtL, dLen, bFinite) ;
       dplNewV.GetDist( dDistNew) ;
       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
+     // mi avvicino cercando di annullare la distanza in un colpo solo
       double dr = - dDistPre / ( dfdU + dfdV) ;
       pSurfBz->GetPointD1D2(( ptSP.x + dr * dfdU) / SBZ_TREG_COEFF, ( ptSP.y + dr * dfdV) / SBZ_TREG_COEFF, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, ptIBz) ;
       DistPointLine dplNew( ptIBz, ptL, vtL, dLen, bFinite) ;
       dplNew.GetDist( dDistNew) ;
-      ++nCount ;
+      ++ nCount ;
    }
 
-   return nCount != 99 ;
+   return ( nCount != 99) ;
 }
 
-////----------------------------------------------------------------------------
-//bool
-//RefineIntersBisec( const Point3d& ptL, const Vector3d& vtL, double dLen, bool bFinite,
-//                   const ISurfBezier* pSurfBz, Point3d& ptSP, Point3d& ptIBz) {
-//
-//}
-
 //----------------------------------------------------------------------------
-void
+static void
 UpdateInfoIntersLineSurfBz( const Point3d& ptL, const Vector3d& vtDir, int nILT, int nT, const Point3d& ptSP, const Point3d& ptIBz, double dCos,
-   const Point3d& ptSP2, const Point3d& ptIBz2, double dCos2, ILSBIVECTOR& vInfo)
+                            const Point3d& ptSP2, const Point3d& ptIBz2, double dCos2, ILSBIVECTOR& vInfo)
 {
    if ( nILT == ILTT_IN  ||  nILT == ILTT_EDGE  ||  nILT == ILTT_VERT) {
       double dU = ( ptIBz - ptL) * vtDir ;
@@ -105,7 +79,7 @@ UpdateInfoIntersLineSurfBz( const Point3d& ptL, const Vector3d& vtDir, int nILT,
 }
 
 //----------------------------------------------------------------------------
-void
+static void
 OrderInfoIntersLineSurfBz( ILSBIVECTOR& vInfo)
 {
    // se non trovati, esco 
@@ -120,7 +94,7 @@ OrderInfoIntersLineSurfBz( ILSBIVECTOR& vInfo)
 }
 
 //----------------------------------------------------------------------------
-//  Intersezione di una linea con una superficie TriMesh
+//  Intersezione di una linea con una superficie di Bezier
 //----------------------------------------------------------------------------
 bool
 IntersLineSurfBz( const Point3d& ptL, const Vector3d& vtL, double dLen, const ISurfBezier* pSurfBz,
@@ -257,4 +231,4 @@ FilterLineSurfBzInters( const ILSBIVECTOR& vInfo, INTDBLVECTOR& vInters)
       ++ j ;
    }
    return true ;
-}
+}

IntersLineSurfBez.h

@@ -1 +0,0 @@
-#pragma once

ProjectCurveSurfTm.cpp

@@ -15,6 +15,8 @@
 #include "stdafx.h"
 #include "DistPointLine.h"
 #include "GeoConst.h"
+#include "/EgtDev/Include/EGkDistPointCurve.h"
+#include "/EgtDev/Include/EGkDistPointSurfTm.h"
 #include "/EgtDev/Include/EGkIntersLineSurfTm.h"
 #include "/EgtDev/Include/EGkProjectCurveSurfTm.h"
 
@@ -33,63 +35,9 @@ PointsInTolerance( const PNT5AXVECTOR& vPt5ax, int nPrec, int nCurr, int nNext,
 }
 
 //----------------------------------------------------------------------------
-bool
-ProjectCurveOnSurfTm( const ICurve& crCrv, const ISurfTriMesh& tmSurf, const Vector3d& vtDir, double dLinTol, double dMaxSegmLen,
-                      PNT5AXVECTOR& vPt5ax)
+static bool
+RemovePointsInExcess( PNT5AXVECTOR& vMyPt5ax, double dLinTol, double dMaxSegmLen)
 {
-  // controllo le tolleranze
-   dLinTol = max( dLinTol, LIN_TOL_MIN) ;
-   dMaxSegmLen = max( dMaxSegmLen, 10 * EPS_SMALL) ;
-  // approssimo la curva con una polilinea entro la metà della tolleranza
-   PolyLine PL ;
-   if ( ! crCrv.ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_STD, PL))
-      return false ;
-   const double MAX_SEG_LEN = min( dMaxSegmLen, 1.) ;
-   if ( ! PL.AdjustForMaxSegmentLen( MAX_SEG_LEN))
-      return false ;
-
-  // Oggetto per calcolo massivo intersezioni tra linee di proiezione e superficie
-   Frame3d frRefLine ;
-   if ( ! frRefLine.Set( ORIG, vtDir))
-      return false ;
-   IntersParLinesSurfTm intPLSTM( frRefLine, tmSurf) ;
-
-  // Vettore locale dei punti risultanti
-   PNT5AXVECTOR vMyPt5ax ;
-   vMyPt5ax.reserve( PL.GetPointNbr()) ;
-
-  // proietto i punti della polilinea sulla superficie
-   double dU ;
-   Point3d ptP ;
-   bool bFound = PL.GetFirstUPoint( &dU, &ptP) ;
-   while ( bFound) {
-      Point3d ptL = GetToLoc( ptP, frRefLine) ;
-      ILSIVECTOR vIntRes ;
-      intPLSTM.GetInters( ptL, 1, vIntRes, false) ;
-      if ( vIntRes.size() > 0) {
-        // calcolo il punto
-         int nI = int( vIntRes.size()) - 1 ;
-         Point3d ptInt ;
-         if ( vIntRes[nI].nILTT == ILTT_SEGM  ||  vIntRes[nI].nILTT == ILTT_SEGM_ON_EDGE)
-            ptInt = vIntRes[nI].ptI2 ;
-         else
-            ptInt = vIntRes[nI].ptI ;
-        // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
-         Triangle3dEx trTria ;
-         if ( ! tmSurf.GetTriangle( vIntRes[nI].nT, trTria))
-            return false ;
-         Vector3d vtN ;
-         double dU, dV, dW ;
-         if ( BarycentricCoord( ptInt, trTria, dU, dV, dW)) 
-            vtN = dU * trTria.GetVertexNorm( 0) + dV * trTria.GetVertexNorm( 1) + dW * trTria.GetVertexNorm( 2) ;
-         if ( ! vtN.Normalize())
-            vtN = trTria.GetN() ;
-        // aggiungo al vettore dei proiettati
-         vMyPt5ax.emplace_back( ptInt, vtN, dU, 1) ;
-      }
-      bFound = PL.GetNextUPoint( &dU, &ptP) ;
-   }
-
   // rimuovo i punti allineati entro la tolleranza e non più lontani tra loro del massimo
    double dSqMaxLen = dMaxSegmLen * dMaxSegmLen ;
    double dSqTol = dLinTol * dLinTol ;
@@ -130,6 +78,309 @@ ProjectCurveOnSurfTm( const ICurve& crCrv, const ISurfTriMesh& tmSurf, const Vec
          ++ nNext ;
       }
    }
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ProjectCurveOnSurfTm( const ICurve& crCrv, const ISurfTriMesh& tmSurf, const Vector3d& vtDir, double dLinTol, double dMaxSegmLen,
+                      PNT5AXVECTOR& vPt5ax)
+{
+  // controllo le tolleranze
+   dLinTol = max( dLinTol, LIN_TOL_MIN) ;
+   dMaxSegmLen = max( dMaxSegmLen, 10 * EPS_SMALL) ;
+  // approssimo la curva con una polilinea entro la metà della tolleranza
+   PolyLine PL ;
+   if ( ! crCrv.ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_STD, PL))
+      return false ;
+   const double MAX_SEG_LEN = min( dMaxSegmLen, 1.) ;
+   if ( ! PL.AdjustForMaxSegmentLen( MAX_SEG_LEN))
+      return false ;
+
+  // Oggetto per calcolo massivo intersezioni tra linee di proiezione e superficie
+   Frame3d frRefLine ;
+   if ( ! frRefLine.Set( ORIG, vtDir))
+      return false ;
+   IntersParLinesSurfTm intPLSTM( frRefLine, tmSurf) ;
+
+  // Vettore locale dei punti risultanti
+   PNT5AXVECTOR vMyPt5ax ;
+   vMyPt5ax.reserve( PL.GetPointNbr()) ;
+
+  // proietto i punti della polilinea sulla superficie
+   double dPar ;
+   Point3d ptP ;
+   bool bFound = PL.GetFirstUPoint( &dPar, &ptP) ;
+   while ( bFound) {
+      Point3d ptL = GetToLoc( ptP, frRefLine) ;
+      ILSIVECTOR vIntRes ;
+      intPLSTM.GetInters( ptL, 1, vIntRes, false) ;
+      if ( ! vIntRes.empty()) {
+        // calcolo il punto
+         int nI = int( vIntRes.size()) - 1 ;
+         Point3d ptInt ;
+         if ( vIntRes[nI].nILTT == ILTT_SEGM  ||  vIntRes[nI].nILTT == ILTT_SEGM_ON_EDGE)
+            ptInt = vIntRes[nI].ptI2 ;
+         else
+            ptInt = vIntRes[nI].ptI ;
+        // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
+         Triangle3dEx trTria ;
+         if ( ! tmSurf.GetTriangle( vIntRes[nI].nT, trTria))
+            return false ;
+         Vector3d vtN ;
+         double dU, dV, dW ;
+         if ( BarycentricCoord( ptInt, trTria, dU, dV, dW)) 
+            vtN = dU * trTria.GetVertexNorm( 0) + dV * trTria.GetVertexNorm( 1) + dW * trTria.GetVertexNorm( 2) ;
+         if ( ! vtN.Normalize())
+            vtN = trTria.GetN() ;
+        // aggiungo al vettore dei proiettati
+         vMyPt5ax.emplace_back( ptInt, vtN, dPar, 1) ;
+      }
+      bFound = PL.GetNextUPoint( &dPar, &ptP) ;
+   }
+
+  // rimuovo i punti allineati entro la tolleranza e non più lontani tra loro del massimo
+   RemovePointsInExcess( vMyPt5ax, dLinTol, dMaxSegmLen) ;
+
+  // copio i punti rimasti nel vettore di ritorno
+   vPt5ax.clear() ;
+   for ( const auto& Pt5ax : vMyPt5ax) {
+      if ( Pt5ax.nFlag != -1)
+         vPt5ax.emplace_back( Pt5ax) ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ProjectCurveOnSurfTm( const ICurve& crCrv, const ISurfTriMesh& tmSurf, const IGeoPoint3d& gpRef,
+                      double dLinTol, double dMaxSegmLen, PNT5AXVECTOR& vPt5ax)
+{
+  // controllo le tolleranze
+   dLinTol = max( dLinTol, LIN_TOL_MIN) ;
+   dMaxSegmLen = max( dMaxSegmLen, 10 * EPS_SMALL) ;
+
+  // approssimo la curva con una polilinea entro la metà della tolleranza
+   PolyLine PL ;
+   if ( ! crCrv.ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_STD, PL))
+      return false ;
+   const double MAX_SEG_LEN = min( dMaxSegmLen, 1.) ;
+   if ( ! PL.AdjustForMaxSegmentLen( MAX_SEG_LEN))
+      return false ;
+
+  // Vettore locale dei punti risultanti
+   PNT5AXVECTOR vMyPt5ax ;
+   vMyPt5ax.reserve( PL.GetPointNbr()) ;
+
+  // proietto i punti della polilinea sulla superficie con direzione data dal punto di riferimento
+   double dPar ;
+   Point3d ptP ;
+   bool bFound = PL.GetFirstUPoint( &dPar, &ptP) ;
+   while ( bFound) {
+     // punto di riferimento
+      Point3d ptMin = gpRef.GetPoint() ;
+     // intersezione della retta di minima distanza con la superficie
+      Vector3d vtLine = ptP - ptMin ;
+      double dLineLen = vtLine.Len() ;
+      if ( dLineLen > EPS_SMALL) {        
+         vtLine /= dLineLen ;
+         ILSIVECTOR vIntRes ;
+         if ( IntersLineSurfTm( ptP, vtLine, dLineLen, tmSurf, vIntRes, false)) {
+            if ( vIntRes.size() > 0) {
+              // calcolo il punto
+               int nI = int( vIntRes.size()) - 1 ;
+               Point3d ptInt ;
+               if ( vIntRes[nI].nILTT == ILTT_SEGM  ||  vIntRes[nI].nILTT == ILTT_SEGM_ON_EDGE)
+                  ptInt = vIntRes[nI].ptI2 ;
+               else
+                  ptInt = vIntRes[nI].ptI ;
+              // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
+               Triangle3dEx trTria ;
+               if ( ! tmSurf.GetTriangle( vIntRes[nI].nT, trTria))
+                  return false ;
+               Vector3d vtN ;
+               double dU, dV, dW ;
+               if ( BarycentricCoord( ptInt, trTria, dU, dV, dW)) 
+                  vtN = dU * trTria.GetVertexNorm( 0) + dV * trTria.GetVertexNorm( 1) + dW * trTria.GetVertexNorm( 2) ;
+               if ( ! vtN.Normalize())
+                  vtN = trTria.GetN() ;
+              // aggiungo al vettore dei proiettati
+               vMyPt5ax.emplace_back( ptInt, vtN, vtLine, dPar, 1) ;
+            }
+         }
+      }
+      bFound = PL.GetNextUPoint( &dPar, &ptP) ;
+   }
+
+  // rimuovo i punti allineati entro la tolleranza e non più lontani tra loro del massimo
+   RemovePointsInExcess( vMyPt5ax, dLinTol, dMaxSegmLen) ;
+
+  // copio i punti rimasti nel vettore di ritorno
+   vPt5ax.clear() ;
+   for ( const auto& Pt5ax : vMyPt5ax) {
+      if ( Pt5ax.nFlag != -1)
+         vPt5ax.emplace_back( Pt5ax) ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ProjectCurveOnSurfTm( const ICurve& crCrv, const ISurfTriMesh& tmSurf, const ICurve& crRef,
+                      double dLinTol, double dMaxSegmLen, PNT5AXVECTOR& vPt5ax)
+{
+  // controllo le tolleranze
+   dLinTol = max( dLinTol, LIN_TOL_MIN) ;
+   dMaxSegmLen = max( dMaxSegmLen, 10 * EPS_SMALL) ;
+
+  // approssimo la curva con una polilinea entro la metà della tolleranza
+   PolyLine PL ;
+   if ( ! crCrv.ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_STD, PL))
+      return false ;
+   const double MAX_SEG_LEN = min( dMaxSegmLen, 1.) ;
+   if ( ! PL.AdjustForMaxSegmentLen( MAX_SEG_LEN))
+      return false ;
+
+  // Vettore locale dei punti risultanti
+   PNT5AXVECTOR vMyPt5ax ;
+   vMyPt5ax.reserve( PL.GetPointNbr()) ;
+
+  // proietto i punti della polilinea sulla superficie con direzione normale alla curva di riferimento
+   double dPar ;
+   Point3d ptP ;
+   bool bFound = PL.GetFirstUPoint( &dPar, &ptP) ;
+   while ( bFound) {
+     // punto sulla curva a minima distanza
+      DistPointCurve dPC( ptP, crRef) ;
+      Point3d ptMin ;
+      int nFlag ;
+      if ( dPC.GetMinDistPoint( 0, ptMin, nFlag)) {
+        // intersezione della retta di minima distanza con la superficie
+         Vector3d vtLine = ptP - ptMin ;
+         double dLineLen = vtLine.Len() ;
+         if ( dLineLen > EPS_SMALL) {        
+            vtLine /= dLineLen ;
+            ILSIVECTOR vIntRes ;
+            if ( IntersLineSurfTm( ptP, vtLine, dLineLen, tmSurf, vIntRes, false)) {
+               if ( vIntRes.size() > 0) {
+                 // calcolo il punto
+                  int nI = int( vIntRes.size()) - 1 ;
+                  Point3d ptInt ;
+                  if ( vIntRes[nI].nILTT == ILTT_SEGM  ||  vIntRes[nI].nILTT == ILTT_SEGM_ON_EDGE)
+                     ptInt = vIntRes[nI].ptI2 ;
+                  else
+                     ptInt = vIntRes[nI].ptI ;
+                 // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
+                  Triangle3dEx trTria ;
+                  if ( ! tmSurf.GetTriangle( vIntRes[nI].nT, trTria))
+                     return false ;
+                  Vector3d vtN ;
+                  double dU, dV, dW ;
+                  if ( BarycentricCoord( ptInt, trTria, dU, dV, dW)) 
+                     vtN = dU * trTria.GetVertexNorm( 0) + dV * trTria.GetVertexNorm( 1) + dW * trTria.GetVertexNorm( 2) ;
+                  if ( ! vtN.Normalize())
+                     vtN = trTria.GetN() ;
+                 // aggiungo al vettore dei proiettati
+                  vMyPt5ax.emplace_back( ptInt, vtN, vtLine, dPar, 1) ;
+               }
+            }
+         }
+      }
+      bFound = PL.GetNextUPoint( &dPar, &ptP) ;
+   }
+
+  // rimuovo i punti allineati entro la tolleranza e non più lontani tra loro del massimo
+   RemovePointsInExcess( vMyPt5ax, dLinTol, dMaxSegmLen) ;
+
+  // copio i punti rimasti nel vettore di ritorno
+   vPt5ax.clear() ;
+   for ( const auto& Pt5ax : vMyPt5ax) {
+      if ( Pt5ax.nFlag != -1)
+         vPt5ax.emplace_back( Pt5ax) ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ProjectCurveOnSurfTm( const ICurve& crCrv, const ISurfTriMesh& tmSurf, const ISurfTriMesh& tmRef,
+                      double dLinTol, double dMaxSegmLen, PNT5AXVECTOR& vPt5ax)
+{
+  // controllo le tolleranze
+   dLinTol = max( dLinTol, LIN_TOL_MIN) ;
+   dMaxSegmLen = max( dMaxSegmLen, 10 * EPS_SMALL) ;
+
+  // approssimo la curva con una polilinea entro la metà della tolleranza
+   PolyLine PL ;
+   if ( ! crCrv.ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_STD, PL))
+      return false ;
+   const double MAX_SEG_LEN = min( dMaxSegmLen, 1.) ;
+   if ( ! PL.AdjustForMaxSegmentLen( MAX_SEG_LEN))
+      return false ;
+
+  // Vettore locale dei punti risultanti
+   PNT5AXVECTOR vMyPt5ax ;
+   vMyPt5ax.reserve( PL.GetPointNbr()) ;
+
+  // proietto i punti della polilinea sulla superficie con direzione normale alla curva di riferimento
+   double dPar ;
+   Point3d ptP ;
+   bool bFound = PL.GetFirstUPoint( &dPar, &ptP) ;
+   while ( bFound) {
+     // punto sulla superficie guida a minima distanza
+      DistPointSurfTm dPS( ptP, tmRef) ;
+      Point3d ptMin ;
+      int nTriaMin ;
+      if ( dPS.GetMinDistPoint( ptMin)  &&  dPS.GetMinDistTriaIndex ( nTriaMin)) {
+        // intersezione della retta di minima distanza con la superficie
+         Vector3d vtLine = ptP - ptMin ;
+         double dLineLen = vtLine.Len() ;
+         if ( dLineLen > EPS_SMALL) {        
+            vtLine /= dLineLen ;
+            ILSIVECTOR vIntRes ;
+            if ( IntersLineSurfTm( ptP, vtLine, dLineLen, tmSurf, vIntRes, false)) {
+               if ( vIntRes.size() > 0) {
+                 // calcolo il punto
+                  int nI = int( vIntRes.size()) - 1 ;
+                  Point3d ptInt ;
+                  if ( vIntRes[nI].nILTT == ILTT_SEGM  ||  vIntRes[nI].nILTT == ILTT_SEGM_ON_EDGE)
+                     ptInt = vIntRes[nI].ptI2 ;
+                  else
+                     ptInt = vIntRes[nI].ptI ;
+                 // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
+                  Triangle3dEx trTria ;
+                  if ( ! tmSurf.GetTriangle( vIntRes[nI].nT, trTria))
+                     return false ;
+                  Vector3d vtN ;
+                  double dU, dV, dW ;
+                  if ( BarycentricCoord( ptInt, trTria, dU, dV, dW)) 
+                     vtN = dU * trTria.GetVertexNorm( 0) + dV * trTria.GetVertexNorm( 1) + dW * trTria.GetVertexNorm( 2) ;
+                  if ( ! vtN.Normalize())
+                     vtN = trTria.GetN() ;
+                 // calcolo la normale della superficie guida
+                  Triangle3dEx trGuide ;
+                  if ( ! tmRef.GetTriangle( nTriaMin, trGuide))
+                     return false ;
+                  Vector3d vtN2 ;
+                  double dU2, dV2, dW2 ;
+                  if ( BarycentricCoord( ptMin, trGuide, dU2, dV2, dW2)) 
+                     vtN2 = dU2 * trGuide.GetVertexNorm( 0) + dV2 * trGuide.GetVertexNorm( 1) + dW2 * trGuide.GetVertexNorm( 2) ;
+                  if ( ! vtN2.Normalize())
+                     vtN2 = trGuide.GetN() ;
+                 // aggiungo al vettore dei proiettati
+                  vMyPt5ax.emplace_back( ptInt, vtN, vtN2, dPar, 1) ;
+               }
+            }
+         }
+      }
+      bFound = PL.GetNextUPoint( &dPar, &ptP) ;
+   }
+
+  // rimuovo i punti allineati entro la tolleranza e non più lontani tra loro del massimo
+   RemovePointsInExcess( vMyPt5ax, dLinTol, dMaxSegmLen) ;
 
   // copio i punti rimasti nel vettore di ritorno
    vPt5ax.clear() ;

Quaternion.cpp

@@ -0,0 +1,158 @@
+//----------------------------------------------------------------------------
+//  EgalTech 2024-2024
+//----------------------------------------------------------------------------
+// File : Quaternion.cpp          Data : 13.04.24          Versione : 2.6d4
+// Contenuto : Funzioni della classe Quaternion.
+//
+//
+//
+// Modifiche : 13.04.24 DS Creazione modulo.
+//
+//
+//----------------------------------------------------------------------------
+
+//--------------------------- Include ----------------------------------------
+#include "stdafx.h"
+#include "\EgtDev\Include\EGkQuaternion.h"
+
+//----------------------------------------------------------------------------
+// Lunghezza o Modulo
+//----------------------------------------------------------------------------
+double 
+Quaternion::Len( void) const
+{
+   if ( abs( x) < EPS_ZERO  &&  abs( y) < EPS_ZERO  &&  abs( z) < EPS_ZERO)
+      return abs( w) ;
+   if ( abs( w) < EPS_ZERO  &&  abs( y) < EPS_ZERO  &&  abs( z) < EPS_ZERO)
+      return abs( x) ;
+   if ( abs( w) < EPS_ZERO  &&  abs( z) < EPS_ZERO  &&  abs( x) < EPS_ZERO)
+      return abs( y) ;
+   if ( abs( w) < EPS_ZERO  &&  abs( x) < EPS_ZERO  &&  abs( y) < EPS_ZERO)
+      return abs( z) ;
+
+   return sqrt( w * w + x * x + y * y + z * z) ; 
+}
+
+//----------------------------------------------------------------------------
+// Normalizzazione
+//----------------------------------------------------------------------------
+bool
+Quaternion::Normalize( double dEps)
+{
+  // se già normalizzato, ok
+   double dSqLen = w * w + x * x + y * y + z * z ;
+   if ( abs( 1.0 - dSqLen) < ( 2 * 1000 * DBL_EPSILON))
+      return true ;
+
+  // se troppo piccolo, errore
+   if ( dSqLen < ( dEps * dEps))
+      return false ;
+
+  // eseguo la normalizzazione
+   double dLen = sqrt( dSqLen) ;
+   *this /= dLen ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+Quaternion
+FromAxisAngle( const Vector3d& vtAx, double dAngDeg)
+{
+   if ( abs( dAngDeg) < EPS_ANG_ZERO)
+      return Q_UNIT ;
+   double dLen = vtAx.Len() ;
+   if ( dLen < EPS_ZERO)
+      return Q_NULL ;
+   double dDenom = 1 / dLen ;
+   double dHCos = cos( dAngDeg / 2 * DEGTORAD) ;
+   double dHSin = sin( dAngDeg / 2 * DEGTORAD) ;
+   return Quaternion( dHCos, dHSin * vtAx.x * dDenom, dHSin * vtAx.y * dDenom, dHSin * vtAx.z * dDenom) ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ToAxisAngle( const Quaternion& qtQ, Vector3d& vtAx, double& dAngDeg)
+{
+   if ( ! qtQ.IsNormalized())
+      return false ;
+   if ( qtQ.IsUnit()  ||  (-qtQ).IsUnit()) {
+      dAngDeg = 0 ;
+      vtAx = Z_AX ;
+   }
+   else {
+      dAngDeg = 2 * acos( qtQ.w) * RADTODEG ;
+      double dDenom = 1 / sqrt( 1 - qtQ.w * qtQ.w) ;
+      vtAx.x = qtQ.x * dDenom ;
+      vtAx.y = qtQ.y * dDenom ;
+      vtAx.z = qtQ.z * dDenom ;
+   }
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+Quaternion
+FromFrame( const Frame3d& frRef)
+{
+  // verifico il riferimento
+   if ( ! frRef.IsValid())
+      return Q_NULL ;
+  // eseguo il calcolo
+   Quaternion qtQ ;
+  // traccia della matrice di rotazione
+   double dTr = frRef.VersX().x + frRef.VersY().y + frRef.VersZ().z ;
+  // se traccia positiva o nulla
+   if ( dTr > - 10 * EPS_ZERO) {
+      double dS = sqrt( 1 + dTr) ;
+      qtQ.w = dS / 2 ;
+      dS = 1 / ( 2 * dS) ;
+      qtQ.x = ( frRef.VersY().z - frRef.VersZ().y) * dS ;
+      qtQ.y = ( frRef.VersZ().x - frRef.VersX().z) * dS ;
+      qtQ.z = ( frRef.VersX().y - frRef.VersY().x) * dS ;
+   }
+  // altrimenti traccia negativa
+   else {
+      if ( frRef.VersX().x > frRef.VersY().y - 10 * EPS_ZERO  &&  frRef.VersX().x > frRef.VersZ().z - 10 * EPS_ZERO) {
+         double dS = sqrt( 1 + frRef.VersX().x - ( frRef.VersY().y + frRef.VersZ().z)) ;
+         qtQ.x = dS / 2 ;
+         dS = 1 / ( 2 * dS) ;
+         qtQ.y = ( frRef.VersY().x + frRef.VersX().y) * dS ;
+         qtQ.z = ( frRef.VersX().z + frRef.VersZ().x) * dS ;
+         qtQ.w = ( frRef.VersY().z - frRef.VersZ().y) * dS ;
+      }
+      else if ( frRef.VersY().y > frRef.VersZ().z - 10 * EPS_ZERO  &&  frRef.VersY().y > frRef.VersX().x - 10 * EPS_ZERO) {
+         double dS = sqrt( 1 + frRef.VersY().y - ( frRef.VersZ().z + frRef.VersX().x)) ;
+         qtQ.y = dS / 2 ;
+         dS = 1 / ( 2 * dS) ;
+         qtQ.z = ( frRef.VersZ().y + frRef.VersY().z) * dS ;
+         qtQ.x = ( frRef.VersY().x + frRef.VersX().y) * dS ;
+         qtQ.w = ( frRef.VersZ().x - frRef.VersX().z) * dS ;
+      }
+      else {
+         double dS = sqrt( 1 + frRef.VersZ().z - ( frRef.VersX().x + frRef.VersY().y)) ;
+         qtQ.z = dS / 2 ;
+         dS = 1 / ( 2 * dS) ;
+         qtQ.x = ( frRef.VersX().z + frRef.VersZ().x) * dS ;
+         qtQ.y = ( frRef.VersZ().y + frRef.VersY().z) * dS ;
+         qtQ.w = ( frRef.VersX().y - frRef.VersY().x) * dS ;
+      }
+   }
+
+   return qtQ ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ToFrame( const Quaternion& qtQ, Frame3d& frRef)
+{
+   double dNrm = qtQ.Len() ;
+   double dS = ( dNrm > EPS_ZERO ? 2 / dNrm : 0) ;
+   double dXs = qtQ.x * dS, dYs = qtQ.y * dS, dZs = qtQ.z * dS ;
+   double dWx = qtQ.w * dXs, dWy = qtQ.w * dYs, dWz = qtQ.w * dZs ;
+   double dXx = qtQ.x * dXs, dXy = qtQ.x * dYs, dXz = qtQ.x * dZs ;
+   double dYy = qtQ.y * dYs, dYz = qtQ.y * dZs, dZz = qtQ.z * dZs ;
+   Vector3d vtX( 1 - ( dYy + dZz), dXy + dWz, dXz - dWy) ;
+   Vector3d vtY( dXy - dWz, 1 - ( dXx + dZz), dYz + dWx) ;
+   Vector3d vtZ( dXz + dWy, dYz - dWx, 1 - ( dXx + dYy)) ;
+   return frRef.Set( ORIG, vtX, vtY, vtZ) ;
+}

RotationMinimizingFrame.cpp

@@ -0,0 +1,231 @@
+//----------------------------------------------------------------------------
+//  EgalTech 2024-2024
+//----------------------------------------------------------------------------
+// File : RotationMinimizingFrame.cpp   Data : 05.03.24      Versione : 2.6d1
+// Contenuto : Classe per RotationMinimizeFrame
+//
+//
+//
+// Modifiche : 05.03.24 RE Creazione modulo.
+//
+//
+//----------------------------------------------------------------------------
+
+//--------------------------- Include ----------------------------------------
+#include "stdafx.h"
+#include "GeoConst.h"
+#include "/EgtDev/Include/EGkRotationMinimizingFrame.h"
+
+using namespace std ;
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::Set( const ICurve* pCrv, const Frame3d& fr_Start)
+{
+  // pulisco
+   Clear() ;
+  // verifico i parametri
+   if ( pCrv == nullptr  ||  ! pCrv->IsValid()  ||
+        ! fr_Start.IsValid())
+      return false ;
+  // assegno i parametri
+   m_pCrv = pCrv->Clone() ;
+   m_Frame0 = fr_Start ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::Clear( void)
+{
+  // pulizia della curva
+   if ( m_pCrv != nullptr)
+      delete m_pCrv ;
+   m_pCrv = nullptr ;
+  // reset del frame di partenza
+   m_Frame0.Reset() ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::IsValid( void)
+{
+  // controllo validità della curva
+   if ( m_pCrv == nullptr  ||  ! m_pCrv->IsValid())
+      return false ;
+
+  // controllo del frame iniziale
+   if ( ! m_Frame0.IsValid())
+      return false ;
+
+  // controllo che l'origine del frame sia sulla curva e che l'asse Z sia tangente alla curva
+   Point3d ptS ;
+   Vector3d vtZ ;
+   if ( ! m_pCrv->GetPointD1D2( 0., ICurve::FROM_MINUS, ptS, &vtZ) ||
+        ! vtZ.Normalize()  ||
+        ! AreSamePointApprox( ptS, m_Frame0.Orig())  ||
+        ! AreSameVectorEpsilon( vtZ, m_Frame0.VersZ(), 5 * EPS_SMALL))
+      return false ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::GetFrameAtParam( const Frame3d& frAct, const double dParNext, Frame3d& frNext)
+{
+  /*
+                                 Double Reflection
+            Computation of Rotation Minimizing Frame in Computer Graphics
+                  Wenping Wang Bert Juttler Dayue Zheng Yang Liu
+  */
+
+  // ricavo i parametri dal frame
+   if ( ! frAct.IsValid())
+      return false ;
+
+  // origine del frame e versori ( ptCurr, [ vtCurrR, vtCurrS, vtCurrT])
+   Point3d ptCurr = frAct.Orig() ;
+   Vector3d vtCurrR = frAct.VersX() ;
+   Vector3d vtCurrT = frAct.VersZ() ;
+
+  // punto i-esimo sulla curva e suo vettore tangente
+   Point3d ptNext ;
+   Vector3d vtNextT ;
+   if ( ! m_pCrv->GetPointD1D2( dParNext, ICurve::FROM_MINUS, ptNext, &vtNextT)  ||
+        ! vtNextT.Normalize())
+      return false ;
+
+  // controllo per casi degeneri
+   if ( AreSamePointEpsilon( ptCurr, ptNext, EPS_ZERO)  ||               // non esiste il piano R1
+        abs( ( ptNext - ptCurr) * ( vtNextT + vtCurrT)) < EPS_ZERO)      // non esiste il piano R2
+      return false ;
+
+  // ricavo il vettore di riflessione rispetto al piano R1
+   Vector3d vR1_norm = ptNext - ptCurr ;
+   if ( ! vR1_norm.IsValid())
+      return false ;
+
+  // parametro di riflessione per R1
+   double dPar1 = vR1_norm * vR1_norm ;
+
+  // riflessione rispetto al piano R1 ( sistema sinistrorso L )
+   Vector3d vt_r_L = vtCurrR - ( 2 / dPar1) * ( vR1_norm * vtCurrR) * vR1_norm ;
+   Vector3d vt_t_L = vtCurrT - ( 2 / dPar1) * ( vR1_norm * vtCurrT) * vR1_norm ;
+
+  // ricavo il vettore di riflessione rispetto al piano R1
+   Vector3d vR2_norm = vtNextT - vt_t_L ;
+   if ( ! vR2_norm.IsValid())
+      return false ;
+
+  // parametro di riflessione per R2
+   double dPar2 = vR2_norm * vR2_norm ;
+  // versore r del nuovo frame
+   Vector3d vt_r_next = vt_r_L - ( 2 / dPar2) * ( vR2_norm * vt_r_L) * vR2_norm ;
+
+  // imposto il nuovo frame
+   return frNext.Set( ptNext, vtNextT, vt_r_next) ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::GetFramesByStep( double dStep, bool bUniform, FRAME3DVECTOR& vRMFrames)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+  // controllo sullo step
+   dStep = max( 10 * EPS_SMALL, dStep) ;
+
+  // lunghezza della curva
+   double dCrvLen = 0. ;
+   if ( ! m_pCrv->GetLength( dCrvLen)  ||  dCrvLen < 10 * EPS_SMALL)
+      return false ;
+
+  // ricavo il numero degli step
+   int nStep = int( ceil( dCrvLen / dStep)) ;
+   double dMyStep = ( bUniform ? dCrvLen / nStep : dStep) ;
+
+  // inserisco il frame iniziale nel vettore dei riferimenti
+   vRMFrames.clear() ;
+   vRMFrames.reserve( nStep + 1) ;
+   vRMFrames.push_back( m_Frame0) ;
+
+  // ciclo sugli step in cui la curva è suddivisa
+   for ( int i = 1 ; i <= nStep ; ++ i) {
+     // ricavo il parametro della curva allo step i-esimo
+      double dParNext ;
+      if ( ! m_pCrv->GetParamAtLength( min( i * dMyStep, dCrvLen - EPS_SMALL), dParNext))
+         return false ;
+     // ricavo il frame alla posizione calcolata
+      Frame3d frNext ;
+      if ( ! GetFrameAtParam( vRMFrames[i-1], dParNext, frNext))
+         return false ;
+     // inserisco nuovo frame nel vettore dei riferimenti
+      vRMFrames.push_back( frNext) ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::GetFramesBySplit( int nIntervals, FRAME3DVECTOR& vRMFrames)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+  // controllo sul numero di intervalli
+   nIntervals = max( 1, nIntervals) ;
+
+  // ricavo lo step
+   double dLen = 0 ;
+   if ( ! m_pCrv->GetLength( dLen))
+      return false ;
+   double dStep = dLen / nIntervals ;
+
+   return GetFramesByStep( dStep, true, vRMFrames) ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationMinimizingFrame::GetFramesByTolerance( double dTol, FRAME3DVECTOR& vRMFrames)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+  // controllo sulla tolleranza
+   dTol = max( EPS_SMALL, dTol) ;
+
+  // ricavo la PolyLine associata alla curva mediante tale tolleranza
+   PolyLine PL ;
+   if ( ! m_pCrv->ApproxWithLines( dTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, PL))
+      return false ;
+   int nStep = PL.GetLineNbr() ;
+
+  // inserisco il frame iniziale nel vettore dei riferimenti
+   vRMFrames.clear() ;
+   vRMFrames.reserve( nStep + 1) ;
+   vRMFrames.push_back( m_Frame0) ;
+
+  // eseguo il calcolo dei frame su ogni punto ricavato dall'approssimazione
+   Point3d ptCurr ;
+   PL.GetFirstPoint( ptCurr) ;
+   double dParNext ;
+   Point3d ptNext ;
+   while ( PL.GetNextUPoint( &dParNext, &ptNext)) {
+     // ricavo il Frame associato a questa posizione
+      Frame3d frNext ;
+      if ( ! GetFrameAtParam( vRMFrames.back(), dParNext, frNext))
+         return false ;
+      vRMFrames.emplace_back( frNext) ;
+   }
+
+   return true ;
+}

RotationXplaneFrame.cpp

@@ -0,0 +1,223 @@
+//----------------------------------------------------------------------------
+//  EgalTech 2024-2024
+//----------------------------------------------------------------------------
+// File : RotationXplaneFrame.cpp      Data : 05.04.24      Versione : 2.6d1
+// Contenuto : Classe per RotationXplaneFrame.
+//
+//
+//
+// Modifiche : 05.04.24 DS Creazione modulo.
+//
+//
+//----------------------------------------------------------------------------
+
+//--------------------------- Include ----------------------------------------
+#include "stdafx.h"
+#include "GeoConst.h"
+#include "/EgtDev/Include/EGkRotationXplaneFrame.h"
+
+using namespace std ;
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::Set( const ICurve* pCrv, const Vector3d& vtNorm, const Vector3d& vtNearX)
+{
+  // pulisco
+   Clear() ;
+  // verifico i parametri
+   if ( pCrv == nullptr  ||  ! pCrv->IsValid()  ||
+        ! vtNorm.IsValid()  ||  vtNorm.IsSmall()  ||  ! vtNearX.IsValid())
+      return false ;
+  // assegno i parametri
+   m_pCrv = pCrv->Clone() ;
+   m_vtNorm = vtNorm ;
+   m_vtNearX = vtNearX ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::Clear( void)
+{
+  // pulizia della curva
+   if ( m_pCrv != nullptr)
+      delete m_pCrv ;
+   m_pCrv = nullptr ;
+  // reset dei vettori
+   m_vtNorm = V_NULL ;
+   m_vtNearX = V_NULL ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::IsValid( void)
+{
+  // controllo validità della curva
+   if ( m_pCrv == nullptr  ||  ! m_pCrv->IsValid())
+      return false ;
+
+  // controllo della normale al piano
+   if ( m_vtNorm.IsSmall())
+      return false ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::GetFrameAtParam( const Frame3d& frAct, const double dParNext, Frame3d& frNext)
+{
+  // verifico il riferimento corrente
+   if ( ! frAct.IsValid())
+      return false ;
+
+  // recupero punto e tangente nel nuovo punto
+   Point3d ptNext ;
+   Vector3d vtNext ;
+   if ( ! m_pCrv->GetPointD1D2( dParNext, ICurve::FROM_MINUS, ptNext, &vtNext)  ||
+        ! vtNext.Normalize())
+      return false ;
+
+  // calcolo il nuovo riferimento
+   Vector3d vtAxX = m_vtNorm ^ vtNext ;
+   if ( vtAxX.IsSmall())
+      vtAxX = frAct.VersX() ;
+   else if ( vtAxX * frAct.VersX() < 0)
+      vtAxX.Invert() ;
+  
+  // lo imposto
+   return frNext.Set( ptNext, vtNext, vtAxX) ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::GetFramesByStep( double dStep, bool bUniform, FRAME3DVECTOR& vRXFrames)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+  // controllo sullo step
+   dStep = max( 10 * EPS_SMALL, dStep) ;
+
+  // lunghezza della curva
+   double dCrvLen = 0. ;
+   if ( ! m_pCrv->GetLength( dCrvLen)  ||  dCrvLen < 10 * EPS_SMALL)
+      return false ;
+
+  // numero e lunghezza effettiva di ogni step
+   int nStep = int( ceil( dCrvLen / dStep)) ;
+   double dMyStep = ( bUniform ? dCrvLen / nStep : dStep) ;
+
+  // calcolo il riferimento iniziale
+   Point3d ptStart ; m_pCrv->GetStartPoint( ptStart) ;
+   Vector3d vtStart ; m_pCrv->GetStartDir( vtStart) ;
+   Vector3d vtAxX = m_vtNorm ^ vtStart ;
+   if ( vtAxX.IsSmall()) {
+      vtAxX = OrthoCompo( m_vtNearX, m_vtNorm) ;
+      if ( vtAxX.IsSmall()) {
+         vtAxX = FromUprightOrtho( m_vtNorm) ;
+         vtAxX.Rotate( m_vtNorm, 0, 1) ;
+      }
+   }
+   Frame3d frStart ;
+   if ( ! frStart.Set( ptStart, vtStart, vtAxX))
+      return false ;
+
+  // inserisco questo frame nel vettore dei riferimenti
+   vRXFrames.clear() ;
+   vRXFrames.reserve( nStep + 1) ;
+   vRXFrames.push_back( frStart) ;
+
+  // ciclo sugli step in cui la curva è suddivisa
+   for ( int i = 1 ; i <= nStep ; ++ i) {
+     // ricavo il parametro della curva allo step i-esimo
+      double dParNext ;
+      if ( ! m_pCrv->GetParamAtLength( min( i * dMyStep, dCrvLen - EPS_SMALL), dParNext))
+         return false ;
+     // ricavo il frame alla posizione calcolata
+      Frame3d frNext ;
+      if ( ! GetFrameAtParam( vRXFrames[i-1], dParNext, frNext))
+         return false ;
+     // inserisco nuovo frame nel vettore dei riferimenti
+      vRXFrames.push_back( frNext) ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::GetFramesBySplit( int nIntervals, FRAME3DVECTOR& vRXFrames)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+  // controllo sul numero di intervalli
+   nIntervals = max( 1, nIntervals) ;
+
+  // ricavo lo step
+   double dLen = 0 ;
+   if ( ! m_pCrv->GetLength( dLen))
+      return false ;
+   double dStep = dLen / nIntervals ;
+
+   return GetFramesByStep( dStep, true, vRXFrames) ;
+}
+
+//----------------------------------------------------------------------------
+bool
+RotationXplaneFrame::GetFramesByTolerance( double dTol, FRAME3DVECTOR& vRXFrames)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+  // controllo sulla tolleranza
+   dTol = max( EPS_SMALL, dTol) ;
+
+  // ricavo la PolyLine associata alla curva mediante tale tolleranza
+   PolyLine PL ;
+   if ( ! m_pCrv->ApproxWithLines( dTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, PL))
+      return false ;
+   int nStep = PL.GetLineNbr() ;
+
+  // calcolo il riferimento iniziale
+   Point3d ptStart ; m_pCrv->GetStartPoint( ptStart) ;
+   Vector3d vtStart ; m_pCrv->GetStartDir( vtStart) ;
+   Vector3d vtAxX = m_vtNorm ^ vtStart ;
+   if ( vtAxX.IsSmall()) {
+      vtAxX = OrthoCompo( m_vtNearX, m_vtNorm) ;
+      if ( vtAxX.IsSmall()) {
+         vtAxX = FromUprightOrtho( m_vtNorm) ;
+         vtAxX.Rotate( m_vtNorm, 0, 1) ;
+      }
+   }
+   Frame3d frStart ;
+   if ( ! frStart.Set( ptStart, vtStart, vtAxX))
+      return false ;
+
+  // inserisco questo frame nel vettore dei riferimenti
+   vRXFrames.clear() ;
+   vRXFrames.reserve( nStep + 1) ;
+   vRXFrames.push_back( frStart) ;
+
+  // eseguo il calcolo dei frame sui punti della polyline approssimante
+   Point3d ptCurr ;
+   PL.GetFirstPoint( ptCurr) ;
+   double dParNext ;
+   Point3d ptNext ;
+   while ( PL.GetNextUPoint( &dParNext, &ptNext)) {
+     // ricavo il Frame associato a questa posizione
+      Frame3d frNext ;
+      if ( ! GetFrameAtParam( vRXFrames.back(), dParNext, frNext))
+         return false ;
+      vRXFrames.emplace_back( frNext) ;
+   }
+
+   return true ;
+}

StmFromCurves.cpp

@@ -1,4 +1,4 @@
-//----------------------------------------------------------------------------
+//----------------------------------------------------------------------------
 //  EgalTech 2015-2015
 //----------------------------------------------------------------------------
 // File : StmFromCurves.cpp    Data : 01.02.15       Versione : 1.6b1
@@ -19,14 +19,15 @@
 #include "CurveComposite.h"
 #include "SurfTriMesh.h"
 #include "Voronoi.h"
+#include "/EgtDev/Include/EGkSfrCreate.h"
 #include "/EgtDev/Include/EGkOffsetCurve.h"
 #include "/EgtDev/Include/EGkStmFromCurves.h"
 #include "/EgtDev/Include/EGkStmFromTriangleSoup.h"
+#include "/EgtDev/Include/EGkRotationMinimizingFrame.h"
+#include "/EgtDev/Include/EGkRotationXplaneFrame.h"
+#include "/EgtDev/Include/EGkIntersCurves.h"
 #include "/EgtDev/Include/EgtPointerOwner.h"
-#include <algorithm>
-#include <thread>
 #include <future>
-#include <EgtDev/Include/EGkIntersCurves.h>
 
 using namespace std ;
 
@@ -110,9 +111,9 @@ GetSurfTriMeshByExtrusion( const ICurve* pCurve, const Vector3d& vtExtr,
    PtrOwner<SurfTriMesh> pSTM( CreateBasicSurfTriMesh()) ;
    if ( IsNull( pSTM)  ||  ! pSTM->CreateByExtrusion( PL, vtExtr))
       return nullptr ;
-  // se da fare, metto i tappi sulle estremità
+  // se da fare, metto i tappi sulle estremità
    if ( bDoCapEnds) {
-     // creo la prima superficie di estremità
+     // creo la prima superficie di estremità
       SurfTriMesh STM1 ;
       if ( ! STM1.CreateByFlatContour( PL))
          return nullptr ;
@@ -156,7 +157,7 @@ GetSurfTriMeshByRegionExtrusion( const CICURVEPVECTOR& vpCurve, const Vector3d&
       for ( int i = 0 ; i < int( vPL.size()) ; ++ i)
          vPL[i].Invert() ;
    }
-  // creo la prima superficie di estremità
+  // creo la prima superficie di estremità
    PtrOwner<SurfTriMesh> pSTM( CreateBasicSurfTriMesh()) ;
    if ( IsNull( pSTM)  ||  ! pSTM->CreateByRegion( vPL))
       return nullptr ;
@@ -288,7 +289,7 @@ GetSurfTriMeshByScrewing( const ICurve* pCurve, const Point3d& ptAx, const Vecto
       return nullptr ;
   // se richiesti caps
    if ( bCapEnds) {
-     // determino se la sezione è chiusa e piatta
+     // determino se la sezione è chiusa e piatta
       Plane3d plPlane ; double dArea ;
       bool bSectClosedFlat = PL.IsClosedAndFlat( plPlane, dArea, 10 * EPS_SMALL) ;
      // determino non sia una semplice rivoluzione
@@ -332,7 +333,7 @@ GetSurfTriMeshSharpRectSwept( double dDimH, double dDimV, const ICurve* pGuide,
    if ( ! pGuide->IsFlat( plGuide, false, 10 * EPS_SMALL))
       return nullptr ;
    Vector3d vtNorm = plGuide.GetVersN() ;
-  // determino se la guida è chiusa
+  // determino se la guida è chiusa
    bool bGuideClosed = pGuide->IsClosed() ;
   // curve di offset
    OffsetCurve OffsCrvR ;
@@ -373,7 +374,7 @@ GetSurfTriMeshSharpRectSwept( double dDimH, double dDimV, const ICurve* pGuide,
    pSTM->SetSmoothAngle( 20) ;
   // se guida aperta e tappi piatti
    if ( ! bGuideClosed  &&  nCapType == RSCAP_FLAT) {
-     // verifico che le due estremità siano chiuse e piatte
+     // verifico che le due estremità siano chiuse e piatte
       POLYLINEVECTOR vPL ;
       if ( ! pSTM->GetLoops( vPL)  ||  vPL.size() != 2)
          return nullptr ;
@@ -452,7 +453,7 @@ GetSurfTriMeshBeveledRectSwept( double dDimH, double dDimV, double dBevelH, doub
    Point3d ptCen ;
    pGuide->GetStartPoint( ptCen) ;
    ptCen -= dDimV / 2 * vtNorm ;
-  // determino se la guida è chiusa
+  // determino se la guida è chiusa
    bool bGuideClosed = pGuide->IsClosed() ;
   // curve di offset
    const int NUM_OFFS = 4 ;
@@ -475,7 +476,7 @@ GetSurfTriMeshBeveledRectSwept( double dDimH, double dDimV, double dBevelH, doub
       }   
    }
    else {
-     // se Voronoi non è possibile calcolare gli offset di una stessa curva in parallelo 
+     // se Voronoi non è possibile calcolare gli offset di una stessa curva in parallelo 
       for ( int i = 0 ; i < NUM_OFFS && bOk ; ++ i)
          bOk = vOffsCrv[i].Make( pGuide, vDist[i], ICurve::OFF_FILLET) ;
    }
@@ -552,7 +553,7 @@ GetSurfTriMeshBeveledRectSwept( double dDimH, double dDimV, double dBevelH, doub
       StmFromTriangleSoup stmCapSoup ;
       if ( ! stmCapSoup.Start( nBuckets))
          return nullptr ;
-     // verifico che le due estremità siano chiuse e piatte
+     // verifico che le due estremità siano chiuse e piatte
       POLYLINEVECTOR vPL ;
       if ( ! pSTM->GetLoops( vPL)  ||  vPL.size() != 2)
          return nullptr ;
@@ -650,47 +651,42 @@ GetSurfTriMeshRectSwept( double dDimH, double dDimV, double dBevelH, double dBev
 }
 
 //-------------------------------------------------------------------------------
-ISurfTriMesh*
-GetSurfTriMeshSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEnds, double dLinTol)
+static ISurfTriMesh*
+GetSurfTriMeshSweptInPlane( const ICurve* pSect, const ICurve* pGuide, const Vector3d& vtNorm, bool bCapEnds, double dLinTol)
 {
-  // verifica parametri
-   if ( pSect == nullptr  ||  pGuide == nullptr)
-      return nullptr ;
-  // calcolo la polilinea che approssima la sezione
-   PolyLine PL ;
-   if ( ! pSect->ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, PL))
-      return nullptr ;
-  // determino se la sezione è chiusa
-   bool bSectClosed = PL.IsClosed() ;
-  // verifico che la linea guida sia piana
-   Plane3d plGuide ;
-   if ( ! pGuide->IsFlat( plGuide, false, 10 * EPS_SMALL))
-      return nullptr ;
-  // determino se la guida è chiusa
+  // determino se la sezione è chiusa
+   bool bSectClosed = pSect->IsClosed() ;
+  // determino se la guida è chiusa
    bool bGuideClosed = pGuide->IsClosed() ;
+
   // riferimento all'inizio della linea guida
    Frame3d frStart ;
    Point3d ptStart ;
    pGuide->GetStartPoint( ptStart) ;
    Vector3d vtStart ;
    pGuide->GetStartDir( vtStart) ;
-   Vector3d vtNorm = plGuide.GetVersN() ;
    frStart.Set( ptStart, -vtStart, vtStart ^ vtNorm) ;
-  // porto la sezione in questo riferimento e ve la appiattisco
+
+   // calcolo la polilinea che approssima la sezione
+   PolyLine PL ;
+   if ( ! pSect->ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, PL))
+      return nullptr ;
+
+   // porto la sezione in questo riferimento e ve la appiattisco
    if ( ! PL.ToLoc( frStart)  ||  ! PL.Flatten())
       return nullptr ;
-  // calcolo la superficie
-   PtrOwner<SurfTriMesh> pSTM( CreateBasicSurfTriMesh()) ;
-   if ( IsNull( pSTM))
+
+  // preparo collettore delle superfici componenti
+   StmFromTriangleSoup StmFts ;
+   if ( ! StmFts.Start())
       return nullptr ;
-  // salvo tolleranza lineare usata
-   pSTM->SetLinearTolerance( dLinTol) ;
+
   // superficie swept
    PtrOwner<ICurve> pPrevCrv ;
    Point3d ptP ;
    bool bPoint = PL.GetFirstPoint( ptP) ;
    while ( bPoint) {
-     // nuova curva
+     // nuova curva ( definita dall'Offset )
       OffsetCurve OffsCrv ;
       if ( ! OffsCrv.Make( pGuide, ptP.x, ICurve::OFF_FILLET)  ||  OffsCrv.GetCurveCount() == 0)
          return nullptr ;
@@ -698,21 +694,32 @@ GetSurfTriMeshSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEnds, d
       if ( IsNull( pCurrCrv))
          return nullptr ;
       pCurrCrv->Translate( ptP.y * frStart.VersY()) ;
-     // se esiste la curva precedente, costruisco la rigata (di tipo minima distanza)
+     // se esiste la curva precedente, costruisco la rigata ( di tipo minima distanza)
       if ( ! IsNull( pPrevCrv)) {
          PtrOwner<ISurfTriMesh> pSr( GetSurfTriMeshRuled( pPrevCrv, pCurrCrv, ISurfTriMesh::RLT_MINDIST, dLinTol)) ;
          if ( IsNull( pSr))
             return nullptr ;
-         pSTM->DoSewing( *pSr) ;
+        // inserisco nel collettore
+         StmFts.AddSurfTriMesh( *pSr) ;
       }
      // salvo la curva come prossima precedente
       pPrevCrv.Set( pCurrCrv) ;
      // prossimo punto
       bPoint = PL.GetNextPoint( ptP) ;
    }
+
+  // recupero la supeficie risultante
+   if ( ! StmFts.End())
+      return nullptr ;
+   PtrOwner<SurfTriMesh> pSTM( GetBasicSurfTriMesh( StmFts.GetSurf())) ;
+   if ( IsNull( pSTM))
+      return nullptr ;
+  // salvo tolleranza lineare usata
+   pSTM->SetLinearTolerance( dLinTol) ;
+
   // se richiesti caps e sezione chiusa e guida aperta
    if ( bCapEnds  &&  bSectClosed  &&  ! bGuideClosed) {
-     // verifico che le due estremità siano chiuse e piatte
+     // verifico che le due estremità siano chiuse e piatte
       POLYLINEVECTOR vPL ;
       if ( ! pSTM->GetLoops( vPL)  ||  vPL.size() != 2)
          return nullptr ;
@@ -726,6 +733,7 @@ GetSurfTriMeshSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEnds, d
       if ( IsNull( pSci)  ||  ! pSci->CreateByFlatContour( PL))
          return nullptr ;
       pSci->ToGlob( frStart) ;
+     // unisco
       pSTM->DoSewing( *pSci) ;
      // riferimento alla fine della linea guida
       Frame3d frEnd ;
@@ -740,6 +748,7 @@ GetSurfTriMeshSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEnds, d
          return nullptr ;
       pSce->Invert() ;
       pSce->ToGlob( frEnd) ;
+     // unisco
       pSTM->DoSewing( *pSce) ;
    }
   // se superficie risultante chiusa, verifico che la normale sia verso l'esterno
@@ -750,6 +759,244 @@ GetSurfTriMeshSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEnds, d
    return Release( pSTM) ;
 }
 
+//-------------------------------------------------------------------------------
+static ISurfTriMesh*
+GetSurfTriMeshSwept3d( const ICurve* pSect, const ICurve* pGuide, const Vector3d& vtAx, bool bCapEnds, double dLinTol)
+{
+  // determino se la sezione è chiusa
+   bool bSectClosed = pSect->IsClosed() ;
+  // determino se la guida è chiusa
+   bool bGuideClosed = pGuide->IsClosed() ;
+  // determino algoritmo da usare per calcolare i riferimenti lungo la curva
+   bool bRMF = vtAx.IsSmall() ;
+
+  // riferimento all'inizio della linea guida
+   Point3d ptStart ;
+   pGuide->GetStartPoint( ptStart) ;
+   Vector3d vtStart ;
+   pGuide->GetStartDir( vtStart) ;
+   Frame3d frStart ;
+   if ( bRMF) {
+      if ( ! frStart.Set( ptStart, vtStart))
+         return nullptr ;
+   }
+   else {
+      Vector3d vtAxX = vtAx ^ vtStart ;
+      if ( vtAxX.IsSmall()) {
+         vtAxX = FromUprightOrtho( vtAx) ;
+         vtAxX.Rotate( vtAx, 0, 1) ;
+      }
+      if ( ! frStart.Set( ptStart, vtStart, vtAxX))
+         return nullptr ;
+   }
+
+  // calcolo la polilinea che approssima la sezione
+   PolyLine PL ;
+   if ( ! pSect->ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, PL))
+      return nullptr ;
+
+  // recupero la sezione dalla PolyLine approssimata come Curva
+   PtrOwner<CurveComposite> pSecLocApprox( CreateBasicCurveComposite()) ;
+   if ( IsNull( pSecLocApprox)  ||
+        ! pSecLocApprox->FromPolyLine( PL)  ||
+        ! pSecLocApprox->IsValid())
+      return nullptr ;
+
+  // porto la PolyLine e la curva della sezione nel riferimento nel punto iniziale della guida
+   PL.ToLoc( frStart) ;
+   pSecLocApprox->ToLoc( frStart) ;
+
+  // calcolo il vettore di Frames campionati lungo la guida mediante la tolleranza definita
+   FRAME3DVECTOR vFrames ;
+   if ( bRMF) {
+      RotationMinimizingFrame RMF ;
+      if ( ! RMF.Set( pGuide, frStart)  ||
+           ! RMF.GetFramesByTolerance( dLinTol, vFrames)  ||  vFrames.empty())
+         return nullptr ;
+   }
+   else {
+      RotationXplaneFrame RXF ;
+      if ( ! RXF.Set( pGuide, vtAx, frStart.VersX())  ||
+           ! RXF.GetFramesByTolerance( dLinTol, vFrames)  ||  vFrames.empty())
+         return nullptr ;
+   }
+
+  // preparo collettore delle superfici componenti
+   StmFromTriangleSoup StmFts ;
+   if ( ! StmFts.Start())
+      return nullptr ;
+
+  // per ogni Frame calcolato, la sezione va roto-traslata lungo la guida
+   for ( int i = 0 ; i < int( vFrames.size()) - 1 ; ++ i) {
+
+     // definisco la sezione allo step corrente
+      PtrOwner<ICurve> pSecCurr( pSecLocApprox->Clone()) ;
+      if ( IsNull( pSecCurr)  ||  ! pSecCurr->IsValid())
+         return nullptr ;
+
+     // considero la sezione ( in locale ) come vista dal globale
+      if ( ! pSecCurr->ToGlob( vFrames[i]))
+         return nullptr ;
+
+     // definisco la sezione allo step successivo
+      PtrOwner<ICurve> pSecSucc( pSecLocApprox->Clone()) ;
+      if ( IsNull( pSecSucc)  ||  ! pSecSucc->IsValid())
+         return nullptr ;
+
+     // considero la sezione ( in locale ) come vista dal globale
+      if ( ! pSecSucc->ToGlob( vFrames[i+1]))
+         return nullptr ;
+
+     // creo la rigata tra queste due sezioni
+      PtrOwner<ISurfTriMesh> pSr( GetSurfTriMeshRuled( pSecCurr, pSecSucc, ISurfTriMesh::RLT_ISOPAR_SMOOTH, dLinTol)) ;
+      if ( IsNull( pSr)  ||  ! pSr->IsValid())
+         return nullptr ;
+     // la inverto
+      pSr->Invert() ;
+     // inserisco la superficie nel collettore
+      StmFts.AddSurfTriMesh( *pSr) ;
+   }
+
+  // se richiesti caps e sezione chiusa e guida aperta
+   if ( bCapEnds  &&  bSectClosed  &&  ! bGuideClosed) {
+
+     // aggiungo il cap sull'inizio ( portandolo nel frame del punto iniziale della guida )
+      PtrOwner<SurfTriMesh> pSci( CreateBasicSurfTriMesh()) ;
+      if ( IsNull( pSci)  ||  ! pSci->CreateByFlatContour( PL))
+         return nullptr ;
+      pSci->ToGlob( vFrames.front()) ;
+     // aggiungo
+      StmFts.AddSurfTriMesh( *pSci) ;
+
+     // aggiungo il cap sulla fine ( portandolo nel frame del punto finale della guida )
+      PtrOwner<SurfTriMesh> pSce( CreateBasicSurfTriMesh()) ;
+      if ( IsNull( pSce)  ||  ! pSce->CreateByFlatContour( PL))
+         return nullptr ;
+      pSce->ToGlob( vFrames.back()) ;
+     // inverto
+      pSce->Invert() ;
+     // aggiungo
+      StmFts.AddSurfTriMesh( *pSce) ;
+   }
+
+  // recupero la supeficie risultante
+   if ( ! StmFts.End())
+      return nullptr ;
+   PtrOwner<SurfTriMesh> pSTM( GetBasicSurfTriMesh( StmFts.GetSurf())) ;
+   if ( IsNull( pSTM))
+      return nullptr ;
+  // salvo tolleranza lineare usata
+   pSTM->SetLinearTolerance( dLinTol) ;
+
+  // se superficie risultante chiusa, verifico che la normale sia verso l'esterno
+   double dVol ;
+   if ( pSTM->GetVolume( dVol)  &&  dVol < 0)
+      pSTM->Invert() ;
+
+  // restituisco la superficie
+   return Release( pSTM) ;
+}
+
+//-------------------------------------------------------------------------------
+ISurfTriMesh*
+GetSurfTriMeshSwept( const ICurve* pSect, const ICurve* pGuide, const Vector3d& vtAx,
+                     bool bCapEnds, double dLinTol)
+{
+  // verifica parametri
+   if ( pSect == nullptr  ||  pGuide == nullptr)
+      return nullptr ;
+     
+   bool bIsLine = false ;
+   if ( pGuide->GetType() == CRV_LINE)
+      bIsLine = true ;
+   else {
+      const CurveComposite* pCompo = GetBasicCurveComposite( pGuide) ;
+      Point3d ptStart, ptEnd ;
+      if ( pCompo != nullptr  &&  pCompo->IsALine( 10 * EPS_SMALL, ptStart, ptEnd))
+         bIsLine = true ;
+   }
+  // se la guida è piana
+   Plane3d plGuide ;
+   if ( pGuide->IsFlat( plGuide, bIsLine, 10 * EPS_SMALL))
+      return GetSurfTriMeshSweptInPlane( pSect, pGuide, plGuide.GetVersN(), bCapEnds, dLinTol) ;
+
+  // altrimenti swept 3d
+   return GetSurfTriMeshSwept3d( pSect, pGuide, vtAx, bCapEnds, dLinTol) ;
+}
+
+//-------------------------------------------------------------------------------
+ISurfTriMesh*
+GetSurfTriMeshSwept( const ISurfFlatRegion* pSfrSect, const ICurve* pGuide, const Vector3d& vtAx,
+                     bool bCapEnds, double dLinTol)
+{
+  // verifica dei parametri
+   if ( pSfrSect == nullptr  ||  pGuide == nullptr)
+      return nullptr ;
+
+  // predispongo collettore superfici componenti
+   StmFromTriangleSoup StmSoup ;
+   StmSoup.Start() ;
+
+  // per ogni loop della superficie, creo una Swept
+   for ( int nC = 0 ; nC < pSfrSect->GetChunkCount() ; ++ nC) {
+      for ( int nL = 0 ; nL < pSfrSect->GetLoopCount( nC) ; ++ nL) {
+        // recupero il loop
+         PtrOwner<ICurve> pCrvLoop( pSfrSect->GetLoop( nC, nL)) ;
+         if ( IsNull( pCrvLoop)  ||  ! pCrvLoop->IsValid())
+            return nullptr ;
+        // creo la Trimesh Swept
+         PtrOwner<ISurfTriMesh> pStmLoopSwept( GetSurfTriMeshSwept( pCrvLoop, pGuide, vtAx, false, dLinTol)) ;
+         if ( IsNull( pStmLoopSwept)  ||  ! pStmLoopSwept->IsValid())
+            return nullptr ;
+        // aggiungo la Swept ricavata al risultato finale ( come triangoli )
+         StmSoup.AddSurfTriMesh( *pStmLoopSwept) ;
+      }
+   }
+
+  // Recupero la superficie
+   if ( ! StmSoup.End())
+      return nullptr ;
+   PtrOwner<ISurfTriMesh> pStmSwept( StmSoup.GetSurf()) ;
+   if ( IsNull( pStmSwept))
+      return nullptr ;
+
+  // se rischiesta chiusura...
+  // Controllo solo che la guida non sia chiusa, la sezione derivando da una Flatregion è sempre chiusa
+   if ( bCapEnds  &&  ! pGuide->IsClosed()) {
+     // recupero i loop all'inizio (dalla regione e apportunamente approssimati)
+      POLYLINEVECTOR vPLi ;
+      for ( int nC = 0 ; nC < pSfrSect->GetChunkCount() ; ++ nC) {
+         for ( int nL = 0 ; nL < pSfrSect->GetLoopCount( nC) ; ++ nL) {
+            vPLi.emplace_back() ;
+            if ( ! pSfrSect->ApproxLoopWithLines( nC, nL, dLinTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, vPLi.back()))
+               return nullptr ;
+         }
+      }
+     // creo il cap sull'inizio e lo attacco alla swept ( è già in posizione giusta)
+      PtrOwner<ISurfTriMesh> pSci( CreateSurfTriMesh()) ;
+      if ( ! pSci->CreateByRegion( vPLi))
+         return nullptr ;
+      pStmSwept->DoSewing( *pSci) ;
+     // recupero i loops alla fine
+      POLYLINEVECTOR vPLe ;
+      if ( ! pStmSwept->GetLoops( vPLe))
+         return nullptr ;
+     // creo la superficie alla fine e la attacco
+      PtrOwner<ISurfTriMesh> pSce( CreateSurfTriMesh()) ;
+      if ( ! pSce->CreateByRegion( vPLe))
+         return nullptr ;
+     // attacco la superficie finale alla swept
+      pSce->Invert() ;
+      pStmSwept->DoSewing( *pSce) ;
+   }
+  // se superficie risultante chiusa, verifico che la normale sia verso l'esterno
+   double dVol ;
+   if ( pStmSwept->GetVolume( dVol)  &&  dVol < 0)
+      pStmSwept->Invert() ;
+
+   return Release( pStmSwept) ;
+}
+
 //-------------------------------------------------------------------------------
 ISurfTriMesh*
 GetSurfTriMeshTransSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEnds, double dLinTol)
@@ -757,7 +1004,7 @@ GetSurfTriMeshTransSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEn
   // verifica parametri
    if ( pSect == nullptr  ||  pGuide == nullptr)
       return nullptr ;
-  // determino se la sezione è chiusa
+  // determino se la sezione è chiusa
    bool bSectClosed = pSect->IsClosed() ;
   // punto iniziale della sezione e vettore a inizio guida
    Point3d ptStart ;
@@ -771,7 +1018,7 @@ GetSurfTriMeshTransSwept( const ICurve* pSect, const ICurve* pGuide, bool bCapEn
    PolyLine PLG ;
    if ( ! pGuide->ApproxWithLines( dLinTol, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL, PLG))
       return nullptr ;
-  // determino se la guida è chiusa
+  // determino se la guida è chiusa
    bool bGuideClosed = PLG.IsClosed() ;
   // calcolo la superficie
    PtrOwner<SurfTriMesh> pSTM( CreateBasicSurfTriMesh()) ;
@@ -937,7 +1184,7 @@ CalcRegionPolyLines( const CICURVEPVECTOR& vpCurve, double dLinTol,
       vCrvCompo[i]->ToLoc( frRef) ;
    }
 
-  // creo una matrice di interi ; ogni riga corrisponde ad un chunk, dove in posizione 0 c'è il loop esterno e nelle
+  // creo una matrice di interi ; ogni riga corrisponde ad un chunk, dove in posizione 0 c'è il loop esterno e nelle
   // successive i loop interni
    INTMATRIX vnPLIndMat ;
 
@@ -970,10 +1217,10 @@ CalcRegionPolyLines( const CICURVEPVECTOR& vpCurve, double dLinTol,
             }
             bFirstCrv = false ;
          }
-        // ... altrimenti verifico se il loop è interno o no
+        // ... altrimenti verifico se il loop è interno o no
          else {
-           // il loop è interno se è sia interno al loop esterno della riga di vnPLIndMat e allo stesso tempo
-           // esterno a tutti i loop già inseriti nella riga attuale.
+           // il loop è interno se è sia interno al loop esterno della riga di vnPLIndMat e allo stesso tempo
+           // esterno a tutti i loop già inseriti nella riga attuale.
            // verifica rispetto loop esterno
             IntersCurveCurve ccInt( *vCrvCompo[vnPLIndMat.back().front()], *vCrvCompo[j]) ;
             CRVCVECTOR ccClass ;

SurfAux.cpp

@@ -20,12 +20,14 @@
 #include "CurveBezier.h"
 #include "CurveComposite.h"
 #include "SurfBezier.h"
-#include "/EgtDev/Include/EgtPointerOwner.h"
 #include "/EgtDev/Include/EGkSurf.h"
 #include "/EgtDev/Include/EGkSurfAux.h"
+#include "/EgtDev/Include/EGkSfrCreate.h"
+#include "/EgtDev/Include/EgtPointerOwner.h"
 
 using namespace std ;
 
+//----------------------------------------------------------------------------
 bool
 NurbsSurfaceCanonicalize( SNurbsSurfData& snData)
 {
@@ -108,26 +110,24 @@ NurbsSurfaceCanonicalize( SNurbsSurfData& snData)
    return true;
 }
 
-
 //----------------------------------------------------------------------------
 ISurf*
 NurbsToBezierSurface(const SNurbsSurfData& snData)
 {
-   // la superficie Nurbs deve essere in forma canonica
+  // la superficie Nurbs deve essere in forma canonica
    if ( snData.bPeriodicU || snData.bPeriodicV || snData.bExtraKnotes )
       return nullptr ;
-   // controllo sul numero dei nodi
+  // controllo sul numero dei nodi
    int nU = snData.nCPU + snData.nDegU - 1 ;
    int nV = snData.nCPV + snData.nDegV - 1 ;
-   // controllo nodi e punti di controllo
-   //if ( nU != int( snData.vU.size())  ||  nV != int( snData.vV.size())  ||  snData.nCPU * snData.nCPV != int( snData.vCP.size()))
+  // controllo nodi e punti di controllo
    if ( nU != int(snData.vU.size()) || nV != int(snData.vV.size())) {
       return nullptr ;
    }
 
-   // verifico le condizioni agli estremi sui nodi (i primi nDeg nodi e gli ultimi nDeg nodi devono essere uguali tra loro)
+  // verifico le condizioni agli estremi sui nodi (i primi nDeg nodi e gli ultimi nDeg nodi devono essere uguali tra loro)
    bool bOk = true ;
-   // direzione U
+  // direzione U
    for ( int i = 1 ; i < snData.nDegU ; ++ i) {
       if ( abs( snData.vU[i] - snData.vU[0]) >= EPS_ZERO)
          bOk = false ;
@@ -136,7 +136,7 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
       if ( abs( snData.vU[nU - 1 - i] - snData.vU[nU - 1]) >= EPS_ZERO)
          bOk = false ;
    }
-   // direzione V
+  // direzione V
    for ( int i = 1 ; i < snData.nDegV ; ++ i) {
       if ( abs( snData.vV[i] - snData.vV[0]) >= EPS_ZERO)
          bOk = false ;
@@ -148,8 +148,8 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    if ( ! bOk)
       return nullptr ;
 
-   // algoritmo 5.7 del libro "The NURBS book"//////////////////////////////////////////////////////////////////////////////////////////////////////////////
-   // creazione delle strips nella direzione U ( trasformo le curve iso con U costante in bezier)
+  // algoritmo 5.7 del libro "The NURBS book"
+  // creazione delle strips nella direzione U ( trasformo le curve iso con U costante in bezier)
    int a = snData.nDegU - 1 ;
    int b = snData.nDegU ;
    int nb = 0 ;  // numero di strisce in U ( lunghezza con U costante)
@@ -164,14 +164,14 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    DBLVECTOR vAlpha ;
    vAlpha.resize( snData.nDegU - 1) ;
    if ( ! snData.bRat ) {
-      for ( int i = 0 ; i <= snData.nDegU ; ++i ) {
+      for ( int i = 0 ; i <= snData.nDegU ; ++i) {
          for ( int row = 0 ; row < snData.nCPV ; ++row ) {
             mBC[i][row] = snData.mCP[i][row] ;
          }
       }
    }
    else {
-      for ( int i = 0 ; i <= snData.nDegU ; ++i ) {
+      for ( int i = 0 ; i <= snData.nDegU ; ++i) {
          for (int row = 0 ; row < snData.nCPV ; ++ row) {
             mW[i][row] = snData.mW[i][row] ;
             mBC[i][row] = snData.mCP[i][row] * snData.mW[i][row] ;
@@ -185,25 +185,25 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
       while ( b < nU - 1  &&  abs( snData.vU[b+1] - snData.vU[b]) < EPS_ZERO)
          ++ b ;
       int mult = b - i + 1 ;
-      if ( mult < snData.nDegU ) {
+      if ( mult < snData.nDegU) {
          bRef = true ;
          // calcolo numeratore e alpha
          double numer = snData.vU[b] - snData.vU[a] ;
          for ( int j = snData.nDegU ; j > mult ; -- j)
             vAlpha[j-mult-1] = numer / ( snData.vU[a+j] - snData.vU[a]) ;
          int r = snData.nDegU - mult ;
-         for ( int j = 1 ; j <= snData.nDegU - mult ; ++j ) {
+         for ( int j = 1 ; j <= snData.nDegU - mult ; ++j) {
             int save = r - j ;
             int s = mult + j ;
             if ( ! snData.bRat ) {
-               for ( int k = snData.nDegU ; k >= s ; --k ) {
+               for ( int k = snData.nDegU ; k >= s ; --k) {
                   for ( int row = 0 ; row < snData.nCPV ; ++row) {
                      mBC[k][row] = vAlpha[k-s] * mBC[k][row] + ( 1 -  vAlpha[k-s]) * mBC[k-1][row] ;
                   }
                }
             }
             else {
-               for ( int k = snData.nDegU ; k >= s ; --k ) {
+               for ( int k = snData.nDegU ; k >= s ; --k) {
                   for ( int row = 0 ; row < snData.nCPV ; ++row) {
                      mBC[k][row] = vAlpha[k-s] * mBC[k][row] + ( 1 -  vAlpha[k-s]) * mBC[k-1][row] ;
                      mW[k][row] = vAlpha[k-s] * mW[k][row] + ( 1 - vAlpha[k-s]) * mW[k-1][row] ;
@@ -222,26 +222,26 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
             }
          }
       }
-      mPC_strip.resize( snData.nDegU * ( nb  + 1) + 1 , vCPV) ;
+      mPC_strip.resize( snData.nDegU * ( nb  + 1) + 1, vCPV) ;
       mW_strip.resize( snData.nDegU * ( nb  + 1) + 1, vV_W) ;
       if ( ! snData.bRat)
          for ( int i = 0 ; i <= snData.nDegU ; ++i) {
-            for ( int row = 0 ; row < snData.nCPV ; ++row ) {
+            for ( int row = 0 ; row < snData.nCPV ; ++row) {
                mPC_strip[i+ nb * snData.nDegU][row] = mBC[i][row] ;
             }
          }
       else {
          for ( int i = 0 ; i <= snData.nDegU ; ++i) {
-            for ( int row = 0 ; row < snData.nCPV ; ++row ) {
+            for ( int row = 0 ; row < snData.nCPV ; ++row) {
                mPC_strip[i+ nb * snData.nDegU][row] = mBC[i][row]/mW[i][row] ;
                mW_strip[i+ nb * snData.nDegU][row] = mW[i][row] ;
             }
          }
       }
+     // ho finito di definire la patch di Bezier attuale e passo alla successiva
       ++ nb ;
-      // ho finito di definire la patch di Bezier attuale e passo alla successiva
 
-      // aggiorno mBC con i valori della prossima pezza di Bezier // corrisponde a nb = nb + 1 
+     // aggiorno mBC con i valori della prossima pezza di Bezier // corrisponde a nb = nb + 1 
       if ( ! snData.bRat){
          for (int i = 0 ; i < snData.nDegU - 1 ; ++ i) {
             for ( int row = 0 ; row < snData.nCPV ; ++row) {
@@ -283,7 +283,7 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
       }
    }
 
-   // se non ho raffinato allora tutti i nodi avevano gi? molteplicit? massima. Converto direttamente in Bezier la dir U
+  // se non ho raffinato allora tutti i nodi avevano gi? molteplicit? massima. Converto direttamente in Bezier la dir U
    int nCPU_ref ; // numero dei punti di controllo in U dopo il raffinamento
    if ( ! bRef ) {
       nCPU_ref = snData.nCPU ;
@@ -308,8 +308,8 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    else
       nCPU_ref = snData.nDegU * nb + 1 ; // numero dei punti di controllo in U dopo il raffinamento
 
-   // ora ho ottenuto le strisce nDegU x nCPV
-   // devo ripetere la procedura, sulla dir V, per ottenere le patch nDegU x nDegV
+  // ora ho ottenuto le strisce nDegU x nCPV
+  // devo ripetere la procedura, sulla dir V, per ottenere le patch nDegU x nDegV
    a = snData.nDegV - 1 ;
    b = snData.nDegV ;
    int nc = 0 ;  // numero di strisce in V ( lunghezza con V costante)
@@ -325,14 +325,14 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    vector<vector<Point3d>> mPC_tot( nCPU_ref, vDegV) ;
    vector<DBLVECTOR> mW_tot( nCPU_ref, vV1_W) ;
    if ( ! snData.bRat ) {
-      for ( int i = 0 ; i < nCPU_ref ; ++i ) {
+      for ( int i = 0 ; i < nCPU_ref ; ++i) {
          for ( int row = 0 ; row <= snData.nDegV ; ++row ) {
             m_BC1[i][row] = mPC_strip[i][row] ;
          }
       }
    }
    else {
-      for ( int i = 0 ; i < nCPU_ref ; ++i ) {
+      for ( int i = 0 ; i < nCPU_ref ; ++i) {
          for (int row = 0 ; row <= snData.nDegV ; ++ row) {
             mW1[i][row] = mW_strip[i][row] ;
             m_BC1[i][row] = mPC_strip[i][row] * mW_strip[i][row] ;
@@ -348,24 +348,24 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
       int mult = b - i + 1 ;
       if ( mult < snData.nDegV ) {
          bRef = true ;
-         // calcolo numeratore e alpha
+        // calcolo numeratore e alpha
          double numer = snData.vV[b] - snData.vV[a] ;
          for ( int j = snData.nDegV ; j > mult ; -- j)
             vAlpha1[j-mult-1] = numer / ( snData.vV[a+j] - snData.vV[a]) ;
          int r = snData.nDegV - mult ;
-         for ( int j = 1 ; j <= snData.nDegV - mult ; ++j ) {
+         for ( int j = 1 ; j <= snData.nDegV - mult ; ++j) {
             int save = r - j ;
             int s = mult + j ;
             if ( ! snData.bRat) {
                for ( int k = 0 ; k < nCPU_ref ; ++k) {
-                  for ( int row = snData.nDegV ; row >= s ; --row ) {
+                  for ( int row = snData.nDegV ; row >= s ; --row) {
                      m_BC1[k][row] = vAlpha1[row-s] * m_BC1[k][row] + ( 1 -  vAlpha1[row-s]) * m_BC1[k][row-1] ;
                   }
                }
             }
             else {
                for ( int k = 0 ; k < nCPU_ref ; ++k) {
-                  for ( int row = snData.nDegV ; row >= s ; --row ) {
+                  for ( int row = snData.nDegV ; row >= s ; --row) {
                      m_BC1[k][row] = vAlpha1[row-s] * m_BC1[k][row] + ( 1 -  vAlpha1[row-s]) * m_BC1[k][row-1] ;
                      mW1[k][row] = vAlpha1[row-s] * mW1[k][row] + ( 1 -  vAlpha1[row-s]) * mW1[k][row-1] ;
                   }
@@ -388,28 +388,28 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
          }
       }
       int nRef = snData.nDegV * ( nc + 1) + 1 ;
-      for ( int k = 0 ; k < nCPU_ref; ++k){
+      for ( int k = 0 ; k < nCPU_ref; ++k) {
          mPC_tot[k].resize( nRef) ;
          mW_tot[k].resize( nRef) ;
       }
       if ( ! snData.bRat)
          for ( int i = 0 ; i < nCPU_ref ; ++i) {
-            for ( int row = 0 ; row <= snData.nDegV ; ++row ) {
+            for ( int row = 0 ; row <= snData.nDegV ; ++row) {
                mPC_tot[i][row + nc * snData.nDegV] = m_BC1[i][row] ;
             }
          }
       else {
          for ( int i = 0 ; i < nCPU_ref ; ++i) {
-            for ( int row = 0 ; row <= snData.nDegV ; ++row ) {
-               mPC_tot[i][row + nc * snData.nDegV] = m_BC1[i][row]/mW1[i][row] ;
+            for ( int row = 0 ; row <= snData.nDegV ; ++row) {
+               mPC_tot[i][row + nc * snData.nDegV] = m_BC1[i][row] / mW1[i][row] ;
                mW_tot[i][row + nc * snData.nDegV] = mW1[i][row] ;
             }
          }
       }
+     // ho finito di definire la patch di Bezier attuale e passo alla successiva
       ++ nc ;
-      // ho finito di definire la patch di Bezier attuale e passo alla successiva
 
-      // aggiorno mBC con i valori della prossima pezza di Bezier // corrisponde a nc = nc + 1 
+     // aggiorno mBC con i valori della prossima pezza di Bezier // corrisponde a nc = nc + 1 
       if ( ! snData.bRat){
          for (int i = 0 ; i < nCPU_ref ; ++ i) {
             for ( int row = 0 ; row < snData.nDegV - 1 ; ++row) {
@@ -450,7 +450,7 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
          ++b ;
       }
    }
-   // se non ho raffinato allora aggiungo direttamente alle matrici della superficie totale
+  // se non ho raffinato allora aggiungo direttamente alle matrici della superficie totale
    int nCPV_ref ; // numero dei punti di controllo in V dopo il raffinamento
    if ( ! bRef) {
       nCPV_ref = snData.nCPV ;
@@ -477,7 +477,7 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    else
       nCPV_ref = snData.nDegV * nc + 1 ;
 
-   // finalmente setto la superficie di bezier totale divisa in nb patch in U e nc patch in V
+  // finalmente setto la superficie di bezier totale divisa in nb patch in U e nc patch in V
    PtrOwner<SurfBezier> pSrfBz( CreateBasicSurfBezier()) ;
    if ( IsNull( pSrfBz))
       return nullptr ;
@@ -498,3 +498,178 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    }
    return Release( pSrfBz) ;
 }
+
+//----------------------------------------------------------------------------
+bool
+MakeUniform( ISurfFlatRegion*& pSfr, bool& bRescaled, const DBLVECTOR& vU0, const DBLVECTOR& vV0,
+             int nDegU, int nDegV, double dScaleU, double dScaleV, bool bRetry) 
+{
+  // la superficie in input arriva già scalata
+   bool bRescaledU = false ;
+   bool bRescaledV = false ;
+   int nSpanU = 1, nSpanV = 1 ;
+   PtrOwner<ISurfFlatRegion> pRescaledSfr( CreateSurfFlatRegion()) ;
+   for ( int nDir = 0 ; nDir <= 1 ; ++ nDir) {
+     // vettore dei nodi
+      DBLVECTOR vU ;
+      int nExtraKnots = 0 ;
+     // controllo in U
+      if ( nDir == 0) {
+         if ( nDegU > 1) {
+            nExtraKnots = nDegU - 1 ;
+         }
+         for ( int i = nExtraKnots ; i < int( vU0.size()) - nExtraKnots ; ++i ) {
+            double dKnot = vU0[i] * SBZ_TREG_COEFF ;
+            // lo aggiungo solo se è diverso dal precedente
+            if ( i == nExtraKnots  ||  dKnot > vU.back() + EPS_SMALL  ||  dKnot < vU.back() - EPS_SMALL)
+               vU.push_back( dKnot) ;
+         }
+         nSpanU = (int)vU.size() - 1 ;
+      }
+     // controllo in V
+      else if ( nDir == 1 ) {
+         if ( nDegV > 1) {
+            nExtraKnots = nDegV - 1 ;
+         }
+         for ( int i = nExtraKnots ; i < int( vV0.size()) - nExtraKnots ; ++i ) {
+            double dKnot = vV0[i] * SBZ_TREG_COEFF ;
+            // lo aggiungo solo se è diverso dal precedente
+            if ( i == nExtraKnots  ||  dKnot > vU.back() + EPS_SMALL  ||  dKnot < vU.back() - EPS_SMALL)
+               vU.push_back( dKnot) ;
+         }
+         nSpanV = (int)vU.size() - 1 ;
+      }
+
+     // controllo se il vettore dei nodi è uniforme
+      int a = 0, b = 1 ;
+      double d0 = abs( vU[b] - vU[a]), d1 = d0 ;
+      // il vettore è uniforme quando la distanza tra nodi consecutivi è sempre zero o un valore costante
+      while ( b < (int)vU.size()  &&  ( ( d1 < d0 + EPS_SMALL  &&  d1 > d0 - EPS_SMALL)  ||  d1 < EPS_SMALL)) {
+         a = b ;
+         ++b ;
+         if ( b < (int)vU.size())
+            d1 = abs( vU[b] - vU[a]) ;
+      }
+      if ( b != (int)vU.size()) {
+         nDir == 0 ? bRescaledU = true : bRescaledV = true ;
+         pRescaledSfr.Set( CreateSurfFlatRegion()) ;
+         if ( IsNull( pRescaledSfr))
+            return false ;
+         for ( int p = 0 ; p < (int)vU.size() - 1 ; ++p) {
+            PtrOwner<ISurfFlatRegion> pSfr_copy( pSfr->Clone()) ;
+            if ( IsNull( pSfr_copy))
+               return false ;
+            double dLenStrip = abs( vU[p+1] - vU[p]) ;
+            if ( dLenStrip < EPS_SMALL)
+               continue ;
+            // creo la maschera per tagliare la superficie originale e ottenere una striscia
+            PtrOwner<ISurfFlatRegion> pSfrTrim( CreateSurfFlatRegion()) ;
+
+            // ricavo la maschera del trim, con cui poi farò l'intersezione con la sfr iniziale
+            Vector3d vtTrim ;
+            if ( nDir == 0) {
+               pSfrTrim.Set( GetSurfFlatRegionRectangle( dLenStrip, dScaleV + 2)) ;
+               vtTrim.Set( abs(vU[p] - vU.front()), - 1, 0) ;
+            }
+            else{
+               pSfrTrim.Set( GetSurfFlatRegionRectangle( dScaleU + 2, dLenStrip)) ;
+               vtTrim.Set( - 1, abs(vU[p] - vU.front()), 0) ;
+            }
+            pSfrTrim->Translate( vtTrim) ;
+
+            if ( ! pSfr_copy->Intersect( *pSfrTrim))
+               return false ;
+
+            // aggiungo la nuova striscia solo se è valida
+            if ( pSfr_copy->IsValid() ) {
+               if ( nDir == 0)
+                  pSfr_copy->Scale( GLOB_FRM, SBZ_TREG_COEFF / dLenStrip, 1, 1) ;
+               else
+                  pSfr_copy->Scale( GLOB_FRM, 1, SBZ_TREG_COEFF / dLenStrip, 1) ;
+
+               // prima di riunire la striscia al resto devo traslarla sul bordo destro della superificie che sto ricostruendo
+
+               Point3d pt ;
+               nDir == 0 ? pt.Set( abs(vU[p] - vU.front()), 0, 0) : pt.Set( 0,abs(vU[p] - vU.front()), 0) ;
+               if ( nDir == 0)
+                  pt.Scale( GLOB_FRM, SBZ_TREG_COEFF / dLenStrip, 1, 1) ;
+               else
+                  pt.Scale( GLOB_FRM, 1, SBZ_TREG_COEFF / dLenStrip, 1) ;
+
+               Vector3d vtJoin ;
+               if ( nDir == 0)
+                  vtJoin.Set( p * SBZ_TREG_COEFF - pt.x, 0, 0) ;
+               else
+                  vtJoin.Set( 0, p  * SBZ_TREG_COEFF - pt.y, 0) ;
+               pSfr_copy->Translate( vtJoin) ;
+               // se sto ritentando MakeUniform, allora faccio anche OFFSET e controOFFSET
+               if ( bRetry)
+                  pSfr_copy->Offset( 10 * EPS_SMALL, ICurve::OFF_FILLET) ; // OFFSET
+               if ( pRescaledSfr->IsValid()) {
+                  if ( ! pRescaledSfr->Add( *pSfr_copy))
+                     return false ;
+               }
+               else
+                  pRescaledSfr.Set( pSfr_copy) ;
+            }
+         }
+         if ( nDir == 0) {
+            dScaleU = ((int)vU.size() - 1) * SBZ_TREG_COEFF ;
+            if ( pRescaledSfr->IsValid()) {
+               if ( bRetry)
+                  pRescaledSfr->Offset( -10 * EPS_SMALL, ICurve::OFF_FILLET) ;  //contro OFFSET
+               pSfr = Release( pRescaledSfr) ;
+            }
+         }
+         else
+            dScaleV = ((int)vU.size() - 1) * SBZ_TREG_COEFF ;
+      }
+   }
+
+   if ( ! IsNull( pRescaledSfr)  &&  pRescaledSfr->IsValid()) {
+      if ( bRetry)
+         pRescaledSfr->Offset( -10 * EPS_SMALL, ICurve::OFF_FILLET) ;  // contro OFFSET
+      pSfr = Release( pRescaledSfr) ;
+   }
+
+   if ( ! bRescaledU  &&  ! bRescaledV)
+      pSfr->Scale( GLOB_FRM, nSpanU / dScaleU * SBZ_TREG_COEFF, nSpanV / dScaleV * SBZ_TREG_COEFF, 1) ;
+   else if ( bRescaledU  &&  ! bRescaledV)
+      pSfr->Scale( GLOB_FRM, 1, nSpanV / dScaleV * SBZ_TREG_COEFF, 1) ;
+   else if ( ! bRescaledU  &&  bRescaledV)
+      pSfr->Scale( GLOB_FRM, nSpanU / dScaleU * SBZ_TREG_COEFF, 1, 1) ;
+
+   bRescaled = ( bRescaledU || bRescaledV) ;
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+OnWhichEdge( double u0, double u1, double v0, double v1, const Point3d& ptToAssign, int& nEdge)
+{
+   Point3d ptTR( u1, v1) ;
+   Point3d ptTl( u0, v1) ;
+   Point3d ptBL( u0, v0) ;
+   Point3d ptBr( u1, v0) ;
+
+   double dEps = 0.1 ;
+   if ( AreSamePointEpsilon( ptToAssign, ptTR, dEps))
+      nEdge = 7 ;
+   else if ( AreSamePointEpsilon( ptToAssign, ptTl, dEps))
+      nEdge = 4 ;
+   else if ( AreSamePointEpsilon( ptToAssign, ptBL, dEps))
+      nEdge = 5 ;
+   else if ( AreSamePointEpsilon( ptToAssign, ptBr, dEps))
+      nEdge = 6 ;
+   else if ( ptToAssign.x > ptBL.x - dEps  &&  ptToAssign.x < ptTR.x + dEps  &&  abs( ptToAssign.y - ptTR.y) < dEps)
+      nEdge = 0 ;
+   else if ( ptToAssign.y > ptBL.y - dEps  &&  ptToAssign.y < ptTR.y + dEps  &&  abs( ptToAssign.x - ptBL.x) < dEps)
+      nEdge = 1 ;
+   else if ( ptToAssign.x > ptBL.x - dEps  &&  ptToAssign.x < ptTR.x + dEps  &&  abs( ptToAssign.y - ptBL.y) < dEps)
+      nEdge = 2 ;
+   else if ( ptToAssign.y > ptBL.y - dEps  &&  ptToAssign.y < ptTR.y + dEps  &&  abs( ptToAssign.x - ptTR.x) < dEps)
+      nEdge = 3 ;
+   else
+      return false ;
+   return true ;
+}

SurfBezier.cpp

@@ -33,12 +33,11 @@
 #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 <limits>
-#include "/EgtDev/Include/EGkGeoObjSave.h"
 #include "/EgtDev/Include/EGkGeoPoint3d.h"
 #include "/EgtDev/Include/EGkIntervals.h"
+#include "/EgtDev/Extern/Eigen/Dense"
+#include <limits>
 
 using namespace std ;
 
@@ -153,13 +152,21 @@ SurfBezier::SetTrimRegion( ISurfFlatRegion& sfrTrimReg, bool bIntersectOrSubtrac
    PtrOwner< ISurfFlatRegion> pSfrTrim( GetSurfFlatRegionRectangle( SBZ_TREG_COEFF * m_nSpanU, SBZ_TREG_COEFF * m_nSpanV)) ;
    // 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 ;
+      if ( IsNull( pSfrTrim)  ||  ! pSfrTrim->Intersect( sfrTrimReg)  ||  ! pSfrTrim->IsValid()) {
+         // provo a offsettare il rettangolo parametrico ( ingrandendolo) per vedere se risolvo problemi di intersezione
+         pSfrTrim->Offset( 10* EPS_SMALL, ICurve::OFF_EXTEND) ;
+         if ( ! 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 ;
+      if ( IsNull( pSfrTrim)  ||  ! pSfrTrim->Subtract( sfrTrimReg)  ||  ! pSfrTrim->IsValid()) {
+         // provo a offsettare il rettangolo parametrico ( ingrandendolo) per vedere se risolvo problemi di sottrazione
+         pSfrTrim->Offset( 10* EPS_SMALL, ICurve::OFF_EXTEND) ;
+         if ( ! pSfrTrim->Subtract( sfrTrimReg)  ||  ! pSfrTrim->IsValid())
+            return false ;
+      }
    }
    ResetAuxSurf() ;
   // assegno la regione di trim
@@ -1487,6 +1494,23 @@ SurfBezier::GetAuxSurf( void) const
   // se già calcolata, la restituisco
    if ( m_pSTM != nullptr)
       return m_pSTM ;
+  // eseguo calcolo
+   m_pSTM = GetApproxSurf( 50 * EPS_SMALL, 100 * EPS_SMALL) ;
+   return m_pSTM ;
+}
+
+//----------------------------------------------------------------------------
+SurfTriMesh*
+SurfBezier::GetApproxSurf( double dTol, double dSideMin) const
+{
+  // la superficie deve essere validata
+   if ( m_nStatus != OK)
+      return nullptr ;
+
+  // se c'è ausiliaria e richiesta con gli stessi parametri, ne restituisco una copia
+   if ( m_pSTM != nullptr  &&
+        abs( dTol - 50 * EPS_SMALL) < EPS_SMALL  &&  abs( dSideMin - 100 * EPS_SMALL) < EPS_SMALL)
+      return m_pSTM->Clone() ;
   
   // costruttore della superficie
    POLYLINEMATRIX vvPL ;
@@ -1497,11 +1521,11 @@ SurfBezier::GetAuxSurf( void) const
    bool bTest = false ; // per debug
    // resetto il vettore degli edge
    m_mCCEdge.clear() ;
-   m_mCCEdge = vector<ICRVCOMPOPVECTOR>(4) ;
+   m_mCCEdge = vector<ICRVCOMPOPOVECTOR>(4) ;
    m_vCCLoop.clear() ;
    for ( int i = 0 ; i < (int) vTrees.size() ; ++ i) {
-      Point3d ptMin = std::get<0>( vTrees[i]) ;
-      Point3d ptMax = std::get<1>( vTrees[i]) ;
+      Point3d ptMin = get<0>( vTrees[i]) ;
+      Point3d ptMax = get<1>( vTrees[i]) ;
       Tree.SetSurf( this, true, ptMin, ptMax) ;
       if ( bTest) {
          Tree.BuildTree_test() ;        // per debug
@@ -1509,17 +1533,19 @@ SurfBezier::GetAuxSurf( void) const
          Tree.SetTestMode() ;
       }
       else {
-         Tree.BuildTree( 5 * LIN_TOL_FINE, 0.1) ;
+         //Tree.BuildTree( 5 * LIN_TOL_FINE, 0.1) ;
+         Tree.BuildTree( dTol, dSideMin) ;
       }
       Tree.GetPolygons( vvPL) ;
       //Tree.GetPolygonsBasic( vPL) ;  // per usare i polygon basic
 
-      // aggiorno la chiusura della superficie
+     // aggiorno la chiusura della superficie
       m_bClosedU = m_bClosedU  ||  Tree.IsClosedU() ;
       m_bClosedV = m_bClosedV  ||  Tree.IsClosedV() ;
-      // salvo i bordi in 3d, che servono in caso si voglia trimmare la superficie DOPO aver costruito la trimesh ausiliaria
+     // salvo i bordi in 3d, che servono in caso si voglia trimmare la superficie DOPO aver costruito la trimesh ausiliaria
       UpdateEdgesFromTree( Tree) ;
    }
+
   //// per usare i polygon basic//////////////////////
    //for (int k = 0 ; k < (int)vPL.size(); ++k) {
    //   vvPL.emplace_back() ;
@@ -1527,9 +1553,8 @@ SurfBezier::GetAuxSurf( void) const
    //}
   //// per usare i polygon basic///////////////////
 
-
-   if ( int(vvPL.size()) == 0)
-      LOG_DBG_ERR( GetEGkLogger(), "ERROR : Bezier Surface couldn't be triangulated, hence wasn't drawn") ;
+   if ( vvPL.empty())
+      LOG_DBG_ERR( GetEGkLogger(), "ERROR : Bezier Surface couldn't be triangulated") ;
    
    StmFromTriangleSoup stmSoup ;
    if ( ! stmSoup.Start())
@@ -1561,12 +1586,11 @@ SurfBezier::GetAuxSurf( void) const
       }
    } 
 
-  // la salvo
+  // termino
    if ( ! stmSoup.End())
       return nullptr ;
-   m_pSTM = GetBasicSurfTriMesh( stmSoup.GetSurf()) ;
-
-   return m_pSTM ;
+  // restituisco
+   return GetBasicSurfTriMesh( stmSoup.GetSurf()) ;
 }
 
 //----------------------------------------------------------------------------
@@ -1837,11 +1861,12 @@ SurfBezier::AddCurveCompoToCuts( ICurveComposite* pCrvCompo, ICRVCOMPOPOVECTOR&
    return true ;
 }
 
+//----------------------------------------------------------------------------
 typedef tuple<int,int,int> TRINT ;
 
 template<>
 struct hash<TRINT> {
-   std::size_t operator()(const TRINT& t) const
+   size_t operator()(const TRINT& t) const
    {
       // Compute individual hash values for first, second and third and combine them using XOR and bit shifting:
       return ((hash<int>()(get<0>(t))) ^ (hash<int>()(get<1>(t)) << 1) >> 1) ^ (hash<int>()(get<2>(t)) << 1) ;
@@ -2047,7 +2072,7 @@ SurfBezier::Cut( const Plane3d& plPlane, bool bSaveOnEq)
 
       // vettore di flag che mi indica quali tagli aperti sono stati aggiunti al nuovo bordo
       BOOLVECTOR vbAdded( vpCCOpen.size()) ;
-      std::fill( vbAdded.begin(), vbAdded.end(), false) ;
+      fill( vbAdded.begin(), vbAdded.end(), false) ;
       PtrOwner<ICurveComposite> pCCNewEdge( CreateCurveComposite()) ;
       PtrOwner<ICurveLine> pCL( CreateCurveLine()) ;
       TRINT tiFirstInters ;
@@ -2480,7 +2505,7 @@ SurfBezier::UnprojectPointFromStm( int nT, const Point3d& ptI, Point3d& ptSP, in
          return false ;
       // trovo quale vertice è sull'edge di polo
       BOOLVECTOR vbOn(3) ;
-      std::fill( vbOn.begin(), vbOn.end(), false) ;
+      fill( vbOn.begin(), vbOn.end(), false) ;
       for ( int p = 0 ; p < 3; ++p ) {
          for ( int c = 0 ; c < 4; ++c) {
             for( int i = 0 ; int( m_mCCEdge[c].size()) ; ++i) {
@@ -3169,4 +3194,4 @@ SurfBezier::GetSingleEdge3D( bool bLineOrBezier, int nEdge) const
    }
    }
    return pCrvCompo ;
-}
+}

SurfBezier.h

@@ -19,7 +19,6 @@
 #include "CurveComposite.h"
 #include "SurfTriMesh.h"
 #include "SurfFlatRegion.h"
-//#include "Tree.h"
 #include "/EgtDev/Include/EGkSurfBezier.h"
 #include "/EgtDev/Include/EGkGeoCollection.h"
 
@@ -114,6 +113,7 @@ class SurfBezier : public ISurfBezier, public IGeoObjRW
       bool GetControlCurveOnU( int nIndV, PolyLine& plCtrlU) const override ;
       bool GetControlCurveOnV( int nIndU, PolyLine& plCtrlV) const override ;
       const SurfTriMesh* GetAuxSurf( void) const override ;
+      SurfTriMesh* GetApproxSurf( double dTol, double dSideMin = 100 * EPS_SMALL) const override ;
       // funzione per ottenere la suddivisione dello spazio parametrico nelle celle utilizzate per la triangolazione.
       bool GetLeaves( std::vector<std::tuple<int, Point3d, Point3d>>& vLeaves) const override ;
       bool GetTriangles2D( std::vector<std::tuple<int,Point3d, Point3d, Point3d>>& vTria2D) const override ;
@@ -212,7 +212,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 vector<ICRVCOMPOPVECTOR> m_mCCEdge ;// vettore dei vettori che contengono le curve compo degli edge della superficie nello spazio 3D
+      mutable vector<ICRVCOMPOPOVECTOR> m_mCCEdge ;// vettore dei vettori che contengono le curve compo degli edge della superficie nello spazio 3D
       mutable ICRVCOMPOPOVECTOR m_vCCLoop ; // vettore dei loop della superficie trimmata
 } ;
 

SurfFlatRegion.cpp

@@ -146,7 +146,7 @@ SurfFlatRegion::AddSimpleExtLoop( ICurve* pCrv)
    else {
      // verifico che il piano della curva coincida con quello XY intrinseco
       plPlane.ToLoc( m_frF) ;
-      if ( ! plPlane.GetVersN().IsZ()  ||  abs( plPlane.GetDist()) > EPS_SMALL)
+      if ( ! plPlane.GetVersN().IsZEpsilon( 2 * EPS_ZERO)  ||  abs( plPlane.GetDist()) > EPS_SMALL)
          return false ;
      // sistemo il senso di rotazione (deve essere CCW -> area > 0)
       if ( ( plPlane.GetVersN().z > 0  &&  dArea < 0)  ||
@@ -288,7 +288,7 @@ SurfFlatRegion::AddSimpleIntLoop( ICurve* pCrv)
    Plane3d plPlane ;
    if ( ! pMyCrv->GetArea( plPlane, dArea))
       return false ;
-   if ( ! plPlane.GetVersN().IsZ()  ||  abs( plPlane.GetDist()) > EPS_SMALL)
+   if ( ! plPlane.GetVersN().IsZEpsilon( 2 * EPS_ZERO)  ||  abs( plPlane.GetDist()) > EPS_SMALL)
       return false ;
   // sistemo il senso di rotazione (deve essere CW -> se N==Z+ area < 0, se N==Z- area > 0)
    if ( ( plPlane.GetVersN().z > 0  &&  dArea > 0)  ||  ( plPlane.GetVersN().z < 0  &&  dArea < 0))

SurfTriMesh.cpp

@@ -3244,7 +3244,6 @@ SurfTriMesh::Scale( const Frame3d& frRef, double dCoeffX, double dCoeffY, double
    bMirror = ( bMirror ? ( dCoeffZ > 0) : ( dCoeffZ < 0)) ; 
 
   // aggiorno le facce
-   bool bRecalc = ( abs( dCoeffX) < EPS_ZERO  ||  abs( dCoeffY) < EPS_ZERO  ||  abs( dCoeffZ) < EPS_ZERO) ;
    for ( int i = 0 ; i < GetTriangleSize() ; ++ i) {
       if ( m_vTria[i].nIdVert[0] != SVT_DEL) {
         // se c'è mirror, devo invertire la faccia
@@ -3254,15 +3253,12 @@ SurfTriMesh::Scale( const Frame3d& frRef, double dCoeffX, double dCoeffY, double
          if ( ! CalcTriangleNormal( i)) {
            // elimino il triangolo
             RemoveTriangle( i) ;
-           // richiedo ricalcolo vertici e topologia
-            bRecalc = true ;
          }
       }
    }
-   if ( bRecalc)
-      return DoCompacting() ;
 
-   return true ;
+   return DoCompacting() ;
+
 }
 
 //----------------------------------------------------------------------------

SurfTriMeshUtilities.cpp

@@ -2,7 +2,7 @@
 //  EgalTech 2021-2021
 //----------------------------------------------------------------------------
 // File : SurfTriMeshUtilities.cpp     Data : 01.11.21       Versione : 2.3k1
-// Contenuto : Implementazione funzioni di utilità di Superfici TriMesh.
+// Contenuto : Implementazione funzioni di utilità di Superfici TriMesh.
 //
 //
 //
@@ -123,7 +123,7 @@ SurfTriMesh::RemoveTJunctions( bool& bModified)
 
   // Ciclo sui triangoli della superficie per determinare gli altri vertici sul loro perimetro
    for ( int nT = 0 ; nT < int( m_vTria.size()) ; ++ nT) {
-     // Se il triangolo non è valido, passo al successivo
+     // Se il triangolo non è valido, passo al successivo
       Triangle3d trTria ;
       if ( ! GetTriangle( nT, trTria)  ||  ! trTria.Validate( true))
          continue ;
@@ -138,7 +138,7 @@ SurfTriMesh::RemoveTJunctions( bool& bModified)
       for ( int nSeg = 0 ; nSeg < 3 ; ++ nSeg) {
         // aggiungo al vettore il vertice iniziale del lato
          vVertOtl.emplace_back( m_vTria[nT].nIdVert[nSeg]) ;
-        // Se in questo lato il triangolo è adiacente a un altro, lo salto.
+        // Se in questo lato il triangolo è adiacente a un altro, lo salto.
          if ( m_vTria[nT].nIdAdjac[nSeg] != SVT_DEL  &&  m_vTria[nT].nIdAdjac[nSeg] != SVT_NULL)
             continue ;
          int nPrevSize = int( vVertOtl.size()) ;
@@ -152,7 +152,7 @@ SurfTriMesh::RemoveTJunctions( bool& bModified)
          vtSeg /= dSegLen ;
         // Ciclo sui triangoli vicini
          for ( int nI = 0 ; nI < int( vNearTria.size()) ; ++ nI) {
-           // Salto il triangolo se è quello di riferimento
+           // Salto il triangolo se è quello di riferimento
             if ( vNearTria[nI] == nT)
                continue ;
            // Cerco i vertici che stanno sul lato del triangolo
@@ -175,11 +175,11 @@ SurfTriMesh::RemoveTJunctions( bool& bModified)
                                  } ;
          sort( vVertOtl.begin() + nPrevSize, vVertOtl.end(), SortVertices) ;
       }
-     // Se ci sono più di 3 vertici
+     // Se ci sono più di 3 vertici
       if ( vVertOtl.size() > 3) {
         // Elimino i vertici ripetuti
          vVertOtl.erase( unique( vVertOtl.begin(), vVertOtl.end()), vVertOtl.end()) ;
-        // Se ci sono ancora più di 3 vertici, inserisco nel Map
+        // Se ci sono ancora più di 3 vertici, inserisco nel Map
          if ( vVertOtl.size() > 3)
             TriaMap.emplace( nT, vVertOtl) ;
       }
@@ -191,7 +191,7 @@ SurfTriMesh::RemoveTJunctions( bool& bModified)
       int nT = itT->first ;
      // Vettore degli altri vertici sul perimetro
       const INTVECTOR& vVertOtl = itT->second ;
-     // Se il triangolo non è valido, passo al successivo
+     // Se il triangolo non è valido, passo al successivo
       Triangle3d trTria ;
       if ( ! GetTriangle( nT, trTria)  ||  ! trTria.Validate( true))
          continue ;
@@ -207,7 +207,7 @@ SurfTriMesh::RemoveTJunctions( bool& bModified)
       int nLastNewTria = SVT_NULL ;
      // Se ci sono 4 vertici, inserisco due triangoli
       if ( vVertOtl.size() == 4) {
-        // se 1-2-3 è triangolo (e quindi 0-1-3)
+        // se 1-2-3 è triangolo (e quindi 0-1-3)
          int nNew1Id[3] = { vVertOtl[1], vVertOtl[2], vVertOtl[3]} ;
          int nNew1Tria = AddTriangle( nNew1Id, nTFlag) ;
          if ( nNew1Tria != SVT_NULL  &&  nNew1Tria != SVT_DEL) {
@@ -285,13 +285,13 @@ IsVertex( PNTULIST& PointList, PNTULIST::const_iterator itCurr)
 static bool
 ChooseGoodStartPoint( PNTULIST& PointList)
 {
-  // se il punto iniziale è un vertice, non devo fare alcunché
+  // se il punto iniziale è un vertice, non devo fare alcunché
    if ( IsVertex( PointList, PointList.begin()))
       return true ;
-  // altrimenti cerco il vertice più vicino
+  // altrimenti cerco il vertice più vicino
    for ( auto it = next( PointList.begin()) ; it != PointList.end() ; ++it) {
       if ( IsVertex( PointList, it)) {
-        // se ultimo punto non devo fare alcunché
+        // se ultimo punto non devo fare alcunché
          if ( next( it) == PointList.end())
             return false ;
         // cancello ultimo punto ( coincide con primo)
@@ -317,11 +317,32 @@ AdjustLoop( PNTULIST& PointList, double dMaxEdgeLen, bool& bModif)
    
      // Se dal punto corrente inizia un segmento adiacente a un'altra faccia
       if ( itLast->second != it->second) {
-        // Elimino i punti interni
-         auto itNextToLast = next( itLast) ;
-         for ( auto itInn = itNextToLast ; itInn != it ; ) {
-            itInn = PointList.erase( itInn) ;
-            bModif = true ;
+        // Raccolgo i punti in una polyline
+         PolyLine PL ;
+         int nPar = -1 ;
+         for ( auto itInn = itLast ; itInn != it ; ) {
+            PL.AddUPoint( ++nPar, itInn->first) ;
+            itInn = next( itInn) ;
+         }
+         PL.AddUPoint( ++nPar, it->first) ;
+        // Provo ad eliminare i punti allineati
+         PL.RemoveAlignedPoints( 50 * EPS_SMALL) ;
+         if ( PL.GetPointNbr() < nPar + 1) {
+           // rimuovo dalla lista dei punti gli eliminati (salto gli estremi)
+            int nUCurr = 1 ; 
+            for ( auto itInn = next( itLast) ; itInn != it ; ) {
+               bool bFound = false ;
+               double dU ;
+               for ( bool bOk = PL.GetFirstU( dU) ; bOk && ! bFound ; bOk = PL.GetNextU( dU))
+                  bFound = abs( dU - nUCurr) < EPS_SMALL ;
+               if ( ! bFound) {
+                  itInn = PointList.erase( itInn) ;
+                  bModif = true ;
+               }
+               else
+                  itInn = next( itInn) ;
+               ++ nUCurr ;
+            }
          }
         // Se la lunghezza del segmento supera il limite imposto
          double dSegLen = Dist( it->first, itLast->first) ;
@@ -354,7 +375,7 @@ AdjustLoop( PNTULIST& PointList, double dMaxEdgeLen, bool& bModif)
          }
         // Se i punti sono allineati
          if ( bAreAligned) {
-           // Verifico se il successivo punto non è più allineato
+           // Verifico se il successivo punto non è più allineato
             auto itNextToCurr = next( it) ;
             if ( itNextToCurr != PointList.end()) {
                for ( auto itInn = itNextToLast ; itInn != itNextToCurr  &&  bAreAligned ; ++ itInn) {
@@ -402,7 +423,7 @@ SurfTriMesh::SimplifyFacets( double dMaxEdgeLen, bool bForced)
   // La trimesh deve essere valida
    if ( ! IsValid())
       return false ;
-  // Se la lunghezza massima del lato del triangolo sul bordo della faccia è nulla e non forzata, non devo fare alcunché
+  // Se la lunghezza massima del lato del triangolo sul bordo della faccia è nulla e non forzata, non devo fare alcunché
    if ( dMaxEdgeLen < EPS_SMALL  &&  ! bForced)
       return true ;
 
@@ -505,7 +526,7 @@ SurfTriMesh::SimplifyFacets( double dMaxEdgeLen, bool bForced)
       }
    }
    
-  // dichiaro necessità ricalcolo della grafica e di hashgrids3d
+  // dichiaro necessità ricalcolo della grafica e di hashgrids3d
    m_OGrMgr.Reset() ;
    ResetHashGrids3d() ;
 
@@ -517,13 +538,13 @@ SurfTriMesh::SimplifyFacets( double dMaxEdgeLen, bool bForced)
 bool
 SurfTriMesh::AddChainToChain( const Chain& ChainToAdd, PNTVECTOR& OrigChain)
 {
-  // Se la catena da aggiungere è vuota, non devo fare alcunchè
+  // Se la catena da aggiungere è vuota, non devo fare alcunchè
    if ( ChainToAdd.size() == 0)
       return true ;
-  // Se la catena originale è vuota, non è possibile aggiungere nulla
+  // Se la catena originale è vuota, non è possibile aggiungere nulla
    if ( OrigChain.size() == 0)
       return false ;
-  // Se la catena originale è chiusa non posso aggiungere nulla
+  // Se la catena originale è chiusa non posso aggiungere nulla
    int nLastOrig = max( int( OrigChain.size()) - 1, 0) ;
    if ( AreSamePointApprox( OrigChain[0], OrigChain[nLastOrig]))
       return false ;
@@ -562,11 +583,11 @@ SurfTriMesh::DistPointFacet( const Point3d& ptP, const POLYLINEVECTOR& vPolyVec,
       return false ;
    double dDistPtPl = DistPointPlane( ptP, plPlane) ;
    Point3d ptProjP = ptP + dDistPtPl * plPlane.GetVersN() ;
-  // Verifico se il punto proiettato è esterno al loop esterno
+  // Verifico se il punto proiettato è esterno al loop esterno
    int nPtOut = -1 ;
    if ( ! IsPointInsidePolyLine( ptProjP, vPolyVec[0], EPS_SMALL))
       nPtOut = 0 ;
-  // Verifico se il punto proiettato è interno ai loop interni (quindi esterno alla faccia)
+  // Verifico se il punto proiettato è interno ai loop interni (quindi esterno alla faccia)
    for ( int nLoop = 1 ; nLoop < int( vPolyVec.size()) && nPtOut < 0 ; ++ nLoop) {
       Plane3d plPlane ;
       double dArea ;
@@ -580,7 +601,7 @@ SurfTriMesh::DistPointFacet( const Point3d& ptP, const POLYLINEVECTOR& vPolyVec,
       dPointFacetDist = abs( dDistPtPl) ;
       return true ;
    }
-  // Altrimenti calcolo la minima distanza del punto dalla polilinea del contorno a cui è esterno
+  // Altrimenti calcolo la minima distanza del punto dalla polilinea del contorno a cui è esterno
    double dDist ;
    if ( DistPointPolyLine( ptP, vPolyVec[nPtOut], dDist)) {
       dPointFacetDist = dDist ; 
@@ -593,7 +614,7 @@ SurfTriMesh::DistPointFacet( const Point3d& ptP, const POLYLINEVECTOR& vPolyVec,
 bool
 SurfTriMesh::ChangeStart( const Point3d& ptNewStart, PNTVECTOR& Loop)
 {
-  // Cerco il tratto del loop chiuso più vicino al punto
+  // Cerco il tratto del loop chiuso più vicino al punto
    int nMinSeg = - 1 ;  
    double dMinSqDinst = DBL_MAX ;
    for ( int nPt = 0 ; nPt < int( Loop.size()) ; ++ nPt) {
@@ -657,7 +678,7 @@ SurfTriMesh::ChangeStart( const Point3d& ptNewStart, PNTVECTOR& Loop)
 bool
 SurfTriMesh::SplitAtPoint( const Point3d& ptStop, const PNTVECTOR& Loop, PNTVECTOR& Loop1, PNTVECTOR& Loop2)
 {
-  // Cerco il tratto del loop chiuso più vicino al punto
+  // Cerco il tratto del loop chiuso più vicino al punto
    int nMinSeg = -1 ;
    double dMinSqDinst = DBL_MAX ;
    for ( int nPt = 0 ; nPt < int( Loop.size()) ; ++ nPt) {
@@ -679,13 +700,13 @@ SurfTriMesh::SplitAtPoint( const Point3d& ptStop, const PNTVECTOR& Loop, PNTVECT
   // Verifico che il punto stia su un vertice, in tal caso non devo aggiungerlo
    bool bFirst = AreSamePointApprox( Loop[nMinSeg], ptStop) ;
    bool bLast = AreSamePointApprox( Loop[( nMinSeg + 1) % int( Loop.size())], ptStop) ;
-  // Se il punto è sul vertice finale del segmento, aggiungo il vertice alla lista da inglobare al primo loop
+  // Se il punto è sul vertice finale del segmento, aggiungo il vertice alla lista da inglobare al primo loop
    if ( bLast)
       ++ nMinSeg ;
   // Inglobo fino a nSeg nel primo loop
    for ( int nPt = 0 ; nPt <= nMinSeg  &&  nPt < int( Loop.size()) ; ++ nPt)
       Loop1.emplace_back( Loop[nPt]) ;
-  // Se il punto è interno al segmento, lo inglobo in entrambi i loop
+  // Se il punto è interno al segmento, lo inglobo in entrambi i loop
    if ( ! ( bFirst  ||  bLast)) {
       Loop1.emplace_back( ptStop) ;
       Loop2.emplace_back( ptStop) ;

Tree.cpp

@@ -24,6 +24,7 @@
 #include "/EgtDev/Include/EGkDistPointCurve.h"
 #include "/EgtDev/Include/EGkCurve.h"
 #include "DistPointCrvComposite.h"
+#include "/EgtDev/Include/EGkSfrCreate.h"
 #include <algorithm>
 
 using namespace std ;
@@ -60,6 +61,38 @@ Tree::Tree( const Point3d ptBl, const Point3d ptTr)
    m_mTree.insert( pair< int, Cell>( -1, cRoot)) ;
 }
 
+//----------------------------------------------------------------------------
+bool
+Tree::LimitLoop( PolyLine& pl, POLYLINEVECTOR& vPl) const
+{
+   // creo la flat region di trim, quella del parametrico e li interseco
+   PtrOwner<ISurfFlatRegion> pSfrTrim( GetSurfFlatRegionFromPolyLine( pl)) ;
+   bool bInverted = false ;
+   if ( ! pSfrTrim->GetNormVersor().IsZplus()) {
+      pSfrTrim->Invert() ;
+      bInverted = true ;
+   }
+   PtrOwner<ISurfFlatRegion> pParamTrim( GetSurfFlatRegionRectangle( SBZ_TREG_COEFF * m_nSpanU, SBZ_TREG_COEFF * m_nSpanV)) ;
+   if ( ! pParamTrim->Intersect( *pSfrTrim)  ||  ! pParamTrim->IsValid())
+      return false ;
+
+   // ricostruisco la curva tenendo solo le parti dentro lo spazio parametrico
+   // devo recuperare la polyline dei bordi dei vari chunk creati
+   for ( int c = 0 ; c < int( pParamTrim->GetChunkCount()) ; ++c) {
+      PtrOwner<ICurve> pCrv ( pParamTrim->GetLoop( c, 0)) ;
+      if ( bInverted)
+         pCrv->Invert() ;
+      PolyLine plApprox ;
+      double dLinTol = 10 * EPS_SMALL, dAngTolDeg = 5 ;
+      int nType = 0 ;
+      pCrv->ApproxWithLines( dLinTol, dAngTolDeg, nType, plApprox) ;
+      // aggiungo la polyline del chunk
+      vPl.push_back( plApprox) ;
+   }
+   
+   return false ;
+}
+
 //----------------------------------------------------------------------------
 void
 Tree::SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches, const Point3d& ptMin, const Point3d& ptMax)
@@ -71,7 +104,6 @@ Tree::SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches, const Point3d& ptMi
    m_vLoop.clear() ;
    m_mChunk.clear() ;
    m_vPlApprox.clear() ;
-   m_vChunk.clear() ;
    m_vPolygons.clear() ;
    m_vPlLoop2D.clear() ;
 
@@ -93,8 +125,6 @@ Tree::SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches, const Point3d& ptMi
   // recupero i loop di trim e li divido per chunk
    if ( m_bTrimmed) {
       int nLoop = 0 ;
-      INTVECTOR vChunk ;
-
      // recupero la superficie di trim per avere accesso diretto ai loop e mantenendo le informazioni sui chunk
       PtrOwner<SurfFlatRegion> pTrimReg( m_pSrfBz->GetTrimRegion()->Clone()) ;
       double dLinTol = 0.01 ; // questo è riferito allo spazio parametrico
@@ -102,7 +132,7 @@ Tree::SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches, const Point3d& ptMi
       for ( int i = 0 ; i < pTrimReg->GetChunkCount() ; ++ i) {
          PtrOwner<SurfFlatRegion> pChunk( pTrimReg->CloneChunk( i)) ;
          for ( int j = 0 ; j < pChunk->GetLoopCount( 0) ; ++ j) {
-            vChunk.push_back( nLoop) ;
+            //vChunk.push_back( nLoop) ;
            // i chunk della falt region sono ancora flat region composte da 1 chunk
             PtrOwner<ICurve> pLoop ( pChunk->GetLoop( 0, j)) ;
             // rimuovo i difetti dei loop prima di salvarli
@@ -127,21 +157,13 @@ Tree::SetSurf( const SurfBezier* pSrfBz, bool bSplitPatches, const Point3d& ptMi
                bCCW = true ;
             else
                bCCW = false ;
-            m_vPlApprox.push_back( tuple<PolyLine,bool>(plApprox,bCCW)) ;
+            // limito il loop allo spazio parametrico  // potrei ottenere più loop a partire da quello originale
+            POLYLINEVECTOR vPlLimited ;
+            LimitLoop( plApprox, vPlLimited) ;
+            for ( int k = 0 ; k < int( vPlLimited.size()) ; ++k)
+               m_vPlApprox.push_back( tuple<PolyLine,bool>(vPlLimited[k], bCCW)) ;
          }
-        // li riordino dal più grande al più piccolo
-         sort(vChunk.begin(), vChunk.end(), [this]( const int& a, const int& b) { double dArea1, dArea2 ;
-                                                                                  m_vLoop[a]->GetAreaXY( dArea1) ;
-                                                                                  m_vLoop[b]->GetAreaXY( dArea2) ;
-                                                                                  return dArea1 > dArea2 ; }) ;
-         m_vChunk.push_back( vChunk) ;
-         vChunk.clear() ;
       }
-     // li riordino dal più grande al più piccolo
-      sort(m_vChunk.begin(), m_vChunk.end(), [this]( const INTVECTOR& a, const INTVECTOR& b) { double dArea1, dArea2 ;
-                                                                               m_vLoop[a[0]]->GetAreaXY( dArea1) ;
-                                                                               m_vLoop[b[0]]->GetAreaXY( dArea2) ;
-                                                                               return dArea1 > dArea2 ;}) ;
    }
   // salvo i vertici 3d della cella root
    Point3d ptTop( nSpanU * SBZ_TREG_COEFF, nSpanV * SBZ_TREG_COEFF) ;
@@ -3509,6 +3531,8 @@ Tree::CheckIfBetween( const Inters& inA, const Inters& inB) const
    }
    else if ( nEdge  == 7)
       nEdge = 0 ;
+   if ( AreSameEdge( inA.nIn, inA.nOut)  &&  ! CheckIfBefore( inA.nIn, inA.vpt[0], inA.vpt.back()))
+      vEdges.push_back( nEdge) ;
    // creo la sequenza di Edges da scorrere per trovare i possibili validNextStart
    while ( ! AreSameEdge( nEdge, inA.nIn) || (int) vEdges.size() == 0) {
       vEdges.push_back( nEdge) ;

Tree.h

@@ -254,6 +254,7 @@ class Tree
       std::vector<bool> GetPoles( void) { return m_vbPole ;} ; // funzione che restituisce i flag che indicano se i lati sono collassati in dei poli
 
    private :
+      bool LimitLoop( PolyLine& pl, POLYLINEVECTOR& vPl) const ;              // funzione che limita i loop di trim allo spazio parametrico
       bool Split( int nId, double dSplitValue) ;                              // funzione di split di una cella al parametro indicato nella direzione data da bVert
       bool Split( int nId) ;                                                  // funzione di split di una cella dell'albero a met� nella direzione data da bVert
       void Balance( void) ;                                                       // creo rami in modo che tutte tutte le foglie abbiano come adiacenti foglie ad una profondit� di +- 1
@@ -290,7 +291,7 @@ class Tree
       const SurfBezier*           m_pSrfBz ;           // superficie di bezier
       DBLVECTOR                   m_vDim ;             // distanze tra i vertici della superficie di bezier in 3d in ordine antiorario a partire da ptP00
       bool                        m_bTrimmed ;         // superficie trimmata
-      INTMATRIX                   m_vChunk ;           // elenco dei loop divisi per chunk
+      //INTMATRIX                   m_vChunk ;           // elenco dei loop divisi per chunk
       std::map<int,int>           m_mChunk ;           // mappa in cui vengono salvati chunk di appartenza per ogni loop di trim
       ICURVEPOVECTOR              m_vLoop ;            // curve di loop
       std::vector<std::tuple<PolyLine,bool>> m_vPlApprox ;  // vettore contenente le approssimazioni dei loop  // il bool indica se la curva è CCW

Triangulate.cpp

@@ -1225,7 +1225,7 @@ Triangulate::GetOuterPntToJoin( const PNTVECTOR& vPt, const Point3d& ptP, int& n
                dMinDist = dZ - ptP.z ;
                nI = ( vPt[i].z >= vPt[j].z) ? i : j ;
                double dY = vPt[i].y + ( vPt[j].y - vPt[i].y) * dCoeff ;
-               ptInt.Set( ptP.y, dY, dZ) ;
+               ptInt.Set( ptP.x, dY, dZ) ;
             }
          }
          break ;

Vector3d.cpp

@@ -169,6 +169,10 @@ Vector3d::ToSpherical( double* pdLen, double* pdAngVertDeg, double* pdAngOrizzDe
 bool
 Vector3d::Normalize( double dEps)
 {
+  // verifico validità
+   if ( ! IsValid())
+      return false ;
+
   // se già normalizzato, ok
    double dSqLen = x * x + y * y + z * z ;
    if ( abs( 1.0 - dSqLen) < ( 2 * 1000 * DBL_EPSILON))

VolZmap.cpp

@@ -21,12 +21,14 @@
 #include "DistPointLine.h"
 #include "GeoConst.h"
 #include "/EgtDev/Include/EGkIntersLinePlane.h"
+#include "/EgtDev/Include/EGkIntersLineSurfTm.h"
 #include "/EgtDev/Include/EGkUiUnits.h"
 #include "/EgtDev/Include/EGkIntervals.h"
 #include "/EgtDev/Include/EGkStringUtils3d.h"
 #include "/EgtDev/Include/EgtNumUtils.h"
 #include <thread>
 #include <future>
+#include "IntersLineBox.h"
 
 using namespace std ;
 
@@ -93,6 +95,54 @@ VolZmap::Clear( void)
    return true ;
 }
 
+//----------------------------------------------------------------------------
+bool
+VolZmap::IsEmpty( void)
+{
+  // controllo validità
+   if ( ! IsValid())
+      return false ;
+
+   bool bIsEmpty = true ; // flag
+
+  // ciclo sul numero di mappe
+   for ( int i = 0 ; bIsEmpty  &&   i < m_nMapNum ; ++ i) {
+     // ciclo sulle dimensioni di ogni mappa
+      for ( int j = 0 ; bIsEmpty  &&  j < m_nDim[i] ; ++ j) {
+        // controllo se esiste un tratto di dexel non vuoto
+         bIsEmpty = ( int( m_Values[i][j].size()) == 1  &&
+                    abs( m_Values[i][j][0].dMin) < EPS_ZERO  &&
+                    m_Values[i][j][0].nToolMin == 0  &&
+                    m_Values[i][j][0].nCompo == 1) ;
+
+         if ( bIsEmpty) {
+            switch ( i) {
+               case 0 :
+                  bIsEmpty = ( m_Values[i][j][0].vtMinN.IsZminus()  &&
+                             abs( m_Values[i][j][0].dMax) < EPS_ZERO  &&
+                             m_Values[i][j][0].vtMaxN.IsZplus()  &&
+                             m_Values[i][j][0].nToolMax == 0) ;
+                  break ;
+               case 1 :
+                  bIsEmpty = ( m_Values[i][j][0].vtMinN.IsZminus()  &&
+                             abs( m_Values[i][j][0].dMax) < EPS_ZERO  &&
+                             m_Values[i][j][0].vtMaxN.IsXplus()  &&
+                             m_Values[i][j][0].nToolMax == 0 ) ;
+                  break ;
+               case 2 :
+                  bIsEmpty = ( m_Values[i][j][0].vtMinN.IsYminus()  &&
+                             abs( m_Values[i][j][0].dMax) < EPS_ZERO  &&
+                             m_Values[i][j][0].vtMaxN.IsYplus()  &&
+                             m_Values[i][j][0].nToolMax == 0 ) ;
+                  break ;
+            }
+         }
+      }
+   }
+
+   return bIsEmpty ;
+}
+
 //----------------------------------------------------------------------------
 VolZmap*
 VolZmap::Clone( void) const
@@ -1657,6 +1707,134 @@ VolZmap::SetToModifyDexelBlocks( int nGrid, int nDex, int nInt)
    return true ;
 }
 
+//----------------------------------------------------------------------------
+bool
+VolZmap::AddSrfTm( const ISurfTriMesh* pStm)
+{
+  // controllo sulla superficie
+   double dVol ;
+   if ( pStm == nullptr  ||  ! pStm->IsValid()  ||  ! pStm->IsClosed()  ||
+      ! pStm->GetVolume( dVol)  ||  dVol < 0)
+      return false ;
+
+  // controllo se il Box3d della superficie si interseca con il Box3d dello Zmap corrente
+   BBox3d BBox_stm ;
+   if ( ! pStm->GetLocalBBox( BBox_stm))
+      return false ;
+   BBox3d BBox_curr ;
+   if ( ! GetLocalBBox( BBox_curr))
+      return false ;
+   BBox3d BBox_inters ;
+   if ( BBox_stm.FindIntersection( BBox_curr, BBox_inters)  &&  BBox_inters.IsEmpty())
+      return true ; // se non ci sono intersezioni, la superficie non influenza lo Zmap
+
+  // determino i punti estremi del bounding box corrente
+   BBox3d BBoxCurr ;
+   if ( ! GetLocalBBox( BBoxCurr))
+      return false ;
+   Point3d ptMapOrig, ptMapEnd ;
+   BBoxCurr.GetMinMax( ptMapOrig, ptMapEnd) ;
+  // determino le dimensioni lineari del BBox
+   Vector3d vtLen = ptMapEnd - ptMapOrig ;
+
+  // creo uno Zmap vuoto per la TriMesh
+   PtrOwner<VolZmap> pZmapStm( CreateBasicVolZmap()) ;
+   if ( IsNull( pZmapStm)  ||
+        ! pZmapStm->CreateEmptyMap( ptMapOrig, BBoxCurr.GetDimX() + 10 * EPS_SMALL,
+                                               BBoxCurr.GetDimY() + 10 * EPS_SMALL,
+                                               BBoxCurr.GetDimZ() + 10 * EPS_SMALL,
+                                    m_dStep, IsTriDexel()))
+      return false ;
+
+
+  // ciclo sulle griglie
+   bool bCompleted = true ;
+   for ( int g = 0 ; g < pZmapStm->m_nMapNum ; ++ g) {
+     // definisco dei sistemi di riferimento ausiliari
+      Frame3d frMapFrame ;
+      if ( g == 0)
+         frMapFrame = m_MapFrame ;
+      else if ( g == 1)
+         frMapFrame.Set( ptMapOrig, Y_AX, Z_AX, X_AX) ;
+      else if ( g == 2)
+         frMapFrame.Set( ptMapOrig, Z_AX, X_AX, Y_AX) ;
+
+     // oggetto per calcolo massivo intersezioni
+      IntersParLinesSurfTm intPLSTM( frMapFrame, *pStm) ;
+
+     // numero massimo di thread
+      int nThreadMax = max( 1, int( thread::hardware_concurrency()) - 1) ;
+      vector<future<bool>> vRes ;
+      vRes.resize( nThreadMax) ;
+      if ( pZmapStm->m_nNx[g] > pZmapStm->m_nNy[g]) {
+         int nDexNum = pZmapStm->m_nNx[g] / nThreadMax ;
+         int nRemainder = pZmapStm->m_nNx[g] % nThreadMax ;
+         int nInfI = 0 ;
+         int nSupI = 0 ;
+         for ( int nThread = 0 ; nThread < nThreadMax ; ++ nThread) {
+            nInfI = nSupI ;
+            nSupI = nInfI + ( nThread < nRemainder ? nDexNum + 1 : nDexNum) ;
+            vRes[nThread] =  async( launch::async, &VolZmap::CreateMapPart, this, g,
+                                    nInfI, nSupI, 0, pZmapStm->m_nNy[g], ref( vtLen), ref( ptMapOrig), ref( *pStm), ref( intPLSTM)) ;
+         }
+      }
+      else {
+         int nDexNum = pZmapStm->m_nNy[g] / nThreadMax ;
+         int nRemainder = pZmapStm->m_nNy[g] % nThreadMax ;
+         int nInfJ = 0 ;
+         int nSupJ = 0 ;
+         for ( int nThread = 0 ; nThread < nThreadMax ; ++ nThread) {
+            nInfJ = nSupJ ;
+            nSupJ = nInfJ + ( nThread < nRemainder ? nDexNum + 1 : nDexNum) ;
+            vRes[nThread] =  async( launch::async, &VolZmap::CreateMapPart, this, g,
+                                    0, pZmapStm->m_nNx[g], nInfJ, nSupJ, ref( vtLen), ref( ptMapOrig), ref( *pStm), ref( intPLSTM)) ;
+         }
+      }
+
+     // ciclo per attendere che tutti gli async abbiano terminato.
+      int nTerminated = 0 ;
+      while ( nTerminated < nThreadMax) {
+         for ( int nL = 0 ; nL < nThreadMax ; ++ nL) {
+            // Async terminato
+            if ( vRes[nL].valid()  &&  vRes[nL].wait_for( chrono::microseconds{ 1}) == future_status::ready) {
+               ++ nTerminated ;
+               bCompleted = bCompleted  &&  vRes[nL].get() ;
+            }
+         }
+      }
+
+     // bho ( ???????? )
+      if ( ! bCompleted)
+         return false ;
+
+     // aggiungo ora gli intervalli ricavati allo Zmap corrente ( *this)
+     // ciclo sui dexel della mappa
+      for ( int nD = 0 ; nD < int( pZmapStm->m_Values[g].size()) ; ++ nD) {
+        // se spillone vuoto, passo al successivo
+         if ( pZmapStm->m_Values[g][nD].empty())
+            continue ;
+        // indici di spillone
+         int nI = nD % pZmapStm->m_nNx[g] ;
+         int nJ = nD / pZmapStm->m_nNx[g] ;
+        // ciclo sui voxel associati
+         for ( int nV = 0 ; nV < int( pZmapStm->m_Values[g][nD].size()) ; ++ nV) {
+           // estremi del voxel
+            double dMin = pZmapStm->m_Values[g][nD][nV].dMin ;
+            double dMax = pZmapStm->m_Values[g][nD][nV].dMax ;
+           // vettori associati agli estremi
+            Vector3d vtMaxN = pZmapStm->m_Values[g][nD][nV].vtMaxN ;
+            Vector3d vtMinN = pZmapStm->m_Values[g][nD][nV].vtMinN ;
+           // aggiungo l'intervallo                           // per ora... (???)
+            AddIntervals( g, nI, nJ, dMin, dMax, vtMinN, vtMaxN, 0) ;
+         }
+      }
+
+   }
+
+   return true ;
+}
+
+
 //----------------------------------------------------------------------------
 bool
 VolZmap::IsMapPartABox( int nMap, int nInfI, int nSupI, int nInfJ, int nSupJ, double& dMinZ, double& dMaxZ)

VolZmap.h

@@ -77,7 +77,9 @@ class VolZmap : public IVolZmap, public IGeoObjRW
    public :                      // IVolZmap
       bool CopyFrom( const IGeoObj* pGObjSrc) override ;
       bool Clear( void) override ;
+      bool IsEmpty( void) override ;
       bool Create( const Point3d& ptO, double dDimX, double dDimY, double dDimZ, double dStep, bool bTriDex) override ;
+      bool CreateEmptyMap( const Point3d& ptO, double dLengthX, double dLengthY, double dLengthZ, double dStep, bool bTriDex) override ;
       bool CreateFromFlatRegion( const ISurfFlatRegion& Surf, double dDimZ, double dStep, bool bTriDex) override ;
       bool CreateFromTriMesh( const ISurfTriMesh& Surf, double dStep, bool bTriDex) override ;
       int  GetBlockCount( void) const override ;
@@ -89,7 +91,7 @@ class VolZmap : public IVolZmap, public IGeoObjRW
               { return m_nMapNum == 3 ; }
       bool GetDexelLines( int nDir, int nPos1, int nPos2, POLYLINELIST& lstPL) const override ;
       int  GetResolution( void) const override
-              { return m_nDexVoxRatio ; }   
+              { return m_nDexVoxRatio ; }
       bool ChangeResolution( int nDexVoxRatio) override ;
       void SetShowEdges( bool bShow) override
               { m_bShowEdges = bShow ;  // qui è necessario far ricreare la grafica
@@ -140,6 +142,7 @@ class VolZmap : public IVolZmap, public IGeoObjRW
       VolZmap* ClonePart( int nPart) const override ;
       bool RemovePart( int nPart) override ;
       int  GetPartMinDistFromPoint( const Point3d& ptP) const override ;
+      bool AddSrfTm( const ISurfTriMesh* pStm) override ;
 
    public :                      // IGeoObjRW
       int GetNgeId( void) const override ;
@@ -411,6 +414,9 @@ class VolZmap : public IVolZmap, public IGeoObjRW
       // Funzione per crezione solido in parallelo
       bool CreateMapPart( int nMap, int nInfI, int nSupI, int nInfJ, int nSupJ, const Vector3d& vtLen, const Point3d& ptMapOrig,
                           const ISurfTriMesh& Surf, IntersParLinesSurfTm& intPLSTM) ;
+      // Funzione per aggiornamento solido in parallelo
+      bool UpdateMapPart( int nMap, int nInfI, int nSupI, int nInfJ, int nSupJ, const Vector3d& vtLen, const Point3d& ptMapOrig,
+                          const ISurfTriMesh& Surf, IntersParLinesSurfTm& intPLSTM) ;
 
    private :
       enum Status { ERR = 0, OK = 1, TO_VERIFY = 2} ;

VolZmapCreation.cpp

@@ -1,4 +1,4 @@
-//----------------------------------------------------------------------------
+//----------------------------------------------------------------------------
 //  EgalTech 2015-2016
 //----------------------------------------------------------------------------
 // File : VolZmap.cpp         Data : 22.01.15           Versione : 1.6a4
@@ -26,14 +26,124 @@ using namespace std ;
 // ------------------------- CREAZIONE MAPPA --------------------------------------------------------------------------------------
 
 //----------------------------------------------------------------------------
-bool 
+bool
 VolZmap::Create( const Point3d& ptO, double dLengthX, double dLengthY, double dLengthZ, double dStep, bool bTriDex)
 {
-  // Controlli sull'ammissibilità delle dimensioni lineari del grezzo e del passo
+
+  // Controlli sull'ammissibilità delle dimensioni lineari del grezzo e del passo
    if ( dStep < EPS_SMALL  ||  dLengthX < EPS_SMALL  ||  dLengthY < EPS_SMALL  ||  dLengthZ < EPS_SMALL)
       return false ;
 
-  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
+  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
+   m_dStep = max( dStep, 100 * EPS_SMALL) ;
+
+  // Aggiorno la dimensione della mappa 1 o 3
+   m_nMapNum = ( bTriDex ? 3 : 1) ;
+
+  // Disponendo i sistemi di riferimento in una successione, le coordinate x,y,z
+  // di uno si ottengono da una permutazione ciclica di quelle del precedente sistema.
+  // es: X(n) = Z(n-1), Y(n) = X(n-1), Z(n) = Y(n-1)
+
+  // Definisco il sistema di riferimento intrinseco
+   m_MapFrame.Set( ptO, X_AX, Y_AX, Z_AX) ;
+
+  // Definisco i vettori dei limiti su indici
+   m_nNx[0] = max( int( ( dLengthX + EPS_SMALL) / m_dStep + 0.5), 1) ;
+   m_nNy[0] = max( int( ( dLengthY + EPS_SMALL) / m_dStep + 0.5), 1) ;
+
+  // Numero di componenti connesse 
+   m_nConnectedCompoCount = 1 ;
+
+  // Se tridexel
+   if ( bTriDex) {
+      m_nNx[1] = m_nNy[0] ;
+      m_nNy[1] = max( int( ( dLengthZ + EPS_SMALL) / m_dStep + 0.5), 1) ;
+      m_nNx[2] = m_nNy[1] ;
+      m_nNy[2] = m_nNx[0] ;
+   }
+
+  // altrimenti mono dexel
+   else {
+      m_nNx[1] = 0 ;
+      m_nNy[1] = 0 ;
+      m_nNx[2] = 0 ;
+      m_nNy[2] = 0 ; 
+   }
+
+  // Definisco il numero di blocchi lungo x,y e z
+   if ( ! CalcBlockNum())
+      return false ;
+
+  // Definizione della mappa
+
+  // Creazione delle mappe
+  // Calcolo del numero di celle per ogni mappa
+   for ( int i = 0 ; i < m_nMapNum ; ++ i)
+      m_nDim[i] = m_nNx[i] * m_nNy[i] ;
+
+  // Creazione delle celle per ogni mappa
+   for ( int i = 0 ; i < m_nMapNum ; ++ i)
+      m_Values[i].resize( m_nDim[i]) ;
+
+  // Riempimento delle celle
+   for ( int i = 0 ; i < m_nMapNum ; ++ i) {
+      for ( int j = 0 ; j < m_nDim[i] ; ++ j) {
+
+        // Aggiungo il tratto al dexel vuoto
+         m_Values[i][j].resize( 1) ;
+         m_Values[i][j][0].dMin = 0 ;
+         m_Values[i][j][0].nToolMin = 0 ;
+         m_Values[i][j][0].nCompo = 1 ;
+
+         switch ( i) {
+         case 0 :
+            m_Values[i][j][0].vtMinN = - Z_AX ;
+            m_Values[i][j][0].dMax = dLengthZ ; 
+            m_Values[i][j][0].vtMaxN = Z_AX ;
+            m_Values[i][j][0].nToolMax = 0 ;
+            break ;
+         case 1 :
+            m_Values[i][j][0].vtMinN = - X_AX ;
+            m_Values[i][j][0].dMax = dLengthX ;
+            m_Values[i][j][0].vtMaxN = X_AX ;
+            m_Values[i][j][0].nToolMax = 0 ;
+            break ;
+         case 2 :
+            m_Values[i][j][0].vtMinN = - Y_AX ;
+            m_Values[i][j][0].dMax = dLengthY ;
+            m_Values[i][j][0].vtMaxN = Y_AX ;
+            m_Values[i][j][0].nToolMax = 0 ;
+            break ; 
+         }
+      }
+   }
+
+  // Definizione delle limitazioni iniziali in Z per ogni mappa
+   m_dMinZ[0] = 0 ;
+   m_dMaxZ[0] = dLengthZ ;
+   m_dMinZ[1] = 0 ;
+   m_dMaxZ[1] = ( bTriDex ? dLengthX : 0) ;
+   m_dMinZ[2] = 0 ;  
+   m_dMaxZ[2] = ( bTriDex ? dLengthY : 0) ;
+
+  // Tipologia
+   m_nShape = BOX ;
+
+  // Aggiornamento dello stato
+   m_nStatus = OK ;
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool 
+VolZmap::CreateEmptyMap( const Point3d& ptO, double dLengthX, double dLengthY, double dLengthZ, double dStep, bool bTriDex)
+{
+  // Controlli sull'ammissibilità delle dimensioni lineari del grezzo e del passo
+   if ( dStep < EPS_SMALL  ||  dLengthX < EPS_SMALL  ||  dLengthY < EPS_SMALL  ||  dLengthZ < EPS_SMALL)
+      return false ;
+
+  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
    m_dStep = max( dStep, 100 * EPS_SMALL) ;
 
   // Aggiorno la dimensione della mappa 1 o 3
@@ -73,8 +183,6 @@ VolZmap::Create( const Point3d& ptO, double dLengthX, double dLengthY, double dL
    if ( ! CalcBlockNum())
       return false ;
 
-  // Definizione della mappa
-
   // Creazione delle mappe
   // Calcolo del numero di celle per ogni mappa
    for ( int i = 0 ; i < m_nMapNum ; ++ i)
@@ -84,38 +192,6 @@ VolZmap::Create( const Point3d& ptO, double dLengthX, double dLengthY, double dL
    for ( int i = 0 ; i < m_nMapNum ; ++ i)
       m_Values[i].resize( m_nDim[i]) ;
 
-  // Riempimento delle celle
-   for ( int i = 0 ; i < m_nMapNum ; ++ i)
-      for ( int j = 0 ; j < m_nDim[i] ; ++ j) {
-
-        // Aggiungo il tratto al dexel vuoto
-         m_Values[i][j].resize( 1) ;
-         m_Values[i][j][0].dMin = 0 ;
-         m_Values[i][j][0].nToolMin = 0 ;
-         m_Values[i][j][0].nCompo = 1 ;
-
-         switch ( i) {
-         case 0 :
-            m_Values[i][j][0].vtMinN = - Z_AX ;
-            m_Values[i][j][0].dMax = dLengthZ ; 
-            m_Values[i][j][0].vtMaxN = Z_AX ;
-            m_Values[i][j][0].nToolMax = 0 ;
-            break ;
-         case 1 :
-            m_Values[i][j][0].vtMinN = - X_AX ;
-            m_Values[i][j][0].dMax = dLengthX ;
-            m_Values[i][j][0].vtMaxN = X_AX ;
-            m_Values[i][j][0].nToolMax = 0 ;
-            break ;
-         case 2 :
-            m_Values[i][j][0].vtMinN = - Y_AX ;
-            m_Values[i][j][0].dMax = dLengthY ;
-            m_Values[i][j][0].vtMaxN = Y_AX ;
-            m_Values[i][j][0].nToolMax = 0 ;
-            break ; 
-         }
-      }
-
   // Definizione delle limitazioni iniziali in Z per ogni mappa
    m_dMinZ[0] = 0 ;
    m_dMaxZ[0] = dLengthZ ;
@@ -124,8 +200,8 @@ VolZmap::Create( const Point3d& ptO, double dLengthX, double dLengthY, double dL
    m_dMinZ[2] = 0 ;  
    m_dMaxZ[2] = ( bTriDex ? dLengthY : 0) ;
 
-   // Tipologia
-   m_nShape = BOX ;
+  // Tipologia
+   m_nShape = GENERIC ;
 
   // Aggiornamento dello stato
    m_nStatus = OK ;
@@ -140,7 +216,7 @@ VolZmap::CreateFromFlatRegion( const ISurfFlatRegion& Surf, double dDimZ, double
   // Aggiorno la dimensione della mappa 1 o 3
    m_nMapNum = ( bTriDex ? 3 : 1) ;
 
-  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
+  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
    m_dStep = max( dStep, 100 * EPS_SMALL) ;
  
   // Determino il bounding box della flat region
@@ -223,11 +299,11 @@ VolZmap::CreateFromFlatRegion( const ISurfFlatRegion& Surf, double dDimZ, double
       CRVCVECTOR IntersectionResults ;
       Surf.GetCurveClassification( GridLine, EPS_SMALL, IntersectionResults) ;
 
-     // Analizzo le parti in cui la retta è stata divisa
+     // Analizzo le parti in cui la retta è stata divisa
       int nPart = int( IntersectionResults.size()) ;
       for ( int k = 0 ; k < nPart ; ++ k) {
 
-        // Se la retta è interna alla regione o coincidente con parte della sua frontiera
+        // Se la retta è interna alla regione o coincidente con parte della sua frontiera
          int nType = IntersectionResults[k].nClass ;
          if ( nType == CRVC_IN  ||  nType == CRVC_ON_P  ||  nType == CRVC_ON_M) {
 
@@ -334,7 +410,7 @@ VolZmap::CreateFromFlatRegion( const ISurfFlatRegion& Surf, double dDimZ, double
          int nPart = int( IntersectionResults.size()) ; 
          for ( int k = 0 ; k < nPart ; ++ k) {
 
-           // Se la retta è interna alla regione o coincidente con parte della sua frontiera
+           // Se la retta è interna alla regione o coincidente con parte della sua frontiera
             int nType = IntersectionResults[k].nClass ;
             if ( nType == CRVC_IN  ||  nType == CRVC_ON_P  ||  nType == CRVC_ON_M) {
 
@@ -471,7 +547,7 @@ VolZmap::CreateMapPart( int nMap, int nInfI, int nSupI, int nInfJ, int nSupJ, co
 
             int nIntType = IntersectionResults[k].nILTT ;
 
-            // Se c'è intersezione
+            // Se c'è intersezione
             if ( nIntType != ILTT_NO) {
 
                double dCos = IntersectionResults[k].dCosDN ;
@@ -523,11 +599,115 @@ VolZmap::CreateMapPart( int nMap, int nInfI, int nSupI, int nInfJ, int nSupJ, co
    return true ;
 }
 
+//----------------------------------------------------------------------------
+bool
+VolZmap::UpdateMapPart( int nMap, int nInfI, int nSupI, int nInfJ, int nSupJ, const Vector3d& vtLen, const Point3d& ptMapOrig,
+                          const ISurfTriMesh& Surf, IntersParLinesSurfTm& intPLSTM)
+{
+      if ( nMap < 0  ||  nMap > 2  ||
+        nInfI < 0  ||  nInfI > m_nNx[nMap]  ||
+        nSupI < 0  ||  nSupI > m_nNx[nMap]  ||
+        nInfJ < 0  ||  nInfJ > m_nNy[nMap]  ||
+        nSupJ < 0  ||  nSupJ > m_nNy[nMap])
+      return false ;
+
+   // Determinazione e ridimensionamento dei dexel interni alla trimesh
+   for ( int i = nInfI ; i < nSupI ; ++ i) {
+      for ( int j = nInfJ ; j < nSupJ ; ++ j) {
+
+         // Definisco la retta da intersecare con la trimesh
+         double dX = ( i + 0.5) * m_dStep ;
+         double dY = ( j + 0.5) * m_dStep ;
+         Point3d ptP0( dX, dY, 0) ;
+
+         // Determino le intersezioni della retta con la TriMesh
+         ILSIVECTOR IntersectionResults ;
+         intPLSTM.GetInters( ptP0, vtLen.v[(nMap+2)%3], IntersectionResults) ;
+
+         for ( int nI = 0 ; nI < int( IntersectionResults.size()) - 3 ; ++ nI) {
+            int nJ = nI + 1 ;  
+            int nK = nJ + 1 ;
+            int nT = nK + 1 ;
+            int nSgnI = IntersectionResults[nI].dCosDN > EPS_SMALL ? 1 : IntersectionResults[nI].dCosDN > -EPS_SMALL ? 0 : - 1 ;
+            int nSgnJ = IntersectionResults[nJ].dCosDN > EPS_SMALL ? 1 : IntersectionResults[nJ].dCosDN > -EPS_SMALL ? 0 : - 1 ;
+            int nSgnK = IntersectionResults[nK].dCosDN > EPS_SMALL ? 1 : IntersectionResults[nK].dCosDN > -EPS_SMALL ? 0 : - 1 ;
+            int nSgnT = IntersectionResults[nT].dCosDN > EPS_SMALL ? 1 : IntersectionResults[nT].dCosDN > -EPS_SMALL ? 0 : - 1 ;
+            double dUJ = IntersectionResults[nJ].dU ;
+            double dUK = IntersectionResults[nK].dU ;
+            if ( nSgnI != 0  &&  nSgnI == nSgnJ  &&  nSgnK != 0  &&  nSgnK == nSgnT  &&  nSgnI == - nSgnT  &&  abs( dUJ - dUK) < EPS_SMALL) {
+               IntersectionResults.erase( IntersectionResults.begin() + nK) ;
+               IntersectionResults.erase( IntersectionResults.begin() + nJ) ;
+            }
+         }
+
+         int nInt = int( IntersectionResults.size()) ;
+
+         int nPos = j * m_nNx[nMap] + i ;
+
+         bool bInside = false ;
+         Point3d ptIn ;
+         Vector3d vtInN ;
+
+         for ( int k = 0 ; k < nInt ; ++ k) {
+
+            int nIntType = IntersectionResults[k].nILTT ;
+
+            // Se c'è intersezione
+            if ( nIntType != ILTT_NO) {
+
+               double dCos = IntersectionResults[k].dCosDN ;
+
+               // entro nella superficie trimesh
+               if ( dCos < - EPS_SMALL) {
+
+                  ptIn = IntersectionResults[k].ptI ;
+
+                  int nT = IntersectionResults[k].nT ;
+                  int nF = Surf.GetFacetFromTria( nT) ;
+
+                  Surf.GetFacetNormal( nF, vtInN) ;
+
+                  bInside = true ;
+               }
+
+               // esco dalla superficie trimesh
+               else if ( dCos > EPS_SMALL  &&  bInside) {
+
+                  Point3d ptOut = IntersectionResults[k].ptI ;
+
+                  int nT = IntersectionResults[k].nT ;
+                  int nF = Surf.GetFacetFromTria( nT) ;
+
+                  Vector3d vtOutN ;
+                  Surf.GetFacetNormal( nF, vtOutN) ; 
+
+                  int nCurrentSize = int( m_Values[nMap][nPos].size()) ;
+
+                  // Aggiungo un tratto al dexel
+                  m_Values[nMap][nPos].resize( nCurrentSize + 1) ;
+
+                  // Aggiorno dati del tratto di dexel 
+                  // Aggiungo un tratto al dexel attuale ( allo spillone )
+                   AddIntervals( nMap, nInfI, nInfJ,
+                                 ptIn.v[(nMap+2)%3] - ptMapOrig.v[(nMap+2)%3],
+                                 ptOut.v[(nMap+2)%3] - ptMapOrig.v[(nMap+2)%3],
+                                 vtInN, vtOutN, 0) ;
+
+                  bInside = false ;
+               }
+            }   
+         }
+      }
+   }
+   return true ;
+
+}
+
 //----------------------------------------------------------------------------
 bool
 VolZmap::CreateFromTriMesh( const ISurfTriMesh& Surf, double dStep, bool bTriDex)
 {
-  // Se la superficie non è chiusa oppure orientata al contrario non ha senso continuare
+  // Se la superficie non è chiusa oppure orientata al contrario non ha senso continuare
    double dVol ;
    if ( ! Surf.IsClosed()  ||  ! Surf.GetVolume( dVol)  ||  dVol < 0)
       return false ;
@@ -543,14 +723,14 @@ VolZmap::CreateFromTriMesh( const ISurfTriMesh& Surf, double dStep, bool bTriDex
    Point3d ptMapOrig, ptMapEnd ; 
    SurfBBox.GetMinMax( ptMapOrig, ptMapEnd) ;
 
-  // Il dexel se parte da un triangolo della trimesh può non trovare l'intersezione,
+  // Il dexel se parte da un triangolo della trimesh può non trovare l'intersezione,
   // quindi espandiamo il bounding box per ovviare al problema.
    SurfBBox.Expand( 100 * EPS_SMALL, 100 * EPS_SMALL, 100 * EPS_SMALL) ;
 
   // Sistema di riferimento intrinseco dello Zmap
    m_MapFrame.Set( ptMapOrig, Frame3d::TOP) ;
 
-  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
+  // Il passo di discretizzazione non può essere inferiore a 100 * EPS_SMALL
    m_dStep = max( dStep, 100 * EPS_SMALL) ;
 
   // Determino le dimensioni lineari del BBox

VolZmapVolume.cpp

@@ -132,10 +132,10 @@ VolZmap::SubtractIntervals( int nGrid, int nI, int nJ,
 
   // Imposto ricalcolo della grafica
    m_OGrMgr.Reset() ;
-  // Imposto forma generica 
-   m_nShape = GENERIC ;    
+  // Imposto forma generica
+   m_nShape = GENERIC ;
   // Imposto ricalcolo numero di componenti connesse
-   m_nConnectedCompoCount = - 1 ; 
+   m_nConnectedCompoCount = - 1 ;
 
   // Passo da indici di dexel a indici di voxel
    nI /= m_nDexVoxRatio ;
@@ -160,7 +160,7 @@ VolZmap::SubtractIntervals( int nGrid, int nI, int nJ,
          ++ nYStop ;
       }
      // Voxel lungo Z
-      int nVoxNumZ = int( m_nNy[1] / m_nDexVoxRatio + ( m_nNy[1] % m_nDexVoxRatio == 0 ? 1 : 2)) ;  
+      int nVoxNumZ = int( m_nNy[1] / m_nDexVoxRatio + ( m_nNy[1] % m_nDexVoxRatio == 0 ? 1 : 2)) ;
       int nMinK = Clamp( int( floor( ( ( dMin - 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) - EPS_SMALL))), 0, nVoxNumZ - 2) ;
       int nMaxK = Clamp( int( floor( ( ( dMax + 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) + EPS_SMALL))), 0, nVoxNumZ - 2) ;
       int nMinZBlock = ( m_nMapNum == 1 ? 0 : Clamp( nMinK / int( m_nVoxNumPerBlock), 0, int( m_nFracLin[2] - 1))) ;
@@ -227,8 +227,8 @@ VolZmap::SubtractIntervals( int nGrid, int nI, int nJ,
          nZBlock[1] = nZBlock[0] - 1 ;
          ++ nZStop ;
       }
-     // Voxel lungo Y 
-      int nVoxNumY = int( m_nNy[0] / m_nDexVoxRatio + ( m_nNy[0] % m_nDexVoxRatio == 0 ? 1 : 2)) ; 
+     // Voxel lungo Y
+      int nVoxNumY = int( m_nNy[0] / m_nDexVoxRatio + ( m_nNy[0] % m_nDexVoxRatio == 0 ? 1 : 2)) ;
       int nMinJ = Clamp( int( floor( ( ( dMin - 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) - EPS_SMALL))), 0, nVoxNumY - 2) ; 
       int nMaxJ = Clamp( int( floor( ( ( dMax + 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) + EPS_SMALL))), 0, nVoxNumY - 2) ; 
       int nMinYBlock = Clamp( nMinJ / int( m_nVoxNumPerBlock), 0, int( m_nFracLin[1] - 1)) ;
@@ -249,232 +249,355 @@ VolZmap::SubtractIntervals( int nGrid, int nI, int nJ,
 
 //----------------------------------------------------------------------------
 bool
-VolZmap::AddIntervals( int nGrid, int nI, int nJ,
-                       double dMin, double dMax, const Vector3d& vtNMin, const Vector3d& vtNMax, int nToolNum)
+VolZmap::AddIntervals( int nGrid, int nI, int nJ, double dMin, double dMax, const Vector3d& vtNMin,
+                       const Vector3d& vtNMax, int nToolNum)
 {
-//   // Controllo che dMin e dMax non siano quasi coincidenti
-//   if ( abs( dMax - dMin) < EPS_SMALL)
-//      return true ;
-//
-//   // Controllo che il numero di griglia sia entro i limiti
-//   if ( nGrid < 0 || nGrid > 2)
-//      return false ;
-//
-//   // Controllo che dMin < dMax
-//   if ( dMin > dMax)
-//      swap( dMin, dMax) ; 
-//   
-//   // Controllo che indici nI, nJ siano entro i limiti
-//   if ( nI < 0 && nI >= m_nNx[nGrid] && 
-//        nJ < 0 && nJ >= m_nNy[nGrid])
-//      return false ;
-//
-//  // Riporto le coordinate cicliche nell'ordine di partenza
-//   Vector3d vtNmi = vtNMin ;
-//   Vector3d vtNma = vtNMax ;
-//   if ( nGrid == 1) {
-//      swap( vtNmi.x, vtNmi.z) ;
-//      swap( vtNmi.y, vtNmi.z) ;
-//      swap( vtNma.x, vtNma.z) ;
-//      swap( vtNma.y, vtNma.z) ;
-//   }
-//   else if ( nGrid == 2) {
-//      swap( vtNmi.y, vtNmi.z) ;
-//      swap( vtNmi.x, vtNmi.z) ;
-//      swap( vtNma.y, vtNma.z) ;
-//      swap( vtNma.x, vtNma.z) ;
-//   }
-//
-//   // Calcolo nPos
-//   unsigned int nPos = nJ * m_nNx[nGrid] + nI ;
-//
-//   
-//   // Se spillone vuoto
-//   if ( m_Values[nGrid][nPos].size() == 0) {
-//   
-//      m_Values[nGrid][nPos].resize( 1) ;
-//   
-//      m_Values[nGrid][nPos][0].dMin = dMin ;
-//      m_Values[nGrid][nPos][0].dMax = dMax ;
-//
-//      m_Values[nGrid][nPos][0].vtMinN = vtNmi ;
-//      m_Values[nGrid][nPos][0].vtMaxN = vtNma ;
-//
-//      if ( dMax > m_dMaxZ[nGrid])
-//         m_dMinZ[nGrid] = dMax ;
-//
-//      if ( dMin < m_dMinZ[nGrid])
-//         m_dMinZ[nGrid] = dMin ;
-//
-//      m_OGrMgr.Reset() ;
-//
-//      return true ;
-//   }
-//
-//   // Ciclo sugli intervalli dello spillone
-//   bool bModified = false ;
-//   unsigned int i = 0 ;
-//   while ( i < m_Values[nGrid][nPos].size()) {
-//
-//      // Eventuale aggiustamento di intervalli sovrapposti
-//      if ( i < m_Values[nGrid][nPos].size() - 1) {
-//         if ( m_Values[nGrid][nPos][i].dMax > m_Values[nGrid][nPos][i + 1].dMin - EPS_SMALL) {
-//
-//            // Se l'intervallo corrente non è contenuto totalmente si esegue l'istruzione successiva
-//            if ( m_Values[nGrid][nPos][i].dMin < m_Values[nGrid][nPos][i + 1].dMin + EPS_SMALL) {
-//         
-//               m_Values[nGrid][nPos][i].dMax = m_Values[nGrid][nPos][i + 1].dMax ;
-//               m_Values[nGrid][nPos][i].vtMaxN = m_Values[nGrid][nPos][i + 1].vtMaxN ;
-//            }
-//            // altrimenti
-//            else {
-//            
-//               m_Values[nGrid][nPos][i].dMin = m_Values[nGrid][nPos][i].dMin ;
-//               m_Values[nGrid][nPos][i].vtMinN = m_Values[nGrid][nPos][i].vtMinN ;
-//
-//               m_Values[nGrid][nPos][i].dMax = m_Values[nGrid][nPos][i].dMax ;
-//               m_Values[nGrid][nPos][i].vtMaxN = m_Values[nGrid][nPos][i].vtMaxN ;
-//            }
-//         
-//            for ( unsigned int j = i + 1 ; j < m_Values[nGrid][nPos].size() - 1 ; ++ j) {
-//
-//               m_Values[nGrid][nPos][j].dMin = m_Values[nGrid][nPos][j + 1].dMin ;
-//               m_Values[nGrid][nPos][j].vtMinN = m_Values[nGrid][nPos][j + 1].vtMinN ;
-//
-//               m_Values[nGrid][nPos][j].dMax = m_Values[nGrid][nPos][j + 1].dMax ;
-//               m_Values[nGrid][nPos][j].vtMaxN = m_Values[nGrid][nPos][j + 1].vtMaxN ;
-//            }
-//         
-//            m_Values[nGrid][nPos].resize( m_Values[nGrid][nPos].size() - 1) ; 
-//            
-//            i = i - 1 ; 
-//         }
-//      }
-//     
-//
-//      // Caso in cui devo aggiungere un intervallo a sinistra dell'intervallo corrente
-//      if ( m_Values[nGrid][nPos][i].dZVal > dMax + EPS_SMALL) {
-//
-//         bModified = true ;   
-//                  
-//         m_Values[nGrid][nPos].resize( m_Values[nGrid][nPos].size() + 2) ;
-//         
-//         for ( size_t j = m_Values[nGrid][nPos].size() - 1 ; j >= i + 2 ; -- j) {
-//         
-//            m_Values[nGrid][nPos][j].dZVal = m_Values[nGrid][nPos][j - 2].dZVal ;
-//            m_Values[nGrid][nPos][j].vtN = m_Values[nGrid][nPos][j - 2].vtN ;
-//         }
-//            
-//         
-//         m_Values[nGrid][nPos][i].dZVal = dMin ;
-//         m_Values[nGrid][nPos][i + 1].dZVal = dMax ;
-//
-//         m_Values[nGrid][nPos][i].vtN = vtNMin ;
-//         m_Values[nGrid][nPos][i + 1].vtN = vtNMax ;
-//
-//         i = i + 2 ;
-//      } 
-//
-//      // Casi d'intersezione:
-//      else if ( m_Values[nGrid][nPos][i + 1].dZVal > dMax - EPS_SMALL) {
-//
-//         // Se l'intervallo da aggiungere sconfina a sinistra modifico il minimo dell'intervalo corrente
-//         if ( m_Values[nGrid][nPos][i].dZVal > dMin - EPS_SMALL) { 
-//
-//            bModified = true ;         
-//            m_Values[nGrid][nPos][i].dZVal = dMin ;
-//            m_Values[nGrid][nPos][i].vtN = vtNmi ;
-//         }
-//      }
-//
-//      else {
-//         // Se l'intervallo corrente è tutto contenuto nell'intervallo da aggungere modifico gli estremi    
-//         if ( m_Values[nGrid][nPos][i].dZVal > dMin + EPS_SMALL) {
-//
-//               bModified = true ;         
-//               m_Values[nGrid][nPos][i].dZVal = dMin ;
-//               m_Values[nGrid][nPos][i + 1].dZVal = dMax ;
-//               m_Values[nGrid][nPos][i].vtN = vtNMin ;
-//               m_Values[nGrid][nPos][i + 1].vtN = vtNma ;
-//         }
-//         // Se l'intervallo da aggiungere sconfina a destra modifico il massimo dell'intervallo corrente 
-//         else if ( m_Values[nGrid][nPos][i + 1].dZVal > dMin - EPS_SMALL) {
-//
-//            bModified = true ;         
-//            m_Values[nGrid][nPos][i + 1].dZVal = dMax ;
-//            m_Values[nGrid][nPos][i + 1].vtN = vtNma ;
-//         }
-//         else {
-//            // Aggiungo intervallo a destra dell'ultimo intervallo
-//            if ( i == m_Values[nGrid][nPos].size() - 2) {
-//
-//               bModified = true ;         
-//               m_Values[nGrid][nPos].resize( m_Values[nGrid][nPos].size() + 2) ;
-//               
-//               m_Values[nGrid][nPos][i + 2].dZVal = dMin ;
-//               m_Values[nGrid][nPos][i + 3].dZVal = dMax ;
-//               m_Values[nGrid][nPos][i + 2].vtN = vtNmi ;
-//               m_Values[nGrid][nPos][i + 3].vtN = vtNma ;
-//
-//               i = i + 2 ;
-//            }
-//         }
-//      }
-//
-//      i = i + 2 ;
-//   }
-//
-//   // se eseguita modifica, imposto ricalcolo della grafica
-//   if ( bModified) {
-//
-//      // Determino quali blocchi sono stati modificati
-//      int nLayerBlock = m_nFracLin[0] * m_nFracLin[1] ;
-//  
-//      if ( nGrid == 0) {
-//
-//         int nXBlock = min( nI / m_nDexNumPBlock, m_nFracLin[0] - 1) ;
-//         int nYBlock = min( nJ / m_nDexNumPBlock, m_nFracLin[1] - 1) ;
-//         int nMinZBlock = max( 0, int( floor( ( dMin / m_dStep))) / int( m_nDexNumPBlock)) ;
-//         int nMaxZBlock = min( int( m_nFracLin[2] - 1), int( floor( ( dMax / m_dStep))) / int( m_nDexNumPBlock)) ;
-//         
-//         for ( int k = nMinZBlock ; k <= nMaxZBlock ; ++ k) 
-//
-//            m_BlockToUpdate[k * nLayerBlock + nYBlock * m_nFracLin[0] + nXBlock] = true ;
-//      }
-//      else if ( nGrid == 1) {
-//
-//         int nYBlock = min( nI / m_nDexNumPBlock, m_nFracLin[1] - 1) ;
-//         int nZBlock = min( nJ / m_nDexNumPBlock, m_nFracLin[2] - 1) ;
-//         int nMinXBlock = max( 0, int( floor( ( dMin / m_dStep))) / int( m_nDexNumPBlock)) ;
-//         int nMaxXBlock = min( int( m_nFracLin[0] - 1), int( floor( ( dMax / m_dStep))) / int( m_nDexNumPBlock)) ;
-//         
-//         for ( int k = nMinXBlock ; k <= nMaxXBlock ; ++ k) 
-//
-//            m_BlockToUpdate[nZBlock * nLayerBlock + nYBlock * m_nFracLin[0] + k] = true ;
-//      }
-//      else if ( nGrid == 2) {
-//
-//         int nXBlock = min( nJ / m_nDexNumPBlock, m_nFracLin[0] - 1) ;
-//         int nZBlock = min( nI / m_nDexNumPBlock, m_nFracLin[2] - 1) ;
-//         int nMinYBlock = max( 0, int( floor( ( dMin / m_dStep))) / int( m_nDexNumPBlock)) ;
-//         int nMaxYBlock = min( int( m_nFracLin[1] - 1), int( floor( ( dMax / m_dStep))) / int( m_nDexNumPBlock)) ;
-//         
-//         for ( int k = nMinYBlock ; k <= nMaxYBlock ; ++ k) 
-//
-//            m_BlockToUpdate[nZBlock * nLayerBlock + k * m_nFracLin[0] + nXBlock] = true ;
-//      }
-//
-//      m_OGrMgr.Reset() ;
-//
-//      // Aggiorno massima e minima Z
-//      // sullo Zmap
-//      if ( dMax > m_dMaxZ[nGrid])
-//         m_dMinZ[nGrid] = dMax ;
-//
-//      if ( dMin < m_dMinZ[nGrid])
-//         m_dMinZ[nGrid] = dMin ;
-//   }
-//
+  // Controllo che il numero di griglia sia entro i limiti
+   if ( nGrid < 0 || nGrid > 2)
+      return false ;
+
+   // Controllo che indici nI, nJ siano entro i limiti
+   if ( nI < 0 && nI >= m_nNx[nGrid]  &&  nJ < 0 && nJ >= m_nNy[nGrid])
+      return false ;
+
+  // Controllo che dMin < dMax ( nel caso swap parametri e vettori )
+   Vector3d vtNmi = vtNMin ;
+   Vector3d vtNma = vtNMax ;
+   if ( dMin > dMax) {
+      swap( dMin, dMax) ;
+      swap( vtNmi, vtNma) ;
+   }
+
+  // Restringo minimo e massimo entro i limiti della mappa
+   if ( dMin < m_dMinZ[nGrid]) {
+      dMin = m_dMinZ[nGrid] ;
+      vtNmi = - Z_AX ;
+   }
+   else if ( dMin > m_dMaxZ[nGrid]) {
+      dMin = m_dMaxZ[nGrid] ;
+      vtNmi = - Z_AX ;
+   }
+   if ( dMax < m_dMinZ[nGrid]) {
+      dMax = m_dMinZ[nGrid] ;
+      vtNma = Z_AX ;
+   }
+   else if ( dMax > m_dMaxZ[nGrid]) {
+      dMax = m_dMaxZ[nGrid] ;
+      vtNma = Z_AX ;
+   }
+
+  // Controllo che dMin e dMax non siano quasi coincidenti ( la limitazione dei range può
+  // far collassare gli intervalli )
+   if ( abs( dMax - dMin) < EPS_SMALL)
+      return true ;
+
+  // Riporto le coordinate cicliche nell'ordine di partenza
+   if ( nGrid == 1) {
+      swap( vtNmi.x, vtNmi.z) ;
+      swap( vtNmi.y, vtNmi.z) ;
+      swap( vtNma.x, vtNma.z) ;
+      swap( vtNma.y, vtNma.z) ;
+   }
+   else if ( nGrid == 2) {
+      swap( vtNmi.y, vtNmi.z) ;
+      swap( vtNmi.x, vtNmi.z) ;
+      swap( vtNma.y, vtNma.z) ;
+      swap( vtNma.x, vtNma.z) ;
+   }
+
+   // Calcolo nPos
+   unsigned int nPos = nJ * m_nNx[nGrid] + nI ;
+   vector<Data>& vDexel = m_Values[nGrid][nPos] ;
+
+   bool bModified = false ;
+
+  // se non esistono segmenti...
+   if ( int( vDexel.size()) == 0) { // ... definisco il primo
+      vDexel.emplace_back() ;
+      vDexel.back().dMin = dMin ;
+      vDexel.back().vtMinN = vtNmi ;
+      vDexel.back().nToolMin = nToolNum ;
+      vDexel.back().dMax = dMax ;
+      vDexel.back().vtMaxN = vtNma ;
+      vDexel.back().nToolMax = nToolNum ;
+      m_OGrMgr.Reset() ;
+      bModified = true ; // modifica effettuata
+   }
+  // se ne esiste almeno uno, devo controllare possibili overlaps
+   else {
+     // Cerco l'ultimo intervallo a sinistra e l'ultimo intervallo a destra di quello da aggiungere,
+     // che non interferiscono con quest'ultimo
+      auto itLastLeft = vDexel.end() ;
+      auto itFirstRight = vDexel.end() ; // per bloccare la prima iterazione...
+      for ( auto it = vDexel.begin() ; it != vDexel.end() ; ++ it) {
+        // se il minimo attuale è maggiore del massimo it-esimo, l'intervallo it-esimo è a sinistra
+         if ( dMin > it->dMax + EPS_SMALL)
+            itLastLeft = it ; // aggiorno
+        // se il massimo attuale è minore del minimo it-esimo, l'interallo it-esimo è a destra
+         if ( dMax < it->dMin - EPS_SMALL  &&  itFirstRight == vDexel.end())
+            itFirstRight = it ; // aggiorno
+      }
+
+     // se esistono intervalli a sinistra... ( inizialmente l'iteratore è inizializzato su end())
+      if ( itLastLeft != vDexel.end()) {
+        // prendo il successivo dell'ultimo più vicino
+         auto itNextToLastLeft = itLastLeft ;
+         ++ itNextToLastLeft ;
+
+        // e se il successivo a sinistra non esiste...
+         if ( itNextToLastLeft == vDexel.end()) {
+           // aggiorno l'ultimo intervallo
+            vDexel.emplace_back() ;
+            vDexel.back().dMin = dMin ;
+            vDexel.back().dMax = dMax ;
+            vDexel.back().vtMinN = vtNmi ;
+            vDexel.back().vtMaxN = vtNma;
+            vDexel.back().nToolMin = nToolNum ;
+            vDexel.back().nToolMax = nToolNum ;
+            bModified = true ;
+         }
+        // se esiste...
+         else {
+           // controllo se coincide con il primo trovato a destra...
+            if ( itNextToLastLeft == itFirstRight) {
+              // in questo caso definisco il nuovo intervallo
+               Data NewSegment ;
+               NewSegment.dMin = dMin ;
+               NewSegment.dMax = dMax ;
+               NewSegment.vtMinN = vtNmi ;
+               NewSegment.vtMaxN = vtNma ;
+               NewSegment.nToolMin = nToolNum ;
+               NewSegment.nToolMax = nToolNum ;
+               vDexel.insert( itFirstRight, NewSegment) ; // inserimento intervallo
+               bModified = true ;
+            }
+            else {
+              // se il successivo non esce a sinistra da quello da aggiungere
+               if ( itNextToLastLeft->dMin > dMin + EPS_SMALL) {
+                  itNextToLastLeft->dMin = dMin ;
+                  itNextToLastLeft->vtMinN = vtNmi ;
+                  itNextToLastLeft->nToolMin = nToolNum ;
+               }
+              // cerco poi l'ultimo segmento che interferisce con quello da aggiungere
+               auto itPrevToFirstRight = vDexel.end() ;
+               for ( auto it = itNextToLastLeft ; it != itFirstRight ; ++ it)
+                  itPrevToFirstRight = it ;
+              // anche in qesto caso, se l'ultimo che interferisce non esce a destra da quello da aggiungere
+               if ( itPrevToFirstRight->dMax < dMax - EPS_SMALL) {
+                  itNextToLastLeft->dMax = dMax ;
+                  itNextToLastLeft->vtMaxN = vtNma ;
+                  itNextToLastLeft->nToolMax = nToolNum ;
+                  bModified = true ;
+               }
+               else {
+                  itNextToLastLeft->dMax = itPrevToFirstRight->dMax ;
+                  itNextToLastLeft->vtMaxN = itPrevToFirstRight->vtMaxN ;
+                  itNextToLastLeft->nToolMax = nToolNum ;
+                  bModified = true ;
+               }
+              // cancello il successivo...
+               auto itFirstToCancel = itNextToLastLeft ;
+               ++ itFirstToCancel ;
+               vDexel.erase( itFirstToCancel, itFirstRight) ;
+            }
+         }
+      }
+     // se non esiste da destra ...
+      else if ( itFirstRight == m_Values[nGrid][nPos].end()) {
+         // e il primo intervallo non sporge a sinistra ...
+         if ( vDexel.begin()->dMin > dMin + EPS_SMALL) {
+            vDexel.begin()->dMin = dMin ;
+            vDexel.begin()->vtMinN = vtNmi ;
+            vDexel.begin()->nToolMin = nToolNum ;
+            bModified = true ;
+         }
+        // se sporge da destra
+         if ( m_Values[nGrid][nPos].back().dMax > dMax + EPS_SMALL) {
+           // allora ci sono più segmenti, inglobo tutti nel primo
+            if ( vDexel.back().dMax > vDexel.begin()->dMax + EPS_SMALL) {
+               vDexel.begin()->dMax = vDexel.back().dMax ;
+               vDexel.begin()->vtMaxN = vDexel.back().vtMaxN ;
+               vDexel.begin()->nToolMax = nToolNum ;
+               bModified = true ;
+            }
+         }
+        // se l'ultimo intervallo non sporge a destra.
+         else {
+            vDexel.begin()->dMax = dMax ;
+            vDexel.begin()->vtMaxN = vtNma ;
+            vDexel.begin()->nToolMax = nToolNum ;
+            bModified = true ;
+         }
+        // cancello quelli inglobati
+         vDexel.erase( vDexel.begin() + 1, vDexel.end()) ;
+      }
+     // Casi Estremi...
+      else {
+        // 1) Tutti i segmenti sono a destra da quello da aggiungere...
+         if ( itFirstRight == vDexel.begin()) {
+           // -> inserisco il nuovo intervallo
+            Data NewSegment ;
+            NewSegment.dMin = dMin ;
+            NewSegment.dMax = dMax ;
+            NewSegment.vtMinN = vtNmi ;
+            NewSegment.vtMaxN = vtNma ;
+            NewSegment.nToolMin = nToolNum ;
+            NewSegment.nToolMax = nToolNum ;
+            vDexel.insert( vDexel.begin(), NewSegment) ; // posizione iniziale
+            bModified = true ;
+         }
+         else {
+           // se il primo segmento non esce a sinistra da quello da aggiungere, cambio l'inizio
+            if ( vDexel.begin()->dMin > dMin + EPS_SMALL) {
+               vDexel.begin()->dMin = dMin ;
+               vDexel.begin()->vtMinN = vtNmi ;
+               vDexel.begin()->nToolMin = nToolNum ;
+               bModified = true ;
+            }
+           // cerco l'ultimo segmento che interferisce con quello da aggiungere
+            auto itPrevToFirstRight = vDexel.begin() ;
+            for ( auto it = m_Values[nGrid][nPos].begin() ; it != itFirstRight ; ++ it)
+               itPrevToFirstRight = it ;
+           // se l'ultimo che interferisce non esce a destra da quello da aggiungere...
+            if ( itPrevToFirstRight->dMax < dMax - EPS_SMALL) {
+               vDexel.begin()->dMax = dMax ;
+               vDexel.begin()->vtMaxN = vtNma ;
+               vDexel.begin()->nToolMax = nToolNum ;
+               bModified = true ;
+            }
+           // ... altrimenti
+            else {
+               vDexel.begin()->dMax = itPrevToFirstRight->dMax ;
+               vDexel.begin()->vtMaxN = itPrevToFirstRight->vtMaxN ;
+               vDexel.begin()->nToolMax = nToolNum ;
+               bModified = true ;
+            }
+           // cancello gli inglobati ( come prima )
+            auto itFirstToCancel = vDexel.begin() ;
+            ++ itFirstToCancel ;
+            vDexel.erase( itFirstToCancel, itFirstRight) ;
+         }
+      }
+   }
+
+  // Se nessuna modifica, esco
+   if ( ! bModified)
+      return true ;
+
+  // Imposto ricalcolo della grafica
+   m_OGrMgr.Reset() ;
+  // Imposto forma generica
+   m_nShape = GENERIC ;
+  // Imposto ricalcolo numero di componenti connesse
+   m_nConnectedCompoCount = - 1 ;
+
+  // Passo da indici di dexel a indici di voxel
+   nI /= m_nDexVoxRatio ;
+   nJ /= m_nDexVoxRatio ;
+
+  // Determino quali blocchi sono stati modificati a seconda della griglia
+   if ( nGrid == 0) {
+      // Voxel lungo X
+      int nXStop = 1 ;
+      int nXBlock[2] ;
+      nXBlock[0] = min( nI / m_nVoxNumPerBlock, m_nFracLin[0] - 1) ;
+      if ( nI % N_VOXBLOCK == 0  &&  nXBlock[0] > 0) {
+         nXBlock[1] = nXBlock[0] - 1 ;
+         ++ nXStop ;
+      }
+     // Voxel lungo Y
+      int nYStop = 1 ;
+      int nYBlock[2] ;
+      nYBlock[0] = min( nJ / m_nVoxNumPerBlock, m_nFracLin[1] - 1) ;
+      if ( nJ % N_VOXBLOCK == 0  &&  nYBlock[0] > 0) {
+         nYBlock[1] = nYBlock[0] - 1 ;
+         ++ nYStop ;
+      }
+     // Voxel lungo Z
+      int nVoxNumZ = int( m_nNy[1] / m_nDexVoxRatio + ( m_nNy[1] % m_nDexVoxRatio == 0 ? 1 : 2)) ;
+      int nMinK = Clamp( int( floor( ( ( dMin - 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) - EPS_SMALL))), 0, nVoxNumZ - 2) ;
+      int nMaxK = Clamp( int( floor( ( ( dMax + 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) + EPS_SMALL))), 0, nVoxNumZ - 2) ;
+      int nMinZBlock = ( m_nMapNum == 1 ? 0 : Clamp( nMinK / int( m_nVoxNumPerBlock), 0, int( m_nFracLin[2] - 1))) ;
+      int nMaxZBlock = min( int( m_nFracLin[2] - 1), nMaxK / int( m_nVoxNumPerBlock)) ;
+     // Assegno flag ai voxel
+      for ( int tI = 0 ; tI < nXStop ; ++ tI) {
+         for ( int tJ = 0 ; tJ < nYStop ; ++ tJ) {
+            for ( int k = nMinZBlock ; k <= nMaxZBlock ; ++ k) {
+               int nBlockNum = k * m_nFracLin[0] * m_nFracLin[1] + nYBlock[tJ] * m_nFracLin[0] + nXBlock[tI] ;
+               m_BlockToUpdate[nBlockNum] = true ;
+            }
+         }
+      }
+   }
+
+   else if ( nGrid == 1) {
+     // Voxel lungo Y
+      int nYStop = 1 ;
+      int nYBlock[2] ;
+      nYBlock[0] = min( nI / m_nVoxNumPerBlock, m_nFracLin[1] - 1) ;
+      if ( nI % N_VOXBLOCK == 0  &&  nYBlock[0] > 0) {
+         nYBlock[1] = nYBlock[0] - 1 ;
+         ++ nYStop ;
+      }
+     // Voxel lungo Z
+      int nZStop = 1 ;
+      int nZBlock[2] ;
+      nZBlock[0] = min( nJ / m_nVoxNumPerBlock, m_nFracLin[2] - 1) ;
+      if ( nJ % N_VOXBLOCK == 0  &&  nZBlock[0] > 0) {
+         nZBlock[1] = nZBlock[0] - 1 ;
+         ++ nZStop ;
+      }
+      // Voxel lungo X
+      int nVoxNumX = int( m_nNx[0] / m_nDexVoxRatio + ( m_nNx[0] % m_nDexVoxRatio == 0 ? 1 : 2)) ;
+      int nMinI = Clamp( int( floor( ( ( dMin - 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) - EPS_SMALL))), 0, nVoxNumX - 2) ;
+      int nMaxI = Clamp( int( floor( ( ( dMax + 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) + EPS_SMALL))), 0, nVoxNumX - 2) ;
+      int nMinXBlock = Clamp( nMinI / int( m_nVoxNumPerBlock), 0, int( m_nFracLin[0] - 1)) ;
+      int nMaxXBlock = min( int( m_nFracLin[0] - 1), nMaxI / int( m_nVoxNumPerBlock)) ;
+      // Assegno flag ai voxel
+      for ( int tI = 0 ; tI < nYStop ; ++ tI) {
+         for ( int tJ = 0 ; tJ < nZStop ; ++ tJ) {
+            for ( int k = nMinXBlock ; k <= nMaxXBlock ; ++ k) {
+               int nBlockNum = nZBlock[tJ] * m_nFracLin[0] * m_nFracLin[1] + nYBlock[tI] * m_nFracLin[0] + k ;
+               m_BlockToUpdate[nBlockNum] = true ;
+            }
+         }
+      }
+   }
+
+   else if ( nGrid == 2) {
+      // Voxel lungo X
+      int nXStop = 1 ;
+      int nXBlock[2] ;
+      nXBlock[0] = min( nJ / m_nVoxNumPerBlock, m_nFracLin[0] - 1) ;
+      if ( nJ % N_VOXBLOCK == 0  &&  nXBlock[0] > 0) {
+         nXBlock[1] = nXBlock[0] - 1 ;
+         ++ nXStop ;
+      }
+      // Voxel lungo Z
+      int nZStop = 1 ;
+      int nZBlock[2] ;
+      nZBlock[0] = min( nI / m_nVoxNumPerBlock, m_nFracLin[2] - 1) ;
+      if ( nI % N_VOXBLOCK == 0  &&  nZBlock[0] > 0) {
+         nZBlock[1] = nZBlock[0] - 1 ;
+         ++ nZStop ;
+      }
+      // Voxel lungo Y 
+      int nVoxNumY = int( m_nNy[0] / m_nDexVoxRatio + ( m_nNy[0] % m_nDexVoxRatio == 0 ? 1 : 2)) ; 
+      int nMinJ = Clamp( int( floor( ( ( dMin - 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) - EPS_SMALL))), 0, nVoxNumY - 2) ; 
+      int nMaxJ = Clamp( int( floor( ( ( dMax + 0.5 * m_dStep) / ( m_nDexVoxRatio * m_dStep) + EPS_SMALL))), 0, nVoxNumY - 2) ; 
+      int nMinYBlock = Clamp( nMinJ / int( m_nVoxNumPerBlock), 0, int( m_nFracLin[1] - 1)) ;
+      int nMaxYBlock = min( int( m_nFracLin[1] - 1), nMaxJ / int( m_nVoxNumPerBlock)) ;
+      // Assegno flag ai voxel
+      for ( int tI = 0 ; tI < nZStop ; ++ tI) {
+         for ( int tJ = 0 ; tJ < nXStop ; ++ tJ) {
+            for ( int k = nMinYBlock ; k <= nMaxYBlock ; ++ k) {
+               int nBlockNum = nZBlock[tI] * m_nFracLin[0] * m_nFracLin[1] + k * m_nFracLin[0] + nXBlock[tJ] ;
+               m_BlockToUpdate[nBlockNum] = true ;
+            }
+         }
+      }
+   }
+
    return true ;
 }
 

Voronoi.cpp

@@ -414,6 +414,11 @@ Voronoi::GetBisectorCurve( int i)
       CurveLine* pLine = CreateBasicCurveLine() ;
       if ( pLine == nullptr)
          return nullptr ;
+     // costruisco il bisettore orientato dal parametro minore al maggiore 
+      if ( dParS > dParE) {
+         swap( ptS, ptE) ;
+         swap( dParS, dParE) ;
+      }
       pLine->Set( ptS, ptE) ;
       pLine->SetTempParam( dParS, 0) ;
       pLine->SetTempParam( dParE, 1) ;
@@ -433,7 +438,7 @@ Voronoi::GetBisectorCurve( int i)
          return nullptr ;
       pArc->SetC2P( ptC, ptS, ptE) ;
       pArc->SetTempParam( dParS, 0) ;
-      pArc->SetTempParam( dParS, 1) ;
+      pArc->SetTempParam( dParS, 1) ; // dParE = dParS
       pArc->ToGlob( m_Frame) ;
       return pArc ;   
    }
@@ -449,15 +454,23 @@ Voronoi::GetBisectorCurve( int i)
       if ( pCompo == nullptr)
          return nullptr ;
    
+     // verifico se devo leggere i punti del bisettore al contrario per averlo orientato dal parametro minore al maggiore
+      bool bInvert = false ;
+      double dPar1, dPar2 ;
+      m_vroni->GetApproxedBisectorParams( i, dPar1, dPar2) ;
+      if ( dPar1 > dPar2 + EPS_SMALL)
+         bInvert = true ;
+
      // punto iniziale
       Point3d pt ;
       double dParS ;
-      m_vroni->GetApproxedBisectorPoint( i, 0, pt.v, dParS) ;
+      m_vroni->GetApproxedBisectorPoint( i, bInvert ? nPoints - 1 : 0, pt.v, dParS) ;
       pCompo->AddPoint( pt) ;
       int nCrvCount = 0 ;
    
       double dParPrev = dParS ;
-      for ( int j = 1 ; j < nPoints ; j ++) {
+      int j = bInvert ? nPoints - 2 : 1 ;
+      while ( ( bInvert  &&  j >= 0)  ||  ( ! bInvert  &&  j < nPoints)) {
          double dPar ;
          m_vroni->GetApproxedBisectorPoint( i, j, pt.v, dPar) ;
          if ( pCompo->AddLine( pt)) {           
@@ -467,7 +480,12 @@ Voronoi::GetBisectorCurve( int i)
            // aggiorno parametro precedente
             dParPrev = dPar ;
             nCrvCount ++ ;
-         }     
+         }
+        // aggiorno per punto successivo
+         if ( bInvert)
+            j -- ;
+         else
+            j ++ ;
       }
    
      // setto parametri sulla curva