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