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