src/GEOMUtils/GEOMUtils.cxx

   1 // Copyright (C) 2007-2016  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, or (at your option) any later version.
  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
  23 #include <Standard_Stream.hxx>
  24
  25 #include <GEOMUtils.hxx>
  26
  27 #include <Basics_OCCTVersion.hxx>
  28
  29 #include <OpUtil.hxx>
  30 #include <Utils_ExceptHandlers.hxx>
  31
  32 // OCCT Includes
  33 #include <BRepMesh_IncrementalMesh.hxx>
  34
  35 #include <BRepExtrema_DistShapeShape.hxx>
  36
  37 #include <BRep_Builder.hxx>
  38 #include <BRep_Tool.hxx>
  39 #include <BRepBndLib.hxx>
  40 #include <BRepGProp.hxx>
  41 #include <BRepTools.hxx>
  42
  43 #include <BRepClass3d_SolidClassifier.hxx>
  44
  45 #include <BRepBuilderAPI_MakeFace.hxx>
  46 #include <BRepBuilderAPI_Sewing.hxx>
  47
  48 #include <BRepCheck_Analyzer.hxx>
  49
  50 #include <Bnd_Box.hxx>
  51
  52 #include <BOPAlgo_ArgumentAnalyzer.hxx>
  53 #include <BOPTools_AlgoTools.hxx>
  54
  55 #include <TopAbs.hxx>
  56 #include <TopExp.hxx>
  57 #include <TopoDS.hxx>
  58 #include <TopoDS_Edge.hxx>
  59 #include <TopoDS_Face.hxx>
  60 #include <TopoDS_Shape.hxx>
  61 #include <TopoDS_Vertex.hxx>
  62 #include <TopoDS_Compound.hxx>
  63 #include <TopoDS_Iterator.hxx>
  64 #include <TopExp_Explorer.hxx>
  65 #include <TopTools_MapOfShape.hxx>
  66 #include <TopTools_ListOfShape.hxx>
  67 #include <TopTools_ListIteratorOfListOfShape.hxx>
  68 #include <TopTools_Array1OfShape.hxx>
  69
  70 #include <Geom_Circle.hxx>
  71 #include <Geom_Surface.hxx>
  72 #include <Geom_Plane.hxx>
  73 #include <Geom_SphericalSurface.hxx>
  74 #include <Geom_ToroidalSurface.hxx>
  75 #include <Geom_RectangularTrimmedSurface.hxx>
  76
  77 #include <GeomLProp_CLProps.hxx>
  78 #include <GeomLProp_SLProps.hxx>
  79
  80 #include <GProp_GProps.hxx>
  81 #include <GProp_PrincipalProps.hxx>
  82
  83 #include <TColStd_Array1OfReal.hxx>
  84
  85 #include <gp_Pln.hxx>
  86 #include <gp_Lin.hxx>
  87
  88 #include <ShapeAnalysis.hxx>
  89 #include <ShapeFix_Shape.hxx>
  90 #include <ShapeFix_ShapeTolerance.hxx>
  91
  92 #include <ProjLib.hxx>
  93 #include <ElSLib.hxx>
  94
  95 #include <vector>
  96 #include <sstream>
  97 #include <algorithm>
  98
  99 #if OCC_VERSION_LARGE >= 0x07010000
 100 #include <V3d_Coordinate.hxx>
 101 #endif
 102
 103 #include <Standard_Failure.hxx>
 104 #include <Standard_NullObject.hxx>
 105 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
 106
 107 #define MAX2(X, Y)    (Abs(X) > Abs(Y) ? Abs(X) : Abs(Y))
 108 #define MAX3(X, Y, Z) (MAX2(MAX2(X,Y), Z))
 109
 110 #define STD_SORT_ALGO 1
 111
 112 #define DEFAULT_TOLERANCE_TOLERANCE     1.e-02
 113 #define DEFAULT_MAX_TOLERANCE_TOLERANCE 1.e-06
 114
 115 // When the following macro is defined, ShapeFix_ShapeTolerance function is used to set max tolerance of curve
 116 // in GEOMUtils::FixShapeCurves function; otherwise less restrictive BRep_Builder::UpdateEdge/UpdateVertex
 117 // approach is used
 118 // VSR (29/12/2014): macro disabled
 119 //#define USE_LIMIT_TOLERANCE
 120
 121 namespace
 122 {
 123   /**
 124    * This function constructs and returns modified shape from the original one
 125    * for singular cases. It is used for the method GetMinDistanceSingular.
 126    *
 127    * \param theShape the original shape
 128    * \param theModifiedShape output parameter. The modified shape.
 129    * \param theAddDist output parameter. The added distance for modified shape.
 130    * \retval true if the shape is modified; false otherwise.
 131    *
 132    * \internal
 133    */
 134   Standard_Boolean ModifyShape(const TopoDS_Shape  &theShape,
 135                                TopoDS_Shape  &theModifiedShape,
 136                                Standard_Real &theAddDist)
 137   {
 138     TopExp_Explorer anExp;
 139     int nbf = 0;
 140
 141     theAddDist = 0.;
 142     theModifiedShape.Nullify();
 143
 144     for ( anExp.Init( theShape, TopAbs_FACE ); anExp.More(); anExp.Next() ) {
 145       nbf++;
 146       theModifiedShape = anExp.Current();
 147     }
 148     if(nbf==1) {
 149       TopoDS_Shape sh = theShape;
 150       while(sh.ShapeType()==TopAbs_COMPOUND) {
 151         TopoDS_Iterator it(sh);
 152         sh = it.Value();
 153       }
 154       Handle(Geom_Surface) S = BRep_Tool::Surface(TopoDS::Face(theModifiedShape));
 155       if( S->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) ||
 156           S->IsKind(STANDARD_TYPE(Geom_ToroidalSurface)) ||
 157           S->IsUPeriodic()) {
 158         const Standard_Boolean isShell =
 159           (sh.ShapeType()==TopAbs_SHELL || sh.ShapeType()==TopAbs_FACE);
 160
 161         if ( !isShell && S->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) ) {
 162           Handle(Geom_SphericalSurface) SS = Handle(Geom_SphericalSurface)::DownCast(S);
 163           gp_Pnt PC = SS->Location();
 164           BRep_Builder B;
 165           TopoDS_Vertex V;
 166           B.MakeVertex(V,PC,1.e-7);
 167           theModifiedShape = V;
 168           theAddDist = SS->Radius();
 169           return Standard_True;
 170         }
 171         if ( !isShell && S->IsKind(STANDARD_TYPE(Geom_ToroidalSurface)) ) {
 172           Handle(Geom_ToroidalSurface) TS = Handle(Geom_ToroidalSurface)::DownCast(S);
 173           gp_Ax3 ax3 = TS->Position();
 174           Handle(Geom_Circle) C = new Geom_Circle(ax3.Ax2(),TS->MajorRadius());
 175           BRep_Builder B;
 176           TopoDS_Edge E;
 177           B.MakeEdge(E,C,1.e-7);
 178           theModifiedShape = E;
 179           theAddDist = TS->MinorRadius();
 180           return Standard_True;
 181         }
 182
 183         // non solid case or any periodic surface (Mantis 22454).
 184         double U1,U2,V1,V2;
 185         // changes for 0020677: EDF 1219 GEOM: MinDistance gives 0 instead of 20.88
 186         //S->Bounds(U1,U2,V1,V2); changed by
 187         ShapeAnalysis::GetFaceUVBounds(TopoDS::Face(theModifiedShape),U1,U2,V1,V2);
 188         // end of changes for 020677 (dmv)
 189         Handle(Geom_RectangularTrimmedSurface) TrS1 =
 190           new Geom_RectangularTrimmedSurface(S,U1,(U1+U2)/2.,V1,V2);
 191         Handle(Geom_RectangularTrimmedSurface) TrS2 =
 192           new Geom_RectangularTrimmedSurface(S,(U1+U2)/2.,U2,V1,V2);
 193         TopoDS_Shape aMShape;
 194
 195         TopoDS_Face F1 = BRepBuilderAPI_MakeFace(TrS1, Precision::Confusion());
 196         TopoDS_Face F2 = BRepBuilderAPI_MakeFace(TrS2, Precision::Confusion());
 197
 198         if (isShell) {
 199           BRep_Builder B;
 200           B.MakeCompound(TopoDS::Compound(aMShape));
 201           B.Add(aMShape, F1);
 202           B.Add(aMShape, F2);
 203         } else {
 204           // The original shape is a solid.
 205           BRepBuilderAPI_Sewing aSewing (Precision::Confusion()*10.0);
 206           aSewing.Add(F1);
 207           aSewing.Add(F2);
 208           aSewing.Perform();
 209           aMShape = aSewing.SewedShape();
 210           BRep_Builder B;
 211           TopoDS_Solid aSolid;
 212           B.MakeSolid(aSolid);
 213           B.Add(aSolid, aMShape);
 214           aMShape = aSolid;
 215         }
 216
 217         Handle(ShapeFix_Shape) sfs = new ShapeFix_Shape;
 218         sfs->Init(aMShape);
 219         sfs->SetPrecision(1.e-6);
 220         sfs->SetMaxTolerance(1.0);
 221         sfs->Perform();
 222         theModifiedShape = sfs->Shape();
 223         return Standard_True;
 224       }
 225     }
 226
 227     theModifiedShape = theShape;
 228     return Standard_False;
 229   }
 230
 231   void parseWard( const GEOMUtils::LevelsList &theLevelList, std::string &treeStr )
 232   {
 233     treeStr.append( "{" );
 234     for( GEOMUtils::LevelsList::const_iterator j = theLevelList.begin();
 235          j != theLevelList.end(); ++j ) {
 236       if ( j != theLevelList.begin() ) {
 237         treeStr.append( ";" );
 238       }
 239       GEOMUtils::LevelInfo level = (*j);
 240       GEOMUtils::LevelInfo::iterator upIter;
 241       for ( upIter = level.begin(); upIter != level.end(); ++upIter ) {
 242         if ( upIter != level.begin() ) {
 243           treeStr.append( "," );
 244         }
 245         treeStr.append( upIter->first );
 246         for ( std::vector<std::string>::iterator k = upIter->second.begin(); k != upIter->second.end(); ++k ) {
 247           treeStr.append( "_" );
 248           treeStr.append( *k );
 249         }
 250       }
 251     }
 252     treeStr.append( "}" );
 253   }
 254
 255   GEOMUtils::LevelsList parseWard( const std::string& theData, std::size_t& theCursor )
 256   {
 257     std::size_t indexStart = theData.find( "{", theCursor ) + 1;
 258     std::size_t indexEnd = theData.find( "}", indexStart );
 259
 260     std::string ward = theData.substr( indexStart, indexEnd - indexStart );
 261     std::stringstream ss(ward);
 262     std::string substr;
 263     std::vector<std::string> levelsListStr;
 264     while ( std::getline( ss, substr, ';' ) ) {
 265       if ( !substr.empty() )
 266         levelsListStr.push_back( substr );
 267     }
 268     GEOMUtils::LevelsList levelsListData;
 269     for( size_t level = 0; level < levelsListStr.size(); level++ ) {
 270       std::vector<std::string> namesListStr;
 271       std::stringstream ss1( levelsListStr[level] );
 272       while ( std::getline( ss1, substr, ',' ) ) {
 273         if ( !substr.empty() )
 274           namesListStr.push_back( substr );
 275       }
 276       GEOMUtils::LevelInfo levelInfoData;
 277       for( size_t node = 0; node < namesListStr.size(); node++ ) {
 278         std::vector<std::string> linksListStr;
 279         std::stringstream ss2( namesListStr[node] );
 280         while ( std::getline( ss2, substr, '_' ) ) {
 281           if ( !substr.empty() )
 282             linksListStr.push_back( substr );
 283         }
 284         std::string nodeItem = linksListStr[0];
 285         if( !nodeItem.empty() ) {
 286           GEOMUtils::NodeLinks linksListData;
 287           for( size_t link = 1; link < linksListStr.size(); link++ ) {
 288             std::string linkItem = linksListStr[link];
 289             linksListData.push_back( linkItem );
 290           }// Links
 291           levelInfoData[nodeItem] = linksListData;
 292         }
 293       }// Level's objects
 294       levelsListData.push_back(levelInfoData);
 295     }// Levels
 296
 297     theCursor = indexEnd + 1;
 298     return levelsListData;
 299   }
 300
 301 }
 302
 303 //=======================================================================
 304 //function : ShapeToDouble
 305 //purpose  : used by CompareShapes::operator()
 306 //=======================================================================
 307 std::pair<double, double> GEOMUtils::ShapeToDouble (const TopoDS_Shape& S, bool isOldSorting)
 308 {
 309   // Computing of CentreOfMass
 310   gp_Pnt GPoint;
 311   double Len;
 312
 313   if (S.ShapeType() == TopAbs_VERTEX) {
 314     GPoint = BRep_Tool::Pnt(TopoDS::Vertex(S));
 315     Len = (double)S.Orientation();
 316   }
 317   else {
 318     GProp_GProps GPr;
 319     // BEGIN: fix for Mantis issue 0020842
 320     if (isOldSorting) {
 321       BRepGProp::LinearProperties(S, GPr);
 322     }
 323     else {
 324       if (S.ShapeType() == TopAbs_EDGE || S.ShapeType() == TopAbs_WIRE) {
 325         BRepGProp::LinearProperties(S, GPr);
 326       }
 327       else if (S.ShapeType() == TopAbs_FACE || S.ShapeType() == TopAbs_SHELL) {
 328         BRepGProp::SurfaceProperties(S, GPr);
 329       }
 330       else {
 331         BRepGProp::VolumeProperties(S, GPr);
 332       }
 333     }
 334     // END: fix for Mantis issue 0020842
 335     GPoint = GPr.CentreOfMass();
 336     Len = GPr.Mass();
 337   }
 338
 339   double dMidXYZ = GPoint.X() * 999.0 + GPoint.Y() * 99.0 + GPoint.Z() * 0.9;
 340   return std::make_pair(dMidXYZ, Len);
 341 }
 342
 343 //=======================================================================
 344 //function : GetPosition
 345 //purpose  :
 346 //=======================================================================
 347 gp_Ax3 GEOMUtils::GetPosition (const TopoDS_Shape& theShape)
 348 {
 349   gp_Ax3 aResult;
 350
 351   if (theShape.IsNull())
 352     return aResult;
 353
 354   // Axes
 355   aResult.Transform(theShape.Location().Transformation());
 356   if (theShape.ShapeType() == TopAbs_FACE) {
 357     Handle(Geom_Surface) aGS = BRep_Tool::Surface(TopoDS::Face(theShape));
 358     if (!aGS.IsNull() && aGS->IsKind(STANDARD_TYPE(Geom_Plane))) {
 359       Handle(Geom_Plane) aGPlane = Handle(Geom_Plane)::DownCast(aGS);
 360       gp_Pln aPln = aGPlane->Pln();
 361       aResult = aPln.Position();
 362       // In case of reverse orinetation of the face invert the plane normal
 363       // (the face's normal does not mathc the plane's normal in this case)
 364       if(theShape.Orientation() == TopAbs_REVERSED)
 365       {
 366         gp_Dir Vx =  aResult.XDirection();
 367         gp_Dir N  =  aResult.Direction().Mirrored(Vx);
 368         gp_Pnt P  =  aResult.Location();
 369         aResult = gp_Ax3(P, N, Vx);
 370       }
 371     }
 372   }
 373
 374   // Origin
 375   gp_Pnt aPnt;
 376
 377   TopAbs_ShapeEnum aShType = theShape.ShapeType();
 378
 379   if (aShType == TopAbs_VERTEX) {
 380     aPnt = BRep_Tool::Pnt(TopoDS::Vertex(theShape));
 381   }
 382   else {
 383     if (aShType == TopAbs_COMPOUND) {
 384       aShType = GetTypeOfSimplePart(theShape);
 385     }
 386
 387     GProp_GProps aSystem;
 388     if (aShType == TopAbs_EDGE || aShType == TopAbs_WIRE)
 389       BRepGProp::LinearProperties(theShape, aSystem);
 390     else if (aShType == TopAbs_FACE || aShType == TopAbs_SHELL)
 391       BRepGProp::SurfaceProperties(theShape, aSystem);
 392     else
 393       BRepGProp::VolumeProperties(theShape, aSystem);
 394
 395     aPnt = aSystem.CentreOfMass();
 396   }
 397
 398   aResult.SetLocation(aPnt);
 399
 400   return aResult;
 401 }
 402
 403 //=======================================================================
 404 //function : GetVector
 405 //purpose  :
 406 //=======================================================================
 407 gp_Vec GEOMUtils::GetVector (const TopoDS_Shape& theShape,
 408                              Standard_Boolean doConsiderOrientation)
 409 {
 410   if (theShape.IsNull())
 411     Standard_NullObject::Raise("Null shape is given for a vector");
 412
 413   if (theShape.ShapeType() != TopAbs_EDGE)
 414     Standard_TypeMismatch::Raise("Invalid shape is given, must be a vector or an edge");
 415
 416   TopoDS_Edge anE = TopoDS::Edge(theShape);
 417
 418   TopoDS_Vertex V1, V2;
 419   TopExp::Vertices(anE, V1, V2, doConsiderOrientation);
 420
 421   if (V1.IsNull() || V2.IsNull())
 422     Standard_NullObject::Raise("Invalid edge is given, it must have two points");
 423
 424   gp_Vec aV (BRep_Tool::Pnt(V1), BRep_Tool::Pnt(V2));
 425   if (aV.Magnitude() < gp::Resolution()) {
 426     Standard_ConstructionError::Raise("Vector of zero length is given");
 427   }
 428
 429   return aV;
 430 }
 431
 432 //=======================================================================
 433 //function : CompareShapes::operator()
 434 //purpose  : used by std::sort(), called from SortShapes()
 435 //=======================================================================
 436 bool GEOMUtils::CompareShapes::operator() (const TopoDS_Shape& theShape1,
 437                                            const TopoDS_Shape& theShape2)
 438 {
 439   if (!myMap.IsBound(theShape1)) {
 440     myMap.Bind(theShape1, ShapeToDouble(theShape1, myIsOldSorting));
 441   }
 442
 443   if (!myMap.IsBound(theShape2)) {
 444     myMap.Bind(theShape2, ShapeToDouble(theShape2, myIsOldSorting));
 445   }
 446
 447   std::pair<double, double> val1 = myMap.Find(theShape1);
 448   std::pair<double, double> val2 = myMap.Find(theShape2);
 449
 450   double tol = Precision::Confusion();
 451   bool exchange = Standard_False;
 452
 453   double dMidXYZ = val1.first - val2.first;
 454   if (dMidXYZ >= tol) {
 455     exchange = Standard_True;
 456   }
 457   else if (Abs(dMidXYZ) < tol) {
 458     double dLength = val1.second - val2.second;
 459     if (dLength >= tol) {
 460       exchange = Standard_True;
 461     }
 462     else if (Abs(dLength) < tol && theShape1.ShapeType() <= TopAbs_FACE) {
 463       // PAL17233
 464       // equal values possible on shapes such as two halves of a sphere and
 465       // a membrane inside the sphere
 466       Bnd_Box box1,box2;
 467       BRepBndLib::Add(theShape1, box1);
 468       if (!box1.IsVoid()) {
 469         BRepBndLib::Add(theShape2, box2);
 470         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
 471         if (dSquareExtent >= tol) {
 472           exchange = Standard_True;
 473         }
 474         else if (Abs(dSquareExtent) < tol) {
 475           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
 476           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 477           val1 = (aXmin+aXmax)*999.0 + (aYmin+aYmax)*99.0 + (aZmin+aZmax)*0.9;
 478           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 479           val2 = (aXmin+aXmax)*999.0 + (aYmin+aYmax)*99.0 + (aZmin+aZmax)*0.9;
 480           if ((val1 - val2) >= tol) {
 481             exchange = Standard_True;
 482           }
 483         }
 484       }
 485     }
 486   }
 487
 488   //return val1 < val2;
 489   return !exchange;
 490 }
 491
 492 //=======================================================================
 493 //function : SortShapes
 494 //purpose  :
 495 //=======================================================================
 496 void GEOMUtils::SortShapes (TopTools_ListOfShape& SL,
 497                             const Standard_Boolean isOldSorting)
 498 {
 499 #ifdef STD_SORT_ALGO
 500   std::vector<TopoDS_Shape> aShapesVec;
 501   aShapesVec.reserve(SL.Extent());
 502
 503   TopTools_ListIteratorOfListOfShape it (SL);
 504   for (; it.More(); it.Next()) {
 505     aShapesVec.push_back(it.Value());
 506   }
 507   SL.Clear();
 508
 509   CompareShapes shComp (isOldSorting);
 510   std::stable_sort(aShapesVec.begin(), aShapesVec.end(), shComp);
 511   //std::sort(aShapesVec.begin(), aShapesVec.end(), shComp);
 512
 513   std::vector<TopoDS_Shape>::const_iterator anIter = aShapesVec.begin();
 514   for (; anIter != aShapesVec.end(); ++anIter) {
 515     SL.Append(*anIter);
 516   }
 517 #else
 518   // old implementation
 519   Standard_Integer MaxShapes = SL.Extent();
 520   TopTools_Array1OfShape  aShapes (1,MaxShapes);
 521   TColStd_Array1OfInteger OrderInd(1,MaxShapes);
 522   TColStd_Array1OfReal    MidXYZ  (1,MaxShapes); //X,Y,Z;
 523   TColStd_Array1OfReal    Length  (1,MaxShapes); //X,Y,Z;
 524
 525   // Computing of CentreOfMass
 526   Standard_Integer Index;
 527   GProp_GProps GPr;
 528   gp_Pnt GPoint;
 529   TopTools_ListIteratorOfListOfShape it(SL);
 530   for (Index=1;  it.More();  Index++)
 531   {
 532     TopoDS_Shape S = it.Value();
 533     SL.Remove( it ); // == it.Next()
 534     aShapes(Index) = S;
 535     OrderInd.SetValue (Index, Index);
 536     if (S.ShapeType() == TopAbs_VERTEX) {
 537       GPoint = BRep_Tool::Pnt( TopoDS::Vertex( S ));
 538       Length.SetValue( Index, (Standard_Real) S.Orientation());
 539     }
 540     else {
 541       // BEGIN: fix for Mantis issue 0020842
 542       if (isOldSorting) {
 543         BRepGProp::LinearProperties (S, GPr);
 544       }
 545       else {
 546         if (S.ShapeType() == TopAbs_EDGE || S.ShapeType() == TopAbs_WIRE) {
 547           BRepGProp::LinearProperties (S, GPr);
 548         }
 549         else if (S.ShapeType() == TopAbs_FACE || S.ShapeType() == TopAbs_SHELL) {
 550           BRepGProp::SurfaceProperties(S, GPr);
 551         }
 552         else {
 553           BRepGProp::VolumeProperties(S, GPr);
 554         }
 555       }
 556       // END: fix for Mantis issue 0020842
 557       GPoint = GPr.CentreOfMass();
 558       Length.SetValue(Index, GPr.Mass());
 559     }
 560     MidXYZ.SetValue(Index, GPoint.X()*999.0 + GPoint.Y()*99.0 + GPoint.Z()*0.9);
 561     //cout << Index << " L: " << Length(Index) << "CG: " << MidXYZ(Index) << endl;
 562   }
 563
 564   // Sorting
 565   Standard_Integer aTemp;
 566   Standard_Boolean exchange, Sort = Standard_True;
 567   Standard_Real    tol = Precision::Confusion();
 568   while (Sort)
 569   {
 570     Sort = Standard_False;
 571     for (Index=1; Index < MaxShapes; Index++)
 572     {
 573       exchange = Standard_False;
 574       Standard_Real dMidXYZ = MidXYZ(OrderInd(Index)) - MidXYZ(OrderInd(Index+1));
 575       Standard_Real dLength = Length(OrderInd(Index)) - Length(OrderInd(Index+1));
 576       if ( dMidXYZ >= tol ) {
 577 //         cout << "MidXYZ: " << MidXYZ(OrderInd(Index))<< " > " <<MidXYZ(OrderInd(Index+1))
 578 //              << " d: " << dMidXYZ << endl;
 579         exchange = Standard_True;
 580       }
 581       else if ( Abs(dMidXYZ) < tol && dLength >= tol ) {
 582 //         cout << "Length: " << Length(OrderInd(Index))<< " > " <<Length(OrderInd(Index+1))
 583 //              << " d: " << dLength << endl;
 584         exchange = Standard_True;
 585       }
 586       else if ( Abs(dMidXYZ) < tol && Abs(dLength) < tol &&
 587                 aShapes(OrderInd(Index)).ShapeType() <= TopAbs_FACE) {
 588         // PAL17233
 589         // equal values possible on shapes such as two halves of a sphere and
 590         // a membrane inside the sphere
 591         Bnd_Box box1,box2;
 592         BRepBndLib::Add( aShapes( OrderInd(Index) ), box1 );
 593         if ( box1.IsVoid() ) continue;
 594         BRepBndLib::Add( aShapes( OrderInd(Index+1) ), box2 );
 595         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
 596         if ( dSquareExtent >= tol ) {
 597 //           cout << "SquareExtent: " << box1.SquareExtent()<<" > "<<box2.SquareExtent() << endl;
 598           exchange = Standard_True;
 599         }
 600         else if ( Abs(dSquareExtent) < tol ) {
 601           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
 602           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 603           val1 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
 604           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 605           val2 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
 606           //exchange = val1 > val2;
 607           if ((val1 - val2) >= tol) {
 608             exchange = Standard_True;
 609           }
 610           //cout << "box: " << val1<<" > "<<val2 << endl;
 611         }
 612       }
 613
 614       if (exchange)
 615       {
 616 //         cout << "exchange " << Index << " & " << Index+1 << endl;
 617         aTemp = OrderInd(Index);
 618         OrderInd(Index) = OrderInd(Index+1);
 619         OrderInd(Index+1) = aTemp;
 620         Sort = Standard_True;
 621       }
 622     }
 623   }
 624
 625   for (Index=1; Index <= MaxShapes; Index++)
 626     SL.Append( aShapes( OrderInd(Index) ));
 627 #endif
 628 }
 629
 630 //=======================================================================
 631 //function : CompsolidToCompound
 632 //purpose  :
 633 //=======================================================================
 634 TopoDS_Shape GEOMUtils::CompsolidToCompound (const TopoDS_Shape& theCompsolid)
 635 {
 636   if (theCompsolid.ShapeType() != TopAbs_COMPSOLID) {
 637     return theCompsolid;
 638   }
 639
 640   TopoDS_Compound aCompound;
 641   BRep_Builder B;
 642   B.MakeCompound(aCompound);
 643
 644   TopTools_MapOfShape mapShape;
 645   TopoDS_Iterator It (theCompsolid, Standard_True, Standard_True);
 646
 647   for (; It.More(); It.Next()) {
 648     TopoDS_Shape aShape_i = It.Value();
 649     if (mapShape.Add(aShape_i)) {
 650       B.Add(aCompound, aShape_i);
 651     }
 652   }
 653
 654   return aCompound;
 655 }
 656
 657 //=======================================================================
 658 //function : AddSimpleShapes
 659 //purpose  :
 660 //=======================================================================
 661 void GEOMUtils::AddSimpleShapes (const TopoDS_Shape& theShape, TopTools_ListOfShape& theList)
 662 {
 663   if (theShape.ShapeType() != TopAbs_COMPOUND &&
 664       theShape.ShapeType() != TopAbs_COMPSOLID) {
 665     theList.Append(theShape);
 666     return;
 667   }
 668
 669   TopTools_MapOfShape mapShape;
 670   TopoDS_Iterator It (theShape, Standard_True, Standard_True);
 671
 672   for (; It.More(); It.Next()) {
 673     TopoDS_Shape aShape_i = It.Value();
 674     if (mapShape.Add(aShape_i)) {
 675       if (aShape_i.ShapeType() == TopAbs_COMPOUND ||
 676           aShape_i.ShapeType() == TopAbs_COMPSOLID) {
 677         AddSimpleShapes(aShape_i, theList);
 678       } else {
 679         theList.Append(aShape_i);
 680       }
 681     }
 682   }
 683 }
 684
 685 //=======================================================================
 686 //function : CheckTriangulation
 687 //purpose  :
 688 //=======================================================================
 689 bool GEOMUtils::CheckTriangulation (const TopoDS_Shape& aShape)
 690 {
 691   bool isTriangulation = true;
 692
 693   TopExp_Explorer exp (aShape, TopAbs_FACE);
 694   if (exp.More())
 695   {
 696     TopLoc_Location aTopLoc;
 697     Handle(Poly_Triangulation) aTRF;
 698     aTRF = BRep_Tool::Triangulation(TopoDS::Face(exp.Current()), aTopLoc);
 699     if (aTRF.IsNull()) {
 700       isTriangulation = false;
 701     }
 702   }
 703   else // no faces, try edges
 704   {
 705     TopExp_Explorer expe (aShape, TopAbs_EDGE);
 706     if (!expe.More()) {
 707       return false;
 708     }
 709     TopLoc_Location aLoc;
 710     Handle(Poly_Polygon3D) aPE = BRep_Tool::Polygon3D(TopoDS::Edge(expe.Current()), aLoc);
 711     if (aPE.IsNull()) {
 712       isTriangulation = false;
 713     }
 714   }
 715
 716   if (!isTriangulation) {
 717     // calculate deflection
 718     Standard_Real aDeviationCoefficient = 0.001;
 719
 720     Bnd_Box B;
 721     BRepBndLib::Add(aShape, B);
 722     Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
 723     B.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 724
 725     Standard_Real dx = aXmax - aXmin, dy = aYmax - aYmin, dz = aZmax - aZmin;
 726     Standard_Real aDeflection = Max(Max(dx, dy), dz) * aDeviationCoefficient * 4;
 727     Standard_Real aHLRAngle = 0.349066;
 728
 729     BRepMesh_IncrementalMesh Inc (aShape, aDeflection, Standard_False, aHLRAngle);
 730   }
 731
 732   return true;
 733 }
 734
 735 //=======================================================================
 736 //function : GetTypeOfSimplePart
 737 //purpose  :
 738 //=======================================================================
 739 TopAbs_ShapeEnum GEOMUtils::GetTypeOfSimplePart (const TopoDS_Shape& theShape)
 740 {
 741   TopAbs_ShapeEnum aType = theShape.ShapeType();
 742   if      (aType == TopAbs_VERTEX)                             return TopAbs_VERTEX;
 743   else if (aType == TopAbs_EDGE  || aType == TopAbs_WIRE)      return TopAbs_EDGE;
 744   else if (aType == TopAbs_FACE  || aType == TopAbs_SHELL)     return TopAbs_FACE;
 745   else if (aType == TopAbs_SOLID || aType == TopAbs_COMPSOLID) return TopAbs_SOLID;
 746   else if (aType == TopAbs_COMPOUND) {
 747     // Only the iType of the first shape in the compound is taken into account
 748     TopoDS_Iterator It (theShape, Standard_False, Standard_False);
 749     if (It.More()) {
 750       return GetTypeOfSimplePart(It.Value());
 751     }
 752   }
 753   return TopAbs_SHAPE;
 754 }
 755
 756 //=======================================================================
 757 //function : GetEdgeNearPoint
 758 //purpose  :
 759 //=======================================================================
 760 TopoDS_Shape GEOMUtils::GetEdgeNearPoint (const TopoDS_Shape& theShape,
 761                                           const TopoDS_Vertex& thePoint)
 762 {
 763   TopoDS_Shape aResult;
 764
 765   // 1. Explode the shape on edges
 766   TopTools_MapOfShape mapShape;
 767   Standard_Integer nbEdges = 0;
 768   TopExp_Explorer exp (theShape, TopAbs_EDGE);
 769   for (; exp.More(); exp.Next()) {
 770     if (mapShape.Add(exp.Current())) {
 771       nbEdges++;
 772     }
 773   }
 774
 775   if (nbEdges == 0)
 776     Standard_NullObject::Raise("Given shape contains no edges");
 777
 778   mapShape.Clear();
 779   Standard_Integer ind = 1;
 780   TopTools_Array1OfShape anEdges (1, nbEdges);
 781   TColStd_Array1OfReal aDistances (1, nbEdges);
 782   for (exp.Init(theShape, TopAbs_EDGE); exp.More(); exp.Next()) {
 783     if (mapShape.Add(exp.Current())) {
 784       TopoDS_Shape anEdge = exp.Current();
 785       anEdges(ind) = anEdge;
 786
 787       // 2. Classify the point relatively each edge
 788       BRepExtrema_DistShapeShape aDistTool (thePoint, anEdges(ind));
 789       if (!aDistTool.IsDone())
 790         Standard_ConstructionError::Raise("Cannot find a distance from the given point to one of edges");
 791
 792       aDistances(ind) = aDistTool.Value();
 793       ind++;
 794     }
 795   }
 796
 797   // 3. Define edge, having minimum distance to the point
 798   Standard_Real nearest = RealLast(), nbFound = 0;
 799   Standard_Real prec = Precision::Confusion();
 800   for (ind = 1; ind <= nbEdges; ind++) {
 801     if (Abs(aDistances(ind) - nearest) < prec) {
 802       nbFound++;
 803     }
 804     else if (aDistances(ind) < nearest) {
 805       nearest = aDistances(ind);
 806       aResult = anEdges(ind);
 807       nbFound = 1;
 808     }
 809     else {
 810     }
 811   }
 812   if (nbFound > 1) {
 813     Standard_ConstructionError::Raise("Multiple edges near the given point are found");
 814   }
 815   else if (nbFound == 0) {
 816     Standard_ConstructionError::Raise("There are no edges near the given point");
 817   }
 818   else {
 819   }
 820
 821   return aResult;
 822 }
 823
 824 //=======================================================================
 825 //function : PreciseBoundingBox
 826 //purpose  :
 827 //=======================================================================
 828 Standard_Boolean GEOMUtils::PreciseBoundingBox
 829                           (const TopoDS_Shape &theShape, Bnd_Box &theBox)
 830 {
 831   if ( theBox.IsVoid() ) BRepBndLib::Add( theShape, theBox );
 832   if ( theBox.IsVoid() ) return Standard_False;
 833
 834   Standard_Real aBound[6];
 835   theBox.Get(aBound[0], aBound[2], aBound[4], aBound[1], aBound[3], aBound[5]);
 836
 837   Standard_Integer i;
 838   const gp_Pnt aMid(0.5*(aBound[1] + aBound[0]),  // XMid
 839                     0.5*(aBound[3] + aBound[2]),  // YMid
 840                     0.5*(aBound[5] + aBound[4])); // ZMid
 841   const gp_XYZ aSize(aBound[1] - aBound[0],       // DX
 842                      aBound[3] - aBound[2],       // DY
 843                      aBound[5] - aBound[4]);      // DZ
 844   const gp_Pnt aPnt[6] =
 845     {
 846       gp_Pnt(aBound[0] - (aBound[1] - aBound[0]), aMid.Y(), aMid.Z()), // XMin
 847       gp_Pnt(aBound[1] + (aBound[1] - aBound[0]), aMid.Y(), aMid.Z()), // XMax
 848       gp_Pnt(aMid.X(), aBound[2] - (aBound[3] - aBound[2]), aMid.Z()), // YMin
 849       gp_Pnt(aMid.X(), aBound[3] + (aBound[3] - aBound[2]), aMid.Z()), // YMax
 850       gp_Pnt(aMid.X(), aMid.Y(), aBound[4] - (aBound[5] - aBound[4])), // ZMin
 851       gp_Pnt(aMid.X(), aMid.Y(), aBound[5] + (aBound[5] - aBound[4]))  // ZMax
 852     };
 853   const gp_Dir aDir[3] = { gp::DX(), gp::DY(), gp::DZ() };
 854   const Standard_Real aPlnSize[3] =
 855     {
 856       0.5*Max(aSize.Y(), aSize.Z()), // XMin, XMax planes
 857       0.5*Max(aSize.X(), aSize.Z()), // YMin, YMax planes
 858       0.5*Max(aSize.X(), aSize.Y())  // ZMin, ZMax planes
 859     };
 860   gp_Pnt aPMin[2];
 861
 862   for (i = 0; i < 6; i++) {
 863     const Standard_Integer iHalf = i/2;
 864     const gp_Pln aPln(aPnt[i], aDir[iHalf]);
 865     BRepBuilderAPI_MakeFace aMkFace(aPln, -aPlnSize[iHalf], aPlnSize[iHalf],
 866                                     -aPlnSize[iHalf], aPlnSize[iHalf]);
 867
 868     if (!aMkFace.IsDone()) {
 869       return Standard_False;
 870     }
 871
 872     TopoDS_Shape aFace = aMkFace.Shape();
 873
 874     // Get minimal distance between planar face and shape.
 875     Standard_Real aMinDist =
 876       GEOMUtils::GetMinDistance(aFace, theShape, aPMin[0], aPMin[1]);
 877
 878     if (aMinDist < 0.) {
 879       return Standard_False;
 880     }
 881
 882     aBound[i] = aPMin[1].Coord(iHalf + 1);
 883   }
 884
 885   // Update Bounding box with the new values.
 886   theBox.SetVoid();
 887   theBox.Update(aBound[0], aBound[2], aBound[4], aBound[1], aBound[3], aBound[5]);
 888
 889   return Standard_True;
 890 }
 891
 892 //=======================================================================
 893 //function : GetMinDistanceSingular
 894 //purpose  :
 895 //=======================================================================
 896 double GEOMUtils::GetMinDistanceSingular(const TopoDS_Shape& aSh1,
 897                                          const TopoDS_Shape& aSh2,
 898                                          gp_Pnt& Ptmp1, gp_Pnt& Ptmp2)
 899 {
 900   TopoDS_Shape     tmpSh1;
 901   TopoDS_Shape     tmpSh2;
 902   Standard_Real    AddDist1 = 0.;
 903   Standard_Real    AddDist2 = 0.;
 904   Standard_Boolean IsChange1 = ModifyShape(aSh1, tmpSh1, AddDist1);
 905   Standard_Boolean IsChange2 = ModifyShape(aSh2, tmpSh2, AddDist2);
 906
 907   if( !IsChange1 && !IsChange2 )
 908     return -2.0;
 909
 910   BRepExtrema_DistShapeShape dst(tmpSh1,tmpSh2);
 911   if (dst.IsDone()) {
 912     double MinDist = 1.e9;
 913     gp_Pnt PMin1, PMin2, P1, P2;
 914     for (int i = 1; i <= dst.NbSolution(); i++) {
 915       P1 = dst.PointOnShape1(i);
 916       P2 = dst.PointOnShape2(i);
 917       Standard_Real Dist = P1.Distance(P2);
 918       if (MinDist > Dist) {
 919         MinDist = Dist;
 920         PMin1 = P1;
 921         PMin2 = P2;
 922       }
 923     }
 924     if(MinDist<1.e-7) {
 925       Ptmp1 = PMin1;
 926       Ptmp2 = PMin2;
 927     }
 928     else {
 929       gp_Dir aDir(gp_Vec(PMin1,PMin2));
 930       if( MinDist > (AddDist1+AddDist2) ) {
 931         Ptmp1 = gp_Pnt( PMin1.X() + aDir.X()*AddDist1,
 932                         PMin1.Y() + aDir.Y()*AddDist1,
 933                         PMin1.Z() + aDir.Z()*AddDist1 );
 934         Ptmp2 = gp_Pnt( PMin2.X() - aDir.X()*AddDist2,
 935                         PMin2.Y() - aDir.Y()*AddDist2,
 936                         PMin2.Z() - aDir.Z()*AddDist2 );
 937         return (MinDist - AddDist1 - AddDist2);
 938       }
 939       else {
 940         if( AddDist1 > 0 ) {
 941           Ptmp1 = gp_Pnt( PMin1.X() + aDir.X()*AddDist1,
 942                           PMin1.Y() + aDir.Y()*AddDist1,
 943                           PMin1.Z() + aDir.Z()*AddDist1 );
 944           Ptmp2 = Ptmp1;
 945         }
 946         else {
 947           Ptmp2 = gp_Pnt( PMin2.X() - aDir.X()*AddDist2,
 948                           PMin2.Y() - aDir.Y()*AddDist2,
 949                           PMin2.Z() - aDir.Z()*AddDist2 );
 950           Ptmp1 = Ptmp2;
 951         }
 952       }
 953     }
 954     double res = MinDist - AddDist1 - AddDist2;
 955     if(res<0.) res = 0.0;
 956     return res;
 957   }
 958   return -2.0;
 959 }
 960
 961 //=======================================================================
 962 //function : GetMinDistance
 963 //purpose  :
 964 //=======================================================================
 965 Standard_Real GEOMUtils::GetMinDistance
 966                                (const TopoDS_Shape& theShape1,
 967                                 const TopoDS_Shape& theShape2,
 968                                 gp_Pnt& thePnt1, gp_Pnt& thePnt2)
 969 {
 970   Standard_Real aResult = 1.e9;
 971
 972   // Issue 0020231: A min distance bug with torus and vertex.
 973   // Make GetMinDistance() return zero if a sole VERTEX is inside any of SOLIDs
 974
 975   // which of shapes consists of only one vertex?
 976   TopExp_Explorer exp1(theShape1,TopAbs_VERTEX), exp2(theShape2,TopAbs_VERTEX);
 977   TopoDS_Shape V1 = exp1.More() ? exp1.Current() : TopoDS_Shape();
 978   TopoDS_Shape V2 = exp2.More() ? exp2.Current() : TopoDS_Shape();
 979   exp1.Next(); exp2.Next();
 980   if ( exp1.More() ) V1.Nullify();
 981   if ( exp2.More() ) V2.Nullify();
 982   // vertex and container of solids
 983   TopoDS_Shape V = V1.IsNull() ? V2 : V1;
 984   TopoDS_Shape S = V1.IsNull() ? theShape1 : theShape2;
 985   if ( !V.IsNull() ) {
 986     // classify vertex against solids
 987     gp_Pnt p = BRep_Tool::Pnt( TopoDS::Vertex( V ) );
 988     for ( exp1.Init( S, TopAbs_SOLID ); exp1.More(); exp1.Next() ) {
 989       BRepClass3d_SolidClassifier classifier( exp1.Current(), p, 1e-6);
 990       if ( classifier.State() == TopAbs_IN ) {
 991         thePnt1 = p;
 992         thePnt2 = p;
 993         return 0.0;
 994       }
 995     }
 996   }
 997   // End Issue 0020231
 998
 999   // skl 30.06.2008
1000   // additional workaround for bugs 19899, 19908 and 19910 from Mantis
1001   double dist = GEOMUtils::GetMinDistanceSingular
1002       (theShape1, theShape2, thePnt1, thePnt2);
1003
1004   if (dist > -1.0) {
1005     return dist;
1006   }
1007
1008   BRepExtrema_DistShapeShape dst (theShape1, theShape2);
1009   if (dst.IsDone()) {
1010     gp_Pnt P1, P2;
1011
1012     for (int i = 1; i <= dst.NbSolution(); i++) {
1013       P1 = dst.PointOnShape1(i);
1014       P2 = dst.PointOnShape2(i);
1015
1016       Standard_Real Dist = P1.Distance(P2);
1017       if (aResult > Dist) {
1018         aResult = Dist;
1019         thePnt1 = P1;
1020         thePnt2 = P2;
1021       }
1022     }
1023   }
1024
1025   return aResult;
1026 }
1027
1028 //=======================================================================
1029 // function : ConvertClickToPoint()
1030 // purpose  : Returns the point clicked in 3D view
1031 //=======================================================================
1032 gp_Pnt GEOMUtils::ConvertClickToPoint( int x, int y, Handle(V3d_View) aView )
1033 {
1034   V3d_Coordinate XEye, YEye, ZEye, XAt, YAt, ZAt;
1035   aView->Eye( XEye, YEye, ZEye );
1036
1037   aView->At( XAt, YAt, ZAt );
1038   gp_Pnt EyePoint( XEye, YEye, ZEye );
1039   gp_Pnt AtPoint( XAt, YAt, ZAt );
1040
1041   gp_Vec EyeVector( EyePoint, AtPoint );
1042   gp_Dir EyeDir( EyeVector );
1043
1044   gp_Pln PlaneOfTheView = gp_Pln( AtPoint, EyeDir );
1045   Standard_Real X, Y, Z;
1046   aView->Convert( x, y, X, Y, Z );
1047   gp_Pnt ConvertedPoint( X, Y, Z );
1048
1049   gp_Pnt2d ConvertedPointOnPlane = ProjLib::Project( PlaneOfTheView, ConvertedPoint );
1050   gp_Pnt ResultPoint = ElSLib::Value( ConvertedPointOnPlane.X(), ConvertedPointOnPlane.Y(), PlaneOfTheView );
1051   return ResultPoint;
1052 }
1053
1054 //=======================================================================
1055 // function : ConvertTreeToString()
1056 // purpose  : Returns the string representation of dependency tree
1057 //=======================================================================
1058 void GEOMUtils::ConvertTreeToString( const TreeModel& tree,
1059                                      std::string& dependencyStr )
1060 {
1061   TreeModel::const_iterator i;
1062   for ( i = tree.begin(); i != tree.end(); ++i ) {
1063     dependencyStr.append( i->first );
1064     dependencyStr.append( "-" );
1065     std::vector<LevelInfo> upLevelList = i->second.first;
1066     dependencyStr.append( "upward" );
1067     parseWard( upLevelList, dependencyStr );
1068     std::vector<LevelInfo> downLevelList = i->second.second;
1069     dependencyStr.append( "downward" );
1070     parseWard( downLevelList, dependencyStr );
1071   }
1072 }
1073
1074 //=======================================================================
1075 // function : ConvertStringToTree()
1076 // purpose  : Returns the dependency tree
1077 //=======================================================================
1078 void GEOMUtils::ConvertStringToTree( const std::string& dependencyStr,
1079                                      TreeModel& tree )
1080 {
1081   std::size_t cursor = 0;
1082
1083   while( dependencyStr.find('-',cursor) != std::string::npos ) //find next selected object
1084   {
1085     std::size_t objectIndex = dependencyStr.find( '-', cursor );
1086     std::string objectEntry = dependencyStr.substr( cursor, objectIndex - cursor );
1087     cursor = objectIndex;
1088
1089     std::size_t upwardIndexBegin = dependencyStr.find("{",cursor) + 1;
1090     //std::size_t upwardIndexFinish = dependencyStr.find("}",upwardIndexBegin);
1091     LevelsList upwardList = parseWard( dependencyStr, cursor );
1092
1093     LevelsList downwardList = parseWard( dependencyStr, cursor );
1094
1095     tree[objectEntry] = std::pair<LevelsList,LevelsList>( upwardList, downwardList );
1096   }
1097 }
1098
1099 bool GEOMUtils::CheckShape( TopoDS_Shape& shape,
1100                             bool checkGeometry )
1101 {
1102   BRepCheck_Analyzer analyzer( shape, checkGeometry );
1103   return analyzer.IsValid();
1104 }
1105
1106 bool GEOMUtils::CheckBOPArguments(const TopoDS_Shape &theShape)
1107 {
1108   BOPAlgo_ArgumentAnalyzer aChecker;
1109
1110   aChecker.SetShape1(theShape);
1111
1112   // Set default options
1113   aChecker.ArgumentTypeMode()   = Standard_True;
1114   aChecker.SelfInterMode()      = Standard_True;
1115   aChecker.SmallEdgeMode()      = Standard_True;
1116   aChecker.RebuildFaceMode()    = Standard_True;
1117   aChecker.ContinuityMode()     = Standard_True;
1118   aChecker.CurveOnSurfaceMode() = Standard_True;
1119
1120   aChecker.StopOnFirstFaulty() = Standard_True;
1121   aChecker.Perform();
1122
1123   // process result of checking
1124   const bool isValid = !aChecker.HasFaulty();
1125
1126   return isValid;
1127 }
1128
1129 bool GEOMUtils::FixShapeTolerance( TopoDS_Shape& shape,
1130                                    TopAbs_ShapeEnum type,
1131                                    Standard_Real tolerance,
1132                                    bool checkGeometry )
1133 {
1134   ShapeFix_ShapeTolerance aSft;
1135   aSft.LimitTolerance( shape, tolerance, tolerance, type );
1136   Handle(ShapeFix_Shape) aSfs = new ShapeFix_Shape( shape );
1137   aSfs->Perform();
1138   shape = aSfs->Shape();
1139   return CheckShape( shape, checkGeometry );
1140 }
1141
1142 bool GEOMUtils::FixShapeTolerance( TopoDS_Shape& shape,
1143                                    Standard_Real tolerance,
1144                                    bool checkGeometry )
1145 {
1146   return FixShapeTolerance( shape, TopAbs_SHAPE, tolerance, checkGeometry );
1147 }
1148
1149 bool GEOMUtils::FixShapeTolerance( TopoDS_Shape& shape,
1150                                    bool checkGeometry )
1151 {
1152   return FixShapeTolerance( shape, Precision::Confusion(), checkGeometry );
1153 }
1154
1155 bool GEOMUtils::FixShapeCurves( TopoDS_Shape& shape )
1156 {
1157   Standard_Real aT, aTolE, aD, aDMax;
1158   TopExp_Explorer aExpF, aExpE;
1159   NCollection_DataMap<TopoDS_Edge, Standard_Real, TopTools_ShapeMapHasher> aDMETol;
1160   aExpF.Init(shape, TopAbs_FACE);
1161   for (; aExpF.More(); aExpF.Next()) {
1162     const TopoDS_Face& aF = *(TopoDS_Face*)&aExpF.Current();
1163     aExpE.Init(aF, TopAbs_EDGE);
1164     for (; aExpE.More(); aExpE.Next()) {
1165       const TopoDS_Edge& aE = *(TopoDS_Edge*)&aExpE.Current();
1166       try {
1167         if (!BOPTools_AlgoTools::ComputeTolerance(aF, aE, aDMax, aT)) {
1168           continue;
1169         }
1170       }
1171       catch(...) {
1172         continue;
1173       }
1174       aTolE = BRep_Tool::Tolerance(aE);
1175       if (aDMax < aTolE) continue;
1176       if (aDMETol.IsBound(aE)) {
1177         aD = aDMETol.Find(aE);
1178         if (aDMax > aD) {
1179           aDMETol.UnBind(aE);
1180           aDMETol.Bind(aE, aDMax);
1181         }
1182       }
1183       else {
1184         aDMETol.Bind(aE, aDMax);
1185       }
1186     }
1187   }
1188   NCollection_DataMap<TopoDS_Edge, Standard_Real, TopTools_ShapeMapHasher>::Iterator aDMETolIt(aDMETol);
1189 #ifdef USE_LIMIT_TOLERANCE
1190   ShapeFix_ShapeTolerance sat;
1191 #else
1192   BRep_Builder b;
1193 #endif
1194   for (; aDMETolIt.More(); aDMETolIt.Next()) {
1195 #ifdef USE_LIMIT_TOLERANCE
1196     sat.LimitTolerance(aDMETolIt.Key(), aDMETolIt.Value()*1.001);
1197 #else
1198     TopoDS_Iterator itv(aDMETolIt.Key());
1199     for (; itv.More(); itv.Next())
1200       b.UpdateVertex(TopoDS::Vertex(itv.Value()), aDMETolIt.Value()*1.001);
1201     b.UpdateEdge(aDMETolIt.Key(), aDMETolIt.Value()*1.001);
1202 #endif
1203   }
1204   return CheckShape( shape );
1205 }
1206
1207 bool GEOMUtils::Write( const TopoDS_Shape& shape, const char* fileName )
1208 {
1209   return BRepTools::Write( shape, fileName );
1210 }
1211
1212 TopoDS_Shape GEOMUtils::ReduceCompound( const TopoDS_Shape& shape )
1213 {
1214   TopoDS_Shape result = shape;
1215
1216   if ( shape.ShapeType() == TopAbs_COMPOUND ||
1217        shape.ShapeType() == TopAbs_COMPSOLID ) {
1218
1219     TopTools_ListOfShape l;
1220
1221     TopoDS_Iterator it ( shape );
1222     for ( ; it.More(); it.Next() )
1223       l.Append( it.Value() );
1224     if ( l.Extent() == 1 && l.First() != shape )
1225       result = ReduceCompound( l.First() );
1226   }
1227
1228   return result;
1229 }
1230
1231 void GEOMUtils::MeshShape( const TopoDS_Shape shape,
1232                            double deflection, bool theForced )
1233 {
1234   Standard_Real aDeflection = ( deflection <= 0 ) ? DefaultDeflection() : deflection;
1235
1236   // Is shape triangulated?
1237   Standard_Boolean alreadyMeshed = true;
1238   TopExp_Explorer ex;
1239   TopLoc_Location aLoc;
1240   for ( ex.Init( shape, TopAbs_FACE ); ex.More() && alreadyMeshed; ex.Next() ) {
1241     const TopoDS_Face& aFace = TopoDS::Face( ex.Current() );
1242     Handle(Poly_Triangulation) aPoly = BRep_Tool::Triangulation( aFace, aLoc );
1243     alreadyMeshed = !aPoly.IsNull();
1244   }
1245
1246   if ( !alreadyMeshed || theForced ) {
1247     // Compute bounding box
1248     Bnd_Box B;
1249     BRepBndLib::Add( shape, B );
1250     if ( B.IsVoid() )
1251       return; // NPAL15983 (Bug when displaying empty groups)
1252     Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
1253     B.Get( aXmin, aYmin, aZmin, aXmax, aYmax, aZmax );
1254
1255     // This magic line comes from Prs3d_ShadedShape.gxx in OCCT
1256     aDeflection = MAX3(aXmax-aXmin, aYmax-aYmin, aZmax-aZmin) * aDeflection * 4;
1257
1258     // Clean triangulation before compute incremental mesh
1259     BRepTools::Clean( shape );
1260
1261     // Compute triangulation
1262     BRepMesh_IncrementalMesh mesh( shape, aDeflection );
1263   }
1264 }
1265
1266 double GEOMUtils::DefaultDeflection()
1267 {
1268   return 0.001;
1269 }
1270
1271 //=======================================================================
1272 //function : IsOpenPath
1273 //purpose  :
1274 //=======================================================================
1275 bool GEOMUtils::IsOpenPath(const TopoDS_Shape &theShape)
1276 {
1277   bool isOpen = true;
1278
1279   if (theShape.IsNull() == Standard_False) {
1280     if (theShape.Closed()) {
1281       // The shape is closed
1282       isOpen = false;
1283     } else {
1284       const TopAbs_ShapeEnum aType = theShape.ShapeType();
1285
1286       if (aType == TopAbs_EDGE || aType == TopAbs_WIRE) {
1287         // Check if path ends are coinsident.
1288         TopoDS_Vertex aV[2];
1289
1290         if (aType == TopAbs_EDGE) {
1291           // Edge
1292           TopExp::Vertices(TopoDS::Edge(theShape), aV[0], aV[1]);
1293         } else {
1294           // Wire
1295           TopExp::Vertices(TopoDS::Wire(theShape), aV[0], aV[1]);
1296         }
1297
1298         if (aV[0].IsNull() == Standard_False &&
1299             aV[1].IsNull() == Standard_False) {
1300           if (aV[0].IsSame(aV[1])) {
1301             // The shape is closed
1302             isOpen = false;
1303           } else {
1304             const Standard_Real aTol1 = BRep_Tool::Tolerance(aV[0]);
1305             const Standard_Real aTol2 = BRep_Tool::Tolerance(aV[1]);
1306             const gp_Pnt        aPnt1 = BRep_Tool::Pnt(aV[0]);
1307             const gp_Pnt        aPnt2 = BRep_Tool::Pnt(aV[1]);
1308
1309             if (aPnt1.Distance(aPnt2) <= aTol1 + aTol2) {
1310               // The shape is closed
1311               isOpen = false;
1312             }
1313           }
1314         }
1315       }
1316     }
1317   }
1318
1319   return isOpen;
1320 }
1321
1322 //=======================================================================
1323 //function : CompareToleranceValues
1324 //purpose  :
1325 //=======================================================================
1326 int GEOMUtils::CompareToleranceValues(const double theTolShape,
1327                                       const double theTolRef)
1328 {
1329   const double aTolTol = Min(DEFAULT_MAX_TOLERANCE_TOLERANCE,
1330                              theTolRef*DEFAULT_TOLERANCE_TOLERANCE);
1331
1332   int aResult = 0;
1333
1334   if (theTolShape < theTolRef - aTolTol) {
1335     aResult = -1;
1336   } else if (theTolShape > theTolRef + aTolTol) {
1337     aResult = 1;
1338   }
1339
1340   return aResult;
1341 }
1342
1343 //=======================================================================
1344 //function : IsFitCondition
1345 //purpose  :
1346 //=======================================================================
1347 bool GEOMUtils::IsFitCondition(const ComparisonCondition theCondition,
1348                                const double              theTolShape,
1349                                const double              theTolRef)
1350 {
1351   const int aCompValue = CompareToleranceValues(theTolShape, theTolRef);
1352   bool      isFit      = false;
1353
1354   switch (theCondition) {
1355     case CC_GT:
1356       isFit = aCompValue == 1;
1357       break;
1358     case GEOMUtils::CC_GE:
1359       isFit = aCompValue != -1;
1360       break;
1361     case GEOMUtils::CC_LT:
1362       isFit = aCompValue == -1;
1363       break;
1364     case GEOMUtils::CC_LE:
1365       isFit = aCompValue != 1;
1366       break;
1367     default:
1368       break;
1369   }
1370
1371   return isFit;
1372 }