EgtGeomKernel :

- raffinamento finale della intersCompoCompo.

CDeClosedSurfTmClosedSurfTm.cpp

@@ -1,12 +1,13 @@
 //----------------------------------------------------------------------------
-//  EgalTech 2020-2024
+//  EgalTech 2020-2026
 //----------------------------------------------------------------------------
-// File : CDeClosedSurfTmClosedSurfTm.h    Data : 24.03.24  Versione : 2.6c2 
+// File : CDeClosedSurfTmClosedSurfTm.h    Data : 23.01.26  Versione : 3.1a3 
 // Contenuto : Implementazione funzione verifica collisione tra
 //             SurfTm e SurfTm.
 //
 // Modifiche : 13.11.20 LM Creazione modulo.
 //             24.03.24 DS Aggiunta TestSurfTmSurfTm.
+//             23.01.26 DS In TestSurfTmSurfTm aggiunto flag per collisione quando una delle due è chiusa e contiene l'altra.
 //
 //----------------------------------------------------------------------------
 
@@ -45,10 +46,10 @@ CDeClosedSurfTmClosedSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& Surf
    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.
+  // 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.
+  // 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
@@ -61,13 +62,13 @@ CDeClosedSurfTmClosedSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& Surf
          continue ;
       BBox3d b3BoxTriaA ;
       trTriaA.GetLocalBBox( b3BoxTriaA) ;
-     // Se è necessario, espando il box di una costante additiva pari alla distanza di sicurezza.
+     // 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.
+     // 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.
@@ -84,14 +85,14 @@ CDeClosedSurfTmClosedSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& Surf
            // 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.
+              // 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.
+              // 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.
+              // 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 ;
@@ -111,10 +112,10 @@ CDeClosedSurfTmClosedSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& Surf
       }
    }
   // Non ho trovato collisioni fra triangoli delle superfici.
-  // Se il BBox della prima superficie non è interno a quello della seconda e viceversa, non c'è collisione
+  // Se il BBox della prima superficie non è interno a quello della seconda e viceversa, non c'è collisione
    if ( ! b3BoxA.Encloses( b3BoxB)  &&  ! b3BoxB.Encloses( b3BoxA))
       return false ;
-  // La collisione c'è se una superficie è dentro l'altra.
+  // La collisione c'è se una superficie è dentro l'altra.
    Point3d ptPointA, ptPointB ;
    pSrfA->GetFirstVertex( ptPointA) ;
    pSrfB->GetFirstVertex( ptPointB) ;
@@ -127,7 +128,7 @@ CDeClosedSurfTmClosedSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& Surf
 // Verifica l'interferenza tra le due superfici : restituisce true in caso di interferenza.
 //----------------------------------------------------------------------------
 bool
-TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double dSafeDist)
+TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double dSafeDist, bool bTestEnclosion)
 {
   // Recupero le superfici base
    const SurfTriMesh* pSrfA = GetBasicSurfTriMesh( &SurfA) ;
@@ -140,10 +141,10 @@ TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double d
    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.
+  // 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.
+  // 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
@@ -156,13 +157,13 @@ TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double d
          continue ;
       BBox3d b3BoxTriaA ;
       trTriaA.GetLocalBBox( b3BoxTriaA) ;
-     // Se è necessario, espando il box di una costante additiva pari alla distanza di sicurezza.
+     // 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.
+     // 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.
@@ -179,14 +180,14 @@ TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double d
            // 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.
+              // 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.
+              // 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.
+              // 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 ;
@@ -205,6 +206,25 @@ TestSurfTmSurfTm( const ISurfTriMesh& SurfA, const ISurfTriMesh& SurfB, double d
          pSrfB->SetTempInt( nTB, 1) ;
       }
    }
-  // Non c'è interferenza
+  // Se non richiesto test di inclusione, non c'è interferenza
+   if ( ! bTestEnclosion)
+      return false ;
+  // Se la prima superficie è chiusa, verifico se include totalmente la seconda
+   if ( pSrfA->IsClosed()  &&  b3BoxA.Encloses( b3BoxB)) {
+      Point3d ptPointB ;
+      pSrfB->GetFirstVertex( ptPointB) ;
+      DistPointSurfTm DistPoinBSrfA( ptPointB, *pSrfA) ;
+      if ( DistPoinBSrfA.IsPointInside()  ||  DistPoinBSrfA.IsEpsilon( dSafeDist))
+         return true ;
+   }
+  // Se la seconda superficie è chiusa, verifico se include totalmente la prima
+   if ( pSrfB->IsClosed()  &&  b3BoxB.Encloses( b3BoxA)) {
+      Point3d ptPointA ;
+      pSrfA->GetFirstVertex( ptPointA) ;
+      DistPointSurfTm DistPoinASrfB( ptPointA, *pSrfB) ;
+      if ( DistPoinASrfB.IsPointInside()  ||  DistPoinASrfB.IsEpsilon( dSafeDist))
+         return true ;
+   }
+  // Non c'è interferenza
    return false ;
 }

CalcDerivate.cpp

@@ -0,0 +1,233 @@
+//----------------------------------------------------------------------------
+//  EgalTech 2026
+//----------------------------------------------------------------------------
+// File : CalcDerivate.cpp   Data : 03.02.26       Versione : 1.5h1
+// Contenuto : Funzioni per calcolo derivate secondo Bessel e Akima.
+//
+//
+//
+// Modifiche : 03.02.26 DB Creazione modulo.
+//
+//
+//----------------------------------------------------------------------------
+
+#pragma once
+
+#include "stdafx.h"
+#include "CalcDerivate.h"
+#include "/EgtDev/Include/EGkPoint3d.h"
+#include "/EgtDev/Include/EgtNumCollection.h"
+#include "/EgtDev/Include/EGkGeoCollection.h"
+
+//----------------------------------------------------------------------------
+bool
+ComputeAkimaTangents( bool bDetectCorner, const DBLVECTOR& vPar, const PNTVECTOR& vPnt, VCT3DVECTOR& vPrevDer, VCT3DVECTOR& vNextDer)
+{
+  // pulisco i vettori dei parametri e delle tangenti
+   vPrevDer.clear() ;
+   vNextDer.clear() ;
+
+  // numero di punti
+   int nSize = int( vPnt.size()) ;
+
+  // sono necessari almeno due punti
+   if ( nSize < 2)
+      return false ;
+
+  // calcolo le derivate
+   vPrevDer.reserve( nSize) ;
+   vNextDer.reserve( nSize) ;
+  // se ci sono solo 2 punti, le tangenti devono essere dirette lungo la linea che li unisce
+   if ( nSize == 2) {
+     // non esiste derivata prima del primo punto
+      vPrevDer.emplace_back( 0, 0, 0) ;
+      vNextDer.push_back( ( vPnt[1] - vPnt[0]) / ( vPar[1] - vPar[0])) ;
+      vPrevDer.push_back( vNextDer[0]) ;
+     // non esiste derivata dopo il secondo e ultimo punto
+      vNextDer.emplace_back( 0, 0, 0) ;
+      return true ;
+   }
+  // verifico se curva chiusa (primo e ultimo punto coincidono)
+   bool bClosed = AreSamePointApprox( vPnt.front(), vPnt.back()) ;
+  // calcolo le derivate
+   for ( int i = 0 ; i < nSize ; ++ i) {
+      Vector3d vtPrevDer ;
+      Vector3d vtNextDer ;
+     // primo punto
+      if ( i == 0) {
+        // se curva chiusa, come precedente uso il penultimo punto
+         if ( bClosed) {
+           // se non ci sono almeno 5 punti
+            if ( nSize < 5) {
+               if ( ! CalcCircleMidDer( vPar[nSize-2] - vPar[nSize-1], vPnt[nSize-2], vPar[i], vPnt[i],
+                                        vPar[i+1], vPnt[i+1], vtNextDer))
+                  return false ;
+               vtPrevDer = vtNextDer ;
+            }
+           // altrimenti
+            else {
+               if ( ! CalcAkimaMidDer( vPar[nSize-3] - vPar[nSize-1], vPnt[nSize-3], vPar[nSize-2] - vPar[nSize-1], vPnt[nSize-2],
+                                       vPar[i], vPnt[i], vPar[i+1], vPnt[i+1],
+                                       vPar[i+2], vPnt[i+2], bDetectCorner,
+                                       vtPrevDer, vtNextDer))
+                  return false ;
+            }
+         }
+        // altrimenti, uso arco sui primi tre punti
+         else {
+            if ( ! CalcCircleStartDer( vPar[i], vPnt[i], vPar[i+1], vPnt[i+1],
+                                       vPar[i+2], vPnt[i+2], vtNextDer))
+               return false ;
+            vtPrevDer = Vector3d( 0, 0, 0) ;
+         }
+      }
+     // ultimo punto
+      else if ( i == nSize - 1) {
+        // se curva chiusa, le tg devono coincidere con quelle del primo
+         if ( bClosed) {
+            vtPrevDer = vPrevDer[0] ;
+            vtNextDer = vNextDer[0] ;
+         }
+        // altrimenti, uso arco sugli ultimi tre punti
+         else {
+            if ( ! CalcCircleEndDer( vPar[i-2], vPnt[i-2], vPar[i-1], vPnt[i-1],
+                                       vPar[i], vPnt[i], vtPrevDer))
+               return false ;
+            vtNextDer = Vector3d( 0, 0, 0) ;
+         }
+      }
+     // punti intermedi
+      else {
+        // se secondo punto
+         if ( i == 1) {
+           // se curva aperta o non ci sono almeno 5 punti
+            if ( ! bClosed  ||  nSize < 5) {
+               if ( ! CalcCircleMidDer( vPar[i-1], vPnt[i-1], vPar[i], vPnt[i],
+                                        vPar[i+1], vPnt[i+1], vtPrevDer))
+                  return false ;
+               vtNextDer = vtPrevDer ;
+            }
+           // altrimenti
+            else {
+               if ( ! CalcAkimaMidDer( vPar[nSize-2] - vPar[nSize-1], vPnt[nSize-2], vPar[i-1], vPnt[i-1],
+                                       vPar[i], vPnt[i], vPar[i+1], vPnt[i+1],
+                                       vPar[i+2], vPnt[i+2], bDetectCorner,
+                                       vtPrevDer, vtNextDer))
+                  return false ;
+            }
+         }
+        // se penultimo punto
+         else if ( i == nSize - 2) {
+           // se curva aperta o non ci sono almeno 5 punti
+            if ( ! bClosed  ||  nSize < 5) {
+               if ( ! CalcCircleMidDer( vPar[i-1], vPnt[i-1], vPar[i], vPnt[i],
+                                        vPar[i+1], vPnt[i+1], vtPrevDer))
+                  return false ;
+               vtNextDer = vtPrevDer ;
+            }
+           // altrimenti
+            else {
+               if ( ! CalcAkimaMidDer( vPar[i-2], vPnt[i-2], vPar[i-1], vPnt[i-1],
+                                       vPar[i], vPnt[i], vPar[i+1], vPnt[i+1],
+                                       vPar[1] + vPar[i+1], vPnt[1], bDetectCorner,
+                                       vtPrevDer, vtNextDer))
+                  return false ;
+            }
+         }
+        // altrimenti
+         else {
+            if ( ! CalcAkimaMidDer( vPar[i-2], vPnt[i-2], vPar[i-1], vPnt[i-1],
+                                    vPar[i], vPnt[i], vPar[i+1], vPnt[i+1],
+                                    vPar[i+2], vPnt[i+2], bDetectCorner,
+                                    vtPrevDer, vtNextDer))
+               return false ;
+         }
+      }
+     // salvo la derivata
+      vPrevDer.push_back( vtPrevDer) ;
+      vNextDer.push_back( vtNextDer) ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+bool
+ComputeBesselTangents( const DBLVECTOR& vPar, const PNTVECTOR& vPnt, VCT3DVECTOR& vPrevDer, VCT3DVECTOR& vNextDer)
+{
+  // pulisco i vettori dei parametri e delle tangenti
+   vPrevDer.clear() ;
+   vNextDer.clear() ;
+
+  // numero di punti
+   int nSize = int( vPnt.size()) ;
+
+  // sono necessari almeno due punti
+   if ( nSize < 2)
+      return false ;
+
+  // calcolo le derivate
+   vPrevDer.reserve( nSize) ;
+   vNextDer.reserve( nSize) ;
+  // se ci sono solo 2 punti, le tangenti devono essere dirette lungo la linea che li unisce
+   if ( nSize == 2) {
+     // non esiste derivata prima del primo punto
+      vPrevDer.emplace_back( 0, 0, 0) ;
+      vNextDer.push_back( ( vPnt[1] - vPnt[0]) / ( vPar[1] - vPar[0])) ;
+      vPrevDer.push_back( vNextDer[0]) ;
+     // non esiste derivata dopo il secondo e ultimo punto
+      vNextDer.emplace_back( 0, 0, 0) ;
+      return true ;
+   }
+  // verifico se curva chiusa (primo e ultimo punto coincidono)
+   bool bClosed = AreSamePointApprox( vPnt.front(), vPnt.back()) ;
+  // calcolo le derivate
+   for ( int i = 0 ; i < nSize ; ++ i) {
+      Vector3d vtPrevDer ;
+      Vector3d vtNextDer ;
+     // primo punto
+      if ( i == 0) {
+        // se curva chiusa, come precedente uso il penultimo punto
+         if ( bClosed) {
+            if ( ! CalcBesselMidDer( vPar[nSize-2] - vPar[nSize-1], vPnt[nSize-2], vPar[i], vPnt[i],
+                                     vPar[i+1], vPnt[i+1], vtNextDer))
+               return false ;
+            vtPrevDer = vtNextDer ;
+         }
+        // altrimenti, uso i primi tre punti
+         else {
+            if ( ! CalcBesselStartDer( vPar[i], vPnt[i], vPar[i+1], vPnt[i+1],
+                                       vPar[i+2], vPnt[i+2], vtNextDer))
+               return false ;
+            vtPrevDer = Vector3d( 0, 0, 0) ;
+         }
+      }
+     // ultimo punto
+      else if ( i == nSize - 1) {
+        // se curva chiusa, le tg devono coincidere con quelle del primo
+         if ( bClosed) {
+            vtPrevDer = vPrevDer[0] ;
+            vtNextDer = vNextDer[0] ;
+         }
+        // altrimenti, uso gli ultimi tre punti
+         else {
+            if ( ! CalcBesselEndDer( vPar[i-2], vPnt[i-2], vPar[i-1], vPnt[i-1],
+                                     vPar[i], vPnt[i], vtPrevDer))
+               return false ;
+            vtNextDer = Vector3d( 0, 0, 0) ;
+         }
+      }
+     // punti intermedi
+      else {
+         if ( ! CalcBesselMidDer( vPar[i-1], vPnt[i-1], vPar[i], vPnt[i],
+                                  vPar[i+1], vPnt[i+1], vtPrevDer))
+            return false ;
+         vtNextDer = vtPrevDer ;
+      }
+     // salvo la derivata
+      vPrevDer.push_back( vtPrevDer) ;
+      vNextDer.push_back( vtNextDer) ;
+   }
+
+   return true ;
+}

CalcDerivate.h

@@ -14,6 +14,8 @@
 #pragma once
 
 #include "/EgtDev/Include/EGkPoint3d.h"
+#include "/EgtDev/Include/EgtNumCollection.h"
+#include "/EgtDev/Include/EGkGeoCollection.h"
 
 
 //----------------------------------------------------------------------------
@@ -179,3 +181,6 @@ CalcAkimaMidDer( double dU0, const Point3d& ptP0, double dU1, const Point3d& ptP
    }
    return ( ! vtPrevDer.IsZero()  &&  ! vtNextDer.IsZero()) ;
 }
+
+bool ComputeAkimaTangents( bool bDetectCorner, const DBLVECTOR& vPar, const PNTVECTOR& vPnt, VCT3DVECTOR& vPrevDer, VCT3DVECTOR& vNextDer) ;
+bool ComputeBesselTangents( const DBLVECTOR& vPar, const PNTVECTOR& vPnt, VCT3DVECTOR& vPrevDer, VCT3DVECTOR& vNextDer) ;

CalcPocketing.cpp

@@ -722,7 +722,7 @@ AssignFeedSpiral( ICurveComposite* pCrv, const ISurfFlatRegion* pSrfRemoved, boo
 
 //----------------------------------------------------------------------------
 static bool
-AssignFeedSpiralOpt( const int nOptType, const PocketParams& PockParams, ICurveComposite* pCrv )
+AssignFeedSpiralOpt( int nOptType, const PocketParams& PockParams, ICurveComposite* pCrv )
 {
   // controllo della curva corrente
    if ( pCrv == nullptr  ||  ! pCrv->IsValid()  ||  pCrv->GetCurveCount() == 0)
@@ -732,24 +732,24 @@ AssignFeedSpiralOpt( const int nOptType, const PocketParams& PockParams, ICurveC
    if ( ! PockParams.bCalcFeed)
       return AssignMaxFeed( pCrv, PockParams) ;
 
-   switch ( PockParams.nType) {
-      case POCKET_SPIRALIN :
-         if ( nOptType == 0) { // Spirale dall'Esterno
-            for ( int u = 0 ; u < pCrv->GetCurveCount() ; ++ u) {
-               if ( u == 0) // prima circonferenza
-                  pCrv->SetCurveTempParam( 0, GetMinFeed( PockParams), 0) ;
-               else // semi cerchi in tangenza
-                  pCrv->SetCurveTempParam( u, GetMaxFeed( PockParams), 0) ;
-            }
+   if ( PockParams.nType == POCKET_SPIRALIN  ||  PockParams.nType == POCKET_CONFORMAL_ZIGZAG  ||
+        PockParams.nType == POCKET_CONFORMAL_ONEWAY) {
+      if ( nOptType == 0) { // Spirale dall'Esterno
+         for ( int u = 0 ; u < pCrv->GetCurveCount() ; ++ u) {
+            if ( u == 0) // prima circonferenza
+               pCrv->SetCurveTempParam( 0, GetMinFeed( PockParams), 0) ;
+            else // semi cerchi in tangenza
+               pCrv->SetCurveTempParam( u, GetMaxFeed( PockParams), 0) ;
          }
-         else if ( nOptType == 1) { // Trapezoidi
-            for ( int u = 0 ; u < pCrv->GetCurveCount() ; ++ u)
-               pCrv->SetCurveTempParam( u, GetMinFeed( PockParams), 0) ;
-         }
-      break ;
+      }
+      else if ( nOptType == 1) { // Trapezoidi
+         for ( int u = 0 ; u < pCrv->GetCurveCount() ; ++ u)
+            pCrv->SetCurveTempParam( u, GetMinFeed( PockParams), 0) ;
+      }
+   }
    /* NB. Essendo la funzione CalcSpiral richiamata sia per lo SpiralIN che per lo SpiralOUT le curve sono sempre
    orientate nello stesso modo, solamente alla fine viene invertita la curva finale per la svuotatura... */
-   case POCKET_SPIRALOUT :
+   else {
       if ( nOptType == 0) { // Spiral verso l'esterno
          for ( int u = 0 ; u < pCrv->GetCurveCount() ; ++ u) {
             if ( u > pCrv->GetCurveCount() - 3) // prime semi circonferenze
@@ -762,9 +762,6 @@ AssignFeedSpiralOpt( const int nOptType, const PocketParams& PockParams, ICurveC
          for ( int u = 0 ; u < pCrv->GetCurveCount() ; ++ u)
             pCrv->SetCurveTempParam( u, GetMinFeed( PockParams), 0) ;
       }
-      break ;
-   default :
-      break ;
    }
 
    return true ;
@@ -1061,9 +1058,9 @@ ExtendPath( ICurveComposite* pCompo, const ISurfFlatRegion* pSfr, const PocketPa
       bool bInside = true ;
       bool bOkOut = true ;
       if ( PockParams.SfrLimit.IsValid()) {
-         bOkOut = ( IsPointInsideSurfFr( ptFall, &PockParams.SfrLimit, dMinDist, bInside)  &&  ! bInside) ;
+         bOkOut = ( IsPointInsideSurfFr( ptFall, &PockParams.SfrLimit, dMinDist - 10. * EPS_SMALL, bInside)  &&  ! bInside) ;
          if ( bOkOut)
-            bOkOut = ( IsPointInsideSurfFr( Media( ptFall, pt), &PockParams.SfrLimit, dMinDist, bInside)  &&  ! bInside) ;
+            bOkOut = ( IsPointInsideSurfFr( Media( ptFall, pt), &PockParams.SfrLimit, dMinDist - 10. * EPS_SMALL, bInside)  &&  ! bInside) ;
       }
       if ( bOkOut)
          bOkOut = ( IsPointInsideSurfFr( ptFall, pSfr, dMinDist, bInside)  &&  ! bInside) ;
@@ -2638,6 +2635,87 @@ PreparareTrapezoidTwoBases( const ICurveComposite* pCrvCompo, const double dDiam
    return true ;
 }
 
+//----------------------------------------------------------------------------
+static bool
+IsForcedStepTrapezoid( const ICurveComposite* pCrvTrap, const PocketParams& PockParam,
+                       int nBase, int nSecondBase, bool& bForced)
+{
+   bForced = false ;
+  // se la curva non è valida, allora non può essere forzato
+   if ( pCrvTrap == nullptr  ||  ! pCrvTrap->IsValid())
+      return false ;
+
+  // scorro la curva e ricavo le TempProps
+   array<int, 4> vnProps ;
+   int nClose = 0 ;
+   for ( int i = 0 ; i < 4 ; ++ i) {
+      if ( ! pCrvTrap->GetCurveTempProp( i, vnProps[i], 0))
+         return false ;
+      if ( vnProps[i] == TEMP_PROP_CLOSE_EDGE)
+         ++ nClose ;
+   }
+
+   double dDiam = 2. * PockParam.dRad ;
+   switch ( nClose) {
+  // se trapezio tutto aperto, allora non è forzato
+   case 0 :
+      bForced = false ;
+   break ;
+  // se ho un lato chiuso, non è forzato
+   case 1 :
+      bForced = false ;
+   break ;
+  // se ho due lati chiusi
+   case 2 : {
+      if ( nBase < 0  ||  nBase > 4  ||  nSecondBase < 0  ||  nSecondBase > 4)
+         return false ;
+     // se entrambe le basi sono chiuse, è forzato
+      if ( vnProps[nBase] == TEMP_PROP_CLOSE_EDGE  &&  vnProps[nSecondBase] == TEMP_PROP_CLOSE_EDGE)
+         bForced = true ;
+     // se entrambe le basi sono aperte
+      else if ( vnProps[nBase] == TEMP_PROP_OPEN_EDGE  &&  vnProps[nSecondBase] == TEMP_PROP_OPEN_EDGE) {
+         const ICurve* pCrvOpenBase = pCrvTrap->GetCurve( nBase) ;
+         const ICurve* pCrvOpenSecondBase = pCrvTrap->GetCurve( nSecondBase) ;
+         if ( pCrvOpenBase == nullptr  ||  ! pCrvOpenBase->IsValid()  ||
+              pCrvOpenSecondBase == nullptr  ||  ! pCrvOpenSecondBase->IsValid())
+            return false ;
+         double dLenOpenBase ; pCrvOpenBase->GetLength( dLenOpenBase) ;
+         double dLenSecondOpenBase ; pCrvOpenSecondBase->GetLength( dLenSecondOpenBase) ;
+         bForced = ( dLenOpenBase < dDiam + 10. * EPS_SMALL  &&
+                     dLenSecondOpenBase < dDiam + 10. * EPS_SMALL) ;
+      }
+     // se alternate, non forzo
+      else
+         bForced = false ;
+   }
+   break ;
+  // se ho tre lati chiusi
+   case 3 : {
+     // diventa forzato se il lato aperto non è grande
+      double dLenOpen = 0. ;
+      for ( int i = 0 ; i < 4 ; ++ i) {
+         if ( vnProps[i] == TEMP_PROP_OPEN_EDGE) {
+            const ICurve* pCrvOpen = pCrvTrap->GetCurve( i) ;
+            if ( pCrvOpen == nullptr  ||  ! pCrvOpen->IsValid())
+               return false ;
+            pCrvOpen->GetLength( dLenOpen) ;
+            break ;
+         }
+      }
+      bForced = ( dLenOpen < dDiam + 10. * EPS_SMALL) ;
+   }
+   break ;
+  // se tutto chiuso, è forzato
+   case 4 :
+      bForced = true ;
+   break ;
+   default :
+      return false ;
+   }
+
+   return true ;
+}
+
 //----------------------------------------------------------------------------
 static bool
 GetTrapezoidFromShape( const ICurveComposite* pCrvCompo, ICurveComposite* pCrvTrap,
@@ -3040,8 +3118,12 @@ GetTrapezoidFromShape( const ICurveComposite* pCrvCompo, ICurveComposite* pCrvTr
       return false ;
 
   // se parametro MaxOptSize non compatibile => non è ottimizzato
-   if ( PockParams.dMaxOptSize > EPS_SMALL  &&  dPocketSize > PockParams.dMaxOptSize)
-      pCrvTrap->Clear() ;
+   bool bForced = false ;
+   IsForcedStepTrapezoid( pCrvTrap, PockParams, nBase, nSecondBase, bForced) ;
+   if ( ! bForced) {
+      if ( PockParams.dMaxOptSize > EPS_SMALL  &&  dPocketSize > PockParams.dMaxOptSize)
+         pCrvTrap->Clear() ;
+   }
 
    return true ;
 }
@@ -3444,18 +3526,18 @@ CalcTrapezoidSpiral( ICurveComposite* pCrvPocket, const Frame3d& frTrap, double
       Vector3d vtX ; pCrvPocket->GetCurve( 1)->GetStartDir( vtX) ;
       Frame3d frDim ; frDim.Set( ptOrig, Z_AX, vtX) ; frDim.Invert() ;
       BBox3d b3Dim ; pCrvPocket->GetBBox( frDim, b3Dim, BBF_EXACT) ;
-      dMaxLarg = ( vnProp[0] != 0 ? b3Dim.GetDimY() : dPocketSize) ;
+      dMaxLarg = ( vnProp[0] != TEMP_PROP_CLOSE_EDGE ? b3Dim.GetDimY() : dPocketSize) ;
    }
 
   // calcolo percorso di svuotatura 
   // se lati obliqui sono entrambi chiusi e dimensione svuotatura è maggiore di diametro fresa e minore del doppio gestione speciale
-   if ( ( bRealTrap  &&  dMaxLarg > PockParams.dRad * 2 + 10 * EPS_SMALL)  &&
+   if ( ( bRealTrap  &&  dMaxLarg > dDiam + 10. * EPS_SMALL)  &&
         ( ( ( vnProp[0] != TEMP_PROP_CLOSE_EDGE  &&  vnProp[2] != TEMP_PROP_CLOSE_EDGE)  &&
             ( vnProp[3] == TEMP_PROP_CLOSE_EDGE  &&  vnProp[1] == TEMP_PROP_CLOSE_EDGE)  &&
-            ( max( dLen0, dLen2) < 2 * dDiam + EPS_SMALL))  ||
+            ( max( dLen0, dLen2) < 2. * dDiam + EPS_SMALL))  ||
           ( ( vnProp[1] != TEMP_PROP_CLOSE_EDGE  &&  vnProp[3] != TEMP_PROP_CLOSE_EDGE)  &&
             ( vnProp[0] == TEMP_PROP_CLOSE_EDGE  &&  vnProp[2] == TEMP_PROP_CLOSE_EDGE)  &&
-            ( max( dLen1, dLen3) < 2 * dDiam + EPS_SMALL)))) {
+            ( max( dLen0, dLen2) < 2. * dDiam + EPS_SMALL)))) {
       if ( ! SpecialAdjustTrapezoidSpiralForAngles( pMCrv, ( vnProp[0] == TEMP_PROP_CLOSE_EDGE), pCrvPocket, PockParams, ptRef)) {
          pMCrv->Clear() ;
          return false ;
@@ -3608,14 +3690,8 @@ CalcTrapezoidSpiral( ICurveComposite* pCrvPocket, const Frame3d& frTrap, double
 
    if ( pMCrv->GetCurveCount() == 0)
       return true ;
-
    pMCrv->ToGlob( frTrap) ;
 
-   if ( PockParams.bInvert) {
-      pMCrv->Invert() ;
-     // inverto le proprietà in modo che nProp3 sia sempre legata al punto iniziale e nProp1 a quello finale
-      swap( vnProp[1], vnProp[3]) ;
-   }
   // segno i lati aperti come temp prop della curva
    int nOpenEdges = vnProp[0] + vnProp[1] * 2 + vnProp[2] * 4 + vnProp[3] * 8 ;
    pMCrv->SetTempProp( nOpenEdges, 0) ;
@@ -3625,6 +3701,112 @@ CalcTrapezoidSpiral( ICurveComposite* pCrvPocket, const Frame3d& frTrap, double
    return true ;
 }
 
+//----------------------------------------------------------------------------
+static bool
+AdjustTrapeziodLeadIn( ICurveComposite* pCrvRes, const PocketParams& PockParam,
+                       const ISurfFlatRegion* pSfrChunk)
+{
+  // recupero la TempProp
+   int nTmpProp = pCrvRes->GetTempProp( 0) ;
+  // se esiste almeno un aperto
+   if ( nTmpProp > 0) {
+     // se solo lato3 aperto
+      bool bCheckHead = ( nTmpProp != 8  &&  nTmpProp != 2) ;
+      if ( nTmpProp == 2)
+         pCrvRes->Invert() ; // entro dall'unico aperto
+      if ( bCheckHead) {
+        // recupero gli estremi della curva corrente e la inverto in base alla Testa
+         Point3d ptS ; pCrvRes->GetStartPoint( ptS) ;
+         Point3d ptE ; pCrvRes->GetEndPoint( ptE) ;
+         Point3d ptSGlob = GetToGlob( ptS, PockParam.frLocXY) ;
+         Point3d ptEGlob = GetToGlob( ptE, PockParam.frLocXY) ;
+         if ( ( PockParam.bAboveHead  &&  ptEGlob.z > ptSGlob.z)  ||
+              ( ! PockParam.bAboveHead  &&  ptEGlob.z < ptSGlob.z))
+            pCrvRes->Invert() ;
+      }
+   }
+
+  // Assegno la Feed
+      AssignFeedSpiralOpt( 1, PockParam, pCrvRes) ;
+  // Se curva da invertire, inverto
+   if ( PockParam.bInvert)
+      pCrvRes->Invert() ;
+
+  // se esiste almeno un aperto, provo ad estendere il percorso
+   if ( nTmpProp > 0) {
+     // Calcolo eventuale entrata da fuori
+      Vector3d vtRef ; pCrvRes->GetStartDir( vtRef) ;
+      vtRef.Invert() ;
+      bool bIsStartExtended = false ;
+      if ( ! ExtendPath( pCrvRes, pSfrChunk, PockParam, vtRef, false, PockParam.dRad + PockParam.dOpenMinSafe, bIsStartExtended))
+         return false ;
+   }
+
+   return true ;
+}
+
+//----------------------------------------------------------------------------
+static bool
+GetZigZagOptimizedCurves( const ISurfFlatRegion* pSfrChunk, const PocketParams& PockParam,
+                          ICurveComposite* pCrvRes)
+{
+  // controllo dei parametri
+   if ( pSfrChunk == nullptr  ||  ! pSfrChunk->IsValid())
+      return false ;
+   pCrvRes->Clear() ;
+
+  // ricavo la curva di bordo del chunk corrente
+   PtrOwner<ICurveComposite> pCrvBorder( ConvertCurveToComposite( pSfrChunk->GetLoop( 0, 0))) ;
+   if ( IsNull( pCrvBorder)  ||  ! pCrvBorder->IsValid())
+      return false ;
+   pCrvBorder->MergeCurves( 10. * EPS_SMALL, 10. * EPS_ANG_SMALL, true, true) ;
+   pCrvBorder->SetExtrusion( pSfrChunk->GetNormVersor()) ;
+
+   /* TRAPEZI
+     - E' richiesto che una dimensione del box della curva sia compatibile con il primo Offset, il
+       quale sarebbe una singola curva aperta
+  */
+   PtrOwner<ICurveComposite> pCrvTrap( CreateCurveComposite()) ;
+   if ( IsNull( pCrvTrap))
+      return false ;
+   Frame3d frTrap ;
+   double dPocketSize ;
+   int nBase, nSecondBase ;
+   bool bOkTrap = GetTrapezoidFromShape( pCrvBorder, pCrvTrap, frTrap, PockParam, dPocketSize, nBase, nSecondBase) ;
+   if ( bOkTrap  &&  pCrvTrap->IsValid()) {
+     // verifico se il trapezio ottenuto deve o meno rispettare il SideStep
+      bool bForcedTrap = false ;
+      IsForcedStepTrapezoid( pCrvTrap, PockParam, nBase, nSecondBase, bForcedTrap) ;
+      if ( ! bForcedTrap)
+         bOkTrap = ( dPocketSize < PockParam.dMaxOptSize + 10. * EPS_SMALL) ;
+   }
+   if ( bOkTrap  &&  pCrvTrap->IsValid()) {
+      pCrvTrap->SetExtrusion( Z_AX) ;
+      CalcTrapezoidSpiral( pCrvTrap, frTrap, dPocketSize, nBase, nSecondBase, PockParam, pCrvRes, bOkTrap) ;
+      if ( bOkTrap) {
+        // verifico che tale curva non interferisca con la regione limite
+         if ( PockParam.SfrLimit.IsValid()) {
+            double dOffsCheck = PockParam.dRad + PockParam.dRadialOffset - 50. * EPS_SMALL ; // restrittivo per sicurezza
+            PtrOwner<ISurfFlatRegion> pSfrToolShape( GetSurfFlatRegionFromFatCurve( pCrvRes->Clone(), dOffsCheck, false, false)) ;
+            bOkTrap = ( ! IsNull( pSfrToolShape)  &&  pSfrToolShape->IsValid()) ;
+            if ( bOkTrap) {
+               bOkTrap = ( pSfrToolShape->Intersect( PockParam.SfrLimit) &&
+                           ! pSfrToolShape->IsValid()) ;
+               if ( ! bOkTrap)
+                  pCrvRes->Clear() ;
+            }
+         }
+         if ( bOkTrap) {
+           AdjustTrapeziodLeadIn( pCrvRes, PockParam, pSfrChunk) ;
+          // imposto il flag di curva singola
+            pCrvRes->SetTempProp( TEMP_PROP_OPT_TRAPEZOID, 0) ;
+         }
+      }
+   }
+
+   return true ;
+}
+
 //----------------------------------------------------------------------------
 static bool
 GetSpiralOptimizedCurves( const ISurfFlatRegion* pSfrChunk, const PocketParams& PockParam,
@@ -3673,7 +3855,7 @@ GetSpiralOptimizedCurves( const ISurfFlatRegion* pSfrChunk, const PocketParams&
                                 ssize( ccClass) == 1  &&  ccClass[0].nClass == CRVC_OUT) ;
                  // NB. una versione più complessa dovrebbe verificare se la sottrazione tra la
                  //     superficie dell'utensile e la regione limite non genera un'altra circonferenza...
-                 //     In questo caso si la sottrazione potrebbe essere trattata come una circonferenza
+                 //     In questo caso la sottrazione potrebbe essere trattata come una circonferenza
                  //     chiusa ed essere ancora svotata a spirale...
                }
             }
@@ -3714,8 +3896,14 @@ GetSpiralOptimizedCurves( const ISurfFlatRegion* pSfrChunk, const PocketParams&
    Frame3d frTrap ;
    double dPocketSize ;
    int nBase, nSecondBase ;
-   bool bOkTrap = ( GetTrapezoidFromShape( pCrvBorder, pCrvTrap, frTrap, PockParam, dPocketSize, nBase, nSecondBase)  &&
-                    dPocketSize < PockParam.dMaxOptSize + 10. * EPS_SMALL) ;
+   bool bOkTrap = GetTrapezoidFromShape( pCrvBorder, pCrvTrap, frTrap, PockParam, dPocketSize, nBase, nSecondBase) ;
+   if ( bOkTrap  &&  pCrvTrap->IsValid()) {
+     // verifico se il trapezio ottenuto deve o meno rispettare il SideStep
+      bool bForcedTrap = false ;
+      IsForcedStepTrapezoid( pCrvTrap, PockParam, nBase, nSecondBase, bForcedTrap) ;
+      if ( ! bForcedTrap)
+         bOkTrap = ( dPocketSize < PockParam.dMaxOptSize + 10. * EPS_SMALL) ;
+   }
    if ( bOkTrap  &&  pCrvTrap->IsValid()) {
       pCrvTrap->SetExtrusion( Z_AX) ;
       CalcTrapezoidSpiral( pCrvTrap, frTrap, dPocketSize, nBase, nSecondBase, PockParam, pCrvRes, bOkTrap) ;
@@ -3733,20 +3921,9 @@ GetSpiralOptimizedCurves( const ISurfFlatRegion* pSfrChunk, const PocketParams&
             }
          }
          if ( bOkTrap) {
-           // calcolo eventuali uscite e ingressi
-            if ( pCrvRes->GetTempProp( 0) > 0) {
-               Vector3d vtRef ; pCrvRes->GetStartDir( vtRef) ;
-               vtRef.Invert() ;
-               bool bIsStartExtended = false ;
-               if ( ! ExtendPath( pCrvRes, pSfrChunk, PockParam, vtRef, false, PockParam.dRad + PockParam.dOpenMinSafe, bIsStartExtended))
-                  return false ;
-               pCrvRes->GetEndDir( vtRef) ;
-               bool bIsEndExtended = false ;
-               if ( ! ExtendPath( pCrvRes, pSfrChunk, PockParam, vtRef, true, PockParam.dRad + PockParam.dOpenMinSafe, bIsEndExtended))
-                  return false ;
-               if ( bIsEndExtended  &&  ! bIsStartExtended)
-                  pCrvRes->Invert() ;
-            }
+            AdjustTrapeziodLeadIn( pCrvRes, PockParam, pSfrChunk) ;
+           // imposto il flag di curva singola
+            pCrvRes->SetTempProp( TEMP_PROP_OPT_TRAPEZOID, 0) ;
          }
       }
    }
@@ -3789,23 +3966,31 @@ GetPocketingOptimizedCurves( ISurfFlatRegion* pSfr, const PocketParams& PockPara
            PockParam.nType == POCKET_CONFORMAL_ZIGZAG  ||  PockParam.nType == POCKET_CONFORMAL_ONEWAY) {
         // curva da resituire
          PtrOwner<ICurveComposite> pCrvOptSpiral( CreateCurveComposite()) ;
-         if (  IsNull( pCrvOptSpiral)  ||
-               ! GetSpiralOptimizedCurves( pSfrChunk, PockParam, pCrvOptSpiral))
+         if ( IsNull( pCrvOptSpiral)  ||
+              ! GetSpiralOptimizedCurves( pSfrChunk, PockParam, pCrvOptSpiral))
             return false ;
         // se ho ricavato una curva ottimizzata
          if ( ! IsNull( pCrvOptSpiral)  &&  pCrvOptSpiral->IsValid()  &&  pCrvOptSpiral->GetCurveCount() > 0) {
-            vCrvOptCurves.emplace_back( Release( pCrvOptSpiral)) ; 
+            vCrvOptCurves.emplace_back( Release( pCrvOptSpiral)) ;
+            pSfr->EraseChunk( nCurrChunk) ;
+         }
+         else
+            ++ nCurrChunk ;
+      }
+      else if ( PockParam.nType == POCKET_ZIGZAG  ||  PockParam.nType == POCKET_ONEWAY) {
+        // curva da restituire
+         PtrOwner<ICurveComposite> pCrvOptZigZag( CreateCurveComposite()) ;
+         if ( IsNull( pCrvOptZigZag)  ||
+              ! GetZigZagOptimizedCurves( pSfrChunk, PockParam, pCrvOptZigZag))
+            return false ;
+        // se ho ricavato una curva ottimizzata
+         if ( ! IsNull( pCrvOptZigZag)  &&  pCrvOptZigZag->IsValid()  &&  pCrvOptZigZag->GetCurveCount() > 0) {
+            vCrvOptCurves.emplace_back( Release( pCrvOptZigZag)) ;
             pSfr->EraseChunk( nCurrChunk) ;
          }
          else
             ++ nCurrChunk ;
       }
-     // else if ( PockParam.nType == POCKET_ZIGZAG)
-     //    ;
-     // else if ( PockParam.nType == POCKET_ONEWAY)
-     //    ;
-     // else if ( PockParam.nType == POCKET_CONFORMAL_ONEWAY  ||  PockParam.nType == POCKET_CONFORMAL_ZIGZAG)
-     //    ;
    }
 
    return true ;
@@ -6044,6 +6229,7 @@ CalcSpiral( const ISurfFlatRegion* pSfrPock, const ISurfFlatRegion* pSfrOrig, co
   // ricavo le regioni progressive
    double dOffsPrec = 0. ;
    int nCrvFirstOffs = 0 ;
+   bool bLastNotValid = false ;
    while ( nIter < MAX_ITER) {
      // Offset della regione attuale
       PtrOwner<ISurfFlatRegion> pSfrOffsVR( pSrfAct->CreateOffsetSurf( - dOffs, ICurve::OFF_FILLET)) ;
@@ -6054,6 +6240,7 @@ CalcSpiral( const ISurfFlatRegion* pSfrPock, const ISurfFlatRegion* pSfrOrig, co
           pSfrOffsVR.Set( pSrfAct->CreateOffsetSurf( - dOffs + 5 * EPS_SMALL, ICurve::OFF_FILLET)) ;
           if ( IsNull( pSfrOffsVR))
              return false ;
+          bLastNotValid = true ;
       }
 
      // se primo Offset
@@ -6085,6 +6272,10 @@ CalcSpiral( const ISurfFlatRegion* pSfrPock, const ISurfFlatRegion* pSfrOrig, co
             bool bInsert = true ;
             if ( ! CheckIfOffsetIsNecessary( pSrfAct, pCrvCompoBorder, dOffs, dOffsPrec, nIter, PockParams, bInsert))
                return false ;
+           // per evitare di allacciare una curva di regione non svuotata alla prima curva di Offset
+           // ( quindi avere un entrata nel pieno del grezzo) controllo di non eliminare il secondo Offset
+            if ( ! bInsert  &&  nChunks == 1  &&  bLastNotValid)
+                  bInsert = true ;
             if ( bInsert) {
                // imposto come secondo TempParam il Side di classificazione
                   pCrvCompoBorder->SetTempParam( j == 0 ? MDS_RIGHT : MDS_LEFT, 1) ;
@@ -8901,14 +9092,14 @@ CalcSpiralPocketing( const ISurfFlatRegion* pSfr, int nType, const PocketParams&
   // il tipo può essere solo SpiralIn o SpiralOut
    if ( nType != POCKET_SPIRALIN  &&  nType != POCKET_SPIRALOUT)
       return false ;
+   PtrOwner<ISurfFlatRegion> pSfrLimit( PockParams.SfrLimit.IsValid() ? PockParams.SfrLimit.Clone() : CreateSurfFlatRegion()) ;
   // calcolo il percorso di svuotatura spiral
    return ( CalcPocketing( pSfr, PockParams.dRad, PockParams.dRadialOffset, PockParams.dSideStep,
                            PockParams.dAngle, PockParams.dOpenMinSafe, nType, PockParams.bSmooth,
                            PockParams.bCalcUnclearedRegs, PockParams.bInvert, PockParams.bAvoidOpt,
                            PockParams.bAllowZigZagOneWayBorders, PockParams.bCalcFeed, PockParams.ptStart,
-                           PockParams.SfrLimit.IsValid() ? PockParams.SfrLimit.Clone() : CreateSurfFlatRegion(),
-                           PockParams.bAvoidOpt, PockParams.dMaxOptSize, PockParams.dLiTang, PockParams.nLiType,
-                           vCrvCompoRes)) ;
+                           pSfrLimit, PockParams.bAvoidOpt, PockParams.dMaxOptSize, PockParams.dLiTang,
+                           PockParams.nLiType, vCrvCompoRes)) ;
 }
 
 //----------------------------------------------------------------------------
@@ -9203,9 +9394,10 @@ SmoothExtensionLinesByIntersection( ICRVCOMPOPOVECTOR& vCrvPaths, const PocketPa
 //----------------------------------------------------------------------------
 bool
 CalcPocketing( const ISurfFlatRegion* pSfr, double dRad, double dRadOffs, double dStep, double dAngle,
-               double dOpenMinSafe, int nType, bool bSmooth, bool bCalcUnclReg, bool bInvert, bool bAvoidOpt, bool bAllowZigZagOneWayBorders,
-               bool bCalcFeed, const Point3d& ptEndPrec, const ISurfFlatRegion* pSfrLimit, bool bAllOffs,
-               double dMaxOptSize, double dLiTang, int nLiType, ICRVCOMPOPOVECTOR& vCrvCompoRes)
+               double dOpenMinSafe, int nType, bool bSmooth, bool bCalcUnclReg, bool bInvert, bool bAvoidOpt,
+               bool bAllowZigZagOneWayBorders, bool bCalcFeed, const Point3d& ptEndPrec,
+               const ISurfFlatRegion* pSfrLimit, bool bAllOffs, double dMaxOptSize, double dLiTang,
+               int nLiType, ICRVCOMPOPOVECTOR& vCrvCompoRes)
 {
   // controllo dei parametri
    if ( pSfr == nullptr  ||  ! pSfr->IsValid()  ||

ChainCurves.cpp

@@ -16,6 +16,7 @@
 #include "/EgtDev/Include/EGkChainCurves.h"
 #include "/EgtDev/Include/EGkGeomDB.h"
 #include "/EgtDev/Include/EGkCurve.h"
+#include "/EgtDev/Include/EGkCurveComposite.h"
 #include <algorithm>
 
 using namespace std ;

CurveAux.cpp

@@ -13,6 +13,8 @@
 
 //--------------------------- Include ----------------------------------------
 #include "stdafx.h"
+#include "CalcDerivate.h"
+#include "Bernstein.h"
 #include "CurveAux.h"
 #include "GeoConst.h"
 #include "CurveLine.h"
@@ -23,12 +25,22 @@
 #include "IntersLineLine.h"
 #include "/EgtDev/Include/EGkDistPointCurve.h"
 #include "/EgtDev/Include/EGkStringUtils3d.h"
+#include "/EgtDev/Include/EgtNumUtils.h"
 #include "/EgtDev/Include/EGkUiUnits.h"
 #include "/EgtDev/Include/EgtPointerOwner.h"
 #include "/EgtDev/Include/EGkCurveByInterp.h"
+#include "/EgtDev/Include/EGkChainCurves.h"
 #define EIGEN_NO_IO
 #include "/EgtDev/Extern/Eigen/Dense"
 
+#define SAVEAPPROX 0
+#define SAVECURVEPASSED 0
+#define SAVELINEARAPPROX 0
+#if SAVEAPPROX  || SAVECURVEPASSED  ||  SAVELINEARAPPROX
+   static int nCrvPassed = 0 ;
+   #include "/EgtDev/Include/EGkGeoObjSave.h"
+#endif
+
 using namespace std ;
 
 static bool FindSpan( double dU, int nDeg, const DBLVECTOR& vKnots, int& nSpan) ;
@@ -522,9 +534,9 @@ LineToBezierCurve( const ICurveLine* pCrvLine, int nDeg, bool bMakeRatOrNot)
       return nullptr ;
 
    PtrOwner<ICurveBezier> pCrvBezier( CreateCurveBezier()) ;
-   // rendo tutte le curve di grado 2 e razionali così posso convertire anche archi e avere tutte curve dello stesso grado e razionali
-   pCrvBezier->Init( nDeg, true) ;
-   pCrvBezier->FromLine( *pCrvLine) ;
+   pCrvBezier->Init( nDeg, false) ;
+   if ( ! pCrvBezier->FromLine( *pCrvLine))
+      return nullptr ;
    if ( bMakeRatOrNot)
       pCrvBezier->MakeRational() ;
    return Release( pCrvBezier) ;
@@ -577,6 +589,13 @@ ArcToBezierCurve( const ICurveArc* pArc, int nDeg, bool bMakeRatOrNot)
          PtrOwner<ICurveBezier> pCrvBez( CreateBasicCurveBezier()) ;
          if ( IsNull( pCrvBez)  ||  ! pCrvBez->FromArc( cArc))
             return nullptr ;
+         if ( ! bMakeRatOrNot) {
+            Point3d ptCen = pArc->GetCenter() ;
+            Vector3d vtN = pArc->GetNormVersor() ;
+            pCrvBez.Set( ApproxArcCurveBezierWithSingleCubic( pCrvBez, ptCen, vtN)) ;
+         }
+         if ( IsNull( pCrvBez))
+             return nullptr ;
          // aumento il grado della curva come richiesto
          while ( pCrvBez->GetDegree() < nDeg)
             pCrvBez.Set( BezierIncreaseDegree( pCrvBez)) ;
@@ -1135,7 +1154,7 @@ FindSpan( double dU, int nDeg, const DBLVECTOR& vKnots, int& nSpan)
       return true ;
    }
    // trovo a quale span appartiene il parametro dU
-   int nKnots = int( vKnots.size()) ;
+   int nKnots = ssize( vKnots) ;
    if ( abs( dU - vKnots[nKnots-1]) < EPS_SMALL) {
       nSpan = nKnots - 1 - nDeg ;
       return true ;
@@ -1161,7 +1180,7 @@ static bool
 CalcBasisFunc( double dU, int nSpan, int nDeg, const DBLVECTOR& vKnots, DBLVECTOR& vBasis)
 {
    // mi aspetto che il vettore vBasis sia di lunghezza  nDeg + 1
-   if ( vBasis.size() != nDeg + 1)
+   if ( ssize( vBasis) != nDeg + 1)
       return false ;
 
    vBasis[0] = 1 ;
@@ -1186,7 +1205,8 @@ CalcBasisFunc( double dU, int nSpan, int nDeg, const DBLVECTOR& vKnots, DBLVECTO
 
 //----------------------------------------------------------------------------
 ICurve*
-InterpolatePointSetWithBezierNoIntermedLines( const PNTVECTOR& vPnt, int nStart, int nEnd, int nDeg, const DBLVECTOR& vLen, double dLenTot)
+InterpolatePointSetWithBezierNoIntermedLines( const PNTVECTOR& vPnt, int nStart, int nEnd, int nDeg, const DBLVECTOR& vLen, double dLenTot,
+                                              const Vector3d& vtStartDir = V_NULL, const Vector3d& vtEndDir = V_NULL)
 {
    PtrOwner<ICurve> pCrvInt ;
 
@@ -1215,6 +1235,8 @@ InterpolatePointSetWithBezierNoIntermedLines( const PNTVECTOR& vPnt, int nStart,
          return nullptr ;
    }
 
+   bool bUseStartEndDir = vtStartDir.IsValid()  &&  vtEndDir.IsValid() ;
+
    DBLVECTOR vPntParam ;
    vPntParam.resize( nPoints) ;
    vPntParam[0] = 0 ;
@@ -1223,13 +1245,15 @@ InterpolatePointSetWithBezierNoIntermedLines( const PNTVECTOR& vPnt, int nStart,
       vPntParam[i] = vPntParam[i-1] + vLen[i-1] / dLenTot ;
 
    DBLVECTOR vKnots ;
-   vKnots.resize( nPoints + nDeg - 1) ;
+   int nKnots = bUseStartEndDir ? nPoints + nDeg - 1 + 2 : nPoints + nDeg - 1 ;
+   vKnots.resize( nKnots) ;
    for ( int i = 0 ; i < nDeg ; ++i) {
       vKnots[i] = 0 ;
       vKnots.end()[-i-1] = 1 ;
    }
 
-   for ( int i = nDeg ; i < nPoints - 1 ; ++i) {
+   int nKnotsToEdit = bUseStartEndDir ? nPoints + 1 : nPoints - 1 ;
+   for ( int i = nDeg ; i < nKnotsToEdit ; ++i) {
       double dKnot = 0 ;
       for ( int j = i + 1 ; j < i + nDeg + 1 ; ++j)
          dKnot += vPntParam[j - nDeg] ;
@@ -1237,13 +1261,14 @@ InterpolatePointSetWithBezierNoIntermedLines( const PNTVECTOR& vPnt, int nStart,
       vKnots[i] = dKnot ;
    }
 
-   Eigen::MatrixXd mA( nPoints, nPoints) ;
+   int nEq = bUseStartEndDir ? nPoints + 2 : nPoints ;
+   Eigen::MatrixXd mA( nEq, nEq) ;
    mA.fill( 0) ;
-   for ( int i = 0 ; i < nPoints ; ++i) {
+   for ( int i = 0 ; i < nEq ; ++i) {
       if ( i == 0)
          mA.row(0).col(0) << 1 ;
-      else if ( i == nPoints - 1)
-         mA.row(i).col(nPoints - 1) << 1 ;
+      else if ( i == nEq - 1)
+         mA.row(i).col( nEq - 1) << 1 ;
       else {
          int nSpan = 0 ; FindSpan( vPntParam[i], nDeg, vKnots, nSpan) ;
          DBLVECTOR vBasis ; vBasis.resize( nDeg + 1) ;
@@ -1294,7 +1319,7 @@ InterpolatePointSetWithBezier( const PNTVECTOR& vPnt, double dLinTol, double dMa
    int nItCount = 0 ;
    while ( dErr > dLinTol  &&  nItCount < 10) {
       pCrvInt->Clear() ;
-      int nPoints = int( vPnt.size()) ;
+      int nPoints = ssize( vPnt) ;
       int nDeg = 3 ;
       if ( nPoints < 2)
          return nullptr ;
@@ -1310,7 +1335,7 @@ InterpolatePointSetWithBezier( const PNTVECTOR& vPnt, double dLinTol, double dMa
       else if ( nPoints == 3) {
          // se ho solo tre punti uso un altro algoritmo
          CurveByInterp cbi ;
-         for ( int i = 0 ; i < int( vPnt.size()) ; ++i)
+         for ( int i = 0 ; i < ssize( vPnt) ; ++i)
             cbi.AddPoint( vPnt[i]) ;
          pCrvInt->AddCurve( cbi.GetCurve( CurveByInterp::AKIMA_CORNER, CurveByInterp::CUBIC_BEZIERS)) ;
          if ( ! IsNull( pCrvInt)  &&  pCrvInt->IsValid())
@@ -1342,10 +1367,16 @@ InterpolatePointSetWithBezier( const PNTVECTOR& vPnt, double dLinTol, double dMa
             }
          }
 
-         if ( vLen.size() != 0) {
+         if ( ! vLen.empty()) {
             if ( nEnd == 0)
                nEnd = nPoints - 1 ;
-            pCrvInt->AddCurve( InterpolatePointSetWithBezierNoIntermedLines( vPnt, nStart, nEnd, nDeg, vLen, dLenTot)) ;
+            Vector3d vtStartDir = V_INVALID ; 
+            Vector3d vtEndDir = V_INVALID ;
+            //if ( nStart != 0  &&  nEnd != nPoints - 1) {
+            //   pCrvInt->GetEndDir( vtStartDir) ;
+            //   vtEndDir = 
+            //}
+            pCrvInt->AddCurve( InterpolatePointSetWithBezierNoIntermedLines( vPnt, nStart, nEnd, nDeg, vLen, dLenTot, vtStartDir, vtEndDir)) ;
          }
 
          if ( bFoundLine) {
@@ -1357,6 +1388,7 @@ InterpolatePointSetWithBezier( const PNTVECTOR& vPnt, double dLinTol, double dMa
          nStart = nEnd ;
       }
 
+      //dErr = 0 ;
       CalcApproxError( pCrvOri, pCrvInt, dErr) ;
       if ( dErr > dLinTol  &&  dMaxLen > 200 * EPS_SMALL)
          dMaxLen /= 2 ;
@@ -1371,42 +1403,318 @@ InterpolatePointSetWithBezier( const PNTVECTOR& vPnt, double dLinTol, double dMa
 }
 
 //----------------------------------------------------------------------------
-ICurve*
-ApproxCurveWithBezier( const ICurve* pCrv , double dTol, int nType)
+static bool
+ParamByLen( const PNTVECTOR& vPnt, DBLVECTOR& vParam, int nFirst, int nLast)
 {
-   PolyLine plApprox ;
-   double dAngTolFine = 2 ;
-   pCrv->ApproxWithLines( dTol, dAngTolFine, ICurve::APL_STD, plApprox) ;
+   int nPoints = nLast - nFirst + 1 ;
+   if ( nPoints < 2)
+      return false ;
+   if ( vParam.empty())
+      vParam.resize( nPoints) ;
+   if ( vParam[nFirst] == 0  &&  vParam[nLast] == 1)
+      return true ;
+   vParam[nFirst] = 0 ;
+   for ( int i = nFirst + 1 ; i <= nLast ; ++i) {
+      double dDist = Dist( vPnt[i], vPnt[i-1]) ;
+      vParam[i] = vParam[i- 1] + dDist ;
+   }
+   for ( int i = nFirst + 1 ; i < nLast ; ++i)
+      vParam[i] /= vParam[nLast] ;
+   vParam[nLast] = 1 ;
 
-   PNTVECTOR vPnt ;
-   Point3d pt ; plApprox.GetFirstPoint( pt) ;
-   do {
-      vPnt.push_back( pt) ;
-   } while ( plApprox.GetNextPoint( pt)) ;
+   return true ;
+}
 
-   // campiono punti lungo la curva e poi li interpolo
+//----------------------------------------------------------------------------
+ICurveBezier*
+ApproxPointSetWithSingleBezier( const PNTVECTOR& vPnt, int nFirst, int nLast,
+                                const Vector3d& vtStartDir, const Vector3d& vtEndDir, const DBLVECTOR& vParam)
+{
+   // cerco di approssimare un set di punti con una sola bezier cubica non razionale
+   int nPoints = nLast - nFirst + 1 ;
+   // se ho solo quattro punti allora costruisco direttamente la curva
+   PtrOwner<ICurveBezier> pCrvBez( CreateCurveBezier()) ;
+   int nDeg = 3 ;
+   bool bRat = false ;
+   pCrvBez->Init( nDeg, bRat) ;
+   const Point3d& pt0 = vPnt[nFirst] ;
+   const Point3d& pt3 = vPnt[nLast] ;
+   pCrvBez->SetControlPoint( 0, pt0) ;
+   pCrvBez->SetControlPoint( 3, pt3) ;
+   Eigen::Vector2d mA ;
+   if ( nPoints > 4) {
+      // risoluzione sistema
+      Eigen::Matrix2d mC ; mC.setZero() ;
+      Eigen::Vector2d mX ; mX.setZero() ;
+      for ( int i = nFirst ; i <= nLast ; ++i) {
+         double dU = vParam[i] ;
+         DBLVECTOR vBern(4) ;
+         GetAllBernstein( dU, 3, vBern) ;
 
-   PtrOwner<ICurve> pCC( InterpolatePointSetWithBezier( vPnt, dTol, 100)) ;
-   if ( ! IsNull( pCC)  &&  pCC->IsValid())
-         return Release( pCC) ;
-      else
-         return nullptr ;
+         Vector3d A1 = vtStartDir * vBern[1] ;
+         Vector3d A2 = vtEndDir * vBern[2] ;
+
+         Vector3d tmp = vPnt[i] - ( pt0 * ( vBern[0] + vBern[1])) - (( pt3 * ( vBern[2] + vBern[3])) - ORIG) ; // ORIG serve solo per trasformare l'ultimo termine in un vettore
+
+         mC(0,0) += A1 * A1 ;
+         mC(0,1) += A1 * A2 ;
+         mC(1,0) += A1 * A2 ;
+         mC(1,1) += A2 * A2 ;
+
+         mX(0) += A1 * tmp ;
+         mX(1) += A2 * tmp ;
+      }
+      mA = mC.fullPivLu().solve(mX) ;
+   }
+   // l'algoritmo è fatto in modo che alpha1 e alpha2 siano positivi ( se tutto va bene)
+   // io invece ho tenuto le tangenti con la direzione originale, quindi il primo dovrebbe essere positivo e il secondo negativo
+   if ( mA(0) < 0  ||  mA(1) > 0  ||  nPoints < 4) {
+      if ( mA(0) < 0  ||  mA(1) > 0)
+         LOG_DBG_ERR( GetEGkLogger(), "valori di alfa sballati, potrebbe essere la spaziatura dismogenea tra punti")
+      double dDistCorr = Dist( pt3, pt0) / 3 ;
+      mA(0) = dDistCorr ;
+      mA(1) = - dDistCorr ;
+   }
+
+   Point3d pt1 = pt0 + vtStartDir * mA(0) ;
+   Point3d pt2 = pt3 + vtEndDir * mA(1) ;
+   pCrvBez->SetControlPoint( 1, pt1) ;
+   pCrvBez->SetControlPoint( 2, pt2) ;
+   return Release( pCrvBez) ;
+}
+
+//----------------------------------------------------------------------------
+bool
+CalcPointSetApproxError( const PNTVECTOR& vPntOrig, const DBLVECTOR& vParam,
+                         int nFirst, int nLast, const ICurve* pCrvNew, double& dErr, int& nPointMaxErr)
+{
+   dErr = 0 ;
+   // calcolo l'errore di approssimazione
+   for ( int i = nFirst ; i <= nLast ; ++i) {
+      Point3d ptBez ; pCrvNew->GetPointD1D2( vParam[i], ICurve::Side::FROM_MINUS, ptBez) ;
+      double dErrTemp = Dist( vPntOrig[i], ptBez) ;
+      if ( dErrTemp > dErr) {
+         dErr = dErrTemp ;
+         nPointMaxErr = i ;
+      }
+   }
+
+   return true ;
 }
 
 //----------------------------------------------------------------------------
 ICurve*
-ApproxPointSetWithBezier( const ICurve* pCrv, double dTol)
+FitWithBezier( const ICurve* pCrvOrig, const PNTVECTOR& vPnt, DBLVECTOR& vParam,
+               int nFirst, int nLast, const VCT3DVECTOR& vPrevDer, const VCT3DVECTOR& vNextDer, double dTol, bool bLimitSplit = false)
 {
-   // campiono punti lungo la curva e poi li interpolo
+   ParamByLen( vPnt, vParam, nFirst, nLast) ;
+   PtrOwner<ICurveComposite> pCrvFit( CreateCurveComposite()) ;
+   PtrOwner<ICurveBezier> pCrvBez( CreateCurveBezier()) ;
+   double dErr = INFINITO ;
+   double dErrPrec = INFINITO ;
+   double dErrSplit = dTol * 25 ;
+   int nIter = 0 ;
+   while ( dErr > dTol  &&  nIter < 10) {
+      if ( dErr < INFINITO) {
+        // riparametrizzo i punti
+         for ( int i = nFirst + 1 ; i < nLast ; ++i) {
+           // questo potrebbe diventare un while appena capisco di quanto si aggiusta il parametro ad ogni iterazione
+            double dCorr = 1 ;
+            do {
+               Vector3d vtDer1, vtDer2 ;
+               Point3d ptBez ; pCrvBez->GetPointD1D2( vParam[i], ICurve::Side::FROM_MINUS, ptBez, &vtDer1, &vtDer2) ;
+               Vector3d vtLink = ptBez - vPnt[i] ;
+               double dNum = vtLink * vtDer1 ;
+               double dDen = vtLink * vtDer2 + vtDer1 * vtDer1 ;
+               dCorr = dDen > EPS_ZERO ? dNum / dDen : 0 ;
+               vParam[i] = vParam[i] - dCorr ;
+            } while ( abs( dCorr) > EPS_ZERO) ;
+            Clamp( vParam[i], 0., 1.) ;
+         }
+      }
 
-   PtrOwner<ICurveComposite> pCC( CreateBasicCurveComposite()) ;
-   return Release( pCC) ;
+     // fit della curva
+      Vector3d vtStartDir, vtEndDir ;
+      if ( bLimitSplit) {
+         vtStartDir = vNextDer[nFirst / 3] ;
+         vtEndDir = vPrevDer[nLast / 3] ;
+      }
+      else {
+         vtStartDir = vNextDer[nFirst] ;
+         vtEndDir = vPrevDer[nLast] ;
+      }
+      pCrvBez.Set( ApproxPointSetWithSingleBezier( vPnt, nFirst, nLast, vtStartDir, vtEndDir, vParam)) ;
+      if ( IsNull( pCrvBez)  ||  ! pCrvBez->IsValid())
+         return nullptr ;
+#if SAVEAPPROX
+      SaveGeoObj( pCrvBez->Clone(), "D:\\Temp\\bezier\\approxWithBezier\\"+ToString(nCrvPassed) + "first_approx.nge") ;
+#endif
+
+      int nPointMaxErr = 0 ;
+      CalcPointSetApproxError( vPnt, vParam, nFirst, nLast, pCrvBez, dErr, nPointMaxErr) ;
+     // se sto unendo due punti consecutivi e l'errore è oltre quello richiesto allora restituisco un segmento che unisce i punti
+      if ( ((nLast - nFirst == 1)  ||  ( bLimitSplit  &&  nLast - nFirst == 3)) &&  dErr > dTol) {
+         CurveLine CL ; CL.Set( vPnt[nFirst], vPnt[nLast]) ;
+         pCrvBez.Set( GetCurveBezier( CurveToBezierCurve( &CL))) ;
+         dErr = 0 ;
+      }
+
+      if ( bLimitSplit  &&  nPointMaxErr % 3 != 0) {
+         nPointMaxErr = 3 * int( round( nPointMaxErr / 3.)) ;
+         if ( nPointMaxErr == nFirst)
+            nPointMaxErr += 3 ;
+         else if( nPointMaxErr == nLast)
+            nPointMaxErr -= 3 ;
+      }
+
+      ++nIter ;
+      bool bSplit = false ;
+      if ( ( nIter == 10  &&  dErr > dTol)  || dErr > dErrPrec  ||  dErrPrec - dErr < dErrPrec / 20)
+         bSplit = true ;
+      dErrPrec = dErr ;
+     // se la curva di approssimazione è ancora molto lontana dalla curva originale allora divido il set di punti in due 
+      if ( dErr > dErrSplit  ||  bSplit) {
+         if ( nLast - nFirst > 1  &&  nPointMaxErr - nFirst > 1  &&  nLast - nPointMaxErr > 1) {
+            if ( ! pCrvFit->AddCurve( FitWithBezier( pCrvOrig, vPnt, vParam, nFirst, nPointMaxErr, vPrevDer, vNextDer,dTol, bLimitSplit)) ||
+                ! pCrvFit->AddCurve( FitWithBezier( pCrvOrig, vPnt, vParam, nPointMaxErr, nLast, vPrevDer, vNextDer, dTol, bLimitSplit)))
+               return nullptr ;
+            break ;
+         }
+         else
+            return nullptr ;
+      }
+   }
+
+   if ( pCrvFit->GetCurveCount() > 0)
+      return Release( pCrvFit) ;
+   else if ( dErr < dTol  &&  ! IsNull( pCrvBez)  &&  pCrvBez->IsValid())
+      return Release( pCrvBez) ;
+   else
+      return nullptr ;
+}
+
+//----------------------------------------------------------------------------
+ICurve*
+ApproxCurveWithBezier( const ICurve* pCrv , double dTol)
+{
+
+#if SAVECURVEPASSED
+   SaveGeoObj( pCrv->Clone(), "D:\\Temp\\bezier\\approxWithBezier\\CurveDaApprossimare\\"+ToString(nCrvPassed) + ".nge") ;
+   ++nCrvPassed ;
+#endif
+
+   //// uso l'algoritmo di Schneider in Grafic Gems I
+   // mi aspetto che non ci siano angoli ( discontinuità della derivata prima) nel risultato desiderato
+   PolyLine plApprox ;
+   double dAngTolFine = 1 ;
+   double dLinTolFine = 0.05 ;
+   pCrv->ApproxWithLines( dLinTolFine, dAngTolFine, ICurve::APL_STD, plApprox) ;
+
+#if SAVELINEARAPPROX
+   CurveComposite CC ; CC.FromPolyLine(plApprox) ;
+   SaveGeoObj( CC.Clone(), "D:\\Temp\\bezier\\approxWithBezier\\approssimazione_lineare.nge") ;
+#endif
+
+   PNTVECTOR vPnt ;
+   PNTVECTOR vPntOverSampling ;
+   Point3d pt ; plApprox.GetFirstPoint( pt) ;
+   do {
+      if ( ! vPntOverSampling.empty()) {
+         vPntOverSampling.push_back( Media( vPnt.back(), pt,1./3.)) ;
+         vPntOverSampling.push_back( Media( vPnt.back(), pt,2./3.)) ;
+      }
+      vPntOverSampling.push_back( pt) ;
+      vPnt.push_back( pt) ;
+   } while ( plApprox.GetNextPoint( pt)) ;
+
+   // calcolo la curvatura nei vari punti per identificare zone a curvatura costante
+   DBLVECTOR vRad ( ssize( vPnt)) ;
+   for ( int i = 1 ; i < ssize( vPnt) - 1 ; ++i) {
+      Vector3d vtA = vPnt[i] - vPnt[i-1] ;
+      Vector3d vtB = vPnt[i+1] - vPnt[i-1] ;
+      double dR = ( vtA.Len() * vtB.Len() * ( vtA - vtB).Len()) / ( 2 * ( vtA ^ vtB).Len()) ;
+      vRad[i] = dR ;
+   }
+   vRad[0] = vRad[1] ;
+   vRad.back() = vRad.end()[-2] ;
+   // identifico le zone a curvatura costante // primo e ultimo punto degli intervalli non devono avere curvatura uguale agli altri perché sono di frontiera
+   INTINTVECTOR vConstCurv ;
+   double dRadPrec = vRad[0] ;
+   int nStart = 0 ;
+   int nEnd = 1 ;
+   double dRatio = 1.5 ;
+   while ( nStart < ssize( vPnt) - 1) {
+      double dRadTol = max( max( vRad[nEnd], dRadPrec) / 10 , 1.) ;
+      if ( dRadPrec > dRatio * vRad[nEnd]  ||  dRatio * dRadPrec < vRad[nEnd])
+         dRadTol = 0 ;
+      while ( nEnd < ssize( vPnt) - 1  &&  abs( vRad[nEnd] - dRadPrec) < dRadTol) {
+         dRadPrec = vRad[nEnd] ;
+         ++nEnd ;
+      }
+      vConstCurv.emplace_back( nStart * 3, nEnd * 3) ;
+      nStart = nEnd ;
+      dRadPrec = vRad[nEnd] ;
+      ++nEnd ;
+   }
+   if ( vConstCurv.empty())
+      vConstCurv.emplace_back( 0, ssize( vPnt) - 1) ;
+
+   int nPoints = ssize( vPnt) ;
+   DBLVECTOR vParam( nPoints) ;
+   int nFirst = 0 ;
+   int nLast = ssize( vPnt) - 1 ;
+   ParamByLen( vPnt, vParam, nFirst, nLast) ;
+
+   VCT3DVECTOR vPrevDer ;
+   VCT3DVECTOR vNextDer ;
+   ComputeAkimaTangents( false, vParam, vPnt, vPrevDer, vNextDer) ;
+
+   int nOverSampling = ssize( vPntOverSampling) ;
+   vParam.resize( nOverSampling) ;
+   nFirst = 0 ;
+   nLast = nOverSampling - 1 ;
+   ParamByLen( vPntOverSampling, vParam, nFirst, nLast) ;
+
+   //normalizzo tutte le derivate
+   for ( int i = 0 ; i < ssize( vPrevDer) ; ++i) {
+      vPrevDer[i].Normalize() ;
+      vNextDer[i].Normalize() ;
+   }
+
+   // potrei verificare prima se un tratto è retto e aggiustare le tangenti del tratto precedente e successivo prima di approssimare
+   PtrOwner<ICurveComposite> pCCApproxTot( CreateCurveComposite()) ;
+   for ( INTINT iiSE : vConstCurv) {
+      nFirst = iiSE.first ;
+      nLast = iiSE.second ;
+      // riconosco se ho un tratto retto
+      int nPnt = nFirst / 3 ;
+      if ( nLast - nFirst == 3  &&  ( nPnt > 0 &&  vRad[nPnt] > dRatio * vRad[nPnt - 1])  &&  ( nPnt < ssize( vRad)  &&  vRad[nPnt]> dRatio * vRad[nPnt + 1])) {
+         CurveLine CL ; CL.Set( vPntOverSampling[nFirst], vPntOverSampling[nLast]) ;
+         PtrOwner<ICurveBezier> pCApprox( LineToBezierCurve( &CL, 3, false)) ;
+         if ( ! pCCApproxTot->AddCurve( Release( pCApprox)))
+            return nullptr ;
+      }
+      else {
+         //definisco la bezier che vado a raffinare iterativamente
+         PtrOwner<ICurve> pCApprox( FitWithBezier( pCrv, vPntOverSampling, vParam, nFirst, nLast, vPrevDer, vNextDer, dTol, true)) ;
+         if ( IsNull( pCApprox)  ||  ! pCApprox->IsValid())
+            return nullptr ;
+         if ( ! pCCApproxTot->AddCurve( Release( pCApprox)))
+            return nullptr ;
+      }
+   }
+
+   return Release( pCCApproxTot) ;
 }
 
 //----------------------------------------------------------------------------
 bool
 CalcApproxError( const ICurve* pCrvOri, const ICurve* pCrvNew, double& dErr, int nPoints)
 {
+   if ( pCrvOri == nullptr  ||  ! pCrvOri->IsValid()  ||  pCrvNew == nullptr  ||  ! pCrvNew->IsValid()){
+      dErr = INFINITO ;
+      return false ;
+   }
    // controllo l'errore effettivo campionando più finemente
    double dLenOri = 0 ; pCrvOri->GetLength( dLenOri) ;
    double dLenNew = 0 ; pCrvNew->GetLength( dLenNew) ;
@@ -1488,7 +1796,7 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
 {
   // se con nodi extra
    if ( cnData.bExtraKnotes) {
-      int nKnotesNbr = int( cnData.vU.size()) ;
+      int nKnotesNbr = ssize( cnData.vU) ;
       if ( nKnotesNbr < 4)
          return false ;
       cnData.bExtraKnotes = false ;
@@ -1502,7 +1810,7 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
       bool bAlreadyChecked = false ;
      // se la curva è peridica verifco che effettivamente ci sia un numero di punti ripetituti uguale al grado della curva
      // wrap della curva su se stessa
-      if ( cnData.bPeriodic  &&  (int(cnData.vU.size()) > int(cnData.vCP.size()) + cnData.nDeg - 1)) {
+      if ( cnData.bPeriodic  &&  ( ssize( cnData.vU) > ssize( cnData.vCP) + cnData.nDeg - 1)) {
          bool bRepeated = true ;
          for ( int i = 0 ; i < cnData.nDeg ; ++i) {
             if ( ! AreSamePointApprox( cnData.vCP[i], cnData.vCP.end()[-cnData.nDeg + i]) ) {
@@ -1511,11 +1819,11 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
             }
          }
          bool bFirstAddedAtEnd = false ;
-         if ( ! bRepeated  ||  (bRepeated  && AreSamePointApprox( cnData.vCP[0],cnData.vCP[cnData.nDeg]))){
+         if ( ! bRepeated  ||  ( bRepeated  &&  AreSamePointApprox( cnData.vCP[0], cnData.vCP[cnData.nDeg]))) {
            // salvo il vettore dei nodi in caso mi accorga di avere tra le mani una curva unclamped
             DBLVECTOR vU = cnData.vU ;
            // se effettivamente ho dei nodi in più da togliere allora li tolgo ed eventualmente aggiungo punti di controllo
-            if ( int(cnData.vU.size()) > int(cnData.vCP.size()) + cnData.nDeg - 1 ) {
+            if ( ssize( cnData.vU) > ssize( cnData.vCP) + cnData.nDeg - 1 ) {
               // 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())) {
                   bFirstAddedAtEnd = true ;
@@ -1527,11 +1835,11 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
                cnData.vU = DBLVECTOR( cnData.vU.begin(), cnData.vU.end() - cnData.nDeg) ;
               // controllo eventualmente 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
+               if ( ssize( cnData.vU) == ssize( cnData.vCP) + 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) {
+               else if ( ssize( cnData.vU) == ssize( cnData.vCP) + cnData.nDeg) {
                   cnData.vU = vector<double>( cnData.vU.begin() + 1, cnData.vU.end()) ;
                }
             }
@@ -1558,7 +1866,7 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
               // recupero il vettore dei nodi
                cnData.vU = vU ;
               // verifico se ho nodi extra
-               if ( cnData.vU.size() == int( cnData.vCP.size()) + cnData.nDeg + 1 ) {
+               if ( ssize( cnData.vU) == ssize( cnData.vCP) + cnData.nDeg + 1 ) {
                  // significa che ci sono due nodi extra:
                  // se la curva ha grado maggiore di 1 e i primi due nodi sono uguali allora tolgo quelli
                   if ( cnData.nDeg > 1 && abs(cnData.vU[1] - cnData.vU[0]) < EPS_SMALL) {
@@ -1569,7 +1877,7 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
                      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)
+               else if ( ssize( cnData.vU) == ssize( cnData.vCP) + cnData.nDeg)
                   cnData.vU = vector<double>( cnData.vU.begin() + 1, cnData.vU.end()) ;
             }
             bAlreadyChecked = true ;
@@ -1594,7 +1902,7 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
      // qui aggiungo un controllo se la curva è collassata in un punto ( ho un polo), lascio stare
       bool bCollapsed = true ;
       Point3d ptFirst = cnData.vCP.front() ;
-      for ( int i = 1 ; i < int( cnData.vCP.size()) ; ++i) {
+      for ( int i = 1 ; i < ssize( cnData.vCP) ; ++i) {
          if ( ! AreSamePointApprox( ptFirst, cnData.vCP[i])) {
             bCollapsed = false ;
             break ;
@@ -1611,7 +1919,7 @@ NurbsCurveCanonicalize( CNurbsData& cnData)
      // agli indici perché uso u_p-1 e u_(m-p+1), anziché u_p e u_m-p
 
      // comincio ad aumentare la molteplictià del nodo u_m-p+1
-      int nCP = int( cnData.vCP.size()) ;
+      int nCP = ssize( cnData.vCP) ;
       int nU = nCP + cnData.nDeg - 1 ;
       int nDeg = cnData.nDeg ;
       PNTVECTOR vBC ;
@@ -1807,9 +2115,9 @@ NurbsToBezierCurve( const CNurbsData& cnData)
    if ( cnData.bPeriodic  ||  cnData.bExtraKnotes)
       return nullptr ;
   // numero dei nodi
-   int nU = int( cnData.vCP.size()) + cnData.nDeg - 1 ;
+   int nU = ssize( cnData.vCP) + cnData.nDeg - 1 ;
   // controllo relazione nodi - punti di controllo
-   if ( nU != int( cnData.vU.size()))
+   if ( nU != ssize( cnData.vU))
       return nullptr ;
   // numero degli intervalli
    int nInt = nU - 2 * cnData.nDeg + 1 ;
@@ -2179,7 +2487,7 @@ CalcCurvesVoronoiDiagram( const CICURVEPVECTOR& vCrvC, ICURVEPOVECTOR& vCrvs, in
    PtrOwner<Voronoi> pVoronoiObj( new( std::nothrow) Voronoi()) ;
    if ( pVoronoiObj == nullptr)
       return false ;
-   for ( int i = 0 ; i < int( vCrvC.size()) ; i ++) {
+   for ( int i = 0 ; i < ssize( vCrvC) ; i ++) {
       if ( ! pVoronoiObj->AddCurve( vCrvC[i]))
          return false ;
    }
@@ -2206,7 +2514,7 @@ CalcCurvesMedialAxis( const CICURVEPVECTOR& vCrvC, ICURVEPOVECTOR& vCrvs, int nS
    PtrOwner<Voronoi> pVoronoiObj( new( std::nothrow) Voronoi()) ;
    if ( pVoronoiObj == nullptr)
       return false ;
-   for ( int i = 0 ; i < int( vCrvC.size()) ; i ++) {
+   for ( int i = 0 ; i < ssize( vCrvC) ; i ++) {
       if ( ! pVoronoiObj->AddCurve( vCrvC[i]))
          return false ;
    }
@@ -2267,7 +2575,7 @@ bool CalcOffsetCurves( const ICURVEPVECTOR& vpCrvs, ICURVEPOVECTOR& vCrvs, doubl
    PtrOwner<Voronoi> pVoronoiObj( new( std::nothrow) Voronoi()) ;
    if ( pVoronoiObj == nullptr)
       return false ;
-   for ( int i = 0 ; i < int( vpCrvs.size()) ; i ++) {
+   for ( int i = 0 ; i < ssize( vpCrvs) ; i ++) {
       if ( ! pVoronoiObj->AddCurve( vpCrvs[i]))
          return false ;
    }
@@ -2297,7 +2605,7 @@ bool CalcFatOffsetCurves( const ICURVEPVECTOR& vpCrvs, ICURVEPOVECTOR& vCrvs, do
    PtrOwner<Voronoi> pVoronoiObj( new( std::nothrow) Voronoi()) ;
    if ( pVoronoiObj == nullptr)
       return false ;
-   for ( int i = 0 ; i < int( vpCrvs.size()) ; i ++) {
+   for ( int i = 0 ; i < ssize( vpCrvs) ; i ++) {
       if ( ! pVoronoiObj->AddCurve( vpCrvs[i]))
          return false ;
    }
@@ -2326,3 +2634,47 @@ ResetCurveVoronoi( const ICurve& crvC)
       break ;
    }
 }
+
+//----------------------------------------------------------------------------
+bool
+GetChainedCurves( ICRVCOMPOPOVECTOR& vCrv, double dChainTol, bool bAllowInvert)
+{
+   if ( ssize( vCrv) == 1)
+      return true ;
+   ChainCurves chainCrv ;
+   // modifico direttamente le curve passate in input
+   chainCrv.Init( bAllowInvert, dChainTol, ssize( vCrv)) ;
+   for ( int c = 0 ; c < ssize( vCrv) ; ++c) {
+      Point3d ptStart, ptEnd ;
+      Vector3d vtStart, vtEnd ;
+      vCrv[c]->GetStartPoint( ptStart) ;
+      vCrv[c]->GetEndPoint( ptEnd) ;
+      vCrv[c]->GetStartDir( vtStart) ;
+      vCrv[c]->GetEndDir( vtEnd) ;
+      chainCrv.AddCurve( 1 + c, ptStart, vtStart, ptEnd, vtEnd) ;
+   }
+   INTVECTOR vIds ;
+   Point3d ptStart = ORIG ;
+   while ( chainCrv.GetChainFromNear( ptStart, false, vIds)) {
+      ICurveComposite* pFirstCrv = vCrv[abs(vIds[0]) - 1] ;
+      for ( int nId : vIds) {
+         bool bInvert = false ;
+         if ( nId < 0)
+            bInvert = true ;
+         nId = abs( nId) - 1 ;
+         if ( bInvert)
+            vCrv[nId]->Invert() ;
+         if ( ! pFirstCrv->AddCurve( Release( vCrv[nId]), true, dChainTol))
+            return false ;
+      }
+      pFirstCrv->GetEndPoint( ptStart) ;
+   }
+   // elimino gli elementi del vettore che non contengono più curve
+   int c = ssize( vCrv) ;
+   while ( c > -1) {
+      if ( IsNull( vCrv[c]))
+         vCrv.erase( vCrv.begin() + c) ;
+      --c ;
+   }
+   return true ;
+}

CurveAux.h

@@ -35,3 +35,4 @@ bool CopyExtrusion( const ICurve* pSouCrv, ICurve* pDestCrv) ;
 bool CopyThickness( const ICurve* pSouCrv, ICurve* pDestCrv) ;
 ICurveBezier* ApproxCurveBezierWithSingleCubic( const ICurve* pCrv) ;
 Voronoi* GetCurveVoronoi( const ICurve& crvC) ;
+bool GetChainedCurves( ICRVCOMPOPOVECTOR& vCrv, double dChainTol, bool bAllowInvert) ;

CurveBezier.cpp

@@ -275,20 +275,19 @@ CurveBezier::FromLine( const ICurveLine& crLine)
 {
    if ( m_nStatus != OK  ||  ! crLine.IsValid())
       return false ;
-   double dWeight = 1 ;
    int nCount = 0 ;
    Point3d ptStart ; crLine.GetStartPoint( ptStart) ;
-   SetControlPoint( nCount, ptStart, dWeight) ;
+   SetControlPoint( nCount, ptStart) ;
    ++nCount ;
    double dPart = 1. / m_nDeg ;
    for ( int i = 1 ; i < m_nDeg ; ++i) {
       double dU = i * dPart ;
       Point3d ptMid ; crLine.GetPointD1D2( dU, ICurve::FROM_MINUS, ptMid) ;
-      SetControlPoint( nCount, ptMid, dWeight) ;
+      SetControlPoint( nCount, ptMid) ;
       ++nCount ;
    }
    Point3d ptEnd ; crLine.GetEndPoint( ptEnd) ;
-   SetControlPoint( nCount, ptEnd, dWeight) ;
+   SetControlPoint( nCount, ptEnd) ;
    ++nCount ;
    return true ;
 }

CurveByApprox.cpp

@@ -202,213 +202,14 @@ CurveByApprox::CalcParameterization( void)
 bool
 CurveByApprox::CalcAkimaTangents( bool bDetectCorner)
 {
-  // pulisco i vettori delle tangenti
-   m_vPrevDer.clear() ;
-   m_vNextDer.clear() ;
-
-  // numero di punti
-   int nSize = int( m_vPnt.size()) ;
-
-  // sono necessari almeno due punti
-   if ( nSize < 2)
-      return false ;
-
-  // calcolo le derivate
-   m_vPrevDer.reserve( nSize) ;
-   m_vNextDer.reserve( nSize) ;
-  // se ci sono solo 2 punti, le tangenti devono essere dirette lungo la linea che li unisce
-   if ( nSize == 2) {
-     // non esiste derivata prima del primo punto
-      m_vPrevDer.emplace_back( 0, 0, 0) ;
-      m_vNextDer.push_back( ( m_vPnt[1] - m_vPnt[0]) / ( m_vPar[1] - m_vPar[0])) ;
-      m_vPrevDer.push_back( m_vNextDer[0]) ;
-     // non esiste derivata dopo il secondo e ultimo punto
-      m_vNextDer.emplace_back( 0, 0, 0) ;
-      return true ;
-   }
-  // verifico se curva chiusa (primo e ultimo punto coincidono)
-   bool bClosed = AreSamePointApprox( m_vPnt.front(), m_vPnt.back()) ;
-  // calcolo le derivate
-   for ( int i = 0 ; i < nSize ; ++ i) {
-      Vector3d vtPrevDer ;
-      Vector3d vtNextDer ;
-     // primo punto
-      if ( i == 0) {
-        // se curva chiusa, come precedente uso il penultimo punto
-         if ( bClosed) {
-           // se non ci sono almeno 5 punti
-            if ( nSize < 5) {
-               if ( ! CalcCircleMidDer( m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2], m_vPar[i], m_vPnt[i],
-                                        m_vPar[i+1], m_vPnt[i+1], vtNextDer))
-                  return false ;
-               vtPrevDer = vtNextDer ;
-            }
-           // altrimenti
-            else {
-               if ( ! CalcAkimaMidDer( m_vPar[nSize-3] - m_vPar[nSize-1], m_vPnt[nSize-3], m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2],
-                                       m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], bDetectCorner,
-                                       vtPrevDer, vtNextDer))
-                  return false ;
-            }
-         }
-        // altrimenti, uso arco sui primi tre punti
-         else {
-            if ( ! CalcCircleStartDer( m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], vtNextDer))
-               return false ;
-            vtPrevDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // ultimo punto
-      else if ( i == nSize - 1) {
-        // se curva chiusa, le tg devono coincidere con quelle del primo
-         if ( bClosed) {
-            vtPrevDer = m_vPrevDer[0] ;
-            vtNextDer = m_vNextDer[0] ;
-         }
-        // altrimenti, uso arco sugli ultimi tre punti
-         else {
-            if ( ! CalcCircleEndDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                       m_vPar[i], m_vPnt[i], vtPrevDer))
-               return false ;
-            vtNextDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // punti intermedi
-      else {
-        // se secondo punto
-         if ( i == 1) {
-           // se curva aperta o non ci sono almeno 5 punti
-            if ( ! bClosed  ||  nSize < 5) {
-               if ( ! CalcCircleMidDer( m_vPar[i-1], m_vPnt[i-1], m_vPar[i], m_vPnt[i],
-                                        m_vPar[i+1], m_vPnt[i+1], vtPrevDer))
-                  return false ;
-               vtNextDer = vtPrevDer ;
-            }
-           // altrimenti
-            else {
-               if ( ! CalcAkimaMidDer( m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2], m_vPar[i-1], m_vPnt[i-1],
-                                       m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], bDetectCorner,
-                                       vtPrevDer, vtNextDer))
-                  return false ;
-            }
-         }
-        // se penultimo punto
-         else if ( i == nSize - 2) {
-           // se curva aperta o non ci sono almeno 5 punti
-            if ( ! bClosed  ||  nSize < 5) {
-               if ( ! CalcCircleMidDer( m_vPar[i-1], m_vPnt[i-1], m_vPar[i], m_vPnt[i],
-                                        m_vPar[i+1], m_vPnt[i+1], vtPrevDer))
-                  return false ;
-               vtNextDer = vtPrevDer ;
-            }
-           // altrimenti
-            else {
-               if ( ! CalcAkimaMidDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                       m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[1] + m_vPar[i+1], m_vPnt[1], bDetectCorner,
-                                       vtPrevDer, vtNextDer))
-                  return false ;
-            }
-         }
-        // altrimenti
-         else {
-            if ( ! CalcAkimaMidDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                    m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                    m_vPar[i+2], m_vPnt[i+2], bDetectCorner,
-                                    vtPrevDer, vtNextDer))
-               return false ;
-         }
-      }
-     // salvo la derivata
-      m_vPrevDer.push_back( vtPrevDer) ;
-      m_vNextDer.push_back( vtNextDer) ;
-   }
-
-   return true ;
+   return ComputeAkimaTangents( bDetectCorner, m_vPar, m_vPnt, m_vPrevDer, m_vNextDer) ;
 }
 
 //----------------------------------------------------------------------------
 bool
 CurveByApprox::CalcBesselTangents( void)
 {
-  // pulisco i vettori delle tangenti
-   m_vPrevDer.clear() ;
-   m_vNextDer.clear() ;
-
-  // numero di punti
-   int nSize = int( m_vPnt.size()) ;
-
-  // sono necessari almeno due punti
-   if ( nSize < 2)
-      return false ;
-
-  // calcolo le derivate
-   m_vPrevDer.reserve( nSize) ;
-   m_vNextDer.reserve( nSize) ;
-  // se ci sono solo 2 punti, le tangenti devono essere dirette lungo la linea che li unisce
-   if ( nSize == 2) {
-     // non esiste derivata prima del primo punto
-      m_vPrevDer.emplace_back( 0, 0, 0) ;
-      m_vNextDer.push_back( ( m_vPnt[1] - m_vPnt[0]) / ( m_vPar[1] - m_vPar[0])) ;
-      m_vPrevDer.push_back( m_vNextDer[0]) ;
-     // non esiste derivata dopo il secondo e ultimo punto
-      m_vNextDer.emplace_back( 0, 0, 0) ;
-      return true ;
-   }
-  // verifico se curva chiusa (primo e ultimo punto coincidono)
-   bool bClosed = AreSamePointApprox( m_vPnt.front(), m_vPnt.back()) ;
-  // calcolo le derivate
-   for ( int i = 0 ; i < nSize ; ++ i) {
-      Vector3d vtPrevDer ;
-      Vector3d vtNextDer ;
-     // primo punto
-      if ( i == 0) {
-        // se curva chiusa, come precedente uso il penultimo punto
-         if ( bClosed) {
-            if ( ! CalcBesselMidDer( m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2], m_vPar[i], m_vPnt[i],
-                                     m_vPar[i+1], m_vPnt[i+1], vtNextDer))
-               return false ;
-            vtPrevDer = vtNextDer ;
-         }
-        // altrimenti, uso i primi tre punti
-         else {
-            if ( ! CalcBesselStartDer( m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], vtNextDer))
-               return false ;
-            vtPrevDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // ultimo punto
-      else if ( i == nSize - 1) {
-        // se curva chiusa, le tg devono coincidere con quelle del primo
-         if ( bClosed) {
-            vtPrevDer = m_vPrevDer[0] ;
-            vtNextDer = m_vNextDer[0] ;
-         }
-        // altrimenti, uso gli ultimi tre punti
-         else {
-            if ( ! CalcBesselEndDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                     m_vPar[i], m_vPnt[i], vtPrevDer))
-               return false ;
-            vtNextDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // punti intermedi
-      else {
-         if ( ! CalcBesselMidDer( m_vPar[i-1], m_vPnt[i-1], m_vPar[i], m_vPnt[i],
-                                  m_vPar[i+1], m_vPnt[i+1], vtPrevDer))
-            return false ;
-         vtNextDer = vtPrevDer ;
-      }
-     // salvo la derivata
-      m_vPrevDer.push_back( vtPrevDer) ;
-      m_vNextDer.push_back( vtNextDer) ;
-   }
-
-   return true ;
+   return ComputeBesselTangents( m_vPar, m_vPnt, m_vPrevDer, m_vNextDer) ;
 }
 
 //----------------------------------------------------------------------------

CurveByInterp.cpp

@@ -50,12 +50,42 @@ CurveByInterp::AddPoint( const Point3d& ptP)
 ICurve*
 CurveByInterp::GetCurve( int nMethod, int nType)
 {
-  // calcolo le tangenti
+   // se richieste curve di Bezier cubiche (ottenute da interpolazione con Nurbs)
+   if ( nType == CUBIC_BEZIERS_LONG) {
+       // creo la curva composita
+       PtrOwner<ICurve> pCrv ;
+       //pCrv.Set( InterpolatePointSetWithBezier( m_vPnt, 50 * EPS_SMALL, 50)) ;
+       //debug
+       pCrv.Set( InterpolatePointSetWithBezier( m_vPnt, 0.1, 100)) ;
+       if ( IsNull(pCrv)  ||  ! pCrv->IsValid())
+          return nullptr ;
+       return Release( pCrv) ;
+   }
+
+  // numero di punti
+   int nSize = int( m_vPnt.size()) ;
+
+  // sono necessari almeno due punti
+   if ( nSize < 2)
+      return nullptr ;
+
+  // calcolo le distanze tra i punti per derivarne i parametri
+   m_vPar.reserve( nSize) ;
+   double dPar = 0 ;
+   m_vPar.push_back( dPar) ;
+
+   for ( int i = 1 ; i < nSize ; ++ i) {
+      double dDist = Dist( m_vPnt[i-1], m_vPnt[i]) ;
+      dPar += dDist ;
+      m_vPar.push_back( dPar) ;
+   }
+   
+   // calcolo le tangenti
    if ( nMethod == BESSEL) {
       if ( ! CalcBesselTangents())
          return nullptr ;
    }
-   else if ( nType != CUBIC_BEZIERS_LONG) {
+   else {
       if ( ! CalcAkimaTangents( nMethod == AKIMA_CORNER))
          return nullptr ;
    }
@@ -103,16 +133,6 @@ CurveByInterp::GetCurve( int nMethod, int nType)
       return ::Release( pCrvCompo) ;
    }
 
-   // se richieste curve di Bezier cubiche (ottenute da interpolazione con Nurbs)
-   if ( nType == CUBIC_BEZIERS_LONG) {
-       // creo la curva composita
-       PtrOwner<ICurve> pCrv ;
-       pCrv.Set( InterpolatePointSetWithBezier( m_vPnt, 50 * EPS_SMALL, 50)) ;
-       if ( IsNull(pCrv)  ||  ! pCrv->IsValid())
-          return nullptr ;
-       return Release( pCrv) ;
-   }
-
    return nullptr ;
 }
 
@@ -120,233 +140,12 @@ CurveByInterp::GetCurve( int nMethod, int nType)
 bool
 CurveByInterp::CalcAkimaTangents( bool bDetectCorner)
 {
-  // pulisco i vettori dei parametri e delle tangenti
-   m_vPar.clear() ;
-   m_vPrevDer.clear() ;
-   m_vNextDer.clear() ;
-
-  // numero di punti
-   int nSize = int( m_vPnt.size()) ;
-
-  // sono necessari almeno due punti
-   if ( nSize < 2)
-      return false ;
-
-  // calcolo le distanze tra i punti per derivarne i parametri
-   m_vPar.reserve( nSize) ;
-   double dPar = 0 ;
-   m_vPar.push_back( dPar) ;
-   for ( int i = 1 ; i < nSize ; ++ i) {
-      double dDist = Dist( m_vPnt[i-1], m_vPnt[i]) ;
-      dPar += dDist ;
-      m_vPar.push_back( dPar) ;
-   }
-
-  // calcolo le derivate
-   m_vPrevDer.reserve( nSize) ;
-   m_vNextDer.reserve( nSize) ;
-  // se ci sono solo 2 punti, le tangenti devono essere dirette lungo la linea che li unisce
-   if ( nSize == 2) {
-     // non esiste derivata prima del primo punto
-      m_vPrevDer.emplace_back( 0, 0, 0) ;
-      m_vNextDer.push_back( ( m_vPnt[1] - m_vPnt[0]) / ( m_vPar[1] - m_vPar[0])) ;
-      m_vPrevDer.push_back( m_vNextDer[0]) ;
-     // non esiste derivata dopo il secondo e ultimo punto
-      m_vNextDer.emplace_back( 0, 0, 0) ;
-      return true ;
-   }
-  // verifico se curva chiusa (primo e ultimo punto coincidono)
-   bool bClosed = AreSamePointApprox( m_vPnt.front(), m_vPnt.back()) ;
-  // calcolo le derivate
-   for ( int i = 0 ; i < nSize ; ++ i) {
-      Vector3d vtPrevDer ;
-      Vector3d vtNextDer ;
-     // primo punto
-      if ( i == 0) {
-        // se curva chiusa, come precedente uso il penultimo punto
-         if ( bClosed) {
-           // se non ci sono almeno 5 punti
-            if ( nSize < 5) {
-               if ( ! CalcCircleMidDer( m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2], m_vPar[i], m_vPnt[i],
-                                        m_vPar[i+1], m_vPnt[i+1], vtNextDer))
-                  return false ;
-               vtPrevDer = vtNextDer ;
-            }
-           // altrimenti
-            else {
-               if ( ! CalcAkimaMidDer( m_vPar[nSize-3] - m_vPar[nSize-1], m_vPnt[nSize-3], m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2],
-                                       m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], bDetectCorner,
-                                       vtPrevDer, vtNextDer))
-                  return false ;
-            }
-         }
-        // altrimenti, uso arco sui primi tre punti
-         else {
-            if ( ! CalcCircleStartDer( m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], vtNextDer))
-               return false ;
-            vtPrevDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // ultimo punto
-      else if ( i == nSize - 1) {
-        // se curva chiusa, le tg devono coincidere con quelle del primo
-         if ( bClosed) {
-            vtPrevDer = m_vPrevDer[0] ;
-            vtNextDer = m_vNextDer[0] ;
-         }
-        // altrimenti, uso arco sugli ultimi tre punti
-         else {
-            if ( ! CalcCircleEndDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                       m_vPar[i], m_vPnt[i], vtPrevDer))
-               return false ;
-            vtNextDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // punti intermedi
-      else {
-        // se secondo punto
-         if ( i == 1) {
-           // se curva aperta o non ci sono almeno 5 punti
-            if ( ! bClosed  ||  nSize < 5) {
-               if ( ! CalcCircleMidDer( m_vPar[i-1], m_vPnt[i-1], m_vPar[i], m_vPnt[i],
-                                        m_vPar[i+1], m_vPnt[i+1], vtPrevDer))
-                  return false ;
-               vtNextDer = vtPrevDer ;
-            }
-           // altrimenti
-            else {
-               if ( ! CalcAkimaMidDer( m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2], m_vPar[i-1], m_vPnt[i-1],
-                                       m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], bDetectCorner,
-                                       vtPrevDer, vtNextDer))
-                  return false ;
-            }
-         }
-        // se penultimo punto
-         else if ( i == nSize - 2) {
-           // se curva aperta o non ci sono almeno 5 punti
-            if ( ! bClosed  ||  nSize < 5) {
-               if ( ! CalcCircleMidDer( m_vPar[i-1], m_vPnt[i-1], m_vPar[i], m_vPnt[i],
-                                        m_vPar[i+1], m_vPnt[i+1], vtPrevDer))
-                  return false ;
-               vtNextDer = vtPrevDer ;
-            }
-           // altrimenti
-            else {
-               if ( ! CalcAkimaMidDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                       m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[1] + m_vPar[i+1], m_vPnt[1], bDetectCorner,
-                                       vtPrevDer, vtNextDer))
-                  return false ;
-            }
-         }
-        // altrimenti
-         else {
-            if ( ! CalcAkimaMidDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                    m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                    m_vPar[i+2], m_vPnt[i+2], bDetectCorner,
-                                    vtPrevDer, vtNextDer))
-               return false ;
-         }
-      }
-     // salvo la derivata
-      m_vPrevDer.push_back( vtPrevDer) ;
-      m_vNextDer.push_back( vtNextDer) ;
-   }
-
-   return true ;
+   return ComputeAkimaTangents( bDetectCorner, m_vPar, m_vPnt, m_vPrevDer, m_vNextDer) ;
 }
 
 //----------------------------------------------------------------------------
 bool
 CurveByInterp::CalcBesselTangents( void)
 {
-  // pulisco i vettori dei parametri e delle tangenti
-   m_vPar.clear() ;
-   m_vPrevDer.clear() ;
-   m_vNextDer.clear() ;
-
-  // numero di punti
-   int nSize = int( m_vPnt.size()) ;
-
-  // sono necessari almeno due punti
-   if ( nSize < 2)
-      return false ;
-
-  // calcolo le distanze tra i punti per derivarne i parametri
-   m_vPar.reserve( nSize) ;
-   double dPar = 0 ;
-   m_vPar.push_back( dPar) ;
-   for ( int i = 1 ; i < nSize ; ++ i) {
-      double dDist = Dist( m_vPnt[i-1], m_vPnt[i]) ;
-      dPar += dDist ;
-      m_vPar.push_back( dPar) ;
-   }
-
-  // calcolo le derivate
-   m_vPrevDer.reserve( nSize) ;
-   m_vNextDer.reserve( nSize) ;
-  // se ci sono solo 2 punti, le tangenti devono essere dirette lungo la linea che li unisce
-   if ( nSize == 2) {
-     // non esiste derivata prima del primo punto
-      m_vPrevDer.emplace_back( 0, 0, 0) ;
-      m_vNextDer.push_back( ( m_vPnt[1] - m_vPnt[0]) / ( m_vPar[1] - m_vPar[0])) ;
-      m_vPrevDer.push_back( m_vNextDer[0]) ;
-     // non esiste derivata dopo il secondo e ultimo punto
-      m_vNextDer.emplace_back( 0, 0, 0) ;
-      return true ;
-   }
-  // verifico se curva chiusa (primo e ultimo punto coincidono)
-   bool bClosed = AreSamePointApprox( m_vPnt.front(), m_vPnt.back()) ;
-  // calcolo le derivate
-   for ( int i = 0 ; i < nSize ; ++ i) {
-      Vector3d vtPrevDer ;
-      Vector3d vtNextDer ;
-     // primo punto
-      if ( i == 0) {
-        // se curva chiusa, come precedente uso il penultimo punto
-         if ( bClosed) {
-            if ( ! CalcBesselMidDer( m_vPar[nSize-2] - m_vPar[nSize-1], m_vPnt[nSize-2], m_vPar[i], m_vPnt[i],
-                                     m_vPar[i+1], m_vPnt[i+1], vtNextDer))
-               return false ;
-            vtPrevDer = vtNextDer ;
-         }
-        // altrimenti, uso i primi tre punti
-         else {
-            if ( ! CalcBesselStartDer( m_vPar[i], m_vPnt[i], m_vPar[i+1], m_vPnt[i+1],
-                                       m_vPar[i+2], m_vPnt[i+2], vtNextDer))
-               return false ;
-            vtPrevDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // ultimo punto
-      else if ( i == nSize - 1) {
-        // se curva chiusa, le tg devono coincidere con quelle del primo
-         if ( bClosed) {
-            vtPrevDer = m_vPrevDer[0] ;
-            vtNextDer = m_vNextDer[0] ;
-         }
-        // altrimenti, uso gli ultimi tre punti
-         else {
-            if ( ! CalcBesselEndDer( m_vPar[i-2], m_vPnt[i-2], m_vPar[i-1], m_vPnt[i-1],
-                                     m_vPar[i], m_vPnt[i], vtPrevDer))
-               return false ;
-            vtNextDer = Vector3d( 0, 0, 0) ;
-         }
-      }
-     // punti intermedi
-      else {
-         if ( ! CalcBesselMidDer( m_vPar[i-1], m_vPnt[i-1], m_vPar[i], m_vPnt[i],
-                                  m_vPar[i+1], m_vPnt[i+1], vtPrevDer))
-            return false ;
-         vtNextDer = vtPrevDer ;
-      }
-     // salvo la derivata
-      m_vPrevDer.push_back( vtPrevDer) ;
-      m_vNextDer.push_back( vtNextDer) ;
-   }
-
-   return true ;
-}
+   return ComputeBesselTangents( m_vPar, m_vPnt, m_vPrevDer, m_vNextDer) ;
+}

DistPointCrvAux.cpp

@@ -116,6 +116,7 @@ PolishMinDistPointCurve( const Point3d& ptP, const ICurve& cCurve,
       vtDiff = ptQ - ptP ;
      // angolo tra vettore e tangente
       dTemp = vtDer1 * vtDiff ;
+      bool bEquiverse = dTemp > 0 ;
       if ( abs( dTemp) > EPS_ZERO) 
          dSqCosA = dTemp * dTemp / ( vtDer1.SqLen() * vtDiff.SqLen()) ;
       else
@@ -123,8 +124,20 @@ PolishMinDistPointCurve( const Point3d& ptP, const ICurve& cCurve,
      // stima prossimo valore del parametro (Newton : Unext = U - F(U) / F'(U))
       dPrevPar = dPar ;
       dTemp = vtDer2 * vtDiff + vtDer1.SqLen() ;
-      if ( abs( dTemp) > EPS_ZERO)
-         dPar = dPrevPar - ( vtDer1 * vtDiff) / dTemp ;
+
+      // se il coseno tra questi due vettori è troppo grande potrei aver avuto una cattiva stima iniziale
+      // provo quindi ad aggiustare a mano, anziché usare il segno suggerito da newton, che con queste premesse potrebbe divergere
+      double dCos75 = 0.2588 ;
+      if ( abs( dTemp) > EPS_ZERO) {
+         double dDelta = ( vtDer1 * vtDiff) / dTemp ;
+         if ( dSqCosA > dCos75) {
+            if ( ( bEquiverse  &&  dDelta > 0)  ||  ( ! bEquiverse  &&  dDelta < 0))
+               dDelta *= -1 ;
+            dPar = dPrevPar + dDelta ;
+         }
+         else
+            dPar = dPrevPar - dDelta ;
+      }
      // clipping parametro
       if ( dPar < approxMin.dParMin) {
          if ( approxMin.bParMinSing  &&  ! bClampedFromSing) {

DistPointCrvBezier.cpp

@@ -26,7 +26,7 @@ DistPointCrvBezier::DistPointCrvBezier( const Point3d& ptP, const ICurveBezier&
   // distanza non calcolata
    m_dDist = - 1 ;
 
-   if ( &CrvBez == nullptr  ||  ! CrvBez.IsValid())
+   if ( ! CrvBez.IsValid())
       return ;
 
   // determino tolleranza di approssimazione in base a ingombro curva
@@ -42,6 +42,17 @@ DistPointCrvBezier::DistPointCrvBezier( const Point3d& ptP, const ICurveBezier&
    if ( ! CrvBez.ApproxWithLines( dLinTol, ANG_TOL_APPROX_DEG, ICurve::APL_STD, PL))
       return ;
 
+   int nDeg = CrvBez.GetDegree() ;
+   if ( PL.GetPointNbr() < nDeg + 1) {
+      // costruisco una polilinea con un numero di curve scelto in base al grado della curva
+      PL.Clear() ; 
+      for ( int i = 0 ; i <= nDeg + 1 ; ++i) {
+         double dU = double(i) / (nDeg + 1) ;
+         Point3d ptBez ;
+         CrvBez.GetPointD1D2( dU, ICurve::Side::FROM_MINUS, ptBez) ;
+         PL.AddUPoint( dU, ptBez) ;
+      }
+   }
   // cerco la minima distanza per la polilinea
    MDCVECTOR vApproxMin ;
    if ( ! CalcMinDistPointPolyLine( ptP, PL, dLinTol, vApproxMin))

DistPointSurfBz.cpp

@@ -23,10 +23,10 @@ using namespace std ;
 
 //----------------------------------------------------------------------------
 DistPointSurfBz::DistPointSurfBz( const Point3d& ptP, const ISurfBezier& pSrfBz)
-   : m_dDist( -1), m_bIsInside( false)
+   : m_dDist( -1), m_bIsInside( false), m_bIsSurfClosed( false)
 {
   // Bezier non valida 
-   if ( &pSrfBz == nullptr  ||  ! pSrfBz.IsValid())
+   if ( ! pSrfBz.IsValid())
       return ;
   // Calcolo la distanza
    Calculate( ptP, pSrfBz) ;
@@ -37,9 +37,9 @@ void
 DistPointSurfBz::Calculate( const Point3d& ptP, const ISurfBezier& srfBz)
 {
   // Inizializzo distanza non calcolata
-   m_dDist = - 1. ;
+   m_dDist = -1 ;
 
-  // Controllo se la superficie è chiusa
+  // Controllo se la superficie è chiusa
    m_bIsSurfClosed = srfBz.IsClosed() ;
 
   // Lavoro con l'oggetto superficie trimesh di base
@@ -49,17 +49,17 @@ DistPointSurfBz::Calculate( const Point3d& ptP, const ISurfBezier& srfBz)
 
    DistPointSurfTm dpst( ptP, *pStmRef) ;
 
-  //recupero il punto a distanza minima sulla trimesh e lo raffino, prima di restituire distanza e punto minimo
+  // recupero il punto a distanza minima sulla trimesh e lo raffino, prima di restituire distanza e punto minimo
    Point3d ptMinTm ; dpst.GetMinDistPoint( ptMinTm) ;
    int nT ; dpst.GetMinDistTriaIndex( nT) ;
-  //salvo il punto corrispondente nel parametrico
+  // salvo il punto corrispondente nel parametrico
    srfBz.UnprojectPointFromStm( nT, ptMinTm, m_ptParam)  ;
   // salvo il punto a minima distanza sulla superficie e la normale alla superficie in quel punto
    srfBz.GetPointNrmD1D2( m_ptParam.x, m_ptParam.y, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS, m_ptMinDistPoint, m_vtN)  ;
 
   // salvo la distanza minima
    m_dDist = Dist( ptP, m_ptMinDistPoint) ;
-  // se il punto è sulla superficie
+  // se il punto è sulla superficie
    if ( m_dDist < EPS_SMALL) {
       m_bIsInside = false ;
       return ;
@@ -96,13 +96,14 @@ DistPointSurfBz::GetMinDistPoint( Point3d& ptMinDistPoint) const
 
 //----------------------------------------------------------------------------
 bool 
-DistPointSurfBz::GetParamPoint( Point3d& ptParamPoint) const
+DistPointSurfBz::GetParamsAtMinDistPoint( double& dU, double& dV) const
 {
   // Distanza non valida
    if ( m_dDist < -EPS_ZERO)
       return false ;
   // Distanza valida
-   ptParamPoint = m_ptParam ;
+   dU = m_ptParam.x ;
+   dV = m_ptParam.y ;
    return true ;
 }
 
@@ -116,4 +117,4 @@ DistPointSurfBz::GetNorm( Vector3d& vtN) const
   // Distanza valida
    vtN = m_vtN ;
    return true ;
-}
+}

EgtGeomKernel.vcxproj

@@ -281,6 +281,7 @@ copy $(TargetPath) \EgtProg\Dll64</Command>
     <ClCompile Include="BBox3d.cpp" />
     <ClCompile Include="BiArcs.cpp" />
     <ClCompile Include="CalcPocketing.cpp" />
+    <ClCompile Include="CalcDerivate.cpp" />
     <ClCompile Include="CAvSilhouetteSurfTm.cpp" />
     <ClCompile Include="CAvSimpleSurfFrMove.cpp" />
     <ClCompile Include="CAvToolSurfTm.cpp" />

EgtGeomKernel.vcxproj.filters

@@ -567,6 +567,9 @@
     <ClCompile Include="Trimming.cpp">
       <Filter>File di origine\GeoStriping</Filter>
     </ClCompile>
+    <ClCompile Include="CalcDerivate.cpp">
+      <Filter>File di origine\Geo</Filter>
+    </ClCompile>
   </ItemGroup>
   <ItemGroup>
     <ClInclude Include="stdafx.h">

IntersCrvCompoCrvCompo.cpp

@@ -34,6 +34,7 @@ static bool CalcATypeFromDisk( const ICurve& CurveA, double dUA, ICurve::Side nS
                                const ICurve& CurveB, double dUB, int& nType) ;
 static bool CalcATypeFromDisk2( const ICurve& CurveA, double dUA, ICurve::Side nSideA,
                                 const ICurve& CurveB, double dUB1, double dUB2, int& nType) ;
+static bool CalcSide( const IntCrvCrvInfo& Icci1, const IntCrvCrvInfo& Icci2, const ICurve* pThisCrv, const ICurve* pOtherCrv, bool bCrvAOrB, int& nType) ;
 
 //----------------------------------------------------------------------------
 IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
@@ -197,7 +198,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
               // se coincidono U e ptInt tra A e B
                if ( abs( m_Info[i].IciA[0].dU - m_Info[j].IciB[0].dU) < EPS_SMALL  &&
                     AreSamePointXYEpsilon( m_Info[i].IciA[0].ptI, m_Info[j].IciB[0].ptI, 10 * EPS_SMALL)) {
-                 // se non è alla fine di curva chiusa 
+                 // se non è alla fine di curva chiusa 
                   if ( ! bCrvAClosed  ||  abs( m_Info[j].IciA[0].dU - dCrvBSpan) > EPS_SMALL)
                     // elimino la seconda
                      EraseOtherInfo( i, j) ;
@@ -264,7 +265,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
             EraseOtherInfo( i, j) ;
             EraseCurrentInfo( i, j) ;
          }
-        // caso DET-NULL -> NULL-DET per prima curva
+        // caso DET-(NULL) -> (NULL)-DET per prima curva
          else if ( m_Info[j].IciA[kj].nPrevTy != ICCT_NULL  &&  m_Info[j].IciA[kj].nNextTy == ICCT_NULL  &&
               m_Info[i].IciA[ki].nPrevTy == ICCT_NULL  &&  m_Info[i].IciA[ki].nNextTy != ICCT_NULL) {
            // per la prima curva tengo i determinati
@@ -272,16 +273,22 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
             m_Info[j].IciA[kj].nNextTy = m_Info[i].IciA[ki].nNextTy ;
            // se overlap equiverso
             if ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq) {
-              // per la seconda curva ogni sottotipo è il duale di quello della prima
+              // per la seconda curva ogni sottotipo è il duale di quello della prima
                m_Info[i].IciB[ki].nPrevTy = GetDualIcct( m_Info[i].IciA[ki].nPrevTy) ;
                m_Info[i].IciB[ki].nNextTy = GetDualIcct( m_Info[i].IciA[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[i].bOverlap  &&  ! m_Info[i].bCBOverEq) {
-              // per la seconda curva ogni sottotipo è come quello della prima ma in posizione invertita
+              // per la seconda curva ogni sottotipo è come quello della prima ma in posizione invertita
                m_Info[i].IciB[ki].nPrevTy = m_Info[i].IciA[ki].nNextTy ;
                m_Info[i].IciB[ki].nNextTy = m_Info[i].IciA[ki].nPrevTy ;
             }
+            else {
+               if ( m_Info[i].IciB[ki].nPrevTy  == ICCT_NULL)
+                  m_Info[i].IciB[ki].nPrevTy = m_Info[j].IciB[kj].nPrevTy ;
+               if ( m_Info[i].IciB[ki].nNextTy  == ICCT_NULL)
+                  m_Info[i].IciB[ki].nNextTy = m_Info[j].IciB[kj].nNextTy ;
+            }
            // se overlap equiverso
             if ( m_Info[j].bOverlap  &&  m_Info[j].bCBOverEq) {
                m_Info[j].IciB[kj].nPrevTy = GetDualIcct( m_Info[i].IciA[ki].nPrevTy) ;
@@ -292,6 +299,13 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
                m_Info[j].IciB[kj].nPrevTy = m_Info[i].IciA[ki].nNextTy ;
                m_Info[j].IciB[kj].nNextTy = m_Info[i].IciA[ki].nPrevTy ;
             }
+            else {
+               if ( m_Info[j].IciB[kj].nPrevTy  == ICCT_NULL)
+                  m_Info[j].IciB[kj].nPrevTy = m_Info[i].IciB[ki].nPrevTy ;
+               if ( m_Info[j].IciB[kj].nNextTy  == ICCT_NULL)
+                  m_Info[j].IciB[kj].nNextTy = m_Info[i].IciB[ki].nNextTy ;
+            }
+
            // medio parametri e punti separatamente per le due curve
             MediaParamPoints( m_Info[i].IciA[ki], m_Info[j].IciA[kj], CCompoA) ;
             MediaParamPoints( m_Info[i].IciB[ki], m_Info[j].IciB[kj], CCompoB) ;
@@ -306,7 +320,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
                EraseCurrentInfo( i, j) ;
             }
          }
-        // caso NULL-DET -> DET-NULL per prima curva (possibile su inizio/fine di curva chiusa)
+        // caso (NULL)-DET -> DET-(NULL) per prima curva (possibile su inizio/fine di curva chiusa)
          else if ( m_Info[j].IciA[kj].nPrevTy == ICCT_NULL  &&  m_Info[j].IciA[kj].nNextTy != ICCT_NULL  &&
                    m_Info[i].IciA[ki].nPrevTy != ICCT_NULL  &&  m_Info[i].IciA[ki].nNextTy == ICCT_NULL) {
            // per la prima curva tengo i determinati
@@ -314,28 +328,41 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
             m_Info[j].IciA[kj].nPrevTy = m_Info[i].IciA[ki].nPrevTy ;
            // se overlap equiverso
             if ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq) {
-              // per la seconda curva ogni sottotipo è il duale di quello della prima
+              // per la seconda curva ogni sottotipo è il duale di quello della prima
                m_Info[i].IciB[ki].nPrevTy = GetDualIcct( m_Info[i].IciA[ki].nPrevTy) ;
                m_Info[i].IciB[ki].nNextTy = GetDualIcct( m_Info[i].IciA[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[i].bOverlap  &&  ! m_Info[i].bCBOverEq) {
-              // per la seconda curva ogni sottotipo è come quello della prima ma in posizione scambiata
+              // per la seconda curva ogni sottotipo è come quello della prima ma in posizione scambiata
                m_Info[i].IciB[ki].nPrevTy = m_Info[i].IciA[ki].nNextTy ;
                m_Info[i].IciB[ki].nNextTy = m_Info[i].IciA[ki].nPrevTy ;
             }
+            else {
+               if ( m_Info[i].IciB[ki].nPrevTy  == ICCT_NULL)
+                  m_Info[i].IciB[ki].nPrevTy = m_Info[j].IciB[kj].nPrevTy ;
+               if ( m_Info[i].IciB[ki].nNextTy  == ICCT_NULL)
+                  m_Info[i].IciB[ki].nNextTy = m_Info[j].IciB[kj].nNextTy ;
+            }
            // se overlap equiverso
             if ( m_Info[j].bOverlap  &&  m_Info[j].bCBOverEq) {
-              // per la seconda curva ogni sottotipo è il duale di quello della prima
+              // per la seconda curva ogni sottotipo è il duale di quello della prima
                m_Info[j].IciB[kj].nPrevTy = GetDualIcct( m_Info[i].IciA[ki].nPrevTy) ;
                m_Info[j].IciB[kj].nNextTy = GetDualIcct( m_Info[i].IciA[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[i].bOverlap  &&  ! m_Info[i].bCBOverEq) {
-              // per la seconda curva ogni sottotipo è come quello della prima ma in posizione scambiata
+              // per la seconda curva ogni sottotipo è come quello della prima ma in posizione scambiata
                m_Info[j].IciB[kj].nPrevTy = m_Info[i].IciA[ki].nNextTy ;
                m_Info[j].IciB[kj].nNextTy = m_Info[i].IciA[ki].nPrevTy ;
             }
+            else {
+               if ( m_Info[j].IciB[kj].nPrevTy  == ICCT_NULL)
+                  m_Info[j].IciB[kj].nPrevTy = m_Info[i].IciB[ki].nPrevTy ;
+               if ( m_Info[j].IciB[kj].nNextTy  == ICCT_NULL)
+                  m_Info[j].IciB[kj].nNextTy = m_Info[i].IciB[ki].nNextTy ;
+            }
+
            // medio parametri e punti separatamente per le due curve
             MediaParamPoints( m_Info[i].IciA[ki], m_Info[j].IciA[kj], CCompoA) ;
             MediaParamPoints( m_Info[i].IciB[ki], m_Info[j].IciB[kj], CCompoB) ;
@@ -370,13 +397,21 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
          }
         // caso NULL-NULL per corrente di prima curva
          else if ( m_Info[i].IciA[ki].nPrevTy == ICCT_NULL  &&  m_Info[i].IciA[ki].nNextTy == ICCT_NULL) {
-            m_Info[j].IciA[kj].nNextTy = m_Info[i].IciA[ki].nNextTy ;
+            m_Info[j].IciA[kj].nNextTy = ICCT_NULL ;
+            if ( m_Info[j].IciB[kj].nNextTy == ICCT_NULL)
+               m_Info[j].IciB[kj].nNextTy = m_Info[i].IciB[ki].nNextTy ;
+            if ( m_Info[j].IciB[kj].nPrevTy == ICCT_NULL)
+               m_Info[j].IciB[kj].nPrevTy = m_Info[i].IciB[ki].nPrevTy ;
            // cancello l'intersezione corrente (non aggiunge nulla rispetto alla precedente)
             EraseCurrentInfo( i, j) ;
          }
         // caso NULL-NULL per precedente di prima curva
          else if ( m_Info[j].IciA[kj].nPrevTy == ICCT_NULL  &&  m_Info[j].IciA[kj].nNextTy == ICCT_NULL) {
-            m_Info[i].IciA[0].nPrevTy = m_Info[j].IciA[0].nPrevTy ;
+            m_Info[i].IciA[ki].nPrevTy = ICCT_NULL ;
+            if ( m_Info[i].IciB[ki].nPrevTy == ICCT_NULL)
+               m_Info[i].IciB[ki].nPrevTy = m_Info[j].IciB[kj].nPrevTy ;
+            if ( m_Info[i].IciB[ki].nNextTy == ICCT_NULL)
+               m_Info[i].IciB[ki].nNextTy = m_Info[j].IciB[kj].nNextTy ;
            // cancello l'intersezione precedente (non aggiunge nulla rispetto alla corrente)
             EraseOtherInfo( i, j) ;
          }
@@ -410,7 +445,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
             EraseOtherInfo( i, j) ;
             EraseCurrentInfo( i, j) ;
          }
-        // caso DET-NULL -> NULL-DET per seconda curva
+        // caso DET-(NULL) -> (NULL)-DET per seconda curva
          else if ( m_Info[j].IciB[kj].nPrevTy != ICCT_NULL  &&  m_Info[j].IciB[kj].nNextTy == ICCT_NULL  &&
                    m_Info[i].IciB[ki].nPrevTy == ICCT_NULL  &&  m_Info[i].IciB[ki].nNextTy != ICCT_NULL) {
            // per la seconda curva tengo i determinati
@@ -418,25 +453,25 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
             m_Info[j].IciB[kj].nNextTy = m_Info[i].IciB[ki].nNextTy ;
            // se overlap equiverso
             if ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq) {
-              // per la prima curva ogni sottotipo è il duale di quello della seconda
+              // per la prima curva ogni sottotipo è il duale di quello della seconda
                m_Info[i].IciA[ki].nPrevTy = GetDualIcct( m_Info[i].IciB[ki].nPrevTy) ;
                m_Info[i].IciA[ki].nNextTy = GetDualIcct( m_Info[i].IciB[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[i].bOverlap  &&  ! m_Info[i].bCBOverEq) {
-              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
+              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
                m_Info[i].IciA[ki].nPrevTy = m_Info[i].IciB[ki].nNextTy ;
                m_Info[i].IciA[ki].nNextTy = m_Info[i].IciB[ki].nPrevTy ;
             }
             // se overlap equiverso
             if ( m_Info[j].bOverlap  &&  m_Info[j].bCBOverEq) {
-              // per la prima curva ogni sottotipo è il duale di quello della seconda
+              // per la prima curva ogni sottotipo è il duale di quello della seconda
                m_Info[j].IciA[kj].nPrevTy = GetDualIcct( m_Info[i].IciB[ki].nPrevTy) ;
                m_Info[j].IciA[kj].nNextTy = GetDualIcct( m_Info[i].IciB[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[j].bOverlap  &&  ! m_Info[j].bCBOverEq) {
-              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
+              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
                m_Info[j].IciA[kj].nPrevTy = m_Info[i].IciB[ki].nNextTy ;
                m_Info[j].IciA[kj].nNextTy = m_Info[i].IciB[ki].nPrevTy ;
             }
@@ -454,7 +489,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
                EraseCurrentInfo( i, j) ;
             }
          }
-        // caso NULL-DET -> DET-NULL per seconda curva (possibile su inizio/fine di curva chiusa)
+        // caso (NULL)-DET -> DET-(NULL) per seconda curva (possibile su inizio/fine di curva chiusa)
          else if ( m_Info[j].IciB[kj].nPrevTy == ICCT_NULL  &&  m_Info[j].IciB[kj].nNextTy != ICCT_NULL  &&
                    m_Info[i].IciB[ki].nPrevTy != ICCT_NULL  &&  m_Info[i].IciB[ki].nNextTy == ICCT_NULL) {
            // per la seconda curva tengo i determinati
@@ -462,25 +497,25 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
             m_Info[j].IciB[kj].nPrevTy = m_Info[i].IciB[ki].nPrevTy ;
            // se overlap equiverso
             if ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq) {
-              // per la prima curva ogni sottotipo è il duale di quello della seconda
+              // per la prima curva ogni sottotipo è il duale di quello della seconda
                m_Info[i].IciA[ki].nPrevTy = GetDualIcct( m_Info[i].IciB[ki].nPrevTy) ;
                m_Info[i].IciA[ki].nNextTy = GetDualIcct( m_Info[i].IciB[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[i].bOverlap  &&  ! m_Info[i].bCBOverEq) {
-              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
+              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
                m_Info[i].IciA[ki].nPrevTy = m_Info[i].IciB[ki].nNextTy ;
                m_Info[i].IciA[ki].nNextTy = m_Info[i].IciB[ki].nPrevTy ;
             }
             // se overlap equiverso
             if ( m_Info[j].bOverlap  &&  m_Info[j].bCBOverEq) {
-              // per la prima curva ogni sottotipo è il duale di quello della seconda
+              // per la prima curva ogni sottotipo è il duale di quello della seconda
                m_Info[j].IciA[kj].nPrevTy = GetDualIcct( m_Info[i].IciB[ki].nPrevTy) ;
                m_Info[j].IciA[kj].nNextTy = GetDualIcct( m_Info[i].IciB[ki].nNextTy) ;
             }
            // se altrimenti overlap controverso
             else if ( m_Info[j].bOverlap  &&  ! m_Info[j].bCBOverEq) {
-              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
+              // per la prima curva ogni sottotipo è come quello della seconda ma in posizione scambiata
                m_Info[j].IciA[kj].nPrevTy = m_Info[i].IciB[ki].nNextTy ;
                m_Info[j].IciA[kj].nNextTy = m_Info[i].IciB[ki].nPrevTy ;
             }
@@ -518,17 +553,21 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
          }
         // caso NULL-NULL per corrente di seconda curva
          else if ( m_Info[i].IciB[ki].nPrevTy == ICCT_NULL  &&  m_Info[i].IciB[ki].nNextTy == ICCT_NULL) {
-            m_Info[j].IciB[kj].nNextTy = m_Info[i].IciB[ki].nNextTy ;
-            if ( m_Info[i].IciA[ki].dU > m_Info[j].IciA[kj].dU + EPS_PARAM)
+            m_Info[j].IciB[kj].nNextTy = ICCT_NULL ;
+            if ( m_Info[j].IciA[kj].nNextTy == ICCT_NULL)
                m_Info[j].IciA[kj].nNextTy = m_Info[i].IciA[ki].nNextTy ;
+            if ( m_Info[j].IciA[kj].nPrevTy == ICCT_NULL)
+               m_Info[j].IciA[kj].nPrevTy = m_Info[i].IciA[ki].nPrevTy ;
            // cancello l'intersezione corrente (non aggiunge nulla rispetto alla precedente)
             EraseCurrentInfo( i, j) ;
          }
         // caso NULL-NULL per precedente di seconda curva
          else if ( m_Info[j].IciB[kj].nPrevTy == ICCT_NULL  &&  m_Info[j].IciB[kj].nNextTy == ICCT_NULL) {
-            m_Info[i].IciB[0].nPrevTy = m_Info[j].IciB[0].nPrevTy ;
-            if ( m_Info[j].IciA[0].dU < m_Info[i].IciA[0].dU - EPS_PARAM)
-               m_Info[i].IciA[0].nPrevTy = m_Info[j].IciA[0].nPrevTy ;
+            m_Info[i].IciB[ki].nPrevTy = ICCT_NULL ;
+            if ( m_Info[i].IciA[ki].nPrevTy == ICCT_NULL)
+               m_Info[i].IciA[ki].nPrevTy = m_Info[j].IciA[kj].nPrevTy ;
+            if ( m_Info[i].IciA[ki].nNextTy == ICCT_NULL)
+               m_Info[i].IciA[ki].nNextTy = m_Info[j].IciA[kj].nNextTy ;
            // cancello l'intersezione precedente (non aggiunge nulla rispetto alla corrente)
             EraseOtherInfo( i, j) ;
          }
@@ -537,16 +576,16 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
   // ripristino ordinamento su prima curva
    stable_sort( m_Info.begin(), m_Info.end(), SortGreaterA) ;
 
-  // verifico se sono rimaste delle intersezioni di tipo non definito sulla curva A e cerco di risolverle per continuità
+  // verifico se sono rimaste delle intersezioni di tipo non definito sulla curva A e cerco di risolverle per continuità
    for ( int i = 0 ; i < m_nNumInters ; ++ i) {
-     // se il tipo di accostamento per la curva A non è definito
+     // se il tipo di accostamento per la curva A non è definito
       if ( m_Info[i].IciA[0].nPrevTy == ICCT_NULL) {
          if ( i > 0  ||  ( bCrvAClosed  &&  ! bAutoInters)) {
             int j = ( i > 0 ? i - 1 : m_nNumInters - 1) ;
             m_Info[i].IciA[0].nPrevTy = ( m_Info[j].bOverlap ? m_Info[j].IciA[1].nNextTy : m_Info[j].IciA[0].nNextTy) ;
          }
       }
-     // se il tipo di allontanamento per la curva A non è definito
+     // se il tipo di allontanamento per la curva A non è definito
       int ki = ( m_Info[i].bOverlap ? 1 : 0) ;
       if ( m_Info[i].IciA[ki].nNextTy == ICCT_NULL) {
          if ( i < m_nNumInters - 1  ||  ( bCrvAClosed  &&  ! bAutoInters)) {
@@ -556,21 +595,22 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
       }
    }
 
-  // verifico se sono rimaste delle intersezioni di tipo non definito sulla curva B e cerco di risolverle per continuità
+  // verifico se sono rimaste delle intersezioni di tipo non definito sulla curva B e cerco di risolverle per continuità
    stable_sort( m_Info.begin(), m_Info.end(), SortGreaterB) ;
    for ( int i = 0 ; i < m_nNumInters ; ++ i) {
-     // se il tipo di accostamento per la curva B non è definito
+     // se il tipo di accostamento per la curva B non è definito
       if ( m_Info[i].IciB[0].nPrevTy == ICCT_NULL) {
          if ( i > 0  ||  ( bCrvBClosed  &&  ! bAutoInters)) {
             int j = ( i > 0 ? i - 1 : m_nNumInters - 1) ;
-            m_Info[i].IciB[0].nPrevTy = ( m_Info[j].bOverlap && ! m_Info[j].bCBOverEq ? m_Info[j].IciB[1].nNextTy : m_Info[j].IciB[0].nNextTy) ;
+            m_Info[i].IciB[0].nPrevTy = ( m_Info[j].bOverlap ? m_Info[j].IciB[1].nNextTy : m_Info[j].IciB[0].nNextTy) ;
          }
       }
-     // se il tipo di allontanamento per la curva B non è definito
-      if ( m_Info[i].IciB[0].nNextTy == ICCT_NULL) {
+     // se il tipo di allontanamento per la curva B non è definito
+      int ki = ( m_Info[i].bOverlap ? 1 : 0) ;
+      if ( m_Info[i].IciB[ki].nNextTy == ICCT_NULL) {
          if ( i < m_nNumInters - 1  ||  ( bCrvBClosed  &&  ! bAutoInters)) {
             int j = ( i < m_nNumInters - 1 ? i + 1 : 0) ;
-            m_Info[i].IciB[0].nNextTy = ( m_Info[j].bOverlap && ! m_Info[j].bCBOverEq ? m_Info[j].IciB[1].nPrevTy : m_Info[j].IciB[0].nPrevTy) ;
+            m_Info[i].IciB[ki].nNextTy = m_Info[j].IciB[0].nPrevTy ;
          }
       }
    }
@@ -581,7 +621,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
    for ( int i = 0 ; i < m_nNumInters ; ++ i) {
      // in assenza di overlap
       if ( ! m_Info[i].bOverlap) {
-        // se il tipo di accostamento per la curva A non è definito
+        // se il tipo di accostamento per la curva A non è definito
          if ( m_Info[i].IciA[0].nPrevTy == ICCT_NULL) {
            // devo studiare un intorno dell'intersezione
             int nType ;
@@ -590,7 +630,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
               // aggiorno il tipo di accostamento della curva A alla curva B
                m_Info[i].IciA[0].nPrevTy = nType ;
          }
-        // se il tipo di accostamento per la curva B non è definito
+        // se il tipo di accostamento per la curva B non è definito
          if ( m_Info[i].IciB[0].nPrevTy == ICCT_NULL) {
            // devo studiare un intorno dell'intersezione
             int nType ;
@@ -599,7 +639,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
               // aggiorno il tipo di accostamento della curva B alla curva A
                m_Info[i].IciB[0].nPrevTy = nType ;
          }
-        // se il tipo di allontanamento per la curva A non è definito
+        // se il tipo di allontanamento per la curva A non è definito
          if ( m_Info[i].IciA[0].nNextTy == ICCT_NULL) {
            // devo studiare un intorno dell'intersezione
             int nType ;
@@ -608,7 +648,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
               // aggiorno il tipo di allontanamento della curva A dalla curva B
                m_Info[i].IciA[0].nNextTy = nType ;
          }
-        // se il tipo di allontanamento per la curva B non è definito
+        // se il tipo di allontanamento per la curva B non è definito
          if ( m_Info[i].IciB[0].nNextTy == ICCT_NULL) {
            // devo studiare un intorno dell'intersezione
             int nType ;
@@ -620,66 +660,68 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
       }
      // in presenza di overlap
       else {
-        // se il tipo di accostamento è non definito per la curva A
+        // se il tipo di accostamento è non definito per la curva A
          if ( m_Info[i].IciA[0].nPrevTy == ICCT_NULL) {
-           // se autointersezione con spike
-            double dDeltaU = abs( m_Info[i].IciA[0].dU - m_Info[i].IciB[0].dU) ;
-            if ( bAutoInters  &&  ( dDeltaU < EPS_PARAM  ||  ( bCrvAClosed  &&  abs( dDeltaU - dCrvASpan) < EPS_PARAM))) {
-              // è punta di spike
-               int nType = ICCT_SPK ;
-              // aggiorno il tipo di allontanamento della curva A dalla curva B
+            // devo studiare un intorno dell'intersezione
+            int nType ;
+            if ( CalcATypeFromDisk( CCompoA, m_Info[i].IciA[0].dU, ICurve::FROM_MINUS,
+                                    CCompoB, m_Info[i].IciB[0].dU, nType)) {
+               // aggiorno il tipo di accostamento della curva A alla curva B
                m_Info[i].IciA[0].nPrevTy = nType ;
-              // aggiorno anche il corrispondente tipo sulla curva B
+               // aggiorno anche il corrispondente tipo sulla curva B
                if ( m_Info[i].bCBOverEq)
                   m_Info[i].IciB[0].nPrevTy = GetDualIcct( nType) ;
                else
                   m_Info[i].IciB[0].nNextTy = nType ;
             }
-           // caso standard
-            else {
-              // devo studiare un intorno dell'intersezione
-               int nType ;
-               if ( CalcATypeFromDisk( CCompoA, m_Info[i].IciA[0].dU, ICurve::FROM_MINUS,
-                                       CCompoB, m_Info[i].IciB[0].dU, nType)) {
-                 // aggiorno il tipo di accostamento della curva A alla curva B
-                  m_Info[i].IciA[0].nPrevTy = nType ;
-                 // aggiorno anche il corrispondente tipo sulla curva B
-                  if ( m_Info[i].bCBOverEq)
-                     m_Info[i].IciB[0].nPrevTy = GetDualIcct( nType) ;
-                  else
-                     m_Info[i].IciB[0].nNextTy = nType ;
-               }
-            }
          }
-        // se il tipo di allontanamento è non definito per la curva A
+        // se il tipo di allontanamento è non definito per la curva A
          if ( m_Info[i].IciA[1].nNextTy == ICCT_NULL) {
-           // se autointersezione con spike
-            double dDeltaU = abs( m_Info[i].IciA[1].dU - m_Info[i].IciB[1].dU) ;
-            if ( bAutoInters  &&  ( dDeltaU < EPS_PARAM  ||  ( bCrvAClosed  &&  abs( dDeltaU - dCrvASpan) < EPS_PARAM))) {
-              // è punta di spike
-               int nType = ICCT_SPK ;
-              // aggiorno il tipo di allontanamento della curva A dalla curva B
+            // devo studiare un intorno dell'intersezione
+            int nType ;
+            if ( CalcATypeFromDisk( CCompoA, m_Info[i].IciA[1].dU, ICurve::FROM_PLUS,
+                                    CCompoB, m_Info[i].IciB[1].dU, nType)) {
+               // aggiorno il tipo di allontanamento della curva A dalla curva B
                m_Info[i].IciA[1].nNextTy = nType ;
-              // aggiorno anche il corrispondente tipo sulla curva B
+               // aggiorno anche il corrispondente tipo sulla curva B
                if ( m_Info[i].bCBOverEq)
                   m_Info[i].IciB[1].nNextTy = GetDualIcct( nType) ;
                else
                   m_Info[i].IciB[1].nPrevTy = nType ;
             }
-           // caso standard
-            else {
-              // devo studiare un intorno dell'intersezione
-               int nType ;
-               if ( CalcATypeFromDisk( CCompoA, m_Info[i].IciA[1].dU, ICurve::FROM_PLUS,
-                                       CCompoB, m_Info[i].IciB[1].dU, nType)) {
-                 // aggiorno il tipo di allontanamento della curva A dalla curva B
-                  m_Info[i].IciA[1].nNextTy = nType ;
-                 // aggiorno anche il corrispondente tipo sulla curva B
-                  if ( m_Info[i].bCBOverEq)
-                     m_Info[i].IciB[1].nNextTy = GetDualIcct( nType) ;
-                  else
-                     m_Info[i].IciB[1].nPrevTy = nType ;
-               }
+         }
+      }
+   }
+
+   if ( bAutoInters) {
+      // controllo l'eventuale presenza di spike
+      for ( int i = 0 ; i < m_nNumInters ; ++ i) {
+         if ( m_Info[i].bOverlap  &&  ! m_Info[i].bCBOverEq) {
+           // se secondo tratto dello spike
+            double dDeltaU = abs( m_Info[i].IciA[0].dU - m_Info[i].IciB[0].dU) ;
+            if ( dDeltaU < EPS_PARAM  ||  ( bCrvAClosed  &&  abs( dDeltaU - dCrvASpan) < EPS_PARAM)) {
+               // è punta di spike
+               int nType = ICCT_SPK ;
+               // aggiorno il tipo di allontanamento della curva A dalla curva B
+               m_Info[i].IciA[0].nPrevTy = nType ;
+               // aggiorno anche il corrispondente tipo sulla curva B
+               if ( m_Info[i].bCBOverEq)
+                  m_Info[i].IciB[0].nPrevTy = GetDualIcct( nType) ;
+               else
+                  m_Info[i].IciB[0].nNextTy = nType ;
+            }
+           // se primo tratto dello spike
+            dDeltaU = abs( m_Info[i].IciA[1].dU - m_Info[i].IciB[1].dU) ;
+            if ( dDeltaU < EPS_PARAM  ||  ( bCrvAClosed  &&  abs( dDeltaU - dCrvASpan) < EPS_PARAM)) {
+               // è punta di spike
+               int nType = ICCT_SPK ;
+               // aggiorno il tipo di allontanamento della curva A dalla curva B
+               m_Info[i].IciA[1].nNextTy = nType ;
+               // aggiorno anche il corrispondente tipo sulla curva B
+               if ( m_Info[i].bCBOverEq)
+                  m_Info[i].IciB[1].nNextTy = GetDualIcct( nType) ;
+               else
+                  m_Info[i].IciB[1].nPrevTy = nType ;
             }
          }
       }
@@ -693,10 +735,10 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
         // se prima curva aperta, salto alla prossima
          if ( ! CCompoA.IsClosed())
             continue ;
-        // è chiusa quindi prendo l'ultima
+        // è chiusa quindi prendo l'ultima
          j = m_nNumInters - 1 ;
       }
-     // se i due indici coincidono, c'è una sola intersezione e posso uscire
+     // se i due indici coincidono, c'è una sola intersezione e posso uscire
       if ( i == j)
           break ;
      // assegno sottoindici (considero solo intersezioni overlap)
@@ -724,7 +766,7 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
       }
    }
 
-  // se è rimasta una sola intersezione con overlap anche ad entrambi gli estremi e le curve sono chiuse
+  // se è rimasta una sola intersezione con overlap anche ad entrambi gli estremi e le curve sono chiuse
    if ( m_nNumInters == 1  &&  CCompoA.IsClosed()  &&  CCompoB.IsClosed()  &&
         m_Info[0].IciA[0].nPrevTy == ICCT_ON  &&  m_Info[0].IciA[1].nNextTy == ICCT_ON  &&
         GetCrvBDirANext( m_Info[0]) == ICCT_ON  &&  GetCrvBDirAPrev( m_Info[0]) == ICCT_ON) {
@@ -738,10 +780,103 @@ IntersCrvCompoCrvCompo::IntersCrvCompoCrvCompo( const ICurveComposite& CCompoA,
       m_Info[0].IciB[1].dU = ( m_Info[0].bCBOverEq ? dEnd : dStart) ;
    }
 
+   // verifico se una curva ha tutte le info e queste info sono coerenti
+   if ( m_nNumInters > 1) {
+      bool bCoherent = true ;
+      INTVECTOR vIncoherenceWithPrev ;
+      // salvo eventuali incoerenze col precedente
+      for ( int i = bCrvAClosed ? 0 : 1 ; i < m_nNumInters ; ++i) {
+         int j = i == 0 ? m_nNumInters - 1 : i - 1 ;
+         int kj = m_Info[j].bOverlap ? 1 : 0 ;
+         if ( m_Info[j].IciA[kj].nNextTy != m_Info[i].IciA[0].nPrevTy &&
+              m_Info[j].IciA[kj].nNextTy != ICCT_SPK  && m_Info[i].IciA[0].nPrevTy != ICCT_SPK) {
+            vIncoherenceWithPrev.push_back( i) ;
+            bCoherent = false ;
+         }
+      }
+      // incoerenze sulla curva A
+      if ( ! bCoherent) {
+         for ( int i : vIncoherenceWithPrev) {
+            int j = i == 0 ? m_nNumInters - 1 : i - 1 ;
+            int kj = m_Info[j].bOverlap ? 1 : 0 ;
+            int nType = 0 ;
+            CalcSide( m_Info[j], m_Info[i], &CCompoA, &CCompoB, true, nType) ;
+            if ( nType == MDS_LEFT)
+               nType = ICCT_IN ;
+            else if ( nType == MDS_RIGHT)
+               nType = ICCT_OUT ;
+            m_Info[i].IciA[0].nPrevTy = nType ;
+            m_Info[j].IciA[kj].nNextTy = nType ;
+         }
+      }
+
+      stable_sort( m_Info.begin(), m_Info.end(), SortGreaterB) ;
+      bCoherent = true ;
+      vIncoherenceWithPrev.clear() ;
+      // salvo eventuali incoerenze col precedente
+      for ( int i = bCrvBClosed ? 0 : 1 ; i < m_nNumInters ; ++i) {
+         int j = i == 0 ? m_nNumInters - 1 : i - 1 ;
+         int kj = m_Info[j].bOverlap ? 1 : 0 ;
+         if ( m_Info[j].IciB[kj].nNextTy != m_Info[i].IciB[0].nPrevTy  &&
+              m_Info[j].IciB[kj].nNextTy != ICCT_SPK  && m_Info[i].IciB[0].nPrevTy != ICCT_SPK) {
+            vIncoherenceWithPrev.push_back( i) ;
+            bCoherent = false ;
+         }
+      }
+      // incoerenze sulla curva B
+      if ( ! bCoherent) {
+         for ( int i : vIncoherenceWithPrev) {
+            int j = i == 0 ? m_nNumInters - 1 : i - 1 ;
+            int kj = m_Info[j].bOverlap ? 1 : 0 ;
+            int nType = 0 ;
+            CalcSide( m_Info[j], m_Info[i], &CCompoB, &CCompoA, false, nType) ;
+            if ( nType == MDS_LEFT)
+               nType = ICCT_IN ;
+            else if ( nType == MDS_RIGHT)
+               nType = ICCT_OUT ;
+            m_Info[i].IciB[0].nPrevTy = nType ;
+            m_Info[j].IciB[kj].nNextTy = nType ;
+         }
+      }
+   }
+   else {
+      // posso completare A guardando B e viceversa
+   }
+
+   stable_sort( m_Info.begin(), m_Info.end(), SortGreaterA) ;
   // verifiche su intersezioni in zone non-manifold
    OrderNonManifoldInters( m_Info, CCompoA, CCompoB) ;
 }
 
+//----------------------------------------------------------------------------
+static bool
+CalcSide( const IntCrvCrvInfo& Icci1, const IntCrvCrvInfo& Icci2,const ICurve* pThisCrv, const ICurve* pOtherCrv, bool bCrvAOrB, int& nType)
+{
+   // calcolo tra l'intersezione 1 e 2 se la curva sta dentro o fuori
+   int kj = Icci1.bOverlap ? 1 : 0 ;
+   double dU = 0 ;
+   if ( bCrvAOrB) {
+      // se precedente minore del successivo faccio la media
+      if ( Icci1.IciA[kj].dU < Icci2.IciA[0].dU)
+         dU = (Icci2.IciA[0].dU + Icci1.IciA[kj].dU) / 2 ;
+      // altrimenti guardo tra lo start e il successivo
+      else
+         dU = (Icci2.IciA[0].dU + 0.) / 2 ;
+   }
+   else {
+      // se precedente minore del successivo faccio la media
+      if ( Icci1.IciB[kj].dU < Icci2.IciB[0].dU)
+         dU = (Icci2.IciB[0].dU + Icci1.IciB[kj].dU) / 2 ;
+      // altrimenti guardi tra lo start e il successivo
+      else 
+         dU = (Icci2.IciB[0].dU + 0.) / 2 ;
+   }
+   Point3d ptTest ; pThisCrv->GetPointD1D2( dU, ICurve::FROM_MINUS, ptTest) ;
+   DistPointCurve dpc( ptTest, *pOtherCrv) ;
+   dpc.GetSideAtMinDistPoint( 0, Z_AX, nType) ;
+   return true ;
+}
+
 //----------------------------------------------------------------------------
 bool
 IntersCrvCompoCrvCompo::IntersSimpleCurves( const ICurve& CurveA, int nA, const ICurve& CurveB, int nB,
@@ -879,7 +1014,7 @@ SortGreaterB( const IntCrvCrvInfo& aInfo1, const IntCrvCrvInfo& aInfo2)
 bool
 OrderNonManifoldInters( ICCIVECTOR& Info, const ICurve& CurveA, const ICurve& CurveB)
 {
-  // questi controlli sono validi solo se la curva B è chiusa
+  // questi controlli sono validi solo se la curva B è chiusa
    if ( ! CurveB.IsClosed())
       return false ;
 
@@ -908,7 +1043,7 @@ OrderNonManifoldInters( ICCIVECTOR& Info, const ICurve& CurveA, const ICurve& Cu
         // se prima curva aperta, salto alla prossima
          if ( ! CurveA.IsClosed())
             continue ;
-        // è chiusa quindi prendo l'ultima
+        // è chiusa quindi prendo l'ultima
          j = nNumInters - 1 ;
       }
      // se i due indici coincidono, salto
@@ -998,7 +1133,7 @@ MediaParamPoints( IntCrvInfo& Ici1, IntCrvInfo& Ici2, const ICurveComposite& crv
    //Point3d ptMed ;
    //if ( ! crvCompo.GetPointD1D2( dUmed, ICurve::FROM_MINUS, ptMed))
    //   return ;
-  // verifico che non sia più lontano dei punti originali dal loro medio
+  // verifico che non sia più lontano dei punti originali dal loro medio
    //if ( SqDist( Ici1.ptI, ptMed) > dSqDist / 4)
    //   return ;
   // medio i parametri
@@ -1012,13 +1147,13 @@ MediaParamPoints( IntCrvInfo& Ici1, IntCrvInfo& Ici2, const ICurveComposite& crv
 static int
 GetCrvBDirAPrev( IntCrvCrvInfo& Icci)
 {
-  // non è overlap, è il prev del primo punto 
+  // non è overlap, è il prev del primo punto 
    if ( ! Icci.bOverlap)
       return Icci.IciB[0].nPrevTy ;
-  // è overlap equiverso, è il prev del primo punto
+  // è overlap equiverso, è il prev del primo punto
    if ( Icci.bCBOverEq)
       return Icci.IciB[0].nPrevTy ;
-  // è overlap controverso, è il next del primo punto
+  // è overlap controverso, è il next del primo punto
    return Icci.IciB[0].nNextTy ;
 }
 
@@ -1028,13 +1163,13 @@ GetCrvBDirAPrev( IntCrvCrvInfo& Icci)
 static int
 GetCrvBDirANext( IntCrvCrvInfo& Icci)
 {
-  // non è overlap, è il next del primo punto 
+  // non è overlap, è il next del primo punto 
    if ( ! Icci.bOverlap)
       return Icci.IciB[0].nNextTy ;
-  // è overlap equiverso, è il next del secondo punto
+  // è overlap equiverso, è il next del secondo punto
    if ( Icci.bCBOverEq)
       return Icci.IciB[1].nNextTy ;
-  // è overlap controverso, è il prev del secondo punto
+  // è overlap controverso, è il prev del secondo punto
    return Icci.IciB[1].nPrevTy ;
 }
 
@@ -1078,7 +1213,7 @@ CalcATypeFromDisk( const ICurve& CurveA, double dUA, ICurve::Side nSideA,
       return false ;
    if ( dAngBpDeg < 0)
       dAngBpDeg += ANG_FULL ;
-  // se l'angolo di DirA è compreso tra DirBn e DirBp (muovendosi in senso CCW) allora è IN, altrimenti OUT
+  // se l'angolo di DirA è compreso tra DirBn e DirBp (muovendosi in senso CCW) allora è IN, altrimenti OUT
    if ( dAngADeg > 0  &&  dAngADeg < dAngBpDeg)
       nType = ICCT_IN ; 
    else
@@ -1151,12 +1286,12 @@ CalcATypeFromDisk2( const ICurve& CurveA, double dUA, ICurve::Side nSideA,
    // li ordino in senso crescente
    if ( dAngB1pDeg > dAngB2pDeg)
       swap( dAngB1pDeg, dAngB2pDeg) ;
-  // se non ci sono variazioni angolari significative, non posso decidere alcunché
+  // se non ci sono variazioni angolari significative, non posso decidere alcunché
    const double MIN_DEV_ANG = 2 ;
    if ( abs( DiffAngle( dAngB1pDeg + ANG_STRAIGHT, 0)) < MIN_DEV_ANG  &&
         abs( DiffAngle( dAngB2pDeg + ANG_STRAIGHT, dAngB2nDeg)) < MIN_DEV_ANG)
       return false ;
-  // IN è tra 0 e dAngB1pDeg e tra dAngB2nDeg e dAngB2pDeg
+  // IN è tra 0 e dAngB1pDeg e tra dAngB2nDeg e dAngB2pDeg
    if ( ( dAngADeg > 0  &&  dAngADeg < dAngB1pDeg)  ||
         ( dAngADeg > dAngB2nDeg  &&  dAngADeg < dAngB2pDeg))
       nType = ICCT_IN ; 

IntersCurveCurve.cpp

@@ -47,7 +47,7 @@ IntersCurveCurve::IntersCurveCurve( const ICurve& CurveA, const ICurve& CurveB,
 
   // ciclo sulle curve per verificare se da approssimare
    for ( int i = 0 ; i < 2 ; ++ i) {
-     // se curva è arco da approssimare oppure è curva di Bezier
+     // se curva è arco da approssimare oppure è curva di Bezier
       if ( ( m_pCurve[i]->GetType() == CRV_ARC  &&  IsArcToApprox( *m_pCurve[i]))  ||
            m_pCurve[i]->GetType() == CRV_BEZIER) {
         // approssimo con rette
@@ -127,7 +127,7 @@ IntersCurveCurve::IsArcToApprox( const ICurve& Curve)
    const CurveArc* pArc = GetBasicCurveArc( &Curve) ;
    if ( pArc == nullptr)
       return false ;
-  // verifico se non è nel piano XY o ha più di un giro al centro
+  // verifico se non è nel piano XY o ha più di un giro al centro
    return ( ( ! pArc->GetNormVersor().IsZplus()  &&  ! pArc->GetNormVersor().IsZminus())  ||
             abs( pArc->GetAngCenter()) > ANG_FULL + EPS_ANG_ZERO) ;
 }
@@ -252,10 +252,10 @@ IntersCurveCurve::CrvCompoCrvCompoCalculate( const ICurve& CurveA, const ICurve&
 bool
 IntersCurveCurve::AdjustIntersParams( bool bAdjCrvA, bool bAdjCrvB)
 {
-  // se non ci sono intersezioni, non va fatto alcunché
+  // se non ci sono intersezioni, non va fatto alcunché
    if ( m_Info.empty())
       return true ;
-  // se le curve originali non sono state approssimate, non va fatto alcunché
+  // se le curve originali non sono state approssimate, non va fatto alcunché
    if ( ! bAdjCrvA  &&  ! bAdjCrvB)
       return true ;
   // procedo ad aggiustare
@@ -296,8 +296,8 @@ int
 IntersCurveCurve::GetInters3DCount( void)
 {
    int nCount = 0 ;
-   for( int i = 0 ; i < m_nIntersCount ; ++i) {
-      if( ! m_Info[i].bOverlap  ||  ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq)) {
+   for ( int i = 0 ; i < m_nIntersCount ; ++i) {
+      if ( ! m_Info[i].bOverlap  ||  ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq)) {
          if ( abs( m_Info[i].IciA[0].ptI.z - m_Info[i].IciB[0].ptI.z) < EPS_SMALL)
             ++nCount ;
       }
@@ -365,8 +365,8 @@ IntersCurveCurve::GetInt3DCrvCrvInfo( int nInd, IntCrvCrvInfo& aInfo)
    if ( nInd < 0  ||  nInd >= GetInters3DCount())
       return false ;
    int nCount = - 1 ;
-   for( int i = 0 ; i < m_nIntersCount ; ++i) {
-      if( ! m_Info[i].bOverlap  ||  ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq)) {
+   for ( int i = 0 ; i < m_nIntersCount ; ++i) {
+      if ( ! m_Info[i].bOverlap  ||  ( m_Info[i].bOverlap  &&  m_Info[i].bCBOverEq)) {
          if ( abs( m_Info[i].IciA[0].ptI.z - m_Info[i].IciB[0].ptI.z) < EPS_SMALL)
             ++nCount ;
       }
@@ -374,7 +374,7 @@ IntersCurveCurve::GetInt3DCrvCrvInfo( int nInd, IntCrvCrvInfo& aInfo)
          if ( abs( m_Info[i].IciA[0].ptI.z - m_Info[i].IciB[1].ptI.z) < EPS_SMALL)
             ++nCount ;
       }
-      if( nCount == nInd) {
+      if ( nCount == nInd) {
          aInfo = m_Info[nInd] ;
          return true ;
       }
@@ -389,11 +389,11 @@ IntersCurveCurve::GetIntersPointNearTo( int nCrv, const Point3d& ptNear, Point3d
    if ( m_nIntersCount == 0  ||  nCrv < 0  ||  nCrv > 1)
       return false ;
 
-  // ricerca del punto più vicino tra le intersezioni singole
+  // ricerca del punto più vicino tra le intersezioni singole
    bool bFound = false ;
    double dMinSqDist = SQ_INFINITO ;
    for ( int i = 0 ; i < m_nIntersCount ; ++ i) {
-     // se è un'intersezione singola
+     // se è un'intersezione singola
       if ( ! m_Info[i].bOverlap) {
         // faccio la verifica sul punto
          Point3d ptP = ( nCrv == 0 ? m_Info[i].IciA[0].ptI : m_Info[i].IciB[0].ptI) ;
@@ -458,7 +458,7 @@ IntersCurveCurve::GetCurveClassification( int nCrv, double dLenMin, CRVCVECTOR&
      // se esiste almeno una intersezione
       if ( m_nIntersCount >= 1)
          return CalcCurveClassification( m_pCurve[0], m_Info, dLenMin, ccClass) ;
-     // altrimenti la curva è completamente interna oppure completamente esterna
+     // altrimenti la curva è completamente interna oppure completamente esterna
       else
          return CalcCurveInOrOut( m_pCurve[0], m_pCurve[1], ccClass) ;
    }
@@ -475,7 +475,7 @@ IntersCurveCurve::GetCurveClassification( int nCrv, double dLenMin, CRVCVECTOR&
      // se esiste almeno una intersezione
       if ( m_nIntersCount >= 1)
          return CalcCurveClassification( m_pCurve[1], InfoTmp, dLenMin, ccClass) ;
-     // altrimenti la curva è completamente interna oppure completamente esterna
+     // altrimenti la curva è completamente interna oppure completamente esterna
       else
          return CalcCurveInOrOut( m_pCurve[1], m_pCurve[0], ccClass) ;
    }
@@ -540,7 +540,7 @@ IntersCurveCurve::CalcCurveClassification( const ICurve* pCurve, const ICCIVECTO
             double dU2 = Info[j].IciA[1].dU ;
             if ( dU2 < dU1  &&  pCurve->IsClosed())
                dU2 += dEndPar ;
-           // se cade nell'intervallo è da saltare
+           // se cade nell'intervallo è da saltare
             if ( Info[i].IciA[0].dU >= dU1  &&  Info[i].IciA[0].dU <= dU2) {
                bToSkip = true ;
                break ;
@@ -559,7 +559,7 @@ IntersCurveCurve::CalcCurveClassification( const ICurve* pCurve, const ICCIVECTO
    double dCurrPar = dStartPar ;
    double dCurrLen = 0 ;
    double dEndLen ; pCurve->GetLength( dEndLen) ;
-   // se è chiusa, recupero come finisce
+   // se è chiusa, recupero come finisce
    if ( pCurve->IsClosed()) {
       if ( ! InfoCorr[nNumInters-1].bOverlap)
          nLastTy = InfoCorr[nNumInters-1].IciA[0].nNextTy ;
@@ -577,7 +577,7 @@ IntersCurveCurve::CalcCurveClassification( const ICurve* pCurve, const ICCIVECTO
    }
   // costruisco il vettore delle classificazioni
    for ( int i = 0 ; i < nNumInters ; ++ i) {
-     // se è definito un tratto precedente
+     // se è definito un tratto precedente
       double dLenU ; pCurve->GetLengthAtParam( InfoCorr[i].IciA[0].dU, dLenU) ;
       if ( InfoCorr[i].IciA[0].dU > dCurrPar + EPS_PARAM  &&  dLenU - dCurrLen > dLenMin) {
         // verifico che la definizione sul tratto sia omogenea e valida
@@ -599,7 +599,7 @@ IntersCurveCurve::CalcCurveClassification( const ICurve* pCurve, const ICCIVECTO
      // altrimenti, salvo il tipo
       else
          nLastTy = InfoCorr[i].IciA[0].nNextTy ;
-     // se è definito un tratto in sovrapposizione
+     // se è definito un tratto in sovrapposizione
       if ( InfoCorr[i].bOverlap) {
         // assegno i dati
          CrvClass segClass ;
@@ -639,7 +639,7 @@ IntersCurveCurve::CalcCurveInOrOut( const ICurve* pCurveA, const ICurve* pCurveB
    double dStartParB, dEndParB ;
    if ( ! pCurveB->GetDomain( dStartParB, dEndParB))
       return false ;
-  // se almeno un punto di ciascuna curva è esterno al box dell'altra, sono sicuramente esterne
+  // se almeno un punto di ciascuna curva è esterno al box dell'altra, sono sicuramente esterne
    BBox3d boxCrvA, boxCrvB ;
    if ( ! pCurveA->GetLocalBBox( boxCrvA)  ||
         ! pCurveB->GetLocalBBox( boxCrvB))
@@ -682,13 +682,37 @@ IntersCurveCurve::CalcCurveInOrOut( const ICurve* pCurveA, const ICurve* pCurveB
    IntersCurveCurve iCC( clLine, *pCurveB) ;
    // dichiaro la classe della curva per default
    int nClass = CRVC_OUT ;
-   // se c'è almeno una intersezione
+   // se c'è almeno una intersezione
    if ( iCC.GetIntersCount() > 0) {
-     // se quanto precede la prima intersezione è interno, allora la curva è interna
+     // se quanto precede la prima intersezione è interno, allora la curva è interna
       IntCrvCrvInfo aInfo ;
       iCC.GetIntCrvCrvInfo( 0, aInfo) ;
       if ( aInfo.IciA[0].nPrevTy == ICCT_IN)
          nClass = CRVC_IN ;
+      else if ( aInfo.IciA[0].nPrevTy == ICCT_OUT)
+         nClass = CRVC_OUT ;
+      else if ( aInfo.IciA[0].nPrevTy == ICCT_NULL) {
+         // se il primo punto scelto non va bene allora ne cerco uno che mi dia informazioni sull'essere interno o esterno
+         CurveLine clLine ;
+         Point3d ptNewChoice ; pCurveA->GetPointD1D2( 0.25, ICurve::FROM_MINUS, ptNewChoice) ;
+         if ( ! clLine.SetPDL( ptNewChoice, 0, dLen))
+            return false ;
+         // calcolo l'intersezione
+         IntersCurveCurve iCC( clLine, *pCurveB) ;
+         if ( iCC.GetIntersCount() > 0) {
+           // se quanto precede la prima intersezione è interno, allora la curva è interna
+            IntCrvCrvInfo aInfo ;
+            iCC.GetIntCrvCrvInfo( 0, aInfo) ;
+            if ( aInfo.IciA[0].nPrevTy == ICCT_IN)
+               nClass = CRVC_IN ;
+            else if ( aInfo.IciA[0].nPrevTy == ICCT_OUT)
+               nClass = CRVC_OUT ;
+            else
+               return false ;  // se arrivo qui potrei ritentare la ricerca
+         }
+         else
+            return false ;
+      }
    }
    // altrimenti sono esterni tra loro
    else {

IntersCurvePlane.cpp

@@ -139,8 +139,6 @@ IntersCurvePlane::CalcIntersLinePlane( const Plane3d& plPlane, const ICurve& Cur
 void
 IntersCurvePlane::OrderAndCompleteIntersections()
 {
-   if ( m_Info.size() < 2)
-      return ;
    // cancello le interesezioni puntuali adiacenti a tratti di sovrapposizione
    // riempio le info PrevTy e NexyTy
    sort( m_Info.begin(), m_Info.end(), []( IntCrvPlnInfo& icpA, IntCrvPlnInfo& icpB) { return icpA.Ici[0].dU < icpA.Ici[0].dU ;}) ;

IntersLineLine.cpp

@@ -14,7 +14,6 @@
 //--------------------------- Include ----------------------------------------
 #include "stdafx.h"
 #include "IntersLineLine.h"
-#include "/EgtDev/Include/EgtNumUtils.h"
 #include <algorithm>
 
 using namespace std ;
@@ -157,34 +156,34 @@ IntersLineLine::IntersFiniteLines( const CurveLine& Line1, const CurveLine& Line
    double dCrossXY = CrossXY( vtDir1, vtDir2) ;
   // flag per linee parallele
    bool bParallel = ( abs( dCrossXY) < SIN_EPS_ANG_ZERO * ( dLen1XY * dLen2XY)) ;
-     // flag per segmenti che si allontanano significativamente
+  // flag per segmenti che si allontanano significativamente
    bool bFarEnds = ( nS1Side != 0  ||  nE1Side != 0  ||  nS2Side != 0  ||  nE2Side != 0) ;
-   
-   // analisi casi speciali di quasi parallelismo
+
+  // analisi casi speciali di quasi parallelismo
    // segmento sovrapposto all'altro
    double dDist1, dDist2 ;
-   if( nS1Side == 0 ||  nE1Side == 0  ||  nS1Side == nE1Side) {
+   if ( nS1Side == 0  ||  nE1Side == 0  ||  nS1Side == nE1Side) {
       dDist1 = CrossXY( ptS1 - ptS2, vtDir2) ;
       dDist2 = CrossXY( ptE1 - ptS2, vtDir2) ;
-      if( abs( dDist1 - dDist2) < EPS_SMALL * dLen2XY) {
+      if ( abs( dDist1 - dDist2) < EPS_SMALL * dLen2XY) {
          bParallel = true ;
          bFarEnds = ! ( (nS1Side == 0  &&  nE1Side == 0)  ||  (nS2Side == 0  &&  nE2Side == 0)) ;
       }
    }
-   else if( nS2Side == 0  ||  nE2Side == 0  ||  nS2Side == nE2Side) {
+   else if ( nS2Side == 0  ||  nE2Side == 0  ||  nS2Side == nE2Side) {
       dDist1 = CrossXY( ptS2 - ptS1, vtDir1) ;
       dDist2 = CrossXY( ptE2 - ptS1, vtDir1) ;
-      if( abs( dDist1 - dDist2) < EPS_SMALL * dLen1XY){
+      if ( abs( dDist1 - dDist2) < EPS_SMALL * dLen1XY) {
          bParallel = true ;
          bFarEnds = ! ( (nS1Side == 0  &&  nE1Side == 0)  ||  (nS2Side == 0  &&  nE2Side == 0)) ;
       }
    }
    // estremità sovrapposte di poco
    if ( ! bParallel  &&  abs( dCrossXY) < ( 0.1 * DEGTORAD) * ( dLen1XY * dLen2XY)) {
-      if (( nS1Side == 0  &&  nS2Side == 0  &&  ! AreSamePointXYEpsilon( ptS1, ptS2, 2 * EPS_SMALL))  ||
-          ( nS1Side == 0  &&  nE2Side == 0  &&  ! AreSamePointXYEpsilon( ptS1, ptE2, 2 * EPS_SMALL))  ||
-          ( nE1Side == 0  &&  nS2Side == 0  &&  ! AreSamePointXYEpsilon( ptE1, ptS2, 2 * EPS_SMALL))  ||
-          ( nE1Side == 0  &&  nE2Side == 0  &&  ! AreSamePointXYEpsilon( ptE1, ptE2, 2 * EPS_SMALL))) {
+      if (( nS1Side == 0  &&  nS2Side == 0  &&  ScalarXY( vtDir1, vtDir2) < 0  &&  ! AreSamePointXYEpsilon( ptS1, ptS2, 2 * EPS_SMALL))  ||
+          ( nS1Side == 0  &&  nE2Side == 0  &&  ScalarXY( vtDir1, vtDir2) > 0  &&  ! AreSamePointXYEpsilon( ptS1, ptE2, 2 * EPS_SMALL))  ||
+          ( nE1Side == 0  &&  nS2Side == 0  &&  ScalarXY( vtDir1, vtDir2) > 0  &&  ! AreSamePointXYEpsilon( ptE1, ptS2, 2 * EPS_SMALL))  ||
+          ( nE1Side == 0  &&  nE2Side == 0  &&  ScalarXY( vtDir1, vtDir2) < 0  &&  ! AreSamePointXYEpsilon( ptE1, ptE2, 2 * EPS_SMALL))) {
          bParallel = true ;
          bFarEnds = false ;
       }
@@ -195,54 +194,27 @@ IntersLineLine::IntersFiniteLines( const CurveLine& Line1, const CurveLine& Line
      // posizioni parametriche dell'intersezione sulle linee
       m_Info.IciA[0].dU = CrossXY( ( ptS2 - ptS1), vtDir2) / dCrossXY ;
       m_Info.IciB[0].dU = CrossXY( ( ptS2 - ptS1), vtDir1) / dCrossXY ;
-      // verifica posizione intersezione su prima linea
+     // verifica posizione intersezione su prima linea
       int nPos1 = ICurve::PP_NULL ;     // fuori
-      if ( nS1Side == 0 ||  nE1Side == 0) {
-         if( nS1Side == 0) {
-            nPos1 = ICurve::PP_START ;
-            m_Info.IciA[0].dU = 0 ;
-            m_Info.IciB[0].dU = vtDir2 * ( ptS1 - ptS2) / Pow( vtDir2.Len(), 2) ;
-         }
-         else {
-            nPos1 = ICurve::PP_END ;
-            m_Info.IciA[0].dU = 1 ;
-            m_Info.IciB[0].dU = vtDir2 * ( ptE1 - ptS2) / Pow( vtDir2.Len(), 2) ;
-         }
-      }
-      else {
-         if ( abs( m_Info.IciA[0].dU * dLen1XY) < EPS_SMALL)
-            nPos1 = ICurve::PP_START ;     // vicino a inizio
-         else if ( abs(( 1 - m_Info.IciA[0].dU) * dLen1XY) < EPS_SMALL)
-            nPos1 = ICurve::PP_END ;       // vicino a fine
-         else if ( m_Info.IciA[0].dU > 0  &&  m_Info.IciA[0].dU < 1)
-            nPos1 = ICurve::PP_MID ;       // nell'interno
-         else
-            return ;
-      }
+      if ( abs( m_Info.IciA[0].dU * dLen1XY) < EPS_SMALL)
+         nPos1 = ICurve::PP_START ;     // vicino a inizio
+      else if ( abs(( 1 - m_Info.IciA[0].dU) * dLen1XY) < EPS_SMALL)
+         nPos1 = ICurve::PP_END ;       // vicino a fine
+      else if ( m_Info.IciA[0].dU > 0  &&  m_Info.IciA[0].dU < 1)
+         nPos1 = ICurve::PP_MID ;       // nell'interno
+      else
+         return ;
      // verifica posizione intersezione su seconda linea
       int nPos2 = ICurve::PP_NULL ;     // fuori
-      if ( nS2Side == 0 ||  nE2Side == 0) {
-         if( nS2Side == 0) {
-            nPos2 = ICurve::PP_START ;
-            m_Info.IciB[0].dU = 0 ;
-            m_Info.IciA[0].dU = vtDir1 * (ptS2 - ptS1) / Pow( vtDir1.Len(), 2) ;
-         }
-         else {
-            nPos2 = ICurve::PP_END ;
-            m_Info.IciB[0].dU = 1 ;
-            m_Info.IciA[0].dU = vtDir1 * (ptE2 - ptS1) / Pow( vtDir1.Len(), 2) ;
-         }
-      }
-      else {
-         if ( abs( m_Info.IciB[0].dU * dLen2XY) < EPS_SMALL)
-            nPos2 = ICurve::PP_START ;     // vicino a inizio
-         else if ( abs(( 1 - m_Info.IciB[0].dU) * dLen2XY) < EPS_SMALL)
-            nPos2 = ICurve::PP_END ;       // vicino a fine
-         else if ( m_Info.IciB[0].dU > 0  &&  m_Info.IciB[0].dU < 1)
-            nPos2 = ICurve::PP_MID ;       // nell'interno
-         else
-            return ;
-      }     // limito i parametri a stare sui segmenti (0...1)
+      if ( abs( m_Info.IciB[0].dU * dLen2XY) < EPS_SMALL)
+         nPos2 = ICurve::PP_START ;     // vicino a inizio
+      else if ( abs(( 1 - m_Info.IciB[0].dU) * dLen2XY) < EPS_SMALL)
+         nPos2 = ICurve::PP_END ;       // vicino a fine
+      else if ( m_Info.IciB[0].dU > 0  &&  m_Info.IciB[0].dU < 1)
+         nPos2 = ICurve::PP_MID ;       // nell'interno
+      else
+         return ;
+     // limito i parametri a stare sui segmenti (0...1)
       m_Info.IciA[0].dU = min( max( m_Info.IciA[0].dU, 0.), 1.) ;
       m_Info.IciB[0].dU = min( max( m_Info.IciB[0].dU, 0.), 1.) ;
      // calcolo i punti sulle due linee (possono differire in Z)

PolygonElevation.cpp

@@ -172,9 +172,9 @@ PolygonElevationInClosedSurfTm( const Polygon3d& pgFacet, const ISurfTriMesh& Cl
                return true ;
             }
          }
-         else if ( ( Inters.nILTT == ILTT_VERT   ||  Inters.nILTT == ILTT_EDGE  ||  Inters.nILTT == ILTT_IN)  &&  Inters.dCosDN > EPS_ZERO)
+         else if ( ( Inters.nILTT == ILTT_VERT  ||  Inters.nILTT == ILTT_EDGE  ||  Inters.nILTT == ILTT_IN)  &&  Inters.dCosDN > EPS_ZERO)
             dElev = max( dElev, Inters.dU) ;
-         else if ( Inters.nILTT == ILTT_SEGM   ||  Inters.nILTT == ILTT_SEGM_ON_EDGE)
+         else if ( Inters.nILTT == ILTT_SEGM  ||  Inters.nILTT == ILTT_SEGM_ON_EDGE)
             dElev = max( dElev, Inters.dU2) ;
       }
    }
@@ -206,22 +206,30 @@ PolygonElevationInClosedSurfTm( const Polygon3d& pgFacet, const ISurfTriMesh& Cl
       if ( vEdges[i].first.z < EPS_SMALL  &&  vEdges[i].second.z < EPS_SMALL)
          continue ;
      // calcolo il segmento di linea
-     CurveLine clLine ;
+      CurveLine clLine ;
       if ( ! clLine.Set( vEdges[i].first, vEdges[i].second))
          return false ;
      // l'elevazione va aggiornata con la massima Z delle eventuali intersezioni dell'edge con il loop
       IntersCurveCurve intLL( clLine, ccLoop) ;
-      IntCrvCrvInfo aInfo ;
-      for ( int j = 0 ; intLL.GetIntCrvCrvInfo( j, aInfo) ; ++ j) {
-         dElev = max( dElev, aInfo.IciA[0].ptI.z) ;
-         if ( aInfo.bOverlap)
-            dElev = max( dElev, aInfo.IciA[1].ptI.z) ;
-        // se prima intersezione va da interno ad esterno allora devo considerare il punto iniziale del segmento (vertice)
-         if ( j == 0  &&  aInfo.IciA[0].nPrevTy == ICCT_IN)
-            dElev = max( dElev, vEdges[i].first.z) ;
-        // se ultima intersezione va da esterno a interno allora devo considerare il punto finale del segmento (vertice)
-         else if ( j == intLL.GetIntersCount() - 1  && aInfo.IciA[ aInfo.bOverlap ? 1 : 0].nNextTy == ICCT_IN)
-            dElev = max( dElev, vEdges[i].second.z) ;
+      if ( intLL.GetIntersCount() == 0) {
+         Point3d ptM = Media( vEdges[i].first, vEdges[i].second) ;
+         ptM.z = 0 ;
+         if ( IsPointInsidePolyLine( ptM, PL, -EPS_SMALL))
+            dElev = max( dElev, max( vEdges[i].first.z, vEdges[i].second.z)) ;
+      }
+      else {
+         IntCrvCrvInfo aInfo ;
+         for ( int j = 0 ; intLL.GetIntCrvCrvInfo( j, aInfo) ; ++ j) {
+            dElev = max( dElev, aInfo.IciA[0].ptI.z) ;
+            if ( aInfo.bOverlap)
+               dElev = max( dElev, aInfo.IciA[1].ptI.z) ;
+           // se prima intersezione va da interno ad esterno allora devo considerare il punto iniziale del segmento (vertice)
+            if ( j == 0  &&  aInfo.IciA[0].nPrevTy == ICCT_IN)
+               dElev = max( dElev, vEdges[i].first.z) ;
+           // se ultima intersezione va da esterno a interno allora devo considerare il punto finale del segmento (vertice)
+            else if ( j == intLL.GetIntersCount() - 1  && aInfo.IciA[ aInfo.bOverlap ? 1 : 0].nNextTy == ICCT_IN)
+               dElev = max( dElev, vEdges[i].second.z) ;
+         }
       }
    }
 

ProjectCurveSurf.cpp

@@ -19,6 +19,7 @@
 #include "/EgtDev/Include/EGkDistPointLine.h"
 #include "/EgtDev/Include/EGkDistPointCurve.h"
 #include "/EgtDev/Include/EGkDistPointSurfTm.h"
+#include "/EgtDev/Include/EGkDistPointSurfBz.h"
 #include "/EgtDev/Include/EGkIntersPlanePlane.h"
 #include "/EgtDev/Include/EGkIntersLinePlane.h"
 #include "/EgtDev/Include/EGkIntersLineSurfTm.h"
@@ -92,6 +93,8 @@ AddPointsOnCorners( PNT5AXVECTOR& vPt5ax)
                      Pt5ax.ptP = ptInt - vtLine1 / dLen1 * 2 * EPS_SMALL ;
                      Pt5ax.vtDir1 = vPt5ax[j].vtDir1 ;
                      Pt5ax.vtDir2 = vPt5ax[j].vtDir2 ;
+                     Pt5ax.vtDirU = vPt5ax[j].vtDirU ;
+                     Pt5ax.vtDirV = vPt5ax[j].vtDirV ;
                      Pt5ax.dPar = ( vPt5ax[i].dPar + vPt5ax[j].dPar) / 2 ;
                      Pt5ax.nFlag = P5AX_CVEX ;
                      vPt5ax.insert( vPt5ax.begin() + i, Pt5ax) ;
@@ -104,6 +107,8 @@ AddPointsOnCorners( PNT5AXVECTOR& vPt5ax)
                      Pt5ax.ptP = ptInt + vtLine2 / dLen2 * 2 * EPS_SMALL ;
                      Pt5ax.vtDir1 = vPt5ax[i].vtDir1 ;
                      Pt5ax.vtDir2 = vPt5ax[i].vtDir2 ;
+                     Pt5ax.vtDirU = vPt5ax[i].vtDirU ;
+                     Pt5ax.vtDirV = vPt5ax[i].vtDirV ;
                      Pt5ax.dPar = ( vPt5ax[i].dPar + vPt5ax[j].dPar) / 2 ;
                      Pt5ax.nFlag = P5AX_CVEX ;
                      vPt5ax.insert( vPt5ax.begin() + i, Pt5ax) ;
@@ -118,6 +123,8 @@ AddPointsOnCorners( PNT5AXVECTOR& vPt5ax)
                   Pt5ax.ptP = ptInt ;
                   Pt5ax.vtDir1 = Media( vPt5ax[i].vtDir1, vPt5ax[j].vtDir1) ; Pt5ax.vtDir1.Normalize() ;
                   Pt5ax.vtDir2 = Media( vPt5ax[i].vtDir2, vPt5ax[j].vtDir2) ; Pt5ax.vtDir2.Normalize() ;
+                  Pt5ax.vtDirU = Media( vPt5ax[i].vtDirU, vPt5ax[j].vtDirU) ; Pt5ax.vtDirU.Normalize() ;
+                  Pt5ax.vtDirV = Media( vPt5ax[i].vtDirV, vPt5ax[j].vtDirV) ; Pt5ax.vtDirV.Normalize() ;
                   Pt5ax.dPar = ( vPt5ax[i].dPar + vPt5ax[j].dPar) / 2 ;
                   Pt5ax.nFlag = P5AX_CONC ;
                   vPt5ax.insert( vPt5ax.begin() + i, Pt5ax) ;
@@ -233,23 +240,55 @@ RemovePointsInExcess( PNT5AXVECTOR& vPt5ax, double dLinTol, double dMaxSegmLen,
 }
 
 //----------------------------------------------------------------------------
-static bool
-ProjectPointOnSurf( const Point3d& ptP, const CISRFTMPVECTOR& vpStm, double dPar, Point5ax& Pt5ax)
+static const SurfTriMesh*
+MyGetAuxSurf( const ISurf* pSrf)
 {
-  // punto sulle supefici a minima distanza
+   if ( pSrf == nullptr)
+      return nullptr ;
+   switch ( pSrf->GetType()) {
+   case SRF_TRIMESH :
+      return GetBasicSurfTriMesh( pSrf) ;
+   case SRF_FLATRGN :
+      return GetBasicSurfFlatRegion( pSrf)->GetAuxSurf() ;
+   case SRF_BEZIER :
+      return GetBasicSurfBezier( pSrf)->GetAuxSurf() ;
+   default :
+      return nullptr ;
+   }
+}
+
+//----------------------------------------------------------------------------
+static bool
+ProjectPointOnSurf( const Point3d& ptP, const CISURFPVECTOR& vpSurf, double dPar, Point5ax& Pt5ax)
+{
+  // punto sulle superfici a minima distanza
    int nSurfMin = -1 ;
-   int nTriaMin ;
+   int nTriaMin = -1 ;
+   double dUMin = -1, dVMin = -1 ;
    Point3d ptMin ;
-   double dMinDist ;
-   for ( int i = 0 ; i < ssize( vpStm) ; ++ i) {
+   double dMinDist = NAN ;
+   for ( int i = 0 ; i < ssize( vpSurf) ; ++ i) {
      // punto sulla superficie a minima distanza
-      DistPointSurfTm dPS( ptP, *vpStm[i]) ;
-      double dDist ;
-      if ( dPS.GetDist( dDist)  &&  ( nSurfMin == -1  ||  dDist < dMinDist)) {
-         nSurfMin = i ;
-         dPS.GetMinDistPoint( ptMin) ;
-         dPS.GetMinDistTriaIndex ( nTriaMin) ;
-         dMinDist = dDist ;
+      int nSrfType = ( vpSurf[i] != nullptr ? vpSurf[i]->GetType() : GEO_NONE) ;
+      if ( nSrfType == SRF_TRIMESH  ||  nSrfType == SRF_FLATRGN) {
+         DistPointSurfTm dPS( ptP, *MyGetAuxSurf( vpSurf[i])) ;
+         double dDist ;
+         if ( dPS.GetDist( dDist)  &&  ( nSurfMin == -1  ||  dDist < dMinDist)) {
+            nSurfMin = i ;
+            dPS.GetMinDistPoint( ptMin) ;
+            dPS.GetMinDistTriaIndex ( nTriaMin) ;
+            dMinDist = dDist ;
+         }
+      }
+      else if ( nSrfType == SRF_BEZIER) {
+         DistPointSurfBz dPS( ptP, *GetBasicSurfBezier( vpSurf[i])) ;
+         double dDist ;
+         if ( dPS.GetDist( dDist)  &&  ( nSurfMin == -1  ||  dDist < dMinDist)) {
+            nSurfMin = i ;
+            dPS.GetMinDistPoint( ptMin) ;
+            dPS.GetParamsAtMinDistPoint( dUMin, dVMin) ;
+            dMinDist = dDist ;
+         }
       }
    }
 
@@ -257,19 +296,44 @@ ProjectPointOnSurf( const Point3d& ptP, const CISRFTMPVECTOR& vpStm, double dPar
    if ( nSurfMin >= 0) {
      // assegno il punto
       Point3d ptInt = ptMin ;
-     // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
-      Triangle3dEx trTria ;
-      if ( ! vpStm[nSurfMin]->GetTriangle( nTriaMin, trTria))
-         return false ;
-      Vector3d vtN ;
-      if ( ! CalcNormal( ptMin, trTria, vtN))
-         vtN = trTria.GetN() ;
-     // assegno valori al punto 5assi
-      Pt5ax.ptP = ptInt ;
-      Pt5ax.vtDir1 = vtN ;
-      Pt5ax.vtDir2 = vtN ;
-      Pt5ax.dPar = dPar ;
-      Pt5ax.nFlag = P5AX_STD ;
+     // calcolo gli altri dati
+      int nSrfType = ( vpSurf[nSurfMin] != nullptr ? vpSurf[nSurfMin]->GetType() : GEO_NONE) ;
+      if ( nSrfType == SRF_TRIMESH  ||  nSrfType == SRF_FLATRGN) {
+        // recupero superficie trimesh
+         const SurfTriMesh* pSurfTm = MyGetAuxSurf( vpSurf[nSurfMin]) ;
+        // calcolo la normale (si calcola smooth, in caso di errore si prende quella del triangolo)
+         Triangle3dEx trTria ;
+         if ( ! pSurfTm->GetTriangle( nTriaMin, trTria))
+            return false ;
+         Vector3d vtN ;
+         if ( ! CalcNormal( ptMin, trTria, vtN))
+            vtN = trTria.GetN() ;
+        // assegno valori al punto 5assi
+         Pt5ax.ptP = ptInt ;
+         Pt5ax.vtDir1 = vtN ;
+         Pt5ax.vtDir2 = vtN ;
+         Pt5ax.vtDirU = V_NULL ;
+         Pt5ax.vtDirV = V_NULL ;
+         Pt5ax.dPar = dPar ;
+         Pt5ax.nFlag = P5AX_STD ;
+      }
+      else if ( nSrfType == SRF_BEZIER) {
+         Point3d ptSB ;
+         Vector3d vtN, vtDerU, vtDerV ;
+         if ( ! GetBasicSurfBezier( vpSurf[nSurfMin])->GetPointNrmD1D2( dUMin, dVMin, ISurfBezier::FROM_MINUS, ISurfBezier::FROM_MINUS,
+                                                                        ptSB, vtN, &vtDerU, &vtDerV))
+            return false ;
+         vtDerU.Normalize() ;
+         vtDerV.Normalize() ;
+        // assegno valori al punto 5assi
+         Pt5ax.ptP = ptInt ;
+         Pt5ax.vtDir1 = vtN ;
+         Pt5ax.vtDir2 = vtN ;
+         Pt5ax.vtDirU = vtDerU ;
+         Pt5ax.vtDirV = vtDerV ;
+         Pt5ax.dPar = dPar ;
+         Pt5ax.nFlag = P5AX_STD ;
+      }
      // ritorno con successo
       return true ;
    }
@@ -282,31 +346,6 @@ bool
 ProjectCurveOnSurf( const ICurve& crCrv, const CISURFPVECTOR& vpSurf,
                     double dLinTol, double dMaxSegmLen, bool bSharpEdges, PNT5AXVECTOR& vPt5ax)
 {
-  // sistemazioni per tipo di superficie
-   CISRFTMPVECTOR vpSurfTm ;
-   for ( int i = 0 ; i < ssize( vpSurf) ; ++ i) {
-      const SurfTriMesh* pSurfTm = nullptr ;
-      switch ( vpSurf[i]->GetType()) {
-      case SRF_TRIMESH :
-         pSurfTm = GetBasicSurfTriMesh( vpSurf[i]) ;
-         break ;
-      case SRF_BEZIER :
-       { double dOldLinTol = GetSurfBezierAuxSurfRefinedTol() ;
-         SetSurfBezierAuxSurfRefinedTol( GetSurfBezierTol( dLinTol)) ;
-         pSurfTm = GetBasicSurfBezier( vpSurf[i])->GetAuxSurfRefined() ;
-         SetSurfBezierAuxSurfRefinedTol( dOldLinTol) ;
-       } break ;
-      case SRF_FLATRGN :
-         pSurfTm = GetBasicSurfFlatRegion( vpSurf[i])->GetAuxSurf() ;
-         break ;
-      default :
-         break ;
-      }
-      if ( pSurfTm == nullptr)
-         return false ;
-      vpSurfTm.emplace_back( pSurfTm) ;
-   }
-
   // controllo le tolleranze
    dLinTol = max( dLinTol, LIN_TOL_MIN) ;
    dMaxSegmLen = max( dMaxSegmLen, 10 * EPS_SMALL) ;
@@ -331,7 +370,7 @@ ProjectCurveOnSurf( const ICurve& crCrv, const CISURFPVECTOR& vpSurf,
    while ( bFound) {
      // se trovo proiezione, la salvo
       Point5ax Pt5ax ; 
-      if ( ProjectPointOnSurf( ptP, vpSurfTm, dPar, Pt5ax))
+      if ( ProjectPointOnSurf( ptP, vpSurf, dPar, Pt5ax))
          vPt5ax.emplace_back( Pt5ax) ;
      // passo al successivo
       bFound = PL.GetNextUPoint( &dPar, &ptP) ;
@@ -403,6 +442,8 @@ ProjectPointOnSurf( const Point3d& ptP, const CISRFTMPVECTOR& vpStm, const Frame
       Pt5ax.ptP = ptInt ;
       Pt5ax.vtDir1 = vtN ;
       Pt5ax.vtDir2 = frRefLine.VersZ() ;
+      Pt5ax.vtDirU = V_NULL ;
+      Pt5ax.vtDirV = V_NULL ;
       Pt5ax.dPar = dPar ;
       Pt5ax.nFlag = P5AX_STD ;
      // ritorno con successo
@@ -557,6 +598,8 @@ ProjectPointOnSurf( const Point3d& ptP, const CISRFTMPVECTOR& vpStm, const IGeoP
          Pt5ax.ptP = ptInt ;
          Pt5ax.vtDir1 = vtN ;
          Pt5ax.vtDir2 = vtLine ;
+         Pt5ax.vtDirU = V_NULL ;
+         Pt5ax.vtDirV = V_NULL ;
          Pt5ax.dPar = dPar ;
          Pt5ax.nFlag = P5AX_STD ;
         // ritorno con successo
@@ -696,6 +739,8 @@ ProjectPointOnSurf( const Point3d& ptP, const CISRFTMPVECTOR& vpStm, const ICurv
             Pt5ax.ptP = ptInt ;
             Pt5ax.vtDir1 = vtN ;
             Pt5ax.vtDir2 = vtLine ;
+            Pt5ax.vtDirU = V_NULL ;
+            Pt5ax.vtDirV = V_NULL ;
             Pt5ax.dPar = dPar ;
             Pt5ax.nFlag = P5AX_STD ;
            // ritorno con successo
@@ -852,6 +897,8 @@ ProjectPointOnSurf( const Point3d& ptP, const CISRFTMPVECTOR& vpStm, const SurfT
          Pt5ax.ptP = ptInt ;
          Pt5ax.vtDir1 = vtN ;
          Pt5ax.vtDir2 = vtN2 ;
+         Pt5ax.vtDirU = V_NULL ;
+         Pt5ax.vtDirV = V_NULL ;
          Pt5ax.dPar = dPar ;
          Pt5ax.nFlag = P5AX_STD ;
         // ritorno con successo

SbzFromCurves.cpp

@@ -735,8 +735,7 @@ GetSurfBezierRuled( const ICurve* pCurve1, const ICurve* pCurve2, int nType, dou
 
 //-------------------------------------------------------------------------------
 ISurfBezier*
-GetSurfBezierRuledGuided( const ICurve* pCurve1, const ICurve* pCurve2, const ICURVEPOVECTOR& vCrv, const ICURVEPOVECTOR& vNewCrv,
-                          const INTVECTOR& vShown, const INTINTVECTOR& vNewOrEdited, double dLinTol)
+GetSurfBezierRuledGuided( const ICurve* pCurve1, const ICurve* pCurve2, const BIPNTVECTOR& vCrv, double dLinTol)
 {
    // verifica parametri
    if ( pCurve1 == nullptr  ||  pCurve2 == nullptr)
@@ -763,7 +762,7 @@ GetSurfBezierRuledGuided( const ICurve* pCurve1, const ICurve* pCurve2, const IC
 
    // creo e setto la superficie trimesh
    PtrOwner<SurfBezier> pSbz( CreateBasicSurfBezier()) ;
-   if ( IsNull( pSbz)  ||  ! pSbz->CreateByIsoParamSet( pCC1, pCC2, vCrv, vNewCrv, vShown, vNewOrEdited))
+   if ( IsNull( pSbz)  ||  ! pSbz->CreateByIsoParamSet( pCC1, pCC2, vCrv))
       return nullptr ;
 
    // restituisco la superficie

StringUtils3d.cpp

@@ -76,11 +76,18 @@ FromString( const string& sVal, Frame3d& frFrame)
 bool
 FromString( const string& sVal, Color& cCol)
 {
-  // devono essere 4 parametri : Red, Green, Blue, Alpha
-   int vnVal[4] ;
-   if ( ! FromString( sVal, vnVal))
-      return false ;
-  // assegno il colore
-   cCol.Set( vnVal[0], vnVal[1], vnVal[2], vnVal[3]) ;
-   return true ;
+  // dovrebbero essere 4 parametri : Red, Green, Blue, Alpha
+   int vnRGBA[4] ;
+   if ( FromString( sVal, vnRGBA)) {
+      cCol.Set( vnRGBA[0], vnRGBA[1], vnRGBA[2], vnRGBA[3]) ;
+      return true ;
+   }
+  // riprovo con 3 parametri : Red, Green, Blue
+   int vnRGB[3] ;
+   if ( FromString( sVal, vnRGB)) {
+      cCol.Set( vnRGB[0], vnRGB[1], vnRGB[2]) ;
+      return true ;
+   }
+  // altrimenti errore
+   return false ;
 }

SurfAux.cpp

@@ -24,6 +24,13 @@
 #include "/EgtDev/Include/EGkSurfAux.h"
 #include "/EgtDev/Include/EGkSfrCreate.h"
 #include "/EgtDev/Include/EgtPointerOwner.h"
+#include "/EgtDev/Include/EGkChainCurves.h"
+
+
+#define SAVEMKUNIF_CRVS 0
+#if SAVEMKUNIF_CRVS
+   #include "/EgtDev/Include/EGkGeoObjSave.h"
+#endif
 
 using namespace std ;
 
@@ -530,6 +537,25 @@ NurbsToBezierSurface(const SNurbsSurfData& snData)
    return Release( pSrfBz) ;
 }
 
+//----------------------------------------------------------------------------
+ICurveComposite*
+GetRectangleCurve( const Point3d& ptStart, double dWidth, double dHeight)
+{
+   PolyLine PL ;
+   PL.AddUPoint( 0, ptStart) ;
+   PL.AddUPoint( 1, Point3d( ptStart.x + dWidth, ptStart.y)) ;
+   PL.AddUPoint( 2, Point3d( ptStart.x + dWidth, ptStart.y + dHeight, 0)) ;
+   PL.AddUPoint( 3, Point3d( ptStart.x , ptStart.y + dHeight, 0)) ;
+   PL.AddUPoint( 4, ptStart) ;
+   // creo la curva e la inserisco nel GDB
+   PtrOwner<ICurveComposite> pCrvCompo( CreateCurveComposite()) ;
+   bool bOk = true ;
+   bOk = bOk && ! IsNull( pCrvCompo) ;
+   // inserisco i segmenti che uniscono i punti
+   bOk = bOk && pCrvCompo->FromPolyLine( PL) ;
+   return Release( pCrvCompo) ;
+}
+
 //----------------------------------------------------------------------------
 bool
 MakeUniform( ISurfFlatRegion*& pSfr, bool& bRescaled, const DBLVECTOR& vU0, const DBLVECTOR& vV0,
@@ -539,36 +565,33 @@ MakeUniform( ISurfFlatRegion*& pSfr, bool& bRescaled, const DBLVECTOR& vU0, cons
    bool bRescaledU = false ;
    bool bRescaledV = false ;
    int nSpanU = 1, nSpanV = 1 ;
-   PtrOwner<ISurfFlatRegion> pRescaledSfr( CreateSurfFlatRegion()) ;
+   ICRVCOMPOPOVECTOR vUniformedCurves ;
+   BOOLVECTOR vbUniform(2) ;
+   fill( vbUniform.begin(), vbUniform.end(), true) ;
+   DBLMATRIX mKnots(2) ;
    for ( int nDir = 0 ; nDir <= 1 ; ++ nDir) {
      // vettore dei nodi
-      DBLVECTOR vU ;
+      DBLVECTOR& vU = mKnots[nDir] ;
       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 ) {
+         for ( int i = nExtraKnots ; i < ssize( vU0) - 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)
+            if ( i == nExtraKnots  ||  dKnot > vU.back() + EPS_SMALL)
                vU.push_back( dKnot) ;
          }
-         nSpanU = (int)vU.size() - 1 ;
+         nSpanU = ssize( vU) - 1 ;
       }
      // controllo in V
       else if ( nDir == 1 ) {
-         if ( nDegV > 1) {
-            nExtraKnots = nDegV - 1 ;
-         }
-         for ( int i = nExtraKnots ; i < int( vV0.size()) - nExtraKnots ; ++i ) {
+         for ( int i = nExtraKnots ; i < ssize( vV0) - 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)
+            if ( i == nExtraKnots  ||  dKnot > vU.back() + EPS_SMALL)
                vU.push_back( dKnot) ;
          }
-         nSpanV = (int)vU.size() - 1 ;
+         nSpanV = ssize( vU) - 1 ;
       }
 
      // controllo se il vettore dei nodi è uniforme
@@ -581,45 +604,97 @@ MakeUniform( ISurfFlatRegion*& pSfr, bool& bRescaled, const DBLVECTOR& vU0, cons
          if ( b < (int)vU.size())
             d1 = abs( vU[b] - vU[a]) ;
       }
-      if ( b != (int)vU.size()) {
+      if ( b != (int)vU.size())
+         vbUniform[nDir] = false ;
+   }
+   // vettore delle curve di loop della regione di trim
+   ICRVCOMPOPOVECTOR vLoop ;
+   if ( ! vbUniform[0]  ||  ! vbUniform[1]) {
+      for ( int c = 0 ;  c < pSfr->GetChunkCount() ; ++c) {
+         for ( int l = 0 ; l < pSfr->GetLoopCount( c); ++l)
+            vLoop.emplace_back( ConvertCurveToComposite( pSfr->GetLoop( c, l))) ;
+      }
+   }
+
+#if SAVEMKUNIF_CRVS
+   //debug
+   vector<IGeoObj*> vGeo ;
+   for( int i = 0 ; i < ssize( vLoop); ++i){
+      vGeo.push_back(vLoop[i]->Clone()) ;
+   }
+   SaveGeoObj( vGeo, "D:\\Temp\\bezier\\import3dm\\trim_error\\failed_trim_crv_loops.nge") ;
+   //debug
+#endif
+
+   for ( int nDir = 0 ; nDir <= 1 ; ++ nDir) {
+      DBLVECTOR& vU = mKnots[nDir] ;
+      if ( ! vbUniform[nDir]) {
          nDir == 0 ? bRescaledU = true : bRescaledV = true ;
-         pRescaledSfr.Set( CreateSurfFlatRegion()) ;
-         if ( IsNull( pRescaledSfr))
-            return false ;
+         // creo il vettore delle curve all'interno di una striscia
+         ICRVCOMPOPOVECTOR vCrvStrip ; 
          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()) ;
+            PtrOwner<ICurveComposite> pTrimMask ;
 
             // 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) ;
+               Point3d ptStart( abs(vU[p] - vU.front()), - 1, 0) ;
+               pTrimMask.Set( GetRectangleCurve( ptStart, dLenStrip, dScaleV + 2)) ;
             }
             else{
-               pSfrTrim.Set( GetSurfFlatRegionRectangle( dScaleU + 2, dLenStrip)) ;
-               vtTrim.Set( - 1, abs(vU[p] - vU.front()), 0) ;
+               Point3d ptStart( - 1, abs(vU[p] - vU.front()), 0) ;
+               pTrimMask.Set( GetRectangleCurve( ptStart, dScaleU + 2, dLenStrip)) ;
             }
-            pSfrTrim->Translate( vtTrim) ;
 
-            if ( ! pSfr_copy->Intersect( *pSfrTrim))
-               return false ;
+            // qui potrei decidere di eliminare tutti i tratti dei loop paralleli alla direzione del parametro che sto analizzando ( e di valore coincidente a quello di un nodo)
 
-            // aggiungo la nuova striscia solo se è valida
-            if ( pSfr_copy->IsValid() ) {
+            for ( int l = 0 ; l < ssize( vLoop); ++l) {
+               #if SAVEMKUNIF_CRVS
+                  //debug
+                  vector<IGeoObj*> vGeo ;
+                  vGeo.push_back(pTrimMask->Clone()) ;
+                  vGeo.push_back(vLoop[l]->Clone()) ;
+                  SaveGeoObj( vGeo, "D:\\Temp\\bezier\\import3dm\\trim_error\\failed_trim_crv_inters.nge") ;
+                  //debug
+               #endif
+               IntersCurveCurve icc( *pTrimMask, *vLoop[l]) ;
+               CRVCVECTOR vCurveClass ;
+               icc.GetCurveClassification( 1, 10 * EPS_SMALL, vCurveClass) ;
+               for( int i = 0 ; i < ssize( vCurveClass); ++i) {
+                  if( vCurveClass[i].nClass == CRVC_IN  ||  vCurveClass[i].nClass == CRVC_ON_P)
+                     vCrvStrip.emplace_back( ConvertCurveToComposite( vLoop[l]->CopyParamRange( vCurveClass[i].dParS, vCurveClass[i].dParE))) ;
+               }
+            }
+
+            #if SAVEMKUNIF_CRVS
+               //debug
+               vector<IGeoObj*> vGeo ;
+               for( int i = 0 ; i < ssize( vCrvStrip); ++i){
+                  vGeo.push_back(vCrvStrip[i]->Clone()) ;
+               }
+               SaveGeoObj( vGeo, "D:\\Temp\\bezier\\import3dm\\trim_error\\failed_trim_crv_strip.nge") ;
+               //debug
+            #endif
+
+            // riscalo le curve nella striscia
+            if ( ! vCrvStrip.empty()) {
+               double dCoeffX = 1, dCoeffY = 1 ;
+               
                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
+                  dCoeffX = SBZ_TREG_COEFF / dLenStrip ;
+               else 
+                  dCoeffY = SBZ_TREG_COEFF / dLenStrip ;
+               for( int i = 0 ; i < ssize( vCrvStrip); ++i) {
+                  if( ! IsNull( vCrvStrip[i]))
+                     vCrvStrip[i]->Scale( GLOB_FRM, dCoeffX, dCoeffY, 1) ;
+                  else
+                     return false ;
+               }
 
+               // prima di riunire le curve al resto devo traslarle 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)
@@ -632,25 +707,78 @@ MakeUniform( ISurfFlatRegion*& pSfr, bool& bRescaled, const DBLVECTOR& vU0, cons
                   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_CHAMFER) ; // OFFSET
-               if ( pRescaledSfr->IsValid()) {
-                  if ( ! pRescaledSfr->Add( *pSfr_copy))
-                     return false ;
+               for( int i = 0 ; i < ssize( vCrvStrip); ++i)
+                  vCrvStrip[i]->Translate( vtJoin) ;
+
+            #if SAVEMKUNIF_CRVS
+               //debug
+               vector<IGeoObj*> vGeo ;
+               for( int i = 0 ; i < ssize( vCrvStrip); ++i){
+                  vGeo.push_back(vCrvStrip[i]->Clone()) ;
+               }
+               SaveGeoObj( vGeo, "D:\\Temp\\bezier\\import3dm\\trim_error\\failed_trim_crv_strip.nge") ;
+               //debug
+            #endif
+               
+               // faccio la chain con le curve delle striscie precedenti
+               
+               if ( ! vUniformedCurves.empty()  ||  ! vCrvStrip.empty()) {
+                  #if SAVEMKUNIF_CRVS
+                     //debug
+                     vector<IGeoObj*> vGeo ;
+                     for( int i = 0 ; i < ssize( vUniformedCurves); ++i){
+                        vGeo.push_back(vUniformedCurves[i]->Clone()) ;
+                     }
+                     SaveGeoObj( vGeo, "D:\\Temp\\bezier\\import3dm\\trim_error\\failed_trim_crv_unif.nge") ;
+                     //debug
+                  #endif
+
+                  ChainCurves chainCrv ;
+                  double dChainTol = 5 * EPS_SMALL ;
+                  chainCrv.Init( false, dChainTol, max(ssize( vUniformedCurves), ssize(vCrvStrip))) ;
+                  for ( int c = 0 ; c < ssize( vUniformedCurves) ; ++c) {
+                     Point3d ptStart, ptEnd ;
+                     Vector3d vtStart, vtEnd ;
+                     vUniformedCurves[c]->GetStartPoint( ptStart) ;
+                     vUniformedCurves[c]->GetEndPoint( ptEnd) ;
+                     vUniformedCurves[c]->GetStartDir( vtStart) ;
+                     vUniformedCurves[c]->GetEndDir( vtEnd) ;
+                     chainCrv.AddCurve( 1 + c, ptStart, vtStart, ptEnd, vtEnd) ;
+                  }
+                  for ( int c = 0 ; c < ssize( vCrvStrip); ++c) {
+                     Point3d ptStart, ptEnd ;
+                     Vector3d vtStart, vtEnd ;
+                     vCrvStrip[c]->GetStartPoint( ptStart) ;
+                     vCrvStrip[c]->GetEndPoint( ptEnd) ;
+                     vCrvStrip[c]->GetStartDir( vtStart) ;
+                     vCrvStrip[c]->GetEndDir( vtEnd) ;
+                     chainCrv.AddCurve( 1 + ssize( vUniformedCurves) + c, ptStart, vtStart, ptEnd, vtEnd) ;
+                  }
+                  INTVECTOR vIds ;
+                  ICRVCOMPOPOVECTOR vNewCrv ;
+                  int nCrvPrec = ssize( vUniformedCurves) ;
+                  while ( chainCrv.GetChainFromNear( ORIG, false, vIds)) {
+                     vNewCrv.emplace_back( CreateBasicCurveComposite()) ;
+                     for ( int nId : vIds) {
+                        nId -= 1 ;
+                        if ( nId < nCrvPrec)
+                           vNewCrv.back()->AddCurve( Release( vUniformedCurves[nId]), true, dChainTol) ;
+                        else
+                           vNewCrv.back()->AddCurve( Release( vCrvStrip[nId - ssize( vUniformedCurves)]), true, dChainTol) ;
+                     }
+                  }
+                  // aggiorno le curve 
+                  vUniformedCurves.clear() ;
+                  vUniformedCurves.swap( vNewCrv)  ;
                }
-               else
-                  pRescaledSfr.Set( pSfr_copy) ;
             }
+            vCrvStrip.clear() ;
          }
          if ( nDir == 0) {
             dScaleU = ((int)vU.size() - 1) * SBZ_TREG_COEFF ;
-            if ( pRescaledSfr->IsValid()) {
-               if ( bRetry)
-                  pRescaledSfr->Offset( -10 * EPS_SMALL, ICurve::OFF_CHAMFER) ;  //contro OFFSET
-               delete pSfr ;
-               pSfr = Release( pRescaledSfr) ;
+            if( ! vbUniform[1]) {
+               vLoop.swap( vUniformedCurves) ;
+               vUniformedCurves.clear() ;
             }
          }
          else
@@ -658,12 +786,38 @@ MakeUniform( ISurfFlatRegion*& pSfr, bool& bRescaled, const DBLVECTOR& vU0, cons
       }
    }
 
-   if ( ! IsNull( pRescaledSfr)  &&  pRescaledSfr->IsValid()) {
-      if ( bRetry)
-         pRescaledSfr->Offset( -10 * EPS_SMALL, ICurve::OFF_CHAMFER) ;  // contro OFFSET
-      delete pSfr ;
-      pSfr = Release( pRescaledSfr) ;
+   if ( ! vUniformedCurves.empty()) {
+      #if SAVEMKUNIF_CRVS
+         //debug
+         vector<IGeoObj*> vGeo ;
+         for( int i = 0 ; i < ssize( vUniformedCurves); ++i){
+            vGeo.push_back(vUniformedCurves[i]->Clone()) ;
+         }
+         SaveGeoObj( vGeo, "D:\\Temp\\bezier\\import3dm\\trim_error\\failed_trim_crv_unif.nge") ;
+         //debug
+      #endif
+
+      // controllo che tutte le curve siano chiuse, sennò vuol dire che ho perso qualche pezzo durante le intersezioni
+      for ( int i = 0 ; i < ssize( vUniformedCurves); ++i) {
+         if ( ! vUniformedCurves[i]->IsClosed())
+            return false ;
+      }
+
+      // creo una regione dalle curve riscalate
+      SurfFlatRegionByContours sfrRescaled ;
+      for ( int i = 0 ; i < ssize( vUniformedCurves); ++i)
+         sfrRescaled.AddCurve( Release( vUniformedCurves[i])) ;
+
+      PtrOwner<ISurfFlatRegion> pRescaledSfr( sfrRescaled.GetSurf()) ;
+      if ( ! IsNull( pRescaledSfr)  &&  pRescaledSfr->IsValid()) {
+         delete pSfr ;
+         pSfr = Release( pRescaledSfr) ;
+      }
+      else
+         return false ;
    }
+   else if( ! vbUniform[0]  ||  ! vbUniform[1])
+      return false ;
 
    if ( ! bRescaledU  &&  ! bRescaledV)
       pSfr->Scale( GLOB_FRM, nSpanU / dScaleU * SBZ_TREG_COEFF, nSpanV / dScaleV * SBZ_TREG_COEFF, 1) ;

SurfBezier.h

@@ -116,9 +116,9 @@ class SurfBezier : public ISurfBezier, public IGeoObjRW
       bool GetControlCurveOnV( int nIndU, PolyLine& plCtrlV) const override ;
       const SurfTriMesh* GetAuxSurf( void) const override ;
       const SurfTriMesh* GetAuxSurfRefined( void) const override ;
-      SurfTriMesh* GetApproxSurf( double dTol, double dSideMin = 100 * EPS_SMALL, bool bUpdateEdges = false) const override ;
+      SurfTriMesh* GetApproxSurf( double dTol, double dSideMin = 10 * EPS_SMALL, bool bUpdateEdges = false) 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 GetLeaves( std::vector<std::tuple<int, Point3d, Point3d>>& vLeaves, bool bRefined = false) const override ;
       bool GetTriangles2D( std::vector<std::tuple<int,Point3d, Point3d, Point3d>>& vTria2D) const override ;
       // funzioni che servono per ricavare l'immagine nel parametrico di un punto appartenente alla trimesh ausiliaria della superficie di Bezier
       bool UnprojectPointFromStm( int nT, const Point3d& ptI, Point3d& ptSP, int nIL = IntLineTriaType::ILTT_IN) const override ;
@@ -143,16 +143,16 @@ class SurfBezier : public ISurfBezier, public IGeoObjRW
       bool IsPlanar( void) const override ;
       bool CreateByFlatContour( const PolyLine& PL) override ;
       bool CreateByRegion( const POLYLINEVECTOR& vPL) override ;
-      bool CreateByExtrusion( const ICurve* pCurve, const Vector3d& vtExtr, bool bDeg3OrDeg2 = false) override ;
+      bool CreateByExtrusion( const ICurve* pCurve, const Vector3d& vtExtr) override ;
       bool CreateByScrewing( const ICurve* pCurve, const Point3d& ptAx, const Vector3d& vtAx, double dAngRotDeg, double dMove) override ;
       bool CreateByPointCurve( const Point3d& pt, const ICurve* pCurve) override ;
       bool CreateByTwoCurves( const ICurve* pCurve1, const ICurve* pCurve2, int nType) override ;
       bool CreateBySetOfCurves( const ICURVEPOVECTOR& vCrvBez, bool bReduceToDeg3) override ;
       PNTVECTOR GetAllControlPoints( void) const ;
       bool GetAllPatchesIsocurves( bool bUorV, ICURVEPOVECTOR& vCrv) const ;
-      bool CreateByIsoParamSet( const ICurve* pCurve0, const ICurve* pCurve1, const ICURVEPOVECTOR& vCrv, const ICURVEPOVECTOR& vNewCrv, const INTVECTOR& vShown ,const INTINTVECTOR& vNewOrEdited) ;
-      bool RemoveCollapsedSpans() override ;
-      bool SwapParameters() ;
+      bool CreateByIsoParamSet( const ICurve* pCurve0, const ICurve* pCurve1, const BIPNTVECTOR& vCrv) ;
+      bool RemoveCollapsedSpans( void) override ;
+      bool SwapParameters( void) ;
 
    public :                      // IGeoObjRW
       int  GetNgeId( void) const override ;

SurfTriMesh.cpp

@@ -4415,6 +4415,35 @@ SurfTriMesh::ClonePart( int nPart) const
    return Release( pSurfTM) ;
 }
 
+//----------------------------------------------------------------------------
+bool 
+SurfTriMesh::GetPartAndShellFromFacet( int nFacet, int& nPart, int& nShell) const
+{
+   // l'indice della faccia deve essere nei limiti
+   if ( nFacet < 0  ||  nFacet >= int( m_vFacet.size()))
+      return false ;
+
+   // mi assicuro che siano calcolate il numero di parti e di shell
+   int nParts = GetPartCount() ;
+
+   int nTria = m_vFacet[nFacet] ;
+   nShell = m_vTria[nTria].nShell ;
+
+   // scopro in quale part è la shell
+   int nPartTemp = 0 ;
+   nPart = - 1 ;
+   while ( nPartTemp < nParts  &&  nPart == - 1) {
+      for ( int i : m_vPart[nPartTemp].vShell) {
+         if ( i == nShell) {
+            nPart = nPartTemp ;
+            break ;
+         }
+      }
+      ++nPartTemp ;
+   }
+   return true ;
+}
+
 //----------------------------------------------------------------------------
 bool
 SurfTriMesh::ResetTFlags( void)

SurfTriMesh.h

@@ -347,6 +347,7 @@ class SurfTriMesh : public ISurfTriMesh, public IGeoObjRW
       bool GetPartVolume( int nPart, double& dVolume) const override ;
       bool GetPartLoops( int nPart, POLYLINEVECTOR& vPL) const override ;
       SurfTriMesh* ClonePart( int nPart) const override ;
+      bool GetPartAndShellFromFacet( int nFacet, int& nPart, int& nShell) const override ;
       bool SetTFlag( int nId, int nTFlag) override ;
       bool GetTFlag( int nId, int& nFlag) const override ;
       int  GetMaxTFlag( void) const override 

Tool.h

@@ -36,25 +36,25 @@ class Tool
       bool SetAdditiveTool( const std::string& sToolName, double dH, double dR, double dRc, int nToolNum) ;
       bool SetToolNum( int nToolNum)
                { m_nCurrentNum = nToolNum ; return true ; }
-      int GetType() const
+      int GetType( void) const
                { return m_nType ; }
-      int GetToolNum() const
+      int GetToolNum( void) const
                { return m_nCurrentNum ; }
-      double GetHeigth() const
+      double GetHeigth( void) const
                { return m_dHeight ; }
-      double GetTipHeigth() const
+      double GetTipHeigth( void) const
                { return m_dTipHeight ; }
-      double GetRadius() const
+      double GetRadius( void) const
                { return m_dRadius ; }
-      double GetTipRadius() const
+      double GetTipRadius( void) const
                { return m_dTipRadius ; } 
-      double GetCornRadius() const
+      double GetCornRadius( void) const
                { return m_dRCorner ; }
-      double GetRefRadius() const
+      double GetRefRadius( void) const
                { return m_dRefRadius ; } 
-      double GetMrtChsWidth() const
+      double GetMrtChsWidth( void) const
                { return m_dMrtChsWidth ; }
-      double GetMrtChsThickness() const
+      double GetMrtChsThickness( void) const
                { return m_dMrtChsThickness ; }
       const CurveComposite& GetOutline( void) const
                { return ( m_Outline) ; }

Tree.cpp

@@ -177,7 +177,7 @@ Tree::SetSurf( const SurfBezier* pSrfBz, const Point3d& ptMin, const Point3d& pt
       int nLoop = 0 ;
      // 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
+      double dLinTol = 0.005 ; // questo è riferito allo spazio parametrico
       double dAngTolDeg = 5 ;
       for ( int i = 0 ; i < pTrimReg->GetChunkCount() ; ++ i) {
          PtrOwner<SurfFlatRegion> pChunk( pTrimReg->CloneChunk( i)) ;
@@ -185,6 +185,8 @@ Tree::SetSurf( const SurfBezier* pSrfBz, const Point3d& ptMin, const Point3d& pt
            // i chunk della falt region sono ancora flat region composte da 1 chunk
             // rimuovo i difetti dei loop prima di salvarli
             PtrOwner<ICurveComposite> pLoop( GetBasicCurveComposite( pChunk->GetLoop( 0, j))) ;
+            if ( IsNull( pLoop))
+               return false ;
             pLoop->MergeCurves( dLinTol, dAngTolDeg) ;
             pLoop->RemoveSmallDefects( dLinTol, dAngTolDeg, true) ;
             pLoop->RemoveSmallParts( dLinTol, dAngTolDeg) ;
@@ -461,8 +463,8 @@ Tree::Split( int nId, double dSplitValue)
                           dSplitValue < cToSplit.GetTopRight().y - 10 * EPS_SMALL ;
    Point3d ptP00, ptP01, ptP10, ptP11 ;
    
-   if( bGoodSplitVert) {
-      if( cToSplit.GetBottomRight().x - dSplitValue > dSplitValue - cToSplit.GetBottomLeft().x) {
+   if ( bGoodSplitVert) {
+      if ( cToSplit.GetBottomRight().x - dSplitValue > dSplitValue - cToSplit.GetBottomLeft().x) {
          GetPoint( cToSplit.GetBottomLeft().x, cToSplit.GetBottomLeft().y, ptP00) ;
          GetPoint( dSplitValue, cToSplit.GetBottomRight().y, ptP10) ;
          GetPoint( cToSplit.GetTopLeft().x, cToSplit.GetTopLeft().y, ptP01) ;
@@ -474,12 +476,13 @@ Tree::Split( int nId, double dSplitValue)
          GetPoint( dSplitValue, cToSplit.GetTopLeft().y, ptP01) ;
          GetPoint( cToSplit.GetTopRight().x, cToSplit.GetTopRight().y, ptP11) ;
       }
-      if( AreSamePointApprox( ptP00, ptP10)  &&  AreSamePointApprox( ptP01, ptP11)  &&
-         ( cToSplit.GetBottomRight().x - dSplitValue < SBZ_TREG_COEFF - EPS_SMALL  ||  dSplitValue - cToSplit.GetBottomLeft().x < SBZ_TREG_COEFF - EPS_SMALL))
+      if ( AreSamePointApprox( ptP00, ptP10)  &&  AreSamePointApprox( ptP01, ptP11)  &&
+          ( cToSplit.GetBottomRight().x - dSplitValue < SBZ_TREG_COEFF - EPS_SMALL  ||
+            dSplitValue - cToSplit.GetBottomLeft().x < SBZ_TREG_COEFF - EPS_SMALL))
          bGoodSplitVert = false ;
    }
-   if( bGoodSplitHoriz) {
-      if( cToSplit.GetTopLeft().y - dSplitValue > dSplitValue - cToSplit.GetBottomLeft().y) {
+   if ( bGoodSplitHoriz) {
+      if ( cToSplit.GetTopLeft().y - dSplitValue > dSplitValue - cToSplit.GetBottomLeft().y) {
          GetPoint( cToSplit.GetBottomLeft().x, cToSplit.GetBottomLeft().y, ptP00) ;
          GetPoint( cToSplit.GetBottomRight().x, cToSplit.GetBottomRight().y, ptP10) ;
          GetPoint( cToSplit.GetTopLeft().x, dSplitValue, ptP01) ;
@@ -491,8 +494,9 @@ Tree::Split( int nId, double dSplitValue)
          GetPoint( cToSplit.GetTopLeft().x, cToSplit.GetTopLeft().y, ptP01) ;
          GetPoint( cToSplit.GetTopRight().x, cToSplit.GetTopRight().y, ptP11) ;
       }
-      if( AreSamePointApprox( ptP00, ptP01)  &&  AreSamePointApprox( ptP10, ptP11)  &&
-          ( cToSplit.GetTopLeft().y - dSplitValue < SBZ_TREG_COEFF - EPS_SMALL  ||  dSplitValue - cToSplit.GetBottomLeft().y < SBZ_TREG_COEFF - EPS_SMALL))
+      if ( AreSamePointApprox( ptP00, ptP01)  &&  AreSamePointApprox( ptP10, ptP11)  &&
+          ( cToSplit.GetTopLeft().y - dSplitValue < SBZ_TREG_COEFF - EPS_SMALL  ||
+            dSplitValue - cToSplit.GetBottomLeft().y < SBZ_TREG_COEFF - EPS_SMALL))
          bGoodSplitHoriz = false ;
    }
    
@@ -662,8 +666,18 @@ Tree::BuildTree( double dLinTol, double dSideMin, double dSideMax)
                }
             }
 
+            // calcolo se la parte di superficie nella cella è piatta
+            PolyLine PL ;
+            PL.AddUPoint( 0, ptP00) ;
+            PL.AddUPoint( 1, ptP10) ;
+            PL.AddUPoint( 2, ptP11) ;
+            PL.AddUPoint( 3, ptP01) ;
+            PL.AddUPoint( 4, ptCen) ;
+            Plane3d plPlane ;
+            bool bIsFlat = PL.IsFlat( plPlane, dLinTol) ;
+
            // su isoparametriche in U e V
-            if ( dSagU < dLinTol  &&  dSagV < dLinTol) {
+            if ( dSagU < dLinTol  &&  dSagV < dLinTol  &&  ! bIsFlat) {
               // step di verifica in U e in V
                int nStepU = ( dLenParU > 1. / m_nDegU ? m_nDegU + 1 : 2) ;
                int nStepV = ( dLenParV > 1. / m_nDegV ? m_nDegV + 1 : 2) ;
@@ -724,6 +738,73 @@ Tree::BuildTree( double dLinTol, double dSideMin, double dSideMax)
                //{ string sLog = " Da Isoparam : FrecciaU=" + ToString( dSagU, 3) + " FrecciaV=" + ToString( dSagV, 3) ;
                //  LOG_DBG_INFO( GetEGkLogger(), sLog.c_str())}
             }
+            else if ( dSagU < dLinTol  &&  dSagV < dLinTol  &&  bIsFlat) {
+               // se la cella è piatta devo verificare che i bordi siano dei tratti retti, altrimenti potrei commettere un errore di approssimazione
+               
+               // bordo inferiore e superiore
+               double dMaxDist = 0 ;
+               for ( int i = 0 ; i < 2 ; ++i) {
+                  CurveLine clU ;
+                  if ( i == 0)
+                     clU.Set( ptP00, ptP10) ;
+                  else if ( i == 1)
+                     clU.Set( ptP01, ptP11) ;
+                  double dV = 0 ;
+                  if ( i == 0)
+                     dV = pcToSplit->GetBottomLeft().y ;
+                  else if ( i == 1)
+                     dV = pcToSplit->GetTopRight().y ;
+
+
+                  int nStepU = 4 ;
+                  for ( int j = 1 ; j < nStepU ; ++ j) {
+                    // parametro U
+                     double dCoeffU = double( j) / nStepU ;
+                     double dU = ( 1 - dCoeffU) * pcToSplit->GetBottomLeft().x + dCoeffU * pcToSplit->GetTopRight().x ;
+                     Point3d ptBez ;
+                     GetPoint( dU, dV, ptBez) ;
+                     DistPointCurve dpc( ptBez, clU) ;
+                     double dDist = 0 ;
+                     dpc.GetDist( dDist) ;
+                     if ( dDist > dMaxDist)
+                        dMaxDist = dDist ;
+                  }
+               }
+               if ( dMaxDist > dLinTol)
+                  dSagU = dMaxDist ;
+
+               // bordo sinistro e destro
+               dMaxDist = 0 ;
+               for ( int i = 0 ; i < 2 ; ++i) {
+                  CurveLine clV ;
+                  if ( i == 0)
+                     clV.Set( ptP00, ptP01) ;
+                  else if ( i == 1)
+                     clV.Set( ptP10, ptP11) ;
+                  double dU = 0 ;
+                  if ( i == 0)
+                     dU = pcToSplit->GetBottomLeft().x ;
+                  else if ( i == 1)
+                     dU = pcToSplit->GetTopRight().x ;
+
+                  int nStepV = 4 ;
+                  for ( int j = 1 ; j < nStepV ; ++ j) {
+                    // parametro in V
+                     double dCoeffV = double( j) / nStepV ;
+                     double dV = ( 1 - dCoeffV) * pcToSplit->GetBottomLeft().y + dCoeffV * pcToSplit->GetTopRight().y ;
+                     Point3d ptBez ;
+                     GetPoint( dU, dV, ptBez) ;
+                     DistPointCurve dpc( ptBez, clV) ;
+                     double dDist = 0 ;
+                     dpc.GetDist( dDist) ;
+                     if ( dDist > dMaxDist)
+                        dMaxDist = dDist ;
+                  }
+               }
+               if ( dMaxDist > dLinTol)
+                  dSagV = dMaxDist ;
+            }
+
 
            // per lo split scelgo la direzione che è più vicina alla superficie originale nel punto di maggior distanza
            // misura approssimativa della curvatura in una direzione
@@ -734,7 +815,8 @@ Tree::BuildTree( double dLinTol, double dSideMin, double dSideMax)
                bVert = false ;
             else
                bVert = ( dSagV <= dSagU) ;
-            pcToSplit->SetSplitDirVert( bVert) ;
+            bool bFirstTry = true ;
+            retry :
             // verifico che la cella sia da splittare e che eventualmente sia abbastanza grande da poterlo fare
             double dSideMinVal = 0 ;
             double dLengMinVal = 0 ;
@@ -776,6 +858,11 @@ Tree::BuildTree( double dLinTol, double dSideMin, double dSideMax)
             }
             else if ( dSagV > dLinTol  ||  dSagU > dLinTol) {
                bSplit = bDimOk ;
+               if ( ! bSplit  &&  bFirstTry) {
+                  bFirstTry = false ;
+                  bVert = ! bVert ;
+                  goto retry ;
+               }
                //if ( bSplit)
                //   LOG_DBG_INFO( GetEGkLogger(), "  Split by SagittaUV")
             }
@@ -1455,7 +1542,7 @@ Tree::GetPolygons( POLYLINEMATRIX& vvPolygons, POLYLINEMATRIX& vvPolygons3d, vec
             ++ nPolyInd ;
             continue ;
          }
-         else if( m_mTree[nId].m_nCollapsed != Cell::Collapsed::NO_COLLAPSE)
+         else if ( m_mTree[nId].m_nCollapsed != Cell::Collapsed::NO_COLLAPSE)
             continue ;
          else {
             // vettore in cui salvo il chunk di appartenenza di ogni loop che attraversa la cella
@@ -1980,7 +2067,7 @@ Tree::FindCell( const Point3d& ptToAssign, const CurveLine& cl, INTVECTOR vCells
             nCells.push_back( nCell) ;
          nEdge = -2 ;
       }
-      if( ssize( nCells) == 1)
+      if ( ssize( nCells) == 1)
          return nCells ;
 
       Vector3d vtDir ;
@@ -2009,14 +2096,14 @@ Tree::FindCell( const Point3d& ptToAssign, const CurveLine& cl, INTVECTOR vCells
          if ( abs(vtDir.x) < 1 - EPS_SMALL/100  &&  abs(vtDir.y) < 1 - EPS_SMALL/100 )
             ptIntersPlus = ptIntersPlus + vtDir * EPS_SMALL ;
          // altrimenti ruoto a destra
-         else if( ( nEdge == 4 && vtDir.x > 1 - EPS_SMALL / 100)  ||  ( nEdge == 6 && vtDir.x < - 1 + EPS_SMALL / 100) ||
+         else if (( nEdge == 4 && vtDir.x > 1 - EPS_SMALL / 100)  ||  ( nEdge == 6 && vtDir.x < - 1 + EPS_SMALL / 100) ||
                   ( nEdge == 5 && vtDir.y > 1 - EPS_SMALL / 100)  ||  ( nEdge == 7 && vtDir.y < - 1 + EPS_SMALL / 100)) {
             Vector3d vtDirDX = vtDir ; vtDirDX.Rotate( Z_AX, -45) ;
             ptIntersPlus = ptIntersPlus + vtDirDX * EPS_SMALL ;
          }
          // altrimenti ruoto a sinistra
-         else /*if( ( nEdge == 4 && vtDir.y < - 1 + EPS_SMALL / 100)  ||  ( nEdge == 6 && vtDir.y < 1 - EPS_SMALL / 100) ||
-                  ( nEdge == 5 && vtDir.x > 1 - EPS_SMALL / 100)  ||  ( nEdge == 7 && vtDir.x < - 1 + EPS_SMALL / 100))    // + tutti gli altri casi */ {
+         else /*if (( nEdge == 4 && vtDir.y < - 1 + EPS_SMALL / 100)  ||  ( nEdge == 6 && vtDir.y < 1 - EPS_SMALL / 100) ||
+                    ( nEdge == 5 && vtDir.x > 1 - EPS_SMALL / 100)  ||  ( nEdge == 7 && vtDir.x < - 1 + EPS_SMALL / 100))    // + tutti gli altri casi */ {
             Vector3d vtDirDX = vtDir ; vtDirDX.Rotate( Z_AX, 45) ;
             ptIntersPlus = ptIntersPlus + vtDirDX * EPS_SMALL ;
          }
@@ -2138,7 +2225,7 @@ Tree::TraceLoopLabelCell( const POLYLINEVECTOR& vplPolygons)
       bool bLoopInside = true ;
       Point3d ptCurr ;
       int nIdPolygon = - 1;
-      if( ! pCell->m_vnPolyId.empty())
+      if ( ! pCell->m_vnPolyId.empty())
          nIdPolygon = pCell->m_vnPolyId[0] ;
       else
          return false ;
@@ -2178,9 +2265,9 @@ Tree::TraceLoopLabelCell( const POLYLINEVECTOR& vplPolygons)
             }
             // se l'intersezione e la stessa della precedente allora potrei essere entrato in un loop infinito
             // se per più volte il punto di intersezione resta più o meno lo stesso allora blocco tutto
-            if( AreSamePointEpsilon( vptInters.back(), ptLastInters, 10 * EPS_SMALL)) {
+            if ( AreSamePointEpsilon( vptInters.back(), ptLastInters, 10 * EPS_SMALL)) {
                ++ nInfiniteLoopCount ;
-               if( nInfiniteLoopCount == 4) {
+               if ( nInfiniteLoopCount == 4) {
                   LOG_ERROR( GetEGkLogger(), "Error Triangulating SurfBezier: infinte while loop occured in Tree::TraceLoopLabelCell") 
                   return false ;
                }
@@ -2189,7 +2276,7 @@ Tree::TraceLoopLabelCell( const POLYLINEVECTOR& vplPolygons)
             // aggiorno il puntatore alla cella
             pCell = &m_mTree[nId] ;
             // recupero l'indice del poligono base associato alla cella
-            if( ! pCell->m_vnPolyId.empty())
+            if ( ! pCell->m_vnPolyId.empty())
                nIdPolygon = pCell->m_vnPolyId[0] ;
             else
                return false ;
@@ -4025,7 +4112,8 @@ Tree::OnWhichEdge( int nId, const Point3d& ptToAssign, int& nEdge) const
    Point3d ptTl ( ptBL.x, ptTR.y) ;
    Point3d ptBr ( ptTR.x, ptBL.y) ;
 
-   if( ptToAssign.y < ptBL.y - EPS_SMALL  ||  ptToAssign.y > ptTR.y + EPS_SMALL  ||  ptToAssign.x < ptBL.x - EPS_SMALL ||  ptToAssign.x > ptTR.x + EPS_SMALL)
+   if ( ptToAssign.y < ptBL.y - EPS_SMALL  ||  ptToAssign.y > ptTR.y + EPS_SMALL  || 
+        ptToAssign.x < ptBL.x - EPS_SMALL  ||  ptToAssign.x > ptTR.x + EPS_SMALL)
       return false ;
    else if ( AreSamePointXYApprox( ptToAssign, ptTR))
       nEdge = 7 ;

Tree.h

@@ -23,7 +23,8 @@
 #include <utility>
 
 struct PairHashInt64 {
-   size_t operator()(const std::pair<int64_t, int64_t>& key) const {
+   size_t
+   operator()( const std::pair<int64_t, int64_t>& key) const {
       size_t h1 = std::hash<int64_t>{}(key.first) ;
       size_t h2 = std::hash<int64_t>{}(key.second) ;
       return h1 ^ (h2 << 1);  // Combine hashes
@@ -33,14 +34,20 @@ struct PairHashInt64 {
 //----------------------------------------------------------------------------
 struct Inters {
    int nIn ;
-   PNTVECTOR vpt ;
    int nOut ;
+   PNTVECTOR vpt ;
    bool bCCW ;
    int nChunk ;
    bool bSortedbyStart ;
-   // riordino le intersezioni per lato in senso antiorario dal top
-   // se ho più intersezioni che entrano in un lato le riordino considerando che percorro i lati in senso antiorario a partire da ptTR
-   bool operator < ( Inters& b)
+
+   // nIn e nOut sono flag che indicano da quale lato ho l'ingresso e l'uscita a partire dal lato top in senso antiorario
+   // oltre il 3 sono le celle adiacenti in diagonale al vertice-> 4 corrisponde al ptTl e da lì in senso antiorario
+   // -1 se la curva è sempre dentro la cella
+
+  // riordino le intersezioni per lato in senso antiorario dal top
+  // se ho più intersezioni che entrano in un lato le riordino considerando che percorro i lati in senso antiorario a partire da ptTR
+   bool
+   operator < ( Inters& b)
    {
      // trovo in che ordine stanno i due start, tenendo conto anche della possibilità che siano vertici
       INTVECTOR vEdges = { 7, 0, 4, 1, 5, 2, 6, 3} ;
@@ -75,7 +82,8 @@ struct Inters {
                ( bEqIn  &&  nEdgeIn == 3  &&  vpt[0].y < b.vpt[0].y)) ;
    }
 
-   static bool FirstEncounter( Inters& a, Inters& b)
+   static bool
+   FirstEncounter( Inters& a, Inters& b)
    {
      // riordino in base al lato toccato, o dall'uscita o  dall'ingresso, che viene prima.
      // ottengo l'ordine che avrei percorrendo il bordo da ptTR e considerando i loop che incontro, indipendentemente se li incontro nel punto di uscita o ingresso
@@ -135,18 +143,18 @@ struct Inters {
               ( nPos1 == 3  &&  a.vpt[nFirstA].y < b.vpt[nFirstB].y) ;
    }
 
-   bool operator == ( Inters& b)
+   bool
+   operator == ( Inters& b)
    {
       return AreSamePointExact( vpt[0], b.vpt[0]) ;
    }
-   bool operator != ( Inters& b)
+
+   bool
+   operator != ( Inters& b)
    {
       return ! AreSamePointExact( vpt[0], b.vpt[0]) ;
    }
 } ;
-// nIn e nOut sono flag che indicano da quale lato ho l'ingresso e l'uscita a partire dal lato top in senso antiorario
-// oltre il 3 sono le celle adiacenti in diagonale al vertice-> 4 corrisponde al ptTl e da lì in senso antiorario
-// -1 se la curva è sempre dentro la cella
 
 //----------------------------------------------------------------------------
 class Cell
@@ -160,11 +168,12 @@ class Cell
    //                 |                 |
    //                 |_________________|
    //   Edge 5 ( SW)     Edge 2 (Bottom)     Edge 6 ( SE)
-   public:
-   enum Collapsed { TO_VERIFY = -1,        // da verificare
-                    NO_COLLAPSE = 0,       // non ho coppie di lati collassati
-                    VERT_EDGES = 1,        // coppia di lati verticali(1-3) sono collassati
-                    HORIZ_EDGES  = 2} ;    // coppia di lati verticali(0-2) sono collassati
+
+   public :
+      enum Collapsed { TO_VERIFY = -1,        // da verificare
+                       NO_COLLAPSE = 0,       // non ho coppie di lati collassati
+                       VERT_EDGES = 1,        // coppia di lati verticali(1-3) sono collassati
+                       HORIZ_EDGES  = 2} ;    // coppia di lati verticali(0-2) sono collassati
 
    public :
       ~Cell( void) {}
@@ -336,4 +345,4 @@ class Tree
       INTVECTOR                   m_vnParents ;        // vettore delle celle ottenute dalla divisione preliminare in singole patch
       ICRVCOMPOPOVECTOR           m_vCCLoop2D ;        // vettore che contiene le CurveCompo che rappresentano i loop di trim tenendo conto della divisione in celle
       std::vector<std::pair<BIPNTVECTOR, ChainCurves>>    m_vCEdge2D ;        // vettore che le chain che rappresentano ciò che resta degli edge originali, tenendo conto dei trim.
-} ;
+} ;

VolZmap.h

@@ -21,9 +21,8 @@
 #include "/EgtDev/Include/EGkSurfTriMesh.h"
 #include <unordered_map>
 #include <stack>
-#include <mutex>
-#include <atomic>
 #include <tuple>
+#include <atomic>
 
 typedef std::pair<Point3d, Vector3d>  PNTVEC3D ;
 typedef std::vector<PNTVEC3D>        PNTVEC3DVECTOR ;  // vettore di intersezioni punto, vettore, tipo superficie
@@ -557,9 +556,9 @@ class VolZmap : public IVolZmap, public IGeoObjRW
       mutable std::vector<InterVoxMatter> m_SliceXY ;
       mutable std::vector<InterVoxMatter> m_SliceXZ ;
       mutable std::vector<InterVoxMatter> m_SliceYZ ;
-      mutable std::mutex                  m_SliceMutex ;
+      mutable std::atomic_flag            m_SliceFlag ;
 
-      std::atomic<bool>                   m_bIsBox ;
+      bool                                m_bIsBox ;
 
       int                                 m_nCurrTool ;
       std::vector<Tool>                   m_vTool ;

VolZmapGraphics.cpp

@@ -1271,13 +1271,15 @@ VolZmap::ExtMarchingCubes( int nBlock, VoxelContainer& vVox) const
                         bDefTopology = true ;
                      }
                      if ( GetBlockNFromIJK( nSlBlockIJK, nSlBlockN)) {
-                        m_SliceMutex.lock() ;
+                        while ( m_SliceFlag.test_and_set( memory_order_acquire))
+                           m_SliceFlag.wait( true, memory_order_relaxed) ;
                         auto it = m_SliceYZ[nSlBlockN].find( nSliceN) ;
                         if ( it != m_SliceYZ[nSlBlockN].end()) {
                            bMatOnSlice = it->second ;
                            bDefTopology = true ;
                         }
-                        m_SliceMutex.unlock() ;
+                        m_SliceFlag.clear( memory_order_release) ;
+                        m_SliceFlag.notify_one() ;
                      }
                   }
                   else if ( abs( nAdjVox3[nCount]) == 2) {
@@ -1287,13 +1289,15 @@ VolZmap::ExtMarchingCubes( int nBlock, VoxelContainer& vVox) const
                         bDefTopology = true ;
                      }
                      if ( GetBlockNFromIJK( nSlBlockIJK, nSlBlockN)) {
-                        m_SliceMutex.lock() ;
+                        while ( m_SliceFlag.test_and_set( memory_order_acquire))
+                           m_SliceFlag.wait( true, memory_order_relaxed) ;
                         auto it = m_SliceXZ[nSlBlockN].find( nSliceN) ;
                         if ( it != m_SliceXZ[nSlBlockN].end()) {
                            bMatOnSlice = it->second ;
                            bDefTopology = true ;
                         }
-                        m_SliceMutex.unlock() ;
+                        m_SliceFlag.clear( memory_order_release) ;
+                        m_SliceFlag.notify_one() ;
                      }
                   }
                   else if ( abs( nAdjVox3[nCount]) == 3) {
@@ -1303,13 +1307,15 @@ VolZmap::ExtMarchingCubes( int nBlock, VoxelContainer& vVox) const
                         bDefTopology = true ;
                      }
                      if ( GetBlockNFromIJK( nSlBlockIJK, nSlBlockN)) {
-                        m_SliceMutex.lock() ;
+                        while ( m_SliceFlag.test_and_set( memory_order_acquire))
+                           m_SliceFlag.wait( true, memory_order_relaxed) ;
                         auto it = m_SliceXY[nSlBlockN].find( nSliceN) ;
                         if ( it != m_SliceXY[nSlBlockN].end()) {
                            bMatOnSlice = it->second ;
                            bDefTopology = true ;
                         }
-                        m_SliceMutex.unlock() ;
+                        m_SliceFlag.clear( memory_order_release) ;
+                        m_SliceFlag.notify_one() ;
                      }
                   }
                }
@@ -1374,27 +1380,33 @@ VolZmap::ExtMarchingCubes( int nBlock, VoxelContainer& vVox) const
                      if ( nSlBlockN == nBlock)
                         SliceYZ.emplace( nSliceN, bMatOnSlice) ;
                      else {
-                        m_SliceMutex.lock() ;
+                        while ( m_SliceFlag.test_and_set( memory_order_acquire))
+                           m_SliceFlag.wait( true, memory_order_relaxed) ;
                         m_SliceYZ[nSlBlockN].emplace( nSliceN, bMatOnSlice) ;
-                        m_SliceMutex.unlock() ;
+                        m_SliceFlag.clear( memory_order_release) ;
+                        m_SliceFlag.notify_one() ;
                      }
                   }
                   else if ( abs(nAdjVox3[nCount]) == 2) {
                      if ( nSlBlockN == nBlock)
                         SliceXZ.emplace( nSliceN, bMatOnSlice) ;
                      else {
-                        m_SliceMutex.lock() ;
+                        while ( m_SliceFlag.test_and_set( memory_order_acquire))
+                           m_SliceFlag.wait( true, memory_order_relaxed) ;
                         m_SliceXZ[nSlBlockN].emplace( nSliceN, bMatOnSlice) ;
-                        m_SliceMutex.unlock() ;
+                        m_SliceFlag.clear( memory_order_release) ;
+                        m_SliceFlag.notify_one() ;
                      }
                   }
                   else if ( abs(nAdjVox3[nCount]) == 3) {
                      if ( nSlBlockN == nBlock)
                         SliceXY.emplace(nSliceN, bMatOnSlice) ;
                      else {
-                        m_SliceMutex.lock() ;
+                        while ( m_SliceFlag.test_and_set( memory_order_acquire))
+                           m_SliceFlag.wait( true, memory_order_relaxed) ;
                         m_SliceXY[nSlBlockN].emplace( nSliceN, bMatOnSlice) ;
-                        m_SliceMutex.unlock() ;
+                        m_SliceFlag.clear( memory_order_release) ;
+                        m_SliceFlag.notify_one() ;
                      }
                   }
                }

VolZmapVolume.cpp

@@ -1585,7 +1585,6 @@ VolZmap::Comp_5AxisMilling( int nGrid, const Point3d& ptS, const Point3d& ptE, c
    Vector3d vtDirTip = ptP2T - ptP1T ;
    bool bTopIsPivot = vtDirTop.IsSmall() ;
    bool bTipIsPivot = vtDirTip.IsSmall() ;
-   bool bTopAndTipAreEquiverse = vtDirTop * vtDirTip > 0 ;
 
    // box dell'intero volume spazzato, nel riferimento object oriented
    BBox3d bbVol ;
@@ -1602,9 +1601,9 @@ VolZmap::Comp_5AxisMilling( int nGrid, const Point3d& ptS, const Point3d& ptE, c
    PNTVECTOR vPntTopStartBack(3) ;
    PNTVECTOR vPntTopEndBack(3) ;
    if ( nSub > 1) {
-      if( bTopIsPivot)
+      if ( bTopIsPivot)
          vtDirTop = vtDirTip ;
-      if( bTipIsPivot)
+      if ( bTipIsPivot)
          vtDirTip = vtDirTop ;
       // determino in che modo collegare il cilindro iniziale con quello finale
       Vector3d vtTopBaseEnd = vtDirTop - (( vtDirTop * vtLe) * vtLe) ;
@@ -1813,7 +1812,7 @@ VolZmap::Comp_5AxisMilling( int nGrid, const Point3d& ptS, const Point3d& ptE, c
                }
             }
             
-            if( ! bTopIsPivot) {
+            if ( ! bTopIsPivot) {
                // superiori
                if ( dSide > 0) {
                   PNTVECTOR vPntTop01 = cBezTopStartF1->GetAllControlPoints();