src/PARTITION/Partition_Loop3d.cxx

   1 //  GEOM PARTITION : partition algorithm
   2 //
   3 //  Copyright (C) 2003  CEA/DEN, EDF R&D
   4 //
   5 //
   6 //
   7 //  File   : Partition_Loop3d.cxx
   8 //  Module : GEOM
   9
  10 using namespace std;
  11 #include "Partition_Loop3d.ixx"
  12
  13 #include <TopExp_Explorer.hxx>
  14 #include <TopExp.hxx>
  15 #include <BRep_Builder.hxx>
  16 #include <TopTools_MapOfShape.hxx>
  17 #include <TopTools_ListIteratorOfListOfShape.hxx>
  18 #include <TopoDS_Shell.hxx>
  19 #include <TopoDS_Iterator.hxx>
  20 #include <TopoDS.hxx>
  21 #include <TopTools_MapIteratorOfMapOfShape.hxx>
  22 #include <gp_Vec.hxx>
  23 #include <gp_Pnt.hxx>
  24 #include <Geom2d_Curve.hxx>
  25 #include <BRep_Tool.hxx>
  26 #include <Geom_Surface.hxx>
  27 #include <gp_Pnt2d.hxx>
  28 #include <gp_Vec2d.hxx>
  29 #include <gp_Dir2d.hxx>
  30 #include <Geom_Curve.hxx>
  31
  32 //=======================================================================
  33 //function : Partition_Loop3d
  34 //purpose  :
  35 //=======================================================================
  36
  37 Partition_Loop3d::Partition_Loop3d()
  38 {
  39 }
  40
  41 //=======================================================================
  42 //function : AddConstFaces
  43 //purpose  : Add faces of <S> as unique faces in the result.
  44 //=======================================================================
  45
  46 void Partition_Loop3d::AddConstFaces(const TopoDS_Shape& S)
  47 {
  48   TopExp_Explorer FaceExp(S, TopAbs_FACE);
  49   for (; FaceExp.More(); FaceExp.Next())
  50     myFaces.Append( FaceExp.Current() );
  51
  52   TopExp::MapShapesAndAncestors(S, TopAbs_EDGE, TopAbs_FACE, myEFMap);
  53 }
  54
  55 //=======================================================================
  56 //function : AddSectionFaces
  57 //purpose  : Add faces of <S> as double faces in the result.
  58 //=======================================================================
  59
  60 void Partition_Loop3d::AddSectionFaces(const TopoDS_Shape& S)
  61 {
  62   AddConstFaces( S );
  63   AddConstFaces( S.Reversed() );
  64 }
  65
  66 //=======================================================================
  67 //function : MakeShells
  68 //purpose  : Make and return shells.
  69 //           <AvoidFacesMap> can contain faces that must not be
  70 //           added to result shells.
  71 //=======================================================================
  72
  73 const TopTools_ListOfShape&
  74   Partition_Loop3d::MakeShells (const TopTools_MapOfOrientedShape& AvoidFacesMap)
  75 {
  76   myNewShells.Clear();
  77
  78   BRep_Builder Builder;
  79   TopTools_MapOfShape CheckedEdgesMap;
  80   TopTools_MapOfOrientedShape AddedFacesMap;
  81
  82   TopTools_ListIteratorOfListOfShape itF (myFaces);
  83   for (; itF.More(); itF.Next())
  84   {
  85     const TopoDS_Shape& FF = itF.Value();
  86     if (AvoidFacesMap.Contains( FF ) ||
  87         ! AddedFacesMap.Add( FF ) )
  88       continue;
  89
  90     // make a new shell
  91     TopoDS_Shell Shell;
  92     Builder.MakeShell(Shell);
  93     Builder.Add(Shell,FF);
  94
  95     // clear the maps from shapes added to previous Shell
  96     TopTools_MapIteratorOfMapOfShape itEM (CheckedEdgesMap);
  97     for (; itEM.More(); itEM.Next()) {
  98       TopTools_ListOfShape& FL = myEFMap.ChangeFromKey( itEM.Key());
  99       TopTools_ListIteratorOfListOfShape it (FL);
 100       while ( it.More()) {
 101         if (AddedFacesMap.Contains( it.Value()))
 102           FL.Remove( it );
 103         else
 104           it.Next();
 105       }
 106     }
 107     CheckedEdgesMap.Clear();
 108
 109
 110     // loop on faces added to Shell; add their neighbor faces to Shell and so on
 111     TopoDS_Iterator itAddedF (Shell);
 112     for (; itAddedF.More(); itAddedF.Next())
 113     {
 114       const TopoDS_Face& F = TopoDS::Face (itAddedF.Value());
 115
 116       // loop on edges of F; find a good neighbor face of F by E
 117       TopExp_Explorer EdgeExp(F, TopAbs_EDGE);
 118       for (; EdgeExp.More(); EdgeExp.Next())
 119       {
 120         const TopoDS_Edge& E = TopoDS::Edge( EdgeExp.Current());
 121         if (! CheckedEdgesMap.Add( E ))
 122           continue;
 123
 124         // candidate faces list
 125         const TopTools_ListOfShape& FL = myEFMap.ChangeFromKey(E);
 126         if (FL.IsEmpty())
 127           continue;
 128         // select one of neighbors
 129         TopoDS_Face SelF;
 130         if (FL.Extent() == 2) {
 131           if (! F.IsSame( FL.First() ))
 132             SelF = TopoDS::Face( FL.First() );
 133           else if (!F.IsSame( FL.Last() ))
 134             SelF = TopoDS::Face( FL.Last() );
 135         }
 136         else {
 137           // check if a face already added to Shell shares E
 138           TopTools_ListIteratorOfListOfShape it (FL);
 139           Standard_Boolean found = Standard_False;
 140           for (; !found && it.More(); it.Next())
 141             if (F != it.Value())
 142               found = AddedFacesMap.Contains( it.Value() );
 143           if (found)
 144             continue;
 145           // select basing on geometrical check
 146           Standard_Boolean GoodOri, inside;
 147           Standard_Real dot, MaxDot = -100;
 148           TopTools_ListOfShape TangFL; // tangent faces
 149           for ( it.Initialize( FL ) ; it.More(); it.Next()) {
 150             const TopoDS_Face& NeighborF = TopoDS::Face( it.Value());
 151             if (NeighborF.IsSame( F ))
 152               continue;
 153             inside = Partition_Loop3d::IsInside( E, F, NeighborF, 1, dot, GoodOri);
 154             if (!GoodOri)
 155               continue;
 156             if (!inside)
 157               dot = -dot - 3;
 158             if (dot < MaxDot)
 159               continue;
 160             if ( IsEqual( dot, MaxDot))
 161               TangFL.Append(SelF);
 162             else
 163               TangFL.Clear();
 164             MaxDot = dot;
 165             SelF = NeighborF;
 166           }
 167           if (!TangFL.IsEmpty()) {
 168             for (it.Initialize( TangFL ); it.More(); it.Next()) {
 169               const TopoDS_Face& NeighborF = TopoDS::Face( it.Value());
 170               if (Partition_Loop3d:: IsInside( E, SelF , NeighborF, 0, dot, GoodOri))
 171                 SelF = NeighborF;
 172             }
 173           }
 174         }
 175         if (!SelF.IsNull() &&
 176             AddedFacesMap.Add( SelF ) &&
 177             !AvoidFacesMap.Contains( SelF ))
 178           Builder.Add( Shell, SelF);
 179
 180       } // loop on edges of F
 181
 182     } // loop on the faces added to Shell
 183
 184     // Shell is complete
 185     myNewShells.Append( Shell );
 186
 187   } // loop on myFaces
 188
 189
 190   // prepare to the next call
 191   myFaces.Clear();
 192   myEFMap.Clear();
 193
 194   return myNewShells;
 195 }
 196
 197
 198
 199 //=======================================================================
 200 //function : Normal
 201 //purpose  :
 202 //=======================================================================
 203
 204 gp_Vec Partition_Loop3d::Normal(const TopoDS_Edge& E,
 205                                 const TopoDS_Face& F)
 206 {
 207   gp_Vec Norm, V1, V2;
 208   Standard_Real First, Last;
 209   gp_Pnt Ps;
 210
 211   Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface (E, F, First, Last);
 212   Handle(Geom_Surface) Sf = BRep_Tool::Surface(F);
 213
 214   gp_Pnt2d p = C2d->Value( 0.5*(First+Last) );
 215   Sf->D1(p.X(), p.Y(), Ps, V1, V2);
 216   Norm = V1.Crossed(V2);
 217
 218   if (F.Orientation() == TopAbs_REVERSED )
 219     Norm.Reverse();
 220
 221   return Norm;
 222 }
 223
 224 //=======================================================================
 225 //function : NextNormal
 226 //purpose  : find normal to F at point a little inside F near the middle of E
 227 //warning  : E must be properly oriented in F.
 228 //=======================================================================
 229
 230 static gp_Vec NextNormal(const TopoDS_Edge& E,
 231                          const TopoDS_Face& F)
 232 {
 233   Standard_Real First, Last;
 234
 235   Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface (E, F, First, Last);
 236   Handle(Geom_Surface) Sf = BRep_Tool::Surface(F);
 237
 238   gp_Pnt2d p;
 239   gp_Vec2d v;
 240   C2d->D1( 0.5*(First+Last), p, v);
 241   if (E.Orientation() != F.Orientation())
 242     v.Reverse();
 243   gp_Dir2d dir( -v.Y(), v.X() ); // dir inside F
 244
 245   Standard_Real duv = 1e-6; // this is not Ok and may give incorrect result if
 246   // resolutionUV of compared faces is very different. To have a good result,
 247   //it is necessary to get normal to faces at points equidistant from E in 3D
 248
 249   p.SetX( p.X() + dir.X()*duv );
 250   p.SetY( p.Y() + dir.Y()*duv );
 251
 252   gp_Pnt Ps;
 253   gp_Vec Norm, V1, V2, VV1, VV2;
 254   Sf->D1( p.X(), p.Y(), Ps, V1, V2);
 255   Norm = V1.Crossed(V2);
 256
 257   if (F.Orientation() == TopAbs_REVERSED )
 258     Norm.Reverse();
 259
 260   return Norm;
 261 }
 262
 263
 264 //=======================================================================
 265 //function : FindEinF
 266 //purpose  : find E in F
 267 //=======================================================================
 268
 269 static TopoDS_Edge FindEinF(const TopoDS_Edge& E,
 270                             const TopoDS_Face& F)
 271 {
 272   TopExp_Explorer expl (F, TopAbs_EDGE);
 273   for (; expl.More(); expl.Next())
 274     if( E.IsSame( expl.Current() ))
 275       return TopoDS::Edge(expl.Current());
 276   TopoDS_Edge nullE;
 277   return nullE;
 278 }
 279
 280 //=======================================================================
 281 //function : IsInside
 282 //purpose  : check if <F2> is inside <F1> by edge <E>.
 283 //           if <CountDot>, compute <Dot>: scalar production of
 284 //           normalized  vectors  pointing  inside  faces,  and
 285 //           check if faces are oriented well for sewing
 286 //=======================================================================
 287
 288 Standard_Boolean Partition_Loop3d::IsInside(const TopoDS_Edge& E,
 289                                             const TopoDS_Face& F1,
 290                                             const TopoDS_Face& F2,
 291                                             const Standard_Boolean CountDot,
 292                                             Standard_Real& Dot,
 293                                             Standard_Boolean& GoodOri)
 294 {
 295   Standard_Real f, l;
 296   gp_Pnt P;
 297   gp_Vec Vc1, Vc2, Vin1, Vin2, Nf1, Nf2;
 298   Handle(Geom_Curve) Curve = BRep_Tool::Curve(E,f,l);
 299   Curve->D1( 0.5*(f + l), P, Vc2);
 300   TopoDS_Edge E1, E2 = FindEinF (E, F2);
 301   if (E2.Orientation() == TopAbs_REVERSED ) Vc2.Reverse();
 302
 303   Nf1 = Normal(E,F1);
 304   Nf2 = Normal(E,F2);
 305
 306   Standard_Real sin =
 307     Nf1.CrossSquareMagnitude(Nf2) / Nf1.SquareMagnitude() / Nf2.SquareMagnitude();
 308   Standard_Boolean tangent = sin < 0.001;
 309
 310   Standard_Boolean inside = 0;
 311   if (tangent) {
 312     E1 = FindEinF (E, F1);
 313     gp_Vec NNf1 = NextNormal(E1,F1);
 314     gp_Vec NNf2 = NextNormal(E2,F2);
 315     Vin2 = NNf2.Crossed(Vc2);
 316     inside = Vin2 * NNf1 < 0;
 317   }
 318   else {
 319     Vin2 = Nf2.Crossed(Vc2);
 320     inside = Vin2 * Nf1 < 0;
 321   }
 322
 323   if (!CountDot) return inside;
 324
 325   if (tangent)
 326     Vin2 = Nf2.Crossed(Vc2);
 327   else
 328     E1 = FindEinF (E, F1);
 329
 330   Vc1 = Vc2;
 331   if (E1.Orientation() != E2.Orientation())
 332     Vc1.Reverse();
 333   Vin1 = Nf1.Crossed(Vc1);
 334
 335   if (tangent) {
 336     Standard_Real N1N2 = Nf1 * Nf2;
 337     GoodOri = (Vin2 * Vin1 < 0) ? N1N2 > 0 : N1N2 < 0;
 338   }
 339   else {
 340     Standard_Real V1N2 = Vin1 * Nf2;
 341     GoodOri = ( inside ? V1N2 <= 0 : V1N2 >= 0);
 342   }
 343
 344   Vin1.Normalize();
 345   Vin2.Normalize();
 346
 347   Dot = Vin2 * Vin1;
 348
 349   return inside;
 350 }
 351