1 // Copyright (C) 2007-2022 CEA/DEN, EDF R&D, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
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.
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.
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
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 #include <Standard_Stream.hxx>
25 #include <GEOMImpl_MeasureDriver.hxx>
26 #include <GEOMImpl_IMeasure.hxx>
27 #include <GEOMImpl_Types.hxx>
29 #include <GEOM_Function.hxx>
31 #include <GEOMUtils.hxx>
33 #include <BRep_Tool.hxx>
34 #include <BRepTools.hxx>
35 #include <BRepGProp.hxx>
36 #include <BRepLProp_SLProps.hxx>
37 #include <BRepBuilderAPI_Copy.hxx>
38 #include <BRepBuilderAPI_MakeEdge.hxx>
39 #include <BRepBuilderAPI_MakeVertex.hxx>
40 #include <BRepPrimAPI_MakeBox.hxx>
41 #include <BRepOffsetAPI_NormalProjection.hxx>
45 #include <TopoDS_Shape.hxx>
46 #include <TopoDS_Edge.hxx>
47 #include <TopoDS_Wire.hxx>
49 #include <TopExp_Explorer.hxx>
50 #include <TopTools_IndexedMapOfShape.hxx>
52 #include <GProp_GProps.hxx>
53 #include <GeomLProp_SLProps.hxx>
54 #include <GeomLProp_CLProps.hxx>
55 #include <Geom_Surface.hxx>
56 #include <GeomAPI_ProjectPointOnCurve.hxx>
58 #include <Bnd_Box.hxx>
59 #include <BRepBndLib.hxx>
60 #include <BRepAdaptor_Surface.hxx>
61 #include <ShapeAnalysis_Surface.hxx>
64 #include <Precision.hxx>
65 #include <Standard_NullObject.hxx>
66 #include <StdFail_NotDone.hxx>
68 //=======================================================================
71 //=======================================================================
72 const Standard_GUID& GEOMImpl_MeasureDriver::GetID()
74 static Standard_GUID aMeasureDriver("FF1BBB65-5D14-4df2-980B-3A668264EA16");
75 return aMeasureDriver;
79 //=======================================================================
80 //function : GEOMImpl_MeasureDriver
82 //=======================================================================
83 GEOMImpl_MeasureDriver::GEOMImpl_MeasureDriver()
87 //! This function is designed to evaluate normal curvature of the surface
88 //! in the given point along the given direction.
89 //! param[in] theFace face of interest.
90 //! param[in] thePoint point of interest.
91 //! param[in] theDir edge, giving the direction of interest.
92 //! return Edge, representing the curvature vector
93 TopoDS_Shape EvaluateAlongCurvature(const TopoDS_Shape& theFace,
94 const TopoDS_Shape& thePoint,
95 const TopoDS_Shape& theDir)
98 if (thePoint.IsNull())
99 Standard_NullObject::Raise("Point for curvature measurement is null");
100 if (thePoint.ShapeType() != TopAbs_VERTEX)
101 Standard_TypeMismatch::Raise("Point for curvature calculation is not a vertex");
102 gp_Pnt aPnt = BRep_Tool::Pnt(TopoDS::Vertex(thePoint));
104 // Point projection on the face
106 aPnt = GEOMUtils::ProjectPointOnFace(aPnt, theFace, U, V);
110 TopoDS_Face aFace = TopoDS::Face(theFace);
111 gp_Vec aV = GEOMUtils::GetVector(theDir, Standard_False);
113 // Calculate differential properties
114 BRepAdaptor_Surface aSurfAdapt (aFace);
115 BRepLProp_SLProps Props (aSurfAdapt, UV.X(), UV.Y(), 2, 1e-7);
116 if (!Props.IsNormalDefined())
117 Standard_ConstructionError::Raise
118 ("Curvature calculation failed: normal direction is not defined");
120 // Get differential properties
121 gp_Vec n = Props.Normal();
122 if (aV.IsParallel(n, Precision::Angular()))
123 Standard_ConstructionError::Raise
124 ("Curvature calculation failed: direction is normal to the face");
126 if (!Props.IsCurvatureDefined())
127 Standard_ConstructionError::Raise
128 ("Curvature calculation failed: curvature is not defined");
131 Standard_Real k = 0.;
133 if (Props.IsUmbilic()) {
134 k = Props.MaxCurvature();
137 gp_Dir aMaxDir, aMinDir;
138 Props.CurvatureDirections(aMaxDir, aMinDir);
140 gp_Vec2d T (aV.Dot(aMinDir), aV.Dot(aMaxDir));
141 if (Abs(T.X()) < Precision::Confusion() &&
142 Abs(T.Y()) < Precision::Confusion())
143 Standard_ConstructionError::Raise
144 ("Curvature calculation failed: direction is normal to the face");
146 Standard_Real cosAlpha = aDirT.X();
148 Standard_Real aMinCurv = Props.MinCurvature();
149 Standard_Real aMaxCurv = Props.MaxCurvature();
151 // Euler formula: k = k1 * cos^2(alpha) + k2 * sin^2(alpha)
152 // k = k2 + (k1 - k2) * cos^2(alpha)
153 k = aMaxCurv + (aMinCurv - aMaxCurv) * cosAlpha * cosAlpha;
156 if (Abs(k) < Precision::Confusion())
157 Standard_Failure::Raise("ZERO_CURVATURE");
159 // Radius or -radius of curvature
160 Standard_Real r = 1.0 / k;
164 gp_Pnt aPntEnd (aPnt.XYZ() + aNormal.XYZ() * r);
165 BRepBuilderAPI_MakeEdge aBuilder (aPnt, aPntEnd);
166 if (!aBuilder.IsDone())
167 Standard_ConstructionError::Raise("Curvature calculation failed: edge is not built");
169 return aBuilder.Shape();
172 //=======================================================================
175 //=======================================================================
176 Standard_Integer GEOMImpl_MeasureDriver::Execute(Handle(TFunction_Logbook)& log) const
178 if (Label().IsNull()) return 0;
179 Handle(GEOM_Function) aFunction = GEOM_Function::GetFunction(Label());
181 GEOMImpl_IMeasure aCI (aFunction);
182 Standard_Integer aType = aFunction->GetType();
186 if (aType == CDG_MEASURE)
188 Handle(GEOM_Function) aRefBase = aCI.GetBase();
189 TopoDS_Shape aShapeBase = aRefBase->GetValue();
190 if (aShapeBase.IsNull())
191 Standard_NullObject::Raise("Shape for centre of mass calculation is null");
193 gp_Ax3 aPos = GEOMUtils::GetPosition(aShapeBase);
194 gp_Pnt aCenterMass = aPos.Location();
195 aShape = BRepBuilderAPI_MakeVertex(aCenterMass).Shape();
197 else if (aType == BND_BOX_MEASURE || aType == BND_BOX_MEASURE_PRECISE)
199 Handle(GEOM_Function) aRefBase = aCI.GetBase();
200 TopoDS_Shape aShapeBase = aRefBase->GetValue();
201 if (aShapeBase.IsNull())
202 Standard_NullObject::Raise("Shape for bounding box calculation is null");
204 BRepBuilderAPI_Copy aCopyTool (aShapeBase);
205 if (!aCopyTool.IsDone())
206 Standard_NullObject::Raise("Shape for bounding box calculation is bad");
208 aShapeBase = aCopyTool.Shape();
210 // remove triangulation to obtain more exact boundaries
211 BRepTools::Clean(aShapeBase);
214 BRepBndLib::Add(aShapeBase, B);
216 if (aType == BND_BOX_MEASURE_PRECISE) {
217 if (!GEOMUtils::PreciseBoundingBox(aShapeBase, B)) {
218 Standard_NullObject::Raise("Bounding box cannot be precised");
222 Standard_Real Xmin, Ymin, Zmin, Xmax, Ymax, Zmax;
223 B.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
225 if (Xmax - Xmin < Precision::Confusion()) Xmax += Precision::Confusion();
226 if (Ymax - Ymin < Precision::Confusion()) Ymax += Precision::Confusion();
227 if (Zmax - Zmin < Precision::Confusion()) Zmax += Precision::Confusion();
229 gp_Pnt P1 (Xmin, Ymin, Zmin);
230 gp_Pnt P2 (Xmax, Ymax, Zmax);
232 BRepPrimAPI_MakeBox MB (P1, P2);
234 if (!MB.IsDone()) StdFail_NotDone::Raise("Bounding box cannot be computed from the given shape");
237 else if (aType == VERTEX_BY_INDEX)
239 Handle(GEOM_Function) aRefBase = aCI.GetBase();
240 TopoDS_Shape aShapeBase = aRefBase->GetValue();
241 if (aShapeBase.IsNull()) {
242 Standard_NullObject::Raise("Shape for centre of mass calculation is null");
245 int index = aCI.GetIndex();
246 bool useOri = aCI.GetUseOri();
249 if (aShapeBase.ShapeType() == TopAbs_VERTEX) {
251 Standard_NullObject::Raise("Vertex index is out of range");
253 aVertex = BRep_Tool::Pnt(TopoDS::Vertex(aShapeBase));
255 else if (aShapeBase.ShapeType() == TopAbs_EDGE) {
256 if (index < 0 || index > 1)
257 Standard_NullObject::Raise("Vertex index is out of range");
259 TopoDS_Edge anEdgeE = TopoDS::Edge(aShapeBase);
260 if ( anEdgeE.Orientation() != TopAbs_FORWARD &&
261 anEdgeE.Orientation() != TopAbs_REVERSED )
262 anEdgeE.Orientation( TopAbs_FORWARD );
265 TopExp::Vertices(anEdgeE, aV[0], aV[1], useOri );
266 aVertex = BRep_Tool::Pnt( aV[ index ]);
268 else if (aShapeBase.ShapeType() == TopAbs_WIRE) {
269 TopTools_IndexedMapOfShape anEdgeShapes;
270 TopTools_IndexedMapOfShape aVertexShapes;
271 TopoDS_Vertex aV1, aV2;
272 TopoDS_Wire aWire = TopoDS::Wire(aShapeBase);
273 TopExp_Explorer exp (aWire, TopAbs_EDGE);
274 for (; exp.More(); exp.Next()) {
275 anEdgeShapes.Add(exp.Current());
276 TopoDS_Edge E = TopoDS::Edge(exp.Current());
277 TopExp::Vertices(E, aV1, aV2);
278 if ( aVertexShapes.Extent() == 0)
279 aVertexShapes.Add(aV1);
280 if ( !aV1.IsSame( aVertexShapes(aVertexShapes.Extent()) ) )
281 aVertexShapes.Add(aV1);
282 if ( !aV2.IsSame( aVertexShapes(aVertexShapes.Extent()) ) )
283 aVertexShapes.Add(aV2);
286 if (index < 0 || index > aVertexShapes.Extent())
287 Standard_NullObject::Raise("Vertex index is out of range");
289 if ( aWire.Orientation() == TopAbs_FORWARD || !useOri )
290 aVertex = BRep_Tool::Pnt(TopoDS::Vertex(aVertexShapes(index+1)));
292 aVertex = BRep_Tool::Pnt(TopoDS::Vertex(aVertexShapes(aVertexShapes.Extent() - index)));
295 Standard_NullObject::Raise("Shape for vertex calculation is not an edge or wire");
298 aShape = BRepBuilderAPI_MakeVertex(aVertex).Shape();
300 else if (aType == VECTOR_FACE_NORMALE)
303 Handle(GEOM_Function) aRefBase = aCI.GetBase();
304 TopoDS_Shape aShapeBase = aRefBase->GetValue();
305 if (aShapeBase.IsNull()) {
306 Standard_NullObject::Raise("Face for normale calculation is null");
308 if (aShapeBase.ShapeType() != TopAbs_FACE) {
309 Standard_NullObject::Raise("Shape for normale calculation is not a face");
311 TopoDS_Face aFace = TopoDS::Face(aShapeBase);
316 Handle(GEOM_Function) aPntFunc = aCI.GetPoint();
317 if (!aPntFunc.IsNull())
319 TopoDS_Shape anOptPnt = aPntFunc->GetValue();
320 if (anOptPnt.IsNull())
321 Standard_NullObject::Raise("Invalid shape given for point argument");
322 p1 = BRep_Tool::Pnt(TopoDS::Vertex(anOptPnt));
326 gp_Ax3 aPos = GEOMUtils::GetPosition(aFace);
327 p1 = aPos.Location();
330 // Point parameters on surface
331 Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace);
332 Handle(ShapeAnalysis_Surface) aSurfAna = new ShapeAnalysis_Surface (aSurf);
333 gp_Pnt2d pUV = aSurfAna->ValueOfUV(p1, Precision::Confusion());
337 BRepAdaptor_Surface SF (aFace);
338 SF.D1(pUV.X(), pUV.Y(), p1, Vec1, Vec2);
339 if (Vec1.Magnitude() < Precision::Confusion()) {
342 SF.D1(pUV.X(), pUV.Y()-0.1, tmpP, Vec1, tmpV);
344 else if (Vec2.Magnitude() < Precision::Confusion()) {
347 SF.D1(pUV.X()-0.1, pUV.Y(), tmpP, tmpV, Vec2);
350 gp_Vec V = Vec1.Crossed(Vec2);
351 Standard_Real mod = V.Magnitude();
352 if (mod < Precision::Confusion())
353 Standard_NullObject::Raise("Normal vector of a face has null magnitude");
355 // Set length of normal vector to average radius of curvature
356 Standard_Real radius = 0.0;
357 GeomLProp_SLProps aProperties (aSurf, pUV.X(), pUV.Y(), 2, Precision::Confusion());
358 if (aProperties.IsCurvatureDefined()) {
359 Standard_Real radius1 = Abs(aProperties.MinCurvature());
360 Standard_Real radius2 = Abs(aProperties.MaxCurvature());
361 if (Abs(radius1) > Precision::Confusion()) {
362 radius = 1.0 / radius1;
363 if (Abs(radius2) > Precision::Confusion()) {
364 radius = (radius + 1.0 / radius2) / 2.0;
368 if (Abs(radius2) > Precision::Confusion()) {
369 radius = 1.0 / radius2;
374 // Set length of normal vector to average dimension of the face
375 // (only if average radius of curvature is not appropriate)
376 if (radius < Precision::Confusion()) {
378 Standard_Real Xmin, Xmax, Ymin, Ymax, Zmin, Zmax;
379 BRepBndLib::Add(aFace, B);
380 B.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
381 radius = ((Xmax - Xmin) + (Ymax - Ymin) + (Zmax - Zmin)) / 3.0;
384 if (radius < Precision::Confusion())
389 // consider the face orientation
390 if (aFace.Orientation() == TopAbs_REVERSED ||
391 aFace.Orientation() == TopAbs_INTERNAL) {
396 gp_Pnt p2 = p1.Translated(V);
397 BRepBuilderAPI_MakeEdge aBuilder (p1, p2);
398 if (!aBuilder.IsDone())
399 Standard_NullObject::Raise("Vector construction failed");
400 aShape = aBuilder.Shape();
402 else if (aType == CURVATURE_VEC_MEASURE) {
403 Handle(GEOM_Function) aSrfFunc = aCI.GetBase();
404 Handle(GEOM_Function) aPntFunc = aCI.GetPoint();
405 Handle(GEOM_Function) aDirFunc = aCI.GetDirection();
407 TopoDS_Shape aFace = aSrfFunc->GetValue();
408 TopoDS_Shape aVertex = aPntFunc->GetValue();
409 TopoDS_Shape anEdge = aDirFunc->GetValue();
411 aShape = EvaluateAlongCurvature(aFace, aVertex, anEdge);
416 if (aShape.IsNull()) return 0;
418 aFunction->SetValue(aShape);
420 log->SetTouched(Label());
425 //================================================================================
427 * \brief Returns a name of creation operation and names and values of creation parameters
429 //================================================================================
431 bool GEOMImpl_MeasureDriver::
432 GetCreationInformation(std::string& theOperationName,
433 std::vector<GEOM_Param>& theParams)
435 if (Label().IsNull()) return 0;
436 Handle(GEOM_Function) function = GEOM_Function::GetFunction(Label());
438 GEOMImpl_IMeasure aCI( function );
439 Standard_Integer aType = function->GetType();
443 theOperationName = "MASS_CENTER";
444 AddParam( theParams, "Object", aCI.GetBase() );
446 case BND_BOX_MEASURE:
447 case BND_BOX_MEASURE_PRECISE:
448 theOperationName = "BND_BOX";
449 AddParam( theParams, "Object", aCI.GetBase() );
451 case VERTEX_BY_INDEX:
452 theOperationName = "GetVertexByIndex";
453 AddParam( theParams, "Object", aCI.GetBase() );
454 AddParam( theParams, "Index", aCI.GetIndex() );
455 AddParam( theParams, "Oriented", aCI.GetUseOri() );
457 case VECTOR_FACE_NORMALE:
458 theOperationName = "NORMALE";
459 AddParam( theParams, "Face", aCI.GetBase() );
460 AddParam( theParams, "Point", aCI.GetPoint(), "face center" );
462 case CURVATURE_VEC_MEASURE:
463 theOperationName = "CURVATURE";
464 AddParam( theParams, "Face", aCI.GetBase() );
465 AddParam( theParams, "Point", aCI.GetPoint(), "point of interest" );
466 AddParam( theParams, "Vector", aCI.GetDirection(), "direction of interest" );
475 IMPLEMENT_STANDARD_RTTIEXT (GEOMImpl_MeasureDriver,GEOM_BaseDriver)