0025929: Make Approx_ComputeLine algorithm adaptive
authoraml <aml@opencascade.com>
Tue, 10 Nov 2015 06:57:03 +0000 (09:57 +0300)
committerbugmaster <bugmaster@opencascade.com>
Fri, 4 Dec 2015 10:12:09 +0000 (13:12 +0300)
Adaptive partition algorithm of WLine is implemented and used in ApproxInt_Approx.gxx file.
Refactoring of ApproxInt_Approx class.
Test cases are updated to the new behaviour.

Filtering algorithm improved.

src/ApproxInt/ApproxInt_Approx.gxx
src/ApproxInt/ApproxInt_KnotTools.cxx [new file with mode: 0644]
src/ApproxInt/ApproxInt_KnotTools.hxx [new file with mode: 0644]
src/ApproxInt/FILES
src/BRepApprox/BRepApprox_Approx.hxx
src/GeomInt/GeomInt_WLApprox.hxx

index 6c40c6d..911e90f 100644 (file)
 #include <gp_Trsf2d.hxx>
 #include <IntSurf_PntOn2S.hxx>
 #include <Precision.hxx>
+#include <ApproxInt_KnotTools.hxx>
 
 const Standard_Integer LimRajout = 5;
-const Standard_Integer NbPntMaxDecoupage = 30 ;
-const Standard_Real RatioTol = 1.5 ;
+const Standard_Integer NbPntMaxDecoupage = 30;
+const Standard_Real RatioTol = 1.5;
 
-static Standard_Real MINABS3(Standard_Real a, Standard_Real b,Standard_Real c) { 
+//=======================================================================
+//function : MINABS3
+//purpose  : Compute minimal absolute distance to 0 from 3 values.
+//=======================================================================
+static Standard_Real MINABS3(Standard_Real a, Standard_Real b,Standard_Real c)
+{
   if(a<0.0) a=-a;
   if(b<0.0) b=-b;
   if(c<0.0) c=-c;
@@ -35,7 +41,12 @@ static Standard_Real MINABS3(Standard_Real a, Standard_Real b,Standard_Real c) {
   return(a);
 }
 
-static Standard_Real MINABS4(Standard_Real a, Standard_Real b,Standard_Real c,Standard_Real d) { 
+//=======================================================================
+//function : MINABS4
+//purpose  : Compute minimal absolute distance to 0 from 4 values.
+//=======================================================================
+static Standard_Real MINABS4(Standard_Real a, Standard_Real b,Standard_Real c,Standard_Real d)
+{
   if(a<0.0) a=-a;
   if(b<0.0) b=-b;
   if(c<0.0) c=-c;
@@ -46,1207 +57,991 @@ static Standard_Real MINABS4(Standard_Real a, Standard_Real b,Standard_Real c,St
   return(a);
 }
 
-static void  ComputeTrsf3d(const Handle(TheWLine)& theline,
-                          Standard_Real& Xo, Standard_Real& Ax,
-                          Standard_Real& Yo, Standard_Real& Ay,
-                          Standard_Real& Zo, Standard_Real& Az) {
-  Standard_Integer nbp = theline->NbPnts();
-  Standard_Real z0,z1,x0,x1,y0,y1;
-  z0=y0=x0=RealLast();
-  z1=y1=x1=RealFirst();
-  for(Standard_Integer i=1;i<=nbp;i++) { 
-    const gp_Pnt& P = theline->Point(i).Value();
-    Standard_Real  X = P.X();
-    Standard_Real  Y = P.Y();
-    Standard_Real  Z = P.Z();
-    if(X<x0) x0=X;
-    if(X>x1) x1=X;
-    if(Y<y0) y0=Y;
-    if(Y>y1) y1=Y;
-    if(Z<z0) z0=Z;
-    if(Z>z1) z1=Z;
-  }
-//-deb-  cout << "ComputeTrsf3d -- NbPnt = " << nbp << endl ;
-//-deb-  cout << "ComputeTrsf3d -- Xm = " << x0 << " Ym = " << y0 << " Zm = " << z0 << endl ;
-//-deb-  cout << "ComputeTrsf3d -- XM = " << x1 << " YM = " << y1 << " ZM = " << z1 << endl ;
-  Standard_Real dx = x1-x0;
-  Standard_Real dy = y1-y0;
-  Standard_Real dz = z1-z0;
-  Standard_Real MaxD = dx; 
-  if(MaxD < dy) MaxD=dy;
-  if(MaxD < dz) MaxD=dz;
-  Standard_Real MaxDF = 0.01*MaxD;
-
-  //-- lbr le 22 fev99 : FPE 
-  if(MaxDF<1e-12) 
-    MaxDF=1.0;
-
-
-  if(dx > MaxDF) { Ax = 1.0 / dx;    Xo = -Ax * x0;  }
-  else {     Ax = 1.0/( MaxDF) ; Xo = -Ax*x0;   }
-  if(dy > MaxDF) { Ay = 1.0 / dy;    Yo = -Ay * y0;  }
-  else {     Ay = 1.0/( MaxDF); Yo = -Ay*y0;   }
-  if(dz > MaxDF) { Az = 1.0 / dz;    Zo = -Az * z0;    }
-  else {     Az = 1.0/(MaxDF); Zo = -Az*z0;   } 
+//=======================================================================
+//function : ComputeTrsf3d
+//purpose  : 
+//=======================================================================
+static void ComputeTrsf3d(const Handle(TheWLine)& theline,
+                          Standard_Real& Xo, Standard_Real& Ax,
+                          Standard_Real& Yo, Standard_Real& Ay,
+                          Standard_Real& Zo, Standard_Real& Az)
+{
+    Standard_Integer nbp = theline->NbPnts();
+    Standard_Real z0,z1,x0,x1,y0,y1;
+    z0=y0=x0=RealLast();
+    z1=y1=x1=RealFirst();
+    for(Standard_Integer i=1;i<=nbp;i++) { 
+      const gp_Pnt& P = theline->Point(i).Value();
+      Standard_Real  X = P.X();
+      Standard_Real  Y = P.Y();
+      Standard_Real  Z = P.Z();
+      if(X<x0) x0=X;
+      if(X>x1) x1=X;
+      if(Y<y0) y0=Y;
+      if(Y>y1) y1=Y;
+      if(Z<z0) z0=Z;
+      if(Z>z1) z1=Z;
+    }
+    Standard_Real dx = x1-x0;
+    Standard_Real dy = y1-y0;
+    Standard_Real dz = z1-z0;
+    Standard_Real MaxD = dx; 
+    if(MaxD < dy) MaxD=dy;
+    if(MaxD < dz) MaxD=dz;
+    Standard_Real MaxDF = 0.01*MaxD;
+
+    if(MaxDF<1e-12) 
+      MaxDF=1.0;
+
+    if(dx > MaxDF) { Ax = 1.0 / dx;    Xo = -Ax * x0;  }
+    else {     Ax = 1.0/( MaxDF) ; Xo = -Ax*x0;   }
+    if(dy > MaxDF) { Ay = 1.0 / dy;    Yo = -Ay * y0;  }
+    else {     Ay = 1.0/( MaxDF); Yo = -Ay*y0;   }
+    if(dz > MaxDF) { Az = 1.0 / dz;    Zo = -Az * z0;  }
+    else {     Az = 1.0/(MaxDF); Zo = -Az*z0;   } 
 }
 
-static void  ComputeTrsf2d(const Handle(TheWLine)& theline,
-                          Standard_Real& Uo, Standard_Real& Au,
-                          Standard_Real& Vo, Standard_Real& Av,
-                          const Standard_Boolean onFirst,
-                          const Standard_Real UVResRatio = 1.) { 
-  Standard_Integer nbp = theline->NbPnts();
-  Standard_Real u0,u1,v0,v1;
-  u0 = v0 = RealLast();
-  u1 = v1 = RealFirst();
-  // pointer to a member-function
-  void (IntSurf_PntOn2S::* pfunc)(Standard_Real&,Standard_Real&) const;
-  if (onFirst)
-    pfunc = &IntSurf_PntOn2S::ParametersOnS1;
-  else
-    pfunc = &IntSurf_PntOn2S::ParametersOnS2;
-  for(Standard_Integer i=1;i<=nbp;i++) { 
-    const IntSurf_PntOn2S&  POn2S = theline->Point(i);
-    Standard_Real  U,V;
-    (POn2S.*pfunc)(U,V);
-    if(U<u0) u0=U;
-    if(U>u1) u1=U;
-    if(V<v0) v0=V;
-    if(V>v1) v1=V;
-  }
+//=======================================================================
+//function : ComputeTrsf2d
+//purpose  :
+//=======================================================================
+static void ComputeTrsf2d(const Handle(TheWLine)& theline,
+                          Standard_Real& Uo, Standard_Real& Au,
+                          Standard_Real& Vo, Standard_Real& Av,
+                          const Standard_Boolean onFirst,
+                          const Standard_Real UVResRatio = 1.0)
+{
+    Standard_Integer nbp = theline->NbPnts();
+    Standard_Real u0,u1,v0,v1;
+    u0 = v0 = RealLast();
+    u1 = v1 = RealFirst();
+    // pointer to a member-function
+    void (IntSurf_PntOn2S::* pfunc)(Standard_Real&,Standard_Real&) const;
+    if (onFirst)
+      pfunc = &IntSurf_PntOn2S::ParametersOnS1;
+    else
+      pfunc = &IntSurf_PntOn2S::ParametersOnS2;
+    for(Standard_Integer i=1;i<=nbp;i++) { 
+      const IntSurf_PntOn2S&  POn2S = theline->Point(i);
+      Standard_Real  U,V;
+      (POn2S.*pfunc)(U,V);
+      if(U<u0) u0=U;
+      if(U>u1) u1=U;
+      if(V<v0) v0=V;
+      if(V>v1) v1=V;
+    }
 
-  Standard_Real du = (u1-u0);
-  Standard_Real dv = (v1-v0);
+    Standard_Real du = (u1-u0);
+    Standard_Real dv = (v1-v0);
 
-  if (UVResRatio > 1.)
-    du *= UVResRatio;
-  else if (UVResRatio < 1.)
-    dv /= UVResRatio;
+    if (UVResRatio > 1.)
+      du *= UVResRatio;
+    else if (UVResRatio < 1.)
+      dv /= UVResRatio;
 
-  Standard_Real MaxUV=du;
-  if(MaxUV<dv) MaxUV=dv;
+    Standard_Real MaxUV=du;
+    if(MaxUV<dv) MaxUV=dv;
 
-  Standard_Real MaxUVF=0.01*MaxUV;
+    Standard_Real MaxUVF=0.01*MaxUV;
 
-  //-- lbr le 22 fev 99 (FPE) 
-  if(MaxUVF<1e-12) 
-    MaxUVF=1.0;
+    //-- lbr le 22 fev 99 (FPE) 
+    if(MaxUVF<1e-12) 
+      MaxUVF=1.0;
 
-  if(du > MaxUVF) { Au = 1.0 / du;    Uo = -Au * u0;  }
-  else {     Au = 1.0/(MaxUVF); Uo = -Au*u0;  }
-  if(dv > MaxUVF) { Av = 1.0 / dv;    Vo = -Av * v0;  }
-  else {     Av = 1.0/(MaxUVF); Vo = -Av*v0;  }
+    if(du > MaxUVF) { Au = 1.0 / du;    Uo = -Au * u0;  }
+    else {     Au = 1.0/(MaxUVF); Uo = -Au*u0;  }
+    if(dv > MaxUVF) { Av = 1.0 / dv;    Vo = -Av * v0;  }
+    else {     Av = 1.0/(MaxUVF); Vo = -Av*v0;  }
 }
 
+//=======================================================================
+//function : Parameters
+//purpose  :
+//=======================================================================
+static void Parameters(const ApproxInt_TheMultiLine& Line,
+                       const Standard_Integer firstP,
+                       const Standard_Integer lastP,
+                       const Approx_ParametrizationType Par,
+                       math_Vector& TheParameters)
+{
+  Standard_Integer i, j, nbP2d, nbP3d;
+  Standard_Real dist;
+  gp_Pnt P1, P2;
+  gp_Pnt2d P12d, P22d;
+
+  if (Par == Approx_ChordLength || Par == Approx_Centripetal) {
+    nbP3d = ApproxInt_TheMultiLineTool::NbP3d(Line);
+    nbP2d = ApproxInt_TheMultiLineTool::NbP2d(Line);
+    Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d;
+    if (nbP3d == 0) mynbP3d = 1;
+    if (nbP2d == 0) mynbP2d = 1;
 
+    TheParameters(firstP) = 0.0;
+    dist = 0.0;
+    TColgp_Array1OfPnt tabP(1, mynbP3d);
+    TColgp_Array1OfPnt tabPP(1, mynbP3d);
+    TColgp_Array1OfPnt2d tabP2d(1, mynbP2d);
+    TColgp_Array1OfPnt2d tabPP2d(1, mynbP2d);
 
-ApproxInt_Approx::ApproxInt_Approx():
-       myComputeLine(4,
-                    8,
-                    0.001,
-                    0.001,
-                    5,
-                    Standard_True),
-       myComputeLineBezier(4,
-                          8,
-                          0.001,
-                          0.001,
-                          5,
-                          Standard_True)
-{ 
+    for (i = firstP+1; i <= lastP; i++) {
+      if (nbP3d != 0 && nbP2d != 0) ApproxInt_TheMultiLineTool::Value(Line, i-1, tabP, tabP2d);
+      else if (nbP2d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i-1, tabP2d);
+      else if (nbP3d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i-1, tabP);
+
+      if (nbP3d != 0 && nbP2d != 0) ApproxInt_TheMultiLineTool::Value(Line, i, tabPP, tabPP2d);
+      else if (nbP2d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i, tabPP2d);
+      else if (nbP3d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i, tabPP);
+      dist = 0;
+      for (j = 1; j <= nbP3d; j++) {
+        P1 = tabP(j);
+        P2 = tabPP(j);
+        dist += P2.Distance(P1);
+      }
+      for (j = 1; j <= nbP2d; j++) {
+        P12d = tabP2d(j);
+        P22d = tabPP2d(j);
+        dist += P22d.Distance(P12d);
+      }
+      if(Par == Approx_ChordLength)
+        TheParameters(i) = TheParameters(i-1) + dist;
+      else {// Par == Approx_Centripetal
+        TheParameters(i) = TheParameters(i-1) + Sqrt(dist);
+      }
+    }
+    for (i = firstP; i <= lastP; i++) TheParameters(i) /= TheParameters(lastP);
+  }
+  else {
+    for (i = firstP; i <= lastP; i++) {
+      TheParameters(i) = (Standard_Real(i)-firstP)/
+        (Standard_Real(lastP)-Standard_Real(firstP));
+    }
+  }
+}
+
+//=======================================================================
+//function : ApproxInt_Approx
+//purpose  : Constructor.
+//=======================================================================
+ApproxInt_Approx::ApproxInt_Approx()
+: myComputeLine(4, 8, 0.001, 0.001, 5, Standard_True),
+  myComputeLineBezier(4, 8, 0.001, 0.001, 5, Standard_True)
+{
   myComputeLine.SetContinuity(2);
-  //-- myComputeLineBezier.SetContinuity(2);
-  myApproxBez = Standard_True;
-  
-  myRelativeTol = Standard_True ;
-  myNbPntMax = NbPntMaxDecoupage ;
-  myMinFactorXYZ = 0.0;
-  myMinFactorUV  = 0.0;
-  myTolReached3d = myTolReached2d = 0.;
+  myData.myBezierApprox = Standard_True;
+
+  myRelativeTol = Standard_True;
+  myNbPntMax = NbPntMaxDecoupage;
+  myData.myMinFactorXYZ = 0.0;
+  myData.myMinFactorUV  = 0.0;
+  myTolReached3d = myTolReached2d = 0.0;
+  myUVRes1 = myUVRes2 = 1.0;
 }
-  
 
+//=======================================================================
+//function : Perform
+//purpose  : Build without surfaces information.
+//=======================================================================
 void ApproxInt_Approx::Perform(const Handle(TheWLine)& theline,
-                              const Standard_Boolean ApproxXYZ,
-                              const Standard_Boolean ApproxU1V1,
-                              const Standard_Boolean ApproxU2V2,
-                              const Standard_Integer indicemin,
-                              const Standard_Integer indicemax) { 
-  
-  myMinFactorXYZ = 0.0;
-  myMinFactorUV  = 0.0;
-  myTolReached3d = myTolReached2d = 0.;
-  
-  
+                               const Standard_Boolean ApproxXYZ,
+                               const Standard_Boolean ApproxU1V1,
+                               const Standard_Boolean ApproxU2V2,
+                               const Standard_Integer indicemin,
+                               const Standard_Integer indicemax)
+{
+  // Prepare DS.
+  prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
+
   Standard_Integer nbpntbez = indicemax-indicemin;
-  Standard_Integer nbpntmax = myNbPntMax;
-  Standard_Boolean OtherInter = Standard_False;
   if(nbpntbez < LimRajout) 
-    myApproxBez = Standard_False;
+    myData.myBezierApprox = Standard_False;
   else 
-    myApproxBez = Standard_True;
-  if(myApproxBez) {
-    myBezToBSpl.Reset();
-    Standard_Integer nbi = (indicemax-indicemin)/nbpntmax;
-    if(nbi>1)  { 
-      nbpntbez = (indicemax-indicemin)/nbi;
-    }
-  }
-  Standard_Integer imin = indicemin;
-  Standard_Integer imax = imin + nbpntbez;
-  myTolReached = Standard_True;
-  
-  Standard_Real Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v;
-  if(ApproxXYZ) { 
-    ComputeTrsf3d(theline,Xo,Ax,Yo,Ay,Zo,Az);
-  }
-  else { 
-    Xo=Yo=Zo=0.0; Ax=Ay=Az=1.0;;
-  }
-  if(ApproxU1V1) { 
-    ComputeTrsf2d(theline,U1o,A1u,V1o,A1v,Standard_True);
-  }
-  else { 
-    U1o=V1o=0.0; A1u=A1v=1.0;
-  }
-  if(ApproxU2V2) { 
-    ComputeTrsf2d(theline,U2o,A2u,V2o,A2v,Standard_False);
-  }
-  else { 
-    U2o=V2o=0.0; A2u=A2v=1.0;
-  }
-  
-  //-deb-  cout << "ApproxInt_Approx -- NbPntMax = " << myNbPntMax    << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Tol3D    = " << myTol3d       << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Tol2D    = " << myTol2d       << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- RelTol   = " << (myRelativeTol ? "RELATIVE" : "ABSOLUTE") << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Xo = " << Xo << " Yo = " << Yo << " Zo = " << Zo << endl ;
-  //-deb-  cout << "ApproxInt_Approx -- Ax = " << Ax << " Ay = " << Ay << " Az = " << Az << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- U1o = " << U1o << " V1o = " << V1o << " A1u = " << A1u << " A1v = " << A1v << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- U2o = " << U2o << " V2o = " << V2o << " A2u = " << A2u << " A2v = " << A2v << endl ; 
-  
-  Standard_Real A3d = MINABS3(Ax,Ay,Az);
-  if((A3d < myMinFactorXYZ) || (myMinFactorXYZ == 0.0)) { 
-    myMinFactorXYZ = A3d;
-  }
-  
-  Standard_Real A2d = MINABS4(A1u,A1v,A2u,A2v);
-  if((A2d < myMinFactorUV) || (myMinFactorUV == 0.0)) { 
-    myMinFactorUV = A2d;
-  }
-  
-  Standard_Boolean cut=Standard_True;
-  Approx_ParametrizationType parametrization;
-  myComputeLineBezier.Parametrization(parametrization);
+    myData.myBezierApprox = Standard_True;
+
+  // Fill data structure.
+  fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
+
+  // Build knots.
+  buildKnots(theline, NULL);
 
-  if(myRelativeTol==Standard_False) { 
-    
+  Standard_Boolean cut = Standard_True;
+  if(myRelativeTol==Standard_False)
+  {
     myComputeLine.Init(myDegMin,
-                      myDegMax,
-                      myTol3d*myMinFactorXYZ,
-                      myTol2d*myMinFactorUV,
-                      myNbIterMax,
-                      cut,
-                      parametrization);
+      myDegMax,
+      myTol3d*myData.myMinFactorXYZ,
+      myTol2d*myData.myMinFactorUV,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
     myComputeLineBezier.Init(myDegMin,
-                            myDegMax,
-                            myTol3d*myMinFactorXYZ,
-                            myTol2d*myMinFactorUV,
-                            myNbIterMax,
-                            cut,
-                            parametrization);
-  }
-  
-  do {
-    ApproxInt_TheMultiLine myMultiLine(theline,
-                                      ((ApproxXYZ)? 1 : 0),
-                                      ((ApproxU1V1)? 1: 0) + ((ApproxU2V2)? 1: 0),
-                                      Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v,
-                                      ApproxU1V1,
-                                      imin,
-                                      imax);
-    
-    if(myApproxBez) { 
-      myComputeLineBezier.Perform(myMultiLine);
-      if (myComputeLineBezier.NbMultiCurves() == 0)
-       return;
-      myTolReached&=myComputeLineBezier.IsToleranceReached();
-    }
-    else { 
-      myComputeLine.Perform(myMultiLine);
-    }
-    UpdateTolReached();
-    
-    Standard_Integer indice3d,indice2d1,indice2d2;
-    indice3d = 1; 
-    indice2d1= 2;
-    indice2d2= 3;
-    if(!ApproxXYZ)  { indice2d1--; indice2d2--; } 
-    if(!ApproxU1V1) { indice2d2--; } 
-    if(ApproxXYZ) { 
-      Standard_Real ax,bx,ay,by,az,bz;
-      ax = 1.0/Ax;   bx = -Xo*ax;
-      ay = 1.0/Ay;   by = -Yo*ay;
-      az = 1.0/Az;   bz = -Zo*az;
-      if(myApproxBez) {
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform(indice3d,bx,ax,by,ay,bz,az);
-       }
-      }
-      else { 
-       myComputeLine.ChangeValue().Transform(indice3d,bx,ax,by,ay,bz,az);
-      }
-    }
-    if(ApproxU1V1) { 
-      Standard_Real ax,bx,ay,by;
-      ax = 1.0/A1u;   bx = -U1o*ax;
-      ay = 1.0/A1v;   by = -V1o*ay;
-      if(myApproxBez) {
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d1,bx,ax,by,ay);
-       }
-      }
-      else { 
-       myComputeLine.ChangeValue().Transform2d(indice2d1,bx,ax,by,ay);
-      }
-    }
-    if(ApproxU2V2) { 
-      Standard_Real ax,bx,ay,by;
-      ax = 1.0/A2u;   bx = -U2o*ax;
-      ay = 1.0/A2v;   by = -V2o*ay;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d2,bx,ax,by,ay);
-       }
-      }
-      else { 
-       myComputeLine.ChangeValue().Transform2d(indice2d2,bx,ax,by,ay);
-      }
-    }
-    
-    OtherInter = Standard_False;
-    if(myApproxBez) {
-      for(Standard_Integer nbmc = 1; 
-         nbmc <= myComputeLineBezier.NbMultiCurves() ;
-         nbmc++) { 
-       myBezToBSpl.Append(myComputeLineBezier.Value(nbmc));
-      }
-      if(imax<indicemax)
-      { 
-        imin = imax;    
-        imax = imin+nbpntbez;
-        OtherInter = Standard_True;
-        if((indicemax-imax)<(nbpntbez/2))
-        {
-          imax = indicemax;
-        }
-        imax = CorrectFinishIdx(imin, imax, theline);
-      }
-    }
-  }
-  while(OtherInter);
-  if(myApproxBez) { 
-    myBezToBSpl.Perform();
+      myDegMax,
+      myTol3d*myData.myMinFactorXYZ,
+      myTol2d*myData.myMinFactorUV,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
   }
+
+  buildCurve(theline, NULL);
 }
 
+//=======================================================================
+//function : Perform
+//purpose  : Param-Param perform.
+//=======================================================================
 void ApproxInt_Approx::Perform(const ThePSurface& Surf1,
-                              const ThePSurface& Surf2,
-                              const Handle(TheWLine)& theline,
-                              const Standard_Boolean ApproxXYZ,
-                              const Standard_Boolean ApproxU1V1,
-                              const Standard_Boolean ApproxU2V2,
-                              const Standard_Integer indicemin,
-                              const Standard_Integer indicemax) {
-  myMinFactorXYZ = 0.0;
-  myMinFactorUV  = 0.0;
-  myTolReached3d = myTolReached2d = 0.;
+                               const ThePSurface& Surf2,
+                               const Handle(TheWLine)& theline,
+                               const Standard_Boolean ApproxXYZ,
+                               const Standard_Boolean ApproxU1V1,
+                               const Standard_Boolean ApproxU2V2,
+                               const Standard_Integer indicemin,
+                               const Standard_Integer indicemax) 
+{
 
   GeomAbs_SurfaceType typeS1 = ThePSurfaceTool::GetType(Surf1);
   GeomAbs_SurfaceType typeS2 = ThePSurfaceTool::GetType(Surf2);
   if ((typeS1 != GeomAbs_Plane    &&
        typeS1 != GeomAbs_Cylinder &&
        typeS1 != GeomAbs_Sphere   &&
-       typeS1 != GeomAbs_Cone) 
-      &&
+       typeS1 != GeomAbs_Cone)    &&
       (typeS2 != GeomAbs_Plane    &&
        typeS2 != GeomAbs_Cylinder &&
        typeS2 != GeomAbs_Sphere   &&
-       typeS2 != GeomAbs_Cone)) { 
-    
-    //------------------------------------------------------------
-    //-- Construction du SvSurfaces
-    //------------------------------------------------------------
+       typeS2 != GeomAbs_Cone))
+  {
+    // Prepare DS.
+    prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
+
+    // Non-analytical case: Param-Param perform.
     ApproxInt_ThePrmPrmSvSurfaces myPrmPrmSvSurfaces(Surf1,Surf2);
-    //------------------------------------------------------------
-    //-- Construction de la MultiLine
-    //------------------------------------------------------------
+
     Standard_Integer nbpntbez = indicemax-indicemin;
-    Standard_Integer nbpntmax = myNbPntMax;
-    Standard_Boolean OtherInter = Standard_False;
-    
     if(nbpntbez < LimRajout) 
-      myApproxBez = Standard_False;
+      myData.myBezierApprox = Standard_False;
     else 
-      myApproxBez = Standard_True;
-    
-    Standard_Address ptrsvsurf = NULL;
+      myData.myBezierApprox = Standard_True;
+
     Standard_Boolean cut = Standard_True;
-    if(nbpntbez < LimRajout) {   
+    if(nbpntbez < LimRajout)
+    {
       cut = Standard_False;
-      //-- cout<<" ApproxInt : Nb de points = "<<nbpntbez<<" Pas de rajout "<<endl;
-    }
-    ptrsvsurf = &myPrmPrmSvSurfaces;
-    
-
-    if(myApproxBez) { 
-      myBezToBSpl.Reset();
-      Standard_Integer nbi = (indicemax-indicemin)/nbpntmax;
-      if(nbi>1)  { 
-       nbpntbez = (indicemax-indicemin)/nbi;
-      }
     }
-    Standard_Integer imin = indicemin;
-    Standard_Integer imax = imin + nbpntbez;
-    myTolReached = Standard_True;
-
 
-    Standard_Real Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v;
-    if(ApproxXYZ) { 
-      ComputeTrsf3d(theline,Xo,Ax,Yo,Ay,Zo,Az);
-    }
-    else { 
-      Xo=Yo=Zo=0.0; Ax=Ay=Az=1.0;;
-    }
-    if(ApproxU1V1) { 
-      Standard_Real UVResRatio = ThePSurfaceTool::UResolution(Surf1,1.)/
-                                 ThePSurfaceTool::VResolution(Surf1,1.);
-      ComputeTrsf2d(theline,U1o,A1u,V1o,A1v,Standard_True,UVResRatio);
-    }
-    else { 
-      U1o=V1o=0.0; A1u=A1v=1.0;
-    }      
-    if(ApproxU2V2) { 
-      Standard_Real UVResRatio = ThePSurfaceTool::UResolution(Surf2,1.)/
-                                 ThePSurfaceTool::VResolution(Surf2,1.);
-      ComputeTrsf2d(theline,U2o,A2u,V2o,A2v,Standard_False,UVResRatio);
-    }
-    else { 
-      U2o=V2o=0.0; A2u=A2v=1.0;
-    }
+    Standard_Real aS1URes = ThePSurfaceTool::UResolution(Surf1, 1.0),
+                  aS1VRes = ThePSurfaceTool::VResolution(Surf1, 1.0),
+                  aS2URes = ThePSurfaceTool::UResolution(Surf2, 1.0),
+                  aS2VRes = ThePSurfaceTool::VResolution(Surf2, 1.0);
+    if(ApproxU1V1)
+      myUVRes1 = aS1URes / aS1VRes;
+    if(ApproxU2V2)
+      myUVRes2 = aS2URes / aS2VRes;
 
-//-deb-  cout << "ApproxInt_Approx -- NbPntMax = " << myNbPntMax    << endl ; 
-//-deb-  cout << "ApproxInt_Approx -- Tol3D    = " << myTol3d       << endl ; 
-//-deb-  cout << "ApproxInt_Approx -- Tol2D    = " << myTol2d       << endl ; 
-//-deb-  cout << "ApproxInt_Approx -- RelTol   = " << (myRelativeTol ? "RELATIVE" : "ABSOLUTE") << endl ; 
-//-deb-  cout << "ApproxInt_Approx -- Xo = " << Xo << " Yo = " << Yo << " Zo = " << Zo << endl ;
-//-deb-  cout << "ApproxInt_Approx -- Ax = " << Ax << " Ay = " << Ay << " Az = " << Az << endl ; 
-//-deb-  cout << "ApproxInt_Approx -- U1o = " << U1o << " V1o = " << V1o << " A1u = " << A1u << " A1v = " << A1v << endl ; 
-//-deb-  cout << "ApproxInt_Approx -- U2o = " << U2o << " V2o = " << V2o << " A2u = " << A2u << " A2v = " << A2v << endl ; 
-    
-    
-    Standard_Real A3d = MINABS3(Ax,Ay,Az);
-    if((A3d < myMinFactorXYZ) || (myMinFactorXYZ == 0.0)) { 
-      myMinFactorXYZ = A3d;
-    }
-    
-    Standard_Real A2d = MINABS4(A1u,A1v,A2u,A2v);
-    if((A2d < myMinFactorUV) || (myMinFactorUV == 0.0)) { 
-      myMinFactorUV = A2d;
-    }
-    
+    // Fill data structure.
+    fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
 
-    Approx_ParametrizationType parametrization;
-    myComputeLineBezier.Parametrization(parametrization);
+    // Build knots.
+    Standard_Address ptrsvsurf = &myPrmPrmSvSurfaces;
+    buildKnots(theline, ptrsvsurf);
 
-    if(myRelativeTol==Standard_False) { 
+    if(myRelativeTol==Standard_False)
+    {
       myComputeLine.Init(myDegMin,
-                        myDegMax,
-                        myTol3d*myMinFactorXYZ,
-                        myTol2d*myMinFactorUV,
-                        myNbIterMax,
-                        cut,
-                        parametrization);
+        myDegMax,
+        myTol3d*myData.myMinFactorXYZ,
+        myTol2d*myData.myMinFactorUV,
+        myNbIterMax,
+        cut,
+        myData.parametrization);
       myComputeLineBezier.Init(myDegMin,
-                              myDegMax,
-                              myTol3d*myMinFactorXYZ,
-                              myTol2d*myMinFactorUV,
-                              myNbIterMax,
-                              cut,
-                              parametrization);
+        myDegMax,
+        myTol3d*myData.myMinFactorXYZ,
+        myTol2d*myData.myMinFactorUV,
+        myNbIterMax,
+        cut,
+        myData.parametrization);
     }
-    else { 
+    else
+    {
       myComputeLine.Init(myDegMin,
-                        myDegMax,
-                        myTol3d,
-                        myTol2d,
-                        myNbIterMax,
-                        cut,
-                        parametrization);
+        myDegMax,
+        myTol3d,
+        myTol2d,
+        myNbIterMax,
+        cut,
+        myData.parametrization);
       myComputeLineBezier.Init(myDegMin,
-                              myDegMax,
-                              myTol3d,
-                              myTol2d,
-                              myNbIterMax,
-                              cut,
-                              parametrization);
-    }     
-
-
-
-    
-    do {
-      ApproxInt_TheMultiLine myMultiLine(theline,
-                                        ptrsvsurf,
-                                        ((ApproxXYZ)? 1 : 0),
-                                        ((ApproxU1V1)? 1: 0) + ((ApproxU2V2)? 1: 0),
-                                        Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v,
-                                        ApproxU1V1,
-                                        imin,
-                                        imax);
-      
-      if(myApproxBez) { 
-       myComputeLineBezier.Perform(myMultiLine);
-       if (myComputeLineBezier.NbMultiCurves() == 0)
-         return;
-       myTolReached&=myComputeLineBezier.IsToleranceReached();
-      }
-      else { 
-       myComputeLine.Perform(myMultiLine);
-      }
-      UpdateTolReached();
-      
-      Standard_Integer indice3d,indice2d1,indice2d2;
-      indice3d = 1; 
-      indice2d1= 2;
-      indice2d2= 3;
-      if(!ApproxXYZ)  { indice2d1--; indice2d2--; } 
-      if(!ApproxU1V1) { indice2d2--; } 
-      if(ApproxXYZ) { 
-       Standard_Real ax,bx,ay,by,az,bz;
-       ax = 1.0/Ax;   bx = -Xo*ax;
-       ay = 1.0/Ay;   by = -Yo*ay;
-       az = 1.0/Az;   bz = -Zo*az;
-       if(myApproxBez) { 
-         for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-           myComputeLineBezier.ChangeValue(nbmc).Transform(indice3d,bx,ax,by,ay,bz,az);
-         }
-       }
-       else { 
-         myComputeLine.ChangeValue().Transform(indice3d,bx,ax,by,ay,bz,az);
-       }
-      }
-      if(ApproxU1V1) { 
-       Standard_Real ax,bx,ay,by;
-       ax = 1.0/A1u;   bx = -U1o*ax;
-       ay = 1.0/A1v;   by = -V1o*ay;
-       if(myApproxBez) { 
-         for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-           myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d1,bx,ax,by,ay);
-         }
-       }
-       else { 
-         myComputeLine.ChangeValue().Transform2d(indice2d1,bx,ax,by,ay);
-       }
-      }
-      if(ApproxU2V2) { 
-       Standard_Real ax,bx,ay,by;
-       ax = 1.0/A2u;   bx = -U2o*ax;
-       ay = 1.0/A2v;   by = -V2o*ay;
-       if(myApproxBez) { 
-         for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-           myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d2,bx,ax,by,ay);
-         }
-       }
-       else { 
-         myComputeLine.ChangeValue().Transform2d(indice2d2,bx,ax,by,ay);
-       }
-      }
-      OtherInter = Standard_False;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = 1; 
-         nbmc <= myComputeLineBezier.NbMultiCurves() ;
-         nbmc++) { 
-         myBezToBSpl.Append(myComputeLineBezier.Value(nbmc));
-       }
-        if(imax<indicemax)
-        {
-          imin = imax;
-          imax = imin+nbpntbez;
-          OtherInter = Standard_True;
-          if((indicemax-imax)<(nbpntbez/2))
-          {
-            imax = indicemax;
-          }
-          imax = CorrectFinishIdx(imin, imax, theline);
-        }
-      }
-    }
-    while(OtherInter);
-    if(myApproxBez) { 
-      myBezToBSpl.Perform();
+        myDegMax,
+        myTol3d,
+        myTol2d,
+        myNbIterMax,
+        cut,
+        myData.parametrization);
     }
 
+    buildCurve(theline, ptrsvsurf);
   }
-  else { 
+  else
+  {
     IntSurf_Quadric Quad;
     Standard_Boolean SecondIsImplicit=Standard_False;
-    switch (typeS1) {
-      
+    switch (typeS1)
+    {
+
     case GeomAbs_Plane:
       Quad.SetValue(ThePSurfaceTool::Plane(Surf1));
       break;
-      
+
     case GeomAbs_Cylinder:
       Quad.SetValue(ThePSurfaceTool::Cylinder(Surf1));
       break;
-      
+
     case GeomAbs_Sphere:
       Quad.SetValue(ThePSurfaceTool::Sphere(Surf1));
       break;
-      
+
     case GeomAbs_Cone:
       Quad.SetValue(ThePSurfaceTool::Cone(Surf1));
       break;
 
     default:
       {
-       SecondIsImplicit = Standard_True;
-       switch (typeS2) {
-       case GeomAbs_Plane:
-         Quad.SetValue(ThePSurfaceTool::Plane(Surf2));
-         break;
-         
-       case GeomAbs_Cylinder:
-         Quad.SetValue(ThePSurfaceTool::Cylinder(Surf2));
-         break;
-         
-       case GeomAbs_Sphere:
-         Quad.SetValue(ThePSurfaceTool::Sphere(Surf2));
-         break;
-         
-       case GeomAbs_Cone:
-         Quad.SetValue(ThePSurfaceTool::Cone(Surf2));
-         break;
-         
-       default:
-         break;
-       }
+        SecondIsImplicit = Standard_True;
+        switch (typeS2)
+        {
+        case GeomAbs_Plane:
+          Quad.SetValue(ThePSurfaceTool::Plane(Surf2));
+          break;
+
+        case GeomAbs_Cylinder:
+          Quad.SetValue(ThePSurfaceTool::Cylinder(Surf2));
+          break;
+
+        case GeomAbs_Sphere:
+          Quad.SetValue(ThePSurfaceTool::Sphere(Surf2));
+          break;
+
+        case GeomAbs_Cone:
+          Quad.SetValue(ThePSurfaceTool::Cone(Surf2));
+          break;
+
+        default:
+          break;
+        }
       }
       break;
     } 
-    if(SecondIsImplicit) {
+    if(SecondIsImplicit)
+    {
       Perform(Surf1,Quad,theline,ApproxXYZ,ApproxU1V1,ApproxU2V2,indicemin,indicemax);
     }
-    else {
+    else
+    {
       Perform(Quad,Surf2,theline,ApproxXYZ,ApproxU1V1,ApproxU2V2,indicemin,indicemax);
     }
   }
 }
-//--------------------------------------------------------------------------------
+
+//=======================================================================
+//function : SetParameters
+//purpose  :
+//=======================================================================
 void ApproxInt_Approx::SetParameters(const Standard_Real     Tol3d,
-                                    const Standard_Real     Tol2d,
-                                    const Standard_Integer  DegMin,
-                                    const Standard_Integer  DegMax,
-                                    const Standard_Integer  NbIterMax,
-                                    const Standard_Boolean  ApproxWithTangency,
-                                    const Approx_ParametrizationType Parametrization) {
-  myWithTangency = ApproxWithTangency;
-  myTol3d        = Tol3d / RatioTol;
-  myTol2d        = Tol2d / RatioTol;
-  myDegMin       = DegMin;
-  myDegMax       = DegMax;
-  myNbIterMax    = NbIterMax;
-  myComputeLine.Init(myDegMin,
-                    myDegMax,
-                    myTol3d,
-                    myTol2d,
-                    myNbIterMax,
-                    Standard_True,
-                    Parametrization);
-
-  if(!ApproxWithTangency) { 
-    myComputeLine.SetConstraints(AppParCurves_PassPoint,AppParCurves_PassPoint);
-  }
-  myComputeLineBezier.Init(myDegMin,
-                          myDegMax,
-                          myTol3d,
-                          myTol2d,
-                          myNbIterMax,
-                          Standard_True,
-                          Parametrization);
-  if(!ApproxWithTangency) { 
-    myComputeLineBezier.SetConstraints(AppParCurves_PassPoint,AppParCurves_PassPoint);
-  }
-  myApproxBez = Standard_True;
+                                     const Standard_Real     Tol2d,
+                                     const Standard_Integer  DegMin,
+                                     const Standard_Integer  DegMax,
+                                     const Standard_Integer  NbIterMax,
+                                     const Standard_Boolean  ApproxWithTangency,
+                                     const Approx_ParametrizationType Parametrization)
+{
+    myWithTangency = ApproxWithTangency;
+    myTol3d        = Tol3d / RatioTol;
+    myTol2d        = Tol2d / RatioTol;
+    myDegMin       = DegMin;
+    myDegMax       = DegMax;
+    myNbIterMax    = NbIterMax;
+    myComputeLine.Init(myDegMin,
+      myDegMax,
+      myTol3d,
+      myTol2d,
+      myNbIterMax,
+      Standard_True,
+      Parametrization);
+
+    if(!ApproxWithTangency) { 
+      myComputeLine.SetConstraints(AppParCurves_PassPoint,AppParCurves_PassPoint);
+    }
+    myComputeLineBezier.Init(myDegMin,
+      myDegMax,
+      myTol3d,
+      myTol2d,
+      myNbIterMax,
+      Standard_True,
+      Parametrization);
+    if(!ApproxWithTangency)
+    {
+      myComputeLineBezier.SetConstraints(AppParCurves_PassPoint,AppParCurves_PassPoint);
+    }
+    myData.myBezierApprox = Standard_True;
 }
 
+//=======================================================================
+//function : SetParameters
+//purpose  : Set parameters with RelativeTol flag and NbPntMax value.
+//=======================================================================
 void ApproxInt_Approx::SetParameters(const Standard_Real     Tol3d,
-                                    const Standard_Real     Tol2d,
-                                    const Standard_Boolean  RelativeTol,
-                                    const Standard_Integer  DegMin,
-                                    const Standard_Integer  DegMax,
-                                    const Standard_Integer  NbIterMax,
-                                    const Standard_Integer  NbPntMax,
-                                    const Standard_Boolean  ApproxWithTangency,
-                                    const Approx_ParametrizationType Parametrization) 
+                                     const Standard_Real     Tol2d,
+                                     const Standard_Boolean  RelativeTol,
+                                     const Standard_Integer  DegMin,
+                                     const Standard_Integer  DegMax,
+                                     const Standard_Integer  NbIterMax,
+                                     const Standard_Integer  NbPntMax,
+                                     const Standard_Boolean  ApproxWithTangency,
+                                     const Approx_ParametrizationType Parametrization) 
 {
   myNbPntMax = NbPntMax ;
   myRelativeTol = RelativeTol ;
   SetParameters (Tol3d, Tol2d, DegMin, DegMax, NbIterMax, ApproxWithTangency, Parametrization) ;
 }
 
-//--------------------------------------------------------------------------------
+//=======================================================================
+//function : Perform
+//purpose  : Param-Analytic perform.
+//=======================================================================
 void ApproxInt_Approx::Perform(const ThePSurface& PSurf,
-                              const TheISurface& ISurf,
-                              const Handle(TheWLine)& theline,
-                              const Standard_Boolean ApproxXYZ,
-                              const Standard_Boolean ApproxU1V1,
-                              const Standard_Boolean ApproxU2V2,
-                              const Standard_Integer indicemin,
-                              const Standard_Integer indicemax)
+                               const TheISurface& ISurf,
+                               const Handle(TheWLine)& theline,
+                               const Standard_Boolean ApproxXYZ,
+                               const Standard_Boolean ApproxU1V1,
+                               const Standard_Boolean ApproxU2V2,
+                               const Standard_Integer indicemin,
+                               const Standard_Integer indicemax)
 {
-  myMinFactorXYZ = 0.0;
-  myMinFactorUV  = 0.0;
-  myTolReached3d = myTolReached2d = 0.;
-  
+  // Prepare DS.
+  prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
+
+  // Non-analytical case: Param-Analytic perform.
   ApproxInt_TheImpPrmSvSurfaces myImpPrmSvSurfaces(PSurf,ISurf);
+
   Standard_Integer nbpntbez = indicemax-indicemin;
-  Standard_Integer nbpntmax = myNbPntMax;
-  Standard_Boolean OtherInter = Standard_False;
   if(nbpntbez < LimRajout) 
-    myApproxBez = Standard_False;
+    myData.myBezierApprox = Standard_False;
   else 
-    myApproxBez = Standard_True;
-  
-  Standard_Address ptrsvsurf = NULL;
+    myData.myBezierApprox = Standard_True;
+
   Standard_Boolean cut = Standard_True;
-  if(nbpntbez < LimRajout) {   
+  if(nbpntbez < LimRajout)
+  {
     cut = Standard_False;
-    //-- cout<<" ApproxInt : Nb de points = "<<nbpntbez<<" Pas de rajout "<<endl;
   }
 
+  Standard_Real aS1URes = ThePSurfaceTool::UResolution(PSurf, 1.0),
+                aS1VRes = ThePSurfaceTool::VResolution(PSurf, 1.0);
+  if(ApproxU1V1)
+    myUVRes1 = aS1URes / aS1VRes;
 
-  Approx_ParametrizationType parametrization;
-  myComputeLineBezier.Parametrization(parametrization);
-
+  // Fill data structure.
+  fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
 
-  ptrsvsurf = &myImpPrmSvSurfaces;  
-  myComputeLine.Init(myDegMin,
-                    myDegMax,
-                    myTol3d,
-                    myTol2d,
-                    myNbIterMax,
-                    cut,
-                    parametrization);
-
-  myComputeLineBezier.Init(myDegMin,
-                          myDegMax,
-                          myTol3d,
-                          myTol2d,
-                          myNbIterMax,
-                          cut,
-                          parametrization);
-  if(myApproxBez) { 
-    myBezToBSpl.Reset();
-    Standard_Integer nbi = (indicemax-indicemin)/nbpntmax;
-    if(nbi>1)  { 
-      nbpntbez = (indicemax-indicemin)/nbi;
-    }
-  }
-  Standard_Integer imin = indicemin;
-  Standard_Integer imax = imin + nbpntbez;
-  myTolReached = Standard_True;
-  Standard_Real Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v;
-  if(ApproxXYZ) { 
-    ComputeTrsf3d(theline,Xo,Ax,Yo,Ay,Zo,Az);
-  }
-  else { 
-    Xo=Yo=Zo=0.0; Ax=Ay=Az=1.0;;
-  }
-  if(ApproxU1V1) { 
-    Standard_Real UVResRatio = ThePSurfaceTool::UResolution(PSurf,1.)/
-                               ThePSurfaceTool::VResolution(PSurf,1.);
-    ComputeTrsf2d(theline,U1o,A1u,V1o,A1v,Standard_True,UVResRatio);
-  }
-  else { 
-    U1o=V1o=0.0; A1u=A1v=1.0;
-  }
-  if(ApproxU2V2) { 
-    ComputeTrsf2d(theline,U2o,A2u,V2o,A2v,Standard_False);
-  }
-  else { 
-    U2o=V2o=0.0; A2u=A2v=1.0;
-  }
-  
-  //-deb-  cout << "ApproxInt_Approx -- NbPntMax = " << myNbPntMax    << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Tol3D    = " << myTol3d       << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Tol2D    = " << myTol2d       << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- RelTol   = " << (myRelativeTol ? "RELATIVE" : "ABSOLUTE") << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Xo = " << Xo << " Yo = " << Yo << " Zo = " << Zo << endl ;
-  //-deb-  cout << "ApproxInt_Approx -- Ax = " << Ax << " Ay = " << Ay << " Az = " << Az << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- U1o = " << U1o << " V1o = " << V1o << " A1u = " << A1u << " A1v = " << A1v << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- U2o = " << U2o << " V2o = " << V2o << " A2u = " << A2u << " A2v = " << A2v << endl ; 
-  
-  
-  Standard_Real A3d = MINABS3(Ax,Ay,Az);
-  if((A3d < myMinFactorXYZ) || (myMinFactorXYZ == 0.0)) { 
-    myMinFactorXYZ = A3d;
-  }
-  
-  Standard_Real A2d = MINABS4(A1u,A1v,A2u,A2v);
-  if((A2d < myMinFactorUV) || (myMinFactorUV == 0.0)) { 
-    myMinFactorUV = A2d;
-  }
-  
-  myComputeLineBezier.Parametrization(parametrization);
+  // Build knots.
+  Standard_Address ptrsvsurf = &myImpPrmSvSurfaces;
+  buildKnots(theline, ptrsvsurf);
 
-  if(myRelativeTol==Standard_False) { 
+  if(myRelativeTol==Standard_False)
+  {
     myComputeLine.Init(myDegMin,
-                      myDegMax,
-                      myTol3d*myMinFactorXYZ,
-                      myTol2d*myMinFactorUV,
-                      myNbIterMax,
-                      cut,
-                      parametrization);
+      myDegMax,
+      myTol3d*myData.myMinFactorXYZ,
+      myTol2d*myData.myMinFactorUV,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
     myComputeLineBezier.Init(myDegMin,
-                            myDegMax,
-                            myTol3d*myMinFactorXYZ,
-                            myTol2d*myMinFactorUV,
-                            myNbIterMax,
-                            cut,
-                            parametrization);
+      myDegMax,
+      myTol3d*myData.myMinFactorXYZ,
+      myTol2d*myData.myMinFactorUV,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
   }
-  else { 
+  else
+  {
     myComputeLine.Init(myDegMin,
-                      myDegMax,
-                      myTol3d,
-                      myTol2d,
-                      myNbIterMax,
-                      cut,
-                      parametrization);
+      myDegMax,
+      myTol3d,
+      myTol2d,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
     myComputeLineBezier.Init(myDegMin,
-                            myDegMax,
-                            myTol3d,
-                            myTol2d,
-                            myNbIterMax,
-                            cut,
-                            parametrization);
-  }     
-  
-  
-  do {
-    
-    ApproxInt_TheMultiLine myMultiLine(theline,
-                                      ptrsvsurf,
-                                      ((ApproxXYZ)? 1 : 0),
-                                      ((ApproxU1V1)? 1: 0) + ((ApproxU2V2)? 1: 0),
-                                      Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v,
-                                      ApproxU1V1,
-                                      imin,
-                                      imax);
-    if(myApproxBez) { 
-      myComputeLineBezier.Perform(myMultiLine);
-      if (myComputeLineBezier.NbMultiCurves() == 0)
-       return;
-      myTolReached&=myComputeLineBezier.IsToleranceReached();
-    }
-    else { 
-      myComputeLine.Perform(myMultiLine);
-    }
-    UpdateTolReached();
-    Standard_Integer indice3d,indice2d1,indice2d2;
-    indice3d = 1; 
-    indice2d1= 2;
-    indice2d2= 3;
-    if(!ApproxXYZ)  { indice2d1--; indice2d2--; } 
-    if(!ApproxU1V1) { indice2d2--; } 
-    if(ApproxXYZ) { 
-      Standard_Real ax,bx,ay,by,az,bz;
-      ax = 1.0/Ax;   bx = -Xo*ax;
-      ay = 1.0/Ay;   by = -Yo*ay;
-      az = 1.0/Az;   bz = -Zo*az;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform(indice3d,bx,ax,by,ay,bz,az);
-       }
-      }
-      else { 
-       myComputeLine.ChangeValue().Transform(indice3d,bx,ax,by,ay,bz,az);
-      }
-    }
-    if(ApproxU1V1) { 
-      Standard_Real ax,bx,ay,by;
-      ax = 1.0/A1u;   bx = -U1o*ax;
-      ay = 1.0/A1v;   by = -V1o*ay;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d1,bx,ax,by,ay);
-       }
-      }
-      else { 
-       myComputeLine.ChangeValue().Transform2d(indice2d1,bx,ax,by,ay);
-      }
-    }
-    if(ApproxU2V2) { 
-      Standard_Real ax,bx,ay,by;
-      ax = 1.0/A2u;   bx = -U2o*ax;
-      ay = 1.0/A2v;   by = -V2o*ay;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d2,bx,ax,by,ay);
-       }
-      }
-      else { 
-       myComputeLine.ChangeValue().Transform2d(indice2d2,bx,ax,by,ay);
-      }
-    }
-    OtherInter = Standard_False;
-    if(myApproxBez) { 
-      for(Standard_Integer nbmc = 1; 
-         nbmc <= myComputeLineBezier.NbMultiCurves() ;
-         nbmc++) { 
-       myBezToBSpl.Append(myComputeLineBezier.Value(nbmc));
-      }
-      if(imax<indicemax)
-      {
-        imin = imax;
-        imax = imin+nbpntbez;
-        OtherInter = Standard_True;
-        if((indicemax-imax)<(nbpntbez/2))
-        {
-          imax = indicemax;
-        }
-        imax = CorrectFinishIdx(imin, imax, theline);
-      }
-    }
-  }
-  while(OtherInter);
-  if(myApproxBez) { 
-    myBezToBSpl.Perform();
+      myDegMax,
+      myTol3d,
+      myTol2d,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
   }
+
+  buildCurve(theline, ptrsvsurf);
 }
-//--------------------------------------------------------------------------------
+
+//=======================================================================
+//function : Perform
+//purpose  : Analytic-Param perform.
+//=======================================================================
 void ApproxInt_Approx::Perform(const TheISurface& ISurf,
-                              const ThePSurface& PSurf,
-                              const Handle(TheWLine)& theline,
-                              const Standard_Boolean ApproxXYZ,
-                              const Standard_Boolean ApproxU1V1,
-                              const Standard_Boolean ApproxU2V2,
-                              const Standard_Integer indicemin,
-                              const Standard_Integer indicemax)
+                               const ThePSurface& PSurf,
+                               const Handle(TheWLine)& theline,
+                               const Standard_Boolean ApproxXYZ,
+                               const Standard_Boolean ApproxU1V1,
+                               const Standard_Boolean ApproxU2V2,
+                               const Standard_Integer indicemin,
+                               const Standard_Integer indicemax)
 {
+  // Prepare DS.
+  prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
+
+  // Non-analytical case: Analytic-Param perform.
+  ApproxInt_TheImpPrmSvSurfaces myImpPrmSvSurfaces(ISurf, PSurf);
 
-  myMinFactorXYZ = 0.0;
-  myMinFactorUV  = 0.0;
-  myTolReached3d = myTolReached2d = 0.;
-  
-  ApproxInt_TheImpPrmSvSurfaces myImpPrmSvSurfaces(ISurf,PSurf);
   Standard_Integer nbpntbez = indicemax-indicemin;
-  
-  Standard_Boolean cut = Standard_True;
   if(nbpntbez < LimRajout) 
-    myApproxBez = Standard_False;
+    myData.myBezierApprox = Standard_False;
   else 
-    myApproxBez = Standard_True;
-  
-  if(nbpntbez < LimRajout) {   
+    myData.myBezierApprox = Standard_True;
+
+  Standard_Boolean cut = Standard_True;
+  if(nbpntbez < LimRajout)
+  {
     cut = Standard_False;
-    //-- cout<<" ApproxInt : Nb de points = "<<nbpntbez<<" Pas de rajout "<<endl;
   }
+
+  Standard_Real aS2URes = ThePSurfaceTool::UResolution(PSurf, 1.0),
+                aS2VRes = ThePSurfaceTool::VResolution(PSurf, 1.0);
+  if(ApproxU2V2)
+    myUVRes2 = aS2URes / aS2VRes;
+
+  // Fill data structure.
+  fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
+
+  // Build knots.
   Standard_Address ptrsvsurf = &myImpPrmSvSurfaces;
-  
-  if(nbpntbez < LimRajout) myApproxBez = Standard_False;
-  
-  Standard_Integer nbpntmax = myNbPntMax;
-  Standard_Boolean OtherInter = Standard_False;
-  if(myApproxBez) { 
-    myBezToBSpl.Reset();
-    Standard_Integer nbi = (indicemax-indicemin)/nbpntmax;
-    if(nbi>1)  { 
-      nbpntbez = (indicemax-indicemin)/nbi;
+  buildKnots(theline, ptrsvsurf);
+
+  if(myRelativeTol==Standard_False)
+  {
+    myComputeLine.Init(myDegMin,
+      myDegMax,
+      myTol3d*myData.myMinFactorXYZ,
+      myTol2d*myData.myMinFactorUV,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
+    myComputeLineBezier.Init(myDegMin,
+      myDegMax,
+      myTol3d*myData.myMinFactorXYZ,
+      myTol2d*myData.myMinFactorUV,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
+  }
+  else
+  {
+    myComputeLine.Init(myDegMin,
+      myDegMax,
+      myTol3d,
+      myTol2d,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
+    myComputeLineBezier.Init(myDegMin,
+      myDegMax,
+      myTol3d,
+      myTol2d,
+      myNbIterMax,
+      cut,
+      myData.parametrization);
+  }
+
+  buildCurve(theline, ptrsvsurf);
+}
+
+//=======================================================================
+//function : NbMultiCurves
+//purpose  :
+//=======================================================================
+Standard_Integer ApproxInt_Approx::NbMultiCurves() const
+{
+  return 1;
+}
+
+//=======================================================================
+//function : UpdateTolReached
+//purpose  :
+//=======================================================================
+void ApproxInt_Approx::UpdateTolReached()
+{
+  if (myData.myBezierApprox)
+  {
+    Standard_Integer ICur;
+    Standard_Integer NbCurves = myComputeLineBezier.NbMultiCurves();
+    for (ICur = 1 ; ICur <= NbCurves ; ICur++)
+    {
+      Standard_Real Tol3D, Tol2D ;
+      myComputeLineBezier.Error (ICur, Tol3D, Tol2D) ;
+      myTolReached3d = Max(myTolReached3d, Tol3D);
+      myTolReached2d = Max(myTolReached2d, Tol2D);
     }
   }
-  Standard_Integer imin = indicemin;
-  Standard_Integer imax = imin + nbpntbez;
-  myTolReached = Standard_True;
-  
-  Standard_Real Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v;
-  if(ApproxXYZ) { 
-    ComputeTrsf3d(theline,Xo,Ax,Yo,Ay,Zo,Az);
+  else
+  {
+    myComputeLine.Error (myTolReached3d, myTolReached2d);
   }
-  else { 
-    Xo=Yo=Zo=0.0; Ax=Ay=Az=1.0;;
+}
+
+//=======================================================================
+//function : TolReached3d
+//purpose  :
+//=======================================================================
+Standard_Real ApproxInt_Approx::TolReached3d() const
+{
+  Standard_Real TheTol3D = RatioTol * myTolReached3d ;
+
+  if (myData.myMinFactorXYZ>1.5e-7)
+    TheTol3D = TheTol3D / myData.myMinFactorXYZ;
+
+  //cout << "Tol 3D: " << TheTol3D << endl;
+  return TheTol3D ;
+}
+
+//=======================================================================
+//function : TolReached2d
+//purpose  :
+//=======================================================================
+Standard_Real ApproxInt_Approx::TolReached2d() const
+{
+  Standard_Real TheTol2D = RatioTol * myTolReached2d ;
+
+  if (myData.myMinFactorUV>1.5e-7)
+    TheTol2D = TheTol2D / myData.myMinFactorUV;
+
+  //cout << "Tol 2D: " << TheTol2D << endl;
+  return TheTol2D ;
+}
+
+//=======================================================================
+//function : IsDone
+//purpose  :
+//=======================================================================
+Standard_Boolean ApproxInt_Approx::IsDone() const
+{
+  if(myData.myBezierApprox)
+  { 
+    return(myComputeLineBezier.NbMultiCurves() > 0);
+    //-- Lorsque la tolerance n est pas atteinte et qu il 
+    //-- faudrait rajouter des points sur la ligne
+    //-- les approx sortent avec la meilleure tolerance
+    //-- atteinte.  ( Pas de rajout de points ds cette version)
+    //-- return(myTolReached);
   }
-  if(ApproxU1V1) { 
-    ComputeTrsf2d(theline,U1o,A1u,V1o,A1v,Standard_True);
+  else
+  {
+    return(myComputeLine.IsToleranceReached());
   }
-  else { 
-    U1o=V1o=0.0; A1u=A1v=1.0;
+}
+
+//=======================================================================
+//function : Value
+//purpose  :
+//=======================================================================
+const AppParCurves_MultiBSpCurve& ApproxInt_Approx::Value(const Standard_Integer ) const
+{
+  if(myData.myBezierApprox)
+  {
+    return(myBezToBSpl.Value());
   }
-  if(ApproxU2V2) { 
-    Standard_Real UVResRatio = ThePSurfaceTool::UResolution(PSurf,1.)/
-                               ThePSurfaceTool::VResolution(PSurf,1.);
-    ComputeTrsf2d(theline,U2o,A2u,V2o,A2v,Standard_False,UVResRatio);
+  else
+  {
+    return(myComputeLine.Value());
   }
-  else { 
-    U2o=V2o=0.0; A2u=A2v=1.0;
+}
+
+//=======================================================================
+//function : fillData
+//purpose  : Fill ApproxInt data structure.
+//=======================================================================
+void ApproxInt_Approx::fillData(const Handle(TheWLine)& theline,
+                                const Standard_Boolean ApproxXYZ,
+                                const Standard_Boolean ApproxU1V1,
+                                const Standard_Boolean ApproxU2V2)
+{
+  if(ApproxXYZ)
+  {
+    ComputeTrsf3d(theline, myData.Xo, myData.Ax, myData.Yo, myData.Ay, myData.Zo, myData.Az);
   }
-  
-  //-deb-  cout << "ApproxInt_Approx -- NbPntMax = " << myNbPntMax    << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Tol3D    = " << myTol3d       << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Tol2D    = " << myTol2d       << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- RelTol   = " << (myRelativeTol ? "RELATIVE" : "ABSOLUTE") << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- Xo = " << Xo << " Yo = " << Yo << " Zo = " << Zo << endl ;
-  //-deb-  cout << "ApproxInt_Approx -- Ax = " << Ax << " Ay = " << Ay << " Az = " << Az << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- U1o = " << U1o << " V1o = " << V1o << " A1u = " << A1u << " A1v = " << A1v << endl ; 
-  //-deb-  cout << "ApproxInt_Approx -- U2o = " << U2o << " V2o = " << V2o << " A2u = " << A2u << " A2v = " << A2v << endl ; 
-  
-  
-  Standard_Real A3d = MINABS3(Ax,Ay,Az);
-  if((A3d < myMinFactorXYZ) || (myMinFactorXYZ == 0.0)) { 
-    myMinFactorXYZ = A3d;
+  else
+  {
+    myData.Xo = myData.Yo = myData.Zo = 0.0;
+    myData.Ax = myData.Ay = myData.Az = 1.0;
   }
-  
-  Standard_Real A2d = MINABS4(A1u,A1v,A2u,A2v);
-  if((A2d < myMinFactorUV) || (myMinFactorUV == 0.0)) { 
-    myMinFactorUV = A2d;
+
+  if(ApproxU1V1)
+  {
+    ComputeTrsf2d(theline, myData.U1o, myData.A1u, myData.V1o, myData.A1v,Standard_True, myUVRes1);
+  }
+  else
+  {
+    myData.U1o = myData.V1o = 0.0;
+    myData.A1u = myData.A1v = 1.0;
+  }
+
+  if(ApproxU2V2)
+  {
+    ComputeTrsf2d(theline, myData.U2o, myData.A2u, myData.V2o, myData.A2v, Standard_False, myUVRes2);
+  }
+  else
+  {
+    myData.U2o = myData.V2o = 0.0;
+    myData.A2u = myData.A2v = 1.0;
+  }
+
+  Standard_Real A3d = MINABS3(myData.Ax, myData.Ay, myData.Az);
+  if((A3d < myData.myMinFactorXYZ) || (myData.myMinFactorXYZ == 0.0))
+  {
+    myData.myMinFactorXYZ = A3d;
+  }
+
+  Standard_Real A2d = MINABS4(myData.A1u, myData.A1v, myData.A2u, myData.A2v);
+  if((A2d < myData.myMinFactorUV) || (myData.myMinFactorUV == 0.0))
+  {
+    myData.myMinFactorUV = A2d;
   }
+}
+
+//=======================================================================
+//function : prepareDS
+//purpose  :
+//=======================================================================
+void ApproxInt_Approx::prepareDS(const Standard_Boolean theApproxXYZ,
+                                 const Standard_Boolean theApproxU1V1,
+                                 const Standard_Boolean theApproxU2V2,
+                                 const Standard_Integer theIndicemin,
+                                 const Standard_Integer theIndicemax)
+{
+  myData.myMinFactorXYZ = 0.0;
+  myData.myMinFactorUV  = 0.0;
+  myTolReached3d = myTolReached2d = 0.0;
+  myUVRes1 = myUVRes2 = 1.0;
+  myData.ApproxU1V1 = theApproxU1V1;
+  myData.ApproxU2V2 = theApproxU2V2;
+  myData.ApproxXYZ = theApproxXYZ;
+  myData.indicemin = theIndicemin;
+  myData.indicemax = theIndicemax;
+  myData.nbpntmax = myNbPntMax;
 
   Approx_ParametrizationType parametrization;
   myComputeLineBezier.Parametrization(parametrization);
-  
-  
-  if(myRelativeTol==Standard_False) { 
-    myComputeLine.Init(myDegMin,
-                      myDegMax,
-                      myTol3d*myMinFactorXYZ,
-                      myTol2d*myMinFactorUV,
-                      myNbIterMax,
-                      cut,
-                      parametrization);
-    myComputeLineBezier.Init(myDegMin,
-                            myDegMax,
-                            myTol3d*myMinFactorXYZ,
-                            myTol2d*myMinFactorUV,
-                            myNbIterMax,
-                            cut,
-                            parametrization);
+  myData.parametrization = parametrization;
+}
+
+//=======================================================================
+//function : buildKnots
+//purpose  :
+//=======================================================================
+void ApproxInt_Approx::buildKnots(const Handle(TheWLine)& theline,
+                                  const Standard_Address thePtrSVSurf)
+{
+  myKnots.Clear();
+  if(myData.myBezierApprox)
+  {
+    ApproxInt_TheMultiLine aTestLine(theline,
+      thePtrSVSurf,
+      ((myData.ApproxXYZ)? 1 : 0),
+      ((myData.ApproxU1V1)? 1: 0) + ((myData.ApproxU2V2)? 1: 0),
+      myData.Xo, myData.Ax, myData.Yo, myData.Ay, myData.Zo, myData.Az,
+      myData.U1o, myData.A1u, 
+      myData.V1o, myData.A1v,
+      myData.U2o, myData.A2u,
+      myData.V2o ,myData.A2v,
+      myData.ApproxU1V1,
+      myData.indicemin,
+      myData.indicemax);
+
+    Standard_Integer nbp3d = aTestLine.NbP3d();
+    Standard_Integer nbp2d = aTestLine.NbP2d();
+    TColgp_Array1OfPnt aTabPnt3d(1, Max(1, nbp3d));
+    TColgp_Array1OfPnt2d aTabPnt2d(1, Max(1, nbp2d));
+    TColgp_Array1OfPnt aPntXYZ(myData.indicemin, myData.indicemax);
+    TColgp_Array1OfPnt2d aPntU1V1(myData.indicemin, myData.indicemax);
+    TColgp_Array1OfPnt2d aPntU2V2(myData.indicemin, myData.indicemax);
+    Standard_Integer i;
+    for(i = myData.indicemin; i <= myData.indicemax; ++i)
+    {
+      if (nbp3d != 0 && nbp2d != 0) aTestLine.Value(i, aTabPnt3d, aTabPnt2d);
+      else if (nbp2d != 0)          aTestLine.Value(i, aTabPnt2d);
+      else if (nbp3d != 0)          aTestLine.Value(i, aTabPnt3d);
+      //
+      if(nbp3d > 0)
+      {
+        aPntXYZ(i) = aTabPnt3d(1);
+      }
+      if(nbp2d > 1)
+      {
+        aPntU1V1(i) = aTabPnt2d(1);
+        aPntU2V2(i) = aTabPnt2d(2);
+      }
+      else if(nbp2d > 0)
+      {
+        if(myData.ApproxU1V1)
+        {
+          aPntU1V1(i) = aTabPnt2d(1);
+        }
+        else
+        {
+          aPntU2V2(i) = aTabPnt2d(1);
+        }
+      }
+    }
+
+    Standard_Integer aMinNbPnts = myData.nbpntmax;
+
+    // Expected parametrization.
+    math_Vector aPars(myData.indicemin, myData.indicemax);
+    Parameters(aTestLine, myData.indicemin, myData.indicemax, myData.parametrization, aPars);
+
+    ApproxInt_KnotTools::BuildKnots(aPntXYZ, aPntU1V1, aPntU2V2, aPars,
+      myData.ApproxXYZ, myData.ApproxU1V1, myData.ApproxU2V2, aMinNbPnts, myKnots);
+  }
+  else
+  {
+    myKnots.Append(myData.indicemin);
+    myKnots.Append(myData.indicemax);
+  }
+}
+
+//=======================================================================
+//function : buildCurve
+//purpose  :
+//=======================================================================
+void ApproxInt_Approx::buildCurve(const Handle(TheWLine)& theline,
+                                  const Standard_Address thePtrSVSurf)
+{
+  if(myData.myBezierApprox)
+  {
+    myBezToBSpl.Reset();
   }
-  else { 
-    myComputeLine.Init(myDegMin,
-                      myDegMax,
-                      myTol3d,
-                      myTol2d,
-                      myNbIterMax,
-                      cut,
-                      parametrization);
-    myComputeLineBezier.Init(myDegMin,
-                            myDegMax,
-                            myTol3d,
-                            myTol2d,
-                            myNbIterMax,
-                            cut,
-                            parametrization);
-  }     
-  
-  
-  
-  do {
-    
-    ApproxInt_TheMultiLine myMultiLine(theline,
-                                      ptrsvsurf,
-                                      ((ApproxXYZ)? 1 : 0),
-                                      ((ApproxU1V1)? 1: 0) + ((ApproxU2V2)? 1: 0),
-                                      Xo,Ax,Yo,Ay,Zo,Az,U1o,A1u,V1o,A1v,U2o,A2u,V2o,A2v,
-                                      ApproxU1V1,
-                                      imin,
-                                      imax);
-    if(myApproxBez) { 
-      myComputeLineBezier.Perform(myMultiLine);
 
-#ifdef OCCT_DEBUG
-      //myMultiLine.Dump();
-#endif
+  Standard_Integer kind = myKnots.Lower();
+  Standard_Integer imin, imax;
+  myTolReached = Standard_True;
+  Standard_Boolean OtherInter = Standard_False;
+  do
+  {
+    // Base cycle: iterate over knots.
+    imin = myKnots(kind);
+    imax = myKnots(kind+1);
+    ApproxInt_TheMultiLine myMultiLine(theline,
+      thePtrSVSurf,
+      ((myData.ApproxXYZ)? 1 : 0),
+      ((myData.ApproxU1V1)? 1: 0) + ((myData.ApproxU2V2)? 1: 0),
+      myData.Xo, myData.Ax, myData.Yo,myData.Ay, myData.Zo,myData.Az,
+      myData.U1o, myData.A1u, myData.V1o, myData.A1v,
+      myData.U2o, myData.A2u, myData.V2o, myData.A2v,
+      myData.ApproxU1V1,
+      imin,
+      imax);
 
+    if(myData.myBezierApprox)
+    {
+      myComputeLineBezier.Perform(myMultiLine);
       if (myComputeLineBezier.NbMultiCurves() == 0)
-       return;
+        return;
       myTolReached&=myComputeLineBezier.IsToleranceReached();
     }
-    else { 
+    else
+    {
       myComputeLine.Perform(myMultiLine);
     }
     UpdateTolReached();
-    
+
     Standard_Integer indice3d,indice2d1,indice2d2;
     indice3d = 1; 
     indice2d1= 2;
     indice2d2= 3;
-    if(!ApproxXYZ)  { indice2d1--; indice2d2--; } 
-    if(!ApproxU1V1) { indice2d2--; } 
-    if(ApproxXYZ) { 
+    if(!myData.ApproxXYZ)  { indice2d1--; indice2d2--; } 
+    if(!myData.ApproxU1V1) { indice2d2--; } 
+    if(myData.ApproxXYZ)
+    { 
       Standard_Real ax,bx,ay,by,az,bz;
-      ax = 1.0/Ax;   bx = -Xo*ax;
-      ay = 1.0/Ay;   by = -Yo*ay;
-      az = 1.0/Az;   bz = -Zo*az;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform(indice3d,bx,ax,by,ay,bz,az);
-       }
+      ax = 1.0/myData.Ax;   bx = -myData.Xo*ax;
+      ay = 1.0/myData.Ay;   by = -myData.Yo*ay;
+      az = 1.0/myData.Az;   bz = -myData.Zo*az;
+      if(myData.myBezierApprox)
+      {
+        for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--)
+        {
+          myComputeLineBezier.ChangeValue(nbmc).Transform(indice3d,bx,ax,by,ay,bz,az);
+        }
       }
-      else { 
-       myComputeLine.ChangeValue().Transform(indice3d,bx,ax,by,ay,bz,az);
+      else
+      {
+        myComputeLine.ChangeValue().Transform(indice3d,bx,ax,by,ay,bz,az);
       }
     }
-    if(ApproxU1V1) { 
+    if(myData.ApproxU1V1)
+    {
       Standard_Real ax,bx,ay,by;
-      ax = 1.0/A1u;   bx = -U1o*ax;
-      ay = 1.0/A1v;   by = -V1o*ay;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d1,bx,ax,by,ay);
-       }
+      ax = 1.0/myData.A1u;   bx = -myData.U1o*ax;
+      ay = 1.0/myData.A1v;   by = -myData.V1o*ay;
+      if(myData.myBezierApprox) {
+        for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--)
+        {
+          myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d1,bx,ax,by,ay);
+        }
       }
-      else { 
-       myComputeLine.ChangeValue().Transform2d(indice2d1,bx,ax,by,ay);
+      else
+      {
+        myComputeLine.ChangeValue().Transform2d(indice2d1,bx,ax,by,ay);
       }
     }
-    if(ApproxU2V2) { 
+    if(myData.ApproxU2V2)
+    {
       Standard_Real ax,bx,ay,by;
-      ax = 1.0/A2u;   bx = -U2o*ax;
-      ay = 1.0/A2v;   by = -V2o*ay;
-      if(myApproxBez) { 
-       for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--) { 
-         myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d2,bx,ax,by,ay);
-       }
+      ax = 1.0/myData.A2u;   bx = -myData.U2o*ax;
+      ay = 1.0/myData.A2v;   by = -myData.V2o*ay;
+      if(myData.myBezierApprox)
+      {
+        for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--)
+        {
+          myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d2,bx,ax,by,ay);
+        }
       }
-      else { 
-       myComputeLine.ChangeValue().Transform2d(indice2d2,bx,ax,by,ay);
+      else
+      {
+        myComputeLine.ChangeValue().Transform2d(indice2d2,bx,ax,by,ay);
       }
     }
+
     OtherInter = Standard_False;
-    if(myApproxBez)
+    if(myData.myBezierApprox)
     {
-      for(Standard_Integer nbmc = 1;
-        nbmc <= myComputeLineBezier.NbMultiCurves() ;
+      for(Standard_Integer nbmc = 1; 
+        nbmc <= myComputeLineBezier.NbMultiCurves();
         nbmc++)
-      { 
+      {
           myBezToBSpl.Append(myComputeLineBezier.Value(nbmc));
       }
-      if(imax<indicemax)
+      kind++;
+      if(kind < myKnots.Upper())
       {
-        imin = imax;
-        imax = imin+nbpntbez;
         OtherInter = Standard_True;
-        if((indicemax-imax)<(nbpntbez/2))
-        {
-          imax = indicemax;
-        }
-        imax = CorrectFinishIdx(imin, imax, theline);
       }
     }
   }
-  while(OtherInter); 
-  if(myApproxBez) { 
-    myBezToBSpl.Perform();
-  }
-}
-//--------------------------------------------------------------------------------
-Standard_Integer ApproxInt_Approx::NbMultiCurves() const { 
-  //  return(myComputeLine.NbMultiCurves());
-  return 1;
-}
-//--------------------------------------------------------------------------------
-void ApproxInt_Approx::UpdateTolReached() {
-
-  if (myApproxBez) {
-    Standard_Integer ICur ;
-    Standard_Integer NbCurves = myComputeLineBezier.NbMultiCurves() ;
-    for (ICur = 1 ; ICur <= NbCurves ; ICur++) {
-      Standard_Real Tol3D, Tol2D ;
-      myComputeLineBezier.Error (ICur, Tol3D, Tol2D) ;
-      myTolReached3d = Max(myTolReached3d, Tol3D);
-      myTolReached2d = Max(myTolReached2d, Tol2D);
-    }
-  } else {
-    myComputeLine.Error (myTolReached3d, myTolReached2d);
-  }
-}
-//--------------------------------------------------------------------------------
-Standard_Real ApproxInt_Approx::TolReached3d() const { 
-
-  Standard_Real TheTol3D = RatioTol * myTolReached3d ;
-  //modified by NIZNHY-PKV Mon Aug 27 14:21:33 2007f
-  //if (myMinFactorXYZ)
-  //TheTol3D = TheTol3D / myMinFactorXYZ ;
-  if (myMinFactorXYZ>1.5e-7) {
-    TheTol3D = TheTol3D / myMinFactorXYZ ;
-  }
-  //modified by NIZNHY-PKV Mon Aug 27 14:21:50 2007t
-  return TheTol3D ;
-}
-//--------------------------------------------------------------------------------
-Standard_Real ApproxInt_Approx::TolReached2d() const {
-
-  Standard_Real TheTol2D = RatioTol * myTolReached2d ;
-  //modified by NIZNHY-PKV Mon Aug 27 14:20:50 2007f
-  //if (myMinFactorUV)
-  //TheTol2D = TheTol2D / myMinFactorUV ;
-  if (myMinFactorUV>1.5e-7) {
-    TheTol2D = TheTol2D / myMinFactorUV ;
-  }
-  //modified by NIZNHY-PKV Mon Aug 27 14:20:55 2007t
-  return TheTol2D ;
-}
-//--------------------------------------------------------------------------------
-Standard_Boolean ApproxInt_Approx::IsDone() const { 
-  if(myApproxBez) { 
-    return(myComputeLineBezier.NbMultiCurves() > 0);
-                           //-- Lorsque la tolerance n est pas atteinte et qu il 
-                           //-- faudrait rajouter des points sur la ligne
-                           //-- les approx sortent avec la meilleure tolerance
-                           //-- atteinte.  ( Pas de rajout de points ds cette version)
-    //-- return(myTolReached);
-  }
-  else { 
-    return(myComputeLine.IsToleranceReached());
-  }
-}
-//--------------------------------------------------------------------------------
-const AppParCurves_MultiBSpCurve& ApproxInt_Approx::Value(const Standard_Integer ) const { 
-  if(myApproxBez) { 
-    return(myBezToBSpl.Value());
-  }
-  else { 
-    return(myComputeLine.Value());
-  }
-}
-//--------------------------------------------------------------------------------
-Standard_Integer ApproxInt_Approx::CorrectFinishIdx(const Standard_Integer theMinIdx,
-                                                    const Standard_Integer theMaxIdx,
-                                                    const Handle(TheWLine)& theline)
-{
-  const Standard_Real aNullCoeff = 1.0e-16;
-  Standard_Real aLimitMaxCoeff = 1.0 / 2500.0;
-  Standard_Real aDist = theline->Point(theMinIdx).Value().SquareDistance(
-                        theline->Point(theMinIdx + 1).Value());
+  while(OtherInter);
 
-  for(Standard_Integer anIdx = theMinIdx + 1; anIdx < theMaxIdx - 1; anIdx++)
+  if(myData.myBezierApprox)
   {
-    Standard_Real aNextDist = theline->Point(anIdx).Value().SquareDistance(
-                              theline->Point(anIdx + 1).Value());
-    Standard_Real aCoeff = Min (aNextDist, aDist) / Max (aNextDist, aDist);
-
-    //
-    if (aCoeff < aLimitMaxCoeff &&        // Base criteria.
-        aNextDist > aDist &&              // Step increasing.
-        aNextDist > aNullCoeff &&         // Avoid separation in case of too small step.
-        aDist > aNullCoeff)               // Usually found when purger not invoked (blend).
-    {
-      return anIdx;
-    }
-    aDist = aNextDist;
+    myBezToBSpl.Perform();
   }
-  return theMaxIdx;
-}
+}
\ No newline at end of file
diff --git a/src/ApproxInt/ApproxInt_KnotTools.cxx b/src/ApproxInt/ApproxInt_KnotTools.cxx
new file mode 100644 (file)
index 0000000..f5e7b07
--- /dev/null
@@ -0,0 +1,613 @@
+// Copyright (c) 1999-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#include <ApproxInt_KnotTools.hxx>
+#include <TColgp_Array1OfPnt2d.hxx>
+#include <TColStd_Array1OfReal.hxx>
+#include <TColStd_HArray1OfReal.hxx>
+#include <TColStd_HArray1OfInteger.hxx>
+#include <math_Vector.hxx>
+#include <Geom_BSplineCurve.hxx>
+#include <Geom2d_BSplineCurve.hxx>
+#include <GeomInt_TheMultiLineOfWLApprox.hxx>
+#include <NCollection_Sequence.hxx>
+#include <NCollection_List.hxx>
+#include <PLib.hxx>
+#include <Precision.hxx>
+#include <NCollection_Vector.hxx>
+#include <TColgp_Array1OfPnt.hxx>
+
+// (Sqrt(5.0) - 1.0) / 4.0
+static const Standard_Real aSinCoeff = 0.30901699437494742410229341718282;
+static const Standard_Integer aMaxPntCoeff = 15;
+
+
+//=======================================================================
+//function : EvalCurv
+//purpose  : Evaluate curvature in dim-dimension point.
+//=======================================================================
+static Standard_Real EvalCurv(const Standard_Real dim, 
+                              const Standard_Real* V1,
+                              const Standard_Real* V2)
+{
+  // Really V1 and V2 are arrays of size dim;
+  // V1 is first derivative, V2 is second derivative
+  // of n-dimension curve
+  // Curvature is curv = |V1^V2|/|V1|^3
+  // V1^V2 is outer product of two vectors:
+  // P(i,j) = V1(i)*V2(j) - V1(j)*V2(i);
+  Standard_Real mp = 0.;
+  Standard_Integer i, j;
+  Standard_Real p;
+  for(i = 1; i < dim; ++i)
+  {
+    for(j = 0; j < i; ++j)
+    {
+      p = V1[i]*V2[j] - V1[j]*V2[i];
+      mp += p*p;
+    }
+  }
+  //mp *= 2.; //P(j,i) = -P(i,j);
+  mp = Sqrt(mp);
+  //
+  Standard_Real q = 0.;
+  for(i = 0; i < dim; ++i)
+  {
+    q += V1[i]*V1[i];
+  }
+  q = Sqrt(q);
+  //
+  Standard_Real curv = mp / (q*q*q);
+
+  return curv;
+}
+
+//=======================================================================
+//function : ComputeKnotInds
+//purpose  : 
+//=======================================================================
+void ApproxInt_KnotTools::ComputeKnotInds(const NCollection_LocalArray<Standard_Real>& theCoords,
+                                          const Standard_Integer theDim, 
+                                          const math_Vector& thePars,
+                                          NCollection_Sequence<Standard_Integer>& theInds)
+{
+  //I: Create discrete curvature.
+  NCollection_Sequence<Standard_Integer> aFeatureInds;
+  TColStd_Array1OfReal aCurv(thePars.Lower(), thePars.Upper());
+  // Arrays are allocated for max theDim = 7: 1 3d curve + 2 2d curves.
+  Standard_Real Val[21], Par[3], Res[21];
+  Standard_Integer i, j, k, l, m, ic;
+  Standard_Real aMaxCurv = 0.;
+  Standard_Integer dim = theDim;
+  //
+  i = aCurv.Lower();
+  for(j = 0; j < 3; ++j)
+  {
+    k = i+j;
+    ic = (k - aCurv.Lower()) * dim;
+    l = dim*j;
+    for(m = 0; m < dim; ++m)
+    {
+      Val[l + m] = theCoords[ic + m];
+    }
+    Par[j] = thePars(k);
+  }
+  PLib::EvalLagrange(Par[0], 2, 2, dim, *Val, *Par, *Res);
+  //
+  aCurv(i) = EvalCurv(dim, &Res[dim], &Res[2*dim]);
+  //
+  if(aCurv(i) > aMaxCurv)
+  {
+    aMaxCurv = aCurv(i);
+  }
+  //
+  for(i = aCurv.Lower()+1; i < aCurv.Upper(); ++i)
+  {
+    for(j = 0; j < 3; ++j)
+    {
+      k = i+j-1;
+      ic = (k - aCurv.Lower()) * dim;
+      l = dim*j;
+      for(m = 0; m < dim; ++m)
+      {
+        Val[l + m] = theCoords[ic + m];
+      }
+      Par[j] = thePars(k);
+    }
+    PLib::EvalLagrange(Par[1], 2, 2, dim, *Val, *Par, *Res);
+    //
+    aCurv(i) = EvalCurv(dim, &Res[dim], &Res[2*dim]);
+    if(aCurv(i) > aMaxCurv)
+    {
+      aMaxCurv = aCurv(i);
+    }
+  }
+  //
+  i = aCurv.Upper();
+  for(j = 0; j < 3; ++j)
+  {
+    k = i+j-2;
+    ic = (k - aCurv.Lower()) * dim;
+    l = dim*j;
+    for(m = 0; m < dim; ++m)
+    {
+      Val[l + m] = theCoords[ic + m];
+    }
+    Par[j] = thePars(k);
+  }
+  PLib::EvalLagrange(Par[2], 2, 2, dim, *Val, *Par, *Res);
+  //
+  aCurv(i) = EvalCurv(dim, &Res[dim], &Res[2*dim]);
+  if(aCurv(i) > aMaxCurv)
+  {
+    aMaxCurv = aCurv(i);
+  }
+
+#if defined(APPROXINT_KNOTTOOLS_DEBUG) || defined(OCCT_DEBUG)
+  cout << "Discrete curvature array is" << endl;
+  for(i = aCurv.Lower(); i <= aCurv.Upper(); ++i)
+  {
+    cout << i << " " << aCurv(i) << endl;
+  }
+#endif
+
+  theInds.Append(aCurv.Lower());
+  if(aMaxCurv <= Precision::Confusion())
+  {
+    // Linear case.
+    theInds.Append(aCurv.Upper());
+    return;
+  }
+
+  // II: Find extremas of curvature.
+  // Not used Precision::PConfusion, by different from "param space" eps nature.
+  Standard_Real eps  = 1.0e-9,
+                eps1 = 1.0e3 * eps;
+  for(i = aCurv.Lower() + 1; i < aCurv.Upper(); ++i)
+  {
+    Standard_Real d1 = aCurv(i) - aCurv(i - 1),
+                  d2 = aCurv(i) - aCurv(i + 1),
+                  ad1 = Abs(d1), ad2 = Abs(d2);
+
+    if(d1*d2 > 0. && ad1 > eps && ad2 > eps)
+    {
+      if(i != theInds.Last())
+      {
+        theInds.Append(i);
+        aFeatureInds.Append(i);
+      }
+    }
+    else if((ad1 < eps && ad2 > eps1) || (ad1 > eps1 && ad2 < eps))
+    {
+      if(i != theInds.Last())
+      {
+        theInds.Append(i);
+        aFeatureInds.Append(i);
+      }
+    }
+    else if(aCurv(i)*aCurv(i + 1) < 0.0)
+    {
+      if(Abs(aCurv(i)) < Abs(aCurv(i + 1)))
+      {
+        if(i != theInds.Last())
+        {
+          theInds.Append(i);
+          aFeatureInds.Append(i);
+        }
+      }
+      else
+      {
+        if(i+1 != theInds.Last())
+        {
+          theInds.Append(i + 1);
+          aFeatureInds.Append(i + 1);
+        }
+      }
+    }
+  }
+  if(aCurv.Upper() != theInds.Last())
+  {
+    theInds.Append(aCurv.Upper());
+  }
+
+#if defined(APPROXINT_KNOTTOOLS_DEBUG) || defined(OCCT_DEBUG)
+  {
+    cout << "Feature indices new: " << endl;
+    i;
+    for(i = theInds.Lower(); i <= theInds.Upper();  ++i)
+    {
+      cout << i << " : " << theInds(i) << endl;
+    }
+  }
+#endif
+
+  //III: Put knots in monotone intervals of curvature.
+  Standard_Boolean Ok;
+  i = 1;
+  do
+  {
+    i++;
+    //
+    Ok = InsKnotBefI(i, aCurv, theCoords, dim, theInds, Standard_True); 
+    if(Ok)
+    {
+      i--;
+    }
+  }
+  while(i < theInds.Length());
+
+  //IV: Cheking feature points.
+  j = 2;
+  for(i = 1; i <= aFeatureInds.Length(); ++i)
+  {
+    Standard_Integer anInd = aFeatureInds(i);
+    for(; j <= theInds.Length() - 1;)
+    {
+      if(theInds(j) == anInd)
+      {
+        Standard_Integer anIndPrev = theInds(j-1);
+        Standard_Integer anIndNext = theInds(j+1);
+        Standard_Real sina;
+        Standard_Integer ici = (anIndPrev - aCurv.Lower()) * theDim,
+          ici1 = (anIndNext - aCurv.Lower()) * theDim,
+          icm = (anInd - aCurv.Lower()) * theDim;
+        NCollection_LocalArray<Standard_Real> V1(theDim), V2(theDim);
+        Standard_Integer k,l;
+        Standard_Real mp = 0., m1 = 0., m2 = 0.;
+        Standard_Real p;
+        for(k = 0; k < theDim; ++k)
+        {
+          V1[k] = theCoords[icm + k] - theCoords[ici + k];
+          m1 += V1[k]*V1[k];
+          V2[k] = theCoords[ici1 + k] - theCoords[icm + k];
+          m2 += V2[k]*V2[k];
+        }
+        for(k = 1; k < theDim; ++k)
+        {
+          for(l = 0; l < k; ++l)
+          {
+            p = V1[k]*V2[l] - V1[l]*V2[k];
+            mp += p*p;
+          }
+        }
+        //mp *= 2.; //P(j,i) = -P(i,j);
+        //
+        sina = mp/(m1*m2);
+        sina = Sqrt(sina);
+
+        if(sina  > aSinCoeff) 
+        {
+          //Insert new knots
+          Standard_Real d1 = Abs(aCurv(anInd) - aCurv(anIndPrev));
+          Standard_Real d2 = Abs(aCurv(anInd) - aCurv(anIndNext));
+          if(d1 > d2)
+          {
+            Ok = InsKnotBefI(j, aCurv, theCoords, dim, theInds, Standard_False);
+            if(Ok)
+            {
+              j++;
+            }
+            else
+            {
+              break;
+            }
+          }
+          else
+          {
+            Ok = InsKnotBefI(j+1, aCurv, theCoords, dim, theInds, Standard_False);
+            if(!Ok)
+            {
+              break;
+            }
+          }
+        }
+        else
+        {
+          j++;
+          break;
+        }
+      }
+      else
+      {
+        j++;
+      }
+    }
+  }
+  //
+}
+
+
+//=======================================================================
+//function : FilterKnots
+//purpose  :
+//=======================================================================
+void ApproxInt_KnotTools::FilterKnots(NCollection_Sequence<Standard_Integer>& theInds, 
+                                      const Standard_Integer theMinNbPnts,
+                                      NCollection_Vector<Standard_Integer>& theLKnots)
+{
+  // Maximum number of points per knot interval.
+  Standard_Integer aMaxNbPnts = aMaxPntCoeff*theMinNbPnts;
+  Standard_Integer i = 1;
+  Standard_Integer aMinNbStep = theMinNbPnts / 2;
+
+  // I: Filter too big number of points per knot interval.
+  while(i < theInds.Length())
+  {
+    Standard_Integer nbint = theInds(i + 1) - theInds(i) + 1;
+    if(nbint <= aMaxNbPnts)
+    {
+      ++i;
+      continue;
+    }
+    else
+    {
+      Standard_Integer ind = theInds(i) + nbint / 2;
+      theInds.InsertAfter(i, ind);     
+    }
+  }
+
+  // II: Filter poins with too small amount of points per knot interval.
+  i = 1;
+  theLKnots.Append(theInds(i));
+  Standard_Integer anIndsPrev = theInds(i);
+  for(i = 2; i <= theInds.Length(); ++i)
+  {
+    if(theInds(i) - anIndsPrev <= theMinNbPnts)
+    {
+      if (i != theInds.Length())
+      {
+        Standard_Integer anIdx = i + 1;
+        for( ; anIdx <= theInds.Length(); ++anIdx)
+        {
+          if (theInds(anIdx) - anIndsPrev > theMinNbPnts)
+            break;
+        }
+        anIdx--;
+
+        Standard_Integer aMidIdx = (theInds(anIdx) + anIndsPrev) / 2;
+        if (aMidIdx - anIndsPrev        < theMinNbPnts &&
+            aMidIdx - theInds(anIdx)    < theMinNbPnts &&
+            theInds(anIdx) - anIndsPrev >= aMinNbStep)
+        {
+          // Bad distribution points merge into one knot interval.
+          theLKnots.Append(theInds(anIdx));
+          anIndsPrev = theInds(anIdx);
+          i = anIdx;
+        }
+        else if (anIdx == theInds.Upper() && // Last point obtained.
+                 theLKnots.Length() >= 2) // It is possible to modify last item.
+        {
+          // Current bad interval from i to last.
+          // Trying to add knot to divide sequence on two parts:
+          // Last good index -> Last index - theMinNbPnts -> Last index
+          Standard_Integer aLastGoodIdx = theLKnots.Value(theLKnots.Upper() - 1);
+          if ( theInds.Last() - 2 * theMinNbPnts >= aLastGoodIdx)
+          {
+            theLKnots(theLKnots.Upper()) = theInds.Last() - theMinNbPnts;
+            theLKnots.Append(theInds.Last());
+            anIndsPrev = theInds(anIdx);
+            i = anIdx;
+          }
+        }
+      } // if (i != theInds.Length())
+      continue;
+    }
+    else
+    {
+      theLKnots.Append(theInds(i));
+      anIndsPrev = theInds(i);
+    }
+  }
+
+  // III: Fill Last Knot.
+  if(theLKnots.Length() < 2)
+  {
+    theLKnots.Append(theInds.Last());
+  }
+  else
+  {
+    if(theLKnots.Last() < theInds.Last())
+    {
+      theLKnots(theLKnots.Upper()) = theInds.Last();
+    }
+  }
+}
+//=======================================================================
+//function : InsKnotBefI
+//purpose  :
+//=======================================================================
+Standard_Boolean ApproxInt_KnotTools::InsKnotBefI(const Standard_Integer theI,
+                                                  const TColStd_Array1OfReal& theCurv,
+                                                  const NCollection_LocalArray<Standard_Real>& theCoords,
+                                                  const Standard_Integer theDim, 
+                                                  NCollection_Sequence<Standard_Integer>& theInds,
+                                                  const Standard_Boolean ChkCurv)
+{
+  Standard_Integer anInd1 = theInds(theI);
+  Standard_Integer anInd = theInds(theI - 1);
+  //
+  if((anInd1-anInd) == 1)
+  {
+    return Standard_False;
+  }
+  //
+  Standard_Real curv = 0.5*(theCurv(anInd) + theCurv(anInd1));
+  Standard_Integer mid = 0, j, jj;
+  const Standard_Real aLimitCurvatureChange = 3.0;
+  for(j = anInd+1; j < anInd1; ++j)
+  {
+    mid = 0;
+
+    // I: Curvature change criteria:
+    // Non-null curvature.
+    if (theCurv(j)     > Precision::Confusion() && 
+        theCurv(anInd) > Precision::Confusion() )
+    {
+      if (theCurv(j) / theCurv(anInd) > aLimitCurvatureChange || 
+          theCurv(j) / theCurv(anInd) < 1.0 / aLimitCurvatureChange)
+      {
+        // Curvature on current interval changed more than 3 times.
+        mid = j;
+        theInds.InsertBefore(theI, mid);
+        return Standard_True;
+      }
+    }
+
+    // II: Angular criteria:
+    Standard_Real ac = theCurv(j - 1), ac1 = theCurv(j);
+    if((curv >= ac && curv <= ac1) || (curv >= ac1 && curv <= ac))
+    {
+      if(Abs(curv - ac) < Abs(curv - ac1))
+      {
+        mid = j - 1;
+      }
+      else
+      {
+        mid = j;
+      }
+    }
+    if(mid == anInd)
+    {
+      mid++;
+    }
+    if(mid == anInd1)
+    {
+      mid--;
+    }
+    if(mid > 0)
+    {
+      if(ChkCurv)
+      {
+        Standard_Real sina;
+        Standard_Integer ici = (anInd - theCurv.Lower()) * theDim,
+          ici1 = (anInd1 - theCurv.Lower()) * theDim,
+          icm = (mid - theCurv.Lower()) * theDim;
+        NCollection_LocalArray<Standard_Real> V1(theDim), V2(theDim);
+        Standard_Integer i;
+        Standard_Real mp = 0., m1 = 0., m2 = 0.;
+        Standard_Real p;
+        for(i = 0; i < theDim; ++i)
+        {
+          V1[i] = theCoords[icm + i] - theCoords[ici + i];
+          m1 += V1[i]*V1[i];
+          V2[i] = theCoords[ici1 + i] - theCoords[icm + i];
+          m2 += V2[i]*V2[i];
+        }
+        for(i = 1; i < theDim; ++i)
+        {
+          for(jj = 0; jj < i; ++jj)
+          {
+            p = V1[i]*V2[jj] - V1[jj]*V2[i];
+            mp += p*p;
+          }
+        }
+        //mp *= 2.; //P(j,i) = -P(i,j);
+        //
+        sina = mp/(m1*m2);
+        sina = Sqrt(sina);
+
+        if(sina > aSinCoeff)
+        {
+          theInds.InsertBefore(theI, mid);
+          return Standard_True;
+        }
+      }
+      else
+      {
+        theInds.InsertBefore(theI, mid);
+        return Standard_True;
+      }
+    }
+  }
+
+  return Standard_False;
+}
+
+//=======================================================================
+//function : BuildKnots
+//purpose  :
+//=======================================================================
+void ApproxInt_KnotTools::BuildKnots(const TColgp_Array1OfPnt& thePntsXYZ,
+                                     const TColgp_Array1OfPnt2d& thePntsU1V1,
+                                     const TColgp_Array1OfPnt2d& thePntsU2V2,
+                                     const math_Vector& thePars,
+                                     const Standard_Boolean theApproxXYZ,
+                                     const Standard_Boolean theApproxU1V1,
+                                     const Standard_Boolean theApproxU2V2,
+                                     const Standard_Integer theMinNbPnts,
+                                     NCollection_Vector<Standard_Integer>& theKnots)
+{
+  NCollection_Sequence<Standard_Integer> aKnots;
+  Standard_Integer aDim = 0;
+
+  // I: Convert input data to the corresponding format.
+  if(theApproxXYZ)
+    aDim += 3;
+  if(theApproxU1V1)
+    aDim += 2;
+  if(theApproxU2V2)
+    aDim += 2;
+
+  NCollection_LocalArray<Standard_Real> aCoords(thePars.Length()*aDim);
+  Standard_Integer i, j;
+  for(i = thePars.Lower(); i <= thePars.Upper(); ++i)
+  {
+    j = (i - thePars.Lower()) * aDim;
+    if(theApproxXYZ)
+    {
+      aCoords[j] = thePntsXYZ.Value(i).X();
+      ++j;
+      aCoords[j] = thePntsXYZ.Value(i).Y();
+      ++j;
+      aCoords[j] = thePntsXYZ.Value(i).Z();
+      ++j;
+    }
+    if(theApproxU1V1)
+    {
+      aCoords[j] = thePntsU1V1.Value(i).X();
+      ++j;
+      aCoords[j] = thePntsU1V1.Value(i).Y();
+      ++j;
+    }
+    if(theApproxU2V2)
+    {
+      aCoords[j] = thePntsU2V2.Value(i).X();
+      ++j;
+      aCoords[j] = thePntsU2V2.Value(i).Y();
+      ++j;
+    }
+  }
+
+  // II: Build draft knot sequence.
+  ComputeKnotInds(aCoords, aDim, thePars, aKnots);
+
+#if defined(APPROXINT_KNOTTOOLS_DEBUG) || defined(OCCT_DEBUG)
+    cout << "Draft knot sequence: " << endl;
+    for(i = aKnots.Lower(); i <= aKnots.Upper();  ++i)
+    {
+      cout << i << " : " << aKnots(i) << endl;
+    }
+#endif
+
+  // III: Build output knot sequence.
+  FilterKnots(aKnots, theMinNbPnts, theKnots);
+
+#if defined(APPROXINT_KNOTTOOLS_DEBUG) || defined(OCCT_DEBUG)
+    cout << "Result knot sequence: " << endl;
+    for(i = theKnots.Lower(); i <= theKnots.Upper();  ++i)
+    {
+      cout << i << " : " << theKnots(i) << endl;
+    }
+#endif
+
+}
diff --git a/src/ApproxInt/ApproxInt_KnotTools.hxx b/src/ApproxInt/ApproxInt_KnotTools.hxx
new file mode 100644 (file)
index 0000000..d80b1aa
--- /dev/null
@@ -0,0 +1,132 @@
+// Copyright (c) 1999-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#ifndef _ApproxInt_KnotTools_HeaderFile
+#define _ApproxInt_KnotTools_HeaderFile
+
+#ifndef _Standard_HeaderFile
+#include <Standard.hxx>
+#endif
+#ifndef _Standard_DefineAlloc_HeaderFile
+#include <Standard_DefineAlloc.hxx>
+#endif
+#ifndef _Standard_Macro_HeaderFile
+#include <Standard_Macro.hxx>
+#endif
+
+#ifndef _Standard_Boolean_HeaderFile
+#include <Standard_Boolean.hxx>
+#endif
+#ifndef _Standard_Real_HeaderFile
+#include <Standard_Real.hxx>
+#endif
+#ifndef _Standard_Integer_HeaderFile
+#include <Standard_Integer.hxx>
+#endif
+
+#include <TColgp_Array1OfPnt2d.hxx>
+#include <TColgp_Array1OfPnt.hxx>
+#include <TColStd_Array1OfReal.hxx>
+#include <NCollection_LocalArray.hxx>
+
+class math_Vector;
+template <class A> class NCollection_Sequence;
+template <class A> class NCollection_List;
+template <class A> class NCollection_Vector;
+
+// Corresponds for debug information output.
+// Debug information is also printed when OCCT_DEBUG defined.
+//#define APPROXINT_KNOTTOOLS_DEBUG
+
+//! This class intended to build knots sequence on discrete set of points for further approximation into bspline curve.
+//!
+//! Short description of algorithm:
+//! 1) Build discrete curvature on points set.
+//! 2) According to special rules build draft knots sequence.
+//! 3) Filter draft sequence to build output sequence.
+//!
+//! For more details look at:
+//! Anshuman Razdan - Knot Placement for B-Spline curve Approximation.
+class ApproxInt_KnotTools
+{
+public:
+
+  DEFINE_STANDARD_ALLOC
+
+  //! Main function to build optimal knot sequence.
+  //! At least one set from (thePntsXYZ, thePntsU1V1, thePntsU2V2) should exist.
+  //! @param thePntsXYZ - Set of 3d points.
+  //! @param thePntsU1V1 - Set of 2d points.
+  //! @param thePntsU2V2 - Set of 2d points.
+  //! @param thePars - Expected parameters assoiated with set.
+  //! @param theApproxXYZ - Flag, existence of 3d set.
+  //! @param theApproxU1V1 - Flag existence of first 2d set.
+  //! @param theApproxU2V2 - Flag existence of second 2d set.
+  //! @param theMinNbPnts - Minimal number of points per knot interval.
+  //! @param theKnots - output knots sequence.
+  Standard_EXPORT static void BuildKnots(const TColgp_Array1OfPnt& thePntsXYZ,
+                                         const TColgp_Array1OfPnt2d& thePntsU1V1,
+                                         const TColgp_Array1OfPnt2d& thePntsU2V2,
+                                         const math_Vector& thePars,
+                                         const Standard_Boolean theApproxXYZ,
+                                         const Standard_Boolean theApproxU1V1,
+                                         const Standard_Boolean theApproxU2V2,
+                                         const Standard_Integer theMinNbPnts,
+                                         NCollection_Vector<Standard_Integer>& theKnots);
+
+private:
+
+  //! Compute indices of knots:
+  //!
+  //! I: Build discrete curvature in points set,
+  //! using outer product of two vectors.
+  //!
+  //! II: Put knots in points which has extremity on discrete curvature.
+  //!
+  //! III: Put knots in monotone intervals of curvature.
+  //!
+  //! IV: Put additional knots near extrema points.
+  static void ComputeKnotInds(const NCollection_LocalArray<Standard_Real>& theCoords,
+                              const Standard_Integer theDim,
+                              const math_Vector& thePars,
+                              NCollection_Sequence<Standard_Integer>& theInds);
+
+  //! Insert knots before index I.
+  //!
+  //! I: Check curvature change:
+  //! if ( maxCurvature / minCurvature ) of current interval greater than 
+  //! threshold value, then stop and use upper index as knot.
+  //!
+  //! II: Check midpoint criteria:
+  //! If exist point between two knot indices with angle greater than
+  //! threshold value, then stop and put this index as knot.
+  static Standard_Boolean InsKnotBefI(const Standard_Integer theI,
+                                      const TColStd_Array1OfReal& theCurv,
+                                      const NCollection_LocalArray<Standard_Real>& theCoords,
+                                      const Standard_Integer theDim, 
+                                      NCollection_Sequence<Standard_Integer>& theInds,
+                                      const Standard_Boolean ChkCurv);
+
+  //! Perform knots filtration.
+  //!
+  //! I: Filter too big number of points per knot interval.
+  //!
+  //! II: Filter poins with too small amount of points per knot interval.
+  //!
+  //! III: Fill Last Knot.
+  static void FilterKnots(NCollection_Sequence<Standard_Integer>& theInds, 
+                          const Standard_Integer theMinNbPnts,
+                          NCollection_Vector<Standard_Integer>& theLKnots);
+};
+
+#endif
index f7d483d..55ab458 100755 (executable)
@@ -2,6 +2,8 @@ ApproxInt_Approx.gxx
 ApproxInt_ImpPrmSvSurfaces.gxx
 ApproxInt_MultiLine.gxx
 ApproxInt_MultiLineTool.gxx
+ApproxInt_KnotTools.hxx
+ApproxInt_KnotTools.cxx
 ApproxInt_MultiLineTool.lxx
 ApproxInt_PrmPrmSvSurfaces.gxx
 ApproxInt_SvSurfaces.cxx
index c3c5bbd..c00efbe 100644 (file)
@@ -21,6 +21,7 @@
 #include <Standard_DefineAlloc.hxx>
 #include <Standard_Handle.hxx>
 
+#include <NCollection_Vector.hxx>
 #include <BRepApprox_TheComputeLineOfApprox.hxx>
 #include <BRepApprox_TheComputeLineBezierOfApprox.hxx>
 #include <Approx_MCurvesToBSpCurve.hxx>
@@ -48,6 +49,22 @@ class BRepApprox_TheComputeLineBezierOfApprox;
 class BRepApprox_MyGradientOfTheComputeLineBezierOfApprox;
 class AppParCurves_MultiBSpCurve;
 
+struct Approx_Data 
+{
+  Approx_Data()
+  {
+    myMinFactorXYZ = 0.0;
+    myMinFactorUV  = 0.0;
+  }
+
+  Standard_Boolean myBezierApprox;
+  Standard_Real  Xo, Ax, Yo, Ay, Zo, Az,
+    U1o, A1u, V1o, A1v, U2o, A2u, V2o, A2v;
+  Standard_Boolean ApproxXYZ, ApproxU1V1, ApproxU2V2;
+  Standard_Integer indicemin, indicemax, nbpntmax;
+  Approx_ParametrizationType parametrization;
+  Standard_Real myMinFactorXYZ, myMinFactorUV;
+};
 
 
 class BRepApprox_Approx 
@@ -100,12 +117,31 @@ private:
   
   Standard_EXPORT void UpdateTolReached();
 
+  //! Fill data structure for intersection approximation.
+  Standard_EXPORT void fillData(const Handle(BRepApprox_ApproxLine)& theLine,
+                                const Standard_Boolean theApproxXYZ,
+                                const Standard_Boolean theApproxU1V1,
+                                const Standard_Boolean theApproxU2V2);
+
+  //! Prepare data structure for further computations.
+  Standard_EXPORT void prepareDS(const Standard_Boolean theApproxXYZ,
+                                 const Standard_Boolean theApproxU1V1,
+                                 const Standard_Boolean theApproxU2V2,
+                                 const Standard_Integer indicemin,
+                                 const Standard_Integer indicemax);
+
+  //! Build knot sequence.
+  Standard_EXPORT void buildKnots(const Handle(BRepApprox_ApproxLine)& theline,
+                                  const Standard_Address thePtrSVSurf);
+
+  //! Build curve.
+  Standard_EXPORT void buildCurve(const Handle(BRepApprox_ApproxLine)& theline,
+                                  const Standard_Address thePtrSVSurf);
 
   BRepApprox_TheComputeLineOfApprox myComputeLine;
   BRepApprox_TheComputeLineBezierOfApprox myComputeLineBezier;
   Approx_MCurvesToBSpCurve myBezToBSpl;
   Standard_Boolean myTolReached;
-  Standard_Boolean myApproxBez;
   Standard_Boolean myWithTangency;
   Standard_Real myTol3d;
   Standard_Real myTol2d;
@@ -114,11 +150,11 @@ private:
   Standard_Integer myDegMax;
   Standard_Integer myNbPntMax;
   Standard_Integer myNbIterMax;
-  Standard_Real myMinFactorXYZ;
-  Standard_Real myMinFactorUV;
   Standard_Real myTolReached3d;
   Standard_Real myTolReached2d;
-
+  Approx_Data myData;
+  Standard_Real myUVRes1, myUVRes2;
+  NCollection_Vector<Standard_Integer> myKnots;
 
 };
 
index b86aadc..ba734a4 100644 (file)
@@ -21,6 +21,7 @@
 #include <Standard_DefineAlloc.hxx>
 #include <Standard_Handle.hxx>
 
+#include <NCollection_Vector.hxx>
 #include <GeomInt_TheComputeLineOfWLApprox.hxx>
 #include <GeomInt_TheComputeLineBezierOfWLApprox.hxx>
 #include <Approx_MCurvesToBSpCurve.hxx>
@@ -48,6 +49,22 @@ class GeomInt_TheComputeLineBezierOfWLApprox;
 class GeomInt_MyGradientOfTheComputeLineBezierOfWLApprox;
 class AppParCurves_MultiBSpCurve;
 
+struct Approx_Data 
+{
+  Approx_Data()
+  {
+    myMinFactorXYZ = 0.0;
+    myMinFactorUV  = 0.0;
+  }
+
+  Standard_Boolean myBezierApprox;
+  Standard_Real  Xo, Ax, Yo, Ay, Zo, Az,
+    U1o, A1u, V1o, A1v, U2o, A2u, V2o, A2v;
+  Standard_Boolean ApproxXYZ, ApproxU1V1, ApproxU2V2;
+  Standard_Integer indicemin, indicemax, nbpntmax;
+  Approx_ParametrizationType parametrization;
+  Standard_Real myMinFactorXYZ, myMinFactorUV;
+};
 
 
 class GeomInt_WLApprox 
@@ -100,12 +117,31 @@ private:
   
   Standard_EXPORT void UpdateTolReached();
 
+  //! Fill data structure for intersection approximation.
+  Standard_EXPORT void fillData(const Handle(IntPatch_WLine)& theLine,
+                                const Standard_Boolean theApproxXYZ,
+                                const Standard_Boolean theApproxU1V1,
+                                const Standard_Boolean theApproxU2V2);
+
+  //! Prepare data structure for further computations.
+  Standard_EXPORT void prepareDS(const Standard_Boolean theApproxXYZ,
+                                 const Standard_Boolean theApproxU1V1,
+                                 const Standard_Boolean theApproxU2V2,
+                                 const Standard_Integer indicemin,
+                                 const Standard_Integer indicemax);
+
+  //! Build knot sequence.
+  Standard_EXPORT void buildKnots(const Handle(IntPatch_WLine)& theline,
+                                  const Standard_Address thePtrSVSurf);
+
+  //! Build curve.
+  Standard_EXPORT void buildCurve(const Handle(IntPatch_WLine)& theline,
+                                  const Standard_Address thePtrSVSurf);
 
   GeomInt_TheComputeLineOfWLApprox myComputeLine;
   GeomInt_TheComputeLineBezierOfWLApprox myComputeLineBezier;
   Approx_MCurvesToBSpCurve myBezToBSpl;
   Standard_Boolean myTolReached;
-  Standard_Boolean myApproxBez;
   Standard_Boolean myWithTangency;
   Standard_Real myTol3d;
   Standard_Real myTol2d;
@@ -114,11 +150,11 @@ private:
   Standard_Integer myDegMax;
   Standard_Integer myNbPntMax;
   Standard_Integer myNbIterMax;
-  Standard_Real myMinFactorXYZ;
-  Standard_Real myMinFactorUV;
   Standard_Real myTolReached3d;
   Standard_Real myTolReached2d;
-
+  Approx_Data myData;
+  Standard_Real myUVRes1, myUVRes2;
+  NCollection_Vector<Standard_Integer> myKnots;
 
 };