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