//----------------------------------------------------------------------------
//  EgalTech 2013-2013
//----------------------------------------------------------------------------
// File : CurveAux.cpp      Data : 22.11.13          Versione : 1.3a1
// Contenuto : Implementazione di alcune funzioni di utilità per le curve.
//
//
//
// Modifiche : 22.11.13 DS Creazione modulo.
//
//
//----------------------------------------------------------------------------

//--------------------------- Include ----------------------------------------
#include "stdafx.h"
#include "CurveAux.h"
#include "GeoConst.h"
#include "CurveLine.h"
#include "CurveArc.h"
#include "CurveBezier.h"
#include "CurveComposite.h"
#include "/EgtDev/Include/EGkDistPointCurve.h"
#include "/EgtDev/Include/EGkStringUtils3d.h"
#include "/EgtDev/Include/EGkUiUnits.h"
#include "/EgtDev/Include/EgtPointerOwner.h"

using namespace std ;

//----------------------------------------------------------------------------
bool
IsClosed( const ICurve& crvC)
{
   Point3d ptStart, ptEnd ;
   return ( crvC.GetStartPoint( ptStart)  &&  crvC.GetEndPoint( ptEnd)  &&  AreSamePointApprox( ptStart, ptEnd)) ;
}

//----------------------------------------------------------------------------
bool
IsValidParam( const ICurve& crvC, double dPar, ICurve::Side nSide)
{
  // recupero l'intervallo di validità del parametro
   double dStart, dEnd ;
   if ( ! crvC.GetDomain( dStart, dEnd))
      return false ;
  // il parametro deve stare nei limiti
   if ( dPar < dStart - EPS_PARAM  ||  dPar > dEnd + EPS_PARAM)
      return false ;
  // se curva chiusa non sono necessari altri controlli
   if ( crvC.IsClosed())
      return true ;
  // per curve aperte non posso studiare la curva prima dell'inizio e dopo la fine
   if ( ( dPar < dStart + EPS_PARAM  &&  nSide == ICurve::FROM_MINUS)  ||
        ( dPar > dEnd - EPS_PARAM  &&  nSide == ICurve::FROM_PLUS))
      return false ;
   return true ;
}

//----------------------------------------------------------------------------
bool
IsStartParam( const ICurve& crvC, double dPar)
{
  // recupero l'intervallo di validità del parametro
   double dStart, dEnd ;
   if ( ! crvC.GetDomain( dStart, dEnd))
      return false ;
  // se il parametro  non è nell'intorno dell'inizio
   if ( abs( dPar - dStart) > EPS_PARAM)
      return false ;
   return true ;
}

//----------------------------------------------------------------------------
bool
IsEndParam( const ICurve& crvC, double dPar)
{
  // recupero l'intervallo di validità del parametro
   double dStart, dEnd ;
   if ( ! crvC.GetDomain( dStart, dEnd))
      return false ;
  // se il parametro  non è nell'intorno della fine
   if ( abs( dPar - dEnd) > EPS_PARAM)
      return false ;
   return true ;
}

//----------------------------------------------------------------------------
bool
GetNearestExtremityToPoint( const Point3d& ptP, const ICurve& Curve, bool& bStart)
{
  // recupero punto iniziale
   Point3d ptStart ;
   if ( ! Curve.GetStartPoint( ptStart))
      return false ;
  // recupero punto finale
   Point3d ptEnd ;
   if ( ! Curve.GetEndPoint( ptEnd))
      return false ;
  // se curva aperta
   if ( ! AreSamePointApprox( ptStart, ptEnd)) {
     // è il più vicino tra inizio e fine
      bStart = ( SqDist( ptP, ptStart) <= SqDist( ptP, ptEnd)) ;
      return true ;
   }
  // altrimenti curva chiusa, devo proiettare il punto sulla curva e misurare la distanza su di essa
   DistPointCurve distPC( ptP, Curve) ;
   int nFlag ;
   double dLen, dU, dULen ;
   if ( Curve.GetLength( dLen)  &&
        distPC.GetParamAtMinDistPoint( 0, dU, nFlag)  &&  Curve.GetLengthAtParam( dU, dULen)) {
      bStart = ( dULen <= ( dLen - dULen)) ;
      return true ; 
   }
   else
      return false ;
}

//----------------------------------------------------------------------------
bool
MoveParamToAvoidTg( double& dU, ICurve::Side nSide, const ICurve& Curve)
{
  // verifico che il parametro sia accettabile
   if ( ! Curve.IsValidParam( dU, nSide))
      return false ;

  // determino un incremento piccolo ma valido del parametro U
   double dDeltaU = 1000 * EPS_PARAM ;
   Point3d ptPos ;
   Vector3d vtDer1 ;
   if ( Curve.GetPointD1D2( dU, nSide, ptPos, &vtDer1)) {
      double dDer1 = vtDer1.LenXY() ;
      if ( dDer1 > EPS_ZERO)
         dDeltaU = 2 * EPS_SMALL / dDer1 ;
   }

  // cerco di spostarmi per evitare eventuali problemi di tangenza
   if ( nSide == ICurve::FROM_MINUS) {
      double dUm = dU - dDeltaU ;
      if ( ! Curve.IsValidParam( dUm, nSide)) {
         if ( Curve.IsClosed()) {
            double dStart, dEnd ;
            Curve.GetDomain( dStart, dEnd) ;
            dUm = ( dEnd - dStart) - dDeltaU ;
         }
         else
            dUm = dU ;
      }
      dU = dUm ;
      return true ;
   }
   else {  // nSide == ICurve::FROM_PLUS
      double dUm = dU + dDeltaU ;
      if ( ! Curve.IsValidParam( dUm, nSide)) {
         if ( Curve.IsClosed()) {
            double dStart, dEnd ;
            Curve.GetDomain( dStart, dEnd) ;
            dUm = dStart + dDeltaU ;
         }
         else
            dUm = dU ;
      }
      dU = dUm ;
      return true ;
   }
}

//----------------------------------------------------------------------------
bool
GetTang( const ICurve& crvC, double dU, ICurve::Side nS, Vector3d& vtTang)
{
   Point3d ptPos ;
   return GetPointTang( crvC, dU, nS, ptPos, vtTang) ;
}

//----------------------------------------------------------------------------
bool
GetPointTang( const ICurve& crvC, double dU, ICurve::Side nS, Point3d& ptPos, Vector3d& vtTang)
{
   if ( ! crvC.GetPointD1D2( dU, nS, ptPos, &vtTang))
      return false ;
   if ( vtTang.Normalize( EPS_ZERO))
      return true ;
  // nel caso la derivata prima sia nulla, utilizziamo la seconda 
   Vector3d vtDummy ;
   if ( ! crvC.GetPointD1D2( dU, nS, ptPos, &vtDummy, &vtTang))
      return false ;
   double dUmod = dU + ( nS == ICurve::FROM_MINUS ? -1 : +1) * 100 * EPS_PARAM ;
   Point3d ptDummy ;
  // se anche la derivata seconda è nulla, provo a spostarmi di poco
   if ( ! vtTang.Normalize(  EPS_ZERO)) {
      if ( ! crvC.GetPointD1D2( dUmod, nS, ptDummy, &vtDummy, &vtTang))
         return false ;
      if ( ! vtTang.Normalize(  EPS_ZERO))
         return false ;
   }
  // per il verso, lo confrontiamo con quello della derivata prima un poco discosta
   Vector3d vtD1near ;
   if ( ! crvC.GetPointD1D2( dUmod, nS, ptDummy, &vtD1near))
      return false ;
   if ( ( vtTang * vtD1near) < 0)
      vtTang.Invert() ;
   return true ;
}

//----------------------------------------------------------------------------
bool
GetPointDiffGeom( const ICurve& crvC, double dU, ICurve::Side nS, CrvPointDiffGeom& oDiffG)
{
  // calcolo punto e derivate
   Point3d  ptPos ;
   Vector3d vtDer1, vtDer2 ;
   if ( ! crvC.GetPointD1D2( dU, nS, ptPos, &vtDer1, &vtDer2))
      return false ;
  // assegno parametro e posizione
   oDiffG.dU = dU ;
   oDiffG.ptP = ptPos ;
  // se esiste la derivata prima non nulla
   double dSqLenD1 = vtDer1.SqLen() ;
   if ( vtDer1.Normalize()) {
      oDiffG.vtT = vtDer1 ;
     // del vettore deriv2^ tengo la sola componente perpendicolare al vettore tangente
      oDiffG.vtN = vtDer2 - ( vtDer2 * vtDer1) * vtDer1 ;
      if ( oDiffG.vtN.Normalize()) {
         oDiffG.dCurv = ( vtDer1 ^ vtDer2).Len() / dSqLenD1 ;
         oDiffG.nStatus = CrvPointDiffGeom::NCRV ;
      }
      else {
         oDiffG.dCurv = 0 ;
         oDiffG.nStatus = CrvPointDiffGeom::TANG ;
      }
   }
  // se esiste almeno derivata seconda non nulla, definisce la tangente
   else if ( vtDer2.Normalize()) {
      oDiffG.vtT = vtDer2 ;
      oDiffG.vtN.Set( 0, 0, 0) ;
      oDiffG.dCurv = 0 ;
     // per il verso, lo confrontiamo con quello della derivata prima un poco discosta
      Point3d ptDummy ;
      Vector3d vtD1near ;
      dU += ( nS == ICurve::FROM_MINUS ? -1 : +1) * 100 * EPS_PARAM ;
      if ( crvC.GetPointD1D2( dU, nS, ptDummy, &vtD1near)) {
         if ( ( oDiffG.vtT * vtD1near) < 0)
            oDiffG.vtT.Invert() ;
         oDiffG.nStatus = CrvPointDiffGeom::TANG ;
      }
      else {
         oDiffG.vtT.Set( 0, 0, 0) ;
         oDiffG.nStatus = CrvPointDiffGeom::POS ;
      }
   }
  // altrimenti non sono definite la tangente e la normale
   else {
      oDiffG.vtT.Set( 0, 0, 0) ;
      oDiffG.vtN.Set( 0, 0, 0) ;
      oDiffG.dCurv = 0 ;
      oDiffG.nStatus = CrvPointDiffGeom::POS ;
   }

   return true ;
}

//----------------------------------------------------------------------------
bool
ImproveCurveParamAtPoint( double& dU, const Point3d& ptP, const ICurve* pCrv)
{
  // Da usare quando il parametro è già molto vicino a quello esatto

  // calcolo il punto della curva in corrispondenza al parametro
   Point3d ptQ ;
   Vector3d vtD ;
   if ( ! pCrv->GetPointD1D2( dU, ICurve::FROM_MINUS, ptQ, &vtD))
      return false ;
  // se sono uguali, è già tutto ok
   if ( AreSamePointExact( ptP, ptQ))
      return true ;
  // se derivata nulla, non posso migliorare
   if ( vtD.IsZero())
      return false ;
  // incremento parametro su linea tangente (uso la derivata)
   double dDeltaU = ( ptP - ptQ) * vtD / vtD.SqLen() ;
  // se migliorato, tengo nuovo valore
   Point3d ptR ;
   if ( pCrv->GetPointD1D2( dU + dDeltaU, ICurve::FROM_MINUS, ptR)  &&
        ( ptP - ptR).SqLen() < ( ptP - ptQ).SqLen())
      dU += dDeltaU ;
   return true ;
}

//----------------------------------------------------------------------------
bool
CurveGetAreaXY( const ICurve& crvC, double& dArea)
{
  // verifico sia chiusa
   if ( ! crvC.IsClosed())
      return false ;
  // approssimo la curva con una polilinea
   PolyLine PL ;
   crvC.ApproxWithLines( LIN_TOL_STD, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL_INT, PL) ;
  // calcolo l'area
   double dAreaXY = 0 ;
   PL.GetAreaXY( dAreaXY) ;
  // restituisco il valore
   if ( &dArea != nullptr)
      dArea = dAreaXY ;
   return true ;
}

//----------------------------------------------------------------------------
bool
CurveGetArea( const ICurve& crvC, Plane3d& plPlane, double& dArea)
{
  // verifico sia chiusa
   if ( ! crvC.IsClosed())
      return false ;
  // approssimo la curva con una polilinea
   PolyLine PL ;
   crvC.ApproxWithLines( LIN_TOL_STD, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL_INT, PL) ;
  // calcolo l'area
   Plane3d plMyPlane ;
   double dMyArea = 0 ;
   if ( ! PL.IsClosedAndFlat( plMyPlane, dMyArea, 50 * EPS_SMALL))
      return false ;
  // restituisco i valori
   if ( &plPlane != nullptr)
      plPlane = plMyPlane ;
   if ( &dArea != nullptr)
      dArea = dMyArea ;
   return true ;
}

//----------------------------------------------------------------------------
bool
CurveDump( const ICurve& crvC, string& sOut, bool bMM, const char* szNewLine)
{
  // verifico validità curva
   if ( ! crvC.IsValid())
      return false ;

  // punti iniziale e finale
   Point3d ptPS ;
   crvC.GetStartPoint( ptPS) ;
   sOut += "PS(" + ToString( GetInUiUnits( ptPS, bMM), 3) + ")" + szNewLine ;
   Point3d ptPE ;
   crvC.GetEndPoint( ptPE) ;
   sOut += "PE(" + ToString( GetInUiUnits( ptPE, bMM), 3) + ")" + szNewLine ;
  // direzioni iniziale e finale
   Vector3d vtDirS ;
   crvC.GetStartDir( vtDirS) ;
   sOut += "VS(" + ToString( vtDirS, 3) + ")" + szNewLine ; 
   Vector3d vtDirE ;
   crvC.GetEndDir( vtDirE) ;
   sOut += "VE(" + ToString( vtDirE, 3) + ")" + szNewLine ;
  // eventuali direzione di estrusione e spessore
   Vector3d vtExtr ;
   crvC.GetExtrusion( vtExtr) ;
   if ( ! vtExtr.IsSmall()) {
     double dThick ;
     crvC.GetThickness( dThick) ;
     sOut += "Extr(" + ToString( vtExtr, 3) + ") Th=" + ToString( GetInUiUnits( dThick, bMM), 3) + szNewLine ;
   }
  // lunghezza
   double dLen = 0 ;
   crvC.GetLength( dLen) ;
   sOut += "Len=" + ToString( GetInUiUnits( dLen, bMM), 3) + szNewLine ;
  // altri dati per curva chiusa
   double dAreaXY ;
   if ( CurveGetAreaXY( crvC, dAreaXY)) {
      bool bCCW = ( dAreaXY > 0) ;
      double dAreaUi = GetAreaInUiUnits( abs( dAreaXY), bMM) ;
      int nDec = ( dAreaUi > 100 ? 1 : ( dAreaUi > 0.1 ? 3 : 6)) ;
      sOut += string( "Closed") + ( bCCW ? " CCW" : " CW") + " AreaXY=" +
              ToString( dAreaUi, nDec) + szNewLine ;
   }

   return true ;
}

//----------------------------------------------------------------------------
bool
CopyExtrusion( const ICurve* pSouCrv, ICurve* pDestCrv)
{
  // verifico puntatori
   if ( pSouCrv == nullptr  ||  pDestCrv == nullptr)
      return false ;
  // eseguo copia
   Vector3d vtExtr ;
   pSouCrv->GetExtrusion( vtExtr) ;
   pDestCrv->SetExtrusion( vtExtr) ;
   return true ;
}

//----------------------------------------------------------------------------
bool
CopyThickness( const ICurve* pSouCrv, ICurve* pDestCrv)
{
  // verifico puntatori
   if ( pSouCrv == nullptr  ||  pDestCrv == nullptr)
      return false ;
  // eseguo copia
   double dThick = 0 ;
   pSouCrv->GetThickness( dThick) ;
   pDestCrv->SetThickness( dThick) ;
   return true ;
}

//----------------------------------------------------------------------------
ICurve*
ArcToBezierCurve( const ICurve* pCrv)
{
  // verifico sia un arco
   const CurveArc* pArc = GetBasicCurveArc( pCrv) ;
   if ( pArc == nullptr)
      return nullptr ;

  // se angolo al centro sotto il limite, basta una curva
   if ( abs( pArc->GetAngCenter()) < BEZARC_ANG_CEN_MAX + EPS_ANG_SMALL) {
     // creo la curva di Bezier equivalente all'arco
      PtrOwner<CurveBezier> pCrvBez( CreateBasicCurveBezier()) ;
      if ( IsNull( pCrvBez)  ||  ! pCrvBez->FromArc( *pArc))
         return nullptr ;
     // restituisco la curva
      return Release( pCrvBez) ;
   }
  // altrimenti curva composita di Bezier
   else {
     // creo la curva composita
      PtrOwner<CurveComposite> pCrvCompo( CreateBasicCurveComposite()) ;
      if ( IsNull( pCrvCompo))
         return nullptr ;
     // inserisco nella CC le curve di Bezier equivalenti alle parti dell'arco
      int nParts = (int) ceil( abs( pArc->GetAngCenter()) / ( BEZARC_ANG_CEN_MAX + EPS_ANG_SMALL)) ;
      nParts = max( nParts, 2) ;
      for ( int i = 0 ; i < nParts ; ++ i) {
        // copio l'arco originale
         CurveArc cArc = *pArc ;
        // lo limito alla parte di interesse
         cArc.TrimStartEndAtParam( i / double( nParts), ( i + 1) / double( nParts)) ;
        // creo la curva di Bezier equivalente
         PtrOwner<CurveBezier> pCrvBez( CreateBasicCurveBezier()) ;
         if ( IsNull( pCrvBez)  ||  ! pCrvBez->FromArc( cArc))
            return nullptr ;
        // aggiungo la curva di Bezier a quella composita
         if ( ! pCrvCompo->AddCurve( Release( pCrvBez)))
            return nullptr ;
      }
     // copio estrusione e spessore
      CopyExtrusion( pArc, pCrvCompo) ;
      CopyThickness( pArc, pCrvCompo) ;
     // restituisco la curva
      return Release( pCrvCompo) ;
   }
}

//----------------------------------------------------------------------------
ICurve*
CurveToNoArcsCurve( const ICurve* pCrv)
{
  // verifico validità curva
   if ( pCrv == nullptr)
      return nullptr ;
  // se arco, devo trasformarlo in curva di Bezier (semplice o composta)
   if ( pCrv->GetType() == CRV_ARC) {
      return ArcToBezierCurve( pCrv) ;
   }
  // se curva composita, devo trasformarla in composita senza archi
   else if ( pCrv->GetType() == CRV_COMPO) {
      PtrOwner<CurveComposite> pCopyCC( GetBasicCurveComposite( pCrv->Clone())) ;
      if ( IsNull( pCopyCC)  ||  ! pCopyCC->ArcsToBezierCurves())
         return nullptr ;
      return Release( pCopyCC) ;
   }
  // altrimenti devo solo copiarla
   else {
      return pCrv->Clone() ;
   }
}

//----------------------------------------------------------------------------
ICurve*
CurveToArcsPerpExtrCurve( const ICurve* pCrv, double dLinTol, double dAngTolDeg)
{
  // verifico validità curva
   if ( pCrv == nullptr)
      return nullptr ;
  // se arco in piano non perpendicolare ad estrusione, curva composita o curva di Bezier trasformo 
   if ( ( pCrv->GetType() == CRV_ARC  &&  ! GetBasicCurveArc(pCrv)->IsInPlanePerpExtr())  ||
        pCrv->GetType() == CRV_COMPO  ||
        pCrv->GetType() == CRV_BEZIER) {
     // creo un poliarco
      PolyArc PA ;
      if ( ! pCrv->ApproxWithArcs( dLinTol, dAngTolDeg, PA))
         return nullptr ;
     // creo una nuova curva composita a partire dal poliarco
      PtrOwner<CurveComposite> pCopyCC( CreateBasicCurveComposite()) ;
      if ( IsNull( pCopyCC))
         return nullptr ;
      if ( ! pCopyCC->FromPolyArc( PA))
         return nullptr ;
     // copio estrusione e spessore
      CopyExtrusion( pCrv, pCopyCC) ;
      CopyThickness( pCrv, pCopyCC) ;
      return Release( pCopyCC) ;
   }
  // altrimenti devo solo copiarla
   else {
      return pCrv->Clone() ;
   }
}

//----------------------------------------------------------------------------
bool
NurbsCurveCanonicalize( CNurbsData& cnData)
{
  // se periodica
   if ( cnData.bPeriodic) {
     // va trasformata in non-periodica (clamped)
     // vedere The NurbsBook di Les Piegl e Tiller
     // mancano esempi per testare
      return false ;
   }
  // se con nodi extra
   if ( cnData.bExtraKnotes) {
      int nKnotesNbr = int( cnData.vU.size()) ;
      if ( nKnotesNbr < 4)
         return false ;
      cnData.bExtraKnotes = false ;
      for ( int i = 0 ; i < nKnotesNbr - 2 ; ++ i)
         cnData.vU[i] = cnData.vU[i+1] ;
      cnData.vU.resize( nKnotesNbr - 2) ;
      return true ;
   }

   return true ;
}

//----------------------------------------------------------------------------
ICurve*
NurbsToBezierCurve( const CNurbsData& cnData)
{
  // la curva Nurbs deve essere in forma canonica
   if ( cnData.bPeriodic  ||  cnData.bExtraKnotes)
      return nullptr ;
  // numero dei nodi
   int nU = int( cnData.vCP.size()) + cnData.nDeg - 1 ;
  // controllo relazione nodi - punti di controllo
   if ( nU != int( cnData.vU.size()))
      return nullptr ;
  // numero degli intervalli
   int nInt = nU - 2 * cnData.nDeg + 1 ;
  // verifico le condizioni agli estremi sui nodi (i primi nDeg nodi e gli ultimi nDeg nodi devono essere uguali tra loro)
   bool bOk = true ;
   for ( int i = 1 ; i < cnData.nDeg ; ++ i) {
      if ( abs( cnData.vU[i] - cnData.vU[0]) >= EPS_ZERO)
         bOk = false ;
   }
   for ( int i = 1 ; i < cnData.nDeg ; ++ i) {
      if ( abs( cnData.vU[nU - 1 - i] - cnData.vU[nU - 1]) >= EPS_ZERO)
         bOk = false ;
   }
   if ( ! bOk)
      return nullptr ;

  // se 1 solo intervallo, la Nurbs è già una curva di Bezier
   if ( nInt == 1) {
     // creo la curva di Bezier
      PtrOwner<ICurveBezier> pCrvBez( CreateCurveBezier()) ;
      if ( IsNull( pCrvBez))
         return nullptr ;
     // la inizializzo
      if ( ! pCrvBez->Init( cnData.nDeg, cnData.bRat))
         return nullptr ;
      for ( int i = 0 ; i <= cnData.nDeg ; ++ i) {
         if ( ! cnData.bRat) {
            if ( ! pCrvBez->SetControlPoint( i, cnData.vCP[i]))
               return nullptr ;
         }
         else {
            if ( ! pCrvBez->SetControlPoint( i, cnData.vCP[i], cnData.vW[i]))
               return nullptr ;
         }
      }
     // se non è una curva ma un punto, la invalido
      if ( pCrvBez->IsAPoint())
         pCrvBez->Init( cnData.nDeg, cnData.bRat) ;
     // restituisco la curva
      return Release( pCrvBez) ;
   }

  // altrimenti è equivalente ad una curva composita, la creo
   PtrOwner<ICurveComposite> pCrvCompo( CreateCurveComposite()) ;
   if ( IsNull( pCrvCompo))
      return nullptr ;

  // vettore dei punti di controllo della curva di Bezier
   PNTVECTOR vBC ;
   vBC.resize( cnData.nDeg + 1) ;
   DBLVECTOR vBW ;
   vBW.resize( cnData.nDeg + 1) ;
   if ( ! cnData.bRat) {
      for ( int i = 0 ; i <= cnData.nDeg ; ++ i)
         vBC[i] = cnData.vCP[i] ;
   }
   else {
      for ( int i = 0 ; i <= cnData.nDeg ; ++ i) {
         vBC[i] = cnData.vCP[i] * cnData.vW[i] ;
         vBW[i] = cnData.vW[i] ;
      }
   }
  // primi coefficienti della successiva
   PNTVECTOR vNextBC ;
   vNextBC.resize( cnData.nDeg - 1) ;
   DBLVECTOR vNextBW ;
   vNextBW.resize( cnData.nDeg - 1) ;
  // ...
   DBLVECTOR vAlfa ;
   vAlfa.resize( cnData.nDeg - 1) ;
   int a = cnData.nDeg - 1 ;
   int b = cnData.nDeg ;
   bool bPrevRejected = false ;
  // ciclo
   while ( b < nU - 1) {
      int i = b ;
      while ( b < nU - 1  &&  abs( cnData.vU[b+1] - cnData.vU[b]) < EPS_ZERO)
         ++ b ;
      int mult = min( b - i + 1, cnData.nDeg) ;
      if ( mult < cnData.nDeg) {
        // numeratore di alfa
         double numer = cnData.vU[b] - cnData.vU[a] ;
        // calcola e salva gli alfa
         for ( int j = cnData.nDeg ; j > mult ; -- j)
            vAlfa[j-mult-1] = numer / ( cnData.vU[a+j] - cnData.vU[a]) ;
        // inserisco il nodo r volte
         int r = cnData.nDeg - mult ;
         for ( int j = 1 ; j <= r ; ++ j) {
            int save = r - j ;
            int s = mult + j ;
            for ( int k = cnData.nDeg ; k >= s ; -- k)
               vBC[k] = vAlfa[k-s] * vBC[k] + ( 1 - vAlfa[k-s]) * vBC[k-1] ;
            if ( cnData.bRat) {
               for ( int k = cnData.nDeg ; k >= s ; -- k)
                  vBW[k] = vAlfa[k-s] * vBW[k] + ( 1 - vAlfa[k-s]) * vBW[k-1] ;
            }
            if ( b < nU - 1) {
               vNextBC[save] = vBC[cnData.nDeg] ;
               if ( cnData.bRat)
                  vNextBW[save] = vBW[cnData.nDeg] ;
            }
         }
      }

     // costruisco la curva di Bezier e la inserisco nella curva composita
      PtrOwner<ICurveBezier> pCrvBez( CreateCurveBezier()) ;
      if ( IsNull( pCrvBez))
         return nullptr ;
     // se precedente saltata
      if ( bPrevRejected) {
        // prendo l'ultimo punto della curva composita per garantire la continuità
         Point3d ptEnd ;
         if ( pCrvCompo->GetEndPoint( ptEnd))
            vBC[0] = ptEnd ;
      }
     // la inizializzo
      if ( ! pCrvBez->Init( cnData.nDeg, cnData.bRat))
         return nullptr ;
      if ( ! cnData.bRat) {
         for ( int i = 0 ; i <= cnData.nDeg ; ++ i) {
            if ( ! pCrvBez->SetControlPoint( i, vBC[i]))
               return nullptr ;
         }
      }
      else {
         for ( int i = 0 ; i <= cnData.nDeg ; ++ i) {
            if ( ! pCrvBez->SetControlPoint( i, vBC[i] / vBW[i], vBW[i]))
               return nullptr ;
         }
      }
     // se è una vera curva, la aggiungo alla curva composita
      if ( ! pCrvBez->IsAPoint()) {
         if ( ! pCrvCompo->AddCurve( Release( pCrvBez)))
            return nullptr ;
         bPrevRejected = false ;
      }
     // altrimenti è un punto, la cancello
      else {
         pCrvBez.Reset() ;
         bPrevRejected = true ;
      }

     // inizializzazioni per la prossima curva di Bezier
      if ( b < nU - 1) {
         if ( ! cnData.bRat) {
            for ( int i = 0 ; i < cnData.nDeg - 1 ; ++ i)
               vBC[i] = vNextBC[i] ;
            for ( int i = cnData.nDeg - mult ; i <= cnData.nDeg ; ++ i)
               vBC[i] = cnData.vCP[b-cnData.nDeg+i+1] ;
         }
         else {
            for ( int i = 0 ; i < cnData.nDeg - 1 ; ++ i) {
               vBC[i] = vNextBC[i] ;
               vBW[i] = vNextBW[i] ;
            }
            for ( int i = cnData.nDeg - mult ; i <= cnData.nDeg ; ++ i) {
               vBC[i] = cnData.vCP[b-cnData.nDeg+i+1] * cnData.vW[b-cnData.nDeg+i+1] ;
               vBW[i] = cnData.vW[b-cnData.nDeg+i+1] ;
            }
         }
         a = b ;
         ++ b ;
      }
   }

  // restituisco la curva composita
   return Release( pCrvCompo) ;
}

//----------------------------------------------------------------------------
ICurve*
FlattenCurve( const ICurve& crCrv, double dToler, double dAngToler, int nFlag)
{
  // Determino il piano medio della curva e verifico se scostamento accettabile
   Plane3d plMid ;
   if ( ! crCrv.IsFlat( plMid, ( nFlag == FLTCRV_USE_EXTR), dToler))
      return nullptr ;
  // Recupero estrusione rispetto alla normale al piano medio
   Vector3d vtExtr ; crCrv.GetExtrusion( vtExtr) ;
   if ( nFlag == FLTCRV_USE_EXTR)
      plMid.Set( plMid.GetPoint(), vtExtr) ;
  // Determino il suo centro geometrico
   Point3d ptCen ;
   if ( ! crCrv.GetCentroid( ptCen))
      return nullptr ;
  // Verifico se curva già piatta
   PolyLine PL ;
   if ( ! crCrv.ApproxWithLines( LIN_TOL_FINE, ANG_TOL_STD_DEG, ICurve::APL_SPECIAL_INT, PL))
      return nullptr ;
   bool bFlat = true ;
   Plane3d plFlat ; plFlat.Set( ptCen, plMid.GetVersN()) ;
   Point3d ptP ;
   bool bPoint = PL.GetFirstPoint( ptP) ;
   while ( bPoint  &&  bFlat) {
      if ( DistPointPlane( ptP, plFlat) > 0.1 * EPS_SMALL)
         bFlat = false ;
      bPoint = PL.GetNextPoint( ptP) ;
   }
  // Se curva già piatta, la copio ed esco
   if ( bFlat) {
      PtrOwner<ICurve> pCrv( crCrv.Clone()) ;
      if ( IsNull( pCrv))
         return nullptr ;
     // assegno estrusione
      if ( nFlag != FLTCRV_SET_EXTR)
         pCrv->SetExtrusion( vtExtr) ;
      else {
         if ( vtExtr * plMid.GetVersN() >= 0)
            pCrv->SetExtrusion( plMid.GetVersN()) ;
         else
            pCrv->SetExtrusion( - plMid.GetVersN()) ;
      }
     // restituisco la copia della curva
      return Release( pCrv) ;
   }
  // altrimenti la appiattisco
   else {
     // mi assicuro che la curva non contenga archi
      PtrOwner<ICurve> pCrv( CurveToNoArcsCurve( &crCrv)) ;
      if ( IsNull( pCrv))
         return nullptr ;
     // calcolo un riferimento con piano XY coincidente con il piano di proiezione
      Frame3d frRef ;
      if ( ! frRef.Set( ptCen, plMid.GetVersN()))
         return nullptr ;
     // eseguo scalatura con fattori X e Y unitari e fattore Z nullo
      if ( ! pCrv->Scale( frRef, 1, 1, 0))
         return nullptr ;
     // assegno estrusione
      if ( nFlag != FLTCRV_SET_EXTR)
         pCrv->SetExtrusion( vtExtr) ;
      else {
         if ( vtExtr * plMid.GetVersN() >= 0)
            pCrv->SetExtrusion( plMid.GetVersN()) ;
         else
            pCrv->SetExtrusion( - plMid.GetVersN()) ;
      }
     // restituisco la nuova curva
      return Release( pCrv) ;
   }
}

//----------------------------------------------------------------------------
ICurve*
ProjectCurveOnPlane( const ICurve& crCrv, const Plane3d& plPlane)
{
  // determino se curva piana e suo eventuale piano
   Plane3d plCrv ;
   if ( crCrv.IsFlat( plCrv, false, EPS_SMALL / 2)) {
     // se il piano della curva è parallelo a quello di proiezione
      if ( AreSameOrOppositeVectorExact( plCrv.GetVersN(), plPlane.GetVersN())) {
        // copio la curva
         PtrOwner<ICurve> pCrv( crCrv.Clone()) ;
         if ( IsNull( pCrv))
            return nullptr ;
        // se non coincidenti, basta eseguire una traslazione
         Point3d ptOC = ORIG + plCrv.GetDist() * plCrv.GetVersN() ;
         Point3d ptOP = ORIG + plPlane.GetDist() * plPlane.GetVersN() ;
         if ( ! AreSamePointApprox( ptOC, ptOP))
            pCrv->Translate( ptOP - ptOC) ;
        // restituisco la nuova curva
         return Release( pCrv) ;
      }
   }

  // mi assicuro che la curva non contenga archi
   PtrOwner<ICurve> pCrv( CurveToNoArcsCurve( &crCrv)) ;
   if ( IsNull( pCrv))
      return nullptr ;
  // calcolo un riferimento con piano XY coincidente con il piano di proiezione
   Frame3d frRef ;
   if ( ! frRef.Set( ORIG + plPlane.GetDist() * plPlane.GetVersN(), plPlane.GetVersN()))
      return nullptr ;
  // eseguo scalatura con fattori X e Y unitari e fattore Z nullo
   if ( ! pCrv->Scale( frRef, 1, 1, 0))
      return nullptr ;
  // restituisco la nuova curva
   return Release( pCrv) ;
}

//----------------------------------------------------------------------------
bool
AdjustCurveSlope( ICurveComposite* pCrv, double dNini, double dNfin)
{
  // verifico curva
   if ( pCrv == nullptr)
      return false ;
  // determino versore estrusione
   Vector3d vtN ;
   if ( ! pCrv->GetExtrusion( vtN)  ||  vtN.IsSmall())
      vtN = Z_AX ;
  // assegno la corretta pendenza
   int i = 0 ;
   double dCurrLen = 0 ;
   double dLen ; pCrv->GetLength( dLen) ;
   const ICurve* pSCrv = pCrv->GetFirstCurve() ;
   while ( pSCrv != nullptr) {
      double dCrvLen ;
      pSCrv->GetLength( dCrvLen) ;
      double dCoeff = dCurrLen / dLen ;
      double dCurrN = dNini * ( 1.0 - dCoeff) + dNfin * dCoeff ;
      Point3d ptJoin ;
      pSCrv->GetStartPoint( ptJoin) ;
      pCrv->ModifyJoint( i, ptJoin + vtN * ( dCurrN - ( ptJoin - ORIG) * vtN)) ;
      // passo al successivo
      dCurrLen += dCrvLen ;
      pSCrv = pCrv->GetNextCurve() ;
      ++ i ;
   }
   Point3d ptFin ; pCrv->GetEndPoint( ptFin) ;
   ptFin += vtN * ( dNfin - ( ptFin - ORIG) * vtN) ;
   pCrv->ModifyEnd( ptFin) ;
   return true ;
}