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