0031035: Coding - uninitialized class fields reported by Visual Studio Code Analysis
[occt.git] / src / BlendFunc / BlendFunc_CSConstRad.cxx
1 // Created on: 1993-12-02
2 // Created by: Jacques GOUSSARD
3 // Copyright (c) 1993-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 // Modified 10/09/1996 PMN Ajout de (Nb)Intervalles, IsRationnal 
18 //                       + Utilisation de GeomFill::GetCircle dans Section.
19
20 #include <Adaptor3d_HCurve.hxx>
21 #include <Adaptor3d_HSurface.hxx>
22 #include <Blend_Point.hxx>
23 #include <BlendFunc.hxx>
24 #include <BlendFunc_CSConstRad.hxx>
25 #include <ElCLib.hxx>
26 #include <GeomFill.hxx>
27 #include <gp.hxx>
28 #include <gp_Circ.hxx>
29 #include <gp_Pnt.hxx>
30 #include <gp_Pnt2d.hxx>
31 #include <gp_Vec.hxx>
32 #include <gp_Vec2d.hxx>
33 #include <math_Gauss.hxx>
34 #include <math_Matrix.hxx>
35 #include <Precision.hxx>
36 #include <Standard_DomainError.hxx>
37 #include <Standard_NotImplemented.hxx>
38
39 //=======================================================================
40 //function : BlendFunc_CSConstRad
41 //purpose  : 
42 //=======================================================================
43 BlendFunc_CSConstRad::BlendFunc_CSConstRad(const Handle(Adaptor3d_HSurface)& S,
44                                            const Handle(Adaptor3d_HCurve)& C,
45                                            const Handle(Adaptor3d_HCurve)& CG) :
46
47        surf(S),curv(C),guide(CG), prmc(0.0), 
48        istangent(Standard_True), ray(0.0),
49        choix(0), normtg(0.0), theD(0.0),
50        maxang(RealFirst()), minang(RealLast()),
51        mySShape(BlendFunc_Rational)
52 {
53 }
54
55
56 //=======================================================================
57 //function : NbEquations
58 //purpose  : 
59 //=======================================================================
60
61 Standard_Integer BlendFunc_CSConstRad::NbEquations () const
62 {
63   return 3;
64 }
65
66
67 //=======================================================================
68 //function : Set
69 //purpose  : 
70 //=======================================================================
71
72 void BlendFunc_CSConstRad::Set(const Standard_Real Radius, const Standard_Integer Choix)
73 {
74   choix = Choix;
75   ray   = -Abs(Radius);
76 }
77
78
79 //=======================================================================
80 //function : Set
81 //purpose  : 
82 //=======================================================================
83
84 void BlendFunc_CSConstRad::Set(const BlendFunc_SectionShape TypeSection)
85 {
86   mySShape = TypeSection;
87 }
88
89
90 //=======================================================================
91 //function : Set
92 //purpose  : 
93 //=======================================================================
94
95 void BlendFunc_CSConstRad::Set(const Standard_Real Param)
96 {
97   guide->D2(Param,ptgui,d1gui,d2gui);
98   normtg = d1gui.Magnitude();
99   nplan  = d1gui.Normalized();
100   theD = - (nplan.XYZ().Dot(ptgui.XYZ()));
101 }
102
103
104 //=======================================================================
105 //function : Set
106 //purpose  : 
107 //=======================================================================
108
109 void BlendFunc_CSConstRad::Set(const Standard_Real, const Standard_Real)
110 {
111   throw Standard_NotImplemented("BlendFunc_CSConstRad::Set");
112 }
113
114
115 //=======================================================================
116 //function : GetTolerance
117 //purpose  : 
118 //=======================================================================
119
120 void BlendFunc_CSConstRad::GetTolerance(math_Vector& Tolerance, const Standard_Real Tol) const
121 {
122   Tolerance(1) = surf->UResolution(Tol);
123   Tolerance(2) = surf->VResolution(Tol);
124   Tolerance(3) = curv->Resolution(Tol);
125 }
126
127
128 //=======================================================================
129 //function : GetBounds
130 //purpose  : 
131 //=======================================================================
132
133 void BlendFunc_CSConstRad::GetBounds(math_Vector& InfBound, math_Vector& SupBound) const
134 {
135   InfBound(1) = surf->FirstUParameter();
136   InfBound(2) = surf->FirstVParameter();
137   InfBound(3) = curv->FirstParameter();
138   SupBound(1) = surf->LastUParameter();
139   SupBound(2) = surf->LastVParameter();
140   SupBound(3) = curv->LastParameter();
141
142   if(!Precision::IsInfinite(InfBound(1)) &&
143      !Precision::IsInfinite(SupBound(1))) {
144     const Standard_Real range = (SupBound(1) - InfBound(1));
145     InfBound(1) -= range;
146     SupBound(1) += range;
147   }
148   if(!Precision::IsInfinite(InfBound(2)) &&
149      !Precision::IsInfinite(SupBound(2))) {
150     const Standard_Real range = (SupBound(2) - InfBound(2));
151     InfBound(2) -= range;
152     SupBound(2) += range;
153   }
154 }
155
156
157 //=======================================================================
158 //function : IsSolution
159 //purpose  : 
160 //=======================================================================
161
162 Standard_Boolean BlendFunc_CSConstRad::IsSolution(const math_Vector& Sol, const Standard_Real Tol)
163 {
164   math_Vector valsol(1,3),secmember(1,3);
165   math_Matrix gradsol(1,3,1,3);
166
167   gp_Vec dnplan,d1u1,d1v1,d1,temp,ns,ns2,ncrossns,resul;
168   Standard_Real norm,ndotns,grosterme;
169   Standard_Real Cosa,Sina,Angle;
170
171   Values(Sol,valsol,gradsol);
172   if (Abs(valsol(1)) <= Tol &&
173       Abs(valsol(2)) <= Tol &&
174       Abs(valsol(3)) <= Tol*Tol) {
175
176     // Calcul des tangentes
177
178     pt2d = gp_Pnt2d(Sol(1),Sol(2));
179     prmc = Sol(3);
180     surf->D1(Sol(1),Sol(2),pts,d1u1,d1v1);
181     curv->D1(Sol(3),ptc,d1);
182     dnplan.SetLinearForm(1./normtg,d2gui,
183                          -1./normtg*(nplan.Dot(d2gui)),nplan);
184
185     temp.SetXYZ(pts.XYZ() - ptgui.XYZ());
186     secmember(1) = normtg - dnplan.Dot(temp);
187
188     temp.SetXYZ(ptc.XYZ() - ptgui.XYZ());
189     secmember(2) = normtg - dnplan.Dot(temp);
190
191     ns = d1u1.Crossed(d1v1);
192     ncrossns = nplan.Crossed(ns);
193     ndotns = nplan.Dot(ns);
194     norm = ncrossns.Magnitude();
195
196     grosterme = ncrossns.Dot(dnplan.Crossed(ns))/norm/norm;
197     temp.SetLinearForm(ray/norm*(dnplan.Dot(ns)-grosterme*ndotns),nplan,
198                        ray*ndotns/norm,dnplan,
199                        ray*grosterme/norm,ns);
200
201     ns.SetLinearForm(ndotns/norm,nplan, -1./norm,ns);
202     resul.SetLinearForm(ray,ns,gp_Vec(ptc,pts));
203     secmember(3) = -2.*(temp.Dot(resul));
204
205     math_Gauss Resol(gradsol);
206     if (Resol.IsDone()) {
207
208       Resol.Solve(secmember);
209
210       tgs.SetLinearForm(secmember(1),d1u1,secmember(2),d1v1);
211       tgc = secmember(3)*d1;
212       tg2d.SetCoord(secmember(1),secmember(2));
213       istangent = Standard_False;
214     }
215     else {
216       istangent = Standard_True;
217     }
218     // mise a jour de maxang
219
220     ns2 = -resul.Normalized();
221
222     Cosa = ns.Dot(ns2);
223     Sina = nplan.Dot(ns.Crossed(ns2));
224     if (choix%2 != 0) {
225       Sina = -Sina;  //nplan est change en -nplan
226     }
227
228     Angle = ACos(Cosa);
229     if (Sina <0.) {
230       Angle = 2.*M_PI - Angle;
231     }
232
233    if (Angle>maxang) {maxang = Angle;}
234    if (Angle<minang) {minang = Angle;}
235
236     return Standard_True;
237   }
238   istangent = Standard_True;
239   return Standard_False;
240 }
241
242
243 //=======================================================================
244 //function : Value
245 //purpose  : 
246 //=======================================================================
247
248 Standard_Boolean BlendFunc_CSConstRad::Value(const math_Vector& X, math_Vector& F)
249 {
250   gp_Vec d1u1,d1v1;
251   surf->D1(X(1),X(2),pts,d1u1,d1v1);
252   ptc = curv->Value(X(3));
253
254   F(1) = nplan.XYZ().Dot(pts.XYZ()) + theD;
255   F(2) = nplan.XYZ().Dot(ptc.XYZ()) + theD;
256
257   gp_Vec vref, ns = d1u1.Crossed(d1v1);
258   const Standard_Real norm = nplan.Crossed(ns).Magnitude();
259   ns.SetLinearForm(nplan.Dot(ns)/norm,nplan, -1./norm,ns);
260   vref.SetLinearForm(ray,ns,gp_Vec(ptc,pts));
261
262   F(3) = vref.SquareMagnitude() - ray*ray;
263
264   pt2d = gp_Pnt2d(X(1),X(2));
265   return Standard_True;
266 }
267
268
269 //=======================================================================
270 //function : Derivatives
271 //purpose  : 
272 //=======================================================================
273
274 Standard_Boolean BlendFunc_CSConstRad::Derivatives(const math_Vector& X, math_Matrix& D)
275 {
276   gp_Vec d1u1,d1v1,d2u1,d2v1,d2uv1,d1;
277   gp_Vec ns,ncrossns,resul,temp, vref;
278
279   Standard_Real norm,ndotns,grosterme;
280
281   surf->D2(X(1),X(2),pts,d1u1,d1v1,d2u1,d2v1,d2uv1);
282   curv->D1(X(3),ptc,d1);
283
284   D(1,1) = nplan.Dot(d1u1);
285   D(1,2) = nplan.Dot(d1v1);
286   D(1,3) = 0.;
287
288   D(2,1) = 0.;
289   D(2,2) = 0.;
290   D(2,3) = nplan.Dot(d1);
291
292
293   ns = d1u1.Crossed(d1v1);
294   ncrossns = nplan.Crossed(ns);
295   norm = ncrossns.Magnitude();
296   ndotns = nplan.Dot(ns);
297
298   vref.SetLinearForm(ndotns,nplan,-1.,ns);
299   vref.Divide(norm);
300   vref.SetLinearForm(ray,vref,gp_Vec(ptc,pts));
301
302    // Derivee par rapport a u1
303   temp = d2u1.Crossed(d1v1).Added(d1u1.Crossed(d2uv1));
304   grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
305   resul.SetLinearForm(-ray/norm*(grosterme*ndotns-nplan.Dot(temp)),nplan,
306                       ray*grosterme/norm,ns,
307                       -ray/norm,temp,
308                       d1u1);
309
310   D(3,1) = 2.*(resul.Dot(vref));
311
312
313   // Derivee par rapport a v1
314   temp = d2uv1.Crossed(d1v1).Added(d1u1.Crossed(d2v1));
315   grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
316   resul.SetLinearForm(-ray/norm*(grosterme*ndotns-nplan.Dot(temp)),nplan,
317                       ray*grosterme/norm,ns,
318                       -ray/norm,temp,
319                       d1v1);
320
321   D(3,2) = 2.*(resul.Dot(vref));
322
323   D(3,3) = -2.*(d1.Dot(vref));
324
325   pt2d = gp_Pnt2d(X(1),X(2));
326   return Standard_True;
327 }
328
329
330 //=======================================================================
331 //function : Values
332 //purpose  : 
333 //=======================================================================
334
335 Standard_Boolean BlendFunc_CSConstRad::Values(const math_Vector& X, math_Vector& F, math_Matrix& D)
336 {
337   gp_Vec d1u1,d1v1,d1;
338   gp_Vec d2u1,d2v1,d2uv1;
339   gp_Vec ns,ncrossns,resul,temp,vref;
340
341   Standard_Real norm,ndotns,grosterme;
342
343   surf->D2(X(1),X(2),pts,d1u1,d1v1,d2u1,d2v1,d2uv1);
344   curv->D1(X(3),ptc,d1);
345
346   F(1) = nplan.XYZ().Dot(pts.XYZ()) + theD;
347   F(2) = nplan.XYZ().Dot(ptc.XYZ()) + theD;
348
349   D(1,1) = nplan.Dot(d1u1);
350   D(1,2) = nplan.Dot(d1v1);
351   D(1,3) = 0.;
352
353   D(2,1) = 0.;
354   D(2,2) = 0.;
355   D(2,3) = nplan.Dot(d1);
356
357
358   ns = d1u1.Crossed(d1v1);
359   ncrossns = nplan.Crossed(ns);
360   norm = ncrossns.Magnitude();
361   ndotns = nplan.Dot(ns);
362
363   vref.SetLinearForm(ndotns,nplan,-1.,ns);
364   vref.Divide(norm);
365   vref.SetLinearForm(ray,vref,gp_Vec(ptc,pts));
366
367   F(3) = vref.SquareMagnitude() - ray*ray;
368
369
370    // Derivee par rapport a u1
371   temp = d2u1.Crossed(d1v1).Added(d1u1.Crossed(d2uv1));
372   grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
373   resul.SetLinearForm(-ray/norm*(grosterme*ndotns-nplan.Dot(temp)),nplan,
374                       ray*grosterme/norm,ns,
375                       -ray/norm,temp,
376                       d1u1);
377
378   D(3,1) = 2.*(resul.Dot(vref));
379
380
381   // Derivee par rapport a v1
382   temp = d2uv1.Crossed(d1v1).Added(d1u1.Crossed(d2v1));
383   grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
384   resul.SetLinearForm(-ray/norm*(grosterme*ndotns-nplan.Dot(temp)),nplan,
385                       ray*grosterme/norm,ns,
386                       -ray/norm,temp,
387                       d1v1);
388
389   D(3,2) = 2.*(resul.Dot(vref));
390
391   D(3,3) = -2.*(d1.Dot(vref));
392
393   pt2d = gp_Pnt2d(X(1),X(2));
394   return Standard_True;
395 }
396
397
398
399 //=======================================================================
400 //function : PointOnS
401 //purpose  : 
402 //=======================================================================
403
404 const gp_Pnt& BlendFunc_CSConstRad::PointOnS () const
405 {
406   return pts;
407 }
408
409 //=======================================================================
410 //function : PointOnC
411 //purpose  : 
412 //=======================================================================
413
414 const gp_Pnt& BlendFunc_CSConstRad::PointOnC () const
415 {
416   return ptc;
417 }
418
419 //=======================================================================
420 //function : Pnt2d
421 //purpose  : 
422 //=======================================================================
423
424 const gp_Pnt2d& BlendFunc_CSConstRad::Pnt2d () const
425 {
426   return pt2d;
427 }
428
429 //=======================================================================
430 //function : ParameterOnC
431 //purpose  : 
432 //=======================================================================
433
434 Standard_Real BlendFunc_CSConstRad::ParameterOnC () const
435 {
436   return prmc;
437 }
438
439
440 //=======================================================================
441 //function : IsTangencyPoint
442 //purpose  : 
443 //=======================================================================
444
445 Standard_Boolean BlendFunc_CSConstRad::IsTangencyPoint () const
446 {
447   return istangent;
448 }
449
450 //=======================================================================
451 //function : TangentOnS
452 //purpose  : 
453 //=======================================================================
454
455 const gp_Vec& BlendFunc_CSConstRad::TangentOnS () const
456 {
457   if (istangent)
458     throw Standard_DomainError("BlendFunc_CSConstRad::TangentOnS");
459   return tgs;
460 }
461
462 //=======================================================================
463 //function : TangentOnC
464 //purpose  : 
465 //=======================================================================
466
467 const gp_Vec& BlendFunc_CSConstRad::TangentOnC () const
468 {
469   if (istangent)
470     throw Standard_DomainError("BlendFunc_CSConstRad::TangentOnC");
471   return tgc;
472 }
473
474 //=======================================================================
475 //function : Tangent2d
476 //purpose  : 
477 //=======================================================================
478
479 const gp_Vec2d& BlendFunc_CSConstRad::Tangent2d () const
480 {
481   if (istangent)
482     throw Standard_DomainError("BlendFunc_CSConstRad::Tangent2d");
483   return tg2d;
484 }
485
486
487 //=======================================================================
488 //function : Tangent
489 //purpose  : 
490 //=======================================================================
491
492 void BlendFunc_CSConstRad::Tangent(const Standard_Real U, const Standard_Real V,
493                                    gp_Vec& TgS, gp_Vec& NmS) const
494 {
495   gp_Pnt bid;
496   gp_Vec d1u,d1v;
497   surf->D1(U,V,bid,d1u,d1v);
498
499   gp_Vec ns;
500   NmS = ns = d1u.Crossed(d1v);
501   
502   const Standard_Real norm = nplan.Crossed(ns).Magnitude();
503   ns.SetLinearForm(nplan.Dot(ns)/norm,nplan, -1./norm,ns);
504
505   gp_Pnt Center(bid.XYZ()+ray*ns.XYZ());
506   TgS = nplan.Crossed(gp_Vec(Center,bid));
507   if (choix%2 == 1)
508     TgS.Reverse();
509 }
510
511
512 //=======================================================================
513 //function : Section
514 //purpose  : 
515 //=======================================================================
516
517 void BlendFunc_CSConstRad::Section(const Standard_Real Param,
518                                    const Standard_Real U,
519                                    const Standard_Real V,
520                                    const Standard_Real W,
521                                    Standard_Real& Pdeb,
522                                    Standard_Real& Pfin,
523                                    gp_Circ& C)
524 {
525   gp_Vec d1u1,d1v1;
526   gp_Vec ns;
527   Standard_Real norm;
528   gp_Pnt Center;
529
530   guide->D1(Param,ptgui,d1gui);
531   nplan  = d1gui.Normalized();
532
533   surf->D1(U,V,pts,d1u1,d1v1);
534   ptc = curv->Value(W);
535
536   ns = d1u1.Crossed(d1v1);
537   
538   norm = nplan.Crossed(ns).Magnitude();
539   ns.SetLinearForm(nplan.Dot(ns)/norm,nplan, -1./norm,ns);
540   Center.SetXYZ(pts.XYZ()+ray*ns.XYZ());
541   C.SetRadius(Abs(ray));
542
543   if (choix%2 == 0) {
544     C.SetPosition(gp_Ax2(Center,nplan,ns));
545   }
546   else {
547     C.SetPosition(gp_Ax2(Center,nplan.Reversed(),ns));
548   }
549   Pdeb = 0.;
550   Pfin = ElCLib::Parameter(C,ptc);
551 }
552
553 Standard_Boolean BlendFunc_CSConstRad::Section(const Blend_Point&, TColgp_Array1OfPnt&, TColgp_Array1OfVec&, TColgp_Array1OfVec&, TColgp_Array1OfPnt2d&, TColgp_Array1OfVec2d&, TColgp_Array1OfVec2d&, TColStd_Array1OfReal&, TColStd_Array1OfReal&, TColStd_Array1OfReal&)
554 {
555  throw Standard_DomainError("BlendFunc_CSConstRad::Section : Not implemented");
556 }
557
558 //=======================================================================
559 //function : GetSection
560 //purpose  : 
561 //=======================================================================
562
563 Standard_Boolean BlendFunc_CSConstRad::GetSection(const Standard_Real Param,
564                                                   const Standard_Real U,
565                                                   const Standard_Real V,
566                                                   const Standard_Real W,
567                                                   TColgp_Array1OfPnt& tabP,
568                                                   TColgp_Array1OfVec& tabV)
569 {
570   Standard_Integer NbPoint=tabP.Length();
571   if (NbPoint != tabV.Length() || NbPoint < 2) {throw Standard_RangeError();}
572
573   Standard_Integer i, lowp = tabP.Lower(), lowv = tabV.Lower();
574
575   gp_Vec d1u1,d1v1,d2u1,d2v1,d2uv1,d1; //,d1u2,d1v2;
576   gp_Vec ns,dnplan,dnw,dn2w,ncrn,dncrn,ns2;
577   gp_Vec ncrossns,resul;
578   gp_Vec resulu,resulv,temp;
579
580   Standard_Real norm,ndotns,grosterme;
581   Standard_Real lambda,Cosa,Sina;
582   Standard_Real Angle = 0.,Dangle = 0.;
583   math_Vector sol(1,3),valsol(1,3),secmember(1,3);
584   math_Matrix gradsol(1,3,1,3);
585
586   guide->D2(Param,ptgui,d1gui,d2gui);
587   normtg = d1gui.Magnitude();
588   nplan  = d1gui.Normalized();
589   dnplan.SetLinearForm(1./normtg,d2gui,
590                        -1./normtg*(nplan.Dot(d2gui)),nplan);
591
592   sol(1) = U; sol(2) = V; sol(3) = W;
593
594   Values(sol,valsol,gradsol);
595
596   surf->D2(U,V,pts,d1u1,d1v1,d2u1,d2v1,d2uv1);
597   curv->D1(W,ptc,d1);
598
599   temp.SetXYZ(pts.XYZ()- ptgui.XYZ());
600   secmember(1) = normtg - dnplan.Dot(temp);
601
602   temp.SetXYZ(ptc.XYZ()- ptgui.XYZ());
603   secmember(2) = normtg - dnplan.Dot(temp);
604
605   ns = d1u1.Crossed(d1v1);
606   ncrossns = nplan.Crossed(ns);
607   ndotns = nplan.Dot(ns);
608   norm = ncrossns.Magnitude();
609
610   // Derivee de n1 par rapport a w
611
612   grosterme = ncrossns.Dot(dnplan.Crossed(ns))/norm/norm;
613   dnw.SetLinearForm((dnplan.Dot(ns)-grosterme*ndotns)/norm,nplan,
614                      ndotns/norm,dnplan,
615                      grosterme/norm,ns);
616
617   temp.SetLinearForm(ndotns/norm,nplan, -1./norm,ns);
618   resul.SetLinearForm(ray,temp,gp_Vec(ptc,pts));
619   secmember(3) = -2.*ray*(dnw.Dot(resul)); // jag 950105 il manquait ray
620
621   math_Gauss Resol(gradsol);
622   if (Resol.IsDone()) {
623     
624     Resol.Solve(secmember);
625     
626     tgs.SetLinearForm(secmember(1),d1u1,secmember(2),d1v1);
627     tgc = secmember(3)*d1;
628
629     // Derivee de n1 par rapport a u1
630     temp = d2u1.Crossed(d1v1).Added(d1u1.Crossed(d2uv1));
631     grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
632     resulu.SetLinearForm(-(grosterme*ndotns-nplan.Dot(temp))/norm,nplan,
633                          grosterme/norm,ns,
634                          -1./norm,temp);
635
636     // Derivee de n1 par rapport a v1
637     temp = d2uv1.Crossed(d1v1).Added(d1u1.Crossed(d2v1));
638     grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
639     resulv.SetLinearForm(-(grosterme*ndotns-nplan.Dot(temp))/norm,nplan,
640                          grosterme/norm,ns,
641                          -1./norm,temp);
642
643
644     dnw.SetLinearForm(secmember(1),resulu,secmember(2),resulv,dnw);
645     ns.SetLinearForm(ndotns/norm,nplan, -1./norm,ns);
646
647     dn2w.SetLinearForm(ray, dnw, -1., tgc, tgs);
648     norm = resul.Magnitude();
649     dn2w.Divide(norm);
650     ns2 = -resul.Normalized();
651     dn2w.SetLinearForm(ns2.Dot(dn2w),ns2,-1.,dn2w);
652
653     if (choix%2 != 0) {
654       nplan.Reverse();
655       dnplan.Reverse();
656     }
657
658
659     tabP(lowp) = pts;
660     tabP(lowp+NbPoint-1)  = ptc;
661
662     tabV(lowv) = tgs;
663     tabV(lowv+NbPoint-1)  = tgc;
664
665     if (NbPoint >2) {
666
667       Cosa = ns.Dot(ns2);
668       Sina = nplan.Dot(ns.Crossed(ns2));
669       Angle = ACos(Cosa);
670       if (Sina <0.) {
671         Angle = 2.*M_PI - Angle;
672       }
673       Dangle = -(dnw.Dot(ns2) + ns.Dot(dn2w))/Sina;
674       ncrn = nplan.Crossed(ns);
675       dncrn = dnplan.Crossed(ns).Added(nplan.Crossed(dnw));
676     }
677
678     for (i=2; i <= NbPoint-1; i++) {
679       lambda = (Standard_Real)(i-1)/(Standard_Real)(NbPoint-1);
680       Cosa = Cos(lambda*Angle);
681       Sina = Sin(lambda*Angle);
682       tabP(lowp+i-1).SetXYZ(pts.XYZ() + 
683                      Abs(ray)*((Cosa-1)*ns.XYZ() + Sina*ncrn.XYZ()));
684
685       temp.SetLinearForm(-Sina,ns,Cosa,ncrn);
686       temp.Multiply(lambda*Dangle);
687       temp.Add(((Cosa-1)*dnw).Added(Sina*dncrn));
688       temp.Multiply(Abs(ray));
689       temp.Add(tgs);
690       tabV(lowv+i-1)= temp;
691     } 
692     return Standard_True;
693   }
694   return Standard_False;
695 }
696
697
698 //=======================================================================
699 //function : IsRational
700 //purpose  : 
701 //=======================================================================
702
703 Standard_Boolean BlendFunc_CSConstRad::IsRational () const
704 {
705   return (mySShape==BlendFunc_Rational || mySShape==BlendFunc_QuasiAngular);
706 }
707
708 //=======================================================================
709 //function : GetSectionSize
710 //purpose  :
711 //=======================================================================
712 Standard_Real BlendFunc_CSConstRad::GetSectionSize() const 
713 {
714   return maxang*Abs(ray);
715 }
716
717 //=======================================================================
718 //function : GetMinimalWeight
719 //purpose  : 
720 //=======================================================================
721 void BlendFunc_CSConstRad::GetMinimalWeight(TColStd_Array1OfReal& Weights) const 
722 {
723   BlendFunc::GetMinimalWeights(mySShape, myTConv, minang, maxang, Weights );
724   // On suppose que cela ne depend pas du Rayon! 
725 }
726
727 //=======================================================================
728 //function : NbIntervals
729 //purpose  : 
730 //=======================================================================
731
732 Standard_Integer BlendFunc_CSConstRad::NbIntervals (const GeomAbs_Shape S) const
733 {
734  return curv->NbIntervals(BlendFunc::NextShape(S));
735 }
736
737
738 //=======================================================================
739 //function : Intervals
740 //purpose  : 
741 //=======================================================================
742
743 void BlendFunc_CSConstRad::Intervals (TColStd_Array1OfReal& T,
744                                     const GeomAbs_Shape S) const
745 {
746   curv->Intervals(T, BlendFunc::NextShape(S));
747 }
748
749 //=======================================================================
750 //function : GetShape
751 //purpose  : 
752 //=======================================================================
753
754 void BlendFunc_CSConstRad::GetShape (Standard_Integer& NbPoles,
755                                      Standard_Integer& NbKnots,
756                                      Standard_Integer& Degree,
757                                      Standard_Integer& NbPoles2d)
758 {
759   NbPoles2d = 1;
760   BlendFunc::GetShape(mySShape,maxang,NbPoles,NbKnots,Degree,myTConv);
761 }
762
763 //=======================================================================
764 //function : GetTolerance
765 //purpose  : Determine les Tolerance a utiliser dans les approximations.
766 //=======================================================================
767 void BlendFunc_CSConstRad::GetTolerance(const Standard_Real BoundTol, 
768                                         const Standard_Real SurfTol, 
769                                         const Standard_Real AngleTol, 
770                                         math_Vector& Tol3d, 
771                                         math_Vector& Tol1d) const
772 {
773  const Standard_Integer low = Tol3d.Lower();
774  const Standard_Integer up = Tol3d.Upper();
775  const Standard_Real Tol = GeomFill::GetTolerance(myTConv, minang, Abs(ray), AngleTol, SurfTol);
776  Tol1d.Init(SurfTol);
777  Tol3d.Init(SurfTol);
778  Tol3d(low+1) = Tol3d(up-1) = Min( Tol, SurfTol);
779  Tol3d(low) = Tol3d(up) = Min( Tol, BoundTol);
780 }
781
782 //=======================================================================
783 //function : Knots
784 //purpose  : 
785 //=======================================================================
786
787 void BlendFunc_CSConstRad::Knots(TColStd_Array1OfReal& TKnots)
788 {
789   GeomFill::Knots(myTConv,TKnots);
790 }
791
792
793 //=======================================================================
794 //function : Mults
795 //purpose  : 
796 //=======================================================================
797
798 void BlendFunc_CSConstRad::Mults(TColStd_Array1OfInteger& TMults)
799 {
800   GeomFill::Mults(myTConv,TMults);
801 }
802
803
804 //=======================================================================
805 //function : Section
806 //purpose  : 
807 //=======================================================================
808
809 void BlendFunc_CSConstRad::Section(const Blend_Point& P,
810                                    TColgp_Array1OfPnt& Poles,
811                                    TColgp_Array1OfPnt2d& Poles2d,
812                                    TColStd_Array1OfReal& Weights)
813 {
814   gp_Vec d1u1,d1v1;//,d1;
815   gp_Vec ns,ns2;//,temp,np2;
816   gp_Pnt Center;
817
818   Standard_Real norm,u1,v1,w;
819
820   Standard_Real prm = P.Parameter();
821   Standard_Integer low = Poles.Lower();
822   Standard_Integer upp = Poles.Upper();
823
824   guide->D1(prm,ptgui,d1gui);
825   nplan  = d1gui.Normalized();
826
827   P.ParametersOnS(u1,v1);
828   w = P.ParameterOnC();
829
830   surf->D1(u1,v1,pts,d1u1,d1v1);
831   ptc = curv->Value(w);
832
833   Poles2d(Poles2d.Lower()).SetCoord(u1,v1);
834
835   // Cas Linear
836   if (mySShape == BlendFunc_Linear) {
837     Poles(low) = pts;
838     Poles(upp) = ptc;
839     Weights(low) = 1.0;
840     Weights(upp) = 1.0;
841     return;
842   }
843
844   ns = d1u1.Crossed(d1v1);
845   norm = nplan.Crossed(ns).Magnitude();
846
847   ns.SetLinearForm(nplan.Dot(ns)/norm,nplan, -1./norm,ns);
848
849   Center.SetXYZ(pts.XYZ()+ray*ns.XYZ());
850
851   ns2 = gp_Vec(Center,ptc).Normalized();
852
853   if (choix%2 != 0) {
854     nplan.Reverse();
855   }
856
857   GeomFill::GetCircle(myTConv,
858                       ns, ns2, 
859                       nplan, pts, ptc,
860                       Abs(ray), Center, 
861                       Poles, Weights);
862 }
863  
864
865 //=======================================================================
866 //function : Section
867 //purpose  : 
868 //=======================================================================
869
870 Standard_Boolean BlendFunc_CSConstRad::Section
871   (const Blend_Point& P,
872    TColgp_Array1OfPnt& Poles,
873    TColgp_Array1OfVec& DPoles,
874    TColgp_Array1OfPnt2d& Poles2d,
875    TColgp_Array1OfVec2d& DPoles2d,
876    TColStd_Array1OfReal& Weights,
877    TColStd_Array1OfReal& DWeights)
878 {
879
880   gp_Vec d1u1,d1v1,d2u1,d2v1,d2uv1,d1;
881   gp_Vec ns,ns2,dnplan,dnw,dn2w; //,np2,dnp2;
882   gp_Vec ncrossns;
883   gp_Vec resulu,resulv,temp,tgct,resul;
884
885   gp_Pnt Center;
886
887   Standard_Real norm,ndotns,grosterme;
888
889   math_Vector sol(1,3),valsol(1,3),secmember(1,3);
890   math_Matrix gradsol(1,3,1,3);
891
892   Standard_Real prm = P.Parameter();
893   Standard_Integer low = Poles.Lower();
894   Standard_Integer upp = Poles.Upper();
895   Standard_Boolean istgt;
896
897   guide->D2(prm,ptgui,d1gui,d2gui);
898   normtg = d1gui.Magnitude();
899   nplan  = d1gui.Normalized();
900   dnplan.SetLinearForm(1./normtg,d2gui,
901                        -1./normtg*(nplan.Dot(d2gui)),nplan);
902
903   P.ParametersOnS(sol(1),sol(2));
904   sol(3) = P.ParameterOnC();
905
906   Values(sol,valsol,gradsol);
907
908   surf->D2(sol(1),sol(2),pts,d1u1,d1v1,d2u1,d2v1,d2uv1);
909   curv->D1(sol(3),ptc,d1);
910
911   temp.SetXYZ(pts.XYZ()- ptgui.XYZ());
912   secmember(1) = normtg - dnplan.Dot(temp);
913
914   temp.SetXYZ(ptc.XYZ()- ptgui.XYZ());
915   secmember(2) = normtg - dnplan.Dot(temp);
916
917   ns = d1u1.Crossed(d1v1);
918   ncrossns = nplan.Crossed(ns);
919   ndotns = nplan.Dot(ns);
920   norm = ncrossns.Magnitude();
921
922   // Derivee de n1 par rapport a w
923
924   grosterme = ncrossns.Dot(dnplan.Crossed(ns))/norm/norm;
925   dnw.SetLinearForm((dnplan.Dot(ns)-grosterme*ndotns)/norm,nplan,
926                      ndotns/norm,dnplan,
927                      grosterme/norm,ns);
928
929   temp.SetLinearForm(ndotns/norm,nplan, -1./norm,ns);
930   resul.SetLinearForm(ray,temp,gp_Vec(ptc,pts));
931   secmember(3) = -2.*ray*(dnw.Dot(resul)); // jag 950105 il manquait ray
932
933   math_Gauss Resol(gradsol);
934
935   if (Resol.IsDone()) {
936     
937     Resol.Solve(secmember);
938     
939     tgs.SetLinearForm(secmember(1),d1u1,secmember(2),d1v1);
940     tgc = secmember(3)*d1;
941
942     // Derivee de n1 par rapport a u1
943     temp = d2u1.Crossed(d1v1).Added(d1u1.Crossed(d2uv1));
944     grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
945     resulu.SetLinearForm(-(grosterme*ndotns-nplan.Dot(temp))/norm,nplan,
946                          grosterme/norm,ns,
947                          -1./norm,temp);
948
949     // Derivee de n1 par rapport a v1
950     temp = d2uv1.Crossed(d1v1).Added(d1u1.Crossed(d2v1));
951     grosterme = ncrossns.Dot(nplan.Crossed(temp))/norm/norm;
952     resulv.SetLinearForm(-(grosterme*ndotns-nplan.Dot(temp))/norm,nplan,
953                          grosterme/norm,ns,
954                          -1./norm,temp);
955
956
957     dnw.SetLinearForm(secmember(1),resulu,secmember(2),resulv,dnw);
958     ns.SetLinearForm(ndotns/norm,nplan, -1./norm,ns);
959
960     dn2w.SetLinearForm(ray, dnw, -1., tgc, tgs);
961     norm = resul.Magnitude();
962     dn2w.Divide(norm);
963     ns2 = -resul.Normalized();
964     dn2w.SetLinearForm(ns2.Dot(dn2w),ns2,-1.,dn2w);
965
966     istgt = Standard_False;
967   }
968   else {
969     ns.SetLinearForm(ndotns/norm,nplan, -1./norm,ns);
970     ns2 = -resul.Normalized();
971     istgt = Standard_True;
972   }
973   
974   // Les poles 2d
975   
976   Poles2d(Poles2d.Lower()).SetCoord(sol(1),sol(2));
977   if (!istgt) {
978     DPoles2d(Poles2d.Lower()).SetCoord(secmember(1),secmember(2));
979   }
980
981   // Cas Linear
982   if (mySShape == BlendFunc_Linear) {
983     Poles(low) = pts;
984     Poles(upp) = ptc;
985     Weights(low) = 1.0;
986     Weights(upp) = 1.0;
987     if (!istgt) {
988       DPoles(low) = tgs;
989       DPoles(upp) = tgc;
990       DWeights(low) = 0.0;
991       DWeights(upp) = 0.0;
992     }
993     return (!istgt);
994   }
995
996   // Cas du cercle
997   Center.SetXYZ(pts.XYZ()+ray*ns.XYZ());
998   if (!istgt) {
999     tgct = tgs.Added(ray*dnw);
1000   }
1001
1002   if (choix%2 != 0) {
1003     nplan.Reverse();
1004     dnplan.Reverse();
1005   }
1006   if (!istgt) {
1007     return GeomFill::GetCircle(myTConv, 
1008                                ns, ns2, 
1009                                dnw, dn2w, 
1010                                nplan, dnplan, 
1011                                pts, ptc, 
1012                                tgs, tgc, 
1013                                Abs(ray), 0, 
1014                                Center, tgct, 
1015                                Poles, 
1016                                DPoles,
1017                                Weights, 
1018                                DWeights); 
1019   }
1020   else {
1021     GeomFill::GetCircle(myTConv,
1022                         ns, ns2, 
1023                         nplan, pts, ptc,
1024                         Abs(ray), Center, 
1025                         Poles, Weights);
1026     return Standard_False;
1027   }
1028
1029 }
1030
1031
1032 void BlendFunc_CSConstRad::Resolution(const Standard_Integer , const Standard_Real Tol,
1033                                       Standard_Real& TolU, Standard_Real& TolV) const
1034 {
1035   TolU = surf->UResolution(Tol);
1036   TolV = surf->VResolution(Tol);
1037 }