src/GEOMUtils/GEOMUtils.cxx

   1 // Copyright (C) 2007-2013  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
  23 #include <Standard_Stream.hxx>
  24
  25 #include <GEOMUtils.hxx>
  26
  27 #include <Basics_OCCTVersion.hxx>
  28
  29 #include <utilities.h>
  30 #include <OpUtil.hxx>
  31 #include <Utils_ExceptHandlers.hxx>
  32
  33 // OCCT Includes
  34 #include <BRepMesh_IncrementalMesh.hxx>
  35
  36 #include <BRepExtrema_DistShapeShape.hxx>
  37
  38 #include <BRep_Builder.hxx>
  39 #include <BRep_Tool.hxx>
  40 #include <BRepBndLib.hxx>
  41 #include <BRepGProp.hxx>
  42 #include <BRepTools.hxx>
  43
  44 #include <Bnd_Box.hxx>
  45
  46 #include <TopAbs.hxx>
  47 #include <TopExp.hxx>
  48 #include <TopoDS.hxx>
  49 #include <TopoDS_Edge.hxx>
  50 #include <TopoDS_Face.hxx>
  51 #include <TopoDS_Shape.hxx>
  52 #include <TopoDS_Vertex.hxx>
  53 #include <TopoDS_Compound.hxx>
  54 #include <TopoDS_Iterator.hxx>
  55 #include <TopExp_Explorer.hxx>
  56 #include <TopTools_MapOfShape.hxx>
  57 #include <TopTools_ListOfShape.hxx>
  58 #include <TopTools_ListIteratorOfListOfShape.hxx>
  59 #include <TopTools_Array1OfShape.hxx>
  60
  61 #include <Geom_Surface.hxx>
  62 #include <Geom_Plane.hxx>
  63
  64 #include <GeomLProp_CLProps.hxx>
  65 #include <GeomLProp_SLProps.hxx>
  66
  67 #include <GProp_GProps.hxx>
  68 #include <GProp_PrincipalProps.hxx>
  69
  70 #include <TColStd_Array1OfReal.hxx>
  71
  72 #include <gp_Pln.hxx>
  73 #include <gp_Lin.hxx>
  74
  75 #include <vector>
  76
  77 #include <Standard_Failure.hxx>
  78 #include <Standard_NullObject.hxx>
  79 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
  80
  81 #define STD_SORT_ALGO 1
  82
  83 //=======================================================================
  84 //function : GetPosition
  85 //purpose  :
  86 //=======================================================================
  87 gp_Ax3 GEOMUtils::GetPosition (const TopoDS_Shape& theShape)
  88 {
  89   gp_Ax3 aResult;
  90
  91   if (theShape.IsNull())
  92     return aResult;
  93
  94   // Axes
  95   aResult.Transform(theShape.Location().Transformation());
  96   if (theShape.ShapeType() == TopAbs_FACE) {
  97     Handle(Geom_Surface) aGS = BRep_Tool::Surface(TopoDS::Face(theShape));
  98     if (!aGS.IsNull() && aGS->IsKind(STANDARD_TYPE(Geom_Plane))) {
  99       Handle(Geom_Plane) aGPlane = Handle(Geom_Plane)::DownCast(aGS);
 100       gp_Pln aPln = aGPlane->Pln();
 101       aResult = aPln.Position();
 102       // In case of reverse orinetation of the face invert the plane normal
 103       // (the face's normal does not mathc the plane's normal in this case)
 104       if(theShape.Orientation() == TopAbs_REVERSED)
 105       {
 106         gp_Dir Vx =  aResult.XDirection();
 107         gp_Dir N  =  aResult.Direction().Mirrored(Vx);
 108         gp_Pnt P  =  aResult.Location();
 109         aResult = gp_Ax3(P, N, Vx);
 110       }
 111     }
 112   }
 113
 114   // Origin
 115   gp_Pnt aPnt;
 116
 117   TopAbs_ShapeEnum aShType = theShape.ShapeType();
 118
 119   if (aShType == TopAbs_VERTEX) {
 120     aPnt = BRep_Tool::Pnt(TopoDS::Vertex(theShape));
 121   }
 122   else {
 123     if (aShType == TopAbs_COMPOUND) {
 124       aShType = GetTypeOfSimplePart(theShape);
 125     }
 126
 127     GProp_GProps aSystem;
 128     if (aShType == TopAbs_EDGE || aShType == TopAbs_WIRE)
 129       BRepGProp::LinearProperties(theShape, aSystem);
 130     else if (aShType == TopAbs_FACE || aShType == TopAbs_SHELL)
 131       BRepGProp::SurfaceProperties(theShape, aSystem);
 132     else
 133       BRepGProp::VolumeProperties(theShape, aSystem);
 134
 135     aPnt = aSystem.CentreOfMass();
 136   }
 137
 138   aResult.SetLocation(aPnt);
 139
 140   return aResult;
 141 }
 142
 143 //=======================================================================
 144 //function : GetVector
 145 //purpose  :
 146 //=======================================================================
 147 gp_Vec GEOMUtils::GetVector (const TopoDS_Shape& theShape,
 148                              Standard_Boolean doConsiderOrientation)
 149 {
 150   if (theShape.IsNull())
 151     Standard_NullObject::Raise("Null shape is given for a vector");
 152
 153   if (theShape.ShapeType() != TopAbs_EDGE)
 154     Standard_TypeMismatch::Raise("Invalid shape is given, must be a vector or an edge");
 155
 156   TopoDS_Edge anE = TopoDS::Edge(theShape);
 157
 158   TopoDS_Vertex V1, V2;
 159   TopExp::Vertices(anE, V1, V2, doConsiderOrientation);
 160
 161   if (V1.IsNull() || V2.IsNull())
 162     Standard_NullObject::Raise("Invalid edge is given, it must have two points");
 163
 164   gp_Vec aV (BRep_Tool::Pnt(V1), BRep_Tool::Pnt(V2));
 165   if (aV.Magnitude() < gp::Resolution()) {
 166     Standard_ConstructionError::Raise("Vector of zero length is given");
 167   }
 168
 169   return aV;
 170 }
 171
 172 //=======================================================================
 173 //function : ShapeToDouble
 174 //purpose  : used by CompareShapes::operator()
 175 //=======================================================================
 176 std::pair<double, double> ShapeToDouble (const TopoDS_Shape& S, bool isOldSorting)
 177 {
 178   // Computing of CentreOfMass
 179   gp_Pnt GPoint;
 180   double Len;
 181
 182   if (S.ShapeType() == TopAbs_VERTEX) {
 183     GPoint = BRep_Tool::Pnt(TopoDS::Vertex(S));
 184     Len = (double)S.Orientation();
 185   }
 186   else {
 187     GProp_GProps GPr;
 188     // BEGIN: fix for Mantis issue 0020842
 189     if (isOldSorting) {
 190       BRepGProp::LinearProperties(S, GPr);
 191     }
 192     else {
 193       if (S.ShapeType() == TopAbs_EDGE || S.ShapeType() == TopAbs_WIRE) {
 194         BRepGProp::LinearProperties(S, GPr);
 195       }
 196       else if (S.ShapeType() == TopAbs_FACE || S.ShapeType() == TopAbs_SHELL) {
 197         BRepGProp::SurfaceProperties(S, GPr);
 198       }
 199       else {
 200         BRepGProp::VolumeProperties(S, GPr);
 201       }
 202     }
 203     // END: fix for Mantis issue 0020842
 204     GPoint = GPr.CentreOfMass();
 205     Len = GPr.Mass();
 206   }
 207
 208   double dMidXYZ = GPoint.X() * 999.0 + GPoint.Y() * 99.0 + GPoint.Z() * 0.9;
 209   return std::make_pair(dMidXYZ, Len);
 210 }
 211
 212 //=======================================================================
 213 //function : CompareShapes::operator()
 214 //purpose  : used by std::sort(), called from SortShapes()
 215 //=======================================================================
 216 bool GEOMUtils::CompareShapes::operator() (const TopoDS_Shape& theShape1,
 217                                            const TopoDS_Shape& theShape2)
 218 {
 219   if (!myMap.IsBound(theShape1)) {
 220     myMap.Bind(theShape1, ShapeToDouble(theShape1, myIsOldSorting));
 221   }
 222
 223   if (!myMap.IsBound(theShape2)) {
 224     myMap.Bind(theShape2, ShapeToDouble(theShape2, myIsOldSorting));
 225   }
 226
 227   std::pair<double, double> val1 = myMap.Find(theShape1);
 228   std::pair<double, double> val2 = myMap.Find(theShape2);
 229
 230   double tol = Precision::Confusion();
 231   bool exchange = Standard_False;
 232
 233   double dMidXYZ = val1.first - val2.first;
 234   if (dMidXYZ >= tol) {
 235     exchange = Standard_True;
 236   }
 237   else if (Abs(dMidXYZ) < tol) {
 238     double dLength = val1.second - val2.second;
 239     if (dLength >= tol) {
 240       exchange = Standard_True;
 241     }
 242     else if (Abs(dLength) < tol && theShape1.ShapeType() <= TopAbs_FACE) {
 243       // PAL17233
 244       // equal values possible on shapes such as two halves of a sphere and
 245       // a membrane inside the sphere
 246       Bnd_Box box1,box2;
 247       BRepBndLib::Add(theShape1, box1);
 248       if (!box1.IsVoid()) {
 249         BRepBndLib::Add(theShape2, box2);
 250         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
 251         if (dSquareExtent >= tol) {
 252           exchange = Standard_True;
 253         }
 254         else if (Abs(dSquareExtent) < tol) {
 255           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
 256           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 257           val1 = (aXmin+aXmax)*999.0 + (aYmin+aYmax)*99.0 + (aZmin+aZmax)*0.9;
 258           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 259           val2 = (aXmin+aXmax)*999.0 + (aYmin+aYmax)*99.0 + (aZmin+aZmax)*0.9;
 260           if ((val1 - val2) >= tol) {
 261             exchange = Standard_True;
 262           }
 263         }
 264       }
 265     }
 266   }
 267
 268   //return val1 < val2;
 269   return !exchange;
 270 }
 271
 272 //=======================================================================
 273 //function : SortShapes
 274 //purpose  :
 275 //=======================================================================
 276 void GEOMUtils::SortShapes (TopTools_ListOfShape& SL,
 277                             const Standard_Boolean isOldSorting)
 278 {
 279 #ifdef STD_SORT_ALGO
 280   std::vector<TopoDS_Shape> aShapesVec;
 281   aShapesVec.reserve(SL.Extent());
 282
 283   TopTools_ListIteratorOfListOfShape it (SL);
 284   for (; it.More(); it.Next()) {
 285     aShapesVec.push_back(it.Value());
 286   }
 287   SL.Clear();
 288
 289   CompareShapes shComp (isOldSorting);
 290   std::stable_sort(aShapesVec.begin(), aShapesVec.end(), shComp);
 291   //std::sort(aShapesVec.begin(), aShapesVec.end(), shComp);
 292
 293   std::vector<TopoDS_Shape>::const_iterator anIter = aShapesVec.begin();
 294   for (; anIter != aShapesVec.end(); ++anIter) {
 295     SL.Append(*anIter);
 296   }
 297 #else
 298   // old implementation
 299   Standard_Integer MaxShapes = SL.Extent();
 300   TopTools_Array1OfShape  aShapes (1,MaxShapes);
 301   TColStd_Array1OfInteger OrderInd(1,MaxShapes);
 302   TColStd_Array1OfReal    MidXYZ  (1,MaxShapes); //X,Y,Z;
 303   TColStd_Array1OfReal    Length  (1,MaxShapes); //X,Y,Z;
 304
 305   // Computing of CentreOfMass
 306   Standard_Integer Index;
 307   GProp_GProps GPr;
 308   gp_Pnt GPoint;
 309   TopTools_ListIteratorOfListOfShape it(SL);
 310   for (Index=1;  it.More();  Index++)
 311   {
 312     TopoDS_Shape S = it.Value();
 313     SL.Remove( it ); // == it.Next()
 314     aShapes(Index) = S;
 315     OrderInd.SetValue (Index, Index);
 316     if (S.ShapeType() == TopAbs_VERTEX) {
 317       GPoint = BRep_Tool::Pnt( TopoDS::Vertex( S ));
 318       Length.SetValue( Index, (Standard_Real) S.Orientation());
 319     }
 320     else {
 321       // BEGIN: fix for Mantis issue 0020842
 322       if (isOldSorting) {
 323         BRepGProp::LinearProperties (S, GPr);
 324       }
 325       else {
 326         if (S.ShapeType() == TopAbs_EDGE || S.ShapeType() == TopAbs_WIRE) {
 327           BRepGProp::LinearProperties (S, GPr);
 328         }
 329         else if (S.ShapeType() == TopAbs_FACE || S.ShapeType() == TopAbs_SHELL) {
 330           BRepGProp::SurfaceProperties(S, GPr);
 331         }
 332         else {
 333           BRepGProp::VolumeProperties(S, GPr);
 334         }
 335       }
 336       // END: fix for Mantis issue 0020842
 337       GPoint = GPr.CentreOfMass();
 338       Length.SetValue(Index, GPr.Mass());
 339     }
 340     MidXYZ.SetValue(Index, GPoint.X()*999.0 + GPoint.Y()*99.0 + GPoint.Z()*0.9);
 341     //cout << Index << " L: " << Length(Index) << "CG: " << MidXYZ(Index) << endl;
 342   }
 343
 344   // Sorting
 345   Standard_Integer aTemp;
 346   Standard_Boolean exchange, Sort = Standard_True;
 347   Standard_Real    tol = Precision::Confusion();
 348   while (Sort)
 349   {
 350     Sort = Standard_False;
 351     for (Index=1; Index < MaxShapes; Index++)
 352     {
 353       exchange = Standard_False;
 354       Standard_Real dMidXYZ = MidXYZ(OrderInd(Index)) - MidXYZ(OrderInd(Index+1));
 355       Standard_Real dLength = Length(OrderInd(Index)) - Length(OrderInd(Index+1));
 356       if ( dMidXYZ >= tol ) {
 357 //         cout << "MidXYZ: " << MidXYZ(OrderInd(Index))<< " > " <<MidXYZ(OrderInd(Index+1))
 358 //              << " d: " << dMidXYZ << endl;
 359         exchange = Standard_True;
 360       }
 361       else if ( Abs(dMidXYZ) < tol && dLength >= tol ) {
 362 //         cout << "Length: " << Length(OrderInd(Index))<< " > " <<Length(OrderInd(Index+1))
 363 //              << " d: " << dLength << endl;
 364         exchange = Standard_True;
 365       }
 366       else if ( Abs(dMidXYZ) < tol && Abs(dLength) < tol &&
 367                 aShapes(OrderInd(Index)).ShapeType() <= TopAbs_FACE) {
 368         // PAL17233
 369         // equal values possible on shapes such as two halves of a sphere and
 370         // a membrane inside the sphere
 371         Bnd_Box box1,box2;
 372         BRepBndLib::Add( aShapes( OrderInd(Index) ), box1 );
 373         if ( box1.IsVoid() ) continue;
 374         BRepBndLib::Add( aShapes( OrderInd(Index+1) ), box2 );
 375         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
 376         if ( dSquareExtent >= tol ) {
 377 //           cout << "SquareExtent: " << box1.SquareExtent()<<" > "<<box2.SquareExtent() << endl;
 378           exchange = Standard_True;
 379         }
 380         else if ( Abs(dSquareExtent) < tol ) {
 381           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
 382           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 383           val1 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
 384           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 385           val2 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
 386           //exchange = val1 > val2;
 387           if ((val1 - val2) >= tol) {
 388             exchange = Standard_True;
 389           }
 390           //cout << "box: " << val1<<" > "<<val2 << endl;
 391         }
 392       }
 393
 394       if (exchange)
 395       {
 396 //         cout << "exchange " << Index << " & " << Index+1 << endl;
 397         aTemp = OrderInd(Index);
 398         OrderInd(Index) = OrderInd(Index+1);
 399         OrderInd(Index+1) = aTemp;
 400         Sort = Standard_True;
 401       }
 402     }
 403   }
 404
 405   for (Index=1; Index <= MaxShapes; Index++)
 406     SL.Append( aShapes( OrderInd(Index) ));
 407 #endif
 408 }
 409
 410 //=======================================================================
 411 //function : CompsolidToCompound
 412 //purpose  :
 413 //=======================================================================
 414 TopoDS_Shape GEOMUtils::CompsolidToCompound (const TopoDS_Shape& theCompsolid)
 415 {
 416   if (theCompsolid.ShapeType() != TopAbs_COMPSOLID) {
 417     return theCompsolid;
 418   }
 419
 420   TopoDS_Compound aCompound;
 421   BRep_Builder B;
 422   B.MakeCompound(aCompound);
 423
 424   TopTools_MapOfShape mapShape;
 425   TopoDS_Iterator It (theCompsolid, Standard_True, Standard_True);
 426
 427   for (; It.More(); It.Next()) {
 428     TopoDS_Shape aShape_i = It.Value();
 429     if (mapShape.Add(aShape_i)) {
 430       B.Add(aCompound, aShape_i);
 431     }
 432   }
 433
 434   return aCompound;
 435 }
 436
 437 //=======================================================================
 438 //function : AddSimpleShapes
 439 //purpose  :
 440 //=======================================================================
 441 void GEOMUtils::AddSimpleShapes (const TopoDS_Shape& theShape, TopTools_ListOfShape& theList)
 442 {
 443   if (theShape.ShapeType() != TopAbs_COMPOUND &&
 444       theShape.ShapeType() != TopAbs_COMPSOLID) {
 445     theList.Append(theShape);
 446     return;
 447   }
 448
 449   TopTools_MapOfShape mapShape;
 450   TopoDS_Iterator It (theShape, Standard_True, Standard_True);
 451
 452   for (; It.More(); It.Next()) {
 453     TopoDS_Shape aShape_i = It.Value();
 454     if (mapShape.Add(aShape_i)) {
 455       if (aShape_i.ShapeType() == TopAbs_COMPOUND ||
 456           aShape_i.ShapeType() == TopAbs_COMPSOLID) {
 457         AddSimpleShapes(aShape_i, theList);
 458       } else {
 459         theList.Append(aShape_i);
 460       }
 461     }
 462   }
 463 }
 464
 465 //=======================================================================
 466 //function : CheckTriangulation
 467 //purpose  :
 468 //=======================================================================
 469 bool GEOMUtils::CheckTriangulation (const TopoDS_Shape& aShape)
 470 {
 471   bool isTriangulation = true;
 472
 473   TopExp_Explorer exp (aShape, TopAbs_FACE);
 474   if (exp.More())
 475   {
 476     TopLoc_Location aTopLoc;
 477     Handle(Poly_Triangulation) aTRF;
 478     aTRF = BRep_Tool::Triangulation(TopoDS::Face(exp.Current()), aTopLoc);
 479     if (aTRF.IsNull()) {
 480       isTriangulation = false;
 481     }
 482   }
 483   else // no faces, try edges
 484   {
 485     TopExp_Explorer expe (aShape, TopAbs_EDGE);
 486     if (!expe.More()) {
 487       return false;
 488     }
 489     TopLoc_Location aLoc;
 490     Handle(Poly_Polygon3D) aPE = BRep_Tool::Polygon3D(TopoDS::Edge(expe.Current()), aLoc);
 491     if (aPE.IsNull()) {
 492       isTriangulation = false;
 493     }
 494   }
 495
 496   if (!isTriangulation) {
 497     // calculate deflection
 498     Standard_Real aDeviationCoefficient = 0.001;
 499
 500     Bnd_Box B;
 501     BRepBndLib::Add(aShape, B);
 502     Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
 503     B.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
 504
 505     Standard_Real dx = aXmax - aXmin, dy = aYmax - aYmin, dz = aZmax - aZmin;
 506     Standard_Real aDeflection = Max(Max(dx, dy), dz) * aDeviationCoefficient * 4;
 507     Standard_Real aHLRAngle = 0.349066;
 508
 509     BRepMesh_IncrementalMesh Inc (aShape, aDeflection, Standard_False, aHLRAngle);
 510   }
 511
 512   return true;
 513 }
 514
 515 //=======================================================================
 516 //function : GetTypeOfSimplePart
 517 //purpose  :
 518 //=======================================================================
 519 TopAbs_ShapeEnum GEOMUtils::GetTypeOfSimplePart (const TopoDS_Shape& theShape)
 520 {
 521   TopAbs_ShapeEnum aType = theShape.ShapeType();
 522   if      (aType == TopAbs_VERTEX)                             return TopAbs_VERTEX;
 523   else if (aType == TopAbs_EDGE  || aType == TopAbs_WIRE)      return TopAbs_EDGE;
 524   else if (aType == TopAbs_FACE  || aType == TopAbs_SHELL)     return TopAbs_FACE;
 525   else if (aType == TopAbs_SOLID || aType == TopAbs_COMPSOLID) return TopAbs_SOLID;
 526   else if (aType == TopAbs_COMPOUND) {
 527     // Only the iType of the first shape in the compound is taken into account
 528     TopoDS_Iterator It (theShape, Standard_False, Standard_False);
 529     if (It.More()) {
 530       return GetTypeOfSimplePart(It.Value());
 531     }
 532   }
 533   return TopAbs_SHAPE;
 534 }
 535
 536 //=======================================================================
 537 //function : GetEdgeNearPoint
 538 //purpose  :
 539 //=======================================================================
 540 TopoDS_Shape GEOMUtils::GetEdgeNearPoint (const TopoDS_Shape& theShape,
 541                                           const TopoDS_Vertex& thePoint)
 542 {
 543   TopoDS_Shape aResult;
 544
 545   // 1. Explode the shape on edges
 546   TopTools_MapOfShape mapShape;
 547   Standard_Integer nbEdges = 0;
 548   TopExp_Explorer exp (theShape, TopAbs_EDGE);
 549   for (; exp.More(); exp.Next()) {
 550     if (mapShape.Add(exp.Current())) {
 551       nbEdges++;
 552     }
 553   }
 554
 555   if (nbEdges == 0)
 556     Standard_NullObject::Raise("Given shape contains no edges");
 557
 558   mapShape.Clear();
 559   Standard_Integer ind = 1;
 560   TopTools_Array1OfShape anEdges (1, nbEdges);
 561   TColStd_Array1OfReal aDistances (1, nbEdges);
 562   for (exp.Init(theShape, TopAbs_EDGE); exp.More(); exp.Next()) {
 563     if (mapShape.Add(exp.Current())) {
 564       TopoDS_Shape anEdge = exp.Current();
 565       anEdges(ind) = anEdge;
 566
 567       // 2. Classify the point relatively each edge
 568       BRepExtrema_DistShapeShape aDistTool (thePoint, anEdges(ind));
 569       if (!aDistTool.IsDone())
 570         Standard_ConstructionError::Raise("Cannot find a distance from the given point to one of edges");
 571
 572       aDistances(ind) = aDistTool.Value();
 573       ind++;
 574     }
 575   }
 576
 577   // 3. Define edge, having minimum distance to the point
 578   Standard_Real nearest = RealLast(), nbFound = 0;
 579   Standard_Real prec = Precision::Confusion();
 580   for (ind = 1; ind <= nbEdges; ind++) {
 581     if (Abs(aDistances(ind) - nearest) < prec) {
 582       nbFound++;
 583     }
 584     else if (aDistances(ind) < nearest) {
 585       nearest = aDistances(ind);
 586       aResult = anEdges(ind);
 587       nbFound = 1;
 588     }
 589     else {
 590     }
 591   }
 592   if (nbFound > 1) {
 593     Standard_ConstructionError::Raise("Multiple edges near the given point are found");
 594   }
 595   else if (nbFound == 0) {
 596     Standard_ConstructionError::Raise("There are no edges near the given point");
 597   }
 598   else {
 599   }
 600
 601   return aResult;
 602 }