Skip to content
Snippets Groups Projects
Select Git revision
  • 14324364f03421b62cd557d8af3c1d0f14a2bec0
  • master default protected
  • alphashapes
  • quadMeshingTools
  • cygwin_conv_path
  • macos_arm64
  • add-transfiniteautomatic-to-geo
  • patch_releases_4_10
  • HierarchicalHDiv
  • isuruf-master-patch-63355
  • hyperbolic
  • hexdom
  • hxt_update
  • jf
  • 1618-pythonocc-and-gmsh-api-integration
  • octreeSizeField
  • hexbl
  • alignIrregularVertices
  • getEdges
  • patch_releases_4_8
  • isuruf-master-patch-51992
  • gmsh_4_11_0
  • gmsh_4_10_5
  • gmsh_4_10_4
  • gmsh_4_10_3
  • gmsh_4_10_2
  • gmsh_4_10_1
  • gmsh_4_10_0
  • gmsh_4_9_5
  • gmsh_4_9_4
  • gmsh_4_9_3
  • gmsh_4_9_2
  • gmsh_4_9_1
  • gmsh_4_9_0
  • gmsh_4_8_4
  • gmsh_4_8_3
  • gmsh_4_8_2
  • gmsh_4_8_1
  • gmsh_4_8_0
  • gmsh_4_7_1
  • gmsh_4_7_0
41 results

Partition_Loop2d.cxx

Blame
  • Forked from gmsh / gmsh
    Source project has a limited visibility.
    Partition_Loop2d.cxx 36.64 KiB
    #ifdef OCCGEOMETRY
    
    //  GEOM PARTITION : partition algorithm
    //
    //  Copyright (C) 2003  CEA/DEN, EDF R& D
    //
    //
    //
    //  File   : Partition_Loop2d.cxx
    //  Author : Benedicte MARTIN
    //  Module : GEOM
    //  $Header: /cvs/netgen/netgen/libsrc/occ/Partition_Loop2d.cxx,v 1.6 2008/03/31 14:20:28 wabro Exp $
    
    //using namespace std;
    #if !defined(HAVE_NO_OCC_CONFIG_H)
    #include <config.h>
    #endif
    #include <Standard_Version.hxx>
    
    //using namespace std;
    #include <climits>
    #include "Partition_Loop2d.ixx"
    
    //#include "utilities.h"
    #include <stdio.h>
    
    #include <BRepAdaptor_Curve2d.hxx>
    #include <BRepAdaptor_Surface.hxx>
    #include <BRepAlgo_AsDes.hxx>
    #include <BRepAlgo_FaceRestrictor.hxx>
    #if (OCC_VERSION_MAJOR == 6) && (OCC_VERSION_MINOR < 5)
    #include <BRepOffset_DataMapOfShapeReal.hxx>
    #else
    #include <TopTools_DataMapOfShapeReal.hxx>
    #endif
    #include <BRepTopAdaptor_FClass2d.hxx>
    #include <BRep_Builder.hxx>
    #include <BRep_Tool.hxx>
    #include <Geom2dInt_GInter.hxx>
    #include <Geom2d_Curve.hxx>
    #include <IntRes2d_IntersectionPoint.hxx>
    #include <Precision.hxx>
    #include <TColStd_MapOfInteger.hxx>
    #include <TColStd_SequenceOfReal.hxx>
    #include <TopExp.hxx>
    #include <TopExp_Explorer.hxx>
    #include <TopTools_DataMapIteratorOfDataMapOfShapeListOfShape.hxx>
    #include <TopTools_DataMapIteratorOfDataMapOfShapeShape.hxx>
    #include <TopTools_DataMapOfShapeInteger.hxx>
    #include <TopTools_DataMapOfShapeShape.hxx>
    #include <TopTools_IndexedMapOfShape.hxx>
    #include <TopTools_ListIteratorOfListOfShape.hxx>
    #include <TopTools_MapIteratorOfMapOfShape.hxx>
    #include <TopTools_MapOfOrientedShape.hxx>
    #include <TopTools_MapOfShape.hxx>
    #include <TopTools_SequenceOfShape.hxx>
    #include <TopoDS.hxx>
    #include <TopoDS_Iterator.hxx>
    #include <TopoDS_Vertex.hxx>
    #include <TopoDS_Wire.hxx>
    #include <gp_Pnt.hxx>
    #include <gp_Pnt2d.hxx>
    
    //=======================================================================
    //function : Partition_Loop2d
    //purpose  :
    //=======================================================================
    
    Partition_Loop2d::Partition_Loop2d()
    {
    }
    
    //=======================================================================
    //function : Init
    //purpose  : Init with <F> the set of edges must have
    //           pcurves on <F>.
    //=======================================================================
    
    void Partition_Loop2d::Init(const TopoDS_Face& F)
    {
      myConstEdges.Clear();
      myNewWires  .Clear();
      myNewFaces  .Clear();
      myFace = F;
      myFaceOri = myFace.Orientation();
      myFace.Orientation( TopAbs_FORWARD );
    }
    
    //=======================================================================
    //function : AddConstEdge
    //purpose  : Add <E> as unique edge in the result.
    //=======================================================================
    
    void Partition_Loop2d::AddConstEdge (const TopoDS_Edge& E)
    {
    #ifdef DEB
      Standard_Real f,l;
      Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( E, myFace, f,l);
      if (pc.IsNull()) {
        INFOS( "AddConstEdge(): EDGE W/O PCURVE on FACE");
      } else
    #endif
      {
        myConstEdges.Append(E);
      }
    }
    
    void Partition_Loop2d::AddSectionEdge (const TopoDS_Edge& E)
    {
    #ifdef DEB
      Standard_Real f,l;
      Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( E, myFace, f,l);
      if (pc.IsNull())
        pc = BRep_Tool::CurveOnSurface( E, myFace, f,l);
      gp_Vec2d Tg1;
      gp_Pnt2d PC;
      pc->D1(0.5*(f+l), PC, Tg1);
      if (Tg1.Magnitude()  <= gp::Resolution()) {
        MESSAGE ("");
      }
      if (pc.IsNull()) {
        INFOS( "AddConstEdge(): EDGE W/O PCURVE on FACE");
      } else
    #endif
      {
        myConstEdges.Append(E);
        myConstEdges.Append(E.Reversed());
        mySectionEdges.Add( E );
      }
    }
    
    //=======================================================================
    //function : preciseU
    //purpose  : find u such that the 3D point on theE is just out of tolerance
    //           of theV
    //=======================================================================
    
    static Standard_Real preciseU (const BRepAdaptor_Surface&  theSurf,
                                   const TopoDS_Edge&          theE,
                                   const TopoDS_Vertex&        theV,
                                   const Handle(Geom2d_Curve)& theC,
                                   const Standard_Boolean      theFirstEnd)
    {
      Standard_Boolean isForward = ( theE.Orientation () == TopAbs_FORWARD );
      if (theFirstEnd) isForward = !isForward;
    
      // find the first point in 2d and 3d
      Standard_Real f,l;
      BRep_Tool::Range( theE, f, l );
      Standard_Real u0 = isForward ? l : f;
      gp_Pnt2d aP2d0 = theC->Value( u0 );
      gp_Pnt aPnt0 = theSurf.Value( aP2d0.X(), aP2d0.Y() );
    
      // shift in 2d and 3d
      Standard_Real du = ( l - f ) / 100, du3d = 0;
      if (isForward)
        du = -du;
    
      // target parameter
      Standard_Real u;
    
      while (du3d < ::RealSmall())
      {
        // u for test
        u = u0 + du;
        du *= 10; // for the next iteration: increase du untill du3d is large enough
    
        // find out how u is far from u0 in 3D
        gp_Pnt2d aP2d  = theC->Value( u );
        gp_Pnt aPnt  = theSurf.Value( aP2d.X(), aP2d.Y() );
        du3d = aPnt0.Distance( aPnt );
      }
    
      // find u such that the 3D point is just out of tolerance of theV
      Standard_Real tolV = BRep_Tool::Tolerance( theV ) + Precision::Confusion();
      u = u0 + du * tolV / du3d;
    
      // check that u is within the range
      if ( isForward ? (u < f) : (u > l) )
        u = u0 + du;
    
      return u;
    }
    
    //=======================================================================
    //function : SelectEdge
    //purpose  : Find in the list <LE> the edge <NE> connected with <CE> by
    //           the vertex <CV>.
    //           <NE> is removed from the list. If <CE> is in <LE>
    //           with the same orientation, it's removed from the list
    //=======================================================================
    
    static Standard_Boolean  SelectEdge(const BRepAdaptor_Surface& Surf,
                                        const TopoDS_Edge&    CE,
                                        const TopoDS_Vertex&  CV,
                                        TopoDS_Edge&          NE,
                                        const TopTools_ListOfShape& LE)
    {
      NE.Nullify();
    
      if (LE.Extent() > 1) {
        //--------------------------------------------------------------
        // Several possible edges.
        // - Test the edges differents of CE
        //--------------------------------------------------------------
        TopoDS_Face FForward = Surf.Face();
        TopoDS_Edge aPrevNE;
    
        gp_Vec2d CTg1, Tg1, CTg2, Tg2;
        gp_Pnt2d PC, P;
    
        Standard_Real f, l;
        Handle(Geom2d_Curve) Cc, C;
        Cc = BRep_Tool::CurveOnSurface(CE,FForward,f,l);
    
        Standard_Boolean isForward = ( CE.Orientation () == TopAbs_FORWARD );
        Standard_Real uc, u, du = Precision::PConfusion();
        uc = isForward ? ( l - du ) : ( f + du );
        Cc->D1(uc, PC, CTg1);
        if (!isForward) CTg1.Reverse();
    
        Standard_Real anglemin = 3 * PI, tolAng = 1.e-8;
    
        // select an edge whose first derivative is most left of CTg1
        // ie an angle between Tg1 and CTg1 is least
        TopTools_ListIteratorOfListOfShape itl;
        for ( itl.Initialize(LE); itl.More(); itl.Next()) {
          const TopoDS_Edge& E = TopoDS::Edge(itl.Value());
          if (E.IsSame(CE))
            continue;
          if (! CV.IsSame( TopExp::FirstVertex( E, Standard_True )))
            continue;
    
          isForward = ( E.Orientation () == TopAbs_FORWARD );
    
          // get E curve
          C = BRep_Tool::CurveOnSurface(E,FForward,f,l);
          // get the first derivative Tg1
          u = isForward ? ( f + du ) : ( l - du );
          C->D1(u, P, Tg1);
          if (!isForward) Tg1.Reverse();
    
          // -PI < angle < PI
          Standard_Real angle = Tg1.Angle(CTg1);
    
          if (PI - Abs(angle) <= tolAng)
          {
            // an angle is too close to PI; assure that an angle sign really
            // reflects an edge position: +PI - an edge is worst,
            // -PI - an edge is best.
            u = preciseU( Surf, CE, CV, Cc, Standard_False);
            gp_Vec2d CTg;
            Cc->D1(u, PC, CTg);
            if (CE.Orientation() == TopAbs_REVERSED) CTg.Reverse();
    
            u = preciseU( Surf, E, CV, C, Standard_True);
            C->D1(u, P, Tg1);
            if (!isForward) Tg1.Reverse();
    
            angle = Tg1.Angle(CTg);
          }
    
          Standard_Boolean isClose = ( Abs( angle - anglemin ) <= tolAng );
          if (angle <= anglemin) {
            if (isClose)
              aPrevNE = NE;
            else
              aPrevNE.Nullify();
            anglemin = angle ;
            NE = E;
          }
          else
            if (isClose)
              aPrevNE = E;
    
        }
        if (!aPrevNE.IsNull()) {
          // select one of close edges, the most left one.
          Cc = BRep_Tool::CurveOnSurface( NE, FForward, f, l );
          uc = preciseU( Surf, NE, CV, Cc, Standard_True);
          Cc->D1(uc, PC, CTg1);
          if (NE.Orientation() != TopAbs_FORWARD) CTg1.Reverse();
          
          u = preciseU( Surf, aPrevNE, CV, C, Standard_True);
          C->D1(u, P, Tg1);
          if (aPrevNE.Orientation() != TopAbs_FORWARD) Tg1.Reverse();
    
          if ( Tg1.Angle(CTg1) < 0)
            NE = aPrevNE;
        }
      }
      else if (LE.Extent() == 1) {
        NE = TopoDS::Edge(LE.First());
      }
      else {
        return Standard_False;
      }
      return !NE.IsNull();
    }
    
    //=======================================================================
    //function : SamePnt2d
    //purpose  :
    //=======================================================================
    
    static Standard_Boolean  SamePnt2d(const TopoDS_Vertex& V1,
                                       const TopoDS_Edge&   E1,
                                       const TopoDS_Vertex& V2,
                                       const TopoDS_Edge&   E2,
                                       const TopoDS_Face&   F)
    {
      Standard_Real   f1,f2,l1,l2;
      Handle(Geom2d_Curve) C1 = BRep_Tool::CurveOnSurface(E1,F,f1,l1);
      Handle(Geom2d_Curve) C2 = BRep_Tool::CurveOnSurface(E2,F,f2,l2);
    
      gp_Pnt2d P1 = C1->Value( BRep_Tool::Parameter(V1,E1));
      gp_Pnt2d P2 = C2->Value( BRep_Tool::Parameter(V2,E2));
    
      Standard_Real Tol  = 100 * BRep_Tool::Tolerance(V1);
      Standard_Real Dist = P1.Distance(P2);
      return Dist < Tol;
    }
    
    
    //=======================================================================
    //function : StoreInMVE
    //purpose  :
    //=======================================================================
    
    static void StoreInMVE (const TopoDS_Face& /*F*/,
                            TopoDS_Edge& E,
                            TopTools_DataMapOfShapeListOfShape& MVE )
    
    {
      TopoDS_Vertex V1, V2;
      TopTools_ListOfShape Empty;
    
      TopExp::Vertices(E,V1,V2);
      if (!MVE.IsBound(V1)) {
        MVE.Bind(V1,Empty);
      }
      MVE(V1).Append(E);
    
      if (!MVE.IsBound(V2)) {
        MVE.Bind(V2,Empty);
      }
      MVE(V2).Append(E);
    }
    
    //=======================================================================
    //function : RemoveFromMVE
    //purpose  :
    //=======================================================================
    
    static void RemoveFromMVE(const TopoDS_Edge& E,
                              TopTools_DataMapOfShapeListOfShape& MVE)
    {
      TopTools_ListIteratorOfListOfShape itl;
      TopoDS_Vertex  V1,V2;
      TopExp::Vertices (E,V1,V2);
      if (MVE.IsBound(V1))
        for ( itl.Initialize(MVE(V1)); itl.More(); itl.Next()) {
          if (itl.Value().IsEqual(E)) {
            MVE(V1).Remove(itl);
            break;
          }
        }
      if (MVE.IsBound(V2))
        for ( itl.Initialize(MVE(V2)); itl.More(); itl.Next()) {
          if (itl.Value().IsEqual(E)) {
            MVE(V2).Remove(itl);
            break;
          }
        }
    }
    //=======================================================================
    //function : addConnected
    //purpose  : add to <EM> all edges reachable from <E>
    //=======================================================================
    
    static void addConnected(const TopoDS_Shape& E,
                             TopTools_MapOfShape& EM,
                             TopTools_MapOfShape& VM,
                             const TopTools_DataMapOfShapeListOfShape& MVE)
    {
      // Loop on vertices of E
      TopoDS_Iterator itV ( E );
      for ( ; itV.More(); itV.Next()) {
    
        if ( ! VM.Add ( itV.Value() )) continue;
    
        // Loop on edges sharing V
        TopTools_ListIteratorOfListOfShape itE( MVE( itV.Value() ) );
        for (; itE.More(); itE.Next()) {
          if ( EM.Add( itE.Value() ))
            addConnected ( itE.Value(), EM, VM, MVE );
        }
      }
    }
    //=======================================================================
    //function : canPassToOld
    //purpose  :
    //=======================================================================
    
    // static Standard_Boolean canPassToOld (const TopoDS_Shape& V,
    //                                    TopTools_MapOfShape& UsedShapesMap,
    //                                    const TopTools_DataMapOfShapeListOfShape& MVE,
    //                                    const TopTools_MapOfShape& SectionEdgesMap)
    // {
    //   TopTools_ListIteratorOfListOfShape itE( MVE(V) );
    //   // Loop on edges sharing V
    //   for (; itE.More(); itE.Next()) {
    //     if ( !UsedShapesMap.Add( itE.Value() ))
    //       continue; // already checked
    
    //     if ( !SectionEdgesMap.Contains( itE.Value() ))
    //       return Standard_True; // WE PASSED
    
    //     TopoDS_Iterator itV( itE.Value() );
    //     // Loop on vertices of an edge
    //     for (; itV.More(); itV.Next()) {
    //       if ( !UsedShapesMap.Add( itV.Value() ))
    //      continue; // already checked
    //       else
    //      return canPassToOld( itV.Value(), UsedShapesMap, MVE, SectionEdgesMap);
    //     }
    //   }
    //   return Standard_False;
    // }
    
    //=======================================================================
    //function : MakeDegenAndSelect
    //purpose  : Find parameter of intersection of <CE> with <DE> and
    //           select an edge with its parameter closest to found one.
    //           Return new degenerated edge trimming <DE> by found parameters
    //=======================================================================
    
    static TopoDS_Edge MakeDegenAndSelect(const TopoDS_Edge& CE,
                                          const TopoDS_Vertex& CV,
                                          TopoDS_Edge& NE,
                                          TopTools_SequenceOfShape& EdgesSeq,
                                          TColStd_SequenceOfReal& USeq,
                                          const TopoDS_Edge& DE)
    {
      if (EdgesSeq.Length() < 3) {
        if (CE == EdgesSeq.First())
          NE = TopoDS::Edge( EdgesSeq.Last() );
        else
          NE = TopoDS::Edge( EdgesSeq.First() );
        return DE;
      }
    
      // find parameter on DE where it intersects CE
    
      Standard_Real U1;
      Standard_Integer i, nb = EdgesSeq.Length();
      for (i=1; i<= nb; ++i) {
        if (CE == EdgesSeq(i)) {
          U1 = USeq(i);
          break;
        }
      }
    
      // select NE with param closest to U1 thus finding U2 for a new degen edge
    
      Standard_Real U2, dU, dUmin = 1.e100;
      Standard_Boolean isReversed = ( DE.Orientation() == TopAbs_REVERSED );
      for (i=1; i<= nb; ++i) {
        dU = USeq(i) - U1;
        if (isReversed ? (dU > 0) : (dU < 0))
            continue;
        dU = Abs( dU );
        if ( dU  > dUmin || IsEqual( dU, 0.))
          continue;
        const TopoDS_Edge& E = TopoDS::Edge ( EdgesSeq(i) );
        if ( ! CV.IsSame( TopExp::FirstVertex( E , Standard_True )))
          continue;
        NE = E;
    #if defined(WIN32) // fix for OCE on win32
        dUmin = dU + (dU+1e-12);
    #else
        dUmin = dU + Epsilon(dU);
    #endif
        U2 = USeq(i);
      }
    
      // make a new degenerated edge
      TopoDS_Edge NewDegen = TopoDS::Edge ( DE.EmptyCopied() );
    
      Standard_Real Tol = BRep_Tool::Tolerance( CV );
      TopoDS_Vertex V = CV;
    
      BRep_Builder B;
      V.Orientation( NewDegen.Orientation() );
      B.UpdateVertex( V, U1, NewDegen, Tol);
      B.Add ( NewDegen , V );
    
      V.Reverse();
      B.UpdateVertex( V, U2, NewDegen, Tol);
      B.Add ( NewDegen , V );
    
      return NewDegen;
    }
    
    //=======================================================================
    //function : prepareDegen
    //purpose  : Intersect <DegEdge> with edges bound to its vertex in <MVE>
    //           and store intersection parameter on <DegEdge> in
    //           <USeq> as well as the edges them-self in <EdgesSeq>.
    //           Bind <DegEdgeIndex> to vertex of <DegEdge> in <MVDEI>
    //=======================================================================
    
    static void prepareDegen (const TopoDS_Edge&                        DegEdge,
                              const TopoDS_Face&                        F,
                              const TopTools_DataMapOfShapeListOfShape& MVE,
                              TopTools_SequenceOfShape&                 EdgesSeq,
                              TColStd_SequenceOfReal&                   USeq,
                              TopTools_DataMapOfShapeInteger&           MVDEI,
                              const Standard_Integer                    DegEdgeIndex)
    {
      const TopoDS_Vertex& V = TopExp::FirstVertex ( DegEdge );
      MVDEI.Bind ( V, DegEdgeIndex );
    
      const TopTools_ListOfShape& EdgesList = MVE ( V );
      // if only 2 edges come to degenerated one, no pb in selection and
      // no need to intersect them, just simulate asked data
      Standard_Boolean doIntersect =  ( EdgesList.Extent() > 2 );
    
      BRepAdaptor_Curve2d DC, C;
      Geom2dInt_GInter InterCC;
      Standard_Real Tol = Precision::PConfusion();
      if ( doIntersect )
        DC.Initialize( DegEdge, F );
    
      // avoid intersecting twice the same edge
    #if (OCC_VERSION_MAJOR == 6) && (OCC_VERSION_MINOR < 5)
      BRepOffset_DataMapOfShapeReal EUMap ( EdgesList.Extent() );
    #else
      TopTools_DataMapOfShapeReal EUMap ( EdgesList.Extent() );
    #endif
    
      Standard_Real U, f, l;
      BRep_Tool::Range (DegEdge, f, l);
    
      TopTools_ListIteratorOfListOfShape itE (EdgesList);
      for (; itE.More(); itE.Next()) {
    
        const TopoDS_Edge& E = TopoDS::Edge ( itE.Value() );
    
        if ( !doIntersect) {
          U = 0.; // it won't be used
        }
        else if ( BRep_Tool::IsClosed( E, F )) {
          // seam edge: select U among f and l
          Standard_Boolean first = Standard_True;
          if ( V.IsSame ( TopExp::FirstVertex( E, Standard_True ) ))
            first = Standard_False;
          if ( DegEdge.Orientation() == TopAbs_REVERSED )
            first = !first;
          U = first ? f : l;
        }
        else if ( EUMap.IsBound( E ) ) {
          // same edge already bound
          U = EUMap( E );
        }
        else {
          // intersect 2d curves
          C.Initialize( E, F );
          InterCC.Perform ( DC, C , Tol, Tol );
          if (! InterCC.IsDone() || InterCC.NbPoints() == 0) {
            //MESSAGE ( "NO 2d INTERSECTION ON DEGENERATED EDGE" );
            continue;
          }
          // hope there is only one point of intersection
          U = InterCC.Point( 1 ).ParamOnFirst();
        }
        USeq.Append ( U );
        EdgesSeq.Append ( E );
      }
    }
    //=======================================================================
    //function : Perform
    //purpose  : Make loops.
    //=======================================================================
    
    void Partition_Loop2d::Perform()
    {
    
      Standard_Integer NbConstEdges = myConstEdges.Extent();
      TopTools_DataMapOfShapeListOfShape MVE(NbConstEdges) , MVE2(NbConstEdges);
      TopTools_DataMapIteratorOfDataMapOfShapeListOfShape Mapit;
      TopTools_ListIteratorOfListOfShape itl;
      TopoDS_Vertex V1,V2;
      BRepAdaptor_Surface Surface ( myFace, Standard_False );
    
      // degenerated edges and parameters of their 2d intersection with other edges
      TopoDS_Edge                    DE [2];
      TopTools_SequenceOfShape       SEID [2]; // seq of edges intersecting degenerated
      TColStd_SequenceOfReal         SeqU [2]; // n-th U corresponds to n-th edge in SEID
      TopTools_DataMapOfShapeInteger MVDEI(2); // map vertex - degenerated edge index
      Standard_Integer               iDeg = 0; // index of degenerated edge [0,1]
    
      //---------------------------------------------------------
      // Construction map vertex => edges, find degenerated edges
      //---------------------------------------------------------
      for (itl.Initialize(myConstEdges); itl.More(); itl.Next()) {
        TopoDS_Edge& E = TopoDS::Edge(itl.Value());
        if ( BRep_Tool::Degenerated( E )) {
          if (DE[0].IsNull()) DE[0] = E;
          else                DE[1] = E;
        }
        else
          StoreInMVE(myFace,E,MVE);
      }
    
      // fill data for degenerated edges
      if ( ! DE[0].IsNull() )
        prepareDegen ( DE[0], myFace, MVE, SEID[0], SeqU[0], MVDEI, 0);
      if ( ! DE[1].IsNull() )
        prepareDegen ( DE[1], myFace, MVE, SEID[1], SeqU[1], MVDEI, 1);
    
    
      // to detect internal wires
      Standard_Boolean isInternCW = 0;
      MVE2 = MVE;
    
    
      //------------------------------
      // Construction of all the wires
      //------------------------------
      // first, we collect wire edges in WEL list looking for same edges that
      // will be then removed possibly exploding a wire into parts;
      // second, build wire(s)
    
      while (!MVE.IsEmpty()) {
    
        TopoDS_Vertex    VF,CV;
        TopoDS_Edge      CE,NE,EF;
        TopoDS_Wire      NW;
        BRep_Builder     B;
        Standard_Boolean End = Standard_False;
        TopTools_ListOfShape WEL;
    
        Mapit.Initialize(MVE);
        if (Mapit.Value().IsEmpty()) {
          MVE.UnBind(Mapit.Key());
          continue;
        }
    
        // EF first edge.
        EF = CE = TopoDS::Edge(Mapit.Value().First());
        // VF first vertex
        VF = TopExp::FirstVertex( CE, Standard_True);
    
        isInternCW = Standard_True;
    
        TopTools_MapOfShape addedEM  (NbConstEdges); // map of edges added to WEL
        TopTools_MapOfShape doubleEM (NbConstEdges); // edges encountered twice in WEL
    
        //-------------------------------
        // Construction of a wire.
        //-------------------------------
        while (!End) {
    
          // only a seam is allowed twice in a wire, the others should be removed
          if (addedEM.Add ( CE ) || BRep_Tool::IsClosed( CE, myFace ) )
            WEL.Append( CE );
          else {
            doubleEM.Add( CE );
            RemoveFromMVE (CE,MVE2);
            TopoDS_Edge CERev = CE;
            CERev.Reverse();
            RemoveFromMVE (CERev,MVE2);
          }
    
          RemoveFromMVE (CE,MVE);
    
          CV = TopExp::LastVertex( CE, Standard_True);
    
          if (isInternCW && !mySectionEdges.Contains(CE))
            // wire is internal if all edges are section ones
            isInternCW = Standard_False;
    
          if (MVDEI.IsBound( CV )) { // CE comes to the degeneration
            iDeg = MVDEI( CV );
            TopoDS_Edge NewDegen;
            NewDegen = MakeDegenAndSelect( CE, CV, NE, SEID[iDeg], SeqU[iDeg], DE[iDeg]);
            WEL.Append( NewDegen );
            CE = NE;
            End = CV.IsSame( VF );
            continue;
          }
    
          //--------------
          // stop test
          //--------------
          if (MVE(CV).IsEmpty()) {
            End=Standard_True;
            MVE.UnBind(CV);
          }
          else if (CV.IsSame(VF) && SamePnt2d(CV,CE, VF,EF, myFace) ) {
            End = Standard_True;
          }
          else {
            //----------------------------
            // select new current edge
            //----------------------------
            if (! SelectEdge (Surface,CE,CV,NE,MVE(CV))) {
              //MESSAGE ( " NOT CLOSED WIRE " );
              End=Standard_True;
            }
            else
              CE = NE;
          }
        } // while ( !End )
    
    
        // WEL is built, built wire(s)
    
    
        itl.Initialize( WEL );
        if ( doubleEM.IsEmpty()) { // no double edges
          B.MakeWire( NW );
          for (; itl.More(); itl.Next())
            B.Add ( NW, itl.Value());
          if (isInternCW) myInternalWL.Append(NW);
          else            myNewWires.Append  (NW);
        }
    
        else {
          // remove double and degenerated edges from WEL
          while (itl.More()) {
            const TopoDS_Edge& E = TopoDS::Edge ( itl.Value() );
            if ( doubleEM.Contains( E ) || BRep_Tool::Degenerated( E ))
              WEL.Remove( itl );
            else
               itl.Next();
          }
          if ( WEL.IsEmpty())
            continue;
          // remove double edges from SEID and SeqU
          Standard_Integer i,j;
          for (j=0; j<2; ++j) {
            for (i=1; i<=SEID[j].Length(); ++i) {
              if (doubleEM.Contains( SEID[j].Value(i))) {
                SEID[j].Remove( i );
                SeqU[j].Remove( i-- );
              }
            }
          }
          // removal of doulbe edges can explode a wire into parts,
          // make new wires of them.
          // A Loop like previous one but without 2d check
          while ( !WEL.IsEmpty() ) {
            CE = TopoDS::Edge( WEL.First() );
            WEL.RemoveFirst();
            B.MakeWire( NW );
            VF = TopExp::FirstVertex ( CE, Standard_True);
    
            End = Standard_False;
            while ( !End) {
              B.Add( NW, CE );
              CV = TopExp::LastVertex  ( CE, Standard_True);
    
              if (MVDEI.IsBound( CV )) {   // CE comes to the degeneration
                iDeg = MVDEI( CV );
                TopoDS_Edge NewDegen;
                NewDegen = MakeDegenAndSelect( CE, CV, NE, SEID[iDeg], SeqU[iDeg], DE[iDeg]);
                B.Add( NW, NewDegen );
                End = CV.IsSame( VF );
                CE = NE;
                if (!NE.IsNull()) { // remove NE from WEL
                  for (itl.Initialize( WEL ); itl.More(); itl.Next())
                    if ( NE == itl.Value()) {
                      WEL.Remove( itl );
                      break;
                    }
                }
              }  // end degeneration
    
              else {
                if (CV.IsSame( VF )) {
                  End = Standard_True;
                  continue;
                }
                // edges in WEL most often are well ordered
                // so try to iterate until the End
                Standard_Boolean add = Standard_False;
                itl.Initialize(WEL);
                while ( itl.More() && !End) {
                  NE = TopoDS::Edge( itl.Value() );
                  if ( CV.IsSame( TopExp::FirstVertex( NE, Standard_True ))) {
                    WEL.Remove( itl );
                    if (add)
                      B.Add( NW, CE );
                    CE = NE;
                    add = Standard_True;
                    CV = TopExp::LastVertex( CE, Standard_True);
                    if (MVDEI.IsBound( CV ) || CV.IsSame( VF ))
                      break;
                  }
                  else
                    itl.Next();
                }
                if (!add)
                  End = Standard_True;
              }
            } // !End
    
            myInternalWL.Append( NW );
          }
        } // end building new wire(s) from WEL
    
      } // end Loop on MVE
    
      // all wires are built
    
    
      // ============================================================
      // select really internal wires i.e. those from which we can`t
      // pass to an old (not section) edge
      // ============================================================
    
      Standard_Integer nbIW = myInternalWL.Extent();
      if (nbIW == 0)
        return;
    
      if ( myNewWires.Extent() != 1 && nbIW > 1) {
        TopTools_MapOfShape outerEM (NbConstEdges); // edges connected to non-section ones
        TopTools_MapOfShape visitedVM (NbConstEdges);
        for ( itl.Initialize( myConstEdges ); itl.More(); itl.Next()) {
          if ( ! mySectionEdges.Contains( itl.Value() ))
            addConnected (itl.Value(), outerEM, visitedVM, MVE2);
        }
        // if an edge of a wire is in <outerEM>, the wire is not internal
        TopExp_Explorer expIWE;
        TopTools_ListIteratorOfListOfShape itIW ( myInternalWL );
        while (itIW.More()) {
          expIWE.Init ( itIW.Value() , TopAbs_EDGE );
          if ( outerEM.Contains( expIWE.Current() )) {
            myNewWires.Append ( itIW.Value() );
            myInternalWL.Remove( itIW ); // == itIW.Next()
          }
          else
            itIW.Next();
        }
      }
    }
    //=======================================================================
    //function : isHole
    //purpose  :
    //=======================================================================
    
    static Standard_Boolean isHole (const TopoDS_Wire& W,
                                    const TopoDS_Face& F)
    {
      BRep_Builder B;
      TopoDS_Shape newFace = F.EmptyCopied();
      B.Add(newFace,W.Oriented(TopAbs_FORWARD));
      BRepTopAdaptor_FClass2d classif (TopoDS::Face(newFace),
                                       Precision::PConfusion());
      return (classif.PerformInfinitePoint() == TopAbs_IN);
    }
    
    //=======================================================================
    //function : IsInside
    //purpose  : check if W1 is inside W2. Suppose W2 is not a hole !!!!
    //=======================================================================
    
    static Standard_Boolean isInside(const TopoDS_Face& F,
                                     const TopoDS_Wire& W1,
                                     const TopoDS_Wire& W2)
    {
      // make a face with wire W2
      BRep_Builder B;
      TopoDS_Shape aLocalShape = F.EmptyCopied();
      TopoDS_Face newFace = TopoDS::Face(aLocalShape);
      B.Add(newFace,W2);
    
      // get any 2d point of W1
      TopExp_Explorer exp(W1,TopAbs_EDGE);
      if (BRep_Tool::Degenerated( TopoDS::Edge( exp.Current() )))
        exp.Next();
      const TopoDS_Edge& e = TopoDS::Edge(exp.Current());
      Standard_Real f,l;
      Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface(e,F,f,l);
      gp_Pnt2d pt2d(C2d->Value( 0.5 * ( f + l )));
    
      BRepTopAdaptor_FClass2d classif(newFace,Precision::PConfusion());
      return (classif.Perform(pt2d) == TopAbs_IN);
    }
    
    //=======================================================================
    //function : NewWires
    //purpose  : Returns the list of wires performed.
    //           can be an empty list.
    //=======================================================================
    
    const TopTools_ListOfShape&  Partition_Loop2d::NewWires() const
    {
      return myNewWires;
    }
    
    //=======================================================================
    //function : NewFaces
    //purpose  : Returns the list of faces.
    //Warning  : The method <WiresToFaces> as to be called before.
    //           can be an empty list.
    //=======================================================================
    
    const TopTools_ListOfShape&  Partition_Loop2d::NewFaces() const
    {
      return myNewFaces;
    }
    
    //=======================================================================
    //function : findEqual
    //purpose  : move wires form <WL> to <EqWL> pairs of wires build of the same edges
    //=======================================================================
    
    static void findEqual (TopTools_ListOfShape& WL,
                           TopTools_DataMapOfShapeShape& EqWM,
                           const TopoDS_Face& F)
    {
      TopTools_ListIteratorOfListOfShape it1, it2;
      Standard_Integer i,j;
      TColStd_MapOfInteger IndMap;
      for (it1.Initialize(WL), i=1;  it1.More();  it1.Next(), i++) {
    
        if (IndMap.Contains(i)) continue;
        const TopoDS_Wire& Wire1 = TopoDS::Wire( it1.Value());
    
        for (it2.Initialize(WL), j=1;  it2.More();  it2.Next(), j++) {
    
          if (j <= i || IndMap.Contains(j)) continue;
    
          TopTools_IndexedMapOfShape EdgesMap;
          TopExp::MapShapes (Wire1, TopAbs_EDGE, EdgesMap);
    
          const TopoDS_Shape& Wire2 = it2.Value();
          TopoDS_Iterator itE ( Wire2);
          for (; itE.More(); itE.Next()) {
            if ( !EdgesMap.Contains( itE.Value()) )
              break;
          }
          if (!itE.More()) { // all edges are same
            if (isHole( Wire1, F)) {
              EqWM.Bind ( Wire1, Wire2 );
            }
            else {
              EqWM.Bind ( Wire2, Wire1 );
            }
            IndMap.Add(i);
            IndMap.Add(j);
            break;
          }
        }
      }
      // clear WL
      it1.Initialize(WL);
      i=1;
      while (it1.More()) {
        if (IndMap.Contains(i))
          WL.Remove(it1); // next node becomes current and with Next() we would miss it
        else
          it1.Next();
        i++;
      }
    }
    
    //=======================================================================
    //function : classify
    //purpose  : bind to a wire a list of internal wires
    //=======================================================================
    
    static void classify(const TopTools_DataMapOfShapeShape& EqWM,
                         BRepAlgo_AsDes& OuterInner,
                         const TopoDS_Face& F)
    {
      TopTools_DataMapIteratorOfDataMapOfShapeShape it1, it2;
    
      for (it1.Initialize(EqWM);  it1.More();  it1.Next()) {
        // find next after it1.Value()
        for (it2.Initialize(EqWM);  it2.More();  it2.Next())
          if (it1.Value().IsSame( it2.Value() ))
          {
            it2.Next();
            break;
          }
        for ( ;  it2.More();  it2.Next()) {
          const TopoDS_Wire& Wire1 = TopoDS::Wire( it1.Value() );
          const TopoDS_Wire& Wire2 = TopoDS::Wire( it2.Value() );
          if (isInside(F, Wire1, Wire2))
            OuterInner.Add (Wire2, Wire1);
          else if (isInside(F, Wire2, Wire1))
            OuterInner.Add (Wire1, Wire2);
        }
      }
    }
    //=======================================================================
    //function : WiresToFaces
    //purpose  : Build faces from the wires result.
    //           <EdgeImage> serves to  find  original edge by new
    //           one. <Section> contains edges resulting from face
    //           intersections
    //=======================================================================
    
    void  Partition_Loop2d::WiresToFaces(const BRepAlgo_Image& )
    {
      Standard_Integer nbW = myNewWires.Extent() + myInternalWL.Extent();
      if (nbW==0)
        return;
    
      BRepAlgo_FaceRestrictor FR;
      FR.Init (myFace,Standard_False);
    
      // FaceRestrictor is instable in rather simple cases
      // (ex. a single face of bellecoque.brep splited by 10 planes:
      // sometimes 1-2 faces are missing ).
      // So we use it as less as possible: no holes -> make faces by hands
    
    
      // are there holes in myFace ?
      Standard_Boolean hasOldHoles = Standard_False;
      TopoDS_Iterator itOldW (myFace);
      if ( itOldW.More()) {
        const TopoDS_Wire& FirstOldWire = TopoDS::Wire( itOldW.Value() );
        itOldW.Next();
        hasOldHoles = itOldW.More() || isHole( FirstOldWire, myFace);
      }
      if (myInternalWL.IsEmpty() && !hasOldHoles) {
        // each wire bounds one face
        BRep_Builder B;
        TopTools_ListIteratorOfListOfShape itNW (myNewWires);
        for (; itNW.More(); itNW.Next()) {
          TopoDS_Face NF = TopoDS::Face ( myFace.EmptyCopied() );
          B.Add ( NF, itNW.Value() );
          NF.Orientation( myFaceOri);
          myNewFaces.Append ( NF );
        }
        return;
      }
    
      // FaceRestrictor can't classify wires build on all the same edges
      // and gives incorrect result in such cases (ex. a plane cut into 2 parts by cylinder)
      // We must make faces of equal wires separately. One of equal wires makes a
      // hole in a face and should come together with outer wires of face.
      // The other of a wires pair bounds a face that may have holes in turn.
    
      // Find equal wires among internal wires
      TopTools_DataMapOfShapeShape EqWM; // key is a hole part of a pair of equal wires
      findEqual (myInternalWL, EqWM, myFace);
    
      if (!EqWM.IsEmpty()) { // there are equal wires
    
        if (hasOldHoles)
          myInternalWL.Append( myNewWires ); // an old wire can be inside an equal wire
    
        // classify equal wire pairs
        BRepAlgo_AsDes OuterInner;
        classify (EqWM,OuterInner,myFace);
    
        // make face of most internal of equal wires and its inner wires
        while ( !EqWM.IsEmpty()) {
    
          TopTools_ListOfShape prevHolesL; // list of hole-part of previous most internal equal wires
    
          // find most internal wires among pairs (key - hole, value - outer part)
          TopTools_DataMapIteratorOfDataMapOfShapeShape it(EqWM);
          Standard_Integer nbEqW = EqWM.Extent(); // protection against infinite loop
          for ( ; it.More(); it.Next()) {
    
            TopoDS_Wire outerW = TopoDS::Wire ( it.Value() );
            if (  OuterInner.HasDescendant( outerW ) && // has internal
                 ! OuterInner.Descendant( outerW ).IsEmpty() )
              continue;
    
            FR.Add( outerW );
    
            // add internal wires that are inside of outerW
            TopTools_ListIteratorOfListOfShape itIW (myInternalWL);
            while ( itIW.More()) {
              TopoDS_Wire IW = TopoDS::Wire ( itIW.Value() );
              if ( isInside (myFace, IW, outerW)) {
                FR.Add (IW);
                myInternalWL.Remove( itIW ); // == itIW.Next() !!!
              }
              else
                itIW.Next();
            }
    
            // the hole-part of current pair of equal wires will be in the next new face
            prevHolesL.Append ( it.Key() );
    
          } // Loop on map of equal pairs searching for innermost wires
    
          // make faces
          FR.Perform();
          if (FR.IsDone()) {
            for (; FR.More(); FR.Next())
              myNewFaces.Append(FR.Current());
          }
    
          FR.Clear();
    
          // add hole-parts to FaceRestrictor,
          // remove them from the EqWM,
          // remove found wires as internal of resting classified wires
          Standard_Boolean clearOuterInner =  ( prevHolesL.Extent() < EqWM.Extent() );
          TopTools_ListIteratorOfListOfShape itPrev (prevHolesL);
          for (; itPrev.More(); itPrev.Next()) {
            TopoDS_Wire& Hole = TopoDS::Wire ( itPrev.Value() );
            FR.Add ( Hole );
            if (clearOuterInner) {
              const TopoDS_Wire& outerW = TopoDS::Wire ( EqWM.Find( Hole ) );
              // Loop on wires including outerW
              TopTools_ListIteratorOfListOfShape itO( OuterInner.Ascendant( outerW ));
              for (; itO.More(); itO.Next()) {
                TopTools_ListOfShape& innerL = OuterInner.ChangeDescendant( itO.Value() );
                TopTools_ListIteratorOfListOfShape itI (innerL);
                // Loop on internal wires of current including wire
                for (; itI.More(); itI.Next())
                  if ( outerW.IsSame( itI.Value() )) {
                    innerL.Remove( itI );   break;
                  }
              }
            }
            EqWM.UnBind ( Hole );
          }
    
          if (nbEqW == EqWM.Extent())
          {
            // error: pb with wires classification
    #ifdef DEB
            MESSAGE("Partition_Loop2d::WiresToFaces(), pb with wires classification");
    #endif
            break;
          }
    
        } // while (!EqWM.IsEmpty)
    
      } //  if !EqWM.IsEmpty()
    
      myNewWires.Append ( myInternalWL );
    
      TopTools_ListIteratorOfListOfShape itW (myNewWires);
      for (; itW.More(); itW.Next()) {
        TopoDS_Wire& W = TopoDS::Wire ( itW.Value() );
        FR.Add(W);
      }
      FR.Perform();
      for (; FR.IsDone() && FR.More(); FR.Next())
        myNewFaces.Append(FR.Current());
    
    
      TopTools_ListIteratorOfListOfShape itNF (myNewFaces);
      for (; itNF.More(); itNF.Next())
        itNF.Value().Orientation( myFaceOri );
    }
    
    #endif