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