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