0025307: Problems resulted from re-triangulation
[occt.git] / src / GeomFill / GeomFill_LocationGuide.cxx
CommitLineData
b311480e 1// Created on: 1998-07-08
2// Created by: Stephanie HUMEAU
3// Copyright (c) 1998-1999 Matra Datavision
973c2be1 4// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 5//
973c2be1 6// This file is part of Open CASCADE Technology software library.
b311480e 7//
d5f74e42 8// This library is free software; you can redistribute it and/or modify it under
9// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 10// by the Free Software Foundation, with special exception defined in the file
11// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12// distribution for complete text of the license and disclaimer of any warranty.
b311480e 13//
973c2be1 14// Alternatively, this file may be used under the terms of Open CASCADE
15// commercial license or contractual agreement.
7fd59977 16
17#include <GeomFill_LocationGuide.ixx>
18#include <gp.hxx>
19#include <gp_Pnt.hxx>
20#include <gp_Vec.hxx>
21#include <gp_Dir.hxx>
22#include <gp_Trsf.hxx>
23#include <gp_GTrsf.hxx>
24#include <gp_XYZ.hxx>
25#include <gp_Ax1.hxx>
26#include <gp_Pnt2d.hxx>
27
28#include <math_Vector.hxx>
29#include <math_Gauss.hxx>
30#include <math_FunctionSetRoot.hxx>
31#include <Precision.hxx>
32
33#include <Geom_SurfaceOfRevolution.hxx>
34#include <Geom_BSplineCurve.hxx>
35#include <Geom_Curve.hxx>
36
37#include <Adaptor3d_SurfaceOfRevolution.hxx>
38#include <Adaptor3d_HSurface.hxx>
39
40#include <IntCurveSurface_IntersectionPoint.hxx>
7fd59977 41#include <Adaptor3d_Surface.hxx>
42#include <GeomAdaptor.hxx>
43#include <GeomAdaptor_HSurface.hxx>
44#include <GeomAdaptor_HCurve.hxx>
45
46
47#include <GeomFill_FunctionGuide.ixx>
48#include <GeomFill_UniformSection.hxx>
49#include <GeomFill_SectionPlacement.hxx>
50#include <Geom_TrimmedCurve.hxx>
51#include <GeomLib.hxx>
52#include <ElCLib.hxx>
53
54#include <TColStd_HArray1OfInteger.hxx>
55#include <TColStd_HArray1OfReal.hxx>
56#include <TColgp_HArray1OfPnt.hxx>
57
f9032cf2 58#include <Extrema_ExtCS.hxx>
59#include <Extrema_POnSurf.hxx>
60
7fd59977 61#if DRAW
62static Standard_Integer Affich = 0;
63#include <Approx_Curve3d.hxx>
64#include <DrawTrSurf.hxx>
65#endif
66
67//=======================================================================
68//function : TraceRevol
69//purpose : Trace la surface de revolution (Debug)
70//=======================================================================
63c629aa 71#if GEOMFILL_DEB
7fd59977 72static void TraceRevol(const Standard_Real t,
73 const Standard_Real s,
74 const Handle(GeomFill_TrihedronWithGuide)& Law,
75 const Handle(GeomFill_SectionLaw)& Section,
76 const Handle(Adaptor3d_HCurve)& Curve,
77 const gp_Mat& Trans)
78
79{
80 gp_Vec T, N, B;
81 gp_Pnt P;
7fd59977 82 gp_Ax3 Rep(gp::Origin(), gp::DZ(), gp::DX());
83
84 Curve->D0(t, P);
96a95605 85 Law->D0(t, T, N, B);
7fd59977 86
87 gp_Mat M(N.XYZ(), B.XYZ(), T.XYZ());
88 M *= Trans;
89
90 gp_Dir D = M.Column(3);
91 gp_Ax1 Ax(P,D); // axe pour la surface de revoltuion
92
93 // calculer transfo entre triedre et Oxyz
94 gp_Dir N2 = N;
95 gp_Ax3 N3(P,D,N2);
96 gp_Trsf Transfo;
97 Transfo.SetTransformation(N3, Rep);
98
99 // transformer la section
100 Standard_Real f, l,e=1.e-7;
101 Handle (Geom_Curve) S, C;
102
103 if (Section->IsConstant(e)) {
104 C = Section->ConstantSection();
105 }
106 else {
107 Standard_Integer NbPoles, NbKnots, Deg;
108 Section->SectionShape(NbPoles, NbKnots, Deg);
109 TColStd_Array1OfInteger Mult(1,NbKnots);
110 Section->Mults( Mult);
111 TColStd_Array1OfReal Knots(1,NbKnots);
112 Section->Knots(Knots);
113 TColgp_Array1OfPnt Poles(1, NbPoles);
114 TColStd_Array1OfReal Weights(1, NbPoles);
115 Section->D0(s, Poles, Weights);
116 if (Section->IsRational())
117 C = new (Geom_BSplineCurve)
118 (Poles, Weights, Knots, Mult ,
119 Deg, Section->IsUPeriodic());
120 else
121 C = new (Geom_BSplineCurve)
122 (Poles, Knots, Mult,
123 Deg, Section->IsUPeriodic());
124
125 }
126
127 f = C->FirstParameter();
128 l = C->LastParameter();
129 S = new (Geom_TrimmedCurve) (C, f, l);
130 S->Transform(Transfo);
131
132 // Surface de revolution
133 Handle (Geom_Surface) Revol = new(Geom_SurfaceOfRevolution) (S, Ax);
134 cout << "Surf Revol at parameter t = " << t << endl;
135
136#if DRAW
137 Standard_CString aName = "TheRevol" ;
138 DrawTrSurf::Set(aName,Revol);
139#endif
140}
141#endif
142
143//==================================================================
144//Function: InGoodPeriod
145//Purpose : Recadre un paramtere
146//==================================================================
147static void InGoodPeriod(const Standard_Real Prec,
148 const Standard_Real Period,
149 Standard_Real& Current)
150{
151 Standard_Real Diff=Current-Prec;
152 Standard_Integer nb = (Standard_Integer ) IntegerPart(Diff/Period);
153 Current -= nb*Period;
154 Diff = Current-Prec;
155 if (Diff > Period/2) Current -= Period;
156 else if (Diff < -Period/2) Current += Period;
157}
158
159//==================================================================
160//Function: GeomFill_LocationGuide
161//Purpose : constructor
162//==================================================================
163 GeomFill_LocationGuide::
164 GeomFill_LocationGuide (const Handle(GeomFill_TrihedronWithGuide)& Triedre)
165 : TolRes(1,3), Inf(1,3,0.), Sup(1,3,0.),
166 X(1,3), R(1,3), myStatus(GeomFill_PipeOk)
167{
168 TolRes.Init(1.e-6);
169 myLaw = Triedre; // loi de triedre
170 mySec.Nullify(); // loi de section
171 myCurve.Nullify();
172 myFirstS = myLastS = -505e77;
173
174 myNbPts = 21; // nb points pour les calculs
175 myGuide = myLaw->Guide(); // courbe guide
176 if (!myGuide->IsPeriodic()) {
177 Standard_Real f, l, delta;
178 f = myGuide->FirstParameter();
179 l = myGuide->LastParameter();
180 delta = (l-f)/100;
181 f-=delta;
182 l+=delta;
183 myGuide = myGuide->Trim(f,l,delta*1.e-7); // courbe guide
184 }// if
185
186 myPoles2d = new (TColgp_HArray2OfPnt2d)(1, 2, 1, myNbPts);
187 rotation = Standard_False; // contact ou non
188 OrigParam1 = 0; // param pour ACR quand trajectoire
189 OrigParam2 = 1; // et guide pas meme sens de parcourt
190 Trans.SetIdentity();
191 WithTrans = Standard_False;
192
193#if DRAW
194 if (Affich) {
195 Approx_Curve3d approx(myGuide, 1.e-4,
196 GeomAbs_C1,
197 15+myGuide->NbIntervals(GeomAbs_CN),
198 14);
199 if (approx.HasResult()) {
200 Standard_CString aName = "TheGuide" ;
201 DrawTrSurf::Set(aName, approx.Curve());
202 }
203 }
204#endif
205}
206
207//==================================================================
208//Function: SetRotation
209//Purpose : init et force la Rotation
210//==================================================================
211 void GeomFill_LocationGuide::SetRotation(const Standard_Real PrecAngle,
212 Standard_Real& LastAngle)
213{
214 if (myCurve.IsNull())
215 Standard_ConstructionError::Raise(
216 "GeomFill_LocationGuide::The path is not setted !!");
217
218 //repere fixe
219 gp_Ax3 Rep(gp::Origin(), gp::DZ(), gp::DX());
220// gp_Pnt P,P1,P2;
221 gp_Pnt P;
222 gp_Vec T,N,B;
223 Standard_Integer ii, Deg;
224 Standard_Boolean isconst, israt=Standard_False;
96a95605 225 Standard_Real t, v,w, OldAngle=0, Angle, DeltaG, Diff;
f9032cf2 226 Standard_Real CurAngle = PrecAngle, a1/*, a2*/;
7fd59977 227 gp_Pnt2d p1,p2;
228 Handle(Geom_SurfaceOfRevolution) Revol; // surface de revolution
229 Handle(GeomAdaptor_HSurface) Pl; // = Revol
230 Handle(Geom_TrimmedCurve) S;
231 IntCurveSurface_IntersectionPoint PInt; // intersection guide/Revol
7fd59977 232 Handle(TColStd_HArray1OfInteger) Mult;
233 Handle(TColStd_HArray1OfReal) Knots, Weights;
234 Handle(TColgp_HArray1OfPnt) Poles;
235
236
237 Standard_Real U=0, UPeriod=0;
238 Standard_Real f = myCurve->FirstParameter();
239 Standard_Real l = myCurve->LastParameter();
240 Standard_Boolean Ok, uperiodic = mySec->IsUPeriodic();
241
242 DeltaG = (myGuide->LastParameter() - myGuide->FirstParameter())/5;
243 Handle(Geom_Curve) mySection;
244 Standard_Real Tol =1.e-9;
245
246 Standard_Integer NbPoles, NbKnots;
247 mySec->SectionShape(NbPoles, NbKnots, Deg);
248
249
250 if (mySec->IsConstant(Tol)) {
251 mySection = mySec->ConstantSection();
252 Uf = mySection->FirstParameter();
253 Ul = mySection->LastParameter();
254
255 isconst = Standard_True;
256 }
257 else {
258 isconst = Standard_False;
259 israt = mySec->IsRational();
260 Mult = new (TColStd_HArray1OfInteger) (1, NbKnots);
261 mySec->Mults( Mult->ChangeArray1());
262 Knots = new (TColStd_HArray1OfReal) (1, NbKnots);
263 mySec->Knots(Knots->ChangeArray1());
264 Poles = new (TColgp_HArray1OfPnt) (1, NbPoles);
265 Weights = new (TColStd_HArray1OfReal) (1, NbPoles);
266 Uf = Knots->Value(1);
267 Ul = Knots->Value(NbKnots);
268 }
269
270 // Bornes de calculs
271 Standard_Real Delta;
272// Standard_Integer bid1, bid2, NbK;
273 Delta = myGuide->LastParameter() - myGuide->FirstParameter();
274 Inf(1) = myGuide->FirstParameter() - Delta/10;
275 Sup(1) = myGuide->LastParameter() + Delta/10;
276
c6541a0c
D
277 Inf(2) = -M_PI;
278 Sup(2) = 3*M_PI;
7fd59977 279
280 Delta = Ul - Uf;
281 Inf(3) = Uf - Delta/10;
282 Sup(3) = Ul + Delta/10;
283
284 // JALONNEMENT
7fd59977 285 if (uperiodic) UPeriod = Ul-Uf;
286
287 for (ii=1; ii<=myNbPts; ii++) {
288 t = Standard_Real(myNbPts - ii)*f + Standard_Real(ii - 1)*l;
289 t /= (myNbPts-1);
290 myCurve->D0(t, P);
291 Ok = myLaw->D0(t, T, N, B);
292 if (!Ok) {
293 myStatus = myLaw->ErrorStatus();
294 return; //Y a rien a faire.
295 }
296 gp_Dir D = T;
297 if (WithTrans) {
298 gp_Mat M(N.XYZ(), B.XYZ(), T.XYZ());
299 M *= Trans;
300 D = M.Column(3);
301 }
302 gp_Ax1 Ax(P,D); // axe pour la surface de revoltuion
303
304 // calculer transfo entre triedre et Oxyz
305 gp_Dir N2 = N;
306 gp_Ax3 N3(P,D,N2);
307 gp_Trsf Transfo;
308 Transfo.SetTransformation(N3, Rep);
309
310 // transformer la section
311 if (! isconst) {
312 U = myFirstS + (t-myCurve->FirstParameter())*ratio;
313 mySec->D0(U, Poles->ChangeArray1(), Weights->ChangeArray1());
314 if (israt)
315 mySection = new (Geom_BSplineCurve)
316 (Poles->Array1(),
317 Weights->Array1(),
318 Knots->Array1(),
319 Mult->Array1(),
320 Deg, mySec->IsUPeriodic());
321 else
322 mySection = new (Geom_BSplineCurve)
323 (Poles->Array1(),
324 Knots->Array1(),
325 Mult->Array1(),
326 Deg, mySec->IsUPeriodic());
327 S = new (Geom_TrimmedCurve) (mySection, Uf, Ul);
328 }
329 else {
330 S = new (Geom_TrimmedCurve)
331 (Handle(Geom_Curve)::DownCast(mySection->Copy()), Uf, Ul);
332 }
333 S->Transform(Transfo);
f9032cf2 334
7fd59977 335 // Surface de revolution
336 Revol = new(Geom_SurfaceOfRevolution) (S, Ax);
337
f9032cf2 338 GeomAdaptor_Surface GArevol(Revol);
339 Extrema_ExtCS DistMini(myGuide->Curve(), GArevol,
340 Precision::Confusion(), Precision::Confusion());
341 Extrema_POnCurv Pc;
342 Extrema_POnSurf Ps;
343 Standard_Real theU = 0., theV = 0.;
344
345 if (!DistMini.IsDone() || DistMini.NbExt() == 0) {
63c629aa 346#if GEOMFILL_DEB
7fd59977 347 cout <<"LocationGuide : Pas d'intersection"<<endl;
348 TraceRevol(t, U, myLaw, mySec, myCurve, Trans);
349#endif
350 Standard_Boolean SOS=Standard_False;
351 if (ii>1) {
352 // Intersection de secour entre surf revol et guide
353 // equation
354 X(1) = myPoles2d->Value(1,ii-1).Y();
355 X(2) = myPoles2d->Value(2,ii-1).X();
356 X(3) = myPoles2d->Value(2,ii-1).Y();
357 GeomFill_FunctionGuide E (mySec, myGuide, U);
358 E.SetParam(U, P, T.XYZ(), N.XYZ());
359 // resolution => angle
360 math_FunctionSetRoot Result(E, X, TolRes,
361 Inf, Sup);
362
363 if (Result.IsDone() &&
364 (Result.FunctionSetErrors().Norm() < TolRes(1)*TolRes(1)) ) {
63c629aa 365#if GEOMFILL_DEB
7fd59977 366 cout << "Ratrappage Reussi !" << endl;
367#endif
368 SOS = Standard_True;
369 math_Vector RR(1,3);
370 Result.Root(RR);
f9032cf2 371 PInt.SetValues(P, RR(2), RR(3), RR(1), IntCurveSurface_Out);
372 theU = PInt.U();
373 theV = PInt.V();
7fd59977 374 }
375 else {
63c629aa 376#if GEOMFILL_DEB
7fd59977 377 cout << "Echec du Ratrappage !" << endl;
378#endif
379 }
380 }
381 if (!SOS) {
382 myStatus = GeomFill_ImpossibleContact;
383 return;
384 }
385 }
386 else { // on prend le point d'intersection
387 // d'angle le plus proche de P
f9032cf2 388
389 Standard_Real MinDist = RealLast();
390 Standard_Integer jref = 0;
391 for (Standard_Integer j = 1; j <= DistMini.NbExt(); j++)
392 {
393 Standard_Real aDist = DistMini.SquareDistance(j);
394 if (aDist < MinDist)
395 {
396 MinDist = aDist;
397 jref = j;
398 }
399 }
400 MinDist = Sqrt(MinDist);
401 DistMini.Points(jref, Pc, Ps);
402
403 Ps.Parameter(theU, theV);
404 a1 = theU;
405
c6541a0c 406 InGoodPeriod (CurAngle, 2*M_PI, a1);
7fd59977 407 }//else
408
409 // Controle de w
f9032cf2 410 w = Pc.Parameter();
411
7fd59977 412 if (ii>1) {
413 Diff = w - myPoles2d->Value(1, ii-1).Y();
414 if (Abs(Diff) > DeltaG) {
415 if (myGuide->IsPeriodic()) {
416 InGoodPeriod (myPoles2d->Value(1, ii-1).Y(),
417 myGuide->Period(), w);
418 Diff = w - myPoles2d->Value(1, ii-1).Y();
419 }
420 }
421
63c629aa 422#if GEOMFILL_DEB
7fd59977 423 if (Abs(Diff) > DeltaG) {
424 cout << "Location :: Diff on Guide : " <<
425 Diff << endl;
426 }
427#endif
428 }
429 //Recadrage de l'angle.
f9032cf2 430 Angle = theU;
431
7fd59977 432 if (ii > 1) {
433 Diff = Angle - OldAngle;
c6541a0c
D
434 if (Abs(Diff) > M_PI) {
435 InGoodPeriod (OldAngle, 2*M_PI, Angle);
7fd59977 436 Diff = Angle - OldAngle;
437 }
63c629aa 438#if GEOMFILL_DEB
c6541a0c 439 if (Abs(Diff) > M_PI/4) {
7fd59977 440 cout << "Diff d'angle trop grand !!" << endl;
441 }
442#endif
443 }
444
445
446 //Recadrage du V
f9032cf2 447 v = theV;
448
7fd59977 449 if (ii > 1) {
450 if (uperiodic) {
451 InGoodPeriod (myPoles2d->Value(2, ii-1).Y(), UPeriod, v);
452 }
453 Diff = v - myPoles2d->Value(2, ii-1).Y();
63c629aa 454#if GEOMFILL_DEB
96a95605 455 if (Abs(Diff) > (Ul-Uf)/(2+NbKnots)) {
7fd59977 456 cout << "Diff sur section trop grand !!" << endl;
457 }
458#endif
459 }
460
461 p1.SetCoord(t, w); // on stocke les parametres
462 p2.SetCoord(Angle , v);
463 CurAngle = Angle;
464 myPoles2d->SetValue(1, ii, p1);
465 myPoles2d->SetValue(2, ii, p2);
466 OldAngle = Angle;
467 }
468
469 LastAngle = CurAngle;
470 rotation = Standard_True; //C'est pret !
471}
472
473
474//==================================================================
475//Function: Set
476//Purpose : init loi de section et force la Rotation
477//==================================================================
478 void GeomFill_LocationGuide::Set(const Handle(GeomFill_SectionLaw)& Section,
479 const Standard_Boolean rotat,
480 const Standard_Real SFirst,
481 const Standard_Real SLast,
482 const Standard_Real PrecAngle,
483 Standard_Real& LastAngle)
484{
485 myStatus = GeomFill_PipeOk;
486 myFirstS = SFirst;
487 myLastS = SLast;
488 LastAngle = PrecAngle;
489 if (myCurve.IsNull())
490 ratio = 0.;
491 else
492 ratio = (SLast-SFirst) / (myCurve->LastParameter() -
493 myCurve->FirstParameter());
494 mySec = Section;
495
496 if (rotat) SetRotation(PrecAngle, LastAngle);
497 else rotation = Standard_False;
498}
499
500//==================================================================
501//Function: EraseRotation
502//Purpose : Supprime la Rotation
503//==================================================================
504 void GeomFill_LocationGuide:: EraseRotation()
505{
506 rotation = Standard_False;
507 if (myStatus == GeomFill_ImpossibleContact) myStatus = GeomFill_PipeOk;
508}
509
510//==================================================================
511//Function: Copy
512//Purpose :
513//==================================================================
514 Handle(GeomFill_LocationLaw) GeomFill_LocationGuide::Copy() const
515{
516 Standard_Real la;
517 Handle(GeomFill_TrihedronWithGuide) L;
518 L = Handle(GeomFill_TrihedronWithGuide)::DownCast(myLaw->Copy());
519 Handle(GeomFill_LocationGuide) copy = new
520 (GeomFill_LocationGuide) (L);
521 copy->SetOrigine(OrigParam1, OrigParam2);
522 copy->Set(mySec, rotation, myFirstS, myLastS,
523 myPoles2d->Value(1,1).X(), la);
524 copy->SetTrsf(Trans);
525
526 return copy;
527}
528
529
530//==================================================================
531//Function: SetCurve
532//Purpose : Calcul des poles sur la surface d'arret (intersection
533// courbe guide / surface de revolution en myNbPts points)
534//==================================================================
535 void GeomFill_LocationGuide::SetCurve(const Handle(Adaptor3d_HCurve)& C)
536{
537 Standard_Real LastAngle;
538 myCurve = C;
539 myTrimmed = C;
540
541 if (!myCurve.IsNull()){
542 myLaw->SetCurve(C);
543 myLaw->Origine(OrigParam1, OrigParam2);
544 myStatus = myLaw->ErrorStatus();
545
546 if (rotation) SetRotation(myPoles2d->Value(1,1).X(), LastAngle);
547 }
548}
549
550//==================================================================
551//Function: GetCurve
552//Purpose : return the trajectoire
553//==================================================================
554 const Handle(Adaptor3d_HCurve)& GeomFill_LocationGuide::GetCurve() const
555{
556 return myCurve;
557}
558
559//==================================================================
560//Function: SetTrsf
561//Purpose :
562//==================================================================
563 void GeomFill_LocationGuide::SetTrsf(const gp_Mat& Transfo)
564{
565 Trans = Transfo;
566 gp_Mat Aux;
567 Aux.SetIdentity();
568 Aux -= Trans;
569 WithTrans = Standard_False; // Au cas ou Trans = I
570 for (Standard_Integer ii=1; ii<=3 && !WithTrans ; ii++)
571 for (Standard_Integer jj=1; jj<=3 && !WithTrans; jj++)
572 if (Abs(Aux.Value(ii, jj)) > 1.e-14) WithTrans = Standard_True;
573}
574
575//==================================================================
576//Function: D0
577//Purpose :
578//==================================================================
579 Standard_Boolean GeomFill_LocationGuide::D0(const Standard_Real Param,
580 gp_Mat& M,
581 gp_Vec& V)
582{
583 Standard_Boolean Ok;
584 gp_Vec T,N,B;
585 gp_Pnt P;
586
587 myCurve->D0(Param, P);
588 V.SetXYZ(P.XYZ());
589 Ok = myLaw->D0(Param, T, N, B);
590 if (!Ok) {
591 myStatus = myLaw->ErrorStatus();
592 return Ok;
593 }
594 M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
595
596 if (WithTrans) {
597 M *= Trans;
598 }
599
600 if(rotation) {
601 Standard_Real U = myFirstS +
602 (Param-myCurve->FirstParameter())*ratio;
603 // initialisations germe
604 InitX(Param);
605
606 Standard_Integer Iter = 100;
607 gp_XYZ t,b,n;
608 t = M.Column(3);
609 b = M.Column(2);
610 n = M.Column(1);
611
612 // Intersection entre surf revol et guide
613 // equation
614 GeomFill_FunctionGuide E (mySec, myGuide, U);
615 E.SetParam(Param, P, t, n);
616 // resolution => angle
617 math_FunctionSetRoot Result(E, X, TolRes,
618 Inf, Sup, Iter);
619
620 if (Result.IsDone()) {
621 // solution
622 Result.Root(R);
623
624 // rotation
625 gp_Mat Rot;
626 Rot.SetRotation(t, R(2));
627 b *= Rot;
628 n *= Rot;
629
630 M.SetCols(n, b, t);
631 }
632 else {
63c629aa 633#if GEOMFILL_DEB
7fd59977 634 cout << "LocationGuide::D0 : No Result !"<<endl;
635 TraceRevol(Param, U, myLaw, mySec, myCurve, Trans);
636#endif
637 myStatus = GeomFill_ImpossibleContact;
638 return Standard_False;
639 }
640 }
641
642 return Standard_True;
643}
644
645//==================================================================
646//Function: D0
647//Purpose : calcul de l'intersection (C0) surface revol / guide
648//==================================================================
649 Standard_Boolean GeomFill_LocationGuide::D0(const Standard_Real Param,
650 gp_Mat& M,
651 gp_Vec& V,
652// TColgp_Array1OfPnt2d& Poles2d)
653 TColgp_Array1OfPnt2d& )
654{
655 gp_Vec T, N, B;
656 gp_Pnt P;
657 Standard_Boolean Ok;
7fd59977 658
659 myCurve->D0(Param, P);
660 V.SetXYZ(P.XYZ());
661 Ok = myLaw->D0(Param, T, N, B);
662 if (!Ok) {
663 myStatus = myLaw->ErrorStatus();
664 return Ok;
665 }
666 M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
667
668 if (WithTrans) {
669 M *= Trans;
670 }
671
672 if (rotation) {
7fd59977 673 //initialisation du germe
674 InitX(Param);
675 Standard_Integer Iter = 100;
676 gp_XYZ b, n, t;
677 t = M.Column(3);
678 b = M.Column(2);
679 n = M.Column(1);
680
681 // equation d'intersection entre surf revol et guide => angle
682 GeomFill_FunctionGuide E (mySec, myGuide, myFirstS +
683 (Param-myCurve->FirstParameter())*ratio);
684 E.SetParam(Param, P, t, n);
685
686 // resolution
687 math_FunctionSetRoot Result(E, X, TolRes,
688 Inf, Sup, Iter);
689
690 if (Result.IsDone()) {
691 // solution
692 Result.Root(R);
693
694 // rotation
695 gp_Mat Rot;
696 Rot.SetRotation(t, R(2));
697
698
699 b *= Rot;
700 n *= Rot;
701
702 M.SetCols(n, b, t);
703 }
704 else {
63c629aa 705#if GEOMFILL_DEB
6e6cd5d9 706 Standard_Real U = myFirstS + ratio*(Param-myCurve->FirstParameter());
7fd59977 707 cout << "LocationGuide::D0 : No Result !"<<endl;
708 TraceRevol(Param, U, myLaw, mySec, myCurve, Trans);
709#endif
710 myStatus = GeomFill_ImpossibleContact;
711 return Standard_False;
712 }
713 }
714
715 return Standard_True;
716}
717
718
719//==================================================================
720//Function: D1
721//Purpose : calcul de l'intersection (C1) surface revol / guide
722//==================================================================
723 Standard_Boolean GeomFill_LocationGuide::D1(const Standard_Real Param,
724 gp_Mat& M,
725 gp_Vec& V,
726 gp_Mat& DM,
727 gp_Vec& DV,
728// TColgp_Array1OfPnt2d& Poles2d,
729 TColgp_Array1OfPnt2d& ,
730// TColgp_Array1OfVec2d& DPoles2d)
731 TColgp_Array1OfVec2d& )
732{
733// gp_Vec T, N, B, DT, DN, DB, T0, N0, B0;
734 gp_Vec T, N, B, DT, DN, DB;
735// gp_Pnt P, P0;
736 gp_Pnt P;
737 Standard_Boolean Ok;
738
739 myCurve->D1(Param, P, DV);
740 V.SetXYZ(P.XYZ());
741 Ok = myLaw->D1(Param, T, DT, N, DN, B, DB);
742 if (!Ok) {
743 myStatus = myLaw->ErrorStatus();
744 return Ok;
745 }
746 M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
747 DM.SetCols(DN.XYZ() , DB.XYZ(), DT.XYZ());
748
749 if (WithTrans) {
750 M *= Trans;
751 DM *= Trans;
752 }
753
754 if (rotation) {
755 return Standard_False;
4e76d93b 756 /*
7fd59977 757#ifdef DEB
758 Standard_Real U = myFirstS + ratio*(Param-myCurve->FirstParameter());
759#else
760 myCurve->FirstParameter() ;
761#endif
762
763 // initialisation du germe
764 InitX(Param);
765
766 Standard_Integer Iter = 100;
767 gp_XYZ t,b,n, dt, db, dn;
768 t = M.Column(3);
769 b = M.Column(2);
770 n = M.Column(1);
771 dt = M.Column(3);
772 db = M.Column(2);
773 dn = M.Column(1);
774
775 // equation d'intersection surf revol / guide => angle
776 GeomFill_FunctionGuide E (mySec, myGuide, myFirstS +
777 (Param-myCurve->FirstParameter())*ratio);
778 E.SetParam(Param, P, t, n);
779
780 // resolution
781 math_FunctionSetRoot Result(E, X, TolRes,
782 Inf, Sup, Iter);
783
784 if (Result.IsDone())
785 {
786 // solution de la fonction
787 Result.Root(R);
788
789 // derivee de la fonction
790 math_Vector DEDT(1,3);
791 E.DerivT(R, DV.XYZ(), dt, DEDT); // dE/dt => DEDT
792
793 math_Vector DSDT (1,3,0);
794 math_Matrix DEDX (1,3,1,3,0);
795 E.Derivatives(R, DEDX); // dE/dx au point R => DEDX
796
797 // resolution du syst. : DEDX*DSDT = -DEDT
798 math_Gauss Ga(DEDX);
799 if (Ga.IsDone())
800 {
801 Ga.Solve (DEDT.Opposite(), DSDT);// resolution du syst.
802 }//if
803 else {
63c629aa 804#if GEOMFILL_DEB
7fd59977 805 cout << "DEDX = " << DEDX << endl;
806 cout << "DEDT = " << DEDT << endl;
807#endif
808 Standard_ConstructionError::Raise(
809 "LocationGuide::D1 : No Result dans la derivee");
810 }
811
812 // transformation = rotation
813 gp_Mat Rot, DRot;
814 Rot.SetRotation(t, R(2));
815
816
817
818 M.SetCols(n*Rot, b*Rot, t);
819
820 // transfo entre triedre (en Q) et Oxyz
821 gp_Ax3 Rep(gp::Origin(),gp::DZ(), gp::DX());
822 gp_Ax3 RepTriedre(gp::Origin(),t,n);
823 gp_Trsf Transfo3;
824 Transfo3.SetTransformation(Rep,RepTriedre);
825 // on se place dans Oxyz
826 Transfo3.Transforms(n);
827 Transfo3.Transforms(b);
828 Transfo3.Transforms(dn);
829 Transfo3.Transforms(db);
830
831 // matrices de rotation et derivees
832 Standard_Real A = R(2);
833 Standard_Real Aprim = DSDT(2);
834
835#ifdef DEB
836 gp_Mat M2 (Cos(A), -Sin(A),0, // rotation autour de T
837 Sin(A), Cos(A),0,
838 0,0,1);
839#endif
840
841 gp_Mat M2prim (-Sin(A), -Cos(A), 0, // derivee rotation autour de T
842 Cos(A), -Sin(A), 0,
843 0, 0, 0);
844 M2prim.Multiply(Aprim);
845
846 // transformations
847
848
849 dn *= Rot;
850 db *= Rot;
851
852 n *= DRot;
853 b *= DRot;
854
855 dn += n;
856 db += b;
857
858 // on repasse dans repere triedre
859 gp_Trsf InvTrsf;
860 InvTrsf = Transfo3.Inverted();
861 InvTrsf.Transforms(dn);
862 InvTrsf.Transforms(db);
863
864 DM.SetCols(dn , db , dt);
865 }//if_Result
866
867 else {
63c629aa 868#if GEOMFILL_DEB
7fd59977 869 cout << "LocationGuide::D1 : No Result !!"<<endl;
870 TraceRevol(Param, U, myLaw, mySec, myCurve, Trans);
871#endif
872 myStatus = GeomFill_ImpossibleContact;
873 return Standard_False;
874 }
4e76d93b 875*/
7fd59977 876 }//if_rotation
877
878
879 return Standard_True;
880
881}
882
883//==================================================================
884//Function: D2
885//Purpose : calcul de l'intersection (C2) surface revol / guide
886//==================================================================
887Standard_Boolean GeomFill_LocationGuide::D2(const Standard_Real Param,
888 gp_Mat& M,
889 gp_Vec& V,
890 gp_Mat& DM,
891 gp_Vec& DV,
892 gp_Mat& D2M,
893 gp_Vec& D2V,
894// TColgp_Array1OfPnt2d& Poles2d,
895 TColgp_Array1OfPnt2d& ,
896// TColgp_Array1OfVec2d& DPoles2d,
897 TColgp_Array1OfVec2d& ,
898// TColgp_Array1OfVec2d& D2Poles2d)
899 TColgp_Array1OfVec2d& )
900{
901 gp_Vec T, N, B, DT, DN, DB, D2T, D2N, D2B;
902// gp_Vec T0, N0, B0, T1, N1, B1;
903// gp_Pnt P, P0, P1;
904 gp_Pnt P;
905 Standard_Boolean Ok;
906
907 myCurve->D2(Param, P, DV, D2V);
908 V.SetXYZ(P.XYZ());
909 Ok = myLaw->D2(Param, T, DT, D2T, N, DN, D2N, B, DB, D2B);
910 if (!Ok) {
911 myStatus = myLaw->ErrorStatus();
912 return Ok;
913 }
914
915 if (WithTrans) {
916 M *= Trans;
917 DM *= Trans;
918 D2M *= Trans;
919 }
920
921 if (rotation)
922 {
923 return Standard_False;
924/*
925 Standard_Real U = myFirstS +
926 (Param-myCurve->FirstParameter())*ratio;
927 // rotation
928 math_Vector X(1,3,0);
929 InitX(Param,X);
930 // tolerance sur X
931
932 TolRes.Init(1.e-6);
933 // tolerance sur E
934// Standard_Real ETol = 1.e-6;
935 Standard_Integer Iter = 100;
936
937
938 // resoudre equation d'intersection entre surf revol et guide => angle
939 GeomFill_FunctionGuide E (mySec, myGuide, myFirstS +
940 (Param-myCurve->FirstParameter())*ratio);
941 E.SetParam(Param, P, T, N);
942
943 // resolution
944 math_FunctionSetRoot Result(E, X, TolRes,
945 Inf, Sup, Iter);
946
947 if (Result.IsDone())
948 {
949 Result.Root(R); // solution
950
951 //gp_Pnt2d p (R(2), R(3)); // point sur la surface (angle, v)
952 //Poles2d.SetValue(1,p);
953
954 // derivee de la fonction
955 math_Vector DEDT(1,3,0);
956 E.DerivT(Param, Param0, R, R0, DEDT); // dE/dt => DEDT
957 math_Vector DSDT (1,3,0);
958 math_Matrix DEDX (1,3,1,3,0);
959 E.Derivatives(R, DEDX); // dE/dx au point R => DEDX
960
961 // resolution du syst. lin. : DEDX*DSDT = -DEDT
962 math_Gauss Ga(DEDX);
963 if (Ga.IsDone())
964 {
965 Ga.Solve (DEDT.Opposite(), DSDT); // resolution du syst. lin.
966 //gp_Vec2d dp (DSDT(2), DSDT(3)); // surface
967 //DPoles2d.SetValue(1, dp);
968 }//if
969 else cout <<"LocationGuide::D2 : No Result dans la derivee premiere"<<endl;
970
971 // deuxieme derivee
972 GeomFill_Tensor D2EDX2(3,3,3);
973 E.Deriv2X(R, D2EDX2); // d2E/dx2
974
975 math_Vector D2EDT2(1,3,0);
976
977 // if(Param1 < Param && Param < Param0)
978 E.Deriv2T(Param1, Param, Param0, R1, R, R0, D2EDT2); // d2E/dt2
979 // else if (Param < Param0 && Param0 < Param1)
980 // E.Deriv2T(Param, Param0, Param1, R, R0, R1, D2EDT2); // d2E/dt2
981 // else
982 // E.Deriv2T(Param0, Param1, Param, R0, R1, R, D2EDT2); // d2E/dt2
983
984 math_Matrix D2EDTDX(1,3,1,3,0);
985 E.DerivTX(Param, Param0, R, R0, D2EDTDX); // d2E/dtdx
986
987 math_Vector D2SDT2(1,3,0); // d2s/dt2
988 math_Matrix M1(1,3,1,3,0);
989 D2EDX2.Multiply(DSDT,M1);
990
991 // resolution du syst. lin.
992 math_Gauss Ga1 (DEDX);
993 if (Ga1.IsDone())
994 {
995 Ga1.Solve ( - M1*DSDT - 2*D2EDTDX*DSDT - D2EDT2 , D2SDT2);
996 //gp_Vec2d d2p (D2SDT2(2), D2SDT2(3)); // surface
997 //D2Poles2d.SetValue(1, d2p);
998 }//if
999 else {
1000 cout <<"LocationGuide::D2 : No Result dans la derivee seconde"<<endl;
1001 myStatus = GeomFill_ImpossibleContact;
1002 }
1003
1004//------------------------------------------
1005// rotation
1006//------------------------------------------
1007
1008 gp_Trsf Tr;
1009 gp_Pnt Q (0, 0 ,0);
1010 gp_Ax1 Axe (Q, D);
1011 Tr.SetRotation(Axe, R(2));
1012
1013 gp_Vec b,b2;
1014 b = b2 = B;
1015 gp_Vec n,n2;
1016 n = n2 = N;
1017
1018 B.Transform(Tr);
1019 N.Transform(Tr);
1020
1021 M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
1022
1023//------------------------------------------
1024// derivees de la rotation
1025// A VERIFIER !!!!
1026//-----------------------------------------
1027 gp_Vec db,dn,db3,dn3;
1028 db = db3 = DB;
1029 dn = dn3 = DN;
1030
1031 gp_Vec db1,dn1,db2,dn2;
1032
1033//transfo entre triedre et Oxyz
1034 gp_Ax3 RepTriedre4(Q,D,B2);
1035 gp_Trsf Transfo3;
1036 Transfo3.SetTransformation(Rep,RepTriedre4);
1037
1038//on passe dans le repere du triedre
1039 n.Transform(Transfo3);
1040 b.Transform(Transfo3);
1041 n2.Transform(Transfo3);
1042 b2.Transform(Transfo3);
1043 dn.Transform(Transfo3);
1044 db.Transform(Transfo3);
1045 dn3.Transform(Transfo3);
1046 db3.Transform(Transfo3);
1047 D2N.Transform(Transfo3);
1048 D2B.Transform(Transfo3);
1049
1050//matrices de rotation et derivees
1051 Standard_Real A = R(2);
1052 Standard_Real Aprim = DSDT(2);
1053 Standard_Real Asec = D2SDT2(2);
1054
1055 gp_Mat M2 (Cos(A),-Sin(A),0, // rotation autour de T
1056 Sin(A), Cos(A),0,
1057 0, 0, 1);
1058
1059 gp_Mat M2prim (-Sin(A),-Cos(A),0, // derivee 1ere rotation autour de T
1060 Cos(A), -Sin(A),0,
1061 0,0,0);
1062
1063 gp_Mat M2sec (-Cos(A), Sin(A), 0, // derivee 2nde rotation autour de T
1064 -Sin(A), -Cos(A), 0,
1065 0,0,0);
1066 M2sec.Multiply(Aprim*Aprim);
1067 gp_Mat M2p = M2prim.Multiplied(Asec);
1068 M2sec.Add(M2p);
1069
1070 M2prim.Multiply(Aprim);
1071
1072// transformation
1073 gp_Trsf Rot;
1074 Rot.SetValues(M2(1,1),M2(1,2),M2(1,3),0,
1075 M2(2,1),M2(2,2),M2(2,3),0,
1076 M2(3,1),M2(3,2),M2(3,3),0,
1077 1.e-8,1.e-8);
1078 gp_Trsf DRot;
1079 DRot.SetValues(M2prim(1,1),M2prim(1,2),M2prim(1,3),0,
1080 M2prim(2,1),M2prim(2,2),M2prim(2,3),0,
1081 M2prim(3,1),M2prim(3,2),M2prim(3,3),0,
1082 1.e-8,1.e-8);
1083
1084 gp_Trsf D2Rot;
1085 D2Rot.SetValues(M2sec(1,1),M2sec(1,2),M2sec(1,3),0,
1086 M2sec(2,1),M2sec(2,2),M2sec(2,3),0,
1087 M2sec(3,1),M2sec(3,2),M2sec(3,3),0,
1088 1.e-8,1.e-8);
1089
1090
1091//derivee premiere
1092 dn.Transform(Rot);
1093 db.Transform(Rot);
1094 n.Transform(DRot);
1095 b.Transform(DRot);
1096 dn1 = dn + n;
1097 db1 = db + b;
1098 dn1.Transform(Transfo3.Inverted());
1099 db1.Transform(Transfo3.Inverted());
1100
1101 DM.SetCols(dn1.XYZ(), db1.XYZ(), DT.XYZ());
1102
1103//derivee seconde
1104 D2N.Transform(Rot);
1105 D2B.Transform(Rot);
1106 dn3.Transform(DRot);
1107 db3.Transform(DRot);
1108 n2.Transform(D2Rot);
1109 b2.Transform(D2Rot);
1110 dn2 = n2 + 2*dn3 + D2N;
1111 db2 = b2 + 2*db3 + D2B;
1112 dn2.Transform(Transfo3.Inverted());
1113 db2.Transform(Transfo3.Inverted());
1114
1115 D2M.SetCols(dn2.XYZ(), db2.XYZ(), D2T.XYZ());
1116
1117 }//if_result
1118 else {
63c629aa 1119#if GEOMFILL_DEB
7fd59977 1120 cout << "LocationGuide::D2 : No Result !!" <<endl;
1121 TraceRevol(Param, U, myLaw, mySec, myCurve, Trans);
1122#endif
1123 return Standard_False;
1124 }*/
1125 }//if_rotation
1126
1127 else
1128 {
1129 M.SetCols(N.XYZ(), B.XYZ(), T.XYZ());
1130 DM.SetCols(DN.XYZ(), DB.XYZ(), DT.XYZ());
1131 D2M.SetCols(D2N.XYZ(), D2B.XYZ(), D2T.XYZ());
1132 }
1133
1134 return Standard_True;
1135// return Standard_False;
1136}
1137
1138//==================================================================
1139//Function : HasFirstRestriction
1140//Purpose :
1141//==================================================================
1142 Standard_Boolean GeomFill_LocationGuide::HasFirstRestriction() const
1143{
1144 return Standard_False;
1145}
1146
1147//==================================================================
1148//Function : HasLastRestriction
1149//Purpose :
1150//==================================================================
1151 Standard_Boolean GeomFill_LocationGuide::HasLastRestriction() const
1152{
1153 return Standard_False;
1154}
1155
1156//==================================================================
1157//Function : TraceNumber
1158//Purpose :
1159//==================================================================
1160 Standard_Integer GeomFill_LocationGuide::TraceNumber() const
1161{
1162 return 0;
1163}
1164
1165//==================================================================
1166//Function : ErrorStatus
1167//Purpose :
1168//==================================================================
1169 GeomFill_PipeError GeomFill_LocationGuide::ErrorStatus() const
1170{
1171 return myStatus;
1172}
1173
1174//==================================================================
1175//Function:NbIntervals
1176//Purpose :
1177//==================================================================
1178 Standard_Integer GeomFill_LocationGuide::NbIntervals
1179 (const GeomAbs_Shape S) const
1180{
1181 Standard_Integer Nb_Sec, Nb_Law;
1182 Nb_Sec = myTrimmed->NbIntervals(S);
1183 Nb_Law = myLaw->NbIntervals(S);
1184
1185 if (Nb_Sec==1) {
1186 return Nb_Law;
1187 }
1188 else if (Nb_Law==1) {
1189 return Nb_Sec;
1190 }
1191
1192 TColStd_Array1OfReal IntC(1, Nb_Sec+1);
1193 TColStd_Array1OfReal IntL(1, Nb_Law+1);
1194 TColStd_SequenceOfReal Inter;
1195 myTrimmed->Intervals(IntC, S);
1196 myLaw->Intervals(IntL, S);
1197
1198 GeomLib::FuseIntervals( IntC, IntL, Inter, Precision::PConfusion()*0.99);
1199 return Inter.Length()-1;
1200
1201}
1202
1203//==================================================================
1204//Function:Intervals
1205//Purpose :
1206//==================================================================
1207 void GeomFill_LocationGuide::Intervals(TColStd_Array1OfReal& T,
1208 const GeomAbs_Shape S) const
1209{
1210 Standard_Integer Nb_Sec, Nb_Law;
1211 Nb_Sec = myTrimmed->NbIntervals(S);
1212 Nb_Law = myLaw->NbIntervals(S);
1213
1214 if (Nb_Sec==1) {
1215 myLaw->Intervals(T, S);
1216 return;
1217 }
1218 else if (Nb_Law==1) {
1219 myTrimmed->Intervals(T, S);
1220 return;
1221 }
1222
1223 TColStd_Array1OfReal IntC(1, Nb_Sec+1);
1224 TColStd_Array1OfReal IntL(1, Nb_Law+1);
1225 TColStd_SequenceOfReal Inter;
1226 myTrimmed->Intervals(IntC, S);
1227 myLaw->Intervals(IntL, S);
1228
1229 GeomLib::FuseIntervals(IntC, IntL, Inter, Precision::PConfusion()*0.99);
1230 for (Standard_Integer ii=1; ii<=Inter.Length(); ii++)
1231 T(ii) = Inter(ii);
1232}
1233
1234//==================================================================
1235//Function:SetInterval
1236//Purpose :
1237//==================================================================
1238 void GeomFill_LocationGuide::SetInterval(const Standard_Real First,
1239 const Standard_Real Last)
1240{
1241 myLaw->SetInterval(First, Last);
1242 myTrimmed = myCurve->Trim(First, Last, 0);
1243}
1244//==================================================================
1245//Function: GetInterval
1246//Purpose :
1247//==================================================================
1248 void GeomFill_LocationGuide::GetInterval(Standard_Real& First,
1249 Standard_Real& Last) const
1250{
1251 First = myTrimmed->FirstParameter();
1252 Last = myTrimmed->LastParameter();
1253}
1254
1255//==================================================================
1256//Function: GetDomain
1257//Purpose :
1258//==================================================================
1259 void GeomFill_LocationGuide::GetDomain(Standard_Real& First,
1260 Standard_Real& Last) const
1261{
1262 First = myCurve->FirstParameter();
1263 Last = myCurve->LastParameter();
1264}
1265
1266//==================================================================
1267//function : SetTolerance
1268//purpose :
1269//==================================================================
1270void GeomFill_LocationGuide::SetTolerance(const Standard_Real Tol3d,
1271 const Standard_Real )
1272{
1273 TolRes(1) = myGuide->Resolution(Tol3d);
1274 Resolution(1, Tol3d, TolRes(2), TolRes(3));
1275
1276}
1277
1278//==================================================================
1279//function : Resolution
1280//purpose : A definir
1281//==================================================================
1282//void GeomFill_LocationGuide::Resolution (const Standard_Integer Index,
1283void GeomFill_LocationGuide::Resolution (const Standard_Integer ,
1284 const Standard_Real Tol,
1285 Standard_Real& TolU,
1286 Standard_Real& TolV) const
1287{
1288 TolU = Tol/100;
1289 TolV = Tol/100;
1290}
1291
1292//==================================================================
1293//Function:GetMaximalNorm
1294//Purpose : On suppose les triedres normes => return 1
1295//==================================================================
1296 Standard_Real GeomFill_LocationGuide::GetMaximalNorm()
1297{
1298 return 1.;
1299}
1300
1301//==================================================================
1302//Function:GetAverageLaw
1303//Purpose :
1304//==================================================================
1305 void GeomFill_LocationGuide::GetAverageLaw(gp_Mat& AM,
1306 gp_Vec& AV)
1307{
1308 Standard_Integer ii;
1309 Standard_Real U, delta;
1310 gp_Vec V, V1, V2, V3;
1311
1312 myLaw->GetAverageLaw(V1, V2, V3);
1313 AM.SetCols(V1.XYZ(), V2.XYZ(), V3.XYZ());
1314
1315 AV.SetCoord(0., 0., 0.);
1316 delta = (myTrimmed->LastParameter() - myTrimmed->FirstParameter())/10;
1317 U = myTrimmed->FirstParameter();
1318 for (ii=0; ii<=myNbPts; ii++, U+=delta) {
1319 V.SetXYZ( myTrimmed->Value(U).XYZ() );
1320 AV += V;
1321 }
1322 AV = AV/(myNbPts+1);
1323}
1324
1325
1326//==================================================================
1327//Function : Section
1328//Purpose :
1329//==================================================================
1330 Handle(Geom_Curve) GeomFill_LocationGuide::Section() const
1331{
1332 return mySec->ConstantSection();
1333}
1334
1335//==================================================================
1336//Function : Guide
1337//Purpose :
1338//==================================================================
1339 Handle(Adaptor3d_HCurve) GeomFill_LocationGuide::Guide() const
1340{
1341 return myGuide;
1342}
1343
1344//==================================================================
1345//Function : IsRotation
1346//Purpose :
1347//==================================================================
1348// Standard_Boolean GeomFill_LocationGuide::IsRotation(Standard_Real& Error) const
1349 Standard_Boolean GeomFill_LocationGuide::IsRotation(Standard_Real& ) const
1350{
1351 return Standard_False;
1352}
1353
1354//==================================================================
1355//Function : Rotation
1356//Purpose :
1357//==================================================================
1358// void GeomFill_LocationGuide::Rotation(gp_Pnt& Centre) const
1359 void GeomFill_LocationGuide::Rotation(gp_Pnt& ) const
1360{
1361 Standard_NotImplemented::Raise("GeomFill_LocationGuide::Rotation");
1362}
1363
1364//==================================================================
1365//Function : IsTranslation
1366//Purpose :
1367//==================================================================
1368// Standard_Boolean GeomFill_LocationGuide::IsTranslation(Standard_Real& Error) const
1369 Standard_Boolean GeomFill_LocationGuide::IsTranslation(Standard_Real& ) const
1370{
1371 return Standard_False;
1372}
1373
1374//==================================================================
1375//Function : InitX
1376//Purpose : recherche par interpolation d'une valeur initiale
1377//==================================================================
1378void GeomFill_LocationGuide::InitX(const Standard_Real Param) const
1379{
1380
1381 Standard_Integer Ideb = 1, Ifin = myPoles2d->RowLength(), Idemi;
1382 Standard_Real Valeur, t1, t2;
1383
1384
1385 Valeur = myPoles2d->Value(1, Ideb).X();
1386 if (Param == Valeur) {
1387 Ifin = Ideb+1;
1388 }
1389
1390 Valeur = myPoles2d->Value(1, Ifin).X();
1391 if (Param == Valeur) {
1392 Ideb = Ifin-1;
1393 }
1394
1395 while ( Ideb+1 != Ifin) {
1396 Idemi = (Ideb+Ifin)/2;
1397 Valeur = myPoles2d->Value(1, Idemi).X();
1398 if (Valeur < Param) {
1399 Ideb = Idemi;
1400 }
1401 else {
1402 if ( Valeur > Param) { Ifin = Idemi;}
1403 else {
1404 Ideb = Idemi;
1405 Ifin = Ideb+1;
1406 }
1407 }
1408 }
1409
1410 t1 = myPoles2d->Value(1,Ideb).X();
1411 t2 = myPoles2d->Value(1,Ifin).X();
1412 Standard_Real diff = t2-t1;
1413
1414 Standard_Real W1, W2;
1415 W1 = myPoles2d->Value(1,Ideb).Coord(2);
1416 W2 = myPoles2d->Value(1,Ifin).Coord(2);
1417 const gp_Pnt2d& P1 = myPoles2d->Value(2, Ideb);
1418 const gp_Pnt2d& P2 = myPoles2d->Value(2, Ifin);
1419
1420 if (diff > 1.e-7) {
1421 Standard_Real b = (Param-t1) / diff,
1422 a = (t2-Param) / diff;
1423 X(1) = a * W1 + b * W2;
1424 X(2) = a * P1.Coord(1) + b * P2.Coord(1); // angle
1425 X(3) = a * P1.Coord(2) + b * P2.Coord(2); // param isov
1426 }
1427 else {
1428 X(1) = (W1+W2) /2;
1429 X(2) = (P1.Coord(1) + P2.Coord(1)) /2;
1430 X(3) = (P1.Coord(2) + P2.Coord(2)) /2;
1431 }
1432
1433 if (myGuide->IsPeriodic()) {
1434 X(1) = ElCLib::InPeriod(X(1), myGuide->FirstParameter(),
1435 myGuide->LastParameter());
1436 }
c6541a0c 1437 X(2) = ElCLib::InPeriod(X(2), 0, 2*M_PI);
7fd59977 1438 if (mySec->IsUPeriodic()) {
1439 X(3) = ElCLib::InPeriod(X(3), Uf, Ul);
1440 }
1441}
1442
1443
1444//==================================================================
1445//Function : SetOrigine
1446//Purpose : utilise pour ACR dans le cas ou la trajectoire est multi-edges
1447//==================================================================
1448void GeomFill_LocationGuide::SetOrigine(const Standard_Real Param1,
1449 const Standard_Real Param2)
1450{
1451 OrigParam1 = Param1;
1452 OrigParam2 = Param2;
1453}
1454
f9032cf2 1455//==================================================================
1456//Function : ComputeAutomaticLaw
1457//Purpose :
1458//==================================================================
1459GeomFill_PipeError GeomFill_LocationGuide::ComputeAutomaticLaw(Handle(TColgp_HArray1OfPnt2d)& ParAndRad) const
1460{
1461 gp_Pnt P;
1462 gp_Vec T,N,B;
1463 Standard_Integer ii;
1464 Standard_Real t;
1465
1466 GeomFill_PipeError theStatus = GeomFill_PipeOk;
1467
1468 Standard_Real f = myCurve->FirstParameter();
1469 Standard_Real l = myCurve->LastParameter();
1470
1471 ParAndRad = new TColgp_HArray1OfPnt2d(1, myNbPts);
1472 for (ii = 1; ii <= myNbPts; ii++)
1473 {
1474 t = Standard_Real(myNbPts - ii)*f + Standard_Real(ii - 1)*l;
1475 t /= (myNbPts-1);
1476 myCurve->D0(t, P);
1477 Standard_Boolean Ok = myLaw->D0(t, T, N, B);
1478 if (!Ok)
1479 {
1480 theStatus = myLaw->ErrorStatus();
1481 return theStatus;
1482 }
1483 gp_Pnt PointOnGuide = myLaw->CurrentPointOnGuide();
1484 Standard_Real CurWidth = P.Distance(PointOnGuide);
1485
1486 gp_Pnt2d aParamWithRadius(t, CurWidth);
1487 ParAndRad->SetValue(ii, aParamWithRadius);
1488 }
1489
1490 return theStatus;
1491}