0025841: Incorrect edge displaying
[occt.git] / src / GCPnts / GCPnts_TangentialDeflection.gxx
index 6670c7c..54c0b3b 100644 (file)
@@ -5,8 +5,8 @@
 //
 // 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 version 2.1 as published
+// 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.
@@ -39,7 +39,7 @@ void GCPnts_TangentialDeflection::EvaluateDu (
     Standard_Real Lc = N.CrossMagnitude (T);
     Standard_Real Ln = Lc/Lt;
     if (Ln > LTol) {
-      Du = sqrt (8.0 * curvatureDeflection / Ln);
+      Du = sqrt (8.0 * Max(curvatureDeflection, myMinLen) / Ln);
       NotDone = Standard_False;
     }
   }
@@ -56,10 +56,11 @@ GCPnts_TangentialDeflection::GCPnts_TangentialDeflection (
  const Standard_Real    AngularDeflection,
  const Standard_Real    CurvatureDeflection,
  const Standard_Integer MinimumOfPoints,
- const Standard_Real    UTol)
+ const Standard_Real    UTol,
+ const Standard_Real    theMinLen)
 
 { 
-  Initialize (C,AngularDeflection,CurvatureDeflection,MinimumOfPoints,UTol); 
+  Initialize (C,AngularDeflection,CurvatureDeflection,MinimumOfPoints,UTol,theMinLen); 
 }
 
 
@@ -75,7 +76,8 @@ GCPnts_TangentialDeflection::GCPnts_TangentialDeflection (
  const Standard_Real    AngularDeflection,
  const Standard_Real    CurvatureDeflection,
  const Standard_Integer MinimumOfPoints,
- const Standard_Real    UTol)
+ const Standard_Real    UTol,
+ const Standard_Real    theMinLen)
 
 { 
   Initialize (C, 
@@ -84,7 +86,7 @@ GCPnts_TangentialDeflection::GCPnts_TangentialDeflection (
         AngularDeflection, 
         CurvatureDeflection, 
         MinimumOfPoints,
-        UTol);
+        UTol, theMinLen);
 }
 
 
@@ -99,7 +101,8 @@ void GCPnts_TangentialDeflection::Initialize (
  const Standard_Real    AngularDeflection, 
  const Standard_Real    CurvatureDeflection,
  const Standard_Integer MinimumOfPoints,
- const Standard_Real    UTol)
+ const Standard_Real    UTol,
+ const Standard_Real    theMinLen)
 
 {
   Initialize (C, 
@@ -108,7 +111,7 @@ void GCPnts_TangentialDeflection::Initialize (
               AngularDeflection, 
               CurvatureDeflection,
               MinimumOfPoints,
-              UTol);
+              UTol, theMinLen);
 }
 
 
@@ -124,7 +127,8 @@ void GCPnts_TangentialDeflection::Initialize (
  const Standard_Real    AngularDeflection, 
  const Standard_Real    CurvatureDeflection,
  const Standard_Integer MinimumOfPoints,
- const Standard_Real    UTol)
+ const Standard_Real    UTol,
+ const Standard_Real    theMinLen)
 
 {
   
@@ -144,6 +148,7 @@ void GCPnts_TangentialDeflection::Initialize (
  angularDeflection   = AngularDeflection;
  curvatureDeflection = CurvatureDeflection;
  minNbPnts           = Max (MinimumOfPoints, 2);
+ myMinLen            = Max(theMinLen, Precision::Confusion());
 
  switch (C.GetType()) {
 
@@ -201,7 +206,6 @@ void GCPnts_TangentialDeflection::PerformLinear (const TheCurve& C) {
   points    .Append (P);
 }
 
-
 //=======================================================================
 //function : PerformCircular
 //purpose  : 
@@ -211,18 +215,18 @@ void GCPnts_TangentialDeflection::PerformCircular (const TheCurve& C)
 {
   // akm 8/01/02 : check the radius before divide by it
   Standard_Real dfR = C.Circle().Radius();
-  Standard_Real Du = 0.;
-  if (Abs(dfR) > Precision::Confusion())
-    Du = Max(1.0e0 - (curvatureDeflection/dfR),0.0e0) ;
-  Du  = acos (Du); Du+=Du;
-  Du               = Min (Du, angularDeflection);
-  Standard_Integer NbPoints = (Standard_Integer )((lastu - firstu) / Du);
-  NbPoints         = Max (NbPoints, minNbPnts-1);
-  Du               = (lastu - firstu) / NbPoints;
+  Standard_Real Du = GCPnts_TangentialDeflection::ArcAngularStep(
+    dfR, curvatureDeflection, angularDeflection, myMinLen);
+    
+  const Standard_Real aDiff = lastu - firstu;
+  Standard_Integer NbPoints = (Standard_Integer)(aDiff / Du);
+  NbPoints = Max(NbPoints, minNbPnts - 1);
+  Du       = aDiff / NbPoints;
 
   gp_Pnt P;
   Standard_Real U = firstu;
-  for (Standard_Integer i = 1; i <= NbPoints; i++) {
+  for (Standard_Integer i = 1; i <= NbPoints; i++)
+  {
     D0 (C, U, P);
     parameters.Append (U);
     points    .Append (P);
@@ -234,17 +238,15 @@ void GCPnts_TangentialDeflection::PerformCircular (const TheCurve& C)
 }
 
 
-
 //=======================================================================
 //function : PerformCurve
 //purpose  : On respecte ll'angle et la fleche, on peut imposer un nombre
 //           minimum de points sur un element lineaire
 //=======================================================================
-
 void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
 
 {
-  Standard_Integer i;       
+  Standard_Integer i, j;       
   gp_XYZ V1, V2;
   gp_Pnt MiddlePoint, CurrentPoint, LastPoint;   
   Standard_Real Du, Dusave, MiddleU, L1, L2;
@@ -264,62 +266,81 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
   parameters.Append (U1);
   points    .Append (CurrentPoint);
 
+  // Used to detect "isLine" current bspline and in Du computation in general handling.
+  Standard_Integer NbInterv = C.NbIntervals(GeomAbs_CN);
+  TColStd_Array1OfReal Intervs(1, NbInterv+1);
+  C.Intervals(Intervs, GeomAbs_CN);
+
   if (NotDone) {
     //C'est soit une droite, soit une singularite :
-    V1 = LastPoint.XYZ ();
-    V1.Subtract (CurrentPoint.XYZ());
+    V1 = (LastPoint.XYZ() - CurrentPoint.XYZ());
     L1 = V1.Modulus ();
-    if (L1 > LTol) {
+    if (L1 > LTol)
+    {
       //Si c'est une droite on verifie en calculant minNbPoints :
       Standard_Boolean IsLine   = Standard_True;
-      Standard_Integer NbPoints = 3;
-      if (minNbPnts > 3) NbPoints = minNbPnts;
-      Du = (lastu-firstu)/NbPoints;
-      MiddleU = firstu + Du;
-      for (i = 2; i < NbPoints; i++) {
-        D0 (C, MiddleU, MiddlePoint);
-        V2 = MiddlePoint.XYZ();
-        V2.Subtract (CurrentPoint.XYZ());
-        L2 = V2.Modulus ();
-        if (L2 > LTol) {
-          if (((V2.CrossMagnitude (V1))/(L1*L2)) >= ATol) {
-            //C'etait une singularite
-            IsLine = Standard_False;
-            break;
-          }
-          if (minNbPnts > 2) {
-            parameters.Append (MiddleU);
-            points    .Append (MiddlePoint);
+      Standard_Integer NbPoints = (minNbPnts > 3) ? minNbPnts : 3;
+      ////
+      Standard_Real param = 0.;
+      for (i = 1; i <= NbInterv && IsLine; ++i)
+      {
+        // Avoid usage intervals out of [firstu, lastu].
+        if ((Intervs(i+1) < firstu) ||
+            (Intervs(i)   > lastu))
+        {
+          continue;
+        }
+        // Fix border points in applicable intervals, to avoid be out of target interval.
+        if ((Intervs(i)   < firstu) &&
+            (Intervs(i+1) > firstu))
+        {
+          Intervs(i) = firstu;
+        }
+        if ((Intervs(i)   < lastu) &&
+            (Intervs(i+1) > lastu))
+        {
+          Intervs(i + 1) = lastu;
+        }
+
+        Standard_Real delta = (Intervs(i+1) - Intervs(i))/NbPoints;
+        for (j = 1; j <= NbPoints && IsLine; ++j)
+        {
+          param = Intervs(i) + j*delta;
+          D0 (C, param, MiddlePoint);
+          V2 = (MiddlePoint.XYZ() - CurrentPoint.XYZ());
+          L2 = V2.Modulus ();
+          if (L2 > LTol)
+          {
+            const Standard_Real aAngle = V2.CrossMagnitude(V1)/(L1*L2);
+            IsLine = (aAngle < ATol);
           }
         }
-        MiddleU += Du;
       }
-      if (IsLine) {
-        //C'etait une droite (plusieurs poles alignes), Calcul termine :
-        parameters.Append (lastu);
-        points    .Append (LastPoint);
+
+      if (IsLine)
+      {
+        parameters.Clear();
+        points    .Clear();
+
+        PerformLinear(C);
         return;
       }
-      else {
+      else
+      {
         //c'etait une singularite on continue :
-        Standard_Integer pointsLength=points.Length ();
-        for (i = 2; i <= pointsLength; i++) {
-          points    .Remove (i);
-          parameters.Remove (i);
-          pointsLength--;
-        }
-        Du = Dusave;
+        //Du = Dusave;
+        EvaluateDu (C, param, MiddlePoint, Du, NotDone);
       }
     }
-    else  {
-      
+    else
+    {
       Du = (lastu-firstu)/2.1;
       MiddleU = firstu + Du;
       D0 (C, MiddleU, MiddlePoint);
-      V1 = MiddlePoint.XYZ ();
-      V1.Subtract (CurrentPoint.XYZ());
+      V1 = (MiddlePoint.XYZ() - CurrentPoint.XYZ());
       L1 = V1.Modulus ();
-      if (L1 < LTol) {
+      if (L1 < LTol)
+      {
         // L1 < LTol C'est une courbe de longueur nulle, calcul termine :
         // on renvoi un segment de 2 points   (protection)
         parameters.Append (lastu);
@@ -345,12 +366,13 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
   Standard_Boolean MorePoints = Standard_True;
   Standard_Real U2            = firstu;   
   Standard_Real AngleMax      = angularDeflection * 0.5;  //car on prend le point milieu
-
+  Standard_Integer aIdx[2] = {Intervs.Lower(), Intervs.Lower()}; // Indexes of intervals of U1 and U2, used to handle non-uniform case.
+  Standard_Boolean isNeedToCheck = Standard_False;
   gp_Pnt aPrevPoint = points.Last();
 
   while (MorePoints) {
-
-    U2 += Du;                             
+    aIdx[0] = getIntervalIdx(U1, Intervs, aIdx[0]);
+    U2 += Du;
 
     if (U2 >= lastu) {                       //Bout de courbe
       U2 = lastu;
@@ -360,37 +382,62 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
     }
     else D0 (C, U2, CurrentPoint);           //Point suivant
 
-    Standard_Real Coef, ACoef = 0., FCoef = 0.;
+    Standard_Real Coef = 0.0, ACoef = 0., FCoef = 0.;
     Standard_Boolean Correction, TooLarge, TooSmall;
     TooLarge   = Standard_False;
-    TooSmall   = Standard_False;
     Correction = Standard_True;
+    TooSmall = Standard_False;
 
-    Standard_Real lastCoef = 0;
-     
     while (Correction) {                     //Ajustement Du
+      if (isNeedToCheck)
+      {
+        aIdx[1] = getIntervalIdx(U2, Intervs, aIdx[0]);
+        if (aIdx[1] > aIdx[0]) // Jump to another polynom.
+        {
+          if (Du > (Intervs(aIdx[0] + 1) - Intervs(aIdx[0]) ) * Us3) // Set Du to the smallest value and check deflection on it.
+          {
+            Du = (Intervs(aIdx[0] + 1) - Intervs(aIdx[0]) ) * Us3;
+            U2 = U1 + Du;
+            if (U2 > lastu)
+              U2 = lastu;
+            D0 (C, U2, CurrentPoint);
+          }
+        }
+      }
       MiddleU = (U1+U2)*0.5;                 //Verif / au point milieu
       D0 (C, MiddleU, MiddlePoint);
 
-      V1 = CurrentPoint.XYZ ();              //Critere de fleche
-      V1.Subtract (aPrevPoint.XYZ());
-      V2 = MiddlePoint.XYZ ();
-      V2.Subtract (aPrevPoint.XYZ());
+      V1 = (CurrentPoint.XYZ() - aPrevPoint.XYZ()); //Critere de fleche
+      V2 = (MiddlePoint.XYZ()  - aPrevPoint.XYZ());
       L1 = V1.Modulus ();
-      if (L1 > LTol) FCoef = V1.CrossMagnitude(V2)/(L1*curvatureDeflection);
-      else           FCoef = 0.0;
 
-      V1 = CurrentPoint.XYZ ();              //Critere d'angle
-      V1.Subtract (MiddlePoint.XYZ ());
+      FCoef = (L1 > myMinLen) ? 
+        V1.CrossMagnitude(V2)/(L1*curvatureDeflection) : 0.0;
+
+      V1 = (CurrentPoint.XYZ() - MiddlePoint.XYZ()); //Critere d'angle
       L1 = V1.Modulus ();
       L2 = V2.Modulus ();
-      Standard_Real angg = V1.CrossMagnitude(V2)/(L1*L2);
-      if (L1 > LTol && L2 > LTol) ACoef = angg/AngleMax;
-      else                        ACoef = 0.0;
+      if (L1 > myMinLen && L2 > myMinLen)
+      {
+        Standard_Real angg = V1.CrossMagnitude(V2) / (L1 * L2);
+        ACoef = angg / AngleMax;
+      }
+      else
+        ACoef = 0.0;
 
-      if (ACoef >= FCoef) Coef = ACoef;      //On retient le plus penalisant
-      else                Coef = FCoef;
+      //On retient le plus penalisant
+      Coef = Max(ACoef, FCoef);
 
+      if (isNeedToCheck && Coef < 0.55)
+      {
+        isNeedToCheck = Standard_False;
+        Du = Dusave;
+        U2 = U1 + Du;
+        if (U2 > lastu)
+          U2 = lastu;
+        D0 (C, U2, CurrentPoint);
+        continue;
+      }
 
       if (Coef <= 1.0) {
         if (Abs (lastu-U2) < uTol) {
@@ -405,6 +452,7 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
             points    .Append (CurrentPoint);
             aPrevPoint = CurrentPoint;
             Correction = Standard_False;
+            isNeedToCheck = Standard_True;
           }
           else if (TooSmall) {
             Correction = Standard_False;
@@ -412,25 +460,24 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
           }
           else {
             TooSmall = Standard_True;
-            lastCoef = Coef;
             //Standard_Real UUU2 = U2;
             Du += Min((U2-U1)*(1.-Coef), Du*Us3);
 
             U2 = U1 + Du;
-            //if (U2 >= lastu) U2 = UUU2;
-            if (U2 >= lastu) {
-              parameters.Append (lastu);
-              points    .Append (LastPoint);
-              MorePoints = Standard_False;
-              Correction = Standard_False;
-            }
-            else D0 (C, U2, CurrentPoint);
+            if (U2 > lastu)
+              U2 = lastu;
+            D0 (C, U2, CurrentPoint);
           }
         }
       }
       else {
 
         if (Coef >= 1.5) {
+          if (!aPrevPoint.IsEqual(points.Last(), Precision::Confusion()))
+          {
+            parameters.Append (U1);
+            points    .Append (aPrevPoint);
+          }
           U2 = MiddleU;
           Du  = U2-U1;
           CurrentPoint = MiddlePoint;
@@ -488,7 +535,6 @@ void GCPnts_TangentialDeflection::PerformCurve (const TheCurve& C)
 //    points.Remove (i+1);
 //    i--;
 //  }
-  
   if (i >= 2) {
     MiddleU = parameters (i-1);
     MiddleU = (lastu + MiddleU)*0.5;