src/PARTITION/Partition_Loop2d.cxx

   1 //  GEOM PARTITION : partition algorithm
   2 //
   3 //  Copyright (C) 2003  CEA/DEN, EDF R& D
   4 //
   5 //
   6 //
   7 //  File   : Partition_Loop2d.cxx
   8 //  Author : Benedicte MARTIN
   9 //  Module : GEOM
  10 //  $Header$
  11
  12 using namespace std;
  13 #include "Partition_Loop2d.ixx"
  14
  15 #include "utilities.h"
  16 #include <stdio.h>
  17
  18 #include <BRepAdaptor_Curve2d.hxx>
  19 #include <BRepAdaptor_Surface.hxx>
  20 #include <BRepAlgo_AsDes.hxx>
  21 #include <BRepAlgo_FaceRestrictor.hxx>
  22 #include <BRepOffset_DataMapOfShapeReal.hxx>
  23 #include <BRepTopAdaptor_FClass2d.hxx>
  24 #include <BRep_Builder.hxx>
  25 #include <BRep_Tool.hxx>
  26 #include <Geom2dInt_GInter.hxx>
  27 #include <Geom2d_Curve.hxx>
  28 #include <IntRes2d_IntersectionPoint.hxx>
  29 #include <Precision.hxx>
  30 #include <TColStd_MapOfInteger.hxx>
  31 #include <TColStd_SequenceOfReal.hxx>
  32 #include <TopExp.hxx>
  33 #include <TopExp_Explorer.hxx>
  34 #include <TopTools_DataMapIteratorOfDataMapOfShapeListOfShape.hxx>
  35 #include <TopTools_DataMapIteratorOfDataMapOfShapeShape.hxx>
  36 #include <TopTools_DataMapOfShapeInteger.hxx>
  37 #include <TopTools_DataMapOfShapeShape.hxx>
  38 #include <TopTools_IndexedMapOfShape.hxx>
  39 #include <TopTools_ListIteratorOfListOfShape.hxx>
  40 #include <TopTools_MapIteratorOfMapOfShape.hxx>
  41 #include <TopTools_MapOfOrientedShape.hxx>
  42 #include <TopTools_MapOfShape.hxx>
  43 #include <TopTools_SequenceOfShape.hxx>
  44 #include <TopoDS.hxx>
  45 #include <TopoDS_Iterator.hxx>
  46 #include <TopoDS_Vertex.hxx>
  47 #include <TopoDS_Wire.hxx>
  48 #include <gp_Pnt.hxx>
  49 #include <gp_Pnt2d.hxx>
  50
  51 //=======================================================================
  52 //function : Partition_Loop2d
  53 //purpose  :
  54 //=======================================================================
  55
  56 Partition_Loop2d::Partition_Loop2d()
  57 {
  58 }
  59
  60 //=======================================================================
  61 //function : Init
  62 //purpose  : Init with <F> the set of edges must have
  63 //           pcurves on <F>.
  64 //=======================================================================
  65
  66 void Partition_Loop2d::Init(const TopoDS_Face& F)
  67 {
  68   myConstEdges.Clear();
  69   myNewWires  .Clear();
  70   myNewFaces  .Clear();
  71   myFace = F;
  72   myFaceOri = myFace.Orientation();
  73   myFace.Orientation( TopAbs_FORWARD );
  74 }
  75
  76 //=======================================================================
  77 //function : AddConstEdge
  78 //purpose  : Add <E> as unique edge in the result.
  79 //=======================================================================
  80
  81 void Partition_Loop2d::AddConstEdge (const TopoDS_Edge& E)
  82 {
  83 #ifdef DEB
  84   Standard_Real f,l;
  85   Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( E, myFace, f,l);
  86   if (pc.IsNull()) {
  87     INFOS( "AddConstEdge(): EDGE W/O PCURVE on FACE");
  88   } else
  89 #endif
  90   {
  91     myConstEdges.Append(E);
  92   }
  93 }
  94
  95 void Partition_Loop2d::AddSectionEdge (const TopoDS_Edge& E)
  96 {
  97 #ifdef DEB
  98   Standard_Real f,l;
  99   Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( E, myFace, f,l);
 100   if (pc.IsNull())
 101     pc = BRep_Tool::CurveOnSurface( E, myFace, f,l);
 102   gp_Vec2d Tg1;
 103   gp_Pnt2d PC;
 104   pc->D1(0.5*(f+l), PC, Tg1);
 105   if (Tg1.Magnitude()  <= gp::Resolution()) {
 106     MESSAGE ("");
 107   }
 108   if (pc.IsNull()) {
 109     INFOS( "AddConstEdge(): EDGE W/O PCURVE on FACE");
 110   } else
 111 #endif
 112   {
 113     myConstEdges.Append(E);
 114     myConstEdges.Append(E.Reversed());
 115     mySectionEdges.Add( E );
 116   }
 117 }
 118
 119 //=======================================================================
 120 //function : preciseU
 121 //purpose  : find u such that the 3D point on theE is just out of tolerance
 122 //           of theV
 123 //=======================================================================
 124
 125 static Standard_Real preciseU (const BRepAdaptor_Surface&  theSurf,
 126                                const TopoDS_Edge&          theE,
 127                                const TopoDS_Vertex&        theV,
 128                                const Handle(Geom2d_Curve)& theC,
 129                                const Standard_Boolean      theFirstEnd)
 130 {
 131   Standard_Boolean isForward = ( theE.Orientation () == TopAbs_FORWARD );
 132   if (theFirstEnd) isForward = !isForward;
 133
 134   // find the first point in 2d and 3d
 135   Standard_Real f,l;
 136   BRep_Tool::Range( theE, f, l );
 137   Standard_Real u0 = isForward ? l : f;
 138   gp_Pnt2d aP2d0 = theC->Value( u0 );
 139   gp_Pnt aPnt0 = theSurf.Value( aP2d0.X(), aP2d0.Y() );
 140
 141   // shift in 2d and 3d
 142   Standard_Real du = ( l - f ) / 100, du3d = 0;
 143   if (isForward)
 144     du = -du;
 145
 146   // target parameter
 147   Standard_Real u;
 148
 149   while (du3d < ::RealSmall())
 150   {
 151     // u for test
 152     u = u0 + du;
 153     du *= 10; // for the next iteration: increase du untill du3d is large enough
 154
 155     // find out how u is far from u0 in 3D
 156     gp_Pnt2d aP2d  = theC->Value( u );
 157     gp_Pnt aPnt  = theSurf.Value( aP2d.X(), aP2d.Y() );
 158     du3d = aPnt0.Distance( aPnt );
 159   }
 160
 161   // find u such that the 3D point is just out of tolerance of theV
 162   Standard_Real tolV = BRep_Tool::Tolerance( theV ) + Precision::Confusion();
 163   u = u0 + du * tolV / du3d;
 164
 165   // check that u is within the range
 166   if ( isForward ? (u < f) : (u > l) )
 167     u = u0 + du;
 168
 169   return u;
 170 }
 171
 172 //=======================================================================
 173 //function : SelectEdge
 174 //purpose  : Find in the list <LE> the edge <NE> connected with <CE> by
 175 //           the vertex <CV>.
 176 //           <NE> is removed from the list. If <CE> is in <LE>
 177 //           with the same orientation, it's removed from the list
 178 //=======================================================================
 179
 180 static Standard_Boolean  SelectEdge(const BRepAdaptor_Surface& Surf,
 181                                     const TopoDS_Edge&    CE,
 182                                     const TopoDS_Vertex&  CV,
 183                                     TopoDS_Edge&          NE,
 184                                     const TopTools_ListOfShape& LE)
 185 {
 186   NE.Nullify();
 187
 188   if (LE.Extent() > 1) {
 189     //--------------------------------------------------------------
 190     // Several possible edges.
 191     // - Test the edges differents of CE
 192     //--------------------------------------------------------------
 193     TopoDS_Face FForward = Surf.Face();
 194     TopoDS_Edge aPrevNE;
 195
 196     gp_Vec2d CTg1, Tg1, CTg2, Tg2;
 197     gp_Pnt2d PC, P;
 198
 199     Standard_Real f, l;
 200     Handle(Geom2d_Curve) Cc, C;
 201     Cc = BRep_Tool::CurveOnSurface(CE,FForward,f,l);
 202
 203     Standard_Boolean isForward = ( CE.Orientation () == TopAbs_FORWARD );
 204     Standard_Real uc, u, du = Precision::PConfusion();
 205     uc = isForward ? ( l - du ) : ( f + du );
 206     Cc->D1(uc, PC, CTg1);
 207     if (!isForward) CTg1.Reverse();
 208
 209     Standard_Real anglemin = 3 * PI, tolAng = 1.e-8;
 210
 211     // select an edge whose first derivative is most left of CTg1
 212     // ie an angle between Tg1 and CTg1 is least
 213     TopTools_ListIteratorOfListOfShape itl;
 214     for ( itl.Initialize(LE); itl.More(); itl.Next()) {
 215       const TopoDS_Edge& E = TopoDS::Edge(itl.Value());
 216       if (E.IsSame(CE))
 217         continue;
 218       if (! CV.IsSame( TopExp::FirstVertex( E, Standard_True )))
 219         continue;
 220
 221       isForward = ( E.Orientation () == TopAbs_FORWARD );
 222
 223       // get E curve
 224       C = BRep_Tool::CurveOnSurface(E,FForward,f,l);
 225       // get the first derivative Tg1
 226       u = isForward ? ( f + du ) : ( l - du );
 227       C->D1(u, P, Tg1);
 228       if (!isForward) Tg1.Reverse();
 229
 230       // -PI < angle < PI
 231       Standard_Real angle = Tg1.Angle(CTg1);
 232
 233       if (PI - Abs(angle) <= tolAng)
 234       {
 235         // an angle is too close to PI; assure that an angle sign really
 236         // reflects an edge position: +PI - an edge is worst,
 237         // -PI - an edge is best.
 238         u = preciseU( Surf, CE, CV, Cc, Standard_False);
 239         gp_Vec2d CTg;
 240         Cc->D1(u, PC, CTg);
 241         if (CE.Orientation() == TopAbs_REVERSED) CTg.Reverse();
 242
 243         u = preciseU( Surf, E, CV, C, Standard_True);
 244         C->D1(u, P, Tg1);
 245         if (!isForward) Tg1.Reverse();
 246
 247         angle = Tg1.Angle(CTg);
 248       }
 249
 250       Standard_Boolean isClose = ( Abs( angle - anglemin ) <= tolAng );
 251       if (angle <= anglemin) {
 252         if (isClose)
 253           aPrevNE = NE;
 254         else
 255           aPrevNE.Nullify();
 256         anglemin = angle ;
 257         NE = E;
 258       }
 259       else
 260         if (isClose)
 261           aPrevNE = E;
 262
 263     }
 264     if (!aPrevNE.IsNull()) {
 265       // select one of close edges, the most left one.
 266       Cc = BRep_Tool::CurveOnSurface( NE, FForward, f, l );
 267       uc = preciseU( Surf, NE, CV, Cc, Standard_True);
 268       Cc->D1(uc, PC, CTg1);
 269       if (NE.Orientation() != TopAbs_FORWARD) CTg1.Reverse();
 270
 271       u = preciseU( Surf, aPrevNE, CV, C, Standard_True);
 272       C->D1(u, P, Tg1);
 273       if (aPrevNE.Orientation() != TopAbs_FORWARD) Tg1.Reverse();
 274
 275       if ( Tg1.Angle(CTg1) < 0)
 276         NE = aPrevNE;
 277     }
 278   }
 279   else if (LE.Extent() == 1) {
 280     NE = TopoDS::Edge(LE.First());
 281   }
 282   else {
 283     return Standard_False;
 284   }
 285   return !NE.IsNull();
 286 }
 287
 288 //=======================================================================
 289 //function : SamePnt2d
 290 //purpose  :
 291 //=======================================================================
 292
 293 static Standard_Boolean  SamePnt2d(const TopoDS_Vertex& V1,
 294                                    const TopoDS_Edge&   E1,
 295                                    const TopoDS_Vertex& V2,
 296                                    const TopoDS_Edge&   E2,
 297                                    const TopoDS_Face&   F)
 298 {
 299   Standard_Real   f1,f2,l1,l2;
 300   Handle(Geom2d_Curve) C1 = BRep_Tool::CurveOnSurface(E1,F,f1,l1);
 301   Handle(Geom2d_Curve) C2 = BRep_Tool::CurveOnSurface(E2,F,f2,l2);
 302
 303   gp_Pnt2d P1 = C1->Value( BRep_Tool::Parameter(V1,E1));
 304   gp_Pnt2d P2 = C2->Value( BRep_Tool::Parameter(V2,E2));
 305
 306   Standard_Real Tol  = 100 * BRep_Tool::Tolerance(V1);
 307   Standard_Real Dist = P1.Distance(P2);
 308   return Dist < Tol;
 309 }
 310
 311
 312 //=======================================================================
 313 //function : StoreInMVE
 314 //purpose  :
 315 //=======================================================================
 316
 317 static void StoreInMVE (const TopoDS_Face& /*F*/,
 318                         TopoDS_Edge& E,
 319                         TopTools_DataMapOfShapeListOfShape& MVE )
 320
 321 {
 322   TopoDS_Vertex V1, V2;
 323   TopTools_ListOfShape Empty;
 324
 325   TopExp::Vertices(E,V1,V2);
 326   if (!MVE.IsBound(V1)) {
 327     MVE.Bind(V1,Empty);
 328   }
 329   MVE(V1).Append(E);
 330
 331   if (!MVE.IsBound(V2)) {
 332     MVE.Bind(V2,Empty);
 333   }
 334   MVE(V2).Append(E);
 335 }
 336
 337 //=======================================================================
 338 //function : RemoveFromMVE
 339 //purpose  :
 340 //=======================================================================
 341
 342 static void RemoveFromMVE(const TopoDS_Edge& E,
 343                           TopTools_DataMapOfShapeListOfShape& MVE)
 344 {
 345   TopTools_ListIteratorOfListOfShape itl;
 346   TopoDS_Vertex  V1,V2;
 347   TopExp::Vertices (E,V1,V2);
 348   if (MVE.IsBound(V1))
 349     for ( itl.Initialize(MVE(V1)); itl.More(); itl.Next()) {
 350       if (itl.Value().IsEqual(E)) {
 351         MVE(V1).Remove(itl);
 352         break;
 353       }
 354     }
 355   if (MVE.IsBound(V2))
 356     for ( itl.Initialize(MVE(V2)); itl.More(); itl.Next()) {
 357       if (itl.Value().IsEqual(E)) {
 358         MVE(V2).Remove(itl);
 359         break;
 360       }
 361     }
 362 }
 363 //=======================================================================
 364 //function : addConnected
 365 //purpose  : add to <EM> all edges reachable from <E>
 366 //=======================================================================
 367
 368 static void addConnected(const TopoDS_Shape& E,
 369                          TopTools_MapOfShape& EM,
 370                          TopTools_MapOfShape& VM,
 371                          const TopTools_DataMapOfShapeListOfShape& MVE)
 372 {
 373   // Loop on vertices of E
 374   TopoDS_Iterator itV ( E );
 375   for ( ; itV.More(); itV.Next()) {
 376
 377     if ( ! VM.Add ( itV.Value() )) continue;
 378
 379     // Loop on edges sharing V
 380     TopTools_ListIteratorOfListOfShape itE( MVE( itV.Value() ) );
 381     for (; itE.More(); itE.Next()) {
 382       if ( EM.Add( itE.Value() ))
 383         addConnected ( itE.Value(), EM, VM, MVE );
 384     }
 385   }
 386 }
 387 //=======================================================================
 388 //function : canPassToOld
 389 //purpose  :
 390 //=======================================================================
 391
 392 // static Standard_Boolean canPassToOld (const TopoDS_Shape& V,
 393 //                                    TopTools_MapOfShape& UsedShapesMap,
 394 //                                    const TopTools_DataMapOfShapeListOfShape& MVE,
 395 //                                    const TopTools_MapOfShape& SectionEdgesMap)
 396 // {
 397 //   TopTools_ListIteratorOfListOfShape itE( MVE(V) );
 398 //   // Loop on edges sharing V
 399 //   for (; itE.More(); itE.Next()) {
 400 //     if ( !UsedShapesMap.Add( itE.Value() ))
 401 //       continue; // already checked
 402
 403 //     if ( !SectionEdgesMap.Contains( itE.Value() ))
 404 //       return Standard_True; // WE PASSED
 405
 406 //     TopoDS_Iterator itV( itE.Value() );
 407 //     // Loop on vertices of an edge
 408 //     for (; itV.More(); itV.Next()) {
 409 //       if ( !UsedShapesMap.Add( itV.Value() ))
 410 //      continue; // already checked
 411 //       else
 412 //      return canPassToOld( itV.Value(), UsedShapesMap, MVE, SectionEdgesMap);
 413 //     }
 414 //   }
 415 //   return Standard_False;
 416 // }
 417
 418 //=======================================================================
 419 //function : MakeDegenAndSelect
 420 //purpose  : Find parameter of intersection of <CE> with <DE> and
 421 //           select an edge with its parameter closest to found one.
 422 //           Return new degenerated edge trimming <DE> by found parameters
 423 //=======================================================================
 424
 425 static TopoDS_Edge MakeDegenAndSelect(const TopoDS_Edge& CE,
 426                                       const TopoDS_Vertex& CV,
 427                                       TopoDS_Edge& NE,
 428                                       TopTools_SequenceOfShape& EdgesSeq,
 429                                       TColStd_SequenceOfReal& USeq,
 430                                       const TopoDS_Edge& DE)
 431 {
 432   if (EdgesSeq.Length() < 3) {
 433     if (CE == EdgesSeq.First())
 434       NE = TopoDS::Edge( EdgesSeq.Last() );
 435     else
 436       NE = TopoDS::Edge( EdgesSeq.First() );
 437     return DE;
 438   }
 439
 440   // find parameter on DE where it intersects CE
 441
 442   Standard_Real U1;
 443   Standard_Integer i, nb = EdgesSeq.Length();
 444   for (i=1; i<= nb; ++i) {
 445     if (CE == EdgesSeq(i)) {
 446       U1 = USeq(i);
 447       break;
 448     }
 449   }
 450
 451   // select NE with param closest to U1 thus finding U2 for a new degen edge
 452
 453   Standard_Real U2, dU, dUmin = 1.e100;
 454   Standard_Boolean isReversed = ( DE.Orientation() == TopAbs_REVERSED );
 455   for (i=1; i<= nb; ++i) {
 456     dU = USeq(i) - U1;
 457     if (isReversed ? (dU > 0) : (dU < 0))
 458         continue;
 459     dU = Abs( dU );
 460     if ( dU  > dUmin || IsEqual( dU, 0.))
 461       continue;
 462     const TopoDS_Edge& E = TopoDS::Edge ( EdgesSeq(i) );
 463     if ( ! CV.IsSame( TopExp::FirstVertex( E , Standard_True )))
 464       continue;
 465     NE = E;
 466     dUmin = dU + Epsilon(dU);
 467     U2 = USeq(i);
 468   }
 469
 470   // make a new degenerated edge
 471   TopoDS_Edge NewDegen = TopoDS::Edge ( DE.EmptyCopied() );
 472
 473   Standard_Real Tol = BRep_Tool::Tolerance( CV );
 474   TopoDS_Vertex V = CV;
 475
 476   BRep_Builder B;
 477   V.Orientation( NewDegen.Orientation() );
 478   B.UpdateVertex( V, U1, NewDegen, Tol);
 479   B.Add ( NewDegen , V );
 480
 481   V.Reverse();
 482   B.UpdateVertex( V, U2, NewDegen, Tol);
 483   B.Add ( NewDegen , V );
 484
 485   return NewDegen;
 486 }
 487
 488 //=======================================================================
 489 //function : prepareDegen
 490 //purpose  : Intersect <DegEdge> with edges bound to its vertex in <MVE>
 491 //           and store intersection parameter on <DegEdge> in
 492 //           <USeq> as well as the edges them-self in <EdgesSeq>.
 493 //           Bind <DegEdgeIndex> to vertex of <DegEdge> in <MVDEI>
 494 //=======================================================================
 495
 496 static void prepareDegen (const TopoDS_Edge&                        DegEdge,
 497                           const TopoDS_Face&                        F,
 498                           const TopTools_DataMapOfShapeListOfShape& MVE,
 499                           TopTools_SequenceOfShape&                 EdgesSeq,
 500                           TColStd_SequenceOfReal&                   USeq,
 501                           TopTools_DataMapOfShapeInteger&           MVDEI,
 502                           const Standard_Integer                    DegEdgeIndex)
 503 {
 504   const TopoDS_Vertex& V = TopExp::FirstVertex ( DegEdge );
 505   MVDEI.Bind ( V, DegEdgeIndex );
 506
 507   const TopTools_ListOfShape& EdgesList = MVE ( V );
 508   // if only 2 edges come to degenerated one, no pb in selection and
 509   // no need to intersect them, just simulate asked data
 510   Standard_Boolean doIntersect =  ( EdgesList.Extent() > 2 );
 511
 512   BRepAdaptor_Curve2d DC, C;
 513   Geom2dInt_GInter InterCC;
 514   Standard_Real Tol = Precision::PConfusion();
 515   if ( doIntersect )
 516     DC.Initialize( DegEdge, F );
 517
 518   // avoid intersecting twice the same edge
 519   BRepOffset_DataMapOfShapeReal EUMap ( EdgesList.Extent() );
 520
 521   Standard_Real U, f, l;
 522   BRep_Tool::Range (DegEdge, f, l);
 523
 524   TopTools_ListIteratorOfListOfShape itE (EdgesList);
 525   for (; itE.More(); itE.Next()) {
 526
 527     const TopoDS_Edge& E = TopoDS::Edge ( itE.Value() );
 528
 529     if ( !doIntersect) {
 530       U = 0.; // it won't be used
 531     }
 532     else if ( BRep_Tool::IsClosed( E, F )) {
 533       // seam edge: select U among f and l
 534       Standard_Boolean first = Standard_True;
 535       if ( V.IsSame ( TopExp::FirstVertex( E, Standard_True ) ))
 536         first = Standard_False;
 537       if ( DegEdge.Orientation() == TopAbs_REVERSED )
 538         first = !first;
 539       U = first ? f : l;
 540     }
 541     else if ( EUMap.IsBound( E ) ) {
 542       // same edge already bound
 543       U = EUMap( E );
 544     }
 545     else {
 546       // intersect 2d curves
 547       C.Initialize( E, F );
 548       InterCC.Perform ( DC, C , Tol, Tol );
 549       if (! InterCC.IsDone() || InterCC.NbPoints() == 0) {
 550         MESSAGE ( "NO 2d INTERSECTION ON DEGENERATED EDGE" );
 551         continue;
 552       }
 553       // hope there is only one point of intersection
 554       U = InterCC.Point( 1 ).ParamOnFirst();
 555     }
 556     USeq.Append ( U );
 557     EdgesSeq.Append ( E );
 558   }
 559 }
 560 //=======================================================================
 561 //function : Perform
 562 //purpose  : Make loops.
 563 //=======================================================================
 564
 565 void Partition_Loop2d::Perform()
 566 {
 567
 568   Standard_Integer NbConstEdges = myConstEdges.Extent();
 569   TopTools_DataMapOfShapeListOfShape MVE(NbConstEdges) , MVE2(NbConstEdges);
 570   TopTools_DataMapIteratorOfDataMapOfShapeListOfShape Mapit;
 571   TopTools_ListIteratorOfListOfShape itl;
 572   TopoDS_Vertex V1,V2;
 573   BRepAdaptor_Surface Surface ( myFace, Standard_False );
 574
 575   // degenerated edges and parameters of their 2d intersection with other edges
 576   TopoDS_Edge                    DE [2];
 577   TopTools_SequenceOfShape       SEID [2]; // seq of edges intersecting degenerated
 578   TColStd_SequenceOfReal         SeqU [2]; // n-th U corresponds to n-th edge in SEID
 579   TopTools_DataMapOfShapeInteger MVDEI(2); // map vertex - degenerated edge index
 580   Standard_Integer               iDeg = 0; // index of degenerated edge [0,1]
 581
 582   //---------------------------------------------------------
 583   // Construction map vertex => edges, find degenerated edges
 584   //---------------------------------------------------------
 585   for (itl.Initialize(myConstEdges); itl.More(); itl.Next()) {
 586     TopoDS_Edge& E = TopoDS::Edge(itl.Value());
 587     if ( BRep_Tool::Degenerated( E )) {
 588       if (DE[0].IsNull()) DE[0] = E;
 589       else                DE[1] = E;
 590     }
 591     else
 592       StoreInMVE(myFace,E,MVE);
 593   }
 594
 595   // fill data for degenerated edges
 596   if ( ! DE[0].IsNull() )
 597     prepareDegen ( DE[0], myFace, MVE, SEID[0], SeqU[0], MVDEI, 0);
 598   if ( ! DE[1].IsNull() )
 599     prepareDegen ( DE[1], myFace, MVE, SEID[1], SeqU[1], MVDEI, 1);
 600
 601
 602   // to detect internal wires
 603   Standard_Boolean isInternCW = 0;
 604   MVE2 = MVE;
 605
 606
 607   //------------------------------
 608   // Construction of all the wires
 609   //------------------------------
 610   // first, we collect wire edges in WEL list looking for same edges that
 611   // will be then removed possibly exploding a wire into parts;
 612   // second, build wire(s)
 613
 614   while (!MVE.IsEmpty()) {
 615
 616     TopoDS_Vertex    VF,CV;
 617     TopoDS_Edge      CE,NE,EF;
 618     TopoDS_Wire      NW;
 619     BRep_Builder     B;
 620     Standard_Boolean End = Standard_False;
 621     TopTools_ListOfShape WEL;
 622
 623     Mapit.Initialize(MVE);
 624     if (Mapit.Value().IsEmpty()) {
 625       MVE.UnBind(Mapit.Key());
 626       continue;
 627     }
 628
 629     // EF first edge.
 630     EF = CE = TopoDS::Edge(Mapit.Value().First());
 631     // VF first vertex
 632     VF = TopExp::FirstVertex( CE, Standard_True);
 633
 634     isInternCW = Standard_True;
 635
 636     TopTools_MapOfShape addedEM  (NbConstEdges); // map of edges added to WEL
 637     TopTools_MapOfShape doubleEM (NbConstEdges); // edges encountered twice in WEL
 638
 639     //-------------------------------
 640     // Construction of a wire.
 641     //-------------------------------
 642     while (!End) {
 643
 644       // only a seam is allowed twice in a wire, the others should be removed
 645       if (addedEM.Add ( CE ) || BRep_Tool::IsClosed( CE, myFace ) )
 646         WEL.Append( CE );
 647       else {
 648         doubleEM.Add( CE );
 649         RemoveFromMVE (CE,MVE2);
 650         TopoDS_Edge CERev = CE;
 651         CERev.Reverse();
 652         RemoveFromMVE (CERev,MVE2);
 653       }
 654
 655       RemoveFromMVE (CE,MVE);
 656
 657       CV = TopExp::LastVertex( CE, Standard_True);
 658
 659       if (isInternCW && !mySectionEdges.Contains(CE))
 660         // wire is internal if all edges are section ones
 661         isInternCW = Standard_False;
 662
 663       if (MVDEI.IsBound( CV )) { // CE comes to the degeneration
 664         iDeg = MVDEI( CV );
 665         TopoDS_Edge NewDegen;
 666         NewDegen = MakeDegenAndSelect( CE, CV, NE, SEID[iDeg], SeqU[iDeg], DE[iDeg]);
 667         WEL.Append( NewDegen );
 668         CE = NE;
 669         End = CV.IsSame( VF );
 670         continue;
 671       }
 672
 673       //--------------
 674       // stop test
 675       //--------------
 676       if (MVE(CV).IsEmpty()) {
 677         End=Standard_True;
 678         MVE.UnBind(CV);
 679       }
 680       else if (CV.IsSame(VF) && SamePnt2d(CV,CE, VF,EF, myFace) ) {
 681         End = Standard_True;
 682       }
 683       else {
 684         //----------------------------
 685         // select new current edge
 686         //----------------------------
 687         if (! SelectEdge (Surface,CE,CV,NE,MVE(CV))) {
 688           MESSAGE ( " NOT CLOSED WIRE " );
 689           End=Standard_True;
 690         }
 691         else
 692           CE = NE;
 693       }
 694     } // while ( !End )
 695
 696
 697     // WEL is built, built wire(s)
 698
 699
 700     itl.Initialize( WEL );
 701     if ( doubleEM.IsEmpty()) { // no double edges
 702       B.MakeWire( NW );
 703       for (; itl.More(); itl.Next())
 704         B.Add ( NW, itl.Value());
 705       if (isInternCW) myInternalWL.Append(NW);
 706       else            myNewWires.Append  (NW);
 707     }
 708
 709     else {
 710       // remove double and degenerated edges from WEL
 711       while (itl.More()) {
 712         const TopoDS_Edge& E = TopoDS::Edge ( itl.Value() );
 713         if ( doubleEM.Contains( E ) || BRep_Tool::Degenerated( E ))
 714           WEL.Remove( itl );
 715         else
 716            itl.Next();
 717       }
 718       if ( WEL.IsEmpty())
 719         continue;
 720       // remove double edges from SEID and SeqU
 721       Standard_Integer i,j;
 722       for (j=0; j<2; ++j) {
 723         for (i=1; i<=SEID[j].Length(); ++i) {
 724           if (doubleEM.Contains( SEID[j].Value(i))) {
 725             SEID[j].Remove( i );
 726             SeqU[j].Remove( i-- );
 727           }
 728         }
 729       }
 730       // removal of doulbe edges can explode a wire into parts,
 731       // make new wires of them.
 732       // A Loop like previous one but without 2d check
 733       while ( !WEL.IsEmpty() ) {
 734         CE = TopoDS::Edge( WEL.First() );
 735         WEL.RemoveFirst();
 736         B.MakeWire( NW );
 737         VF = TopExp::FirstVertex ( CE, Standard_True);
 738
 739         End = Standard_False;
 740         while ( !End) {
 741           B.Add( NW, CE );
 742           CV = TopExp::LastVertex  ( CE, Standard_True);
 743
 744           if (MVDEI.IsBound( CV )) {   // CE comes to the degeneration
 745             iDeg = MVDEI( CV );
 746             TopoDS_Edge NewDegen;
 747             NewDegen = MakeDegenAndSelect( CE, CV, NE, SEID[iDeg], SeqU[iDeg], DE[iDeg]);
 748             B.Add( NW, NewDegen );
 749             End = CV.IsSame( VF );
 750             CE = NE;
 751             if (!NE.IsNull()) { // remove NE from WEL
 752               for (itl.Initialize( WEL ); itl.More(); itl.Next())
 753                 if ( NE == itl.Value()) {
 754                   WEL.Remove( itl );
 755                   break;
 756                 }
 757             }
 758           }  // end degeneration
 759
 760           else {
 761             if (CV.IsSame( VF )) {
 762               End = Standard_True;
 763               continue;
 764             }
 765             // edges in WEL most often are well ordered
 766             // so try to iterate until the End
 767             Standard_Boolean add = Standard_False;
 768             itl.Initialize(WEL);
 769             while ( itl.More() && !End) {
 770               NE = TopoDS::Edge( itl.Value() );
 771               if ( CV.IsSame( TopExp::FirstVertex( NE, Standard_True ))) {
 772                 WEL.Remove( itl );
 773                 if (add)
 774                   B.Add( NW, CE );
 775                 CE = NE;
 776                 add = Standard_True;
 777                 CV = TopExp::LastVertex( CE, Standard_True);
 778                 if (MVDEI.IsBound( CV ) || CV.IsSame( VF ))
 779                   break;
 780               }
 781               else
 782                 itl.Next();
 783             }
 784             if (!add)
 785               End = Standard_True;
 786           }
 787         } // !End
 788
 789         myInternalWL.Append( NW );
 790       }
 791     } // end building new wire(s) from WEL
 792
 793   } // end Loop on MVE
 794
 795   // all wires are built
 796
 797
 798   // ============================================================
 799   // select really internal wires i.e. those from which we can`t
 800   // pass to an old (not section) edge
 801   // ============================================================
 802
 803   Standard_Integer nbIW = myInternalWL.Extent();
 804   if (nbIW == 0)
 805     return;
 806
 807   if ( myNewWires.Extent() != 1 && nbIW > 1) {
 808     TopTools_MapOfShape outerEM (NbConstEdges); // edges connected to non-section ones
 809     TopTools_MapOfShape visitedVM (NbConstEdges);
 810     for ( itl.Initialize( myConstEdges ); itl.More(); itl.Next()) {
 811       if ( ! mySectionEdges.Contains( itl.Value() ))
 812         addConnected (itl.Value(), outerEM, visitedVM, MVE2);
 813     }
 814     // if an edge of a wire is in <outerEM>, the wire is not internal
 815     TopExp_Explorer expIWE;
 816     TopTools_ListIteratorOfListOfShape itIW ( myInternalWL );
 817     while (itIW.More()) {
 818       expIWE.Init ( itIW.Value() , TopAbs_EDGE );
 819       if ( outerEM.Contains( expIWE.Current() )) {
 820         myNewWires.Append ( itIW.Value() );
 821         myInternalWL.Remove( itIW ); // == itIW.Next()
 822       }
 823       else
 824         itIW.Next();
 825     }
 826   }
 827 }
 828 //=======================================================================
 829 //function : isHole
 830 //purpose  :
 831 //=======================================================================
 832
 833 static Standard_Boolean isHole (const TopoDS_Wire& W,
 834                                 const TopoDS_Face& F)
 835 {
 836   BRep_Builder B;
 837   TopoDS_Shape newFace = F.EmptyCopied();
 838   B.Add(newFace,W.Oriented(TopAbs_FORWARD));
 839   BRepTopAdaptor_FClass2d classif (TopoDS::Face(newFace),
 840                                    Precision::PConfusion());
 841   return (classif.PerformInfinitePoint() == TopAbs_IN);
 842 }
 843
 844 //=======================================================================
 845 //function : IsInside
 846 //purpose  : check if W1 is inside W2. Suppose W2 is not a hole !!!!
 847 //=======================================================================
 848
 849 static Standard_Boolean isInside(const TopoDS_Face& F,
 850                                  const TopoDS_Wire& W1,
 851                                  const TopoDS_Wire& W2)
 852 {
 853   // make a face with wire W2
 854   BRep_Builder B;
 855   TopoDS_Shape aLocalShape = F.EmptyCopied();
 856   TopoDS_Face newFace = TopoDS::Face(aLocalShape);
 857   B.Add(newFace,W2);
 858
 859   // get any 2d point of W1
 860   TopExp_Explorer exp(W1,TopAbs_EDGE);
 861   if (BRep_Tool::Degenerated( TopoDS::Edge( exp.Current() )))
 862     exp.Next();
 863   const TopoDS_Edge& e = TopoDS::Edge(exp.Current());
 864   Standard_Real f,l;
 865   Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface(e,F,f,l);
 866   gp_Pnt2d pt2d(C2d->Value( 0.5 * ( f + l )));
 867
 868   BRepTopAdaptor_FClass2d classif(newFace,Precision::PConfusion());
 869   return (classif.Perform(pt2d) == TopAbs_IN);
 870 }
 871
 872 //=======================================================================
 873 //function : NewWires
 874 //purpose  : Returns the list of wires performed.
 875 //           can be an empty list.
 876 //=======================================================================
 877
 878 const TopTools_ListOfShape&  Partition_Loop2d::NewWires() const
 879 {
 880   return myNewWires;
 881 }
 882
 883 //=======================================================================
 884 //function : NewFaces
 885 //purpose  : Returns the list of faces.
 886 //Warning  : The method <WiresToFaces> as to be called before.
 887 //           can be an empty list.
 888 //=======================================================================
 889
 890 const TopTools_ListOfShape&  Partition_Loop2d::NewFaces() const
 891 {
 892   return myNewFaces;
 893 }
 894
 895 //=======================================================================
 896 //function : findEqual
 897 //purpose  : move wires form <WL> to <EqWL> pairs of wires build of the same edges
 898 //=======================================================================
 899
 900 static void findEqual (TopTools_ListOfShape& WL,
 901                        TopTools_DataMapOfShapeShape& EqWM,
 902                        const TopoDS_Face& F)
 903 {
 904   TopTools_ListIteratorOfListOfShape it1, it2;
 905   Standard_Integer i,j;
 906   TColStd_MapOfInteger IndMap;
 907   for (it1.Initialize(WL), i=1;  it1.More();  it1.Next(), i++) {
 908
 909     if (IndMap.Contains(i)) continue;
 910     const TopoDS_Wire& Wire1 = TopoDS::Wire( it1.Value());
 911
 912     for (it2.Initialize(WL), j=1;  it2.More();  it2.Next(), j++) {
 913
 914       if (j <= i || IndMap.Contains(j)) continue;
 915
 916       TopTools_IndexedMapOfShape EdgesMap;
 917       TopExp::MapShapes (Wire1, TopAbs_EDGE, EdgesMap);
 918
 919       const TopoDS_Shape& Wire2 = it2.Value();
 920       TopoDS_Iterator itE ( Wire2);
 921       for (; itE.More(); itE.Next()) {
 922         if ( !EdgesMap.Contains( itE.Value()) )
 923           break;
 924       }
 925       if (!itE.More()) { // all edges are same
 926         if (isHole( Wire1, F)) {
 927           EqWM.Bind ( Wire1, Wire2 );
 928         }
 929         else {
 930           EqWM.Bind ( Wire2, Wire1 );
 931         }
 932         IndMap.Add(i);
 933         IndMap.Add(j);
 934         break;
 935       }
 936     }
 937   }
 938   // clear WL
 939   it1.Initialize(WL);
 940   i=1;
 941   while (it1.More()) {
 942     if (IndMap.Contains(i))
 943       WL.Remove(it1); // next node becomes current and with Next() we would miss it
 944     else
 945       it1.Next();
 946     i++;
 947   }
 948 }
 949
 950 //=======================================================================
 951 //function : classify
 952 //purpose  : bind to a wire a list of internal wires
 953 //=======================================================================
 954
 955 static void classify(const TopTools_DataMapOfShapeShape& EqWM,
 956                      BRepAlgo_AsDes& OuterInner,
 957                      const TopoDS_Face& F)
 958 {
 959   TopTools_DataMapIteratorOfDataMapOfShapeShape it1, it2;
 960
 961   for (it1.Initialize(EqWM);  it1.More();  it1.Next()) {
 962     // find next after it1.Value()
 963     for (it2.Initialize(EqWM);  it2.More();  it2.Next())
 964       if (it1.Value().IsSame( it2.Value() ))
 965       {
 966         it2.Next();
 967         break;
 968       }
 969     for ( ;  it2.More();  it2.Next()) {
 970       const TopoDS_Wire& Wire1 = TopoDS::Wire( it1.Value() );
 971       const TopoDS_Wire& Wire2 = TopoDS::Wire( it2.Value() );
 972       if (isInside(F, Wire1, Wire2))
 973         OuterInner.Add (Wire2, Wire1);
 974       else if (isInside(F, Wire2, Wire1))
 975         OuterInner.Add (Wire1, Wire2);
 976     }
 977   }
 978 }
 979 //=======================================================================
 980 //function : WiresToFaces
 981 //purpose  : Build faces from the wires result.
 982 //           <EdgeImage> serves to  find  original edge by new
 983 //           one. <Section> contains edges resulting from face
 984 //           intersections
 985 //=======================================================================
 986
 987 void  Partition_Loop2d::WiresToFaces(const BRepAlgo_Image& )
 988 {
 989   Standard_Integer nbW = myNewWires.Extent() + myInternalWL.Extent();
 990   if (nbW==0)
 991     return;
 992
 993   BRepAlgo_FaceRestrictor FR;
 994   FR.Init (myFace,Standard_False);
 995
 996   // FaceRestrictor is instable in rather simple cases
 997   // (ex. a single face of bellecoque.brep splited by 10 planes:
 998   // sometimes 1-2 faces are missing ).
 999   // So we use it as less as possible: no holes -> make faces by hands
1000
1001
1002   // are there holes in myFace ?
1003   Standard_Boolean hasOldHoles = Standard_False;
1004   TopoDS_Iterator itOldW (myFace);
1005   if ( itOldW.More()) {
1006     const TopoDS_Wire& FirstOldWire = TopoDS::Wire( itOldW.Value() );
1007     itOldW.Next();
1008     hasOldHoles = itOldW.More() || isHole( FirstOldWire, myFace);
1009   }
1010   if (myInternalWL.IsEmpty() && !hasOldHoles) {
1011     // each wire bounds one face
1012     BRep_Builder B;
1013     TopTools_ListIteratorOfListOfShape itNW (myNewWires);
1014     for (; itNW.More(); itNW.Next()) {
1015       TopoDS_Face NF = TopoDS::Face ( myFace.EmptyCopied() );
1016       B.Add ( NF, itNW.Value() );
1017       NF.Orientation( myFaceOri);
1018       myNewFaces.Append ( NF );
1019     }
1020     return;
1021   }
1022
1023   // FaceRestrictor can't classify wires build on all the same edges
1024   // and gives incorrect result in such cases (ex. a plane cut into 2 parts by cylinder)
1025   // We must make faces of equal wires separately. One of equal wires makes a
1026   // hole in a face and should come together with outer wires of face.
1027   // The other of a wires pair bounds a face that may have holes in turn.
1028
1029   // Find equal wires among internal wires
1030   TopTools_DataMapOfShapeShape EqWM; // key is a hole part of a pair of equal wires
1031   findEqual (myInternalWL, EqWM, myFace);
1032
1033   if (!EqWM.IsEmpty()) { // there are equal wires
1034
1035     if (hasOldHoles)
1036       myInternalWL.Append( myNewWires ); // an old wire can be inside an equal wire
1037
1038     // classify equal wire pairs
1039     BRepAlgo_AsDes OuterInner;
1040     classify (EqWM,OuterInner,myFace);
1041
1042     // make face of most internal of equal wires and its inner wires
1043     while ( !EqWM.IsEmpty()) {
1044
1045       TopTools_ListOfShape prevHolesL; // list of hole-part of previous most internal equal wires
1046
1047       // find most internal wires among pairs (key - hole, value - outer part)
1048       TopTools_DataMapIteratorOfDataMapOfShapeShape it(EqWM);
1049       Standard_Integer nbEqW = EqWM.Extent(); // protection against infinite loop
1050       for ( ; it.More(); it.Next()) {
1051
1052         TopoDS_Wire outerW = TopoDS::Wire ( it.Value() );
1053         if (  OuterInner.HasDescendant( outerW ) && // has internal
1054              ! OuterInner.Descendant( outerW ).IsEmpty() )
1055           continue;
1056
1057         FR.Add( outerW );
1058
1059         // add internal wires that are inside of outerW
1060         TopTools_ListIteratorOfListOfShape itIW (myInternalWL);
1061         while ( itIW.More()) {
1062           TopoDS_Wire IW = TopoDS::Wire ( itIW.Value() );
1063           if ( isInside (myFace, IW, outerW)) {
1064             FR.Add (IW);
1065             myInternalWL.Remove( itIW ); // == itIW.Next() !!!
1066           }
1067           else
1068             itIW.Next();
1069         }
1070
1071         // the hole-part of current pair of equal wires will be in the next new face
1072         prevHolesL.Append ( it.Key() );
1073
1074       } // Loop on map of equal pairs searching for innermost wires
1075
1076       // make faces
1077       FR.Perform();
1078       if (FR.IsDone()) {
1079         for (; FR.More(); FR.Next())
1080           myNewFaces.Append(FR.Current());
1081       }
1082
1083       FR.Clear();
1084
1085       // add hole-parts to FaceRestrictor,
1086       // remove them from the EqWM,
1087       // remove found wires as internal of resting classified wires
1088       Standard_Boolean clearOuterInner =  ( prevHolesL.Extent() < EqWM.Extent() );
1089       TopTools_ListIteratorOfListOfShape itPrev (prevHolesL);
1090       for (; itPrev.More(); itPrev.Next()) {
1091         TopoDS_Wire& Hole = TopoDS::Wire ( itPrev.Value() );
1092         FR.Add ( Hole );
1093         if (clearOuterInner) {
1094           const TopoDS_Wire& outerW = TopoDS::Wire ( EqWM.Find( Hole ) );
1095           // Loop on wires including outerW
1096           TopTools_ListIteratorOfListOfShape itO( OuterInner.Ascendant( outerW ));
1097           for (; itO.More(); itO.Next()) {
1098             TopTools_ListOfShape& innerL = OuterInner.ChangeDescendant( itO.Value() );
1099             TopTools_ListIteratorOfListOfShape itI (innerL);
1100             // Loop on internal wires of current including wire
1101             for (; itI.More(); itI.Next())
1102               if ( outerW.IsSame( itI.Value() )) {
1103                 innerL.Remove( itI );   break;
1104               }
1105           }
1106         }
1107         EqWM.UnBind ( Hole );
1108       }
1109
1110       if (nbEqW == EqWM.Extent())
1111       {
1112         // error: pb with wires classification
1113 #ifdef DEB
1114         cout << "Partition_Loop2d::WiresToFaces(), pb with wires classification" << endl;
1115 #endif
1116         break;
1117       }
1118
1119     } // while (!EqWM.IsEmpty)
1120
1121   } //  if !EqWM.IsEmpty()
1122
1123   myNewWires.Append ( myInternalWL );
1124
1125   TopTools_ListIteratorOfListOfShape itW (myNewWires);
1126   for (; itW.More(); itW.Next()) {
1127     TopoDS_Wire& W = TopoDS::Wire ( itW.Value() );
1128     FR.Add(W);
1129   }
1130   FR.Perform();
1131   for (; FR.IsDone() && FR.More(); FR.Next())
1132     myNewFaces.Append(FR.Current());
1133
1134
1135   TopTools_ListIteratorOfListOfShape itNF (myNewFaces);
1136   for (; itNF.More(); itNF.Next())
1137     itNF.Value().Orientation( myFaceOri );
1138 }