Replacing french comments by english one
[occt.git] / src / BRepBlend / BRepBlend_RstRstConstRad.cxx
1 // File:        BRepBlend_RstRstConstRad.cxx
2 // Created:     Mon Feb 10 10:32:10 1997
3 // Author:      Laurent BOURESCHE
4 // Author:      Jacques GOUSSARD
5 //              <lbo@pomalox.paris1.matra-dtv.fr>
6 // Modif : jlr le 28/07/97 modif de section pour Edge/Face
7 //         modif de set (courbe bornee)
8 #include <BRepBlend_RstRstConstRad.ixx>
9 #include <math_Gauss.hxx>
10 #include <math_SVD.hxx>
11
12 #include <ElCLib.hxx>
13 #include <gp.hxx>
14 #include <BlendFunc.hxx>
15 #include <GeomFill.hxx>
16 #include <Standard_DomainError.hxx>
17 #include <Standard_NotImplemented.hxx>
18 #include <Precision.hxx>
19
20 #define Eps 1.e-15
21
22 static void t3dto2d(Standard_Real& a,
23                     Standard_Real& b,
24                     const gp_Vec& A,
25                     const gp_Vec& B,
26                     const gp_Vec& C)
27 {
28   Standard_Real AB = A.Dot(B);
29   Standard_Real AC = A.Dot(C);
30   Standard_Real BC = B.Dot(C);
31   Standard_Real BB = B.Dot(B);
32   Standard_Real CC = C.Dot(C);
33   Standard_Real deno = (BB*CC-BC*BC);
34   a = (AB*CC-AC*BC)/deno;
35   b = (AC*BB-AB*BC)/deno;
36 }
37
38 //=======================================================================
39 //function : BRepBlend_RstRstConstRad
40 //purpose  : 
41 //=======================================================================
42
43 BRepBlend_RstRstConstRad::BRepBlend_RstRstConstRad
44 (const Handle(Adaptor3d_HSurface)& Surf1,
45  const Handle(Adaptor2d_HCurve2d)& Rst1,
46  const Handle(Adaptor3d_HSurface)& Surf2,
47  const Handle(Adaptor2d_HCurve2d)& Rst2,
48  const Handle(Adaptor3d_HCurve)&   CGuide):
49  surf1(Surf1), surf2(Surf2), rst1(Rst1), rst2(Rst2),
50  cons1(Rst1, Surf1), cons2(Rst2, Surf2), 
51  guide(CGuide), tguide(CGuide),
52  istangent(Standard_True), maxang(RealFirst()), minang(RealLast()),
53  distmin(RealLast()),
54  mySShape(BlendFunc_Rational)
55 {}
56
57 //=======================================================================
58 //function : NbVariables
59 //purpose  : 
60 //=======================================================================
61
62 Standard_Integer BRepBlend_RstRstConstRad::NbVariables() const
63 {
64   return 2;
65 }
66
67 //=======================================================================
68 //function : NbEquations
69 //purpose  : 
70 //=======================================================================
71
72 Standard_Integer BRepBlend_RstRstConstRad::NbEquations() const
73 {
74   return 2;
75 }
76
77 //=======================================================================
78 //function : Value
79 //purpose  : 
80 //=======================================================================
81
82 Standard_Boolean BRepBlend_RstRstConstRad::Value(const math_Vector& X,
83                                                  math_Vector&       F)
84 {
85   ptrst1 = cons1.Value(X(1));
86   ptrst2 = cons2.Value(X(2));
87   
88   F(1)   = nplan.XYZ().Dot(ptrst1.XYZ()) + theD;  
89   F(2)   = nplan.XYZ().Dot(ptrst2.XYZ()) + theD;
90   
91   return Standard_True;
92 }
93
94 //=======================================================================
95 //function : Derivatives
96 //purpose  : 
97 //=======================================================================
98
99 Standard_Boolean BRepBlend_RstRstConstRad::Derivatives(const math_Vector& X,
100                                                        math_Matrix&       D)
101 {
102   gp_Vec d11, d21;
103
104   cons1.D1(X(1), ptrst1, d11);
105   cons2.D1(X(2), ptrst2, d21);
106   
107   D(1,1) = nplan.Dot(d11);
108   D(1,2) = 0.;
109   
110   D(2,1) = 0.;
111   D(2,2) = nplan.Dot(d21);
112     
113   return Standard_True;
114 }
115
116 //=======================================================================
117 //function : Values
118 //purpose  : 
119 //=======================================================================
120
121 Standard_Boolean BRepBlend_RstRstConstRad::Values(const math_Vector& X,
122                                                   math_Vector&       F,
123                                                   math_Matrix&       D)
124 {
125   Standard_Boolean Error;
126
127   Error = Value(X, F);  
128   Error = Derivatives(X, D);
129  
130   return Standard_True;
131 }
132
133 //=======================================================================
134 //function : Set
135 //purpose  : 
136 //=======================================================================
137
138 void BRepBlend_RstRstConstRad::Set(const Handle(Adaptor3d_HSurface)& SurfRef1,
139                                    const Handle(Adaptor2d_HCurve2d)& RstRef1,
140                                    const Handle(Adaptor3d_HSurface)& SurfRef2,
141                                    const Handle(Adaptor2d_HCurve2d)& RstRef2)
142 {
143   surfref1 = SurfRef1;
144   surfref2 = SurfRef2;
145   rstref1  = RstRef1;
146   rstref2  = RstRef2;
147 }
148
149 //=======================================================================
150 //function : Set
151 //purpose  : 
152 //=======================================================================
153
154 void BRepBlend_RstRstConstRad::Set(const Standard_Real Param)
155 {
156   tguide->D2(Param, ptgui, d1gui, d2gui);
157   normtg = d1gui.Magnitude();
158   nplan  = d1gui.Normalized();
159   theD   = - (nplan.XYZ().Dot(ptgui.XYZ()));
160 }
161
162 //=======================================================================
163 //function : Set
164 //purpose  : 
165 //=======================================================================
166
167 void BRepBlend_RstRstConstRad::Set(const Standard_Real First,
168                                    const Standard_Real Last)
169
170  tguide = guide->Trim(First, Last, 1.e-12);
171 }
172
173 //=======================================================================
174 //function : GetTolerance
175 //purpose  : 
176 //=======================================================================
177
178 void BRepBlend_RstRstConstRad::GetTolerance(math_Vector&        Tolerance,
179                                             const Standard_Real Tol) const
180 {
181   Tolerance(1) = cons1.Resolution(Tol);
182   Tolerance(2) = cons2.Resolution(Tol);
183 }
184
185 //=======================================================================
186 //function : GetBounds
187 //purpose  : 
188 //=======================================================================
189
190 void BRepBlend_RstRstConstRad::GetBounds(math_Vector& InfBound,
191                                          math_Vector& SupBound) const
192 {
193   InfBound(1) = cons1.FirstParameter();
194   InfBound(2) = cons2.FirstParameter();
195   SupBound(1) = cons1.LastParameter();
196   SupBound(2) = cons2.LastParameter();
197   
198 }
199
200 //=======================================================================
201 //function : IsSolution
202 //purpose  : 
203 //=======================================================================
204
205 Standard_Boolean BRepBlend_RstRstConstRad::IsSolution(const math_Vector&  Sol,
206                                                       const Standard_Real Tol)
207      
208      
209 {
210   math_Vector valsol(1, 2), secmember(1, 2);
211   math_Matrix gradsol(1, 2, 1, 2);
212   
213   gp_Vec dnplan, d1urst1, d1vrst1, d1urst2, d1vrst2, d11, d21, temp;
214   gp_Pnt bid;
215
216   Standard_Real Cosa, Sina, Angle;
217   
218   Values(Sol, valsol, gradsol);
219
220   if (Abs(valsol(1)) <= Tol &&
221       Abs(valsol(2)) <= Tol ) {
222     
223     // Calculation of tangents
224     prmrst1  = Sol(1);    
225     pt2drst1 = rst1->Value(prmrst1);
226     prmrst2  = Sol(2);
227     pt2drst2 = rst2->Value(prmrst2);
228
229     cons1.D1(Sol(1), ptrst1, d11);
230     cons2.D1(Sol(2), ptrst2, d21);
231
232     dnplan.SetLinearForm(1./normtg, d2gui,
233                          -1./normtg * (nplan.Dot(d2gui)), nplan);
234     
235     temp.SetXYZ(ptrst1.XYZ() - ptgui.XYZ());
236     secmember(1) = normtg - dnplan.Dot(temp);
237     
238     temp.SetXYZ(ptrst2.XYZ() - ptgui.XYZ());
239     secmember(2) = normtg - dnplan.Dot(temp);
240     
241     math_Gauss Resol(gradsol);
242
243     if (Resol.IsDone()) {    
244       Resol.Solve(secmember);
245       istangent = Standard_False;
246     }
247     else {
248       math_SVD SingRS (gradsol);
249       if (SingRS.IsDone()) {
250         math_Vector DEDT(1,3);
251         DEDT = secmember;
252         SingRS.Solve(DEDT, secmember, 1.e-6);
253         istangent = Standard_False;
254       }
255       else istangent = Standard_True;
256     }
257
258
259     if (!istangent) {      
260       tgrst1 = secmember(1) * d11;
261       tgrst2 = secmember(2) * d21;
262
263       Standard_Real a, b;
264       surf1->D1(pt2drst1.X(), pt2drst1.Y(), bid, d1urst1, d1vrst1);
265       t3dto2d(a, b, tgrst1, d1urst1, d1vrst1);
266       tg2drst1.SetCoord(a, b);
267       surf2->D1(pt2drst2.X(), pt2drst2.Y(), bid, d1urst2, d1vrst2);
268       t3dto2d(a, b, tgrst1, d1urst2, d1vrst2);
269       tg2drst2.SetCoord(a, b);
270     }
271  
272     gp_Pnt Center;
273     gp_Vec NotUsed;
274     Standard_Boolean IsCenter;
275
276     IsCenter = CenterCircleRst1Rst2(ptrst1, ptrst2, nplan, Center, NotUsed);
277
278     if (!IsCenter) return Standard_False;
279
280     gp_Vec n1(Center, ptrst1) , n2(Center, ptrst2);
281
282     n1.Normalize();
283     n2.Normalize();
284     
285     Cosa = n1.Dot(n2);
286     Sina = nplan.Dot(n1.Crossed(n2));
287
288     if (choix%2 != 0) {
289       Sina = -Sina;  //nplan is changed into -nplan
290     }
291     
292     Angle = ACos(Cosa);
293     if (Sina < 0.) {
294       Angle = 2.*PI - Angle;
295     }
296     
297     if (Angle > maxang) {maxang = Angle;}
298     if (Angle < minang) {minang = Angle;}
299     distmin = Min( distmin, ptrst1.Distance(ptrst2));
300
301     return Standard_True;
302   }
303   istangent = Standard_True;
304   return Standard_False;
305 }
306
307 //=======================================================================
308 //function : GetMinimalDistance
309 //purpose  : 
310 //=======================================================================
311
312 Standard_Real BRepBlend_RstRstConstRad::GetMinimalDistance() const
313 {
314   return distmin;
315 }
316
317 //=======================================================================
318 //function : PointOnRst1
319 //purpose  : 
320 //=======================================================================
321
322 const gp_Pnt& BRepBlend_RstRstConstRad::PointOnRst1() const
323 {
324   return ptrst1;
325 }
326
327 //=======================================================================
328 //function : PointOnRst2
329 //purpose  : 
330 //=======================================================================
331
332 const gp_Pnt& BRepBlend_RstRstConstRad::PointOnRst2() const
333 {
334   return ptrst2;
335 }
336
337 //=======================================================================
338 //function : Pnt2dOnRst1
339 //purpose  : 
340 //=======================================================================
341
342 const gp_Pnt2d& BRepBlend_RstRstConstRad::Pnt2dOnRst1() const
343 {
344   return pt2drst1;
345 }
346
347 //=======================================================================
348 //function : Pnt2dOnRst2
349 //purpose  : 
350 //=======================================================================
351
352 const gp_Pnt2d& BRepBlend_RstRstConstRad::Pnt2dOnRst2() const
353 {
354   return pt2drst2;
355 }
356
357 //=======================================================================
358 //function : ParameterOnRst1
359 //purpose  : 
360 //=======================================================================
361
362 Standard_Real BRepBlend_RstRstConstRad::ParameterOnRst1() const
363 {
364   return prmrst1;
365 }
366
367 //=======================================================================
368 //function : ParameterOnRst2
369 //purpose  : 
370 //=======================================================================
371
372 Standard_Real BRepBlend_RstRstConstRad::ParameterOnRst2() const
373 {
374   return prmrst2;
375 }
376 //=======================================================================
377 //function : IsTangencyPoint
378 //purpose  : 
379 //=======================================================================
380
381 Standard_Boolean BRepBlend_RstRstConstRad::IsTangencyPoint() const
382 {
383   return istangent;
384 }
385
386 //=======================================================================
387 //function : TangentOnRst1
388 //purpose  : 
389 //=======================================================================
390
391 const gp_Vec& BRepBlend_RstRstConstRad::TangentOnRst1() const
392 {
393   if (istangent) {Standard_DomainError::Raise();}
394   return tgrst1;
395 }
396
397 //=======================================================================
398 //function : Tangent2dOnRst1
399 //purpose  : 
400 //=======================================================================
401
402 const gp_Vec2d& BRepBlend_RstRstConstRad::Tangent2dOnRst1() const
403 {
404   if (istangent) {Standard_DomainError::Raise();}
405   return tg2drst1;
406 }
407
408 //=======================================================================
409 //function : TangentOnRst2
410 //purpose  : 
411 //=======================================================================
412
413 const gp_Vec& BRepBlend_RstRstConstRad::TangentOnRst2() const
414 {
415   if (istangent) {Standard_DomainError::Raise();}
416   return tgrst2;
417 }
418
419 //=======================================================================
420 //function : Tangent2dOnRst2
421 //purpose  : 
422 //=======================================================================
423
424 const gp_Vec2d& BRepBlend_RstRstConstRad::Tangent2dOnRst2() const
425 {
426   if (istangent) {Standard_DomainError::Raise();}
427   return tg2drst2;
428 }
429
430 //=======================================================================
431 //function : Decroch
432 //purpose  : 
433 //=======================================================================
434
435 Blend_DecrochStatus BRepBlend_RstRstConstRad::Decroch(const math_Vector& Sol,
436                                                       gp_Vec&            NRst1,
437                                                       gp_Vec&            TgRst1,
438                                                       gp_Vec&            NRst2,
439                                                       gp_Vec&            TgRst2)const
440 {
441   gp_Vec NRst1InPlane, NRst2InPlane;
442   gp_Pnt PtTmp1, PtTmp2, Center;
443   gp_Vec d1u, d1v, centptrst, NotUsed;
444   Standard_Real norm, unsurnorm;
445   Standard_Real u,v;
446
447   rstref1->Value(Sol(1)).Coord(u, v);
448   surfref1->D1(u, v,PtTmp1,d1u,d1v);
449   // Normal to the reference surface 1
450   NRst1     = d1u.Crossed(d1v);  
451   rstref2->Value(Sol(2)).Coord(u, v);
452   surfref2->D1(u, v, PtTmp2, d1u, d1v);
453   // Normal to the reference surface 2
454   NRst2     = d1u.Crossed(d1v);
455
456   Standard_Boolean IsCenter;
457
458   IsCenter = CenterCircleRst1Rst2(PtTmp1, PtTmp2, nplan, Center, NotUsed);
459
460   norm      = nplan.Crossed(NRst1).Magnitude();
461   unsurnorm = 1. / norm;
462
463   NRst1InPlane.SetLinearForm(nplan.Dot(NRst1) * unsurnorm, nplan, -unsurnorm, NRst1);
464
465   centptrst.SetXYZ(PtTmp1.XYZ() - Center.XYZ());
466
467   if (centptrst.Dot(NRst1InPlane) < 0.) NRst1InPlane.Reverse();
468
469   TgRst1    = nplan.Crossed(centptrst);
470
471   norm      = nplan.Crossed(NRst2).Magnitude();
472   unsurnorm = 1./ norm;
473   NRst2InPlane.SetLinearForm(nplan.Dot(NRst2) * unsurnorm, nplan, -unsurnorm, NRst2);
474   centptrst.SetXYZ(PtTmp2.XYZ() - Center.XYZ());
475
476
477   if (centptrst.Dot(NRst2InPlane) < 0.) NRst2InPlane.Reverse();
478
479   TgRst2 = nplan.Crossed(centptrst);
480
481   if (choix %2 != 0) {
482     TgRst1.Reverse();
483     TgRst2.Reverse();
484   }
485
486   // The vectors are returned 
487   if (NRst1InPlane.Dot(TgRst1) > -1.e-10) {
488     if (NRst2InPlane.Dot(TgRst2) < 1.e-10) {
489       return Blend_DecrochBoth;
490     }        
491     else {
492       return Blend_DecrochRst1;  
493     }
494   }
495   else {
496     if (NRst2InPlane.Dot(TgRst2) < 1.e-10) {
497       return Blend_DecrochRst2;
498     }        
499     else {
500       return Blend_NoDecroch;
501     }
502   }
503   
504 }
505
506 //=======================================================================
507 //function : Set
508 //purpose  : 
509 //=======================================================================
510
511 void BRepBlend_RstRstConstRad::Set(const Standard_Real    Radius,
512                                    const Standard_Integer Choix)
513 {
514   choix = Choix;
515   ray   = Abs(Radius);
516
517 }
518
519 //=======================================================================
520 //function : Set
521 //purpose  : 
522 //=======================================================================
523
524 void BRepBlend_RstRstConstRad::Set(const BlendFunc_SectionShape TypeSection)
525 {
526   mySShape = TypeSection;
527 }
528
529
530
531 //=======================================================================
532 //function : CenterCircleRst1Rst2
533 //purpose  : Calculate the center of the circle passing by two points of restrictions
534 //=======================================================================
535 Standard_Boolean  BRepBlend_RstRstConstRad::CenterCircleRst1Rst2(const gp_Pnt&       PtRst1,
536                                                                  const gp_Pnt&       PtRst2,
537                                                                  const gp_Vec&       np,
538                                                                  gp_Pnt&             Center,
539                                                                  gp_Vec&             VdMed) const
540 {  
541   
542   gp_Vec rst1rst2(PtRst1, PtRst2);
543   gp_Vec   vdmedNor; //,NRst1;    vdmedNor directong vector of the perpendicular bisector  
544   Standard_Real norm2;
545   Standard_Real Dist;// distance between the middle of PtRst1,PtRst2 and Center
546
547   // Calculate the center of the circle 
548   VdMed = rst1rst2.Crossed(np); 
549   norm2 = rst1rst2.SquareMagnitude();
550   Dist  = ray * ray - 0.25 * norm2;
551
552   if (choix > 2) { 
553     VdMed.Reverse();
554   }
555
556   if (Dist < - 1.E-07) return Standard_False;
557
558   if (Dist > 1.E-07) {
559     Dist = sqrt(Dist); 
560     vdmedNor = VdMed.Normalized();
561     Center.SetXYZ(0.5 * rst1rst2.XYZ() + PtRst1.XYZ() + Dist * vdmedNor.XYZ());
562   }
563   else
564   {
565     Center.SetXYZ(0.5 * rst1rst2.XYZ() + PtRst1.XYZ());    
566   }
567
568   return Standard_True;
569
570 }
571
572
573
574
575
576
577 //=======================================================================
578 //function : Section
579 //purpose  : 
580 //=======================================================================
581
582 void BRepBlend_RstRstConstRad::Section(const Standard_Real Param,
583                                        const Standard_Real U,
584                                        const Standard_Real V,
585                                        Standard_Real&      Pdeb,
586                                        Standard_Real&      Pfin,
587                                        gp_Circ&               C)
588 {
589   gp_Vec ns, np, NotUsed;
590   gp_Pnt Center;
591   
592   tguide->D1(Param, ptgui, d1gui);
593   np       = d1gui.Normalized();
594   ptrst1   = cons1.Value(U);
595   ptrst2   = cons2.Value(V);
596
597   Standard_Boolean IsCenter;
598
599   IsCenter = CenterCircleRst1Rst2(ptrst1, ptrst2, np, Center, NotUsed);
600
601   C.SetRadius(Abs(ray));
602   ns = gp_Vec(Center, ptrst1).Normalized(); 
603  
604   if (choix%2 != 0) {
605     np.Reverse();
606   }
607
608   C.SetPosition(gp_Ax2(Center, np, ns));
609   Pdeb = 0; //ElCLib::Parameter(C, pts);
610   Pfin = ElCLib::Parameter(C, ptrst2);
611
612   // Test of angles negative and almost null : Special Case
613   if (Pfin > 1.5 * PI) {
614     np.Reverse();
615     C.SetPosition(gp_Ax2(Center, np, ns));
616     Pfin = ElCLib::Parameter(C, ptrst2);
617   }
618   if (Pfin < Precision::PConfusion()) Pfin += Precision::PConfusion();
619 }
620
621 //=======================================================================
622 //function : IsRational
623 //purpose  : 
624 //=======================================================================
625
626 Standard_Boolean BRepBlend_RstRstConstRad::IsRational () const
627 {
628   return  (mySShape==BlendFunc_Rational || mySShape==BlendFunc_QuasiAngular);
629 }
630
631 //=======================================================================
632 //function : GetSectionSize
633 //purpose  :
634 //=======================================================================
635
636 Standard_Real BRepBlend_RstRstConstRad::GetSectionSize() const 
637 {
638   return maxang * Abs(ray);
639 }
640
641 //=======================================================================
642 //function : GetMinimalWeight
643 //purpose  : 
644 //=======================================================================
645
646 void BRepBlend_RstRstConstRad::GetMinimalWeight(TColStd_Array1OfReal& Weights) const 
647 {
648   BlendFunc::GetMinimalWeights(mySShape, myTConv, minang, maxang, Weights );
649   // It is supposed that it does not depend on the Radius! 
650 }
651
652 //=======================================================================
653 //function : NbIntervals
654 //purpose  : 
655 //=======================================================================
656
657 Standard_Integer BRepBlend_RstRstConstRad::NbIntervals (const GeomAbs_Shape S) const
658 {
659   return guide->NbIntervals(BlendFunc::NextShape(S));
660 }
661
662 //=======================================================================
663 //function : Intervals
664 //purpose  : 
665 //=======================================================================
666
667 void BRepBlend_RstRstConstRad::Intervals (TColStd_Array1OfReal& T,
668                                           const GeomAbs_Shape S) const
669 {
670   guide->Intervals(T, BlendFunc::NextShape(S));
671 }
672
673 //=======================================================================
674 //function : GetShape
675 //purpose  : 
676 //=======================================================================
677
678 void BRepBlend_RstRstConstRad::GetShape (Standard_Integer& NbPoles,
679                                          Standard_Integer& NbKnots,
680                                          Standard_Integer& Degree,
681                                          Standard_Integer& NbPoles2d)
682 {
683   NbPoles2d = 2;
684   BlendFunc::GetShape(mySShape, maxang, NbPoles, NbKnots, Degree, myTConv);
685 }
686
687 //=======================================================================
688 //function : GetTolerance
689 //purpose  : Determine Tolerances to be used in approximations.
690 //=======================================================================
691
692 void BRepBlend_RstRstConstRad::GetTolerance(const Standard_Real BoundTol, 
693                                             const Standard_Real SurfTol, 
694                                             const Standard_Real AngleTol, 
695                                             math_Vector& Tol3d, 
696                                             math_Vector& Tol1d) const
697 {
698   Standard_Integer low = Tol3d.Lower(), up = Tol3d.Upper();
699   Standard_Real Tol;
700   Tol= GeomFill::GetTolerance(myTConv, minang, Abs(ray), 
701                                AngleTol, SurfTol);
702   Tol1d.Init(SurfTol);
703   Tol3d.Init(SurfTol);
704   Tol3d(low+1) = Tol3d(up-1) = Min(Tol, SurfTol);
705   Tol3d(low)   = Tol3d(up)   = Min(Tol, BoundTol);
706 }
707
708 //=======================================================================
709 //function : Knots
710 //purpose  : 
711 //=======================================================================
712
713 void BRepBlend_RstRstConstRad::Knots(TColStd_Array1OfReal& TKnots)
714 {
715   GeomFill::Knots(myTConv, TKnots);
716 }
717
718 //=======================================================================
719 //function : Mults
720 //purpose  : 
721 //=======================================================================
722
723 void BRepBlend_RstRstConstRad::Mults(TColStd_Array1OfInteger& TMults)
724 {
725   GeomFill::Mults(myTConv, TMults);
726 }
727
728 //=======================================================================
729 //function : Section
730 //purpose  : 
731 //=======================================================================
732
733 void BRepBlend_RstRstConstRad::Section(const Blend_Point& P,
734                                        TColgp_Array1OfPnt& Poles,
735                                        TColgp_Array1OfPnt2d& Poles2d,
736                                        TColStd_Array1OfReal& Weights)
737 {
738   gp_Vec ns, ns2, NotUsed;
739   gp_Pnt Center;
740   Standard_Real u, v;
741   
742   Standard_Real prm    = P.Parameter();
743   Standard_Integer low = Poles.Lower();
744   Standard_Integer upp = Poles.Upper();
745   
746   tguide->D1(prm,ptgui, d1gui);
747   nplan  = d1gui.Normalized();
748   
749   u = P.ParameterOnC1(); 
750   v = P.ParameterOnC2();
751
752   gp_Pnt2d  pt2d1 = rst1->Value(u);
753   gp_Pnt2d  pt2d2 = rst2->Value(v);
754
755   ptrst1 = cons1.Value(u); 
756   ptrst2 = cons2.Value(v);
757   distmin = Min (distmin, ptrst1.Distance(ptrst2)); 
758
759   Poles2d(Poles2d.Lower()).SetCoord(pt2d1.X(),pt2d1.Y());
760   Poles2d(Poles2d.Upper()).SetCoord(pt2d2.X(),pt2d2.Y());
761   
762   // Linear case
763   if (mySShape == BlendFunc_Linear) {
764     Poles(low)   = ptrst1;
765     Poles(upp)   = ptrst2;
766     Weights(low) = 1.0;
767     Weights(upp) = 1.0;
768     return;
769   }
770
771   // Calculate the center of the circle
772   Standard_Boolean IsCenter;
773   IsCenter = CenterCircleRst1Rst2(ptrst1, ptrst2, nplan, Center, NotUsed);
774
775   // normals to the section with points 
776   ns  = gp_Vec(Center, ptrst1).Normalized();  
777   ns2 = gp_Vec(Center, ptrst2).Normalized();
778
779   if (choix%2 != 0) {
780     nplan.Reverse();
781   }
782   
783   GeomFill::GetCircle(myTConv,
784                       ns, ns2, 
785                       nplan, ptrst1, ptrst2,
786                       Abs(ray), Center, 
787                       Poles, Weights);
788 }
789
790 //=======================================================================
791 //function : Section
792 //purpose  : 
793 //=======================================================================
794
795 Standard_Boolean BRepBlend_RstRstConstRad::Section(const Blend_Point& P,
796                                                    TColgp_Array1OfPnt& Poles,
797                                                    TColgp_Array1OfVec& DPoles,
798                                                    TColgp_Array1OfPnt2d& Poles2d,
799                                                    TColgp_Array1OfVec2d& DPoles2d,
800                                                    TColStd_Array1OfReal& Weights,
801                                                    TColStd_Array1OfReal& DWeights)
802 {
803   
804   gp_Vec d11, d21;
805   gp_Vec  dnplan, d1n1, d1n2;//,np2, dnp2;
806   gp_Vec temp, tgct;
807   gp_Vec d1urst, d1vrst;
808   gp_Pnt Center, NotUsed;
809   
810   Standard_Real norm2, normmed, Dist;
811   
812   math_Vector sol(1, 2), valsol(1, 2), secmember(1, 2);
813   math_Matrix gradsol(1, 2, 1, 2);
814   
815   Standard_Real prm       = P.Parameter();
816 #ifdef DEB
817   Standard_Integer NbSpan = (Poles.Length() - 1) / 2;
818 #endif
819   Standard_Integer low    = Poles.Lower();
820   Standard_Integer upp    = Poles.Upper();
821   Standard_Boolean istgt;
822   
823   tguide->D2(prm, ptgui, d1gui, d2gui);
824   normtg = d1gui.Magnitude();
825   nplan  = d1gui.Normalized();
826   dnplan.SetLinearForm(1./normtg, d2gui,
827                        -1./normtg * (nplan.Dot(d2gui)), nplan);
828   
829   sol(1)   = prmrst1 = P.ParameterOnC1();
830   sol(2)   = prmrst2 = P.ParameterOnC2();
831   pt2drst1 = rst1->Value(prmrst1);
832   pt2drst2 = rst2->Value(prmrst2);
833   
834   Values(sol, valsol, gradsol);
835   
836   cons1.D1(sol(1), ptrst1, d11);
837   cons2.D1(sol(2), ptrst2, d21);
838   
839   temp.SetXYZ(ptrst1.XYZ() - ptgui.XYZ());
840   secmember(1) = normtg - dnplan.Dot(temp);
841   
842   temp.SetXYZ(ptrst2.XYZ() - ptgui.XYZ());
843   secmember(2) = normtg - dnplan.Dot(temp);
844   
845   math_Gauss Resol(gradsol, 1.e-9);
846   
847   if (Resol.IsDone()) {
848     istgt = Standard_False;
849     Resol.Solve(secmember);
850   }
851   else {
852     math_SVD SingRS (gradsol);
853     if (SingRS.IsDone()) {
854       math_Vector DEDT(1,2);
855       DEDT = secmember;
856       SingRS.Solve(DEDT, secmember, 1.e-6);
857       istgt = Standard_False;
858     }
859     else istgt = Standard_True;
860   }
861
862   gp_Vec med;
863   gp_Vec rst1rst2(ptrst1, ptrst2);
864   Standard_Boolean IsCenter;
865
866   IsCenter = CenterCircleRst1Rst2(ptrst1, ptrst2, nplan, Center, med);
867   if (!IsCenter) return Standard_False;
868     
869   normmed = med.Magnitude();
870   med.Normalize();
871   gp_Vec n1(Center, ptrst1), n2(Center, ptrst2);
872   n1.Normalize();
873   n2.Normalize();
874
875   if (!istgt) {
876     // secmember contains derivatives of parameters on curves
877     // compared to t  
878     tgrst1 = secmember(1) * d11;
879     tgrst2 = secmember(2) * d21;
880
881     gp_Vec d1rst1rst2;
882
883     norm2      = rst1rst2.SquareMagnitude();
884     d1rst1rst2 = tgrst2 - tgrst1;        
885     Dist       = ray * ray - 0.25 * norm2;
886
887     if (Dist >  1.E-07) { 
888       gp_Vec d1P1P2CrosNp, dmed;
889       d1P1P2CrosNp = d1rst1rst2.Crossed(nplan) + rst1rst2.Crossed(dnplan);
890       // derivative of the perpendicular bisector
891       dmed = d1P1P2CrosNp - med.Dot(d1P1P2CrosNp) * med;
892       dmed /= normmed; 
893       Dist = sqrt(Dist);
894       Standard_Real d1Dist = - (0.25 / Dist) * rst1rst2.Dot(d1rst1rst2);
895
896       if  (choix > 2) {
897         dmed.Reverse();
898       }
899
900       // the derivative of coefficient Dist is located in dmed
901       dmed.SetLinearForm(Dist, dmed, d1Dist, med);
902       d1rst1rst2 *= 0.5;   
903       // derivative of the Normal to the curve in P1    
904       d1n1 = - (dmed + d1rst1rst2) / ray;
905
906       // derivative of the Normal to the curve in P2
907       d1n2 = (d1rst1rst2 - dmed) / ray; 
908     }
909     else {
910       d1rst1rst2 *= 0.5;
911       // Normal to the curve in P1    
912       d1n1 = - d1rst1rst2 / ray;
913
914       // Normal to the curve in P2
915       d1n2 = d1rst1rst2 / ray;       
916     }   
917   }
918   
919   // Tops 2d
920   
921   Poles2d(Poles2d.Lower()).SetCoord(pt2drst1.X(), pt2drst1.Y());
922   Poles2d(Poles2d.Upper()).SetCoord(pt2drst2.X(), pt2drst2.Y());
923   if (!istgt) {
924     Standard_Real a, b;
925     surf1->D1(pt2drst1.X(), pt2drst1.Y(), NotUsed, d1urst, d1vrst);
926     t3dto2d(a,b,tgrst1, d1urst, d1vrst);
927     DPoles2d(Poles2d.Lower()).SetCoord(a, b);
928
929     surf2->D1(pt2drst2.X(), pt2drst2.Y(), NotUsed, d1urst, d1vrst);
930     t3dto2d(a, b, tgrst2, d1urst, d1vrst);
931     DPoles2d(Poles2d.Upper()).SetCoord(a, b);
932   }
933   
934   // Linear case
935   if (mySShape == BlendFunc_Linear) {
936     Poles(low)   = ptrst1;
937     Poles(upp)   = ptrst2;
938     Weights(low) = 1.0;
939     Weights(upp) = 1.0;
940     if (!istgt) {
941       DPoles(low)   = tgrst1;
942       DPoles(upp)   = tgrst2;
943       DWeights(low) = 0.0;
944       DWeights(upp) = 0.0;
945     }
946     return (!istgt);
947   }
948   
949   // Case of the circle
950   // tangent to the center of the circle
951   if (!istgt) {
952     tgct.SetLinearForm(-ray, d1n1, tgrst1);
953   }
954
955   
956   if (choix%2 != 0) {
957     nplan.Reverse();
958     dnplan.Reverse();
959   }
960
961   if (!istgt) {
962     return GeomFill::GetCircle(myTConv, 
963                               n1, n2, 
964                               d1n1, d1n2, 
965                               nplan, dnplan, 
966                               ptrst1, ptrst2, 
967                               tgrst1, tgrst2, 
968                               Abs(ray), 0, 
969                               Center, tgct, 
970                               Poles, 
971                               DPoles,
972                               Weights, 
973                               DWeights); 
974   }
975   else {
976     GeomFill::GetCircle(myTConv,
977                        n1, n2, 
978                        nplan, ptrst1, ptrst2,
979                        Abs(ray), Center, 
980                        Poles, Weights);
981     return Standard_False;
982   }
983 }
984
985 //=======================================================================
986 //function : Section
987 //purpose  : 
988 //=======================================================================
989
990 Standard_Boolean BRepBlend_RstRstConstRad::Section
991 (const Blend_Point&,
992  TColgp_Array1OfPnt&,
993  TColgp_Array1OfVec&,
994  TColgp_Array1OfVec&,
995  TColgp_Array1OfPnt2d&,
996  TColgp_Array1OfVec2d&,
997  TColgp_Array1OfVec2d&,
998  TColStd_Array1OfReal&,
999  TColStd_Array1OfReal&,
1000  TColStd_Array1OfReal&)
1001 {
1002   return Standard_False;
1003 }
1004
1005
1006 void BRepBlend_RstRstConstRad::Resolution(const Standard_Integer IC2d,
1007                                           const Standard_Real Tol,
1008                                           Standard_Real& TolU,
1009                                           Standard_Real& TolV) const
1010 {
1011   if(IC2d == 1){
1012     TolU = surf1->UResolution(Tol);
1013     TolV = surf1->VResolution(Tol);
1014   }
1015   else {
1016     TolU = surf2->UResolution(Tol);
1017     TolV = surf2->VResolution(Tol);
1018   }
1019 }
1020