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