projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0022792: Globally defined symbol PI conflicts with VTK definition (Intel compiler)
[occt.git] / src / GeomFill / GeomFill_SectionPlacement.cxx
CommitLineData
7fd59977 1// File: GeomFill_SectionPlacement.cxx
2// Created: Mon Dec 15 17:09:47 1997
3// Author: Philippe MANGIN
4// <pmn@sgi29>
5
6
7#include <GeomFill_SectionPlacement.ixx>
8
9#include <GeomLib.hxx>
10#include <Geom_Plane.hxx>
11#include <GeomLProp_CLProps.hxx>
12#include <GeomAbs_CurveType.hxx>
13#include <GeomAdaptor_HSurface.hxx>
14
15#include <gp_Ax3.hxx>
16#include <gp_Ax2.hxx>
17#include <gp_Pnt.hxx>
18#include <gp_Dir.hxx>
19#include <gp_Trsf.hxx>
20#include <TColgp_HArray1OfPnt.hxx>
21#include <Bnd_Box.hxx>
22#include <BndLib_Add3dCurve.hxx>
23
24#include <Precision.hxx>
25#include <gp.hxx>
26#include <Extrema_ExtCC.hxx>
27#include <Extrema_POnCurv.hxx>
28#include <IntCurveSurface_HInter.hxx>
29#include <IntCurveSurface_IntersectionPoint.hxx>
30
31#include <Geom_BSplineCurve.hxx>
32#include <Geom_TrimmedCurve.hxx>
33#include <Geom_Line.hxx>
34#include <Geom_Conic.hxx>
35#include <Geom_BezierCurve.hxx>
36
37#include <Geom_CartesianPoint.hxx>
38
39
40//===============================================================
41// Function :
42// Purpose :
43//===============================================================
44static void Tangente(const Adaptor3d_Curve& Path,
45 const Standard_Real Param,
46 gp_Pnt& P,
47 gp_Vec& Tang)
48{
49 Path.D1(Param, P, Tang);
50 Standard_Real Norm = Tang.Magnitude();
51
52 for (Standard_Integer ii=2; (ii<12) && (Norm < Precision::Confusion());
53 ii++) {
54 Tang = Path.DN(Param, ii);
55 Norm = Tang.Magnitude();
56 }
57
58 if (Norm > 100*gp::Resolution()) Tang /= Norm;
59}
60
61
62static Standard_Real Penalite(const Standard_Real angle,
63 const Standard_Real dist)
64{
65 Standard_Real penal;
66
67 if (dist < 1)
68 penal = Sqrt(dist);
69 else if (dist <2)
70 penal = Pow(dist, 2);
71 else
72 penal = dist + 2;
73
74 if (angle > 1.e-3) {
c6541a0c 75 penal += 1./angle -2./M_PI;
7fd59977 76 }
77 else {
78 penal += 1.e3;
79 }
80
81 return penal;
82}
83
84static Standard_Real EvalAngle(const gp_Vec& V1,
85 const gp_Vec& V2)
86{
87 Standard_Real angle;
88 angle = V1.Angle(V2);
c6541a0c 89 if (angle > M_PI/2) angle = M_PI - angle;
7fd59977 90 return angle;
91}
92
93static Standard_Integer NbSamples(const Handle(Geom_Curve)& aCurve)
94{
95 Standard_Real nbs = 100.; //on default
96
97 Handle(Geom_Curve) theCurve = aCurve;
98 if (aCurve->IsInstance(STANDARD_TYPE(Geom_TrimmedCurve)))
99 theCurve = (Handle(Geom_TrimmedCurve)::DownCast(aCurve))->BasisCurve();
100
101 if (theCurve->IsInstance(STANDARD_TYPE(Geom_Line)))
102 nbs = 1;
103 else if (theCurve->IsKind(STANDARD_TYPE(Geom_Conic)))
104 nbs = 4;
105 else if (theCurve->IsInstance(STANDARD_TYPE(Geom_BezierCurve)))
106 {
107 Handle(Geom_BezierCurve) BC = Handle(Geom_BezierCurve)::DownCast(theCurve);
108 nbs = 3 + BC->NbPoles();
109 }
110 else if (theCurve->IsInstance(STANDARD_TYPE(Geom_BSplineCurve)))
111 {
112 Handle(Geom_BSplineCurve) BC = Handle(Geom_BSplineCurve)::DownCast(theCurve);
113 nbs = BC->NbKnots();
114 nbs *= BC->Degree();
115 Standard_Real ratio =
116 (aCurve->LastParameter() - aCurve->FirstParameter())/(BC->LastParameter() - BC->FirstParameter());
117 nbs *= ratio;
118 if(nbs < 4.0)
119 nbs = 4;
120 }
121
122 if (nbs > 300.)
123 nbs = 300;
124
125 return ((Standard_Integer)nbs);
126}
127
128//===============================================================
129// Function :DistMini
130// Purpose : Examine un extrema pour updater <Dist> & <Param>
131//===============================================================
132static void DistMini(const Extrema_ExtPC& Ext,
133 const Adaptor3d_Curve& C,
134 Standard_Real& Dist,
135 Standard_Real& Param)
136{
137 Standard_Real dist1, dist2;
138 Standard_Integer ii;
139 gp_Pnt P1, P2;
140 Standard_Real Dist2 = RealLast();
141
142 Ext.TrimmedSquareDistances(dist1, dist2, P1, P2);
143 if ( (dist1<Dist2) || (dist2<Dist2) ) {
144 if (dist1 < dist2) {
145 Dist2 = dist1;
146 Param = C.FirstParameter();
147 }
148 else {
149 Dist2 = dist2;
150 Param = C.LastParameter();
151 }
152 }
153
154 if (Ext.IsDone())
155 {
156 for (ii=1; ii<= Ext.NbExt(); ii++) {
157 if (Ext.SquareDistance(ii) < Dist2) {
158 Dist2 = Ext.SquareDistance(ii);
159 Param = Ext.Point(ii).Parameter();
160 }
161 }
162 }
163 Dist = sqrt (Dist2);
164}
165
166
167//===============================================================
168// Function : Constructeur
169// Purpose :
170//===============================================================
171GeomFill_SectionPlacement::
172GeomFill_SectionPlacement(const Handle(GeomFill_LocationLaw)& L,
173 const Handle(Geom_Geometry)& Section) :
174 myLaw(L), /* myAdpSection(Section), mySection(Section), */
175 Dist(RealLast()), AngleMax(0.)
176{
177
178 done = Standard_False;
179 isplan = Standard_False;
180 myIsPoint = Standard_False;
181
182 if (Section->IsInstance(STANDARD_TYPE(Geom_CartesianPoint)))
183 {
184 myIsPoint = Standard_True;
185 myPoint = Handle(Geom_CartesianPoint)::DownCast(Section)->Pnt();
186 isplan = Standard_True;
187 }
188 else
189 {
190 Handle(Geom_Curve) CurveSection = Handle(Geom_Curve)::DownCast(Section);
191 myAdpSection.Load(CurveSection);
192 mySection = CurveSection;
193 }
194
195 Standard_Integer i, j, NbPoles=0;
196 Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
197
198 // Boite d'encombrement de la section pour en deduire le gabarit
199 Bnd_Box box;
200 if (myIsPoint)
201 box.Add(myPoint);
202 else
203 BndLib_Add3dCurve::Add(myAdpSection, 1.e-4, box);
204 box.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
205
206 Standard_Real DX = aXmax - aXmin ;
207 Standard_Real DY = aYmax - aYmin ;
208 Standard_Real DZ = aZmax - aZmin ;
209 Gabarit = Sqrt( DX*DX+ DY*DY + DZ*DZ )/2. ;
210
211 Gabarit += Precision::Confusion(); // Cas des toute petite
212
213
214 // Initialisation de TheAxe pour les cas singulier
215 if (!myIsPoint)
216 {
217 gp_Pnt P;
218 gp_Vec V;
219 Tangente(myAdpSection,
220 (myAdpSection.FirstParameter()+myAdpSection.LastParameter())/2,
221 P, V);
222 TheAxe.SetLocation(P);
223 TheAxe.SetDirection(V);
224
225 // y a t'il un Plan moyen ?
226 GeomAbs_CurveType TheType = myAdpSection.GetType();
227 switch (TheType) {
228 case GeomAbs_Circle:
229 {
230 isplan = Standard_True;
231 TheAxe = myAdpSection.Circle().Axis();
232 break;
233 }
234 case GeomAbs_Ellipse:
235 {
236 isplan = Standard_True;
237 TheAxe = myAdpSection.Ellipse().Axis();
238 break;
239 }
240 case GeomAbs_Hyperbola:
241 {
242 isplan = Standard_True;
243 TheAxe = myAdpSection.Hyperbola().Axis();
244 break;
245 }
246 case GeomAbs_Parabola:
247 {
248 isplan = Standard_True;
249 TheAxe = myAdpSection.Parabola().Axis();
250 break;
251 }
252 case GeomAbs_Line:
253 {
254 NbPoles = 0; // Pas de Plan !!
255 break;
256 }
257 case GeomAbs_BezierCurve:
258 case GeomAbs_BSplineCurve:
259 {
260 NbPoles = myAdpSection.NbPoles();
261 break;
262 }
263 default:
264 NbPoles = 21;
265 }
266
267
268 if (!isplan && NbPoles>2)
269 {
270 // Calcul d'un plan moyen.
271 Handle(TColgp_HArray1OfPnt) Pnts;
272 Standard_Real first = myAdpSection.FirstParameter();
273 Standard_Real last = myAdpSection.LastParameter();
274 Standard_Real t, delta;
275 if (myAdpSection.GetType() == GeomAbs_BSplineCurve)
276 {
277 Handle(Geom_BSplineCurve) BC =
278 Handle(Geom_BSplineCurve)::DownCast(myAdpSection.Curve());
279 Standard_Integer I1, I2, I3, I4;
280 BC->LocateU( first, Precision::Confusion(), I1, I2 );
281 BC->LocateU( last, Precision::Confusion(), I3, I4 );
282 Standard_Integer NbKnots = I3 - I2 + 1;
283
284 Standard_Integer NbLocalPnts = 10;
285 Standard_Integer NbPnts = (NbKnots-1) * NbLocalPnts;
286 if (I1 != I2)
287 NbPnts += NbLocalPnts;
288 if (I3 != I4 && first < BC->Knot(I3))
289 NbPnts += NbLocalPnts;
290 if (!myAdpSection.IsClosed())
291 NbPnts++;
292 Pnts = new TColgp_HArray1OfPnt(1, NbPnts);
293 Standard_Integer nb = 1;
294 if (I1 != I2)
295 {
296 Standard_Real locallast = (BC->Knot(I2) < last)? BC->Knot(I2) : last;
297 delta = (locallast - first) / NbLocalPnts;
298 for (j = 0; j < NbLocalPnts; j++)
299 {
300 t = first + j*delta;
301 Pnts->SetValue( nb++, myAdpSection.Value(t) );
302 }
303 }
304 for (i = I2; i < I3; i++)
305 {
306 t = BC->Knot(i);
307 delta = (BC->Knot(i+1) - t) / NbLocalPnts;
308 for (j = 0; j < NbLocalPnts; j++)
309 {
310 t += delta;
311 Pnts->SetValue( nb++, myAdpSection.Value(t) );
312 }
313 }
314 if (I3 != I4 && first < BC->Knot(I3))
315 {
316 t = BC->Knot(I3);
317 delta = (last - t) / NbLocalPnts;
318 for (j = 0; j < NbLocalPnts; j++)
319 {
320 t += delta;
321 Pnts->SetValue( nb++, myAdpSection.Value(t) );
322 }
323 }
324 if (!myAdpSection.IsClosed())
325 Pnts->SetValue( nb, myAdpSection.Value(last) );
326 }
327 else // other type
328 {
329 Standard_Integer NbPnts = NbPoles-1;
330 if (!myAdpSection.IsClosed())
331 NbPnts++;
332 Pnts = new TColgp_HArray1OfPnt(1, NbPnts);
333 delta = (last - first) / (NbPoles-1);
334 for (i = 0; i < NbPoles-1; i++)
335 {
336 t = first + i*delta;
337 Pnts->SetValue( i+1, myAdpSection.Value(t) );
338 }
339 if (!myAdpSection.IsClosed())
340 Pnts->SetValue( NbPnts, myAdpSection.Value(last) );
341 }
342
343 Standard_Boolean issing;
344 gp_Ax2 axe;
345 GeomLib::AxeOfInertia(Pnts->Array1(), axe, issing, Precision::Confusion());
346 if (!issing) {
347 isplan = Standard_True;
348 TheAxe.SetLocation(axe.Location());
349 TheAxe.SetDirection(axe.Direction());
350 }
351 }
352
353
354 myExt.Initialize(myAdpSection,
355 myAdpSection.FirstParameter(),
356 myAdpSection.LastParameter(),
357 Precision::Confusion());
358 }
359}
360
361//===============================================================
362// Function :SetLocation
363// Purpose :
364//===============================================================
365void GeomFill_SectionPlacement::
366SetLocation(const Handle(GeomFill_LocationLaw)& L)
367{
368 myLaw = L;
369}
370
371//===============================================================
372// Function : Perform
373// Purpose : Le plus simple
374//===============================================================
375void GeomFill_SectionPlacement::Perform(const Standard_Real Tol)
376{
377 Handle(Adaptor3d_HCurve) Path;
378 Path = myLaw->GetCurve();
379 Perform(Path, Tol);
380}
381
382//===============================================================
383// Function :Perform
384// Purpose : Recherche automatique
385//===============================================================
386void GeomFill_SectionPlacement::Perform(const Handle(Adaptor3d_HCurve)& Path,
387 const Standard_Real Tol)
388{
389 Standard_Real IntTol = 1.e-5;
390 Standard_Real DistCenter = Precision::Infinite();
391
392 if (myIsPoint)
393 {
394 Extrema_ExtPC Projector(myPoint, Path->Curve(), Precision::Confusion());
395 DistMini( Projector, Path->Curve(), Dist, PathParam );
c6541a0c 396 AngleMax = M_PI/2;
7fd59977 397 }
398 else
399 {
400 PathParam = Path->FirstParameter();
401 SecParam = myAdpSection.FirstParameter();
402
403 Standard_Real distaux, taux, alpha;
404 gp_Pnt PonPath, PonSec, P;
405 gp_Vec VRef, dp1;
406 VRef.SetXYZ(TheAxe.Direction().XYZ());
407
408 Tangente( Path->Curve(), PathParam, PonPath, dp1);
409 PonSec = myAdpSection.Value(SecParam);
410 Dist = PonPath.Distance(PonSec);
411 if (Dist > Tol) { // On Cherche un meilleur point sur la section
412 myExt.Perform(PonPath);
413 if ( myExt.IsDone() ) {
414 DistMini(myExt, myAdpSection, Dist, SecParam);
415 PonSec = myAdpSection.Value(SecParam);
416 }
417 }
418 AngleMax = EvalAngle(VRef, dp1);
c6541a0c 419 if (isplan) AngleMax = M_PI/2 - AngleMax;
7fd59977 420
421 Standard_Boolean Trouve = Standard_False;
422 Standard_Integer ii;
423
424 if (isplan) {
425 // (1.1) Distances Point-Plan
426 Standard_Real DistPlan;
427 gp_Vec V1 (PonPath, TheAxe.Location());
428 DistPlan = Abs(V1.Dot(VRef));
429 if (DistPlan <= IntTol)
430 DistCenter = V1.Magnitude();
431
432 gp_Pnt Plast = Path->Value( Path->LastParameter() );
433 V1.SetXYZ( TheAxe.Location().XYZ() - Plast.XYZ() );
434 DistPlan = Abs(V1.Dot(VRef));
435 if (DistPlan <= IntTol)
436 {
437 Standard_Real aDist = V1.Magnitude();
438 if (aDist < DistCenter)
439 {
440 DistCenter = aDist;
441 PonPath = Plast;
442 PathParam = Path->LastParameter();
443 }
444 }
445
446 // (1.2) Intersection Plan-courbe
447 gp_Ax3 axe (TheAxe.Location(), TheAxe.Direction());
448 Handle(Geom_Plane) plan = new (Geom_Plane)(axe);
449 Handle(GeomAdaptor_HSurface) adplan =
450 new (GeomAdaptor_HSurface)(plan);
451 IntCurveSurface_HInter Intersector;
452 Intersector.Perform(Path, adplan);
453 if (Intersector.IsDone())
454 {
455 Standard_Real w;
456 gp_Pnt P;
457 Standard_Real aDist;
458 for (ii=1; ii<=Intersector.NbPoints(); ii++)
459 {
460 w = Intersector.Point(ii).W();
461 P = Path->Value( w );
462 aDist = P.Distance( TheAxe.Location() );
463 if (aDist < DistCenter)
464 {
465 DistCenter = aDist;
466 PonPath = P;
467 PathParam = w;
468 }
469 }
470 }
471 if (!Intersector.IsDone() || Intersector.NbPoints() == 0)
472 {
473 Standard_Integer NbPnts = NbSamples( mySection );
474 TColgp_Array1OfPnt Pnts( 1, NbPnts+1 );
475 Standard_Real delta = (mySection->LastParameter()-mySection->FirstParameter())/NbPnts;
476 for (ii = 0; ii <= NbPnts; ii++)
477 Pnts(ii+1) = mySection->Value( mySection->FirstParameter() + ii*delta );
478
479 gp_Pnt BaryCenter;
480 gp_Dir Xdir, Ydir;
481 Standard_Real Xgap, Ygap, Zgap;
482 GeomLib::Inertia( Pnts, BaryCenter, Xdir, Ydir, Xgap, Ygap, Zgap );
483
484 gp_Pnt Pfirst = Path->Value( Path->FirstParameter() );
485 if (Pfirst.Distance(BaryCenter) < Plast.Distance(BaryCenter))
486 PathParam = Path->FirstParameter();
487 else
488 {
489 PathParam = Path->LastParameter();
490 Tangente( Path->Curve(), PathParam, PonPath, dp1);
491 PonSec = myAdpSection.Value(SecParam);
492 Dist = PonPath.Distance(PonSec);
493 if (Dist > Tol) { // On Cherche un meilleur point sur la section
494 myExt.Perform(PonPath);
495 if ( myExt.IsDone() ) {
496 DistMini(myExt, myAdpSection, Dist, SecParam);
497 PonSec = myAdpSection.Value(SecParam);
498 }
499 }
500 AngleMax = EvalAngle(VRef, dp1);
c6541a0c 501 AngleMax = M_PI/2 - AngleMax;
7fd59977 502 }
503 }
504
505 /*
506 // (1.1) Distances Point-Plan
507 Standard_Real DistPlan;
508 gp_Vec V1 (PonPath, TheAxe.Location());
509 DistPlan = Abs(V1.Dot(VRef));
510
511 // On examine l'autre extremite
512 gp_Pnt P;
513 Tangente(Path->Curve(), Path->LastParameter(), P, dp1);
514 V1.SetXYZ(TheAxe.Location().XYZ()-P.XYZ());
515 if (Abs(V1.Dot(VRef)) <= DistPlan ) { // On prend l'autre extremite
c6541a0c 516 alpha = M_PI/2 - EvalAngle(VRef, dp1);
7fd59977 517 distaux = PonPath.Distance(PonSec);
518 if (distaux > Tol) {
519 myExt.Perform(P);
520 if ( myExt.IsDone() )
521 DistMini(myExt, myAdpSection, distaux, taux);
522 }
523 else
524 taux = SecParam;
525
526 if (Choix(distaux, alpha)) {
527 Dist = distaux;
528 AngleMax = alpha;
529 PonPath = P;
530 PathParam = Path->LastParameter();
531 }
532 }
533
534 // (1.2) Intersection Plan-courbe
535 gp_Ax3 axe (TheAxe.Location(), TheAxe.Direction());
536 Handle(Geom_Plane) plan = new (Geom_Plane)(axe);
537 Handle(GeomAdaptor_HSurface) adplan =
538 new (GeomAdaptor_HSurface)(plan);
539 IntCurveSurface_HInter Intersector;
540 Intersector.Perform(Path, adplan);
541 if (Intersector.IsDone()) {
542 Standard_Real w;
543 gp_Vec V;
544 for (ii=1; ii<=Intersector.NbPoints(); ii++) {
545 w = Intersector.Point(ii).W();
546 //(1.3) test d'angle
547 Tangente( Path->Curve(), w, P, V);
c6541a0c 548 alpha = M_PI/2 - EvalAngle(V, VRef);
7fd59977 549 //(1.4) Test de distance Point-Courbe
550 myExt.Perform(P);
551 if ( myExt.IsDone() ) {
552 DistMini(myExt, myAdpSection, distaux, taux);
553 if (Choix(distaux, alpha)) {
554 Dist = distaux;
555 SecParam = taux;
556 AngleMax = alpha;
557 PonPath = P;
558 PathParam = w;
559 PonSec = myAdpSection.Value(SecParam);
560 }
561 }
562 else {
563 distaux = P.Distance(PonSec);
564 if (Choix(distaux, alpha)) {
565 Dist = distaux;
566 AngleMax = alpha;
567 PonPath = P;
568 PathParam = w;
569 }
570 }
571 }
572 }
573 */
574#if DEB
575 if (Intersector.NbPoints() == 0) {
576 Intersector.Dump();
577 }
578#endif
579 }
580
581 // Cas General
582 if (! isplan) {
583 // (2.1) Distance avec les extremites ...
584 myExt.Perform(PonPath);
585 if ( myExt.IsDone() ) {
586 DistMini(myExt, myAdpSection, distaux, taux);
587 if (distaux < Dist) {
588 Dist = distaux;
589 SecParam = taux;
590 }
591 }
592 Trouve = (Dist <= Tol);
593 if (!Trouve) {
594 Tangente( Path->Curve(), Path->LastParameter(), P, dp1);
595 alpha = EvalAngle(VRef, dp1);
596 myExt.Perform(P);
597 if ( myExt.IsDone() ) {
598 if ( myExt.IsDone() ) {
599 DistMini(myExt, myAdpSection, distaux, taux);
600 if (Choix(distaux, alpha)) {
601 Dist = distaux;
602 SecParam = taux;
603 AngleMax = alpha;
604 PonPath = P;
605 PathParam = Path->LastParameter();
606 }
607 }
608 }
609 Trouve = (Dist <= Tol);
610 }
611
612 // (2.2) Distance courbe-courbe
613 if (!Trouve) {
614 Extrema_ExtCC Ext(Path->Curve(), myAdpSection,
615 Path->FirstParameter(), Path->LastParameter(),
616 myAdpSection.FirstParameter(),
617 myAdpSection.LastParameter(),
618 Path->Resolution(Tol/100),
619 myAdpSection.Resolution(Tol/100));
620 if (Ext.IsDone()) {
621 Extrema_POnCurv P1, P2;
622 for (ii=1; ii<=Ext.NbExt(); ii++) {
623 distaux = sqrt (Ext.SquareDistance(ii));
624 Ext.Points(ii, P1, P2);
625 Tangente(Path->Curve(), P1.Parameter(), P, dp1);
626 alpha = EvalAngle(VRef, dp1);
627 if (Choix(distaux, alpha)) {
628 Trouve = Standard_True;
629 Dist = distaux;
630 PathParam = P1.Parameter();
631 SecParam = P2.Parameter();
632 PonSec = P2.Value();
633 PonPath = P;
634 AngleMax = alpha;
635 }
636 }
637 }
638 if (!Trouve){
639 // Si l'on a toujours rien, on essai une distance point/path
640 // c'est la derniere chance.
641 Extrema_ExtPC PExt;
642 PExt.Initialize(Path->Curve(),
643 Path->FirstParameter(),
644 Path->LastParameter(),
645 Precision::Confusion());
646 PExt.Perform(PonSec);
647 if ( PExt.IsDone() ) {
648 // modified for OCC13595 Tue Oct 17 15:00:08 2006.BEGIN
649 // DistMini(PExt, myAdpSection, distaux, taux);
650 DistMini(PExt, Path->Curve(), distaux, taux);
651 // modified for OCC13595 Tue Oct 17 15:00:11 2006.END
652 Tangente(Path->Curve(), taux, P, dp1);
653 alpha = EvalAngle(VRef, dp1);
654 if (Choix(distaux, alpha)) {
655 Dist = distaux;
656 PonPath = P;
657 AngleMax = alpha;
658 PathParam = taux;
659 }
660 }
661 }
662 }
663 }
664 }
665
666 done = Standard_True;
667}
668
669
670//===============================================================
671// Function : Perform
672// Purpose : Calcul le placement pour un parametre donne.
673//===============================================================
674void GeomFill_SectionPlacement::Perform(const Standard_Real Param,
675 const Standard_Real Tol)
676{
677 done = Standard_True;
678 Handle(Adaptor3d_HCurve) Path;
679 Path = myLaw->GetCurve();
680
681 PathParam = Param;
682 if (myIsPoint)
683 {
684 gp_Pnt PonPath = Path->Value( PathParam );
685 Dist = PonPath.Distance(myPoint);
c6541a0c 686 AngleMax = M_PI/2;
7fd59977 687 }
688 else
689 {
690 SecParam = myAdpSection.FirstParameter();
691
692 // Standard_Real distaux, taux, alpha;
693 // gp_Pnt PonPath, PonSec, P;
694 gp_Pnt PonPath, PonSec;
695 gp_Vec VRef, dp1;
696 VRef.SetXYZ(TheAxe.Direction().XYZ());
697
698 Tangente( Path->Curve(), PathParam, PonPath, dp1);
699 PonSec = myAdpSection.Value(SecParam);
700 Dist = PonPath.Distance(PonSec);
701 if (Dist > Tol) { // On Cherche un meilleur point sur la section
702 myExt.Perform(PonPath);
703 if ( myExt.IsDone() ) {
704 DistMini(myExt, myAdpSection, Dist, SecParam);
705 PonSec = myAdpSection.Value(SecParam);
706 }
707 }
708 AngleMax = EvalAngle(VRef, dp1);
c6541a0c 709 if (isplan) AngleMax = M_PI/2 - AngleMax;
7fd59977 710 }
711
712 done = Standard_True;
713}
714
715//===============================================================
716// Function :IsDone
717// Purpose :
718//===============================================================
719 Standard_Boolean GeomFill_SectionPlacement::IsDone() const
720{
721 return done;
722}
723
724//===============================================================
725// Function : ParameterOnPath
726// Purpose :
727//===============================================================
728 Standard_Real GeomFill_SectionPlacement::ParameterOnPath() const
729{
730 return PathParam;
731}
732
733//===============================================================
734// Function : ParameterOnSection
735// Purpose :
736//===============================================================
737 Standard_Real GeomFill_SectionPlacement::ParameterOnSection() const
738{
739 return SecParam;
740}
741
742//===============================================================
743// Function : Distance
744// Purpose :
745//===============================================================
746 Standard_Real GeomFill_SectionPlacement::Distance() const
747{
748 return Dist;
749}
750
751//===============================================================
752// Function : Angle
753// Purpose :
754//===============================================================
755 Standard_Real GeomFill_SectionPlacement::Angle() const
756{
757 return AngleMax;
758}
759
760//===============================================================
761// Function : Transformation
762// Purpose :
763//===============================================================
764 gp_Trsf GeomFill_SectionPlacement::Transformation(
765 const Standard_Boolean WithTranslation,
766 const Standard_Boolean WithCorrection) const
767{
768 gp_Vec V;
769 gp_Mat M;
770 gp_Dir DT, DN, D;
771// modified by NIZHNY-MKK Fri Oct 17 15:27:07 2003
772 gp_Pnt P(0., 0., 0.), PSection(0., 0., 0.);
773
774 // Calcul des reperes
775 myLaw->D0(PathParam, M, V);
776
777 P.SetXYZ(V.XYZ());
778 D.SetXYZ (M.Column(3));
779 DN.SetXYZ(M.Column(1));
780 gp_Ax3 Paxe(P, D, DN);
781
782
783 if (WithTranslation || WithCorrection) {
784 if (myIsPoint)
785 PSection = myPoint;
786 else
787 PSection = mySection->Value(SecParam);
788 }
789 // modified by NIZHNY-MKK Wed Oct 8 15:03:19 2003.BEGIN
790 gp_Trsf Rot;
791
792 if (WithCorrection && !myIsPoint) {
793 Standard_Real angle;
794
795 if (!isplan)
796 Standard_Failure::Raise("Illegal usage: can't rotate non-planar profile");
797
798 gp_Dir ProfileNormal = TheAxe.Direction();
799 gp_Dir SpineStartDir = Paxe.Direction();
800
801 if (! ProfileNormal.IsParallel( SpineStartDir, Precision::Angular() )) {
802 gp_Dir DirAxeOfRotation = ProfileNormal ^ SpineStartDir;
803 angle = ProfileNormal.AngleWithRef( SpineStartDir, DirAxeOfRotation );
804 gp_Ax1 AxeOfRotation( TheAxe.Location(), DirAxeOfRotation );
805 Rot.SetRotation( AxeOfRotation, angle );
806 }
807 PSection.Transform(Rot);
808 }
809 // modified by NIZHNY-MKK Wed Oct 8 15:03:21 2003.END
810
811 if (WithTranslation) {
812 P.ChangeCoord().SetLinearForm(-1,PSection.XYZ(),
813 V.XYZ());
814 }
815 else {
816 P.SetCoord(0., 0., 0.);
817 }
818
819 gp_Ax3 Saxe(P, gp::DZ(), gp::DX());
820
821 // Transfo
822 gp_Trsf Tf;
823 Tf.SetTransformation(Saxe, Paxe);
824
825 if (WithCorrection) {
826 // modified by NIZHNY-MKK Fri Oct 17 15:26:36 2003.BEGIN
827 // Standard_Real angle;
828 // gp_Trsf Rot;
829
830 // if (!isplan)
831 // Standard_Failure::Raise("Illegal usage: can't rotate non-planar profile");
832
833 // gp_Dir ProfileNormal = TheAxe.Direction();
834 // gp_Dir SpineStartDir = Paxe.Direction();
835 // if (! ProfileNormal.IsParallel( SpineStartDir, Precision::Angular() ))
836 // {
837 // gp_Dir DirAxeOfRotation = ProfileNormal ^ SpineStartDir;
838 // angle = ProfileNormal.AngleWithRef( SpineStartDir, DirAxeOfRotation );
839 // gp_Ax1 AxeOfRotation( TheAxe.Location(), DirAxeOfRotation );
840 // Rot.SetRotation( AxeOfRotation, angle );
841 // }
842 // modified by NIZHNY-MKK Fri Oct 17 15:26:42 2003.END
843
844 Tf *= Rot;
845 }
846
847 return Tf;
848}
849
850//===============================================================
851// Function : Section
852// Purpose :
853//===============================================================
854 Handle(Geom_Curve) GeomFill_SectionPlacement::
855 Section(const Standard_Boolean WithTranslation) const
856{
857 Handle(Geom_Curve) TheSection =
858 Handle(Geom_Curve)::DownCast(mySection->Copy());
859 TheSection->Transform(Transformation(WithTranslation, Standard_False));
860 return TheSection;
861}
862
863
864//===============================================================
865// Function :
866// Purpose :
867//===============================================================
868 Handle(Geom_Curve) GeomFill_SectionPlacement::
869 ModifiedSection(const Standard_Boolean WithTranslation) const
870{
871 Handle(Geom_Curve) TheSection =
872 Handle(Geom_Curve)::DownCast(mySection->Copy());
873 TheSection->Transform(Transformation(WithTranslation, Standard_True));
874 return TheSection;
875}
876
877//===============================================================
878// Function :SectionAxis
879// Purpose :
880//===============================================================
881 void GeomFill_SectionPlacement::SectionAxis(const gp_Mat& M,
882 gp_Vec& T,
883 gp_Vec& N,
884 gp_Vec& BN) const
885{
886 Standard_Real Eps = 1.e-10;
887 gp_Dir D;
888 gp_Vec PathNormal;
889 GeomLProp_CLProps CP(mySection, SecParam, 2, Eps);
890 if (CP.IsTangentDefined()) {
891 CP.Tangent(D);
892 T.SetXYZ(D.XYZ());
893 T.Normalize();
894 if (CP.Curvature() > Eps) {
895 CP.Normal(D);
896 N.SetXYZ(D.XYZ());
897 }
898 else {
899 // Cas ambigu, on essai de recuperer la normal a la trajectoire
900 // Reste le probleme des points d'inflexions, qui n'est pas
901 // bien traiter par LProp (pas de recuperation de la normal) :
902 // A voir...
903 PathNormal.SetXYZ(M.Column(1));
904 PathNormal.Normalize();
905 BN = T ^ PathNormal;
906 if (BN.Magnitude() > Eps) {
907 BN.Normalize();
908 //N = BN ^ T;
909 }
910 N = BN ^ T;
911 }
912 }
913 else { // Cas indefinie, on prend le triedre
914 // complet sur la trajectoire
915 T.SetXYZ(M.Column(3));
916 N.SetXYZ(M.Column(2));
917 }
918 BN = T ^ N;
919}
920
921
922//===============================================================
923// Function :Choix
924// Purpose : Decide si le couple (dist, angle) est "meilleur" que
925// couple courrant.
926//===============================================================
927 Standard_Boolean GeomFill_SectionPlacement::Choix(const Standard_Real dist,
928 const Standard_Real angle) const
929{
930 Standard_Real evoldist, evolangle;
931 evoldist = dist - Dist;
932 evolangle = angle - AngleMax;
933 // (1) Si la gain en distance est > que le gabarit, on prend
934 if (evoldist < - Gabarit)
935 return Standard_True;
936
937 // (2) si l'ecart en distance est de l'ordre du gabarit
938 if (Abs(evoldist) < Gabarit) {
939// (2.1) si le gain en angle est important on garde
940 if (evolangle > 0.5)
941 return Standard_True;
942 // (2.2) si la variation d'angle est moderee on evalue
943 // une fonction de penalite
944 if (Penalite(angle, dist/Gabarit) < Penalite(AngleMax, Dist/Gabarit) )
945 return Standard_True;
946 }
947
948 return Standard_False;
949}
950
951
952
953
954
955
956
957
Open CASCADE Technology repository