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