1 // Copyright (c) 1995-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
4 // This file is part of Open CASCADE Technology software library.
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
17 #include <Extrema_ExtPElS.hxx>
18 #include <Extrema_POnSurf.hxx>
19 #include <gp_Cone.hxx>
20 #include <gp_Cylinder.hxx>
23 #include <gp_Sphere.hxx>
24 #include <gp_Torus.hxx>
25 #include <Standard_NotImplemented.hxx>
26 #include <Standard_OutOfRange.hxx>
27 #include <StdFail_NotDone.hxx>
29 static const Standard_Real ExtPElS_MyEps = Epsilon(2. * M_PI);
30 //=============================================================================
32 Extrema_ExtPElS::Extrema_ExtPElS()
34 myDone = Standard_False;
36 for (Standard_Integer i = 0; i < 4; i++)
38 mySqDist[i] = RealLast();
41 //=============================================================================
43 Extrema_ExtPElS::Extrema_ExtPElS (const gp_Pnt& P,
45 const Standard_Real Tol)
50 /*-----------------------------------------------------------------------------
52 Find 2 extreme distances between point P and cylinder S.
55 Let Pp be the projection of P in plane XOY of the cylinder;
56 2 cases are considered:
57 1- distance(Pp,O) < Tol:
58 There are infinite solutions; IsDone() = Standard_False.
59 2- distance(Pp,O) > Tol:
61 U1 = angle(OX,OPp) with 0 < U1 < 2.*M_PI
62 U2 = U1 + M_PI with 0 < U2 < 2.*M_PI;
63 then (U1,V) corresponds to the min distance.
64 and (U2,V) corresponds to the max distance.
65 -----------------------------------------------------------------------------*/
67 void Extrema_ExtPElS::Perform(const gp_Pnt& P,
69 const Standard_Real Tol)
71 myDone = Standard_False;
74 // Projection of point P in plane XOY of the cylinder ...
75 gp_Ax3 Pos = S.Position();
76 gp_Pnt O = Pos.Location();
77 gp_Vec OZ (Pos.Direction());
78 Standard_Real V = gp_Vec(O,P).Dot(OZ);
79 gp_Pnt Pp = P.Translated(OZ.Multiplied(-V));
81 // Calculation of extrema
83 if (OPp.Magnitude() < Tol) { return; }
84 gp_Vec myZ = Pos.XDirection()^Pos.YDirection();
85 Standard_Real U1 = gp_Vec(Pos.XDirection()).AngleWithRef(OPp,myZ); //-M_PI<U1<M_PI
86 if (U1 > -ExtPElS_MyEps && U1 < ExtPElS_MyEps) { U1 = 0.; }
87 Standard_Real U2 = U1 + M_PI;
88 if (U1 < 0.) { U1 += 2. * M_PI; }
91 Ps = ElSLib::Value(U1,V,S);
92 mySqDist[0] = Ps.SquareDistance(P);
93 myPoint[0] = Extrema_POnSurf(U1,V,Ps);
94 Ps = ElSLib::Value(U2,V,S);
95 mySqDist[1] = Ps.SquareDistance(P);
96 myPoint[1] = Extrema_POnSurf(U2,V,Ps);
99 myDone = Standard_True;
101 //=============================================================================
103 Extrema_ExtPElS::Extrema_ExtPElS (const gp_Pnt& P,
105 const Standard_Real Tol)
109 /*-----------------------------------------------------------------------------
111 Find 2 extreme distances between point P and cone S.
114 Let M the top of the cone.
115 2 cases are considered:
116 1- distance(P,M) < Tol:
117 there is a minimum in M.
118 2- distance(P,M) > Tol:
119 Let Pp the projection of P in the plane XOY of the cone;
120 2 cases are considered:
121 1- distance(Pp,O) < Tol:
122 There is an infinite number of solutions; IsDone() = Standard_False.
123 2- distance(Pp,O) > Tol:
124 There exist 2 extrema:
125 Let Vm = value of v for point M,
126 Vp = value of v for point P,
127 U1 = angle(OX,OPp) if Vp > Vm )
128 -angle(OX,OPp) otherwise ) with 0. < U1 < 2*M_PI,
129 U2 = U1 + M_PI with 0. < U2 < 2*M_PI;
130 We are in plane PpOZ.
131 Let A the angle of the cone,
132 B = angle(MP,MO) with 0. < B < M_PI,
134 V1 = (L * cos(B-A)) + Vm,
135 V2 = (L * cos(B+A)) + Vm;
136 then (U1,V1) and (U2,V2) correspond to min distances.
137 -----------------------------------------------------------------------------*/
139 void Extrema_ExtPElS::Perform(const gp_Pnt& P,
141 const Standard_Real Tol)
143 myDone = Standard_False;
147 gp_Ax3 Pos = S.Position();
148 gp_Pnt O = Pos.Location();
149 Standard_Real A = S.SemiAngle();
150 gp_Vec OZ (Pos.Direction());
151 gp_Vec myZ = Pos.XDirection()^Pos.YDirection();
154 Standard_Real L2 = MP.SquareMagnitude();
155 Standard_Real Vm = -(S.RefRadius() / Sin(A));
157 // Case when P is mixed with S ...
158 if (L2 < Tol * Tol) {
160 myPoint[0] = Extrema_POnSurf(0.,Vm,M);
162 myDone = Standard_True;
166 if (M.SquareDistance(O) < Tol * Tol)
168 DirZ = (A < 0 ? -OZ : OZ);
173 // Projection of P in the reference plane of the cone ...
174 Standard_Real Zp = gp_Vec(O, P).Dot(OZ);
176 gp_Pnt Pp = P.Translated(OZ.Multiplied(-Zp));
178 if (OPp.SquareMagnitude() < Tol * Tol) return;
179 Standard_Real B, U1, V1, U2, V2;
180 Standard_Boolean Same = DirZ.Dot(MP) >= 0.0;
181 U1 = gp_Vec(Pos.XDirection()).AngleWithRef(OPp,myZ); //-M_PI<U1<M_PI
182 if (U1 > -ExtPElS_MyEps && U1 < ExtPElS_MyEps) { U1 = 0.; }
184 if (!Same) { U1 += M_PI; }
186 if (U1 < 0.) { U1 += 2. * M_PI; }
187 if (U2 > 2.*M_PI) { U2 -= 2. * M_PI; }
190 Standard_Real L = sqrt(L2);
200 Standard_Real Sense = OZ.Dot(gp_Dir(DirZ));
201 V1 *= Sense; V2 *= Sense;
205 Ps = ElSLib::Value(U1,V1,S);
206 mySqDist[0] = Ps.SquareDistance(P);
207 myPoint[0] = Extrema_POnSurf(U1,V1,Ps);
208 Ps = ElSLib::Value(U2,V2,S);
209 mySqDist[1] = Ps.SquareDistance(P);
210 myPoint[1] = Extrema_POnSurf(U2,V2,Ps);
213 myDone = Standard_True;
215 //=============================================================================
217 Extrema_ExtPElS::Extrema_ExtPElS (const gp_Pnt& P,
219 const Standard_Real Tol)
223 /*-----------------------------------------------------------------------------
225 Find 2 extreme distances between point P and sphere S.
228 Let O be the origin of the sphere.
229 2 cases are considered:
230 1- distance(P,O) < Tol:
231 There is an infinite number of solutions; IsDone() = Standard_False
232 2- distance(P,O) > Tol:
233 Let Pp be the projection of point P in the plane XOY of the sphere;
234 2 cases are considered:
235 1- distance(Pp,O) < Tol:
236 2 solutions are: (0,-M_PI/2.) and (0.,M_PI/2.)
237 2- distance(Pp,O) > Tol:
238 Let U1 = angle(OX,OPp) with 0. < U1 < 2.*M_PI,
239 U2 = U1 + M_PI avec 0. < U2 < 2*M_PI,
240 V1 = angle(OPp,OP) with -M_PI/2. < V1 < M_PI/2. ,
241 then (U1, V1) corresponds to the min distance
242 and (U2,-V1) corresponds to the max distance.
243 -----------------------------------------------------------------------------*/
245 void Extrema_ExtPElS::Perform(const gp_Pnt& P,
247 const Standard_Real Tol)
249 myDone = Standard_False;
252 gp_Ax3 Pos = S.Position();
253 gp_Vec OP (Pos.Location(),P);
255 // Case when P is mixed with O ...
256 if (OP.SquareMagnitude() < Tol * Tol) { return; }
258 // Projection if P in plane XOY of the sphere ...
259 gp_Pnt O = Pos.Location();
260 gp_Vec OZ (Pos.Direction());
261 Standard_Real Zp = OP.Dot(OZ);
262 gp_Pnt Pp = P.Translated(OZ.Multiplied(-Zp));
264 // Calculation of extrema ...
266 Standard_Real U1, U2, V;
267 if (OPp.SquareMagnitude() < Tol * Tol) {
270 if (Zp < 0.) { V = -M_PI / 2.; }
271 else { V = M_PI / 2.; }
274 gp_Vec myZ = Pos.XDirection()^Pos.YDirection();
275 U1 = gp_Vec(Pos.XDirection()).AngleWithRef(OPp,myZ);
276 if (U1 > -ExtPElS_MyEps && U1 < ExtPElS_MyEps) { U1 = 0.; }
278 if (U1 < 0.) { U1 += 2. * M_PI; }
280 if (Zp < 0.) { V = -V; }
284 Ps = ElSLib::Value(U1,V,S);
285 mySqDist[0] = Ps.SquareDistance(P);
286 myPoint[0] = Extrema_POnSurf(U1,V,Ps);
287 Ps = ElSLib::Value(U2,-V,S);
288 mySqDist[1] = Ps.SquareDistance(P);
289 myPoint[1] = Extrema_POnSurf(U2,-V,Ps);
292 myDone = Standard_True;
294 //=============================================================================
296 Extrema_ExtPElS::Extrema_ExtPElS (const gp_Pnt& P,
298 const Standard_Real Tol)
302 /*-----------------------------------------------------------------------------
304 Find 2 extreme distances between point P and torus S.
307 Let Pp be the projection of point P in plane XOY of the torus;
308 2 cases are consideres:
309 1- distance(Pp,O) < Tol:
310 There is an infinite number of solutions; IsDone() = Standard_False.
311 2- distance(Pp,O) > Tol:
312 One is located in plane PpOZ;
313 Let V1 = angle(OX,OPp) with 0. < V1 < 2.*M_PI,
314 V2 = V1 + M_PI with 0. < V2 < 2.*M_PI,
315 O1 and O2 centers of circles (O1 on coord. posit.)
318 then (U1,V1) corresponds to the min distance
319 and (U2,V2) corresponds to the max distance.
320 -----------------------------------------------------------------------------*/
321 void Extrema_ExtPElS::Perform(const gp_Pnt& P,
323 const Standard_Real Tol)
325 const Standard_Real tol2 = Tol*Tol;
326 myDone = Standard_False;
329 // Projection of P in plane XOY ...
330 gp_Ax3 Pos = S.Position();
331 gp_Pnt O = Pos.Location();
332 gp_Vec OZ (Pos.Direction());
333 gp_Pnt Pp = P.Translated(OZ.Multiplied(-(gp_Vec(O,P).Dot(Pos.Direction()))));
335 // Calculation of extrema ...
337 Standard_Real R2 = OPp.SquareMagnitude();
338 if (R2 < tol2) { return; }
340 gp_Vec myZ = Pos.XDirection()^Pos.YDirection();
341 Standard_Real U1 = gp_Vec(Pos.XDirection()).AngleWithRef(OPp,myZ);
342 if (U1 > -ExtPElS_MyEps && U1 < ExtPElS_MyEps) { U1 = 0.; }
343 Standard_Real U2 = U1 + M_PI;
344 if (U1 < 0.) { U1 += 2. * M_PI; }
345 Standard_Real R = sqrt(R2);
346 gp_Vec OO1 = OPp.Divided(R).Multiplied(S.MajorRadius());
347 gp_Vec OO2 = OO1.Multiplied(-1.);
348 gp_Pnt O1 = O.Translated(OO1);
349 gp_Pnt O2 = O.Translated(OO2);
351 if(O1.SquareDistance(P) < tol2) { return; }
352 if(O2.SquareDistance(P) < tol2) { return; }
354 Standard_Real V1 = OPp.AngleWithRef(gp_Vec(O1,P),OPp.Crossed(OZ));
355 if (V1 > -ExtPElS_MyEps && V1 < ExtPElS_MyEps) { V1 = 0.; }
357 Standard_Real V2 = OPp.AngleWithRef(gp_Vec(P,O2),OPp.Crossed(OZ));
358 if (V2 > -ExtPElS_MyEps && V2 < ExtPElS_MyEps) { V2 = 0.; }
360 if (V1 < 0.) { V1 += 2. * M_PI; }
361 if (V2 < 0.) { V2 += 2. * M_PI; }
364 Ps = ElSLib::Value(U1,V1,S);
365 mySqDist[0] = Ps.SquareDistance(P);
366 myPoint[0] = Extrema_POnSurf(U1,V1,Ps);
368 Ps = ElSLib::Value(U1,V1+M_PI,S);
369 mySqDist[1] = Ps.SquareDistance(P);
370 myPoint[1] = Extrema_POnSurf(U1,V1+M_PI,Ps);
372 Ps = ElSLib::Value(U2,V2,S);
373 mySqDist[2] = Ps.SquareDistance(P);
374 myPoint[2] = Extrema_POnSurf(U2,V2,Ps);
376 Ps = ElSLib::Value(U2,V2+M_PI,S);
377 mySqDist[3] = Ps.SquareDistance(P);
378 myPoint[3] = Extrema_POnSurf(U2,V2+M_PI,Ps);
381 myDone = Standard_True;
385 Extrema_ExtPElS::Extrema_ExtPElS (const gp_Pnt& P,
387 const Standard_Real Tol)
392 void Extrema_ExtPElS::Perform (const gp_Pnt& P,
394 // const Standard_Real Tol)
395 const Standard_Real )
397 myDone = Standard_False;
400 // Projection of point P in plane XOY of the cylinder ...
401 gp_Pnt O = S.Location();
402 gp_Vec OZ (S.Axis().Direction());
403 Standard_Real U, V = gp_Vec(O,P).Dot(OZ);
404 gp_Pnt Pp = P.Translated(OZ.Multiplied(-V));
406 ElSLib::Parameters(S, P, U, V);
407 mySqDist[0] = Pp.SquareDistance(P);
408 myPoint[0] = Extrema_POnSurf(U,V,Pp);
410 myDone = Standard_True;
414 //=============================================================================
416 Standard_Boolean Extrema_ExtPElS::IsDone () const { return myDone; }
417 //=============================================================================
419 Standard_Integer Extrema_ExtPElS::NbExt () const
421 if (!IsDone()) { throw StdFail_NotDone(); }
424 //=============================================================================
426 Standard_Real Extrema_ExtPElS::SquareDistance (const Standard_Integer N) const
428 if ((N < 1) || (N > NbExt()))
430 throw Standard_OutOfRange();
432 return mySqDist[N-1];
434 //=============================================================================
436 const Extrema_POnSurf& Extrema_ExtPElS::Point (const Standard_Integer N) const
438 if ((N < 1) || (N > NbExt()))
440 throw Standard_OutOfRange();
444 //=============================================================================