1 // Created on: 1992-06-30
2 // Created by: Laurent BUCHARD
3 // Copyright (c) 1992-1999 Matra Datavision
4 // Copyright (c) 1999-2012 OPEN CASCADE SAS
6 // The content of this file is subject to the Open CASCADE Technology Public
7 // License Version 6.5 (the "License"). You may not use the content of this file
8 // except in compliance with the License. Please obtain a copy of the License
9 // at http://www.opencascade.org and read it completely before using this file.
11 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
12 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
14 // The Original Code and all software distributed under the License is
15 // distributed on an "AS IS" basis, without warranty of any kind, and the
16 // Initial Developer hereby disclaims all such warranties, including without
17 // limitation, any warranties of merchantability, fitness for a particular
18 // purpose or non-infringement. Please see the License for the specific terms
19 // and conditions governing the rights and limitations under the License.
23 #define No_Standard_RangeError
24 #define No_Standard_OutOfRange
27 //----------------------------------------------------------------------
28 //-- Differents constructeurs sont proposes qui correspondent aux
30 //-- A(Sin(Theta),Cos(Theta)) Z**2
31 //-- + B(Sin(Theta),Cos(Theta)) Z
32 //-- + C(Sin(Theta),Cos(Theta))
34 //-- Une Courbe est definie sur un domaine
36 //-- Value retourne le point de parametre U(Theta),V(Theta)
37 //-- ou V est la solution du polynome A V**2 + B V + C
38 //-- (Selon les cas, on prend V+ ou V-)
40 //-- D1u calcule le vecteur tangent a la courbe
41 //-- et retourne le booleen Standard_False si ce calcul ne peut
42 //-- pas etre mene a bien.
43 //----------------------------------------------------------------------
45 #include <Precision.hxx>
47 #include <IntAna_Curve.ixx>
49 #include <Standard_DomainError.hxx>
50 #include <math_DirectPolynomialRoots.hxx>
54 //=======================================================================
55 //function : IntAna_Curve
57 //=======================================================================
58 IntAna_Curve::IntAna_Curve()
60 typequadric=GeomAbs_OtherSurface;
61 firstbounded=Standard_False;
62 lastbounded=Standard_False;
64 //=======================================================================
65 //function : SetConeQuadValues
66 //purpose : Description de l intersection Cone Quadrique
67 //=======================================================================
68 void IntAna_Curve::SetConeQuadValues(const gp_Cone& Cone,
69 const Standard_Real Qxx,
70 const Standard_Real Qyy,
71 const Standard_Real Qzz,
72 const Standard_Real Qxy,
73 const Standard_Real Qxz,
74 const Standard_Real Qyz,
75 const Standard_Real Qx,
76 const Standard_Real Qy,
77 const Standard_Real Qz,
78 const Standard_Real Q1,
79 const Standard_Real TOL,
80 const Standard_Real DomInf,
81 const Standard_Real DomSup,
82 const Standard_Boolean twocurves,
83 const Standard_Boolean takezpositive)
86 Ax3 = Cone.Position();
87 RCyl = Cone.RefRadius();
89 Angle = Cone.SemiAngle();
90 Standard_Real UnSurTgAngle = 1.0/(Tan(Cone.SemiAngle()));
92 typequadric= GeomAbs_Cone;
94 TwoCurves = twocurves; //-- deux Z pour un meme parametre
95 TakeZPositive = takezpositive; //-- Prendre sur la courbe le Z Positif
96 //-- ( -B + Sqrt()) et non (-B - Sqrt())
99 Z0Cte = Q1; //-- Attention On a Z?Cos Cos(t)
100 Z0Sin = 0.0; //-- et Non 2 Z?Cos Cos(t) !!!
101 Z0Cos = 0.0; //-- Ce pour tous les Parametres
102 Z0CosCos = 0.0; //-- ie pas de Coefficient 2
103 Z0SinSin = 0.0; //-- devant les termes CS C S
106 Z1Cte = 2.0*(UnSurTgAngle)*Qz;
113 Z2Cte = Qzz * UnSurTgAngle*UnSurTgAngle;
114 Z2Sin = (UnSurTgAngle+UnSurTgAngle)*Qyz;
115 Z2Cos = (UnSurTgAngle+UnSurTgAngle)*Qxz;
124 RestrictedInf = RestrictedSup = Standard_True; //-- Le Domaine est Borne
125 firstbounded = lastbounded = Standard_False;
128 //=======================================================================
129 //function : SetCylinderQuadValues
130 //purpose : Description de l intersection Cylindre Quadrique
131 //=======================================================================
132 void IntAna_Curve::SetCylinderQuadValues(const gp_Cylinder& Cyl,
133 const Standard_Real Qxx,
134 const Standard_Real Qyy,
135 const Standard_Real Qzz,
136 const Standard_Real Qxy,
137 const Standard_Real Qxz,
138 const Standard_Real Qyz,
139 const Standard_Real Qx,
140 const Standard_Real Qy,
141 const Standard_Real Qz,
142 const Standard_Real Q1,
143 const Standard_Real TOL,
144 const Standard_Real DomInf,
145 const Standard_Real DomSup,
146 const Standard_Boolean twocurves,
147 const Standard_Boolean takezpositive)
150 Ax3 = Cyl.Position();
152 typequadric= GeomAbs_Cylinder;
154 TwoCurves = twocurves; //-- deux Z pour un meme parametre
155 TakeZPositive = takezpositive; //-- Prendre sur la courbe le Z Positif
156 Standard_Real RCylmul2 = RCyl+RCyl; //-- ( -B + Sqrt())
161 Z0CosCos = Qxx*RCyl*RCyl;
162 Z0SinSin = Qyy*RCyl*RCyl;
163 Z0CosSin = RCylmul2*RCyl*Qxy;
166 Z1Sin = RCylmul2*Qyz;
167 Z1Cos = RCylmul2*Qxz;
183 RestrictedInf = RestrictedSup = Standard_True;
184 firstbounded = lastbounded = Standard_False;
187 //=======================================================================
190 //=======================================================================
191 Standard_Boolean IntAna_Curve::IsOpen() const
193 return(RestrictedInf && RestrictedSup);
196 //=======================================================================
199 //=======================================================================
200 void IntAna_Curve::Domain(Standard_Real& DInf,
201 Standard_Real& DSup) const
203 if(RestrictedInf && RestrictedSup) {
211 Standard_DomainError::Raise("IntAna_Curve::Domain");
214 //=======================================================================
215 //function : IsConstant
217 //=======================================================================
218 Standard_Boolean IntAna_Curve::IsConstant() const
220 //-- ??? Pas facile de decider a la seule vue des Param.
221 return(Standard_False);
224 //=======================================================================
225 //function : IsFirstOpen
227 //=======================================================================
228 Standard_Boolean IntAna_Curve::IsFirstOpen() const
230 return(firstbounded);
233 //=======================================================================
234 //function : IsLastOpen
236 //=======================================================================
237 Standard_Boolean IntAna_Curve::IsLastOpen() const
241 //=======================================================================
242 //function : SetIsFirstOpen
244 //=======================================================================
245 void IntAna_Curve::SetIsFirstOpen(const Standard_Boolean Flag)
250 //=======================================================================
251 //function : SetIsLastOpen
253 //=======================================================================
254 void IntAna_Curve::SetIsLastOpen(const Standard_Boolean Flag)
259 //=======================================================================
260 //function : InternalUVValue
262 //=======================================================================
263 void IntAna_Curve::InternalUVValue(const Standard_Real theta,
264 Standard_Real& Param1,
265 Standard_Real& Param2,
271 Standard_Real& SigneSqrtDis) const
273 Standard_Real Theta=theta;
274 Standard_Boolean SecondSolution=Standard_False;
276 if((Theta<DomainInf) ||
277 ((Theta>DomainSup) && (!TwoCurves)) ||
278 (Theta>(DomainSup+DomainSup-DomainInf+0.00000000000001))) {
279 Standard_DomainError::Raise("IntAna_Curve::Domain");
282 if(Theta>DomainSup) {
283 Theta=DomainSup+DomainSup-Theta;
284 SecondSolution=Standard_True;
290 SecondSolution=TakeZPositive;
295 Standard_Real costsint = cost*sint;
297 A=Z2Cte+sint*(Z2Sin+sint*Z2SinSin)+cost*(Z2Cos+cost*Z2CosCos)
300 B=Z1Cte+sint*(Z1Sin+sint*Z1SinSin)+cost*(Z1Cos+cost*Z1CosCos)
303 C=Z0Cte+sint*(Z0Sin+sint*Z0SinSin)+cost*(Z0Cos+cost*Z0CosCos)
307 Standard_Real Discriminant = B*B-4.0*A*C;
309 if(Abs(A)<=0.000000001) {
310 //-- cout<<" IntAna_Curve:: Internal UV Value : A="<<A<<" -> Abs(A)="<<Abs(A)<<endl;
311 if(Abs(B)<=0.0000000001) {
312 //-- cout<<" Probleme : Pas de solutions "<<endl;
316 //modified by NIZNHY-PKV Fri Dec 2 16:02:46 2005f
319 if(!SecondSolution) {
320 //-- Cas Param2 = (-B+Sqrt(Discriminant))/(A+A);
321 //-- = (-B+Sqrt(B**2 - Eps)) / 2A
326 //-- Cas Param2 = (-B-Sqrt(Discriminant))/(A+A);
327 //-- = (-B-Sqrt(B**2 - Eps)) / 2A
332 Param2 = -B*10000000.0;
336 //modified by NIZNHY-PKV Fri Dec 2 16:02:54 2005t
340 if(Discriminant<=0.0000000001 ||
341 Abs(Discriminant/(4*A))<=0.0000000001) Discriminant=0.0;
342 SigneSqrtDis = (SecondSolution)? Sqrt(Discriminant)
343 : -Sqrt(Discriminant);
344 Param2=(-B+SigneSqrtDis)/(A+A);
348 //=======================================================================
351 //=======================================================================
352 gp_Pnt IntAna_Curve::Value(const Standard_Real theta)
354 Standard_Real A, B, C, U, V, sint, cost, SigneSqrtDis;
360 InternalUVValue(theta,U,V,A,B,C,cost,sint,SigneSqrtDis);
361 //-- checked the parameter U and Raises Domain Error if Error
362 return(InternalValue(U,V));
364 //=======================================================================
367 //=======================================================================
368 Standard_Boolean IntAna_Curve::D1u(const Standard_Real theta,
372 //-- Pour detecter le cas ou le calcul est impossible
373 Standard_Real A, B, C, U, V, sint, cost, SigneSqrtDis;
378 InternalUVValue(theta,U,V,A,B,C,cost,sint,SigneSqrtDis);
381 if(Abs(SigneSqrtDis)<0.0000000001 || Abs(A)<0.0000001) return(Standard_False);
384 //-- Approximation de la derivee (mieux que le calcul mathematique!)
385 Standard_Real dtheta = (DomainSup-DomainInf)*0.000001;
386 Standard_Real theta2 = theta+dtheta;
387 if((theta2<DomainInf) || ((theta2>DomainSup) && (!TwoCurves))
388 || (theta2>(DomainSup+DomainSup-DomainInf+0.00000000000001))) {
390 theta2 = theta+dtheta;
392 gp_Pnt P2 = Value(theta2);
394 Vec.SetCoord((P2.X()-Pt.X())*dtheta,
395 (P2.Y()-Pt.Y())*dtheta,
396 (P2.Z()-Pt.Z())*dtheta);
398 return(Standard_True);
400 //=======================================================================
401 //function : FindParameter
402 //purpose : Para est en sortie le parametre sur la courbe
403 //=======================================================================
404 Standard_Boolean IntAna_Curve::FindParameter (const gp_Pnt& P,
405 Standard_Real& Para) const
407 Standard_Real theta,z, aTolPrecision=0.0001;
408 Standard_Real PIpPI = M_PI + M_PI;
410 switch (typequadric) {
412 case GeomAbs_Cylinder:
414 ElSLib::CylinderParameters(Ax3,RCyl,P,theta,z);
420 ElSLib::ConeParameters(Ax3,RCyl,Angle,P,theta,z);
428 Standard_Real epsAng = 1.e-8;
429 Standard_Real tmin = DomainInf;
430 Standard_Real tmax = DomainSup;
431 Standard_Real U,V,A,B,C,sint,cost,SigneSqrtDis;
438 //U=V=A=B=C=sint=cost=SigneSqrtDis=0.0;
440 if (!firstbounded && tmin > theta && (tmin-theta) <= epsAng) {
443 else if (!lastbounded && theta > tmax && (theta-tmax) <= epsAng) {
448 theta = theta + PIpPI;
450 else if (theta > tmax) {
451 theta = theta - PIpPI;
453 if (theta < tmin || theta > tmax) {
455 InternalUVValue(tmax,U,V,A,B,C,cost,sint,SigneSqrtDis);
456 gp_Pnt PMax(InternalValue(U,V));
457 if(PMax.Distance(P) < aTolPrecision) {
459 return(Standard_True);
463 InternalUVValue(tmin,U,V,A,B,C,cost,sint,SigneSqrtDis);
464 gp_Pnt PMin(InternalValue(U,V));
465 if(PMin.Distance(P) < aTolPrecision) {
467 return(Standard_True);
470 //-- lbr le 14 Fev 96 : On teste malgre tout si le point n est pas le
471 //-- point de debut ou de fin
472 //-- cout<<"False 1 "<<endl;
473 // theta = tmin; le 25 Nov 96
481 InternalUVValue(theta,U,z1,A,B,C,cost,sint,SigneSqrtDis);
482 A = B = C = sint = cost = SigneSqrtDis = 0.0;
483 InternalUVValue(tmax+tmax - theta,U,z2,A,B,C,cost,sint,SigneSqrtDis);
485 if (Abs(z-z1) <= Abs(z-z2)) {
489 Para = tmax+tmax - theta;
496 if((Para<DomainInf) || ((Para>DomainSup) && (!TwoCurves))
497 || (Para>(DomainSup+DomainSup-DomainInf+0.00000000000001))) {
498 return(Standard_False);
501 InternalUVValue(Para,U,V,A,B,C,cost,sint,SigneSqrtDis);
502 gp_Pnt PPara = InternalValue(U,V);
503 Standard_Real Dist = PPara.Distance(P);
504 if(Dist > aTolPrecision) {
505 //-- Il y a eu un probleme
506 //-- On teste si le point est un point double
507 InternalUVValue(tmin,U,V,A,B,C,cost,sint,SigneSqrtDis);
508 PPara = InternalValue(U,V);
509 Dist = PPara.Distance(P);
510 if(Dist <= aTolPrecision) {
512 return(Standard_True);
515 InternalUVValue(tmax,U,V,A,B,C,cost,sint,SigneSqrtDis);
516 PPara = InternalValue(U,V);
517 Dist = PPara.Distance(P);
518 if(Dist <= aTolPrecision) {
520 return(Standard_True);
523 Standard_Real Theta = DomainSup+DomainSup-DomainInf;
524 InternalUVValue(Theta,U,V,A,B,C,cost,sint,SigneSqrtDis);
525 PPara = InternalValue(U,V);
526 Dist = PPara.Distance(P);
527 if(Dist <= aTolPrecision) {
529 return(Standard_True);
532 return(Standard_False);
534 return(Standard_True);
536 //=======================================================================
537 //function : InternalValue
539 //=======================================================================
540 gp_Pnt IntAna_Curve::InternalValue(const Standard_Real U,
541 const Standard_Real _V) const
543 //-- cout<<" ["<<U<<","<<V<<"]";
544 Standard_Real V = _V;
545 if(V > 100000.0 ) { V= 100000.0; }
546 if(V < -100000.0 ) { V=-100000.0; }
548 switch(typequadric) {
551 //------------------------------------------------
552 //-- Parametrage : X = V * Cos(U) ---
553 //-- Y = V * Sin(U) ---
554 //-- Z = (V-RCyl) / Tan(SemiAngle)--
555 //------------------------------------------------
556 //-- Angle Vaut Cone.SemiAngle()
557 return(ElSLib::ConeValue(U,(V-RCyl)/Sin(Angle),Ax3,RCyl,Angle));
561 case GeomAbs_Cylinder:
562 return(ElSLib::CylinderValue(U,V,Ax3,RCyl));
564 return(ElSLib::SphereValue(U,V,Ax3,RCyl));
566 return(gp_Pnt(0.0,0.0,0.0));
568 return(gp_Pnt(0.0,0.0,0.0));
572 //=======================================================================
573 //function : SetDomain
575 //=======================================================================
576 void IntAna_Curve::SetDomain(const Standard_Real DInf,
577 const Standard_Real DSup)
580 Standard_DomainError::Raise("IntAna_Curve::Domain");