src/GEOMImpl/GEOMImpl_IShapesOperations.cxx

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