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