src/GEOMImpl/GEOMImpl_IMeasureOperations.cxx

   1 // Copyright (C) 2007-2023  CEA, EDF, 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 <GEOMImpl_IMeasureOperations.hxx>
  24 #include <GEOMImpl_IMeasure.hxx>
  25 #include <GEOMImpl_MeasureDriver.hxx>
  26
  27 #include <GEOMImpl_IPatchFace.hxx>
  28 #include <GEOMImpl_IProximity.hxx>
  29 #include <GEOMImpl_PatchFaceDriver.hxx>
  30 #include <GEOMImpl_ShapeProximityDriver.hxx>
  31
  32 #include <GEOMImpl_Types.hxx>
  33
  34 #include <GEOMImpl_IConformity.hxx>
  35 #include <GEOMImpl_ConformityDriver.hxx>
  36
  37 #include <GEOMUtils.hxx>
  38
  39 #include <GEOMAlgo_AlgoTools.hxx>
  40 #include <GEOMAlgo_KindOfName.hxx>
  41 #include <GEOMAlgo_ShapeInfoFiller.hxx>
  42
  43 #include <GEOM_PythonDump.hxx>
  44
  45 #include <utilities.h>
  46
  47 // OCCT Includes
  48 #include <Bnd_Box.hxx>
  49 #include <BOPAlgo_CheckerSI.hxx>
  50 #include <BOPDS_DS.hxx>
  51 #include <BOPDS_MapOfPair.hxx>
  52 #include <BOPDS_Pair.hxx>
  53 #include <BOPTools_AlgoTools.hxx>
  54 #include <BRepBndLib.hxx>
  55 #include <BRepBuilderAPI_Copy.hxx>
  56 #include <BRepCheck_ListIteratorOfListOfStatus.hxx>
  57 #include <BRepCheck_Shell.hxx>
  58 #include <BRepClass3d_SolidClassifier.hxx>
  59 #include <BRepClass_FaceClassifier.hxx>
  60 #include <BRepExtrema_DistShapeShape.hxx>
  61 #include <BRepExtrema_ShapeProximity.hxx>
  62 #include <BRepExtrema_SelfIntersection.hxx>
  63 #include <BRepExtrema_MapOfIntegerPackedMapOfInteger.hxx>
  64 #include <BRepGProp.hxx>
  65 #include <BRepTools.hxx>
  66 #include <BRep_Tool.hxx>
  67
  68 #include <Geom_Line.hxx>
  69 #include <Geom_Plane.hxx>
  70 #include <GeomAPI_ProjectPointOnCurve.hxx>
  71 #include <GeomAPI_ProjectPointOnSurf.hxx>
  72 #include <GeomLProp_CLProps.hxx>
  73 #include <GeomLProp_SLProps.hxx>
  74 #include <GProp_GProps.hxx>
  75 #include <GProp_PrincipalProps.hxx>
  76
  77 #include <ShapeAnalysis_ShapeTolerance.hxx>
  78 #include <ShapeAnalysis_Surface.hxx>
  79
  80 #include <TColStd_MapIteratorOfPackedMapOfInteger.hxx>
  81 #include <TopExp.hxx>
  82 #include <TopExp_Explorer.hxx>
  83 #include <TopoDS.hxx>
  84 #include <TopoDS_Edge.hxx>
  85 #include <TopTools_DataMapIteratorOfDataMapOfIntegerListOfShape.hxx>
  86 #include <TopTools_IndexedMapOfShape.hxx>
  87 #include <TopTools_ListIteratorOfListOfShape.hxx>
  88 #include <TopTools_ListOfShape.hxx>
  89
  90 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
  91
  92 #include <set>
  93
  94 //=============================================================================
  95 /*!
  96  *  Constructor
  97  */
  98 //=============================================================================
  99 GEOMImpl_IMeasureOperations::GEOMImpl_IMeasureOperations (GEOM_Engine* theEngine)
 100 : GEOM_IOperations(theEngine)
 101 {
 102   MESSAGE("GEOMImpl_IMeasureOperations::GEOMImpl_IMeasureOperations");
 103 }
 104
 105 //=============================================================================
 106 /*!
 107  *  Destructor
 108  */
 109 //=============================================================================
 110 GEOMImpl_IMeasureOperations::~GEOMImpl_IMeasureOperations()
 111 {
 112   MESSAGE("GEOMImpl_IMeasureOperations::~GEOMImpl_IMeasureOperations");
 113 }
 114
 115 //=============================================================================
 116 /*! Get kind and parameters of the given shape.
 117  */
 118 //=============================================================================
 119 GEOMImpl_IMeasureOperations::ShapeKind GEOMImpl_IMeasureOperations::KindOfShape
 120                              (Handle(GEOM_Object) theShape,
 121                               Handle(TColStd_HSequenceOfInteger)& theIntegers,
 122                               Handle(TColStd_HSequenceOfReal)&    theDoubles)
 123 {
 124   SetErrorCode(KO);
 125   ShapeKind aKind = SK_NO_SHAPE;
 126
 127   if (theIntegers.IsNull()) theIntegers = new TColStd_HSequenceOfInteger;
 128   else                      theIntegers->Clear();
 129
 130   if (theDoubles.IsNull()) theDoubles = new TColStd_HSequenceOfReal;
 131   else                     theDoubles->Clear();
 132
 133   if (theShape.IsNull())
 134     return aKind;
 135
 136   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 137   if (aRefShape.IsNull()) return aKind;
 138
 139   TopoDS_Shape aShape = aRefShape->GetValue();
 140   if (aShape.IsNull()) return aKind;
 141
 142   int geom_type = theShape->GetType();
 143
 144   // check if it's advanced shape
 145   if ( geom_type > USER_TYPE ) {
 146     SetErrorCode(OK);
 147     return SK_ADVANCED;
 148   }
 149
 150   // Call algorithm
 151   GEOMAlgo_ShapeInfoFiller aSF;
 152   aSF.SetShape(aShape);
 153   aSF.Perform();
 154
 155   Standard_Integer iErr = aSF.ErrorStatus();
 156
 157   if (iErr) {
 158     SetErrorCode("Error in GEOMAlgo_ShapeInfoFiller");
 159     return SK_NO_SHAPE;
 160   }
 161   const GEOMAlgo_ShapeInfo& anInfo = aSF.Info();
 162
 163   // specific processing for some "advanced" objects
 164   switch ( geom_type ) {
 165   case GEOM_MARKER:
 166     // local coordinate system
 167     // (+) geompy.kind.LCS  xc yc zc xx xy xz yx yy yz zx zy zz
 168
 169     TopoDS_Face aFace = TopoDS::Face( aShape );
 170     Handle(Geom_Plane) aPlane = Handle(Geom_Plane)::DownCast( BRep_Tool::Surface( aFace ) );
 171     gp_Pnt aC = aPlane->Pln().Location();
 172     gp_Ax3 anAx3 = aPlane->Pln().Position();
 173
 174     theDoubles->Append(aC.X());
 175     theDoubles->Append(aC.Y());
 176     theDoubles->Append(aC.Z());
 177
 178     gp_Dir aD = anAx3.XDirection();
 179     theDoubles->Append(aD.X());
 180     theDoubles->Append(aD.Y());
 181     theDoubles->Append(aD.Z());
 182     aD = anAx3.YDirection();
 183     theDoubles->Append(aD.X());
 184     theDoubles->Append(aD.Y());
 185     theDoubles->Append(aD.Z());
 186     aD = anAx3.Direction();
 187     theDoubles->Append(aD.X());
 188     theDoubles->Append(aD.Y());
 189     theDoubles->Append(aD.Z());
 190
 191     SetErrorCode(OK);
 192     return SK_LCS;
 193   }
 194
 195   // Interpret results
 196   TopAbs_ShapeEnum aType = anInfo.Type();
 197   switch (aType)
 198   {
 199   case TopAbs_COMPOUND:
 200   case TopAbs_COMPSOLID:
 201     {
 202       // (+) geompy.kind.COMPOUND     nb_solids nb_faces nb_edges nb_vertices
 203       // (+) geompy.kind.COMPSOLID    nb_solids nb_faces nb_edges nb_vertices
 204       // ??? "nb_faces" - all faces or only 'standalone' faces?
 205       if (aType == TopAbs_COMPOUND)
 206         aKind = SK_COMPOUND;
 207       else
 208         aKind = SK_COMPSOLID;
 209
 210       //theIntegers->Append(anInfo.NbSubShapes(TopAbs_COMPOUND));
 211       //theIntegers->Append(anInfo.NbSubShapes(TopAbs_COMPSOLID));
 212       theIntegers->Append(anInfo.NbSubShapes(TopAbs_SOLID));
 213       theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
 214       theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 215       theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 216     }
 217     break;
 218
 219   case TopAbs_SHELL:
 220     {
 221       // (+) geompy.kind.SHELL  geompy.info.closed   nb_faces nb_edges nb_vertices
 222       // (+) geompy.kind.SHELL  geompy.info.unclosed nb_faces nb_edges nb_vertices
 223       aKind = SK_SHELL;
 224
 225       theIntegers->Append((int)anInfo.KindOfClosed());
 226
 227       theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
 228       theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 229       theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 230     }
 231     break;
 232
 233   case TopAbs_WIRE:
 234     {
 235       // (+) geompy.kind.WIRE  geompy.info.closed   nb_edges nb_vertices
 236       // (+) geompy.kind.WIRE  geompy.info.unclosed nb_edges nb_vertices
 237       aKind = SK_WIRE;
 238
 239       theIntegers->Append((int)anInfo.KindOfClosed());
 240
 241       theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 242       theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 243     }
 244     break;
 245
 246   case TopAbs_SOLID:
 247     {
 248       aKind = SK_SOLID;
 249
 250       GEOMAlgo_KindOfName aKN = anInfo.KindOfName();
 251       switch (aKN)
 252       {
 253       case GEOMAlgo_KN_SPHERE:
 254         // (+) geompy.kind.SPHERE  xc yc zc  R
 255         {
 256           aKind = SK_SPHERE;
 257
 258           gp_Pnt aC = anInfo.Location();
 259           theDoubles->Append(aC.X());
 260           theDoubles->Append(aC.Y());
 261           theDoubles->Append(aC.Z());
 262
 263           theDoubles->Append(anInfo.Radius1());
 264         }
 265         break;
 266       case GEOMAlgo_KN_CYLINDER:
 267         // (+) geompy.kind.CYLINDER  xb yb zb  dx dy dz  R  H
 268         {
 269           aKind = SK_CYLINDER;
 270
 271           gp_Pnt aC = anInfo.Location();
 272           theDoubles->Append(aC.X());
 273           theDoubles->Append(aC.Y());
 274           theDoubles->Append(aC.Z());
 275
 276           gp_Ax3 anAx3 = anInfo.Position();
 277           gp_Dir aD = anAx3.Direction();
 278           theDoubles->Append(aD.X());
 279           theDoubles->Append(aD.Y());
 280           theDoubles->Append(aD.Z());
 281
 282           theDoubles->Append(anInfo.Radius1());
 283           theDoubles->Append(anInfo.Height());
 284         }
 285         break;
 286       case GEOMAlgo_KN_BOX:
 287         // (+) geompy.kind.BOX  xc yc zc  ax ay az
 288         {
 289           aKind = SK_BOX;
 290
 291           gp_Pnt aC = anInfo.Location();
 292           theDoubles->Append(aC.X());
 293           theDoubles->Append(aC.Y());
 294           theDoubles->Append(aC.Z());
 295
 296           gp_Ax3 anAx3 = anInfo.Position();
 297           gp_Dir aD = anAx3.Direction();
 298           gp_Dir aX = anAx3.XDirection();
 299
 300           // ax ay az
 301           if (aD.IsParallel(gp::DZ(), Precision::Angular()) &&
 302               aX.IsParallel(gp::DX(), Precision::Angular())) {
 303             theDoubles->Append(anInfo.Length()); // ax'
 304             theDoubles->Append(anInfo.Width());  // ay'
 305             theDoubles->Append(anInfo.Height()); // az'
 306           }
 307           else if (aD.IsParallel(gp::DZ(), Precision::Angular()) &&
 308                    aX.IsParallel(gp::DY(), Precision::Angular())) {
 309             theDoubles->Append(anInfo.Width());  // ay'
 310             theDoubles->Append(anInfo.Length()); // ax'
 311             theDoubles->Append(anInfo.Height()); // az'
 312           }
 313           else if (aD.IsParallel(gp::DX(), Precision::Angular()) &&
 314                    aX.IsParallel(gp::DZ(), Precision::Angular())) {
 315             theDoubles->Append(anInfo.Height()); // az'
 316             theDoubles->Append(anInfo.Width());  // ay'
 317             theDoubles->Append(anInfo.Length()); // ax'
 318           }
 319           else if (aD.IsParallel(gp::DX(), Precision::Angular()) &&
 320                    aX.IsParallel(gp::DY(), Precision::Angular())) {
 321             theDoubles->Append(anInfo.Height()); // az'
 322             theDoubles->Append(anInfo.Length()); // ax'
 323             theDoubles->Append(anInfo.Width());  // ay'
 324           }
 325           else if (aD.IsParallel(gp::DY(), Precision::Angular()) &&
 326                    aX.IsParallel(gp::DZ(), Precision::Angular())) {
 327             theDoubles->Append(anInfo.Width());  // ay'
 328             theDoubles->Append(anInfo.Height()); // az'
 329             theDoubles->Append(anInfo.Length()); // ax'
 330           }
 331           else if (aD.IsParallel(gp::DY(), Precision::Angular()) &&
 332                    aX.IsParallel(gp::DX(), Precision::Angular())) {
 333             theDoubles->Append(anInfo.Length()); // ax'
 334             theDoubles->Append(anInfo.Height()); // az'
 335             theDoubles->Append(anInfo.Width());  // ay'
 336           }
 337           else {
 338             // (+) geompy.kind.ROTATED_BOX  xo yo zo  zx zy zz  xx xy xz  ax ay az
 339             aKind = SK_ROTATED_BOX;
 340
 341             // Direction and XDirection
 342             theDoubles->Append(aD.X());
 343             theDoubles->Append(aD.Y());
 344             theDoubles->Append(aD.Z());
 345
 346             theDoubles->Append(aX.X());
 347             theDoubles->Append(aX.Y());
 348             theDoubles->Append(aX.Z());
 349
 350             // ax ay az
 351             theDoubles->Append(anInfo.Length());
 352             theDoubles->Append(anInfo.Width());
 353             theDoubles->Append(anInfo.Height());
 354           }
 355         }
 356         break;
 357       case GEOMAlgo_KN_TORUS:
 358         // (+) geompy.kind.TORUS  xc yc zc  dx dy dz  R_1 R_2
 359         {
 360           aKind = SK_TORUS;
 361
 362           gp_Pnt aO = anInfo.Location();
 363           theDoubles->Append(aO.X());
 364           theDoubles->Append(aO.Y());
 365           theDoubles->Append(aO.Z());
 366
 367           gp_Ax3 anAx3 = anInfo.Position();
 368           gp_Dir aD = anAx3.Direction();
 369           theDoubles->Append(aD.X());
 370           theDoubles->Append(aD.Y());
 371           theDoubles->Append(aD.Z());
 372
 373           theDoubles->Append(anInfo.Radius1());
 374           theDoubles->Append(anInfo.Radius2());
 375         }
 376         break;
 377       case GEOMAlgo_KN_CONE:
 378         // (+) geompy.kind.CONE  xb yb zb  dx dy dz  R_1 R_2  H
 379         {
 380           aKind = SK_CONE;
 381
 382           gp_Pnt aO = anInfo.Location();
 383           theDoubles->Append(aO.X());
 384           theDoubles->Append(aO.Y());
 385           theDoubles->Append(aO.Z());
 386
 387           gp_Ax3 anAx3 = anInfo.Position();
 388           gp_Dir aD = anAx3.Direction();
 389           theDoubles->Append(aD.X());
 390           theDoubles->Append(aD.Y());
 391           theDoubles->Append(aD.Z());
 392
 393           theDoubles->Append(anInfo.Radius1());
 394           theDoubles->Append(anInfo.Radius2());
 395           theDoubles->Append(anInfo.Height());
 396         }
 397         break;
 398       case GEOMAlgo_KN_POLYHEDRON:
 399         // (+) geompy.kind.POLYHEDRON  nb_faces nb_edges nb_vertices
 400         {
 401           aKind = SK_POLYHEDRON;
 402
 403           theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
 404           theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 405           theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 406         }
 407         break;
 408       default:
 409         // (+) geompy.kind.SOLID  nb_faces nb_edges nb_vertices
 410         {
 411           theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
 412           theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 413           theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 414         }
 415       }
 416     }
 417     break;
 418
 419   case TopAbs_FACE:
 420     {
 421       aKind = SK_FACE;
 422
 423       GEOMAlgo_KindOfName aKN = anInfo.KindOfName();
 424       switch (aKN) {
 425       case GEOMAlgo_KN_SPHERE:
 426         // (+) geompy.kind.SPHERE2D  xc yc zc  R
 427         {
 428           aKind = SK_SPHERE2D;
 429
 430           gp_Pnt aC = anInfo.Location();
 431           theDoubles->Append(aC.X());
 432           theDoubles->Append(aC.Y());
 433           theDoubles->Append(aC.Z());
 434
 435           theDoubles->Append(anInfo.Radius1());
 436         }
 437         break;
 438       case GEOMAlgo_KN_CYLINDER:
 439         // (+) geompy.kind.CYLINDER2D  xb yb zb  dx dy dz  R  H
 440         {
 441           aKind = SK_CYLINDER2D;
 442
 443           gp_Pnt aO = anInfo.Location();
 444           theDoubles->Append(aO.X());
 445           theDoubles->Append(aO.Y());
 446           theDoubles->Append(aO.Z());
 447
 448           gp_Ax3 anAx3 = anInfo.Position();
 449           gp_Dir aD = anAx3.Direction();
 450           theDoubles->Append(aD.X());
 451           theDoubles->Append(aD.Y());
 452           theDoubles->Append(aD.Z());
 453
 454           theDoubles->Append(anInfo.Radius1());
 455           theDoubles->Append(anInfo.Height());
 456         }
 457         break;
 458       case GEOMAlgo_KN_TORUS:
 459         // (+) geompy.kind.TORUS2D  xc yc zc  dx dy dz  R_1 R_2
 460         {
 461           aKind = SK_TORUS2D;
 462
 463           gp_Pnt aO = anInfo.Location();
 464           theDoubles->Append(aO.X());
 465           theDoubles->Append(aO.Y());
 466           theDoubles->Append(aO.Z());
 467
 468           gp_Ax3 anAx3 = anInfo.Position();
 469           gp_Dir aD = anAx3.Direction();
 470           theDoubles->Append(aD.X());
 471           theDoubles->Append(aD.Y());
 472           theDoubles->Append(aD.Z());
 473
 474           theDoubles->Append(anInfo.Radius1());
 475           theDoubles->Append(anInfo.Radius2());
 476         }
 477         break;
 478       case GEOMAlgo_KN_CONE:
 479         // (+) geompy.kind.CONE2D  xc yc zc  dx dy dz  R_1 R_2  H
 480         {
 481           aKind = SK_CONE2D;
 482
 483           gp_Pnt aO = anInfo.Location();
 484           theDoubles->Append(aO.X());
 485           theDoubles->Append(aO.Y());
 486           theDoubles->Append(aO.Z());
 487
 488           gp_Ax3 anAx3 = anInfo.Position();
 489           gp_Dir aD = anAx3.Direction();
 490           theDoubles->Append(aD.X());
 491           theDoubles->Append(aD.Y());
 492           theDoubles->Append(aD.Z());
 493
 494           theDoubles->Append(anInfo.Radius1());
 495           theDoubles->Append(anInfo.Radius2());
 496           theDoubles->Append(anInfo.Height());
 497         }
 498         break;
 499       case GEOMAlgo_KN_DISKCIRCLE:
 500         // (+) geompy.kind.DISK_CIRCLE  xc yc zc  dx dy dz  R
 501         {
 502           aKind = SK_DISK_CIRCLE;
 503
 504           gp_Pnt aC = anInfo.Location();
 505           theDoubles->Append(aC.X());
 506           theDoubles->Append(aC.Y());
 507           theDoubles->Append(aC.Z());
 508
 509           gp_Ax3 anAx3 = anInfo.Position();
 510           gp_Dir aD = anAx3.Direction();
 511           theDoubles->Append(aD.X());
 512           theDoubles->Append(aD.Y());
 513           theDoubles->Append(aD.Z());
 514
 515           theDoubles->Append(anInfo.Radius1());
 516         }
 517         break;
 518       case GEOMAlgo_KN_DISKELLIPSE:
 519         // (+) geompy.kind.DISK_ELLIPSE  xc yc zc  dx dy dz  R_1 R_2
 520         {
 521           aKind = SK_DISK_ELLIPSE;
 522
 523           gp_Pnt aC = anInfo.Location();
 524           theDoubles->Append(aC.X());
 525           theDoubles->Append(aC.Y());
 526           theDoubles->Append(aC.Z());
 527
 528           gp_Ax3 anAx3 = anInfo.Position();
 529           gp_Dir aD = anAx3.Direction();
 530           theDoubles->Append(aD.X());
 531           theDoubles->Append(aD.Y());
 532           theDoubles->Append(aD.Z());
 533
 534           theDoubles->Append(anInfo.Radius1());
 535           theDoubles->Append(anInfo.Radius2());
 536         }
 537         break;
 538       case GEOMAlgo_KN_RECTANGLE:
 539       case GEOMAlgo_KN_TRIANGLE:
 540       case GEOMAlgo_KN_QUADRANGLE:
 541       case GEOMAlgo_KN_POLYGON:
 542         // (+) geompy.kind.POLYGON  xo yo zo  dx dy dz  nb_edges nb_vertices
 543         {
 544           aKind = SK_POLYGON;
 545
 546           gp_Pnt aO = anInfo.Location();
 547           theDoubles->Append(aO.X());
 548           theDoubles->Append(aO.Y());
 549           theDoubles->Append(aO.Z());
 550
 551           gp_Ax3 anAx3 = anInfo.Position();
 552           gp_Dir aD = anAx3.Direction();
 553           theDoubles->Append(aD.X());
 554           theDoubles->Append(aD.Y());
 555           theDoubles->Append(aD.Z());
 556
 557           theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 558           theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 559         }
 560         break;
 561       case GEOMAlgo_KN_PLANE: // infinite
 562         // (+) geompy.kind.PLANE  xo yo zo  dx dy dz
 563         {
 564           aKind = SK_PLANE;
 565
 566           gp_Pnt aC = anInfo.Location();
 567           theDoubles->Append(aC.X());
 568           theDoubles->Append(aC.Y());
 569           theDoubles->Append(aC.Z());
 570
 571           gp_Ax3 anAx3 = anInfo.Position();
 572           gp_Dir aD = anAx3.Direction();
 573           theDoubles->Append(aD.X());
 574           theDoubles->Append(aD.Y());
 575           theDoubles->Append(aD.Z());
 576
 577           if (anInfo.KindOfBounds() != GEOMAlgo_KB_INFINITE)
 578           {
 579             // (+) geompy.kind.PLANAR  xo yo zo  dx dy dz  nb_edges nb_vertices
 580             aKind = SK_PLANAR;
 581
 582             theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 583             theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 584           }
 585         }
 586         break;
 587       default:
 588         // ??? geompy.kind.FACE  nb_edges nb_vertices _surface_type_id_
 589         // (+) geompy.kind.FACE  nb_edges nb_vertices
 590         theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
 591         theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 592       }
 593     }
 594     break;
 595
 596   case TopAbs_EDGE:
 597     {
 598       aKind = SK_EDGE;
 599
 600       GEOMAlgo_KindOfName aKN = anInfo.KindOfName();
 601       switch (aKN) {
 602       case GEOMAlgo_KN_CIRCLE:
 603         {
 604           // (+) geompy.kind.CIRCLE  xc yc zc  dx dy dz  R
 605           aKind = SK_CIRCLE;
 606
 607           gp_Pnt aC = anInfo.Location();
 608           theDoubles->Append(aC.X());
 609           theDoubles->Append(aC.Y());
 610           theDoubles->Append(aC.Z());
 611
 612           gp_Ax3 anAx3 = anInfo.Position();
 613           gp_Dir aD = anAx3.Direction();
 614           theDoubles->Append(aD.X());
 615           theDoubles->Append(aD.Y());
 616           theDoubles->Append(aD.Z());
 617
 618           theDoubles->Append(anInfo.Radius1());
 619         }
 620         break;
 621       case GEOMAlgo_KN_ARCCIRCLE:
 622         {
 623           // (+) geompy.kind.ARC_CIRCLE  xc yc zc  dx dy dz  R  x1 y1 z1  x2 y2 z2
 624           aKind = SK_ARC_CIRCLE;
 625
 626           gp_Pnt aC = anInfo.Location();
 627           theDoubles->Append(aC.X());
 628           theDoubles->Append(aC.Y());
 629           theDoubles->Append(aC.Z());
 630
 631           gp_Ax3 anAx3 = anInfo.Position();
 632           gp_Dir aD = anAx3.Direction();
 633           theDoubles->Append(aD.X());
 634           theDoubles->Append(aD.Y());
 635           theDoubles->Append(aD.Z());
 636
 637           theDoubles->Append(anInfo.Radius1());
 638
 639           gp_Pnt aP1 = anInfo.Pnt1();
 640           theDoubles->Append(aP1.X());
 641           theDoubles->Append(aP1.Y());
 642           theDoubles->Append(aP1.Z());
 643
 644           gp_Pnt aP2 = anInfo.Pnt2();
 645           theDoubles->Append(aP2.X());
 646           theDoubles->Append(aP2.Y());
 647           theDoubles->Append(aP2.Z());
 648         }
 649         break;
 650       case GEOMAlgo_KN_ELLIPSE:
 651         {
 652           // (+) geompy.kind.ELLIPSE  xc yc zc  dx dy dz  R_1 R_2  xVx yVx zVx  xVy yVy zVy
 653           aKind = SK_ELLIPSE;
 654
 655           gp_Pnt aC = anInfo.Location();
 656           theDoubles->Append(aC.X());
 657           theDoubles->Append(aC.Y());
 658           theDoubles->Append(aC.Z());
 659
 660           gp_Ax3 anAx3 = anInfo.Position();
 661           gp_Dir aD = anAx3.Direction();
 662           theDoubles->Append(aD.X());
 663           theDoubles->Append(aD.Y());
 664           theDoubles->Append(aD.Z());
 665
 666           theDoubles->Append(anInfo.Radius1());
 667           theDoubles->Append(anInfo.Radius2());
 668
 669           gp_Dir aXD = anAx3.XDirection();
 670           theDoubles->Append(aXD.X());
 671           theDoubles->Append(aXD.Y());
 672           theDoubles->Append(aXD.Z());
 673
 674           gp_Dir aYD = anAx3.YDirection();
 675           theDoubles->Append(aYD.X());
 676           theDoubles->Append(aYD.Y());
 677           theDoubles->Append(aYD.Z());
 678         }
 679         break;
 680       case GEOMAlgo_KN_ARCELLIPSE:
 681         {
 682           // (+) geompy.kind.ARC_ELLIPSE  xc yc zc  dx dy dz  R_1 R_2  x1 y1 z1  x2 y2 z2  xVx yVx zVx  xVy yVy zVy
 683           aKind = SK_ARC_ELLIPSE;
 684
 685           gp_Pnt aC = anInfo.Location();
 686           theDoubles->Append(aC.X());
 687           theDoubles->Append(aC.Y());
 688           theDoubles->Append(aC.Z());
 689
 690           gp_Ax3 anAx3 = anInfo.Position();
 691           gp_Dir aD = anAx3.Direction();
 692           theDoubles->Append(aD.X());
 693           theDoubles->Append(aD.Y());
 694           theDoubles->Append(aD.Z());
 695
 696           theDoubles->Append(anInfo.Radius1());
 697           theDoubles->Append(anInfo.Radius2());
 698
 699           gp_Pnt aP1 = anInfo.Pnt1();
 700           theDoubles->Append(aP1.X());
 701           theDoubles->Append(aP1.Y());
 702           theDoubles->Append(aP1.Z());
 703
 704           gp_Pnt aP2 = anInfo.Pnt2();
 705           theDoubles->Append(aP2.X());
 706           theDoubles->Append(aP2.Y());
 707           theDoubles->Append(aP2.Z());
 708
 709           gp_Dir aXD = anAx3.XDirection();
 710           theDoubles->Append(aXD.X());
 711           theDoubles->Append(aXD.Y());
 712           theDoubles->Append(aXD.Z());
 713
 714           gp_Dir aYD = anAx3.YDirection();
 715           theDoubles->Append(aYD.X());
 716           theDoubles->Append(aYD.Y());
 717           theDoubles->Append(aYD.Z());
 718         }
 719         break;
 720       case GEOMAlgo_KN_LINE:
 721         {
 722           // ??? geompy.kind.LINE  x1 y1 z1  x2 y2 z2
 723           // (+) geompy.kind.LINE  x1 y1 z1  dx dy dz
 724           aKind = SK_LINE;
 725
 726           gp_Pnt aO = anInfo.Location();
 727           theDoubles->Append(aO.X());
 728           theDoubles->Append(aO.Y());
 729           theDoubles->Append(aO.Z());
 730
 731           gp_Dir aD = anInfo.Direction();
 732           theDoubles->Append(aD.X());
 733           theDoubles->Append(aD.Y());
 734           theDoubles->Append(aD.Z());
 735         }
 736         break;
 737       case GEOMAlgo_KN_SEGMENT:
 738         {
 739           // (+) geompy.kind.SEGMENT  x1 y1 z1  x2 y2 z2
 740           aKind = SK_SEGMENT;
 741
 742           gp_Pnt aP1 = anInfo.Pnt1();
 743           theDoubles->Append(aP1.X());
 744           theDoubles->Append(aP1.Y());
 745           theDoubles->Append(aP1.Z());
 746
 747           gp_Pnt aP2 = anInfo.Pnt2();
 748           theDoubles->Append(aP2.X());
 749           theDoubles->Append(aP2.Y());
 750           theDoubles->Append(aP2.Z());
 751         }
 752         break;
 753       case GEOMAlgo_KN_CURVEBSPLINE:
 754         {
 755           // (+) geompy.kind.CRV_BSPLINE  p d np nk nw nm  x1 y1 z1 ... xnp ynp znp  k1 ... knk  w1 ... wnw  m1 ... mnm
 756           aKind = SK_CRV_BSPLINE;
 757           Standard_Integer aNbPoles = anInfo.NbPoles();
 758           Standard_Integer aNbKnots = anInfo.NbKnots();
 759           Standard_Integer aNbWeights = anInfo.NbWeights();
 760           Standard_Integer aNbMultiplicities = anInfo.NbMultiplicities();
 761
 762           theIntegers->Append(anInfo.KindOfPeriod() == GEOMAlgo_KP_PERIODIC ? 1 : 0);
 763           theIntegers->Append(anInfo.Degree());
 764           theIntegers->Append(aNbPoles);
 765           theIntegers->Append(aNbKnots);
 766           theIntegers->Append(aNbWeights);
 767           theIntegers->Append(aNbMultiplicities);
 768           //
 769           Standard_Integer i;
 770           if (aNbPoles > 0) {
 771             Handle(TColgp_HArray1OfPnt) aPoles = anInfo.Poles();
 772             if (aPoles.IsNull() || aPoles->Length() != aNbPoles) {
 773               SetErrorCode("B-Spline Curve: no or wrong number of poles given");
 774               return aKind;
 775             }
 776             for (i=1; i<=aNbPoles; i++) {
 777               const gp_Pnt &aP = aPoles->Value(i);
 778               theDoubles->Append(aP.X());
 779               theDoubles->Append(aP.Y());
 780               theDoubles->Append(aP.Z());
 781             }
 782           }
 783           //
 784           if (aNbKnots > 0) {
 785             Handle(TColStd_HArray1OfReal) aKnots = anInfo.Knots();
 786             if (aKnots.IsNull() || aKnots->Length() != aNbKnots) {
 787               SetErrorCode("B-Spline Curve: no or wrong number of knots given");
 788               return aKind;
 789             }
 790             for (i=1; i<=aNbKnots; i++)
 791               theDoubles->Append(aKnots->Value(i));
 792           }
 793           //
 794           if (aNbWeights > 0) {
 795             Handle(TColStd_HArray1OfReal) aWeights = anInfo.Weights();
 796             if (aNbWeights > 0 && (aWeights.IsNull() || aWeights->Length() != aNbWeights)) {
 797               SetErrorCode("B-Spline Curve: no or wrong number of weights given");
 798               return aKind;
 799             }
 800             for (i=1; i<=aNbWeights; i++)
 801               theDoubles->Append(aWeights->Value(i));
 802           }
 803           //
 804           if (aNbMultiplicities > 0) {
 805             Handle(TColStd_HArray1OfInteger) aMults = anInfo.Multiplicities();
 806             if (aMults.IsNull() || aMults->Length() != aNbMultiplicities) {
 807               SetErrorCode("B-Spline Curve: no or wrong number of multiplicities given");
 808               return aKind;
 809             }
 810             for (i=1; i<=aNbMultiplicities; i++)
 811               theIntegers->Append(aMults->Value(i));
 812           }
 813         }
 814         break;
 815       case GEOMAlgo_KN_CURVEBEZIER:
 816         {
 817           // (+) geompy.kind.CRV_BEZIER  np nw  x1 y1 z1 ... xnp ynp znp  w1 ... wnw
 818           aKind = SK_CRV_BEZIER;
 819           Standard_Integer aNbPoles = anInfo.NbPoles();
 820           Standard_Integer aNbWeights = anInfo.NbWeights();
 821
 822           theIntegers->Append(aNbPoles);
 823           theIntegers->Append(aNbWeights);
 824           //
 825           Standard_Integer i;
 826           if (aNbPoles > 0) {
 827             Handle(TColgp_HArray1OfPnt) aPoles = anInfo.Poles();
 828             if (aPoles.IsNull() || aPoles->Length() != aNbPoles) {
 829               SetErrorCode("Bezier Curve: no or wrong number of poles given");
 830               return aKind;
 831             }
 832             for (i=1; i<=aNbPoles; i++) {
 833               const gp_Pnt &aP = aPoles->Value(i);
 834               theDoubles->Append(aP.X());
 835               theDoubles->Append(aP.Y());
 836               theDoubles->Append(aP.Z());
 837             }
 838           }
 839           //
 840           if (aNbWeights > 0) {
 841             Handle(TColStd_HArray1OfReal) aWeights = anInfo.Weights();
 842             if (aNbWeights > 0 && (aWeights.IsNull() || aWeights->Length() != aNbWeights)) {
 843               SetErrorCode("B-Spline Curve: no or wrong number of weights given");
 844               return aKind;
 845             }
 846             for (i=1; i<=aNbWeights; i++)
 847               theDoubles->Append(aWeights->Value(i));
 848           }
 849         }
 850         break;
 851       case GEOMAlgo_KN_HYPERBOLA:
 852         {
 853           // (+) geompy.kind.HYPERBOLA  xc yc zc  dx dy dz  R_1 R_2  xVx yVx zVx  xVy yVy zVy
 854           aKind = SK_HYPERBOLA;
 855
 856           gp_Pnt aC = anInfo.Location();
 857           theDoubles->Append(aC.X());
 858           theDoubles->Append(aC.Y());
 859           theDoubles->Append(aC.Z());
 860
 861           gp_Ax3 anAx3 = anInfo.Position();
 862           gp_Dir aD = anAx3.Direction();
 863           theDoubles->Append(aD.X());
 864           theDoubles->Append(aD.Y());
 865           theDoubles->Append(aD.Z());
 866
 867           theDoubles->Append(anInfo.Radius1());
 868           theDoubles->Append(anInfo.Radius2());
 869
 870           gp_Dir aXD = anAx3.XDirection();
 871           theDoubles->Append(aXD.X());
 872           theDoubles->Append(aXD.Y());
 873           theDoubles->Append(aXD.Z());
 874
 875           gp_Dir aYD = anAx3.YDirection();
 876           theDoubles->Append(aYD.X());
 877           theDoubles->Append(aYD.Y());
 878           theDoubles->Append(aYD.Z());
 879         }
 880         break;
 881       case GEOMAlgo_KN_PARABOLA:
 882         {
 883           // (+) geompy.kind.PARABOLA  xc yc zc  dx dy dz  F  xVx yVx zVx  xVy yVy zVy
 884           aKind = SK_PARABOLA;
 885
 886           gp_Pnt aC = anInfo.Location();
 887           theDoubles->Append(aC.X());
 888           theDoubles->Append(aC.Y());
 889           theDoubles->Append(aC.Z());
 890
 891           gp_Ax3 anAx3 = anInfo.Position();
 892           gp_Dir aD = anAx3.Direction();
 893           theDoubles->Append(aD.X());
 894           theDoubles->Append(aD.Y());
 895           theDoubles->Append(aD.Z());
 896
 897           theDoubles->Append(anInfo.Radius1());
 898
 899           gp_Dir aXD = anAx3.XDirection();
 900           theDoubles->Append(aXD.X());
 901           theDoubles->Append(aXD.Y());
 902           theDoubles->Append(aXD.Z());
 903
 904           gp_Dir aYD = anAx3.YDirection();
 905           theDoubles->Append(aYD.X());
 906           theDoubles->Append(aYD.Y());
 907           theDoubles->Append(aYD.Z());
 908         }
 909         break;
 910       default:
 911         // ??? geompy.kind.EDGE  nb_vertices _curve_type_id_
 912         // (+) geompy.kind.EDGE  nb_vertices
 913         theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
 914       }
 915     }
 916     break;
 917
 918   case TopAbs_VERTEX:
 919     {
 920       // (+) geompy.kind.VERTEX  x y z
 921       aKind = SK_VERTEX;
 922
 923       gp_Pnt aP = anInfo.Location();
 924       theDoubles->Append(aP.X());
 925       theDoubles->Append(aP.Y());
 926       theDoubles->Append(aP.Z());
 927     }
 928     break;
 929   default:;
 930   }
 931
 932   SetErrorCode(OK);
 933   return aKind;
 934 }
 935
 936 //=============================================================================
 937 /*!
 938  *  GetPosition
 939  */
 940 //=============================================================================
 941 void GEOMImpl_IMeasureOperations::GetPosition
 942                    (Handle(GEOM_Object) theShape,
 943                     Standard_Real& Ox, Standard_Real& Oy, Standard_Real& Oz,
 944                     Standard_Real& Zx, Standard_Real& Zy, Standard_Real& Zz,
 945                     Standard_Real& Xx, Standard_Real& Xy, Standard_Real& Xz)
 946 {
 947   SetErrorCode(KO);
 948
 949   //Set default values: global CS
 950   Ox = Oy = Oz = Zx = Zy = Xy = Xz = 0.;
 951   Zz = Xx = 1.;
 952
 953   if (theShape.IsNull()) return;
 954
 955   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 956   if (aRefShape.IsNull()) return;
 957
 958   TopoDS_Shape aShape = aRefShape->GetValue();
 959   if (aShape.IsNull()) {
 960     SetErrorCode("The Objects has NULL Shape");
 961     return;
 962   }
 963
 964   try {
 965     OCC_CATCH_SIGNALS;
 966
 967     gp_Ax3 anAx3 = GEOMUtils::GetPosition(aShape);
 968
 969     gp_Pnt anOri = anAx3.Location();
 970     gp_Dir aDirZ = anAx3.Direction();
 971     gp_Dir aDirX = anAx3.XDirection();
 972
 973     // Output values
 974     anOri.Coord(Ox, Oy, Oz);
 975     aDirZ.Coord(Zx, Zy, Zz);
 976     aDirX.Coord(Xx, Xy, Xz);
 977   }
 978   catch (Standard_Failure& aFail) {
 979     SetErrorCode(aFail.GetMessageString());
 980     return;
 981   }
 982
 983   SetErrorCode(OK);
 984 }
 985
 986 //=============================================================================
 987 /*!
 988  *  GetCentreOfMass
 989  */
 990 //=============================================================================
 991 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetCentreOfMass
 992                                                 (Handle(GEOM_Object) theShape)
 993 {
 994   SetErrorCode(KO);
 995
 996   if (theShape.IsNull()) return NULL;
 997
 998   //Add a new CentreOfMass object
 999   Handle(GEOM_Object) aCDG = GetEngine()->AddObject(GEOM_CDG);
1000
1001   //Add a new CentreOfMass function
1002   Handle(GEOM_Function) aFunction =
1003     aCDG->AddFunction(GEOMImpl_MeasureDriver::GetID(), CDG_MEASURE);
1004   if (aFunction.IsNull()) return NULL;
1005
1006   //Check if the function is set correctly
1007   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
1008
1009   GEOMImpl_IMeasure aCI (aFunction);
1010
1011   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1012   if (aRefShape.IsNull()) return NULL;
1013
1014   aCI.SetBase(aRefShape);
1015
1016   //Compute the CentreOfMass value
1017   try {
1018     OCC_CATCH_SIGNALS;
1019     if (!GetSolver()->ComputeFunction(aFunction)) {
1020       SetErrorCode("Measure driver failed to compute centre of mass");
1021       return NULL;
1022     }
1023   }
1024   catch (Standard_Failure& aFail) {
1025     SetErrorCode(aFail.GetMessageString());
1026     return NULL;
1027   }
1028
1029   //Make a Python command
1030   GEOM::TPythonDump(aFunction) << aCDG << " = geompy.MakeCDG(" << theShape << ")";
1031
1032   SetErrorCode(OK);
1033   return aCDG;
1034 }
1035
1036 //=============================================================================
1037 /*!
1038  *  GetVertexByIndex
1039  */
1040 //=============================================================================
1041 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetVertexByIndex
1042                                                 (Handle(GEOM_Object) theShape,
1043                                                  Standard_Integer theIndex,
1044                                                  Standard_Boolean theUseOri)
1045 {
1046   SetErrorCode(KO);
1047
1048   if (theShape.IsNull()) return NULL;
1049
1050   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1051   if (aRefShape.IsNull()) return NULL;
1052
1053   //Add a new Vertex object
1054   Handle(GEOM_Object) aVertex = GetEngine()->AddObject(GEOM_POINT);
1055
1056   //Add a function
1057   Handle(GEOM_Function) aFunction =
1058     aVertex->AddFunction(GEOMImpl_MeasureDriver::GetID(), VERTEX_BY_INDEX);
1059   if (aFunction.IsNull()) return NULL;
1060
1061   //Check if the function is set correctly
1062   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
1063
1064   GEOMImpl_IMeasure aCI (aFunction);
1065   aCI.SetBase(aRefShape);
1066   aCI.SetIndex(theIndex);
1067   aCI.SetUseOri(theUseOri);
1068
1069   //Compute
1070   try {
1071     OCC_CATCH_SIGNALS;
1072     if (!GetSolver()->ComputeFunction(aFunction)) {
1073       SetErrorCode("Vertex by index driver failed.");
1074       return NULL;
1075     }
1076   }
1077   catch (Standard_Failure& aFail) {
1078     SetErrorCode(aFail.GetMessageString());
1079     return NULL;
1080   }
1081
1082   //Make a Python command
1083   GEOM::TPythonDump(aFunction) << aVertex << " = geompy.GetVertexByIndex("
1084                                << theShape << ", "
1085                                << theIndex << ", "
1086                                << theUseOri << ")";
1087
1088   SetErrorCode(OK);
1089   return aVertex;
1090 }
1091
1092 //=============================================================================
1093 /*!
1094  *  GetNormal
1095  */
1096 //=============================================================================
1097 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetNormal
1098                                          (Handle(GEOM_Object) theFace,
1099                                           Handle(GEOM_Object) theOptionalPoint)
1100 {
1101   SetErrorCode(KO);
1102
1103   if (theFace.IsNull()) return NULL;
1104
1105   //Add a new Normale object
1106   Handle(GEOM_Object) aNorm = GetEngine()->AddObject(GEOM_VECTOR);
1107
1108   //Add a new Normale function
1109   Handle(GEOM_Function) aFunction =
1110     aNorm->AddFunction(GEOMImpl_MeasureDriver::GetID(), VECTOR_FACE_NORMALE);
1111   if (aFunction.IsNull()) return NULL;
1112
1113   //Check if the function is set correctly
1114   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
1115
1116   GEOMImpl_IMeasure aCI (aFunction);
1117
1118   Handle(GEOM_Function) aFace = theFace->GetLastFunction();
1119   if (aFace.IsNull()) return NULL;
1120
1121   aCI.SetBase(aFace);
1122
1123   if (!theOptionalPoint.IsNull()) {
1124     Handle(GEOM_Function) anOptPnt = theOptionalPoint->GetLastFunction();
1125     aCI.SetPoint(anOptPnt);
1126   }
1127
1128   //Compute the Normale value
1129   try {
1130     OCC_CATCH_SIGNALS;
1131     if (!GetSolver()->ComputeFunction(aFunction)) {
1132       SetErrorCode("Measure driver failed to compute normake of face");
1133       return NULL;
1134     }
1135   }
1136   catch (Standard_Failure& aFail) {
1137     SetErrorCode(aFail.GetMessageString());
1138     return NULL;
1139   }
1140
1141   //Make a Python command
1142   GEOM::TPythonDump pd (aFunction);
1143   pd << aNorm << " = geompy.GetNormal(" << theFace;
1144   if (!theOptionalPoint.IsNull()) {
1145     pd << ", " << theOptionalPoint;
1146   }
1147   pd << ")";
1148
1149   SetErrorCode(OK);
1150   return aNorm;
1151 }
1152
1153 //=============================================================================
1154 /*!
1155  *  GetBasicProperties
1156  */
1157 //=============================================================================
1158 void GEOMImpl_IMeasureOperations::GetBasicProperties (Handle(GEOM_Object) theShape,
1159                                                       const Standard_Real theTolerance,
1160                                                       Standard_Real& theLength,
1161                                                       Standard_Real& theSurfArea,
1162                                                       Standard_Real& theVolume)
1163 {
1164   SetErrorCode(KO);
1165
1166   if (theShape.IsNull()) return;
1167
1168   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1169   if (aRefShape.IsNull()) return;
1170
1171   TopoDS_Shape aShape = aRefShape->GetValue();
1172   if (aShape.IsNull()) {
1173     SetErrorCode("The Objects has NULL Shape");
1174     return;
1175   }
1176
1177   //Compute the parameters
1178   GProp_GProps LProps, SProps;
1179   Standard_Real anEps = theTolerance >= 0 ? theTolerance : 1.e-6;
1180   try {
1181     OCC_CATCH_SIGNALS;
1182     BRepGProp::LinearProperties(aShape, LProps, Standard_True);
1183     theLength = LProps.Mass();
1184
1185     BRepGProp::SurfaceProperties(aShape, SProps, anEps, Standard_True);
1186     theSurfArea = SProps.Mass();
1187
1188     theVolume = 0.0;
1189     if (aShape.ShapeType() < TopAbs_SHELL) {
1190       for (TopExp_Explorer Exp (aShape, TopAbs_SOLID); Exp.More(); Exp.Next()) {
1191         GProp_GProps VProps;
1192         BRepGProp::VolumeProperties(Exp.Current(), VProps, anEps, Standard_True);
1193         theVolume += VProps.Mass();
1194       }
1195     }
1196   }
1197   catch (Standard_Failure& aFail) {
1198     SetErrorCode(aFail.GetMessageString());
1199     return;
1200   }
1201
1202   SetErrorCode(OK);
1203 }
1204
1205 //=============================================================================
1206 /*!
1207  *  GetInertia
1208  */
1209 //=============================================================================
1210 void GEOMImpl_IMeasureOperations::GetInertia
1211                    (Handle(GEOM_Object) theShape,
1212                     Standard_Real& I11, Standard_Real& I12, Standard_Real& I13,
1213                     Standard_Real& I21, Standard_Real& I22, Standard_Real& I23,
1214                     Standard_Real& I31, Standard_Real& I32, Standard_Real& I33,
1215                     Standard_Real& Ix , Standard_Real& Iy , Standard_Real& Iz)
1216 {
1217   SetErrorCode(KO);
1218
1219   if (theShape.IsNull()) return;
1220
1221   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1222   if (aRefShape.IsNull()) return;
1223
1224   TopoDS_Shape aShape = aRefShape->GetValue();
1225   if (aShape.IsNull()) {
1226     SetErrorCode("The Objects has NULL Shape");
1227     return;
1228   }
1229
1230   //Compute the parameters
1231   GProp_GProps System;
1232
1233   try {
1234     OCC_CATCH_SIGNALS;
1235     if (aShape.ShapeType() == TopAbs_VERTEX ||
1236         aShape.ShapeType() == TopAbs_EDGE ||
1237         aShape.ShapeType() == TopAbs_WIRE) {
1238       BRepGProp::LinearProperties(aShape, System, Standard_True);
1239     } else if (aShape.ShapeType() == TopAbs_FACE ||
1240                aShape.ShapeType() == TopAbs_SHELL) {
1241       BRepGProp::SurfaceProperties(aShape, System, Standard_True);
1242     } else {
1243       BRepGProp::VolumeProperties(aShape, System, Standard_True);
1244     }
1245     gp_Mat I = System.MatrixOfInertia();
1246
1247     I11 = I(1,1);
1248     I12 = I(1,2);
1249     I13 = I(1,3);
1250
1251     I21 = I(2,1);
1252     I22 = I(2,2);
1253     I23 = I(2,3);
1254
1255     I31 = I(3,1);
1256     I32 = I(3,2);
1257     I33 = I(3,3);
1258
1259     GProp_PrincipalProps Pr = System.PrincipalProperties();
1260     Pr.Moments(Ix,Iy,Iz);
1261   }
1262   catch (Standard_Failure& aFail) {
1263     SetErrorCode(aFail.GetMessageString());
1264     return;
1265   }
1266
1267   SetErrorCode(OK);
1268 }
1269
1270 //=============================================================================
1271 /*!
1272  *  GetBoundingBox
1273  */
1274 //=============================================================================
1275 void GEOMImpl_IMeasureOperations::GetBoundingBox
1276                                      (Handle(GEOM_Object) theShape,
1277                                       const Standard_Boolean precise,
1278                                       Standard_Real& Xmin, Standard_Real& Xmax,
1279                                       Standard_Real& Ymin, Standard_Real& Ymax,
1280                                       Standard_Real& Zmin, Standard_Real& Zmax)
1281 {
1282   SetErrorCode(KO);
1283
1284   if (theShape.IsNull()) return;
1285
1286   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1287   if (aRefShape.IsNull()) return;
1288
1289   TopoDS_Shape aShape = aRefShape->GetValue();
1290   if (aShape.IsNull()) {
1291     SetErrorCode("The Objects has NULL Shape");
1292     return;
1293   }
1294
1295   //Compute the parameters
1296   Bnd_Box B;
1297
1298   try {
1299     OCC_CATCH_SIGNALS;
1300     BRepBuilderAPI_Copy aCopyTool (aShape);
1301     if (!aCopyTool.IsDone()) {
1302       SetErrorCode("GetBoundingBox Error: Bad shape detected");
1303       return;
1304     }
1305
1306     aShape = aCopyTool.Shape();
1307
1308     // remove triangulation to obtain more exact boundaries
1309     BRepTools::Clean(aShape);
1310
1311     BRepBndLib::Add(aShape, B);
1312
1313     if (precise) {
1314       if (!GEOMUtils::PreciseBoundingBox(aShape, B)) {
1315         SetErrorCode("GetBoundingBox Error: Bounding box cannot be precised");
1316         return;
1317       }
1318     }
1319
1320     B.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
1321   }
1322   catch (Standard_Failure& aFail) {
1323     SetErrorCode(aFail.GetMessageString());
1324     return;
1325   }
1326
1327   SetErrorCode(OK);
1328 }
1329
1330 //=============================================================================
1331 /*!
1332  *  GetBoundingBox
1333  */
1334 //=============================================================================
1335 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetBoundingBox
1336                                                 (Handle(GEOM_Object) theShape,
1337                                                  const Standard_Boolean precise)
1338 {
1339   SetErrorCode(KO);
1340
1341   if (theShape.IsNull()) return NULL;
1342
1343   //Add a new BoundingBox object
1344   Handle(GEOM_Object) aBnd = GetEngine()->AddObject(GEOM_BOX);
1345
1346   //Add a new BoundingBox function
1347   const int aType = (precise ? BND_BOX_MEASURE_PRECISE : BND_BOX_MEASURE);
1348   Handle(GEOM_Function) aFunction =
1349     aBnd->AddFunction(GEOMImpl_MeasureDriver::GetID(), aType);
1350   if (aFunction.IsNull()) return NULL;
1351
1352   //Check if the function is set correctly
1353   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
1354
1355   GEOMImpl_IMeasure aCI (aFunction);
1356
1357   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1358   if (aRefShape.IsNull()) return NULL;
1359
1360   aCI.SetBase(aRefShape);
1361
1362   //Compute the BoundingBox value
1363   try {
1364     OCC_CATCH_SIGNALS;
1365     if (!GetSolver()->ComputeFunction(aFunction)) {
1366       SetErrorCode("Measure driver failed to compute a bounding box");
1367       return NULL;
1368     }
1369   }
1370   catch (Standard_Failure& aFail) {
1371     SetErrorCode(aFail.GetMessageString());
1372     return NULL;
1373   }
1374
1375   //Make a Python command
1376   GEOM::TPythonDump aPd(aFunction);
1377
1378   aPd << aBnd << " = geompy.MakeBoundingBox(" << theShape;
1379
1380   if (precise) {
1381     aPd << ", True";
1382   }
1383
1384   aPd << ")";
1385
1386   SetErrorCode(OK);
1387   return aBnd;
1388 }
1389
1390 //=============================================================================
1391 /*!
1392  *  GetTolerance
1393  */
1394 //=============================================================================
1395 void GEOMImpl_IMeasureOperations::GetTolerance
1396                                (Handle(GEOM_Object) theShape,
1397                                 Standard_Real& FaceMin, Standard_Real& FaceMax,
1398                                 Standard_Real& EdgeMin, Standard_Real& EdgeMax,
1399                                 Standard_Real& VertMin, Standard_Real& VertMax)
1400 {
1401   SetErrorCode(KO);
1402
1403   if (theShape.IsNull()) return;
1404
1405   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1406   if (aRefShape.IsNull()) return;
1407
1408   TopoDS_Shape aShape = aRefShape->GetValue();
1409   if (aShape.IsNull()) {
1410     SetErrorCode("The Objects has NULL Shape");
1411     return;
1412   }
1413
1414   //Compute the parameters
1415   Standard_Real T;
1416   FaceMin = EdgeMin = VertMin = RealLast();
1417   FaceMax = EdgeMax = VertMax = -RealLast();
1418
1419   try {
1420     OCC_CATCH_SIGNALS;
1421     for (TopExp_Explorer ExF (aShape, TopAbs_FACE); ExF.More(); ExF.Next()) {
1422       TopoDS_Face Face = TopoDS::Face(ExF.Current());
1423       T = BRep_Tool::Tolerance(Face);
1424       if (T > FaceMax)
1425         FaceMax = T;
1426       if (T < FaceMin)
1427         FaceMin = T;
1428     }
1429     for (TopExp_Explorer ExE (aShape, TopAbs_EDGE); ExE.More(); ExE.Next()) {
1430       TopoDS_Edge Edge = TopoDS::Edge(ExE.Current());
1431       T = BRep_Tool::Tolerance(Edge);
1432       if (T > EdgeMax)
1433         EdgeMax = T;
1434       if (T < EdgeMin)
1435         EdgeMin = T;
1436     }
1437     for (TopExp_Explorer ExV (aShape, TopAbs_VERTEX); ExV.More(); ExV.Next()) {
1438       TopoDS_Vertex Vertex = TopoDS::Vertex(ExV.Current());
1439       T = BRep_Tool::Tolerance(Vertex);
1440       if (T > VertMax)
1441         VertMax = T;
1442       if (T < VertMin)
1443         VertMin = T;
1444     }
1445   }
1446   catch (Standard_Failure& aFail) {
1447     SetErrorCode(aFail.GetMessageString());
1448     return;
1449   }
1450
1451   SetErrorCode(OK);
1452 }
1453
1454 //=============================================================================
1455 /*!
1456  *  CheckShape
1457  */
1458 //=============================================================================
1459 bool GEOMImpl_IMeasureOperations::CheckShape (Handle(GEOM_Object)     theShape,
1460                                               const Standard_Boolean  theIsCheckGeom,
1461                                               std::list<ShapeError>  &theErrors)
1462 {
1463   SetErrorCode(KO);
1464   theErrors.clear();
1465
1466   if (theShape.IsNull()) return false;
1467
1468   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1469   if (aRefShape.IsNull()) return false;
1470
1471   TopoDS_Shape aShape = aRefShape->GetValue();
1472   if (aShape.IsNull()) {
1473     SetErrorCode("The Objects has NULL Shape");
1474     return false;
1475   }
1476
1477   //Compute the parameters
1478   bool isValid = false;
1479   try {
1480     OCC_CATCH_SIGNALS;
1481     BRepCheck_Analyzer ana (aShape, theIsCheckGeom);
1482     if (ana.IsValid()) {
1483       isValid = true;
1484     } else {
1485       FillErrors(ana, aShape, theErrors);
1486     }
1487   }
1488   catch (Standard_Failure& aFail) {
1489     SetErrorCode(aFail.GetMessageString());
1490     return false;
1491   }
1492
1493   SetErrorCode(OK);
1494   return isValid;
1495 }
1496
1497 //=============================================================================
1498 /*!
1499  *  PrintShapeErrors
1500  */
1501 //=============================================================================
1502 TCollection_AsciiString GEOMImpl_IMeasureOperations::PrintShapeErrors
1503                                  (Handle(GEOM_Object)          theShape,
1504                                   const std::list<ShapeError> &theErrors)
1505 {
1506   TCollection_AsciiString aDump;
1507
1508   if (theShape.IsNull()) {
1509     return aDump;
1510   }
1511
1512   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1513
1514   if (aRefShape.IsNull()) {
1515     return aDump;
1516   }
1517
1518   TopoDS_Shape aShape = aRefShape->GetValue();
1519
1520   if (aShape.IsNull()) {
1521     SetErrorCode("The Objects has NULL Shape");
1522     return aDump;
1523   }
1524
1525   if (!theErrors.empty()) {
1526     // The shape is not valid. Prepare errors for dump.
1527     TopTools_IndexedMapOfShape anIndices;
1528     std::list<ShapeError>::const_iterator anIter = theErrors.begin();
1529     Standard_Integer nbv, nbe, nbw, nbf, nbs, nbo;
1530     nbv = nbe = nbw = nbf = nbs = nbo = 0;
1531
1532     // Map sub-shapes and their indices
1533     TopExp::MapShapes(aShape, anIndices);
1534
1535     const Standard_Integer aNbSubShapes = anIndices.Extent();
1536     TColStd_MapOfInteger   aMapPbInd;
1537
1538     aDump += " -- The Shape has problems :\n";
1539     aDump += "  Check                                    Count\n";
1540     aDump += " ------------------------------------------------\n";
1541
1542     for (; anIter != theErrors.end(); anIter++) {
1543       Standard_Integer aNbShapes = anIter->incriminated.size();
1544
1545       switch(anIter->error) {
1546       case BRepCheck_InvalidPointOnCurve:
1547         aDump += "  Invalid Point on Curve ................... ";
1548         break;
1549       case BRepCheck_InvalidPointOnCurveOnSurface:
1550         aDump += "  Invalid Point on CurveOnSurface .......... ";
1551         break;
1552       case BRepCheck_InvalidPointOnSurface:
1553         aDump += "  Invalid Point on Surface ................. ";
1554         break;
1555       case BRepCheck_No3DCurve:
1556         aDump += "  No 3D Curve .............................. ";
1557         break;
1558       case BRepCheck_Multiple3DCurve:
1559         aDump += "  Multiple 3D Curve ........................ ";
1560         break;
1561       case BRepCheck_Invalid3DCurve:
1562         aDump += "  Invalid 3D Curve ......................... ";
1563         break;
1564       case BRepCheck_NoCurveOnSurface:
1565         aDump += "  No Curve on Surface ...................... ";
1566         break;
1567       case BRepCheck_InvalidCurveOnSurface:
1568         aDump += "  Invalid Curve on Surface ................. ";
1569         break;
1570       case BRepCheck_InvalidCurveOnClosedSurface:
1571         aDump += "  Invalid Curve on closed Surface .......... ";
1572         break;
1573       case BRepCheck_InvalidSameRangeFlag:
1574         aDump += "  Invalid SameRange Flag ................... ";
1575         break;
1576       case BRepCheck_InvalidSameParameterFlag:
1577         aDump += "  Invalid SameParameter Flag ............... ";
1578         break;
1579       case BRepCheck_InvalidDegeneratedFlag:
1580         aDump += "  Invalid Degenerated Flag ................. ";
1581         break;
1582       case BRepCheck_FreeEdge:
1583         aDump += "  Free Edge ................................ ";
1584         break;
1585       case BRepCheck_InvalidMultiConnexity:
1586         aDump += "  Invalid MultiConnexity ................... ";
1587         break;
1588       case BRepCheck_InvalidRange:
1589         aDump += "  Invalid Range ............................ ";
1590         break;
1591       case BRepCheck_EmptyWire:
1592         aDump += "  Empty Wire ............................... ";
1593         break;
1594       case BRepCheck_RedundantEdge:
1595         aDump += "  Redundant Edge ........................... ";
1596         break;
1597       case BRepCheck_SelfIntersectingWire:
1598         aDump += "  Self Intersecting Wire ................... ";
1599         break;
1600       case BRepCheck_NoSurface:
1601         aDump += "  No Surface ............................... ";
1602         break;
1603       case BRepCheck_InvalidWire:
1604         aDump += "  Invalid Wire ............................. ";
1605         break;
1606       case BRepCheck_RedundantWire:
1607         aDump += "  Redundant Wire ........................... ";
1608         break;
1609       case BRepCheck_IntersectingWires:
1610         aDump += "  Intersecting Wires ....................... ";
1611         break;
1612       case BRepCheck_InvalidImbricationOfWires:
1613         aDump += "  Invalid Imbrication of Wires ............. ";
1614         break;
1615       case BRepCheck_EmptyShell:
1616         aDump += "  Empty Shell .............................. ";
1617         break;
1618       case BRepCheck_RedundantFace:
1619         aDump += "  Redundant Face ........................... ";
1620         break;
1621       case BRepCheck_UnorientableShape:
1622         aDump += "  Unorientable Shape ....................... ";
1623         break;
1624       case BRepCheck_NotClosed:
1625         aDump += "  Not Closed ............................... ";
1626         break;
1627       case BRepCheck_NotConnected:
1628         aDump += "  Not Connected ............................ ";
1629         break;
1630       case BRepCheck_SubshapeNotInShape:
1631         aDump += "  Sub-shape not in Shape ................... ";
1632         break;
1633       case BRepCheck_BadOrientation:
1634         aDump += "  Bad Orientation .......................... ";
1635         break;
1636       case BRepCheck_BadOrientationOfSubshape:
1637         aDump += "  Bad Orientation of Sub-shape ............. ";
1638         break;
1639       case BRepCheck_InvalidToleranceValue:
1640         aDump += "  Invalid Tolerance Value .................. ";
1641         break;
1642       case BRepCheck_CheckFail:
1643         aDump += "  Check Shape Failure ...................... ";
1644         break;
1645       default:
1646         break;
1647       }
1648
1649       aDump += TCollection_AsciiString(aNbShapes) + "\n";
1650
1651       // Count types of shape.
1652       std::list<int>::const_iterator aShIter = anIter->incriminated.begin();
1653
1654       for (; aShIter != anIter->incriminated.end(); aShIter++) {
1655         const Standard_Integer anIndex = *aShIter;
1656
1657         if (anIndex > 0 && anIndex <= aNbSubShapes && aMapPbInd.Add(anIndex)) {
1658           const TopoDS_Shape     &aSubShape = anIndices.FindKey(anIndex);
1659           const TopAbs_ShapeEnum  aType     = aSubShape.ShapeType();
1660
1661           switch (aType) {
1662             case TopAbs_VERTEX : nbv++; break;
1663             case TopAbs_EDGE   : nbe++; break;
1664             case TopAbs_WIRE   : nbw++; break;
1665             case TopAbs_FACE   : nbf++; break;
1666             case TopAbs_SHELL  : nbs++; break;
1667             case TopAbs_SOLID  : nbo++; break;
1668             default            : break;
1669           }
1670         }
1671       }
1672     }
1673
1674     const Standard_Integer aNbFaultyShapes = nbv + nbe + nbw + nbf + nbs + nbo;
1675     aDump += " ------------------------------------------------\n";
1676     aDump += "*** Shapes with problems : ";
1677     aDump += TCollection_AsciiString(aNbFaultyShapes) + "\n";
1678
1679     if (nbv > 0) {
1680       aDump += "VERTEX : ";
1681       if (nbv < 10) aDump += " ";
1682       aDump += TCollection_AsciiString(nbv) + "\n";
1683     }
1684     if (nbe > 0) {
1685       aDump += "EDGE   : ";
1686       if (nbe < 10) aDump += " ";
1687       aDump += TCollection_AsciiString(nbe) + "\n";
1688     }
1689     if (nbw > 0) {
1690       aDump += "WIRE   : ";
1691       if (nbw < 10) aDump += " ";
1692       aDump += TCollection_AsciiString(nbw) + "\n";
1693     }
1694     if (nbf > 0) {
1695       aDump += "FACE   : ";
1696       if (nbf < 10) aDump += " ";
1697       aDump += TCollection_AsciiString(nbf) + "\n";
1698     }
1699     if (nbs > 0) {
1700       aDump += "SHELL  : ";
1701       if (nbs < 10) aDump += " ";
1702       aDump += TCollection_AsciiString(nbs) + "\n";
1703     }
1704     if (nbo > 0) {
1705       aDump += "SOLID  : ";
1706       if (nbo < 10) aDump += " ";
1707       aDump += TCollection_AsciiString(nbo) + "\n";
1708     }
1709   }
1710
1711   return aDump;
1712 }
1713
1714 //=============================================================================
1715 /*!
1716  *  CheckSelfIntersections
1717  */
1718 //=============================================================================
1719 bool GEOMImpl_IMeasureOperations::CheckSelfIntersections
1720                          (Handle(GEOM_Object)                 theShape,
1721                           const SICheckLevel                  theCheckLevel,
1722                           Handle(TColStd_HSequenceOfInteger)& theIntersections)
1723 {
1724   SetErrorCode(KO);
1725
1726   if (theIntersections.IsNull())
1727     theIntersections = new TColStd_HSequenceOfInteger;
1728   else
1729     theIntersections->Clear();
1730
1731   if (theShape.IsNull())
1732     return false;
1733
1734   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1735   if (aRefShape.IsNull()) return false;
1736
1737   TopoDS_Shape aShape = aRefShape->GetValue();
1738   if (aShape.IsNull()) return false;
1739
1740   // 0. Prepare data
1741   TopoDS_Shape aScopy;
1742   //
1743   GEOMAlgo_AlgoTools::CopyShape(aShape, aScopy);
1744   //
1745   // Map sub-shapes and their indices
1746   TopTools_IndexedMapOfShape anIndices;
1747   TopExp::MapShapes(aScopy, anIndices);
1748
1749   TopTools_ListOfShape aLCS;
1750   aLCS.Append(aScopy);
1751   //
1752   BOPAlgo_CheckerSI aCSI; // checker of self-interferences
1753   aCSI.SetArguments(aLCS);
1754   aCSI.SetLevelOfCheck(theCheckLevel);
1755
1756   // 1. Launch the checker
1757   aCSI.Perform();
1758   Standard_Boolean iErr = aCSI.HasErrors();
1759
1760   //
1761   Standard_Integer aNbS, n1, n2;
1762   BOPDS_MapIteratorOfMapOfPair aItMPK;
1763   //
1764   // 2. Take the shapes from DS
1765   const BOPDS_DS& aDS = aCSI.DS();
1766   aNbS=aDS.NbShapes();
1767   //
1768   // 3. Get the pairs of interfered shapes
1769   const BOPDS_MapOfPair& aMPK=aDS.Interferences();
1770   aItMPK.Initialize(aMPK);
1771   for (; aItMPK.More(); aItMPK.Next()) {
1772     const BOPDS_Pair& aPK=aItMPK.Value();
1773     aPK.Indices(n1, n2);
1774     //
1775     if (n1 > aNbS || n2 > aNbS){
1776       return false; // Error
1777     }
1778     const TopoDS_Shape& aS1 = aDS.Shape(n1);
1779     const TopoDS_Shape& aS2 = aDS.Shape(n2);
1780
1781     theIntersections->Append(anIndices.FindIndex(aS1));
1782     theIntersections->Append(anIndices.FindIndex(aS2));
1783   }
1784
1785   if (!iErr) {
1786     SetErrorCode(OK);
1787   }
1788
1789   return theIntersections->IsEmpty();
1790 }
1791
1792 //=============================================================================
1793 /*!
1794  *  CheckSelfIntersectionsFast
1795  */
1796 //=============================================================================
1797 bool GEOMImpl_IMeasureOperations::CheckSelfIntersectionsFast
1798                          (Handle(GEOM_Object) theShape,
1799                           float theDeflection, double theTolerance,
1800                           Handle(TColStd_HSequenceOfInteger)& theIntersections)
1801 {
1802   SetErrorCode(KO);
1803
1804   if (theIntersections.IsNull())
1805     theIntersections = new TColStd_HSequenceOfInteger;
1806   else
1807     theIntersections->Clear();
1808
1809   if (theShape.IsNull())
1810     return false;
1811
1812   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1813   if (aRefShape.IsNull()) return false;
1814
1815   TopoDS_Shape aShape = aRefShape->GetValue();
1816   if (aShape.IsNull()) return false;
1817
1818   // Prepare data
1819   TopoDS_Shape aScopy;
1820
1821   GEOMAlgo_AlgoTools::CopyShape(aShape, aScopy);
1822   GEOMUtils::MeshShape(aScopy, theDeflection);
1823
1824   // Map sub-shapes and their indices
1825   TopTools_IndexedMapOfShape anIndices;
1826   TopExp::MapShapes(aScopy, anIndices);
1827
1828   // Checker of fast interferences
1829   BRepExtrema_SelfIntersection aTool(aScopy, (theTolerance <= 0.) ? 0.0 : theTolerance);
1830
1831   // Launch the checker
1832   aTool.Perform();
1833
1834   const BRepExtrema_MapOfIntegerPackedMapOfInteger& intersections = aTool.OverlapElements();
1835
1836   std::set<Standard_Integer> processed;
1837
1838   for (BRepExtrema_MapOfIntegerPackedMapOfInteger::Iterator it(intersections); it.More(); it.Next()) {
1839     Standard_Integer idxLeft = it.Key();
1840     if (processed.count(idxLeft) > 0) continue; // already added
1841     processed.insert(idxLeft);
1842     const TColStd_PackedMapOfInteger& overlaps = it.Value();
1843     for (TColStd_MapIteratorOfPackedMapOfInteger subit(overlaps); subit.More(); subit.Next()) {
1844       Standard_Integer idxRight = subit.Key();
1845       if (processed.count(idxRight) > 0) continue; // already added
1846       const TopoDS_Shape& aS1 = aTool.GetSubShape(idxLeft);
1847       const TopoDS_Shape& aS2 = aTool.GetSubShape(idxRight);
1848       theIntersections->Append(anIndices.FindIndex(aS1));
1849       theIntersections->Append(anIndices.FindIndex(aS2));
1850     }
1851   }
1852
1853   if (aTool.IsDone())
1854     SetErrorCode(OK);
1855
1856   return theIntersections->IsEmpty();
1857 }
1858
1859 //=============================================================================
1860 /*!
1861  *  CheckBOPArguments
1862  */
1863 //=============================================================================
1864 bool GEOMImpl_IMeasureOperations::CheckBOPArguments
1865                                       (const Handle(GEOM_Object) &theShape)
1866 {
1867   SetErrorCode(KO);
1868
1869   if (theShape.IsNull()) {
1870     return false;
1871   }
1872
1873   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1874
1875   if (aRefShape.IsNull()) {
1876     return false;
1877   }
1878
1879   TopoDS_Shape aShape = aRefShape->GetValue();
1880
1881   if (aShape.IsNull()) {
1882     return false;
1883   }
1884
1885   //Compute the parameters
1886   bool isValid = GEOMUtils::CheckBOPArguments(aShape);
1887
1888   SetErrorCode(OK);
1889
1890   return isValid;
1891 }
1892
1893 //=============================================================================
1894 /*!
1895  *  FastIntersect
1896  */
1897 //=============================================================================
1898 bool GEOMImpl_IMeasureOperations::FastIntersect (Handle(GEOM_Object) theShape1, Handle(GEOM_Object) theShape2,
1899                                                  double theTolerance, float theDeflection,
1900                                                  Handle(TColStd_HSequenceOfInteger)& theIntersections1,
1901                                                  Handle(TColStd_HSequenceOfInteger)& theIntersections2)
1902 {
1903   SetErrorCode(KO);
1904   bool isGood = false;
1905
1906   if (theIntersections1.IsNull())
1907     theIntersections1 = new TColStd_HSequenceOfInteger;
1908   else
1909     theIntersections1->Clear();
1910
1911   if (theIntersections2.IsNull())
1912     theIntersections2 = new TColStd_HSequenceOfInteger;
1913   else
1914     theIntersections2->Clear();
1915
1916   if (theShape1.IsNull() || theShape2.IsNull()) {
1917     SetErrorCode("Objects have NULL Shape");
1918     return isGood;
1919   }
1920
1921   if (theShape1 == theShape2) {
1922     SetErrorCode("Objects are equal");
1923     return isGood;
1924   }
1925   Handle(GEOM_Function) aRefShape1 = theShape1->GetLastFunction();
1926   Handle(GEOM_Function) aRefShape2 = theShape2->GetLastFunction();
1927   if (aRefShape1.IsNull() || aRefShape2.IsNull()) return isGood;
1928
1929   TopoDS_Shape aShape1 = aRefShape1->GetValue();
1930   TopoDS_Shape aShape2 = aRefShape2->GetValue();
1931   if (aShape1.IsNull() || aShape2.IsNull()) return isGood;
1932
1933   // 0. Prepare data
1934   TopoDS_Shape aScopy1, aScopy2;
1935   GEOMAlgo_AlgoTools::CopyShape(aShape1, aScopy1);
1936   GEOMAlgo_AlgoTools::CopyShape(aShape2, aScopy2);
1937
1938   GEOMUtils::MeshShape(aScopy1, theDeflection);
1939   GEOMUtils::MeshShape(aScopy2, theDeflection);
1940   //
1941   // Map sub-shapes and their indices
1942   TopTools_IndexedMapOfShape anIndices1, anIndices2;
1943   TopExp::MapShapes(aScopy1, anIndices1);
1944   TopExp::MapShapes(aScopy2, anIndices2);
1945
1946   TopTools_ListOfShape aLCS1, aLCS2;
1947   aLCS1.Append(aScopy1); aLCS2.Append(aScopy2);
1948   //
1949   BRepExtrema_ShapeProximity aBSP; // checker of fast interferences
1950   aBSP.LoadShape1(aScopy1); aBSP.LoadShape2(aScopy2);
1951   aBSP.SetTolerance((theTolerance <= 0.) ? 0.0 : theTolerance);
1952
1953   // 1. Launch the checker
1954   aBSP.Perform();
1955
1956   // 2. Get sets of IDs of overlapped faces
1957   for (BRepExtrema_MapOfIntegerPackedMapOfInteger::Iterator anIt1 (aBSP.OverlapSubShapes1()); anIt1.More(); anIt1.Next())
1958   {
1959     const TopoDS_Shape& aS1 = aBSP.GetSubShape1(anIt1.Key());
1960     theIntersections1->Append(anIndices1.FindIndex(aS1));
1961   }
1962
1963   for (BRepExtrema_MapOfIntegerPackedMapOfInteger::Iterator anIt2 (aBSP.OverlapSubShapes2()); anIt2.More(); anIt2.Next())
1964   {
1965     const TopoDS_Shape& aS2 = aBSP.GetSubShape2(anIt2.Key());
1966     theIntersections2->Append(anIndices2.FindIndex(aS2));
1967   }
1968
1969   isGood = !theIntersections1->IsEmpty() && !theIntersections1->IsEmpty();
1970
1971   if (aBSP.IsDone())
1972     SetErrorCode(OK);
1973
1974   return isGood;
1975 }
1976
1977 //=============================================================================
1978 /*!
1979  *  IsGoodForSolid
1980  */
1981 //=============================================================================
1982 TCollection_AsciiString GEOMImpl_IMeasureOperations::IsGoodForSolid (Handle(GEOM_Object) theShape)
1983 {
1984   SetErrorCode(KO);
1985
1986   TCollection_AsciiString aRes = "";
1987
1988   if (theShape.IsNull()) {
1989     aRes = "WRN_NULL_OBJECT_OR_SHAPE";
1990   }
1991   else {
1992     Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1993     if (aRefShape.IsNull()) {
1994       aRes = "WRN_NULL_OBJECT_OR_SHAPE";
1995     }
1996     else {
1997       TopoDS_Shape aShape = aRefShape->GetValue();
1998       if (aShape.IsNull()) {
1999         aRes = "WRN_NULL_OBJECT_OR_SHAPE";
2000       }
2001       else {
2002         if (aShape.ShapeType() == TopAbs_COMPOUND) {
2003           TopoDS_Iterator It (aShape, Standard_True, Standard_True);
2004           if (It.More()) aShape = It.Value();
2005         }
2006         if (aShape.ShapeType() == TopAbs_SHELL) {
2007           BRepCheck_Shell chkShell (TopoDS::Shell(aShape));
2008           if (chkShell.Closed() == BRepCheck_NotClosed) {
2009             aRes = "WRN_SHAPE_UNCLOSED";
2010           }
2011         }
2012         else {
2013           aRes = "WRN_SHAPE_NOT_SHELL";
2014         }
2015       }
2016     }
2017   }
2018
2019   if (aRes.IsEmpty())
2020     SetErrorCode(OK);
2021
2022   return aRes;
2023 }
2024
2025 //=============================================================================
2026 /*!
2027  *  WhatIs
2028  */
2029 //=============================================================================
2030 TCollection_AsciiString GEOMImpl_IMeasureOperations::WhatIs (Handle(GEOM_Object) theShape)
2031 {
2032   SetErrorCode(KO);
2033
2034   TCollection_AsciiString Astr;
2035
2036   if (theShape.IsNull()) return Astr;
2037
2038   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
2039   if (aRefShape.IsNull()) return Astr;
2040
2041   TopoDS_Shape aShape = aRefShape->GetValue();
2042   if (aShape.IsNull()) {
2043     SetErrorCode("The Objects has NULL Shape");
2044     return Astr;
2045   }
2046
2047   //Compute the parameters
2048   if (aShape.ShapeType() == TopAbs_EDGE) {
2049     if (BRep_Tool::Degenerated(TopoDS::Edge(aShape))) {
2050       Astr = Astr + " It is a degenerated edge \n";
2051     }
2052   }
2053
2054   Astr = Astr + " Number of sub-shapes : \n";
2055
2056   try {
2057     OCC_CATCH_SIGNALS;
2058     int iType, nbTypes [TopAbs_SHAPE], nbFlatType [TopAbs_SHAPE];
2059     for (iType = 0; iType < TopAbs_SHAPE; ++iType) {
2060       nbTypes[iType] = 0;
2061       nbFlatType[iType] = 0;
2062     }
2063     nbTypes[aShape.ShapeType()]++;
2064
2065     TopTools_MapOfShape aMapOfShape;
2066     aMapOfShape.Add(aShape);
2067     TopTools_ListOfShape aListOfShape;
2068     aListOfShape.Append(aShape);
2069
2070     TopTools_ListIteratorOfListOfShape itL (aListOfShape);
2071     for (; itL.More(); itL.Next()) {
2072       TopoDS_Shape sp = itL.Value();
2073       TopoDS_Iterator it (sp);
2074       for (; it.More(); it.Next()) {
2075         TopoDS_Shape s = it.Value();
2076         if (aMapOfShape.Add(s)) {
2077           aListOfShape.Append(s);
2078           nbTypes[s.ShapeType()]++;
2079           if ((sp.ShapeType() == TopAbs_COMPOUND) || (sp.ShapeType() == TopAbs_COMPSOLID)) {
2080             nbFlatType[s.ShapeType()]++;
2081           }
2082         }
2083       }
2084     }
2085
2086     Astr = Astr + " VERTEX : " + TCollection_AsciiString(nbTypes[TopAbs_VERTEX]) + "\n";
2087     Astr = Astr + " EDGE : " + TCollection_AsciiString(nbTypes[TopAbs_EDGE]) + "\n";
2088     Astr = Astr + " WIRE : " + TCollection_AsciiString(nbTypes[TopAbs_WIRE]) + "\n";
2089     Astr = Astr + " FACE : " + TCollection_AsciiString(nbTypes[TopAbs_FACE]) + "\n";
2090     Astr = Astr + " SHELL : " + TCollection_AsciiString(nbTypes[TopAbs_SHELL]) + "\n";
2091     Astr = Astr + " SOLID : " + TCollection_AsciiString(nbTypes[TopAbs_SOLID]) + "\n";
2092     Astr = Astr + " COMPSOLID : " + TCollection_AsciiString(nbTypes[TopAbs_COMPSOLID]) + "\n";
2093     Astr = Astr + " COMPOUND : " + TCollection_AsciiString(nbTypes[TopAbs_COMPOUND]) + "\n";
2094     Astr = Astr + " SHAPE : " + TCollection_AsciiString(aMapOfShape.Extent()) + "\n";
2095
2096     if ((aShape.ShapeType() == TopAbs_COMPOUND) || (aShape.ShapeType() == TopAbs_COMPSOLID)){
2097       Astr = Astr + " --------------------- \n Flat content : \n";
2098       if (nbFlatType[TopAbs_VERTEX] > 0)
2099         Astr = Astr + " VERTEX : " + TCollection_AsciiString(nbFlatType[TopAbs_VERTEX]) + "\n";
2100       if (nbFlatType[TopAbs_EDGE] > 0)
2101         Astr = Astr + " EDGE : " + TCollection_AsciiString(nbFlatType[TopAbs_EDGE]) + "\n";
2102       if (nbFlatType[TopAbs_WIRE] > 0)
2103         Astr = Astr + " WIRE : " + TCollection_AsciiString(nbFlatType[TopAbs_WIRE]) + "\n";
2104       if (nbFlatType[TopAbs_FACE] > 0)
2105         Astr = Astr + " FACE : " + TCollection_AsciiString(nbFlatType[TopAbs_FACE]) + "\n";
2106       if (nbFlatType[TopAbs_SHELL] > 0)
2107         Astr = Astr + " SHELL : " + TCollection_AsciiString(nbFlatType[TopAbs_SHELL]) + "\n";
2108       if (nbFlatType[TopAbs_SOLID] > 0)
2109         Astr = Astr + " SOLID : " + TCollection_AsciiString(nbFlatType[TopAbs_SOLID]) + "\n";
2110     }
2111   }
2112   catch (Standard_Failure& aFail) {
2113     SetErrorCode(aFail.GetMessageString());
2114     return Astr;
2115   }
2116
2117   SetErrorCode(OK);
2118   return Astr;
2119 }
2120
2121 //=============================================================================
2122 /*!
2123  *  AreCoordsInside
2124  */
2125 //=============================================================================
2126 std::vector<bool>
2127 GEOMImpl_IMeasureOperations::AreCoordsInside(Handle(GEOM_Object)        theShape,
2128                                              const std::vector<double>& coords,
2129                                              double                     tolerance)
2130 {
2131   std::vector<bool> isInsideRes;
2132   if (!theShape.IsNull()) {
2133     Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
2134     if (!aRefShape.IsNull()) {
2135       TopoDS_Shape aShape = aRefShape->GetValue();
2136       if (!aShape.IsNull())
2137       {
2138         TopTools_IndexedMapOfShape mapShape;
2139         {
2140           TopExp_Explorer anExp;
2141           for ( anExp.Init( aShape, TopAbs_SOLID ); anExp.More(); anExp.Next() )
2142             mapShape.Add( anExp.Current() );
2143           for ( anExp.Init( aShape, TopAbs_FACE, TopAbs_SOLID ); anExp.More(); anExp.Next() )
2144             mapShape.Add( anExp.Current() );
2145           for ( anExp.Init( aShape, TopAbs_EDGE, TopAbs_FACE ); anExp.More(); anExp.Next() )
2146             mapShape.Add( anExp.Current() );
2147           for ( anExp.Init( aShape, TopAbs_VERTEX, TopAbs_EDGE ); anExp.More(); anExp.Next() )
2148             mapShape.Add( anExp.Current() ); //// ?????????
2149         }
2150         size_t nb_points = coords.size()/3, nb_points_inside = 0;
2151         isInsideRes.resize( nb_points, false );
2152
2153         for ( int iS = 1; iS <= mapShape.Extent(); ++iS )
2154         {
2155           if ( nb_points_inside == nb_points )
2156             break;
2157           aShape = mapShape( iS );
2158           switch ( aShape.ShapeType() ) {
2159           case TopAbs_SOLID:
2160           {
2161             BRepClass3d_SolidClassifier SC( TopoDS::Solid( aShape ));
2162             for ( size_t i = 0; i < nb_points; i++)
2163             {
2164               if ( isInsideRes[ i ]) continue;
2165               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
2166               SC.Perform( aPnt, tolerance );
2167               isInsideRes[ i ] = (( SC.State() == TopAbs_IN ) || ( SC.State() == TopAbs_ON ));
2168               nb_points_inside += isInsideRes[ i ];
2169             }
2170             break;
2171           }
2172           case TopAbs_FACE:
2173           {
2174             Standard_Real u1,u2,v1,v2;
2175             const TopoDS_Face&   face = TopoDS::Face( aShape );
2176             Handle(Geom_Surface) surf = BRep_Tool::Surface( face );
2177             surf->Bounds( u1,u2,v1,v2 );
2178             GeomAPI_ProjectPointOnSurf project;
2179             project.Init(surf, u1,u2, v1,v2, tolerance );
2180             for ( size_t i = 0; i < nb_points; i++)
2181             {
2182               if ( isInsideRes[ i ]) continue;
2183               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
2184               project.Perform( aPnt );
2185               if ( project.IsDone() &&
2186                    project.NbPoints() > 0 &&
2187                    project.LowerDistance() <= tolerance )
2188               {
2189                 Standard_Real u, v;
2190                 project.LowerDistanceParameters(u, v);
2191                 gp_Pnt2d uv( u, v );
2192                 BRepClass_FaceClassifier FC ( face, uv, tolerance );
2193                 isInsideRes[ i ] = (( FC.State() == TopAbs_IN ) || ( FC.State() == TopAbs_ON ));
2194                 nb_points_inside += isInsideRes[ i ];
2195               }
2196             }
2197             break;
2198           }
2199           case TopAbs_EDGE:
2200           {
2201             Standard_Real f,l;
2202             const TopoDS_Edge&  edge = TopoDS::Edge( aShape );
2203             Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, f, l );
2204             GeomAPI_ProjectPointOnCurve project;
2205             project.Init( curve, f, l );
2206             for ( size_t i = 0; i < nb_points; i++)
2207             {
2208               if ( isInsideRes[ i ]) continue;
2209               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
2210               project.Perform( aPnt );
2211               isInsideRes[ i ] = ( project.NbPoints() > 0 &&
2212                                    project.LowerDistance() <= tolerance );
2213               nb_points_inside += isInsideRes[ i ];
2214             }
2215             break;
2216           }
2217           case TopAbs_VERTEX:
2218           {
2219             gp_Pnt aVPnt = BRep_Tool::Pnt( TopoDS::Vertex( aShape ));
2220             for ( size_t i = 0; i < nb_points; i++)
2221             {
2222               if ( isInsideRes[ i ]) continue;
2223               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
2224               isInsideRes[ i ] = ( aPnt.SquareDistance( aVPnt ) <= tolerance * tolerance );
2225               nb_points_inside += isInsideRes[ i ];
2226             }
2227             break;
2228           }
2229           default:;
2230           } // switch ( aShape.ShapeType() )
2231         }
2232       }
2233     }
2234   }
2235   return isInsideRes;
2236 }
2237
2238 //=============================================================================
2239 /*!
2240  *  GetMinDistance
2241  */
2242 //=============================================================================
2243 Standard_Real
2244 GEOMImpl_IMeasureOperations::GetMinDistance (Handle(GEOM_Object) theShape1,
2245                                              Handle(GEOM_Object) theShape2,
2246                                              Standard_Real& X1,
2247                                              Standard_Real& Y1,
2248                                              Standard_Real& Z1,
2249                                              Standard_Real& X2,
2250                                              Standard_Real& Y2,
2251                                              Standard_Real& Z2)
2252 {
2253   SetErrorCode(KO);
2254   Standard_Real MinDist = 1.e9;
2255
2256   if (theShape1.IsNull() || theShape2.IsNull()) return MinDist;
2257
2258   Handle(GEOM_Function) aRefShape1 = theShape1->GetLastFunction();
2259   Handle(GEOM_Function) aRefShape2 = theShape2->GetLastFunction();
2260   if (aRefShape1.IsNull() || aRefShape2.IsNull()) return MinDist;
2261
2262   TopoDS_Shape aShape1 = aRefShape1->GetValue();
2263   TopoDS_Shape aShape2 = aRefShape2->GetValue();
2264   if (aShape1.IsNull() || aShape2.IsNull()) {
2265     SetErrorCode("One of Objects has NULL Shape");
2266     return MinDist;
2267   }
2268
2269   //Compute the parameters
2270   try {
2271     OCC_CATCH_SIGNALS;
2272
2273     gp_Pnt aPnt1, aPnt2;
2274
2275     MinDist = GEOMUtils::GetMinDistance(aShape1, aShape2, aPnt1, aPnt2);
2276
2277     if (MinDist >= 0.0) {
2278       aPnt1.Coord(X1, Y1, Z1);
2279       aPnt2.Coord(X2, Y2, Z2);
2280     } else {
2281       return MinDist;
2282     }
2283   }
2284   catch (Standard_Failure& aFail) {
2285     SetErrorCode(aFail.GetMessageString());
2286     return MinDist;
2287   }
2288
2289   SetErrorCode(OK);
2290   return MinDist;
2291 }
2292
2293 //=======================================================================
2294 /*!
2295  *  Get coordinates of closest points of two shapes
2296  */
2297 //=======================================================================
2298 Standard_Integer GEOMImpl_IMeasureOperations::ClosestPoints (Handle(GEOM_Object) theShape1,
2299                                                              Handle(GEOM_Object) theShape2,
2300                                                              Handle(TColStd_HSequenceOfReal)& theDoubles)
2301 {
2302   SetErrorCode(KO);
2303   Standard_Integer nbSolutions = 0;
2304
2305   if (theShape1.IsNull() || theShape2.IsNull()) return nbSolutions;
2306
2307   Handle(GEOM_Function) aRefShape1 = theShape1->GetLastFunction();
2308   Handle(GEOM_Function) aRefShape2 = theShape2->GetLastFunction();
2309   if (aRefShape1.IsNull() || aRefShape2.IsNull()) return nbSolutions;
2310
2311   TopoDS_Shape aShape1 = aRefShape1->GetValue();
2312   TopoDS_Shape aShape2 = aRefShape2->GetValue();
2313   if (aShape1.IsNull() || aShape2.IsNull()) {
2314     SetErrorCode("One of Objects has NULL Shape");
2315     return nbSolutions;
2316   }
2317
2318   // Compute the extremities
2319   try {
2320     OCC_CATCH_SIGNALS;
2321
2322     // skl 30.06.2008
2323     // additional workaround for bugs 19899, 19908 and 19910 from Mantis
2324     gp_Pnt P1s, P2s;
2325     double dist = GEOMUtils::GetMinDistanceSingular(aShape1, aShape2, P1s, P2s);
2326     bool singularBetter = dist >= 0;
2327
2328     BRepExtrema_DistShapeShape dst (aShape1, aShape2);
2329     if (dst.IsDone()) {
2330       nbSolutions = dst.NbSolution();
2331       if (theDoubles.IsNull()) theDoubles = new TColStd_HSequenceOfReal;
2332
2333       gp_Pnt P1, P2;
2334       for (int i = 1; i <= nbSolutions; i++) {
2335         P1 = dst.PointOnShape1(i);
2336         P2 = dst.PointOnShape2(i);
2337
2338         theDoubles->Append(P1.X());
2339         theDoubles->Append(P1.Y());
2340         theDoubles->Append(P1.Z());
2341         theDoubles->Append(P2.X());
2342         theDoubles->Append(P2.Y());
2343         theDoubles->Append(P2.Z());
2344
2345         Standard_Real Dist = P1.Distance(P2);
2346         singularBetter = singularBetter && dist < Dist;
2347       }
2348     }
2349
2350     if (singularBetter) {
2351       if (theDoubles.IsNull()) theDoubles = new TColStd_HSequenceOfReal;
2352       else theDoubles->Clear();
2353
2354       nbSolutions = 1;
2355
2356       theDoubles->Append(P1s.X());
2357       theDoubles->Append(P1s.Y());
2358       theDoubles->Append(P1s.Z());
2359       theDoubles->Append(P2s.X());
2360       theDoubles->Append(P2s.Y());
2361       theDoubles->Append(P2s.Z());
2362     }
2363   }
2364   catch (Standard_Failure& aFail) {
2365     SetErrorCode(aFail.GetMessageString());
2366     return nbSolutions;
2367   }
2368
2369   SetErrorCode(OK);
2370   return nbSolutions;
2371 }
2372
2373 //=======================================================================
2374 /*!
2375  *  Get coordinates of point
2376  */
2377 //=======================================================================
2378 void GEOMImpl_IMeasureOperations::PointCoordinates (Handle(GEOM_Object) theShape,
2379                         Standard_Real& theX, Standard_Real& theY, Standard_Real& theZ)
2380 {
2381   SetErrorCode(KO);
2382
2383   if (theShape.IsNull())
2384     return;
2385
2386   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
2387   if (aRefShape.IsNull())
2388     return;
2389
2390   TopoDS_Shape aShape = aRefShape->GetValue();
2391   if (aShape.IsNull() || aShape.ShapeType() != TopAbs_VERTEX)
2392   {
2393     SetErrorCode( "Shape must be a vertex" );
2394     return;
2395   }
2396
2397   try {
2398     OCC_CATCH_SIGNALS;
2399     gp_Pnt aPnt = BRep_Tool::Pnt( TopoDS::Vertex( aShape ) );
2400     theX = aPnt.X();
2401     theY = aPnt.Y();
2402     theZ = aPnt.Z();
2403
2404     SetErrorCode(OK);
2405   }
2406   catch (Standard_Failure& aFail)
2407   {
2408     SetErrorCode( aFail.GetMessageString() );
2409   }
2410 }
2411
2412 //=======================================================================
2413 /*!
2414  *  Compute angle (in degrees) between two lines
2415  */
2416 //=======================================================================
2417 Standard_Real GEOMImpl_IMeasureOperations::GetAngle (Handle(GEOM_Object) theLine1,
2418                                                      Handle(GEOM_Object) theLine2)
2419 {
2420   if (theLine1->GetType() == GEOM_VECTOR &&
2421       theLine2->GetType() == GEOM_VECTOR)
2422     return GetAngleBtwVectors(theLine1, theLine2);
2423
2424   SetErrorCode(KO);
2425
2426   Standard_Real anAngle = -1.0;
2427
2428   if (theLine1.IsNull() || theLine2.IsNull())
2429     return anAngle;
2430
2431   Handle(GEOM_Function) aRefLine1 = theLine1->GetLastFunction();
2432   Handle(GEOM_Function) aRefLine2 = theLine2->GetLastFunction();
2433   if (aRefLine1.IsNull() || aRefLine2.IsNull())
2434     return anAngle;
2435
2436   TopoDS_Shape aLine1 = aRefLine1->GetValue();
2437   TopoDS_Shape aLine2 = aRefLine2->GetValue();
2438   if (aLine1.IsNull() || aLine2.IsNull() ||
2439       aLine1.ShapeType() != TopAbs_EDGE ||
2440       aLine2.ShapeType() != TopAbs_EDGE)
2441   {
2442     SetErrorCode("Two edges must be given");
2443     return anAngle;
2444   }
2445
2446   try {
2447     OCC_CATCH_SIGNALS;
2448     TopoDS_Edge E1 = TopoDS::Edge(aLine1);
2449     TopoDS_Edge E2 = TopoDS::Edge(aLine2);
2450
2451     double fp,lp;
2452     Handle(Geom_Curve) C1 = BRep_Tool::Curve(E1,fp,lp);
2453     Handle(Geom_Curve) C2 = BRep_Tool::Curve(E2,fp,lp);
2454
2455     if ( C1.IsNull() || C2.IsNull() ||
2456         !C1->IsKind(STANDARD_TYPE(Geom_Line)) ||
2457         !C2->IsKind(STANDARD_TYPE(Geom_Line)))
2458     {
2459       SetErrorCode("The edges must be linear");
2460       return anAngle;
2461     }
2462
2463     Handle(Geom_Line) L1 = Handle(Geom_Line)::DownCast(C1);
2464     Handle(Geom_Line) L2 = Handle(Geom_Line)::DownCast(C2);
2465
2466     gp_Lin aLin1 = L1->Lin();
2467     gp_Lin aLin2 = L2->Lin();
2468
2469     anAngle = aLin1.Angle(aLin2);
2470     anAngle *= 180. / M_PI; // convert radians into degrees
2471
2472     if (anAngle > 90.0) {
2473       anAngle = 180.0 - anAngle;
2474     }
2475
2476     SetErrorCode(OK);
2477   }
2478   catch (Standard_Failure& aFail)
2479   {
2480     SetErrorCode(aFail.GetMessageString());
2481   }
2482
2483   return anAngle;
2484 }
2485
2486 //=======================================================================
2487 /*!
2488  *  Compute angle (in degrees) between two vectors
2489  */
2490 //=======================================================================
2491 Standard_Real GEOMImpl_IMeasureOperations::GetAngleBtwVectors (Handle(GEOM_Object) theVec1,
2492                                                                Handle(GEOM_Object) theVec2)
2493 {
2494   SetErrorCode(KO);
2495
2496   Standard_Real anAngle = -1.0;
2497
2498   if (theVec1.IsNull() || theVec2.IsNull())
2499     return anAngle;
2500
2501   if (theVec1->GetType() != GEOM_VECTOR || theVec2->GetType() != GEOM_VECTOR) {
2502     SetErrorCode("Two vectors must be given");
2503     return anAngle;
2504   }
2505
2506   Handle(GEOM_Function) aRefVec1 = theVec1->GetLastFunction();
2507   Handle(GEOM_Function) aRefVec2 = theVec2->GetLastFunction();
2508   if (aRefVec1.IsNull() || aRefVec2.IsNull())
2509     return anAngle;
2510
2511   TopoDS_Shape aVec1 = aRefVec1->GetValue();
2512   TopoDS_Shape aVec2 = aRefVec2->GetValue();
2513   if (aVec1.IsNull() || aVec2.IsNull() ||
2514       aVec1.ShapeType() != TopAbs_EDGE ||
2515       aVec2.ShapeType() != TopAbs_EDGE)
2516   {
2517     SetErrorCode("Two edges must be given");
2518     return anAngle;
2519   }
2520
2521   try {
2522     OCC_CATCH_SIGNALS;
2523     TopoDS_Edge aE1 = TopoDS::Edge(aVec1);
2524     TopoDS_Edge aE2 = TopoDS::Edge(aVec2);
2525
2526     TopoDS_Vertex aP11, aP12, aP21, aP22;
2527     TopExp::Vertices(aE1, aP11, aP12, Standard_True);
2528     TopExp::Vertices(aE2, aP21, aP22, Standard_True);
2529     if (aP11.IsNull() || aP12.IsNull() || aP21.IsNull() || aP22.IsNull()) {
2530       SetErrorCode("Bad edge given");
2531       return anAngle;
2532     }
2533
2534     gp_Vec aV1 (BRep_Tool::Pnt(aP11), BRep_Tool::Pnt(aP12));
2535     gp_Vec aV2 (BRep_Tool::Pnt(aP21), BRep_Tool::Pnt(aP22)) ;
2536
2537     anAngle = aV1.Angle(aV2);
2538     anAngle *= 180. / M_PI; // convert radians into degrees
2539
2540     SetErrorCode(OK);
2541   }
2542   catch (Standard_Failure& aFail)
2543   {
2544     SetErrorCode(aFail.GetMessageString());
2545   }
2546
2547   return anAngle;
2548 }
2549
2550
2551 //=============================================================================
2552 /*!
2553  *  PatchFace
2554  */
2555  //=============================================================================
2556 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IMeasureOperations::PatchFace(Handle(GEOM_Object) theShape)
2557 {
2558   SetErrorCode(KO);
2559
2560   if (theShape.IsNull()) return NULL;
2561
2562   Handle(GEOM_Object) aPatchFace = GetEngine()->AddObject(GEOM_PATCH_FACE);
2563   Handle(GEOM_Function) aFunction = aPatchFace->AddFunction(GEOMImpl_PatchFaceDriver::GetID(), 1);
2564   if (aFunction.IsNull()) return NULL;
2565
2566   //Check if the function is set correctly
2567   if (aFunction->GetDriverGUID() != GEOMImpl_PatchFaceDriver::GetID()) return NULL;
2568
2569   GEOMImpl_IPatchFace aPI(aFunction);
2570   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
2571   if (aRefShape.IsNull()) return NULL;
2572
2573   aPI.SetShape(aRefShape);
2574   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2575
2576   // Perform
2577   try
2578   {
2579     OCC_CATCH_SIGNALS;
2580     if (!GetSolver()->ComputeFunction(aFunction))
2581     {
2582       SetErrorCode("patch face driver failed");
2583       return NULL;
2584     }
2585
2586     // Get result compound and collect all faces into result sequence
2587     TopoDS_Shape aResCompound = aFunction->GetValue();
2588     TopTools_IndexedMapOfShape anIndices;
2589     TopExp::MapShapes(aResCompound, anIndices);
2590
2591     Handle(TColStd_HArray1OfInteger) anArray;
2592     for (TopExp_Explorer anExpW(aResCompound, TopAbs_FACE); anExpW.More(); anExpW.Next())
2593     {
2594       TopoDS_Shape aValue = anExpW.Value();
2595       anArray = new TColStd_HArray1OfInteger(1, 1);
2596       anArray->SetValue(1, anIndices.FindIndex(aValue));
2597
2598       Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(aPatchFace, anArray);
2599       if (!anObj.IsNull())
2600       {
2601         aSeq->Append(anObj);
2602       }
2603     }
2604   }
2605   catch (Standard_Failure& aFail)
2606   {
2607     SetErrorCode(aFail.GetMessageString());
2608     return aSeq;
2609   }
2610
2611   //Make a Python command
2612   GEOM::TPythonDump(aFunction, true)
2613     << "[" << aSeq << "] = geompy.PatchFace(" << theShape << ")";
2614
2615   SetErrorCode(OK);
2616   return aSeq;
2617 }
2618
2619 //=============================================================================
2620 /*!
2621  *  CurveCurvatureByParam
2622  */
2623 //=============================================================================
2624 Standard_Real GEOMImpl_IMeasureOperations::CurveCurvatureByParam
2625                         (Handle(GEOM_Object) theCurve, Standard_Real& theParam)
2626 {
2627   SetErrorCode(KO);
2628   Standard_Real aRes = -1.0;
2629
2630   if(theCurve.IsNull()) return aRes;
2631
2632   Handle(GEOM_Function) aRefShape = theCurve->GetLastFunction();
2633   if(aRefShape.IsNull()) return aRes;
2634
2635   TopoDS_Shape aShape = aRefShape->GetValue();
2636   if(aShape.IsNull()) {
2637     SetErrorCode("One of Objects has NULL Shape");
2638     return aRes;
2639   }
2640
2641   Standard_Real aFP, aLP, aP;
2642   Handle(Geom_Curve) aCurve = BRep_Tool::Curve(TopoDS::Edge(aShape), aFP, aLP);
2643   aP = aFP + (aLP - aFP) * theParam;
2644
2645   if(aCurve.IsNull()) return aRes;
2646
2647   //Compute curvature
2648   try {
2649     OCC_CATCH_SIGNALS;
2650     GeomLProp_CLProps Prop = GeomLProp_CLProps
2651       (aCurve, aP, 2, Precision::Confusion());
2652     aRes = fabs(Prop.Curvature());
2653     SetErrorCode(OK);
2654   }
2655   catch (Standard_Failure& aFail) {
2656     SetErrorCode(aFail.GetMessageString());
2657     return aRes;
2658   }
2659
2660   if( aRes > Precision::Confusion() )
2661     aRes = 1/aRes;
2662   else
2663     aRes = RealLast();
2664
2665   return aRes;
2666 }
2667
2668
2669 //=============================================================================
2670 /*!
2671  *  CurveCurvatureByPoint
2672  */
2673 //=============================================================================
2674 Standard_Real GEOMImpl_IMeasureOperations::CurveCurvatureByPoint
2675                    (Handle(GEOM_Object) theCurve, Handle(GEOM_Object) thePoint)
2676 {
2677   SetErrorCode(KO);
2678   Standard_Real aRes = -1.0;
2679
2680   if( theCurve.IsNull() || thePoint.IsNull() ) return aRes;
2681
2682   Handle(GEOM_Function) aRefCurve = theCurve->GetLastFunction();
2683   Handle(GEOM_Function) aRefPoint = thePoint->GetLastFunction();
2684   if( aRefCurve.IsNull() || aRefPoint.IsNull() ) return aRes;
2685
2686   TopoDS_Edge anEdge = TopoDS::Edge(aRefCurve->GetValue());
2687   TopoDS_Vertex aPnt = TopoDS::Vertex(aRefPoint->GetValue());
2688   if( anEdge.IsNull() || aPnt.IsNull() ) {
2689     SetErrorCode("One of Objects has NULL Shape");
2690     return aRes;
2691   }
2692
2693   Standard_Real aFP, aLP;
2694   Handle(Geom_Curve) aCurve = BRep_Tool::Curve(anEdge, aFP, aLP);
2695   if(aCurve.IsNull()) return aRes;
2696   gp_Pnt aPoint = BRep_Tool::Pnt(aPnt);
2697
2698   //Compute curvature
2699   try {
2700     OCC_CATCH_SIGNALS;
2701     GeomAPI_ProjectPointOnCurve PPCurve(aPoint, aCurve, aFP, aLP);
2702     if(PPCurve.NbPoints()>0) {
2703       GeomLProp_CLProps Prop = GeomLProp_CLProps
2704         (aCurve, PPCurve.LowerDistanceParameter(), 2, Precision::Confusion());
2705       aRes = fabs(Prop.Curvature());
2706       SetErrorCode(OK);
2707     }
2708   }
2709   catch (Standard_Failure& aFail) {
2710     SetErrorCode(aFail.GetMessageString());
2711     return aRes;
2712   }
2713
2714   if( aRes > Precision::Confusion() )
2715     aRes = 1/aRes;
2716   else
2717     aRes = RealLast();
2718
2719   return aRes;
2720 }
2721
2722
2723 //=============================================================================
2724 /*!
2725  *  getSurfaceCurvatures
2726  */
2727 //=============================================================================
2728 Standard_Real GEOMImpl_IMeasureOperations::getSurfaceCurvatures
2729                                           (const Handle(Geom_Surface)& aSurf,
2730                                            Standard_Real theUParam,
2731                                            Standard_Real theVParam,
2732                                            Standard_Boolean theNeedMaxCurv)
2733 {
2734   SetErrorCode(KO);
2735   Standard_Real aRes = 1.0;
2736
2737   if (aSurf.IsNull()) return aRes;
2738
2739   try {
2740     OCC_CATCH_SIGNALS;
2741     GeomLProp_SLProps Prop = GeomLProp_SLProps
2742       (aSurf, theUParam, theVParam, 2, Precision::Confusion());
2743     if(Prop.IsCurvatureDefined()) {
2744       if(Prop.IsUmbilic()) {
2745         //cout<<"is umbilic"<<endl;
2746         aRes = fabs(Prop.MeanCurvature());
2747       }
2748       else {
2749         //cout<<"is not umbilic"<<endl;
2750         double c1 = fabs(Prop.MaxCurvature());
2751         double c2 = fabs(Prop.MinCurvature());
2752         if(theNeedMaxCurv)
2753           aRes = Max(c1,c2);
2754         else
2755           aRes = Min(c1,c2);
2756       }
2757       SetErrorCode(OK);
2758     }
2759   }
2760   catch (Standard_Failure& aFail) {
2761     SetErrorCode(aFail.GetMessageString());
2762     return aRes;
2763   }
2764
2765   if( fabs(aRes) > Precision::Confusion() )
2766     aRes = 1/aRes;
2767   else
2768     aRes = RealLast();
2769
2770   return aRes;
2771 }
2772
2773
2774 //=============================================================================
2775 /*!
2776  *  MaxSurfaceCurvatureByParam
2777  */
2778 //=============================================================================
2779 Standard_Real GEOMImpl_IMeasureOperations::MaxSurfaceCurvatureByParam
2780                                                   (Handle(GEOM_Object) theSurf,
2781                                                    Standard_Real& theUParam,
2782                                                    Standard_Real& theVParam)
2783 {
2784   SetErrorCode(KO);
2785   Standard_Real aRes = -1.0;
2786
2787   if (theSurf.IsNull()) return aRes;
2788
2789   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2790   if(aRefShape.IsNull()) return aRes;
2791
2792   TopoDS_Shape aShape = aRefShape->GetValue();
2793   if(aShape.IsNull()) {
2794     SetErrorCode("One of Objects has NULL Shape");
2795     return aRes;
2796   }
2797
2798   TopoDS_Face F = TopoDS::Face(aShape);
2799   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(F);
2800
2801   //Compute the parameters
2802   Standard_Real U1,U2,V1,V2;
2803   BRepTools::UVBounds(F,U1,U2,V1,V2);
2804   Standard_Real U = U1 + (U2-U1)*theUParam;
2805   Standard_Real V = V1 + (V2-V1)*theVParam;
2806
2807   return getSurfaceCurvatures(aSurf, U, V, true);
2808 }
2809
2810
2811 //=============================================================================
2812 /*!
2813  *  MaxSurfaceCurvatureByPoint
2814  */
2815 //=============================================================================
2816 Standard_Real GEOMImpl_IMeasureOperations::MaxSurfaceCurvatureByPoint
2817                     (Handle(GEOM_Object) theSurf, Handle(GEOM_Object) thePoint)
2818 {
2819   SetErrorCode(KO);
2820   Standard_Real aRes = -1.0;
2821
2822   if( theSurf.IsNull() || thePoint.IsNull() ) return aRes;
2823
2824   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2825   Handle(GEOM_Function) aRefPoint = thePoint->GetLastFunction();
2826   if( aRefShape.IsNull() || aRefPoint.IsNull() ) return aRes;
2827
2828   TopoDS_Face aFace = TopoDS::Face(aRefShape->GetValue());
2829   TopoDS_Vertex aPnt = TopoDS::Vertex(aRefPoint->GetValue());
2830   if( aFace.IsNull() || aPnt.IsNull() ) {
2831     SetErrorCode("One of Objects has NULL Shape");
2832     return 0;
2833   }
2834
2835   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace);
2836   if(aSurf.IsNull()) return aRes;
2837   gp_Pnt aPoint = BRep_Tool::Pnt(aPnt);
2838
2839   //Compute the parameters
2840   ShapeAnalysis_Surface sas(aSurf);
2841   gp_Pnt2d UV = sas.ValueOfUV(aPoint,Precision::Confusion());
2842
2843   return getSurfaceCurvatures(aSurf, UV.X(), UV.Y(), true);
2844 }
2845
2846
2847 //=============================================================================
2848 /*!
2849  *  MinSurfaceCurvatureByParam
2850  */
2851 //=============================================================================
2852 Standard_Real GEOMImpl_IMeasureOperations::MinSurfaceCurvatureByParam
2853                                                   (Handle(GEOM_Object) theSurf,
2854                                                    Standard_Real& theUParam,
2855                                                    Standard_Real& theVParam)
2856 {
2857   SetErrorCode(KO);
2858   Standard_Real aRes = -1.0;
2859
2860   if (theSurf.IsNull()) return aRes;
2861
2862   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2863   if(aRefShape.IsNull()) return aRes;
2864
2865   TopoDS_Shape aShape = aRefShape->GetValue();
2866   if(aShape.IsNull()) {
2867     SetErrorCode("One of Objects has NULL Shape");
2868     return aRes;
2869   }
2870
2871   TopoDS_Face F = TopoDS::Face(aShape);
2872   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(F);
2873
2874   //Compute the parameters
2875   Standard_Real U1,U2,V1,V2;
2876   BRepTools::UVBounds(F,U1,U2,V1,V2);
2877   Standard_Real U = U1 + (U2-U1)*theUParam;
2878   Standard_Real V = V1 + (V2-V1)*theVParam;
2879
2880   return getSurfaceCurvatures(aSurf, U, V, false);
2881 }
2882
2883
2884 //=============================================================================
2885 /*!
2886  *  MinSurfaceCurvatureByPoint
2887  */
2888 //=============================================================================
2889 Standard_Real GEOMImpl_IMeasureOperations::MinSurfaceCurvatureByPoint
2890                     (Handle(GEOM_Object) theSurf, Handle(GEOM_Object) thePoint)
2891 {
2892   SetErrorCode(KO);
2893   Standard_Real aRes = -1.0;
2894
2895   if( theSurf.IsNull() || thePoint.IsNull() ) return aRes;
2896
2897   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2898   Handle(GEOM_Function) aRefPoint = thePoint->GetLastFunction();
2899   if( aRefShape.IsNull() || aRefPoint.IsNull() ) return aRes;
2900
2901   TopoDS_Face aFace = TopoDS::Face(aRefShape->GetValue());
2902   TopoDS_Vertex aPnt = TopoDS::Vertex(aRefPoint->GetValue());
2903   if( aFace.IsNull() || aPnt.IsNull() ) {
2904     SetErrorCode("One of Objects has NULL Shape");
2905     return 0;
2906   }
2907
2908   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace);
2909   if(aSurf.IsNull()) return aRes;
2910   gp_Pnt aPoint = BRep_Tool::Pnt(aPnt);
2911
2912   //Compute the parameters
2913   ShapeAnalysis_Surface sas(aSurf);
2914   gp_Pnt2d UV = sas.ValueOfUV(aPoint,Precision::Confusion());
2915
2916   return getSurfaceCurvatures(aSurf, UV.X(), UV.Y(), false);
2917 }
2918
2919 //=============================================================================
2920 /*!
2921  *  SurfaceCurvatureByPointAndDirection
2922  */
2923 //=============================================================================
2924 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::SurfaceCurvatureByPointAndDirection
2925                                                  (Handle(GEOM_Object) theSurf,
2926                                                   Handle(GEOM_Object) thePoint,
2927                                                   Handle(GEOM_Object) theDirection)
2928 {
2929   SetErrorCode(KO);
2930
2931   if (theSurf.IsNull() || thePoint.IsNull() || theDirection.IsNull()) return NULL;
2932
2933   Handle(GEOM_Function) aSurf = theSurf->GetLastFunction();
2934   Handle(GEOM_Function) aPoint = thePoint->GetLastFunction();
2935   Handle(GEOM_Function) aDirection = theDirection->GetLastFunction();
2936   if (aSurf.IsNull() || aPoint.IsNull() || aDirection.IsNull()) return NULL;
2937
2938   //Add a new CurvatureVector object
2939   //Handle(GEOM_Object) aCV = GetEngine()->AddObject(GEOM_CURVATURE_VEC);
2940   Handle(GEOM_Object) aCV = GetEngine()->AddObject(GEOM_VECTOR);
2941
2942   //Add a new CurvatureVector function
2943   Handle(GEOM_Function) aFunction =
2944     aCV->AddFunction(GEOMImpl_MeasureDriver::GetID(), CURVATURE_VEC_MEASURE);
2945   if (aFunction.IsNull()) return NULL;
2946
2947   //Check if the function is set correctly
2948   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
2949
2950   GEOMImpl_IMeasure aCI (aFunction);
2951   aCI.SetBase(aSurf);
2952   aCI.SetPoint(aPoint);
2953   aCI.SetDirection(aDirection);
2954
2955   //Compute the CurvatureVector
2956   try {
2957     OCC_CATCH_SIGNALS;
2958     if (!GetSolver()->ComputeFunction(aFunction)) {
2959       SetErrorCode("Measure driver failed to compute a surface curvature");
2960       return NULL;
2961     }
2962   }
2963   catch (Standard_Failure& aFail) {
2964     SetErrorCode(aFail.GetMessageString());
2965     return NULL;
2966   }
2967
2968   //Make a Python command
2969   GEOM::TPythonDump(aFunction) << aCV << " = geompy.CurvatureOnFace(" << theSurf
2970                                << ", " << thePoint << ", " << theDirection << ")";
2971
2972   SetErrorCode(OK);
2973   return aCV;
2974 }
2975
2976 //=============================================================================
2977 /*!
2978  *  XYZtoUV
2979  */
2980  //=============================================================================
2981 Handle(TColStd_HArray1OfReal) GEOMImpl_IMeasureOperations::XYZtoUV
2982                                                   (Handle(GEOM_Object) theSurf,
2983                                                    const Handle(TColStd_HArray1OfReal)& theXYZlist,
2984                                                    bool isNormalized)
2985 {
2986   SetErrorCode(KO);
2987
2988   Handle(TColStd_HArray1OfReal) aRet;
2989
2990   // Check list of coordinates
2991   int nbC = theXYZlist->Length();
2992   int nbP = nbC / 3;
2993   if (nbP * 3 != nbC) {
2994     SetErrorCode("Coordinates list length is not divisible by 3");
2995     return aRet;
2996   }
2997
2998   // Check face
2999   if (theSurf.IsNull()) {
3000     SetErrorCode("The shape is NULL");
3001     return aRet;
3002   }
3003
3004   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
3005   if (aRefShape.IsNull()) {
3006     SetErrorCode("The shape is NULL");
3007     return aRet;
3008   }
3009
3010   TopoDS_Shape aShape = aRefShape->GetValue();
3011   if (aShape.IsNull()) {
3012     SetErrorCode("The shape is NULL");
3013     return aRet;
3014   }
3015
3016   // The shape can be a face, a shell of one face or a compound with one face
3017   TopoDS_Face F;
3018   if (aShape.ShapeType() == TopAbs_FACE) {
3019     F = TopoDS::Face(aShape);
3020   }
3021   else if (aShape.ShapeType() < TopAbs_FACE) {
3022     TopExp_Explorer Exp (aShape, TopAbs_FACE);
3023     if (Exp.More()) {
3024       F = TopoDS::Face(Exp.Current());
3025       Exp.Next();
3026       if (Exp.More()) {
3027         SetErrorCode("There should be only one face");
3028         return aRet;
3029       }
3030     }
3031   }
3032   if (F.IsNull()) {
3033     SetErrorCode("There are no faces");
3034     return aRet;
3035   }
3036
3037   // Face tolerance
3038   ShapeAnalysis_ShapeTolerance sat;
3039   sat.InitTolerance();
3040   sat.AddTolerance(F);
3041   Standard_Real aTol = sat.GlobalTolerance(1);
3042   Standard_Real squareTolerance = aTol * aTol;
3043
3044   // Compute parameters
3045   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(F);
3046   aRet = new TColStd_HArray1OfReal (0, nbP * 2 - 1);
3047
3048   Standard_Real U1,U2, V1,V2;
3049   BRepTools::UVBounds(F, U1, U2, V1, V2);
3050   Standard_Real dU = U2 - U1;
3051   Standard_Real dV = V2 - V1;
3052
3053   int iCLower = theXYZlist->Lower();
3054   for (int iP = 0; iP < nbP; iP++) {
3055     gp_Pnt aP (theXYZlist->Value(iCLower + iP * 3),
3056                theXYZlist->Value(iCLower + iP * 3 + 1),
3057                theXYZlist->Value(iCLower + iP * 3 + 2));
3058     try { // as GEOMUtils::ProjectPointOnFace can throw exceptions
3059       Standard_Real U, V;
3060       gp_Pnt aPonF = GEOMUtils::ProjectPointOnFace(aP, F, U, V, aTol);
3061       if (aP.SquareDistance(aPonF) < squareTolerance) {
3062         if (isNormalized) {
3063           // Normalize parameters to be in [0, 1]
3064           U = (U - U1) / dU;
3065           V = (V - V1) / dV;
3066         }
3067         aRet->SetValue(iP * 2    , U);
3068         aRet->SetValue(iP * 2 + 1, V);
3069       }
3070       else {
3071         SetErrorCode("Point too far from face");
3072         return aRet;
3073       }
3074     }
3075     catch (Standard_Failure& aFail) {
3076       SetErrorCode(aFail.GetMessageString());
3077       return aRet;
3078     }
3079   }
3080
3081   SetErrorCode(OK);
3082   return aRet;
3083 }
3084
3085 //=============================================================================
3086 /*!
3087  *  UVtoXYZ
3088  */
3089  //=============================================================================
3090 Handle(TColStd_HArray1OfReal) GEOMImpl_IMeasureOperations::UVtoXYZ
3091                                                   (Handle(GEOM_Object) theSurf,
3092                                                    const Handle(TColStd_HArray1OfReal)& theUVlist,
3093                                                    bool isNormalized)
3094 {
3095   SetErrorCode(KO);
3096
3097   Handle(TColStd_HArray1OfReal) aRet;
3098
3099   // Check list of parameters
3100   int nbC = theUVlist->Length();
3101   int nbP = nbC / 2;
3102   if (nbP * 2 != nbC) {
3103     SetErrorCode("Parameters list length is not divisible by 2");
3104     return aRet;
3105   }
3106
3107   // Check face
3108   if (theSurf.IsNull()) {
3109     SetErrorCode("The shape is NULL");
3110     return aRet;
3111   }
3112
3113   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
3114   if (aRefShape.IsNull()) {
3115     SetErrorCode("The shape is NULL");
3116     return aRet;
3117   }
3118
3119   TopoDS_Shape aShape = aRefShape->GetValue();
3120   if (aShape.IsNull()) {
3121     SetErrorCode("The shape is NULL");
3122     return aRet;
3123   }
3124
3125   // The shape can be a face, a shell of one face or a compound with one face
3126   TopoDS_Face F;
3127   if (aShape.ShapeType() == TopAbs_FACE) {
3128     F = TopoDS::Face(aShape);
3129   }
3130   else if (aShape.ShapeType() < TopAbs_FACE) {
3131     TopExp_Explorer Exp (aShape, TopAbs_FACE);
3132     if (Exp.More()) {
3133       F = TopoDS::Face(Exp.Current());
3134       Exp.Next();
3135       if (Exp.More()) {
3136         SetErrorCode("There should be only one face");
3137         return aRet;
3138       }
3139     }
3140   }
3141   if (F.IsNull()) {
3142     SetErrorCode("There are no faces");
3143     return aRet;
3144   }
3145
3146   // Face tolerance
3147   Standard_Real squareTolerance = BRep_Tool::Tolerance(F);
3148   squareTolerance = squareTolerance * squareTolerance;
3149
3150   // Compute coordinates
3151   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(F);
3152   aRet = new TColStd_HArray1OfReal (0, nbP * 3 - 1);
3153
3154   Standard_Real U1,U2, V1,V2;
3155   BRepTools::UVBounds(F, U1, U2, V1, V2);
3156   Standard_Real dU = U2 - U1;
3157   Standard_Real dV = V2 - V1;
3158
3159   Standard_Real tol = 1.e-4;
3160   Standard_Real pc = Precision::Confusion();
3161
3162   int iCLower = theUVlist->Lower();
3163   for (int iP = 0; iP < nbP; iP++) {
3164     Standard_Real U = theUVlist->Value(iCLower + iP * 2);
3165     Standard_Real V = theUVlist->Value(iCLower + iP * 2 + 1);
3166
3167     if (isNormalized) {
3168       // Get real parameters from given normalized ones in [0, 1]
3169       if (!(-pc < U && U < 1 + pc) || !(-pc < V && V < 1 + pc)) {
3170         SetErrorCode("Normalized parameter is out of range [0,1]");
3171         return aRet;
3172       }
3173       U = U1 + dU * U;
3174       V = V1 + dV * V;
3175     }
3176
3177     gp_Pnt2d aP2d (U, V);
3178
3179     BRepClass_FaceClassifier aClsf (F, aP2d, tol);
3180     if (aClsf.State() != TopAbs_IN && aClsf.State() != TopAbs_ON) {
3181       SetErrorCode("Given parameters are out of face");
3182       return aRet;
3183     }
3184     gp_Pnt surfPnt = aSurf->Value(U, V);
3185
3186     aRet->SetValue(iP * 3    , surfPnt.X());
3187     aRet->SetValue(iP * 3 + 1, surfPnt.Y());
3188     aRet->SetValue(iP * 3 + 2, surfPnt.Z());
3189   }
3190
3191   SetErrorCode(OK);
3192   return aRet;
3193 }
3194
3195 //=============================================================================
3196 /*!
3197  *  SelfIntersected2D
3198  *  Find all self-intersected 2D curves.
3199  *  \param theChecks list of failed checks, contains type of check and failed shapes
3200  */
3201  //=============================================================================
3202 std::list<GEOMImpl_IMeasureOperations::CoupleOfObjects>
3203   GEOMImpl_IMeasureOperations::SelfIntersected2D(const std::list<FailedChecks>& theChecks)
3204 {
3205   SetErrorCode(KO);
3206   MESSAGE("GEOMImpl_IMeasureOperations::selfIntersected2D");
3207
3208   std::list<GEOMImpl_IMeasureOperations::CoupleOfObjects> aSelfInters2D;
3209   try
3210   {
3211     OCC_CATCH_SIGNALS;
3212     for (std::list<FailedChecks>::const_iterator anIter(theChecks.begin());
3213          anIter != theChecks.end(); ++anIter)
3214     {
3215       if (anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_InvalidCurveOnSurface)
3216         aSelfInters2D.push_back(anIter->FailedShapes);
3217     }
3218   }
3219   catch (Standard_Failure& aFail)
3220   {
3221     SetErrorCode(aFail.GetMessageString());
3222     return aSelfInters2D;
3223   }
3224
3225   SetErrorCode(OK);
3226   return aSelfInters2D;
3227 }
3228
3229 namespace
3230 {
3231   static bool checkTypes(const GEOMImpl_IMeasureOperations::CoupleOfObjects& theShapes,
3232                          const int theShapeType1,
3233                          const int theShapeType2)
3234   {
3235     if (theShapeType1 == -1 && theShapeType2 == -1)
3236       return true;
3237
3238     TopAbs_ShapeEnum aShapeType1 = theShapes.first.IsNull()
3239       ? TopAbs_SHAPE
3240       : theShapes.first->GetValue().ShapeType();
3241     TopAbs_ShapeEnum aShapeType2 = theShapes.second.IsNull()
3242       ? TopAbs_SHAPE
3243       : theShapes.second->GetValue().ShapeType();
3244
3245     if (theShapeType1 == -1)
3246       return aShapeType1 == theShapeType2 || aShapeType2 == theShapeType2;
3247     else if (theShapeType2 == -1)
3248       return aShapeType1 == theShapeType1 || aShapeType2 == theShapeType1;
3249     return (aShapeType1 == theShapeType1 && aShapeType2 == theShapeType2) ||
3250       (aShapeType1 == theShapeType2 && aShapeType2 == theShapeType1);
3251   }
3252 } // namespace
3253
3254 //=============================================================================
3255 /*!
3256  *  InterferingSubshapes
3257  *  Find pairs of interfering sub-shapes, by default all pairs of interfering shapes are returned.
3258  *  \param theChecks list of failed checks, contains type of check and failed shapes
3259  *  \param theShapeType1 Type of shape.
3260  *  \param theShapeType2 Type of shape.
3261  */
3262  //=============================================================================
3263 std::list<GEOMImpl_IMeasureOperations::CoupleOfObjects>
3264   GEOMImpl_IMeasureOperations::InterferingSubshapes
3265     (const std::list<FailedChecks>& theChecks,
3266      const int                      theShapeType1,
3267      const int                      theShapeType2)
3268 {
3269   SetErrorCode(KO);
3270   MESSAGE("GEOMImpl_IMeasureOperations::interferingSubshapes");
3271
3272   std::list<GEOMImpl_IMeasureOperations::CoupleOfObjects> anInterfer;
3273   try
3274   {
3275     OCC_CATCH_SIGNALS;
3276     for (std::list<FailedChecks>::const_iterator anIter(theChecks.begin());
3277          anIter != theChecks.end(); ++anIter)
3278     {
3279       if (anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_SelfIntersect &&
3280           checkTypes(anIter->FailedShapes, theShapeType1, theShapeType2))
3281         anInterfer.push_back(anIter->FailedShapes);
3282     }
3283   }
3284   catch (Standard_Failure& aFail)
3285   {
3286     SetErrorCode(aFail.GetMessageString());
3287     return anInterfer;
3288   }
3289
3290   SetErrorCode(OK);
3291   return anInterfer;
3292 }
3293
3294 //=============================================================================
3295 /*!
3296  *  SmallEdges
3297  *  Find edges, which are fully covered by tolerances of vertices.
3298  *  \param theChecks list of failed checks, contains type of check and failed shapes
3299  */
3300  //=============================================================================
3301 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IMeasureOperations::SmallEdges(
3302     const std::list<FailedChecks>& theChecks)
3303 {
3304   SetErrorCode(KO);
3305   MESSAGE("GEOMImpl_IMeasureOperations::smallEdges");
3306
3307   Handle(TColStd_HSequenceOfTransient) aSmallEdges = new TColStd_HSequenceOfTransient;
3308   try
3309   {
3310     OCC_CATCH_SIGNALS;
3311     for (std::list<FailedChecks>::const_iterator anIter(theChecks.begin());
3312          anIter != theChecks.end(); ++anIter)
3313     {
3314       if (anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_TooSmallEdge)
3315         aSmallEdges->Append(anIter->FailedShapes.first);
3316     }
3317   }
3318   catch (Standard_Failure& aFail)
3319   {
3320     SetErrorCode(aFail.GetMessageString());
3321     return NULL;
3322   }
3323
3324   SetErrorCode(OK);
3325   return aSmallEdges;
3326 }
3327
3328 //=============================================================================
3329 /*!
3330  *  DistantShapes
3331  *  find remote objects (sub-shape on a shape).
3332  *  \param theShape Shape for check.
3333  *  \param theShapeType Type of shape.
3334  *  \param theSubShapeType Type of sub-shape.
3335  *  \param theTolerance tolerance.
3336  */
3337  //=============================================================================
3338 std::list<GEOMImpl_IMeasureOperations::CoupleOfObjects>
3339   GEOMImpl_IMeasureOperations::DistantShapes
3340     (const std::list<FailedChecks>& theChecks,
3341      const int                      theShapeType,
3342      const int                      theSubShapeType,
3343      double                         theTolerance)
3344 {
3345   SetErrorCode(KO);
3346   MESSAGE("GEOMImpl_IMeasureOperations::distantShapes");
3347
3348   std::list<GEOMImpl_IMeasureOperations::CoupleOfObjects> aDistShapes;
3349   try
3350   {
3351     OCC_CATCH_SIGNALS;
3352     for (std::list<FailedChecks>::const_iterator anIter(theChecks.begin());
3353          anIter != theChecks.end(); ++anIter)
3354     {
3355       Handle(GEOM_Object) aSubShape = anIter->FailedShapes.first;
3356       Handle(GEOM_Object) aShape = anIter->FailedShapes.second;
3357       if ((anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_InvalidCurveOnSurface ||
3358            anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_IncompatibilityOfVertex ||
3359            anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_IncompatibilityOfEdge ||
3360            anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_IncompatibilityOfFace) &&
3361           aShape && (theShapeType == -1 || aShape->GetValue().ShapeType() == theShapeType) &&
3362           aSubShape && (theSubShapeType == -1 || aSubShape->GetValue().ShapeType() == theSubShapeType))
3363       {
3364         gp_XYZ aP1, aP2;
3365         Standard_Real aDist = Precision::Infinite();
3366         if (anIter->TypeOfCheck == BOPAlgo_CheckStatus::BOPAlgo_InvalidCurveOnSurface)
3367           aDist = ComputeTolerance(aSubShape, aShape);
3368         if (aDist > theTolerance)
3369           aDistShapes.push_back(anIter->FailedShapes);
3370       }
3371     }
3372   }
3373   catch (Standard_Failure& aFail)
3374   {
3375     SetErrorCode(aFail.GetMessageString());
3376     return aDistShapes;
3377   }
3378
3379   SetErrorCode(OK);
3380   return aDistShapes;
3381 }
3382
3383 //=============================================================================
3384 /*!
3385  *  CheckConformityShape
3386  *  Perform analyse of shape and find imperfections in the shape.
3387  *  \param theShape Shape for analyse.
3388  */
3389  //=============================================================================
3390 void GEOMImpl_IMeasureOperations::CheckConformityShape(Handle(GEOM_Object) theShape, std::list<FailedChecks>& theChecks)
3391 {
3392   SetErrorCode(KO);
3393   MESSAGE("GEOMImpl_IMeasureOperations::checkShape");
3394
3395   Handle(GEOM_Object) aConformity = GetEngine()->AddObject(GEOM_CHECKCONFORMITY);
3396   Handle(GEOM_Function) aFunction = aConformity->AddFunction(GEOMImpl_ConformityDriver::GetID(), CONFORMITY_CHECK_SHAPE);
3397   if (aFunction.IsNull()) return;
3398
3399   //Check if the function is set correctly
3400   if (aFunction->GetDriverGUID() != GEOMImpl_ConformityDriver::GetID()) return;
3401
3402   GEOMImpl_IConformity aCI(aFunction);
3403
3404   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
3405   if (aRefShape.IsNull()) return;
3406
3407   aCI.SetShape(aRefShape);
3408
3409   try
3410   {
3411     OCC_CATCH_SIGNALS;
3412     if (!GetSolver()->ComputeFunction(aFunction))
3413     {
3414       SetErrorCode("Failed: checkShape");
3415       return;
3416     }
3417     Handle(TColStd_HArray1OfInteger) aTypesChecks = aFunction->GetIntegerArray(CHECKCONFORMITY_RET_TYPES_CHECKS);
3418     Handle(TColStd_HArray2OfInteger) aRes = aCI.GetListOfShapesIndices();
3419     if (aRes.IsNull())
3420       return;
3421
3422     for (Standard_Integer anIndex = 1; anIndex <= aRes->NbRows(); ++anIndex)
3423     {
3424       std::pair<Handle(GEOM_Object), Handle(GEOM_Object)> aPair;
3425       Handle(TColStd_HArray1OfInteger) anArray;
3426       anArray = new TColStd_HArray1OfInteger(1, 1);
3427       anArray->SetValue(1, aRes->Value(anIndex, 1));
3428
3429       Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theShape, anArray);
3430       if (!anObj.IsNull())
3431         aPair.first = anObj;
3432
3433       anArray = new TColStd_HArray1OfInteger(1, 1);
3434       anArray->SetValue(1, aRes->Value(anIndex, 2));
3435
3436       anObj = GetEngine()->AddSubShape(theShape, anArray);
3437       if (!anObj.IsNull())
3438         aPair.second = anObj;
3439       theChecks.push_back({ aTypesChecks->Value(anIndex), aPair });
3440     }
3441   }
3442   catch (Standard_Failure& aFail)
3443   {
3444     SetErrorCode(aFail.GetMessageString());
3445     return;
3446   }
3447
3448   SetErrorCode(OK);
3449   return;
3450 }
3451
3452 //=============================================================================
3453 /*!
3454  *  UpdateTolerance
3455  *  Compute possible tolerance for the shape, minimize tolerance of shape as well
3456  *  as tolerance of sub-shapes as much as possible
3457  *  \param theShape Shape for compute tolerance.
3458  */
3459  //=============================================================================
3460 double GEOMImpl_IMeasureOperations::UpdateTolerance(Handle(GEOM_Object) theShape)
3461 {
3462   SetErrorCode(KO);
3463   MESSAGE("GEOMImpl_IMeasureOperations::updateTolerance");
3464
3465   double aResTol = -1;
3466   Handle(GEOM_Object) aConformity = GetEngine()->AddObject(GEOM_CHECKCONFORMITY);
3467   Handle(GEOM_Function) aFunction = aConformity->AddFunction(GEOMImpl_ConformityDriver::GetID(), CONFORMITY_UPDATE_TOL);
3468   if (aFunction.IsNull()) return aResTol;
3469
3470   //Check if the function is set correctly
3471   if (aFunction->GetDriverGUID() != GEOMImpl_ConformityDriver::GetID()) return aResTol;
3472
3473   GEOMImpl_IConformity aCI(aFunction);
3474
3475   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
3476   if (aRefShape.IsNull()) return aResTol;
3477
3478   aCI.SetShape(aRefShape);
3479
3480   try
3481   {
3482     OCC_CATCH_SIGNALS;
3483     if (!GetSolver()->ComputeFunction(aFunction))
3484     {
3485       SetErrorCode("Failed: updateTolerance");
3486       return aResTol;
3487     }
3488     aResTol = aFunction->GetReal(CHECKCONFORMITY_RET_TOLERANCE);
3489   }
3490   catch (Standard_Failure& aFail)
3491   {
3492     SetErrorCode(aFail.GetMessageString());
3493     return aResTol;
3494   }
3495
3496   SetErrorCode(OK);
3497   return aResTol;
3498 }
3499
3500 //=============================================================================
3501 /*!
3502  *  ComputeTolerance
3503  *  Compute distance from the edge to the face.
3504  */
3505  //=============================================================================
3506 double GEOMImpl_IMeasureOperations::ComputeTolerance(Handle(GEOM_Object) theEdge,
3507                                                      Handle(GEOM_Object) theFace)
3508 {
3509   double aMaxDist = Precision::Infinite();
3510   if (theEdge.IsNull() || theFace.IsNull())
3511     return aMaxDist;
3512
3513   Handle(GEOM_Function) aRefEdge = theEdge->GetLastFunction();
3514   Handle(GEOM_Function) aRefFace = theFace->GetLastFunction();
3515   if (aRefEdge.IsNull() || aRefFace.IsNull())
3516     return aMaxDist;
3517
3518   TopoDS_Edge aEdge = TopoDS::Edge(aRefEdge->GetValue());
3519   TopoDS_Face aFace = TopoDS::Face(aRefFace->GetValue());
3520   if (aEdge.IsNull() || aFace.IsNull())
3521     return aMaxDist;
3522
3523   double aParam = 0.0;
3524   BOPTools_AlgoTools::ComputeTolerance(aFace, aEdge, aMaxDist, aParam);
3525   return aMaxDist;
3526 }
3527
3528 //=======================================================================
3529 //function : FillErrorsSub
3530 //purpose  : Fill the errors list of subshapes on shape.
3531 //=======================================================================
3532 void GEOMImpl_IMeasureOperations::FillErrorsSub
3533            (const BRepCheck_Analyzer                   &theAna,
3534             const TopoDS_Shape                         &theShape,
3535             const TopAbs_ShapeEnum                     theSubType,
3536                   TopTools_DataMapOfIntegerListOfShape &theMapErrors) const
3537 {
3538   TopExp_Explorer anExp(theShape, theSubType);
3539   TopTools_MapOfShape aMapSubShapes;
3540
3541   for (; anExp.More(); anExp.Next()) {
3542     const TopoDS_Shape &aSubShape = anExp.Current();
3543
3544     if (aMapSubShapes.Add(aSubShape)) {
3545       const Handle(BRepCheck_Result) &aRes = theAna.Result(aSubShape);
3546
3547       for (aRes->InitContextIterator();
3548            aRes->MoreShapeInContext();
3549            aRes->NextShapeInContext()) {
3550         if (aRes->ContextualShape().IsSame(theShape)) {
3551           BRepCheck_ListIteratorOfListOfStatus itl(aRes->StatusOnShape());
3552
3553           if (itl.Value() != BRepCheck_NoError) {
3554             // Add all errors for theShape and its sub-shape.
3555             for (;itl.More(); itl.Next()) {
3556               const Standard_Integer aStat = (Standard_Integer)itl.Value();
3557
3558               if (!theMapErrors.IsBound(aStat)) {
3559                 TopTools_ListOfShape anEmpty;
3560
3561                 theMapErrors.Bind(aStat, anEmpty);
3562               }
3563
3564               TopTools_ListOfShape &theShapes = theMapErrors.ChangeFind(aStat);
3565
3566               theShapes.Append(aSubShape);
3567               theShapes.Append(theShape);
3568             }
3569           }
3570         }
3571
3572         break;
3573       }
3574     }
3575   }
3576 }
3577
3578 //=======================================================================
3579 //function : FillErrors
3580 //purpose  : Fill the errors list.
3581 //=======================================================================
3582 void GEOMImpl_IMeasureOperations::FillErrors
3583              (const BRepCheck_Analyzer                   &theAna,
3584               const TopoDS_Shape                         &theShape,
3585                     TopTools_DataMapOfIntegerListOfShape &theMapErrors,
3586                     TopTools_MapOfShape                  &theMapShapes) const
3587 {
3588   if (theMapShapes.Add(theShape)) {
3589     // Fill errors of child shapes.
3590     for (TopoDS_Iterator iter(theShape); iter.More(); iter.Next()) {
3591       FillErrors(theAna, iter.Value(), theMapErrors, theMapShapes);
3592     }
3593
3594     // Fill errors of theShape.
3595     const Handle(BRepCheck_Result) &aRes = theAna.Result(theShape);
3596
3597     if (!aRes.IsNull()) {
3598       BRepCheck_ListIteratorOfListOfStatus itl(aRes->Status());
3599
3600       if (itl.Value() != BRepCheck_NoError) {
3601         // Add all errors for theShape.
3602         for (;itl.More(); itl.Next()) {
3603           const Standard_Integer aStat = (Standard_Integer)itl.Value();
3604
3605           if (!theMapErrors.IsBound(aStat)) {
3606             TopTools_ListOfShape anEmpty;
3607
3608             theMapErrors.Bind(aStat, anEmpty);
3609           }
3610
3611           theMapErrors.ChangeFind(aStat).Append(theShape);
3612         }
3613       }
3614     }
3615
3616     // Add errors of subshapes on theShape.
3617     const TopAbs_ShapeEnum aType = theShape.ShapeType();
3618
3619     switch (aType) {
3620     case TopAbs_EDGE:
3621       FillErrorsSub(theAna, theShape, TopAbs_VERTEX, theMapErrors);
3622       break;
3623     case TopAbs_FACE:
3624       FillErrorsSub(theAna, theShape, TopAbs_WIRE, theMapErrors);
3625       FillErrorsSub(theAna, theShape, TopAbs_EDGE, theMapErrors);
3626       FillErrorsSub(theAna, theShape, TopAbs_VERTEX, theMapErrors);
3627       break;
3628     case TopAbs_SOLID:
3629       FillErrorsSub(theAna, theShape, TopAbs_SHELL, theMapErrors);
3630       break;
3631     default:
3632       break;
3633     }
3634   }
3635 }
3636
3637 //=======================================================================
3638 //function : FillErrors
3639 //purpose  : Fill the errors list.
3640 //=======================================================================
3641 void GEOMImpl_IMeasureOperations::FillErrors
3642                   (const BRepCheck_Analyzer    &theAna,
3643                    const TopoDS_Shape          &theShape,
3644                          std::list<ShapeError> &theErrors) const
3645 {
3646   // Fill the errors map.
3647   TopTools_DataMapOfIntegerListOfShape aMapErrors;
3648   TopTools_MapOfShape                  aMapShapes;
3649
3650   FillErrors(theAna, theShape, aMapErrors, aMapShapes);
3651
3652   // Map sub-shapes and their indices
3653   TopTools_IndexedMapOfShape anIndices;
3654
3655   TopExp::MapShapes(theShape, anIndices);
3656
3657   TopTools_DataMapIteratorOfDataMapOfIntegerListOfShape aMapIter(aMapErrors);
3658
3659   for (; aMapIter.More(); aMapIter.Next()) {
3660     ShapeError anError;
3661
3662     anError.error = (BRepCheck_Status)aMapIter.Key();
3663
3664     TopTools_ListIteratorOfListOfShape aListIter(aMapIter.Value());
3665     TopTools_MapOfShape                aMapUnique;
3666
3667     for (; aListIter.More(); aListIter.Next()) {
3668       const TopoDS_Shape &aShape = aListIter.Value();
3669
3670       if (aMapUnique.Add(aShape)) {
3671         const Standard_Integer anIndex = anIndices.FindIndex(aShape);
3672
3673         anError.incriminated.push_back(anIndex);
3674       }
3675     }
3676
3677     if (!anError.incriminated.empty()) {
3678       theErrors.push_back(anError);
3679     }
3680   }
3681 }
3682
3683 //=======================================================================
3684 //function : ShapeProximityCalculator
3685 //purpose  : returns an object to compute the proximity value
3686 //=======================================================================
3687 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::ShapeProximityCalculator
3688                                           (Handle(GEOM_Object) theShape1,
3689                                            Handle(GEOM_Object) theShape2)
3690 {
3691   SetErrorCode(KO);
3692
3693   if (theShape1.IsNull() || theShape2.IsNull())
3694     return NULL;
3695
3696   Handle(GEOM_Function) aShapeFunc1 = theShape1->GetLastFunction();
3697   Handle(GEOM_Function) aShapeFunc2 = theShape2->GetLastFunction();
3698   if (aShapeFunc1.IsNull() || aShapeFunc2.IsNull())
3699     return NULL;
3700
3701   Handle(GEOM_Object) aProximityCalc = GetEngine()->AddObject(GEOM_SHAPE_PROXIMITY);
3702   if (aProximityCalc.IsNull())
3703     return NULL;
3704
3705   Handle(GEOM_Function) aProximityFuncCoarse =
3706       aProximityCalc->AddFunction(GEOMImpl_ShapeProximityDriver::GetID(), PROXIMITY_COARSE);
3707   //Check if the function is set correctly
3708   if (aProximityFuncCoarse.IsNull() ||
3709       aProximityFuncCoarse->GetDriverGUID() != GEOMImpl_ShapeProximityDriver::GetID())
3710     return NULL;
3711
3712   GEOMImpl_IProximity aProximity (aProximityFuncCoarse);
3713   aProximity.SetShapes(aShapeFunc1, aShapeFunc2);
3714
3715   //Make a Python command
3716   GEOM::TPythonDump pd (aProximityFuncCoarse);
3717   pd << "p = geompy.ShapeProximity()\n";
3718   pd << "p.setShapes(" << theShape1 << ", " << theShape2 << ")";
3719
3720   SetErrorCode(OK);
3721   return aProximityCalc;
3722 }
3723
3724 //=======================================================================
3725 //function : SetShapeSampling
3726 //purpose  : set number sample points to compute the coarse proximity
3727 //=======================================================================
3728 void GEOMImpl_IMeasureOperations::SetShapeSampling(Handle(GEOM_Object) theCalculator,
3729                                                    Handle(GEOM_Object) theShape,
3730                                                    const Standard_Integer theNbSamples)
3731 {
3732   SetErrorCode(KO);
3733   if (theShape.IsNull() ||
3734       theCalculator.IsNull() ||
3735       theCalculator->GetNbFunctions() <= 0 ||
3736       theNbSamples <= 0)
3737     return ;
3738
3739   Handle(GEOM_Function) aProximityFuncCoarse = theCalculator->GetFunction(1);
3740   if (aProximityFuncCoarse.IsNull() ||
3741       aProximityFuncCoarse->GetDriverGUID() != GEOMImpl_ShapeProximityDriver::GetID())
3742     return ;
3743
3744   Handle(GEOM_Function) aShapeFunc = theShape->GetLastFunction();
3745   if (aShapeFunc.IsNull())
3746     return ;
3747
3748   GEOMImpl_IProximity aProximity(aProximityFuncCoarse);
3749   Handle(GEOM_Function) aShape1, aShape2;
3750   aProximity.GetShapes(aShape1, aShape2);
3751   if (aShape1->GetValue() == aShapeFunc->GetValue())
3752     aProximity.SetNbSamples(PROXIMITY_ARG_SAMPLES1, theNbSamples);
3753   else if (aShape2->GetValue() == aShapeFunc->GetValue())
3754     aProximity.SetNbSamples(PROXIMITY_ARG_SAMPLES2, theNbSamples);
3755
3756   //Make a Python command
3757   GEOM::TPythonDump(aProximityFuncCoarse, /*append=*/true) <<
3758     "p.setSampling(" << theShape << ", " << theNbSamples << ")";
3759
3760   SetErrorCode(OK);
3761 }
3762
3763 //=======================================================================
3764 //function : GetCoarseProximity
3765 //purpose  : compute coarse proximity
3766 //=======================================================================
3767 Standard_Real GEOMImpl_IMeasureOperations::GetCoarseProximity(Handle(GEOM_Object) theCalculator,
3768                                                               bool doPythonDump)
3769 {
3770   SetErrorCode(KO);
3771   if (theCalculator.IsNull())
3772     return -1;
3773
3774   Handle(GEOM_Function) aProximityFuncCoarse = theCalculator->GetFunction(1);
3775   if (aProximityFuncCoarse.IsNull() ||
3776       aProximityFuncCoarse->GetDriverGUID() != GEOMImpl_ShapeProximityDriver::GetID() ||
3777       aProximityFuncCoarse->GetType() != PROXIMITY_COARSE)
3778     return -1;
3779
3780   // Perform
3781   // We have to recompute the function each time,
3782   // because the number of samples can be changed
3783   try {
3784     OCC_CATCH_SIGNALS;
3785     if (!GetSolver()->ComputeFunction(aProximityFuncCoarse)) {
3786       SetErrorCode("shape proximity driver failed");
3787       return -1;
3788     }
3789   }
3790   catch (Standard_Failure& aFail) {
3791     SetErrorCode(aFail.GetMessageString());
3792     return -1;
3793   }
3794
3795   //Make a Python command
3796   if (doPythonDump)
3797     GEOM::TPythonDump(aProximityFuncCoarse, /*append=*/true) << "value = p.coarseProximity()";
3798
3799   SetErrorCode(OK);
3800   GEOMImpl_IProximity aProximity (aProximityFuncCoarse);
3801   return aProximity.GetValue();
3802 }
3803
3804 //=======================================================================
3805 //function : GetPreciseProximity
3806 //purpose  : compute precise proximity
3807 //=======================================================================
3808 Standard_Real GEOMImpl_IMeasureOperations::GetPreciseProximity(Handle(GEOM_Object) theCalculator)
3809 {
3810   SetErrorCode(KO);
3811   if (theCalculator.IsNull())
3812     return -1;
3813
3814   Handle(GEOM_Function) aProximityFuncCoarse = theCalculator->GetFunction(1);
3815   Handle(GEOM_Function) aProximityFuncFine = theCalculator->GetFunction(2);
3816   if (aProximityFuncFine.IsNull())
3817     aProximityFuncFine = theCalculator->AddFunction
3818       (GEOMImpl_ShapeProximityDriver::GetID(), PROXIMITY_PRECISE);
3819
3820   //Check if the functions are set correctly
3821   if (aProximityFuncCoarse.IsNull() ||
3822       aProximityFuncCoarse->GetDriverGUID() != GEOMImpl_ShapeProximityDriver::GetID() ||
3823       aProximityFuncFine.IsNull() ||
3824       aProximityFuncFine->GetDriverGUID() != GEOMImpl_ShapeProximityDriver::GetID())
3825     return -1;
3826
3827   // perform coarse computation beforehand
3828   GetCoarseProximity(theCalculator, /*doPythonDump=*/false);
3829
3830   // transfer parameters from the coarse to precise calculator
3831   GEOMImpl_IProximity aCoarseProximity (aProximityFuncCoarse);
3832   Handle(GEOM_Function) aShape1, aShape2;
3833   aCoarseProximity.GetShapes(aShape1, aShape2);
3834   if (aShape1.IsNull() || aShape2.IsNull())
3835     return -1;
3836   gp_Pnt aProxPnt1, aProxPnt2;
3837   Standard_Integer intStatus1, intStatus2;
3838   aCoarseProximity.GetProximityPoints(aProxPnt1, aProxPnt2);
3839   aCoarseProximity.GetStatusOfPoints(intStatus1, intStatus2);
3840   Standard_Real aResultValue = aCoarseProximity.GetValue();
3841
3842   GEOMImpl_IProximity aFineProximity (aProximityFuncFine);
3843   aFineProximity.SetShapes(aShape1, aShape2);
3844   aFineProximity.SetProximityPoints(aProxPnt1, aProxPnt2);
3845   aFineProximity.SetStatusOfPoints(intStatus1, intStatus2);
3846   aFineProximity.SetValue(aResultValue); // in some cases this value cannot be precised
3847
3848   // Perform
3849   try {
3850     OCC_CATCH_SIGNALS;
3851     if (!GetSolver()->ComputeFunction(aProximityFuncFine)) {
3852       SetErrorCode("shape proximity driver failed");
3853       return -1;
3854     }
3855   }
3856   catch (Standard_Failure& aFail) {
3857     SetErrorCode(aFail.GetMessageString());
3858     return -1;
3859   }
3860
3861   aResultValue = aFineProximity.GetValue();
3862   aFineProximity.GetProximityPoints(aProxPnt1, aProxPnt2);
3863
3864   //Make a Python command
3865   GEOM::TPythonDump(aProximityFuncCoarse, /*append=*/true) << "value = p.preciseProximity()";
3866
3867   SetErrorCode(OK);
3868   return aResultValue;
3869 }