0031524: Unify same domain corrupts shape representing a cylindrical tube

- Added test case bug31524;
- Added test data shapee bug31524.brep;
- Modified ShapeUpgrade_UnifySameDomain.cxx.

src/ShapeUpgrade/ShapeUpgrade_UnifySameDomain.cxx

@@ -21,7 +21,11 @@
 #include <BRepLib_MakeFace.hxx>
 #include <BRepTopAdaptor_TopolTool.hxx>
 #include <GC_MakeCircle.hxx>
+#include <Geom2d_Circle.hxx>
+#include <Geom2d_Ellipse.hxx>
+#include <Geom2d_Hyperbola.hxx>
 #include <Geom2d_Line.hxx>
+#include <Geom2d_Parabola.hxx>
 #include <GCE2d_MakeLine.hxx>
 #include <Geom2d_TrimmedCurve.hxx>
 #include <Geom2dConvert.hxx>
@@ -67,6 +71,7 @@
 #include <TColGeom_SequenceOfSurface.hxx>
 #include <TColStd_Array1OfReal.hxx>
 #include <TColStd_MapOfInteger.hxx>
+#include <TColStd_SequenceOfBoolean.hxx>
 #include <TopExp.hxx>
 #include <TopExp_Explorer.hxx>
 #include <TopoDS.hxx>
@@ -1423,6 +1428,267 @@ static TopoDS_Edge GlueEdgesWithPCurves(const TopTools_SequenceOfShape& aChain,
   return ResEdge;
 }
 
+//=======================================================================
+//function : Concat2dCurves
+//purpose  : Concatenates 2d curves
+//=======================================================================
+static void Concat2dCurves(
+  TopoDS_Edge& theEdge,
+  const TopTools_SequenceOfShape& theChain)
+{
+  TColGeom_SequenceOfSurface aSurfSeq;
+  NCollection_Sequence<TopLoc_Location> aLocSeq;
+  TopoDS_Edge aFirstEdge = TopoDS::Edge(theChain(1));
+  for (int aCurveIndex = 1; ; aCurveIndex++)
+  { //Searh surfaces of edge
+    Handle(Geom2d_Curve) aCurve;
+    Handle(Geom_Surface) aSurface;
+    TopLoc_Location aLocation;
+    Standard_Real aFirst, aLast;
+    BRep_Tool::CurveOnSurface(aFirstEdge, aCurve, aSurface, aLocation, aFirst, aLast, aCurveIndex);
+    if (aCurve.IsNull())
+      break;
+    aSurfSeq.Append(aSurface);
+    aLocSeq.Append(aLocation);
+  } //End of search
+
+  Standard_Integer i, j;
+  Standard_Integer aNbCurves = theChain.Length();
+  TColGeom2d_SequenceOfBoundedCurve ResPCurves;
+  TopoDS_Vertex aPrevVertex, aFirstVOfEdge, aLastVOfEdge;
+  TopExp::Vertices(theEdge, aFirstVOfEdge, aLastVOfEdge);
+  Standard_Real fpar, lpar;
+  TColStd_Array1OfReal tabtolvertex(0, aNbCurves - 1);
+  Standard_Real MaxTol = Precision::Confusion();
+
+  for (j = 1; j <= aSurfSeq.Length(); j++)
+  { //Cycle of surfaces
+    TColGeom2d_SequenceOfCurve CurveSeq;
+    TColStd_SequenceOfReal FparSeq;
+    TColStd_SequenceOfReal LparSeq;
+    TColStd_SequenceOfReal TolSeq;
+    TColStd_SequenceOfBoolean IsFwdSeq;
+    GeomAbs_CurveType CurType = GeomAbs_OtherCurve;
+    TopoDS_Vertex FirstVertex, LastVertex;
+
+    aPrevVertex = aFirstVOfEdge;
+    for (i = 1; i <= aNbCurves; i++)
+    { //Cycle of curves
+      TopoDS_Edge anEdge = TopoDS::Edge(theChain(i));
+      TopoDS_Vertex VF, VL;
+      TopExp::Vertices(anEdge, VF, VL);
+
+      Handle(Geom2d_Curve) aPCurve =
+        BRep_Tool::CurveOnSurface(anEdge, aSurfSeq(j), aLocSeq(j), fpar, lpar);
+      if (aPCurve.IsNull())
+        continue;
+
+      Geom2dAdaptor_Curve anAdaptor(aPCurve);
+      GeomAbs_CurveType aType = anAdaptor.GetType();
+      Handle(Geom2d_Curve) aBasisCurve = anAdaptor.Curve();
+      Standard_Boolean isForward = (anEdge.Orientation() != TopAbs_REVERSED);
+
+      if (aBasisCurve->IsPeriodic()) {
+        ElCLib::AdjustPeriodic
+        (aBasisCurve->FirstParameter(), aBasisCurve->LastParameter(),
+          Precision::PConfusion(), fpar, lpar);
+      }
+
+      if (CurveSeq.IsEmpty()) {
+        CurveSeq.Append(aPCurve);
+        FparSeq.Append(fpar);
+        LparSeq.Append(lpar);
+        IsFwdSeq.Append(isForward);
+        CurType = aType;
+        FirstVertex = VF;
+      }
+      else {
+        Standard_Real aNewFpar = RealFirst(), aNewLpar = RealLast();
+        Standard_Boolean isSameCurve = Standard_False;
+        Geom2dAdaptor_Curve aPrevAdaptor(CurveSeq.Last());
+
+        if (aPCurve == CurveSeq.Last()) {
+          aNewFpar = fpar;
+          aNewLpar = lpar;
+          isSameCurve = Standard_True;
+        }
+        else if (aType == CurType) {
+          switch (aType)
+          {
+          case GeomAbs_Line:
+          {
+            Geom2d_Line aLine = anAdaptor.Line();
+            Geom2d_Line aPrevLine = aPrevAdaptor.Line();
+            if (aLine.Location().IsEqual(aPrevLine.Location(), Precision::Confusion()) &&
+              aLine.Direction().IsParallel(aPrevLine.Direction(), Precision::Angular())) {
+              gp_Pnt2d P1 = ElCLib::Value(fpar, aLine.Lin2d());
+              gp_Pnt2d P2 = ElCLib::Value(lpar, aLine.Lin2d());
+
+              aNewFpar = ElCLib::Parameter(aPrevLine.Lin2d(), P1);
+              aNewLpar = ElCLib::Parameter(aPrevLine.Lin2d(), P2);
+              isSameCurve = Standard_True;
+            }
+          }
+            break;
+          case GeomAbs_Circle:
+            break;
+          case GeomAbs_Ellipse:
+            break;
+          case GeomAbs_Hyperbola:
+            break;
+          case GeomAbs_Parabola:
+            break;
+          case GeomAbs_BezierCurve:
+            break;
+          case GeomAbs_BSplineCurve:
+            break;
+          case GeomAbs_OffsetCurve:
+            break;
+          case GeomAbs_OtherCurve:
+            break;
+          default:
+            break;
+          } //end of switch
+
+          if (isSameCurve) {
+            const Standard_Boolean isSameDir = (isForward == IsFwdSeq.Last());
+
+            if (aBasisCurve->IsPeriodic()) {
+              // Treat periodic curves.
+              const Standard_Real aPeriod = aBasisCurve->Period();
+
+              if (isSameDir) {
+                // Check if first parameter is greater then the last one.
+                while (aNewFpar > aNewLpar) {
+                  aNewFpar -= aPeriod;
+                }
+              }
+              else { // !isSameDir
+                     // Check if last parameter is greater then the first one.
+                while (aNewLpar > aNewFpar) {
+                  aNewLpar -= aPeriod;
+                }
+
+                // Change parameters
+                const Standard_Real aTmpPar = aNewLpar;
+
+                aNewLpar = aNewFpar;
+                aNewFpar = aTmpPar;
+              }
+
+              // Udjust parameters on periodic curves.
+              if (IsFwdSeq.Last()) {
+                // The current curve should be after the previous one.
+                ElCLib::AdjustPeriodic(LparSeq.Last(), LparSeq.Last() + aPeriod,
+                  Precision::PConfusion(), aNewFpar, aNewLpar);
+              }
+              else {
+                // The current curve should be before the previous one.
+                ElCLib::AdjustPeriodic(FparSeq.Last() - aPeriod, FparSeq.Last(),
+                  Precision::PConfusion(), aNewFpar, aNewLpar);
+              }
+            }
+            else if (!isSameDir) {
+              // Not periodic curves. Opposite dirs.
+              const Standard_Real aTmpPar = aNewLpar;
+
+              aNewLpar = aNewFpar;
+              aNewFpar = aTmpPar;
+            }
+
+            if (IsFwdSeq.Last()) {
+              // Update last parameter 
+              LparSeq(LparSeq.Length()) = aNewLpar;
+            }
+            else {
+              // Update first parameter 
+              FparSeq(FparSeq.Length()) = aNewFpar;
+            }
+          }
+          else {
+            // Add new curve.
+            CurveSeq.Append(aPCurve);
+            FparSeq.Append(fpar);
+            LparSeq.Append(lpar);
+            IsFwdSeq.Append(isForward);
+            TolSeq.Append(BRep_Tool::Tolerance(VF));
+            CurType = aType;
+          }
+        } //end of else
+      }
+    } //end of cycle searching and forming curves
+
+    LastVertex = aLastVOfEdge;
+    TolSeq.Append(BRep_Tool::Tolerance(LastVertex));
+
+    Standard_Boolean isReverse = Standard_False;
+    if (!IsFwdSeq.IsEmpty()) {
+      isReverse = !IsFwdSeq(1);
+    }
+    if (isReverse)
+    {
+      aFirstVOfEdge = LastVertex;
+      aLastVOfEdge = FirstVertex;
+    }
+    else
+    {
+      aFirstVOfEdge = FirstVertex;
+      aLastVOfEdge = LastVertex;
+    }
+    aFirstVOfEdge.Orientation(TopAbs_FORWARD);
+    aLastVOfEdge.Orientation(TopAbs_REVERSED);
+
+    TColGeom2d_Array1OfBSplineCurve tab_c2d(0, CurveSeq.Length() - 1); //array of the pcurves
+    for (i = 1; i <= CurveSeq.Length(); i++) {
+      Handle(Geom2d_TrimmedCurve) aTrPCurve = 
+        new Geom2d_TrimmedCurve(CurveSeq(i), FparSeq(i), LparSeq(i));
+      tab_c2d(i - 1) = Geom2dConvert::CurveToBSplineCurve(aTrPCurve);
+      Geom2dConvert::C0BSplineToC1BSplineCurve(tab_c2d(i - 1), Precision::Confusion());
+      if (!IsFwdSeq(i))
+        tab_c2d(i - 1)->Reverse();
+    }
+    
+    if (tab_c2d.Length() > 1) {
+      //Prepare C1 concatenation
+      Handle(TColGeom2d_HArray1OfBSplineCurve)  concatc2d;     //array of the concatenated curves
+      Handle(TColStd_HArray1OfInteger)        ArrayOfInd2d;  //array of the remining Vertex
+      Standard_Boolean closed_flag = Standard_False;
+      Geom2dConvert::ConcatC1(tab_c2d,
+        tabtolvertex,
+        ArrayOfInd2d,
+        concatc2d,
+        closed_flag,
+        Precision::Confusion());   //C1 concatenation
+
+      if (concatc2d->Length() > 1) //concatenation
+      {
+        Geom2dConvert_CompCurveToBSplineCurve Concat2d(concatc2d->Value(concatc2d->Lower()));
+
+        for (i = concatc2d->Lower() + 1; i <= concatc2d->Upper(); i++)
+          Concat2d.Add(concatc2d->Value(i), MaxTol, Standard_True);
+
+        concatc2d->SetValue(concatc2d->Lower(), Concat2d.BSplineCurve());
+      } //end of concatenation
+      Handle(Geom2d_BSplineCurve) aResPCurve = concatc2d->Value(concatc2d->Lower());
+      ResPCurves.Append(aResPCurve);
+    }
+    else {
+      ResPCurves.Append(tab_c2d(0));
+    }
+  } //end of cycle of surfaces
+
+  //Update result edge
+  BRep_Builder aBuilder;
+  for (j = 1; j <= ResPCurves.Length(); j++)
+  {
+    aBuilder.UpdateEdge(theEdge, ResPCurves(j), aSurfSeq(j), aLocSeq(j), MaxTol);
+    aBuilder.Range(theEdge, aSurfSeq(j), aLocSeq(j), ResPCurves(j)->FirstParameter(), ResPCurves(j)->LastParameter());
+  }
+
+  BRepLib::SameParameter(theEdge, MaxTol, Standard_True);
+
+} //end of function
+
 //=======================================================================
 //function : MergeSubSeq
 //purpose  : Merges a sequence of edges into one edge if possible
@@ -1500,12 +1766,12 @@ static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
     if(c3d1.IsNull() || c3d2.IsNull()) 
       return Standard_False;
 
-    while(c3d1->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
+    if (c3d1->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
       Handle(Geom_TrimmedCurve) tc =
         Handle(Geom_TrimmedCurve)::DownCast(c3d1);
       c3d1 = tc->BasisCurve();
     }
-    while(c3d2->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
+    if (c3d2->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
       Handle(Geom_TrimmedCurve) tc =
         Handle(Geom_TrimmedCurve)::DownCast(c3d2);
       c3d2 = tc->BasisCurve();
@@ -1562,6 +1828,7 @@ static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
     B.Add (E,V[0]);  B.Add (E,V[1]);
     B.UpdateVertex(V[0], 0., E, 0.);
     B.UpdateVertex(V[1], dist, E, 0.);
+    Concat2dCurves(E, theChain);
     OutEdge = E;
     return Standard_True;
   }
@@ -1572,7 +1839,7 @@ static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
     TopoDS_Edge FE = TopoDS::Edge(theChain.First());
     Handle(Geom_Curve) c3d = BRep_Tool::Curve(FE,f,l);
 
-    while(c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
+    if (c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
       Handle(Geom_TrimmedCurve) tc =
         Handle(Geom_TrimmedCurve)::DownCast(c3d);
       c3d = tc->BasisCurve();
@@ -1617,7 +1884,12 @@ static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
         B.Add(E,V[1]);
       }
     }
-    else {
+    else //open chain
+    {
+      Standard_Real ParamFirst = BRep_Tool::Parameter(V[0], FE);
+      TopoDS_Vertex VertexLastOnFE = sae.LastVertex(FE);
+      Standard_Real ParamLast  = BRep_Tool::Parameter(VertexLastOnFE, FE);
+      
       if (isSafeInputMode) {
         for (int k = 0; k < 2; k++) {
           if (!theContext->IsRecorded(V[k])) {
@@ -1629,39 +1901,32 @@ static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
             V[k] = TopoDS::Vertex(theContext->Apply(V[k]));
         }
       }
-      gp_Pnt PV1 = BRep_Tool::Pnt(V[0]);
-      gp_Pnt PV2 = BRep_Tool::Pnt(V[1]);
-      TopoDS_Vertex VM = sae.LastVertex(FE);
-      gp_Pnt PVM = BRep_Tool::Pnt(VM);
-      GC_MakeCircle MC (PV1,PVM,PV2);
-      Handle(Geom_Circle) C = MC.Value();
-      gp_Pnt P0 = C->Location();
-      gp_Dir D1(gp_Vec(P0,PV1));
-      gp_Dir D2(gp_Vec(P0,PV2));
-      Standard_Real fpar = C->XAxis().Direction().Angle(D1);
-      if(fabs(fpar)>Precision::Confusion()) {
-        // check orientation
-        gp_Dir ND =  C->XAxis().Direction().Crossed(D1);
-        if(ND.IsOpposite(C->Axis().Direction(),Precision::Confusion())) {
-          fpar = -fpar;
-        }
-      }
-      Standard_Real lpar = C->XAxis().Direction().Angle(D2);
-      if(fabs(lpar)>Precision::Confusion()) {
-        // check orientation
-        gp_Dir ND =  C->XAxis().Direction().Crossed(D2);
-        if(ND.IsOpposite(C->Axis().Direction(),Precision::Confusion())) {
-          lpar = -lpar;
-        }
-      }
-      if (lpar < fpar) lpar += 2*M_PI;
-      Handle(Geom_TrimmedCurve) tc = new Geom_TrimmedCurve(C,fpar,lpar);
+      
+      gp_Pnt PointFirst = BRep_Tool::Pnt(V[0]);
+      while (Abs(ParamLast - ParamFirst) > 7*M_PI/8)
+        ParamLast = (ParamFirst + ParamLast)/2;
+      BRepAdaptor_Curve BAcurveFE(FE);
+      gp_Pnt PointLast = BAcurveFE.Value(ParamLast);
+      gp_Pnt Origin = Cir->Circ().Location();
+      gp_Dir Dir1 = gp_Vec(Origin, PointFirst);
+      gp_Dir Dir2 = gp_Vec(Origin, PointLast);
+      gp_Dir Vdir = Dir1 ^ Dir2;
+      gp_Ax2 anAx2(Origin, Vdir, Dir1);
+      Handle(Geom_Circle) aNewCircle = new Geom_Circle(anAx2, Cir->Radius());
+      gp_Pnt PointLastInChain = BRep_Tool::Pnt(V[1]);
+      gp_Dir DirLastInChain = gp_Vec(Origin, PointLastInChain);
+      Standard_Real lpar = Dir1.AngleWithRef(DirLastInChain, Vdir);
+      if (lpar < 0.)
+        lpar += 2*M_PI;
+      
+      Handle(Geom_TrimmedCurve) tc = new Geom_TrimmedCurve(aNewCircle,0.,lpar);
       B.MakeEdge (E,tc,Precision::Confusion());
       B.Add(E,V[0]);
       B.Add(E,V[1]);
-      B.UpdateVertex(V[0], fpar, E, 0.);
+      B.UpdateVertex(V[0], 0., E, 0.);
       B.UpdateVertex(V[1], lpar, E, 0.);
     }
+    Concat2dCurves(E, theChain);
     OutEdge = E;
     return Standard_True;
   }
@@ -1678,7 +1943,7 @@ static Standard_Boolean MergeSubSeq(const TopTools_SequenceOfShape& theChain,
       TopLoc_Location Loc;
       Handle(Geom_Curve) c3d = BRep_Tool::Curve(edge,Loc,fp1,lp1);
       if(c3d.IsNull()) continue;
-      while(c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
+      if (c3d->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
         Handle(Geom_TrimmedCurve) tc =
           Handle(Geom_TrimmedCurve)::DownCast(c3d);
         c3d = tc->BasisCurve();
@@ -3211,4 +3476,4 @@ void ShapeUpgrade_UnifySameDomain::FillHistory()
 
   // Merge the history of the operation into global history
   myHistory->Merge(aUSDHistory);
-}
+}

tests/bugs/heal/bug31524

@@ -0,0 +1,9 @@
+puts "=========================================================================="
+puts "OCC31524: Unify same domain corrupts shape representing a cylindrical tube"
+puts "=========================================================================="
+puts ""
+
+pload ALL
+restore [locate_data_file bug31524.brep] sh
+unifysamedom res sh
+checkshape sh

tests/bugs/heal/bug31778

@@ -0,0 +1,28 @@
+puts "============================================="
+puts "OCC31778: UnifySameDomain fails in Debug mode"
+puts "============================================="
+puts ""
+
+brestore [locate_data_file bug31778.brep] s
+explode s
+bclearobjects
+bcleartools
+baddobjects s_1
+baddtools s_2 s_3
+bfillds
+bbop q 1
+explode q
+
+unifysamedom result q_1
+
+checkshape result
+
+checknbshapes result -solid 1 -shell 1 -face 19 -wire 21 -edge 51 -vertex 34
+
+set tolres [checkmaxtol result]
+
+if { ${tolres} > 5.e-5} {
+   puts "Error: bad tolerance of result"
+}
+
+checkprops result -v 15173.9