0027383: Modeling - improve handling of regularity on edges

1. There has been implemented calculation of all possible types of continuity for shared edges:
  * G1 is set if tangential planes are the same for connected faces in each control points through the edge;
  * C1 is set in addition to G1 conditions if derivatives, orthogonal to the edge on each face, are equal vectors;
  * G2 is set in addition to G1 if the centers of principal curvatures are the same for connected faces in each control points through the edge;
  * C2 is set in addition to C1 and G2 if directions of principal curvatures are equal;
  * CN continuity is set only if both connected faces are based on elementary surfaces (the conditions for this case are similar to C2 continuity).

2. ShapeFix::EncodeRegularity() is merged into BRepLib::EncodeRegularity().
3. Implemented several test cases to check correct handling of regularity.
4. Fix incorrect usage of BRepLib::EncodeRegularity() in BRepBuilderAPI_Sewing.
5. Implement a method for calculation of regularity on the given list of edges.
6. Documentation updates

dox/user_guides/modeling_algos/modeling_algos.md

@@ -2937,7 +2937,7 @@ The algorithm of shape triangulation is provided by the functionality of *BRepMe
 ~~~~~
 const Standard_Real aRadius = 10.0; 
 const Standard_Real aHeight = 25.0; 
-BRepBuilderAPI_MakeCylinder aCylinder(aRadius, aHeight); 
+BRepPrimAPI_MakeCylinder aCylinder(aRadius, aHeight); 
 TopoDS_Shape aShape = aCylinder.Shape();
  
 const Standard_Real aLinearDeflection   = 0.01;

dox/user_guides/modeling_data/modeling_data.md

@@ -603,6 +603,50 @@ To check the concavity of a surface, proceed as follows:
   2. If the value of the curvature changes of sign, the surface is concave or convex depending on the point of view.
   3. To compute a Gaussian curvature, use the class <i> SLprops</i> from <i> GeomLProp</i>, which instantiates the generic class <i> SLProps </i>from <i> LProp</i> and use the method <i> GaussianCurvature</i>.
 
+@subsection occt_modat_4_2a Continuity of Curves and Surfaces
+
+Types of supported continuities for curves and surfaces are described in *GeomAbs_Shape* enumeration.
+
+In respect of curves, the following types of continuity are supported (see the figure below):
+  * C0 (*GeomAbs_C0*) - parametric continuity. It is the same as G0 (geometric continuity), so the last one is not represented by separate variable.
+  * G1 (*GeomAbs_G1*) - tangent vectors on left and on right are parallel.
+  * C1 (*GeomAbs_C1*) - indicates the continuity of the first derivative.
+  * G2 (*GeomAbs_G2*) - in addition to G1 continuity, the centers of curvature on left and on right are the same.
+  * C2 (*GeomAbs_C2*) - continuity of all derivatives till the second order.
+  * C3 (*GeomAbs_C3*) - continuity of all derivatives till the third order.
+  * CN (*GeomAbs_CN*) - continuity of all derivatives till the N-th order (infinite order of continuity).
+
+*Note:* Geometric continuity (G1, G2) means that the curve can be reparametrized to have parametric (C1, C2) continuity.
+
+@image html /user_guides/modeling_data/images/modeling_data_continuity_curves.svg "Continuity of Curves"
+@image latex /user_guides/modeling_data/images/modeling_data_continuity_curves.svg "Continuity of Curves" width=\\textwidth
+
+The following types of surface continuity are supported:
+  * C0 (*GeomAbs_C0*) - parametric continuity (the surface has no points or curves of discontinuity).
+  * G1 (*GeomAbs_G1*) - surface has single tangent plane in each point.
+  * C1 (*GeomAbs_C1*) - indicates the continuity of the first derivatives.
+  * G2 (*GeomAbs_G2*) - in addition to G1 continuity, principal curvatures and directions are continuous.
+  * C2 (*GeomAbs_C2*) - continuity of all derivatives till the second order.
+  * C3 (*GeomAbs_C3*) - continuity of all derivatives till the third order.
+  * CN (*GeomAbs_CN*) - continuity of all derivatives till the N-th order (infinite order of continuity).
+
+@image html /user_guides/modeling_data/images/modeling_data_continuity_surfaces.svg "Continuity of Surfaces"
+@image latex /user_guides/modeling_data/images/modeling_data_continuity_surfaces.svg "Continuity of Surfaces" width=\\textwidth
+
+Against single surface, the connection of two surfaces (see the figure above) defines its continuity in each intersection point only. Smoothness of connection is a minimal value of continuities on the intersection curve.
+
+
+@subsection occt_modat_4_2b Regularity of Shared Edges
+
+Regularity of an edge is a smoothness of connection of two faces sharing this edge. In other words, regularity is a minimal continuity between connected faces in each point on edge.
+
+Edge's regularity can be set by *BRep_Builder::Continuity* method. To get the regularity use *BRep_Tool::Continuity* method.
+
+Some algorithms like @ref occt_modalg_6 "Fillet" set regularity of produced edges by their own algorithms. On the other hand, some other algorithms (like @ref occt_user_guides__boolean_operations "Boolean Operations", @ref occt_user_guides__shape_healing "Shape Healing", etc.) do not set regularity. If the regularity is needed to be set correctly on a shape, the method *BRepLib::EncodeRegularity* can be used. It calculates and sets correct values for all edges of the shape.
+
+The regularity flag is extensively used by the following high level algorithms: @ref occt_modalg_6_1_2 "Chamfer", @ref occt_modalg_7_3 "Draft Angle", @ref occt_modalg_10 "Hidden Line Removal", @ref occt_modalg_9_2_3 "Gluer".
+
+  
 @subsection occt_modat_4_3 Global Properties of Shapes
 
 The Global Properties of Shapes component provides algorithms for computing the global 
@@ -1233,7 +1277,7 @@ Below is the auxiliary function, which copies the element of rank *i* from the m
 For example, in the wire in the image we want to recuperate the edges in the order {e1, e2, e3,e4, e5} :
 
 @image html /user_guides/modeling_data/images/modeling_data_image014.png "A wire composed of 6 edges."
-@image latex /user_guides/modeling_data/images/modeling_data_image014.png "A wire composed of 6 edges.
+@image latex /user_guides/modeling_data/images/modeling_data_image014.png "A wire composed of 6 edges."
 
 *TopExp_Explorer*, however, recuperates the lines in any order.
  

src/BRepBuilderAPI/BRepBuilderAPI_Sewing.cxx

@@ -890,7 +890,6 @@ TopoDS_Edge BRepBuilderAPI_Sewing::SameParameterEdge(const TopoDS_Edge& edgeFirs
       else {
         aBuilder.UpdateEdge(edge, c2d2, surf2, loc2, Precision::Confusion());
       }
-
     }
   }
   Standard_Real tolReached = Precision::Infinite();
@@ -976,7 +975,6 @@ TopoDS_Edge BRepBuilderAPI_Sewing::SameParameterEdge(const TopoDS_Edge& edgeFirs
     }
   }
 
-  BRepLib::EncodeRegularity(edge,0.01);
   Standard_Real tolEdge1 = BRep_Tool::Tolerance(edge);
   if (tolEdge1 > MaxTolerance()) edge.Nullify();
   return edge;
@@ -1923,6 +1921,9 @@ void BRepBuilderAPI_Sewing::Perform(const Handle(Message_ProgressIndicator)& the
         mySewedShape.Nullify();
         return;
       }
+
+      EdgeRegularity (thePI);
+
       if (mySameParameterMode && myFaceMode)
         SameParameterShape();
       if (!aPS.More())
@@ -4030,6 +4031,28 @@ void BRepBuilderAPI_Sewing::EdgeProcessing(const Handle(Message_ProgressIndicato
       }
     }
   }
+}
+
+//=======================================================================
+//function : EdgeRegularity
+//purpose  : update Continuity flag on newly created edges
+//=======================================================================
+
+void BRepBuilderAPI_Sewing::EdgeRegularity(const Handle(Message_ProgressIndicator)& thePI)
+{
+  TopTools_IndexedDataMapOfShapeListOfShape aMapEF;
+  TopExp::MapShapesAndAncestors(mySewedShape, TopAbs_EDGE, TopAbs_FACE, aMapEF);
+
+  Message_ProgressSentry aPS(thePI, "Encode edge regularity", 0, myMergedEdges.Extent(), 1);
+  for (TopTools_MapIteratorOfMapOfShape aMEIt(myMergedEdges); aMEIt.More() && aPS.More(); aMEIt.Next(), aPS.Next())
+  {
+    TopoDS_Edge anEdge = TopoDS::Edge(myReShape->Apply(aMEIt.Value()));
+    const TopTools_ListOfShape* aFaces = aMapEF.Seek(anEdge);
+    // encode regularity if and only if edges is shared by two faces
+    if (aFaces && aFaces->Extent() == 2)
+      BRepLib::EncodeRegularity(anEdge, TopoDS::Face(aFaces->First()), TopoDS::Face(aFaces->Last()));
+  }
+
   myMergedEdges.Clear();
 }
 

src/BRepBuilderAPI/BRepBuilderAPI_Sewing.hxx

@@ -264,6 +264,9 @@ protected:
   
   Standard_EXPORT void EdgeProcessing (const Handle(Message_ProgressIndicator)& thePI = 0);
 
+  //! Recompute regularity on merged edges
+  Standard_EXPORT void EdgeRegularity (const Handle(Message_ProgressIndicator)& thePI = 0);
+  
   Standard_EXPORT void CreateOutputInformations();
   
   //! Defines if surface is U closed.

src/BRepLib/BRepLib.cxx

@@ -71,6 +71,7 @@
 #include <Poly_Triangulation.hxx>
 #include <Precision.hxx>
 #include <ProjLib_ProjectedCurve.hxx>
+#include <Standard_ErrorHandler.hxx>
 #include <Standard_Real.hxx>
 #include <TColgp_Array1OfPnt.hxx>
 #include <TColgp_Array1OfPnt2d.hxx>
@@ -1656,144 +1657,316 @@ Standard_Boolean BRepLib::OrientClosedSolid(TopoDS_Solid& solid)
 
   return Standard_True;
 }
+
+// Structure for calculation of properties, necessary for decision about continuity
+class SurfaceProperties
+{
+public:
+  SurfaceProperties(const Handle(Geom_Surface)& theSurface,
+                    const gp_Trsf&              theSurfaceTrsf,
+                    const Handle(Geom2d_Curve)& theCurve2D,
+                    const Standard_Boolean      theReversed)
+    : mySurfaceProps(theSurface, 2, Precision::Confusion()),
+      mySurfaceTrsf(theSurfaceTrsf),
+      myCurve2d(theCurve2D),
+      myIsReversed(theReversed)
+  {}
+
+  // Calculate derivatives on surface related to the point on curve
+  void Calculate(const Standard_Real theParamOnCurve)
+  {
+    gp_Pnt2d aUV;
+    myCurve2d->D1(theParamOnCurve, aUV, myCurveTangent);
+    mySurfaceProps.SetParameters(aUV.X(), aUV.Y());
+  }
+
+  // Returns point just calculated
+  gp_Pnt Value() 
+  { return mySurfaceProps.Value().Transformed(mySurfaceTrsf); }
+
+  // Calculate a derivative orthogonal to curve's tangent vector
+  gp_Vec Derivative()
+  {
+    gp_Vec aDeriv;
+    // direction orthogonal to tangent vector of the curve
+    gp_Vec2d anOrtho(-myCurveTangent.Y(), myCurveTangent.X());
+    Standard_Real aLen = anOrtho.Magnitude();
+    if (aLen < Precision::Confusion())
+      return aDeriv;
+    anOrtho /= aLen;
+    if (myIsReversed)
+      anOrtho.Reverse();
+
+    aDeriv.SetLinearForm(anOrtho.X(), mySurfaceProps.D1U(),
+                         anOrtho.Y(), mySurfaceProps.D1V());
+    return aDeriv.Transformed(mySurfaceTrsf);
+  }
+
+  // Calculate principal curvatures, which consist of minimal and maximal normal curvatures and
+  // the directions on the tangent plane (principal direction) where the extremums are reached
+  void Curvature(gp_Dir& thePrincipalDir1, Standard_Real& theCurvature1,
+                 gp_Dir& thePrincipalDir2, Standard_Real& theCurvature2)
+  {
+    mySurfaceProps.CurvatureDirections(thePrincipalDir1, thePrincipalDir2);
+    theCurvature1 = mySurfaceProps.MaxCurvature();
+    theCurvature2 = mySurfaceProps.MinCurvature();
+    if (myIsReversed)
+    {
+      theCurvature1 = -theCurvature1;
+      theCurvature2 = -theCurvature2;
+    }
+    thePrincipalDir1.Transform(mySurfaceTrsf);
+    thePrincipalDir2.Transform(mySurfaceTrsf);
+  }
+
+private:
+  GeomLProp_SLProps    mySurfaceProps; // properties calculator
+  gp_Trsf              mySurfaceTrsf;
+  Handle(Geom2d_Curve) myCurve2d;
+  Standard_Boolean     myIsReversed; // the face based on the surface is reversed
+
+  // tangent vector to Pcurve in UV
+  gp_Vec2d myCurveTangent;
+};
+
 //=======================================================================
 //function : tgtfaces
 //purpose  : check the angle at the border between two squares.
 //           Two shares should have a shared front edge.
 //=======================================================================
 static GeomAbs_Shape tgtfaces(const TopoDS_Edge& Ed,
-  const TopoDS_Face& F1,
+                              const TopoDS_Face& F1,
-  const TopoDS_Face& F2,
+                              const TopoDS_Face& F2,
-  const Standard_Real theAngleTol,
+                              const Standard_Real theAngleTol)
-  const Standard_Boolean couture)
 {
+  Standard_Boolean isSeam = F1.IsEqual(F2);
+
+  TopoDS_Edge E = Ed;
+
   // Check if pcurves exist on both faces of edge
   Standard_Real aFirst,aLast;
-  Handle(Geom2d_Curve) aCurve;
+  E.Orientation(TopAbs_FORWARD);
-  aCurve = BRep_Tool::CurveOnSurface(Ed,F1,aFirst,aLast);
+  Handle(Geom2d_Curve) aCurve1 = BRep_Tool::CurveOnSurface(E, F1, aFirst, aLast);
-  if(aCurve.IsNull())
+  if(aCurve1.IsNull())
-    return GeomAbs_C0;
-  aCurve = BRep_Tool::CurveOnSurface(Ed,F2,aFirst,aLast);
-  if(aCurve.IsNull())
     return GeomAbs_C0;
 
-  Standard_Real u;
+  if (isSeam)
-  TopoDS_Edge E = Ed;
+    E.Orientation(TopAbs_REVERSED);
-  BRepAdaptor_Surface aBAS1(F1,Standard_False);
+  Handle(Geom2d_Curve) aCurve2 = BRep_Tool::CurveOnSurface(E, F2, aFirst, aLast);
-  BRepAdaptor_Surface aBAS2(F2,Standard_False);
+  if(aCurve2.IsNull())
+    return GeomAbs_C0;
+
+  TopLoc_Location aLoc1, aLoc2;
+  Handle(Geom_Surface) aSurface1 = BRep_Tool::Surface(F1, aLoc1);
+  const gp_Trsf& aSurf1Trsf = aLoc1.Transformation();
+  Handle(Geom_Surface) aSurface2 = BRep_Tool::Surface(F2, aLoc2);
+  const gp_Trsf& aSurf2Trsf = aLoc2.Transformation();
+
+  if (aSurface1->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)))
+    aSurface1 = Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface1)->BasisSurface();
+  if (aSurface2->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)))
+    aSurface2 = Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface2)->BasisSurface();
 
   // seam edge on elementary surface is always CN
   Standard_Boolean isElementary =
-    (aBAS1.Surface().Surface()->IsKind(STANDARD_TYPE(Geom_ElementarySurface)) &&
+    (aSurface1->IsKind(STANDARD_TYPE(Geom_ElementarySurface)) &&
-     aBAS1.Surface().Surface()->IsKind(STANDARD_TYPE(Geom_ElementarySurface)));
+     aSurface2->IsKind(STANDARD_TYPE(Geom_ElementarySurface)));
-  if (couture && isElementary)
+  if (isSeam && isElementary)
   {
     return GeomAbs_CN;
   }
 
-  Handle(BRepAdaptor_HSurface) HS1 = new BRepAdaptor_HSurface (aBAS1);
+  SurfaceProperties aSP1(aSurface1, aSurf1Trsf, aCurve1, F1.Orientation() == TopAbs_REVERSED);
-  Handle(BRepAdaptor_HSurface) HS2;
+  SurfaceProperties aSP2(aSurface2, aSurf2Trsf, aCurve2, F2.Orientation() == TopAbs_REVERSED);
-  if(couture) HS2 = HS1;
-  else HS2 = new BRepAdaptor_HSurface(aBAS2);
-  //case when edge lies on the one face
 
-  E.Orientation(TopAbs_FORWARD);
+  Standard_Real f, l, eps;
-  Handle(BRepAdaptor_HCurve2d) HC2d1 = new BRepAdaptor_HCurve2d();
-  HC2d1->ChangeCurve2d().Initialize(E,F1);
-  if(couture) E.Orientation(TopAbs_REVERSED);
-  Handle(BRepAdaptor_HCurve2d) HC2d2 = new BRepAdaptor_HCurve2d();
-  HC2d2->ChangeCurve2d().Initialize(E,F2);
-  Adaptor3d_CurveOnSurface C1(HC2d1,HS1);
-  Adaptor3d_CurveOnSurface C2(HC2d2,HS2);
-
-  Standard_Boolean rev1 = (F1.Orientation() == TopAbs_REVERSED);
-  Standard_Boolean rev2 = (F2.Orientation() == TopAbs_REVERSED);
-  Standard_Real f,l,eps;
   BRep_Tool::Range(E,f,l);
   Extrema_LocateExtPC ext;
-  Standard_Boolean IsInitialized = Standard_False;
+  Handle(BRepAdaptor_HCurve) aHC2;
 
   eps = (l - f)/100.;
   f += eps; // to avoid calculations on  
   l -= eps; // points of pointed squares.
-  gp_Pnt2d p;
-  gp_Pnt pp1,pp2;//,PP;
-  gp_Vec du1, dv1, d2u1, d2v1, d2uv1;
-  gp_Vec du2, dv2, d2u2, d2v2, d2uv2;
-  gp_Vec d1,d2;
-  Standard_Real uu, vv, norm;
 
-  Standard_Integer i;
+  const Standard_Real anAngleTol2 = theAngleTol * theAngleTol;
+
+  gp_Vec aDer1, aDer2;
+  gp_Vec aNorm1;
+  Standard_Real aSqLen1, aSqLen2;
+  gp_Dir aCrvDir1[2], aCrvDir2[2];
+  Standard_Real aCrvLen1[2], aCrvLen2[2];
+
   GeomAbs_Shape aCont = (isElementary ? GeomAbs_CN : GeomAbs_C2);
-  for(i = 0; i<= 20 && aCont > GeomAbs_C0; i++)
+  GeomAbs_Shape aCurCont;
+  Standard_Real u;
+  for (Standard_Integer i = 0; i <= 20 && aCont > GeomAbs_C0; i++)
   {
     // First suppose that this is sameParameter
     u = f + (l-f)*i/20;
 
-    // take derivatives of surfaces at the same u, and compute normals
+    // Check conditions for G1 and C1 continuity:
-    HC2d1->D0(u,p);
+    // * calculate a derivative in tangent plane of each surface
-    HS1->D2 (p.X(), p.Y(), pp1, du1, dv1, d2u1, d2v1, d2uv1);
+    //   orthogonal to curve's tangent vector
-    d1 = (du1.Crossed(dv1));
+    // * continuity is C1 if the vectors are equal
-    norm = d1.Magnitude(); 
+    // * continuity is G1 if the vectors are just parallel
-    if (norm > 1.e-12) d1 /= norm;
+    aCurCont = GeomAbs_C0;
-    else continue; // skip degenerated point
-    if(rev1) d1.Reverse();
 
-    HC2d2->D0(u,p);
+    aSP1.Calculate(u);
-    HS2->D2 (p.X(), p.Y(), pp2, du2, dv2, d2u2, d2v2, d2uv2);
+    aSP2.Calculate(u);
-    d2 = (du2.Crossed(dv2));
-    norm = d2.Magnitude();
-    if (norm > 1.e-12) d2 /= norm;
-    else continue; // skip degenerated point
-    if(rev2) d2.Reverse();
 
-    // check 
+    aDer1 = aSP1.Derivative();
-    Standard_Real ang = d1.Angle(d2);
+    aSqLen1 = aDer1.SquareMagnitude();
-
+    aDer2 = aSP2.Derivative();
-    // check special case of precise equality of derivatives,
+    aSqLen2 = aDer2.SquareMagnitude();
-    // occurring when edge connects two faces built on equally 
+    Standard_Boolean isSmoothSuspect = (aDer1.CrossSquareMagnitude(aDer2) <= anAngleTol2 * aSqLen1 * aSqLen2);
-    // defined surfaces (e.g. seam-like edges on periodic surfaces, 
+    if (!isSmoothSuspect)
-    // or planar faces on the same plane)
-    if (aCont >= GeomAbs_C2 && ang < Precision::Angular() &&
-        d2u1 .IsEqual (d2u2,  Precision::PConfusion(), Precision::Angular()) &&
-        d2v1 .IsEqual (d2v2,  Precision::PConfusion(), Precision::Angular()) &&
-        d2uv1.IsEqual (d2uv2, Precision::PConfusion(), Precision::Angular()))
     {
-      continue;
+      // Refine by projection
+      if (aHC2.IsNull())
+      {
+        // adaptor for pcurve on the second surface
+        aHC2 = new BRepAdaptor_HCurve(BRepAdaptor_Curve(E, F2));
+        ext.Initialize(aHC2->Curve(), f, l, Precision::PConfusion());
+      }
+      ext.Perform(aSP1.Value(), u);
+      if (ext.IsDone() && ext.IsMin())
+      {
+        const Extrema_POnCurv& poc = ext.Point();
+        aSP2.Calculate(poc.Parameter());
+        aDer2 = aSP2.Derivative();
+        aSqLen2 = aDer2.SquareMagnitude();
+      }
+      isSmoothSuspect = (aDer1.CrossSquareMagnitude(aDer2) <= anAngleTol2 * aSqLen1 * aSqLen2);
+    }
+    if (isSmoothSuspect)
+    {
+      aCurCont = GeomAbs_G1;
+      if (Abs(Sqrt(aSqLen1) - Sqrt(aSqLen2)) < Precision::Confusion() &&
+          aDer1.Dot(aDer2) > Precision::SquareConfusion()) // <= check vectors are codirectional
+        aCurCont = GeomAbs_C1;
+    }
+    else
+      return GeomAbs_C0;
+
+    if (aCont < GeomAbs_G2)
+      continue; // no need further processing, because maximal continuity is less than G2
+
+    // Check conditions for G2 and C2 continuity:
+    // * calculate principal curvatures on each surface
+    // * continuity is C2 if directions of principal curvatures are equal on differenct surfaces
+    // * continuity is G2 if directions of principal curvatures are just parallel
+    //   and values of curvatures are the same
+    aSP1.Curvature(aCrvDir1[0], aCrvLen1[0], aCrvDir1[1], aCrvLen1[1]);
+    aSP2.Curvature(aCrvDir2[0], aCrvLen2[0], aCrvDir2[1], aCrvLen2[1]);
+    for (Standard_Integer aStep = 0; aStep <= 1; ++aStep)
+    {
+      if (aCrvDir1[0].XYZ().CrossSquareMagnitude(aCrvDir2[aStep].XYZ()) <= Precision::SquareConfusion() &&
+          Abs(aCrvLen1[0] - aCrvLen2[aStep]) < Precision::Confusion() &&
+          aCrvDir1[1].XYZ().CrossSquareMagnitude(aCrvDir2[1 - aStep].XYZ()) <= Precision::SquareConfusion() &&
+          Abs(aCrvLen1[1] - aCrvLen2[1 - aStep]) < Precision::Confusion())
+      {
+        if (aCurCont == GeomAbs_C1 &&
+            aCrvDir1[0].Dot(aCrvDir2[aStep]) > Precision::Confusion() &&
+            aCrvDir1[1].Dot(aCrvDir2[1 - aStep]) > Precision::Confusion())
+          aCurCont = GeomAbs_C2;
+        else
+          aCurCont = GeomAbs_G2;
+        break;
+      }
     }
 
-    aCont = GeomAbs_G1;
+    if (aCurCont < aCont)
-
+      aCont = aCurCont;
-    // Refine by projection
-    if (ang > theAngleTol)
-    {
-      if (! IsInitialized ) {
-        ext.Initialize(C2,f,l,Precision::PConfusion());
-        IsInitialized = Standard_True;
-      }      
-      ext.Perform(pp1,u);
-      if(ext.IsDone() && ext.IsMin()){
-        Extrema_POnCurv poc = ext.Point();
-        Standard_Real v = poc.Parameter();
-
-        HC2d2->D0(v,p);
-        p.Coord(uu,vv);
-        HS2->D1(p.X(), p.Y(), pp2, du2, dv2);
-        d2 = (du2.Crossed(dv2));
-        norm = d2.Magnitude();
-        if (norm> 1.e-12) d2 /= norm;
-        else continue; // degenerated point
-        if(rev2) d2.Reverse();
-        ang = d1.Angle(d2);
-      }
-      if (ang > theAngleTol)
-        return GeomAbs_C0;
-    }
   }
 
+  // according to the list of supported elementary surfaces,
+  // if the continuity is C2, than it is totally CN
+  if (isElementary && aCont == GeomAbs_C2)
+    aCont = GeomAbs_CN;
   return aCont;
 }
 
+//=======================================================================
+// function : EncodeRegularity
+// purpose  : Code the regularities on all edges of the shape, boundary of 
+//            two faces that do not have it.
+//            Takes into account that compound may consists of same solid
+//            placed with different transformations
+//=======================================================================
+static void EncodeRegularity(const TopoDS_Shape&        theShape,
+                             const Standard_Real        theTolAng,
+                             TopTools_MapOfShape&       theMap,
+                             const TopTools_MapOfShape& theEdgesToEncode = TopTools_MapOfShape())
+{
+  TopoDS_Shape aShape = theShape;
+  TopLoc_Location aNullLoc;
+  aShape.Location(aNullLoc); // nullify location
+  if (!theMap.Add(aShape))
+    return; // do not need to process shape twice
+
+  if (aShape.ShapeType() == TopAbs_COMPOUND ||
+      aShape.ShapeType() == TopAbs_COMPSOLID)
+  {
+    for (TopoDS_Iterator it(aShape); it.More(); it.Next())
+      EncodeRegularity(it.Value(), theTolAng, theMap, theEdgesToEncode);
+    return;
+  }
+
+  try {
+    OCC_CATCH_SIGNALS
+
+    TopTools_IndexedDataMapOfShapeListOfShape M;
+    TopExp::MapShapesAndAncestors(aShape, TopAbs_EDGE, TopAbs_FACE, M);
+    TopTools_ListIteratorOfListOfShape It;
+    TopExp_Explorer Ex;
+    TopoDS_Face F1,F2;
+    Standard_Boolean found;
+    for (Standard_Integer i = 1; i <= M.Extent(); i++){
+      TopoDS_Edge E = TopoDS::Edge(M.FindKey(i));
+      if (!theEdgesToEncode.IsEmpty())
+      {
+        // process only the edges from the list to update their regularity
+        TopoDS_Shape aPureEdge = E.Located(aNullLoc);
+        aPureEdge.Orientation(TopAbs_FORWARD);
+        if (!theEdgesToEncode.Contains(aPureEdge))
+          continue;
+      }
+
+      found = Standard_False;                                     
+      F1.Nullify();
+      for (It.Initialize(M.FindFromIndex(i)); It.More() && !found; It.Next()){
+        if (F1.IsNull()) { F1 = TopoDS::Face(It.Value()); }
+        else {
+          const TopoDS_Face& aTmpF2 = TopoDS::Face(It.Value());
+          if (!F1.IsSame(aTmpF2)){
+            found = Standard_True;
+            F2 = aTmpF2;
+          }
+        }
+      }
+      if (!found && !F1.IsNull()){//is it a sewing edge?
+        TopAbs_Orientation orE = E.Orientation();
+        TopoDS_Edge curE;
+        for (Ex.Init(F1, TopAbs_EDGE); Ex.More() && !found; Ex.Next()){
+          curE = TopoDS::Edge(Ex.Current());
+          if (E.IsSame(curE) && orE != curE.Orientation()) {
+            found = Standard_True;
+            F2 = F1;
+          }
+        }
+      }
+      if (found)
+        BRepLib::EncodeRegularity(E, F1, F2, theTolAng);
+    }
+  }
+  catch (Standard_Failure) {
+#ifdef OCCT_DEBUG
+    cout << "Warning: Exception in BRepLib::EncodeRegularity(): ";
+    Standard_Failure::Caught()->Print(cout);
+    cout << endl;
+#endif
+  }
+}
 
 //=======================================================================
 // function : EncodeRegularity
@@ -1804,56 +1977,39 @@ static GeomAbs_Shape tgtfaces(const TopoDS_Edge& Ed,
 void BRepLib::EncodeRegularity(const TopoDS_Shape& S,
   const Standard_Real TolAng)
 {
-  BRep_Builder B;
+  TopTools_MapOfShape aMap;
-  TopTools_IndexedDataMapOfShapeListOfShape M;
+  ::EncodeRegularity(S, TolAng, aMap);
-  TopExp::MapShapesAndAncestors(S,TopAbs_EDGE,TopAbs_FACE,M);
-  TopTools_ListIteratorOfListOfShape It;
-  TopExp_Explorer Ex;
-  TopoDS_Face F1,F2;
-  Standard_Boolean found, couture;
-  for(Standard_Integer i = 1; i <= M.Extent(); i++){
-    TopoDS_Edge E = TopoDS::Edge(M.FindKey(i));
-    found = Standard_False; couture = Standard_False;
-    F1.Nullify();
-    for(It.Initialize(M.FindFromIndex(i));It.More() && !found;It.Next()){
-      if(F1.IsNull()) { F1 = TopoDS::Face(It.Value()); }
-      else {
-        if(!F1.IsSame(TopoDS::Face(It.Value()))){
-          found = Standard_True;
-          F2 = TopoDS::Face(It.Value());
-        }
-      }
-    }
-    if (!found && !F1.IsNull()){//is it a sewing edge?
-      TopAbs_Orientation orE = E.Orientation();
-      TopoDS_Edge curE;
-      for(Ex.Init(F1,TopAbs_EDGE);Ex.More() && !found;Ex.Next()){
-        curE= TopoDS::Edge(Ex.Current());
-        if(E.IsSame(curE) && orE != curE.Orientation()) {
-          found = Standard_True;
-          couture = Standard_True;
-          F2 = F1;
-        }
-      }
-    }
-    if(found){
-      if(BRep_Tool::Continuity(E,F1,F2)<=GeomAbs_C0){
-
-        try {
-          GeomAbs_Shape aCont = tgtfaces(E, F1, F2, TolAng, couture);
-          B.Continuity(E,F1,F2,aCont);
-        }
-        catch(Standard_Failure)
-        {
-        }
-      }
-    }
-  }
 }
 
 //=======================================================================
 // function : EncodeRegularity
-// purpose  : code the regularity between 2 faces on an edge 
+// purpose  : code the regularities on all edges in the list that do not 
+//            have it, and which are boundary of two faces on the shape.
+//=======================================================================
+
+void BRepLib::EncodeRegularity(const TopoDS_Shape& S,
+  const TopTools_ListOfShape& LE,
+  const Standard_Real TolAng)
+{
+  // Collect edges without location and orientation
+  TopTools_MapOfShape aPureEdges;
+  TopLoc_Location aNullLoc;
+  TopTools_ListIteratorOfListOfShape anEdgeIt(LE);
+  for (; anEdgeIt.More(); anEdgeIt.Next())
+  {
+    TopoDS_Shape anEdge = anEdgeIt.Value();
+    anEdge.Location(aNullLoc);
+    anEdge.Orientation(TopAbs_FORWARD);
+    aPureEdges.Add(anEdge);
+  }
+
+  TopTools_MapOfShape aMap;
+  ::EncodeRegularity(S, TolAng, aMap, aPureEdges);
+}
+
+//=======================================================================
+// function : EncodeRegularity
+// purpose  : code the regularity between 2 faces connected by edge 
 //=======================================================================
 
 void BRepLib::EncodeRegularity(TopoDS_Edge& E,
@@ -1864,12 +2020,15 @@ void BRepLib::EncodeRegularity(TopoDS_Edge& E,
   BRep_Builder B;
   if(BRep_Tool::Continuity(E,F1,F2)<=GeomAbs_C0){
     try {
-      GeomAbs_Shape aCont = tgtfaces(E, F1, F2, TolAng, F1.IsEqual(F2));
+      GeomAbs_Shape aCont = tgtfaces(E, F1, F2, TolAng);
       B.Continuity(E,F1,F2,aCont);
       
     }
     catch(Standard_Failure)
     {
+#ifdef OCCT_DEBUG
+      cout << "Failure: Exception in BRepLib::EncodeRegularity" << endl;
+#endif
     }
   }
 }

src/BRepLib/BRepLib.hxx

@@ -163,11 +163,17 @@ public:
   //! is done.
   Standard_EXPORT static void EncodeRegularity (const TopoDS_Shape& S, const Standard_Real TolAng = 1.0e-10);
 
+  //! Encodes the Regularity of edges in list <LE> on the shape <S>
+  //! Warning: <TolAng> is an angular tolerance, expressed in Rad.
+  //! Warning: If the edges's regularity are coded before, nothing
+  //! is done.
+  Standard_EXPORT static void EncodeRegularity(const TopoDS_Shape& S, const TopTools_ListOfShape& LE, const Standard_Real TolAng = 1.0e-10);
+  
   //! Encodes the Regularity beetween <F1> and <F2> by <E>
   //! Warning: <TolAng> is an angular tolerance, expressed in Rad.
   //! Warning: If the edge's regularity is coded before, nothing
   //! is done.
-  Standard_EXPORT static void EncodeRegularity (TopoDS_Edge& S, const TopoDS_Face& F1, const TopoDS_Face& F2, const Standard_Real TolAng = 1.0e-10);
+  Standard_EXPORT static void EncodeRegularity (TopoDS_Edge& E, const TopoDS_Face& F1, const TopoDS_Face& F2, const Standard_Real TolAng = 1.0e-10);
   
   //! Sorts in  LF the Faces of   S on the  complexity of
   //! their                  surfaces

src/ShapeFix/ShapeFix.cxx

@@ -29,7 +29,6 @@
 
 #include <Precision.hxx>
 
-#include <Standard_ErrorHandler.hxx>
 #include <Standard_Failure.hxx>
 
 #include <TopExp_Explorer.hxx>
@@ -255,41 +254,10 @@ Standard_Boolean ShapeFix::SameParameter(const TopoDS_Shape& shape,
 //purpose  : 
 //=======================================================================
 
-static void EncodeRegularity (const TopoDS_Shape& shape,
-                              const Standard_Real tolang,
-                              TopTools_MapOfShape &aMap)
-{
-  TopoDS_Shape S = shape;
-  TopLoc_Location L;
-  S.Location ( L );
-  if ( ! aMap.Add ( S ) ) return;
-  
-  if ( S.ShapeType() == TopAbs_COMPOUND || 
-       S.ShapeType() == TopAbs_COMPSOLID ) {
-    for ( TopoDS_Iterator it(S); it.More(); it.Next() ) {
-      EncodeRegularity ( it.Value(), tolang, aMap );
-    }
-    return;
-  }
-  
-  try {
-    OCC_CATCH_SIGNALS
-    BRepLib::EncodeRegularity ( S, tolang );
-  }
-  catch(Standard_Failure) {
-#ifdef OCCT_DEBUG
-    cout << "Warning: Exception in ShapeFix::EncodeRegularity(): ";
-    Standard_Failure::Caught()->Print ( cout );
-    cout << endl;
-#endif
-  }
-}
-
 void ShapeFix::EncodeRegularity (const TopoDS_Shape& shape,
                                  const Standard_Real tolang)
 {
-  TopTools_MapOfShape aMap;
+  BRepLib::EncodeRegularity(shape, tolang);
-  ::EncodeRegularity ( shape, tolang, aMap );
 }
 
 

tests/bugs/modalg_6/bug27341_313

@@ -23,6 +23,6 @@ build3d result
 fit
 
 checkprops result -l 9662.5
-checknbshapes result -vertex 4423 -edge 2218
+checknbshapes result -vertex 4419 -edge 2216
 
 checkview -screenshot -2d -path ${imagedir}/${test_image}.png

tests/bugs/modalg_6/bug27341_314

@@ -22,7 +22,7 @@ build3d result
 
 fit
 
-checkprops result -l 6141.2
+checkprops result -l 5934.34
-checknbshapes result -vertex 1444 -edge 722
+checknbshapes result -vertex 1414 -edge 707
 
 checkview -screenshot -2d -path ${imagedir}/${test_image}.png

tests/bugs/modalg_6/bug27383_1

@@ -2,7 +2,9 @@ puts "========"
 puts "0027383: Modeling - improve handling of regularity on edges"
 puts "========"
 puts ""
-puts "Check regularity setting on edges between same-domain faces after fuse"
+#################################################################
+# Check regularity setting on edges between same-domain faces after fuse
+#################################################################
 
 # fuse of two boxes
 pload MODELING
@@ -14,9 +16,6 @@ explode r_7 e
 getedgeregularity r_7_3 r_3 r_7 ;# returns C0: regularity is not set at all
 encoderegularity r
 
-explode res f
-explode res_3 e
-explode res_4 e
 if { ! [regexp "CN" [getedgeregularity r_7_3 r_3 r_7]] } {
   puts "Error: Invalid regularity of the edge, expected CN"
 }

tests/bugs/modalg_6/bug27383_2

@@ -2,7 +2,9 @@ puts "========"
 puts "0027383: Modeling - improve handling of regularity on edges"
 puts "========"
 puts ""
-puts "Check regularity setting on edges of cylinder split by angle"
+#################################################################
+# Check regularity setting on edges of cylinder split by angle
+#################################################################
 
 # split of cylinder
 pload MODELING

tests/bugs/modalg_6/bug27383_3

@@ -0,0 +1,44 @@
+puts "========"
+puts "0027383: Modeling - improve handling of regularity on edges"
+puts "========"
+puts ""
+#################################################################
+# Check regularity setting on edges between rotated or mirrored half spheres
+#################################################################
+
+pload MODELING
+psphere s1 10 180
+psphere s2 10 180
+explode s1 f
+explode s2 f
+tmirror s2_1 0 0 0 0 1 0
+trotate s2_1 0 0 0 0 1 0 45
+sewing r s1_1 s2_1
+encoderegularity r
+
+explode r F
+set edges [explode r_1 E]
+set len 0
+foreach e ${edges} {
+  # check for degenerated edge
+  set e_props [lprops $e 1.e-4]
+  regexp {Mass : +([0-9.+-eE]+)} $e_props full len
+
+  if { [expr abs($len) > 1.e-7] } {
+    set cont "G2"
+  } else {
+    # degenerated edges have always C0 continuity
+    set cont "C0"
+  }
+
+  if { ! [regexp "${cont}" [getedgeregularity $e r_1 r_2]] } {
+    puts "Error: Invalid regularity of the edge, expected ${cont}"
+  }
+}
+
+# make image in HLR mode as illustration
+pload VISUALIZATION
+vdisplay r
+vfit
+vhlr on
+vdump ${imagedir}/${test_image}.png

tests/bugs/modalg_6/bug27383_4

@@ -0,0 +1,30 @@
+puts "========"
+puts "0027383: Modeling - improve handling of regularity on edges"
+puts "========"
+puts ""
+#################################################################
+# Check regularity setting on edges between patches of B-Spline surface obtained by conversion of cylinder
+#################################################################
+
+pload MODELING
+pcylinder p 2 10
+explode p f
+nurbsconvert q p_1
+DT_ClosedSplit r q
+encoderegularity r
+
+explode r F
+explode r_1 E
+if { ! [regexp "G2" [getedgeregularity r_1_2 r_1 r_2]] } {
+  puts "Error: Invalid regularity of the edge, expected G2"
+}
+if { ! [regexp "G2" [getedgeregularity r_1_4 r_1 r_2]] } {
+  puts "Error: Invalid regularity of the edge, expected G2"
+}
+
+# make image in HLR mode as illustration
+pload VISUALIZATION
+vdisplay r
+vfit
+vhlr on
+vdump ${imagedir}/${test_image}.png

tests/bugs/modalg_6/bug27383_5

@@ -0,0 +1,28 @@
+puts "========"
+puts "0027383: Modeling - improve handling of regularity on edges"
+puts "========"
+puts ""
+#################################################################
+# Check regularity setting on generatrix shared by two cylinders
+#################################################################
+
+pload MODELING
+pcylinder cyl 2 10 90
+explode cyl F
+tcopy cyl_1 c2
+trotate c2 0 2 0 0 0 1 180
+sewing r cyl_1 c2
+encoderegularity r
+
+explode r F
+explode r_1 E
+if { ! [regexp "C1" [getedgeregularity r_1_2 r_1 r_2]] } {
+  puts "Error: Invalid regularity of the edge, expected C1"
+}
+
+# make image in HLR mode as illustration
+pload VISUALIZATION
+vdisplay r
+vfit
+vhlr on
+vdump ${imagedir}/${test_image}.png

tests/bugs/modalg_6/bug27383_6

@@ -0,0 +1,40 @@
+puts "========"
+puts "0027383: Modeling - improve handling of regularity on edges"
+puts "========"
+puts ""
+#################################################################
+# Check regularity setting on edges between cylinder and a plane
+#################################################################
+
+pload MODELING
+plane p 0 2 0 0 1 0
+mkface pln p 0 10 -2 0
+pcylinder cyl 2 10 90
+explode cyl F
+sewing r cyl_1 pln
+encoderegularity r
+
+set cont "G1"
+
+explode r F
+explode r_1 E
+if { ! [regexp "${cont}" [getedgeregularity r_1_2 r_1 r_2]] } {
+  puts "Error: Invalid regularity of the edge, expected ${cont}"
+}
+
+nurbsconvert ncyl cyl_1
+sewing r ncyl pln
+encoderegularity r
+
+explode r F
+explode r_1 E
+if { ! [regexp "${cont}" [getedgeregularity r_1_2 r_1 r_2]] } {
+  puts "Error: Invalid regularity of the edge, expected ${cont}"
+}
+
+# make image in HLR mode as illustration
+pload VISUALIZATION
+vdisplay r
+vfit
+vhlr on
+vdump ${imagedir}/${test_image}.png

tests/bugs/modalg_6/bug27383_7

@@ -0,0 +1,68 @@
+puts "========"
+puts "0027383: Modeling - improve handling of regularity on edges"
+puts "========"
+puts ""
+#################################################################
+# Check regularity setting on edges between parts of cylinder and spheres
+#################################################################
+
+pload MODELING
+pcylinder p1 2 10 180
+pcylinder p2 2 10 180
+tmirror p2 0 0 0 0 1 0
+nurbsconvert q p2
+
+psphere s1 2
+ttranslate s1 0 0 10
+psphere s2 2
+trotate s2 0 0 0 0 1 0 -90
+
+bclearobjects
+bcleartools
+baddobjects p1 p2
+baddtools s1 s2
+bfillds
+# fuse all solids
+bbop r 1
+
+encoderegularity r
+
+# check that minimal continuity is G1
+set rfaces [explode r F]
+set len [llength $rfaces]
+set nbshared 0
+set nbexpected 6
+for {set i 0} {$i < [expr $len-1]} {incr i} {
+  set f1 [lindex $rfaces $i]
+  set redges1 [explode $f1 E]
+
+  for {set j [expr $i+1]} {$j < $len} {incr j} {
+    set f2 [lindex $rfaces $j]
+    if {$f1 == $f2} continue
+
+    set redges2 [explode $f2 E]
+    foreach e1 $redges1 {
+      foreach e2 $redges2 {
+        if {[regexp "not" [compare $e1 $e2]]} continue
+
+        incr nbshared
+        if { [regexp "C0" [getedgeregularity $e1 $f1 $f2]] } {
+          puts "Error: Invalid regularity of the edge, expected at least G1"
+        }
+      }
+    }
+  }
+}
+
+if {${nbshared} == ${nbexpected} } {
+  puts "OK: Number of processed edges: ${nbshared}"
+} else {
+  puts "Error: incorrect number of processed edges (${nbshared} instead of ${nbexpected})"
+}
+
+# make image in HLR mode as illustration
+pload VISUALIZATION
+vdisplay r
+vfit
+vhlr on
+vdump ${imagedir}/${test_image}.png