src/NMTAlgo/NMTAlgo_Loop3d.cxx

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