src/GEOMImpl/GEOMImpl_IShapesOperations.cxx

   1 #include <Standard_Stream.hxx>
   2
   3 #include <GEOMImpl_IShapesOperations.hxx>
   4
   5 #include <GEOMImpl_Types.hxx>
   6
   7 #include <GEOMImpl_VectorDriver.hxx>
   8 #include <GEOMImpl_ShapeDriver.hxx>
   9 #include <GEOMImpl_CopyDriver.hxx>
  10 #include <GEOMImpl_GlueDriver.hxx>
  11
  12 #include <GEOMImpl_IVector.hxx>
  13 #include <GEOMImpl_IShapes.hxx>
  14 #include <GEOMImpl_IGlue.hxx>
  15
  16 #include <GEOMImpl_Block6Explorer.hxx>
  17
  18 #include <GEOM_Function.hxx>
  19 #include <GEOM_PythonDump.hxx>
  20
  21 #include <GEOMAlgo_FinderShapeOn1.hxx>
  22
  23 #include "utilities.h"
  24 #include <OpUtil.hxx>
  25 #include <Utils_ExceptHandlers.hxx>
  26
  27 #include <TFunction_DriverTable.hxx>
  28 #include <TFunction_Driver.hxx>
  29 #include <TFunction_Logbook.hxx>
  30 #include <TDataStd_Integer.hxx>
  31 #include <TDataStd_IntegerArray.hxx>
  32 #include <TDF_Tool.hxx>
  33
  34 #include <BRepExtrema_ExtCF.hxx>
  35
  36 #include <BRep_Tool.hxx>
  37 #include <BRepGProp.hxx>
  38 #include <BRepAdaptor_Curve.hxx>
  39 #include <BRepBndLib.hxx>
  40 #include <BRepBuilderAPI_MakeFace.hxx>
  41 #include <BRepMesh_IncrementalMesh.hxx>
  42
  43 #include <TopAbs.hxx>
  44 #include <TopExp.hxx>
  45 #include <TopoDS.hxx>
  46 #include <TopoDS_Shape.hxx>
  47 #include <TopoDS_Face.hxx>
  48 #include <TopoDS_Edge.hxx>
  49 #include <TopoDS_Vertex.hxx>
  50 #include <TopoDS_Iterator.hxx>
  51 #include <TopExp_Explorer.hxx>
  52 #include <TopLoc_Location.hxx>
  53 #include <TopTools_MapOfShape.hxx>
  54 #include <TopTools_Array1OfShape.hxx>
  55 #include <TopTools_ListIteratorOfListOfShape.hxx>
  56 #include <TopTools_IndexedMapOfShape.hxx>
  57
  58 #include <Geom_Surface.hxx>
  59 #include <Geom_Plane.hxx>
  60 #include <Geom_SphericalSurface.hxx>
  61 #include <Geom_CylindricalSurface.hxx>
  62 #include <GeomAdaptor_Surface.hxx>
  63
  64 #include <Geom2d_Curve.hxx>
  65
  66 #include <Bnd_Box.hxx>
  67 #include <GProp_GProps.hxx>
  68 #include <gp_Pnt.hxx>
  69 #include <gp_Lin.hxx>
  70 #include <TColStd_Array1OfReal.hxx>
  71 #include <TColStd_HArray1OfInteger.hxx>
  72 #include <TColStd_ListOfInteger.hxx>
  73 #include <TColStd_ListIteratorOfListOfInteger.hxx>
  74 #include <TColStd_DataMapOfIntegerInteger.hxx>
  75 #include <TColStd_DataMapIteratorOfDataMapOfIntegerInteger.hxx>
  76
  77 #include <vector>
  78 //#include <iostream>
  79
  80 //#include <OSD_Timer.hxx>
  81
  82 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
  83
  84 //=============================================================================
  85 /*!
  86  *   constructor:
  87  */
  88 //=============================================================================
  89 GEOMImpl_IShapesOperations::GEOMImpl_IShapesOperations (GEOM_Engine* theEngine, int theDocID)
  90 : GEOM_IOperations(theEngine, theDocID)
  91 {
  92   MESSAGE("GEOMImpl_IShapesOperations::GEOMImpl_IShapesOperations");
  93 }
  94
  95 //=============================================================================
  96 /*!
  97  *  destructor
  98  */
  99 //=============================================================================
 100 GEOMImpl_IShapesOperations::~GEOMImpl_IShapesOperations()
 101 {
 102   MESSAGE("GEOMImpl_IShapesOperations::~GEOMImpl_IShapesOperations");
 103 }
 104
 105
 106 //=============================================================================
 107 /*!
 108  *  MakeEdge
 109  */
 110 //=============================================================================
 111 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeEdge
 112                      (Handle(GEOM_Object) thePnt1, Handle(GEOM_Object) thePnt2)
 113 {
 114   SetErrorCode(KO);
 115
 116   if (thePnt1.IsNull() || thePnt2.IsNull()) return NULL;
 117
 118   //Add a new Edge object
 119   Handle(GEOM_Object) anEdge = GetEngine()->AddObject(GetDocID(), GEOM_EDGE);
 120
 121   //Add a new Vector function
 122   Handle(GEOM_Function) aFunction =
 123     anEdge->AddFunction(GEOMImpl_VectorDriver::GetID(), VECTOR_TWO_PNT);
 124
 125   //Check if the function is set correctly
 126   if (aFunction->GetDriverGUID() != GEOMImpl_VectorDriver::GetID()) return NULL;
 127
 128   GEOMImpl_IVector aPI (aFunction);
 129
 130   Handle(GEOM_Function) aRef1 = thePnt1->GetLastFunction();
 131   Handle(GEOM_Function) aRef2 = thePnt2->GetLastFunction();
 132   if (aRef1.IsNull() || aRef2.IsNull()) return NULL;
 133
 134   aPI.SetPoint1(aRef1);
 135   aPI.SetPoint2(aRef2);
 136
 137   //Compute the Edge value
 138   try {
 139     if (!GetSolver()->ComputeFunction(aFunction)) {
 140       SetErrorCode("Vector driver failed");
 141       return NULL;
 142     }
 143   }
 144   catch (Standard_Failure) {
 145     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 146     SetErrorCode(aFail->GetMessageString());
 147     return NULL;
 148   }
 149
 150   //Make a Python command
 151   GEOM::TPythonDump(aFunction) << anEdge << " = geompy.MakeEdge("
 152                                << thePnt1 << ", " << thePnt2 << ")";
 153
 154   SetErrorCode(OK);
 155   return anEdge;
 156 }
 157
 158 //=============================================================================
 159 /*!
 160  *  MakeWire
 161  */
 162 //=============================================================================
 163 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeWire
 164                              (list<Handle(GEOM_Object)> theShapes)
 165 {
 166   return MakeShape(theShapes, GEOM_WIRE, WIRE_EDGES, "MakeWire");
 167 }
 168
 169 //=============================================================================
 170 /*!
 171  *  MakeFace
 172  */
 173 //=============================================================================
 174 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeFace (Handle(GEOM_Object) theWire,
 175                                                           const bool isPlanarWanted)
 176 {
 177   SetErrorCode(KO);
 178
 179   if (theWire.IsNull()) return NULL;
 180
 181   //Add a new Face object
 182   Handle(GEOM_Object) aFace = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 183
 184   //Add a new Shape function for creation of a face from a wire
 185   Handle(GEOM_Function) aFunction =
 186     aFace->AddFunction(GEOMImpl_ShapeDriver::GetID(), FACE_WIRE);
 187   if (aFunction.IsNull()) return NULL;
 188
 189   //Check if the function is set correctly
 190   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 191
 192   GEOMImpl_IShapes aCI (aFunction);
 193
 194   Handle(GEOM_Function) aRefWire = theWire->GetLastFunction();
 195
 196   if (aRefWire.IsNull()) return NULL;
 197
 198   aCI.SetBase(aRefWire);
 199   aCI.SetIsPlanar(isPlanarWanted);
 200
 201   //Compute the Face value
 202   try {
 203     if (!GetSolver()->ComputeFunction(aFunction)) {
 204       SetErrorCode("Shape driver failed to compute a face");
 205       return NULL;
 206     }
 207   }
 208   catch (Standard_Failure) {
 209     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 210     SetErrorCode(aFail->GetMessageString());
 211     return NULL;
 212   }
 213
 214   //Make a Python command
 215   GEOM::TPythonDump(aFunction) << aFace << " = geompy.MakeFace("
 216     << theWire << ", " << (int)isPlanarWanted << ")";
 217
 218   SetErrorCode(OK);
 219   return aFace;
 220 }
 221
 222 //=============================================================================
 223 /*!
 224  *  MakeFaceWires
 225  */
 226 //=============================================================================
 227 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeFaceWires
 228                              (list<Handle(GEOM_Object)> theShapes,
 229                               const bool isPlanarWanted)
 230 {
 231   SetErrorCode(KO);
 232
 233   //Add a new object
 234   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 235
 236   //Add a new function
 237   Handle(GEOM_Function) aFunction =
 238     aShape->AddFunction(GEOMImpl_ShapeDriver::GetID(), FACE_WIRES);
 239   if (aFunction.IsNull()) return NULL;
 240
 241   //Check if the function is set correctly
 242   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 243
 244   GEOMImpl_IShapes aCI (aFunction);
 245
 246   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 247
 248   // Shapes
 249   list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 250   for (; it != theShapes.end(); it++) {
 251     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 252     if (aRefSh.IsNull()) {
 253       SetErrorCode("NULL argument shape for the face construction");
 254       return NULL;
 255     }
 256     aShapesSeq->Append(aRefSh);
 257   }
 258   aCI.SetShapes(aShapesSeq);
 259
 260   aCI.SetIsPlanar(isPlanarWanted);
 261
 262   //Compute the shape
 263   try {
 264     if (!GetSolver()->ComputeFunction(aFunction)) {
 265       SetErrorCode("Shape driver failed");
 266       return NULL;
 267     }
 268   }
 269   catch (Standard_Failure) {
 270     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 271     SetErrorCode(aFail->GetMessageString());
 272     return NULL;
 273   }
 274
 275   //Make a Python command
 276   GEOM::TPythonDump pd (aFunction);
 277   pd << aShape << " = geompy.MakeFaceWires([";
 278
 279   // Shapes
 280   it = theShapes.begin();
 281   if (it != theShapes.end()) {
 282     pd << (*it++);
 283     while (it != theShapes.end()) {
 284       pd << ", " << (*it++);
 285     }
 286   }
 287   pd << "], " << (int)isPlanarWanted << ")";
 288
 289   SetErrorCode(OK);
 290   return aShape;
 291 }
 292
 293 //=============================================================================
 294 /*!
 295  *  MakeShell
 296  */
 297 //=============================================================================
 298 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeShell
 299                              (list<Handle(GEOM_Object)> theShapes)
 300 {
 301   return MakeShape(theShapes, GEOM_SHELL, SHELL_FACES, "MakeShell");
 302 }
 303
 304 //=============================================================================
 305 /*!
 306  *  MakeSolidShells
 307  */
 308 //=============================================================================
 309 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeSolidShells
 310                              (list<Handle(GEOM_Object)> theShapes)
 311 {
 312   return MakeShape(theShapes, GEOM_SOLID, SOLID_SHELLS, "MakeSolidShells");
 313 }
 314
 315 //=============================================================================
 316 /*!
 317  *  MakeSolidShell
 318  */
 319 //=============================================================================
 320 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeSolidShell (Handle(GEOM_Object) theShell)
 321 {
 322   SetErrorCode(KO);
 323
 324   if (theShell.IsNull()) return NULL;
 325
 326   //Add a new Solid object
 327   Handle(GEOM_Object) aSolid = GetEngine()->AddObject(GetDocID(), GEOM_SOLID);
 328
 329   //Add a new Solid function for creation of a solid from a shell
 330   Handle(GEOM_Function) aFunction =
 331     aSolid->AddFunction(GEOMImpl_ShapeDriver::GetID(), SOLID_SHELL);
 332   if (aFunction.IsNull()) return NULL;
 333
 334   //Check if the function is set correctly
 335   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 336
 337   GEOMImpl_IShapes aCI (aFunction);
 338
 339   Handle(GEOM_Function) aRefShell = theShell->GetLastFunction();
 340
 341   if (aRefShell.IsNull()) return NULL;
 342
 343   aCI.SetBase(aRefShell);
 344
 345   //Compute the Solid value
 346   try {
 347     if (!GetSolver()->ComputeFunction(aFunction)) {
 348       SetErrorCode("Solid driver failed");
 349       return NULL;
 350     }
 351   }
 352   catch (Standard_Failure) {
 353     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 354     SetErrorCode(aFail->GetMessageString());
 355     return NULL;
 356   }
 357
 358   //Make a Python command
 359   GEOM::TPythonDump(aFunction) << aSolid
 360     << " = geompy.MakeSolid(" << theShell << ")";
 361
 362   SetErrorCode(OK);
 363   return aSolid;
 364 }
 365
 366 //=============================================================================
 367 /*!
 368  *  MakeCompound
 369  */
 370 //=============================================================================
 371 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeCompound
 372                              (list<Handle(GEOM_Object)> theShapes)
 373 {
 374   return MakeShape(theShapes, GEOM_COMPOUND, COMPOUND_SHAPES, "MakeCompound");
 375 }
 376
 377 //=============================================================================
 378 /*!
 379  *  MakeShape
 380  */
 381 //=============================================================================
 382 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeShape
 383                              (list<Handle(GEOM_Object)>      theShapes,
 384                               const Standard_Integer         theObjectType,
 385                               const Standard_Integer         theFunctionType,
 386                               const TCollection_AsciiString& theMethodName)
 387 {
 388   SetErrorCode(KO);
 389
 390   //Add a new object
 391   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), theObjectType);
 392
 393   //Add a new function
 394   Handle(GEOM_Function) aFunction =
 395     aShape->AddFunction(GEOMImpl_ShapeDriver::GetID(), theFunctionType);
 396   if (aFunction.IsNull()) return NULL;
 397
 398   //Check if the function is set correctly
 399   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 400
 401   GEOMImpl_IShapes aCI (aFunction);
 402
 403   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 404
 405   // Shapes
 406   list<Handle(GEOM_Object)>::iterator it = theShapes.begin();
 407   for (; it != theShapes.end(); it++) {
 408     Handle(GEOM_Function) aRefSh = (*it)->GetLastFunction();
 409     if (aRefSh.IsNull()) {
 410       SetErrorCode("NULL argument shape for the shape construction");
 411       return NULL;
 412     }
 413     aShapesSeq->Append(aRefSh);
 414   }
 415   aCI.SetShapes(aShapesSeq);
 416
 417   //Compute the shape
 418   try {
 419     if (!GetSolver()->ComputeFunction(aFunction)) {
 420       SetErrorCode("Shape driver failed");
 421       return NULL;
 422     }
 423   }
 424   catch (Standard_Failure) {
 425     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 426     SetErrorCode(aFail->GetMessageString());
 427     return NULL;
 428   }
 429
 430   //Make a Python command
 431   GEOM::TPythonDump pd (aFunction);
 432   pd << aShape << " = geompy." << theMethodName.ToCString() << "([";
 433
 434   // Shapes
 435   it = theShapes.begin();
 436   if (it != theShapes.end()) {
 437     pd << (*it++);
 438     while (it != theShapes.end()) {
 439       pd << ", " << (*it++);
 440     }
 441   }
 442   pd << "])";
 443
 444   SetErrorCode(OK);
 445   return aShape;
 446 }
 447
 448 //=============================================================================
 449 /*!
 450  *  MakeGlueFaces
 451  */
 452 //=============================================================================
 453 Handle(GEOM_Object) GEOMImpl_IShapesOperations::MakeGlueFaces
 454                                                 (Handle(GEOM_Object) theShape,
 455                                                  const Standard_Real theTolerance)
 456 {
 457   SetErrorCode(KO);
 458
 459   if (theShape.IsNull()) return NULL;
 460
 461   //Add a new Glued object
 462   Handle(GEOM_Object) aGlued = GetEngine()->AddObject(GetDocID(), GEOM_GLUED);
 463
 464   //Add a new Glue function
 465   Handle(GEOM_Function) aFunction;
 466   aFunction = aGlued->AddFunction(GEOMImpl_GlueDriver::GetID(), GLUE_FACES);
 467   if (aFunction.IsNull()) return NULL;
 468
 469   //Check if the function is set correctly
 470   if (aFunction->GetDriverGUID() != GEOMImpl_GlueDriver::GetID()) return NULL;
 471
 472   GEOMImpl_IGlue aCI (aFunction);
 473
 474   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 475   if (aRefShape.IsNull()) return NULL;
 476
 477   aCI.SetBase(aRefShape);
 478   aCI.SetTolerance(theTolerance);
 479
 480   //Compute the sub-shape value
 481   Standard_Boolean isWarning = Standard_False;
 482   try {
 483     if (!GetSolver()->ComputeFunction(aFunction)) {
 484       SetErrorCode("Shape driver failed to glue faces");
 485       return NULL;
 486     }
 487   }
 488   catch (Standard_Failure) {
 489     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 490     SetErrorCode(aFail->GetMessageString());
 491     // to provide warning
 492     if (!aFunction->GetValue().IsNull()) {
 493       isWarning = Standard_True;
 494     } else {
 495       return NULL;
 496     }
 497   }
 498
 499   //Make a Python command
 500   GEOM::TPythonDump(aFunction) << aGlued << " = geompy.MakeGlueFaces("
 501     << theShape << ", " << theTolerance << ")";
 502
 503   // to provide warning
 504   if (!isWarning) SetErrorCode(OK);
 505   return aGlued;
 506 }
 507
 508 //=============================================================================
 509 /*!
 510  *  MakeExplode
 511  */
 512 //=============================================================================
 513 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::MakeExplode
 514                                           (Handle(GEOM_Object)    theShape,
 515                                            const Standard_Integer theShapeType,
 516                                            const Standard_Boolean isSorted)
 517 {
 518 //  OSD_Timer timer1, timer2, timer3, timer4;
 519 //  timer1.Start();
 520
 521   SetErrorCode(KO);
 522
 523   if (theShape.IsNull()) return NULL;
 524   TopoDS_Shape aShape = theShape->GetValue();
 525   if (aShape.IsNull()) return NULL;
 526
 527   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
 528   Handle(GEOM_Object) anObj;
 529   Handle(GEOM_Function) aFunction;
 530   TopTools_MapOfShape mapShape;
 531   TopTools_ListOfShape listShape;
 532
 533   if (aShape.ShapeType() == TopAbs_COMPOUND &&
 534       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 535        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
 536        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND)) {
 537     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
 538     for (; It.More(); It.Next()) {
 539       if (mapShape.Add(It.Value())) {
 540         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 541             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
 542           listShape.Append(It.Value());
 543         }
 544       }
 545     }
 546   } else {
 547     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
 548     for (; exp.More(); exp.Next())
 549       if (mapShape.Add(exp.Current()))
 550         listShape.Append(exp.Current());
 551   }
 552
 553   if (listShape.IsEmpty()) {
 554     SetErrorCode("The given shape has no sub-shapes of the requested type");
 555     return aSeq;
 556   }
 557
 558 //  timer1.Stop();
 559 //  timer2.Start();
 560
 561   if (isSorted)
 562     SortShapes(listShape);
 563
 564 //  timer2.Stop();
 565 //  timer3.Start();
 566
 567   TopTools_IndexedMapOfShape anIndices;
 568   TopExp::MapShapes(aShape, anIndices);
 569   Handle(TColStd_HArray1OfInteger) anArray;
 570
 571   TopTools_ListIteratorOfListOfShape itSub (listShape);
 572   TCollection_AsciiString anAsciiList, anEntry;
 573   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
 574     TopoDS_Shape aValue = itSub.Value();
 575     anArray = new TColStd_HArray1OfInteger(1,1);
 576     anArray->SetValue(1, anIndices.FindIndex(aValue));
 577     anObj = GetEngine()->AddSubShape(theShape, anArray);
 578     aSeq->Append(anObj);
 579
 580     // for python command
 581     TDF_Tool::Entry(anObj->GetEntry(), anEntry);
 582     anAsciiList += anEntry;
 583     anAsciiList += ",";
 584   }
 585
 586   //Make a Python command
 587   anAsciiList.Trunc(anAsciiList.Length() - 1);
 588
 589   aFunction = theShape->GetLastFunction();
 590   TCollection_AsciiString anOldDescr = aFunction->GetDescription();
 591
 592   GEOM::TPythonDump pd (aFunction);
 593   pd << anOldDescr.ToCString() << "\n\t[" << anAsciiList.ToCString();
 594   pd << "] = geompy.SubShapeAll" << (isSorted ? "Sorted(" : "(");
 595   pd << theShape << ", " << theShapeType << ")";
 596
 597   SetErrorCode(OK);
 598
 599 //  timer4.Stop();
 600
 601 //  cout << "Explosure takes:" << endl;
 602 //  timer1.Show();
 603 //  cout << "Sorting takes:" << endl;
 604 //  timer2.Show();
 605 //  cout << "Sub-shapes addition takes:" << endl;
 606 //  timer3.Show();
 607 //  cout << "Update Description takes:" << endl;
 608 //  timer4.Show();
 609
 610   return aSeq;
 611 }
 612
 613 //=============================================================================
 614 /*!
 615  *  GetSubShapeAllIDs
 616  */
 617 //=============================================================================
 618 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::SubShapeAllIDs
 619                                           (Handle(GEOM_Object)    theShape,
 620                                            const Standard_Integer theShapeType,
 621                                            const Standard_Boolean isSorted)
 622 {
 623   SetErrorCode(KO);
 624
 625   if (theShape.IsNull()) return NULL;
 626   TopoDS_Shape aShape = theShape->GetValue();
 627   if (aShape.IsNull()) return NULL;
 628
 629   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
 630   TopTools_MapOfShape mapShape;
 631   TopTools_ListOfShape listShape;
 632
 633   if (aShape.ShapeType() == TopAbs_COMPOUND &&
 634       (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 635        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPSOLID ||
 636        TopAbs_ShapeEnum(theShapeType) == TopAbs_COMPOUND)) {
 637     TopoDS_Iterator It (aShape, Standard_True, Standard_True);
 638     for (; It.More(); It.Next()) {
 639       if (mapShape.Add(It.Value())) {
 640         if (TopAbs_ShapeEnum(theShapeType) == TopAbs_SHAPE ||
 641             TopAbs_ShapeEnum(theShapeType) == It.Value().ShapeType()) {
 642           listShape.Append(It.Value());
 643         }
 644       }
 645     }
 646   } else {
 647     TopExp_Explorer exp (aShape, TopAbs_ShapeEnum(theShapeType));
 648     for (; exp.More(); exp.Next())
 649       if (mapShape.Add(exp.Current()))
 650         listShape.Append(exp.Current());
 651   }
 652
 653   if (listShape.IsEmpty()) {
 654     SetErrorCode("The given shape has no sub-shapes of the requested type");
 655     return aSeq;
 656   }
 657
 658   if (isSorted)
 659     SortShapes(listShape);
 660
 661   TopTools_IndexedMapOfShape anIndices;
 662   TopExp::MapShapes(aShape, anIndices);
 663   Handle(TColStd_HArray1OfInteger) anArray;
 664
 665   TopTools_ListIteratorOfListOfShape itSub (listShape);
 666   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
 667     TopoDS_Shape aValue = itSub.Value();
 668     aSeq->Append(anIndices.FindIndex(aValue));
 669   }
 670
 671   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 672   TCollection_AsciiString anOldDescr = aFunction->GetDescription();
 673
 674   //Make a Python command
 675   GEOM::TPythonDump pd (aFunction);
 676   pd << anOldDescr.ToCString() << "\n\tlistSubShapeIDs = geompy.SubShapeAll";
 677   pd << (isSorted ? "SortedIDs(" : "IDs(");
 678   pd << theShape << ", " << theShapeType << ")";
 679
 680   SetErrorCode(OK);
 681   return aSeq;
 682 }
 683
 684 //=============================================================================
 685 /*!
 686  *  GetSubShape
 687  */
 688 //=============================================================================
 689 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetSubShape
 690                                           (Handle(GEOM_Object)    theMainShape,
 691                                            const Standard_Integer theID)
 692 {
 693   SetErrorCode(KO);
 694
 695   if (theMainShape.IsNull()) return NULL;
 696
 697   Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
 698   anArray->SetValue(1, theID);
 699   Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theMainShape, anArray,true);
 700   if (anObj.IsNull()) {
 701     SetErrorCode("Can not get a sub-shape with the given ID");
 702     return NULL;
 703   }
 704
 705   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
 706
 707   //Make a Python command
 708   GEOM::TPythonDump(aFunction) << anObj << " = geompy.GetSubShape("
 709                                << theMainShape << ", [" << theID << "])";
 710
 711   SetErrorCode(OK);
 712   return anObj;
 713 }
 714
 715
 716 //=============================================================================
 717 /*!
 718  *  NumberOfFaces
 719  */
 720 //=============================================================================
 721 Standard_Integer GEOMImpl_IShapesOperations::NumberOfFaces (Handle(GEOM_Object) theShape)
 722 {
 723   SetErrorCode(KO);
 724
 725   Standard_Integer nb = 0;
 726
 727   if (theShape.IsNull()) return -1;
 728   TopoDS_Shape aShape = theShape->GetValue();
 729   if (aShape.IsNull()) return -1;
 730
 731   TopTools_MapOfShape mapShape;
 732
 733   TopExp_Explorer exp (aShape, TopAbs_FACE);
 734   for (; exp.More(); exp.Next())
 735     if (mapShape.Add(exp.Current()))
 736       nb++;
 737
 738   SetErrorCode(OK);
 739   return nb;
 740 }
 741
 742 //=============================================================================
 743 /*!
 744  *  NumberOfEdges
 745  */
 746 //=============================================================================
 747 Standard_Integer GEOMImpl_IShapesOperations::NumberOfEdges (Handle(GEOM_Object) theShape)
 748 {
 749   SetErrorCode(KO);
 750
 751   Standard_Integer nb = 0;
 752
 753   if (theShape.IsNull()) return -1;
 754   TopoDS_Shape aShape = theShape->GetValue();
 755   if (aShape.IsNull()) return -1;
 756
 757   TopTools_MapOfShape mapShape;
 758
 759   TopExp_Explorer exp (aShape, TopAbs_EDGE);
 760   for (; exp.More(); exp.Next())
 761     if (mapShape.Add(exp.Current()))
 762       nb++;
 763
 764   SetErrorCode(OK);
 765   return nb;
 766 }
 767
 768 //=============================================================================
 769 /*!
 770  *  ReverseShape
 771  */
 772 //=============================================================================
 773 Handle(GEOM_Object) GEOMImpl_IShapesOperations::ReverseShape(Handle(GEOM_Object) theShape)
 774 {
 775   SetErrorCode(KO);
 776
 777   if (theShape.IsNull()) return NULL;
 778
 779   //Add a new reversed object
 780   Handle(GEOM_Object) aReversed = GetEngine()->AddObject(GetDocID(), theShape->GetType());
 781
 782   //Add a new Revese function
 783   Handle(GEOM_Function) aFunction;
 784   aFunction = aReversed->AddFunction(GEOMImpl_ShapeDriver::GetID(), REVERSE_ORIENTATION);
 785   if (aFunction.IsNull()) return NULL;
 786
 787   //Check if the function is set correctly
 788   if (aFunction->GetDriverGUID() != GEOMImpl_ShapeDriver::GetID()) return NULL;
 789
 790   GEOMImpl_IShapes aSI (aFunction);
 791
 792   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
 793   if (aRefShape.IsNull()) return NULL;
 794
 795   aSI.SetBase(aRefShape);
 796
 797   //Compute the sub-shape value
 798   try {
 799     if (!GetSolver()->ComputeFunction(aFunction)) {
 800       SetErrorCode("Shape driver failed to reverse shape");
 801       return NULL;
 802     }
 803   }
 804   catch (Standard_Failure) {
 805     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 806     SetErrorCode(aFail->GetMessageString());
 807     return NULL;
 808   }
 809
 810   //Make a Python command
 811   GEOM::TPythonDump(aFunction) << aReversed
 812     << " = geompy.ChangeOrientation(" << theShape << ")";
 813
 814   SetErrorCode(OK);
 815   return aReversed;
 816 }
 817
 818 //=============================================================================
 819 /*!
 820  *  GetFreeFacesIDs
 821  */
 822 //=============================================================================
 823 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetFreeFacesIDs
 824                                                  (Handle(GEOM_Object) theShape)
 825 {
 826   SetErrorCode(KO);
 827
 828   if (theShape.IsNull()) return NULL;
 829   TopoDS_Shape aShape = theShape->GetValue();
 830   if (aShape.IsNull()) return NULL;
 831
 832   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
 833
 834   TopTools_IndexedDataMapOfShapeListOfShape mapFaceBlocks;
 835   GEOMImpl_Block6Explorer::MapShapesAndAncestors
 836     (aShape, TopAbs_FACE, TopAbs_SOLID, mapFaceBlocks);
 837
 838   Standard_Integer ind = 1, nbFaces = mapFaceBlocks.Extent();
 839
 840   if (nbFaces == 0) {
 841     SetErrorCode("The given shape has no faces");
 842     return aSeq;
 843   }
 844
 845   TopTools_IndexedMapOfShape anIndices;
 846   TopExp::MapShapes(aShape, anIndices);
 847
 848   Standard_Integer id;
 849   for (; ind <= nbFaces; ind++) {
 850     if (mapFaceBlocks.FindFromIndex(ind).Extent() != 2) {
 851       id = anIndices.FindIndex(mapFaceBlocks.FindKey(ind));
 852       aSeq->Append(id);
 853     }
 854   }
 855
 856   //The explode doesn't change object so no new function is required.
 857   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 858   TCollection_AsciiString anOldDescr = aFunction->GetDescription();
 859
 860   //Make a Python command
 861   GEOM::TPythonDump(aFunction) << anOldDescr.ToCString()
 862     << "\n\tlistFreeFacesIDs = geompy.GetFreeFacesIDs(" << theShape << ")";
 863
 864   SetErrorCode(OK);
 865   return aSeq;
 866 }
 867
 868 //=======================================================================
 869 //function : GetSharedShapes
 870 //purpose  :
 871 //=======================================================================
 872
 873 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetSharedShapes
 874                                                 (Handle(GEOM_Object)    theShape1,
 875                                                  Handle(GEOM_Object)    theShape2,
 876                                                  const Standard_Integer theShapeType)
 877 {
 878   SetErrorCode(KO);
 879
 880   if (theShape1.IsNull() || theShape2.IsNull()) return NULL;
 881
 882   TopoDS_Shape aShape1 = theShape1->GetValue();
 883   TopoDS_Shape aShape2 = theShape2->GetValue();
 884
 885   if (aShape1.IsNull() || aShape2.IsNull()) return NULL;
 886
 887   TopTools_IndexedMapOfShape anIndices;
 888   TopExp::MapShapes(aShape1, anIndices);
 889   Handle(TColStd_HArray1OfInteger) anArray;
 890
 891   TopTools_IndexedMapOfShape mapShape1;
 892   TopExp::MapShapes(aShape1, TopAbs_ShapeEnum(theShapeType), mapShape1);
 893
 894   Handle(GEOM_Object) anObj;
 895   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
 896   TCollection_AsciiString anAsciiList, anEntry;
 897
 898   TopTools_MapOfShape mapShape2;
 899   TopExp_Explorer exp (aShape2, TopAbs_ShapeEnum(theShapeType));
 900   for (; exp.More(); exp.Next()) {
 901     TopoDS_Shape aSS = exp.Current();
 902     if (mapShape2.Add(aSS) && mapShape1.Contains(aSS)) {
 903       anArray = new TColStd_HArray1OfInteger(1,1);
 904       anArray->SetValue(1, anIndices.FindIndex(aSS));
 905       anObj = GetEngine()->AddSubShape(theShape1, anArray);
 906       aSeq->Append(anObj);
 907
 908       // for python command
 909       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
 910       anAsciiList += anEntry;
 911       anAsciiList += ",";
 912     }
 913   }
 914
 915   if (aSeq->IsEmpty()) {
 916     SetErrorCode("The given shapes have no shared sub-shapes of the requested type");
 917     return aSeq;
 918   }
 919
 920   //Make a Python command
 921   anAsciiList.Trunc(anAsciiList.Length() - 1);
 922
 923   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
 924
 925   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
 926     << "] = geompy.GetSharedShapes(" << theShape1 << ", "
 927       << theShape2 << ", " << theShapeType << ")";
 928
 929   SetErrorCode(OK);
 930   return aSeq;
 931 }
 932
 933 //=============================================================================
 934 /*!
 935  *
 936  */
 937 //=============================================================================
 938 static GEOM::TPythonDump& operator<< (GEOM::TPythonDump&   theDump,
 939                                       const GEOMAlgo_State theState)
 940 {
 941   switch (theState) {
 942   case GEOMAlgo_ST_IN:
 943     theDump << "geompy.GEOM.ST_IN";
 944     break;
 945   case GEOMAlgo_ST_OUT:
 946     theDump << "geompy.GEOM.ST_OUT";
 947     break;
 948   case GEOMAlgo_ST_ON:
 949     theDump << "geompy.GEOM.ST_ON";
 950     break;
 951   case GEOMAlgo_ST_ONIN:
 952     theDump << "geompy.GEOM.ST_ONIN";
 953     break;
 954   case GEOMAlgo_ST_ONOUT:
 955     theDump << "geompy.GEOM.ST_ONOUT";
 956     break;
 957   default:
 958     theDump << "geompy.GEOM.ST_UNKNOWN";
 959     break;
 960   }
 961   return theDump;
 962 }
 963
 964 //=======================================================================
 965 //function : checkTypeShapesOn
 966 /*!
 967  * \brief Checks if theShapeType parameter of GetShapesOnXXX() is OK
 968  * \param theShapeType - the shape type to check
 969  * \retval bool  - result of the check
 970  */
 971 //=======================================================================
 972
 973 bool GEOMImpl_IShapesOperations::checkTypeShapesOn(const Standard_Integer theShapeType)
 974 {
 975   if (theShapeType != TopAbs_VERTEX &&
 976       theShapeType != TopAbs_EDGE &&
 977       theShapeType != TopAbs_FACE &&
 978       theShapeType != TopAbs_SOLID) {
 979     SetErrorCode("Only solids, vertices, edges or faces can be found by this method");
 980     return false;
 981   }
 982   return true;
 983 }
 984
 985 //=======================================================================
 986 //function : makePlane
 987   /*!
 988    * \brief Creates Geom_Plane
 989     * \param theAx1 - shape object defining plane parameters
 990     * \retval Handle(Geom_Surface) - resulting surface
 991    */
 992 //=======================================================================
 993
 994 Handle(Geom_Surface) GEOMImpl_IShapesOperations::makePlane(const TopoDS_Shape& anAx1)
 995 {
 996   if (anAx1.ShapeType() != TopAbs_EDGE) return NULL;
 997   TopoDS_Edge anEdge = TopoDS::Edge(anAx1);
 998   TopoDS_Vertex V1, V2;
 999   TopExp::Vertices(anEdge, V1, V2, Standard_True);
1000   if (V1.IsNull() || V2.IsNull()) {
1001     SetErrorCode("Bad edge given for the plane normal vector");
1002     return NULL;
1003   }
1004   gp_Pnt aLoc = BRep_Tool::Pnt(V1);
1005   gp_Vec aVec (aLoc, BRep_Tool::Pnt(V2));
1006   if (aVec.Magnitude() < Precision::Confusion()) {
1007     SetErrorCode("Vector with null magnitude given");
1008     return NULL;
1009   }
1010   return new Geom_Plane(aLoc, aVec);
1011 }
1012
1013 //=======================================================================
1014 //function : makeCylinder
1015   /*!
1016    * \brief Creates Geom_CylindricalSurface
1017     * \param theAx1 - edge defining cylinder axis
1018     * \param theRadius - cylinder radius
1019     * \retval Handle(Geom_Surface) - resulting surface
1020    */
1021 //=======================================================================
1022
1023 Handle(Geom_Surface) GEOMImpl_IShapesOperations::makeCylinder(const TopoDS_Shape& anAxis,
1024                                                               const Standard_Real theRadius)
1025 {
1026   //Axis of the cylinder
1027   if (anAxis.ShapeType() != TopAbs_EDGE) {
1028     SetErrorCode("Not an edge given for the axis");
1029     return NULL;
1030   }
1031   TopoDS_Edge anEdge = TopoDS::Edge(anAxis);
1032   TopoDS_Vertex V1, V2;
1033   TopExp::Vertices(anEdge, V1, V2, Standard_True);
1034   if (V1.IsNull() || V2.IsNull()) {
1035     SetErrorCode("Bad edge given for the axis");
1036     return NULL;
1037   }
1038   gp_Pnt aLoc = BRep_Tool::Pnt(V1);
1039   gp_Vec aVec (aLoc, BRep_Tool::Pnt(V2));
1040   if (aVec.Magnitude() < Precision::Confusion()) {
1041     SetErrorCode("Vector with null magnitude given");
1042     return NULL;
1043   }
1044
1045   gp_Ax3 anAx3 (aLoc, aVec);
1046   return new Geom_CylindricalSurface(anAx3, theRadius);
1047 }
1048
1049
1050 //=======================================================================
1051 //function : getShapesOnSurfaceIDs
1052   /*!
1053    * \brief Find IDs of subshapes complying with given status about surface
1054     * \param theSurface - the surface to check state of subshapes against
1055     * \param theShape - the shape to explore
1056     * \param theShapeType - type of subshape of theShape
1057     * \param theState - required state
1058     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1059    */
1060 //=======================================================================
1061
1062 Handle(TColStd_HSequenceOfInteger)
1063   GEOMImpl_IShapesOperations::getShapesOnSurfaceIDs(const Handle(Geom_Surface)& theSurface,
1064                                                     const TopoDS_Shape&         theShape,
1065                                                     TopAbs_ShapeEnum            theShapeType,
1066                                                     GEOMAlgo_State              theState)
1067 {
1068   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs;
1069 //  MESSAGE("--------------------------- GetShapesOnPlane phase 1 takes:");
1070 //  OSD_Timer timer1;
1071 //  timer1.Start();
1072
1073   // Check presence of triangulation, build if need
1074   if (!CheckTriangulation(theShape))
1075     return aSeqOfIDs;
1076
1077   // Call algo
1078   GEOMAlgo_FinderShapeOn1 aFinder;
1079   Standard_Real aTol = 0.0001; // default value
1080
1081   aFinder.SetShape(theShape);
1082   aFinder.SetTolerance(aTol);
1083   aFinder.SetSurface(theSurface);
1084   aFinder.SetShapeType(theShapeType);
1085   aFinder.SetState(theState);
1086
1087   // Sets the minimal number of inner points for the faces that do not have own
1088   // inner points at all (for e.g. rectangular planar faces have just 2 triangles).
1089   // Default value=3
1090   aFinder.SetNbPntsMin(3);
1091   // Sets the maximal number of inner points for edges or faces.
1092   // It is usefull for the cases when this number is very big (e.g =2000) to improve
1093   // the performance. If this value =0, all inner points will be taken into account.
1094   // Default value=0
1095   aFinder.SetNbPntsMax(100);
1096
1097 //  timer1.Stop();
1098 //  timer1.Show();
1099
1100 //  MESSAGE("--------------------------- Perform on Plane takes:");
1101 //  timer1.Reset();
1102 //  timer1.Start();
1103   aFinder.Perform();
1104 //  timer1.Stop();
1105 //  timer1.Show();
1106
1107 //  MESSAGE("--------------------------- GetShapesOnPlane phase 3 takes:");
1108 //  timer1.Reset();
1109 //  timer1.Start();
1110
1111   // Interprete results
1112   Standard_Integer iErr = aFinder.ErrorStatus();
1113   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
1114   if (iErr) {
1115     MESSAGE(" iErr : " << iErr);
1116     TCollection_AsciiString aMsg (" iErr : ");
1117     aMsg += TCollection_AsciiString(iErr);
1118     SetErrorCode(aMsg);
1119     return aSeqOfIDs;
1120   }
1121   Standard_Integer iWrn = aFinder.WarningStatus();
1122   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
1123   if (iWrn) {
1124     MESSAGE(" *** iWrn : " << iWrn);
1125   }
1126
1127   const TopTools_ListOfShape& listSS = aFinder.Shapes(); // the result
1128
1129   if (listSS.Extent() < 1) {
1130     SetErrorCode("Not a single sub-shape of the requested type found on the given surface");
1131     return aSeqOfIDs;
1132   }
1133
1134 //  timer1.Stop();
1135 //  timer1.Show();
1136
1137 //  MESSAGE("--------------------------- GetShapesOnPlane phase 4 takes:");
1138 //  timer1.Reset();
1139 //  timer1.Start();
1140
1141   // Fill sequence of object IDs
1142   aSeqOfIDs = new TColStd_HSequenceOfInteger;
1143
1144   TopTools_IndexedMapOfShape anIndices;
1145   TopExp::MapShapes(theShape, anIndices);
1146
1147   TopTools_ListIteratorOfListOfShape itSub (listSS);
1148   for (int index = 1; itSub.More(); itSub.Next(), ++index) {
1149     int id = anIndices.FindIndex(itSub.Value());
1150     aSeqOfIDs->Append(id);
1151   }
1152 //  timer1.Stop();
1153 //  timer1.Show();
1154   return aSeqOfIDs;
1155 }
1156
1157 //=======================================================================
1158 //function : getObjectsShapesOn
1159 /*!
1160  * \brief Find shape objects and their entries by their ids
1161  * \param theShapeIDs - incoming shape ids
1162  * \param theShapeEntries - outgoing entries like "entry1, entry2, ..."
1163  * \retval Handle(TColStd_HSequenceOfTransient) - found shape objects
1164  */
1165 //=======================================================================
1166
1167 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::
1168  getObjectsShapesOn(const Handle(GEOM_Object)&                theShape,
1169                     const Handle(TColStd_HSequenceOfInteger)& theShapeIDs,
1170                     TCollection_AsciiString &                 theShapeEntries)
1171 {
1172   Handle(TColStd_HSequenceOfTransient) aSeq;
1173
1174   if ( !theShapeIDs.IsNull() && theShapeIDs->Length() > 0 )
1175   {
1176     aSeq = new TColStd_HSequenceOfTransient;
1177     Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1178     TCollection_AsciiString anEntry;
1179     for ( int i = 1; i <= theShapeIDs->Length(); ++i )
1180     {
1181       anArray->SetValue(1, theShapeIDs->Value( i ));
1182       Handle(GEOM_Object) anObj = GetEngine()->AddSubShape(theShape, anArray);
1183       aSeq->Append( anObj );
1184
1185       TDF_Tool::Entry(anObj->GetEntry(), anEntry);
1186       if ( i != 1 ) theShapeEntries += ",";
1187       theShapeEntries += anEntry;
1188     }
1189   }
1190   return aSeq;
1191 }
1192
1193 //=======================================================================
1194 //function : getShapesOnSurface
1195 /*!
1196    * \brief Find subshapes complying with given status about surface
1197     * \param theSurface - the surface to check state of subshapes against
1198     * \param theShape - the shape to explore
1199     * \param theShapeType - type of subshape of theShape
1200     * \param theState - required state
1201     * \param theShapeEntries - outgoing entries like "entry1, entry2, ..."
1202     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1203  */
1204 //=======================================================================
1205
1206 Handle(TColStd_HSequenceOfTransient)
1207     GEOMImpl_IShapesOperations::getShapesOnSurface(const Handle(Geom_Surface)& theSurface,
1208                                                    const Handle(GEOM_Object)&  theShape,
1209                                                    TopAbs_ShapeEnum            theShapeType,
1210                                                    GEOMAlgo_State              theState,
1211                                                    TCollection_AsciiString &   theShapeEntries)
1212 {
1213   // Find subshapes ids
1214   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1215     getShapesOnSurfaceIDs (theSurface, theShape->GetValue(), theShapeType, theState);
1216   if ( aSeqOfIDs.IsNull()  || aSeqOfIDs->Length() == 0 )
1217     return NULL;
1218
1219   return getObjectsShapesOn( theShape, aSeqOfIDs, theShapeEntries );
1220 }
1221
1222 //=============================================================================
1223 /*!
1224  *  GetShapesOnPlane
1225  */
1226 //=============================================================================
1227 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnPlane
1228                                         (const Handle(GEOM_Object)& theShape,
1229                                          const Standard_Integer     theShapeType,
1230                                          const Handle(GEOM_Object)& theAx1,
1231                                          const GEOMAlgo_State       theState)
1232 {
1233   SetErrorCode(KO);
1234
1235 //  MESSAGE("--------------------------- GetShapesOnPlane phase 1 takes:");
1236 //  OSD_Timer timer1;
1237 //  timer1.Start();
1238
1239   if (theShape.IsNull() || theAx1.IsNull()) return NULL;
1240
1241   TopoDS_Shape aShape = theShape->GetValue();
1242   TopoDS_Shape anAx1  = theAx1->GetValue();
1243
1244   if (aShape.IsNull() || anAx1.IsNull()) return NULL;
1245
1246   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1247   if ( !checkTypeShapesOn( theShapeType ))
1248     return NULL;
1249
1250   // Create plane
1251   Handle(Geom_Surface) aPlane = makePlane( anAx1 );
1252   if ( aPlane.IsNull() )
1253     return NULL;
1254
1255   // Find objects
1256   TCollection_AsciiString anAsciiList;
1257   Handle(TColStd_HSequenceOfTransient) aSeq;
1258   aSeq = getShapesOnSurface( aPlane, theShape, aShapeType, theState, anAsciiList );
1259   if ( aSeq.IsNull() || aSeq->Length() == 0 )
1260     return NULL;
1261
1262 //  timer1.Stop();
1263 //  timer1.Show();
1264
1265 //  MESSAGE("--------------------------- GetShapesOnPlane phase 5 takes:");
1266 //  timer1.Reset();
1267 //  timer1.Start();
1268
1269   // Make a Python command
1270
1271   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
1272   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1273
1274   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
1275     << "] = geompy.GetShapesOnPlane(" << theShape << ", "
1276       << theShapeType << ", " << theAx1 << ", " << theState << ")";
1277
1278   SetErrorCode(OK);
1279   return aSeq;
1280 }
1281
1282 //=============================================================================
1283 /*!
1284  *  GetShapesOnCylinder
1285  */
1286 //=============================================================================
1287 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnCylinder
1288                                           (const Handle(GEOM_Object)& theShape,
1289                                            const Standard_Integer     theShapeType,
1290                                            const Handle(GEOM_Object)& theAxis,
1291                                            const Standard_Real        theRadius,
1292                                            const GEOMAlgo_State       theState)
1293 {
1294   SetErrorCode(KO);
1295
1296   if (theShape.IsNull() || theAxis.IsNull()) return NULL;
1297
1298   TopoDS_Shape aShape = theShape->GetValue();
1299   TopoDS_Shape anAxis = theAxis->GetValue();
1300
1301   if (aShape.IsNull() || anAxis.IsNull()) return NULL;
1302
1303   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1304   if ( !checkTypeShapesOn( aShapeType ))
1305     return NULL;
1306
1307   // Create a cylinder surface
1308   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
1309   if ( aCylinder.IsNull() )
1310     return NULL;
1311
1312   // Find objects
1313   TCollection_AsciiString anAsciiList;
1314   Handle(TColStd_HSequenceOfTransient) aSeq;
1315   aSeq = getShapesOnSurface( aCylinder, theShape, aShapeType, theState, anAsciiList );
1316   if ( aSeq.IsNull() || aSeq->Length() == 0 )
1317     return NULL;
1318
1319   // Make a Python command
1320
1321   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
1322   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1323
1324   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
1325     << "] = geompy.GetShapesOnCylinder(" << theShape << ", " << theShapeType
1326       << ", " << theAxis << ", " << theRadius << ", " << theState << ")";
1327
1328   SetErrorCode(OK);
1329   return aSeq;
1330 }
1331
1332 //=============================================================================
1333 /*!
1334  *  GetShapesOnSphere
1335  */
1336 //=============================================================================
1337 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IShapesOperations::GetShapesOnSphere
1338                                           (const Handle(GEOM_Object)& theShape,
1339                                            const Standard_Integer     theShapeType,
1340                                            const Handle(GEOM_Object)& theCenter,
1341                                            const Standard_Real        theRadius,
1342                                            const GEOMAlgo_State       theState)
1343 {
1344   SetErrorCode(KO);
1345
1346   if (theShape.IsNull() || theCenter.IsNull()) return NULL;
1347
1348   TopoDS_Shape aShape  = theShape->GetValue();
1349   TopoDS_Shape aCenter = theCenter->GetValue();
1350
1351   if (aShape.IsNull() || aCenter.IsNull()) return NULL;
1352
1353   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1354   if ( !checkTypeShapesOn( aShapeType ))
1355     return NULL;
1356
1357   // Center of the sphere
1358   if (aCenter.ShapeType() != TopAbs_VERTEX) return NULL;
1359   gp_Pnt aLoc = BRep_Tool::Pnt(TopoDS::Vertex(aCenter));
1360
1361   gp_Ax3 anAx3 (aLoc, gp::DZ());
1362   Handle(Geom_SphericalSurface) aSphere =
1363     new Geom_SphericalSurface(anAx3, theRadius);
1364
1365   // Find objects
1366   TCollection_AsciiString anAsciiList;
1367   Handle(TColStd_HSequenceOfTransient) aSeq;
1368   aSeq = getShapesOnSurface( aSphere, theShape, aShapeType, theState, anAsciiList );
1369   if ( aSeq.IsNull() || aSeq->Length() == 0 )
1370     return NULL;
1371
1372   // Make a Python command
1373
1374   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
1375   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1376
1377   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
1378     << "] = geompy.GetShapesOnSphere(" << theShape << ", " << theShapeType
1379       << ", " << theCenter << ", " << theRadius << ", " << theState << ")";
1380
1381   SetErrorCode(OK);
1382   return aSeq;
1383 }
1384 //=============================================================================
1385 /*!
1386  *  GetShapesOnPlaneIDs
1387  */
1388 //=============================================================================
1389 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnPlaneIDs
1390                                         (const Handle(GEOM_Object)& theShape,
1391                                          const Standard_Integer     theShapeType,
1392                                          const Handle(GEOM_Object)& theAx1,
1393                                          const GEOMAlgo_State       theState)
1394 {
1395   SetErrorCode(KO);
1396
1397   if (theShape.IsNull() || theAx1.IsNull()) return NULL;
1398
1399   TopoDS_Shape aShape = theShape->GetValue();
1400   TopoDS_Shape anAx1  = theAx1->GetValue();
1401
1402   if (aShape.IsNull() || anAx1.IsNull()) return NULL;
1403
1404   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1405   if ( !checkTypeShapesOn( aShapeType ))
1406     return NULL;
1407
1408   // Create plane
1409   Handle(Geom_Surface) aPlane = makePlane( anAx1 );
1410   if ( aPlane.IsNull() )
1411     return NULL;
1412
1413   // Find object IDs
1414   Handle(TColStd_HSequenceOfInteger) aSeq;
1415   aSeq = getShapesOnSurfaceIDs( aPlane, aShape, aShapeType, theState );
1416
1417   // The GetShapesOnPlaneIDs() doesn't change object so no new function is required.
1418   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1419
1420   // Make a Python command
1421   const bool append = true;
1422   GEOM::TPythonDump(aFunction,append)
1423     << "listShapesOnPlane = IShapesOperations.GetShapesOnPlaneIDs"
1424     << "(" << theShape << "," << theShapeType << "," << theAx1 << "," << theState << ")";
1425
1426   SetErrorCode(OK);
1427   return aSeq;
1428 }
1429
1430 //=============================================================================
1431 /*!
1432  *  GetShapesOnCylinderIDs
1433  */
1434 //=============================================================================
1435 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnCylinderIDs
1436                                           (const Handle(GEOM_Object)& theShape,
1437                                            const Standard_Integer     theShapeType,
1438                                            const Handle(GEOM_Object)& theAxis,
1439                                            const Standard_Real        theRadius,
1440                                            const GEOMAlgo_State       theState)
1441 {
1442   SetErrorCode(KO);
1443
1444   if (theShape.IsNull() || theAxis.IsNull()) return NULL;
1445
1446   TopoDS_Shape aShape = theShape->GetValue();
1447   TopoDS_Shape anAxis = theAxis->GetValue();
1448
1449   if (aShape.IsNull() || anAxis.IsNull()) return NULL;
1450
1451   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1452   if ( !checkTypeShapesOn( aShapeType ))
1453     return NULL;
1454
1455   // Create a cylinder surface
1456   Handle(Geom_Surface) aCylinder = makeCylinder( anAxis, theRadius );
1457   if ( aCylinder.IsNull() )
1458     return NULL;
1459
1460   // Find object IDs
1461   Handle(TColStd_HSequenceOfInteger) aSeq;
1462   aSeq = getShapesOnSurfaceIDs( aCylinder, aShape, aShapeType, theState );
1463
1464   // The GetShapesOnCylinder() doesn't change object so no new function is required.
1465   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1466
1467   // Make a Python command
1468   const bool append = true;
1469   GEOM::TPythonDump(aFunction,append)
1470     << "listShapesOnCylinder = IShapesOperations.GetShapesOnCylinderIDs"
1471     << "(" << theShape << ", " << theShapeType << ", " << theAxis << ", "
1472     << theRadius << ", " << theState << ")";
1473
1474   SetErrorCode(OK);
1475   return aSeq;
1476 }
1477
1478 //=============================================================================
1479 /*!
1480  *  GetShapesOnSphereIDs
1481  */
1482 //=============================================================================
1483 Handle(TColStd_HSequenceOfInteger) GEOMImpl_IShapesOperations::GetShapesOnSphereIDs
1484                                           (const Handle(GEOM_Object)& theShape,
1485                                            const Standard_Integer     theShapeType,
1486                                            const Handle(GEOM_Object)& theCenter,
1487                                            const Standard_Real        theRadius,
1488                                            const GEOMAlgo_State       theState)
1489 {
1490   SetErrorCode(KO);
1491
1492   if (theShape.IsNull() || theCenter.IsNull()) return NULL;
1493
1494   TopoDS_Shape aShape  = theShape->GetValue();
1495   TopoDS_Shape aCenter = theCenter->GetValue();
1496
1497   if (aShape.IsNull() || aCenter.IsNull()) return NULL;
1498
1499   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1500   if ( !checkTypeShapesOn( aShapeType ))
1501     return NULL;
1502
1503   // Center of the sphere
1504   if (aCenter.ShapeType() != TopAbs_VERTEX) return NULL;
1505   gp_Pnt aLoc = BRep_Tool::Pnt(TopoDS::Vertex(aCenter));
1506
1507   gp_Ax3 anAx3 (aLoc, gp::DZ());
1508   Handle(Geom_SphericalSurface) aSphere =
1509     new Geom_SphericalSurface(anAx3, theRadius);
1510
1511   // Find object IDs
1512   Handle(TColStd_HSequenceOfInteger) aSeq;
1513   aSeq = getShapesOnSurfaceIDs( aSphere, aShape, aShapeType, theState );
1514
1515   // The GetShapesOnSphere() doesn't change object so no new function is required.
1516   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1517
1518   // Make a Python command
1519   const bool append = true;
1520   GEOM::TPythonDump(aFunction,append)
1521     << "listShapesOnCylinder = IShapesOperations.GetShapesOnCylinderIDs"
1522     << "(" << theShape << ", " << theShapeType << ", " << theCenter << ", "
1523     << theRadius << ", " << theState << ")";
1524
1525   SetErrorCode(OK);
1526   return aSeq;
1527 }
1528
1529 //=======================================================================
1530 //function : getShapesOnQuadrangleIDs
1531   /*!
1532    * \brief Find IDs of subshapes complying with given status about quadrangle
1533     * \param theShape - the shape to explore
1534     * \param theShapeType - type of subshape of theShape
1535     * \param theTopLeftPoint - top left quadrangle corner
1536     * \param theTopRigthPoint - top right quadrangle corner
1537     * \param theBottomLeftPoint - bottom left quadrangle corner
1538     * \param theBottomRigthPoint - bottom right quadrangle corner
1539     * \param theState - required state
1540     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1541    */
1542 //=======================================================================
1543
1544 Handle(TColStd_HSequenceOfInteger)
1545   GEOMImpl_IShapesOperations::getShapesOnQuadrangleIDs (const Handle(GEOM_Object)& theShape,
1546                                                         const Standard_Integer     theShapeType,
1547                                                         const Handle(GEOM_Object)& theTopLeftPoint,
1548                                                         const Handle(GEOM_Object)& theTopRigthPoint,
1549                                                         const Handle(GEOM_Object)& theBottomLeftPoint,
1550                                                         const Handle(GEOM_Object)& theBottomRigthPoint,
1551                                                         const GEOMAlgo_State       theState)
1552 {
1553   SetErrorCode(KO);
1554
1555   if ( theShape.IsNull() ||
1556        theTopLeftPoint.IsNull() ||
1557        theTopRigthPoint.IsNull() ||
1558        theBottomLeftPoint.IsNull() ||
1559        theBottomRigthPoint.IsNull() )
1560     return NULL;
1561
1562   TopoDS_Shape aShape = theShape->GetValue();
1563   TopoDS_Shape aTL = theTopLeftPoint->GetValue();
1564   TopoDS_Shape aTR = theTopRigthPoint->GetValue();
1565   TopoDS_Shape aBL = theBottomLeftPoint->GetValue();
1566   TopoDS_Shape aBR = theBottomRigthPoint->GetValue();
1567
1568   if (aShape.IsNull() ||
1569       aTL.IsNull() ||
1570       aTR.IsNull() ||
1571       aBL.IsNull() ||
1572       aBR.IsNull() ||
1573       aTL.ShapeType() != TopAbs_VERTEX ||
1574       aTR.ShapeType() != TopAbs_VERTEX ||
1575       aBL.ShapeType() != TopAbs_VERTEX ||
1576       aBR.ShapeType() != TopAbs_VERTEX )
1577     return NULL;
1578
1579   TopAbs_ShapeEnum aShapeType = TopAbs_ShapeEnum(theShapeType);
1580   if ( !checkTypeShapesOn( aShapeType ))
1581     return NULL;
1582
1583   vector< gp_Pnt > points(4);
1584   points[0] = BRep_Tool::Pnt(TopoDS::Vertex(aTL));
1585   points[1] = BRep_Tool::Pnt(TopoDS::Vertex(aTR));
1586   points[2] = BRep_Tool::Pnt(TopoDS::Vertex(aBR));
1587   points[3] = BRep_Tool::Pnt(TopoDS::Vertex(aBL));
1588
1589   // Algo: for each pair of neighboring point of a quadrangle, make a plane
1590   // and find subshape indices having theState. Then
1591   // - for IN state, find indices common for all pairs.
1592   // - else, keep IDs that are OK for any plane
1593   const bool keepOkOnAllPlanes = ( theState == GEOMAlgo_ST_IN ||
1594                                    theState == GEOMAlgo_ST_ONIN ||
1595                                    theState == GEOMAlgo_ST_INOUT );
1596
1597   // Find plane normal defined by corner points, it will be used to define a plane
1598   // for each pair of points. For that, find non coincident corner points
1599   vector< gp_Pnt > farPoints;
1600   farPoints.reserve( 5 );
1601   farPoints.push_back( points[0] );
1602   for ( int i = 1; i < 4; ++i )
1603   {
1604     // check if i-th point is far from all farPoints
1605     bool tooClose = false;
1606     vector< gp_Pnt >::iterator p = farPoints.begin();
1607     for ( ; p != farPoints.end(); ++p )
1608       if ( p->SquareDistance( points[ i ]) <= DBL_MIN )
1609         tooClose = true;
1610     if ( !tooClose )
1611       farPoints.push_back( points[ i ]);
1612   }
1613   if ( farPoints.size() < 3 ) {
1614     SetErrorCode("Coincident input points");
1615     return NULL;
1616   }
1617   gp_Vec aVecX =
1618     gp_Vec( farPoints[0], farPoints[1] ) ^ gp_Vec( farPoints[0], farPoints[2] );
1619   //std::cout << " X Vec : " << aVecX.X() << " " <<aVecX.Y() << " " <<aVecX.Z() << " " << endl;
1620
1621   // Use datamap to find IDs which have good state with all planes:
1622   // count nb of OK states for each ID
1623   TColStd_DataMapOfIntegerInteger nbOkStatesOfID;
1624
1625   // loop on point pairs
1626   int nbPlanes = 0;
1627   farPoints[ farPoints.size() ] = farPoints[ 0 ];
1628   for ( int i = 0; i < farPoints.size(); ++i )
1629   {
1630     // point1 -> point2 vector
1631     gp_Vec aVecY( farPoints[ i ], farPoints[ i + 1 ]);
1632     //std::cout << " Y Vec : " << aVecY.X() << " " <<aVecY.Y() << " " <<aVecY.Z() << " " << endl;
1633
1634   // plane normal
1635     gp_Vec aVecZ = aVecX ^ aVecY;
1636     //std::cout << " Z Vec : " << aVecZ.X() << " " <<aVecZ.Y() << " " <<aVecZ.Z() << " " << endl;
1637     if ( aVecZ.SquareMagnitude() <= DBL_MIN )
1638       continue;
1639
1640     // Check that normal direction is outside a quadrangle
1641     // (Suppose there are no concave corners in a quadrangle)
1642     int iPrev = i ? i - 1 : farPoints.size() - 1;
1643     if ( aVecZ * gp_Vec( farPoints[ i ], farPoints[ iPrev ]) >= 0. )
1644       aVecZ.Reverse();
1645
1646     ++nbPlanes;
1647     Handle(Geom_Plane) aPlane = new Geom_Plane( farPoints[ i ], aVecZ );
1648
1649     // Find subshape indices
1650     Handle(TColStd_HSequenceOfInteger) aSeq;
1651     aSeq = getShapesOnSurfaceIDs( aPlane, aShape, aShapeType, theState );
1652
1653     if ( aSeq.IsNull() || aSeq->Length() == 0 )
1654     {
1655       if ( keepOkOnAllPlanes )
1656         return NULL;
1657     }
1658     else
1659     {
1660       // put IDs to the datamap
1661       //std::cout << " IDS in plane " << nbPlanes << " : ";
1662       for ( int iID = 1; iID <= aSeq->Length(); ++iID )
1663       {
1664         int id = aSeq->Value( iID );
1665         //std::cout << id << " ";
1666         if ( nbOkStatesOfID.IsBound( id ))
1667           nbOkStatesOfID( id )++;
1668         else
1669           nbOkStatesOfID.Bind( id, 1 );
1670       }
1671       //std::cout << endl;
1672     }
1673   }
1674
1675   // select IDs that are OK with all planes
1676   Handle(TColStd_HSequenceOfInteger) aSeq = new TColStd_HSequenceOfInteger;
1677   TColStd_DataMapIteratorOfDataMapOfIntegerInteger id_nb;
1678   for ( id_nb.Initialize( nbOkStatesOfID ); id_nb.More(); id_nb.Next() )
1679   {
1680     //std::cout << id_nb.Key() << " in " << id_nb.Value() << " planes " << endl;
1681     if ( !keepOkOnAllPlanes || id_nb.Value() == nbPlanes )
1682       aSeq->Append( id_nb.Key() );
1683   }
1684   return aSeq;
1685 }
1686
1687 //=======================================================================
1688 //function : GetShapesOnQuadrangle
1689   /*!
1690    * \brief Find subshapes complying with given status about quadrangle
1691     * \param theShape - the shape to explore
1692     * \param theShapeType - type of subshape of theShape
1693     * \param theTopLeftPoint - top left quadrangle corner
1694     * \param theTopRigthPoint - top right quadrangle corner
1695     * \param theBottomLeftPoint - bottom left quadrangle corner
1696     * \param theBottomRigthPoint - bottom right quadrangle corner
1697     * \param theState - required state
1698     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1699    */
1700 //=======================================================================
1701
1702 Handle(TColStd_HSequenceOfTransient)
1703     GEOMImpl_IShapesOperations::GetShapesOnQuadrangle (const Handle(GEOM_Object)& theShape,
1704                                                        const Standard_Integer     theShapeType,
1705                                                        const Handle(GEOM_Object)& theTopLeftPoint,
1706                                                        const Handle(GEOM_Object)& theTopRigthPoint,
1707                                                        const Handle(GEOM_Object)& theBottomLeftPoint,
1708                                                        const Handle(GEOM_Object)& theBottomRigthPoint,
1709                                                        const GEOMAlgo_State       theState)
1710 {
1711   // Find indices
1712   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1713     getShapesOnQuadrangleIDs( theShape,
1714                               theShapeType,
1715                               theTopLeftPoint,
1716                               theTopRigthPoint,
1717                               theBottomLeftPoint,
1718                               theBottomRigthPoint,
1719                               theState);
1720   if ( aSeqOfIDs.IsNull() || aSeqOfIDs->IsEmpty() )
1721     return NULL;
1722
1723   // Find objects by indices
1724   TCollection_AsciiString anAsciiList;
1725   Handle(TColStd_HSequenceOfTransient) aSeq;
1726   aSeq = getObjectsShapesOn( theShape, aSeqOfIDs, anAsciiList );
1727   if ( aSeq.IsNull() || aSeq->IsEmpty() )
1728     return NULL;
1729
1730   // Make a Python command
1731
1732   Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast( aSeq->Value( 1 ));
1733   Handle(GEOM_Function) aFunction = anObj->GetLastFunction();
1734
1735   GEOM::TPythonDump(aFunction)
1736     << "[" << anAsciiList.ToCString() << "] = geompy.GetShapesOnQuadrangle("
1737     << theShape << ", "
1738     << theShapeType << ", "
1739     << theTopLeftPoint << ", "
1740     << theTopRigthPoint << ", "
1741     << theBottomLeftPoint << ", "
1742     << theBottomRigthPoint << ", "
1743     << theState << ")";
1744
1745   SetErrorCode(OK);
1746   return aSeq;
1747 }
1748
1749 //=======================================================================
1750 //function : GetShapesOnQuadrangleIDs
1751   /*!
1752    * \brief Find IDs of subshapes complying with given status about quadrangle
1753     * \param theShape - the shape to explore
1754     * \param theShapeType - type of subshape of theShape
1755     * \param theTopLeftPoint - top left quadrangle corner
1756     * \param theTopRigthPoint - top right quadrangle corner
1757     * \param theBottomLeftPoint - bottom left quadrangle corner
1758     * \param theBottomRigthPoint - bottom right quadrangle corner
1759     * \param theState - required state
1760     * \retval Handle(TColStd_HSequenceOfInteger) - IDs of found subshapes
1761    */
1762 //=======================================================================
1763
1764 Handle(TColStd_HSequenceOfInteger)
1765   GEOMImpl_IShapesOperations::GetShapesOnQuadrangleIDs (const Handle(GEOM_Object)& theShape,
1766                                                         const Standard_Integer     theShapeType,
1767                                                         const Handle(GEOM_Object)& theTopLeftPoint,
1768                                                         const Handle(GEOM_Object)& theTopRigthPoint,
1769                                                         const Handle(GEOM_Object)& theBottomLeftPoint,
1770                                                         const Handle(GEOM_Object)& theBottomRigthPoint,
1771                                                         const GEOMAlgo_State       theState)
1772 {
1773   // Find indices
1774   Handle(TColStd_HSequenceOfInteger) aSeqOfIDs =
1775     getShapesOnQuadrangleIDs( theShape,
1776                               theShapeType,
1777                               theTopLeftPoint,
1778                               theTopRigthPoint,
1779                               theBottomLeftPoint,
1780                               theBottomRigthPoint,
1781                               theState);
1782   if ( aSeqOfIDs.IsNull() || aSeqOfIDs->IsEmpty() )
1783     return NULL;
1784
1785   // Make a Python command
1786
1787   // The GetShapesOnCylinder() doesn't change object so no new function is required.
1788   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1789
1790   const bool append = true;
1791   GEOM::TPythonDump(aFunction,append)
1792     << "listShapesOnQuadrangle = geompy.GetShapesOnQuadrangleIDs("
1793     << theShape << ", "
1794     << theShapeType << ", "
1795     << theTopLeftPoint << ", "
1796     << theTopRigthPoint << ", "
1797     << theBottomLeftPoint << ", "
1798     << theBottomRigthPoint << ", "
1799     << theState << ")";
1800
1801   SetErrorCode(OK);
1802   return aSeqOfIDs;
1803 }
1804
1805
1806 //=============================================================================
1807 /*!
1808  *  GetInPlace
1809  */
1810 //=============================================================================
1811 Handle(GEOM_Object) GEOMImpl_IShapesOperations::GetInPlace
1812                                           (Handle(GEOM_Object) theShapeWhere,
1813                                            Handle(GEOM_Object) theShapeWhat)
1814 {
1815   SetErrorCode(KO);
1816
1817   if (theShapeWhere.IsNull() || theShapeWhat.IsNull()) return NULL;
1818
1819   TopoDS_Shape aWhere = theShapeWhere->GetValue();
1820   TopoDS_Shape aWhat  = theShapeWhat->GetValue();
1821
1822   if (aWhere.IsNull() || aWhat.IsNull()) return NULL;
1823
1824   //Fill array of indices
1825   Handle(TColStd_HArray1OfInteger) aModifiedArray;
1826
1827   Handle(GEOM_Function) aWhereFunction = theShapeWhere->GetLastFunction();
1828
1829   TopTools_IndexedMapOfShape aWhereIndices;
1830   TopExp::MapShapes(aWhere, aWhereIndices);
1831
1832   if (aWhereIndices.Contains(aWhat)) {
1833
1834     // entity was not changed by the operation
1835     Standard_Integer aWhatIndex = aWhereIndices.FindIndex(aWhat);
1836     aModifiedArray = new TColStd_HArray1OfInteger(1,1);
1837     aModifiedArray->SetValue(1, aWhatIndex);
1838
1839   } else {
1840
1841     TDF_Label aHistoryLabel = aWhereFunction->GetHistoryEntry(Standard_False);
1842     if (aHistoryLabel.IsNull()) {
1843       SetErrorCode("Modifications history does not exist for the shape under consideration.");
1844       return NULL;
1845     }
1846
1847     // search in history for all argument shapes
1848     Standard_Boolean isFound = Standard_False;
1849
1850     TDF_LabelSequence aLabelSeq;
1851     aWhereFunction->GetDependency(aLabelSeq);
1852     Standard_Integer nbArg = aLabelSeq.Length();
1853
1854     for (Standard_Integer iarg = 1; iarg <= nbArg && !isFound; iarg++) {
1855
1856       TDF_Label anArgumentRefLabel = aLabelSeq.Value(iarg);
1857
1858       Handle(GEOM_Object) anArgumentObject = GEOM_Object::GetReferencedObject(anArgumentRefLabel);
1859       TopoDS_Shape anArgumentShape = anArgumentObject->GetValue();
1860
1861       TopTools_IndexedMapOfShape anArgumentIndices;
1862       TopExp::MapShapes(anArgumentShape, anArgumentIndices);
1863
1864       if (anArgumentIndices.Contains(aWhat)) {
1865         isFound = Standard_True;
1866         Standard_Integer aWhatIndex = anArgumentIndices.FindIndex(aWhat);
1867
1868         // Find corresponding label in history
1869         TDF_Label anArgumentHistoryLabel =
1870           aWhereFunction->GetArgumentHistoryEntry(anArgumentRefLabel, Standard_False);
1871         if (anArgumentHistoryLabel.IsNull()) {
1872           // Lost History of operation argument. Possibly, all its entities was removed.
1873           SetErrorCode(OK);
1874           return NULL;
1875         }
1876
1877         TDF_Label aWhatHistoryLabel = anArgumentHistoryLabel.FindChild(aWhatIndex, Standard_False);
1878         if (aWhatHistoryLabel.IsNull()) {
1879           // Removed entity
1880           SetErrorCode(OK);
1881           return NULL;
1882         }
1883
1884         Handle(TDataStd_IntegerArray) anIntegerArray;
1885         if (!aWhatHistoryLabel.FindAttribute(TDataStd_IntegerArray::GetID(), anIntegerArray)) {
1886           SetErrorCode("Error: Empty modifications history for the sought shape.");
1887           return NULL;
1888         }
1889
1890         aModifiedArray = anIntegerArray->Array();
1891         if (aModifiedArray->Length() == 0) {
1892           SetErrorCode("Error: Empty modifications history for the sought shape.");
1893           return NULL;
1894         }
1895       }
1896     }
1897
1898     if (!isFound) {
1899       SetErrorCode("The sought shape does not belong to any operation argument.");
1900       return NULL;
1901     }
1902   }
1903
1904   //Add a new object
1905   Handle(GEOM_Object) aResult = GetEngine()->AddSubShape(theShapeWhere, aModifiedArray);
1906
1907   if (aModifiedArray->Length() > 1) {
1908     //Set a GROUP type
1909     aResult->SetType(GEOM_GROUP);
1910
1911     //Set a sub shape type
1912     TDF_Label aFreeLabel = aResult->GetFreeLabel();
1913     TopAbs_ShapeEnum aShapeType = aWhat.ShapeType();
1914     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)aShapeType);
1915   }
1916
1917   //Make a Python command
1918   Handle(GEOM_Function) aFunction = aResult->GetFunction(1);
1919
1920   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetInPlace("
1921     << theShapeWhere << ", " << theShapeWhat << ")";
1922
1923   SetErrorCode(OK);
1924   return aResult;
1925 }
1926
1927 //=======================================================================
1928 //function : SortShapes
1929 //purpose  :
1930 //=======================================================================
1931 void GEOMImpl_IShapesOperations::SortShapes(TopTools_ListOfShape& SL)
1932 {
1933   Standard_Integer MaxShapes = SL.Extent();
1934   TopTools_Array1OfShape  aShapes (1,MaxShapes);
1935   TColStd_Array1OfInteger OrderInd(1,MaxShapes);
1936   TColStd_Array1OfReal    MidXYZ  (1,MaxShapes); //X,Y,Z;
1937   TColStd_Array1OfReal    Length  (1,MaxShapes); //X,Y,Z;
1938
1939   // Computing of CentreOfMass
1940   Standard_Integer Index;
1941   GProp_GProps GPr;
1942   gp_Pnt GPoint;
1943   TopTools_ListIteratorOfListOfShape it(SL);
1944   for (Index=1;  it.More();  Index++)
1945   {
1946     TopoDS_Shape S = it.Value();
1947     SL.Remove( it ); // == it.Next()
1948     aShapes(Index) = S;
1949     OrderInd.SetValue (Index, Index);
1950     if (S.ShapeType() == TopAbs_VERTEX)
1951     {
1952       GPoint = BRep_Tool::Pnt( TopoDS::Vertex( S ));
1953       Length.SetValue( Index, (Standard_Real) S.Orientation());
1954     }
1955     else
1956     {
1957       BRepGProp::LinearProperties (S, GPr);
1958       GPoint = GPr.CentreOfMass();
1959       Length.SetValue( Index, GPr.Mass() );
1960     }
1961     MidXYZ.SetValue(Index,
1962                     GPoint.X()*999 + GPoint.Y()*99 + GPoint.Z()*0.9);
1963   }
1964   // Sorting
1965   Standard_Integer aTemp;
1966   Standard_Boolean exchange, Sort = Standard_True;
1967   while (Sort)
1968   {
1969     Sort = Standard_False;
1970     for (Index=1; Index < MaxShapes; Index++)
1971     {
1972       if (MidXYZ(OrderInd(Index)) > MidXYZ(OrderInd(Index+1)))
1973         exchange = Standard_True;
1974       else if (MidXYZ(OrderInd(Index)) == MidXYZ(OrderInd(Index+1)) &&
1975                Length(OrderInd(Index)) >  Length(OrderInd(Index+1)) )
1976         exchange = Standard_True;
1977       else
1978         exchange = Standard_False;
1979       if (exchange)
1980       {
1981         aTemp = OrderInd(Index);
1982         OrderInd(Index) = OrderInd(Index+1);
1983         OrderInd(Index+1) = aTemp;
1984         Sort = Standard_True;
1985       }
1986     }
1987   }
1988   for (Index=1; Index <= MaxShapes; Index++)
1989     SL.Append( aShapes( OrderInd(Index) ));
1990 }
1991
1992 //=======================================================================
1993 //function : CheckTriangulation
1994 //purpose  :
1995 //=======================================================================
1996 bool GEOMImpl_IShapesOperations::CheckTriangulation (const TopoDS_Shape& aShape)
1997 {
1998 //  MESSAGE("CheckTriangulation");
1999 //
2000 //  OSD_Timer timer1;
2001 //  timer1.Start();
2002
2003   TopExp_Explorer exp (aShape, TopAbs_FACE);
2004   if (!exp.More()) {
2005     SetErrorCode("Shape without faces given");
2006     return false;
2007   }
2008
2009   TopLoc_Location aTopLoc;
2010   Handle(Poly_Triangulation) aTRF;
2011   aTRF = BRep_Tool::Triangulation(TopoDS::Face(exp.Current()), aTopLoc);
2012   if (aTRF.IsNull()) {
2013     // calculate deflection
2014     Standard_Real aDeviationCoefficient = 0.001;
2015
2016     Bnd_Box B;
2017     BRepBndLib::Add(aShape, B);
2018     Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
2019     B.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
2020
2021     Standard_Real dx = aXmax - aXmin, dy = aYmax - aYmin, dz = aZmax - aZmin;
2022     Standard_Real aDeflection = Max(Max(dx, dy), dz) * aDeviationCoefficient * 4;
2023
2024 //    MESSAGE("Deflection = " << aDeflection);
2025
2026     Standard_Real aHLRAngle = 0.349066;
2027
2028     BRepMesh_IncrementalMesh Inc (aShape, aDeflection, Standard_False, aHLRAngle);
2029   }
2030 //  timer1.Stop();
2031 //  timer1.Show();
2032
2033   return true;
2034 }