src/GEOMImpl/GEOMImpl_IMeasureOperations.cxx

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