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