1 // Copyright (C) 2014-2022 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "GeomAlgoAPI_BoundingBox.h"
22 #include <BRepBuilderAPI_Copy.hxx>
23 #include <Bnd_Box.hxx>
24 #include <BRepTools.hxx>
25 #include <BRep_Tool.hxx>
26 #include <BRepBndLib.hxx>
27 #include <BRep_Builder.hxx>
28 #include <Geom_Circle.hxx>
29 #include <ShapeAnalysis.hxx>
30 #include <TopoDS_Shape.hxx>
33 #include <BRepBuilderAPI_MakeFace.hxx>
34 #include <Geom_RectangularTrimmedSurface.hxx>
35 #include <BRepExtrema_DistShapeShape.hxx>
36 #include <ShapeFix_Shape.hxx>
37 #include <BRepBuilderAPI_Sewing.hxx>
38 #include <Standard_ErrorHandler.hxx>
39 #include <TopExp_Explorer.hxx>
40 #include <BRepClass3d_SolidClassifier.hxx>
41 #include <Geom_SphericalSurface.hxx>
42 #include <Geom_ToroidalSurface.hxx>
45 * This function constructs and returns modified shape from the original one
46 * for singular cases. It is used for the method GetMinDistanceSingular.
48 * \param theShape the original shape
49 * \param theModifiedShape output parameter. The modified shape.
50 * \param theAddDist output parameter. The added distance for modified shape.
51 * \retval true if the shape is modified; false otherwise.
55 Standard_Boolean ModifyShape(const TopoDS_Shape &theShape,
56 TopoDS_Shape &theModifiedShape,
57 Standard_Real &theAddDist)
59 TopExp_Explorer anExp;
63 theModifiedShape.Nullify();
65 for ( anExp.Init( theShape, TopAbs_FACE ); anExp.More(); anExp.Next() ) {
67 theModifiedShape = anExp.Current();
70 TopoDS_Shape sh = theShape;
71 while(sh.ShapeType()==TopAbs_COMPOUND) {
72 TopoDS_Iterator it(sh);
75 Handle(Geom_Surface) S = BRep_Tool::Surface(TopoDS::Face(theModifiedShape));
76 if(S->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) ||
77 S->IsKind(STANDARD_TYPE(Geom_ToroidalSurface)) ||
79 const Standard_Boolean isShell =
80 (sh.ShapeType()==TopAbs_SHELL || sh.ShapeType()==TopAbs_FACE);
82 if (!isShell && S->IsKind(STANDARD_TYPE(Geom_SphericalSurface))) {
83 Handle(Geom_SphericalSurface) SS = Handle(Geom_SphericalSurface)::DownCast(S);
84 gp_Pnt PC = SS->Location();
87 B.MakeVertex(V,PC,1.e-7);
89 theAddDist = SS->Radius();
92 if (!isShell && S->IsKind(STANDARD_TYPE(Geom_ToroidalSurface))) {
93 Handle(Geom_ToroidalSurface) TS = Handle(Geom_ToroidalSurface)::DownCast(S);
94 gp_Ax3 ax3 = TS->Position();
95 Handle(Geom_Circle) C = new Geom_Circle(ax3.Ax2(),TS->MajorRadius());
98 B.MakeEdge(E,C,1.e-7);
100 theAddDist = TS->MinorRadius();
101 return Standard_True;
104 // non solid case or any periodic surface (Mantis 22454).
106 // changes for 0020677: EDF 1219 GEOM: MinDistance gives 0 instead of 20.88
107 //S->Bounds(U1,U2,V1,V2); changed by
108 ShapeAnalysis::GetFaceUVBounds(TopoDS::Face(theModifiedShape),U1,U2,V1,V2);
109 // end of changes for 020677 (dmv)
110 Handle(Geom_RectangularTrimmedSurface) TrS1 =
111 new Geom_RectangularTrimmedSurface(S,U1,(U1+U2)/2.,V1,V2);
112 Handle(Geom_RectangularTrimmedSurface) TrS2 =
113 new Geom_RectangularTrimmedSurface(S,(U1+U2)/2.,U2,V1,V2);
114 TopoDS_Shape aMShape;
116 TopoDS_Face F1 = BRepBuilderAPI_MakeFace(TrS1, Precision::Confusion());
117 TopoDS_Face F2 = BRepBuilderAPI_MakeFace(TrS2, Precision::Confusion());
121 B.MakeCompound(TopoDS::Compound(aMShape));
125 // The original shape is a solid.
126 BRepBuilderAPI_Sewing aSewing (Precision::Confusion()*10.0);
130 aMShape = aSewing.SewedShape();
134 B.Add(aSolid, aMShape);
138 Handle(ShapeFix_Shape) sfs = new ShapeFix_Shape;
140 sfs->SetPrecision(1.e-6);
141 sfs->SetMaxTolerance(1.0);
143 theModifiedShape = sfs->Shape();
144 return Standard_True;
148 theModifiedShape = theShape;
149 return Standard_False;
152 //=======================================================================
153 // function : GetMinDistanceSingular
154 //=======================================================================
155 double GetMinDistanceSingular(const TopoDS_Shape& aSh1,
156 const TopoDS_Shape& aSh2,
157 gp_Pnt& Ptmp1, gp_Pnt& Ptmp2)
161 Standard_Real AddDist1 = 0.;
162 Standard_Real AddDist2 = 0.;
163 Standard_Boolean IsChange1 = ModifyShape(aSh1, tmpSh1, AddDist1);
164 Standard_Boolean IsChange2 = ModifyShape(aSh2, tmpSh2, AddDist2);
166 if( !IsChange1 && !IsChange2 )
169 BRepExtrema_DistShapeShape dst(tmpSh1,tmpSh2);
171 double MinDist = 1.e9;
172 gp_Pnt PMin1, PMin2, P1, P2;
173 for (int i = 1; i <= dst.NbSolution(); i++) {
174 P1 = dst.PointOnShape1(i);
175 P2 = dst.PointOnShape2(i);
176 Standard_Real Dist = P1.Distance(P2);
177 if (MinDist > Dist) {
188 gp_Dir aDir(gp_Vec(PMin1,PMin2));
189 if( MinDist > (AddDist1+AddDist2) ) {
190 Ptmp1 = gp_Pnt(PMin1.X() + aDir.X()*AddDist1,
191 PMin1.Y() + aDir.Y()*AddDist1,
192 PMin1.Z() + aDir.Z()*AddDist1);
193 Ptmp2 = gp_Pnt(PMin2.X() - aDir.X()*AddDist2,
194 PMin2.Y() - aDir.Y()*AddDist2,
195 PMin2.Z() - aDir.Z()*AddDist2);
196 return (MinDist - AddDist1 - AddDist2);
200 Ptmp1 = gp_Pnt(PMin1.X() + aDir.X()*AddDist1,
201 PMin1.Y() + aDir.Y()*AddDist1,
202 PMin1.Z() + aDir.Z()*AddDist1);
206 Ptmp2 = gp_Pnt(PMin2.X() - aDir.X()*AddDist2,
207 PMin2.Y() - aDir.Y()*AddDist2,
208 PMin2.Z() - aDir.Z()*AddDist2);
213 double res = MinDist - AddDist1 - AddDist2;
214 if(res<0.) res = 0.0;
220 //=======================================================================
221 // function : GetMinDistance
222 //=======================================================================
223 Standard_Real GetMinDistance(const TopoDS_Shape& theShape1,
224 const TopoDS_Shape& theShape2,
225 gp_Pnt& thePnt1, gp_Pnt& thePnt2)
227 Standard_Real aResult = 1.e9;
229 // Issue 0020231: A min distance bug with torus and vertex.
230 // Make GetMinDistance() return zero if a sole VERTEX is inside any of SOLIDs
232 // which of shapes consists of only one vertex?
233 TopExp_Explorer exp1(theShape1,TopAbs_VERTEX), exp2(theShape2,TopAbs_VERTEX);
234 TopoDS_Shape V1 = exp1.More() ? exp1.Current() : TopoDS_Shape();
235 TopoDS_Shape V2 = exp2.More() ? exp2.Current() : TopoDS_Shape();
236 exp1.Next(); exp2.Next();
237 if ( exp1.More() ) V1.Nullify();
238 if ( exp2.More() ) V2.Nullify();
239 // vertex and container of solids
240 TopoDS_Shape V = V1.IsNull() ? V2 : V1;
241 TopoDS_Shape S = V1.IsNull() ? theShape1 : theShape2;
243 // classify vertex against solids
244 gp_Pnt p = BRep_Tool::Pnt( TopoDS::Vertex( V ) );
245 for ( exp1.Init( S, TopAbs_SOLID ); exp1.More(); exp1.Next() ) {
246 BRepClass3d_SolidClassifier classifier( exp1.Current(), p, 1e-6);
247 if ( classifier.State() == TopAbs_IN ) {
255 aResult = GetMinDistanceSingular(theShape1, theShape2, thePnt1, thePnt2);
258 BRepExtrema_DistShapeShape dst (theShape1, theShape2);
262 for (int i = 1; i <= dst.NbSolution(); i++) {
263 P1 = dst.PointOnShape1(i);
264 P2 = dst.PointOnShape2(i);
266 Standard_Real Dist = P1.Distance(P2);
267 if (aResult < 0 || aResult > Dist) {
278 //=======================================================================
279 // function : PreciseBoundingBox
280 //=======================================================================
281 Standard_Boolean PreciseBoundingBox(const TopoDS_Shape &theShape, Bnd_Box &theBox)
283 if (theBox.IsVoid()) BRepBndLib::Add( theShape, theBox );
284 if (theBox.IsVoid()) return Standard_False;
286 Standard_Real aBound[6];
287 theBox.Get(aBound[0], aBound[2], aBound[4], aBound[1], aBound[3], aBound[5]);
290 const gp_Pnt aMid(0.5*(aBound[1] + aBound[0]), // XMid
291 0.5*(aBound[3] + aBound[2]), // YMid
292 0.5*(aBound[5] + aBound[4])); // ZMid
293 const gp_XYZ aSize(aBound[1] - aBound[0], // DX
294 aBound[3] - aBound[2], // DY
295 aBound[5] - aBound[4]); // DZ
296 const gp_Pnt aPnt[6] =
298 gp_Pnt(aBound[0] - (aBound[1] - aBound[0]), aMid.Y(), aMid.Z()), // XMin
299 gp_Pnt(aBound[1] + (aBound[1] - aBound[0]), aMid.Y(), aMid.Z()), // XMax
300 gp_Pnt(aMid.X(), aBound[2] - (aBound[3] - aBound[2]), aMid.Z()), // YMin
301 gp_Pnt(aMid.X(), aBound[3] + (aBound[3] - aBound[2]), aMid.Z()), // YMax
302 gp_Pnt(aMid.X(), aMid.Y(), aBound[4] - (aBound[5] - aBound[4])), // ZMin
303 gp_Pnt(aMid.X(), aMid.Y(), aBound[5] + (aBound[5] - aBound[4])) // ZMax
305 const gp_Dir aDir[3] = { gp::DX(), gp::DY(), gp::DZ() };
306 const Standard_Real aPlnSize[3] =
308 0.5*Max(aSize.Y(), aSize.Z()), // XMin, XMax planes
309 0.5*Max(aSize.X(), aSize.Z()), // YMin, YMax planes
310 0.5*Max(aSize.X(), aSize.Y()) // ZMin, ZMax planes
314 for (i = 0; i < 6; i++) {
315 const Standard_Integer iHalf = i/2;
316 const gp_Pln aPln(aPnt[i], aDir[iHalf]);
317 BRepBuilderAPI_MakeFace aMkFace(aPln, -aPlnSize[iHalf], aPlnSize[iHalf],
318 -aPlnSize[iHalf], aPlnSize[iHalf]);
320 if (!aMkFace.IsDone()) {
321 return Standard_False;
324 TopoDS_Shape aFace = aMkFace.Shape();
326 // Get minimal distance between planar face and shape.
327 Standard_Real aMinDist = GetMinDistance(aFace, theShape, aPMin[0], aPMin[1]);
330 return Standard_False;
333 aBound[i] = aPMin[1].Coord(iHalf + 1);
336 // Update Bounding box with the new values.
338 theBox.Update(aBound[0], aBound[2], aBound[4], aBound[1], aBound[3], aBound[5]);
340 return Standard_True;
343 //=================================================================================================
344 bool GetBoundingBox(const std::shared_ptr<GeomAPI_Shape>& theShape,
345 double& theXmin, double& theXmax,
346 double& theYmin, double& theYmax,
347 double& theZmin, double& theZmax,
348 std::string& theError)
351 std::cout << "GetBoundingBox " << std::endl;
354 if (!theShape.get()) {
355 theError = "GetBoundingBox : An invalid argument";
359 TopoDS_Shape aShape = theShape->impl<TopoDS_Shape>();
361 //Compute the parameters
365 BRepBuilderAPI_Copy aCopyTool (aShape);
366 if (!aCopyTool.IsDone()) {
367 theError = "GetBoundingBox Error: Bad shape detected";
371 aShape = aCopyTool.Shape();
373 // remove triangulation to obtain more exact boundaries
374 BRepTools::Clean(aShape);
376 BRepBndLib::Add(aShape, B);
378 if (!PreciseBoundingBox(aShape, B)) {
379 theError = "GetBoundingBox Error: Bounding box cannot be precised";
383 B.Get(theXmin, theYmin, theZmin, theXmax, theYmax, theZmax);
385 catch (Standard_Failure& aFail) {
386 theError = aFail.GetMessageString();