src/GEOMImpl/GEOMImpl_IShapesOperations.cxx

   1 // Copyright (C) 2007-2011  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21
  22 //  File      : GEOMImpl_IShapesOperations.cxx
  23 //  Created   :
  24 //  Author    : modified by Lioka RAZAFINDRAZAKA (CEA) 22/06/2007
  25 //  Project   : SALOME
  26
  27 #include <Standard_Stream.hxx>
  28
  29 #include "GEOMImpl_IShapesOperations.hxx"
  30
  31 #include "GEOMImpl_Types.hxx"
  32
  33 #include "GEOMImpl_VectorDriver.hxx"
  34 #include "GEOMImpl_ShapeDriver.hxx"
  35 #include "GEOMImpl_CopyDriver.hxx"
  36 #include "GEOMImpl_GlueDriver.hxx"
  37
  38 #include "GEOMImpl_IVector.hxx"
  39 #include "GEOMImpl_IShapes.hxx"
  40 #include "GEOMImpl_IGlue.hxx"
  41
  42 #include "GEOMImpl_Block6Explorer.hxx"
  43
  44 #include "GEOM_Function.hxx"
  45 #include "GEOM_ISubShape.hxx"
  46 #include "GEOM_PythonDump.hxx"
  47
  48 #include "GEOMAlgo_ClsfBox.hxx"
  49 #include "GEOMAlgo_ClsfSolid.hxx"
  50 #include "GEOMAlgo_CoupleOfShapes.hxx"
  51 #include "GEOMAlgo_FinderShapeOn1.hxx"
  52 #include "GEOMAlgo_FinderShapeOnQuad.hxx"
  53 #include "GEOMAlgo_FinderShapeOn2.hxx"
  54 #include "GEOMAlgo_GetInPlace.hxx"
  55 #include "GEOMAlgo_GlueDetector.hxx"
  56 #include "GEOMAlgo_ListIteratorOfListOfCoupleOfShapes.hxx"
  57 #include "GEOMAlgo_ListOfCoupleOfShapes.hxx"
  58
  59 #include <Basics_OCCTVersion.hxx>
  60
  61 #include "utilities.h"
  62 #include "OpUtil.hxx"
  63 #include "Utils_ExceptHandlers.hxx"
  64
  65 #include <TFunction_DriverTable.hxx>
  66 #include <TFunction_Driver.hxx>
  67 #include <TFunction_Logbook.hxx>
  68 #include <TDataStd_Integer.hxx>
  69 #include <TDataStd_IntegerArray.hxx>
  70 #include <TDataStd_ListIteratorOfListOfExtendedString.hxx>
  71 #include <TDF_Tool.hxx>
  72
  73 #include <BRepExtrema_ExtCF.hxx>
  74 #include <BRepExtrema_DistShapeShape.hxx>
  75
  76 #include <BRep_Tool.hxx>
  77 #include <BRep_Builder.hxx>
  78 #include <BRepTools.hxx>
  79 #include <BRepGProp.hxx>
  80 #include <BRepAdaptor_Curve.hxx>
  81 #include <BRepAdaptor_Surface.hxx>
  82 #include <BRepBndLib.hxx>
  83 #include <BRepBuilderAPI_MakeFace.hxx>
  84 #include <BRepMesh_IncrementalMesh.hxx>
  85
  86 #include <TopAbs.hxx>
  87 #include <TopExp.hxx>
  88 #include <TopExp_Explorer.hxx>
  89 #include <TopLoc_Location.hxx>
  90 #include <TopoDS.hxx>
  91 #include <TopoDS_Shape.hxx>
  92 #include <TopoDS_Solid.hxx>
  93 #include <TopoDS_Face.hxx>
  94 #include <TopoDS_Edge.hxx>
  95 #include <TopoDS_Vertex.hxx>
  96 #include <TopoDS_Compound.hxx>
  97 #include <TopoDS_Iterator.hxx>
  98 #include <TopTools_Array1OfShape.hxx>
  99 #include <TopTools_DataMapIteratorOfDataMapOfShapeListOfShape.hxx>
 100 #include <TopTools_IndexedMapOfShape.hxx>
 101 #include <TopTools_ListIteratorOfListOfShape.hxx>
 102 #include <TopTools_MapOfShape.hxx>
 103 #include <TopTools_MapOfOrientedShape.hxx>
 104
 105 #include <Geom_Surface.hxx>
 106 #include <Geom_Plane.hxx>
 107 #include <Geom_SphericalSurface.hxx>
 108 #include <Geom_CylindricalSurface.hxx>
 109 #include <GeomAdaptor_Surface.hxx>
 110
 111 #include <GeomLib_Tool.hxx>
 112 #include <Geom2d_Curve.hxx>
 113
 114 #include <Bnd_Box.hxx>
 115 #include <GProp_GProps.hxx>
 116 #include <TColStd_Array1OfReal.hxx>
 117 #include <TColStd_HArray1OfInteger.hxx>
 118 #include <TColStd_ListIteratorOfListOfInteger.hxx>
 119 #include <TColStd_ListOfInteger.hxx>
 120 #include <gp_Cylinder.hxx>
 121 #include <gp_Lin.hxx>
 122 #include <gp_Pnt.hxx>
 123
 124 #include <vector>
 125 #include <algorithm>
 126 #include <functional>
 127
 128 #include <Standard_NullObject.hxx>
 129 #include <Standard_Failure.hxx>
 130 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
 131
 132 // Includes added for GetInPlace algorithm improvement
 133
 134 #include <GEOMImpl_MeasureDriver.hxx>
 135 #include <GEOMImpl_IMeasure.hxx>
 136 #include <BRepBuilderAPI_MakeVertex.hxx>
 137
 138 #include <BRepClass_FaceClassifier.hxx>
 139 #include <BRepClass3d_SolidClassifier.hxx>
 140 #include <Precision.hxx>
 141
 142 #define STD_SORT_ALGO 1
 143
 144 //=============================================================================
 145 /*!
 146  *   constructor:
 147  */
 148 //=============================================================================
 149 GEOMImpl_IShapesOperations::GEOMImpl_IShapesOperations (GEOM_Engine* theEngine, int theDocID)
 150 : GEOM_IOperations(theEngine, theDocID)
 151 {
 152   MESSAGE("GEOMImpl_IShapesOperations::GEOMImpl_IShapesOperations");
 153 }
 154
 155 //=============================================================================
 156 /*!
 157  *  destructor
 158  */
 159 //=============================================================================
 160 GEOMImpl_IShapesOperations::~GEOMImpl_IShapesOperations()
 161 {
 162   MESSAGE("GEOMImpl_IShapesOperations::~GEOMImpl_IShapesOperations");
 163 }
 164
 165 //=============================================================================
 166 /*!
 167  *  MakeEdge
 168  */
 169 //=============================================================================
 170 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeEdge
 171                      (Handle(GEOM_Object) thePnt1, Handle(GEOM_Object) thePnt2)
 172 {
 173   SetErrorCode(KO);
 174
 175   if (thePnt1.IsNull() || thePnt2.IsNull()) return NULL;
 176
 177   //Add a new Edge object
 178   Handle(GEOM_Object) anEdge = GetEngine()->AddObject(GetDocID(), GEOM_EDGE);
 179
 180   //Add a new Vector function
 181   Handle(GEOM_Function) aFunction =
 182     anEdge->AddFunction(GEOMImpl_VectorDriver::GetID(), VECTOR_TWO_PNT);
 183
 184   //Check if the function is set correctly
 185   if (aFunction->GetDriverGUID() != GEOMImpl_VectorDriver::GetID()) return NULL;
 186
 187   GEOMImpl_IVector aPI (aFunction);
 188
 189   Handle(GEOM_Function) aRef1 = thePnt1->GetLastFunction();
 190   Handle(GEOM_Function) aRef2 = thePnt2->GetLastFunction();
 191   if (aRef1.IsNull() || aRef2.IsNull()) return NULL;
 192
 193   aPI.SetPoint1(aRef1);
 194   aPI.SetPoint2(aRef2);
 195
 196   //Compute the Edge value
 197   try {
 198 #if OCC_VERSION_LARGE > 0x06010000
 199     OCC_CATCH_SIGNALS;
 200 #endif
 201     if (!GetSolver()->ComputeFunction(aFunction)) {
 202       SetErrorCode("Vector driver failed");
 203       return NULL;
 204     }
 205   }
 206   catch (Standard_Failure) {
 207     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 208     SetErrorCode(aFail->GetMessageString());
 209     return NULL;
 210   }
 211
 212   //Make a Python command
 213   GEOM::TPythonDump(aFunction) << anEdge << " = geompy.MakeEdge("
 214                                << thePnt1 << ", " << thePnt2 << ")";
 215
 216   SetErrorCode(OK);
 217   return anEdge;
 218 }
 219
 220 //=============================================================================
 221 /*!
 222  *  MakeEdgeOnCurveByLength
 223  */
 224 //=============================================================================
 225 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeEdgeOnCurveByLength
 226                      (Handle(GEOM_Object) theRefCurve,
 227                       const Standard_Real theLength,
 228                       Handle(GEOM_Object) theStartPoint)
 229 {
 230   SetErrorCode(KO);
 231
 232   if (theRefCurve.IsNull()) return NULL;
 233
 234   //Add a new Edge object
 235   Handle(GEOM_Object) anEdge = GetEngine()->AddObject(GetDocID(), GEOM_EDGE);
 236
 237   //Add a new Vector function
 238   Handle(GEOM_Function) aFunction =
 239     anEdge->AddFunction(GEOMImpl_ShapeDriver::GetID(), EDGE_CURVE_LENGTH);
 240
 241   //Check if the function is set correctly
 242   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 243
 244   GEOMImpl_IVector aPI (aFunction);
 245
 246   Handle(GEOM_Function) aRef1 = theRefCurve->GetLastFunction();
 247   if (aRef1.IsNull()) return NULL;
 248   aPI.SetPoint1(aRef1);
 249
 250   if (!theStartPoint.IsNull()) {
 251     Handle(GEOM_Function) aRef2 = theStartPoint->GetLastFunction();
 252     aPI.SetPoint2(aRef2);
 253   }
 254
 255   aPI.SetParameter(theLength);
 256
 257   //Compute the Edge value
 258   try {
 259 #if OCC_VERSION_LARGE > 0x06010000
 260     OCC_CATCH_SIGNALS;
 261 #endif
 262     if (!GetSolver()->ComputeFunction(aFunction)) {
 263       SetErrorCode("Vector driver failed");
 264       return NULL;
 265     }
 266   }
 267   catch (Standard_Failure) {
 268     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 269     SetErrorCode(aFail->GetMessageString());
 270     return NULL;
 271   }
 272
 273   //Make a Python command
 274   GEOM::TPythonDump(aFunction) << anEdge << " = geompy.MakeEdgeOnCurveByLength("
 275                                << theRefCurve << ", " << theLength << ", " << theStartPoint << ")";
 276
 277   SetErrorCode(OK);
 278   return anEdge;
 279 }
 280
 281 //=============================================================================
 282 /*!
 283  *  MakeEdgeWire
 284  */
 285 //=============================================================================
 286 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeEdgeWire
 287                     (Handle(GEOM_Object) theWire,
 288                      const Standard_Real theLinearTolerance,
 289                      const Standard_Real theAngularTolerance)
 290 {
 291   SetErrorCode(KO);
 292
 293   if (theWire.IsNull()) return NULL;
 294
 295   //Add a new Edge object
 296   Handle(GEOM_Object) anEdge = GetEngine()->AddObject(GetDocID(), GEOM_EDGE);
 297
 298   //Add a new Vector function
 299   Handle(GEOM_Function) aFunction =
 300     anEdge->AddFunction(GEOMImpl_ShapeDriver::GetID(), EDGE_WIRE);
 301
 302   //Check if the function is set correctly
 303   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 304
 305   GEOMImpl_IShapes aCI (aFunction);
 306
 307   Handle(GEOM_Function) aWire = theWire->GetLastFunction();
 308
 309   if (aWire.IsNull()) return NULL;
 310
 311   aCI.SetBase(aWire);
 312   aCI.SetTolerance(theLinearTolerance);
 313   aCI.SetAngularTolerance(theAngularTolerance);
 314
 315   //Compute the Edge value
 316   try {
 317 #if OCC_VERSION_LARGE > 0x06010000
 318     OCC_CATCH_SIGNALS;
 319 #endif
 320     if (!GetSolver()->ComputeFunction(aFunction)) {
 321       SetErrorCode("Shape driver failed");
 322       return NULL;
 323     }
 324   }
 325   catch (Standard_Failure) {
 326     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 327     SetErrorCode(aFail->GetMessageString());
 328     return NULL;
 329   }
 330
 331   const double DEF_LIN_TOL = Precision::Confusion();
 332   const double DEF_ANG_TOL = Precision::Angular();
 333   //Make a Python command
 334   if ( theAngularTolerance == DEF_ANG_TOL ) {
 335     if ( theLinearTolerance == DEF_LIN_TOL )
 336       GEOM::TPythonDump(aFunction) << anEdge  << " = geompy.MakeEdgeWire("
 337                                    << theWire << ")";
 338     else
 339       GEOM::TPythonDump(aFunction) << anEdge  << " = geompy.MakeEdgeWire("
 340                                    << theWire << ", " << theLinearTolerance << ")";
 341   }
 342   else {
 343     GEOM::TPythonDump(aFunction) << anEdge  << " = geompy.MakeEdgeWire("
 344                                  << theWire << ", " << theLinearTolerance << ", "
 345                                  << theAngularTolerance << ")";
 346   }
 347
 348   SetErrorCode(OK);
 349   return anEdge;
 350 }
 351
 352 //=============================================================================
 353 /*!
 354  *  MakeWire
 355  */
 356 //=============================================================================
 357 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeWire
 358                              (std::list<Handle(GEOM_Object)> theShapes,
 359                               const Standard_Real            theTolerance)
 360 {
 361   SetErrorCode(KO);
 362
 363   //Add a new object
 364   Handle(GEOM_Object) aWire = GetEngine()->AddObject(GetDocID(), GEOM_WIRE);
 365
 366   //Add a new function
 367   Handle(GEOM_Function) aFunction =
 368     aWire->AddFunction(GEOMImpl_ShapeDriver::GetID(), WIRE_EDGES);
 369   if (aFunction.IsNull()) return NULL;
 370
 371   //Check if the function is set correctly
 372   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 373
 374   GEOMImpl_IShapes aCI (aFunction);
 375   aCI.SetTolerance(theTolerance);
 376
 377   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 378
 379   // Shapes
 380   std::list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 381   for (; it != theShapes.end(); it++) {
 382     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 383     if (aRefSh.IsNull()) {
 384       SetErrorCode("NULL argument shape for the shape construction");
 385       return NULL;
 386     }
 387     aShapesSeq->Append(aRefSh);
 388   }
 389   aCI.SetShapes(aShapesSeq);
 390
 391   //Compute the shape
 392   try {
 393 #if OCC_VERSION_LARGE > 0x06010000
 394     OCC_CATCH_SIGNALS;
 395 #endif
 396     if (!GetSolver()->ComputeFunction(aFunction)) {
 397       SetErrorCode("Shape driver failed");
 398       return NULL;
 399     }
 400   }
 401   catch (Standard_Failure) {
 402     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 403     SetErrorCode(aFail->GetMessageString());
 404     return NULL;
 405   }
 406
 407   //Make a Python command
 408   GEOM::TPythonDump pd (aFunction);
 409   pd << aWire << " = geompy.MakeWire([";
 410
 411   // Shapes
 412   it = theShapes.begin();
 413   if (it != theShapes.end()) {
 414     pd << (*it++);
 415     while (it != theShapes.end()) {
 416       pd << ", " << (*it++);
 417     }
 418   }
 419   pd << "], " << theTolerance << ")";
 420
 421   SetErrorCode(OK);
 422   return aWire;
 423 }
 424
 425 //=============================================================================
 426 /*!
 427  *  MakeFace
 428  */
 429 //=============================================================================
 430 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeFace (Handle(GEOM_Object) theWire,
 431                                                           const bool isPlanarWanted)
 432 {
 433   SetErrorCode(KO);
 434
 435   if (theWire.IsNull()) return NULL;
 436
 437   //Add a new Face object
 438   Handle(GEOM_Object) aFace = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 439
 440   //Add a new Shape function for creation of a face from a wire
 441   Handle(GEOM_Function) aFunction =
 442     aFace->AddFunction(GEOMImpl_ShapeDriver::GetID(), FACE_WIRE);
 443   if (aFunction.IsNull()) return NULL;
 444
 445   //Check if the function is set correctly
 446   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 447
 448   GEOMImpl_IShapes aCI (aFunction);
 449
 450   Handle(GEOM_Function) aRefWire = theWire->GetLastFunction();
 451
 452   if (aRefWire.IsNull()) return NULL;
 453
 454   aCI.SetBase(aRefWire);
 455   aCI.SetIsPlanar(isPlanarWanted);
 456
 457   //Compute the Face value
 458   try {
 459 #if OCC_VERSION_LARGE > 0x06010000
 460     OCC_CATCH_SIGNALS;
 461 #endif
 462     if (!GetSolver()->ComputeFunction(aFunction)) {
 463       SetErrorCode("Shape driver failed to compute a face");
 464       return NULL;
 465     }
 466   }
 467   catch (Standard_Failure) {
 468     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 469     SetErrorCode(aFail->GetMessageString());
 470     return NULL;
 471   }
 472
 473   //Make a Python command
 474   GEOM::TPythonDump(aFunction) << aFace << " = geompy.MakeFace("
 475     << theWire << ", " << (int)isPlanarWanted << ")";
 476
 477   SetErrorCode(OK);
 478   return aFace;
 479 }
 480
 481 //=============================================================================
 482 /*!
 483  *  MakeFaceWires
 484  */
 485 //=============================================================================
 486 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeFaceWires
 487                              (std::list<Handle(GEOM_Object)> theShapes,
 488                               const bool isPlanarWanted)
 489 {
 490   SetErrorCode(KO);
 491
 492   //Add a new object
 493   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 494
 495   //Add a new function
 496   Handle(GEOM_Function) aFunction =
 497     aShape->AddFunction(GEOMImpl_ShapeDriver::GetID(), FACE_WIRES);
 498   if (aFunction.IsNull()) return NULL;
 499
 500   //Check if the function is set correctly
 501   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 502
 503   GEOMImpl_IShapes aCI (aFunction);
 504
 505   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 506
 507   // Shapes
 508   std::list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 509   for (; it != theShapes.end(); it++) {
 510     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 511     if (aRefSh.IsNull()) {
 512       SetErrorCode("NULL argument shape for the face construction");
 513       return NULL;
 514     }
 515     aShapesSeq->Append(aRefSh);
 516   }
 517   aCI.SetShapes(aShapesSeq);
 518
 519   aCI.SetIsPlanar(isPlanarWanted);
 520
 521   //Compute the shape
 522   try {
 523 #if OCC_VERSION_LARGE > 0x06010000
 524     OCC_CATCH_SIGNALS;
 525 #endif
 526     if (!GetSolver()->ComputeFunction(aFunction)) {
 527       SetErrorCode("Shape driver failed");
 528       return NULL;
 529     }
 530   }
 531   catch (Standard_Failure) {
 532     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 533     SetErrorCode(aFail->GetMessageString());
 534     return NULL;
 535   }
 536
 537   //Make a Python command
 538   GEOM::TPythonDump pd (aFunction);
 539   pd << aShape << " = geompy.MakeFaceWires([";
 540
 541   // Shapes
 542   it = theShapes.begin();
 543   if (it != theShapes.end()) {
 544     pd << (*it++);
 545     while (it != theShapes.end()) {
 546       pd << ", " << (*it++);
 547     }
 548   }
 549   pd << "], " << (int)isPlanarWanted << ")";
 550
 551   SetErrorCode(OK);
 552   return aShape;
 553 }
 554
 555 //=============================================================================
 556 /*!
 557  *  MakeShell
 558  */
 559 //=============================================================================
 560 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeShell
 561                              (std::list<Handle(GEOM_Object)> theShapes)
 562 {
 563   return MakeShape(theShapes, GEOM_SHELL, SHELL_FACES, "MakeShell");
 564 }
 565
 566 //=============================================================================
 567 /*!
 568  *  MakeSolidShells
 569  */
 570 //=============================================================================
 571 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeSolidShells
 572                              (std::list<Handle(GEOM_Object)> theShapes)
 573 {
 574   return MakeShape(theShapes, GEOM_SOLID, SOLID_SHELLS, "MakeSolid");
 575 }
 576
 577 //=============================================================================
 578 /*!
 579  *  MakeCompound
 580  */
 581 //=============================================================================
 582 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeCompound
 583                              (std::list<Handle(GEOM_Object)> theShapes)
 584 {
 585   return MakeShape(theShapes, GEOM_COMPOUND, COMPOUND_SHAPES, "MakeCompound");
 586 }
 587
 588 //=============================================================================
 589 /*!
 590  *  MakeShape
 591  */
 592 //=============================================================================
 593 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeShape
 594                              (std::list<Handle(GEOM_Object)> theShapes,
 595                               const Standard_Integer         theObjectType,
 596                               const Standard_Integer         theFunctionType,
 597                               const TCollection_AsciiString& theMethodName)
 598 {
 599   SetErrorCode(KO);
 600
 601   //Add a new object
 602   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), theObjectType);
 603
 604   //Add a new function
 605   Handle(GEOM_Function) aFunction =
 606     aShape->AddFunction(GEOMImpl_ShapeDriver::GetID(), theFunctionType);
 607   if (aFunction.IsNull()) return NULL;
 608
 609   //Check if the function is set correctly
 610   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 611
 612   GEOMImpl_IShapes aCI (aFunction);
 613
 614   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 615
 616   // Shapes
 617   std::list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 618   for (; it != theShapes.end(); it++) {
 619     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 620     if (aRefSh.IsNull()) {
 621       SetErrorCode("NULL argument shape for the shape construction");
 622       return NULL;
 623     }
 624     aShapesSeq->Append(aRefSh);
 625   }
 626   aCI.SetShapes(aShapesSeq);
 627
 628   //Compute the shape
 629   try {
 630 #if OCC_VERSION_LARGE > 0x06010000
 631     OCC_CATCH_SIGNALS;
 632 #endif
 633     if (!GetSolver()->ComputeFunction(aFunction)) {
 634       SetErrorCode("Shape driver failed");
 635       return NULL;
 636     }
 637   }
 638   catch (Standard_Failure) {
 639     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 640     SetErrorCode(aFail->GetMessageString());
 641     return NULL;
 642   }
 643
 644   //Make a Python command
 645   GEOM::TPythonDump pd (aFunction);
 646   pd << aShape << " = geompy." << theMethodName.ToCString() << "([";
 647
 648   // Shapes
 649   it = theShapes.begin();
 650   if (it != theShapes.end()) {
 651     pd << (*it++);
 652     while (it != theShapes.end()) {
 653       pd << ", " << (*it++);
 654     }
 655   }
 656   pd << "])";
 657
 658   SetErrorCode(OK);
 659   return aShape;
 660 }
 661
 662 //=============================================================================
 663 /*!
 664  *  MakeGlueFaces
 665  */
 666 //=============================================================================
 667 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueFaces
 668                                                 (Handle(GEOM_Object) theShape,
 669                                                  const Standard_Real theTolerance,
 670                                                  const Standard_Boolean doKeepNonSolids)
 671 {
 672   SetErrorCode(KO);
 673
 674   if (theShape.IsNull()) return NULL;
 675
 676   //Add a new Glued object
 677   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
 678
 679   //Add a new Glue function
 680   Handle(GEOM_Function) aFunction;
 681   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_FACES);
 682   if (aFunction.IsNull()) return NULL;
 683
 684   //Check if the function is set correctly
 685   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
 686
 687   GEOMImpl_IGlue aCI (aFunction);
 688
 689   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 690   if (aRefShape.IsNull()) return NULL;
 691
 692   aCI.SetBase(aRefShape);
 693   aCI.SetTolerance(theTolerance);
 694   aCI.SetKeepNonSolids(doKeepNonSolids);
 695
 696   //Compute the sub-shape value
 697   Standard_Boolean isWarning = Standard_False;
 698   try {
 699 #if OCC_VERSION_LARGE > 0x06010000
 700     OCC_CATCH_SIGNALS;
 701 #endif
 702     if (!GetSolver()->ComputeFunction(aFunction)) {
 703       SetErrorCode("Shape driver failed to glue faces");
 704       return NULL;
 705     }
 706   }
 707   catch (Standard_Failure) {
 708     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 709     SetErrorCode(aFail->GetMessageString());
 710     // to provide warning
 711     if (!aFunction->GetValue().IsNull()) {
 712       isWarning = Standard_True;
 713     } else {
 714       return NULL;
 715     }
 716   }
 717
 718   //Make a Python command
 719   GEOM::TPythonDump(aFunction) << aGlued << " = geompy.MakeGlueFaces("
 720     << theShape << ", " << theTolerance << ")";
 721
 722   // to provide warning
 723   if (!isWarning) SetErrorCode(OK);
 724   return aGlued;
 725 }
 726
 727 //=============================================================================
 728 /*!
 729  *  GetGlueFaces
 730  */
 731 //=============================================================================
 732 /*
 733 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetGlueFaces
 734                                                 (Handle(GEOM_Object) theShape,
 735                                                  const Standard_Real theTolerance)
 736 {
 737   SetErrorCode(KO);
 738
 739   if (theShape.IsNull()) return NULL;
 740   TopoDS_Shape aShape = theShape->GetValue();
 741   if (aShape.IsNull()) return NULL;
 742
 743   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
 744
 745   Standard_Integer iErr;
 746   TopoDS_Shape aS;
 747   GEOMAlgo_Gluer1 aGluer;
 748   GEOMAlgo_ListIteratorOfListOfCoupleOfShapes aItCS;
 749   GEOMAlgo_CoupleOfShapes aCS;
 750   GEOMAlgo_ListOfCoupleOfShapes aLCS;
 751
 752   //aGluer = new GEOMAlgo_Gluer1;
 753   aGluer.SetShape(aShape);
 754   aGluer.SetTolerance(theTolerance);
 755   aGluer.Perform();
 756   iErr = aGluer.ErrorStatus();
 757   if (iErr) return NULL;
 758
 759   TopTools_ListOfShape listShape;
 760   const GEOMAlgo_ListOfCoupleOfShapes& aLCSG = aGluer.GluedFaces();
 761   // Access to faces
 762   aItCS.Initialize(aLCSG);
 763   for (; aItCS.More(); aItCS.Next()) {
 764     const GEOMAlgo_CoupleOfShapes& aCSG = aItCS.Value();
 765     listShape.Append(aCSG.Shape1());
 766   }
 767
 768   TopTools_ListIteratorOfListOfShape itSub (listShape);
 769   TCollection_AsciiString anAsciiList, anEntry;
 770   TopTools_IndexedMapOfShape anIndices;
 771   TopExp::MapShapes(aShape, anIndices);
 772   Handle(TColStd_HArray1OfInteger) anArray;
 773   Handle(GEOM_Object) anObj;
 774   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
 775     TopoDS_Shape aValue = itSub.Value();
 776     anArray = new TColStd_HArray1OfInteger(1,1);
 777     anArray->SetValue(1, anIndices.FindIndex(aValue));
 778     anObj = GetEngine()->AddSubShape(theShape, anArray);
 779     if (!anObj.IsNull()) {
 780       aSeq->Append(anObj);
 781
 782       // for python command
 783       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
 784       anAsciiList += anEntry;
 785       anAsciiList += ",";
 786     }
 787   }
 788
 789   //Make a Python command
 790   if( anAsciiList.Length() > 0 ) {
 791     anAsciiList.Trunc(anAsciiList.Length() - 1);
 792     Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 793     GEOM::TPythonDump pd (aFunction, true);
 794     pd << "[" << anAsciiList.ToCString();
 795     pd << "] = geompy.GetGlueFaces(" << theShape << ", " << theTolerance << ")";
 796   }
 797
 798   SetErrorCode(OK);
 799
 800   return aSeq;
 801 }
 802 */
 803
 804 //=============================================================================
 805 /*!
 806  *  MakeGlueFacesByList
 807  */
 808 //=============================================================================
 809 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueFacesByList
 810                                                 (Handle(GEOM_Object) theShape,
 811                                                  const Standard_Real theTolerance,
 812                                                  std::list<Handle(GEOM_Object)> theFaces,
 813                                                  const Standard_Boolean doKeepNonSolids,
 814                                                  const Standard_Boolean doGlueAllEdges)
 815 {
 816   SetErrorCode(KO);
 817
 818   if (theShape.IsNull()) return NULL;
 819
 820   //Add a new Glued object
 821   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
 822
 823   //Add a new Glue function
 824   Handle(GEOM_Function) aFunction;
 825   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_FACES_BY_LIST);
 826   if (aFunction.IsNull()) return NULL;
 827
 828   //Check if the function is set correctly
 829   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
 830
 831   GEOMImpl_IGlue aCI (aFunction);
 832
 833   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 834   if (aRefShape.IsNull()) return NULL;
 835
 836   aCI.SetBase(aRefShape);
 837   aCI.SetTolerance(theTolerance);
 838   aCI.SetKeepNonSolids(doKeepNonSolids);
 839   aCI.SetGlueAllEdges(doGlueAllEdges);
 840
 841   Handle(TColStd_HSequenceOfTransient) aFaces = new TColStd_HSequenceOfTransient;
 842   std::list<Handle(GEOM_Object)>::iterator it = theFaces.begin();
 843   for (; it != theFaces.end(); it++) {
 844     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 845     if (aRefSh.IsNull()) {
 846       SetErrorCode("NULL argument shape for the shape construction");
 847       return NULL;
 848     }
 849     aFaces->Append(aRefSh);
 850   }
 851   aCI.SetFaces(aFaces);
 852
 853   //Compute the sub-shape value
 854   Standard_Boolean isWarning = Standard_False;
 855   try {
 856 #if OCC_VERSION_LARGE > 0x06010000
 857     OCC_CATCH_SIGNALS;
 858 #endif
 859     if (!GetSolver()->ComputeFunction(aFunction)) {
 860       SetErrorCode("Shape driver failed to glue faces");
 861       return NULL;
 862     }
 863   }
 864   catch (Standard_Failure) {
 865     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 866     SetErrorCode(aFail->GetMessageString());
 867     // to provide warning
 868     if (!aFunction->GetValue().IsNull()) {
 869       isWarning = Standard_True;
 870     } else {
 871       return NULL;
 872     }
 873   }
 874
 875   //Make a Python command
 876
 877   GEOM::TPythonDump pd(aFunction);
 878   pd << aGlued << " = geompy.MakeGlueFacesByList("
 879      << theShape << ", " << theTolerance << ", [";
 880   // Faces
 881   it = theFaces.begin();
 882   if (it != theFaces.end()) {
 883     pd << (*it++);
 884     while (it != theFaces.end()) {
 885       pd << ", " << (*it++);
 886     }
 887   }
 888   pd << "], " << (bool)doKeepNonSolids << ", " << (bool)doGlueAllEdges << ")";
 889
 890   // to provide warning
 891   if (!isWarning) SetErrorCode(OK);
 892   return aGlued;
 893 }
 894
 895 //=============================================================================
 896 /*!
 897  *  MakeGlueEdges
 898  */
 899 //=============================================================================
 900 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueEdges
 901                                                 (Handle(GEOM_Object) theShape,
 902                                                  const Standard_Real theTolerance)
 903 {
 904   SetErrorCode(KO);
 905
 906   if (theShape.IsNull()) return NULL;
 907
 908   //Add a new Glued object
 909   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
 910
 911   //Add a new Glue function
 912   Handle(GEOM_Function) aFunction;
 913   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_EDGES);
 914   if (aFunction.IsNull()) return NULL;
 915
 916   //Check if the function is set correctly
 917   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
 918
 919   GEOMImpl_IGlue aCI (aFunction);
 920
 921   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 922   if (aRefShape.IsNull()) return NULL;
 923
 924   aCI.SetBase(aRefShape);
 925   aCI.SetTolerance(theTolerance);
 926   aCI.SetKeepNonSolids(true);
 927
 928   //Compute the sub-shape value
 929   Standard_Boolean isWarning = Standard_False;
 930   try {
 931 #if OCC_VERSION_LARGE > 0x06010000
 932     OCC_CATCH_SIGNALS;
 933 #endif
 934     if (!GetSolver()->ComputeFunction(aFunction)) {
 935       SetErrorCode("Shape driver failed to glue edges");
 936       return NULL;
 937     }
 938   }
 939   catch (Standard_Failure) {
 940     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 941     SetErrorCode(aFail->GetMessageString());
 942     // to provide warning
 943     if (!aFunction->GetValue().IsNull()) {
 944       isWarning = Standard_True;
 945     } else {
 946       return NULL;
 947     }
 948   }
 949
 950   //Make a Python command
 951   GEOM::TPythonDump(aFunction) << aGlued << " = geompy.MakeGlueEdges("
 952     << theShape << ", " << theTolerance << ")";
 953
 954   // to provide warning
 955   if (!isWarning) SetErrorCode(OK);
 956   return aGlued;
 957 }
 958
 959 //=============================================================================
 960 /*!
 961  *  GetGlueShapes
 962  */
 963 //=============================================================================
 964 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetGlueShapes
 965                                                 (Handle(GEOM_Object)    theShape,
 966                                                  const Standard_Real    theTolerance,
 967                                                  const TopAbs_ShapeEnum theType)
 968 {
 969   SetErrorCode(KO);
 970
 971   if (theShape.IsNull()) return NULL;
 972   TopoDS_Shape aShape = theShape->GetValue();
 973   if (aShape.IsNull()) return NULL;
 974
 975   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
 976
 977   GEOMAlgo_GlueDetector aGluer;
 978   aGluer.SetArgument(aShape);
 979   aGluer.SetTolerance(theTolerance);
 980   aGluer.Perform();
 981   Standard_Integer iErr = aGluer.ErrorStatus();
 982   if (iErr) return NULL;
 983
 984   TCollection_AsciiString anAsciiList, anEntry;
 985   TopTools_IndexedMapOfShape anIndices;
 986   TopExp::MapShapes(aShape, anIndices);
 987   Handle(TColStd_HArray1OfInteger) anArray;
 988   Handle(GEOM_Object) anObj;
 989
 990   TopTools_ListOfShape listOnePerSet;
 991
 992   const TopTools_DataMapOfShapeListOfShape& aImages = aGluer.Images();
 993   TopTools_DataMapIteratorOfDataMapOfShapeListOfShape aItDMSLS (aImages);
 994   for (int index = 1; aItDMSLS.More(); aItDMSLS.Next(), ++index) {
 995     // some key shape
 996     //const TopoDS_Shape& aSkey = aItDMSLS.Key();
 997
 998     // list of shapes of the argument that can be glued
 999     const TopTools_ListOfShape& aLSD = aItDMSLS.Value();
1000
1001     //listShape.Append(aLSD.First());
1002     TopoDS_Shape aValue = aLSD.First();
1003
1004     if (aValue.ShapeType() == theType) {
1005       listOnePerSet.Append(aValue);
1006     }
1007   }
1008
1009   // for stable order of returned entities
1010   GEOMImpl_IShapesOperations::SortShapes(listOnePerSet, Standard_False);
1011
1012   TopTools_ListIteratorOfListOfShape aListIt (listOnePerSet);
1013   for (; aListIt.More(); aListIt.Next()) {
1014     TopoDS_Shape aValue = aListIt.Value();
1015     anArray = new TColStd_HArray1OfInteger(1,1);
1016     anArray->SetValue(1, anIndices.FindIndex(aValue));
1017     anObj = GetEngine()->AddSubShape(theShape, anArray);
1018     if (!anObj.IsNull()) {
1019       aSeq->Append(anObj);
1020
1021       // for python command
1022       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1023       anAsciiList += anEntry;
1024       anAsciiList += ",";
1025     }
1026   }
1027
1028   // Make a Python command
1029   if (anAsciiList.Length() > 0) {
1030     anAsciiList.Trunc(anAsciiList.Length() - 1);
1031     Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1032     GEOM::TPythonDump pd (aFunction, /*append=*/true);
1033     pd << "[" << anAsciiList.ToCString();
1034     if (theType == TopAbs_FACE)
1035       pd << "] = geompy.GetGlueFaces(" << theShape << ", " << theTolerance << ")";
1036     else if (theType == TopAbs_EDGE)
1037       pd << "] = geompy.GetGlueEdges(" << theShape << ", " << theTolerance << ")";
1038   }
1039
1040   SetErrorCode(OK);
1041
1042   return aSeq;
1043 }
1044
1045 //=============================================================================
1046 /*!
1047  *  MakeGlueEdgesByList
1048  */
1049 //=============================================================================
1050 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueEdgesByList
1051                                                 (Handle(GEOM_Object) theShape,
1052                                                  const Standard_Real theTolerance,
1053                                                  std::list<Handle(GEOM_Object)> theEdges)
1054 {
1055   SetErrorCode(KO);
1056
1057   if (theShape.IsNull()) return NULL;
1058
1059   //Add a new Glued object
1060   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
1061
1062   //Add a new Glue function
1063   Handle(GEOM_Function) aFunction;
1064   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_EDGES_BY_LIST);
1065   if (aFunction.IsNull()) return NULL;
1066
1067   //Check if the function is set correctly
1068   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
1069
1070   GEOMImpl_IGlue aCI (aFunction);
1071
1072   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1073   if (aRefShape.IsNull()) return NULL;
1074
1075   aCI.SetBase(aRefShape);
1076   aCI.SetTolerance(theTolerance);
1077   aCI.SetKeepNonSolids(true);
1078
1079   Handle(TColStd_HSequenceOfTransient) anEdges = new TColStd_HSequenceOfTransient;
1080   std::list<Handle(GEOM_Object)>::iterator it = theEdges.begin();
1081   for (; it != theEdges.end(); it++) {
1082     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
1083     if (aRefSh.IsNull()) {
1084       SetErrorCode("NULL argument shape for the shape construction");
1085       return NULL;
1086     }
1087     anEdges->Append(aRefSh);
1088   }
1089   aCI.SetFaces(anEdges);
1090
1091   //Compute the sub-shape value
1092   Standard_Boolean isWarning = Standard_False;
1093   try {
1094 #if OCC_VERSION_LARGE > 0x06010000
1095     OCC_CATCH_SIGNALS;
1096 #endif
1097     if (!GetSolver()->ComputeFunction(aFunction)) {
1098       SetErrorCode("Shape driver failed to glue edges");
1099       return NULL;
1100     }
1101   }
1102   catch (Standard_Failure) {
1103     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1104     SetErrorCode(aFail->GetMessageString());
1105     // to provide warning
1106     if (!aFunction->GetValue().IsNull()) {
1107       isWarning = Standard_True;
1108     } else {
1109       return NULL;
1110     }
1111   }
1112
1113   //Make a Python command
1114
1115   GEOM::TPythonDump pd (aFunction);
1116   pd << aGlued << " = geompy.MakeGlueEdgesByList("
1117      << theShape << ", " << theTolerance << ", [";
1118   // Edges
1119   it = theEdges.begin();
1120   if (it != theEdges.end()) {
1121     pd << (*it++);
1122     while (it != theEdges.end()) {
1123       pd << ", " << (*it++);
1124     }
1125   }
1126   pd << "])";
1127
1128   // to provide warning
1129   if (!isWarning) SetErrorCode(OK);
1130   return aGlued;
1131 }
1132
1133 //=============================================================================
1134 /*!
1135  *  GetExistingSubObjects
1136  */
1137 //=============================================================================
1138 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetExistingSubObjects
1139                                           (Handle(GEOM_Object)    theShape,
1140                                            const Standard_Boolean theGroupsOnly)
1141 {
1142   SetErrorCode(KO);
1143
1144   if (theShape.IsNull()) return NULL;
1145
1146   Handle(GEOM_Function) aMainShape = theShape->GetLastFunction();
1147   if (aMainShape.IsNull()) return NULL;
1148
1149   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1150   SetErrorCode(NOT_FOUND_ANY);
1151
1152   if (!aMainShape->HasSubShapeReferences()) return aSeq;
1153   const TDataStd_ListOfExtendedString& aListEntries = aMainShape->GetSubShapeReferences();
1154   if (aListEntries.IsEmpty()) return aSeq;
1155
1156   SetErrorCode(KO);
1157
1158   TCollection_AsciiString anAsciiList;
1159
1160   TDataStd_ListIteratorOfListOfExtendedString anIt (aListEntries);
1161   for (; anIt.More(); anIt.Next()) {
1162     TCollection_ExtendedString anEntry = anIt.Value();
1163     Standard_Integer aStrLen = anEntry.LengthOfCString();
1164     char* anEntryStr = new char[aStrLen];
1165     anEntry.ToUTF8CString(anEntryStr);
1166     Handle(GEOM_Object) anObj = GetEngine()->GetObject(GetDocID(), anEntryStr, false);
1167     if (!anObj.IsNull()) {
1168       if (!theGroupsOnly || anObj->GetType() == GEOM_GROUP) {
1169         aSeq->Append(anObj);
1170
1171         // for python command
1172         anAsciiList += anEntryStr;
1173         anAsciiList += ",";
1174       }
1175     }
1176     delete [] anEntryStr;
1177   }
1178
1179   if (aSeq->Length() == 0) {
1180     SetErrorCode(NOT_FOUND_ANY);
1181     return aSeq;
1182   }
1183
1184   //Make a Python command
1185   anAsciiList.Trunc(anAsciiList.Length() - 1);
1186
1187   GEOM::TPythonDump pd (aMainShape, /*append=*/true);
1188   pd << "[" << anAsciiList.ToCString();
1189   pd << "] = geompy.GetExistingSubObjects(";
1190   pd << theShape << ", " << (bool)theGroupsOnly << ")";
1191
1192   SetErrorCode(OK);
1193
1194   return aSeq;
1195 }
1196
1197 //=============================================================================
1198 /*!
1199  *  MakeExplode
1200  */
1201 //=============================================================================
1202 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::MakeExplode
1203                                           (Handle(GEOM_Object)    theShape,
1204                                            const Standard_Integer theShapeType,
1205                                            const Standard_Boolean isSorted,
1206                                            const ExplodeType      theExplodeType)
1207 {
1208   SetErrorCode(KO);
1209
1210   if (theShape.IsNull()) return NULL;
1211   TopoDS_Shape aShape = theShape->GetValue();
1212   if (aShape.IsNull()) return NULL;
1213
1214   Handle(GEOM_Function) aMainShape = theShape->GetLastFunction();
1215
1216   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1217   Handle(GEOM_Object) anObj;
1218   TopTools_MapOfShape mapShape;
1219   TopTools_ListOfShape listShape;
1220
1221   if (aShape.ShapeType() == TopAbs_COMPOUND &&
1222       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
1223        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
1224        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND))
1225   {
1226     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
1227     for (; It.More(); It.Next()) {
1228       if (mapShape.Add(It.Value())) {
1229         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
1230             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
1231           listShape.Append(It.Value());
1232         }
1233       }
1234     }
1235   }
1236   else if (theExplodeType != EXPLODE_NEW_EXCLUDE_MAIN || aShape.ShapeType() != theShapeType) // issue 0021079
1237   {
1238     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
1239     for (; exp.More(); exp.Next())
1240       if (mapShape.Add(exp.Current()))
1241         listShape.Append(exp.Current());
1242   }
1243
1244   if (listShape.IsEmpty()) {
1245     //SetErrorCode("The given shape has no sub-shapes of the requested type");
1246     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
1247     return aSeq;
1248   }
1249
1250   if (isSorted) {
1251     bool isOldSorting = false;
1252     if (theExplodeType == EXPLODE_OLD_INCLUDE_MAIN)
1253       isOldSorting = true;
1254     SortShapes(listShape, isOldSorting);
1255   }
1256
1257   TopTools_IndexedMapOfShape anIndices;
1258   TopExp::MapShapes(aShape, anIndices);
1259   Handle(TColStd_HArray1OfInteger) anArray;
1260
1261   TopTools_ListIteratorOfListOfShape itSub (listShape);
1262   TCollection_AsciiString anAsciiList, anEntry;
1263   for (int index = 1; itSub.More(); itSub.Next(), ++index)
1264   {
1265     TopoDS_Shape aValue = itSub.Value();
1266     anArray = new TColStd_HArray1OfInteger(1,1);
1267     anArray->SetValue(1, anIndices.FindIndex(aValue));
1268
1269     //anObj = GetEngine()->AddSubShape(theShape, anArray);
1270     {
1271       anObj = GetEngine()->AddObject(GetDocID(), GEOM_SUBSHAPE);
1272       Handle(GEOM_Function) aFunction = anObj->AddFunction(GEOM_Object::GetSubShapeID(), 1);
1273       if (aFunction.IsNull()) return aSeq;
1274
1275       GEOM_ISubShape aSSI (aFunction);
1276       aSSI.SetMainShape(aMainShape);
1277       aSSI.SetIndices(anArray);
1278
1279       // Set function value directly, as we know it.
1280       // Usage of Solver here would lead to significant loss of time,
1281       // because GEOM_SubShapeDriver will build TopTools_IndexedMapOfShape
1282       // on the main shape for each being calculated sub-shape separately.
1283       aFunction->SetValue(aValue);
1284
1285       // Put this subshape in the list of subshapes of theMainShape
1286       aMainShape->AddSubShapeReference(aFunction);
1287     }
1288
1289     if (!anObj.IsNull()) {
1290       aSeq->Append(anObj);
1291
1292       // for python command
1293       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1294       anAsciiList += anEntry;
1295       anAsciiList += ",";
1296     }
1297   }
1298
1299   //Make a Python command
1300   anAsciiList.Trunc(anAsciiList.Length() - 1);
1301
1302   GEOM::TPythonDump pd (aMainShape, /*append=*/true);
1303   pd << "[" << anAsciiList.ToCString() << "] = geompy.";
1304   switch (theExplodeType) {
1305   case EXPLODE_NEW_EXCLUDE_MAIN:
1306     pd << "ExtractShapes(" << theShape << ", "
1307        << TopAbs_ShapeEnum(theShapeType) << ", " << (isSorted ? "True" : "False") << ")";
1308     break;
1309   case EXPLODE_NEW_INCLUDE_MAIN:
1310     pd << "SubShapeAll" << (isSorted ? "SortedCentres(" : "(")
1311        << theShape << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
1312     break;
1313   case EXPLODE_OLD_INCLUDE_MAIN:
1314     pd << "SubShapeAll" << (isSorted ? "Sorted(" : "(")
1315        << theShape << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
1316     break;
1317   default: ;
1318   }
1319   SetErrorCode(OK);
1320
1321   return aSeq;
1322 }
1323
1324 //=============================================================================
1325 /*!
1326  *  SubShapeAllIDs
1327  */
1328 //=============================================================================
1329 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::SubShapeAllIDs
1330                                           (Handle(GEOM_Object)    theShape,
1331                                            const Standard_Integer theShapeType,
1332                                            const Standard_Boolean isSorted,
1333                                            const ExplodeType      theExplodeType)
1334 {
1335   SetErrorCode(KO);
1336
1337   if (theShape.IsNull()) return NULL;
1338   TopoDS_Shape aShape = theShape->GetValue();
1339   if (aShape.IsNull()) return NULL;
1340
1341   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
1342   TopTools_MapOfShape mapShape;
1343   TopTools_ListOfShape listShape;
1344
1345   if (aShape.ShapeType() == TopAbs_COMPOUND &&
1346       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
1347        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
1348        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND))
1349   {
1350     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
1351     for (; It.More(); It.Next()) {
1352       if (mapShape.Add(It.Value())) {
1353         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
1354             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
1355           listShape.Append(It.Value());
1356         }
1357       }
1358     }
1359   }
1360   else if (theExplodeType != EXPLODE_NEW_EXCLUDE_MAIN || aShape.ShapeType() != theShapeType) // issue 0021079
1361   {
1362     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
1363     for (; exp.More(); exp.Next())
1364       if (mapShape.Add(exp.Current()))
1365         listShape.Append(exp.Current());
1366   }
1367
1368   if (listShape.IsEmpty()) {
1369     //SetErrorCode("The given shape has no sub-shapes of the requested type");
1370     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
1371     return aSeq;
1372   }
1373
1374   if (isSorted) {
1375     bool isOldSorting = false;
1376     if (theExplodeType == EXPLODE_OLD_INCLUDE_MAIN)
1377       isOldSorting = true;
1378     SortShapes(listShape, isOldSorting);
1379   }
1380
1381   TopTools_IndexedMapOfShape anIndices;
1382   TopExp::MapShapes(aShape, anIndices);
1383   Handle(TColStd_HArray1OfInteger) anArray;
1384
1385   TopTools_ListIteratorOfListOfShape itSub (listShape);
1386   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
1387     TopoDS_Shape aValue = itSub.Value();
1388     aSeq->Append(anIndices.FindIndex(aValue));
1389   }
1390
1391   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1392
1393   //Make a Python command
1394   GEOM::TPythonDump pd (aFunction, /*append=*/true);
1395   pd << "listSubShapeIDs = geompy.SubShapeAll";
1396   switch (theExplodeType) {
1397   case EXPLODE_NEW_EXCLUDE_MAIN:
1398     break;
1399   case EXPLODE_NEW_INCLUDE_MAIN:
1400     pd << (isSorted ? "SortedCentresIDs(" : "IDs(")
1401        << theShape << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
1402     break;
1403   case EXPLODE_OLD_INCLUDE_MAIN:
1404     pd << (isSorted ? "SortedIDs(" : "IDs(")
1405        << theShape << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
1406     break;
1407   default: ;
1408   }
1409
1410   SetErrorCode(OK);
1411   return aSeq;
1412 }
1413
1414 //=============================================================================
1415 /*!
1416  *  GetSubShape
1417  */
1418 //=============================================================================
1419 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetSubShape
1420                                           (Handle(GEOM_Object)    theMainShape,
1421                                            const Standard_Integer theID)
1422 {
1423   SetErrorCode(KO);
1424
1425   if (theMainShape.IsNull()) return NULL;
1426
1427   Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1428   anArray->SetValue(1, theID);
1429   Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theMainShape, anArray,true);
1430   if (anObj.IsNull()) {
1431     SetErrorCode("Can not get a sub-shape with the given ID");
1432     return NULL;
1433   }
1434
1435   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1436
1437   //Make a Python command
1438   GEOM::TPythonDump(aFunction) << anObj << " = geompy.GetSubShape("
1439                                << theMainShape << ", [" << theID << "])";
1440
1441   SetErrorCode(OK);
1442   return anObj;
1443 }
1444
1445 //=============================================================================
1446 /*!
1447  *  MakeSubShapes
1448  */
1449 //=============================================================================
1450 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::MakeSubShapes
1451                                 (Handle(GEOM_Object)              theMainShape,
1452                                  Handle(TColStd_HArray1OfInteger) theIndices)
1453 {
1454   SetErrorCode(KO);
1455
1456   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1457
1458   if (!theIndices->Length()) {
1459     SetErrorCode(NOT_FOUND_ANY);
1460     return aSeq;
1461   }
1462
1463   if (theMainShape.IsNull()) return NULL;
1464   TopoDS_Shape aShape = theMainShape->GetValue();
1465   if (aShape.IsNull()) return NULL;
1466
1467   Handle(GEOM_Function) aMainShape = theMainShape->GetLastFunction();
1468
1469   TopTools_IndexedMapOfShape anIndices;
1470   TopExp::MapShapes(aShape, anIndices);
1471
1472   Handle(TColStd_HArray1OfInteger) anArray;
1473   Handle(GEOM_Object) anObj;
1474
1475   TCollection_AsciiString anAsciiList, anEntry;
1476   Standard_Integer i, low = theIndices->Lower(), up = theIndices->Upper();
1477   for (i = low; i <= up; i++) {
1478     int id = theIndices->Value(i);
1479     if (1 <= id && id <= anIndices.Extent()) {
1480       TopoDS_Shape aValue = anIndices.FindKey(id);
1481       anArray = new TColStd_HArray1OfInteger(1,1);
1482       anArray->SetValue(1, id);
1483
1484       anObj = GetEngine()->AddObject(GetDocID(), GEOM_SUBSHAPE);
1485       if (!anObj.IsNull()) {
1486         Handle(GEOM_Function) aFunction = anObj->AddFunction(GEOM_Object::GetSubShapeID(), 1);
1487         if (aFunction.IsNull()) return aSeq;
1488
1489         GEOM_ISubShape aSSI (aFunction);
1490         aSSI.SetMainShape(aMainShape);
1491         aSSI.SetIndices(anArray);
1492
1493         // Set function value directly, as we know it.
1494         // Usage of Solver here would lead to significant loss of time,
1495         // because GEOM_SubShapeDriver will build TopTools_IndexedMapOfShape
1496         // on the main shape for each being calculated sub-shape separately.
1497         aFunction->SetValue(aValue);
1498
1499         // Put this subshape in the list of subshapes of theMainShape
1500         aMainShape->AddSubShapeReference(aFunction);
1501
1502         aSeq->Append(anObj);
1503
1504         // for python command
1505         TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1506         anAsciiList += anEntry;
1507         anAsciiList += ",";
1508       }
1509     }
1510   }
1511
1512   //Make a Python command
1513   anAsciiList.Trunc(anAsciiList.Length() - 1);
1514
1515   GEOM::TPythonDump pd (aMainShape, /*append=*/true);
1516   pd << "[" << anAsciiList.ToCString() << "] = geompy.SubShapes("
1517      << theMainShape << ", [" ;
1518   for (i = low; i <= up - 1; i++) {
1519     pd << theIndices->Value(i) << ", ";
1520   }
1521   pd << theIndices->Value(up) << "])";
1522
1523   SetErrorCode(OK);
1524
1525   return aSeq;
1526 }
1527
1528 //=============================================================================
1529 /*!
1530  *  GetSubShapeIndex
1531  */
1532 //=============================================================================
1533 Standard_Integer GEOMImpl_IShapesOperations::GetSubShapeIndex (Handle(GEOM_Object) theMainShape,
1534                                                                Handle(GEOM_Object) theSubShape)
1535 {
1536   SetErrorCode(KO);
1537
1538   TopoDS_Shape aMainShape = theMainShape->GetValue();
1539   TopoDS_Shape aSubShape = theSubShape->GetValue();
1540
1541   if (aMainShape.IsNull() || aSubShape.IsNull()) return -1;
1542
1543   TopTools_IndexedMapOfShape anIndices;
1544   TopExp::MapShapes(aMainShape, anIndices);
1545   if (anIndices.Contains(aSubShape)) {
1546     SetErrorCode(OK);
1547     return anIndices.FindIndex(aSubShape);
1548   }
1549
1550   return -1;
1551 }
1552
1553 //=============================================================================
1554 /*!
1555  *  GetTopologyIndex
1556  */
1557 //=============================================================================
1558 Standard_Integer GEOMImpl_IShapesOperations::GetTopologyIndex (Handle(GEOM_Object) theMainShape,
1559                                                                Handle(GEOM_Object) theSubShape)
1560 {
1561   SetErrorCode(OK);
1562
1563   TopoDS_Shape aMainShape = theMainShape->GetValue();
1564   TopoDS_Shape aSubShape = theSubShape->GetValue();
1565
1566   if (aMainShape.IsNull() || aSubShape.IsNull()) {
1567     SetErrorCode("Null argument shape given");
1568     return -1;
1569   }
1570
1571   int index = 1;
1572   if (aSubShape.ShapeType() == TopAbs_COMPOUND) {
1573     TopoDS_Iterator it;
1574     TopTools_ListOfShape CL;
1575     CL.Append(aMainShape);
1576     TopTools_ListIteratorOfListOfShape itC;
1577     for (itC.Initialize(CL); itC.More(); itC.Next()) {
1578       for (it.Initialize(itC.Value()); it.More(); it.Next()) {
1579         if (it.Value().ShapeType() == TopAbs_COMPOUND) {
1580           if (it.Value().IsSame(aSubShape))
1581             return index;
1582           else
1583             index++;
1584           CL.Append(it.Value());
1585         }
1586       }
1587     }
1588   } else {
1589     TopExp_Explorer anExp (aMainShape, aSubShape.ShapeType());
1590     TopTools_MapOfShape M;
1591     for (; anExp.More(); anExp.Next()) {
1592       if (M.Add(anExp.Current())) {
1593         if (anExp.Current().IsSame(aSubShape))
1594           return index;
1595         index++;
1596       }
1597     }
1598   }
1599
1600   SetErrorCode("The sub-shape does not belong to the main shape");
1601   return -1;
1602 }
1603
1604 //=============================================================================
1605 /*!
1606  *  GetShapeTypeString
1607  */
1608 //=============================================================================
1609 TCollection_AsciiString GEOMImpl_IShapesOperations::GetShapeTypeString (Handle(GEOM_Object) theShape)
1610 {
1611   SetErrorCode(KO);
1612
1613   TCollection_AsciiString aTypeName ("Null Shape");
1614
1615   TopoDS_Shape aShape = theShape->GetValue();
1616   if (aShape.IsNull())
1617     return aTypeName;
1618
1619   switch (aShape.ShapeType() )
1620   {
1621   case TopAbs_COMPOUND:
1622     aTypeName = "Compound";
1623     break;
1624   case  TopAbs_COMPSOLID:
1625     aTypeName = "Compound Solid";
1626     break;
1627   case TopAbs_SOLID:
1628     aTypeName = "Solid";
1629     break;
1630   case TopAbs_SHELL:
1631     aTypeName = "Shell";
1632     break;
1633   case TopAbs_FACE:
1634     {
1635       BRepAdaptor_Surface surf (TopoDS::Face(aShape));
1636       if (surf.GetType() == GeomAbs_Plane)
1637         aTypeName = "Plane";
1638       else if (surf.GetType() == GeomAbs_Cylinder)
1639         aTypeName = "Cylindrical Face";
1640       else if (surf.GetType() == GeomAbs_Sphere)
1641         aTypeName = "Spherical Face";
1642       else if (surf.GetType() == GeomAbs_Torus)
1643         aTypeName = "Toroidal Face";
1644       else if (surf.GetType() == GeomAbs_Cone)
1645         aTypeName = "Conical Face";
1646       else
1647         aTypeName = "GEOM::FACE";
1648     }
1649     break;
1650   case TopAbs_WIRE:
1651     aTypeName = "Wire";
1652     break;
1653   case TopAbs_EDGE:
1654     {
1655       BRepAdaptor_Curve curv (TopoDS::Edge(aShape));
1656       if (curv.GetType() == GeomAbs_Line) {
1657         if ((Abs(curv.FirstParameter()) >= 1E6) ||
1658             (Abs(curv.LastParameter()) >= 1E6))
1659           aTypeName = "Line";
1660         else
1661           aTypeName = "Edge";
1662       } else if (curv.GetType() == GeomAbs_Circle) {
1663         if (curv.IsClosed())
1664           aTypeName = "Circle";
1665         else
1666           aTypeName = "Arc";
1667       } else {
1668         aTypeName = "Edge";
1669       }
1670     }
1671     break;
1672   case TopAbs_VERTEX:
1673     aTypeName = "Vertex";
1674     break;
1675   case TopAbs_SHAPE:
1676     aTypeName = "Shape";
1677     break;
1678   default:
1679     aTypeName = "Shape of unknown type";
1680   }
1681
1682   return aTypeName;
1683 }
1684
1685 //=============================================================================
1686 /*!
1687  *  NumberOfSubShapes
1688  */
1689 //=============================================================================
1690 Standard_Integer GEOMImpl_IShapesOperations::NumberOfSubShapes
1691                                           (Handle(GEOM_Object)    theShape,
1692                                            const Standard_Integer theShapeType)
1693 {
1694   SetErrorCode(KO);
1695   Standard_Integer nbShapes = 0;
1696
1697   if (theShape.IsNull()) return -1;
1698   TopoDS_Shape aShape = theShape->GetValue();
1699   if (aShape.IsNull()) return -1;
1700
1701   /*
1702   TopTools_MapOfShape mapShape;
1703
1704   if (aShape.ShapeType() == TopAbs_COMPOUND &&
1705       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
1706        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
1707        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND)) {
1708     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
1709     for (; It.More(); It.Next()) {
1710       if (mapShape.Add(It.Value())) {
1711         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
1712             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
1713           nbShapes++;
1714         }
1715       }
1716     }
1717   } else {
1718     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
1719     for (; exp.More(); exp.Next())
1720       if (mapShape.Add(exp.Current()))
1721         nbShapes++;
1722   }
1723   */
1724
1725   try {
1726 #if OCC_VERSION_LARGE > 0x06010000
1727     OCC_CATCH_SIGNALS;
1728 #endif
1729     int iType, nbTypes [TopAbs_SHAPE];
1730     for (iType = 0; iType < TopAbs_SHAPE; ++iType)
1731       nbTypes[iType] = 0;
1732     nbTypes[aShape.ShapeType()]++;
1733
1734     TopTools_MapOfShape aMapOfShape;
1735     aMapOfShape.Add(aShape);
1736     TopTools_ListOfShape aListOfShape;
1737     aListOfShape.Append(aShape);
1738
1739     TopTools_ListIteratorOfListOfShape itL (aListOfShape);
1740     for (; itL.More(); itL.Next()) {
1741       TopoDS_Iterator it (itL.Value());
1742       for (; it.More(); it.Next()) {
1743         TopoDS_Shape s = it.Value();
1744         if (aMapOfShape.Add(s)) {
1745           aListOfShape.Append(s);
1746           nbTypes[s.ShapeType()]++;
1747         }
1748       }
1749     }
1750
1751     if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE)
1752       nbShapes = aMapOfShape.Extent();
1753     else
1754       nbShapes = nbTypes[theShapeType];
1755   }
1756   catch (Standard_Failure) {
1757     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1758     SetErrorCode(aFail->GetMessageString());
1759     return -1;
1760   }
1761
1762   SetErrorCode(OK);
1763   return nbShapes;
1764 }
1765
1766 //=============================================================================
1767 /*!
1768  *  ReverseShape
1769  */
1770 //=============================================================================
1771 Handle(GEOM_Object) GEOMImpl_IShapesOperations::ReverseShape(Handle(GEOM_Object) theShape)
1772 {
1773   SetErrorCode(KO);
1774
1775   if (theShape.IsNull()) return NULL;
1776
1777   //Add a new reversed object
1778   Handle(GEOM_Object) aReversed = GetEngine()->AddObject(GetDocID(), theShape->GetType());
1779
1780   //Add a new Revese function
1781   Handle(GEOM_Function) aFunction;
1782   aFunction = aReversed->AddFunction(GEOMImpl_ShapeDriver::GetID(), REVERSE_ORIENTATION);
1783   if (aFunction.IsNull()) return NULL;
1784
1785   //Check if the function is set correctly
1786   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
1787
1788   GEOMImpl_IShapes aSI (aFunction);
1789
1790   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1791   if (aRefShape.IsNull()) return NULL;
1792
1793   aSI.SetBase(aRefShape);
1794
1795   //Compute the sub-shape value
1796   try {
1797 #if OCC_VERSION_LARGE > 0x06010000
1798     OCC_CATCH_SIGNALS;
1799 #endif
1800     if (!GetSolver()->ComputeFunction(aFunction)) {
1801       SetErrorCode("Shape driver failed to reverse shape");
1802       return NULL;
1803     }
1804   }
1805   catch (Standard_Failure) {
1806     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1807     SetErrorCode(aFail->GetMessageString());
1808     return NULL;
1809   }
1810
1811   //Make a Python command
1812   GEOM::TPythonDump(aFunction) << aReversed
1813     << " = geompy.ChangeOrientation(" << theShape << ")";
1814
1815   SetErrorCode(OK);
1816   return aReversed;
1817 }
1818
1819 //=============================================================================
1820 /*!
1821  *  GetFreeFacesIDs
1822  */
1823 //=============================================================================
1824 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetFreeFacesIDs
1825                                                  (Handle(GEOM_Object) theShape)
1826 {
1827   SetErrorCode(KO);
1828
1829   if (theShape.IsNull()) return NULL;
1830   TopoDS_Shape aShape = theShape->GetValue();
1831   if (aShape.IsNull()) return NULL;
1832
1833   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
1834
1835   TopTools_IndexedDataMapOfShapeListOfShape mapFaceBlocks;
1836   GEOMImpl_Block6Explorer::MapShapesAndAncestors
1837     (aShape, TopAbs_FACE, TopAbs_SOLID, mapFaceBlocks);
1838
1839   Standard_Integer ind = 1, nbFaces = mapFaceBlocks.Extent();
1840
1841   if (nbFaces == 0) {
1842     SetErrorCode("The given shape has no faces");
1843     return aSeq;
1844   }
1845
1846   TopTools_IndexedMapOfShape anIndices;
1847   TopExp::MapShapes(aShape, anIndices);
1848
1849   Standard_Integer id;
1850   for (; ind <= nbFaces; ind++) {
1851     if (mapFaceBlocks.FindFromIndex(ind).Extent() != 2) {
1852       id = anIndices.FindIndex(mapFaceBlocks.FindKey(ind));
1853       aSeq->Append(id);
1854     }
1855   }
1856
1857   //The explode doesn't change object so no new function is required.
1858   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1859
1860   //Make a Python command
1861   GEOM::TPythonDump(aFunction, /*append=*/true)
1862     << "listFreeFacesIDs = geompy.GetFreeFacesIDs(" << theShape << ")";
1863
1864   SetErrorCode(OK);
1865   return aSeq;
1866 }
1867
1868 //=======================================================================
1869 //function : GetSharedShapes
1870 //purpose  :
1871 //=======================================================================
1872 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetSharedShapes
1873                                                 (Handle(GEOM_Object)    theShape1,
1874                                                  Handle(GEOM_Object)    theShape2,
1875                                                  const Standard_Integer theShapeType)
1876 {
1877   SetErrorCode(KO);
1878
1879   if (theShape1.IsNull() || theShape2.IsNull()) return NULL;
1880
1881   TopoDS_Shape aShape1 = theShape1->GetValue();
1882   TopoDS_Shape aShape2 = theShape2->GetValue();
1883
1884   if (aShape1.IsNull() || aShape2.IsNull()) return NULL;
1885
1886   TopTools_IndexedMapOfShape anIndices;
1887   TopExp::MapShapes(aShape1, anIndices);
1888   Handle(TColStd_HArray1OfInteger) anArray;
1889
1890   TopTools_IndexedMapOfShape mapShape1;
1891   TopExp::MapShapes(aShape1, TopAbs_ShapeEnum(theShapeType), mapShape1);
1892
1893   Handle(GEOM_Object) anObj;
1894   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1895   TCollection_AsciiString anAsciiList, anEntry;
1896
1897   TopTools_MapOfShape mapShape2;
1898   TopExp_Explorer exp (aShape2, TopAbs_ShapeEnum(theShapeType));
1899   for (; exp.More(); exp.Next()) {
1900     TopoDS_Shape aSS = exp.Current();
1901     if (mapShape2.Add(aSS) && mapShape1.Contains(aSS)) {
1902       anArray = new TColStd_HArray1OfInteger(1,1);
1903       anArray->SetValue(1, anIndices.FindIndex(aSS));
1904       anObj = GetEngine()->AddSubShape(theShape1, anArray);
1905       aSeq->Append(anObj);
1906
1907       // for python command
1908       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1909       anAsciiList += anEntry;
1910       anAsciiList += ",";
1911     }
1912   }
1913
1914   if (aSeq->IsEmpty()) {
1915     SetErrorCode("The given shapes have no shared sub-shapes of the requested type");
1916     return aSeq;
1917   }
1918
1919   //Make a Python command
1920   anAsciiList.Trunc(anAsciiList.Length() - 1);
1921
1922   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1923
1924   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
1925     << "] = geompy.GetSharedShapes(" << theShape1 << ", "
1926       << theShape2 << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
1927
1928   SetErrorCode(OK);
1929   return aSeq;
1930 }
1931
1932 //=======================================================================
1933 //function : GetSharedShapes
1934 //purpose  :
1935 //=======================================================================
1936 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetSharedShapes
1937                                      (std::list<Handle(GEOM_Object)> theShapes,
1938                                       const Standard_Integer         theShapeType)
1939 {
1940   SetErrorCode(KO);
1941
1942   int aLen = theShapes.size();
1943   if (aLen < 1) return NULL;
1944
1945   int ind = 1;
1946   std::list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
1947
1948   Handle(GEOM_Object) aMainObj = (*it++);
1949   Handle(GEOM_Function) aMainShape = aMainObj->GetLastFunction();
1950   if (aMainShape.IsNull()) {
1951     SetErrorCode("NULL shape for GetSharedShapes");
1952     return NULL;
1953   }
1954
1955   TopoDS_Shape aShape1 = aMainShape->GetValue();
1956   if (aShape1.IsNull()) return NULL;
1957
1958   TopTools_IndexedMapOfShape anIndices;
1959   TopExp::MapShapes(aShape1, anIndices);
1960
1961   TopTools_IndexedMapOfShape mapSelected;
1962   TopExp::MapShapes(aShape1, TopAbs_ShapeEnum(theShapeType), mapSelected);
1963
1964   // Find shared shapes
1965   BRep_Builder B;
1966   TopoDS_Compound aCurrSelection;
1967
1968   for (; it != theShapes.end(); it++, ind++) {
1969     Handle(GEOM_Function) aRefShape = (*it)->GetLastFunction();
1970     if (aRefShape.IsNull()) {
1971       SetErrorCode("NULL shape for GetSharedShapes");
1972       return NULL;
1973     }
1974
1975     TopoDS_Compound aCompound;
1976     B.MakeCompound(aCompound);
1977
1978     TopoDS_Shape aShape2 = aRefShape->GetValue();
1979     if (aShape2.IsNull()) return NULL;
1980
1981     TopTools_MapOfShape mapShape2;
1982     TopExp_Explorer exp (aShape2, TopAbs_ShapeEnum(theShapeType));
1983     for (; exp.More(); exp.Next()) {
1984       TopoDS_Shape aSS = exp.Current();
1985       if (mapShape2.Add(aSS) && mapSelected.Contains(aSS)) {
1986         B.Add(aCompound, aSS);
1987       }
1988     }
1989
1990     mapSelected.Clear();
1991     TopExp::MapShapes(aCompound, TopAbs_ShapeEnum(theShapeType), mapSelected);
1992     aCurrSelection = aCompound;
1993   }
1994
1995   // Create GEOM_Object for each found shared shape (collected in aCurrSelection)
1996   Handle(GEOM_Object) anObj;
1997   Handle(TColStd_HArray1OfInteger) anArray;
1998   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1999   TCollection_AsciiString anAsciiList, anEntry;
2000
2001   TopoDS_Iterator itSel (aCurrSelection, Standard_True, Standard_True);
2002   for (; itSel.More(); itSel.Next()) {
2003     anArray = new TColStd_HArray1OfInteger(1,1);
2004     anArray->SetValue(1, anIndices.FindIndex(itSel.Value()));
2005     anObj = GetEngine()->AddSubShape(aMainObj, anArray);
2006     aSeq->Append(anObj);
2007
2008     // for python command
2009     TDF_Tool::Entry(anObj->GetEntry(), anEntry);
2010     anAsciiList += anEntry;
2011     anAsciiList += ",";
2012   }
2013
2014   if (aSeq->IsEmpty()) {
2015     SetErrorCode("The given shapes have no shared sub-shapes of the requested type");
2016     return aSeq;
2017   }
2018
2019   // Make a Python command
2020   anAsciiList.Trunc(anAsciiList.Length() - 1);
2021
2022   GEOM::TPythonDump pd (aMainShape, /*append=*/true);
2023   pd << "[" << anAsciiList.ToCString()
2024      << "] = geompy.GetSharedShapesMulti([";
2025
2026   it = theShapes.begin();
2027   pd << (*it++);
2028   while (it != theShapes.end()) {
2029     pd << ", " << (*it++);
2030   }
2031
2032   pd << "], " << TopAbs_ShapeEnum(theShapeType) << ")";
2033
2034   SetErrorCode(OK);
2035   return aSeq;
2036 }
2037
2038 //=============================================================================
2039 /*!
2040  *
2041  */
2042 //=============================================================================
2043 static GEOM::TPythonDump& operator<< (GEOM::TPythonDump&   theDump,
2044                                       const GEOMAlgo_State theState)
2045 {
2046   switch (theState) {
2047   case GEOMAlgo_ST_IN:
2048     theDump << "geompy.GEOM.ST_IN";
2049     break;
2050   case GEOMAlgo_ST_OUT:
2051     theDump << "geompy.GEOM.ST_OUT";
2052     break;
2053   case GEOMAlgo_ST_ON:
2054     theDump << "geompy.GEOM.ST_ON";
2055     break;
2056   case GEOMAlgo_ST_ONIN:
2057     theDump << "geompy.GEOM.ST_ONIN";
2058     break;
2059   case GEOMAlgo_ST_ONOUT:
2060     theDump << "geompy.GEOM.ST_ONOUT";
2061     break;
2062   default:
2063     theDump << "geompy.GEOM.ST_UNKNOWN";
2064     break;
2065   }
2066   return theDump;
2067 }
2068
2069 //=======================================================================
2070 //function : checkTypeShapesOn
2071 /*!
2072  * \brief Checks if theShapeType parameter of GetShapesOnXXX() is OK
2073  * \param theShapeType - the shape type to check
2074  * \retval bool  - result of the check
2075  */
2076 //=======================================================================
2077 bool GEOMImpl_IShapesOperations::checkTypeShapesOn(const Standard_Integer theShapeType)
2078 {
2079   if (theShapeType != TopAbs_VERTEX &&
2080       theShapeType != TopAbs_EDGE &&
2081       theShapeType != TopAbs_FACE &&
2082       theShapeType != TopAbs_SOLID) {
2083     SetErrorCode("Only solids, vertices, edges or faces can be found by this method");
2084     return false;
2085   }
2086   return true;
2087 }
2088
2089 //=======================================================================
2090 //function : makePlane
2091   /*!
2092    * \brief Creates Geom_Plane
2093     * \param theAx1 - shape object defining plane parameters
2094     * \retval Handle(Geom_Surface) - resulting surface
2095    */
2096 //=======================================================================
2097 Handle(Geom_Surface) GEOMImpl_IShapesOperations::makePlane(const TopoDS_Shape& anAx1)
2098 {
2099   if (anAx1.ShapeType() != TopAbs_EDGE) return NULL;
2100   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
2101   TopoDS_Vertex V1, V2;
2102   TopExp::Vertices(anEdge, V1, V2, Standard_True);
2103   if (V1.IsNull() || V2.IsNull()) {
2104     SetErrorCode("Bad edge given for the plane normal vector");
2105     return NULL;
2106   }
2107   gp_Pnt aLoc = BRep_Tool::Pnt(V1);
2108   gp_Vec aVec (aLoc, BRep_Tool::Pnt(V2));
2109   if (aVec.Magnitude() < Precision::Confusion()) {
2110     SetErrorCode("Vector with null magnitude given");
2111     return NULL;
2112   }
2113   return new Geom_Plane(aLoc, aVec);
2114 }
2115
2116 //=======================================================================
2117 //function : makeCylinder
2118   /*!
2119    * \brief Creates Geom_CylindricalSurface
2120     * \param theAx1 - edge defining cylinder axis
2121     * \param theRadius - cylinder radius
2122     * \retval Handle(Geom_Surface) - resulting surface
2123    */
2124 //=======================================================================
2125 Handle(Geom_Surface) GEOMImpl_IShapesOperations::makeCylinder(const TopoDS_Shape& anAxis,
2126                                                               const Standard_Real theRadius)
2127 {
2128   //Axis of the cylinder
2129   if (anAxis.ShapeType() != TopAbs_EDGE) {
2130     SetErrorCode("Not an edge given for the axis");
2131     return NULL;
2132   }
2133   TopoDS_Edge anEdge = TopoDS::Edge(anAxis);
2134   TopoDS_Vertex V1, V2;
2135   TopExp::Vertices(anEdge, V1, V2, Standard_True);
2136   if (V1.IsNull() || V2.IsNull()) {
2137     SetErrorCode("Bad edge given for the axis");
2138     return NULL;
2139   }
2140   gp_Pnt aLoc = BRep_Tool::Pnt(V1);
2141   gp_Vec aVec (aLoc, BRep_Tool::Pnt(V2));
2142   if (aVec.Magnitude() < Precision::Confusion()) {
2143     SetErrorCode("Vector with null magnitude given");
2144     return NULL;
2145   }
2146
2147   gp_Ax3 anAx3 (aLoc, aVec);
2148   return new Geom_CylindricalSurface(anAx3, theRadius);
2149 }
2150
2151 //=======================================================================
2152 //function : getShapesOnBoxIDs
2153   /*!
2154    * \brief Find IDs of subshapes complying with given status about surface
2155     * \param theBox - the box to check state of subshapes against
2156     * \param theShape - the shape to explore
2157     * \param theShapeType - type of subshape of theShape
2158     * \param theState - required state
2159     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2160    */
2161 //=======================================================================
2162 Handle(TColStd_HSequenceOfInteger)
2163   GEOMImpl_IShapesOperations::getShapesOnBoxIDs(const Handle(GEOM_Object)& theBox,
2164                                                 const Handle(GEOM_Object)& theShape,
2165                                                 const Standard_Integer theShapeType,
2166                                                 GEOMAlgo_State theState)
2167 {
2168   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
2169
2170   TopoDS_Shape aBox = theBox->GetValue();
2171   TopoDS_Shape aShape = theShape->GetValue();
2172
2173   // Check presence of triangulation, build if need
2174   if (!CheckTriangulation(aShape)) {
2175     SetErrorCode("Cannot build triangulation on the shape");
2176     return aSeqOfIDs;
2177   }
2178
2179   // Call algo
2180   GEOMAlgo_FinderShapeOn2 aFinder;
2181   Standard_Real aTol = 0.0001; // default value
2182
2183   Handle(GEOMAlgo_ClsfBox) aClsfBox = new GEOMAlgo_ClsfBox;
2184   aClsfBox->SetBox(aBox);
2185
2186   aFinder.SetShape(aShape);
2187   aFinder.SetTolerance(aTol);
2188   aFinder.SetClsf(aClsfBox);
2189   aFinder.SetShapeType( (TopAbs_ShapeEnum)theShapeType );
2190   aFinder.SetState(theState);
2191   aFinder.Perform();
2192
2193   // Interprete results
2194   Standard_Integer iErr = aFinder.ErrorStatus();
2195   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
2196   if (iErr) {
2197     MESSAGE(" iErr : " << iErr);
2198     TCollection_AsciiString aMsg (" iErr : ");
2199     aMsg += TCollection_AsciiString(iErr);
2200     SetErrorCode(aMsg);
2201     return aSeqOfIDs;
2202   }
2203   Standard_Integer iWrn = aFinder.WarningStatus();
2204   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
2205   if (iWrn) {
2206     MESSAGE(" *** iWrn : " << iWrn);
2207   }
2208
2209   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
2210
2211   if (listSS.Extent() < 1) {
2212     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
2213     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
2214     return aSeqOfIDs;
2215   }
2216
2217   // Fill sequence of object IDs
2218   aSeqOfIDs = new TColStd_HSequenceOfInteger;
2219
2220   TopTools_IndexedMapOfShape anIndices;
2221   TopExp::MapShapes(aShape, anIndices);
2222
2223   TopTools_ListIteratorOfListOfShape itSub (listSS);
2224   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
2225     int id = anIndices.FindIndex(itSub.Value());
2226     aSeqOfIDs->Append(id);
2227   }
2228
2229   return aSeqOfIDs;
2230 }
2231
2232 //=======================================================================
2233 //function : GetShapesOnBoxIDs
2234 /*!
2235    * \brief Find subshapes complying with given status about surface
2236     * \param theBox - the box to check state of subshapes against
2237     * \param theShape - the shape to explore
2238     * \param theShapeType - type of subshape of theShape
2239     * \param theState - required state
2240     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2241  */
2242 //=======================================================================
2243 Handle(TColStd_HSequenceOfInteger)
2244     GEOMImpl_IShapesOperations::GetShapesOnBoxIDs(const Handle(GEOM_Object)& theBox,
2245                                                   const Handle(GEOM_Object)& theShape,
2246                                                   const Standard_Integer theShapeType,
2247                                                   GEOMAlgo_State theState)
2248 {
2249   // Find subshapes ids
2250   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2251     getShapesOnBoxIDs (theBox, theShape, theShapeType, theState);
2252   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
2253     return NULL;
2254
2255   // The GetShapesOnBox() doesn't change object so no new function is required.
2256   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theBox)->GetLastFunction();
2257
2258   // Make a Python command
2259   GEOM::TPythonDump(aFunction)
2260     << "listShapesOnBoxIDs = geompy.GetShapesOnBoxIDs("
2261     << theBox << ", "
2262     << theShape << ", "
2263     << TopAbs_ShapeEnum(theShapeType) << ", "
2264     << theState << ")";
2265
2266   SetErrorCode(OK);
2267   return aSeqOfIDs;
2268 }
2269
2270 //=======================================================================
2271 //function : GetShapesOnBox
2272 /*!
2273    * \brief Find subshapes complying with given status about surface
2274     * \param theBox - the box to check state of subshapes against
2275     * \param theShape - the shape to explore
2276     * \param theShapeType - type of subshape of theShape
2277     * \param theState - required state
2278     * \retval Handle(TColStd_HSequenceOfTransient) - found subshapes
2279  */
2280 //=======================================================================
2281 Handle(TColStd_HSequenceOfTransient)
2282     GEOMImpl_IShapesOperations::GetShapesOnBox(const Handle(GEOM_Object)& theBox,
2283                                                const Handle(GEOM_Object)&  theShape,
2284                                                const Standard_Integer theShapeType,
2285                                                GEOMAlgo_State theState)
2286 {
2287   // Find subshapes ids
2288   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2289     getShapesOnBoxIDs (theBox, theShape, theShapeType, theState);
2290   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
2291     return NULL;
2292
2293   // Find objects by indices
2294   TCollection_AsciiString anAsciiList;
2295   Handle(TColStd_HSequenceOfTransient) aSeq;
2296   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
2297   if ( aSeq.IsNull() || aSeq->IsEmpty() )
2298     return NULL;
2299
2300   // Make a Python command
2301
2302   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2303   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2304
2305   GEOM::TPythonDump(aFunction)
2306     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnBox("
2307     << theBox << ", "
2308     << theShape << ", "
2309     << TopAbs_ShapeEnum(theShapeType) << ", "
2310     << theState << ")";
2311
2312   SetErrorCode(OK);
2313   return aSeq;
2314 }
2315
2316 //=======================================================================
2317 //function : getShapesOnShapeIDs
2318 /*!
2319  * \brief Find IDs of subshapes complying with given status about surface
2320  * \param theCheckShape - the shape to check state of subshapes against
2321  * \param theShape - the shape to explore
2322  * \param theShapeType - type of subshape of theShape
2323  * \param theState - required state
2324  * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2325  */
2326 //=======================================================================
2327 Handle(TColStd_HSequenceOfInteger)
2328   GEOMImpl_IShapesOperations::getShapesOnShapeIDs
2329                                  (const Handle(GEOM_Object)& theCheckShape,
2330                                   const Handle(GEOM_Object)& theShape,
2331                                   const Standard_Integer theShapeType,
2332                                   GEOMAlgo_State theState)
2333 {
2334   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
2335
2336   TopoDS_Shape aCheckShape = theCheckShape->GetValue();
2337   TopoDS_Shape aShape = theShape->GetValue();
2338   TopTools_ListOfShape res;
2339
2340   // Check presence of triangulation, build if need
2341   if (!CheckTriangulation(aShape)) {
2342     SetErrorCode("Cannot build triangulation on the shape");
2343     return aSeqOfIDs;
2344   }
2345
2346   // Call algo
2347   GEOMAlgo_FinderShapeOn2 aFinder;
2348   Standard_Real aTol = 0.0001; // default value
2349
2350   Handle(GEOMAlgo_ClsfSolid) aClsfSolid = new GEOMAlgo_ClsfSolid;
2351   aClsfSolid->SetShape(aCheckShape);
2352
2353   aFinder.SetShape(aShape);
2354   aFinder.SetTolerance(aTol);
2355   aFinder.SetClsf(aClsfSolid);
2356   aFinder.SetShapeType( (TopAbs_ShapeEnum)theShapeType );
2357   aFinder.SetState(theState);
2358   aFinder.Perform();
2359
2360   // Interprete results
2361   Standard_Integer iErr = aFinder.ErrorStatus();
2362   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
2363   if (iErr) {
2364     if (iErr == 41) {
2365       SetErrorCode("theCheckShape must be a solid");
2366     }
2367     else {
2368       MESSAGE(" iErr : " << iErr);
2369       TCollection_AsciiString aMsg (" iErr : ");
2370       aMsg += TCollection_AsciiString(iErr);
2371       SetErrorCode(aMsg);
2372     }
2373     return aSeqOfIDs;
2374   }
2375   Standard_Integer iWrn = aFinder.WarningStatus();
2376   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
2377   if (iWrn) {
2378     MESSAGE(" *** iWrn : " << iWrn);
2379   }
2380
2381   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
2382
2383   if (listSS.Extent() < 1) {
2384     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
2385     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
2386   }
2387
2388   // Fill sequence of object IDs
2389   aSeqOfIDs = new TColStd_HSequenceOfInteger;
2390
2391   TopTools_IndexedMapOfShape anIndices;
2392   TopExp::MapShapes(aShape, anIndices);
2393
2394   TopTools_ListIteratorOfListOfShape itSub (listSS);
2395   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
2396     int id = anIndices.FindIndex(itSub.Value());
2397     aSeqOfIDs->Append(id);
2398   }
2399
2400   return aSeqOfIDs;
2401 }
2402
2403 //=======================================================================
2404 //function : GetShapesOnShapeIDs
2405 /*!
2406  * \brief Find subshapes complying with given status about surface
2407  * \param theCheckShape - the shape to check state of subshapes against
2408  * \param theShape - the shape to explore
2409  * \param theShapeType - type of subshape of theShape
2410  * \param theState - required state
2411  * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2412  */
2413 //=======================================================================
2414 Handle(TColStd_HSequenceOfInteger)
2415     GEOMImpl_IShapesOperations::GetShapesOnShapeIDs
2416                             (const Handle(GEOM_Object)& theCheckShape,
2417                              const Handle(GEOM_Object)& theShape,
2418                              const Standard_Integer theShapeType,
2419                              GEOMAlgo_State theState)
2420 {
2421   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2422     getShapesOnShapeIDs (theCheckShape, theShape, theShapeType, theState);
2423
2424   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
2425     return NULL;
2426
2427   // The GetShapesOnShape() doesn't change object so no new function is required.
2428   Handle(GEOM_Function) aFunction =
2429     GEOM::GetCreatedLast(theShape,theCheckShape)->GetLastFunction();
2430
2431   // Make a Python command
2432   GEOM::TPythonDump(aFunction)
2433     << "listShapesOnBoxIDs = geompy.GetShapesOnShapeIDs("
2434     << theCheckShape << ", "
2435     << theShape << ", "
2436     << TopAbs_ShapeEnum(theShapeType) << ", "
2437     << theState << ")";
2438
2439   SetErrorCode(OK);
2440   return aSeqOfIDs;
2441 }
2442
2443 //=======================================================================
2444 //function : GetShapesOnShape
2445 /*!
2446  * \brief Find subshapes complying with given status about surface
2447  * \param theCheckShape - the shape to check state of subshapes against
2448  * \param theShape - the shape to explore
2449  * \param theShapeType - type of subshape of theShape
2450  * \param theState - required state
2451  * \retval Handle(TColStd_HSequenceOfTransient) - found subshapes
2452  */
2453 //=======================================================================
2454 Handle(TColStd_HSequenceOfTransient)
2455   GEOMImpl_IShapesOperations::GetShapesOnShape
2456                              (const Handle(GEOM_Object)& theCheckShape,
2457                               const Handle(GEOM_Object)&  theShape,
2458                               const Standard_Integer theShapeType,
2459                               GEOMAlgo_State theState)
2460 {
2461   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2462     getShapesOnShapeIDs (theCheckShape, theShape, theShapeType, theState);
2463   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
2464     return NULL;
2465
2466   // Find objects by indices
2467   TCollection_AsciiString anAsciiList;
2468   Handle(TColStd_HSequenceOfTransient) aSeq;
2469   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
2470
2471   if ( aSeq.IsNull() || aSeq->IsEmpty() )
2472     return NULL;
2473
2474   // Make a Python command
2475
2476   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2477   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2478
2479   GEOM::TPythonDump(aFunction)
2480     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnShape("
2481     << theCheckShape << ", "
2482     << theShape << ", "
2483     << TopAbs_ShapeEnum(theShapeType) << ", "
2484     << theState << ")";
2485
2486   SetErrorCode(OK);
2487   return aSeq;
2488 }
2489
2490 //=======================================================================
2491 //function : GetShapesOnShapeAsCompound
2492 //=======================================================================
2493 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetShapesOnShapeAsCompound
2494                              (const Handle(GEOM_Object)& theCheckShape,
2495                               const Handle(GEOM_Object)&  theShape,
2496                               const Standard_Integer theShapeType,
2497                               GEOMAlgo_State theState)
2498 {
2499   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2500     getShapesOnShapeIDs (theCheckShape, theShape, theShapeType, theState);
2501
2502   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
2503     return NULL;
2504
2505   // Find objects by indices
2506   TCollection_AsciiString anAsciiList;
2507   Handle(TColStd_HSequenceOfTransient) aSeq;
2508   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
2509
2510   if ( aSeq.IsNull() || aSeq->IsEmpty() )
2511     return NULL;
2512
2513   TopoDS_Compound aCompound;
2514   BRep_Builder B;
2515   B.MakeCompound(aCompound);
2516   int i = 1;
2517   for(; i<=aSeq->Length(); i++) {
2518     Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(aSeq->Value(i));
2519     TopoDS_Shape aShape_i = anObj->GetValue();
2520     B.Add(aCompound,aShape_i);
2521   }
2522
2523   //Add a new result object
2524   Handle(GEOM_Object) aRes = GetEngine()->AddObject(GetDocID(), GEOM_SHAPES_ON_SHAPE);
2525   Handle(GEOM_Function) aFunction =
2526     aRes->AddFunction(GEOMImpl_ShapeDriver::GetID(), SHAPES_ON_SHAPE);
2527   aFunction->SetValue(aCompound);
2528
2529   GEOM::TPythonDump(aFunction)
2530     << aRes << " = geompy.GetShapesOnShapeAsCompound("
2531     << theCheckShape << ", "
2532     << theShape << ", "
2533     << TopAbs_ShapeEnum(theShapeType) << ", "
2534     << theState << ")";
2535
2536   SetErrorCode(OK);
2537
2538   return aRes;
2539 }
2540
2541 //=======================================================================
2542 //function : getShapesOnSurfaceIDs
2543   /*!
2544    * \brief Find IDs of subshapes complying with given status about surface
2545     * \param theSurface - the surface to check state of subshapes against
2546     * \param theShape - the shape to explore
2547     * \param theShapeType - type of subshape of theShape
2548     * \param theState - required state
2549     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2550    */
2551 //=======================================================================
2552 Handle(TColStd_HSequenceOfInteger)
2553   GEOMImpl_IShapesOperations::getShapesOnSurfaceIDs(const Handle(Geom_Surface)& theSurface,
2554                                                     const TopoDS_Shape&         theShape,
2555                                                     TopAbs_ShapeEnum            theShapeType,
2556                                                     GEOMAlgo_State              theState)
2557 {
2558   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
2559
2560   // Check presence of triangulation, build if need
2561   if (!CheckTriangulation(theShape)) {
2562     SetErrorCode("Cannot build triangulation on the shape");
2563     return aSeqOfIDs;
2564   }
2565
2566   // BEGIN: Mantis issue 0020961: Error on a pipe T-Shape
2567   // Compute tolerance
2568   Standard_Real T, VertMax = -RealLast();
2569   try {
2570 #if OCC_VERSION_LARGE > 0x06010000
2571     OCC_CATCH_SIGNALS;
2572 #endif
2573     for (TopExp_Explorer ExV (theShape, TopAbs_VERTEX); ExV.More(); ExV.Next()) {
2574       TopoDS_Vertex Vertex = TopoDS::Vertex(ExV.Current());
2575       T = BRep_Tool::Tolerance(Vertex);
2576       if (T > VertMax)
2577         VertMax = T;
2578     }
2579   }
2580   catch (Standard_Failure) {
2581     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2582     SetErrorCode(aFail->GetMessageString());
2583     return aSeqOfIDs;
2584   }
2585   // END: Mantis issue 0020961
2586
2587   // Call algo
2588   GEOMAlgo_FinderShapeOn1 aFinder;
2589   //Standard_Real aTol = 0.0001; // default value
2590   Standard_Real aTol = VertMax; // Mantis issue 0020961
2591
2592   aFinder.SetShape(theShape);
2593   aFinder.SetTolerance(aTol);
2594   aFinder.SetSurface(theSurface);
2595   aFinder.SetShapeType(theShapeType);
2596   aFinder.SetState(theState);
2597
2598   // Sets the minimal number of inner points for the faces that do not have own
2599   // inner points at all (for e.g. rectangular planar faces have just 2 triangles).
2600   // Default value=3
2601   aFinder.SetNbPntsMin(3);
2602   // Sets the maximal number of inner points for edges or faces.
2603   // It is usefull for the cases when this number is very big (e.g =2000) to improve
2604   // the performance. If this value =0, all inner points will be taken into account.
2605   // Default value=0
2606   aFinder.SetNbPntsMax(100);
2607
2608   aFinder.Perform();
2609
2610   // Interprete results
2611   Standard_Integer iErr = aFinder.ErrorStatus();
2612   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
2613   if (iErr) {
2614     MESSAGE(" iErr : " << iErr);
2615     TCollection_AsciiString aMsg (" iErr : ");
2616     aMsg += TCollection_AsciiString(iErr);
2617     SetErrorCode(aMsg);
2618     return aSeqOfIDs;
2619   }
2620   Standard_Integer iWrn = aFinder.WarningStatus();
2621   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
2622   if (iWrn) {
2623     MESSAGE(" *** iWrn : " << iWrn);
2624   }
2625
2626   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
2627
2628   if (listSS.Extent() < 1) {
2629     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
2630     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
2631     return aSeqOfIDs;
2632   }
2633
2634   // Fill sequence of object IDs
2635   aSeqOfIDs = new TColStd_HSequenceOfInteger;
2636
2637   TopTools_IndexedMapOfShape anIndices;
2638   TopExp::MapShapes(theShape, anIndices);
2639
2640   TopTools_ListIteratorOfListOfShape itSub (listSS);
2641   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
2642     int id = anIndices.FindIndex(itSub.Value());
2643     aSeqOfIDs->Append(id);
2644   }
2645
2646   return aSeqOfIDs;
2647 }
2648
2649 //=======================================================================
2650 //function : getObjectsShapesOn
2651 /*!
2652  * \brief Find shape objects and their entries by their ids
2653  * \param theShapeIDs - incoming shape ids
2654  * \param theShapeEntries - outgoing entries like "entry1, entry2, ..."
2655  * \retval Handle(TColStd_HSequenceOfTransient) - found shape objects
2656  */
2657 //=======================================================================
2658 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::
2659  getObjectsShapesOn(const Handle(GEOM_Object)&                theShape,
2660                     const Handle(TColStd_HSequenceOfInteger)& theShapeIDs,
2661                     TCollection_AsciiString &                 theShapeEntries)
2662 {
2663   Handle(TColStd_HSequenceOfTransient) aSeq;
2664
2665   if ( !theShapeIDs.IsNull() && theShapeIDs->Length() > 0 )
2666   {
2667     aSeq = new TColStd_HSequenceOfTransient;
2668     Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
2669     TCollection_AsciiString anEntry;
2670     for ( int i = 1; i <= theShapeIDs->Length(); ++i )
2671     {
2672       anArray->SetValue(1, theShapeIDs->Value( i ));
2673       Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theShape, anArray);
2674       aSeq->Append( anObj );
2675
2676       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
2677       if ( i != 1 ) theShapeEntries += ",";
2678       theShapeEntries += anEntry;
2679     }
2680   }
2681   return aSeq;
2682 }
2683
2684 //=======================================================================
2685 //function : getShapesOnSurface
2686 /*!
2687    * \brief Find subshapes complying with given status about surface
2688     * \param theSurface - the surface to check state of subshapes against
2689     * \param theShape - the shape to explore
2690     * \param theShapeType - type of subshape of theShape
2691     * \param theState - required state
2692     * \param theShapeEntries - outgoing entries like "entry1, entry2, ..."
2693     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2694  */
2695 //=======================================================================
2696 Handle(TColStd_HSequenceOfTransient)
2697     GEOMImpl_IShapesOperations::getShapesOnSurface(const Handle(Geom_Surface)& theSurface,
2698                                                    const Handle(GEOM_Object)&  theShape,
2699                                                    TopAbs_ShapeEnum            theShapeType,
2700                                                    GEOMAlgo_State              theState,
2701                                                    TCollection_AsciiString &   theShapeEntries)
2702 {
2703   // Find subshapes ids
2704   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2705     getShapesOnSurfaceIDs (theSurface, theShape->GetValue(), theShapeType, theState);
2706   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
2707     return NULL;
2708
2709   return getObjectsShapesOn( theShape, aSeqOfIDs, theShapeEntries );
2710 }
2711
2712 //=============================================================================
2713 /*!
2714  *  GetShapesOnPlane
2715  */
2716 //=============================================================================
2717 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnPlane
2718                                         (const Handle(GEOM_Object)& theShape,
2719                                          const Standard_Integer     theShapeType,
2720                                          const Handle(GEOM_Object)& theAx1,
2721                                          const GEOMAlgo_State       theState)
2722 {
2723   SetErrorCode(KO);
2724
2725   if (theShape.IsNull() || theAx1.IsNull()) return NULL;
2726
2727   TopoDS_Shape aShape = theShape->GetValue();
2728   TopoDS_Shape anAx1  = theAx1->GetValue();
2729
2730   if (aShape.IsNull() || anAx1.IsNull()) return NULL;
2731
2732   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2733   if ( !checkTypeShapesOn( theShapeType ))
2734     return NULL;
2735
2736   // Create plane
2737   Handle(Geom_Surface) aPlane = makePlane( anAx1 );
2738   if ( aPlane.IsNull() )
2739     return NULL;
2740
2741   // Find objects
2742   TCollection_AsciiString anAsciiList;
2743   Handle(TColStd_HSequenceOfTransient) aSeq;
2744   aSeq = getShapesOnSurface( aPlane, theShape, aShapeType, theState, anAsciiList );
2745   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2746     return NULL;
2747
2748   // Make a Python command
2749
2750   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2751   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2752
2753   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2754     << "] = geompy.GetShapesOnPlane(" << theShape << ", "
2755       << aShapeType << ", " << theAx1 << ", " << theState << ")";
2756
2757   SetErrorCode(OK);
2758   return aSeq;
2759 }
2760
2761 //=============================================================================
2762 /*!
2763  *  GetShapesOnPlaneWithLocation
2764  */
2765 //=============================================================================
2766 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnPlaneWithLocation
2767                                         (const Handle(GEOM_Object)& theShape,
2768                                          const Standard_Integer     theShapeType,
2769                                          const Handle(GEOM_Object)& theAx1,
2770                                          const Handle(GEOM_Object)& thePnt,
2771                                          const GEOMAlgo_State       theState)
2772 {
2773   SetErrorCode(KO);
2774
2775   if (theShape.IsNull() || theAx1.IsNull() || thePnt.IsNull()) return NULL;
2776
2777   TopoDS_Shape aShape = theShape->GetValue();
2778   TopoDS_Shape anAx1  = theAx1->GetValue();
2779   TopoDS_Shape anPnt = thePnt->GetValue();
2780
2781   if (aShape.IsNull() || anAx1.IsNull() || anPnt.IsNull()) return NULL;
2782
2783   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2784   if ( !checkTypeShapesOn( theShapeType ))
2785     return NULL;
2786
2787   // Create plane
2788   if ( anAx1.ShapeType() != TopAbs_EDGE || anPnt.ShapeType() != TopAbs_VERTEX ) return NULL;
2789   TopoDS_Vertex V1, V2, V3;
2790   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
2791   TopExp::Vertices(anEdge, V1, V2, Standard_True);
2792
2793   if (V1.IsNull() || V2.IsNull()) {
2794     SetErrorCode("Bad edge given for the plane normal vector");
2795     return NULL;
2796   }
2797   V3 = TopoDS::Vertex(anPnt);
2798
2799   if(V3.IsNull()) {
2800     SetErrorCode("Bad vertex given for the plane location");
2801       return NULL;
2802   }
2803   gp_Pnt aLoc = BRep_Tool::Pnt(V3);
2804   gp_Vec aVec(BRep_Tool::Pnt(V1),BRep_Tool::Pnt(V2));
2805
2806   if (aVec.Magnitude() < Precision::Confusion()) {
2807     SetErrorCode("Vector with null magnitude given");
2808     return NULL;
2809   }
2810   Handle(Geom_Surface) aPlane = new Geom_Plane(aLoc, aVec);
2811
2812   if ( aPlane.IsNull() )
2813     return NULL;
2814
2815   // Find objects
2816   TCollection_AsciiString anAsciiList;
2817   Handle(TColStd_HSequenceOfTransient) aSeq;
2818   aSeq = getShapesOnSurface( aPlane, theShape, aShapeType, theState, anAsciiList );
2819   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2820     return NULL;
2821
2822   // Make a Python command
2823
2824   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2825   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2826
2827   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2828     << "] = geompy.GetShapesOnPlaneWithLocation(" << theShape << ", "
2829     << aShapeType << ", " << theAx1 << ", "<< thePnt <<", " << theState << ")";
2830
2831   SetErrorCode(OK);
2832   return aSeq;
2833 }
2834
2835 //=============================================================================
2836 /*!
2837  *  GetShapesOnCylinder
2838  */
2839 //=============================================================================
2840 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnCylinder
2841                                           (const Handle(GEOM_Object)& theShape,
2842                                            const Standard_Integer     theShapeType,
2843                                            const Handle(GEOM_Object)& theAxis,
2844                                            const Standard_Real        theRadius,
2845                                            const GEOMAlgo_State       theState)
2846 {
2847   SetErrorCode(KO);
2848
2849   if (theShape.IsNull() || theAxis.IsNull()) return NULL;
2850
2851   TopoDS_Shape aShape = theShape->GetValue();
2852   TopoDS_Shape anAxis = theAxis->GetValue();
2853
2854   if (aShape.IsNull() || anAxis.IsNull()) return NULL;
2855
2856   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2857   if ( !checkTypeShapesOn( aShapeType ))
2858     return NULL;
2859
2860   // Create a cylinder surface
2861   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
2862   if ( aCylinder.IsNull() )
2863     return NULL;
2864
2865   // Find objects
2866   TCollection_AsciiString anAsciiList;
2867   Handle(TColStd_HSequenceOfTransient) aSeq;
2868   aSeq = getShapesOnSurface( aCylinder, theShape, aShapeType, theState, anAsciiList );
2869   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2870     return NULL;
2871
2872   // Make a Python command
2873
2874   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2875   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2876
2877   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2878     << "] = geompy.GetShapesOnCylinder(" << theShape << ", " << aShapeType
2879       << ", " << theAxis << ", " << theRadius << ", " << theState << ")";
2880
2881   SetErrorCode(OK);
2882   return aSeq;
2883 }
2884
2885 //=============================================================================
2886 /*!
2887  *  GetShapesOnCylinderWithLocation
2888  */
2889 //=============================================================================
2890 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnCylinderWithLocation
2891                                           (const Handle(GEOM_Object)& theShape,
2892                                            const Standard_Integer     theShapeType,
2893                                            const Handle(GEOM_Object)& theAxis,
2894                                            const Handle(GEOM_Object)& thePnt,
2895                                            const Standard_Real        theRadius,
2896                                            const GEOMAlgo_State       theState)
2897 {
2898   SetErrorCode(KO);
2899
2900   if (theShape.IsNull() || theAxis.IsNull() || thePnt.IsNull()) return NULL;
2901
2902   TopoDS_Shape aShape = theShape->GetValue();
2903   TopoDS_Shape anAxis = theAxis->GetValue();
2904   TopoDS_Shape aPnt   = thePnt->GetValue();
2905
2906   if (aShape.IsNull() || anAxis.IsNull() || aPnt.IsNull()) return NULL;
2907
2908   if (aPnt.ShapeType() != TopAbs_VERTEX )
2909   {
2910     SetErrorCode("Bottom location point must be vertex");
2911     return NULL;
2912   }
2913
2914   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2915   if ( !checkTypeShapesOn( aShapeType ))
2916     return NULL;
2917
2918   // Create a cylinder surface
2919   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
2920   if ( aCylinder.IsNull() )
2921     return NULL;
2922
2923   // translate the surface
2924   Handle(Geom_CylindricalSurface) aCylSurface =
2925     Handle(Geom_CylindricalSurface)::DownCast( aCylinder );
2926   if ( aCylSurface.IsNull() )
2927   {
2928     SetErrorCode("Unexpected surface type instead of Geom_CylindricalSurface");
2929     return NULL;
2930   }
2931   gp_Pnt fromLoc = aCylSurface->Cylinder().Location();
2932   gp_Pnt toLoc   = BRep_Tool::Pnt( TopoDS::Vertex( aPnt ));
2933   aCylinder->Translate( fromLoc, toLoc );
2934
2935   // Find objects
2936   TCollection_AsciiString anAsciiList;
2937   Handle(TColStd_HSequenceOfTransient) aSeq;
2938   aSeq = getShapesOnSurface( aCylinder, theShape, aShapeType, theState, anAsciiList );
2939   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2940     return NULL;
2941
2942   // Make a Python command
2943
2944   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2945   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2946
2947   GEOM::TPythonDump(aFunction)
2948     << "[" << anAsciiList.ToCString()
2949     << "] = geompy.GetShapesOnCylinderWithLocation(" << theShape << ", " << aShapeType << ", "
2950     << theAxis << ", " << thePnt << ", " << theRadius << ", " << theState << ")";
2951
2952   SetErrorCode(OK);
2953   return aSeq;
2954 }
2955
2956 //=============================================================================
2957 /*!
2958  *  GetShapesOnSphere
2959  */
2960 //=============================================================================
2961 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnSphere
2962                                           (const Handle(GEOM_Object)& theShape,
2963                                            const Standard_Integer     theShapeType,
2964                                            const Handle(GEOM_Object)& theCenter,
2965                                            const Standard_Real        theRadius,
2966                                            const GEOMAlgo_State       theState)
2967 {
2968   SetErrorCode(KO);
2969
2970   if (theShape.IsNull() || theCenter.IsNull()) return NULL;
2971
2972   TopoDS_Shape aShape  = theShape->GetValue();
2973   TopoDS_Shape aCenter = theCenter->GetValue();
2974
2975   if (aShape.IsNull() || aCenter.IsNull()) return NULL;
2976
2977   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2978   if ( !checkTypeShapesOn( aShapeType ))
2979     return NULL;
2980
2981   // Center of the sphere
2982   if (aCenter.ShapeType() != TopAbs_VERTEX) return NULL;
2983   gp_Pnt aLoc = BRep_Tool::Pnt(TopoDS::Vertex(aCenter));
2984
2985   gp_Ax3 anAx3 (aLoc, gp::DZ());
2986   Handle(Geom_SphericalSurface) aSphere =
2987     new Geom_SphericalSurface(anAx3, theRadius);
2988
2989   // Find objects
2990   TCollection_AsciiString anAsciiList;
2991   Handle(TColStd_HSequenceOfTransient) aSeq;
2992   aSeq = getShapesOnSurface( aSphere, theShape, aShapeType, theState, anAsciiList );
2993   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2994     return NULL;
2995
2996   // Make a Python command
2997
2998   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2999   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
3000
3001   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
3002     << "] = geompy.GetShapesOnSphere(" << theShape << ", " << aShapeType
3003       << ", " << theCenter << ", " << theRadius << ", " << theState << ")";
3004
3005   SetErrorCode(OK);
3006   return aSeq;
3007 }
3008
3009 //=============================================================================
3010 /*!
3011  *  GetShapesOnPlaneIDs
3012  */
3013 //=============================================================================
3014 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnPlaneIDs
3015                                         (const Handle(GEOM_Object)& theShape,
3016                                          const Standard_Integer     theShapeType,
3017                                          const Handle(GEOM_Object)& theAx1,
3018                                          const GEOMAlgo_State       theState)
3019 {
3020   SetErrorCode(KO);
3021
3022   if (theShape.IsNull() || theAx1.IsNull()) return NULL;
3023
3024   TopoDS_Shape aShape = theShape->GetValue();
3025   TopoDS_Shape anAx1  = theAx1->GetValue();
3026
3027   if (aShape.IsNull() || anAx1.IsNull()) return NULL;
3028
3029   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
3030   if ( !checkTypeShapesOn( aShapeType ))
3031     return NULL;
3032
3033   // Create plane
3034   Handle(Geom_Surface) aPlane = makePlane( anAx1 );
3035   if ( aPlane.IsNull() )
3036     return NULL;
3037
3038   // Find object IDs
3039   Handle(TColStd_HSequenceOfInteger) aSeq;
3040   aSeq = getShapesOnSurfaceIDs( aPlane, aShape, aShapeType, theState );
3041
3042   // The GetShapesOnPlaneIDs() doesn't change object so no new function is required.
3043   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theAx1)->GetLastFunction();
3044
3045   // Make a Python command
3046   GEOM::TPythonDump(aFunction, /*append=*/true)
3047     << "listShapesOnPlane = geompy.GetShapesOnPlaneIDs"
3048     << "(" << theShape << "," << aShapeType << "," << theAx1 << "," << theState << ")";
3049
3050   SetErrorCode(OK);
3051   return aSeq;
3052 }
3053
3054 //=============================================================================
3055 /*!
3056  *  GetShapesOnPlaneWithLocationIDs
3057  */
3058 //=============================================================================
3059 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnPlaneWithLocationIDs
3060                                         (const Handle(GEOM_Object)& theShape,
3061                                          const Standard_Integer     theShapeType,
3062                                          const Handle(GEOM_Object)& theAx1,
3063                                          const Handle(GEOM_Object)& thePnt,
3064                                          const GEOMAlgo_State       theState)
3065 {
3066   SetErrorCode(KO);
3067
3068   if (theShape.IsNull() || theAx1.IsNull() || thePnt.IsNull()) return NULL;
3069
3070   TopoDS_Shape aShape = theShape->GetValue();
3071   TopoDS_Shape anAx1  = theAx1->GetValue();
3072   TopoDS_Shape anPnt  = thePnt->GetValue();
3073
3074   if (aShape.IsNull() || anAx1.IsNull() || anPnt.IsNull()) return NULL;
3075
3076   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
3077   if ( !checkTypeShapesOn( aShapeType ))
3078     return NULL;
3079
3080   // Create plane
3081   if (anAx1.ShapeType() != TopAbs_EDGE || anPnt.ShapeType() != TopAbs_VERTEX) return NULL;
3082   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
3083   TopoDS_Vertex V1, V2, V3;
3084   TopExp::Vertices(anEdge, V1, V2, Standard_True);
3085   if (V1.IsNull() || V2.IsNull()) {
3086     SetErrorCode("Bad edge given for the plane normal vector");
3087     return NULL;
3088   }
3089   V3 = TopoDS::Vertex(anPnt);
3090   if(V3.IsNull()) {
3091     SetErrorCode("Bad vertex given for the plane location");
3092       return NULL;
3093   }
3094   gp_Pnt aLoc = BRep_Tool::Pnt(V3);
3095   gp_Vec aVec(BRep_Tool::Pnt(V1),BRep_Tool::Pnt(V2));
3096   if (aVec.Magnitude() < Precision::Confusion()) {
3097     SetErrorCode("Vector with null magnitude given");
3098     return NULL;
3099   }
3100
3101   Handle(Geom_Surface) aPlane = new Geom_Plane(aLoc, aVec);
3102   if ( aPlane.IsNull() )
3103     return NULL;
3104
3105   // Find object IDs
3106   Handle(TColStd_HSequenceOfInteger) aSeq;
3107   aSeq = getShapesOnSurfaceIDs( aPlane, aShape, aShapeType, theState );
3108
3109   // The GetShapesOnPlaneIDs() doesn't change object so no new function is required.
3110   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theAx1)->GetLastFunction();
3111
3112   // Make a Python command
3113   GEOM::TPythonDump(aFunction, /*append=*/true)
3114     << "listShapesOnPlane = geompy.GetShapesOnPlaneWithLocationIDs"
3115     << "(" << theShape << ", " << aShapeType << ", " << theAx1 << ", "<< thePnt << ", "  << theState << ")";
3116
3117   SetErrorCode(OK);
3118   return aSeq;
3119 }
3120
3121 //=============================================================================
3122 /*!
3123  *  GetShapesOnCylinderIDs
3124  */
3125 //=============================================================================
3126 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnCylinderIDs
3127                                           (const Handle(GEOM_Object)& theShape,
3128                                            const Standard_Integer     theShapeType,
3129                                            const Handle(GEOM_Object)& theAxis,
3130                                            const Standard_Real        theRadius,
3131                                            const GEOMAlgo_State       theState)
3132 {
3133   SetErrorCode(KO);
3134
3135   if (theShape.IsNull() || theAxis.IsNull()) return NULL;
3136
3137   TopoDS_Shape aShape = theShape->GetValue();
3138   TopoDS_Shape anAxis = theAxis->GetValue();
3139
3140   if (aShape.IsNull() || anAxis.IsNull()) return NULL;
3141
3142   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
3143   if ( !checkTypeShapesOn( aShapeType ))
3144     return NULL;
3145
3146   // Create a cylinder surface
3147   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
3148   if ( aCylinder.IsNull() )
3149     return NULL;
3150
3151   // Find object IDs
3152   Handle(TColStd_HSequenceOfInteger) aSeq;
3153   aSeq = getShapesOnSurfaceIDs( aCylinder, aShape, aShapeType, theState );
3154
3155   // The GetShapesOnCylinder() doesn't change object so no new function is required.
3156   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theAxis)->GetLastFunction();
3157
3158   // Make a Python command
3159   GEOM::TPythonDump(aFunction, /*append=*/true)
3160     << "listShapesOnCylinder = geompy.GetShapesOnCylinderIDs"
3161     << "(" << theShape << ", " << aShapeType << ", " << theAxis << ", "
3162     << theRadius << ", " << theState << ")";
3163
3164   SetErrorCode(OK);
3165   return aSeq;
3166 }
3167
3168 //=============================================================================
3169 /*!
3170  *  GetShapesOnCylinderWithLocationIDs
3171  */
3172 //=============================================================================
3173 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnCylinderWithLocationIDs
3174                                           (const Handle(GEOM_Object)& theShape,
3175                                            const Standard_Integer     theShapeType,
3176                                            const Handle(GEOM_Object)& theAxis,
3177                                            const Handle(GEOM_Object)& thePnt,
3178                                            const Standard_Real        theRadius,
3179                                            const GEOMAlgo_State       theState)
3180 {
3181   SetErrorCode(KO);
3182
3183   if (theShape.IsNull() || theAxis.IsNull() || thePnt.IsNull()) return NULL;
3184
3185   TopoDS_Shape aShape = theShape->GetValue();
3186   TopoDS_Shape anAxis = theAxis->GetValue();
3187   TopoDS_Shape aPnt   = thePnt->GetValue();
3188
3189   if (aShape.IsNull() || anAxis.IsNull() || aPnt.IsNull()) return NULL;
3190
3191   if (aPnt.ShapeType() != TopAbs_VERTEX )
3192   {
3193     SetErrorCode("Bottom location point must be vertex");
3194     return NULL;
3195   }
3196
3197   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
3198   if ( !checkTypeShapesOn( aShapeType ))
3199     return NULL;
3200
3201   // Create a cylinder surface
3202   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
3203   if ( aCylinder.IsNull() )
3204     return NULL;
3205
3206   // translate the surface
3207   Handle(Geom_CylindricalSurface) aCylSurface =
3208     Handle(Geom_CylindricalSurface)::DownCast( aCylinder );
3209   if ( aCylSurface.IsNull() )
3210   {
3211     SetErrorCode("Unexpected surface type instead of Geom_CylindricalSurface");
3212     return NULL;
3213   }
3214   gp_Pnt fromLoc = aCylSurface->Cylinder().Location();
3215   gp_Pnt toLoc   = BRep_Tool::Pnt( TopoDS::Vertex( aPnt ));
3216   aCylinder->Translate( fromLoc, toLoc );
3217
3218   // Find object IDs
3219   Handle(TColStd_HSequenceOfInteger) aSeq;
3220   aSeq = getShapesOnSurfaceIDs( aCylinder, aShape, aShapeType, theState );
3221
3222   // The GetShapesOnCylinder() doesn't change object so no new function is required.
3223   Handle(GEOM_Function) aFunction =
3224     GEOM::GetCreatedLast(theShape, GEOM::GetCreatedLast(thePnt,theAxis))->GetLastFunction();
3225
3226   // Make a Python command
3227   GEOM::TPythonDump(aFunction, /*append=*/true)
3228     << "listShapesOnCylinder = geompy.GetShapesOnCylinderWithLocationIDs"
3229     << "(" << theShape << ", " << aShapeType << ", " << theAxis << ", "
3230     << thePnt << ", " << theRadius << ", " << theState << ")";
3231
3232   SetErrorCode(OK);
3233   return aSeq;
3234 }
3235
3236 //=============================================================================
3237 /*!
3238  *  GetShapesOnSphereIDs
3239  */
3240 //=============================================================================
3241 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnSphereIDs
3242                                           (const Handle(GEOM_Object)& theShape,
3243                                            const Standard_Integer     theShapeType,
3244                                            const Handle(GEOM_Object)& theCenter,
3245                                            const Standard_Real        theRadius,
3246                                            const GEOMAlgo_State       theState)
3247 {
3248   SetErrorCode(KO);
3249
3250   if (theShape.IsNull() || theCenter.IsNull()) return NULL;
3251
3252   TopoDS_Shape aShape  = theShape->GetValue();
3253   TopoDS_Shape aCenter = theCenter->GetValue();
3254
3255   if (aShape.IsNull() || aCenter.IsNull()) return NULL;
3256
3257   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
3258   if ( !checkTypeShapesOn( aShapeType ))
3259     return NULL;
3260
3261   // Center of the sphere
3262   if (aCenter.ShapeType() != TopAbs_VERTEX) return NULL;
3263   gp_Pnt aLoc = BRep_Tool::Pnt(TopoDS::Vertex(aCenter));
3264
3265   gp_Ax3 anAx3 (aLoc, gp::DZ());
3266   Handle(Geom_SphericalSurface) aSphere =
3267     new Geom_SphericalSurface(anAx3, theRadius);
3268
3269   // Find object IDs
3270   Handle(TColStd_HSequenceOfInteger) aSeq;
3271   aSeq = getShapesOnSurfaceIDs( aSphere, aShape, aShapeType, theState );
3272
3273   // The GetShapesOnSphere() doesn't change object so no new function is required.
3274   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theCenter)->GetLastFunction();
3275
3276   // Make a Python command
3277   GEOM::TPythonDump(aFunction, /*append=*/true)
3278     << "listShapesOnCylinder = geompy.GetShapesOnCylinderIDs"
3279     << "(" << theShape << ", " << aShapeType << ", " << theCenter << ", "
3280     << theRadius << ", " << theState << ")";
3281
3282   SetErrorCode(OK);
3283   return aSeq;
3284 }
3285
3286 //=======================================================================
3287 //function : getShapesOnQuadrangleIDs
3288   /*!
3289    * \brief Find IDs of subshapes complying with given status about quadrangle
3290     * \param theShape - the shape to explore
3291     * \param theShapeType - type of subshape of theShape
3292     * \param theTopLeftPoint - top left quadrangle corner
3293     * \param theTopRigthPoint - top right quadrangle corner
3294     * \param theBottomLeftPoint - bottom left quadrangle corner
3295     * \param theBottomRigthPoint - bottom right quadrangle corner
3296     * \param theState - required state
3297     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
3298    */
3299 //=======================================================================
3300 Handle(TColStd_HSequenceOfInteger)
3301   GEOMImpl_IShapesOperations::getShapesOnQuadrangleIDs (const Handle(GEOM_Object)& theShape,
3302                                                         const Standard_Integer     theShapeType,
3303                                                         const Handle(GEOM_Object)& theTopLeftPoint,
3304                                                         const Handle(GEOM_Object)& theTopRigthPoint,
3305                                                         const Handle(GEOM_Object)& theBottomLeftPoint,
3306                                                         const Handle(GEOM_Object)& theBottomRigthPoint,
3307                                                         const GEOMAlgo_State       theState)
3308 {
3309   SetErrorCode(KO);
3310
3311   if ( theShape.IsNull() ||
3312        theTopLeftPoint.IsNull() ||
3313        theTopRigthPoint.IsNull() ||
3314        theBottomLeftPoint.IsNull() ||
3315        theBottomRigthPoint.IsNull() )
3316     return NULL;
3317
3318   TopoDS_Shape aShape = theShape->GetValue();
3319   TopoDS_Shape aTL = theTopLeftPoint->GetValue();
3320   TopoDS_Shape aTR = theTopRigthPoint->GetValue();
3321   TopoDS_Shape aBL = theBottomLeftPoint->GetValue();
3322   TopoDS_Shape aBR = theBottomRigthPoint->GetValue();
3323
3324   if (aShape.IsNull() ||
3325       aTL.IsNull() ||
3326       aTR.IsNull() ||
3327       aBL.IsNull() ||
3328       aBR.IsNull() ||
3329       aTL.ShapeType() != TopAbs_VERTEX ||
3330       aTR.ShapeType() != TopAbs_VERTEX ||
3331       aBL.ShapeType() != TopAbs_VERTEX ||
3332       aBR.ShapeType() != TopAbs_VERTEX )
3333     return NULL;
3334
3335   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
3336   if ( !checkTypeShapesOn( aShapeType ))
3337     return NULL;
3338
3339   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
3340
3341   // Check presence of triangulation, build if need
3342   if (!CheckTriangulation(aShape)) {
3343     SetErrorCode("Cannot build triangulation on the shape");
3344     return aSeqOfIDs;
3345   }
3346
3347   // Call algo
3348   gp_Pnt aPntTL = BRep_Tool::Pnt(TopoDS::Vertex(aTL));
3349   gp_Pnt aPntTR = BRep_Tool::Pnt(TopoDS::Vertex(aTR));
3350   gp_Pnt aPntBL = BRep_Tool::Pnt(TopoDS::Vertex(aBL));
3351   gp_Pnt aPntBR = BRep_Tool::Pnt(TopoDS::Vertex(aBR));
3352
3353   GEOMAlgo_FinderShapeOnQuad aFinder( aPntTL, aPntTR, aPntBL, aPntBR );
3354   Standard_Real aTol = 0.0001; // default value
3355
3356   aFinder.SetShape(aShape);
3357   aFinder.SetTolerance(aTol);
3358   //aFinder.SetSurface(theSurface);
3359   aFinder.SetShapeType(aShapeType);
3360   aFinder.SetState(theState);
3361
3362   // Sets the minimal number of inner points for the faces that do not have own
3363   // inner points at all (for e.g. rectangular planar faces have just 2 triangles).
3364   // Default value=3
3365   aFinder.SetNbPntsMin(3);
3366   // Sets the maximal number of inner points for edges or faces.
3367   // It is usefull for the cases when this number is very big (e.g =2000) to improve
3368   // the performance. If this value =0, all inner points will be taken into account.
3369   // Default value=0
3370   aFinder.SetNbPntsMax(100);
3371
3372   aFinder.Perform();
3373
3374   // Interprete results
3375   Standard_Integer iErr = aFinder.ErrorStatus();
3376   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
3377   if (iErr) {
3378     MESSAGE(" iErr : " << iErr);
3379     TCollection_AsciiString aMsg (" iErr : ");
3380     aMsg += TCollection_AsciiString(iErr);
3381     SetErrorCode(aMsg);
3382     return aSeqOfIDs;
3383   }
3384   Standard_Integer iWrn = aFinder.WarningStatus();
3385   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
3386   if (iWrn) {
3387     MESSAGE(" *** iWrn : " << iWrn);
3388   }
3389
3390   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
3391
3392   if (listSS.Extent() < 1) {
3393     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
3394     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
3395     return aSeqOfIDs;
3396   }
3397
3398   // Fill sequence of object IDs
3399   aSeqOfIDs = new TColStd_HSequenceOfInteger;
3400
3401   TopTools_IndexedMapOfShape anIndices;
3402   TopExp::MapShapes(aShape, anIndices);
3403
3404   TopTools_ListIteratorOfListOfShape itSub (listSS);
3405   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
3406     int id = anIndices.FindIndex(itSub.Value());
3407     aSeqOfIDs->Append(id);
3408   }
3409   return aSeqOfIDs;
3410 }
3411
3412 //=======================================================================
3413 //function : GetShapesOnQuadrangle
3414   /*!
3415    * \brief Find subshapes complying with given status about quadrangle
3416     * \param theShape - the shape to explore
3417     * \param theShapeType - type of subshape of theShape
3418     * \param theTopLeftPoint - top left quadrangle corner
3419     * \param theTopRigthPoint - top right quadrangle corner
3420     * \param theBottomLeftPoint - bottom left quadrangle corner
3421     * \param theBottomRigthPoint - bottom right quadrangle corner
3422     * \param theState - required state
3423     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
3424    */
3425 //=======================================================================
3426 Handle(TColStd_HSequenceOfTransient)
3427     GEOMImpl_IShapesOperations::GetShapesOnQuadrangle (const Handle(GEOM_Object)& theShape,
3428                                                        const Standard_Integer     theShapeType,
3429                                                        const Handle(GEOM_Object)& theTopLeftPoint,
3430                                                        const Handle(GEOM_Object)& theTopRigthPoint,
3431                                                        const Handle(GEOM_Object)& theBottomLeftPoint,
3432                                                        const Handle(GEOM_Object)& theBottomRigthPoint,
3433                                                        const GEOMAlgo_State       theState)
3434 {
3435   // Find indices
3436   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
3437     getShapesOnQuadrangleIDs( theShape,
3438                               theShapeType,
3439                               theTopLeftPoint,
3440                               theTopRigthPoint,
3441                               theBottomLeftPoint,
3442                               theBottomRigthPoint,
3443                               theState);
3444   if ( aSeqOfIDs.IsNull() || aSeqOfIDs->IsEmpty() )
3445     return NULL;
3446
3447   // Find objects by indices
3448   TCollection_AsciiString anAsciiList;
3449   Handle(TColStd_HSequenceOfTransient) aSeq;
3450   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
3451   if ( aSeq.IsNull() || aSeq->IsEmpty() )
3452     return NULL;
3453
3454   // Make a Python command
3455
3456   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
3457   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
3458
3459   GEOM::TPythonDump(aFunction)
3460     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnQuadrangle("
3461     << theShape << ", "
3462     << TopAbs_ShapeEnum(theShapeType) << ", "
3463     << theTopLeftPoint << ", "
3464     << theTopRigthPoint << ", "
3465     << theBottomLeftPoint << ", "
3466     << theBottomRigthPoint << ", "
3467     << theState << ")";
3468
3469   SetErrorCode(OK);
3470   return aSeq;
3471 }
3472
3473 //=======================================================================
3474 //function : GetShapesOnQuadrangleIDs
3475   /*!
3476    * \brief Find IDs of subshapes complying with given status about quadrangle
3477     * \param theShape - the shape to explore
3478     * \param theShapeType - type of subshape of theShape
3479     * \param theTopLeftPoint - top left quadrangle corner
3480     * \param theTopRigthPoint - top right quadrangle corner
3481     * \param theBottomLeftPoint - bottom left quadrangle corner
3482     * \param theBottomRigthPoint - bottom right quadrangle corner
3483     * \param theState - required state
3484     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
3485    */
3486 //=======================================================================
3487 Handle(TColStd_HSequenceOfInteger)
3488   GEOMImpl_IShapesOperations::GetShapesOnQuadrangleIDs (const Handle(GEOM_Object)& theShape,
3489                                                         const Standard_Integer     theShapeType,
3490                                                         const Handle(GEOM_Object)& theTopLeftPoint,
3491                                                         const Handle(GEOM_Object)& theTopRigthPoint,
3492                                                         const Handle(GEOM_Object)& theBottomLeftPoint,
3493                                                         const Handle(GEOM_Object)& theBottomRigthPoint,
3494                                                         const GEOMAlgo_State       theState)
3495 {
3496   // Find indices
3497   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
3498     getShapesOnQuadrangleIDs( theShape,
3499                               theShapeType,
3500                               theTopLeftPoint,
3501                               theTopRigthPoint,
3502                               theBottomLeftPoint,
3503                               theBottomRigthPoint,
3504                               theState);
3505   if ( aSeqOfIDs.IsNull() || aSeqOfIDs->IsEmpty() )
3506     return NULL;
3507
3508   // Make a Python command
3509
3510   // The GetShapesOnCylinder() doesn't change object so no new function is required.
3511   Handle(GEOM_Object) lastObj = GEOM::GetCreatedLast(theShape,theTopLeftPoint);
3512   lastObj = GEOM::GetCreatedLast(lastObj,theTopRigthPoint);
3513   lastObj = GEOM::GetCreatedLast(lastObj,theBottomRigthPoint);
3514   lastObj = GEOM::GetCreatedLast(lastObj,theBottomLeftPoint);
3515   Handle(GEOM_Function) aFunction = lastObj->GetLastFunction();
3516
3517   GEOM::TPythonDump(aFunction, /*append=*/true)
3518     << "listShapesOnQuadrangle = geompy.GetShapesOnQuadrangleIDs("
3519     << theShape << ", "
3520     << TopAbs_ShapeEnum(theShapeType) << ", "
3521     << theTopLeftPoint << ", "
3522     << theTopRigthPoint << ", "
3523     << theBottomLeftPoint << ", "
3524     << theBottomRigthPoint << ", "
3525     << theState << ")";
3526
3527   SetErrorCode(OK);
3528   return aSeqOfIDs;
3529 }
3530
3531 //=============================================================================
3532 /*!
3533  *  GetInPlaceOfShape
3534  */
3535 //=============================================================================
3536 static bool GetInPlaceOfShape (const Handle(GEOM_Function)& theWhereFunction,
3537                                const TopTools_IndexedMapOfShape& theWhereIndices,
3538                                const TopoDS_Shape& theWhat,
3539                                TColStd_ListOfInteger& theModifiedList)
3540 {
3541   if (theWhereFunction.IsNull() || theWhat.IsNull()) return false;
3542
3543   if (theWhereIndices.Contains(theWhat)) {
3544     // entity was not changed by the operation
3545     Standard_Integer aWhatIndex = theWhereIndices.FindIndex(theWhat);
3546     theModifiedList.Append(aWhatIndex);
3547     return true;
3548   }
3549
3550   // try to find in history
3551   TDF_Label aHistoryLabel = theWhereFunction->GetHistoryEntry(Standard_False);
3552
3553   // search in history for all argument shapes
3554   Standard_Boolean isFound = Standard_False;
3555   Standard_Boolean isGood = Standard_False;
3556
3557   TDF_LabelSequence aLabelSeq;
3558   theWhereFunction->GetDependency(aLabelSeq);
3559   Standard_Integer nbArg = aLabelSeq.Length();
3560
3561   for (Standard_Integer iarg = 1; iarg <= nbArg && !isFound; iarg++) {
3562
3563     TDF_Label anArgumentRefLabel = aLabelSeq.Value(iarg);
3564
3565     Handle(GEOM_Object) anArgumentObject = GEOM_Object::GetReferencedObject(anArgumentRefLabel);
3566     TopoDS_Shape anArgumentShape = anArgumentObject->GetValue();
3567
3568     TopTools_IndexedMapOfShape anArgumentIndices;
3569     TopExp::MapShapes(anArgumentShape, anArgumentIndices);
3570
3571     if (anArgumentIndices.Contains(theWhat)) {
3572       isFound = Standard_True;
3573       Standard_Integer aWhatIndex = anArgumentIndices.FindIndex(theWhat);
3574
3575       // Find corresponding label in history
3576       TDF_Label anArgumentHistoryLabel =
3577         theWhereFunction->GetArgumentHistoryEntry(anArgumentRefLabel, Standard_False);
3578       if (anArgumentHistoryLabel.IsNull()) {
3579         // Lost History of operation argument. Possibly, all its entities was removed.
3580         isGood = Standard_True;
3581       }
3582       else {
3583         TDF_Label aWhatHistoryLabel = anArgumentHistoryLabel.FindChild(aWhatIndex, Standard_False);
3584
3585         if (aWhatHistoryLabel.IsNull()) {
3586           // Removed entity ? Compound ? Compsolid ? Shell ? Wire
3587           isGood = Standard_False;
3588         } else {
3589           Handle(TDataStd_IntegerArray) anIntegerArray;
3590           if (!aWhatHistoryLabel.FindAttribute(TDataStd_IntegerArray::GetID(), anIntegerArray)) {
3591             //Error: Empty modifications history for the sought shape.
3592             isGood = Standard_False;
3593           }
3594           else {
3595             isGood = Standard_True;
3596             Standard_Integer imod, aModifLen = anIntegerArray->Array()->Length();
3597             for (imod = 1; imod <= aModifLen; imod++) {
3598               theModifiedList.Append(anIntegerArray->Array()->Value(imod));
3599             }
3600           }
3601         }
3602       }
3603     }
3604   }
3605
3606   isFound = isGood;
3607
3608   if (!isFound) {
3609     // try compound/compsolid/shell/wire element by element
3610     bool isFoundAny = false;
3611     TopTools_MapOfShape mapShape;
3612
3613     if (theWhat.ShapeType() == TopAbs_COMPOUND ||
3614         theWhat.ShapeType() == TopAbs_COMPSOLID) {
3615       // recursive processing of compound/compsolid
3616       TopoDS_Iterator anIt (theWhat, Standard_True, Standard_True);
3617       for (; anIt.More(); anIt.Next()) {
3618         if (mapShape.Add(anIt.Value())) {
3619           TopoDS_Shape curWhat = anIt.Value();
3620           isFoundAny = GetInPlaceOfShape(theWhereFunction, theWhereIndices, curWhat, theModifiedList);
3621           if (isFoundAny) isFound = Standard_True;
3622         }
3623       }
3624     }
3625     else if (theWhat.ShapeType() == TopAbs_SHELL) {
3626       // try to replace a shell by its faces images
3627       TopExp_Explorer anExp (theWhat, TopAbs_FACE);
3628       for (; anExp.More(); anExp.Next()) {
3629         if (mapShape.Add(anExp.Current())) {
3630           TopoDS_Shape curWhat = anExp.Current();
3631           isFoundAny = GetInPlaceOfShape(theWhereFunction, theWhereIndices, curWhat, theModifiedList);
3632           if (isFoundAny) isFound = Standard_True;
3633         }
3634       }
3635     }
3636     else if (theWhat.ShapeType() == TopAbs_WIRE) {
3637       // try to replace a wire by its edges images
3638       TopExp_Explorer anExp (theWhat, TopAbs_EDGE);
3639       for (; anExp.More(); anExp.Next()) {
3640         if (mapShape.Add(anExp.Current())) {
3641           TopoDS_Shape curWhat = anExp.Current();
3642           isFoundAny = GetInPlaceOfShape(theWhereFunction, theWhereIndices, curWhat, theModifiedList);
3643           if (isFoundAny) isFound = Standard_True;
3644         }
3645       }
3646     }
3647     else {
3648       // Removed entity
3649     }
3650   }
3651
3652   return isFound;
3653 }
3654
3655 //=============================================================================
3656 /*!
3657  *  GetShapeProperties
3658  */
3659 //=============================================================================
3660 void GEOMImpl_IShapesOperations::GetShapeProperties( const TopoDS_Shape aShape, Standard_Real tab[],
3661                                                      gp_Pnt & aVertex )
3662 {
3663   GProp_GProps theProps;
3664   gp_Pnt aCenterMass;
3665   //TopoDS_Shape aPntShape;
3666   Standard_Real aShapeSize;
3667
3668   if    (aShape.ShapeType() == TopAbs_VERTEX) aCenterMass = BRep_Tool::Pnt( TopoDS::Vertex( aShape ) );
3669   else if (aShape.ShapeType() == TopAbs_EDGE) BRepGProp::LinearProperties(aShape,  theProps);
3670   else if (aShape.ShapeType() == TopAbs_FACE) BRepGProp::SurfaceProperties(aShape, theProps);
3671   else                                        BRepGProp::VolumeProperties(aShape,  theProps);
3672
3673   if (aShape.ShapeType() == TopAbs_VERTEX)
3674     aShapeSize = 1;
3675   else {
3676     aCenterMass = theProps.CentreOfMass();
3677     aShapeSize  = theProps.Mass();
3678   }
3679
3680 //   aPntShape = BRepBuilderAPI_MakeVertex(aCenterMass).Shape();
3681 //   aVertex   = BRep_Tool::Pnt( TopoDS::Vertex( aPntShape ) );
3682   aVertex = aCenterMass;
3683   tab[0] = aVertex.X();
3684   tab[1] = aVertex.Y();
3685   tab[2] = aVertex.Z();
3686   tab[3] = aShapeSize;
3687   return;
3688 }
3689
3690 namespace {
3691
3692   //================================================================================
3693   /*!
3694    * \brief Return normal to face at extrema point
3695    */
3696   //================================================================================
3697
3698   gp_Vec GetNormal (const TopoDS_Face& face, const BRepExtrema_DistShapeShape& extrema)
3699   {
3700     gp_Vec defaultNorm(1,0,0); // to have same normals on different faces
3701     try {
3702       // get UV at extrema point
3703       Standard_Real u,v, f,l;
3704       switch ( extrema.SupportTypeShape2(1) ) {
3705       case BRepExtrema_IsInFace: {
3706         extrema.ParOnFaceS2(1, u, v );
3707         break;
3708       }
3709       case BRepExtrema_IsOnEdge: {
3710         TopoDS_Edge edge = TopoDS::Edge( extrema.SupportOnShape2(1));
3711         Handle(Geom2d_Curve) pcurve = BRep_Tool::CurveOnSurface( edge, face, f,l );
3712         extrema.ParOnEdgeS2( 1, u );
3713         gp_Pnt2d uv = pcurve->Value( u );
3714         u = uv.Coord(1);
3715         v = uv.Coord(2);
3716         break;
3717       }
3718       case BRepExtrema_IsVertex: return defaultNorm;
3719       }
3720       // get derivatives
3721       BRepAdaptor_Surface surface( face, false );
3722       gp_Vec du, dv; gp_Pnt p;
3723       surface.D1( u, v, p, du, dv );
3724
3725       return du ^ dv;
3726
3727     } catch (Standard_Failure ) {
3728     }
3729     return defaultNorm;
3730   }
3731 }
3732
3733 //================================================================================
3734 /*!
3735  * \brief Return type of shape for explode. In case of compound it will be a type of sub shape.
3736  */
3737 //================================================================================
3738 TopAbs_ShapeEnum GEOMImpl_IShapesOperations::GetTypeOfSimplePart (const TopoDS_Shape& theShape)
3739 {
3740   TopAbs_ShapeEnum aType = theShape.ShapeType();
3741   if      (aType == TopAbs_VERTEX)                             return TopAbs_VERTEX;
3742   else if (aType == TopAbs_EDGE  || aType == TopAbs_WIRE)      return TopAbs_EDGE;
3743   else if (aType == TopAbs_FACE  || aType == TopAbs_SHELL)     return TopAbs_FACE;
3744   else if (aType == TopAbs_SOLID || aType == TopAbs_COMPSOLID) return TopAbs_SOLID;
3745   else if (aType == TopAbs_COMPOUND) {
3746     // Only the iType of the first shape in the compound is taken into account
3747     TopoDS_Iterator It (theShape, Standard_False, Standard_False);
3748     if (It.More()) {
3749       return GetTypeOfSimplePart(It.Value());
3750     }
3751   }
3752   return TopAbs_SHAPE;
3753 }
3754
3755 //=============================================================================
3756 /*!
3757  *  case GetInPlace:
3758  *  default:
3759  */
3760 //=============================================================================
3761 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetInPlace (Handle(GEOM_Object) theShapeWhere,
3762                                                             Handle(GEOM_Object) theShapeWhat)
3763 {
3764   SetErrorCode(KO);
3765
3766   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
3767
3768   TopoDS_Shape aWhere = theShapeWhere->GetValue();
3769   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
3770   TopoDS_Shape aPntShape;
3771   TopoDS_Vertex aVertex;
3772
3773   if (aWhere.IsNull() || aWhat.IsNull()) {
3774     SetErrorCode("Error: aWhere and aWhat TopoDS_Shape are Null.");
3775     return NULL;
3776   }
3777
3778   Handle(GEOM_Function) aWhereFunction = theShapeWhere->GetLastFunction();
3779   if (aWhereFunction.IsNull()) {
3780     SetErrorCode("Error: aWhereFunction is Null.");
3781     return NULL;
3782   }
3783
3784   TopTools_IndexedMapOfShape aWhereIndices;
3785   TopExp::MapShapes(aWhere, aWhereIndices);
3786
3787   TopAbs_ShapeEnum iType = TopAbs_SOLID;
3788   Standard_Real    dl_l = 1e-3;
3789   Standard_Real    min_l, Tol_0D, Tol_1D, Tol_2D, Tol_3D, Tol_Mass;
3790   Standard_Real    aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
3791   Bnd_Box          BoundingBox;
3792   gp_Pnt           aPnt, aPnt_aWhat, tab_Pnt[2];
3793   GProp_GProps     aProps;
3794
3795   // Find the iType of the aWhat shape
3796   iType = GetTypeOfSimplePart(aWhat);
3797   if (iType == TopAbs_SHAPE) {
3798     SetErrorCode("Error: An attempt to extract a shape of not supported type.");
3799     return NULL;
3800   }
3801
3802   TopExp_Explorer Exp_aWhat  ( aWhat,  iType );
3803   TopExp_Explorer Exp_aWhere ( aWhere, iType );
3804   TopExp_Explorer Exp_Edge   ( aWhere, TopAbs_EDGE );
3805
3806   // Find the shortest edge in theShapeWhere shape
3807   BRepBndLib::Add(aWhere, BoundingBox);
3808   BoundingBox.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
3809   min_l = fabs(aXmax - aXmin);
3810   if( min_l < fabs(aYmax - aYmin) ) min_l = fabs(aYmax - aYmin);
3811   if( min_l < fabs(aZmax - aZmin) ) min_l = fabs(aZmax - aZmin);
3812   min_l /= dl_l;
3813   // Mantis issue 0020908 BEGIN
3814   if (!Exp_Edge.More()) {
3815     min_l = Precision::Confusion();
3816   }
3817   // Mantis issue 0020908 END
3818   for ( Standard_Integer nbEdge = 0; Exp_Edge.More(); Exp_Edge.Next(), nbEdge++ ) {
3819     TopExp_Explorer Exp_Vertex( Exp_Edge.Current(), TopAbs_VERTEX);
3820     for ( Standard_Integer nbVertex = 0; Exp_Vertex.More(); Exp_Vertex.Next(), nbVertex++ ) {
3821       aPnt = BRep_Tool::Pnt( TopoDS::Vertex( Exp_Vertex.Current() ) );
3822       tab_Pnt[nbVertex] = aPnt;
3823     }
3824     if ( ! tab_Pnt[0].IsEqual(tab_Pnt[1], dl_l) ) {
3825       BRepGProp::LinearProperties(Exp_Edge.Current(), aProps);
3826       if ( aProps.Mass() < min_l ) min_l = aProps.Mass();
3827     }
3828   }
3829
3830   // Compute tolerances
3831   Tol_0D = dl_l;
3832   Tol_1D = dl_l * min_l;
3833   Tol_2D = dl_l * ( min_l * min_l) * ( 2. + dl_l);
3834   Tol_3D = dl_l * ( min_l * min_l * min_l ) * ( 3. + (3 * dl_l) + (dl_l * dl_l) );
3835
3836   if (Tol_0D < Precision::Confusion()) Tol_0D = Precision::Confusion();
3837   if (Tol_1D < Precision::Confusion()) Tol_1D = Precision::Confusion();
3838   if (Tol_2D < Precision::Confusion()) Tol_2D = Precision::Confusion();
3839   if (Tol_3D < Precision::Confusion()) Tol_3D = Precision::Confusion();
3840
3841   Tol_Mass = Tol_3D;
3842   if ( iType == TopAbs_VERTEX )    Tol_Mass = Tol_0D;
3843   else if ( iType == TopAbs_EDGE ) Tol_Mass = Tol_1D;
3844   else if ( iType == TopAbs_FACE ) Tol_Mass = Tol_2D;
3845
3846   // Searching for the sub-shapes inside the ShapeWhere shape
3847   GEOMAlgo_GetInPlace aGIP;
3848   aGIP.SetTolerance(Tol_1D);
3849   aGIP.SetTolMass(Tol_Mass);
3850   aGIP.SetTolCG(Tol_1D);
3851
3852   aGIP.SetArgument(aWhat);
3853   aGIP.SetShapeWhere(aWhere);
3854
3855   aGIP.Perform();
3856   int iErr = aGIP.ErrorStatus();
3857   if (iErr) {
3858     SetErrorCode("Error in GEOMAlgo_GetInPlace");
3859     return NULL;
3860   }
3861
3862   // aGIP.IsFound() returns true only when the whole theShapeWhat
3863   // is found (as one shape or several parts). But we are also interested
3864   // in the partial result, that is why this check is commented.
3865   //if (!aGIP.IsFound()) {
3866   //  SetErrorCode(NOT_FOUND_ANY);
3867   //  return NULL;
3868   //}
3869
3870   const TopTools_DataMapOfShapeListOfShape& aDMSLS = aGIP.Images();
3871   if (!aDMSLS.IsBound(aWhat)) {
3872     SetErrorCode(NOT_FOUND_ANY);
3873     return NULL;
3874   }
3875
3876   // the list of shapes aLSA contains the shapes
3877   // of the Shape For Search that corresponds
3878   // to the Argument aWhat
3879   const TopTools_ListOfShape& aLSA = aDMSLS.Find(aWhat);
3880   if (aLSA.Extent() == 0) {
3881     SetErrorCode(NOT_FOUND_ANY); // Not found any Results
3882     return NULL;
3883   }
3884
3885   Handle(TColStd_HArray1OfInteger) aModifiedArray = new TColStd_HArray1OfInteger (1, aLSA.Extent());
3886   TopTools_ListIteratorOfListOfShape anIterModif (aLSA);
3887   for (Standard_Integer imod = 1; anIterModif.More(); anIterModif.Next(), imod++) {
3888     if (aWhereIndices.Contains(anIterModif.Value())) {
3889       aModifiedArray->SetValue(imod, aWhereIndices.FindIndex(anIterModif.Value()));
3890     }
3891     else {
3892       SetErrorCode("Error: wrong sub shape returned");
3893       return NULL;
3894     }
3895   }
3896
3897   //Add a new object
3898   Handle(GEOM_Object) aResult = GetEngine()->AddSubShape(theShapeWhere, aModifiedArray);
3899   if (aResult.IsNull()) {
3900     SetErrorCode("Error in algorithm: result found, but cannot be returned.");
3901     return NULL;
3902   }
3903
3904   if (aModifiedArray->Length() > 1) {
3905     //Set a GROUP type
3906     aResult->SetType(GEOM_GROUP);
3907
3908     //Set a sub shape type
3909     TopoDS_Shape aFirstFound = aLSA.First();
3910     TopAbs_ShapeEnum aShapeType = aFirstFound.ShapeType();
3911
3912     TDF_Label aFreeLabel = aResult->GetFreeLabel();
3913     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)aShapeType);
3914   }
3915
3916   //Make a Python command
3917   Handle(GEOM_Function) aFunction = aResult->GetFunction(1);
3918
3919   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetInPlace("
3920     << theShapeWhere << ", " << theShapeWhat << ", True)";
3921
3922   SetErrorCode(OK);
3923   return aResult;
3924 }
3925
3926 //=============================================================================
3927 /*!
3928  *  case GetInPlaceOld:
3929  *  default:
3930  */
3931 //=============================================================================
3932 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetInPlaceOld (Handle(GEOM_Object) theShapeWhere,
3933                                                                Handle(GEOM_Object) theShapeWhat)
3934 {
3935   SetErrorCode(KO);
3936
3937   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
3938
3939   TopoDS_Shape aWhere = theShapeWhere->GetValue();
3940   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
3941   TopoDS_Shape aPntShape;
3942   TopoDS_Vertex aVertex;
3943
3944   if (aWhere.IsNull() || aWhat.IsNull()) {
3945     SetErrorCode("Error: aWhere and aWhat TopoDS_Shape are Null.");
3946     return NULL;
3947   }
3948
3949   Handle(GEOM_Function) aWhereFunction = theShapeWhere->GetLastFunction();
3950   if (aWhereFunction.IsNull()) {
3951     SetErrorCode("Error: aWhereFunction is Null.");
3952     return NULL;
3953   }
3954
3955   TopTools_IndexedMapOfShape aWhereIndices;
3956   TopExp::MapShapes(aWhere, aWhereIndices);
3957
3958   TColStd_ListOfInteger aModifiedList;
3959   Standard_Integer aWhereIndex;
3960   Handle(TColStd_HArray1OfInteger) aModifiedArray;
3961   Handle(GEOM_Object) aResult;
3962
3963   bool isFound = false;
3964   TopAbs_ShapeEnum iType = TopAbs_SOLID;
3965   //Standard_Real    aWhat_Mass = 0., aWhere_Mass = 0.;
3966   Standard_Real    tab_aWhat[4],    tab_aWhere[4];
3967   Standard_Real    dl_l = 1e-3;
3968   Standard_Real    min_l, Tol_0D, Tol_1D, Tol_2D, Tol_3D, Tol_Mass;
3969   Standard_Real    aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
3970   Bnd_Box          BoundingBox;
3971   gp_Pnt           aPnt, aPnt_aWhat, tab_Pnt[2];
3972   GProp_GProps     aProps;
3973
3974   // Find the iType of the aWhat shape
3975   /*
3976   if      ( aWhat.ShapeType() == TopAbs_VERTEX )                                         iType = TopAbs_VERTEX;
3977   else if ( aWhat.ShapeType() == TopAbs_EDGE  || aWhat.ShapeType() == TopAbs_WIRE )      iType = TopAbs_EDGE;
3978   else if ( aWhat.ShapeType() == TopAbs_FACE  || aWhat.ShapeType() == TopAbs_SHELL )     iType = TopAbs_FACE;
3979   else if ( aWhat.ShapeType() == TopAbs_SOLID || aWhat.ShapeType() == TopAbs_COMPSOLID ) iType = TopAbs_SOLID;
3980   else if ( aWhat.ShapeType() == TopAbs_COMPOUND ) {
3981     // Only the iType of the first shape in the compound is taken into account
3982     TopoDS_Iterator It (aWhat, Standard_False, Standard_False);
3983     if ( !It.More() ) {
3984       SetErrorCode("Error: theShapeWhat is an empty COMPOUND.");
3985       return NULL;
3986     }
3987     TopAbs_ShapeEnum compType = It.Value().ShapeType();
3988     if      ( compType == TopAbs_VERTEX )                               iType = TopAbs_VERTEX;
3989     else if ( compType == TopAbs_EDGE  || compType == TopAbs_WIRE )     iType = TopAbs_EDGE;
3990     else if ( compType == TopAbs_FACE  || compType == TopAbs_SHELL)     iType = TopAbs_FACE;
3991     else if ( compType == TopAbs_SOLID || compType == TopAbs_COMPSOLID) iType = TopAbs_SOLID;
3992   }
3993   else {
3994     SetErrorCode("Error: An attempt to extract a shape of not supported type.");
3995     return NULL;
3996   }
3997   */
3998   iType = GetTypeOfSimplePart(aWhat);
3999   if (iType == TopAbs_SHAPE) {
4000     SetErrorCode("Error: An attempt to extract a shape of not supported type.");
4001     return NULL;
4002   }
4003
4004   TopExp_Explorer Exp_aWhat  ( aWhat,  iType );
4005   TopExp_Explorer Exp_aWhere ( aWhere, iType );
4006   TopExp_Explorer Exp_Edge   ( aWhere, TopAbs_EDGE );
4007
4008   // Find the shortest edge in theShapeWhere shape
4009   BRepBndLib::Add(aWhere, BoundingBox);
4010   BoundingBox.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
4011   min_l = fabs(aXmax - aXmin);
4012   if( min_l < fabs(aYmax - aYmin) ) min_l = fabs(aYmax - aYmin);
4013   if( min_l < fabs(aZmax - aZmin) ) min_l = fabs(aZmax - aZmin);
4014   min_l /= dl_l;
4015   // Mantis issue 0020908 BEGIN
4016   if (!Exp_Edge.More()) {
4017     min_l = Precision::Confusion();
4018   }
4019   // Mantis issue 0020908 END
4020   for ( Standard_Integer nbEdge = 0; Exp_Edge.More(); Exp_Edge.Next(), nbEdge++ ) {
4021     TopExp_Explorer Exp_Vertex( Exp_Edge.Current(), TopAbs_VERTEX);
4022     for ( Standard_Integer nbVertex = 0; Exp_Vertex.More(); Exp_Vertex.Next(), nbVertex++ ) {
4023       aPnt = BRep_Tool::Pnt( TopoDS::Vertex( Exp_Vertex.Current() ) );
4024       tab_Pnt[nbVertex] = aPnt;
4025     }
4026     if ( ! tab_Pnt[0].IsEqual(tab_Pnt[1], dl_l) ) {
4027       BRepGProp::LinearProperties(Exp_Edge.Current(), aProps);
4028       if ( aProps.Mass() < min_l ) min_l = aProps.Mass();
4029     }
4030   }
4031
4032   // Compute tolerances
4033   Tol_0D = dl_l;
4034   Tol_1D = dl_l * min_l;
4035   Tol_2D = dl_l * ( min_l * min_l) * ( 2. + dl_l);
4036   Tol_3D = dl_l * ( min_l * min_l * min_l ) * ( 3. + (3 * dl_l) + (dl_l * dl_l) );
4037
4038   if (Tol_0D < Precision::Confusion()) Tol_0D = Precision::Confusion();
4039   if (Tol_1D < Precision::Confusion()) Tol_1D = Precision::Confusion();
4040   if (Tol_2D < Precision::Confusion()) Tol_2D = Precision::Confusion();
4041   if (Tol_3D < Precision::Confusion()) Tol_3D = Precision::Confusion();
4042
4043   //if (Tol_1D > 1.0) Tol_1D = 1.0;
4044   //if (Tol_2D > 1.0) Tol_2D = 1.0;
4045   //if (Tol_3D > 1.0) Tol_3D = 1.0;
4046
4047   Tol_Mass = Tol_3D;
4048   if ( iType == TopAbs_VERTEX )    Tol_Mass = Tol_0D;
4049   else if ( iType == TopAbs_EDGE ) Tol_Mass = Tol_1D;
4050   else if ( iType == TopAbs_FACE ) Tol_Mass = Tol_2D;
4051
4052   // Compute the ShapeWhat Mass
4053   /*
4054   for ( ; Exp_aWhat.More(); Exp_aWhat.Next() ) {
4055     if ( iType == TopAbs_VERTEX ) {
4056       aWhat_Mass += 1;
4057       continue;
4058     }
4059     else if ( iType == TopAbs_EDGE ) BRepGProp::LinearProperties(Exp_aWhat.Current(),  aProps);
4060     else if ( iType == TopAbs_FACE ) BRepGProp::SurfaceProperties(Exp_aWhat.Current(), aProps);
4061     else                             BRepGProp::VolumeProperties(Exp_aWhat.Current(),  aProps);
4062     aWhat_Mass += aProps.Mass();
4063   }
4064   */
4065
4066   // Searching for the sub-shapes inside the ShapeWhere shape
4067   TopTools_MapOfShape map_aWhere;
4068   for ( Exp_aWhere.ReInit(); Exp_aWhere.More(); Exp_aWhere.Next() ) {
4069     if (!map_aWhere.Add(Exp_aWhere.Current()))
4070       continue; // skip repeated shape to avoid mass addition
4071     GetShapeProperties( Exp_aWhere.Current(), tab_aWhere, aPnt );
4072     for ( Exp_aWhat.ReInit(); Exp_aWhat.More(); Exp_aWhat.Next() ) {
4073       GetShapeProperties( Exp_aWhat.Current(), tab_aWhat, aPnt_aWhat );
4074       if ( fabs(tab_aWhat[3] - tab_aWhere[3]) <= Tol_Mass && aPnt_aWhat.Distance(aPnt) <= Tol_1D )
4075         isFound = true;
4076       else {
4077         if ( (tab_aWhat[3] - tab_aWhere[3]) > Tol_Mass ) {
4078           aPntShape = BRepBuilderAPI_MakeVertex( aPnt ).Shape();
4079           aVertex   = TopoDS::Vertex( aPntShape );
4080           BRepExtrema_DistShapeShape aWhereDistance ( aVertex, Exp_aWhere.Current() );
4081           BRepExtrema_DistShapeShape aWhatDistance  ( aVertex, Exp_aWhat.Current() );
4082           if ( aWhereDistance.IsDone() && aWhatDistance.IsDone() &&
4083                fabs(aWhereDistance.Value() - aWhatDistance.Value()) <= Tol_1D )
4084           {
4085             // 0020162: "EDF 961 GEOM : Getinplace is getting additionnal orthogonal faces"
4086             // aVertex must be projected to the same point on Where and on What
4087             gp_Pnt pOnWhat  = aWhatDistance.PointOnShape2(1);
4088             gp_Pnt pOnWhere = aWhereDistance.PointOnShape2(1);
4089             isFound = ( pOnWhat.Distance(pOnWhere) <= Tol_1D );
4090             if ( isFound && iType == TopAbs_FACE )
4091             {
4092               // check normals at pOnWhat and pOnWhere
4093               const double angleTol = PI/180.;
4094               gp_Vec normToWhat  = GetNormal( TopoDS::Face(Exp_aWhat.Current()), aWhatDistance);
4095               gp_Vec normToWhere = GetNormal( TopoDS::Face(Exp_aWhere.Current()), aWhereDistance);
4096               if ( normToWhat * normToWhere < 0 )
4097                 normToWhat.Reverse();
4098               isFound = ( normToWhat.Angle( normToWhere ) < angleTol );
4099             }
4100           }
4101         }
4102       }
4103       if ( isFound ) {
4104         aWhereIndex = aWhereIndices.FindIndex(Exp_aWhere.Current());
4105         aModifiedList.Append(aWhereIndex);
4106         //aWhere_Mass += tab_aWhere[3];
4107         isFound = false;
4108         break;
4109       }
4110     }
4111     //if ( fabs( aWhat_Mass - aWhere_Mass ) <= Tol_Mass )
4112       //break;
4113   }
4114
4115   if (aModifiedList.Extent() == 0) { // Not found any Results
4116     SetErrorCode(NOT_FOUND_ANY);
4117     return NULL;
4118   }
4119
4120   aModifiedArray = new TColStd_HArray1OfInteger (1, aModifiedList.Extent());
4121   TColStd_ListIteratorOfListOfInteger anIterModif (aModifiedList);
4122   for (Standard_Integer imod = 1; anIterModif.More(); anIterModif.Next(), imod++)
4123     aModifiedArray->SetValue(imod, anIterModif.Value());
4124
4125   //Add a new object
4126   aResult = GetEngine()->AddSubShape(theShapeWhere, aModifiedArray);
4127   if (aResult.IsNull()) {
4128     SetErrorCode("Error in algorithm: result found, but cannot be returned.");
4129     return NULL;
4130   }
4131
4132   if (aModifiedArray->Length() > 1) {
4133     //Set a GROUP type
4134     aResult->SetType(GEOM_GROUP);
4135
4136     //Set a sub shape type
4137     TopoDS_Shape aFirstFound = aWhereIndices.FindKey(aModifiedArray->Value(1));
4138     TopAbs_ShapeEnum aShapeType = aFirstFound.ShapeType();
4139
4140     TDF_Label aFreeLabel = aResult->GetFreeLabel();
4141     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)aShapeType);
4142   }
4143
4144   //Make a Python command
4145   Handle(GEOM_Function) aFunction = aResult->GetFunction(1);
4146
4147   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetInPlace("
4148     << theShapeWhere << ", " << theShapeWhat << ", False)";
4149
4150   SetErrorCode(OK);
4151   return aResult;
4152 }
4153
4154 //=======================================================================
4155 //function : GetInPlaceByHistory
4156 //purpose  :
4157 //=======================================================================
4158 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetInPlaceByHistory
4159                                           (Handle(GEOM_Object) theShapeWhere,
4160                                            Handle(GEOM_Object) theShapeWhat)
4161 {
4162   SetErrorCode(KO);
4163
4164   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
4165
4166   TopoDS_Shape aWhere = theShapeWhere->GetValue();
4167   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
4168
4169   if (aWhere.IsNull() || aWhat.IsNull()) return NULL;
4170
4171   Handle(GEOM_Function) aWhereFunction = theShapeWhere->GetLastFunction();
4172   if (aWhereFunction.IsNull()) return NULL;
4173
4174   //Fill array of indices
4175   TopTools_IndexedMapOfShape aWhereIndices;
4176   TopExp::MapShapes(aWhere, aWhereIndices);
4177
4178   // process shape
4179   TColStd_ListOfInteger aModifiedList;
4180   bool isFound = GetInPlaceOfShape(aWhereFunction, aWhereIndices, aWhat, aModifiedList);
4181
4182   if (!isFound || aModifiedList.Extent() < 1) {
4183     SetErrorCode("Error: No history found for the sought shape or its sub-shapes.");
4184     return NULL;
4185   }
4186
4187   Handle(TColStd_HArray1OfInteger) aModifiedArray =
4188     new TColStd_HArray1OfInteger (1, aModifiedList.Extent());
4189   TColStd_ListIteratorOfListOfInteger anIterModif (aModifiedList);
4190   for (Standard_Integer imod = 1; anIterModif.More(); anIterModif.Next(), imod++) {
4191     aModifiedArray->SetValue(imod, anIterModif.Value());
4192   }
4193
4194   //Add a new object
4195   Handle(GEOM_Object) aResult = GetEngine()->AddSubShape(theShapeWhere, aModifiedArray);
4196   if (aResult.IsNull()) {
4197     SetErrorCode("Error in algorithm: result found, but cannot be returned.");
4198     return NULL;
4199   }
4200
4201   if (aModifiedArray->Length() > 1) {
4202     //Set a GROUP type
4203     aResult->SetType(GEOM_GROUP);
4204
4205     //Set a sub shape type
4206     TopoDS_Shape aFirstFound = aWhereIndices.FindKey(aModifiedArray->Value(1));
4207     TopAbs_ShapeEnum aShapeType = aFirstFound.ShapeType();
4208
4209     TDF_Label aFreeLabel = aResult->GetFreeLabel();
4210     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)aShapeType);
4211   }
4212
4213   //Make a Python command
4214   Handle(GEOM_Function) aFunction = aResult->GetFunction(1);
4215
4216   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetInPlaceByHistory("
4217     << theShapeWhere << ", " << theShapeWhat << ")";
4218
4219   SetErrorCode(OK);
4220   return aResult;
4221 }
4222
4223 //=======================================================================
4224 //function : ShapeToDouble
4225 //purpose  : used by CompareShapes::operator()
4226 //=======================================================================
4227 std::pair<double, double> ShapeToDouble (const TopoDS_Shape& S, bool isOldSorting)
4228 {
4229   // Computing of CentreOfMass
4230   gp_Pnt GPoint;
4231   double Len;
4232
4233   if (S.ShapeType() == TopAbs_VERTEX) {
4234     GPoint = BRep_Tool::Pnt(TopoDS::Vertex(S));
4235     Len = (double)S.Orientation();
4236   }
4237   else {
4238     GProp_GProps GPr;
4239     // BEGIN: fix for Mantis issue 0020842
4240     if (isOldSorting) {
4241       BRepGProp::LinearProperties(S, GPr);
4242     }
4243     else {
4244       if (S.ShapeType() == TopAbs_EDGE || S.ShapeType() == TopAbs_WIRE) {
4245         BRepGProp::LinearProperties(S, GPr);
4246       }
4247       else if (S.ShapeType() == TopAbs_FACE || S.ShapeType() == TopAbs_SHELL) {
4248         BRepGProp::SurfaceProperties(S, GPr);
4249       }
4250       else {
4251         BRepGProp::VolumeProperties(S, GPr);
4252       }
4253     }
4254     // END: fix for Mantis issue 0020842
4255     GPoint = GPr.CentreOfMass();
4256     Len = GPr.Mass();
4257   }
4258
4259   double dMidXYZ = GPoint.X() * 999.0 + GPoint.Y() * 99.0 + GPoint.Z() * 0.9;
4260   return std::make_pair(dMidXYZ, Len);
4261 }
4262
4263 //=======================================================================
4264 //function : CompareShapes::operator()
4265 //purpose  : used by std::sort(), called from SortShapes()
4266 //=======================================================================
4267 bool GEOMImpl_IShapesOperations::CompareShapes::operator()(const TopoDS_Shape& theShape1,
4268                                                            const TopoDS_Shape& theShape2)
4269 {
4270   if (!myMap.IsBound(theShape1)) {
4271     myMap.Bind(theShape1, ShapeToDouble(theShape1, myIsOldSorting));
4272   }
4273
4274   if (!myMap.IsBound(theShape2)) {
4275     myMap.Bind(theShape2, ShapeToDouble(theShape2, myIsOldSorting));
4276   }
4277
4278   std::pair<double, double> val1 = myMap.Find(theShape1);
4279   std::pair<double, double> val2 = myMap.Find(theShape2);
4280
4281   double tol = Precision::Confusion();
4282   bool exchange = Standard_False;
4283
4284   double dMidXYZ = val1.first - val2.first;
4285   if (dMidXYZ >= tol) {
4286     exchange = Standard_True;
4287   }
4288   else if (Abs(dMidXYZ) < tol) {
4289     double dLength = val1.second - val2.second;
4290     if (dLength >= tol) {
4291       exchange = Standard_True;
4292     }
4293     else if (Abs(dLength) < tol && theShape1.ShapeType() <= TopAbs_FACE) {
4294       // PAL17233
4295       // equal values possible on shapes such as two halves of a sphere and
4296       // a membrane inside the sphere
4297       Bnd_Box box1,box2;
4298       BRepBndLib::Add(theShape1, box1);
4299       if (!box1.IsVoid()) {
4300         BRepBndLib::Add(theShape2, box2);
4301         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
4302         if (dSquareExtent >= tol) {
4303           exchange = Standard_True;
4304         }
4305         else if (Abs(dSquareExtent) < tol) {
4306           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
4307           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
4308           val1 = (aXmin+aXmax)*999.0 + (aYmin+aYmax)*99.0 + (aZmin+aZmax)*0.9;
4309           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
4310           val2 = (aXmin+aXmax)*999.0 + (aYmin+aYmax)*99.0 + (aZmin+aZmax)*0.9;
4311           if ((val1 - val2) >= tol) {
4312             exchange = Standard_True;
4313           }
4314         }
4315       }
4316     }
4317   }
4318
4319   //return val1 < val2;
4320   return !exchange;
4321 }
4322
4323 //=======================================================================
4324 //function : SortShapes
4325 //purpose  :
4326 //=======================================================================
4327 void GEOMImpl_IShapesOperations::SortShapes(TopTools_ListOfShape& SL,
4328                                             const Standard_Boolean isOldSorting)
4329 {
4330 #ifdef STD_SORT_ALGO
4331   std::vector<TopoDS_Shape> aShapesVec;
4332   aShapesVec.reserve(SL.Extent());
4333
4334   TopTools_ListIteratorOfListOfShape it (SL);
4335   for (; it.More(); it.Next()) {
4336     aShapesVec.push_back(it.Value());
4337   }
4338   SL.Clear();
4339
4340   CompareShapes shComp (isOldSorting);
4341   std::stable_sort(aShapesVec.begin(), aShapesVec.end(), shComp);
4342   //std::sort(aShapesVec.begin(), aShapesVec.end(), shComp);
4343
4344   std::vector<TopoDS_Shape>::const_iterator anIter = aShapesVec.begin();
4345   for (; anIter != aShapesVec.end(); ++anIter) {
4346     SL.Append(*anIter);
4347   }
4348 #else
4349   // old implementation
4350   Standard_Integer MaxShapes = SL.Extent();
4351   TopTools_Array1OfShape  aShapes (1,MaxShapes);
4352   TColStd_Array1OfInteger OrderInd(1,MaxShapes);
4353   TColStd_Array1OfReal    MidXYZ  (1,MaxShapes); //X,Y,Z;
4354   TColStd_Array1OfReal    Length  (1,MaxShapes); //X,Y,Z;
4355
4356   // Computing of CentreOfMass
4357   Standard_Integer Index;
4358   GProp_GProps GPr;
4359   gp_Pnt GPoint;
4360   TopTools_ListIteratorOfListOfShape it(SL);
4361   for (Index=1;  it.More();  Index++)
4362   {
4363     TopoDS_Shape S = it.Value();
4364     SL.Remove( it ); // == it.Next()
4365     aShapes(Index) = S;
4366     OrderInd.SetValue (Index, Index);
4367     if (S.ShapeType() == TopAbs_VERTEX) {
4368       GPoint = BRep_Tool::Pnt( TopoDS::Vertex( S ));
4369       Length.SetValue( Index, (Standard_Real) S.Orientation());
4370     }
4371     else {
4372       // BEGIN: fix for Mantis issue 0020842
4373       if (isOldSorting) {
4374         BRepGProp::LinearProperties (S, GPr);
4375       }
4376       else {
4377         if (S.ShapeType() == TopAbs_EDGE || S.ShapeType() == TopAbs_WIRE) {
4378           BRepGProp::LinearProperties (S, GPr);
4379         }
4380         else if (S.ShapeType() == TopAbs_FACE || S.ShapeType() == TopAbs_SHELL) {
4381           BRepGProp::SurfaceProperties(S, GPr);
4382         }
4383         else {
4384           BRepGProp::VolumeProperties(S, GPr);
4385         }
4386       }
4387       // END: fix for Mantis issue 0020842
4388       GPoint = GPr.CentreOfMass();
4389       Length.SetValue(Index, GPr.Mass());
4390     }
4391     MidXYZ.SetValue(Index, GPoint.X()*999.0 + GPoint.Y()*99.0 + GPoint.Z()*0.9);
4392     //cout << Index << " L: " << Length(Index) << "CG: " << MidXYZ(Index) << endl;
4393   }
4394
4395   // Sorting
4396   Standard_Integer aTemp;
4397   Standard_Boolean exchange, Sort = Standard_True;
4398   Standard_Real    tol = Precision::Confusion();
4399   while (Sort)
4400   {
4401     Sort = Standard_False;
4402     for (Index=1; Index < MaxShapes; Index++)
4403     {
4404       exchange = Standard_False;
4405       Standard_Real dMidXYZ = MidXYZ(OrderInd(Index)) - MidXYZ(OrderInd(Index+1));
4406       Standard_Real dLength = Length(OrderInd(Index)) - Length(OrderInd(Index+1));
4407       if ( dMidXYZ >= tol ) {
4408 //         cout << "MidXYZ: " << MidXYZ(OrderInd(Index))<< " > " <<MidXYZ(OrderInd(Index+1))
4409 //              << " d: " << dMidXYZ << endl;
4410         exchange = Standard_True;
4411       }
4412       else if ( Abs(dMidXYZ) < tol && dLength >= tol ) {
4413 //         cout << "Length: " << Length(OrderInd(Index))<< " > " <<Length(OrderInd(Index+1))
4414 //              << " d: " << dLength << endl;
4415         exchange = Standard_True;
4416       }
4417       else if ( Abs(dMidXYZ) < tol && Abs(dLength) < tol &&
4418                 aShapes(OrderInd(Index)).ShapeType() <= TopAbs_FACE) {
4419         // PAL17233
4420         // equal values possible on shapes such as two halves of a sphere and
4421         // a membrane inside the sphere
4422         Bnd_Box box1,box2;
4423         BRepBndLib::Add( aShapes( OrderInd(Index) ), box1 );
4424         if ( box1.IsVoid() ) continue;
4425         BRepBndLib::Add( aShapes( OrderInd(Index+1) ), box2 );
4426         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
4427         if ( dSquareExtent >= tol ) {
4428 //           cout << "SquareExtent: " << box1.SquareExtent()<<" > "<<box2.SquareExtent() << endl;
4429           exchange = Standard_True;
4430         }
4431         else if ( Abs(dSquareExtent) < tol ) {
4432           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
4433           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
4434           val1 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
4435           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
4436           val2 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
4437           //exchange = val1 > val2;
4438           if ((val1 - val2) >= tol) {
4439             exchange = Standard_True;
4440           }
4441           //cout << "box: " << val1<<" > "<<val2 << endl;
4442         }
4443       }
4444
4445       if (exchange)
4446       {
4447 //         cout << "exchange " << Index << " & " << Index+1 << endl;
4448         aTemp = OrderInd(Index);
4449         OrderInd(Index) = OrderInd(Index+1);
4450         OrderInd(Index+1) = aTemp;
4451         Sort = Standard_True;
4452       }
4453     }
4454   }
4455
4456   for (Index=1; Index <= MaxShapes; Index++)
4457     SL.Append( aShapes( OrderInd(Index) ));
4458 #endif
4459 }
4460
4461 //=======================================================================
4462 //function : CompsolidToCompound
4463 //purpose  :
4464 //=======================================================================
4465 TopoDS_Shape GEOMImpl_IShapesOperations::CompsolidToCompound (const TopoDS_Shape& theCompsolid)
4466 {
4467   if (theCompsolid.ShapeType() != TopAbs_COMPSOLID) {
4468     return theCompsolid;
4469   }
4470
4471   TopoDS_Compound aCompound;
4472   BRep_Builder B;
4473   B.MakeCompound(aCompound);
4474
4475   TopTools_MapOfShape mapShape;
4476   TopoDS_Iterator It (theCompsolid, Standard_True, Standard_True);
4477
4478   for (; It.More(); It.Next()) {
4479     TopoDS_Shape aShape_i = It.Value();
4480     if (mapShape.Add(aShape_i)) {
4481       B.Add(aCompound, aShape_i);
4482     }
4483   }
4484
4485   return aCompound;
4486 }
4487
4488 //=======================================================================
4489 //function : CheckTriangulation
4490 //purpose  :
4491 //=======================================================================
4492 bool GEOMImpl_IShapesOperations::CheckTriangulation (const TopoDS_Shape& aShape)
4493 {
4494   bool isTriangulation = true;
4495
4496   TopExp_Explorer exp (aShape, TopAbs_FACE);
4497   if (exp.More())
4498   {
4499     TopLoc_Location aTopLoc;
4500     Handle(Poly_Triangulation) aTRF;
4501     aTRF = BRep_Tool::Triangulation(TopoDS::Face(exp.Current()), aTopLoc);
4502     if (aTRF.IsNull()) {
4503       isTriangulation = false;
4504     }
4505   }
4506   else // no faces, try edges
4507   {
4508     TopExp_Explorer expe (aShape, TopAbs_EDGE);
4509     if (!expe.More()) {
4510       return false;
4511     }
4512     TopLoc_Location aLoc;
4513     Handle(Poly_Polygon3D) aPE = BRep_Tool::Polygon3D(TopoDS::Edge(expe.Current()), aLoc);
4514     if (aPE.IsNull()) {
4515       isTriangulation = false;
4516     }
4517   }
4518
4519   if (!isTriangulation) {
4520     // calculate deflection
4521     Standard_Real aDeviationCoefficient = 0.001;
4522
4523     Bnd_Box B;
4524     BRepBndLib::Add(aShape, B);
4525     Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
4526     B.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
4527
4528     Standard_Real dx = aXmax - aXmin, dy = aYmax - aYmin, dz = aZmax - aZmin;
4529     Standard_Real aDeflection = Max(Max(dx, dy), dz) * aDeviationCoefficient * 4;
4530     Standard_Real aHLRAngle = 0.349066;
4531
4532     BRepMesh_IncrementalMesh Inc (aShape, aDeflection, Standard_False, aHLRAngle);
4533   }
4534
4535   return true;
4536 }
4537
4538 #define MAX_TOLERANCE 1.e-7
4539
4540 //=======================================================================
4541 //function : isSameEdge
4542 //purpose  : Returns True if two edges coincide
4543 //=======================================================================
4544 static bool isSameEdge(const TopoDS_Edge& theEdge1, const TopoDS_Edge& theEdge2)
4545 {
4546   TopoDS_Vertex V11, V12, V21, V22;
4547   TopExp::Vertices(theEdge1, V11, V12);
4548   TopExp::Vertices(theEdge2, V21, V22);
4549   gp_Pnt P11 = BRep_Tool::Pnt(V11);
4550   gp_Pnt P12 = BRep_Tool::Pnt(V12);
4551   gp_Pnt P21 = BRep_Tool::Pnt(V21);
4552   gp_Pnt P22 = BRep_Tool::Pnt(V22);
4553   bool coincide = false;
4554
4555   //Check that ends of edges coincide
4556   if(P11.Distance(P21) <= MAX_TOLERANCE) {
4557     if(P12.Distance(P22) <= MAX_TOLERANCE) coincide =  true;
4558   }
4559   else if(P11.Distance(P22) <= MAX_TOLERANCE) {
4560     if(P12.Distance(P21) <= MAX_TOLERANCE) coincide = true;
4561   }
4562
4563   if(!coincide) return false;
4564
4565   if (BRep_Tool::Degenerated(theEdge1))
4566     if (BRep_Tool::Degenerated(theEdge2)) return true;
4567     else return false;
4568   else
4569     if (BRep_Tool::Degenerated(theEdge2)) return false;
4570
4571   double U11, U12, U21, U22;
4572   Handle(Geom_Curve) C1 = BRep_Tool::Curve(theEdge1, U11, U12);
4573   Handle(Geom_Curve) C2 = BRep_Tool::Curve(theEdge2, U21, U22);
4574   if(C1->DynamicType() == C2->DynamicType()) return true;
4575
4576   //Check that both edges has the same geometry
4577   double range = U12-U11;
4578   double U = U11+ range/3.0;
4579   gp_Pnt P1 = C1->Value(U);     //Compute a point on one third of the edge's length
4580   U = U11+range*2.0/3.0;
4581   gp_Pnt P2 = C1->Value(U);     //Compute a point on two thirds of the edge's length
4582
4583   if(!GeomLib_Tool::Parameter(C2, P1, MAX_TOLERANCE, U) ||  U < U21 || U > U22)
4584     return false;
4585
4586   if(P1.Distance(C2->Value(U)) > MAX_TOLERANCE) return false;
4587
4588   if(!GeomLib_Tool::Parameter(C2, P2, MAX_TOLERANCE, U) || U < U21 || U > U22)
4589     return false;
4590
4591   if(P2.Distance(C2->Value(U)) > MAX_TOLERANCE) return false;
4592
4593   return true;
4594 }
4595
4596 #include <TopoDS_TShape.hxx>
4597 //=======================================================================
4598 //function : isSameFace
4599 //purpose  : Returns True if two faces coincide
4600 //=======================================================================
4601 static bool isSameFace(const TopoDS_Face& theFace1, const TopoDS_Face& theFace2)
4602 {
4603   TopExp_Explorer E(theFace1, TopAbs_EDGE);
4604   TopTools_ListOfShape LS1, LS2;
4605   for(; E.More(); E.Next()) LS1.Append(E.Current());
4606
4607   E.Init(theFace2, TopAbs_EDGE);
4608   for(; E.More(); E.Next()) LS2.Append(E.Current());
4609
4610   //Compare the number of edges in the faces
4611   if(LS1.Extent() != LS2.Extent()) return false;
4612
4613   double aMin = RealFirst(), aMax = RealLast();
4614   double xminB1=aMax, yminB1=aMax, zminB1=aMax, xminB2=aMax, yminB2=aMax, zminB2=aMax;
4615   double xmaxB1=aMin, ymaxB1=aMin, zmaxB1=aMin, xmaxB2=aMin, ymaxB2=aMin, zmaxB2=aMin;
4616
4617   for(E.Init(theFace1, TopAbs_VERTEX); E.More(); E.Next()) {
4618     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
4619     if(P.X() < xminB1) xminB1 = P.X();
4620     if(P.Y() < yminB1) yminB1 = P.Y();
4621     if(P.Z() < zminB1) zminB1 = P.Z();
4622     if(P.X() > xmaxB1) xmaxB1 = P.X();
4623     if(P.Y() > ymaxB1) ymaxB1 = P.Y();
4624     if(P.Z() > zmaxB1) zmaxB1 = P.Z();
4625   }
4626
4627   for(E.Init(theFace2, TopAbs_VERTEX); E.More(); E.Next()) {
4628     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
4629     if(P.X() < xminB2) xminB2 = P.X();
4630     if(P.Y() < yminB2) yminB2 = P.Y();
4631     if(P.Z() < zminB2) zminB2 = P.Z();
4632     if(P.X() > xmaxB2) xmaxB2 = P.X();
4633     if(P.Y() > ymaxB2) ymaxB2 = P.Y();
4634     if(P.Z() > zmaxB2) zmaxB2 = P.Z();
4635   }
4636
4637   //Compare the bounding boxes of both faces
4638   if(gp_Pnt(xminB1, yminB1, zminB1).Distance(gp_Pnt(xminB2, yminB2, zminB2)) > MAX_TOLERANCE)
4639     return false;
4640
4641   if(gp_Pnt(xmaxB1, ymaxB1, zmaxB1).Distance(gp_Pnt(xmaxB2, ymaxB2, zmaxB2)) > MAX_TOLERANCE)
4642     return false;
4643
4644   Handle(Geom_Surface) S1 = BRep_Tool::Surface(theFace1);
4645   Handle(Geom_Surface) S2 = BRep_Tool::Surface(theFace2);
4646
4647   //Check if there a coincidence of two surfaces at least in two points
4648   double U11, U12, V11, V12, U21, U22, V21, V22;
4649   BRepTools::UVBounds(theFace1, U11, U12, V11, V12);
4650   BRepTools::UVBounds(theFace2, U21, U22, V21, V22);
4651
4652   double rangeU = U12-U11;
4653   double rangeV = V12-V11;
4654   double U = U11 + rangeU/3.0;
4655   double V = V11 + rangeV/3.0;
4656   gp_Pnt P1 = S1->Value(U, V);
4657   U = U11+rangeU*2.0/3.0;
4658   V = V11+rangeV*2.0/3.0;
4659   gp_Pnt P2 = S1->Value(U, V);
4660
4661   if (!GeomLib_Tool::Parameters(S2, P1, MAX_TOLERANCE, U, V) || U < U21 || U > U22 || V < V21 || V > V22)
4662     return false;
4663
4664   if (P1.Distance(S2->Value(U,V)) > MAX_TOLERANCE) return false;
4665
4666   if (!GeomLib_Tool::Parameters(S2, P2, MAX_TOLERANCE, U, V) || U < U21 || U > U22 || V < V21 || V > V22)
4667     return false;
4668
4669   if (P2.Distance(S2->Value(U, V)) > MAX_TOLERANCE) return false;
4670
4671   //Check that each edge of the Face1 has a counterpart in the Face2
4672   TopTools_MapOfOrientedShape aMap;
4673   TopTools_ListIteratorOfListOfShape LSI1(LS1);
4674   for(; LSI1.More(); LSI1.Next()) {
4675     TopoDS_Edge E = TopoDS::Edge(LSI1.Value());
4676     bool isFound = false;
4677     TopTools_ListIteratorOfListOfShape LSI2(LS2);
4678     for(; LSI2.More(); LSI2.Next()) {
4679       TopoDS_Shape aValue = LSI2.Value();
4680       if(aMap.Contains(aValue)) continue; //To avoid checking already found edge several times
4681       if(isSameEdge(E, TopoDS::Edge(aValue))) {
4682         aMap.Add(aValue);
4683         isFound = true;
4684         break;
4685       }
4686     }
4687     if(!isFound) return false;
4688   }
4689
4690   return true;
4691 }
4692
4693 //=======================================================================
4694 //function : isSameSolid
4695 //purpose  : Returns True if two solids coincide
4696 //=======================================================================
4697 bool isSameSolid(const TopoDS_Solid& theSolid1, const TopoDS_Solid& theSolid2)
4698 {
4699   TopExp_Explorer E(theSolid1, TopAbs_FACE);
4700   TopTools_ListOfShape LS1, LS2;
4701   for(; E.More(); E.Next()) LS1.Append(E.Current());
4702   E.Init(theSolid2, TopAbs_FACE);
4703   for(; E.More(); E.Next()) LS2.Append(E.Current());
4704
4705   if(LS1.Extent() != LS2.Extent()) return false;
4706
4707   double aMin = RealFirst(), aMax = RealLast();
4708   double xminB1=aMax, yminB1=aMax, zminB1=aMax, xminB2=aMax, yminB2=aMax, zminB2=aMax;
4709   double xmaxB1=aMin, ymaxB1=aMin, zmaxB1=aMin, xmaxB2=aMin, ymaxB2=aMin, zmaxB2=aMin;
4710
4711   for(E.Init(theSolid1, TopAbs_VERTEX); E.More(); E.Next()) {
4712     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
4713     if(P.X() < xminB1) xminB1 = P.X();
4714     if(P.Y() < yminB1) yminB1 = P.Y();
4715     if(P.Z() < zminB1) zminB1 = P.Z();
4716     if(P.X() > xmaxB1) xmaxB1 = P.X();
4717     if(P.Y() > ymaxB1) ymaxB1 = P.Y();
4718     if(P.Z() > zmaxB1) zmaxB1 = P.Z();
4719   }
4720
4721   for(E.Init(theSolid2, TopAbs_VERTEX); E.More(); E.Next()) {
4722     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
4723     if(P.X() < xminB2) xminB2 = P.X();
4724     if(P.Y() < yminB2) yminB2 = P.Y();
4725     if(P.Z() < zminB2) zminB2 = P.Z();
4726     if(P.X() > xmaxB2) xmaxB2 = P.X();
4727     if(P.Y() > ymaxB2) ymaxB2 = P.Y();
4728     if(P.Z() > zmaxB2) zmaxB2 = P.Z();
4729   }
4730
4731   //Compare the bounding boxes of both solids
4732   if(gp_Pnt(xminB1, yminB1, zminB1).Distance(gp_Pnt(xminB2, yminB2, zminB2)) > MAX_TOLERANCE)
4733     return false;
4734
4735   if(gp_Pnt(xmaxB1, ymaxB1, zmaxB1).Distance(gp_Pnt(xmaxB2, ymaxB2, zmaxB2)) > MAX_TOLERANCE)
4736     return false;
4737
4738   //Check that each face of the Solid1 has a counterpart in the Solid2
4739   TopTools_MapOfOrientedShape aMap;
4740   TopTools_ListIteratorOfListOfShape LSI1(LS1);
4741   for(; LSI1.More(); LSI1.Next()) {
4742     TopoDS_Face F = TopoDS::Face(LSI1.Value());
4743     bool isFound = false;
4744     TopTools_ListIteratorOfListOfShape LSI2(LS2);
4745     for(; LSI2.More(); LSI2.Next()) {
4746       if(aMap.Contains(LSI2.Value())) continue; //To avoid checking already found faces several times
4747       if(isSameFace(F, TopoDS::Face(LSI2.Value()))) {
4748         aMap.Add(LSI2.Value());
4749         isFound = true;
4750         break;
4751       }
4752     }
4753     if(!isFound) return false;
4754   }
4755
4756   return true;
4757 }
4758
4759 //=======================================================================
4760 //function : GetSame
4761 //purpose  :
4762 //=======================================================================
4763 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetSame(const Handle(GEOM_Object)& theShapeWhere,
4764                                                         const Handle(GEOM_Object)& theShapeWhat)
4765 {
4766   SetErrorCode(KO);
4767   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
4768
4769   TopoDS_Shape aWhere = theShapeWhere->GetValue();
4770   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
4771
4772   if (aWhere.IsNull() || aWhat.IsNull()) return NULL;
4773
4774   int anIndex = -1;
4775   bool isFound = false;
4776   TopoDS_Shape aSubShape;
4777   TopTools_MapOfShape aMap;
4778
4779   if (aWhat.ShapeType() == TopAbs_COMPOUND || aWhat.ShapeType() == TopAbs_COMPSOLID) {
4780     TopoDS_Iterator It (aWhat, Standard_True, Standard_True);
4781     if (It.More()) aWhat = It.Value();
4782     It.Next();
4783     if (It.More()) {
4784       SetErrorCode("Compounds of two or more shapes are not allowed for aWhat argument");
4785       return NULL;
4786     }
4787   }
4788
4789   switch (aWhat.ShapeType()) {
4790     case TopAbs_VERTEX: {
4791       gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(aWhat));
4792       TopExp_Explorer E(aWhere, TopAbs_VERTEX);
4793       for(; E.More(); E.Next()) {
4794         if(!aMap.Add(E.Current())) continue;
4795         gp_Pnt P2 = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
4796         if(P.Distance(P2) <= MAX_TOLERANCE) {
4797           isFound = true;
4798           aSubShape = E.Current();
4799           break;
4800         }
4801       }
4802       break;
4803                         }
4804     case TopAbs_EDGE: {
4805       TopoDS_Edge anEdge = TopoDS::Edge(aWhat);
4806       TopExp_Explorer E(aWhere, TopAbs_EDGE);
4807       for(; E.More(); E.Next()) {
4808         if(!aMap.Add(E.Current())) continue;
4809         if(isSameEdge(anEdge, TopoDS::Edge(E.Current()))) {
4810           aSubShape = E.Current();
4811           isFound = true;
4812           break;
4813         }
4814       }
4815       break;
4816                       }
4817     case TopAbs_FACE: {
4818       TopoDS_Face aFace = TopoDS::Face(aWhat);
4819       TopExp_Explorer E(aWhere, TopAbs_FACE);
4820       for(; E.More(); E.Next()) {
4821         if(!aMap.Add(E.Current())) continue;
4822         if(isSameFace(aFace, TopoDS::Face(E.Current()))) {
4823           aSubShape = E.Current();
4824           isFound = true;
4825           break;
4826         }
4827       }
4828       break;
4829                       }
4830     case TopAbs_SOLID: {
4831       TopoDS_Solid aSolid = TopoDS::Solid(aWhat);
4832       TopExp_Explorer E(aWhere, TopAbs_SOLID);
4833       for(; E.More(); E.Next()) {
4834         if(!aMap.Add(E.Current())) continue;
4835         if(isSameSolid(aSolid, TopoDS::Solid(E.Current()))) {
4836           aSubShape = E.Current();
4837           isFound = true;
4838           break;
4839         }
4840       }
4841       break;
4842                        }
4843     default:
4844       return NULL;
4845   }
4846
4847   if (isFound) {
4848     TopTools_IndexedMapOfShape anIndices;
4849     TopExp::MapShapes(aWhere, anIndices);
4850     if (anIndices.Contains(aSubShape))
4851       anIndex = anIndices.FindIndex(aSubShape);
4852   }
4853
4854   if (anIndex < 0) return NULL;
4855
4856   Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
4857
4858   anArray->SetValue(1, anIndex);
4859
4860   Handle(GEOM_Object) aResult = GetEngine()->AddSubShape(theShapeWhere, anArray);
4861   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
4862
4863   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetSame("
4864     << theShapeWhere << ", " << theShapeWhat << ")";
4865
4866   SetErrorCode(OK);
4867
4868   return aResult;
4869 }