projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0025418: Debug output to be limited to OCC development environment
[occt.git] / src / IntAna / IntAna_Curve.cxx
CommitLineData
b311480e 1// Created on: 1992-06-30
2// Created by: Laurent BUCHARD
3// Copyright (c) 1992-1999 Matra Datavision
973c2be1 4// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 5//
973c2be1 6// This file is part of Open CASCADE Technology software library.
b311480e 7//
d5f74e42 8// This library is free software; you can redistribute it and/or modify it under
9// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 10// by the Free Software Foundation, with special exception defined in the file
11// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12// distribution for complete text of the license and disclaimer of any warranty.
b311480e 13//
973c2be1 14// Alternatively, this file may be used under the terms of Open CASCADE
15// commercial license or contractual agreement.
7fd59977 16
0797d9d3 17#ifndef OCCT_DEBUG
7fd59977 18#define No_Standard_RangeError
19#define No_Standard_OutOfRange
20#endif
21
22//----------------------------------------------------------------------
23//-- Differents constructeurs sont proposes qui correspondent aux
24//-- polynomes en Z :
25//-- A(Sin(Theta),Cos(Theta)) Z**2
26//-- + B(Sin(Theta),Cos(Theta)) Z
27//-- + C(Sin(Theta),Cos(Theta))
28//--
29//-- Une Courbe est definie sur un domaine
30//--
31//-- Value retourne le point de parametre U(Theta),V(Theta)
32//-- ou V est la solution du polynome A V**2 + B V + C
33//-- (Selon les cas, on prend V+ ou V-)
34//--
35//-- D1u calcule le vecteur tangent a la courbe
36//-- et retourne le booleen Standard_False si ce calcul ne peut
37//-- pas etre mene a bien.
38//----------------------------------------------------------------------
39
40#include <Precision.hxx>
41
42#include <IntAna_Curve.ixx>
43
44#include <Standard_DomainError.hxx>
45#include <math_DirectPolynomialRoots.hxx>
46#include <ElSLib.hxx>
47#include <gp_XYZ.hxx>
48
49//=======================================================================
50//function : IntAna_Curve
51//purpose :
52//=======================================================================
53 IntAna_Curve::IntAna_Curve()
54{
55 typequadric=GeomAbs_OtherSurface;
56 firstbounded=Standard_False;
57 lastbounded=Standard_False;
58}
59//=======================================================================
60//function : SetConeQuadValues
61//purpose : Description de l intersection Cone Quadrique
62//=======================================================================
63 void IntAna_Curve::SetConeQuadValues(const gp_Cone& Cone,
64 const Standard_Real Qxx,
65 const Standard_Real Qyy,
66 const Standard_Real Qzz,
67 const Standard_Real Qxy,
68 const Standard_Real Qxz,
69 const Standard_Real Qyz,
70 const Standard_Real Qx,
71 const Standard_Real Qy,
72 const Standard_Real Qz,
73 const Standard_Real Q1,
74 const Standard_Real TOL,
75 const Standard_Real DomInf,
76 const Standard_Real DomSup,
77 const Standard_Boolean twocurves,
78 const Standard_Boolean takezpositive)
79{
80
81 Ax3 = Cone.Position();
82 RCyl = Cone.RefRadius();
83
84 Angle = Cone.SemiAngle();
85 Standard_Real UnSurTgAngle = 1.0/(Tan(Cone.SemiAngle()));
86
87 typequadric= GeomAbs_Cone;
88
89 TwoCurves = twocurves; //-- deux Z pour un meme parametre
90 TakeZPositive = takezpositive; //-- Prendre sur la courbe le Z Positif
91 //-- ( -B + Sqrt()) et non (-B - Sqrt())
92
93
94 Z0Cte = Q1; //-- Attention On a Z?Cos Cos(t)
95 Z0Sin = 0.0; //-- et Non 2 Z?Cos Cos(t) !!!
96 Z0Cos = 0.0; //-- Ce pour tous les Parametres
97 Z0CosCos = 0.0; //-- ie pas de Coefficient 2
98 Z0SinSin = 0.0; //-- devant les termes CS C S
99 Z0CosSin = 0.0;
100
101 Z1Cte = 2.0*(UnSurTgAngle)*Qz;
102 Z1Sin = Qy+Qy;
103 Z1Cos = Qx+Qx;
104 Z1CosCos = 0.0;
105 Z1SinSin = 0.0;
106 Z1CosSin = 0.0;
107
108 Z2Cte = Qzz * UnSurTgAngle*UnSurTgAngle;
109 Z2Sin = (UnSurTgAngle+UnSurTgAngle)*Qyz;
110 Z2Cos = (UnSurTgAngle+UnSurTgAngle)*Qxz;
111 Z2CosCos = Qxx;
112 Z2SinSin = Qyy;
113 Z2CosSin = Qxy+Qxy;
114
115 Tolerance = TOL;
116 DomainInf = DomInf;
117 DomainSup = DomSup;
118
119 RestrictedInf = RestrictedSup = Standard_True; //-- Le Domaine est Borne
120 firstbounded = lastbounded = Standard_False;
121}
122
123//=======================================================================
124//function : SetCylinderQuadValues
125//purpose : Description de l intersection Cylindre Quadrique
126//=======================================================================
127 void IntAna_Curve::SetCylinderQuadValues(const gp_Cylinder& Cyl,
128 const Standard_Real Qxx,
129 const Standard_Real Qyy,
130 const Standard_Real Qzz,
131 const Standard_Real Qxy,
132 const Standard_Real Qxz,
133 const Standard_Real Qyz,
134 const Standard_Real Qx,
135 const Standard_Real Qy,
136 const Standard_Real Qz,
137 const Standard_Real Q1,
138 const Standard_Real TOL,
139 const Standard_Real DomInf,
140 const Standard_Real DomSup,
141 const Standard_Boolean twocurves,
142 const Standard_Boolean takezpositive)
143{
144
145 Ax3 = Cyl.Position();
146 RCyl = Cyl.Radius();
147 typequadric= GeomAbs_Cylinder;
148
149 TwoCurves = twocurves; //-- deux Z pour un meme parametre
150 TakeZPositive = takezpositive; //-- Prendre sur la courbe le Z Positif
151 Standard_Real RCylmul2 = RCyl+RCyl; //-- ( -B + Sqrt())
152
153 Z0Cte = Q1;
154 Z0Sin = RCylmul2*Qy;
155 Z0Cos = RCylmul2*Qx;
156 Z0CosCos = Qxx*RCyl*RCyl;
157 Z0SinSin = Qyy*RCyl*RCyl;
158 Z0CosSin = RCylmul2*RCyl*Qxy;
159
160 Z1Cte = Qz+Qz;
161 Z1Sin = RCylmul2*Qyz;
162 Z1Cos = RCylmul2*Qxz;
163 Z1CosCos = 0.0;
164 Z1SinSin = 0.0;
165 Z1CosSin = 0.0;
166
167 Z2Cte = Qzz;
168 Z2Sin = 0.0;
169 Z2Cos = 0.0;
170 Z2CosCos = 0.0;
171 Z2SinSin = 0.0;
172 Z2CosSin = 0.0;
173
174 Tolerance = TOL;
175 DomainInf = DomInf;
176 DomainSup = DomSup;
177
178 RestrictedInf = RestrictedSup = Standard_True;
179 firstbounded = lastbounded = Standard_False;
180}
181
182//=======================================================================
183//function : IsOpen
184//purpose :
185//=======================================================================
186 Standard_Boolean IntAna_Curve::IsOpen() const
187{
188 return(RestrictedInf && RestrictedSup);
189}
190
191//=======================================================================
192//function : Domain
193//purpose :
194//=======================================================================
195 void IntAna_Curve::Domain(Standard_Real& DInf,
196 Standard_Real& DSup) const
197{
198 if(RestrictedInf && RestrictedSup) {
199 DInf=DomainInf;
200 DSup=DomainSup;
201 if(TwoCurves) {
202 DSup+=DSup-DInf;
203 }
204 }
205 else {
206 Standard_DomainError::Raise("IntAna_Curve::Domain");
207 }
208}
209//=======================================================================
210//function : IsConstant
211//purpose :
212//=======================================================================
213 Standard_Boolean IntAna_Curve::IsConstant() const
214{
215 //-- ??? Pas facile de decider a la seule vue des Param.
216 return(Standard_False);
217}
218
219//=======================================================================
220//function : IsFirstOpen
221//purpose :
222//=======================================================================
223 Standard_Boolean IntAna_Curve::IsFirstOpen() const
224{
225 return(firstbounded);
226}
227
228//=======================================================================
229//function : IsLastOpen
230//purpose :
231//=======================================================================
232 Standard_Boolean IntAna_Curve::IsLastOpen() const
233{
234 return(lastbounded);
235}
236//=======================================================================
237//function : SetIsFirstOpen
238//purpose :
239//=======================================================================
240 void IntAna_Curve::SetIsFirstOpen(const Standard_Boolean Flag)
241{
242 firstbounded = Flag;
243}
244
245//=======================================================================
246//function : SetIsLastOpen
247//purpose :
248//=======================================================================
249 void IntAna_Curve::SetIsLastOpen(const Standard_Boolean Flag)
250{
251 lastbounded = Flag;
252}
253
254//=======================================================================
255//function : InternalUVValue
256//purpose :
257//=======================================================================
258 void IntAna_Curve::InternalUVValue(const Standard_Real theta,
259 Standard_Real& Param1,
260 Standard_Real& Param2,
261 Standard_Real& A,
262 Standard_Real& B,
263 Standard_Real& C,
264 Standard_Real& cost,
265 Standard_Real& sint,
266 Standard_Real& SigneSqrtDis) const
267{
268 Standard_Real Theta=theta;
269 Standard_Boolean SecondSolution=Standard_False;
270
271 if((Theta<DomainInf) ||
272 ((Theta>DomainSup) && (!TwoCurves)) ||
273 (Theta>(DomainSup+DomainSup-DomainInf+0.00000000000001))) {
274 Standard_DomainError::Raise("IntAna_Curve::Domain");
275 }
276
277 if(Theta>DomainSup) {
278 Theta=DomainSup+DomainSup-Theta;
279 SecondSolution=Standard_True;
280 }
281
282 Param1=Theta;
283
284 if(!TwoCurves) {
285 SecondSolution=TakeZPositive;
286 }
287 //
288 cost = Cos(Theta);
289 sint = Sin(Theta);
290 Standard_Real costsint = cost*sint;
291
292 A=Z2Cte+sint*(Z2Sin+sint*Z2SinSin)+cost*(Z2Cos+cost*Z2CosCos)
293 +Z2CosSin*costsint;
294
295 B=Z1Cte+sint*(Z1Sin+sint*Z1SinSin)+cost*(Z1Cos+cost*Z1CosCos)
296 +Z1CosSin*costsint;
297
298 C=Z0Cte+sint*(Z0Sin+sint*Z0SinSin)+cost*(Z0Cos+cost*Z0CosCos)
299 +Z0CosSin*costsint;
300
301
302 Standard_Real Discriminant = B*B-4.0*A*C;
303
304 if(Abs(A)<=0.000000001) {
305 //-- cout<<" IntAna_Curve:: Internal UV Value : A="<<A<<" -> Abs(A)="<<Abs(A)<<endl;
306 if(Abs(B)<=0.0000000001) {
307 //-- cout<<" Probleme : Pas de solutions "<<endl;
308 Param2=0.0;
309 }
310 else {
311 //modified by NIZNHY-PKV Fri Dec 2 16:02:46 2005f
312 Param2 = -C/B;
313 /*
314 if(!SecondSolution) {
315 //-- Cas Param2 = (-B+Sqrt(Discriminant))/(A+A);
316 //-- = (-B+Sqrt(B**2 - Eps)) / 2A
317 //-- = -C / B
318 Param2 = -C/B;
319 }
320 else {
321 //-- Cas Param2 = (-B-Sqrt(Discriminant))/(A+A);
322 //-- = (-B-Sqrt(B**2 - Eps)) / 2A
323 if(A) {
324 Param2 = -B/A;
325 }
326 else {
327 Param2 = -B*10000000.0;
328 }
329 }
330 */
331 //modified by NIZNHY-PKV Fri Dec 2 16:02:54 2005t
332 }
333 }
334 else {
335 if(Discriminant<=0.0000000001 ||
336 Abs(Discriminant/(4*A))<=0.0000000001) Discriminant=0.0;
337 SigneSqrtDis = (SecondSolution)? Sqrt(Discriminant)
338 : -Sqrt(Discriminant);
339 Param2=(-B+SigneSqrtDis)/(A+A);
340 }
341}
342
343//=======================================================================
344//function : Value
345//purpose :
346//=======================================================================
347 gp_Pnt IntAna_Curve::Value(const Standard_Real theta)
348{
349 Standard_Real A, B, C, U, V, sint, cost, SigneSqrtDis;
350 //
351 A=0.0; B=0.0; C=0.0;
352 U=0.0; V=0.0;
353 sint=0.0; cost=0.0;
354 SigneSqrtDis=0.0;
355 InternalUVValue(theta,U,V,A,B,C,cost,sint,SigneSqrtDis);
356 //-- checked the parameter U and Raises Domain Error if Error
357 return(InternalValue(U,V));
358}
359//=======================================================================
360//function : D1u
361//purpose :
362//=======================================================================
363 Standard_Boolean IntAna_Curve::D1u(const Standard_Real theta,
364 gp_Pnt& Pt,
365 gp_Vec& Vec)
366{
367 //-- Pour detecter le cas ou le calcul est impossible
368 Standard_Real A, B, C, U, V, sint, cost, SigneSqrtDis;
369 A=0.0; B=0.0; C=0.0;
370 U=0.0; V=0.0;
371 sint=0.0; cost=0.0;
372 //
373 InternalUVValue(theta,U,V,A,B,C,cost,sint,SigneSqrtDis);
374 //
375 Pt = Value(theta);
1ef32e96 376 if(Abs(A)<0.0000001 || Abs(SigneSqrtDis)<0.0000000001) return(Standard_False);
7fd59977 377
378
379 //-- Approximation de la derivee (mieux que le calcul mathematique!)
380 Standard_Real dtheta = (DomainSup-DomainInf)*0.000001;
381 Standard_Real theta2 = theta+dtheta;
382 if((theta2<DomainInf) || ((theta2>DomainSup) && (!TwoCurves))
383 || (theta2>(DomainSup+DomainSup-DomainInf+0.00000000000001))) {
384 dtheta = -dtheta;
385 theta2 = theta+dtheta;
386 }
387 gp_Pnt P2 = Value(theta2);
388 dtheta = 1.0/dtheta;
389 Vec.SetCoord((P2.X()-Pt.X())*dtheta,
390 (P2.Y()-Pt.Y())*dtheta,
391 (P2.Z()-Pt.Z())*dtheta);
392
393 return(Standard_True);
394}
395//=======================================================================
396//function : FindParameter
397//purpose : Para est en sortie le parametre sur la courbe
398//=======================================================================
399 Standard_Boolean IntAna_Curve::FindParameter (const gp_Pnt& P,
400 Standard_Real& Para) const
401{
402 Standard_Real theta,z, aTolPrecision=0.0001;
c6541a0c 403 Standard_Real PIpPI = M_PI + M_PI;
7fd59977 404 //
405 switch (typequadric) {
406
407 case GeomAbs_Cylinder:
408 {
409 ElSLib::CylinderParameters(Ax3,RCyl,P,theta,z);
410 }
411 break;
412
413 case GeomAbs_Cone :
414 {
415 ElSLib::ConeParameters(Ax3,RCyl,Angle,P,theta,z);
416 }
417 break;
418
419 default:
841d6672 420 return Standard_False;
7fd59977 421 break;
422 }
423 //
424 Standard_Real epsAng = 1.e-8;
425 Standard_Real tmin = DomainInf;
426 Standard_Real tmax = DomainSup;
427 Standard_Real U,V,A,B,C,sint,cost,SigneSqrtDis;
428 Standard_Real z1,z2;
429
430 A=0.0; B=0.0; C=0.0;
431 U=0.0; V=0.0;
432 sint=0.0; cost=0.0;
433 SigneSqrtDis=0.0;
434 //U=V=A=B=C=sint=cost=SigneSqrtDis=0.0;
435 //
436 if (!firstbounded && tmin > theta && (tmin-theta) <= epsAng) {
437 theta = tmin;
438 }
439 else if (!lastbounded && theta > tmax && (theta-tmax) <= epsAng) {
440 theta = tmax;
441 }
442 //
443 if (theta < tmin ) {
444 theta = theta + PIpPI;
445 }
446 else if (theta > tmax) {
447 theta = theta - PIpPI;
448 }
449 if (theta < tmin || theta > tmax) {
450 if(theta>tmax) {
451 InternalUVValue(tmax,U,V,A,B,C,cost,sint,SigneSqrtDis);
452 gp_Pnt PMax(InternalValue(U,V));
453 if(PMax.Distance(P) < aTolPrecision) {
454 Para = tmax;
455 return(Standard_True);
456 }
457 }
458 if(theta<tmin) {
459 InternalUVValue(tmin,U,V,A,B,C,cost,sint,SigneSqrtDis);
460 gp_Pnt PMin(InternalValue(U,V));
461 if(PMin.Distance(P) < aTolPrecision) {
462 Para = tmin;
463 return(Standard_True);
464 }
465 }
466 //-- lbr le 14 Fev 96 : On teste malgre tout si le point n est pas le
467 //-- point de debut ou de fin
468 //-- cout<<"False 1 "<<endl;
469 // theta = tmin; le 25 Nov 96
470 }
471
472 if (TwoCurves) {
473 if(theta > tmax)
474 theta = tmax;
475 if(theta < tmin)
476 theta = tmin;
477 InternalUVValue(theta,U,z1,A,B,C,cost,sint,SigneSqrtDis);
478 A = B = C = sint = cost = SigneSqrtDis = 0.0;
479 InternalUVValue(tmax+tmax - theta,U,z2,A,B,C,cost,sint,SigneSqrtDis);
480
481 if (Abs(z-z1) <= Abs(z-z2)) {
482 Para = theta;
483 }
484 else {
485 Para = tmax+tmax - theta;
486 }
487 }
488 else {
489 Para = theta;
490 }
491
492 if((Para<DomainInf) || ((Para>DomainSup) && (!TwoCurves))
493 || (Para>(DomainSup+DomainSup-DomainInf+0.00000000000001))) {
494 return(Standard_False);
495 }
496
497 InternalUVValue(Para,U,V,A,B,C,cost,sint,SigneSqrtDis);
498 gp_Pnt PPara = InternalValue(U,V);
499 Standard_Real Dist = PPara.Distance(P);
500 if(Dist > aTolPrecision) {
501 //-- Il y a eu un probleme
502 //-- On teste si le point est un point double
503 InternalUVValue(tmin,U,V,A,B,C,cost,sint,SigneSqrtDis);
504 PPara = InternalValue(U,V);
505 Dist = PPara.Distance(P);
506 if(Dist <= aTolPrecision) {
507 Para = tmin;
508 return(Standard_True);
509 }
510
511 InternalUVValue(tmax,U,V,A,B,C,cost,sint,SigneSqrtDis);
512 PPara = InternalValue(U,V);
513 Dist = PPara.Distance(P);
514 if(Dist <= aTolPrecision) {
515 Para = tmax;
516 return(Standard_True);
517 }
518 if (TwoCurves) {
519 Standard_Real Theta = DomainSup+DomainSup-DomainInf;
520 InternalUVValue(Theta,U,V,A,B,C,cost,sint,SigneSqrtDis);
521 PPara = InternalValue(U,V);
522 Dist = PPara.Distance(P);
523 if(Dist <= aTolPrecision) {
524 Para = Theta;
525 return(Standard_True);
526 }
527 }
528 return(Standard_False);
529 }
530 return(Standard_True);
531}
532//=======================================================================
533//function : InternalValue
534//purpose :
535//=======================================================================
536 gp_Pnt IntAna_Curve::InternalValue(const Standard_Real U,
537 const Standard_Real _V) const
538{
539 //-- cout<<" ["<<U<<","<<V<<"]";
540 Standard_Real V = _V;
541 if(V > 100000.0 ) { V= 100000.0; }
542 if(V < -100000.0 ) { V=-100000.0; }
543
544 switch(typequadric) {
545 case GeomAbs_Cone:
546 {
547 //------------------------------------------------
548 //-- Parametrage : X = V * Cos(U) ---
549 //-- Y = V * Sin(U) ---
550 //-- Z = (V-RCyl) / Tan(SemiAngle)--
551 //------------------------------------------------
552 //-- Angle Vaut Cone.SemiAngle()
553 return(ElSLib::ConeValue(U,(V-RCyl)/Sin(Angle),Ax3,RCyl,Angle));
554 }
555 break;
556
557 case GeomAbs_Cylinder:
558 return(ElSLib::CylinderValue(U,V,Ax3,RCyl));
559 case GeomAbs_Sphere:
560 return(ElSLib::SphereValue(U,V,Ax3,RCyl));
561 default:
562 return(gp_Pnt(0.0,0.0,0.0));
563 }
7fd59977 564}
565
566//
567//=======================================================================
568//function : SetDomain
569//purpose :
570//=======================================================================
571 void IntAna_Curve::SetDomain(const Standard_Real DInf,
572 const Standard_Real DSup)
573{
574 if(DInf>=DSup) {
575 Standard_DomainError::Raise("IntAna_Curve::Domain");
576 }
577 //
578 DomainInf=DInf;
579 DomainSup=DSup;
580}
Open CASCADE Technology repository