src/PARTITION/Partition_Loop2d.cxx

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