src/GEOMImpl/GEOMImpl_IShapesOperations.cxx

   1 // Copyright (C) 2005  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   2 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   3 //
   4 // This library is free software; you can redistribute it and/or
   5 // modify it under the terms of the GNU Lesser General Public
   6 // License as published by the Free Software Foundation; either
   7 // version 2.1 of the License.
   8 //
   9 // This library is distributed in the hope that it will be useful
  10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  12 // Lesser General Public License for more details.
  13 //
  14 // You should have received a copy of the GNU Lesser General Public
  15 // License along with this library; if not, write to the Free Software
  16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  17 //
  18 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  19 //
  20 // File      : GEOMImpl_IShapesOperations.cxx
  21 // Created   :
  22 // Author    : modified by Lioka RAZAFINDRAZAKA (CEA) 22/06/2007
  23 // Project   : SALOME
  24 // $Header$
  25
  26
  27 #include <Standard_Stream.hxx>
  28
  29 #include "GEOMImpl_IShapesOperations.hxx"
  30
  31 #include "GEOMImpl_Types.hxx"
  32
  33 #include "GEOMImpl_VectorDriver.hxx"
  34 #include "GEOMImpl_ShapeDriver.hxx"
  35 #include "GEOMImpl_CopyDriver.hxx"
  36 #include "GEOMImpl_GlueDriver.hxx"
  37
  38 #include "GEOMImpl_IVector.hxx"
  39 #include "GEOMImpl_IShapes.hxx"
  40 #include "GEOMImpl_IGlue.hxx"
  41
  42 #include "GEOMImpl_Block6Explorer.hxx"
  43
  44 #include "GEOM_Function.hxx"
  45 #include "GEOM_ISubShape.hxx"
  46 #include "GEOM_PythonDump.hxx"
  47
  48 #include "GEOMAlgo_FinderShapeOn1.hxx"
  49 #include "GEOMAlgo_FinderShapeOnQuad.hxx"
  50 #include "GEOMAlgo_FinderShapeOn2.hxx"
  51 #include "GEOMAlgo_ClsfBox.hxx"
  52 #include "GEOMAlgo_ClsfSolid.hxx"
  53 #include "GEOMAlgo_Gluer1.hxx"
  54 #include "GEOMAlgo_ListIteratorOfListOfCoupleOfShapes.hxx"
  55 #include "GEOMAlgo_CoupleOfShapes.hxx"
  56 #include "GEOMAlgo_ListOfCoupleOfShapes.hxx"
  57
  58 #include "utilities.h"
  59 #include "OpUtil.hxx"
  60 #include "Utils_ExceptHandlers.hxx"
  61
  62 #include <TFunction_DriverTable.hxx>
  63 #include <TFunction_Driver.hxx>
  64 #include <TFunction_Logbook.hxx>
  65 #include <TDataStd_Integer.hxx>
  66 #include <TDataStd_IntegerArray.hxx>
  67 #include <TDF_Tool.hxx>
  68
  69 #include <BRepExtrema_ExtCF.hxx>
  70 #include <BRepExtrema_DistShapeShape.hxx>
  71
  72 #include <BRep_Tool.hxx>
  73 #include <BRep_Builder.hxx>
  74 #include <BRepTools.hxx>
  75 #include <BRepGProp.hxx>
  76 #include <BRepAdaptor_Curve.hxx>
  77 #include <BRepAdaptor_Surface.hxx>
  78 #include <BRepBndLib.hxx>
  79 #include <BRepBuilderAPI_MakeFace.hxx>
  80 #include <BRepMesh_IncrementalMesh.hxx>
  81
  82 #include <TopAbs.hxx>
  83 #include <TopExp.hxx>
  84 #include <TopoDS.hxx>
  85 #include <TopoDS_Shape.hxx>
  86 #include <TopoDS_Solid.hxx>
  87 #include <TopoDS_Face.hxx>
  88 #include <TopoDS_Edge.hxx>
  89 #include <TopoDS_Vertex.hxx>
  90 #include <TopoDS_Compound.hxx>
  91 #include <TopoDS_Iterator.hxx>
  92 #include <TopExp_Explorer.hxx>
  93 #include <TopLoc_Location.hxx>
  94 #include <TopTools_MapOfShape.hxx>
  95 #include <TopTools_MapOfOrientedShape.hxx>
  96 #include <TopTools_Array1OfShape.hxx>
  97 #include <TopTools_ListIteratorOfListOfShape.hxx>
  98 #include <TopTools_IndexedMapOfShape.hxx>
  99
 100 #include <Geom_Surface.hxx>
 101 #include <Geom_Plane.hxx>
 102 #include <Geom_SphericalSurface.hxx>
 103 #include <Geom_CylindricalSurface.hxx>
 104 #include <GeomAdaptor_Surface.hxx>
 105
 106 #include <GeomLib_Tool.hxx>
 107 #include <Geom2d_Curve.hxx>
 108
 109 #include <Bnd_Box.hxx>
 110 #include <GProp_GProps.hxx>
 111 #include <gp_Pnt.hxx>
 112 #include <gp_Lin.hxx>
 113 #include <TColStd_ListOfInteger.hxx>
 114 #include <TColStd_ListIteratorOfListOfInteger.hxx>
 115 #include <TColStd_Array1OfReal.hxx>
 116 #include <TColStd_HArray1OfInteger.hxx>
 117
 118 #include <vector>
 119
 120 #include <Standard_NullObject.hxx>
 121 #include <Standard_Failure.hxx>
 122 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
 123
 124 // Includes added for GetInPlace algorithm improvement
 125
 126 #include <GEOMImpl_MeasureDriver.hxx>
 127 #include <GEOMImpl_IMeasure.hxx>
 128 #include <BRepBuilderAPI_MakeVertex.hxx>
 129
 130 #include <BRepClass_FaceClassifier.hxx>
 131 #include <BRepClass3d_SolidClassifier.hxx>
 132 #include <Precision.hxx>
 133
 134 using namespace std;
 135
 136
 137 //=============================================================================
 138 /*!
 139  *   constructor:
 140  */
 141 //=============================================================================
 142 GEOMImpl_IShapesOperations::GEOMImpl_IShapesOperations (GEOM_Engine* theEngine, int theDocID)
 143 : GEOM_IOperations(theEngine, theDocID)
 144 {
 145   MESSAGE("GEOMImpl_IShapesOperations::GEOMImpl_IShapesOperations");
 146 }
 147
 148 //=============================================================================
 149 /*!
 150  *  destructor
 151  */
 152 //=============================================================================
 153 GEOMImpl_IShapesOperations::~GEOMImpl_IShapesOperations()
 154 {
 155   MESSAGE("GEOMImpl_IShapesOperations::~GEOMImpl_IShapesOperations");
 156 }
 157
 158
 159 //=============================================================================
 160 /*!
 161  *  MakeEdge
 162  */
 163 //=============================================================================
 164 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeEdge
 165                      (Handle(GEOM_Object) thePnt1, Handle(GEOM_Object) thePnt2)
 166 {
 167   SetErrorCode(KO);
 168
 169   if (thePnt1.IsNull() || thePnt2.IsNull()) return NULL;
 170
 171   //Add a new Edge object
 172   Handle(GEOM_Object) anEdge = GetEngine()->AddObject(GetDocID(), GEOM_EDGE);
 173
 174   //Add a new Vector function
 175   Handle(GEOM_Function) aFunction =
 176     anEdge->AddFunction(GEOMImpl_VectorDriver::GetID(), VECTOR_TWO_PNT);
 177
 178   //Check if the function is set correctly
 179   if (aFunction->GetDriverGUID() != GEOMImpl_VectorDriver::GetID()) return NULL;
 180
 181   GEOMImpl_IVector aPI (aFunction);
 182
 183   Handle(GEOM_Function) aRef1 = thePnt1->GetLastFunction();
 184   Handle(GEOM_Function) aRef2 = thePnt2->GetLastFunction();
 185   if (aRef1.IsNull() || aRef2.IsNull()) return NULL;
 186
 187   aPI.SetPoint1(aRef1);
 188   aPI.SetPoint2(aRef2);
 189
 190   //Compute the Edge value
 191   try {
 192 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 193     OCC_CATCH_SIGNALS;
 194 #endif
 195     if (!GetSolver()->ComputeFunction(aFunction)) {
 196       SetErrorCode("Vector driver failed");
 197       return NULL;
 198     }
 199   }
 200   catch (Standard_Failure) {
 201     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 202     SetErrorCode(aFail->GetMessageString());
 203     return NULL;
 204   }
 205
 206   //Make a Python command
 207   GEOM::TPythonDump(aFunction) << anEdge << " = geompy.MakeEdge("
 208                                << thePnt1 << ", " << thePnt2 << ")";
 209
 210   SetErrorCode(OK);
 211   return anEdge;
 212 }
 213
 214 //=============================================================================
 215 /*!
 216  *  MakeWire
 217  */
 218 //=============================================================================
 219 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeWire
 220                              (list<Handle(GEOM_Object)> theShapes)
 221 {
 222   return MakeShape(theShapes, GEOM_WIRE, WIRE_EDGES, "MakeWire");
 223 }
 224
 225 //=============================================================================
 226 /*!
 227  *  MakeFace
 228  */
 229 //=============================================================================
 230 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeFace (Handle(GEOM_Object) theWire,
 231                                                           const bool isPlanarWanted)
 232 {
 233   SetErrorCode(KO);
 234
 235   if (theWire.IsNull()) return NULL;
 236
 237   //Add a new Face object
 238   Handle(GEOM_Object) aFace = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 239
 240   //Add a new Shape function for creation of a face from a wire
 241   Handle(GEOM_Function) aFunction =
 242     aFace->AddFunction(GEOMImpl_ShapeDriver::GetID(), FACE_WIRE);
 243   if (aFunction.IsNull()) return NULL;
 244
 245   //Check if the function is set correctly
 246   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 247
 248   GEOMImpl_IShapes aCI (aFunction);
 249
 250   Handle(GEOM_Function) aRefWire = theWire->GetLastFunction();
 251
 252   if (aRefWire.IsNull()) return NULL;
 253
 254   aCI.SetBase(aRefWire);
 255   aCI.SetIsPlanar(isPlanarWanted);
 256
 257   //Compute the Face value
 258   try {
 259 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 260     OCC_CATCH_SIGNALS;
 261 #endif
 262     if (!GetSolver()->ComputeFunction(aFunction)) {
 263       SetErrorCode("Shape driver failed to compute a face");
 264       return NULL;
 265     }
 266   }
 267   catch (Standard_Failure) {
 268     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 269     SetErrorCode(aFail->GetMessageString());
 270     return NULL;
 271   }
 272
 273   //Make a Python command
 274   GEOM::TPythonDump(aFunction) << aFace << " = geompy.MakeFace("
 275     << theWire << ", " << (int)isPlanarWanted << ")";
 276
 277   SetErrorCode(OK);
 278   return aFace;
 279 }
 280
 281 //=============================================================================
 282 /*!
 283  *  MakeFaceWires
 284  */
 285 //=============================================================================
 286 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeFaceWires
 287                              (list<Handle(GEOM_Object)> theShapes,
 288                               const bool isPlanarWanted)
 289 {
 290   SetErrorCode(KO);
 291
 292   //Add a new object
 293   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 294
 295   //Add a new function
 296   Handle(GEOM_Function) aFunction =
 297     aShape->AddFunction(GEOMImpl_ShapeDriver::GetID(), FACE_WIRES);
 298   if (aFunction.IsNull()) return NULL;
 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(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 306
 307   // Shapes
 308   list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 309   for (; it != theShapes.end(); it++) {
 310     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 311     if (aRefSh.IsNull()) {
 312       SetErrorCode("NULL argument shape for the face construction");
 313       return NULL;
 314     }
 315     aShapesSeq->Append(aRefSh);
 316   }
 317   aCI.SetShapes(aShapesSeq);
 318
 319   aCI.SetIsPlanar(isPlanarWanted);
 320
 321   //Compute the shape
 322   try {
 323 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 324     OCC_CATCH_SIGNALS;
 325 #endif
 326     if (!GetSolver()->ComputeFunction(aFunction)) {
 327       SetErrorCode("Shape driver failed");
 328       return NULL;
 329     }
 330   }
 331   catch (Standard_Failure) {
 332     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 333     SetErrorCode(aFail->GetMessageString());
 334     return NULL;
 335   }
 336
 337   //Make a Python command
 338   GEOM::TPythonDump pd (aFunction);
 339   pd << aShape << " = geompy.MakeFaceWires([";
 340
 341   // Shapes
 342   it = theShapes.begin();
 343   if (it != theShapes.end()) {
 344     pd << (*it++);
 345     while (it != theShapes.end()) {
 346       pd << ", " << (*it++);
 347     }
 348   }
 349   pd << "], " << (int)isPlanarWanted << ")";
 350
 351   SetErrorCode(OK);
 352   return aShape;
 353 }
 354
 355 //=============================================================================
 356 /*!
 357  *  MakeShell
 358  */
 359 //=============================================================================
 360 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeShell
 361                              (list<Handle(GEOM_Object)> theShapes)
 362 {
 363   return MakeShape(theShapes, GEOM_SHELL, SHELL_FACES, "MakeShell");
 364 }
 365
 366 //=============================================================================
 367 /*!
 368  *  MakeSolidShells
 369  */
 370 //=============================================================================
 371 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeSolidShells
 372                              (list<Handle(GEOM_Object)> theShapes)
 373 {
 374   return MakeShape(theShapes, GEOM_SOLID, SOLID_SHELLS, "MakeSolid");
 375 }
 376
 377 //=============================================================================
 378 /*!
 379  *  MakeSolidShell
 380  */
 381 //=============================================================================
 382 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeSolidShell (Handle(GEOM_Object) theShell)
 383 {
 384   SetErrorCode(KO);
 385
 386   if (theShell.IsNull()) return NULL;
 387
 388   //Add a new Solid object
 389   Handle(GEOM_Object) aSolid = GetEngine()->AddObject(GetDocID(), GEOM_SOLID);
 390
 391   //Add a new Solid function for creation of a solid from a shell
 392   Handle(GEOM_Function) aFunction =
 393     aSolid->AddFunction(GEOMImpl_ShapeDriver::GetID(), SOLID_SHELL);
 394   if (aFunction.IsNull()) return NULL;
 395
 396   //Check if the function is set correctly
 397   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 398
 399   GEOMImpl_IShapes aCI (aFunction);
 400
 401   Handle(GEOM_Function) aRefShell = theShell->GetLastFunction();
 402
 403   if (aRefShell.IsNull()) return NULL;
 404
 405   aCI.SetBase(aRefShell);
 406
 407   //Compute the Solid value
 408   try {
 409 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 410     OCC_CATCH_SIGNALS;
 411 #endif
 412     if (!GetSolver()->ComputeFunction(aFunction)) {
 413       SetErrorCode("Solid driver failed");
 414       return NULL;
 415     }
 416   }
 417   catch (Standard_Failure) {
 418     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 419     SetErrorCode(aFail->GetMessageString());
 420     return NULL;
 421   }
 422
 423   //Make a Python command
 424   GEOM::TPythonDump(aFunction) << aSolid
 425     << " = geompy.MakeSolid(" << theShell << ")";
 426
 427   SetErrorCode(OK);
 428   return aSolid;
 429 }
 430
 431 //=============================================================================
 432 /*!
 433  *  MakeCompound
 434  */
 435 //=============================================================================
 436 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeCompound
 437                              (list<Handle(GEOM_Object)> theShapes)
 438 {
 439   return MakeShape(theShapes, GEOM_COMPOUND, COMPOUND_SHAPES, "MakeCompound");
 440 }
 441
 442 //=============================================================================
 443 /*!
 444  *  MakeShape
 445  */
 446 //=============================================================================
 447 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeShape
 448                              (list<Handle(GEOM_Object)>      theShapes,
 449                               const Standard_Integer         theObjectType,
 450                               const Standard_Integer         theFunctionType,
 451                               const TCollection_AsciiString& theMethodName)
 452 {
 453   SetErrorCode(KO);
 454
 455   //Add a new object
 456   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), theObjectType);
 457
 458   //Add a new function
 459   Handle(GEOM_Function) aFunction =
 460     aShape->AddFunction(GEOMImpl_ShapeDriver::GetID(), theFunctionType);
 461   if (aFunction.IsNull()) return NULL;
 462
 463   //Check if the function is set correctly
 464   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 465
 466   GEOMImpl_IShapes aCI (aFunction);
 467
 468   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 469
 470   // Shapes
 471   list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 472   for (; it != theShapes.end(); it++) {
 473     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 474     if (aRefSh.IsNull()) {
 475       SetErrorCode("NULL argument shape for the shape construction");
 476       return NULL;
 477     }
 478     aShapesSeq->Append(aRefSh);
 479   }
 480   aCI.SetShapes(aShapesSeq);
 481
 482   //Compute the shape
 483   try {
 484 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 485     OCC_CATCH_SIGNALS;
 486 #endif
 487     if (!GetSolver()->ComputeFunction(aFunction)) {
 488       SetErrorCode("Shape driver failed");
 489       return NULL;
 490     }
 491   }
 492   catch (Standard_Failure) {
 493     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 494     SetErrorCode(aFail->GetMessageString());
 495     return NULL;
 496   }
 497
 498   //Make a Python command
 499   GEOM::TPythonDump pd (aFunction);
 500   pd << aShape << " = geompy." << theMethodName.ToCString() << "([";
 501
 502   // Shapes
 503   it = theShapes.begin();
 504   if (it != theShapes.end()) {
 505     pd << (*it++);
 506     while (it != theShapes.end()) {
 507       pd << ", " << (*it++);
 508     }
 509   }
 510   pd << "])";
 511
 512   SetErrorCode(OK);
 513   return aShape;
 514 }
 515
 516 //=============================================================================
 517 /*!
 518  *  MakeGlueFaces
 519  */
 520 //=============================================================================
 521 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueFaces
 522                                                 (Handle(GEOM_Object) theShape,
 523                                                  const Standard_Real theTolerance,
 524                                                  const Standard_Boolean doKeepNonSolids)
 525 {
 526   SetErrorCode(KO);
 527
 528   if (theShape.IsNull()) return NULL;
 529
 530   //Add a new Glued object
 531   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
 532
 533   //Add a new Glue function
 534   Handle(GEOM_Function) aFunction;
 535   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_FACES);
 536   if (aFunction.IsNull()) return NULL;
 537
 538   //Check if the function is set correctly
 539   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
 540
 541   GEOMImpl_IGlue aCI (aFunction);
 542
 543   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 544   if (aRefShape.IsNull()) return NULL;
 545
 546   aCI.SetBase(aRefShape);
 547   aCI.SetTolerance(theTolerance);
 548   aCI.SetKeepNonSolids(doKeepNonSolids);
 549
 550   //Compute the sub-shape value
 551   Standard_Boolean isWarning = Standard_False;
 552   try {
 553 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 554     OCC_CATCH_SIGNALS;
 555 #endif
 556     if (!GetSolver()->ComputeFunction(aFunction)) {
 557       SetErrorCode("Shape driver failed to glue faces");
 558       return NULL;
 559     }
 560   }
 561   catch (Standard_Failure) {
 562     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 563     SetErrorCode(aFail->GetMessageString());
 564     // to provide warning
 565     if (!aFunction->GetValue().IsNull()) {
 566       isWarning = Standard_True;
 567     } else {
 568       return NULL;
 569     }
 570   }
 571
 572   //Make a Python command
 573   GEOM::TPythonDump(aFunction) << aGlued << " = geompy.MakeGlueFaces("
 574     << theShape << ", " << theTolerance << ")";
 575
 576   // to provide warning
 577   if (!isWarning) SetErrorCode(OK);
 578   return aGlued;
 579 }
 580
 581 //=============================================================================
 582 /*!
 583  *  GetGlueFaces
 584  */
 585 //=============================================================================
 586 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetGlueFaces
 587                                                 (Handle(GEOM_Object) theShape,
 588                                                  const Standard_Real theTolerance)
 589 {
 590   SetErrorCode(KO);
 591
 592   if (theShape.IsNull()) return NULL;
 593   TopoDS_Shape aShape = theShape->GetValue();
 594   if (aShape.IsNull()) return NULL;
 595
 596   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
 597
 598   Standard_Integer iErr;
 599   TopoDS_Shape aS;
 600   GEOMAlgo_Gluer1 aGluer;
 601   GEOMAlgo_ListIteratorOfListOfCoupleOfShapes aItCS;
 602   GEOMAlgo_CoupleOfShapes aCS;
 603   GEOMAlgo_ListOfCoupleOfShapes aLCS;
 604
 605   //aGluer = new GEOMAlgo_Gluer1;
 606   aGluer.SetShape(aShape);
 607   aGluer.SetTolerance(theTolerance);
 608   aGluer.Perform();
 609   iErr = aGluer.ErrorStatus();
 610   if (iErr) return NULL;
 611
 612   TopTools_ListOfShape listShape;
 613   const GEOMAlgo_ListOfCoupleOfShapes& aLCSG = aGluer.GluedFaces();
 614   // Access to faces
 615   aItCS.Initialize(aLCSG);
 616   for (; aItCS.More(); aItCS.Next()) {
 617     const GEOMAlgo_CoupleOfShapes& aCSG = aItCS.Value();
 618     listShape.Append(aCSG.Shape1());
 619   }
 620
 621   TopTools_ListIteratorOfListOfShape itSub (listShape);
 622   TCollection_AsciiString anAsciiList, anEntry;
 623   TopTools_IndexedMapOfShape anIndices;
 624   TopExp::MapShapes(aShape, anIndices);
 625   Handle(TColStd_HArray1OfInteger) anArray;
 626   Handle(GEOM_Object) anObj;
 627   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
 628     TopoDS_Shape aValue = itSub.Value();
 629     anArray = new TColStd_HArray1OfInteger(1,1);
 630     anArray->SetValue(1, anIndices.FindIndex(aValue));
 631     anObj = GetEngine()->AddSubShape(theShape, anArray);
 632     if (!anObj.IsNull()) {
 633       aSeq->Append(anObj);
 634
 635       // for python command
 636       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
 637       anAsciiList += anEntry;
 638       anAsciiList += ",";
 639     }
 640   }
 641
 642   //Make a Python command
 643   if(anAsciiList.Length()>0)
 644     anAsciiList.Trunc(anAsciiList.Length() - 1);
 645   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 646   GEOM::TPythonDump pd (aFunction, /*append=*/true);
 647   pd << "[" << anAsciiList.ToCString();
 648   pd << "] = geompy.GetGlueFaces(" << theShape << ", " << theTolerance << ")";
 649
 650   SetErrorCode(OK);
 651
 652   return aSeq;
 653 }
 654
 655
 656 //=============================================================================
 657 /*!
 658  *  MakeGlueFacesByList
 659  */
 660 //=============================================================================
 661 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueFacesByList
 662                                                 (Handle(GEOM_Object) theShape,
 663                                                  const Standard_Real theTolerance,
 664                                                  list<Handle(GEOM_Object)> theFaces,
 665                                                  const Standard_Boolean doKeepNonSolids)
 666 {
 667   SetErrorCode(KO);
 668
 669   if (theShape.IsNull()) return NULL;
 670
 671   //Add a new Glued object
 672   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
 673
 674   //Add a new Glue function
 675   Handle(GEOM_Function) aFunction;
 676   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_FACES_BY_LIST);
 677   if (aFunction.IsNull()) return NULL;
 678
 679   //Check if the function is set correctly
 680   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
 681
 682   GEOMImpl_IGlue aCI (aFunction);
 683
 684   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 685   if (aRefShape.IsNull()) return NULL;
 686
 687   aCI.SetBase(aRefShape);
 688   aCI.SetTolerance(theTolerance);
 689   aCI.SetKeepNonSolids(doKeepNonSolids);
 690
 691   Handle(TColStd_HSequenceOfTransient) aFaces = new TColStd_HSequenceOfTransient;
 692   list<Handle(GEOM_Object)>::iterator it = theFaces.begin();
 693   for (; it != theFaces.end(); it++) {
 694     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 695     if (aRefSh.IsNull()) {
 696       SetErrorCode("NULL argument shape for the shape construction");
 697       return NULL;
 698     }
 699     aFaces->Append(aRefSh);
 700   }
 701   aCI.SetFaces(aFaces);
 702
 703   //Compute the sub-shape value
 704   Standard_Boolean isWarning = Standard_False;
 705   try {
 706 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
 707     OCC_CATCH_SIGNALS;
 708 #endif
 709     if (!GetSolver()->ComputeFunction(aFunction)) {
 710       SetErrorCode("Shape driver failed to glue faces");
 711       return NULL;
 712     }
 713   }
 714   catch (Standard_Failure) {
 715     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 716     SetErrorCode(aFail->GetMessageString());
 717     // to provide warning
 718     if (!aFunction->GetValue().IsNull()) {
 719       isWarning = Standard_True;
 720     } else {
 721       return NULL;
 722     }
 723   }
 724
 725   //Make a Python command
 726
 727   GEOM::TPythonDump pd(aFunction);
 728   pd << aGlued << " = geompy.MakeGlueFacesByList("
 729      << theShape << ", " << theTolerance << ", [";
 730   // Faces
 731   it = theFaces.begin();
 732   if (it != theFaces.end()) {
 733     pd << (*it++);
 734     while (it != theFaces.end()) {
 735       pd << ", " << (*it++);
 736     }
 737   }
 738   pd << "])";
 739
 740
 741   // to provide warning
 742   if (!isWarning) SetErrorCode(OK);
 743   return aGlued;
 744 }
 745
 746
 747
 748 //=============================================================================
 749 /*!
 750  *  MakeExplode
 751  */
 752 //=============================================================================
 753 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::MakeExplode
 754                                           (Handle(GEOM_Object)    theShape,
 755                                            const Standard_Integer theShapeType,
 756                                            const Standard_Boolean isSorted)
 757 {
 758   SetErrorCode(KO);
 759
 760   if (theShape.IsNull()) return NULL;
 761   TopoDS_Shape aShape = theShape->GetValue();
 762   if (aShape.IsNull()) return NULL;
 763
 764   Handle(GEOM_Function) aMainShape = theShape->GetLastFunction();
 765
 766   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
 767   Handle(GEOM_Object) anObj;
 768   TopTools_MapOfShape mapShape;
 769   TopTools_ListOfShape listShape;
 770
 771   if (aShape.ShapeType() == TopAbs_COMPOUND &&
 772       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 773        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
 774        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND)) {
 775     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
 776     for (; It.More(); It.Next()) {
 777       if (mapShape.Add(It.Value())) {
 778         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 779             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
 780           listShape.Append(It.Value());
 781         }
 782       }
 783     }
 784   } else {
 785     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
 786     for (; exp.More(); exp.Next())
 787       if (mapShape.Add(exp.Current()))
 788         listShape.Append(exp.Current());
 789   }
 790
 791   if (listShape.IsEmpty()) {
 792     //SetErrorCode("The given shape has no sub-shapes of the requested type");
 793     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
 794     return aSeq;
 795   }
 796
 797   if (isSorted)
 798     SortShapes(listShape);
 799
 800   TopTools_IndexedMapOfShape anIndices;
 801   TopExp::MapShapes(aShape, anIndices);
 802   Handle(TColStd_HArray1OfInteger) anArray;
 803
 804   TopTools_ListIteratorOfListOfShape itSub (listShape);
 805   TCollection_AsciiString anAsciiList, anEntry;
 806   for (int index = 1; itSub.More(); itSub.Next(), ++index)
 807   {
 808     TopoDS_Shape aValue = itSub.Value();
 809     anArray = new TColStd_HArray1OfInteger(1,1);
 810     anArray->SetValue(1, anIndices.FindIndex(aValue));
 811
 812     //anObj = GetEngine()->AddSubShape(theShape, anArray);
 813     {
 814       anObj = GetEngine()->AddObject(GetDocID(), GEOM_SUBSHAPE);
 815       Handle(GEOM_Function) aFunction = anObj->AddFunction(GEOM_Object::GetSubShapeID(), 1);
 816       if (aFunction.IsNull()) return aSeq;
 817
 818       GEOM_ISubShape aSSI (aFunction);
 819       aSSI.SetMainShape(aMainShape);
 820       aSSI.SetIndices(anArray);
 821
 822       // Set function value directly, as we know it.
 823       // Usage of Solver here would lead to significant loss of time,
 824       // because GEOM_SubShapeDriver will build TopTools_IndexedMapOfShape
 825       // on the main shape for each being calculated sub-shape separately.
 826       aFunction->SetValue(aValue);
 827     }
 828
 829     if (!anObj.IsNull()) {
 830       aSeq->Append(anObj);
 831
 832       // for python command
 833       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
 834       anAsciiList += anEntry;
 835       anAsciiList += ",";
 836     }
 837   }
 838
 839   //Make a Python command
 840   anAsciiList.Trunc(anAsciiList.Length() - 1);
 841
 842   GEOM::TPythonDump pd (aMainShape, /*append=*/true);
 843   pd << "[" << anAsciiList.ToCString();
 844   pd << "] = geompy.SubShapeAll" << (isSorted ? "Sorted(" : "(");
 845   pd << theShape << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
 846
 847   SetErrorCode(OK);
 848
 849   return aSeq;
 850 }
 851
 852 //=============================================================================
 853 /*!
 854  *  SubShapeAllIDs
 855  */
 856 //=============================================================================
 857 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::SubShapeAllIDs
 858                                           (Handle(GEOM_Object)    theShape,
 859                                            const Standard_Integer theShapeType,
 860                                            const Standard_Boolean isSorted)
 861 {
 862   SetErrorCode(KO);
 863
 864   if (theShape.IsNull()) return NULL;
 865   TopoDS_Shape aShape = theShape->GetValue();
 866   if (aShape.IsNull()) return NULL;
 867
 868   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
 869   TopTools_MapOfShape mapShape;
 870   TopTools_ListOfShape listShape;
 871
 872   if (aShape.ShapeType() == TopAbs_COMPOUND &&
 873       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 874        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
 875        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND)) {
 876     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
 877     for (; It.More(); It.Next()) {
 878       if (mapShape.Add(It.Value())) {
 879         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 880             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
 881           listShape.Append(It.Value());
 882         }
 883       }
 884     }
 885   } else {
 886     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
 887     for (; exp.More(); exp.Next())
 888       if (mapShape.Add(exp.Current()))
 889         listShape.Append(exp.Current());
 890   }
 891
 892   if (listShape.IsEmpty()) {
 893     //SetErrorCode("The given shape has no sub-shapes of the requested type");
 894     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
 895     return aSeq;
 896   }
 897
 898   if (isSorted)
 899     SortShapes(listShape);
 900
 901   TopTools_IndexedMapOfShape anIndices;
 902   TopExp::MapShapes(aShape, anIndices);
 903   Handle(TColStd_HArray1OfInteger) anArray;
 904
 905   TopTools_ListIteratorOfListOfShape itSub (listShape);
 906   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
 907     TopoDS_Shape aValue = itSub.Value();
 908     aSeq->Append(anIndices.FindIndex(aValue));
 909   }
 910
 911   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 912
 913   //Make a Python command
 914   GEOM::TPythonDump pd (aFunction, /*append=*/true);
 915   pd << "listSubShapeIDs = geompy.SubShapeAll";
 916   pd << (isSorted ? "SortedIDs(" : "IDs(");
 917   pd << theShape << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
 918
 919   SetErrorCode(OK);
 920   return aSeq;
 921 }
 922
 923 //=============================================================================
 924 /*!
 925  *  GetSubShape
 926  */
 927 //=============================================================================
 928 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetSubShape
 929                                           (Handle(GEOM_Object)    theMainShape,
 930                                            const Standard_Integer theID)
 931 {
 932   SetErrorCode(KO);
 933
 934   if (theMainShape.IsNull()) return NULL;
 935
 936   Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
 937   anArray->SetValue(1, theID);
 938   Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theMainShape, anArray,true);
 939   if (anObj.IsNull()) {
 940     SetErrorCode("Can not get a sub-shape with the given ID");
 941     return NULL;
 942   }
 943
 944   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
 945
 946   //Make a Python command
 947   GEOM::TPythonDump(aFunction) << anObj << " = geompy.GetSubShape("
 948                                << theMainShape << ", [" << theID << "])";
 949
 950   SetErrorCode(OK);
 951   return anObj;
 952 }
 953
 954 //=============================================================================
 955 /*!
 956  *  GetSubShapeIndex
 957  */
 958 //=============================================================================
 959 Standard_Integer GEOMImpl_IShapesOperations::GetSubShapeIndex (Handle(GEOM_Object) theMainShape,
 960                                                                Handle(GEOM_Object) theSubShape)
 961 {
 962   SetErrorCode(KO);
 963
 964   TopoDS_Shape aMainShape = theMainShape->GetValue();
 965   TopoDS_Shape aSubShape = theSubShape->GetValue();
 966
 967   if (aMainShape.IsNull() || aSubShape.IsNull()) return -1;
 968
 969   TopTools_IndexedMapOfShape anIndices;
 970   TopExp::MapShapes(aMainShape, anIndices);
 971   if (anIndices.Contains(aSubShape)) {
 972     SetErrorCode(OK);
 973     return anIndices.FindIndex(aSubShape);
 974   }
 975
 976   return -1;
 977 }
 978
 979 //=============================================================================
 980 /*!
 981  *  GetTopologyIndex
 982  */
 983 //=============================================================================
 984 Standard_Integer GEOMImpl_IShapesOperations::GetTopologyIndex (Handle(GEOM_Object) theMainShape,
 985                                                                Handle(GEOM_Object) theSubShape)
 986 {
 987   SetErrorCode(OK);
 988
 989   TopoDS_Shape aMainShape = theMainShape->GetValue();
 990   TopoDS_Shape aSubShape = theSubShape->GetValue();
 991
 992   if (aMainShape.IsNull() || aSubShape.IsNull()) {
 993     SetErrorCode("Null argument shape given");
 994     return -1;
 995   }
 996
 997   int index = 1;
 998   if (aSubShape.ShapeType() == TopAbs_COMPOUND) {
 999     TopoDS_Iterator it;
1000     TopTools_ListOfShape CL;
1001     CL.Append(aMainShape);
1002     TopTools_ListIteratorOfListOfShape itC;
1003     for (itC.Initialize(CL); itC.More(); itC.Next()) {
1004       for (it.Initialize(itC.Value()); it.More(); it.Next()) {
1005         if (it.Value().ShapeType() == TopAbs_COMPOUND) {
1006           if (it.Value().IsSame(aSubShape))
1007             return index;
1008           else
1009             index++;
1010           CL.Append(it.Value());
1011         }
1012       }
1013     }
1014   } else {
1015     TopExp_Explorer anExp (aMainShape, aSubShape.ShapeType());
1016     TopTools_MapOfShape M;
1017     for (; anExp.More(); anExp.Next()) {
1018       if (M.Add(anExp.Current())) {
1019         if (anExp.Current().IsSame(aSubShape))
1020           return index;
1021         index++;
1022       }
1023     }
1024   }
1025
1026   SetErrorCode("The sub-shape does not belong to the main shape");
1027   return -1;
1028 }
1029
1030 //=============================================================================
1031 /*!
1032  *  GetShapeTypeString
1033  */
1034 //=============================================================================
1035 TCollection_AsciiString GEOMImpl_IShapesOperations::GetShapeTypeString (Handle(GEOM_Object) theShape)
1036 {
1037   SetErrorCode(KO);
1038
1039   TCollection_AsciiString aTypeName ("Null Shape");
1040
1041   TopoDS_Shape aShape = theShape->GetValue();
1042   if (aShape.IsNull())
1043     return aTypeName;
1044
1045   switch (aShape.ShapeType() )
1046   {
1047   case TopAbs_COMPOUND:
1048     aTypeName = "Compound";
1049     break;
1050   case  TopAbs_COMPSOLID:
1051     aTypeName = "Compound Solid";
1052     break;
1053   case TopAbs_SOLID:
1054     aTypeName = "Solid";
1055     break;
1056   case TopAbs_SHELL:
1057     aTypeName = "Shell";
1058     break;
1059   case TopAbs_FACE:
1060     {
1061       BRepAdaptor_Surface surf (TopoDS::Face(aShape));
1062       if (surf.GetType() == GeomAbs_Plane)
1063         aTypeName = "Plane";
1064       else if (surf.GetType() == GeomAbs_Cylinder)
1065         aTypeName = "Cylindrical Face";
1066       else if (surf.GetType() == GeomAbs_Sphere)
1067         aTypeName = "Spherical Face";
1068       else if (surf.GetType() == GeomAbs_Torus)
1069         aTypeName = "Toroidal Face";
1070       else if (surf.GetType() == GeomAbs_Cone)
1071         aTypeName = "Conical Face";
1072       else
1073         aTypeName = "GEOM::FACE";
1074     }
1075     break;
1076   case TopAbs_WIRE:
1077     aTypeName = "Wire";
1078     break;
1079   case TopAbs_EDGE:
1080     {
1081       BRepAdaptor_Curve curv (TopoDS::Edge(aShape));
1082       if (curv.GetType() == GeomAbs_Line) {
1083         if ((Abs(curv.FirstParameter()) >= 1E6) ||
1084             (Abs(curv.LastParameter()) >= 1E6))
1085           aTypeName = "Line";
1086         else
1087           aTypeName = "Edge" ;
1088       } else if (curv.GetType() == GeomAbs_Circle) {
1089         if (curv.IsClosed())
1090           aTypeName = "Circle";
1091         else
1092           aTypeName = "Arc";
1093       } else {
1094         aTypeName = "Edge";
1095       }
1096     }
1097     break;
1098   case TopAbs_VERTEX:
1099     aTypeName = "Vertex";
1100     break;
1101   case TopAbs_SHAPE:
1102     aTypeName = "Shape";
1103     break;
1104   default:
1105     aTypeName = "Shape of unknown type";
1106   }
1107
1108   return aTypeName;
1109 }
1110
1111
1112 //=============================================================================
1113 /*!
1114  *  NumberOfFaces
1115  */
1116 //=============================================================================
1117 Standard_Integer GEOMImpl_IShapesOperations::NumberOfFaces (Handle(GEOM_Object) theShape)
1118 {
1119   SetErrorCode(KO);
1120
1121   Standard_Integer nb = 0;
1122
1123   if (theShape.IsNull()) return -1;
1124   TopoDS_Shape aShape = theShape->GetValue();
1125   if (aShape.IsNull()) return -1;
1126
1127   TopTools_MapOfShape mapShape;
1128
1129   TopExp_Explorer exp (aShape, TopAbs_FACE);
1130   for (; exp.More(); exp.Next())
1131     if (mapShape.Add(exp.Current()))
1132       nb++;
1133
1134   SetErrorCode(OK);
1135   return nb;
1136 }
1137
1138 //=============================================================================
1139 /*!
1140  *  NumberOfEdges
1141  */
1142 //=============================================================================
1143 Standard_Integer GEOMImpl_IShapesOperations::NumberOfEdges (Handle(GEOM_Object) theShape)
1144 {
1145   SetErrorCode(KO);
1146
1147   Standard_Integer nb = 0;
1148
1149   if (theShape.IsNull()) return -1;
1150   TopoDS_Shape aShape = theShape->GetValue();
1151   if (aShape.IsNull()) return -1;
1152
1153   TopTools_MapOfShape mapShape;
1154
1155   TopExp_Explorer exp (aShape, TopAbs_EDGE);
1156   for (; exp.More(); exp.Next())
1157     if (mapShape.Add(exp.Current()))
1158       nb++;
1159
1160   SetErrorCode(OK);
1161   return nb;
1162 }
1163
1164 //=============================================================================
1165 /*!
1166  *  ReverseShape
1167  */
1168 //=============================================================================
1169 Handle(GEOM_Object) GEOMImpl_IShapesOperations::ReverseShape(Handle(GEOM_Object) theShape)
1170 {
1171   SetErrorCode(KO);
1172
1173   if (theShape.IsNull()) return NULL;
1174
1175   //Add a new reversed object
1176   Handle(GEOM_Object) aReversed = GetEngine()->AddObject(GetDocID(), theShape->GetType());
1177
1178   //Add a new Revese function
1179   Handle(GEOM_Function) aFunction;
1180   aFunction = aReversed->AddFunction(GEOMImpl_ShapeDriver::GetID(), REVERSE_ORIENTATION);
1181   if (aFunction.IsNull()) return NULL;
1182
1183   //Check if the function is set correctly
1184   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
1185
1186   GEOMImpl_IShapes aSI (aFunction);
1187
1188   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1189   if (aRefShape.IsNull()) return NULL;
1190
1191   aSI.SetBase(aRefShape);
1192
1193   //Compute the sub-shape value
1194   try {
1195 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
1196     OCC_CATCH_SIGNALS;
1197 #endif
1198     if (!GetSolver()->ComputeFunction(aFunction)) {
1199       SetErrorCode("Shape driver failed to reverse shape");
1200       return NULL;
1201     }
1202   }
1203   catch (Standard_Failure) {
1204     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1205     SetErrorCode(aFail->GetMessageString());
1206     return NULL;
1207   }
1208
1209   //Make a Python command
1210   GEOM::TPythonDump(aFunction) << aReversed
1211     << " = geompy.ChangeOrientation(" << theShape << ")";
1212
1213   SetErrorCode(OK);
1214   return aReversed;
1215 }
1216
1217 //=============================================================================
1218 /*!
1219  *  GetFreeFacesIDs
1220  */
1221 //=============================================================================
1222 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetFreeFacesIDs
1223                                                  (Handle(GEOM_Object) theShape)
1224 {
1225   SetErrorCode(KO);
1226
1227   if (theShape.IsNull()) return NULL;
1228   TopoDS_Shape aShape = theShape->GetValue();
1229   if (aShape.IsNull()) return NULL;
1230
1231   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
1232
1233   TopTools_IndexedDataMapOfShapeListOfShape mapFaceBlocks;
1234   GEOMImpl_Block6Explorer::MapShapesAndAncestors
1235     (aShape, TopAbs_FACE, TopAbs_SOLID, mapFaceBlocks);
1236
1237   Standard_Integer ind = 1, nbFaces = mapFaceBlocks.Extent();
1238
1239   if (nbFaces == 0) {
1240     SetErrorCode("The given shape has no faces");
1241     return aSeq;
1242   }
1243
1244   TopTools_IndexedMapOfShape anIndices;
1245   TopExp::MapShapes(aShape, anIndices);
1246
1247   Standard_Integer id;
1248   for (; ind <= nbFaces; ind++) {
1249     if (mapFaceBlocks.FindFromIndex(ind).Extent() != 2) {
1250       id = anIndices.FindIndex(mapFaceBlocks.FindKey(ind));
1251       aSeq->Append(id);
1252     }
1253   }
1254
1255   //The explode doesn't change object so no new function is required.
1256   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1257
1258   //Make a Python command
1259   GEOM::TPythonDump(aFunction, /*append=*/true)
1260     << "listFreeFacesIDs = geompy.GetFreeFacesIDs(" << theShape << ")";
1261
1262   SetErrorCode(OK);
1263   return aSeq;
1264 }
1265
1266 //=======================================================================
1267 //function : GetSharedShapes
1268 //purpose  :
1269 //=======================================================================
1270
1271 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetSharedShapes
1272                                                 (Handle(GEOM_Object)    theShape1,
1273                                                  Handle(GEOM_Object)    theShape2,
1274                                                  const Standard_Integer theShapeType)
1275 {
1276   SetErrorCode(KO);
1277
1278   if (theShape1.IsNull() || theShape2.IsNull()) return NULL;
1279
1280   TopoDS_Shape aShape1 = theShape1->GetValue();
1281   TopoDS_Shape aShape2 = theShape2->GetValue();
1282
1283   if (aShape1.IsNull() || aShape2.IsNull()) return NULL;
1284
1285   TopTools_IndexedMapOfShape anIndices;
1286   TopExp::MapShapes(aShape1, anIndices);
1287   Handle(TColStd_HArray1OfInteger) anArray;
1288
1289   TopTools_IndexedMapOfShape mapShape1;
1290   TopExp::MapShapes(aShape1, TopAbs_ShapeEnum(theShapeType), mapShape1);
1291
1292   Handle(GEOM_Object) anObj;
1293   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1294   TCollection_AsciiString anAsciiList, anEntry;
1295
1296   TopTools_MapOfShape mapShape2;
1297   TopExp_Explorer exp (aShape2, TopAbs_ShapeEnum(theShapeType));
1298   for (; exp.More(); exp.Next()) {
1299     TopoDS_Shape aSS = exp.Current();
1300     if (mapShape2.Add(aSS) && mapShape1.Contains(aSS)) {
1301       anArray = new TColStd_HArray1OfInteger(1,1);
1302       anArray->SetValue(1, anIndices.FindIndex(aSS));
1303       anObj = GetEngine()->AddSubShape(theShape1, anArray);
1304       aSeq->Append(anObj);
1305
1306       // for python command
1307       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1308       anAsciiList += anEntry;
1309       anAsciiList += ",";
1310     }
1311   }
1312
1313   if (aSeq->IsEmpty()) {
1314     SetErrorCode("The given shapes have no shared sub-shapes of the requested type");
1315     return aSeq;
1316   }
1317
1318   //Make a Python command
1319   anAsciiList.Trunc(anAsciiList.Length() - 1);
1320
1321   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1322
1323   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
1324     << "] = geompy.GetSharedShapes(" << theShape1 << ", "
1325       << theShape2 << ", " << TopAbs_ShapeEnum(theShapeType) << ")";
1326
1327   SetErrorCode(OK);
1328   return aSeq;
1329 }
1330
1331 //=============================================================================
1332 /*!
1333  *
1334  */
1335 //=============================================================================
1336 static GEOM::TPythonDump& operator<< (GEOM::TPythonDump&   theDump,
1337                                       const GEOMAlgo_State theState)
1338 {
1339   switch (theState) {
1340   case GEOMAlgo_ST_IN:
1341     theDump << "geompy.GEOM.ST_IN";
1342     break;
1343   case GEOMAlgo_ST_OUT:
1344     theDump << "geompy.GEOM.ST_OUT";
1345     break;
1346   case GEOMAlgo_ST_ON:
1347     theDump << "geompy.GEOM.ST_ON";
1348     break;
1349   case GEOMAlgo_ST_ONIN:
1350     theDump << "geompy.GEOM.ST_ONIN";
1351     break;
1352   case GEOMAlgo_ST_ONOUT:
1353     theDump << "geompy.GEOM.ST_ONOUT";
1354     break;
1355   default:
1356     theDump << "geompy.GEOM.ST_UNKNOWN";
1357     break;
1358   }
1359   return theDump;
1360 }
1361
1362 //=======================================================================
1363 //function : checkTypeShapesOn
1364 /*!
1365  * \brief Checks if theShapeType parameter of GetShapesOnXXX() is OK
1366  * \param theShapeType - the shape type to check
1367  * \retval bool  - result of the check
1368  */
1369 //=======================================================================
1370
1371 bool GEOMImpl_IShapesOperations::checkTypeShapesOn(const Standard_Integer theShapeType)
1372 {
1373   if (theShapeType != TopAbs_VERTEX &&
1374       theShapeType != TopAbs_EDGE &&
1375       theShapeType != TopAbs_FACE &&
1376       theShapeType != TopAbs_SOLID) {
1377     SetErrorCode("Only solids, vertices, edges or faces can be found by this method");
1378     return false;
1379   }
1380   return true;
1381 }
1382
1383 //=======================================================================
1384 //function : makePlane
1385   /*!
1386    * \brief Creates Geom_Plane
1387     * \param theAx1 - shape object defining plane parameters
1388     * \retval Handle(Geom_Surface) - resulting surface
1389    */
1390 //=======================================================================
1391
1392 Handle(Geom_Surface) GEOMImpl_IShapesOperations::makePlane(const TopoDS_Shape& anAx1)
1393 {
1394   if (anAx1.ShapeType() != TopAbs_EDGE) return NULL;
1395   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
1396   TopoDS_Vertex V1, V2;
1397   TopExp::Vertices(anEdge, V1, V2, Standard_True);
1398   if (V1.IsNull() || V2.IsNull()) {
1399     SetErrorCode("Bad edge given for the plane normal vector");
1400     return NULL;
1401   }
1402   gp_Pnt aLoc = BRep_Tool::Pnt(V1);
1403   gp_Vec aVec (aLoc, BRep_Tool::Pnt(V2));
1404   if (aVec.Magnitude() < Precision::Confusion()) {
1405     SetErrorCode("Vector with null magnitude given");
1406     return NULL;
1407   }
1408   return new Geom_Plane(aLoc, aVec);
1409 }
1410
1411 //=======================================================================
1412 //function : makeCylinder
1413   /*!
1414    * \brief Creates Geom_CylindricalSurface
1415     * \param theAx1 - edge defining cylinder axis
1416     * \param theRadius - cylinder radius
1417     * \retval Handle(Geom_Surface) - resulting surface
1418    */
1419 //=======================================================================
1420
1421 Handle(Geom_Surface) GEOMImpl_IShapesOperations::makeCylinder(const TopoDS_Shape& anAxis,
1422                                                               const Standard_Real theRadius)
1423 {
1424   //Axis of the cylinder
1425   if (anAxis.ShapeType() != TopAbs_EDGE) {
1426     SetErrorCode("Not an edge given for the axis");
1427     return NULL;
1428   }
1429   TopoDS_Edge anEdge = TopoDS::Edge(anAxis);
1430   TopoDS_Vertex V1, V2;
1431   TopExp::Vertices(anEdge, V1, V2, Standard_True);
1432   if (V1.IsNull() || V2.IsNull()) {
1433     SetErrorCode("Bad edge given for the axis");
1434     return NULL;
1435   }
1436   gp_Pnt aLoc = BRep_Tool::Pnt(V1);
1437   gp_Vec aVec (aLoc, BRep_Tool::Pnt(V2));
1438   if (aVec.Magnitude() < Precision::Confusion()) {
1439     SetErrorCode("Vector with null magnitude given");
1440     return NULL;
1441   }
1442
1443   gp_Ax3 anAx3 (aLoc, aVec);
1444   return new Geom_CylindricalSurface(anAx3, theRadius);
1445 }
1446
1447
1448 //=======================================================================
1449 //function : getShapesOnBoxIDs
1450   /*!
1451    * \brief Find IDs of subshapes complying with given status about surface
1452     * \param theBox - the box to check state of subshapes against
1453     * \param theShape - the shape to explore
1454     * \param theShapeType - type of subshape of theShape
1455     * \param theState - required state
1456     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1457    */
1458 //=======================================================================
1459
1460 Handle(TColStd_HSequenceOfInteger)
1461   GEOMImpl_IShapesOperations::getShapesOnBoxIDs(const Handle(GEOM_Object)& theBox,
1462                                                 const Handle(GEOM_Object)& theShape,
1463                                                 const Standard_Integer theShapeType,
1464                                                 GEOMAlgo_State theState)
1465 {
1466   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
1467
1468   TopoDS_Shape aBox = theBox->GetValue();
1469   TopoDS_Shape aShape = theShape->GetValue();
1470
1471   // Check presence of triangulation, build if need
1472   if (!CheckTriangulation(aShape)) {
1473     SetErrorCode("Cannot build triangulation on the shape");
1474     return aSeqOfIDs;
1475   }
1476
1477   // Call algo
1478   GEOMAlgo_FinderShapeOn2 aFinder;
1479   Standard_Real aTol = 0.0001; // default value
1480
1481   Handle(GEOMAlgo_ClsfBox) aClsfBox = new GEOMAlgo_ClsfBox;
1482   aClsfBox->SetBox(aBox);
1483
1484   aFinder.SetShape(aShape);
1485   aFinder.SetTolerance(aTol);
1486   aFinder.SetClsf(aClsfBox);
1487   aFinder.SetShapeType( (TopAbs_ShapeEnum)theShapeType );
1488   aFinder.SetState(theState);
1489   aFinder.Perform();
1490
1491   // Interprete results
1492   Standard_Integer iErr = aFinder.ErrorStatus();
1493   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
1494   if (iErr) {
1495     MESSAGE(" iErr : " << iErr);
1496     TCollection_AsciiString aMsg (" iErr : ");
1497     aMsg += TCollection_AsciiString(iErr);
1498     SetErrorCode(aMsg);
1499     return aSeqOfIDs;
1500   }
1501   Standard_Integer iWrn = aFinder.WarningStatus();
1502   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
1503   if (iWrn) {
1504     MESSAGE(" *** iWrn : " << iWrn);
1505   }
1506
1507   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
1508
1509   if (listSS.Extent() < 1) {
1510     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
1511     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
1512     return aSeqOfIDs;
1513   }
1514
1515   // Fill sequence of object IDs
1516   aSeqOfIDs = new TColStd_HSequenceOfInteger;
1517
1518   TopTools_IndexedMapOfShape anIndices;
1519   TopExp::MapShapes(aShape, anIndices);
1520
1521   TopTools_ListIteratorOfListOfShape itSub (listSS);
1522   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
1523     int id = anIndices.FindIndex(itSub.Value());
1524     aSeqOfIDs->Append(id);
1525   }
1526
1527   return aSeqOfIDs;
1528 }
1529
1530
1531 //=======================================================================
1532 //function : GetShapesOnBoxIDs
1533 /*!
1534    * \brief Find subshapes complying with given status about surface
1535     * \param theBox - the box to check state of subshapes against
1536     * \param theShape - the shape to explore
1537     * \param theShapeType - type of subshape of theShape
1538     * \param theState - required state
1539     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1540  */
1541 //=======================================================================
1542
1543 Handle(TColStd_HSequenceOfInteger)
1544     GEOMImpl_IShapesOperations::GetShapesOnBoxIDs(const Handle(GEOM_Object)& theBox,
1545                                                   const Handle(GEOM_Object)& theShape,
1546                                                   const Standard_Integer theShapeType,
1547                                                   GEOMAlgo_State theState)
1548 {
1549   // Find subshapes ids
1550   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1551     getShapesOnBoxIDs (theBox, theShape, theShapeType, theState);
1552   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1553     return NULL;
1554
1555   // The GetShapesOnBox() doesn't change object so no new function is required.
1556   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theBox)->GetLastFunction();
1557
1558   // Make a Python command
1559   GEOM::TPythonDump(aFunction)
1560     << "listShapesOnBoxIDs = geompy.GetShapesOnBoxIDs("
1561     << theBox << ", "
1562     << theShape << ", "
1563     << TopAbs_ShapeEnum(theShapeType) << ", "
1564     << theState << ")";
1565
1566   SetErrorCode(OK);
1567   return aSeqOfIDs;
1568 }
1569
1570 //=======================================================================
1571 //function : GetShapesOnBox
1572 /*!
1573    * \brief Find subshapes complying with given status about surface
1574     * \param theBox - the box to check state of subshapes against
1575     * \param theShape - the shape to explore
1576     * \param theShapeType - type of subshape of theShape
1577     * \param theState - required state
1578     * \retval Handle(TColStd_HSequenceOfTransient) - found subshapes
1579  */
1580 //=======================================================================
1581
1582 Handle(TColStd_HSequenceOfTransient)
1583     GEOMImpl_IShapesOperations::GetShapesOnBox(const Handle(GEOM_Object)& theBox,
1584                                                const Handle(GEOM_Object)&  theShape,
1585                                                const Standard_Integer theShapeType,
1586                                                GEOMAlgo_State theState)
1587 {
1588   // Find subshapes ids
1589   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1590     getShapesOnBoxIDs (theBox, theShape, theShapeType, theState);
1591   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1592     return NULL;
1593
1594   // Find objects by indices
1595   TCollection_AsciiString anAsciiList;
1596   Handle(TColStd_HSequenceOfTransient) aSeq;
1597   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
1598   if ( aSeq.IsNull() || aSeq->IsEmpty() )
1599     return NULL;
1600
1601   // Make a Python command
1602
1603   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
1604   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1605
1606   GEOM::TPythonDump(aFunction)
1607     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnBox("
1608     << theBox << ", "
1609     << theShape << ", "
1610     << TopAbs_ShapeEnum(theShapeType) << ", "
1611     << theState << ")";
1612
1613   SetErrorCode(OK);
1614   return aSeq;
1615 }
1616
1617
1618 //=======================================================================
1619 //function : getShapesOnShapeIDs
1620 /*!
1621  * \brief Find IDs of subshapes complying with given status about surface
1622  * \param theCheckShape - the shape to check state of subshapes against
1623  * \param theShape - the shape to explore
1624  * \param theShapeType - type of subshape of theShape
1625  * \param theState - required state
1626  * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1627  */
1628 //=======================================================================
1629
1630 Handle(TColStd_HSequenceOfInteger)
1631   GEOMImpl_IShapesOperations::getShapesOnShapeIDs
1632                                  (const Handle(GEOM_Object)& theCheckShape,
1633                                   const Handle(GEOM_Object)& theShape,
1634                                   const Standard_Integer theShapeType,
1635                                   GEOMAlgo_State theState)
1636 {
1637   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
1638
1639   TopoDS_Shape aCheckShape = theCheckShape->GetValue();
1640   TopoDS_Shape aShape = theShape->GetValue();
1641   TopTools_ListOfShape res;
1642
1643   // Check presence of triangulation, build if need
1644   if (!CheckTriangulation(aShape)) {
1645     SetErrorCode("Cannot build triangulation on the shape");
1646     return aSeqOfIDs;
1647   }
1648
1649   // Call algo
1650   GEOMAlgo_FinderShapeOn2 aFinder;
1651   Standard_Real aTol = 0.0001; // default value
1652
1653   Handle(GEOMAlgo_ClsfSolid) aClsfSolid = new GEOMAlgo_ClsfSolid;
1654   aClsfSolid->SetShape(aCheckShape);
1655
1656   aFinder.SetShape(aShape);
1657   aFinder.SetTolerance(aTol);
1658   aFinder.SetClsf(aClsfSolid);
1659   aFinder.SetShapeType( (TopAbs_ShapeEnum)theShapeType );
1660   aFinder.SetState(theState);
1661   aFinder.Perform();
1662
1663   // Interprete results
1664   Standard_Integer iErr = aFinder.ErrorStatus();
1665   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
1666   if (iErr) {
1667     MESSAGE(" iErr : " << iErr);
1668     TCollection_AsciiString aMsg (" iErr : ");
1669     aMsg += TCollection_AsciiString(iErr);
1670     SetErrorCode(aMsg);
1671     return aSeqOfIDs;
1672   }
1673   Standard_Integer iWrn = aFinder.WarningStatus();
1674   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
1675   if (iWrn) {
1676     MESSAGE(" *** iWrn : " << iWrn);
1677   }
1678
1679   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
1680
1681   if (listSS.Extent() < 1) {
1682     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
1683     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
1684   }
1685
1686   // Fill sequence of object IDs
1687   aSeqOfIDs = new TColStd_HSequenceOfInteger;
1688
1689   TopTools_IndexedMapOfShape anIndices;
1690   TopExp::MapShapes(aShape, anIndices);
1691
1692   TopTools_ListIteratorOfListOfShape itSub (listSS);
1693   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
1694     int id = anIndices.FindIndex(itSub.Value());
1695     aSeqOfIDs->Append(id);
1696   }
1697
1698   return aSeqOfIDs;
1699 }
1700
1701
1702 //=======================================================================
1703 //function : GetShapesOnShapeIDs
1704 /*!
1705  * \brief Find subshapes complying with given status about surface
1706  * \param theCheckShape - the shape to check state of subshapes against
1707  * \param theShape - the shape to explore
1708  * \param theShapeType - type of subshape of theShape
1709  * \param theState - required state
1710  * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1711  */
1712 //=======================================================================
1713
1714 Handle(TColStd_HSequenceOfInteger)
1715     GEOMImpl_IShapesOperations::GetShapesOnShapeIDs
1716                             (const Handle(GEOM_Object)& theCheckShape,
1717                              const Handle(GEOM_Object)& theShape,
1718                              const Standard_Integer theShapeType,
1719                              GEOMAlgo_State theState)
1720 {
1721   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1722     getShapesOnShapeIDs (theCheckShape, theShape, theShapeType, theState);
1723
1724   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1725     return NULL;
1726
1727   // The GetShapesOnShape() doesn't change object so no new function is required.
1728   Handle(GEOM_Function) aFunction =
1729     GEOM::GetCreatedLast(theShape,theCheckShape)->GetLastFunction();
1730
1731   // Make a Python command
1732   GEOM::TPythonDump(aFunction)
1733     << "listShapesOnBoxIDs = geompy.GetShapesOnShapeIDs("
1734     << theCheckShape << ", "
1735     << theShape << ", "
1736     << TopAbs_ShapeEnum(theShapeType) << ", "
1737     << theState << ")";
1738
1739   SetErrorCode(OK);
1740   return aSeqOfIDs;
1741 }
1742
1743
1744 //=======================================================================
1745 //function : GetShapesOnShape
1746 /*!
1747  * \brief Find subshapes complying with given status about surface
1748  * \param theCheckShape - the shape to check state of subshapes against
1749  * \param theShape - the shape to explore
1750  * \param theShapeType - type of subshape of theShape
1751  * \param theState - required state
1752  * \retval Handle(TColStd_HSequenceOfTransient) - found subshapes
1753  */
1754 //=======================================================================
1755
1756 Handle(TColStd_HSequenceOfTransient)
1757   GEOMImpl_IShapesOperations::GetShapesOnShape
1758                              (const Handle(GEOM_Object)& theCheckShape,
1759                               const Handle(GEOM_Object)&  theShape,
1760                               const Standard_Integer theShapeType,
1761                               GEOMAlgo_State theState)
1762 {
1763   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1764     getShapesOnShapeIDs (theCheckShape, theShape, theShapeType, theState);
1765   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1766     return NULL;
1767
1768   // Find objects by indices
1769   TCollection_AsciiString anAsciiList;
1770   Handle(TColStd_HSequenceOfTransient) aSeq;
1771   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
1772
1773   if ( aSeq.IsNull() || aSeq->IsEmpty() )
1774     return NULL;
1775
1776   // Make a Python command
1777
1778   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
1779   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1780
1781   GEOM::TPythonDump(aFunction)
1782     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnShape("
1783     << theCheckShape << ", "
1784     << theShape << ", "
1785     << TopAbs_ShapeEnum(theShapeType) << ", "
1786     << theState << ")";
1787
1788
1789   SetErrorCode(OK);
1790   return aSeq;
1791 }
1792
1793
1794 //=======================================================================
1795 //function : GetShapesOnShapeAsCompound
1796 //=======================================================================
1797
1798 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetShapesOnShapeAsCompound
1799                              (const Handle(GEOM_Object)& theCheckShape,
1800                               const Handle(GEOM_Object)&  theShape,
1801                               const Standard_Integer theShapeType,
1802                               GEOMAlgo_State theState)
1803 {
1804   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1805     getShapesOnShapeIDs (theCheckShape, theShape, theShapeType, theState);
1806
1807   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1808     return NULL;
1809
1810   // Find objects by indices
1811   TCollection_AsciiString anAsciiList;
1812   Handle(TColStd_HSequenceOfTransient) aSeq;
1813   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
1814
1815   if ( aSeq.IsNull() || aSeq->IsEmpty() )
1816     return NULL;
1817
1818   TopoDS_Compound aCompound;
1819   BRep_Builder B;
1820   B.MakeCompound(aCompound);
1821   int i = 1;
1822   for(; i<=aSeq->Length(); i++) {
1823     Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(aSeq->Value(i));
1824     TopoDS_Shape aShape_i = anObj->GetValue();
1825     B.Add(aCompound,aShape_i);
1826   }
1827
1828   //Add a new result object
1829   Handle(GEOM_Object) aRes = GetEngine()->AddObject(GetDocID(), GEOM_SHAPES_ON_SHAPE);
1830   Handle(GEOM_Function) aFunction =
1831     aRes->AddFunction(GEOMImpl_ShapeDriver::GetID(), SHAPES_ON_SHAPE);
1832   aFunction->SetValue(aCompound);
1833
1834   GEOM::TPythonDump(aFunction)
1835     << aRes << " = geompy.GetShapesOnShapeAsCompound("
1836     << theCheckShape << ", "
1837     << theShape << ", "
1838     << TopAbs_ShapeEnum(theShapeType) << ", "
1839     << theState << ")";
1840
1841   SetErrorCode(OK);
1842
1843   return aRes;
1844 }
1845
1846
1847 //=======================================================================
1848 //function : getShapesOnSurfaceIDs
1849   /*!
1850    * \brief Find IDs of subshapes complying with given status about surface
1851     * \param theSurface - the surface to check state of subshapes against
1852     * \param theShape - the shape to explore
1853     * \param theShapeType - type of subshape of theShape
1854     * \param theState - required state
1855     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1856    */
1857 //=======================================================================
1858
1859 Handle(TColStd_HSequenceOfInteger)
1860   GEOMImpl_IShapesOperations::getShapesOnSurfaceIDs(const Handle(Geom_Surface)& theSurface,
1861                                                     const TopoDS_Shape&         theShape,
1862                                                     TopAbs_ShapeEnum            theShapeType,
1863                                                     GEOMAlgo_State              theState)
1864 {
1865   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
1866
1867   // Check presence of triangulation, build if need
1868   if (!CheckTriangulation(theShape)) {
1869     SetErrorCode("Cannot build triangulation on the shape");
1870     return aSeqOfIDs;
1871   }
1872
1873   // Call algo
1874   GEOMAlgo_FinderShapeOn1 aFinder;
1875   Standard_Real aTol = 0.0001; // default value
1876
1877   aFinder.SetShape(theShape);
1878   aFinder.SetTolerance(aTol);
1879   aFinder.SetSurface(theSurface);
1880   aFinder.SetShapeType(theShapeType);
1881   aFinder.SetState(theState);
1882
1883   // Sets the minimal number of inner points for the faces that do not have own
1884   // inner points at all (for e.g. rectangular planar faces have just 2 triangles).
1885   // Default value=3
1886   aFinder.SetNbPntsMin(3);
1887   // Sets the maximal number of inner points for edges or faces.
1888   // It is usefull for the cases when this number is very big (e.g =2000) to improve
1889   // the performance. If this value =0, all inner points will be taken into account.
1890   // Default value=0
1891   aFinder.SetNbPntsMax(100);
1892
1893   aFinder.Perform();
1894
1895   // Interprete results
1896   Standard_Integer iErr = aFinder.ErrorStatus();
1897   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
1898   if (iErr) {
1899     MESSAGE(" iErr : " << iErr);
1900     TCollection_AsciiString aMsg (" iErr : ");
1901     aMsg += TCollection_AsciiString(iErr);
1902     SetErrorCode(aMsg);
1903     return aSeqOfIDs;
1904   }
1905   Standard_Integer iWrn = aFinder.WarningStatus();
1906   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
1907   if (iWrn) {
1908     MESSAGE(" *** iWrn : " << iWrn);
1909   }
1910
1911   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
1912
1913   if (listSS.Extent() < 1) {
1914     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
1915     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
1916     return aSeqOfIDs;
1917   }
1918
1919   // Fill sequence of object IDs
1920   aSeqOfIDs = new TColStd_HSequenceOfInteger;
1921
1922   TopTools_IndexedMapOfShape anIndices;
1923   TopExp::MapShapes(theShape, anIndices);
1924
1925   TopTools_ListIteratorOfListOfShape itSub (listSS);
1926   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
1927     int id = anIndices.FindIndex(itSub.Value());
1928     aSeqOfIDs->Append(id);
1929   }
1930
1931   return aSeqOfIDs;
1932 }
1933
1934 //=======================================================================
1935 //function : getObjectsShapesOn
1936 /*!
1937  * \brief Find shape objects and their entries by their ids
1938  * \param theShapeIDs - incoming shape ids
1939  * \param theShapeEntries - outgoing entries like "entry1, entry2, ..."
1940  * \retval Handle(TColStd_HSequenceOfTransient) - found shape objects
1941  */
1942 //=======================================================================
1943
1944 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::
1945  getObjectsShapesOn(const Handle(GEOM_Object)&                theShape,
1946                     const Handle(TColStd_HSequenceOfInteger)& theShapeIDs,
1947                     TCollection_AsciiString &                 theShapeEntries)
1948 {
1949   Handle(TColStd_HSequenceOfTransient) aSeq;
1950
1951   if ( !theShapeIDs.IsNull() && theShapeIDs->Length() > 0 )
1952   {
1953     aSeq = new TColStd_HSequenceOfTransient;
1954     Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1955     TCollection_AsciiString anEntry;
1956     for ( int i = 1; i <= theShapeIDs->Length(); ++i )
1957     {
1958       anArray->SetValue(1, theShapeIDs->Value( i ));
1959       Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theShape, anArray);
1960       aSeq->Append( anObj );
1961
1962       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1963       if ( i != 1 ) theShapeEntries += ",";
1964       theShapeEntries += anEntry;
1965     }
1966   }
1967   return aSeq;
1968 }
1969
1970 //=======================================================================
1971 //function : getShapesOnSurface
1972 /*!
1973    * \brief Find subshapes complying with given status about surface
1974     * \param theSurface - the surface to check state of subshapes against
1975     * \param theShape - the shape to explore
1976     * \param theShapeType - type of subshape of theShape
1977     * \param theState - required state
1978     * \param theShapeEntries - outgoing entries like "entry1, entry2, ..."
1979     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1980  */
1981 //=======================================================================
1982
1983 Handle(TColStd_HSequenceOfTransient)
1984     GEOMImpl_IShapesOperations::getShapesOnSurface(const Handle(Geom_Surface)& theSurface,
1985                                                    const Handle(GEOM_Object)&  theShape,
1986                                                    TopAbs_ShapeEnum            theShapeType,
1987                                                    GEOMAlgo_State              theState,
1988                                                    TCollection_AsciiString &   theShapeEntries)
1989 {
1990   // Find subshapes ids
1991   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1992     getShapesOnSurfaceIDs (theSurface, theShape->GetValue(), theShapeType, theState);
1993   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1994     return NULL;
1995
1996   return getObjectsShapesOn( theShape, aSeqOfIDs, theShapeEntries );
1997 }
1998
1999 //=============================================================================
2000 /*!
2001  *  GetShapesOnPlane
2002  */
2003 //=============================================================================
2004 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnPlane
2005                                         (const Handle(GEOM_Object)& theShape,
2006                                          const Standard_Integer     theShapeType,
2007                                          const Handle(GEOM_Object)& theAx1,
2008                                          const GEOMAlgo_State       theState)
2009 {
2010   SetErrorCode(KO);
2011
2012   if (theShape.IsNull() || theAx1.IsNull()) return NULL;
2013
2014   TopoDS_Shape aShape = theShape->GetValue();
2015   TopoDS_Shape anAx1  = theAx1->GetValue();
2016
2017   if (aShape.IsNull() || anAx1.IsNull()) return NULL;
2018
2019   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2020   if ( !checkTypeShapesOn( theShapeType ))
2021     return NULL;
2022
2023   // Create plane
2024   Handle(Geom_Surface) aPlane = makePlane( anAx1 );
2025   if ( aPlane.IsNull() )
2026     return NULL;
2027
2028   // Find objects
2029   TCollection_AsciiString anAsciiList;
2030   Handle(TColStd_HSequenceOfTransient) aSeq;
2031   aSeq = getShapesOnSurface( aPlane, theShape, aShapeType, theState, anAsciiList );
2032   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2033     return NULL;
2034
2035   // Make a Python command
2036
2037   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2038   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2039
2040   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2041     << "] = geompy.GetShapesOnPlane(" << theShape << ", "
2042       << aShapeType << ", " << theAx1 << ", " << theState << ")";
2043
2044   SetErrorCode(OK);
2045   return aSeq;
2046 }
2047
2048 //=============================================================================
2049 /*!
2050  *  GetShapesOnPlaneWithLocation
2051  */
2052 //=============================================================================
2053 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnPlaneWithLocation
2054                                         (const Handle(GEOM_Object)& theShape,
2055                                          const Standard_Integer     theShapeType,
2056                                          const Handle(GEOM_Object)& theAx1,
2057                                          const Handle(GEOM_Object)& thePnt,
2058                                          const GEOMAlgo_State       theState)
2059 {
2060   SetErrorCode(KO);
2061
2062   if (theShape.IsNull() || theAx1.IsNull() || thePnt.IsNull()) return NULL;
2063
2064   TopoDS_Shape aShape = theShape->GetValue();
2065   TopoDS_Shape anAx1  = theAx1->GetValue();
2066   TopoDS_Shape anPnt = thePnt->GetValue();
2067
2068   if (aShape.IsNull() || anAx1.IsNull() || anPnt.IsNull()) return NULL;
2069
2070   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2071   if ( !checkTypeShapesOn( theShapeType ))
2072     return NULL;
2073
2074   // Create plane
2075   if ( anAx1.ShapeType() != TopAbs_EDGE || anPnt.ShapeType() != TopAbs_VERTEX ) return NULL;
2076   TopoDS_Vertex V1, V2, V3;
2077   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
2078   TopExp::Vertices(anEdge, V1, V2, Standard_True);
2079
2080   if (V1.IsNull() || V2.IsNull()) {
2081     SetErrorCode("Bad edge given for the plane normal vector");
2082     return NULL;
2083   }
2084   V3 = TopoDS::Vertex(anPnt);
2085
2086   if(V3.IsNull()) {
2087     SetErrorCode("Bad vertex given for the plane location");
2088       return NULL;
2089   }
2090   gp_Pnt aLoc = BRep_Tool::Pnt(V3);
2091   gp_Vec aVec(BRep_Tool::Pnt(V1),BRep_Tool::Pnt(V2));
2092
2093   if (aVec.Magnitude() < Precision::Confusion()) {
2094      SetErrorCode("Vector with null magnitude given");
2095     return NULL;
2096   }
2097   Handle(Geom_Surface) aPlane = new Geom_Plane(aLoc, aVec);
2098
2099   if ( aPlane.IsNull() )
2100     return NULL;
2101
2102   // Find objects
2103   TCollection_AsciiString anAsciiList;
2104   Handle(TColStd_HSequenceOfTransient) aSeq;
2105   aSeq = getShapesOnSurface( aPlane, theShape, aShapeType, theState, anAsciiList );
2106   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2107     return NULL;
2108
2109   // Make a Python command
2110
2111   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2112   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2113
2114   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2115     << "] = geompy.GetShapesOnPlaneWithLocation(" << theShape << ", "
2116     << aShapeType << ", " << theAx1 << ", "<< thePnt <<", " << theState << ")";
2117
2118   SetErrorCode(OK);
2119   return aSeq;
2120 }
2121
2122 //=============================================================================
2123 /*!
2124  *  GetShapesOnCylinder
2125  */
2126 //=============================================================================
2127 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnCylinder
2128                                           (const Handle(GEOM_Object)& theShape,
2129                                            const Standard_Integer     theShapeType,
2130                                            const Handle(GEOM_Object)& theAxis,
2131                                            const Standard_Real        theRadius,
2132                                            const GEOMAlgo_State       theState)
2133 {
2134   SetErrorCode(KO);
2135
2136   if (theShape.IsNull() || theAxis.IsNull()) return NULL;
2137
2138   TopoDS_Shape aShape = theShape->GetValue();
2139   TopoDS_Shape anAxis = theAxis->GetValue();
2140
2141   if (aShape.IsNull() || anAxis.IsNull()) return NULL;
2142
2143   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2144   if ( !checkTypeShapesOn( aShapeType ))
2145     return NULL;
2146
2147   // Create a cylinder surface
2148   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
2149   if ( aCylinder.IsNull() )
2150     return NULL;
2151
2152   // Find objects
2153   TCollection_AsciiString anAsciiList;
2154   Handle(TColStd_HSequenceOfTransient) aSeq;
2155   aSeq = getShapesOnSurface( aCylinder, theShape, aShapeType, theState, anAsciiList );
2156   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2157     return NULL;
2158
2159   // Make a Python command
2160
2161   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2162   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2163
2164   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2165     << "] = geompy.GetShapesOnCylinder(" << theShape << ", " << aShapeType
2166       << ", " << theAxis << ", " << theRadius << ", " << theState << ")";
2167
2168   SetErrorCode(OK);
2169   return aSeq;
2170 }
2171
2172 //=============================================================================
2173 /*!
2174  *  GetShapesOnSphere
2175  */
2176 //=============================================================================
2177 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnSphere
2178                                           (const Handle(GEOM_Object)& theShape,
2179                                            const Standard_Integer     theShapeType,
2180                                            const Handle(GEOM_Object)& theCenter,
2181                                            const Standard_Real        theRadius,
2182                                            const GEOMAlgo_State       theState)
2183 {
2184   SetErrorCode(KO);
2185
2186   if (theShape.IsNull() || theCenter.IsNull()) return NULL;
2187
2188   TopoDS_Shape aShape  = theShape->GetValue();
2189   TopoDS_Shape aCenter = theCenter->GetValue();
2190
2191   if (aShape.IsNull() || aCenter.IsNull()) return NULL;
2192
2193   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2194   if ( !checkTypeShapesOn( aShapeType ))
2195     return NULL;
2196
2197   // Center of the sphere
2198   if (aCenter.ShapeType() != TopAbs_VERTEX) return NULL;
2199   gp_Pnt aLoc = BRep_Tool::Pnt(TopoDS::Vertex(aCenter));
2200
2201   gp_Ax3 anAx3 (aLoc, gp::DZ());
2202   Handle(Geom_SphericalSurface) aSphere =
2203     new Geom_SphericalSurface(anAx3, theRadius);
2204
2205   // Find objects
2206   TCollection_AsciiString anAsciiList;
2207   Handle(TColStd_HSequenceOfTransient) aSeq;
2208   aSeq = getShapesOnSurface( aSphere, theShape, aShapeType, theState, anAsciiList );
2209   if ( aSeq.IsNull() || aSeq->Length() == 0 )
2210     return NULL;
2211
2212   // Make a Python command
2213
2214   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2215   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2216
2217   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2218     << "] = geompy.GetShapesOnSphere(" << theShape << ", " << aShapeType
2219       << ", " << theCenter << ", " << theRadius << ", " << theState << ")";
2220
2221   SetErrorCode(OK);
2222   return aSeq;
2223 }
2224
2225 //=============================================================================
2226 /*!
2227  *  GetShapesOnPlaneIDs
2228  */
2229 //=============================================================================
2230 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnPlaneIDs
2231                                         (const Handle(GEOM_Object)& theShape,
2232                                          const Standard_Integer     theShapeType,
2233                                          const Handle(GEOM_Object)& theAx1,
2234                                          const GEOMAlgo_State       theState)
2235 {
2236   SetErrorCode(KO);
2237
2238   if (theShape.IsNull() || theAx1.IsNull()) return NULL;
2239
2240   TopoDS_Shape aShape = theShape->GetValue();
2241   TopoDS_Shape anAx1  = theAx1->GetValue();
2242
2243   if (aShape.IsNull() || anAx1.IsNull()) return NULL;
2244
2245   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2246   if ( !checkTypeShapesOn( aShapeType ))
2247     return NULL;
2248
2249   // Create plane
2250   Handle(Geom_Surface) aPlane = makePlane( anAx1 );
2251   if ( aPlane.IsNull() )
2252     return NULL;
2253
2254   // Find object IDs
2255   Handle(TColStd_HSequenceOfInteger) aSeq;
2256   aSeq = getShapesOnSurfaceIDs( aPlane, aShape, aShapeType, theState );
2257
2258   // The GetShapesOnPlaneIDs() doesn't change object so no new function is required.
2259   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theAx1)->GetLastFunction();
2260
2261   // Make a Python command
2262   GEOM::TPythonDump(aFunction, /*append=*/true)
2263     << "listShapesOnPlane = geompy.GetShapesOnPlaneIDs"
2264     << "(" << theShape << "," << aShapeType << "," << theAx1 << "," << theState << ")";
2265
2266   SetErrorCode(OK);
2267   return aSeq;
2268 }
2269
2270 //=============================================================================
2271 /*!
2272  *  GetShapesOnPlaneWithLocationIDs
2273  */
2274 //=============================================================================
2275 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnPlaneWithLocationIDs
2276                                         (const Handle(GEOM_Object)& theShape,
2277                                          const Standard_Integer     theShapeType,
2278                                          const Handle(GEOM_Object)& theAx1,
2279                                          const Handle(GEOM_Object)& thePnt,
2280                                          const GEOMAlgo_State       theState)
2281 {
2282   SetErrorCode(KO);
2283
2284   if (theShape.IsNull() || theAx1.IsNull() || thePnt.IsNull()) return NULL;
2285
2286   TopoDS_Shape aShape = theShape->GetValue();
2287   TopoDS_Shape anAx1  = theAx1->GetValue();
2288   TopoDS_Shape anPnt  = thePnt->GetValue();
2289
2290   if (aShape.IsNull() || anAx1.IsNull() || anPnt.IsNull()) return NULL;
2291
2292   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2293   if ( !checkTypeShapesOn( aShapeType ))
2294     return NULL;
2295
2296   // Create plane
2297   if (anAx1.ShapeType() != TopAbs_EDGE || anPnt.ShapeType() != TopAbs_VERTEX) return NULL;
2298   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
2299   TopoDS_Vertex V1, V2, V3;
2300   TopExp::Vertices(anEdge, V1, V2, Standard_True);
2301   if (V1.IsNull() || V2.IsNull()) {
2302     SetErrorCode("Bad edge given for the plane normal vector");
2303     return NULL;
2304   }
2305   V3 = TopoDS::Vertex(anPnt);
2306   if(V3.IsNull()) {
2307     SetErrorCode("Bad vertex given for the plane location");
2308       return NULL;
2309   }
2310   gp_Pnt aLoc = BRep_Tool::Pnt(V3);
2311   gp_Vec aVec(BRep_Tool::Pnt(V1),BRep_Tool::Pnt(V2));
2312   if (aVec.Magnitude() < Precision::Confusion()) {
2313     SetErrorCode("Vector with null magnitude given");
2314     return NULL;
2315   }
2316
2317   Handle(Geom_Surface) aPlane = new Geom_Plane(aLoc, aVec);
2318   if ( aPlane.IsNull() )
2319     return NULL;
2320
2321   // Find object IDs
2322   Handle(TColStd_HSequenceOfInteger) aSeq;
2323   aSeq = getShapesOnSurfaceIDs( aPlane, aShape, aShapeType, theState );
2324
2325   // The GetShapesOnPlaneIDs() doesn't change object so no new function is required.
2326   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theAx1)->GetLastFunction();
2327
2328   // Make a Python command
2329   GEOM::TPythonDump(aFunction, /*append=*/true)
2330     << "listShapesOnPlane = geompy.GetShapesOnPlaneWithLocationIDs"
2331     << "(" << theShape << ", " << aShapeType << ", " << theAx1 << ", "<< thePnt << ", "  << theState << ")";
2332
2333   SetErrorCode(OK);
2334   return aSeq;
2335 }
2336
2337 //=============================================================================
2338 /*!
2339  *  GetShapesOnCylinderIDs
2340  */
2341 //=============================================================================
2342 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnCylinderIDs
2343                                           (const Handle(GEOM_Object)& theShape,
2344                                            const Standard_Integer     theShapeType,
2345                                            const Handle(GEOM_Object)& theAxis,
2346                                            const Standard_Real        theRadius,
2347                                            const GEOMAlgo_State       theState)
2348 {
2349   SetErrorCode(KO);
2350
2351   if (theShape.IsNull() || theAxis.IsNull()) return NULL;
2352
2353   TopoDS_Shape aShape = theShape->GetValue();
2354   TopoDS_Shape anAxis = theAxis->GetValue();
2355
2356   if (aShape.IsNull() || anAxis.IsNull()) return NULL;
2357
2358   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2359   if ( !checkTypeShapesOn( aShapeType ))
2360     return NULL;
2361
2362   // Create a cylinder surface
2363   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
2364   if ( aCylinder.IsNull() )
2365     return NULL;
2366
2367   // Find object IDs
2368   Handle(TColStd_HSequenceOfInteger) aSeq;
2369   aSeq = getShapesOnSurfaceIDs( aCylinder, aShape, aShapeType, theState );
2370
2371   // The GetShapesOnCylinder() doesn't change object so no new function is required.
2372   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theAxis)->GetLastFunction();
2373
2374   // Make a Python command
2375   GEOM::TPythonDump(aFunction, /*append=*/true)
2376     << "listShapesOnCylinder = geompy.GetShapesOnCylinderIDs"
2377     << "(" << theShape << ", " << aShapeType << ", " << theAxis << ", "
2378     << theRadius << ", " << theState << ")";
2379
2380   SetErrorCode(OK);
2381   return aSeq;
2382 }
2383
2384 //=============================================================================
2385 /*!
2386  *  GetShapesOnSphereIDs
2387  */
2388 //=============================================================================
2389 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnSphereIDs
2390                                           (const Handle(GEOM_Object)& theShape,
2391                                            const Standard_Integer     theShapeType,
2392                                            const Handle(GEOM_Object)& theCenter,
2393                                            const Standard_Real        theRadius,
2394                                            const GEOMAlgo_State       theState)
2395 {
2396   SetErrorCode(KO);
2397
2398   if (theShape.IsNull() || theCenter.IsNull()) return NULL;
2399
2400   TopoDS_Shape aShape  = theShape->GetValue();
2401   TopoDS_Shape aCenter = theCenter->GetValue();
2402
2403   if (aShape.IsNull() || aCenter.IsNull()) return NULL;
2404
2405   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2406   if ( !checkTypeShapesOn( aShapeType ))
2407     return NULL;
2408
2409   // Center of the sphere
2410   if (aCenter.ShapeType() != TopAbs_VERTEX) return NULL;
2411   gp_Pnt aLoc = BRep_Tool::Pnt(TopoDS::Vertex(aCenter));
2412
2413   gp_Ax3 anAx3 (aLoc, gp::DZ());
2414   Handle(Geom_SphericalSurface) aSphere =
2415     new Geom_SphericalSurface(anAx3, theRadius);
2416
2417   // Find object IDs
2418   Handle(TColStd_HSequenceOfInteger) aSeq;
2419   aSeq = getShapesOnSurfaceIDs( aSphere, aShape, aShapeType, theState );
2420
2421   // The GetShapesOnSphere() doesn't change object so no new function is required.
2422   Handle(GEOM_Function) aFunction = GEOM::GetCreatedLast(theShape,theCenter)->GetLastFunction();
2423
2424   // Make a Python command
2425   GEOM::TPythonDump(aFunction, /*append=*/true)
2426     << "listShapesOnCylinder = geompy.GetShapesOnCylinderIDs"
2427     << "(" << theShape << ", " << aShapeType << ", " << theCenter << ", "
2428     << theRadius << ", " << theState << ")";
2429
2430   SetErrorCode(OK);
2431   return aSeq;
2432 }
2433
2434 //=======================================================================
2435 //function : getShapesOnQuadrangleIDs
2436   /*!
2437    * \brief Find IDs of subshapes complying with given status about quadrangle
2438     * \param theShape - the shape to explore
2439     * \param theShapeType - type of subshape of theShape
2440     * \param theTopLeftPoint - top left quadrangle corner
2441     * \param theTopRigthPoint - top right quadrangle corner
2442     * \param theBottomLeftPoint - bottom left quadrangle corner
2443     * \param theBottomRigthPoint - bottom right quadrangle corner
2444     * \param theState - required state
2445     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2446    */
2447 //=======================================================================
2448
2449 Handle(TColStd_HSequenceOfInteger)
2450   GEOMImpl_IShapesOperations::getShapesOnQuadrangleIDs (const Handle(GEOM_Object)& theShape,
2451                                                         const Standard_Integer     theShapeType,
2452                                                         const Handle(GEOM_Object)& theTopLeftPoint,
2453                                                         const Handle(GEOM_Object)& theTopRigthPoint,
2454                                                         const Handle(GEOM_Object)& theBottomLeftPoint,
2455                                                         const Handle(GEOM_Object)& theBottomRigthPoint,
2456                                                         const GEOMAlgo_State       theState)
2457 {
2458   SetErrorCode(KO);
2459
2460   if ( theShape.IsNull() ||
2461        theTopLeftPoint.IsNull() ||
2462        theTopRigthPoint.IsNull() ||
2463        theBottomLeftPoint.IsNull() ||
2464        theBottomRigthPoint.IsNull() )
2465     return NULL;
2466
2467   TopoDS_Shape aShape = theShape->GetValue();
2468   TopoDS_Shape aTL = theTopLeftPoint->GetValue();
2469   TopoDS_Shape aTR = theTopRigthPoint->GetValue();
2470   TopoDS_Shape aBL = theBottomLeftPoint->GetValue();
2471   TopoDS_Shape aBR = theBottomRigthPoint->GetValue();
2472
2473   if (aShape.IsNull() ||
2474       aTL.IsNull() ||
2475       aTR.IsNull() ||
2476       aBL.IsNull() ||
2477       aBR.IsNull() ||
2478       aTL.ShapeType() != TopAbs_VERTEX ||
2479       aTR.ShapeType() != TopAbs_VERTEX ||
2480       aBL.ShapeType() != TopAbs_VERTEX ||
2481       aBR.ShapeType() != TopAbs_VERTEX )
2482     return NULL;
2483
2484   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
2485   if ( !checkTypeShapesOn( aShapeType ))
2486     return NULL;
2487
2488   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
2489
2490   // Check presence of triangulation, build if need
2491   if (!CheckTriangulation(aShape)) {
2492     SetErrorCode("Cannot build triangulation on the shape");
2493     return aSeqOfIDs;
2494   }
2495
2496   // Call algo
2497   gp_Pnt aPntTL = BRep_Tool::Pnt(TopoDS::Vertex(aTL));
2498   gp_Pnt aPntTR = BRep_Tool::Pnt(TopoDS::Vertex(aTR));
2499   gp_Pnt aPntBL = BRep_Tool::Pnt(TopoDS::Vertex(aBL));
2500   gp_Pnt aPntBR = BRep_Tool::Pnt(TopoDS::Vertex(aBR));
2501
2502   GEOMAlgo_FinderShapeOnQuad aFinder( aPntTL, aPntTR, aPntBL, aPntBR );
2503   Standard_Real aTol = 0.0001; // default value
2504
2505   aFinder.SetShape(aShape);
2506   aFinder.SetTolerance(aTol);
2507   //aFinder.SetSurface(theSurface);
2508   aFinder.SetShapeType(aShapeType);
2509   aFinder.SetState(theState);
2510
2511   // Sets the minimal number of inner points for the faces that do not have own
2512   // inner points at all (for e.g. rectangular planar faces have just 2 triangles).
2513   // Default value=3
2514   aFinder.SetNbPntsMin(3);
2515   // Sets the maximal number of inner points for edges or faces.
2516   // It is usefull for the cases when this number is very big (e.g =2000) to improve
2517   // the performance. If this value =0, all inner points will be taken into account.
2518   // Default value=0
2519   aFinder.SetNbPntsMax(100);
2520
2521   aFinder.Perform();
2522
2523   // Interprete results
2524   Standard_Integer iErr = aFinder.ErrorStatus();
2525   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
2526   if (iErr) {
2527     MESSAGE(" iErr : " << iErr);
2528     TCollection_AsciiString aMsg (" iErr : ");
2529     aMsg += TCollection_AsciiString(iErr);
2530     SetErrorCode(aMsg);
2531     return aSeqOfIDs;
2532   }
2533   Standard_Integer iWrn = aFinder.WarningStatus();
2534   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
2535   if (iWrn) {
2536     MESSAGE(" *** iWrn : " << iWrn);
2537   }
2538
2539   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
2540
2541   if (listSS.Extent() < 1) {
2542     //SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
2543     SetErrorCode(NOT_FOUND_ANY); // NPAL18017
2544     return aSeqOfIDs;
2545   }
2546
2547   // Fill sequence of object IDs
2548   aSeqOfIDs = new TColStd_HSequenceOfInteger;
2549
2550   TopTools_IndexedMapOfShape anIndices;
2551   TopExp::MapShapes(aShape, anIndices);
2552
2553   TopTools_ListIteratorOfListOfShape itSub (listSS);
2554   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
2555     int id = anIndices.FindIndex(itSub.Value());
2556     aSeqOfIDs->Append(id);
2557   }
2558   return aSeqOfIDs;
2559 }
2560
2561 //=======================================================================
2562 //function : GetShapesOnQuadrangle
2563   /*!
2564    * \brief Find subshapes complying with given status about quadrangle
2565     * \param theShape - the shape to explore
2566     * \param theShapeType - type of subshape of theShape
2567     * \param theTopLeftPoint - top left quadrangle corner
2568     * \param theTopRigthPoint - top right quadrangle corner
2569     * \param theBottomLeftPoint - bottom left quadrangle corner
2570     * \param theBottomRigthPoint - bottom right quadrangle corner
2571     * \param theState - required state
2572     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2573    */
2574 //=======================================================================
2575
2576 Handle(TColStd_HSequenceOfTransient)
2577     GEOMImpl_IShapesOperations::GetShapesOnQuadrangle (const Handle(GEOM_Object)& theShape,
2578                                                        const Standard_Integer     theShapeType,
2579                                                        const Handle(GEOM_Object)& theTopLeftPoint,
2580                                                        const Handle(GEOM_Object)& theTopRigthPoint,
2581                                                        const Handle(GEOM_Object)& theBottomLeftPoint,
2582                                                        const Handle(GEOM_Object)& theBottomRigthPoint,
2583                                                        const GEOMAlgo_State       theState)
2584 {
2585   // Find indices
2586   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2587     getShapesOnQuadrangleIDs( theShape,
2588                               theShapeType,
2589                               theTopLeftPoint,
2590                               theTopRigthPoint,
2591                               theBottomLeftPoint,
2592                               theBottomRigthPoint,
2593                               theState);
2594   if ( aSeqOfIDs.IsNull() || aSeqOfIDs->IsEmpty() )
2595     return NULL;
2596
2597   // Find objects by indices
2598   TCollection_AsciiString anAsciiList;
2599   Handle(TColStd_HSequenceOfTransient) aSeq;
2600   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
2601   if ( aSeq.IsNull() || aSeq->IsEmpty() )
2602     return NULL;
2603
2604   // Make a Python command
2605
2606   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
2607   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
2608
2609   GEOM::TPythonDump(aFunction)
2610     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnQuadrangle("
2611     << theShape << ", "
2612     << TopAbs_ShapeEnum(theShapeType) << ", "
2613     << theTopLeftPoint << ", "
2614     << theTopRigthPoint << ", "
2615     << theBottomLeftPoint << ", "
2616     << theBottomRigthPoint << ", "
2617     << theState << ")";
2618
2619   SetErrorCode(OK);
2620   return aSeq;
2621 }
2622
2623 //=======================================================================
2624 //function : GetShapesOnQuadrangleIDs
2625   /*!
2626    * \brief Find IDs of subshapes complying with given status about quadrangle
2627     * \param theShape - the shape to explore
2628     * \param theShapeType - type of subshape of theShape
2629     * \param theTopLeftPoint - top left quadrangle corner
2630     * \param theTopRigthPoint - top right quadrangle corner
2631     * \param theBottomLeftPoint - bottom left quadrangle corner
2632     * \param theBottomRigthPoint - bottom right quadrangle corner
2633     * \param theState - required state
2634     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
2635    */
2636 //=======================================================================
2637
2638 Handle(TColStd_HSequenceOfInteger)
2639   GEOMImpl_IShapesOperations::GetShapesOnQuadrangleIDs (const Handle(GEOM_Object)& theShape,
2640                                                         const Standard_Integer     theShapeType,
2641                                                         const Handle(GEOM_Object)& theTopLeftPoint,
2642                                                         const Handle(GEOM_Object)& theTopRigthPoint,
2643                                                         const Handle(GEOM_Object)& theBottomLeftPoint,
2644                                                         const Handle(GEOM_Object)& theBottomRigthPoint,
2645                                                         const GEOMAlgo_State       theState)
2646 {
2647   // Find indices
2648   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
2649     getShapesOnQuadrangleIDs( theShape,
2650                               theShapeType,
2651                               theTopLeftPoint,
2652                               theTopRigthPoint,
2653                               theBottomLeftPoint,
2654                               theBottomRigthPoint,
2655                               theState);
2656   if ( aSeqOfIDs.IsNull() || aSeqOfIDs->IsEmpty() )
2657     return NULL;
2658
2659   // Make a Python command
2660
2661   // The GetShapesOnCylinder() doesn't change object so no new function is required.
2662   Handle(GEOM_Object) lastObj = GEOM::GetCreatedLast(theShape,theTopLeftPoint);
2663   lastObj = GEOM::GetCreatedLast(lastObj,theTopRigthPoint);
2664   lastObj = GEOM::GetCreatedLast(lastObj,theBottomRigthPoint);
2665   lastObj = GEOM::GetCreatedLast(lastObj,theBottomLeftPoint);
2666   Handle(GEOM_Function) aFunction = lastObj->GetLastFunction();
2667
2668   GEOM::TPythonDump(aFunction, /*append=*/true)
2669     << "listShapesOnQuadrangle = geompy.GetShapesOnQuadrangleIDs("
2670     << theShape << ", "
2671     << TopAbs_ShapeEnum(theShapeType) << ", "
2672     << theTopLeftPoint << ", "
2673     << theTopRigthPoint << ", "
2674     << theBottomLeftPoint << ", "
2675     << theBottomRigthPoint << ", "
2676     << theState << ")";
2677
2678   SetErrorCode(OK);
2679   return aSeqOfIDs;
2680 }
2681
2682
2683 //=============================================================================
2684 /*!
2685  *  GetInPlaceOfShape
2686  */
2687 //=============================================================================
2688 static bool GetInPlaceOfShape (const Handle(GEOM_Function)& theWhereFunction,
2689                                const TopTools_IndexedMapOfShape& theWhereIndices,
2690                                const TopoDS_Shape& theWhat,
2691                                TColStd_ListOfInteger& theModifiedList)
2692 {
2693   if (theWhereFunction.IsNull() || theWhat.IsNull()) return false;
2694
2695   if (theWhereIndices.Contains(theWhat)) {
2696     // entity was not changed by the operation
2697     Standard_Integer aWhatIndex = theWhereIndices.FindIndex(theWhat);
2698     theModifiedList.Append(aWhatIndex);
2699     return true;
2700   }
2701
2702   // try to find in history
2703   TDF_Label aHistoryLabel = theWhereFunction->GetHistoryEntry(Standard_False);
2704
2705   // search in history for all argument shapes
2706   Standard_Boolean isFound = Standard_False;
2707   Standard_Boolean isGood = Standard_False;
2708
2709   TDF_LabelSequence aLabelSeq;
2710   theWhereFunction->GetDependency(aLabelSeq);
2711   Standard_Integer nbArg = aLabelSeq.Length();
2712
2713   for (Standard_Integer iarg = 1; iarg <= nbArg && !isFound; iarg++) {
2714
2715     TDF_Label anArgumentRefLabel = aLabelSeq.Value(iarg);
2716
2717     Handle(GEOM_Object) anArgumentObject = GEOM_Object::GetReferencedObject(anArgumentRefLabel);
2718     TopoDS_Shape anArgumentShape = anArgumentObject->GetValue();
2719
2720     TopTools_IndexedMapOfShape anArgumentIndices;
2721     TopExp::MapShapes(anArgumentShape, anArgumentIndices);
2722
2723     if (anArgumentIndices.Contains(theWhat)) {
2724       isFound = Standard_True;
2725       Standard_Integer aWhatIndex = anArgumentIndices.FindIndex(theWhat);
2726
2727       // Find corresponding label in history
2728       TDF_Label anArgumentHistoryLabel =
2729         theWhereFunction->GetArgumentHistoryEntry(anArgumentRefLabel, Standard_False);
2730       if (anArgumentHistoryLabel.IsNull()) {
2731         // Lost History of operation argument. Possibly, all its entities was removed.
2732         isGood = Standard_True;
2733       }
2734       else {
2735         TDF_Label aWhatHistoryLabel = anArgumentHistoryLabel.FindChild(aWhatIndex, Standard_False);
2736
2737         if (aWhatHistoryLabel.IsNull()) {
2738           // Removed entity ? Compound ? Compsolid ? Shell ? Wire
2739           isGood = Standard_False;
2740         } else {
2741           Handle(TDataStd_IntegerArray) anIntegerArray;
2742           if (!aWhatHistoryLabel.FindAttribute(TDataStd_IntegerArray::GetID(), anIntegerArray)) {
2743             //Error: Empty modifications history for the sought shape.
2744             isGood = Standard_False;
2745           }
2746           else {
2747             isGood = Standard_True;
2748             Standard_Integer imod, aModifLen = anIntegerArray->Array()->Length();
2749             for (imod = 1; imod <= aModifLen; imod++) {
2750               theModifiedList.Append(anIntegerArray->Array()->Value(imod));
2751             }
2752           }
2753         }
2754       }
2755     }
2756   }
2757
2758   isFound = isGood;
2759
2760   if (!isFound) {
2761     // try compound/compsolid/shell/wire element by element
2762     bool isFoundAny = false;
2763     TopTools_MapOfShape mapShape;
2764
2765     if (theWhat.ShapeType() == TopAbs_COMPOUND ||
2766         theWhat.ShapeType() == TopAbs_COMPSOLID) {
2767       // recursive processing of compound/compsolid
2768       TopoDS_Iterator anIt (theWhat, Standard_True, Standard_True);
2769       for (; anIt.More(); anIt.Next()) {
2770         if (mapShape.Add(anIt.Value())) {
2771           TopoDS_Shape curWhat = anIt.Value();
2772           isFoundAny = GetInPlaceOfShape(theWhereFunction, theWhereIndices, curWhat, theModifiedList);
2773           if (isFoundAny) isFound = Standard_True;
2774         }
2775       }
2776     }
2777     else if (theWhat.ShapeType() == TopAbs_SHELL) {
2778       // try to replace a shell by its faces images
2779       TopExp_Explorer anExp (theWhat, TopAbs_FACE);
2780       for (; anExp.More(); anExp.Next()) {
2781         if (mapShape.Add(anExp.Current())) {
2782           TopoDS_Shape curWhat = anExp.Current();
2783           isFoundAny = GetInPlaceOfShape(theWhereFunction, theWhereIndices, curWhat, theModifiedList);
2784           if (isFoundAny) isFound = Standard_True;
2785         }
2786       }
2787     }
2788     else if (theWhat.ShapeType() == TopAbs_WIRE) {
2789       // try to replace a wire by its edges images
2790       TopExp_Explorer anExp (theWhat, TopAbs_EDGE);
2791       for (; anExp.More(); anExp.Next()) {
2792         if (mapShape.Add(anExp.Current())) {
2793           TopoDS_Shape curWhat = anExp.Current();
2794           isFoundAny = GetInPlaceOfShape(theWhereFunction, theWhereIndices, curWhat, theModifiedList);
2795           if (isFoundAny) isFound = Standard_True;
2796         }
2797       }
2798     }
2799     else {
2800       // Removed entity
2801     }
2802   }
2803
2804   return isFound;
2805 }
2806
2807 //=============================================================================
2808 /*!
2809  *  GetShapeProperties
2810  */
2811 //=============================================================================
2812
2813 void GEOMImpl_IShapesOperations::GetShapeProperties( const TopoDS_Shape aShape, Standard_Real tab[],
2814                                                      gp_Pnt & aVertex )
2815 {
2816   GProp_GProps theProps;
2817   gp_Pnt aCenterMass;
2818   TopoDS_Shape aPntShape;
2819   Standard_Real aShapeSize;
2820
2821   if      (aShape.ShapeType() == TopAbs_EDGE) BRepGProp::LinearProperties(aShape,  theProps);
2822   else if (aShape.ShapeType() == TopAbs_FACE) BRepGProp::SurfaceProperties(aShape, theProps);
2823   else                                        BRepGProp::VolumeProperties(aShape,  theProps);
2824
2825   aCenterMass = theProps.CentreOfMass();
2826   aShapeSize  = theProps.Mass();
2827
2828   aPntShape = BRepBuilderAPI_MakeVertex(aCenterMass).Shape();
2829   aVertex   = BRep_Tool::Pnt( TopoDS::Vertex( aPntShape ) );
2830   tab[0] = aVertex.X();
2831   tab[1] = aVertex.Y();
2832   tab[2] = aVertex.Z();
2833   tab[3] = aShapeSize;
2834   return;
2835 }
2836
2837 //=============================================================================
2838 /*!
2839  *  GetInPlace
2840  */
2841 //=============================================================================
2842 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetInPlace (Handle(GEOM_Object) theShapeWhere,
2843                                                             Handle(GEOM_Object) theShapeWhat)
2844 {
2845   SetErrorCode(KO);
2846
2847   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
2848
2849   TopoDS_Shape aWhere = theShapeWhere->GetValue();
2850   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
2851   TopoDS_Shape aPntShape;
2852   TopoDS_Vertex aVertex;
2853
2854   if (aWhere.IsNull() || aWhat.IsNull()) {
2855     SetErrorCode("Error: aWhere and aWhat TopoDS_Shape are Null.");
2856     return NULL;
2857   }
2858
2859   Handle(GEOM_Function) aWhereFunction = theShapeWhere->GetLastFunction();
2860   if (aWhereFunction.IsNull()) {
2861     SetErrorCode("Error: aWhereFunction is Null.");
2862     return NULL;
2863   }
2864
2865   TopTools_IndexedMapOfShape aWhereIndices;
2866   TopExp::MapShapes(aWhere, aWhereIndices);
2867
2868   TColStd_ListOfInteger aModifiedList;
2869   Standard_Integer aWhereIndex;
2870   Handle(TColStd_HArray1OfInteger) aModifiedArray;
2871   Handle(GEOM_Object) aResult;
2872
2873   bool isFound = false;
2874   Standard_Integer iType = TopAbs_SOLID;
2875   Standard_Integer compType = TopAbs_SOLID;
2876   Standard_Real    aWhat_Mass = 0., aWhere_Mass = 0.;
2877   Standard_Real    tab_aWhat[4],    tab_aWhere[4];
2878   Standard_Real    dl_l = 1e-3;
2879   Standard_Real    min_l, Tol_1D, Tol_2D, Tol_3D, Tol_Mass;
2880   Standard_Real    aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
2881   Bnd_Box          BoundingBox;
2882   gp_Pnt           aPnt, aPnt_aWhat, tab_Pnt[2];
2883   GProp_GProps     aProps;
2884
2885   // Find the iType of the aWhat shape
2886   if      ( aWhat.ShapeType() == TopAbs_EDGE  || aWhat.ShapeType() == TopAbs_WIRE )      iType = TopAbs_EDGE;
2887   else if ( aWhat.ShapeType() == TopAbs_FACE  || aWhat.ShapeType() == TopAbs_SHELL )     iType = TopAbs_FACE;
2888   else if ( aWhat.ShapeType() == TopAbs_SOLID || aWhat.ShapeType() == TopAbs_COMPSOLID ) iType = TopAbs_SOLID;
2889   else if ( aWhat.ShapeType() == TopAbs_COMPOUND ) {
2890     // Only the iType of the first shape in the compound is taken into account
2891     TopoDS_Iterator It (aWhat, Standard_True, Standard_True);
2892     compType = It.Value().ShapeType();
2893     if      ( compType == TopAbs_EDGE  || compType == TopAbs_WIRE )     iType = TopAbs_EDGE;
2894     else if ( compType == TopAbs_FACE  || compType == TopAbs_SHELL)     iType = TopAbs_FACE;
2895     else if ( compType == TopAbs_SOLID || compType == TopAbs_COMPSOLID) iType = TopAbs_SOLID;
2896   }
2897   else {
2898     SetErrorCode("Error: An attempt to extract a shape of not supported type.");
2899     return NULL;
2900   }
2901
2902   TopExp_Explorer Exp_aWhat( aWhat,   TopAbs_ShapeEnum( iType ) );
2903   TopExp_Explorer Exp_aWhere( aWhere, TopAbs_ShapeEnum( iType ) );
2904   TopExp_Explorer Exp_Edge( aWhere,   TopAbs_EDGE );
2905
2906   // Find the shortest edge in theShapeWhere shape
2907   BRepBndLib::Add(aWhere, BoundingBox);
2908   BoundingBox.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
2909   min_l = fabs(aXmax - aXmin);
2910   if( min_l < fabs(aYmax - aYmin) ) min_l = fabs(aYmax - aYmin);
2911   if( min_l < fabs(aZmax - aZmin) ) min_l = fabs(aZmax - aZmin);
2912   min_l /= dl_l;
2913   for ( Standard_Integer nbEdge = 0; Exp_Edge.More(); Exp_Edge.Next(), nbEdge++ ) {
2914     TopExp_Explorer Exp_Vertex( Exp_Edge.Current(), TopAbs_VERTEX);
2915     for ( Standard_Integer nbVertex = 0; Exp_Vertex.More(); Exp_Vertex.Next(), nbVertex++ ) {
2916       aPnt = BRep_Tool::Pnt( TopoDS::Vertex( Exp_Vertex.Current() ) );
2917       tab_Pnt[nbVertex] = aPnt;
2918     }
2919     if ( ! tab_Pnt[0].IsEqual(tab_Pnt[1], dl_l) ) {
2920       BRepGProp::LinearProperties(Exp_Edge.Current(), aProps);
2921       if ( aProps.Mass() < min_l ) min_l = aProps.Mass();
2922     }
2923   }
2924
2925   // Compute tolerances
2926   Tol_1D = dl_l * min_l;
2927   Tol_2D = dl_l * ( min_l * min_l) * ( 2. + dl_l);
2928   Tol_3D = dl_l * ( min_l * min_l * min_l ) * ( 3. + (3 * dl_l) + (dl_l * dl_l) );
2929
2930   if (Tol_1D < Precision::Confusion()) Tol_1D = Precision::Confusion();
2931   if (Tol_2D < Precision::Confusion()) Tol_2D = Precision::Confusion();
2932   if (Tol_3D < Precision::Confusion()) Tol_3D = Precision::Confusion();
2933
2934   Tol_Mass = Tol_3D;
2935   if      ( iType == TopAbs_EDGE ) Tol_Mass = Tol_1D;
2936   else if ( iType == TopAbs_FACE ) Tol_Mass = Tol_2D;
2937
2938   // Compute the ShapeWhat Mass
2939   for ( ; Exp_aWhat.More(); Exp_aWhat.Next() ) {
2940     if      ( iType == TopAbs_EDGE ) BRepGProp::LinearProperties(Exp_aWhat.Current(),  aProps);
2941     else if ( iType == TopAbs_FACE ) BRepGProp::SurfaceProperties(Exp_aWhat.Current(), aProps);
2942     else                             BRepGProp::VolumeProperties(Exp_aWhat.Current(),  aProps);
2943     aWhat_Mass += aProps.Mass();
2944   }
2945
2946   // Searching for the sub-shapes inside the ShapeWhere shape
2947   TopTools_MapOfShape map_aWhere;
2948   for ( Exp_aWhere.ReInit(); Exp_aWhere.More(); Exp_aWhere.Next() ) {
2949     if (!map_aWhere.Add(Exp_aWhere.Current()))
2950       continue; // skip repeated shape to avoid mass addition
2951     GetShapeProperties( Exp_aWhere.Current(), tab_aWhere, aPnt );
2952     for ( Exp_aWhat.ReInit(); Exp_aWhat.More(); Exp_aWhat.Next() ) {
2953       GetShapeProperties( Exp_aWhat.Current(), tab_aWhat, aPnt_aWhat );
2954       if ( fabs(tab_aWhat[3] - tab_aWhere[3]) <= Tol_Mass && aPnt_aWhat.Distance(aPnt) <= Tol_1D )
2955         isFound = true;
2956       else {
2957         if ( (tab_aWhat[3] - tab_aWhere[3]) > Tol_Mass ) {
2958           aPntShape = BRepBuilderAPI_MakeVertex( aPnt ).Shape();
2959           aVertex   = TopoDS::Vertex( aPntShape );
2960           BRepExtrema_DistShapeShape aWhereDistance ( aVertex, Exp_aWhere.Current() );
2961           BRepExtrema_DistShapeShape aWhatDistance  ( aVertex, Exp_aWhat.Current() );
2962           if ( fabs(aWhereDistance.Value() - aWhatDistance.Value()) <= Tol_1D )
2963             isFound = true;
2964         }
2965       }
2966       if ( isFound ) {
2967         aWhereIndex = aWhereIndices.FindIndex(Exp_aWhere.Current());
2968         aModifiedList.Append(aWhereIndex);
2969         aWhere_Mass += tab_aWhere[3];
2970         isFound = false;
2971         break;
2972       }
2973     }
2974     if ( fabs( aWhat_Mass - aWhere_Mass ) <= Tol_Mass ) break;
2975   }
2976
2977   if (aModifiedList.Extent() == 0) { // Not found any Results
2978     SetErrorCode(NOT_FOUND_ANY);
2979     return NULL;
2980   }
2981
2982   aModifiedArray = new TColStd_HArray1OfInteger (1, aModifiedList.Extent());
2983   TColStd_ListIteratorOfListOfInteger anIterModif (aModifiedList);
2984   for (Standard_Integer imod = 1; anIterModif.More(); anIterModif.Next(), imod++)
2985     aModifiedArray->SetValue(imod, anIterModif.Value());
2986
2987   //Add a new object
2988   aResult = GetEngine()->AddSubShape(theShapeWhere, aModifiedArray);
2989   if (aResult.IsNull()) {
2990     SetErrorCode("Error in algorithm: result found, but cannot be returned.");
2991     return NULL;
2992   }
2993
2994   if (aModifiedArray->Length() > 1) {
2995     //Set a GROUP type
2996     aResult->SetType(GEOM_GROUP);
2997
2998     //Set a sub shape type
2999     TopoDS_Shape aFirstFound = aWhereIndices.FindKey(aModifiedArray->Value(1));
3000     TopAbs_ShapeEnum aShapeType = aFirstFound.ShapeType();
3001
3002     TDF_Label aFreeLabel = aResult->GetFreeLabel();
3003     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)aShapeType);
3004   }
3005
3006   //Make a Python command
3007   Handle(GEOM_Function) aFunction = aResult->GetFunction(1);
3008
3009   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetInPlace("
3010     << theShapeWhere << ", " << theShapeWhat << ")";
3011
3012   SetErrorCode(OK);
3013   return aResult;
3014 }
3015
3016 //=======================================================================
3017 //function : GetInPlaceByHistory
3018 //purpose  :
3019 //=======================================================================
3020 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetInPlaceByHistory
3021                                           (Handle(GEOM_Object) theShapeWhere,
3022                                            Handle(GEOM_Object) theShapeWhat)
3023 {
3024   SetErrorCode(KO);
3025
3026   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
3027
3028   TopoDS_Shape aWhere = theShapeWhere->GetValue();
3029   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
3030
3031   if (aWhere.IsNull() || aWhat.IsNull()) return NULL;
3032
3033   Handle(GEOM_Function) aWhereFunction = theShapeWhere->GetLastFunction();
3034   if (aWhereFunction.IsNull()) return NULL;
3035
3036   //Fill array of indices
3037   TopTools_IndexedMapOfShape aWhereIndices;
3038   TopExp::MapShapes(aWhere, aWhereIndices);
3039
3040   // process shape
3041   TColStd_ListOfInteger aModifiedList;
3042   bool isFound = GetInPlaceOfShape(aWhereFunction, aWhereIndices, aWhat, aModifiedList);
3043
3044   if (!isFound || aModifiedList.Extent() < 1) {
3045     SetErrorCode("Error: No history found for the sought shape or its sub-shapes.");
3046     return NULL;
3047   }
3048
3049   Handle(TColStd_HArray1OfInteger) aModifiedArray =
3050     new TColStd_HArray1OfInteger (1, aModifiedList.Extent());
3051   TColStd_ListIteratorOfListOfInteger anIterModif (aModifiedList);
3052   for (Standard_Integer imod = 1; anIterModif.More(); anIterModif.Next(), imod++) {
3053     aModifiedArray->SetValue(imod, anIterModif.Value());
3054   }
3055
3056   //Add a new object
3057   Handle(GEOM_Object) aResult = GetEngine()->AddSubShape(theShapeWhere, aModifiedArray);
3058   if (aResult.IsNull()) {
3059     SetErrorCode("Error in algorithm: result found, but cannot be returned.");
3060     return NULL;
3061   }
3062
3063   if (aModifiedArray->Length() > 1) {
3064     //Set a GROUP type
3065     aResult->SetType(GEOM_GROUP);
3066
3067     //Set a sub shape type
3068     TopoDS_Shape aFirstFound = aWhereIndices.FindKey(aModifiedArray->Value(1));
3069     TopAbs_ShapeEnum aShapeType = aFirstFound.ShapeType();
3070
3071     TDF_Label aFreeLabel = aResult->GetFreeLabel();
3072     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)aShapeType);
3073   }
3074
3075   //Make a Python command
3076   Handle(GEOM_Function) aFunction = aResult->GetFunction(1);
3077
3078   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetInPlaceByHistory("
3079     << theShapeWhere << ", " << theShapeWhat << ")";
3080
3081   SetErrorCode(OK);
3082   return aResult;
3083 }
3084
3085 //=======================================================================
3086 //function : SortShapes
3087 //purpose  :
3088 //=======================================================================
3089 void GEOMImpl_IShapesOperations::SortShapes(TopTools_ListOfShape& SL)
3090 {
3091   Standard_Integer MaxShapes = SL.Extent();
3092   TopTools_Array1OfShape  aShapes (1,MaxShapes);
3093   TColStd_Array1OfInteger OrderInd(1,MaxShapes);
3094   TColStd_Array1OfReal    MidXYZ  (1,MaxShapes); //X,Y,Z;
3095   TColStd_Array1OfReal    Length  (1,MaxShapes); //X,Y,Z;
3096
3097   // Computing of CentreOfMass
3098   Standard_Integer Index;
3099   GProp_GProps GPr;
3100   gp_Pnt GPoint;
3101   TopTools_ListIteratorOfListOfShape it(SL);
3102   for (Index=1;  it.More();  Index++)
3103   {
3104     TopoDS_Shape S = it.Value();
3105     SL.Remove( it ); // == it.Next()
3106     aShapes(Index) = S;
3107     OrderInd.SetValue (Index, Index);
3108     if (S.ShapeType() == TopAbs_VERTEX)
3109     {
3110       GPoint = BRep_Tool::Pnt( TopoDS::Vertex( S ));
3111       Length.SetValue( Index, (Standard_Real) S.Orientation());
3112     }
3113     else
3114     {
3115       BRepGProp::LinearProperties (S, GPr);
3116       GPoint = GPr.CentreOfMass();
3117       Length.SetValue( Index, GPr.Mass() );
3118     }
3119     MidXYZ.SetValue(Index,
3120                     GPoint.X()*999 + GPoint.Y()*99 + GPoint.Z()*0.9);
3121     //cout << Index << " L: " << Length(Index) << "CG: " << MidXYZ(Index) << endl;
3122   }
3123
3124   // Sorting
3125   Standard_Integer aTemp;
3126   Standard_Boolean exchange, Sort = Standard_True;
3127   Standard_Real    tol = Precision::Confusion();
3128   while (Sort)
3129   {
3130     Sort = Standard_False;
3131     for (Index=1; Index < MaxShapes; Index++)
3132     {
3133       exchange = Standard_False;
3134       Standard_Real dMidXYZ = MidXYZ(OrderInd(Index)) - MidXYZ(OrderInd(Index+1));
3135       Standard_Real dLength = Length(OrderInd(Index)) - Length(OrderInd(Index+1));
3136       if ( dMidXYZ >= tol ) {
3137 //         cout << "MidXYZ: " << MidXYZ(OrderInd(Index))<< " > " <<MidXYZ(OrderInd(Index+1))
3138 //              << " d: " << dMidXYZ << endl;
3139         exchange = Standard_True;
3140       }
3141       else if ( Abs(dMidXYZ) < tol && dLength >= tol ) {
3142 //         cout << "Length: " << Length(OrderInd(Index))<< " > " <<Length(OrderInd(Index+1))
3143 //              << " d: " << dLength << endl;
3144         exchange = Standard_True;
3145       }
3146       else if ( Abs(dMidXYZ) < tol && Abs(dLength) < tol &&
3147                 aShapes(OrderInd(Index)).ShapeType() <= TopAbs_FACE) {
3148         // PAL17233
3149         // equal values possible on shapes such as two halves of a sphere and
3150         // a membrane inside the sphere
3151         Bnd_Box box1,box2;
3152         BRepBndLib::Add( aShapes( OrderInd(Index) ), box1 );
3153         if ( box1.IsVoid() ) continue;
3154         BRepBndLib::Add( aShapes( OrderInd(Index+1) ), box2 );
3155         Standard_Real dSquareExtent = box1.SquareExtent() - box2.SquareExtent();
3156         if ( dSquareExtent >= tol ) {
3157 //           cout << "SquareExtent: " << box1.SquareExtent()<<" > "<<box2.SquareExtent() << endl;
3158           exchange = Standard_True;
3159         }
3160         else if ( Abs(dSquareExtent) < tol ) {
3161           Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax, val1, val2;
3162           box1.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
3163           val1 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
3164           box2.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
3165           val2 = (aXmin+aXmax)*999 + (aYmin+aYmax)*99 + (aZmin+aZmax)*0.9;
3166           //exchange = val1 > val2;
3167           if ((val1 - val2) >= tol) {
3168             exchange = Standard_True;
3169           }
3170           //cout << "box: " << val1<<" > "<<val2 << endl;
3171         }
3172       }
3173
3174       if (exchange)
3175       {
3176 //         cout << "exchange " << Index << " & " << Index+1 << endl;
3177         aTemp = OrderInd(Index);
3178         OrderInd(Index) = OrderInd(Index+1);
3179         OrderInd(Index+1) = aTemp;
3180         Sort = Standard_True;
3181       }
3182     }
3183   }
3184
3185   for (Index=1; Index <= MaxShapes; Index++)
3186     SL.Append( aShapes( OrderInd(Index) ));
3187 }
3188
3189 //=======================================================================
3190 //function : CompsolidToCompound
3191 //purpose  :
3192 //=======================================================================
3193 TopoDS_Shape GEOMImpl_IShapesOperations::CompsolidToCompound (const TopoDS_Shape& theCompsolid)
3194 {
3195   if (theCompsolid.ShapeType() != TopAbs_COMPSOLID) {
3196     return theCompsolid;
3197   }
3198
3199   TopoDS_Compound aCompound;
3200   BRep_Builder B;
3201   B.MakeCompound(aCompound);
3202
3203   TopTools_MapOfShape mapShape;
3204   TopoDS_Iterator It (theCompsolid, Standard_True, Standard_True);
3205
3206   for (; It.More(); It.Next()) {
3207     TopoDS_Shape aShape_i = It.Value();
3208     if (mapShape.Add(aShape_i)) {
3209       B.Add(aCompound, aShape_i);
3210     }
3211   }
3212
3213   return aCompound;
3214 }
3215
3216 //=======================================================================
3217 //function : CheckTriangulation
3218 //purpose  :
3219 //=======================================================================
3220 bool GEOMImpl_IShapesOperations::CheckTriangulation (const TopoDS_Shape& aShape)
3221 {
3222   bool isTriangulation = true;
3223
3224   TopExp_Explorer exp (aShape, TopAbs_FACE);
3225   if (exp.More())
3226   {
3227     TopLoc_Location aTopLoc;
3228     Handle(Poly_Triangulation) aTRF;
3229     aTRF = BRep_Tool::Triangulation(TopoDS::Face(exp.Current()), aTopLoc);
3230     if (aTRF.IsNull()) {
3231       isTriangulation = false;
3232     }
3233   }
3234   else // no faces, try edges
3235   {
3236     TopExp_Explorer expe (aShape, TopAbs_EDGE);
3237     if (!expe.More()) {
3238       return false;
3239     }
3240     TopLoc_Location aLoc;
3241     Handle(Poly_Polygon3D) aPE = BRep_Tool::Polygon3D(TopoDS::Edge(expe.Current()), aLoc);
3242     if (aPE.IsNull()) {
3243       isTriangulation = false;
3244     }
3245   }
3246
3247   if (!isTriangulation) {
3248     // calculate deflection
3249     Standard_Real aDeviationCoefficient = 0.001;
3250
3251     Bnd_Box B;
3252     BRepBndLib::Add(aShape, B);
3253     Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
3254     B.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
3255
3256     Standard_Real dx = aXmax - aXmin, dy = aYmax - aYmin, dz = aZmax - aZmin;
3257     Standard_Real aDeflection = Max(Max(dx, dy), dz) * aDeviationCoefficient * 4;
3258     Standard_Real aHLRAngle = 0.349066;
3259
3260     BRepMesh_IncrementalMesh Inc (aShape, aDeflection, Standard_False, aHLRAngle);
3261   }
3262
3263   return true;
3264 }
3265
3266 #define MAX_TOLERANCE 1.e-7
3267
3268
3269 //=======================================================================
3270 //function : isSameEdge
3271 //purpose  : Returns True if two edges coincide
3272 //=======================================================================
3273 static bool isSameEdge(const TopoDS_Edge& theEdge1, const TopoDS_Edge& theEdge2)
3274 {
3275   TopoDS_Vertex V11, V12, V21, V22;
3276   TopExp::Vertices(theEdge1, V11, V12);
3277   TopExp::Vertices(theEdge2, V21, V22);
3278   gp_Pnt P11 = BRep_Tool::Pnt(V11);
3279   gp_Pnt P12 = BRep_Tool::Pnt(V12);
3280   gp_Pnt P21 = BRep_Tool::Pnt(V21);
3281   gp_Pnt P22 = BRep_Tool::Pnt(V22);
3282   bool coincide = false;
3283
3284   //Check that ends of edges coincide
3285   if(P11.Distance(P21) <= MAX_TOLERANCE) {
3286     if(P12.Distance(P22) <= MAX_TOLERANCE) coincide =  true;
3287   }
3288   else if(P11.Distance(P22) <= MAX_TOLERANCE) {
3289     if(P12.Distance(P21) <= MAX_TOLERANCE) coincide = true;
3290   }
3291
3292   if(!coincide) return false;
3293
3294   if (BRep_Tool::Degenerated(theEdge1))
3295     if (BRep_Tool::Degenerated(theEdge2)) return true;
3296     else return false;
3297   else
3298     if (BRep_Tool::Degenerated(theEdge2)) return false;
3299
3300   double U11, U12, U21, U22;
3301   Handle(Geom_Curve) C1 = BRep_Tool::Curve(theEdge1, U11, U12);
3302   Handle(Geom_Curve) C2 = BRep_Tool::Curve(theEdge2, U21, U22);
3303   if(C1->DynamicType() == C2->DynamicType()) return true;
3304
3305   //Check that both edges has the same geometry
3306   double range = U12-U11;
3307   double U = U11+ range/3.0;
3308   gp_Pnt P1 = C1->Value(U);     //Compute a point on one third of the edge's length
3309   U = U11+range*2.0/3.0;
3310   gp_Pnt P2 = C1->Value(U);     //Compute a point on two thirds of the edge's length
3311
3312   if(!GeomLib_Tool::Parameter(C2, P1, MAX_TOLERANCE, U) ||  U < U21 || U > U22)
3313     return false;
3314
3315   if(P1.Distance(C2->Value(U)) > MAX_TOLERANCE) return false;
3316
3317   if(!GeomLib_Tool::Parameter(C2, P2, MAX_TOLERANCE, U) || U < U21 || U > U22)
3318     return false;
3319
3320   if(P2.Distance(C2->Value(U)) > MAX_TOLERANCE) return false;
3321
3322   return true;
3323 }
3324
3325 #include <TopoDS_TShape.hxx>
3326 //=======================================================================
3327 //function : isSameFace
3328 //purpose  : Returns True if two faces coincide
3329 //=======================================================================
3330 static bool isSameFace(const TopoDS_Face& theFace1, const TopoDS_Face& theFace2)
3331 {
3332   TopExp_Explorer E(theFace1, TopAbs_EDGE);
3333   TopTools_ListOfShape LS1, LS2;
3334   for(; E.More(); E.Next()) LS1.Append(E.Current());
3335
3336   E.Init(theFace2, TopAbs_EDGE);
3337   for(; E.More(); E.Next()) LS2.Append(E.Current());
3338
3339   //Compare the number of edges in the faces
3340   if(LS1.Extent() != LS2.Extent()) return false;
3341
3342   double aMin = RealFirst(), aMax = RealLast();
3343   double xminB1=aMax, yminB1=aMax, zminB1=aMax, xminB2=aMax, yminB2=aMax, zminB2=aMax;
3344   double xmaxB1=aMin, ymaxB1=aMin, zmaxB1=aMin, xmaxB2=aMin, ymaxB2=aMin, zmaxB2=aMin;
3345
3346   for(E.Init(theFace1, TopAbs_VERTEX); E.More(); E.Next()) {
3347     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
3348     if(P.X() < xminB1) xminB1 = P.X();
3349     if(P.Y() < yminB1) yminB1 = P.Y();
3350     if(P.Z() < zminB1) zminB1 = P.Z();
3351     if(P.X() > xmaxB1) xmaxB1 = P.X();
3352     if(P.Y() > ymaxB1) ymaxB1 = P.Y();
3353     if(P.Z() > zmaxB1) zmaxB1 = P.Z();
3354   }
3355
3356   for(E.Init(theFace2, TopAbs_VERTEX); E.More(); E.Next()) {
3357     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
3358     if(P.X() < xminB2) xminB2 = P.X();
3359     if(P.Y() < yminB2) yminB2 = P.Y();
3360     if(P.Z() < zminB2) zminB2 = P.Z();
3361     if(P.X() > xmaxB2) xmaxB2 = P.X();
3362     if(P.Y() > ymaxB2) ymaxB2 = P.Y();
3363     if(P.Z() > zmaxB2) zmaxB2 = P.Z();
3364   }
3365
3366
3367   //Compare the bounding boxes of both faces
3368   if(gp_Pnt(xminB1, yminB1, zminB1).Distance(gp_Pnt(xminB2, yminB2, zminB2)) > MAX_TOLERANCE)
3369     return false;
3370
3371   if(gp_Pnt(xmaxB1, ymaxB1, zmaxB1).Distance(gp_Pnt(xmaxB2, ymaxB2, zmaxB2)) > MAX_TOLERANCE)
3372     return false;
3373
3374   Handle(Geom_Surface) S1 = BRep_Tool::Surface(theFace1);
3375   Handle(Geom_Surface) S2 = BRep_Tool::Surface(theFace2);
3376
3377   //Check if there a coincidence of two surfaces at least in two points
3378   double U11, U12, V11, V12, U21, U22, V21, V22;
3379   BRepTools::UVBounds(theFace1, U11, U12, V11, V12);
3380   BRepTools::UVBounds(theFace2, U21, U22, V21, V22);
3381
3382   double rangeU = U12-U11;
3383   double rangeV = V12-V11;
3384   double U = U11 + rangeU/3.0;
3385   double V = V11 + rangeV/3.0;
3386   gp_Pnt P1 = S1->Value(U, V);
3387   U = U11+rangeU*2.0/3.0;
3388   V = V11+rangeV*2.0/3.0;
3389   gp_Pnt P2 = S1->Value(U, V);
3390
3391   if(!GeomLib_Tool::Parameters(S2, P1, MAX_TOLERANCE, U, V) || U < U21 || U > U22 || V < V21 || V > V22)
3392     return false;
3393
3394   if(P1.Distance(S2->Value(U,V)) > MAX_TOLERANCE) return false;
3395
3396   if(!GeomLib_Tool::Parameters(S2, P2, MAX_TOLERANCE, U, V) || U < U21 || U > U22 || V < V21 || V > V22)
3397     return false;
3398
3399   if(P2.Distance(S2->Value(U, V)) > MAX_TOLERANCE) return false;
3400
3401   //Check that each edge of the Face1 has a counterpart in the Face2
3402   TopTools_MapOfOrientedShape aMap;
3403   TopTools_ListIteratorOfListOfShape LSI1(LS1);
3404   for(; LSI1.More(); LSI1.Next()) {
3405     TopoDS_Edge E = TopoDS::Edge(LSI1.Value());
3406     bool isFound = false;
3407     TopTools_ListIteratorOfListOfShape LSI2(LS2);
3408     for(; LSI2.More(); LSI2.Next()) {
3409       TopoDS_Shape aValue = LSI2.Value();
3410       if(aMap.Contains(aValue)) continue; //To avoid checking already found edge several times
3411       if(isSameEdge(E, TopoDS::Edge(aValue))) {
3412         aMap.Add(aValue);
3413         isFound = true;
3414         break;
3415       }
3416     }
3417     if(!isFound) return false;
3418   }
3419
3420   return true;
3421 }
3422
3423 //=======================================================================
3424 //function : isSameSolid
3425 //purpose  : Returns True if two solids coincide
3426 //=======================================================================
3427 bool isSameSolid(const TopoDS_Solid& theSolid1, const TopoDS_Solid& theSolid2)
3428 {
3429   TopExp_Explorer E(theSolid1, TopAbs_FACE);
3430   TopTools_ListOfShape LS1, LS2;
3431   for(; E.More(); E.Next()) LS1.Append(E.Current());
3432   E.Init(theSolid2, TopAbs_FACE);
3433   for(; E.More(); E.Next()) LS2.Append(E.Current());
3434
3435   if(LS1.Extent() != LS2.Extent()) return false;
3436
3437   double aMin = RealFirst(), aMax = RealLast();
3438   double xminB1=aMax, yminB1=aMax, zminB1=aMax, xminB2=aMax, yminB2=aMax, zminB2=aMax;
3439   double xmaxB1=aMin, ymaxB1=aMin, zmaxB1=aMin, xmaxB2=aMin, ymaxB2=aMin, zmaxB2=aMin;
3440
3441   for(E.Init(theSolid1, TopAbs_VERTEX); E.More(); E.Next()) {
3442     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
3443     if(P.X() < xminB1) xminB1 = P.X();
3444     if(P.Y() < yminB1) yminB1 = P.Y();
3445     if(P.Z() < zminB1) zminB1 = P.Z();
3446     if(P.X() > xmaxB1) xmaxB1 = P.X();
3447     if(P.Y() > ymaxB1) ymaxB1 = P.Y();
3448     if(P.Z() > zmaxB1) zmaxB1 = P.Z();
3449   }
3450
3451   for(E.Init(theSolid2, TopAbs_VERTEX); E.More(); E.Next()) {
3452     gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
3453     if(P.X() < xminB2) xminB2 = P.X();
3454     if(P.Y() < yminB2) yminB2 = P.Y();
3455     if(P.Z() < zminB2) zminB2 = P.Z();
3456     if(P.X() > xmaxB2) xmaxB2 = P.X();
3457     if(P.Y() > ymaxB2) ymaxB2 = P.Y();
3458     if(P.Z() > zmaxB2) zmaxB2 = P.Z();
3459   }
3460
3461   //Compare the bounding boxes of both solids
3462   if(gp_Pnt(xminB1, yminB1, zminB1).Distance(gp_Pnt(xminB2, yminB2, zminB2)) > MAX_TOLERANCE)
3463     return false;
3464
3465   if(gp_Pnt(xmaxB1, ymaxB1, zmaxB1).Distance(gp_Pnt(xmaxB2, ymaxB2, zmaxB2)) > MAX_TOLERANCE)
3466     return false;
3467
3468   //Check that each face of the Solid1 has a counterpart in the Solid2
3469   TopTools_MapOfOrientedShape aMap;
3470   TopTools_ListIteratorOfListOfShape LSI1(LS1);
3471   for(; LSI1.More(); LSI1.Next()) {
3472     TopoDS_Face F = TopoDS::Face(LSI1.Value());
3473     bool isFound = false;
3474     TopTools_ListIteratorOfListOfShape LSI2(LS2);
3475     for(; LSI2.More(); LSI2.Next()) {
3476       if(aMap.Contains(LSI2.Value())) continue; //To avoid checking already found faces several times
3477       if(isSameFace(F, TopoDS::Face(LSI2.Value()))) {
3478         aMap.Add(LSI2.Value());
3479         isFound = true;
3480         break;
3481       }
3482     }
3483     if(!isFound) return false;
3484   }
3485
3486   return true;
3487 }
3488
3489 //=======================================================================
3490 //function : GetSame
3491 //purpose  :
3492 //=======================================================================
3493 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetSame(const Handle(GEOM_Object)& theShapeWhere,
3494                                                         const Handle(GEOM_Object)& theShapeWhat)
3495 {
3496   SetErrorCode(KO);
3497   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
3498
3499   TopoDS_Shape aWhere = theShapeWhere->GetValue();
3500   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
3501
3502   if (aWhere.IsNull() || aWhat.IsNull()) return NULL;
3503
3504   int anIndex = -1;
3505   bool isFound = false;
3506   TopoDS_Shape aSubShape;
3507   TopTools_MapOfShape aMap;
3508
3509   switch(aWhat.ShapeType()) {
3510     case TopAbs_VERTEX: {
3511       gp_Pnt P = BRep_Tool::Pnt(TopoDS::Vertex(aWhat));
3512       TopExp_Explorer E(aWhere, TopAbs_VERTEX);
3513       for(; E.More(); E.Next()) {
3514         if(!aMap.Add(E.Current())) continue;
3515         gp_Pnt P2 = BRep_Tool::Pnt(TopoDS::Vertex(E.Current()));
3516         if(P.Distance(P2) <= MAX_TOLERANCE) {
3517           isFound = true;
3518           aSubShape = E.Current();
3519           break;
3520         }
3521       }
3522       break;
3523                         }
3524     case TopAbs_FACE: {
3525       TopoDS_Face aFace = TopoDS::Face(aWhat);
3526       TopExp_Explorer E(aWhere, TopAbs_FACE);
3527       for(; E.More(); E.Next()) {
3528         if(!aMap.Add(E.Current())) continue;
3529         if(isSameFace(aFace, TopoDS::Face(E.Current()))) {
3530           aSubShape = E.Current();
3531           isFound = true;
3532           break;
3533         }
3534       }
3535       break;
3536                       }
3537     case TopAbs_EDGE: {
3538       TopoDS_Edge anEdge = TopoDS::Edge(aWhat);
3539       TopExp_Explorer E(aWhere, TopAbs_EDGE);
3540       for(; E.More(); E.Next()) {
3541         if(!aMap.Add(E.Current())) continue;
3542         if(isSameEdge(anEdge, TopoDS::Edge(E.Current()))) {
3543           aSubShape = E.Current();
3544           isFound = true;
3545           break;
3546         }
3547       }
3548       break;
3549                       }
3550     case TopAbs_SOLID: {
3551       TopoDS_Solid aSolid = TopoDS::Solid(aWhat);
3552       TopExp_Explorer E(aWhere, TopAbs_SOLID);
3553       for(; E.More(); E.Next()) {
3554         if(!aMap.Add(E.Current())) continue;
3555         if(isSameSolid(aSolid, TopoDS::Solid(E.Current()))) {
3556           aSubShape = E.Current();
3557           isFound = true;
3558           break;
3559         }
3560       }
3561       break;
3562                        }
3563     default:
3564       return NULL;
3565   }
3566
3567   if(isFound) {
3568     TopTools_IndexedMapOfShape anIndices;
3569     TopExp::MapShapes(aWhere, anIndices);
3570     if (anIndices.Contains(aSubShape))
3571       anIndex = anIndices.FindIndex(aSubShape);
3572   }
3573
3574   if(anIndex < 0) return NULL;
3575
3576   Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
3577
3578   anArray->SetValue(1, anIndex);
3579
3580   Handle(GEOM_Object) aResult = GetEngine()->AddSubShape(theShapeWhere, anArray);
3581   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
3582
3583   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetSame("
3584     << theShapeWhere << ", " << theShapeWhat << ")";
3585
3586   SetErrorCode(OK);
3587
3588   return aResult;
3589 }