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