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