]> SALOME platform Git repositories - modules/geom.git/blob - src/GEOMImpl/GEOMImpl_IMeasureOperations.cxx
Salome HOME
0023247: [CEA 1771] Problem in computing dimensions of objects in GEOM module
[modules/geom.git] / src / GEOMImpl / GEOMImpl_IMeasureOperations.cxx
1 // Copyright (C) 2007-2016  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
5 //
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
10 //
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14 // Lesser General Public License for more details.
15 //
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
19 //
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 //
22
23 #include <Basics_OCCTVersion.hxx>
24
25 #include <GEOMImpl_IMeasureOperations.hxx>
26 #include <GEOMImpl_IMeasure.hxx>
27 #include <GEOMImpl_MeasureDriver.hxx>
28 #include <GEOMImpl_Types.hxx>
29
30 #include <GEOMUtils.hxx>
31
32 #include <GEOMAlgo_AlgoTools.hxx>
33 #include <GEOMAlgo_KindOfName.hxx>
34 #include <GEOMAlgo_ShapeInfoFiller.hxx>
35
36 #include <GEOM_PythonDump.hxx>
37
38 #include <utilities.h>
39
40 // OCCT Includes
41 #include <Bnd_Box.hxx>
42 #include <BOPAlgo_CheckerSI.hxx>
43 #include <BOPCol_ListOfShape.hxx>
44 #include <BOPDS_DS.hxx>
45 #include <BOPDS_MapOfPassKey.hxx>
46 #include <BRepBndLib.hxx>
47 #include <BRepBuilderAPI_Copy.hxx>
48 #include <BRepCheck_ListIteratorOfListOfStatus.hxx>
49 #include <BRepCheck_Shell.hxx>
50 #include <BRepClass3d_SolidClassifier.hxx>
51 #include <BRepClass_FaceClassifier.hxx>
52 #include <BRepExtrema_DistShapeShape.hxx>
53 #include <BRepExtrema_ShapeProximity.hxx>
54 #if OCC_VERSION_LARGE > 0x06090000
55 #include <BRepExtrema_SelfIntersection.hxx>
56 #include <BRepExtrema_MapOfIntegerPackedMapOfInteger.hxx>
57 #endif
58 #include <BRepGProp.hxx>
59 #include <BRepTools.hxx>
60 #include <BRep_Tool.hxx>
61 #include <Geom_Line.hxx>
62 #include <GeomAPI_ProjectPointOnCurve.hxx>
63 #include <GeomAPI_ProjectPointOnSurf.hxx>
64 #include <GeomLProp_CLProps.hxx>
65 #include <GeomLProp_SLProps.hxx>
66 #include <Geom_Plane.hxx>
67 #include <GProp_GProps.hxx>
68 #include <GProp_PrincipalProps.hxx>
69 #include <ShapeAnalysis.hxx>
70 #include <ShapeAnalysis_Surface.hxx>
71 #include <TopExp.hxx>
72 #include <TopExp_Explorer.hxx>
73 #include <TopoDS.hxx>
74 #include <TopoDS_Edge.hxx>
75 #include <TopTools_IndexedMapOfShape.hxx>
76 #include <TopTools_DataMapIteratorOfDataMapOfIntegerListOfShape.hxx>
77 #include <TColStd_MapIteratorOfPackedMapOfInteger.hxx>
78 #include <TopTools_ListIteratorOfListOfShape.hxx>
79 #include <TopTools_ListOfShape.hxx>
80 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
81
82 #include <set>
83
84 //=============================================================================
85 /*!
86  *  Constructor
87  */
88 //=============================================================================
89 GEOMImpl_IMeasureOperations::GEOMImpl_IMeasureOperations (GEOM_Engine* theEngine, int theDocID)
90 : GEOM_IOperations(theEngine, theDocID)
91 {
92   MESSAGE("GEOMImpl_IMeasureOperations::GEOMImpl_IMeasureOperations");
93 }
94
95 //=============================================================================
96 /*!
97  *  Destructor
98  */
99 //=============================================================================
100 GEOMImpl_IMeasureOperations::~GEOMImpl_IMeasureOperations()
101 {
102   MESSAGE("GEOMImpl_IMeasureOperations::~GEOMImpl_IMeasureOperations");
103 }
104
105 //=============================================================================
106 /*! Get kind and parameters of the given shape.
107  */
108 //=============================================================================
109 GEOMImpl_IMeasureOperations::ShapeKind GEOMImpl_IMeasureOperations::KindOfShape
110                              (Handle(GEOM_Object) theShape,
111                               Handle(TColStd_HSequenceOfInteger)& theIntegers,
112                               Handle(TColStd_HSequenceOfReal)&    theDoubles)
113 {
114   SetErrorCode(KO);
115   ShapeKind aKind = SK_NO_SHAPE;
116
117   if (theIntegers.IsNull()) theIntegers = new TColStd_HSequenceOfInteger;
118   else                      theIntegers->Clear();
119
120   if (theDoubles.IsNull()) theDoubles = new TColStd_HSequenceOfReal;
121   else                     theDoubles->Clear();
122
123   if (theShape.IsNull())
124     return aKind;
125
126   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
127   if (aRefShape.IsNull()) return aKind;
128
129   TopoDS_Shape aShape = aRefShape->GetValue();
130   if (aShape.IsNull()) return aKind;
131
132   int geom_type = theShape->GetType();
133
134   // check if it's advanced shape
135   if ( geom_type > USER_TYPE ) {
136     SetErrorCode(OK);
137     return SK_ADVANCED;
138   }
139
140   // Call algorithm
141   GEOMAlgo_ShapeInfoFiller aSF;
142   aSF.SetShape(aShape);
143   aSF.Perform();
144   Standard_Integer iErr = aSF.ErrorStatus();
145   if (iErr) {
146     SetErrorCode("Error in GEOMAlgo_ShapeInfoFiller");
147     return SK_NO_SHAPE;
148   }
149   const GEOMAlgo_ShapeInfo& anInfo = aSF.Info();
150
151   // specific processing for some "advandced" objects
152   switch ( geom_type ) {
153   case GEOM_MARKER:
154     // local coordinate systen
155     // (+) geompy.kind.LCS  xc yc zc xx xy xz yx yy yz zx zy zz
156
157     TopoDS_Face aFace = TopoDS::Face( aShape );
158     Handle(Geom_Plane) aPlane = Handle(Geom_Plane)::DownCast( BRep_Tool::Surface( aFace ) );
159     gp_Pnt aC = aPlane->Pln().Location();
160     gp_Ax3 anAx3 = aPlane->Pln().Position();
161
162     theDoubles->Append(aC.X());
163     theDoubles->Append(aC.Y());
164     theDoubles->Append(aC.Z());
165     
166     gp_Dir aD = anAx3.XDirection();
167     theDoubles->Append(aD.X());
168     theDoubles->Append(aD.Y());
169     theDoubles->Append(aD.Z());
170     aD = anAx3.YDirection();
171     theDoubles->Append(aD.X());
172     theDoubles->Append(aD.Y());
173     theDoubles->Append(aD.Z());
174     aD = anAx3.Direction();
175     theDoubles->Append(aD.X());
176     theDoubles->Append(aD.Y());
177     theDoubles->Append(aD.Z());
178
179     SetErrorCode(OK);
180     return SK_LCS;
181   }
182
183   // Interprete results
184   TopAbs_ShapeEnum aType = anInfo.Type();
185   switch (aType)
186   {
187   case TopAbs_COMPOUND:
188   case TopAbs_COMPSOLID:
189     {
190       // (+) geompy.kind.COMPOUND     nb_solids nb_faces nb_edges nb_vertices
191       // (+) geompy.kind.COMPSOLID    nb_solids nb_faces nb_edges nb_vertices
192       // ??? "nb_faces" - all faces or only 'standalone' faces?
193       if (aType == TopAbs_COMPOUND)
194         aKind = SK_COMPOUND;
195       else
196         aKind = SK_COMPSOLID;
197
198       //theIntegers->Append(anInfo.NbSubShapes(TopAbs_COMPOUND));
199       //theIntegers->Append(anInfo.NbSubShapes(TopAbs_COMPSOLID));
200       theIntegers->Append(anInfo.NbSubShapes(TopAbs_SOLID));
201       theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
202       theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
203       theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
204     }
205     break;
206
207   case TopAbs_SHELL:
208     {
209       // (+) geompy.kind.SHELL  geompy.info.closed   nb_faces nb_edges nb_vertices
210       // (+) geompy.kind.SHELL  geompy.info.unclosed nb_faces nb_edges nb_vertices
211       aKind = SK_SHELL;
212
213       theIntegers->Append((int)anInfo.KindOfClosed());
214
215       theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
216       theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
217       theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
218     }
219     break;
220
221   case TopAbs_WIRE:
222     {
223       // (+) geompy.kind.WIRE  geompy.info.closed   nb_edges nb_vertices
224       // (+) geompy.kind.WIRE  geompy.info.unclosed nb_edges nb_vertices
225       aKind = SK_WIRE;
226
227       theIntegers->Append((int)anInfo.KindOfClosed());
228
229       theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
230       theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
231     }
232     break;
233
234   case TopAbs_SOLID:
235     {
236       aKind = SK_SOLID;
237
238       GEOMAlgo_KindOfName aKN = anInfo.KindOfName();
239       switch (aKN)
240       {
241       case GEOMAlgo_KN_SPHERE:
242         // (+) geompy.kind.SPHERE  xc yc zc  R
243         {
244           aKind = SK_SPHERE;
245
246           gp_Pnt aC = anInfo.Location();
247           theDoubles->Append(aC.X());
248           theDoubles->Append(aC.Y());
249           theDoubles->Append(aC.Z());
250
251           theDoubles->Append(anInfo.Radius1());
252         }
253         break;
254       case GEOMAlgo_KN_CYLINDER:
255         // (+) geompy.kind.CYLINDER  xb yb zb  dx dy dz  R  H
256         {
257           aKind = SK_CYLINDER;
258
259           gp_Pnt aC = anInfo.Location();
260           theDoubles->Append(aC.X());
261           theDoubles->Append(aC.Y());
262           theDoubles->Append(aC.Z());
263
264           gp_Ax3 anAx3 = anInfo.Position();
265           gp_Dir aD = anAx3.Direction();
266           theDoubles->Append(aD.X());
267           theDoubles->Append(aD.Y());
268           theDoubles->Append(aD.Z());
269
270           theDoubles->Append(anInfo.Radius1());
271           theDoubles->Append(anInfo.Height());
272         }
273         break;
274       case GEOMAlgo_KN_BOX:
275         // (+) geompy.kind.BOX  xc yc zc  ax ay az
276         {
277           aKind = SK_BOX;
278
279           gp_Pnt aC = anInfo.Location();
280           theDoubles->Append(aC.X());
281           theDoubles->Append(aC.Y());
282           theDoubles->Append(aC.Z());
283
284           gp_Ax3 anAx3 = anInfo.Position();
285           gp_Dir aD = anAx3.Direction();
286           gp_Dir aX = anAx3.XDirection();
287
288           // ax ay az
289           if (aD.IsParallel(gp::DZ(), Precision::Angular()) &&
290               aX.IsParallel(gp::DX(), Precision::Angular())) {
291             theDoubles->Append(anInfo.Length()); // ax'
292             theDoubles->Append(anInfo.Width());  // ay'
293             theDoubles->Append(anInfo.Height()); // az'
294           }
295           else if (aD.IsParallel(gp::DZ(), Precision::Angular()) &&
296                    aX.IsParallel(gp::DY(), Precision::Angular())) {
297             theDoubles->Append(anInfo.Width());  // ay'
298             theDoubles->Append(anInfo.Length()); // ax'
299             theDoubles->Append(anInfo.Height()); // az'
300           }
301           else if (aD.IsParallel(gp::DX(), Precision::Angular()) &&
302                    aX.IsParallel(gp::DZ(), Precision::Angular())) {
303             theDoubles->Append(anInfo.Height()); // az'
304             theDoubles->Append(anInfo.Width());  // ay'
305             theDoubles->Append(anInfo.Length()); // ax'
306           }
307           else if (aD.IsParallel(gp::DX(), Precision::Angular()) &&
308                    aX.IsParallel(gp::DY(), Precision::Angular())) {
309             theDoubles->Append(anInfo.Height()); // az'
310             theDoubles->Append(anInfo.Length()); // ax'
311             theDoubles->Append(anInfo.Width());  // ay'
312           }
313           else if (aD.IsParallel(gp::DY(), Precision::Angular()) &&
314                    aX.IsParallel(gp::DZ(), Precision::Angular())) {
315             theDoubles->Append(anInfo.Width());  // ay'
316             theDoubles->Append(anInfo.Height()); // az'
317             theDoubles->Append(anInfo.Length()); // ax'
318           }
319           else if (aD.IsParallel(gp::DY(), Precision::Angular()) &&
320                    aX.IsParallel(gp::DX(), Precision::Angular())) {
321             theDoubles->Append(anInfo.Length()); // ax'
322             theDoubles->Append(anInfo.Height()); // az'
323             theDoubles->Append(anInfo.Width());  // ay'
324           }
325           else {
326             // (+) geompy.kind.ROTATED_BOX  xo yo zo  zx zy zz  xx xy xz  ax ay az
327             aKind = SK_ROTATED_BOX;
328
329             // Direction and XDirection
330             theDoubles->Append(aD.X());
331             theDoubles->Append(aD.Y());
332             theDoubles->Append(aD.Z());
333
334             theDoubles->Append(aX.X());
335             theDoubles->Append(aX.Y());
336             theDoubles->Append(aX.Z());
337
338             // ax ay az
339             theDoubles->Append(anInfo.Length());
340             theDoubles->Append(anInfo.Width());
341             theDoubles->Append(anInfo.Height());
342           }
343         }
344         break;
345       case GEOMAlgo_KN_TORUS:
346         // (+) geompy.kind.TORUS  xc yc zc  dx dy dz  R_1 R_2
347         {
348           aKind = SK_TORUS;
349
350           gp_Pnt aO = anInfo.Location();
351           theDoubles->Append(aO.X());
352           theDoubles->Append(aO.Y());
353           theDoubles->Append(aO.Z());
354
355           gp_Ax3 anAx3 = anInfo.Position();
356           gp_Dir aD = anAx3.Direction();
357           theDoubles->Append(aD.X());
358           theDoubles->Append(aD.Y());
359           theDoubles->Append(aD.Z());
360
361           theDoubles->Append(anInfo.Radius1());
362           theDoubles->Append(anInfo.Radius2());
363         }
364         break;
365       case GEOMAlgo_KN_CONE:
366         // (+) geompy.kind.CONE  xb yb zb  dx dy dz  R_1 R_2  H
367         {
368           aKind = SK_CONE;
369
370           gp_Pnt aO = anInfo.Location();
371           theDoubles->Append(aO.X());
372           theDoubles->Append(aO.Y());
373           theDoubles->Append(aO.Z());
374
375           gp_Ax3 anAx3 = anInfo.Position();
376           gp_Dir aD = anAx3.Direction();
377           theDoubles->Append(aD.X());
378           theDoubles->Append(aD.Y());
379           theDoubles->Append(aD.Z());
380
381           theDoubles->Append(anInfo.Radius1());
382           theDoubles->Append(anInfo.Radius2());
383           theDoubles->Append(anInfo.Height());
384         }
385         break;
386       case GEOMAlgo_KN_POLYHEDRON:
387         // (+) geompy.kind.POLYHEDRON  nb_faces nb_edges nb_vertices
388         {
389           aKind = SK_POLYHEDRON;
390
391           theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
392           theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
393           theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
394         }
395         break;
396       default:
397         // (+) geompy.kind.SOLID  nb_faces nb_edges nb_vertices
398         {
399           theIntegers->Append(anInfo.NbSubShapes(TopAbs_FACE));
400           theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
401           theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
402         }
403       }
404     }
405     break;
406
407   case TopAbs_FACE:
408     {
409       aKind = SK_FACE;
410
411       GEOMAlgo_KindOfName aKN = anInfo.KindOfName();
412       switch (aKN) {
413       case GEOMAlgo_KN_SPHERE:
414         // (+) geompy.kind.SPHERE2D  xc yc zc  R
415         {
416           aKind = SK_SPHERE2D;
417
418           gp_Pnt aC = anInfo.Location();
419           theDoubles->Append(aC.X());
420           theDoubles->Append(aC.Y());
421           theDoubles->Append(aC.Z());
422
423           theDoubles->Append(anInfo.Radius1());
424         }
425         break;
426       case GEOMAlgo_KN_CYLINDER:
427         // (+) geompy.kind.CYLINDER2D  xb yb zb  dx dy dz  R  H
428         {
429           aKind = SK_CYLINDER2D;
430
431           gp_Pnt aO = anInfo.Location();
432           theDoubles->Append(aO.X());
433           theDoubles->Append(aO.Y());
434           theDoubles->Append(aO.Z());
435
436           gp_Ax3 anAx3 = anInfo.Position();
437           gp_Dir aD = anAx3.Direction();
438           theDoubles->Append(aD.X());
439           theDoubles->Append(aD.Y());
440           theDoubles->Append(aD.Z());
441
442           theDoubles->Append(anInfo.Radius1());
443           theDoubles->Append(anInfo.Height());
444         }
445         break;
446       case GEOMAlgo_KN_TORUS:
447         // (+) geompy.kind.TORUS2D  xc yc zc  dx dy dz  R_1 R_2
448         {
449           aKind = SK_TORUS2D;
450
451           gp_Pnt aO = anInfo.Location();
452           theDoubles->Append(aO.X());
453           theDoubles->Append(aO.Y());
454           theDoubles->Append(aO.Z());
455
456           gp_Ax3 anAx3 = anInfo.Position();
457           gp_Dir aD = anAx3.Direction();
458           theDoubles->Append(aD.X());
459           theDoubles->Append(aD.Y());
460           theDoubles->Append(aD.Z());
461
462           theDoubles->Append(anInfo.Radius1());
463           theDoubles->Append(anInfo.Radius2());
464         }
465         break;
466       case GEOMAlgo_KN_CONE:
467         // (+) geompy.kind.CONE2D  xc yc zc  dx dy dz  R_1 R_2  H
468         {
469           aKind = SK_CONE2D;
470
471           gp_Pnt aO = anInfo.Location();
472           theDoubles->Append(aO.X());
473           theDoubles->Append(aO.Y());
474           theDoubles->Append(aO.Z());
475
476           gp_Ax3 anAx3 = anInfo.Position();
477           gp_Dir aD = anAx3.Direction();
478           theDoubles->Append(aD.X());
479           theDoubles->Append(aD.Y());
480           theDoubles->Append(aD.Z());
481
482           theDoubles->Append(anInfo.Radius1());
483           theDoubles->Append(anInfo.Radius2());
484           theDoubles->Append(anInfo.Height());
485         }
486         break;
487       case GEOMAlgo_KN_DISKCIRCLE:
488         // (+) geompy.kind.DISK_CIRCLE  xc yc zc  dx dy dz  R
489         {
490           aKind = SK_DISK_CIRCLE;
491
492           gp_Pnt aC = anInfo.Location();
493           theDoubles->Append(aC.X());
494           theDoubles->Append(aC.Y());
495           theDoubles->Append(aC.Z());
496
497           gp_Ax3 anAx3 = anInfo.Position();
498           gp_Dir aD = anAx3.Direction();
499           theDoubles->Append(aD.X());
500           theDoubles->Append(aD.Y());
501           theDoubles->Append(aD.Z());
502
503           theDoubles->Append(anInfo.Radius1());
504         }
505         break;
506       case GEOMAlgo_KN_DISKELLIPSE:
507         // (+) geompy.kind.DISK_ELLIPSE  xc yc zc  dx dy dz  R_1 R_2
508         {
509           aKind = SK_DISK_ELLIPSE;
510
511           gp_Pnt aC = anInfo.Location();
512           theDoubles->Append(aC.X());
513           theDoubles->Append(aC.Y());
514           theDoubles->Append(aC.Z());
515
516           gp_Ax3 anAx3 = anInfo.Position();
517           gp_Dir aD = anAx3.Direction();
518           theDoubles->Append(aD.X());
519           theDoubles->Append(aD.Y());
520           theDoubles->Append(aD.Z());
521
522           theDoubles->Append(anInfo.Radius1());
523           theDoubles->Append(anInfo.Radius2());
524         }
525         break;
526       case GEOMAlgo_KN_RECTANGLE:
527       case GEOMAlgo_KN_TRIANGLE:
528       case GEOMAlgo_KN_QUADRANGLE:
529       case GEOMAlgo_KN_POLYGON:
530         // (+) geompy.kind.POLYGON  xo yo zo  dx dy dz  nb_edges nb_vertices
531         {
532           aKind = SK_POLYGON;
533
534           gp_Pnt aO = anInfo.Location();
535           theDoubles->Append(aO.X());
536           theDoubles->Append(aO.Y());
537           theDoubles->Append(aO.Z());
538
539           gp_Ax3 anAx3 = anInfo.Position();
540           gp_Dir aD = anAx3.Direction();
541           theDoubles->Append(aD.X());
542           theDoubles->Append(aD.Y());
543           theDoubles->Append(aD.Z());
544
545           theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
546           theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
547         }
548         break;
549       case GEOMAlgo_KN_PLANE: // infinite
550         // (+) geompy.kind.PLANE  xo yo zo  dx dy dz
551         {
552           aKind = SK_PLANE;
553
554           gp_Pnt aC = anInfo.Location();
555           theDoubles->Append(aC.X());
556           theDoubles->Append(aC.Y());
557           theDoubles->Append(aC.Z());
558
559           gp_Ax3 anAx3 = anInfo.Position();
560           gp_Dir aD = anAx3.Direction();
561           theDoubles->Append(aD.X());
562           theDoubles->Append(aD.Y());
563           theDoubles->Append(aD.Z());
564
565           if (anInfo.KindOfBounds() != GEOMAlgo_KB_INFINITE)
566           {
567             // (+) geompy.kind.PLANAR  xo yo zo  dx dy dz  nb_edges nb_vertices
568             aKind = SK_PLANAR;
569             
570             theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
571             theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
572           }
573         }
574         break;
575       default:
576         // ??? geompy.kind.FACE  nb_edges nb_vertices _surface_type_id_
577         // (+) geompy.kind.FACE  nb_edges nb_vertices
578         theIntegers->Append(anInfo.NbSubShapes(TopAbs_EDGE));
579         theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
580       }
581     }
582     break;
583
584   case TopAbs_EDGE:
585     {
586       aKind = SK_EDGE;
587
588       GEOMAlgo_KindOfName aKN = anInfo.KindOfName();
589       switch (aKN) {
590       case GEOMAlgo_KN_CIRCLE:
591         {
592           // (+) geompy.kind.CIRCLE  xc yc zc  dx dy dz  R
593           aKind = SK_CIRCLE;
594
595           gp_Pnt aC = anInfo.Location();
596           theDoubles->Append(aC.X());
597           theDoubles->Append(aC.Y());
598           theDoubles->Append(aC.Z());
599
600           gp_Ax3 anAx3 = anInfo.Position();
601           gp_Dir aD = anAx3.Direction();
602           theDoubles->Append(aD.X());
603           theDoubles->Append(aD.Y());
604           theDoubles->Append(aD.Z());
605
606           theDoubles->Append(anInfo.Radius1());
607         }
608         break;
609       case GEOMAlgo_KN_ARCCIRCLE:
610         {
611           // (+) geompy.kind.ARC_CIRCLE  xc yc zc  dx dy dz  R  x1 y1 z1  x2 y2 z2
612           aKind = SK_ARC_CIRCLE;
613
614           gp_Pnt aC = anInfo.Location();
615           theDoubles->Append(aC.X());
616           theDoubles->Append(aC.Y());
617           theDoubles->Append(aC.Z());
618
619           gp_Ax3 anAx3 = anInfo.Position();
620           gp_Dir aD = anAx3.Direction();
621           theDoubles->Append(aD.X());
622           theDoubles->Append(aD.Y());
623           theDoubles->Append(aD.Z());
624
625           theDoubles->Append(anInfo.Radius1());
626
627           gp_Pnt aP1 = anInfo.Pnt1();
628           theDoubles->Append(aP1.X());
629           theDoubles->Append(aP1.Y());
630           theDoubles->Append(aP1.Z());
631
632           gp_Pnt aP2 = anInfo.Pnt2();
633           theDoubles->Append(aP2.X());
634           theDoubles->Append(aP2.Y());
635           theDoubles->Append(aP2.Z());
636         }
637         break;
638       case GEOMAlgo_KN_ELLIPSE:
639         {
640           // (+) geompy.kind.ELLIPSE  xc yc zc  dx dy dz  R_1 R_2
641           aKind = SK_ELLIPSE;
642
643           gp_Pnt aC = anInfo.Location();
644           theDoubles->Append(aC.X());
645           theDoubles->Append(aC.Y());
646           theDoubles->Append(aC.Z());
647
648           gp_Ax3 anAx3 = anInfo.Position();
649           gp_Dir aD = anAx3.Direction();
650           theDoubles->Append(aD.X());
651           theDoubles->Append(aD.Y());
652           theDoubles->Append(aD.Z());
653
654           theDoubles->Append(anInfo.Radius1());
655           theDoubles->Append(anInfo.Radius2());
656         }
657         break;
658       case GEOMAlgo_KN_ARCELLIPSE:
659         {
660           // (+) geompy.kind.ARC_ELLIPSE  xc yc zc  dx dy dz  R_1 R_2  x1 y1 z1  x2 y2 z2
661           aKind = SK_ARC_ELLIPSE;
662
663           gp_Pnt aC = anInfo.Location();
664           theDoubles->Append(aC.X());
665           theDoubles->Append(aC.Y());
666           theDoubles->Append(aC.Z());
667
668           gp_Ax3 anAx3 = anInfo.Position();
669           gp_Dir aD = anAx3.Direction();
670           theDoubles->Append(aD.X());
671           theDoubles->Append(aD.Y());
672           theDoubles->Append(aD.Z());
673
674           theDoubles->Append(anInfo.Radius1());
675           theDoubles->Append(anInfo.Radius2());
676
677           gp_Pnt aP1 = anInfo.Pnt1();
678           theDoubles->Append(aP1.X());
679           theDoubles->Append(aP1.Y());
680           theDoubles->Append(aP1.Z());
681
682           gp_Pnt aP2 = anInfo.Pnt2();
683           theDoubles->Append(aP2.X());
684           theDoubles->Append(aP2.Y());
685           theDoubles->Append(aP2.Z());
686         }
687         break;
688       case GEOMAlgo_KN_LINE:
689         {
690           // ??? geompy.kind.LINE  x1 y1 z1  x2 y2 z2
691           // (+) geompy.kind.LINE  x1 y1 z1  dx dy dz
692           aKind = SK_LINE;
693
694           gp_Pnt aO = anInfo.Location();
695           theDoubles->Append(aO.X());
696           theDoubles->Append(aO.Y());
697           theDoubles->Append(aO.Z());
698
699           gp_Dir aD = anInfo.Direction();
700           theDoubles->Append(aD.X());
701           theDoubles->Append(aD.Y());
702           theDoubles->Append(aD.Z());
703         }
704         break;
705       case GEOMAlgo_KN_SEGMENT:
706         {
707           // (+) geompy.kind.SEGMENT  x1 y1 z1  x2 y2 z2
708           aKind = SK_SEGMENT;
709
710           gp_Pnt aP1 = anInfo.Pnt1();
711           theDoubles->Append(aP1.X());
712           theDoubles->Append(aP1.Y());
713           theDoubles->Append(aP1.Z());
714
715           gp_Pnt aP2 = anInfo.Pnt2();
716           theDoubles->Append(aP2.X());
717           theDoubles->Append(aP2.Y());
718           theDoubles->Append(aP2.Z());
719         }
720         break;
721       default:
722         // ??? geompy.kind.EDGE  nb_vertices _curve_type_id_
723         // (+) geompy.kind.EDGE  nb_vertices
724         theIntegers->Append(anInfo.NbSubShapes(TopAbs_VERTEX));
725       }
726     }
727     break;
728
729   case TopAbs_VERTEX:
730     {
731       // (+) geompy.kind.VERTEX  x y z
732       aKind = SK_VERTEX;
733
734       gp_Pnt aP = anInfo.Location();
735       theDoubles->Append(aP.X());
736       theDoubles->Append(aP.Y());
737       theDoubles->Append(aP.Z());
738     }
739     break;
740   }
741
742   SetErrorCode(OK);
743   return aKind;
744 }
745
746 //=============================================================================
747 /*!
748  *  GetPosition
749  */
750 //=============================================================================
751 void GEOMImpl_IMeasureOperations::GetPosition
752                    (Handle(GEOM_Object) theShape,
753                     Standard_Real& Ox, Standard_Real& Oy, Standard_Real& Oz,
754                     Standard_Real& Zx, Standard_Real& Zy, Standard_Real& Zz,
755                     Standard_Real& Xx, Standard_Real& Xy, Standard_Real& Xz)
756 {
757   SetErrorCode(KO);
758
759   //Set default values: global CS
760   Ox = Oy = Oz = Zx = Zy = Xy = Xz = 0.;
761   Zz = Xx = 1.;
762
763   if (theShape.IsNull()) return;
764
765   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
766   if (aRefShape.IsNull()) return;
767
768   TopoDS_Shape aShape = aRefShape->GetValue();
769   if (aShape.IsNull()) {
770     SetErrorCode("The Objects has NULL Shape");
771     return;
772   }
773
774   try {
775     OCC_CATCH_SIGNALS;
776
777     gp_Ax3 anAx3 = GEOMUtils::GetPosition(aShape);
778
779     gp_Pnt anOri = anAx3.Location();
780     gp_Dir aDirZ = anAx3.Direction();
781     gp_Dir aDirX = anAx3.XDirection();
782
783     // Output values
784     anOri.Coord(Ox, Oy, Oz);
785     aDirZ.Coord(Zx, Zy, Zz);
786     aDirX.Coord(Xx, Xy, Xz);
787   }
788   catch (Standard_Failure) {
789     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
790     SetErrorCode(aFail->GetMessageString());
791     return;
792   }
793
794   SetErrorCode(OK);
795 }
796
797 //=============================================================================
798 /*!
799  *  GetCentreOfMass
800  */
801 //=============================================================================
802 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetCentreOfMass
803                                                 (Handle(GEOM_Object) theShape)
804 {
805   SetErrorCode(KO);
806
807   if (theShape.IsNull()) return NULL;
808
809   //Add a new CentreOfMass object
810   Handle(GEOM_Object) aCDG = GetEngine()->AddObject(GetDocID(), GEOM_CDG);
811
812   //Add a new CentreOfMass function
813   Handle(GEOM_Function) aFunction =
814     aCDG->AddFunction(GEOMImpl_MeasureDriver::GetID(), CDG_MEASURE);
815   if (aFunction.IsNull()) return NULL;
816
817   //Check if the function is set correctly
818   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
819
820   GEOMImpl_IMeasure aCI (aFunction);
821
822   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
823   if (aRefShape.IsNull()) return NULL;
824
825   aCI.SetBase(aRefShape);
826
827   //Compute the CentreOfMass value
828   try {
829     OCC_CATCH_SIGNALS;
830     if (!GetSolver()->ComputeFunction(aFunction)) {
831       SetErrorCode("Measure driver failed to compute centre of mass");
832       return NULL;
833     }
834   }
835   catch (Standard_Failure) {
836     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
837     SetErrorCode(aFail->GetMessageString());
838     return NULL;
839   }
840
841   //Make a Python command
842   GEOM::TPythonDump(aFunction) << aCDG << " = geompy.MakeCDG(" << theShape << ")";
843
844   SetErrorCode(OK);
845   return aCDG;
846 }
847
848 //=============================================================================
849 /*!
850  *  GetVertexByIndex
851  */
852 //=============================================================================
853 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetVertexByIndex
854                                                 (Handle(GEOM_Object) theShape,
855                                                  Standard_Integer theIndex)
856 {
857   SetErrorCode(KO);
858
859   if (theShape.IsNull()) return NULL;
860
861   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
862   if (aRefShape.IsNull()) return NULL;
863
864   //Add a new Vertex object
865   Handle(GEOM_Object) aVertex = GetEngine()->AddObject(GetDocID(), GEOM_POINT);
866
867   //Add a function
868   Handle(GEOM_Function) aFunction =
869     aVertex->AddFunction(GEOMImpl_MeasureDriver::GetID(), VERTEX_BY_INDEX);
870   if (aFunction.IsNull()) return NULL;
871
872   //Check if the function is set correctly
873   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
874
875   GEOMImpl_IMeasure aCI (aFunction);
876   aCI.SetBase(aRefShape);
877   aCI.SetIndex(theIndex);
878
879   //Compute
880   try {
881     OCC_CATCH_SIGNALS;
882     if (!GetSolver()->ComputeFunction(aFunction)) {
883       SetErrorCode("Vertex by index driver failed.");
884       return NULL;
885     }
886   }
887   catch (Standard_Failure) {
888     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
889     SetErrorCode(aFail->GetMessageString());
890     return NULL;
891   }
892
893   //Make a Python command
894   GEOM::TPythonDump(aFunction) << aVertex << " = geompy.GetVertexByIndex(" << theShape << ", " << theIndex << ")";
895
896   SetErrorCode(OK);
897   return aVertex;
898 }
899
900 //=============================================================================
901 /*!
902  *  GetNormal
903  */
904 //=============================================================================
905 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetNormal
906                                          (Handle(GEOM_Object) theFace,
907                                           Handle(GEOM_Object) theOptionalPoint)
908 {
909   SetErrorCode(KO);
910
911   if (theFace.IsNull()) return NULL;
912
913   //Add a new Normale object
914   Handle(GEOM_Object) aNorm = GetEngine()->AddObject(GetDocID(), GEOM_VECTOR);
915
916   //Add a new Normale function
917   Handle(GEOM_Function) aFunction =
918     aNorm->AddFunction(GEOMImpl_MeasureDriver::GetID(), VECTOR_FACE_NORMALE);
919   if (aFunction.IsNull()) return NULL;
920
921   //Check if the function is set correctly
922   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
923
924   GEOMImpl_IMeasure aCI (aFunction);
925
926   Handle(GEOM_Function) aFace = theFace->GetLastFunction();
927   if (aFace.IsNull()) return NULL;
928
929   aCI.SetBase(aFace);
930
931   if (!theOptionalPoint.IsNull()) {
932     Handle(GEOM_Function) anOptPnt = theOptionalPoint->GetLastFunction();
933     aCI.SetPoint(anOptPnt);
934   }
935
936   //Compute the Normale value
937   try {
938     OCC_CATCH_SIGNALS;
939     if (!GetSolver()->ComputeFunction(aFunction)) {
940       SetErrorCode("Measure driver failed to compute normake of face");
941       return NULL;
942     }
943   }
944   catch (Standard_Failure) {
945     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
946     SetErrorCode(aFail->GetMessageString());
947     return NULL;
948   }
949
950   //Make a Python command
951   GEOM::TPythonDump pd (aFunction);
952   pd << aNorm << " = geompy.GetNormal(" << theFace;
953   if (!theOptionalPoint.IsNull()) {
954     pd << ", " << theOptionalPoint;
955   }
956   pd << ")";
957
958   SetErrorCode(OK);
959   return aNorm;
960 }
961
962 //=============================================================================
963 /*!
964  *  GetBasicProperties
965  */
966 //=============================================================================
967 void GEOMImpl_IMeasureOperations::GetBasicProperties (Handle(GEOM_Object) theShape,
968                                                       Standard_Real& theLength,
969                                                       Standard_Real& theSurfArea,
970                                                       Standard_Real& theVolume)
971 {
972   SetErrorCode(KO);
973
974   if (theShape.IsNull()) return;
975
976   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
977   if (aRefShape.IsNull()) return;
978
979   TopoDS_Shape aShape = aRefShape->GetValue();
980   if (aShape.IsNull()) {
981     SetErrorCode("The Objects has NULL Shape");
982     return;
983   }
984
985   //Compute the parameters
986   GProp_GProps LProps, SProps;
987   Standard_Real anEps = 1.e-6;
988   try {
989     OCC_CATCH_SIGNALS;
990     BRepGProp::LinearProperties(aShape, LProps);
991     theLength = LProps.Mass();
992
993     BRepGProp::SurfaceProperties(aShape, SProps, anEps);
994     theSurfArea = SProps.Mass();
995
996     theVolume = 0.0;
997     if (aShape.ShapeType() < TopAbs_SHELL) {
998       for (TopExp_Explorer Exp (aShape, TopAbs_SOLID); Exp.More(); Exp.Next()) {
999         GProp_GProps VProps;
1000         BRepGProp::VolumeProperties(Exp.Current(), VProps, anEps);
1001         theVolume += VProps.Mass();
1002       }
1003     }
1004   }
1005   catch (Standard_Failure) {
1006     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1007     SetErrorCode(aFail->GetMessageString());
1008     return;
1009   }
1010
1011   SetErrorCode(OK);
1012 }
1013
1014 //=============================================================================
1015 /*!
1016  *  GetInertia
1017  */
1018 //=============================================================================
1019 void GEOMImpl_IMeasureOperations::GetInertia
1020                    (Handle(GEOM_Object) theShape,
1021                     Standard_Real& I11, Standard_Real& I12, Standard_Real& I13,
1022                     Standard_Real& I21, Standard_Real& I22, Standard_Real& I23,
1023                     Standard_Real& I31, Standard_Real& I32, Standard_Real& I33,
1024                     Standard_Real& Ix , Standard_Real& Iy , Standard_Real& Iz)
1025 {
1026   SetErrorCode(KO);
1027
1028   if (theShape.IsNull()) return;
1029
1030   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1031   if (aRefShape.IsNull()) return;
1032
1033   TopoDS_Shape aShape = aRefShape->GetValue();
1034   if (aShape.IsNull()) {
1035     SetErrorCode("The Objects has NULL Shape");
1036     return;
1037   }
1038
1039   //Compute the parameters
1040   GProp_GProps System;
1041
1042   try {
1043     OCC_CATCH_SIGNALS;
1044     if (aShape.ShapeType() == TopAbs_VERTEX ||
1045         aShape.ShapeType() == TopAbs_EDGE ||
1046         aShape.ShapeType() == TopAbs_WIRE) {
1047       BRepGProp::LinearProperties(aShape, System);
1048     } else if (aShape.ShapeType() == TopAbs_FACE ||
1049                aShape.ShapeType() == TopAbs_SHELL) {
1050       BRepGProp::SurfaceProperties(aShape, System);
1051     } else {
1052       BRepGProp::VolumeProperties(aShape, System);
1053     }
1054     gp_Mat I = System.MatrixOfInertia();
1055
1056     I11 = I(1,1);
1057     I12 = I(1,2);
1058     I13 = I(1,3);
1059
1060     I21 = I(2,1);
1061     I22 = I(2,2);
1062     I23 = I(2,3);
1063
1064     I31 = I(3,1);
1065     I32 = I(3,2);
1066     I33 = I(3,3);
1067
1068     GProp_PrincipalProps Pr = System.PrincipalProperties();
1069     Pr.Moments(Ix,Iy,Iz);
1070   }
1071   catch (Standard_Failure) {
1072     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1073     SetErrorCode(aFail->GetMessageString());
1074     return;
1075   }
1076
1077   SetErrorCode(OK);
1078 }
1079
1080 //=============================================================================
1081 /*!
1082  *  GetBoundingBox
1083  */
1084 //=============================================================================
1085 void GEOMImpl_IMeasureOperations::GetBoundingBox
1086                                      (Handle(GEOM_Object) theShape,
1087                                       const Standard_Boolean precise,
1088                                       Standard_Real& Xmin, Standard_Real& Xmax,
1089                                       Standard_Real& Ymin, Standard_Real& Ymax,
1090                                       Standard_Real& Zmin, Standard_Real& Zmax)
1091 {
1092   SetErrorCode(KO);
1093
1094   if (theShape.IsNull()) return;
1095
1096   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1097   if (aRefShape.IsNull()) return;
1098
1099   TopoDS_Shape aShape = aRefShape->GetValue();
1100   if (aShape.IsNull()) {
1101     SetErrorCode("The Objects has NULL Shape");
1102     return;
1103   }
1104
1105   //Compute the parameters
1106   Bnd_Box B;
1107
1108   try {
1109     OCC_CATCH_SIGNALS;
1110     BRepBuilderAPI_Copy aCopyTool (aShape);
1111     if (!aCopyTool.IsDone()) {
1112       SetErrorCode("GetBoundingBox Error: Bad shape detected");
1113       return;
1114     }
1115
1116     aShape = aCopyTool.Shape();
1117
1118     // remove triangulation to obtain more exact boundaries
1119     BRepTools::Clean(aShape);
1120
1121     BRepBndLib::Add(aShape, B);
1122
1123     if (precise) {
1124       if (!GEOMUtils::PreciseBoundingBox(aShape, B)) {
1125         SetErrorCode("GetBoundingBox Error: Bounding box cannot be precised");
1126         return;
1127       }
1128     }
1129
1130     B.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
1131   }
1132   catch (Standard_Failure) {
1133     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1134     SetErrorCode(aFail->GetMessageString());
1135     return;
1136   }
1137
1138   SetErrorCode(OK);
1139 }
1140
1141 //=============================================================================
1142 /*!
1143  *  GetBoundingBox
1144  */
1145 //=============================================================================
1146 Handle(GEOM_Object) GEOMImpl_IMeasureOperations::GetBoundingBox
1147                                                 (Handle(GEOM_Object) theShape,
1148                                                  const Standard_Boolean precise)
1149 {
1150   SetErrorCode(KO);
1151
1152   if (theShape.IsNull()) return NULL;
1153
1154   //Add a new BoundingBox object
1155   Handle(GEOM_Object) aBnd = GetEngine()->AddObject(GetDocID(), GEOM_BOX);
1156
1157   //Add a new BoundingBox function
1158   const int aType = (precise ? BND_BOX_MEASURE_PRECISE : BND_BOX_MEASURE);
1159   Handle(GEOM_Function) aFunction =
1160     aBnd->AddFunction(GEOMImpl_MeasureDriver::GetID(), aType);
1161   if (aFunction.IsNull()) return NULL;
1162
1163   //Check if the function is set correctly
1164   if (aFunction->GetDriverGUID() != GEOMImpl_MeasureDriver::GetID()) return NULL;
1165
1166   GEOMImpl_IMeasure aCI (aFunction);
1167
1168   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1169   if (aRefShape.IsNull()) return NULL;
1170
1171   aCI.SetBase(aRefShape);
1172
1173   //Compute the BoundingBox value
1174   try {
1175     OCC_CATCH_SIGNALS;
1176     if (!GetSolver()->ComputeFunction(aFunction)) {
1177       SetErrorCode("Measure driver failed to compute a bounding box");
1178       return NULL;
1179     }
1180   }
1181   catch (Standard_Failure) {
1182     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1183     SetErrorCode(aFail->GetMessageString());
1184     return NULL;
1185   }
1186
1187   //Make a Python command
1188   GEOM::TPythonDump aPd(aFunction);
1189   
1190   aPd << aBnd << " = geompy.MakeBoundingBox(" << theShape;
1191
1192   if (precise) {
1193     aPd << ", True";
1194   }
1195
1196   aPd << ")";
1197
1198   SetErrorCode(OK);
1199   return aBnd;
1200 }
1201
1202 //=============================================================================
1203 /*!
1204  *  GetTolerance
1205  */
1206 //=============================================================================
1207 void GEOMImpl_IMeasureOperations::GetTolerance
1208                                (Handle(GEOM_Object) theShape,
1209                                 Standard_Real& FaceMin, Standard_Real& FaceMax,
1210                                 Standard_Real& EdgeMin, Standard_Real& EdgeMax,
1211                                 Standard_Real& VertMin, Standard_Real& VertMax)
1212 {
1213   SetErrorCode(KO);
1214
1215   if (theShape.IsNull()) return;
1216
1217   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1218   if (aRefShape.IsNull()) return;
1219
1220   TopoDS_Shape aShape = aRefShape->GetValue();
1221   if (aShape.IsNull()) {
1222     SetErrorCode("The Objects has NULL Shape");
1223     return;
1224   }
1225
1226   //Compute the parameters
1227   Standard_Real T;
1228   FaceMin = EdgeMin = VertMin = RealLast();
1229   FaceMax = EdgeMax = VertMax = -RealLast();
1230
1231   try {
1232     OCC_CATCH_SIGNALS;
1233     for (TopExp_Explorer ExF (aShape, TopAbs_FACE); ExF.More(); ExF.Next()) {
1234       TopoDS_Face Face = TopoDS::Face(ExF.Current());
1235       T = BRep_Tool::Tolerance(Face);
1236       if (T > FaceMax)
1237         FaceMax = T;
1238       if (T < FaceMin)
1239         FaceMin = T;
1240     }
1241     for (TopExp_Explorer ExE (aShape, TopAbs_EDGE); ExE.More(); ExE.Next()) {
1242       TopoDS_Edge Edge = TopoDS::Edge(ExE.Current());
1243       T = BRep_Tool::Tolerance(Edge);
1244       if (T > EdgeMax)
1245         EdgeMax = T;
1246       if (T < EdgeMin)
1247         EdgeMin = T;
1248     }
1249     for (TopExp_Explorer ExV (aShape, TopAbs_VERTEX); ExV.More(); ExV.Next()) {
1250       TopoDS_Vertex Vertex = TopoDS::Vertex(ExV.Current());
1251       T = BRep_Tool::Tolerance(Vertex);
1252       if (T > VertMax)
1253         VertMax = T;
1254       if (T < VertMin)
1255         VertMin = T;
1256     }
1257   }
1258   catch (Standard_Failure) {
1259     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1260     SetErrorCode(aFail->GetMessageString());
1261     return;
1262   }
1263
1264   SetErrorCode(OK);
1265 }
1266
1267 //=============================================================================
1268 /*!
1269  *  CheckShape
1270  */
1271 //=============================================================================
1272 bool GEOMImpl_IMeasureOperations::CheckShape (Handle(GEOM_Object)     theShape,
1273                                               const Standard_Boolean  theIsCheckGeom,
1274                                               std::list<ShapeError>  &theErrors)
1275 {
1276   SetErrorCode(KO);
1277   theErrors.clear();
1278
1279   if (theShape.IsNull()) return false;
1280
1281   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1282   if (aRefShape.IsNull()) return false;
1283
1284   TopoDS_Shape aShape = aRefShape->GetValue();
1285   if (aShape.IsNull()) {
1286     SetErrorCode("The Objects has NULL Shape");
1287     return false;
1288   }
1289
1290   //Compute the parameters
1291   bool isValid = false;
1292   try {
1293     OCC_CATCH_SIGNALS;
1294     BRepCheck_Analyzer ana (aShape, theIsCheckGeom);
1295     if (ana.IsValid()) {
1296       isValid = true;
1297     } else {
1298       FillErrors(ana, aShape, theErrors);
1299     }
1300   }
1301   catch (Standard_Failure) {
1302     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1303     SetErrorCode(aFail->GetMessageString());
1304     return false;
1305   }
1306
1307   SetErrorCode(OK);
1308   return isValid;
1309 }
1310
1311 //=============================================================================
1312 /*!
1313  *  PrintShapeErrors
1314  */
1315 //=============================================================================
1316 TCollection_AsciiString GEOMImpl_IMeasureOperations::PrintShapeErrors
1317                                  (Handle(GEOM_Object)          theShape,
1318                                   const std::list<ShapeError> &theErrors)
1319 {
1320   TCollection_AsciiString aDump;
1321
1322   if (theShape.IsNull()) {
1323     return aDump;
1324   }
1325
1326   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1327
1328   if (aRefShape.IsNull()) {
1329     return aDump;
1330   }
1331
1332   TopoDS_Shape aShape = aRefShape->GetValue();
1333
1334   if (aShape.IsNull()) {
1335     SetErrorCode("The Objects has NULL Shape");
1336     return aDump;
1337   }
1338
1339   if (!theErrors.empty()) {
1340     // The shape is not valid. Prepare errors for dump.
1341     TopTools_IndexedMapOfShape anIndices;
1342     std::list<ShapeError>::const_iterator anIter = theErrors.begin();
1343     Standard_Integer nbv, nbe, nbw, nbf, nbs, nbo;
1344     nbv = nbe = nbw = nbf = nbs = nbo = 0;
1345
1346     // Map sub-shapes and their indices
1347     TopExp::MapShapes(aShape, anIndices);
1348
1349     const Standard_Integer aNbSubShapes = anIndices.Extent();
1350     TColStd_MapOfInteger   aMapPbInd;
1351
1352     aDump += " -- The Shape has problems :\n";
1353     aDump += "  Check                                    Count\n";
1354     aDump += " ------------------------------------------------\n";
1355
1356     for (; anIter != theErrors.end(); anIter++) {
1357       Standard_Integer aNbShapes = anIter->incriminated.size();
1358
1359       switch(anIter->error) {
1360       case BRepCheck_InvalidPointOnCurve:
1361         aDump += "  Invalid Point on Curve ................... ";
1362         break;
1363       case BRepCheck_InvalidPointOnCurveOnSurface:
1364         aDump += "  Invalid Point on CurveOnSurface .......... ";
1365         break;
1366       case BRepCheck_InvalidPointOnSurface:
1367         aDump += "  Invalid Point on Surface ................. ";
1368         break;
1369       case BRepCheck_No3DCurve:
1370         aDump += "  No 3D Curve .............................. ";
1371         break;
1372       case BRepCheck_Multiple3DCurve:
1373         aDump += "  Multiple 3D Curve ........................ ";
1374         break;
1375       case BRepCheck_Invalid3DCurve:
1376         aDump += "  Invalid 3D Curve ......................... ";
1377         break;
1378       case BRepCheck_NoCurveOnSurface:
1379         aDump += "  No Curve on Surface ...................... ";
1380         break;
1381       case BRepCheck_InvalidCurveOnSurface:
1382         aDump += "  Invalid Curve on Surface ................. ";
1383         break;
1384       case BRepCheck_InvalidCurveOnClosedSurface:
1385         aDump += "  Invalid Curve on closed Surface .......... ";
1386         break;
1387       case BRepCheck_InvalidSameRangeFlag:
1388         aDump += "  Invalid SameRange Flag ................... ";
1389         break;
1390       case BRepCheck_InvalidSameParameterFlag:
1391         aDump += "  Invalid SameParameter Flag ............... ";
1392         break;
1393       case BRepCheck_InvalidDegeneratedFlag:
1394         aDump += "  Invalid Degenerated Flag ................. ";
1395         break;
1396       case BRepCheck_FreeEdge:
1397         aDump += "  Free Edge ................................ ";
1398         break;
1399       case BRepCheck_InvalidMultiConnexity:
1400         aDump += "  Invalid MultiConnexity ................... ";
1401         break;
1402       case BRepCheck_InvalidRange:
1403         aDump += "  Invalid Range ............................ ";
1404         break;
1405       case BRepCheck_EmptyWire:
1406         aDump += "  Empty Wire ............................... ";
1407         break;
1408       case BRepCheck_RedundantEdge:
1409         aDump += "  Redundant Edge ........................... ";
1410         break;
1411       case BRepCheck_SelfIntersectingWire:
1412         aDump += "  Self Intersecting Wire ................... ";
1413         break;
1414       case BRepCheck_NoSurface:
1415         aDump += "  No Surface ............................... ";
1416         break;
1417       case BRepCheck_InvalidWire:
1418         aDump += "  Invalid Wire ............................. ";
1419         break;
1420       case BRepCheck_RedundantWire:
1421         aDump += "  Redundant Wire ........................... ";
1422         break;
1423       case BRepCheck_IntersectingWires:
1424         aDump += "  Intersecting Wires ....................... ";
1425         break;
1426       case BRepCheck_InvalidImbricationOfWires:
1427         aDump += "  Invalid Imbrication of Wires ............. ";
1428         break;
1429       case BRepCheck_EmptyShell:
1430         aDump += "  Empty Shell .............................. ";
1431         break;
1432       case BRepCheck_RedundantFace:
1433         aDump += "  Redundant Face ........................... ";
1434         break;
1435       case BRepCheck_UnorientableShape:
1436         aDump += "  Unorientable Shape ....................... ";
1437         break;
1438       case BRepCheck_NotClosed:
1439         aDump += "  Not Closed ............................... ";
1440         break;
1441       case BRepCheck_NotConnected:
1442         aDump += "  Not Connected ............................ ";
1443         break;
1444       case BRepCheck_SubshapeNotInShape:
1445         aDump += "  Sub-shape not in Shape ................... ";
1446         break;
1447       case BRepCheck_BadOrientation:
1448         aDump += "  Bad Orientation .......................... ";
1449         break;
1450       case BRepCheck_BadOrientationOfSubshape:
1451         aDump += "  Bad Orientation of Sub-shape ............. ";
1452         break;
1453       case BRepCheck_InvalidToleranceValue:
1454         aDump += "  Invalid Tolerance Value .................. ";
1455         break;
1456       case BRepCheck_CheckFail:
1457         aDump += "  Check Shape Failure ...................... ";
1458         break;
1459       default:
1460         break;
1461       }
1462
1463       aDump += TCollection_AsciiString(aNbShapes) + "\n";
1464
1465       // Count types of shape.
1466       std::list<int>::const_iterator aShIter = anIter->incriminated.begin();
1467
1468       for (; aShIter != anIter->incriminated.end(); aShIter++) {
1469         const Standard_Integer anIndex = *aShIter;
1470
1471         if (anIndex > 0 && anIndex <= aNbSubShapes && aMapPbInd.Add(anIndex)) {
1472           const TopoDS_Shape     &aSubShape = anIndices.FindKey(anIndex);
1473           const TopAbs_ShapeEnum  aType     = aSubShape.ShapeType();
1474
1475           switch (aType) {
1476             case TopAbs_VERTEX : nbv++; break;
1477             case TopAbs_EDGE   : nbe++; break;
1478             case TopAbs_WIRE   : nbw++; break;
1479             case TopAbs_FACE   : nbf++; break;
1480             case TopAbs_SHELL  : nbs++; break;
1481             case TopAbs_SOLID  : nbo++; break;
1482             default            : break;
1483           }
1484         }
1485       }
1486     }
1487
1488     const Standard_Integer aNbFaultyShapes = nbv + nbe + nbw + nbf + nbs + nbo;
1489     aDump += " ------------------------------------------------\n";
1490     aDump += "*** Shapes with problems : ";
1491     aDump += TCollection_AsciiString(aNbFaultyShapes) + "\n";
1492
1493     if (nbv > 0) {
1494       aDump += "VERTEX : ";
1495       if (nbv < 10) aDump += " ";
1496       aDump += TCollection_AsciiString(nbv) + "\n";
1497     }
1498     if (nbe > 0) {
1499       aDump += "EDGE   : ";
1500       if (nbe < 10) aDump += " ";
1501       aDump += TCollection_AsciiString(nbe) + "\n";
1502     }
1503     if (nbw > 0) {
1504       aDump += "WIRE   : ";
1505       if (nbw < 10) aDump += " ";
1506       aDump += TCollection_AsciiString(nbw) + "\n";
1507     }
1508     if (nbf > 0) {
1509       aDump += "FACE   : ";
1510       if (nbf < 10) aDump += " ";
1511       aDump += TCollection_AsciiString(nbf) + "\n";
1512     }
1513     if (nbs > 0) {
1514       aDump += "SHELL  : ";
1515       if (nbs < 10) aDump += " ";
1516       aDump += TCollection_AsciiString(nbs) + "\n";
1517     }
1518     if (nbo > 0) {
1519       aDump += "SOLID  : ";
1520       if (nbo < 10) aDump += " ";
1521       aDump += TCollection_AsciiString(nbo) + "\n";
1522     }
1523   }
1524
1525   return aDump;
1526 }
1527
1528 //=============================================================================
1529 /*!
1530  *  CheckSelfIntersections
1531  */
1532 //=============================================================================
1533 bool GEOMImpl_IMeasureOperations::CheckSelfIntersections
1534                          (Handle(GEOM_Object)                 theShape,
1535                           const SICheckLevel                  theCheckLevel,
1536                           Handle(TColStd_HSequenceOfInteger)& theIntersections)
1537 {
1538   SetErrorCode(KO);
1539
1540   if (theIntersections.IsNull())
1541     theIntersections = new TColStd_HSequenceOfInteger;
1542   else
1543     theIntersections->Clear();
1544
1545   if (theShape.IsNull())
1546     return false;
1547
1548   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1549   if (aRefShape.IsNull()) return false;
1550
1551   TopoDS_Shape aShape = aRefShape->GetValue();
1552   if (aShape.IsNull()) return false;
1553
1554   // 0. Prepare data
1555   TopoDS_Shape aScopy;
1556   //
1557   GEOMAlgo_AlgoTools::CopyShape(aShape, aScopy);
1558   //
1559   // Map sub-shapes and their indices
1560   TopTools_IndexedMapOfShape anIndices;
1561   TopExp::MapShapes(aScopy, anIndices);
1562
1563   BOPCol_ListOfShape aLCS;
1564   aLCS.Append(aScopy);
1565   //
1566   BOPAlgo_CheckerSI aCSI; // checker of self-interferences
1567   aCSI.SetArguments(aLCS);
1568   aCSI.SetLevelOfCheck(theCheckLevel);
1569
1570   // 1. Launch the checker
1571   aCSI.Perform();
1572   Standard_Integer iErr = aCSI.ErrorStatus();
1573
1574   //
1575   Standard_Integer aNbS, n1, n2;
1576   BOPDS_MapIteratorMapOfPassKey aItMPK;
1577   //
1578   // 2. Take the shapes from DS
1579   const BOPDS_DS& aDS = aCSI.DS();
1580   aNbS=aDS.NbShapes();
1581   //
1582   // 3. Get the pairs of interfered shapes
1583   const BOPDS_MapOfPassKey& aMPK=aDS.Interferences();
1584   aItMPK.Initialize(aMPK);
1585   for (; aItMPK.More(); aItMPK.Next()) {
1586     const BOPDS_PassKey& aPK=aItMPK.Value();
1587     aPK.Ids(n1, n2);
1588     //
1589     if (n1 > aNbS || n2 > aNbS){
1590       return false; // Error
1591     }
1592     const TopoDS_Shape& aS1 = aDS.Shape(n1);
1593     const TopoDS_Shape& aS2 = aDS.Shape(n2);
1594
1595     theIntersections->Append(anIndices.FindIndex(aS1));
1596     theIntersections->Append(anIndices.FindIndex(aS2));
1597   }
1598
1599   if (!iErr) {
1600     SetErrorCode(OK);
1601   }
1602
1603   return theIntersections->IsEmpty();
1604 }
1605
1606 //=============================================================================
1607 /*!
1608  *  CheckSelfIntersectionsFast
1609  */
1610 //=============================================================================
1611 bool GEOMImpl_IMeasureOperations::CheckSelfIntersectionsFast
1612                          (Handle(GEOM_Object) theShape,
1613                           float theDeflection, double theTolerance,
1614                           Handle(TColStd_HSequenceOfInteger)& theIntersections)
1615 {
1616   SetErrorCode(KO);
1617
1618   if (theIntersections.IsNull())
1619     theIntersections = new TColStd_HSequenceOfInteger;
1620   else
1621     theIntersections->Clear();
1622
1623   if (theShape.IsNull())
1624     return false;
1625
1626   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1627   if (aRefShape.IsNull()) return false;
1628
1629   TopoDS_Shape aShape = aRefShape->GetValue();
1630   if (aShape.IsNull()) return false;
1631
1632   // Prepare data
1633   TopoDS_Shape aScopy;
1634
1635   GEOMAlgo_AlgoTools::CopyShape(aShape, aScopy);
1636   GEOMUtils::MeshShape(aScopy, theDeflection);
1637
1638   // Map sub-shapes and their indices
1639   TopTools_IndexedMapOfShape anIndices;
1640   TopExp::MapShapes(aScopy, anIndices);
1641
1642 #if OCC_VERSION_LARGE > 0x06090000
1643   // Checker of fast interferences
1644   BRepExtrema_SelfIntersection aTool(aScopy, (theTolerance <= 0.) ? 0.0 : theTolerance);
1645
1646   // Launch the checker
1647   aTool.Perform();
1648   
1649   const BRepExtrema_MapOfIntegerPackedMapOfInteger& intersections = aTool.OverlapElements();
1650   
1651   std::set<Standard_Integer> processed;
1652   
1653   for (BRepExtrema_MapOfIntegerPackedMapOfInteger::Iterator it(intersections); it.More(); it.Next()) {
1654     Standard_Integer idxLeft = it.Key();
1655     if (processed.count(idxLeft) > 0) continue; // already added
1656     processed.insert(idxLeft);
1657     const TColStd_PackedMapOfInteger& overlaps = it.Value();
1658     for (TColStd_MapIteratorOfPackedMapOfInteger subit(overlaps); subit.More(); subit.Next()) {
1659       Standard_Integer idxRight = subit.Key();
1660       if (processed.count(idxRight) > 0) continue; // already added
1661       const TopoDS_Shape& aS1 = aTool.GetSubShape(idxLeft);
1662       const TopoDS_Shape& aS2 = aTool.GetSubShape(idxRight);
1663       theIntersections->Append(anIndices.FindIndex(aS1));
1664       theIntersections->Append(anIndices.FindIndex(aS2));
1665     }
1666   }
1667
1668   if (aTool.IsDone())
1669     SetErrorCode(OK);
1670 #endif
1671
1672   return theIntersections->IsEmpty();
1673 }
1674
1675 //=============================================================================
1676 /*!
1677  *  FastIntersect
1678  */
1679 //=============================================================================
1680 bool GEOMImpl_IMeasureOperations::FastIntersect (Handle(GEOM_Object) theShape1, Handle(GEOM_Object) theShape2,
1681                                                  double theTolerance, float theDeflection,
1682                                                  Handle(TColStd_HSequenceOfInteger)& theIntersections1,
1683                                                  Handle(TColStd_HSequenceOfInteger)& theIntersections2)
1684 {
1685   SetErrorCode(KO);
1686   bool isGood = false;
1687
1688   if (theIntersections1.IsNull())
1689     theIntersections1 = new TColStd_HSequenceOfInteger;
1690   else
1691     theIntersections1->Clear();
1692
1693   if (theIntersections2.IsNull())
1694     theIntersections2 = new TColStd_HSequenceOfInteger;
1695   else
1696     theIntersections2->Clear();
1697
1698   if (theShape1.IsNull() || theShape2.IsNull()) {
1699     SetErrorCode("Objects have NULL Shape");
1700     return isGood;
1701   }
1702
1703   if (theShape1 == theShape2) {
1704     SetErrorCode("Objects are equal");
1705     return isGood;
1706   }
1707   Handle(GEOM_Function) aRefShape1 = theShape1->GetLastFunction();
1708   Handle(GEOM_Function) aRefShape2 = theShape2->GetLastFunction();
1709   if (aRefShape1.IsNull() || aRefShape2.IsNull()) return isGood;
1710
1711   TopoDS_Shape aShape1 = aRefShape1->GetValue();
1712   TopoDS_Shape aShape2 = aRefShape2->GetValue();
1713   if (aShape1.IsNull() || aShape2.IsNull()) return isGood;
1714
1715   // 0. Prepare data
1716   TopoDS_Shape aScopy1, aScopy2;
1717   GEOMAlgo_AlgoTools::CopyShape(aShape1, aScopy1);
1718   GEOMAlgo_AlgoTools::CopyShape(aShape2, aScopy2);
1719
1720   GEOMUtils::MeshShape(aScopy1, theDeflection);
1721   GEOMUtils::MeshShape(aScopy2, theDeflection);
1722   //
1723   // Map sub-shapes and their indices
1724   TopTools_IndexedMapOfShape anIndices1, anIndices2;
1725   TopExp::MapShapes(aScopy1, anIndices1);
1726   TopExp::MapShapes(aScopy2, anIndices2);
1727
1728   BOPCol_ListOfShape aLCS1, aLCS2;
1729   aLCS1.Append(aScopy1); aLCS2.Append(aScopy2);
1730   //
1731   BRepExtrema_ShapeProximity aBSP; // checker of fast interferences
1732   aBSP.LoadShape1(aScopy1); aBSP.LoadShape2(aScopy2);
1733   aBSP.SetTolerance((theTolerance <= 0.) ? 0.0 : theTolerance);
1734
1735   // 1. Launch the checker
1736   aBSP.Perform();
1737  
1738   // 2. Get sets of IDs of overlapped faces
1739 #if OCC_VERSION_LARGE > 0x06090000
1740   for (BRepExtrema_MapOfIntegerPackedMapOfInteger::Iterator anIt1 (aBSP.OverlapSubShapes1()); anIt1.More(); anIt1.Next())
1741 #else
1742   for (BRepExtrema_OverlappedSubShapes::Iterator anIt1 (aBSP.OverlapSubShapes1()); anIt1.More(); anIt1.Next())
1743 #endif
1744   {
1745     const TopoDS_Shape& aS1 = aBSP.GetSubShape1(anIt1.Key());
1746     theIntersections1->Append(anIndices1.FindIndex(aS1));
1747   }
1748   
1749 #if OCC_VERSION_LARGE > 0x06090000
1750   for (BRepExtrema_MapOfIntegerPackedMapOfInteger::Iterator anIt2 (aBSP.OverlapSubShapes2()); anIt2.More(); anIt2.Next())
1751 #else
1752   for (BRepExtrema_OverlappedSubShapes::Iterator anIt2 (aBSP.OverlapSubShapes2()); anIt2.More(); anIt2.Next())
1753 #endif
1754   {
1755     const TopoDS_Shape& aS2 = aBSP.GetSubShape2(anIt2.Key());
1756     theIntersections2->Append(anIndices2.FindIndex(aS2));
1757   }
1758
1759   isGood = !theIntersections1->IsEmpty() && !theIntersections1->IsEmpty();
1760
1761   if (aBSP.IsDone())
1762     SetErrorCode(OK);
1763
1764   return isGood;
1765 }
1766
1767 //=============================================================================
1768 /*!
1769  *  IsGoodForSolid
1770  */
1771 //=============================================================================
1772 TCollection_AsciiString GEOMImpl_IMeasureOperations::IsGoodForSolid (Handle(GEOM_Object) theShape)
1773 {
1774   SetErrorCode(KO);
1775
1776   TCollection_AsciiString aRes = "";
1777
1778   if (theShape.IsNull()) {
1779     aRes = "WRN_NULL_OBJECT_OR_SHAPE";
1780   }
1781   else {
1782     Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1783     if (aRefShape.IsNull()) {
1784       aRes = "WRN_NULL_OBJECT_OR_SHAPE";
1785     }
1786     else {
1787       TopoDS_Shape aShape = aRefShape->GetValue();
1788       if (aShape.IsNull()) {
1789         aRes = "WRN_NULL_OBJECT_OR_SHAPE";
1790       }
1791       else {
1792         if (aShape.ShapeType() == TopAbs_COMPOUND) {
1793           TopoDS_Iterator It (aShape, Standard_True, Standard_True);
1794           if (It.More()) aShape = It.Value();
1795         }
1796         if (aShape.ShapeType() == TopAbs_SHELL) {
1797           BRepCheck_Shell chkShell (TopoDS::Shell(aShape));
1798           if (chkShell.Closed() == BRepCheck_NotClosed) {
1799             aRes = "WRN_SHAPE_UNCLOSED";
1800           }
1801         }
1802         else {
1803           aRes = "WRN_SHAPE_NOT_SHELL";
1804         }
1805       }
1806     }
1807   }
1808
1809   if (aRes.IsEmpty())
1810     SetErrorCode(OK);
1811
1812   return aRes;
1813 }
1814
1815 //=============================================================================
1816 /*!
1817  *  WhatIs
1818  */
1819 //=============================================================================
1820 TCollection_AsciiString GEOMImpl_IMeasureOperations::WhatIs (Handle(GEOM_Object) theShape)
1821 {
1822   SetErrorCode(KO);
1823
1824   TCollection_AsciiString Astr;
1825
1826   if (theShape.IsNull()) return Astr;
1827
1828   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1829   if (aRefShape.IsNull()) return Astr;
1830
1831   TopoDS_Shape aShape = aRefShape->GetValue();
1832   if (aShape.IsNull()) {
1833     SetErrorCode("The Objects has NULL Shape");
1834     return Astr;
1835   }
1836
1837   //Compute the parameters
1838   if (aShape.ShapeType() == TopAbs_EDGE) {
1839     if (BRep_Tool::Degenerated(TopoDS::Edge(aShape))) {
1840       Astr = Astr + " It is a degenerated edge \n";
1841     }
1842   }
1843
1844   Astr = Astr + " Number of sub-shapes : \n";
1845
1846   try {
1847     OCC_CATCH_SIGNALS;
1848     int iType, nbTypes [TopAbs_SHAPE], nbFlatType [TopAbs_SHAPE];
1849     for (iType = 0; iType < TopAbs_SHAPE; ++iType) {
1850       nbTypes[iType] = 0;
1851       nbFlatType[iType] = 0;
1852     }
1853     nbTypes[aShape.ShapeType()]++;
1854
1855     TopTools_MapOfShape aMapOfShape;
1856     aMapOfShape.Add(aShape);
1857     TopTools_ListOfShape aListOfShape;
1858     aListOfShape.Append(aShape);
1859
1860     TopTools_ListIteratorOfListOfShape itL (aListOfShape);
1861     for (; itL.More(); itL.Next()) {
1862       TopoDS_Shape sp = itL.Value();
1863       TopoDS_Iterator it (sp);
1864       for (; it.More(); it.Next()) {
1865         TopoDS_Shape s = it.Value();
1866         if (aMapOfShape.Add(s)) {
1867           aListOfShape.Append(s);
1868           nbTypes[s.ShapeType()]++;
1869           if ((sp.ShapeType() == TopAbs_COMPOUND) || (sp.ShapeType() == TopAbs_COMPSOLID)) {
1870             nbFlatType[s.ShapeType()]++;
1871           }
1872         }
1873       }
1874     }
1875
1876     Astr = Astr + " VERTEX : " + TCollection_AsciiString(nbTypes[TopAbs_VERTEX]) + "\n";
1877     Astr = Astr + " EDGE : " + TCollection_AsciiString(nbTypes[TopAbs_EDGE]) + "\n";
1878     Astr = Astr + " WIRE : " + TCollection_AsciiString(nbTypes[TopAbs_WIRE]) + "\n";
1879     Astr = Astr + " FACE : " + TCollection_AsciiString(nbTypes[TopAbs_FACE]) + "\n";
1880     Astr = Astr + " SHELL : " + TCollection_AsciiString(nbTypes[TopAbs_SHELL]) + "\n";
1881     Astr = Astr + " SOLID : " + TCollection_AsciiString(nbTypes[TopAbs_SOLID]) + "\n";
1882     Astr = Astr + " COMPSOLID : " + TCollection_AsciiString(nbTypes[TopAbs_COMPSOLID]) + "\n";
1883     Astr = Astr + " COMPOUND : " + TCollection_AsciiString(nbTypes[TopAbs_COMPOUND]) + "\n";
1884     Astr = Astr + " SHAPE : " + TCollection_AsciiString(aMapOfShape.Extent()) + "\n";
1885
1886     if ((aShape.ShapeType() == TopAbs_COMPOUND) || (aShape.ShapeType() == TopAbs_COMPSOLID)){
1887       Astr = Astr + " --------------------- \n Flat content : \n";
1888       if (nbFlatType[TopAbs_VERTEX] > 0)
1889         Astr = Astr + " VERTEX : " + TCollection_AsciiString(nbFlatType[TopAbs_VERTEX]) + "\n";
1890       if (nbFlatType[TopAbs_EDGE] > 0)
1891         Astr = Astr + " EDGE : " + TCollection_AsciiString(nbFlatType[TopAbs_EDGE]) + "\n";
1892       if (nbFlatType[TopAbs_WIRE] > 0)
1893         Astr = Astr + " WIRE : " + TCollection_AsciiString(nbFlatType[TopAbs_WIRE]) + "\n";
1894       if (nbFlatType[TopAbs_FACE] > 0)
1895         Astr = Astr + " FACE : " + TCollection_AsciiString(nbFlatType[TopAbs_FACE]) + "\n";
1896       if (nbFlatType[TopAbs_SHELL] > 0)
1897         Astr = Astr + " SHELL : " + TCollection_AsciiString(nbFlatType[TopAbs_SHELL]) + "\n";
1898       if (nbFlatType[TopAbs_SOLID] > 0)
1899         Astr = Astr + " SOLID : " + TCollection_AsciiString(nbFlatType[TopAbs_SOLID]) + "\n";
1900     }
1901   }
1902   catch (Standard_Failure) {
1903     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1904     SetErrorCode(aFail->GetMessageString());
1905     return Astr;
1906   }
1907
1908   SetErrorCode(OK);
1909   return Astr;
1910 }
1911
1912 //=============================================================================
1913 /*!
1914  *  AreCoordsInside
1915  */
1916 //=============================================================================
1917 std::vector<bool>
1918 GEOMImpl_IMeasureOperations::AreCoordsInside(Handle(GEOM_Object)        theShape,
1919                                              const std::vector<double>& coords,
1920                                              double                     tolerance)
1921 {
1922   std::vector<bool> isInsideRes;
1923   if (!theShape.IsNull()) {
1924     Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
1925     if (!aRefShape.IsNull()) {
1926       TopoDS_Shape aShape = aRefShape->GetValue();
1927       if (!aShape.IsNull())
1928       {
1929         TopTools_IndexedMapOfShape mapShape;
1930         {
1931           TopExp_Explorer anExp;
1932           for ( anExp.Init( aShape, TopAbs_SOLID ); anExp.More(); anExp.Next() )
1933             mapShape.Add( anExp.Current() );
1934           for ( anExp.Init( aShape, TopAbs_FACE, TopAbs_SOLID ); anExp.More(); anExp.Next() )
1935             mapShape.Add( anExp.Current() );
1936           for ( anExp.Init( aShape, TopAbs_EDGE, TopAbs_FACE ); anExp.More(); anExp.Next() )
1937             mapShape.Add( anExp.Current() );
1938           for ( anExp.Init( aShape, TopAbs_VERTEX, TopAbs_EDGE ); anExp.More(); anExp.Next() )
1939             mapShape.Add( anExp.Current() ); //// ?????????
1940         }
1941         size_t nb_points = coords.size()/3, nb_points_inside = 0;
1942         isInsideRes.resize( nb_points, false );
1943
1944         for ( int iS = 1; iS <= mapShape.Extent(); ++iS )
1945         {
1946           if ( nb_points_inside == nb_points )
1947             break;
1948           aShape = mapShape( iS );
1949           switch ( aShape.ShapeType() ) {
1950           case TopAbs_SOLID:
1951           {
1952             BRepClass3d_SolidClassifier SC( TopoDS::Solid( aShape ));
1953             for ( size_t i = 0; i < nb_points; i++)
1954             {
1955               if ( isInsideRes[ i ]) continue;
1956               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
1957               SC.Perform( aPnt, tolerance );
1958               isInsideRes[ i ] = (( SC.State() == TopAbs_IN ) || ( SC.State() == TopAbs_ON ));
1959               nb_points_inside += isInsideRes[ i ];
1960             }
1961             break;
1962           }
1963           case TopAbs_FACE:
1964           {
1965             Standard_Real u1,u2,v1,v2;
1966             const TopoDS_Face&   face = TopoDS::Face( aShape );
1967             Handle(Geom_Surface) surf = BRep_Tool::Surface( face );
1968             surf->Bounds( u1,u2,v1,v2 );
1969             GeomAPI_ProjectPointOnSurf project;
1970             project.Init(surf, u1,u2, v1,v2, tolerance );
1971             for ( size_t i = 0; i < nb_points; i++)
1972             {
1973               if ( isInsideRes[ i ]) continue;
1974               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
1975               project.Perform( aPnt );
1976               if ( project.IsDone() &&
1977                    project.NbPoints() > 0 &&
1978                    project.LowerDistance() <= tolerance )
1979               {
1980                 Quantity_Parameter u, v;
1981                 project.LowerDistanceParameters(u, v);
1982                 gp_Pnt2d uv( u, v );
1983                 BRepClass_FaceClassifier FC ( face, uv, tolerance );
1984                 isInsideRes[ i ] = (( FC.State() == TopAbs_IN ) || ( FC.State() == TopAbs_ON ));
1985                 nb_points_inside += isInsideRes[ i ];
1986               }
1987             }
1988             break;
1989           }
1990           case TopAbs_EDGE:
1991           {
1992             Standard_Real f,l;
1993             const TopoDS_Edge&  edge = TopoDS::Edge( aShape );
1994             Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, f, l );
1995             GeomAPI_ProjectPointOnCurve project;
1996             project.Init( curve, f, l );
1997             for ( size_t i = 0; i < nb_points; i++)
1998             {
1999               if ( isInsideRes[ i ]) continue;
2000               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
2001               project.Perform( aPnt );
2002               isInsideRes[ i ] = ( project.NbPoints() > 0 &&
2003                                    project.LowerDistance() <= tolerance );
2004               nb_points_inside += isInsideRes[ i ];
2005             }
2006             break;
2007           }
2008           case TopAbs_VERTEX:
2009           {
2010             gp_Pnt aVPnt = BRep_Tool::Pnt( TopoDS::Vertex( aShape ));
2011             for ( size_t i = 0; i < nb_points; i++)
2012             {
2013               if ( isInsideRes[ i ]) continue;
2014               gp_Pnt aPnt( coords[3*i], coords[3*i+1], coords[3*i+2] );
2015               isInsideRes[ i ] = ( aPnt.SquareDistance( aVPnt ) <= tolerance * tolerance );
2016               nb_points_inside += isInsideRes[ i ];
2017             }
2018             break;
2019           }
2020           default:;
2021           } // switch ( aShape.ShapeType() )
2022         }
2023       }
2024     }
2025   }
2026   return isInsideRes;
2027 }
2028
2029 //=============================================================================
2030 /*!
2031  *  GetMinDistance
2032  */
2033 //=============================================================================
2034 Standard_Real
2035 GEOMImpl_IMeasureOperations::GetMinDistance (Handle(GEOM_Object) theShape1,
2036                                              Handle(GEOM_Object) theShape2,
2037                                              Standard_Real& X1,
2038                                              Standard_Real& Y1,
2039                                              Standard_Real& Z1,
2040                                              Standard_Real& X2,
2041                                              Standard_Real& Y2,
2042                                              Standard_Real& Z2)
2043 {
2044   SetErrorCode(KO);
2045   Standard_Real MinDist = 1.e9;
2046
2047   if (theShape1.IsNull() || theShape2.IsNull()) return MinDist;
2048
2049   Handle(GEOM_Function) aRefShape1 = theShape1->GetLastFunction();
2050   Handle(GEOM_Function) aRefShape2 = theShape2->GetLastFunction();
2051   if (aRefShape1.IsNull() || aRefShape2.IsNull()) return MinDist;
2052
2053   TopoDS_Shape aShape1 = aRefShape1->GetValue();
2054   TopoDS_Shape aShape2 = aRefShape2->GetValue();
2055   if (aShape1.IsNull() || aShape2.IsNull()) {
2056     SetErrorCode("One of Objects has NULL Shape");
2057     return MinDist;
2058   }
2059
2060   //Compute the parameters
2061   try {
2062     OCC_CATCH_SIGNALS;
2063
2064     gp_Pnt aPnt1, aPnt2;
2065
2066     MinDist = GEOMUtils::GetMinDistance(aShape1, aShape2, aPnt1, aPnt2);
2067
2068     if (MinDist >= 0.0) {
2069       aPnt1.Coord(X1, Y1, Z1);
2070       aPnt2.Coord(X2, Y2, Z2);
2071     } else {
2072       return MinDist;
2073     }
2074   }
2075   catch (Standard_Failure) {
2076     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2077     SetErrorCode(aFail->GetMessageString());
2078     return MinDist;
2079   }
2080
2081   SetErrorCode(OK);
2082   return MinDist;
2083 }
2084
2085 //=======================================================================
2086 /*!
2087  *  Get coordinates of closest points of two shapes
2088  */
2089 //=======================================================================
2090 Standard_Integer GEOMImpl_IMeasureOperations::ClosestPoints (Handle(GEOM_Object) theShape1,
2091                                                              Handle(GEOM_Object) theShape2,
2092                                                              Handle(TColStd_HSequenceOfReal)& theDoubles)
2093 {
2094   SetErrorCode(KO);
2095   Standard_Integer nbSolutions = 0;
2096
2097   if (theShape1.IsNull() || theShape2.IsNull()) return nbSolutions;
2098
2099   Handle(GEOM_Function) aRefShape1 = theShape1->GetLastFunction();
2100   Handle(GEOM_Function) aRefShape2 = theShape2->GetLastFunction();
2101   if (aRefShape1.IsNull() || aRefShape2.IsNull()) return nbSolutions;
2102
2103   TopoDS_Shape aShape1 = aRefShape1->GetValue();
2104   TopoDS_Shape aShape2 = aRefShape2->GetValue();
2105   if (aShape1.IsNull() || aShape2.IsNull()) {
2106     SetErrorCode("One of Objects has NULL Shape");
2107     return nbSolutions;
2108   }
2109
2110   // Compute the extremities
2111   try {
2112     OCC_CATCH_SIGNALS;
2113
2114     // skl 30.06.2008
2115     // additional workaround for bugs 19899, 19908 and 19910 from Mantis
2116     gp_Pnt P1, P2;
2117     double dist = GEOMUtils::GetMinDistanceSingular(aShape1, aShape2, P1, P2);
2118     if (dist > -1.0) {
2119       nbSolutions = 1;
2120
2121       theDoubles->Append(P1.X());
2122       theDoubles->Append(P1.Y());
2123       theDoubles->Append(P1.Z());
2124       theDoubles->Append(P2.X());
2125       theDoubles->Append(P2.Y());
2126       theDoubles->Append(P2.Z());
2127
2128       SetErrorCode(OK);
2129       return nbSolutions;
2130     }
2131
2132     BRepExtrema_DistShapeShape dst (aShape1, aShape2);
2133     if (dst.IsDone()) {
2134       nbSolutions = dst.NbSolution();
2135       if (theDoubles.IsNull()) theDoubles = new TColStd_HSequenceOfReal;
2136
2137       gp_Pnt P1, P2;
2138       for (int i = 1; i <= nbSolutions; i++) {
2139         P1 = dst.PointOnShape1(i);
2140         P2 = dst.PointOnShape2(i);
2141
2142         theDoubles->Append(P1.X());
2143         theDoubles->Append(P1.Y());
2144         theDoubles->Append(P1.Z());
2145         theDoubles->Append(P2.X());
2146         theDoubles->Append(P2.Y());
2147         theDoubles->Append(P2.Z());
2148       }
2149     }
2150   }
2151   catch (Standard_Failure) {
2152     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2153     SetErrorCode(aFail->GetMessageString());
2154     return nbSolutions;
2155   }
2156
2157   SetErrorCode(OK);
2158   return nbSolutions;
2159 }
2160
2161 //=======================================================================
2162 /*!
2163  *  Get coordinates of point
2164  */
2165 //=======================================================================
2166 void GEOMImpl_IMeasureOperations::PointCoordinates (Handle(GEOM_Object) theShape,
2167                         Standard_Real& theX, Standard_Real& theY, Standard_Real& theZ)
2168 {
2169   SetErrorCode(KO);
2170
2171   if (theShape.IsNull())
2172     return;
2173
2174   Handle(GEOM_Function) aRefShape = theShape->GetLastFunction();
2175   if (aRefShape.IsNull())
2176     return;
2177
2178   TopoDS_Shape aShape = aRefShape->GetValue();
2179   if (aShape.IsNull() || aShape.ShapeType() != TopAbs_VERTEX)
2180   {
2181     SetErrorCode( "Shape must be a vertex" );
2182     return;
2183   }
2184
2185   try {
2186     OCC_CATCH_SIGNALS;
2187     gp_Pnt aPnt = BRep_Tool::Pnt( TopoDS::Vertex( aShape ) );
2188     theX = aPnt.X();
2189     theY = aPnt.Y();
2190     theZ = aPnt.Z();
2191
2192     SetErrorCode(OK);
2193   }
2194   catch (Standard_Failure)
2195   {
2196     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2197     SetErrorCode( aFail->GetMessageString() );
2198   }
2199 }
2200
2201 //=======================================================================
2202 /*!
2203  *  Compute angle (in degrees) between two lines
2204  */
2205 //=======================================================================
2206 Standard_Real GEOMImpl_IMeasureOperations::GetAngle (Handle(GEOM_Object) theLine1,
2207                                                      Handle(GEOM_Object) theLine2)
2208 {
2209   if (theLine1->GetType() == GEOM_VECTOR &&
2210       theLine2->GetType() == GEOM_VECTOR)
2211     return GetAngleBtwVectors(theLine1, theLine2);
2212
2213   SetErrorCode(KO);
2214
2215   Standard_Real anAngle = -1.0;
2216
2217   if (theLine1.IsNull() || theLine2.IsNull())
2218     return anAngle;
2219
2220   Handle(GEOM_Function) aRefLine1 = theLine1->GetLastFunction();
2221   Handle(GEOM_Function) aRefLine2 = theLine2->GetLastFunction();
2222   if (aRefLine1.IsNull() || aRefLine2.IsNull())
2223     return anAngle;
2224
2225   TopoDS_Shape aLine1 = aRefLine1->GetValue();
2226   TopoDS_Shape aLine2 = aRefLine2->GetValue();
2227   if (aLine1.IsNull() || aLine2.IsNull() ||
2228       aLine1.ShapeType() != TopAbs_EDGE ||
2229       aLine2.ShapeType() != TopAbs_EDGE)
2230   {
2231     SetErrorCode("Two edges must be given");
2232     return anAngle;
2233   }
2234
2235   try {
2236     OCC_CATCH_SIGNALS;
2237     TopoDS_Edge E1 = TopoDS::Edge(aLine1);
2238     TopoDS_Edge E2 = TopoDS::Edge(aLine2);
2239
2240     double fp,lp;
2241     Handle(Geom_Curve) C1 = BRep_Tool::Curve(E1,fp,lp);
2242     Handle(Geom_Curve) C2 = BRep_Tool::Curve(E2,fp,lp);
2243
2244     if ( C1.IsNull() || C2.IsNull() ||
2245         !C1->IsKind(STANDARD_TYPE(Geom_Line)) ||
2246         !C2->IsKind(STANDARD_TYPE(Geom_Line)))
2247     {
2248       SetErrorCode("The edges must be linear");
2249       return anAngle;
2250     }
2251
2252     Handle(Geom_Line) L1 = Handle(Geom_Line)::DownCast(C1);
2253     Handle(Geom_Line) L2 = Handle(Geom_Line)::DownCast(C2);
2254
2255     gp_Lin aLin1 = L1->Lin();
2256     gp_Lin aLin2 = L2->Lin();
2257
2258     anAngle = aLin1.Angle(aLin2);
2259     anAngle *= 180. / M_PI; // convert radians into degrees
2260
2261     if (anAngle > 90.0) {
2262       anAngle = 180.0 - anAngle;
2263     }
2264
2265     SetErrorCode(OK);
2266   }
2267   catch (Standard_Failure)
2268   {
2269     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2270     SetErrorCode(aFail->GetMessageString());
2271   }
2272
2273   return anAngle;
2274 }
2275
2276 //=======================================================================
2277 /*!
2278  *  Compute angle (in degrees) between two vectors
2279  */
2280 //=======================================================================
2281 Standard_Real GEOMImpl_IMeasureOperations::GetAngleBtwVectors (Handle(GEOM_Object) theVec1,
2282                                                                Handle(GEOM_Object) theVec2)
2283 {
2284   SetErrorCode(KO);
2285
2286   Standard_Real anAngle = -1.0;
2287
2288   if (theVec1.IsNull() || theVec2.IsNull())
2289     return anAngle;
2290
2291   if (theVec1->GetType() != GEOM_VECTOR || theVec2->GetType() != GEOM_VECTOR) {
2292     SetErrorCode("Two vectors must be given");
2293     return anAngle;
2294   }
2295
2296   Handle(GEOM_Function) aRefVec1 = theVec1->GetLastFunction();
2297   Handle(GEOM_Function) aRefVec2 = theVec2->GetLastFunction();
2298   if (aRefVec1.IsNull() || aRefVec2.IsNull())
2299     return anAngle;
2300
2301   TopoDS_Shape aVec1 = aRefVec1->GetValue();
2302   TopoDS_Shape aVec2 = aRefVec2->GetValue();
2303   if (aVec1.IsNull() || aVec2.IsNull() ||
2304       aVec1.ShapeType() != TopAbs_EDGE ||
2305       aVec2.ShapeType() != TopAbs_EDGE)
2306   {
2307     SetErrorCode("Two edges must be given");
2308     return anAngle;
2309   }
2310
2311   try {
2312     OCC_CATCH_SIGNALS;
2313     TopoDS_Edge aE1 = TopoDS::Edge(aVec1);
2314     TopoDS_Edge aE2 = TopoDS::Edge(aVec2);
2315
2316     TopoDS_Vertex aP11, aP12, aP21, aP22;
2317     TopExp::Vertices(aE1, aP11, aP12, Standard_True);
2318     TopExp::Vertices(aE2, aP21, aP22, Standard_True);
2319     if (aP11.IsNull() || aP12.IsNull() || aP21.IsNull() || aP22.IsNull()) {
2320       SetErrorCode("Bad edge given");
2321       return anAngle;
2322     }
2323
2324     gp_Vec aV1 (BRep_Tool::Pnt(aP11), BRep_Tool::Pnt(aP12));
2325     gp_Vec aV2 (BRep_Tool::Pnt(aP21), BRep_Tool::Pnt(aP22)) ;
2326
2327     anAngle = aV1.Angle(aV2);
2328     anAngle *= 180. / M_PI; // convert radians into degrees
2329
2330     SetErrorCode(OK);
2331   }
2332   catch (Standard_Failure)
2333   {
2334     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2335     SetErrorCode(aFail->GetMessageString());
2336   }
2337
2338   return anAngle;
2339 }
2340
2341
2342 //=============================================================================
2343 /*!
2344  *  CurveCurvatureByParam
2345  */
2346 //=============================================================================
2347 Standard_Real GEOMImpl_IMeasureOperations::CurveCurvatureByParam
2348                         (Handle(GEOM_Object) theCurve, Standard_Real& theParam)
2349 {
2350   SetErrorCode(KO);
2351   Standard_Real aRes = -1.0;
2352
2353   if(theCurve.IsNull()) return aRes;
2354
2355   Handle(GEOM_Function) aRefShape = theCurve->GetLastFunction();
2356   if(aRefShape.IsNull()) return aRes;
2357
2358   TopoDS_Shape aShape = aRefShape->GetValue();
2359   if(aShape.IsNull()) {
2360     SetErrorCode("One of Objects has NULL Shape");
2361     return aRes;
2362   }
2363
2364   Standard_Real aFP, aLP, aP;
2365   Handle(Geom_Curve) aCurve = BRep_Tool::Curve(TopoDS::Edge(aShape), aFP, aLP);
2366   aP = aFP + (aLP - aFP) * theParam;
2367
2368   if(aCurve.IsNull()) return aRes;
2369
2370   //Compute curvature
2371   try {
2372     OCC_CATCH_SIGNALS;
2373     GeomLProp_CLProps Prop = GeomLProp_CLProps
2374       (aCurve, aP, 2, Precision::Confusion());
2375     aRes = fabs(Prop.Curvature());
2376     SetErrorCode(OK);
2377   }
2378   catch (Standard_Failure) {
2379     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2380     SetErrorCode(aFail->GetMessageString());
2381     return aRes;
2382   }
2383
2384   if( aRes > Precision::Confusion() )
2385     aRes = 1/aRes;
2386   else
2387     aRes = RealLast();
2388
2389   return aRes;
2390 }
2391
2392
2393 //=============================================================================
2394 /*!
2395  *  CurveCurvatureByPoint
2396  */
2397 //=============================================================================
2398 Standard_Real GEOMImpl_IMeasureOperations::CurveCurvatureByPoint
2399                    (Handle(GEOM_Object) theCurve, Handle(GEOM_Object) thePoint)
2400 {
2401   SetErrorCode(KO);
2402   Standard_Real aRes = -1.0;
2403
2404   if( theCurve.IsNull() || thePoint.IsNull() ) return aRes;
2405
2406   Handle(GEOM_Function) aRefCurve = theCurve->GetLastFunction();
2407   Handle(GEOM_Function) aRefPoint = thePoint->GetLastFunction();
2408   if( aRefCurve.IsNull() || aRefPoint.IsNull() ) return aRes;
2409
2410   TopoDS_Edge anEdge = TopoDS::Edge(aRefCurve->GetValue());
2411   TopoDS_Vertex aPnt = TopoDS::Vertex(aRefPoint->GetValue());
2412   if( anEdge.IsNull() || aPnt.IsNull() ) {
2413     SetErrorCode("One of Objects has NULL Shape");
2414     return aRes;
2415   }
2416
2417   Standard_Real aFP, aLP;
2418   Handle(Geom_Curve) aCurve = BRep_Tool::Curve(anEdge, aFP, aLP);
2419   if(aCurve.IsNull()) return aRes;
2420   gp_Pnt aPoint = BRep_Tool::Pnt(aPnt);
2421
2422   //Compute curvature
2423   try {
2424     OCC_CATCH_SIGNALS;
2425     GeomAPI_ProjectPointOnCurve PPCurve(aPoint, aCurve, aFP, aLP);
2426     if(PPCurve.NbPoints()>0) {
2427       GeomLProp_CLProps Prop = GeomLProp_CLProps
2428         (aCurve, PPCurve.LowerDistanceParameter(), 2, Precision::Confusion());
2429       aRes = fabs(Prop.Curvature());
2430       SetErrorCode(OK);
2431     }
2432   }
2433   catch (Standard_Failure) {
2434     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2435     SetErrorCode(aFail->GetMessageString());
2436     return aRes;
2437   }
2438
2439   if( aRes > Precision::Confusion() )
2440     aRes = 1/aRes;
2441   else
2442     aRes = RealLast();
2443
2444   return aRes;
2445 }
2446
2447
2448 //=============================================================================
2449 /*!
2450  *  getSurfaceCurvatures
2451  */
2452 //=============================================================================
2453 Standard_Real GEOMImpl_IMeasureOperations::getSurfaceCurvatures
2454                                           (const Handle(Geom_Surface)& aSurf,
2455                                            Standard_Real theUParam,
2456                                            Standard_Real theVParam,
2457                                            Standard_Boolean theNeedMaxCurv)
2458 {
2459   SetErrorCode(KO);
2460   Standard_Real aRes = 1.0;
2461
2462   if (aSurf.IsNull()) return aRes;
2463
2464   try {
2465     OCC_CATCH_SIGNALS;
2466     GeomLProp_SLProps Prop = GeomLProp_SLProps
2467       (aSurf, theUParam, theVParam, 2, Precision::Confusion());
2468     if(Prop.IsCurvatureDefined()) {
2469       if(Prop.IsUmbilic()) {
2470         //cout<<"is umbilic"<<endl;
2471         aRes = fabs(Prop.MeanCurvature());
2472       }
2473       else {
2474         //cout<<"is not umbilic"<<endl;
2475         double c1 = fabs(Prop.MaxCurvature());
2476         double c2 = fabs(Prop.MinCurvature());
2477         if(theNeedMaxCurv)
2478           aRes = Max(c1,c2);
2479         else
2480           aRes = Min(c1,c2);
2481       }
2482       SetErrorCode(OK);
2483     }
2484   }
2485   catch (Standard_Failure) {
2486     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2487     SetErrorCode(aFail->GetMessageString());
2488     return aRes;
2489   }
2490
2491   if( fabs(aRes) > Precision::Confusion() )
2492     aRes = 1/aRes;
2493   else
2494     aRes = RealLast();
2495
2496   return aRes;
2497 }
2498
2499
2500 //=============================================================================
2501 /*!
2502  *  MaxSurfaceCurvatureByParam
2503  */
2504 //=============================================================================
2505 Standard_Real GEOMImpl_IMeasureOperations::MaxSurfaceCurvatureByParam
2506                                                   (Handle(GEOM_Object) theSurf,
2507                                                    Standard_Real& theUParam,
2508                                                    Standard_Real& theVParam)
2509 {
2510   SetErrorCode(KO);
2511   Standard_Real aRes = -1.0;
2512
2513   if (theSurf.IsNull()) return aRes;
2514
2515   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2516   if(aRefShape.IsNull()) return aRes;
2517
2518   TopoDS_Shape aShape = aRefShape->GetValue();
2519   if(aShape.IsNull()) {
2520     SetErrorCode("One of Objects has NULL Shape");
2521     return aRes;
2522   }
2523
2524   TopoDS_Face F = TopoDS::Face(aShape);
2525   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(F);
2526
2527   //Compute the parameters
2528   Standard_Real U1,U2,V1,V2;
2529   ShapeAnalysis::GetFaceUVBounds(F,U1,U2,V1,V2);
2530   Standard_Real U = U1 + (U2-U1)*theUParam;
2531   Standard_Real V = V1 + (V2-V1)*theVParam;
2532
2533   return getSurfaceCurvatures(aSurf, U, V, true);
2534 }
2535
2536
2537 //=============================================================================
2538 /*!
2539  *  MaxSurfaceCurvatureByPoint
2540  */
2541 //=============================================================================
2542 Standard_Real GEOMImpl_IMeasureOperations::MaxSurfaceCurvatureByPoint
2543                     (Handle(GEOM_Object) theSurf, Handle(GEOM_Object) thePoint)
2544 {
2545   SetErrorCode(KO);
2546   Standard_Real aRes = -1.0;
2547
2548   if( theSurf.IsNull() || thePoint.IsNull() ) return aRes;
2549
2550   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2551   Handle(GEOM_Function) aRefPoint = thePoint->GetLastFunction();
2552   if( aRefShape.IsNull() || aRefPoint.IsNull() ) return aRes;
2553
2554   TopoDS_Face aFace = TopoDS::Face(aRefShape->GetValue());
2555   TopoDS_Vertex aPnt = TopoDS::Vertex(aRefPoint->GetValue());
2556   if( aFace.IsNull() || aPnt.IsNull() ) {
2557     SetErrorCode("One of Objects has NULL Shape");
2558     return 0;
2559   }
2560
2561   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace);
2562   if(aSurf.IsNull()) return aRes;
2563   gp_Pnt aPoint = BRep_Tool::Pnt(aPnt);
2564
2565   //Compute the parameters
2566   ShapeAnalysis_Surface sas(aSurf);
2567   gp_Pnt2d UV = sas.ValueOfUV(aPoint,Precision::Confusion());
2568
2569   return getSurfaceCurvatures(aSurf, UV.X(), UV.Y(), true);
2570 }
2571
2572
2573 //=============================================================================
2574 /*!
2575  *  MinSurfaceCurvatureByParam
2576  */
2577 //=============================================================================
2578 Standard_Real GEOMImpl_IMeasureOperations::MinSurfaceCurvatureByParam
2579                                                   (Handle(GEOM_Object) theSurf,
2580                                                    Standard_Real& theUParam,
2581                                                    Standard_Real& theVParam)
2582 {
2583   SetErrorCode(KO);
2584   Standard_Real aRes = -1.0;
2585
2586   if (theSurf.IsNull()) return aRes;
2587
2588   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2589   if(aRefShape.IsNull()) return aRes;
2590
2591   TopoDS_Shape aShape = aRefShape->GetValue();
2592   if(aShape.IsNull()) {
2593     SetErrorCode("One of Objects has NULL Shape");
2594     return aRes;
2595   }
2596
2597   TopoDS_Face F = TopoDS::Face(aShape);
2598   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(F);
2599
2600   //Compute the parameters
2601   Standard_Real U1,U2,V1,V2;
2602   ShapeAnalysis::GetFaceUVBounds(F,U1,U2,V1,V2);
2603   Standard_Real U = U1 + (U2-U1)*theUParam;
2604   Standard_Real V = V1 + (V2-V1)*theVParam;
2605
2606   return getSurfaceCurvatures(aSurf, U, V, false);
2607 }
2608
2609
2610 //=============================================================================
2611 /*!
2612  *  MinSurfaceCurvatureByPoint
2613  */
2614 //=============================================================================
2615 Standard_Real GEOMImpl_IMeasureOperations::MinSurfaceCurvatureByPoint
2616                     (Handle(GEOM_Object) theSurf, Handle(GEOM_Object) thePoint)
2617 {
2618   SetErrorCode(KO);
2619   Standard_Real aRes = -1.0;
2620
2621   if( theSurf.IsNull() || thePoint.IsNull() ) return aRes;
2622
2623   Handle(GEOM_Function) aRefShape = theSurf->GetLastFunction();
2624   Handle(GEOM_Function) aRefPoint = thePoint->GetLastFunction();
2625   if( aRefShape.IsNull() || aRefPoint.IsNull() ) return aRes;
2626
2627   TopoDS_Face aFace = TopoDS::Face(aRefShape->GetValue());
2628   TopoDS_Vertex aPnt = TopoDS::Vertex(aRefPoint->GetValue());
2629   if( aFace.IsNull() || aPnt.IsNull() ) {
2630     SetErrorCode("One of Objects has NULL Shape");
2631     return 0;
2632   }
2633
2634   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace);
2635   if(aSurf.IsNull()) return aRes;
2636   gp_Pnt aPoint = BRep_Tool::Pnt(aPnt);
2637
2638   //Compute the parameters
2639   ShapeAnalysis_Surface sas(aSurf);
2640   gp_Pnt2d UV = sas.ValueOfUV(aPoint,Precision::Confusion());
2641
2642   return getSurfaceCurvatures(aSurf, UV.X(), UV.Y(), false);
2643 }
2644
2645 //=======================================================================
2646 //function : FillErrorsSub
2647 //purpose  : Fill the errors list of subshapes on shape.
2648 //=======================================================================
2649 void GEOMImpl_IMeasureOperations::FillErrorsSub
2650            (const BRepCheck_Analyzer                   &theAna,
2651             const TopoDS_Shape                         &theShape,
2652             const TopAbs_ShapeEnum                     theSubType,
2653                   TopTools_DataMapOfIntegerListOfShape &theMapErrors) const
2654 {
2655   TopExp_Explorer anExp(theShape, theSubType);
2656   TopTools_MapOfShape aMapSubShapes;
2657
2658   for (; anExp.More(); anExp.Next()) {
2659     const TopoDS_Shape &aSubShape = anExp.Current();
2660
2661     if (aMapSubShapes.Add(aSubShape)) {
2662       const Handle(BRepCheck_Result) &aRes = theAna.Result(aSubShape);
2663
2664       for (aRes->InitContextIterator();
2665            aRes->MoreShapeInContext(); 
2666            aRes->NextShapeInContext()) {
2667         if (aRes->ContextualShape().IsSame(theShape)) {
2668           BRepCheck_ListIteratorOfListOfStatus itl(aRes->StatusOnShape());
2669
2670           if (itl.Value() != BRepCheck_NoError) {
2671             // Add all errors for theShape and its sub-shape.
2672             for (;itl.More(); itl.Next()) {
2673               const Standard_Integer aStat = (Standard_Integer)itl.Value();
2674
2675               if (!theMapErrors.IsBound(aStat)) {
2676                 TopTools_ListOfShape anEmpty;
2677
2678                 theMapErrors.Bind(aStat, anEmpty);
2679               }
2680
2681               TopTools_ListOfShape &theShapes = theMapErrors.ChangeFind(aStat);
2682
2683               theShapes.Append(aSubShape);
2684               theShapes.Append(theShape);
2685             }
2686           }
2687         }
2688
2689         break;
2690       }
2691     }
2692   }
2693 }
2694
2695 //=======================================================================
2696 //function : FillErrors
2697 //purpose  : Fill the errors list.
2698 //=======================================================================
2699 void GEOMImpl_IMeasureOperations::FillErrors
2700              (const BRepCheck_Analyzer                   &theAna,
2701               const TopoDS_Shape                         &theShape,
2702                     TopTools_DataMapOfIntegerListOfShape &theMapErrors,
2703                     TopTools_MapOfShape                  &theMapShapes) const
2704 {
2705   if (theMapShapes.Add(theShape)) {
2706     // Fill errors of child shapes.
2707     for (TopoDS_Iterator iter(theShape); iter.More(); iter.Next()) {
2708       FillErrors(theAna, iter.Value(), theMapErrors, theMapShapes);
2709     }
2710
2711     // Fill errors of theShape.
2712     const Handle(BRepCheck_Result) &aRes = theAna.Result(theShape);
2713
2714     if (!aRes.IsNull()) {
2715       BRepCheck_ListIteratorOfListOfStatus itl(aRes->Status());
2716
2717       if (itl.Value() != BRepCheck_NoError) {
2718         // Add all errors for theShape.
2719         for (;itl.More(); itl.Next()) {
2720           const Standard_Integer aStat = (Standard_Integer)itl.Value();
2721
2722           if (!theMapErrors.IsBound(aStat)) {
2723             TopTools_ListOfShape anEmpty;
2724
2725             theMapErrors.Bind(aStat, anEmpty);
2726           }
2727
2728           theMapErrors.ChangeFind(aStat).Append(theShape);
2729         }
2730       }
2731     }
2732
2733     // Add errors of subshapes on theShape.
2734     const TopAbs_ShapeEnum aType = theShape.ShapeType();
2735
2736     switch (aType) {
2737     case TopAbs_EDGE:
2738       FillErrorsSub(theAna, theShape, TopAbs_VERTEX, theMapErrors);
2739       break;
2740     case TopAbs_FACE:
2741       FillErrorsSub(theAna, theShape, TopAbs_WIRE, theMapErrors);
2742       FillErrorsSub(theAna, theShape, TopAbs_EDGE, theMapErrors);
2743       FillErrorsSub(theAna, theShape, TopAbs_VERTEX, theMapErrors);
2744       break;
2745     case TopAbs_SOLID:
2746       FillErrorsSub(theAna, theShape, TopAbs_SHELL, theMapErrors);
2747       break;
2748     default:
2749       break;
2750     }
2751   }
2752 }
2753
2754 //=======================================================================
2755 //function : FillErrors
2756 //purpose  : Fill the errors list.
2757 //=======================================================================
2758 void GEOMImpl_IMeasureOperations::FillErrors
2759                   (const BRepCheck_Analyzer    &theAna,
2760                    const TopoDS_Shape          &theShape,
2761                          std::list<ShapeError> &theErrors) const
2762 {
2763   // Fill the errors map.
2764   TopTools_DataMapOfIntegerListOfShape aMapErrors;
2765   TopTools_MapOfShape                  aMapShapes;
2766
2767   FillErrors(theAna, theShape, aMapErrors, aMapShapes);
2768
2769   // Map sub-shapes and their indices
2770   TopTools_IndexedMapOfShape anIndices;
2771
2772   TopExp::MapShapes(theShape, anIndices);
2773
2774   TopTools_DataMapIteratorOfDataMapOfIntegerListOfShape aMapIter(aMapErrors);
2775
2776   for (; aMapIter.More(); aMapIter.Next()) {
2777     ShapeError anError;
2778
2779     anError.error = (BRepCheck_Status)aMapIter.Key();
2780
2781     TopTools_ListIteratorOfListOfShape aListIter(aMapIter.Value());
2782     TopTools_MapOfShape                aMapUnique;
2783
2784     for (; aListIter.More(); aListIter.Next()) {
2785       const TopoDS_Shape &aShape = aListIter.Value();
2786
2787       if (aMapUnique.Add(aShape)) {
2788         const Standard_Integer anIndex = anIndices.FindIndex(aShape);
2789
2790         anError.incriminated.push_back(anIndex);
2791       }
2792     }
2793
2794     if (!anError.incriminated.empty()) {
2795       theErrors.push_back(anError);
2796     }
2797   }
2798 }