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