Integration of OCCT 6.5.0 from SVN
[occt.git] / src / Bisector / Bisector_BisecCC.cxx
CommitLineData
7fd59977 1// File: Bisector_BisecCC.cxx
2// Created: Thu Mar 10 17:54:52 1994
3// Author: Yves FRICAUD
4// <yfr@phylox>
5
6#include <Bisector_BisecCC.ixx>
7#include <Bisector_BisecPC.hxx>
8#include <Bisector.hxx>
9#include <Bisector_Curve.hxx>
10#include <Bisector_FunctionH.hxx>
11#include <Bisector_PointOnBis.hxx>
12#include <Geom2dAdaptor_Curve.hxx>
13#include <Geom2d_Curve.hxx>
14#include <Geom2dLProp_CLProps2d.hxx>
15#include <Geom2dGcc.hxx>
16#include <Geom2dGcc_Circ2d2TanRad.hxx>
17#include <Geom2dGcc_QualifiedCurve.hxx>
18#include <Geom2d_TrimmedCurve.hxx>
19#include <Geom2d_Circle.hxx>
20#include <Geom2d_Line.hxx>
21#include <Geom2dInt_GInter.hxx>
22#include <Geom2dAPI_ProjectPointOnCurve.hxx>
23#include <gp_Pnt2d.hxx>
24#include <gp_Vec2d.hxx>
25#include <gp.hxx>
26#include <IntRes2d_IntersectionPoint.hxx>
27#include <Precision.hxx>
28#include <math_FunctionRoot.hxx>
29#include <math_FunctionRoots.hxx>
30#include <math_BissecNewton.hxx>
31
32#include <Standard_OutOfRange.hxx>
33#include <Standard_DivideByZero.hxx>
34#include <Standard_NotImplemented.hxx>
35
36
37static Standard_Real ProjOnCurve (const gp_Pnt2d& P,
38 const Handle(Geom2d_Curve)& C);
39
40static Standard_Real Curvature (const Handle(Geom2d_Curve)& C,
41 Standard_Real U,
42 Standard_Real Tol) ;
43
44static Standard_Boolean TestExtension (const Handle(Geom2d_Curve)& C1,
45 const Handle(Geom2d_Curve)& C2,
46 const Standard_Integer Start_End);
47
48static Standard_Boolean DiscretPar(const Standard_Real DU,
49 const Standard_Real EpsMin,
50 const Standard_Real EpsMax,
51 const Standard_Integer NbMin,
52 const Standard_Integer NbMax,
53 Standard_Real& Eps,
54 Standard_Integer& Nb);
55
56//=============================================================================
57//function :
58//purpose :
59//=============================================================================
60Bisector_BisecCC::Bisector_BisecCC()
61{
62 shiftParameter = 0;
63 isEmpty = Standard_False;
64}
65
66//=============================================================================
67//function :
68//purpose :
69//=============================================================================
70Bisector_BisecCC::Bisector_BisecCC(const Handle(Geom2d_Curve)& Cu1,
71 const Handle(Geom2d_Curve)& Cu2,
72 const Standard_Real Side1,
73 const Standard_Real Side2,
74 const gp_Pnt2d& Origin,
75 const Standard_Real DistMax)
76{
77 Perform (Cu1,Cu2,Side1,Side2,Origin,DistMax);
78}
79
80//=============================================================================
81//function : Perform
82//purpose :
83//=============================================================================
84void Bisector_BisecCC::Perform(const Handle(Geom2d_Curve)& Cu1,
85 const Handle(Geom2d_Curve)& Cu2,
86 const Standard_Real Side1,
87 const Standard_Real Side2,
88 const gp_Pnt2d& Origin,
89 const Standard_Real DistMax)
90{
91 isEmpty = Standard_False;
92 distMax = DistMax;
93
94 curve1 = Handle (Geom2d_Curve)::DownCast(Cu1->Copy());
95 curve2 = Handle (Geom2d_Curve)::DownCast(Cu2->Copy());
96
97 sign1 = Side1;
98 sign2 = Side2;
99
100 isConvex1 = Bisector::IsConvex(curve1,sign1);
101 isConvex2 = Bisector::IsConvex(curve2,sign2);
102
103 Standard_Real U,UC1,UC2,Dist,dU,USol;
104 gp_Pnt2d P;
105 Standard_Integer NbPnts = 21;
106 Standard_Real EpsMin = 10*Precision::Confusion();
107 Standard_Boolean YaPoly = Standard_True;
108 Standard_Boolean OriInPoly = Standard_False;
109 //---------------------------------------------
110 // Calcul premier point du polygone.
111 //---------------------------------------------
112 U = ProjOnCurve (Origin,curve1);
113 P = ValueByInt (U,UC1,UC2,Dist);
114
115 if (Dist < Precision::Infinite()) {
116 //----------------------------------------------------
117 // le parametre du point origine donne un point sur le
118 // polygone.
119 //----------------------------------------------------
120 myPolygon.Append(Bisector_PointOnBis(UC1,UC2,U,Dist,P));
121 startIntervals.Append(U);
122 if (P.IsEqual(Origin,Precision::Confusion())) {
123 //----------------------------------------
124 // test si le prenier point est l origine.
125 //----------------------------------------
126 OriInPoly = Standard_True;
127 }
128 }
129 else {
130 //-------------------------------------------------------
131 // Le point origine est sur un prolongement.
132 // Recherche du premier point du polygone par dichotomie.
133 //-------------------------------------------------------
134 dU = (curve1->LastParameter() - U)/(NbPnts - 1);
135 U += dU;
136 for (Standard_Integer i = 1; i <= NbPnts - 1; i++) {
137 P = ValueByInt(U,UC1,UC2,Dist);
138 if (Dist < Precision::Infinite()) {
139 USol = SearchBound(U - dU,U);
140 P = ValueByInt(USol,UC1,UC2,Dist);
141 startIntervals.Append(USol);
142 myPolygon.Append(Bisector_PointOnBis(UC1,UC2,USol,Dist,P));
143 break;
144 }
145 U += dU;
146 }
147 }
148
149 if ( !myPolygon.Length() == 0) {
150 SupLastParameter();
151 //----------------------------------------------
152 // Construction du polygone de la bissectrice.
153 //---------------------------------------------
154 U = FirstParameter();
155 Standard_Real DU = LastParameter() - U;
156
157 if (DU < EpsMin) {NbPnts = 3;}
158 dU = DU/(NbPnts - 1);
159
160 U += dU;
161// modified by NIZHNY-EAP Fri Jan 21 09:33:20 2000 ___BEGIN___
162// prevent addition of the same point
163 gp_Pnt2d prevPnt = P;
164 for (Standard_Integer i = 1; i <= NbPnts - 1; i++) {
165 P = ValueByInt(U,UC1,UC2,Dist);
166 if (Dist < Precision::Infinite()) {
167 if (P.Distance (prevPnt) > Precision::Confusion())
168 myPolygon.Append(Bisector_PointOnBis(UC1,UC2,U,Dist,P));
169 }
170 else {
171 USol = SearchBound(U - dU,U);
172 P = ValueByInt(USol,UC1,UC2,Dist);
173 endIntervals.SetValue(1,USol);
174 if (P.Distance (prevPnt) > Precision::Confusion())
175 myPolygon.Append(Bisector_PointOnBis(UC1,UC2,USol,Dist,P));
176 break;
177 }
178 U += dU;
179 prevPnt=P;
180// modified by NIZHNY-EAP Fri Jan 21 09:33:24 2000 ___END___
181 }
182 }
183 else {
184 //----------------
185 // Polygone vide.
186 //----------------
187 YaPoly = Standard_False;
188 }
189
190 extensionStart = Standard_False;
191 extensionEnd = Standard_False;
192 pointStart = Origin;
193
194 if (isConvex1 && isConvex2) {
195 if (YaPoly) pointEnd = myPolygon.Last().Point();
196 }
197 else {
198 //-----------------------------------------------------------------------------
199 // Prolongement : La courbe est prolongee au debut ou/et a la fin si
200 // - une des deux courbes est concave.
201 // - Les courbes ont un point commun au debut ou/et a la fin
202 // - l angle d ouverture au point commun entre les deux courbes
203 // vaut PI.
204 // le prolongemt au debut est pris en compte si l origine se trouve dessus.
205 // ie : l origine n est pas dans le polygone.
206 //-----------------------------------------------------------------------------
207
208 //---------------------------------
209 // Existent ils des prolongemnets ?
210 //---------------------------------
211 if (OriInPoly) {
212 extensionStart = Standard_False;
213 }
214 else {
215 extensionStart = TestExtension(curve1,curve2,1);
216 }
217 extensionEnd = TestExtension(curve1,curve2,2);
218
219 //-----------------
220 // Calcul pointEnd.
221 //-----------------
222 if (extensionEnd) {
223 pointEnd = curve1->Value(curve1->LastParameter());
224 }
225 else if (YaPoly) {
226 pointEnd = myPolygon.Last().Point();
227 }
228 else {
229 ComputePointEnd();
230 }
231 //------------------------------------------------------
232 // Mise a jour des Bornes des intervalles de definition.
233 //------------------------------------------------------
234 if (YaPoly) {
235 if (extensionStart) {
236 gp_Pnt2d P1 = myPolygon.First().Point();
237 Standard_Real UFirst = startIntervals.First() - pointStart.Distance(P1);
238 startIntervals.InsertBefore(1,UFirst);
239 endIntervals .InsertBefore(1,startIntervals.Value(2));
240 }
241 if (extensionEnd) {
242 gp_Pnt2d P1;
243 Standard_Real UFirst,ULast;
244 P1 = myPolygon.Last().Point();
245 UFirst = endIntervals.Last();
246 ULast = UFirst + pointEnd.Distance(P1);
247 startIntervals.Append(UFirst);
248 endIntervals .Append(ULast );
249 }
250 }
251 else {
252 //--------------------------------------------------
253 // Pas de polygone => la bissectrise est un segment.
254 //--------------------------------------------------
255 startIntervals.Append(0.);
256 endIntervals .Append(pointEnd.Distance(pointStart));
257 }
258 }
259 if (!YaPoly && !extensionStart && !extensionEnd)
260 isEmpty = Standard_True;
261// modified by NIZHNY-EAP Mon Jan 17 17:32:40 2000 ___BEGIN___
262 if (myPolygon.Length() <= 2)
263 isEmpty = Standard_True;
264// modified by NIZHNY-EAP Mon Jan 17 17:32:42 2000 ___END___
265}
266
267//=============================================================================
268//function : IsExtendAtStart
269//purpose :
270//=============================================================================
271Standard_Boolean Bisector_BisecCC::IsExtendAtStart() const
272{
273 return extensionStart;
274}
275
276//=============================================================================
277//function : IsExtendAtEnd
278//purpose :
279//=============================================================================
280Standard_Boolean Bisector_BisecCC::IsExtendAtEnd() const
281{
282 return extensionEnd;
283}
284
285//=============================================================================
286//function : IsEmpty
287//purpose :
288//=============================================================================
289Standard_Boolean Bisector_BisecCC::IsEmpty() const
290{
291 return isEmpty;
292}
293
294//=============================================================================
295//function : Reverse
296//purpose :
297//=============================================================================
298void Bisector_BisecCC::Reverse()
299{
300 Standard_NotImplemented::Raise();
301}
302
303//=============================================================================
304//function : ReversedParameter
305// purpose :
306//=============================================================================
307Standard_Real Bisector_BisecCC::ReversedParameter(const Standard_Real U) const
308{
309 return LastParameter() + FirstParameter() - U;
310}
311
312//=============================================================================
313//function : Copy
314//purpose :
315//=============================================================================
316Handle(Geom2d_Geometry) Bisector_BisecCC::Copy() const
317{
318 Handle(Geom2d_Curve) CopyCurve1
319 = Handle(Geom2d_Curve)::DownCast(curve1->Copy());
320 Handle(Geom2d_Curve) CopyCurve2
321 = Handle(Geom2d_Curve)::DownCast(curve2->Copy());
322
323 Handle(Bisector_BisecCC) C = new Bisector_BisecCC();
324
325 C -> Curve (1, CopyCurve1) ; C -> Curve (2, CopyCurve2);
326 C -> Sign (1, sign1 ) ; C -> Sign (2, sign2 );
327 C -> IsConvex (1, isConvex1) ; C -> IsConvex (2, isConvex2);
328 C -> Polygon (myPolygon);
329 C -> IsEmpty (isEmpty) ;
330 C -> DistMax (distMax) ;
331 C -> StartIntervals (startIntervals); C -> EndIntervals (endIntervals);
332 C -> ExtensionStart (extensionStart); C -> ExtensionEnd (extensionEnd);
333 C -> PointStart (pointStart) ; C -> PointEnd (pointEnd) ;
334
335 return C;
336}
337
338//=============================================================================
339//function : ChangeGuide
340//purpose : Changement de la ligne guide pour le parametrage de la bissectrice
341// ATTENTION : - Ceci peut inverser le sens de parametrage.
342// - Ceci ne concerne que la partie de la courbe
343// correspondante au polygone.
344//=============================================================================
345Handle(Bisector_BisecCC) Bisector_BisecCC::ChangeGuide() const
346{
347 Handle(Bisector_BisecCC) C = new Bisector_BisecCC();
348
349 C -> Curve (1, curve2) ; C -> Curve (2, curve1);
350 C -> Sign (1, sign2 ) ; C -> Sign (2, sign1 );
351 C -> IsConvex (1, isConvex2); C -> IsConvex (2, isConvex1);
352
353 //-------------------------------------------------------------------------
354 // Construction du nouveau polygone a partir de celui d origine.
355 // inversion des PointOnBis et Calcul nouveau parametre sur la bissectrice.
356 //-------------------------------------------------------------------------
357 Bisector_PolyBis Poly;
358 if (sign1 == sign2 ) {
359 //---------------------------------------------------------------
360 // les elements du nouveau polygone sont ranges dans l autre sens.
361 //---------------------------------------------------------------
362 for (Standard_Integer i = myPolygon.Length(); i >=1; i--) {
363 Bisector_PointOnBis P = myPolygon.Value(i);
364 Bisector_PointOnBis NewP (P.ParamOnC2(), P.ParamOnC1(),
365 P.ParamOnC2(), P.Distance (),
366 P.Point());
367 Poly.Append (NewP);
368 }
369 }
370 else {
371 for (Standard_Integer i = 1; i <= myPolygon.Length(); i ++) {
372 Bisector_PointOnBis P = myPolygon.Value(i);
373 Bisector_PointOnBis NewP (P.ParamOnC2(), P.ParamOnC1(),
374 P.ParamOnC2(), P.Distance (),
375 P.Point());
376 Poly.Append (NewP);
377 }
378 }
379 C -> Polygon (Poly);
380 C -> FirstParameter (Poly.First().ParamOnBis());
381 C -> LastParameter (Poly.Last() .ParamOnBis());
382
383 return C;
384}
385
386//=============================================================================
387//function : Transform
388//purpose :
389//=============================================================================
390void Bisector_BisecCC::Transform (const gp_Trsf2d& T)
391{
392 curve1 ->Transform(T);
393 curve2 ->Transform(T);
394 myPolygon . Transform(T);
395 pointStart. Transform(T);
396 pointEnd . Transform(T);
397}
398
399//=============================================================================
400//function : IsCN
401//purpose :
402//=============================================================================
403Standard_Boolean Bisector_BisecCC::IsCN (const Standard_Integer N) const
404{
405 return (curve1->IsCN(N+1) && curve2->IsCN(N+1));
406}
407
408//=============================================================================
409//function : FirstParameter
410//purpose :
411//=============================================================================
412Standard_Real Bisector_BisecCC::FirstParameter() const
413{
414 return startIntervals.First();
415}
416
417//=============================================================================
418//function : LastParameter
419//purpose :
420//=============================================================================
421Standard_Real Bisector_BisecCC::LastParameter() const
422{
423 return endIntervals.Last();
424}
425
426//=============================================================================
427//function : Continuity
428//purpose :
429//=============================================================================
430GeomAbs_Shape Bisector_BisecCC::Continuity() const
431{
432 GeomAbs_Shape Cont = curve1->Continuity();
433 switch (Cont) {
434 case GeomAbs_C1 : return GeomAbs_C0;
435 case GeomAbs_C2 : return GeomAbs_C1;
436 case GeomAbs_C3 : return GeomAbs_C2;
437 case GeomAbs_CN : return GeomAbs_CN;
438#ifndef DEB
439 default: break;
440#endif
441 }
442 return GeomAbs_C0;
443}
444
445//=============================================================================
446//function : NbIntervals
447//purpose :
448//=============================================================================
449Standard_Integer Bisector_BisecCC::NbIntervals() const
450{
451 return startIntervals.Length();
452}
453
454//=============================================================================
455//function : IntervalFirst
456//purpose :
457//=============================================================================
458Standard_Real Bisector_BisecCC::IntervalFirst(const Standard_Integer Index) const
459{
460 return startIntervals.Value(Index);
461}
462
463//=============================================================================
464//function : IntervalLast
465//purpose :
466//=============================================================================
467Standard_Real Bisector_BisecCC::IntervalLast(const Standard_Integer Index) const
468{
469 return endIntervals.Value(Index);
470}
471
472//=============================================================================
473//function : IntervalContinuity
474//purpose :
475//=============================================================================
476GeomAbs_Shape Bisector_BisecCC::IntervalContinuity() const
477{
478 GeomAbs_Shape Cont = curve1->Continuity();
479 switch (Cont) {
480 case GeomAbs_C1 : return GeomAbs_C0;
481 case GeomAbs_C2 : return GeomAbs_C1;
482 case GeomAbs_C3 : return GeomAbs_C2;
483 case GeomAbs_CN : return GeomAbs_CN;
484#ifndef DEB
485 default: break;
486#endif
487 }
488 return GeomAbs_C0;
489}
490
491//=============================================================================
492//function : IsClosed
493//purpose :
494//=============================================================================
495Standard_Boolean Bisector_BisecCC::IsClosed() const
496{
497 if (curve1->IsClosed()) {
498 if (startIntervals.First() == curve1->FirstParameter() &&
499 endIntervals .Last () == curve1->LastParameter () )
500 return Standard_True;
501 }
502 return Standard_False;
503}
504
505//=============================================================================
506//function : IsPeriodic
507//purpose :
508//=============================================================================
509Standard_Boolean Bisector_BisecCC::IsPeriodic() const
510{
511 return Standard_False;
512}
513
514
515//=============================================================================
516//function : Curvature
517//purpose :
518//=============================================================================
519static Standard_Real Curvature (const Handle(Geom2d_Curve)& C,
520 Standard_Real U,
521 Standard_Real Tol)
522{
523 Standard_Real K1;
524 gp_Vec2d D1,D2;
525 gp_Pnt2d P;
526 C->D2(U,P,D1,D2);
527 Standard_Real Norm2 = D1.SquareMagnitude();;
528 if (Norm2 < Tol) {
529 K1 = 0.0;
530 }
531 else {
532 K1 = (D1^D2)/(Norm2*sqrt(Norm2));
533 }
534 return K1;
535}
536
537//=============================================================================
538//function : Value
539//purpose : CALCUL DU POINT COURANT PAR METHODE ITERATIVE.
540// ----------------------------------------------
541// Calcul du point courant, de la distance du point courant aux deux
542// courbes, des parametres sur chaque courbe de la projection du
543// point courrant.
544//
545//method : - Recheche parametre de depart en utilisant <myPolygon>.
546// - Calcul du parametre U2 sur la courbe C2 solution de H(U,V)= 0
547// - P(U) = F(U,U2)
548//
549// ou :
550// ||P2(v0)P1(u)||**2
551// F(u,v) = P1(u) - 1/2 *----------------* N(u)
552// (N(u).P2(v0)P1(u))
553//
554// H(u,v) = (Tu.P1(u)P2(v))**2||Tv||**2 - (Tv.P1(u)P2(v))**2||Tu||**2
555//
556//=============================================================================
557gp_Pnt2d Bisector_BisecCC::ValueAndDist (const Standard_Real U,
558 Standard_Real& U1,
559 Standard_Real& U2,
560 Standard_Real& Dist) const
561{
562 gp_Vec2d T;
563
564 //-----------------------------------------------
565 // le polygone est il reduit a un point ou vide?
566 //-----------------------------------------------
567 if (myPolygon.Length() <= 1) {
568 return Extension(U,U1,U2,Dist,T);
569 }
570
571 //-----------------------------------------------
572 // test U en dehors des bornes du polygone.
573 //-----------------------------------------------
574 if (U < myPolygon.First().ParamOnBis()) {
575 return Extension(U,U1,U2,Dist,T);
576 }
577 if (U > myPolygon.Last().ParamOnBis()) {
578 return Extension(U,U1,U2,Dist,T);
579 }
580
581 //-------------------------------------------------------
582 // Recheche parametre de depart en utilisant <myPolygon>.
583 //-------------------------------------------------------
584 Standard_Integer IntervalIndex = myPolygon.Interval(U);
585 Standard_Real UMin = myPolygon.Value(IntervalIndex ).ParamOnBis();
586 Standard_Real UMax = myPolygon.Value(IntervalIndex + 1).ParamOnBis();
587 Standard_Real VMin = myPolygon.Value(IntervalIndex ).ParamOnC2();
588 Standard_Real VMax = myPolygon.Value(IntervalIndex + 1).ParamOnC2();
589 Standard_Real Alpha,VInit;
590
591 if (Abs(UMax - UMin) < gp::Resolution()) {
592 VInit = VMin;
593 }
594 else {
595 Alpha = (U - UMin)/(UMax - UMin);
596 VInit = VMin + Alpha*(VMax - VMin);
597 }
598
599 U1 = LinkBisCurve(U);
600 Standard_Real VTemp = Min(VMin,VMax);
601 VMax = Max(VMin,VMax); VMin = VTemp;
602 Standard_Boolean Valid = Standard_True;
603 //---------------------------------------------------------------
604 // Calcul du parametre U2 sur la courbe C2 solution de H(u,v)= 0
605 //---------------------------------------------------------------
606 gp_Pnt2d P1;
607 gp_Vec2d T1;
608 Standard_Real EpsH = 1.E-8;
609 Standard_Real EpsH100 = 1.E-6;
610 curve1->D1 (U1,P1,T1);
611 gp_Vec2d N1(T1.Y(), - T1.X());
612
613 if ((VMax - VMin) < Precision::PConfusion()) {
614 U2 = VInit;
615 }
616 else {
617 Bisector_FunctionH H (curve2,P1,sign1*sign2*T1);
618 Standard_Real FInit;
619 H.Value(VInit,FInit);
620 if (Abs(FInit) < EpsH) {
621 U2 = VInit;
622 }
623 else {
624 math_BissecNewton SolNew (H,VMin - EpsH100,VMax + EpsH100,EpsH,10);
625 if (SolNew.IsDone()) {
626 U2 = SolNew.Root();
627 }
628 else {
629 math_FunctionRoot SolRoot (H,VInit,EpsH,VMin - EpsH100,VMax + EpsH100);
630 if (SolRoot.IsDone()) {
631 U2 = SolRoot.Root();
632 }
633 else { Valid = Standard_False;}
634 }
635 }
636 }
637
638 gp_Pnt2d PBis = pointStart;
639 //----------------
640 // P(U) = F(U1,U2)
641 //----------------
642 if (Valid) {
643 gp_Pnt2d P2 = curve2->Value(U2);
644 gp_Vec2d P2P1(P1.X() - P2.X(),P1.Y() - P2.Y());
645 Standard_Real SquareP2P1 = P2P1.SquareMagnitude();
646 Standard_Real N1P2P1 = N1.Dot(P2P1);
647
648 if (P1.IsEqual(P2,Precision::Confusion())) {
649 PBis = P1 ;
650 Dist = 0.0;
651 }
652 else if (N1P2P1*sign1 < 0) {
653 Valid = Standard_False;
654 }
655 else {
656 PBis = P1.Translated(- (0.5*SquareP2P1/N1P2P1)*N1);
657 Dist = P1.SquareDistance(PBis);
658 }
659 }
660
661 //----------------------------------------------------------------
662 // Si le point n est pas valide
663 // calcul par intersection.
664 //----------------------------------------------------------------
665 if (!Valid) {
666 //--------------------------------------------------------------------
667 // Construction de la bisectrice point courbe et de la droite passant
668 // par P1 et portee par la normale. la curve2 est restreinte par VMin et
669 // VMax.
670 //--------------------------------------------------------------------
671 Standard_Real DMin = Precision::Infinite();
672 gp_Pnt2d P;
673 Handle(Bisector_BisecPC) BisPC
674 = new Bisector_BisecPC(curve2, P1, sign2, VMin, VMax);
675 Handle(Geom2d_Line) NorLi = new Geom2d_Line (P1,N1);
676
677 Geom2dAdaptor_Curve ABisPC(BisPC);
678 Geom2dAdaptor_Curve ANorLi(NorLi);
679 //-------------------------------------------------------------------------
680 Geom2dInt_GInter Intersect(ABisPC,ANorLi,
681 Precision::Confusion(),Precision::Confusion());
682 //-------------------------------------------------------------------------
683
684 if (Intersect.IsDone() && !Intersect.IsEmpty()) {
685 for (Standard_Integer i = 1; i <= Intersect.NbPoints(); i++) {
686 if (Intersect.Point(i).ParamOnSecond()*sign1 < Precision::PConfusion()) {
687 P = Intersect.Point(i).Value();
688 if (P.SquareDistance(P1) < DMin) {
689 DMin = P.SquareDistance(P1);
690 PBis = P;
691 U2 = BisPC->LinkBisCurve(Intersect.Point(i).ParamOnFirst());
692 Dist = DMin;
693 }
694 }
695 }
696 }
697 }
698 return PBis;
699}
700
701//=============================================================================
702//function : ValueByInt
703//purpose : CALCUL DU POINT COURANT PAR INTERSECTION.
704// -----------------------------------------
705// Calcul du point courant, de la distance du point courant aux deux
706// courbes, des parametres sur chaque courbe de la projection du
707// point courrant.
708// le point courrant au parametre U est l intersection de la
709// bissectrice point courbe (P1,curve2) et de la droite passant par
710// P1 de vecteur directeur N1.
711// P1 est le point courrant de parametre U sur la curve1 et N1 la
712// normale en ce point.
713//=============================================================================
714gp_Pnt2d Bisector_BisecCC::ValueByInt (const Standard_Real U,
715 Standard_Real& U1,
716 Standard_Real& U2,
717 Standard_Real& Dist) const
718{
719 //------------------------------------------------------------------
720 // Recuperation des point,tangente,normale sur C1 au parametre U.
721 //-------------------------------------------------------------------
722 U1 = LinkBisCurve(U);
723
724 gp_Pnt2d P1,P2,P,PSol;
725 gp_Vec2d Tan1,Tan2;
726 curve1->D1(U1,P1,Tan1);
727 gp_Vec2d N1( Tan1.Y(), - Tan1.X());
728
729 //--------------------------------------------------------------------------
730 // test de confusion de P1 avec extremite de curve2.
731 //--------------------------------------------------------------------------
732 if (P1.Distance(curve2->Value(curve2->FirstParameter())) < Precision::Confusion()) {
733 U2 = curve2->FirstParameter();
734 curve2->D1(U2,P2,Tan2);
735 if ( isConvex1 && isConvex2 ) {
736 Dist = 0.;
737 return P1;
738 }
739 if (! Tan1.IsParallel(Tan2,Precision::Angular())) {
740 Dist = 0.;
741 return P1;
742 }
743 }
744 if (P1.Distance(curve2->Value(curve2->LastParameter())) < Precision::Confusion()) {
745 U2 = curve2->LastParameter();
746 curve2->D1(U2,P2,Tan2);
747 if ( isConvex1 && isConvex2 ) {
748 Dist = 0.;
749 return P1;
750 }
751 if (! Tan1.IsParallel(Tan2,Precision::Angular())) {
752 Dist = 0.;
753 return P1;
754 }
755 }
756
757 Standard_Boolean YaSol = Standard_False;
758 Standard_Real DMin = Precision::Infinite();
759 Standard_Real USol;
760 Standard_Real EpsMax = 1.E-6;
761 Standard_Real EpsX;
762 Standard_Real EpsH = 1.E-8;
763 Standard_Real DistPP1;
764 Standard_Integer NbSamples =20;
765 Standard_Real UFirstOnC2 = curve2->FirstParameter();
766 Standard_Real ULastOnC2 = curve2->LastParameter();
767
768 if (!myPolygon.IsEmpty()){
769 if (sign1 == sign2) { ULastOnC2 = myPolygon.Last().ParamOnC2();}
770 else { UFirstOnC2 = myPolygon.Last().ParamOnC2();}
771 }
772
773 if (Abs(ULastOnC2 - UFirstOnC2) < Precision::PConfusion()/100.) {
774 Dist = Precision::Infinite();
775 return P1;
776 }
777
778 DiscretPar(Abs(ULastOnC2 - UFirstOnC2),EpsH,EpsMax,2,20,EpsX,NbSamples);
779
780 Bisector_FunctionH H (curve2,P1,sign1*sign2*Tan1);
781 math_FunctionRoots SolRoot (H,
782 UFirstOnC2,
783 ULastOnC2 ,
784 NbSamples,
785 EpsX,EpsH,EpsH);
786 if (SolRoot.IsDone()) {
787 for (Standard_Integer j = 1; j <= SolRoot.NbSolutions(); j++) {
788 USol = SolRoot.Value(j);
789 gp_Pnt2d P2 = curve2->Value(USol);
790 gp_Vec2d P2P1(P1.X() - P2.X(),P1.Y() - P2.Y());
791 Standard_Real SquareP2P1 = P2P1.SquareMagnitude();
792 Standard_Real N1P2P1 = N1.Dot(P2P1);
793
794 // Test si la solution est du bon cote des courbes.
795 if (N1P2P1*sign1 > 0 ) {
796 P = P1.Translated(- (0.5*SquareP2P1/N1P2P1)*N1);
797 DistPP1 = P1.SquareDistance(P);
798 if (DistPP1 < DMin) {
799 DMin = DistPP1;
800 PSol = P;
801 U2 = USol;
802 YaSol = Standard_True;
803 }
804 }
805 }
806 }
807
808/*
809 if (!YaSol) {
810 //--------------------------------------------------------------------
811 // Construction de la bisectrice point courbe et de la droite passant
812 // par P1 et portee par la normale.
813 //--------------------------------------------------------------------
814 Handle(Bisector_BisecPC) BisPC
815 = new Bisector_BisecPC(curve2,P1,sign2,2*distMax);
816 //-------------------------------
817 // Test si la bissectrice existe.
818 //-------------------------------
819 if (BisPC->IsEmpty()) {
820 Dist = Precision::Infinite();
821 PSol = P1;
822 return PSol;
823 }
824
825 Handle(Geom2d_Line) NorLi = new Geom2d_Line (P1,N1);
826 Geom2dAdaptor_Curve NorLiAd;
827 if (sign1 < 0.) {NorLiAd.Load(NorLi,0. ,distMax);}
828 else {NorLiAd.Load(NorLi,- distMax,0. );}
829
830 //-------------------------------------------------------------------------
831 Geom2dInt_GInter Intersect(BisPC,NorLiAd,
832 Precision::Confusion(),Precision::Confusion());
833 //-------------------------------------------------------------------------
834 if (Intersect.IsDone() && !Intersect.IsEmpty()) {
835 for (Standard_Integer i = 1; i <= Intersect.NbPoints(); i++) {
836 if (Intersect.Point(i).ParamOnSecond()*sign1< Precision::PConfusion()) {
837 P = Intersect.Point(i).Value();
838 DistPP1 = P.SquareDistance(P1);
839 if (DistPP1 < DMin) {
840 DMin = DistPP1;
841 PSol = P;
842 U2 = Intersect.Point(i).ParamOnFirst();
843 YaSol = Standard_True;
844 }
845 }
846 }
847 }
848 }
849*/
850
851 if (YaSol) {
852 Dist = DMin;
853 //--------------------------------------------------------------
854 // Point trouve => Test distance courbure + Test angulaire
855 //---------------------------------------------------------------
856 P2 = curve2->Value(U2);
857 gp_Vec2d PP1(P1.X() - PSol.X(),P1.Y() - PSol.Y());
858 gp_Vec2d PP2(P2.X() - PSol.X(),P2.Y() - PSol.Y());
859
860 //-----------------------------------------------
861 // Dist = produit des normes = distance au carre.
862 //-----------------------------------------------
863 if (PP1.Dot(PP2) > (1. - Precision::Angular())*Dist) {
864 YaSol = Standard_False;
865 }
866 else {
867 if ( !isConvex1 ) {
868 Standard_Real K1 = Curvature(curve1,U1,Precision::Confusion());
869 if (K1 != 0.) {
870 if (Dist > 1/(K1*K1)) YaSol = Standard_False;
871 }
872 }
873 if (YaSol) {
874 if ( !isConvex2 ) {
875 Standard_Real K2 = Curvature(curve2,U2,Precision::Confusion());
876 if (K2 != 0.) {
877 if (Dist > 1/(K2*K2)) YaSol = Standard_False;
878 }
879 }
880 }
881 }
882 }
883 if (!YaSol) {
884 Dist = Precision::Infinite();
885 PSol = P1;
886 }
887 return PSol;
888}
889
890//=============================================================================
891//function : D0
892//purpose :
893//=============================================================================
894void Bisector_BisecCC::D0(const Standard_Real U,
895 gp_Pnt2d& P) const
896{
897 Standard_Real U1,U2,Dist;
898
899 P = ValueAndDist(U,U1,U2,Dist);
900}
901
902//=============================================================================
903//function : D1
904//purpose :
905//=============================================================================
906void Bisector_BisecCC::D1(const Standard_Real U,
907 gp_Pnt2d& P,
908 gp_Vec2d& V ) const
909{
910 V.SetCoord(0.,0.);
911 gp_Vec2d V2,V3;
912 Values(U,1,P,V,V2,V3);
913}
914
915//=============================================================================
916//function : D2
917//purpose :
918//=============================================================================
919void Bisector_BisecCC::D2(const Standard_Real U,
920 gp_Pnt2d& P,
921 gp_Vec2d& V1,
922 gp_Vec2d& V2) const
923{
924 V1.SetCoord(0.,0.);
925 V2.SetCoord(0.,0.);
926 gp_Vec2d V3;
927 Values(U,2,P,V1,V2,V3);
928}
929
930//=============================================================================
931//function : D3
932//purpose :
933//=============================================================================
934void Bisector_BisecCC::D3(const Standard_Real U,
935 gp_Pnt2d& P,
936 gp_Vec2d& V1,
937 gp_Vec2d& V2,
938 gp_Vec2d& V3) const
939{
940 V1.SetCoord(0.,0.);
941 V2.SetCoord(0.,0.);
942 V3.SetCoord(0.,0.);
943 Values(U,3,P,V1,V2,V3);
944}
945
946//=============================================================================
947//function : DN
948//purpose :
949//=============================================================================
950gp_Vec2d Bisector_BisecCC::DN(const Standard_Real U,
951 const Standard_Integer N) const
952{
953 gp_Pnt2d P;
954 gp_Vec2d V1(0.,0.);
955 gp_Vec2d V2(0.,0.);
956 gp_Vec2d V3(0.,0.);
957 Values (U,N,P,V1,V2,V3);
958 switch (N) {
959 case 1 : return V1;
960 case 2 : return V2;
961 case 3 : return V3;
962 default: {
963 Standard_NotImplemented::Raise();
964 }
965 }
966 return V1;
967}
968
969//=============================================================================
970//function : Values
971// purpose : la courbe peut etre decrite par les equations suivantes:
972//
973// B(u) = F(u,v0)
974// ou v0 = Phi(u) est donne par H (u,v) = 0.
975//
976// avec :
977// ||P2(v0)P1(u)||**2
978// F(u,v) = P1(u) - 1/2 *----------------* N(u)
979// (N(u).P2(v0)P1(u))
980//
981// H(u,v) = (Tu.P1(u)P2(v))**2||Tv||**2 - (Tv.P1(u)P2(v))**2||Tu||**2
982//
983// => dB(u)/du = dF/du + dF/dv(- dH/du:dH/dv)
984//
985// Remarque : la tangente a la bisectrice est bissectrice aux
986// tangentes T1(u) et T2(v0)
987//
988//=============================================================================
989void Bisector_BisecCC::Values (const Standard_Real U,
990 const Standard_Integer N,
991 gp_Pnt2d& P,
992 gp_Vec2d& V1,
993 gp_Vec2d& V2,
994 gp_Vec2d& V3) const
995{
996 V1 = gp_Vec2d(0.,0.);
997 V2 = gp_Vec2d(0.,0.);
998 V3 = gp_Vec2d(0.,0.);
999 //-------------------------------------------------------------------------
1000 // Calcul du point courant sur la bisectrice et des parametres sur chaque
1001 // courbe.
1002 //-------------------------------------------------------------------------
1003 Standard_Real U0,V0,Dist;
1004
1005 //-----------------------------------------------
1006 // le polygone est il reduit a un point ou vide?
1007 //-----------------------------------------------
1008 if (myPolygon.Length() <= 1) {
1009 P = Extension(U,U0,V0,Dist,V1);
1010 }
1011 if (U < myPolygon.First().ParamOnBis()) {
1012 P = Extension(U,U0,V0,Dist,V1);
1013 return;
1014 }
1015 if (U > myPolygon.Last().ParamOnBis()) {
1016 P = Extension(U,U0,V0,Dist,V1);
1017 return;
1018 }
1019 P = ValueAndDist(U,U0,V0,Dist);
1020
1021 if (N == 0) return;
1022 //------------------------------------------------------------------
1023 // Recuperation des point,tangente,normale sur C1 au parametre U0.
1024 //-------------------------------------------------------------------
1025 gp_Pnt2d P1 ; // point sur C1.
1026 gp_Vec2d Tu ; // tangente a C1 en U0.
1027 gp_Vec2d Tuu ; // derivee seconde a C1 en U0.
1028 curve1->D2(U0,P1,Tu,Tuu);
1029 gp_Vec2d Nor( - Tu .Y() , Tu .X()); // Normale en U0.
1030 gp_Vec2d Nu ( - Tuu.Y() , Tuu.X()); // derivee de la normale en U0.
1031
1032 //-------------------------------------------------------------------
1033 // Recuperation des point,tangente,normale sur C2 au parametre V0.
1034 //-------------------------------------------------------------------
1035 gp_Pnt2d P2 ; // point sur C2.
1036 gp_Vec2d Tv ; // tangente a C2 en V.
1037 gp_Vec2d Tvv ; // derivee seconde a C2 en V.
1038 curve2->D2(V0,P2,Tv,Tvv);
1039
1040 gp_Vec2d PuPv(P2.X() - P1.X(), P2.Y() - P1.Y());
1041
1042 //-----------------------------
1043 // Calcul de dH/du et de dH/dv.
1044 //-----------------------------
1045 Standard_Real TuTu,TvTv,TuTv;
1046 Standard_Real TuPuPv,TvPuPv ;
1047 Standard_Real TuuPuPv,TuTuu ;
1048 Standard_Real TvvPuPv,TvTvv ;
1049
1050 TuTu = Tu.Dot(Tu) ; TvTv = Tv.Dot(Tv) ; TuTv = Tu.Dot(Tv);
1051 TuPuPv = Tu.Dot(PuPv) ; TvPuPv = Tv.Dot(PuPv);
1052 TuuPuPv = Tuu.Dot(PuPv) ; TuTuu = Tu.Dot(Tuu) ;
1053 TvvPuPv = Tvv.Dot(PuPv) ; TvTvv = Tv.Dot(Tvv) ;
1054
1055 Standard_Real dHdu = 2*(TuPuPv*(TuuPuPv - TuTu)*TvTv +
1056 TvPuPv*TuTv*TuTu -TuTuu*TvPuPv*TvPuPv);
1057 Standard_Real dHdv = 2*(TuPuPv*TuTv*TvTv + TvTvv*TuPuPv*TuPuPv -
1058 TvPuPv*(TvvPuPv + TvTv)*TuTu);
1059
1060 //-----------------------------
1061 // Calcul de dF/du et de dF/dv.
1062 //-----------------------------
1063 Standard_Real NorPuPv,NuPuPv,NorTv;
1064 Standard_Real A,B,dAdu,dAdv,dBdu,dBdv,BB;
1065
1066 NorPuPv = Nor.Dot(PuPv);
1067 NuPuPv = Nu .Dot(PuPv);
1068 NorTv = Nor.Dot(Tv) ;
1069
1070 A = 0.5*PuPv.SquareMagnitude();
1071 B = - NorPuPv;
1072 BB = B*B;
1073 dAdu = - TuPuPv;
1074 dBdu = - NuPuPv ;
1075 dAdv = TvPuPv;
1076 dBdv = - NorTv;
1077
1078 //---------------------------------------
1079 // F(u,v) = Pu - (A(u,v)/B(u,v))*Nor(u)
1080 //----------------------------------------
1081 if (BB < gp::Resolution()) {
1082 V1 = Tu.Normalized() + Tv.Normalized();
1083 V1 = 0.5*Tu.SquareMagnitude()*V1;
1084 }
1085 else {
1086 gp_Vec2d dFdu = Tu - (dAdu/B - dBdu*A/BB)*Nor - (A/B)*Nu;
1087 gp_Vec2d dFdv = ( - dAdv/B + dBdv*A/BB)*Nor ;
1088
1089 if (Abs(dHdv) > gp::Resolution()) {
1090 V1 = dFdu + dFdv*( - dHdu / dHdv );
1091 }
1092 else {
1093 V1 = Tu;
1094 }
1095 }
1096 if (N == 1) return;
1097}
1098
1099//=============================================================================
1100//function : Extension
1101// purpose : Calcul du point courant sur les extensions ou prolongement en
1102// tangence de la courbe.
1103//============================================================================
1104gp_Pnt2d Bisector_BisecCC::Extension (const Standard_Real U,
1105 Standard_Real& U1,
1106 Standard_Real& U2,
1107 Standard_Real& Dist,
1108 gp_Vec2d& T ) const
1109{
1110 Bisector_PointOnBis PRef;
1111 gp_Pnt2d P,P1,P2,PBis;
1112 gp_Vec2d T1,Tang;
1113#ifndef DEB
1114 Standard_Real dU = 0.;
1115#else
1116 Standard_Real dU;
1117#endif
1118 Standard_Boolean ExtensionTangent = Standard_False;
1119
1120 if (myPolygon.Length() == 0) {
1121 //---------------------------------------------
1122 // Polygone vide => segment (pointStart,pointEnd)
1123 //---------------------------------------------
1124 dU = U - startIntervals.First();
1125 P = pointStart;
1126 P1 = pointEnd;
1127 U1 = curve1->LastParameter();
1128 if (sign1 == sign2) { U2 = curve2->FirstParameter();}
1129 else { U2 = curve2->LastParameter() ;}
1130 Tang.SetCoord(P1.X() - P.X(),P1.Y() - P.Y());
1131 }
1132 else if (U < myPolygon.First().ParamOnBis()) {
1133 PRef = myPolygon.First();
1134 P = PRef.Point();
1135 dU = U - PRef.ParamOnBis();
1136 if (extensionStart) {
1137 //------------------------------------------------------------
1138 // extension = segment (pointstart,premier point du polygone.)
1139 //------------------------------------------------------------
1140 P1 = pointStart;
1141 U1 = curve1->FirstParameter();
1142 if (sign1 == sign2) { U2 = curve2->LastParameter();}
1143 else { U2 = curve2->FirstParameter();}
1144 Tang.SetCoord(P.X() - P1.X(),P.Y() - P1.Y());
1145 }
1146 else {
1147 ExtensionTangent = Standard_True;
1148 }
1149 }
1150 else if (U > myPolygon.Last().ParamOnBis()) {
1151 PRef = myPolygon.Last();
1152 P = PRef.Point();
1153 dU = U - PRef.ParamOnBis();
1154 if (extensionEnd) {
1155 //------------------------------------------------------------
1156 // extension = segment (dernier point du polygone.pointEnd)
1157 //------------------------------------------------------------
1158 P1 = pointEnd;
1159 U1 = curve1->LastParameter();
1160 if (sign1 == sign2) { U2 = curve2->LastParameter();}
1161 else { U2 = curve2->FirstParameter();}
1162 Tang.SetCoord(P1.X() - P.X(),P1.Y() - P.Y());
1163 }
1164 else {
1165 ExtensionTangent = Standard_True;
1166 }
1167 }
1168
1169 if (ExtensionTangent) {
1170 //-----------------------------------------------------------
1171 // Si la courbe n a pas d extension, celle ci est prolonge
1172 // en tangence.
1173 //------------------------------------------------------------
1174 U1 = PRef.ParamOnC1();
1175 U2 = PRef.ParamOnC2();
1176 P2 = curve2->Value(U2);
1177 curve1->D1(U1,P1,T1);
1178 Tang.SetCoord(2*P.X() - P1.X() - P2.X(), 2*P.Y() - P1.Y() - P2.Y());
1179 if (Tang.Magnitude() < Precision::Confusion()) {
1180 Tang = T1;
1181 }
1182 if (T1.Dot(Tang) < 0.) Tang = - Tang;
1183 }
1184
1185 T = Tang.Normalized();
1186 PBis.SetCoord(P.X() + dU*T.X(),P.Y() + dU*T.Y());
1187 Dist = P1.Distance(PBis);
1188 return PBis;
1189}
1190
1191//=============================================================================
1192//function : PointByInt
1193// purpose :
1194//=============================================================================
1195static Standard_Boolean PointByInt(const Handle(Geom2d_Curve)& CA,
1196 const Handle(Geom2d_Curve)& CB,
1197 const Standard_Real SignA,
1198 const Standard_Real SignB,
1199 const Standard_Real UOnA,
1200 Standard_Real& UOnB,
1201 Standard_Real& Dist)
1202{
1203 //------------------------------------------------------------------
1204 // Recuperation des point,tangente,normale sur CA au parametre UOnA.
1205 //-------------------------------------------------------------------
1206 gp_Pnt2d P1,P2,P,PSol;
1207 gp_Vec2d Tan1,Tan2;
1208 Standard_Boolean IsConvexA = Bisector::IsConvex(CA,SignA);
1209 Standard_Boolean IsConvexB = Bisector::IsConvex(CB,SignB);
1210
1211 CA->D1(UOnA,P1,Tan1);
1212 gp_Vec2d N1(Tan1.Y(), - Tan1.X());
1213
1214 //--------------------------------------------------------------------------
1215 // test de confusion de P1 avec extremite de curve2.
1216 //--------------------------------------------------------------------------
1217 if (P1.Distance(CB->Value(CB->FirstParameter())) < Precision::Confusion()) {
1218 UOnB = CB->FirstParameter();
1219 CB->D1(UOnB,P2,Tan2);
1220 if ( IsConvexA && IsConvexB ) {
1221 Dist = 0.;
1222 return Standard_True;
1223 }
1224 if (! Tan1.IsParallel(Tan2,Precision::Angular())) {
1225 Dist = 0.;
1226 return Standard_False;
1227 }
1228 }
1229 if (P1.Distance(CB->Value(CB->LastParameter())) < Precision::Confusion()) {
1230 UOnB = CB->LastParameter();
1231 CB->D1(UOnB,P2,Tan2);
1232 if ( IsConvexA && IsConvexB ) {
1233 Dist = 0.;
1234 return Standard_True;
1235 }
1236 if (! Tan1.IsParallel(Tan2,Precision::Angular())) {
1237 Dist = 0.;
1238 return Standard_False;
1239 }
1240 }
1241
1242 Standard_Real DMin = Precision::Infinite();
1243 Standard_Real UPC;
1244 Standard_Boolean YaSol = Standard_False;
1245 //--------------------------------------------------------------------
1246 // Construction de la bisectrice point courbe et de la droite passant
1247 // par P1 et portee par la normale.
1248 //--------------------------------------------------------------------
1249 Handle(Bisector_BisecPC) BisPC
1250 = new Bisector_BisecPC(CB,P1,SignB );
1251 //-------------------------------
1252 // Test si la bissectrice existe.
1253 //-------------------------------
1254 if (BisPC->IsEmpty()) {
1255 Dist = Precision::Infinite();
1256 PSol = P1;
1257 return Standard_False;
1258 }
1259
1260 Handle(Geom2d_Line) NorLi = new Geom2d_Line (P1,N1);
1261
1262 Geom2dAdaptor_Curve ABisPC(BisPC);
1263 Geom2dAdaptor_Curve ANorLi(NorLi);
1264 //-------------------------------------------------------------------------
1265 Geom2dInt_GInter Intersect(ABisPC,ANorLi,
1266 Precision::Confusion(),Precision::Confusion());
1267 //-------------------------------------------------------------------------
1268
1269 if (Intersect.IsDone() && !Intersect.IsEmpty()) {
1270 for (Standard_Integer i = 1; i <= Intersect.NbPoints(); i++) {
1271 if (Intersect.Point(i).ParamOnSecond()*SignA < Precision::PConfusion()) {
1272 P = Intersect.Point(i).Value();
1273 if (P.SquareDistance(P1) < DMin) {
1274 DMin = P.SquareDistance(P1);
1275 PSol = P;
1276 UPC = Intersect.Point(i).ParamOnFirst();
1277 UOnB = BisPC->LinkBisCurve(UPC);
1278 Dist = DMin;
1279 YaSol = Standard_True;
1280 }
1281 }
1282 }
1283 }
1284 if (YaSol) {
1285 //--------------------------------------------------------------
1286 // Point trouve => Test distance courbure + Test angulaire
1287 //---------------------------------------------------------------
1288 P2 = CB->Value(UOnB);
1289 gp_Dir2d PP1Unit(P1.X() - PSol.X(),P1.Y() - PSol.Y());
1290 gp_Dir2d PP2Unit(P2.X() - PSol.X(),P2.Y() - PSol.Y());
1291
1292 if (PP1Unit*PP2Unit > 1. - Precision::Angular()) {
1293 YaSol = Standard_False;
1294 }
1295 else {
1296 Dist = sqrt(Dist);
1297 if ( !IsConvexA ) {
1298 Standard_Real K1 = Curvature(CA,UOnA,Precision::Confusion());
1299 if (K1 != 0.) {
1300 if (Dist > Abs(1/K1)) YaSol = Standard_False;
1301 }
1302 }
1303 if (YaSol) {
1304 if ( !IsConvexB ) {
1305 Standard_Real K2 = Curvature(CB,UOnB,Precision::Confusion());
1306 if (K2 != 0.) {
1307 if (Dist > Abs(1/K2)) YaSol = Standard_False;
1308 }
1309 }
1310 }
1311 }
1312 }
1313 return YaSol;
1314}
1315
1316//=============================================================================
1317//function : SupLastParameter
1318// purpose :
1319//=============================================================================
1320void Bisector_BisecCC::SupLastParameter()
1321{
1322 endIntervals.Append(curve1->LastParameter());
1323 // ----------------------------------------------------------------------
1324 // Calcul du parametre sur curve1 associees a l une ou lautre des extremites
1325 // de curve2 suivant les valeurs de sign1 et sign2.
1326 // la bissectrice est restreinte par les parametres obtenus.
1327 //------------------------------------------------------------------------
1328 Standard_Real UOnC1,UOnC2,Dist;
1329 if (sign1 == sign2) {
1330 UOnC2 = curve2->FirstParameter();
1331 }
1332 else {
1333 UOnC2 = curve2->LastParameter();
1334 }
1335 Standard_Boolean YaSol = PointByInt(curve2,curve1,sign2,sign1,UOnC2,UOnC1,Dist);
1336 if (YaSol) {
1337 if (UOnC1 > startIntervals.First() && UOnC1 < endIntervals.Last()) {
1338 endIntervals.SetValue(1,UOnC1);
1339 }
1340 }
1341}
1342
1343//=============================================================================
1344//function : Curve
1345// purpose :
1346//=============================================================================
1347Handle(Geom2d_Curve) Bisector_BisecCC::Curve(const Standard_Integer I) const
1348{
1349 if (I == 1) return curve1;
1350 else if (I == 2) return curve2;
1351 else Standard_OutOfRange::Raise();
1352 return curve1;
1353}
1354
1355//=============================================================================
1356//function : LinkBisCurve
1357//purpose :
1358//=============================================================================
1359Standard_Real Bisector_BisecCC::LinkBisCurve(const Standard_Real U) const
1360{
1361 return (U - shiftParameter);
1362}
1363
1364//=============================================================================
1365//function : LinkCurveBis
1366//purpose :
1367//=============================================================================
1368Standard_Real Bisector_BisecCC::LinkCurveBis(const Standard_Real U) const
1369{
1370 return (U + shiftParameter);
1371}
1372
1373//=============================================================================
1374//function : Indent
1375//purpose :
1376//=============================================================================
1377static void Indent(const Standard_Integer Offset) {
1378 if (Offset > 0) {
1379 for (Standard_Integer i = 0; i < Offset; i++) {cout << " ";}
1380 }
1381}
1382
1383//=============================================================================
1384//function : Polygon
1385// purpose :
1386//=============================================================================
1387const Bisector_PolyBis& Bisector_BisecCC::Polygon() const
1388{
1389 return myPolygon;
1390}
1391
1392//==========================================================================
1393//function : Parameter
1394//purpose :
1395//==========================================================================
1396Standard_Real Bisector_BisecCC::Parameter(const gp_Pnt2d& P) const
1397{
1398 Standard_Real UOnCurve;
1399
1400 if (P.IsEqual(Value(FirstParameter()),Precision::Confusion())) {
1401 UOnCurve = FirstParameter();
1402 }
1403 else if (P.IsEqual(Value(LastParameter()),Precision::Confusion())) {
1404 UOnCurve = LastParameter();
1405 }
1406 else {
1407 UOnCurve = ProjOnCurve(P,curve1);
1408 }
1409 return UOnCurve;
1410}
1411
1412
1413//=============================================================================
1414//function : Dump
1415// purpose :
1416//=============================================================================
1417//void Bisector_BisecCC::Dump(const Standard_Integer Deep,
1418void Bisector_BisecCC::Dump(const Standard_Integer ,
1419 const Standard_Integer Offset) const
1420{
1421 Indent (Offset);
1422 cout <<"Bisector_BisecCC :"<<endl;
1423 Indent (Offset);
1424// cout <<"Curve1 :"<<curve1<<endl;
1425// cout <<"Curve2 :"<<curve2<<endl;
1426 cout <<"Sign1 :"<<sign1<<endl;
1427 cout <<"Sign2 :"<<sign2<<endl;
1428
1429 cout <<"Number Of Intervals :"<<startIntervals.Length()<<endl;
1430 for (Standard_Integer i = 1; i <= startIntervals.Length(); i++) {
1431 cout <<"Interval number :"<<i<<"Start :"<<startIntervals.Value(i)
1432 <<" end :"<< endIntervals.Value(i)<<endl ;
1433 }
1434 cout <<"Index Current Interval :"<<currentInterval<<endl;
1435}
1436
1437//=============================================================================
1438//function : Curve
1439// purpose :
1440//=============================================================================
1441void Bisector_BisecCC::Curve(const Standard_Integer I,
1442 const Handle(Geom2d_Curve)& C)
1443{
1444 if (I == 1) curve1 = C;
1445 else if (I == 2) curve2 = C;
1446 else Standard_OutOfRange::Raise();
1447}
1448
1449//=============================================================================
1450//function : Sign
1451// purpose :
1452//=============================================================================
1453void Bisector_BisecCC::Sign(const Standard_Integer I,
1454 const Standard_Real S)
1455{
1456 if (I == 1) sign1 = S;
1457 else if (I == 2) sign2 = S;
1458 else Standard_OutOfRange::Raise();
1459}
1460
1461//=============================================================================
1462//function : Polygon
1463// purpose :
1464//=============================================================================
1465void Bisector_BisecCC::Polygon(const Bisector_PolyBis& P)
1466{
1467 myPolygon = P;
1468}
1469
1470//=============================================================================
1471//function : DistMax
1472// purpose :
1473//=============================================================================
1474void Bisector_BisecCC::DistMax(const Standard_Real D)
1475{
1476 distMax = D;
1477}
1478
1479//=============================================================================
1480//function : IsConvex
1481// purpose :
1482//=============================================================================
1483void Bisector_BisecCC::IsConvex(const Standard_Integer I,
1484 const Standard_Boolean IsConvex)
1485{
1486 if (I == 1) isConvex1 = IsConvex;
1487 else if (I == 2) isConvex2 = IsConvex;
1488 else Standard_OutOfRange::Raise();
1489}
1490
1491//=============================================================================
1492//function : IsEmpty
1493// purpose :
1494//=============================================================================
1495void Bisector_BisecCC::IsEmpty ( const Standard_Boolean IsEmpty)
1496{
1497 isEmpty = IsEmpty;
1498}
1499
1500//=============================================================================
1501//function : ExtensionStart
1502// purpose :
1503//=============================================================================
1504void Bisector_BisecCC::ExtensionStart( const Standard_Boolean ExtensionStart)
1505{
1506 extensionStart = ExtensionStart;
1507}
1508
1509//=============================================================================
1510//function : ExtensionEnd
1511// purpose :
1512//=============================================================================
1513void Bisector_BisecCC::ExtensionEnd( const Standard_Boolean ExtensionEnd)
1514{
1515 extensionEnd = ExtensionEnd;
1516}
1517
1518//=============================================================================
1519//function : PointStart
1520// purpose :
1521//=============================================================================
1522void Bisector_BisecCC::PointStart( const gp_Pnt2d& Point)
1523{
1524 pointStart = Point;
1525}
1526
1527//=============================================================================
1528//function : PointEnd
1529// purpose :
1530//=============================================================================
1531void Bisector_BisecCC::PointEnd( const gp_Pnt2d& Point)
1532{
1533 pointEnd = Point;
1534}
1535
1536//=============================================================================
1537//function : StartIntervals
1538// purpose :
1539//=============================================================================
1540void Bisector_BisecCC::StartIntervals
1541 (const TColStd_SequenceOfReal& StartIntervals)
1542{
1543 startIntervals = StartIntervals;
1544}
1545
1546//=============================================================================
1547//function : EndIntervals
1548// purpose :
1549//=============================================================================
1550void Bisector_BisecCC::EndIntervals
1551 (const TColStd_SequenceOfReal& EndIntervals)
1552{
1553 endIntervals = EndIntervals;
1554}
1555
1556//=============================================================================
1557//function : FirstParameter
1558// purpose :
1559//=============================================================================
1560void Bisector_BisecCC::FirstParameter (const Standard_Real U)
1561{
1562 startIntervals.Append(U);
1563}
1564
1565//=============================================================================
1566//function : LastParameter
1567// purpose :
1568//=============================================================================
1569void Bisector_BisecCC::LastParameter (const Standard_Real U)
1570{
1571 endIntervals.Append(U);
1572}
1573
1574//=============================================================================
1575//function : SearchBound
1576// purpose :
1577//=============================================================================
1578Standard_Real Bisector_BisecCC::SearchBound (const Standard_Real U1,
1579 const Standard_Real U2) const
1580{
1581 Standard_Real UMid,Dist1,Dist2,DistMid,U11,U22;
1582 Standard_Real UC1,UC2;
1583 gp_Pnt2d PBis,PBisPrec;
1584 Standard_Real TolPnt = Precision::Confusion();
1585 Standard_Real TolPar = Precision::PConfusion();
1586 U11 = U1; U22 = U2;
1587 PBisPrec = ValueByInt(U11,UC1,UC2,Dist1);
1588 PBis = ValueByInt(U22,UC1,UC2,Dist2);
1589
1590 while ((U22 - U11) > TolPar ||
1591 ((Dist1 < Precision::Infinite() &&
1592 Dist2 < Precision::Infinite() &&
1593 !PBis.IsEqual(PBisPrec,TolPnt)))) {
1594 PBisPrec = PBis;
1595 UMid = 0.5*( U22 + U11);
1596 PBis = ValueByInt(UMid,UC1,UC2,DistMid);
1597 if ((Dist1 < Precision::Infinite()) == (DistMid < Precision::Infinite())) {
1598 U11 = UMid;
1599 Dist1 = DistMid;
1600 }
1601 else {
1602 U22 = UMid;
1603 Dist2 = DistMid;
1604 }
1605 }
1606 PBis = ValueByInt(U11,UC1,UC2,Dist1);
1607 if (Dist1 < Precision::Infinite()) {
1608 UMid = U11;
1609 }
1610 else {
1611 UMid = U22;
1612 }
1613 return UMid;
1614}
1615
1616//=============================================================================
1617//function : ProjOnCurve
1618// purpose :
1619//=============================================================================
1620static Standard_Real ProjOnCurve (const gp_Pnt2d& P,
1621 const Handle(Geom2d_Curve)& C)
1622{
1623#ifndef DEB
1624 Standard_Real UOnCurve =0.;
1625#else
1626 Standard_Real UOnCurve;
1627#endif
1628 gp_Pnt2d PF,PL;
1629 gp_Vec2d TF,TL;
1630
1631 C->D1(C->FirstParameter(),PF,TF);
1632 C->D1(C->LastParameter() ,PL,TL);
1633
1634 if (P.IsEqual(PF ,Precision::Confusion())) {
1635 return C->FirstParameter();
1636 }
1637 if (P.IsEqual(PL ,Precision::Confusion())) {
1638 return C->LastParameter();
1639 }
1640 gp_Vec2d PPF(PF.X() - P.X(), PF.Y() - P.Y());
1641 TF.Normalize();
1642 if ( Abs (PPF.Dot(TF)) < Precision::Confusion()) {
1643 return C->FirstParameter();
1644 }
1645 gp_Vec2d PPL (PL.X() - P.X(), PL.Y() - P.Y());
1646 TL.Normalize();
1647 if ( Abs (PPL.Dot(TL)) < Precision::Confusion()) {
1648 return C->LastParameter();
1649 }
1650 Geom2dAPI_ProjectPointOnCurve Proj(P,C,
1651 C->FirstParameter(),
1652 C->LastParameter());
1653 if (Proj.NbPoints() > 0) {
1654 UOnCurve = Proj.LowerDistanceParameter();
1655 }
1656 else {
1657 Standard_OutOfRange::Raise();
1658 }
1659 return UOnCurve;
1660}
1661
1662//=============================================================================
1663//function : TestExtension
1664// purpose :
1665//=============================================================================
1666static Standard_Boolean TestExtension (const Handle(Geom2d_Curve)& C1,
1667 const Handle(Geom2d_Curve)& C2,
1668 const Standard_Integer Start_End)
1669{
1670 gp_Pnt2d P1,P2;
1671 gp_Vec2d T1,T2;
1672 Standard_Boolean Test = Standard_False;
1673 if (Start_End == 1) {
1674 C1->D1(C1->FirstParameter(),P1,T1);
1675 }
1676 else {
1677 C1->D1(C1->LastParameter(),P1,T1);
1678 }
1679 C2->D1(C2->FirstParameter(),P2,T2);
1680 if (P1.IsEqual(P2,Precision::Confusion())) {
1681 T1.Normalize(); T2.Normalize();
1682 if (T1.Dot(T2) > 1.0 - Precision::Confusion()) {
1683 Test = Standard_True;
1684 }
1685 }
1686 else {
1687 C2->D1(C2->LastParameter(),P2,T2);
1688 if (P1.IsEqual(P2,Precision::Confusion())) {
1689 T2.Normalize();
1690 if (T1.Dot(T2) > 1.0 - Precision::Confusion()) {
1691 Test = Standard_True;
1692 }
1693 }
1694 }
1695 return Test;
1696}
1697
1698//=============================================================================
1699//function : ComputePointEnd
1700// purpose :
1701//=============================================================================
1702void Bisector_BisecCC::ComputePointEnd ()
1703{
1704 Standard_Real U1,U2;
1705 Standard_Real KC,RC;
1706 U1 = curve1->FirstParameter();
1707 if (sign1 == sign2) {
1708 U2 = curve2->LastParameter();
1709 }
1710 else {
1711 U2 = curve2->FirstParameter();
1712 }
1713 Standard_Real K1 = Curvature(curve1,U1,Precision::Confusion());
1714 Standard_Real K2 = Curvature(curve2,U2,Precision::Confusion());
1715 if (!isConvex1 && !isConvex2) {
1716 if (K1 < K2) {KC = K1;} else {KC = K2;}
1717 }
1718 else if (!isConvex1) {KC = K1;}
1719 else {KC = K2;}
1720
1721 gp_Pnt2d PF;
1722 gp_Vec2d TF;
1723 curve1->D1(U1,PF,TF);
1724 TF.Normalize();
1725 if (KC != 0.) { RC = Abs(1/KC);}
1726 else { RC = Precision::Infinite();}
1727 pointEnd.SetCoord(PF.X() - sign1*RC*TF.Y(), PF.Y() + sign1*RC*TF.X());
1728
1729}
1730
1731//=============================================================================
1732//function : DiscretPar
1733// purpose :
1734//=============================================================================
1735static Standard_Boolean DiscretPar(const Standard_Real DU,
1736 const Standard_Real EpsMin,
1737 const Standard_Real EpsMax,
1738 const Standard_Integer NbMin,
1739 const Standard_Integer NbMax,
1740 Standard_Real& Eps,
1741 Standard_Integer& Nb)
1742{
1743 if (DU <= NbMin*EpsMin) {
1744 Eps = DU/(NbMin + 1) ;
1745 Nb = NbMin;
1746 return Standard_False;
1747 }
1748
1749 Eps = Min (EpsMax,DU/NbMax);
1750
1751 if (Eps < EpsMin) {
1752 Eps = EpsMin;
1753 Nb = Standard_Integer(DU/EpsMin);
1754 }
1755 else { Nb = NbMax;}
1756
1757 return Standard_True;
1758}
1759
1760