Commit | Line | Data |
---|---|---|
b311480e | 1 | // Created on: 1991-06-25 |
2 | // Created by: JCV | |
3 | // Copyright (c) 1991-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 | // modified by Edward AGAPOV (eap) Jan 28 2002 --- DN(), occ143(BUC60654) | |
18 | ||
19 | ||
20 | ||
21 | #include <Geom2d_OffsetCurve.ixx> | |
22 | #include <gp.hxx> | |
23 | #include <Standard_ConstructionError.hxx> | |
24 | #include <Standard_RangeError.hxx> | |
25 | #include <Standard_NotImplemented.hxx> | |
26 | #include <Geom2d_UndefinedDerivative.hxx> | |
27 | #include <Geom2d_UndefinedValue.hxx> | |
28 | #include <Geom2d_Line.hxx> | |
29 | #include <Geom2d_Circle.hxx> | |
30 | #include <Geom2d_Ellipse.hxx> | |
31 | #include <Geom2d_Hyperbola.hxx> | |
32 | #include <Geom2d_Parabola.hxx> | |
33 | #include <Geom2d_BezierCurve.hxx> | |
34 | #include <Geom2d_BSplineCurve.hxx> | |
35 | #include <Geom2d_TrimmedCurve.hxx> | |
36 | #include <gp_XY.hxx> | |
37 | ||
38 | typedef Handle(Geom2d_OffsetCurve) Handle(OffsetCurve); | |
39 | typedef Geom2d_OffsetCurve OffsetCurve; | |
40 | typedef Handle(Geom2d_Geometry) Handle(Geometry); | |
41 | typedef Handle(Geom2d_Curve) Handle(Curve); | |
42 | typedef Geom2d_Curve Curve; | |
43 | typedef gp_Dir2d Dir2d; | |
44 | typedef gp_Pnt2d Pnt2d; | |
45 | typedef gp_Vec2d Vec2d; | |
46 | typedef gp_Trsf2d Trsf2d; | |
47 | typedef gp_XY XY; | |
48 | ||
49 | ||
7fd59977 | 50 | //ordre de derivation maximum pour la recherche de la premiere |
51 | //derivee non nulle | |
32ca7a51 | 52 | static const int maxDerivOrder = 3; |
53 | static const Standard_Real MinStep = 1e-7; | |
7fd59977 | 54 | |
55 | ||
56 | //======================================================================= | |
57 | //function : Copy | |
58 | //purpose : | |
59 | //======================================================================= | |
60 | ||
61 | Handle(Geom2d_Geometry) Geom2d_OffsetCurve::Copy () const | |
62 | { | |
63 | Handle(OffsetCurve) C; | |
64 | C = new OffsetCurve (basisCurve, offsetValue); | |
65 | return C; | |
66 | } | |
67 | ||
68 | ||
69 | //======================================================================= | |
70 | //function : Geom2d_OffsetCurve | |
71 | //purpose : | |
72 | //======================================================================= | |
73 | ||
74 | Geom2d_OffsetCurve::Geom2d_OffsetCurve (const Handle(Curve)& C, | |
75 | const Standard_Real Offset) | |
76 | : offsetValue (Offset) | |
77 | { | |
78 | if (C->DynamicType() == STANDARD_TYPE(Geom2d_OffsetCurve)) { | |
e26b06c3 | 79 | Handle(OffsetCurve) OC = Handle(OffsetCurve)::DownCast(C); |
06be28a4 | 80 | SetBasisCurve (OC->BasisCurve()); |
7fd59977 | 81 | offsetValue += OC->Offset(); |
06be28a4 RL |
82 | } |
83 | else { | |
84 | SetBasisCurve(C); | |
7fd59977 | 85 | } |
86 | } | |
87 | ||
88 | //======================================================================= | |
89 | //function : Reverse | |
90 | //purpose : | |
91 | //======================================================================= | |
92 | ||
93 | void Geom2d_OffsetCurve::Reverse () | |
94 | { | |
95 | basisCurve->Reverse(); | |
96 | offsetValue = -offsetValue; | |
97 | } | |
98 | ||
99 | //======================================================================= | |
100 | //function : ReversedParameter | |
101 | //purpose : | |
102 | //======================================================================= | |
103 | ||
104 | Standard_Real Geom2d_OffsetCurve::ReversedParameter( const Standard_Real U) const | |
105 | { | |
106 | return basisCurve->ReversedParameter( U); | |
107 | } | |
108 | ||
109 | //======================================================================= | |
110 | //function : SetBasisCurve | |
111 | //purpose : | |
112 | //======================================================================= | |
113 | ||
114 | void Geom2d_OffsetCurve::SetBasisCurve (const Handle(Curve)& C) | |
115 | { | |
06be28a4 RL |
116 | Handle(Geom2d_Curve) aBasisCurve = Handle(Geom2d_Curve)::DownCast(C->Copy()); |
117 | ||
118 | // Basis curve must be at least C1 | |
119 | if (aBasisCurve->Continuity() == GeomAbs_C0) | |
120 | { | |
121 | // For B-splines it is sometimes possible to increase continuity by removing | |
122 | // unnecessarily duplicated knots | |
123 | if (aBasisCurve->IsKind(STANDARD_TYPE(Geom2d_BSplineCurve))) | |
124 | { | |
125 | Handle(Geom2d_BSplineCurve) aBCurve = Handle(Geom2d_BSplineCurve)::DownCast(aBasisCurve); | |
126 | Standard_Integer degree = aBCurve->Degree(); | |
127 | Standard_Real Toler = 1e-7; | |
128 | Standard_Integer start = aBCurve->IsPeriodic() ? 1 : aBCurve->FirstUKnotIndex(), | |
129 | finish = aBCurve->IsPeriodic() ? aBCurve->NbKnots() : aBCurve->LastUKnotIndex(); | |
130 | for (Standard_Integer i = start; i <= finish; i++) | |
131 | { | |
132 | Standard_Integer mult = aBCurve->Multiplicity(i); | |
133 | if ( mult == degree ) | |
134 | aBCurve->RemoveKnot(i,degree - 1, Toler); | |
135 | } | |
136 | } | |
137 | ||
138 | // Raise exception if still C0 | |
139 | if (aBasisCurve->Continuity() == GeomAbs_C0) | |
140 | Standard_ConstructionError::Raise("Offset on C0 curve"); | |
141 | } | |
142 | ||
143 | basisCurve = aBasisCurve; | |
7fd59977 | 144 | } |
145 | ||
146 | //======================================================================= | |
147 | //function : SetOffsetValue | |
148 | //purpose : | |
149 | //======================================================================= | |
150 | ||
151 | void Geom2d_OffsetCurve::SetOffsetValue (const Standard_Real D) { offsetValue = D; } | |
152 | ||
153 | //======================================================================= | |
154 | //function : BasisCurve | |
155 | //purpose : | |
156 | //======================================================================= | |
157 | ||
158 | Handle(Curve) Geom2d_OffsetCurve::BasisCurve () const | |
159 | { | |
160 | return basisCurve; | |
161 | } | |
162 | ||
163 | //======================================================================= | |
164 | //function : Continuity | |
165 | //purpose : | |
166 | //======================================================================= | |
167 | ||
168 | GeomAbs_Shape Geom2d_OffsetCurve::Continuity () const | |
169 | { | |
170 | GeomAbs_Shape OffsetShape=GeomAbs_C0; | |
171 | switch (basisCurve->Continuity()) { | |
172 | case GeomAbs_C0 : OffsetShape = GeomAbs_C0; break; | |
173 | case GeomAbs_C1 : OffsetShape = GeomAbs_C0; break; | |
174 | case GeomAbs_C2 : OffsetShape = GeomAbs_C1; break; | |
175 | case GeomAbs_C3 : OffsetShape = GeomAbs_C2; break; | |
176 | case GeomAbs_CN : OffsetShape = GeomAbs_CN; break; | |
177 | case GeomAbs_G1 : OffsetShape = GeomAbs_G1; break; | |
178 | case GeomAbs_G2 : OffsetShape = GeomAbs_G2; break; | |
179 | } | |
180 | return OffsetShape; | |
181 | } | |
182 | ||
183 | //======================================================================= | |
184 | //function : D0 | |
185 | //purpose : | |
186 | //======================================================================= | |
187 | ||
32ca7a51 | 188 | void Geom2d_OffsetCurve::D0 (const Standard_Real theU, |
189 | Pnt2d& theP ) const | |
190 | { | |
191 | const Standard_Real aTol = gp::Resolution(); | |
7fd59977 | 192 | |
32ca7a51 | 193 | Vec2d vD1; |
194 | ||
195 | basisCurve->D1 (theU, theP, vD1); | |
196 | Standard_Real Ndu = vD1.Magnitude(); | |
197 | ||
198 | if(Ndu <= aTol) | |
199 | { | |
200 | const Standard_Real anUinfium = basisCurve->FirstParameter(); | |
201 | const Standard_Real anUsupremum = basisCurve->LastParameter(); | |
202 | ||
203 | const Standard_Real DivisionFactor = 1.e-3; | |
204 | Standard_Real du; | |
205 | if((anUsupremum >= RealLast()) || (anUinfium <= RealFirst())) | |
206 | du = 0.0; | |
207 | else | |
208 | du = anUsupremum-anUinfium; | |
209 | ||
210 | const Standard_Real aDelta = Max(du*DivisionFactor,MinStep); | |
211 | //Derivative is approximated by Taylor-series | |
212 | ||
213 | Standard_Integer anIndex = 1; //Derivative order | |
214 | Vec2d V; | |
215 | ||
216 | do | |
217 | { | |
218 | V = basisCurve->DN(theU,++anIndex); | |
219 | Ndu = V.Magnitude(); | |
220 | } | |
221 | while((Ndu <= aTol) && anIndex < maxDerivOrder); | |
222 | ||
223 | Standard_Real u; | |
224 | ||
225 | if(theU-anUinfium < aDelta) | |
226 | u = theU+aDelta; | |
227 | else | |
228 | u = theU-aDelta; | |
229 | ||
230 | Pnt2d P1, P2; | |
231 | basisCurve->D0(Min(theU, u),P1); | |
232 | basisCurve->D0(Max(theU, u),P2); | |
233 | ||
234 | Vec2d V1(P1,P2); | |
235 | Standard_Real aDirFactor = V.Dot(V1); | |
236 | ||
237 | if(aDirFactor < 0.0) | |
238 | vD1 = -V; | |
239 | else | |
240 | vD1 = V; | |
241 | ||
242 | Ndu = vD1.Magnitude(); | |
243 | }//if(Ndu <= aTol) | |
244 | ||
245 | if (Ndu <= aTol) | |
246 | Geom2d_UndefinedValue::Raise("Exception: Undefined normal vector " | |
247 | "because tangent vector has zero-magnitude!"); | |
248 | ||
249 | Standard_Real A = vD1.Y(); | |
250 | Standard_Real B = - vD1.X(); | |
251 | A = A * offsetValue/Ndu; | |
252 | B = B * offsetValue/Ndu; | |
253 | theP.SetCoord(theP.X() + A, theP.Y() + B); | |
7fd59977 | 254 | } |
7fd59977 | 255 | |
256 | //======================================================================= | |
257 | //function : D1 | |
258 | //purpose : | |
259 | //======================================================================= | |
32ca7a51 | 260 | void Geom2d_OffsetCurve::D1 (const Standard_Real theU, Pnt2d& P, Vec2d& theV1) const |
261 | { | |
7fd59977 | 262 | // P(u) = p(u) + Offset * Ndir / R |
263 | // with R = || p' ^ Z|| and Ndir = P' ^ Z | |
264 | ||
265 | // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R - Ndir * (DR/R)) | |
266 | ||
32ca7a51 | 267 | const Standard_Real aTol = gp::Resolution(); |
268 | ||
7fd59977 | 269 | Vec2d V2; |
32ca7a51 | 270 | basisCurve->D2 (theU, P, theV1, V2); |
271 | ||
272 | if(theV1.Magnitude() <= aTol) | |
273 | { | |
274 | const Standard_Real anUinfium = basisCurve->FirstParameter(); | |
275 | const Standard_Real anUsupremum = basisCurve->LastParameter(); | |
276 | ||
277 | const Standard_Real DivisionFactor = 1.e-3; | |
278 | Standard_Real du; | |
279 | if((anUsupremum >= RealLast()) || (anUinfium <= RealFirst())) | |
280 | du = 0.0; | |
281 | else | |
282 | du = anUsupremum-anUinfium; | |
283 | ||
284 | const Standard_Real aDelta = Max(du*DivisionFactor,MinStep); | |
285 | //Derivative is approximated by Taylor-series | |
286 | ||
287 | Standard_Integer anIndex = 1; //Derivative order | |
288 | Vec2d V; | |
289 | ||
290 | do | |
291 | { | |
292 | V = basisCurve->DN(theU,++anIndex); | |
293 | } | |
294 | while((V.Magnitude() <= aTol) && anIndex < maxDerivOrder); | |
295 | ||
296 | Standard_Real u; | |
297 | ||
298 | if(theU-anUinfium < aDelta) | |
299 | u = theU+aDelta; | |
300 | else | |
301 | u = theU-aDelta; | |
302 | ||
303 | Pnt2d P1, P2; | |
304 | basisCurve->D0(Min(theU, u),P1); | |
305 | basisCurve->D0(Max(theU, u),P2); | |
306 | ||
307 | Vec2d V1(P1,P2); | |
308 | Standard_Real aDirFactor = V.Dot(V1); | |
309 | ||
310 | if(aDirFactor < 0.0) | |
311 | { | |
312 | theV1 = -V; | |
313 | V2 = - basisCurve->DN (theU, anIndex+1); | |
314 | } | |
315 | else | |
316 | { | |
317 | theV1 = V; | |
318 | V2 = basisCurve->DN (theU, anIndex+1); | |
319 | } | |
320 | }//if(theV1.Magnitude() <= aTol) | |
321 | ||
322 | XY Ndir (theV1.Y(), -theV1.X()); | |
7fd59977 | 323 | XY DNdir (V2.Y(), -V2.X()); |
324 | Standard_Real R2 = Ndir.SquareModulus(); | |
325 | Standard_Real R = Sqrt (R2); | |
326 | Standard_Real R3 = R * R2; | |
327 | Standard_Real Dr = Ndir.Dot (DNdir); | |
328 | if (R3 <= gp::Resolution()) { | |
329 | //We try another computation but the stability is not very good. | |
330 | if (R2 <= gp::Resolution()) Geom2d_UndefinedDerivative::Raise(); | |
331 | DNdir.Multiply(R); | |
332 | DNdir.Subtract (Ndir.Multiplied (Dr/R)); | |
333 | DNdir.Multiply (offsetValue/R2); | |
32ca7a51 | 334 | theV1.Add (Vec2d(DNdir)); |
7fd59977 | 335 | } |
336 | else { | |
337 | // Same computation as IICURV in EUCLID-IS because the stability is | |
338 | // better | |
339 | DNdir.Multiply (offsetValue/R); | |
340 | DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3)); | |
32ca7a51 | 341 | theV1.Add (Vec2d(DNdir)); |
7fd59977 | 342 | } |
32ca7a51 | 343 | |
344 | D0(theU, P); | |
7fd59977 | 345 | } |
346 | ||
347 | //======================================================================= | |
348 | //function : D2 | |
349 | //purpose : | |
350 | //======================================================================= | |
351 | ||
32ca7a51 | 352 | void Geom2d_OffsetCurve::D2 (const Standard_Real theU, |
7fd59977 | 353 | Pnt2d& P, |
32ca7a51 | 354 | Vec2d& theV1, Vec2d& V2) const |
7fd59977 | 355 | { |
356 | // P(u) = p(u) + Offset * Ndir / R | |
357 | // with R = || p' ^ Z|| and Ndir = P' ^ Z | |
358 | ||
359 | // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R - Ndir * (DR/R)) | |
360 | ||
361 | // P"(u) = p"(u) + (Offset / R) * (D2Ndir/DU - DNdir * (2.0 * Dr/ R**2) + | |
362 | // Ndir * ( (3.0 * Dr**2 / R**4) - (D2r / R**2))) | |
363 | ||
7fd59977 | 364 | Vec2d V3; |
32ca7a51 | 365 | basisCurve->D3 (theU, P, theV1, V2, V3); |
366 | ||
367 | const Standard_Real aTol = gp::Resolution(); | |
368 | ||
369 | Standard_Boolean IsDirectionChange = Standard_False; | |
370 | ||
371 | if(theV1.Magnitude() <= aTol) | |
372 | { | |
373 | const Standard_Real anUinfium = basisCurve->FirstParameter(); | |
374 | const Standard_Real anUsupremum = basisCurve->LastParameter(); | |
375 | ||
376 | const Standard_Real DivisionFactor = 1.e-3; | |
377 | Standard_Real du; | |
378 | if((anUsupremum >= RealLast()) || (anUinfium <= RealFirst())) | |
379 | du = 0.0; | |
380 | else | |
381 | du = anUsupremum-anUinfium; | |
382 | ||
383 | const Standard_Real aDelta = Max(du*DivisionFactor,MinStep); | |
384 | //Derivative is approximated by Taylor-series | |
385 | ||
386 | Standard_Integer anIndex = 1; //Derivative order | |
387 | Vec2d V; | |
388 | ||
389 | do | |
390 | { | |
391 | V = basisCurve->DN(theU,++anIndex); | |
392 | } | |
393 | while((V.Magnitude() <= aTol) && anIndex < maxDerivOrder); | |
394 | ||
395 | Standard_Real u; | |
396 | ||
397 | if(theU-anUinfium < aDelta) | |
398 | u = theU+aDelta; | |
399 | else | |
400 | u = theU-aDelta; | |
401 | ||
402 | Pnt2d P1, P2; | |
403 | basisCurve->D0(Min(theU, u),P1); | |
404 | basisCurve->D0(Max(theU, u),P2); | |
405 | ||
406 | Vec2d V1(P1,P2); | |
407 | Standard_Real aDirFactor = V.Dot(V1); | |
408 | ||
409 | if(aDirFactor < 0.0) | |
410 | { | |
411 | theV1 = -V; | |
412 | V2 = -basisCurve->DN (theU, anIndex+1); | |
413 | V3 = -basisCurve->DN (theU, anIndex + 2); | |
414 | ||
415 | IsDirectionChange = Standard_True; | |
416 | } | |
417 | else | |
418 | { | |
419 | theV1 = V; | |
420 | V2 = basisCurve->DN (theU, anIndex+1); | |
421 | V3 = basisCurve->DN (theU, anIndex + 2); | |
422 | } | |
423 | }//if(V1.Magnitude() <= aTol) | |
424 | ||
425 | XY Ndir (theV1.Y(), -theV1.X()); | |
7fd59977 | 426 | XY DNdir (V2.Y(), -V2.X()); |
427 | XY D2Ndir (V3.Y(), -V3.X()); | |
428 | Standard_Real R2 = Ndir.SquareModulus(); | |
429 | Standard_Real R = Sqrt (R2); | |
430 | Standard_Real R3 = R2 * R; | |
431 | Standard_Real R4 = R2 * R2; | |
432 | Standard_Real R5 = R3 * R2; | |
433 | Standard_Real Dr = Ndir.Dot (DNdir); | |
434 | Standard_Real D2r = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir); | |
32ca7a51 | 435 | if (R5 <= gp::Resolution()) |
436 | { | |
7fd59977 | 437 | //We try another computation but the stability is not very good |
438 | //dixit ISG. | |
32ca7a51 | 439 | if (R4 <= gp::Resolution()) |
440 | { | |
441 | Geom2d_UndefinedDerivative::Raise(); | |
442 | } | |
443 | // V2 = P" (U) : | |
444 | Standard_Real R4 = R2 * R2; | |
445 | D2Ndir.Subtract (DNdir.Multiplied (2.0 * Dr / R2)); | |
446 | D2Ndir.Add (Ndir.Multiplied (((3.0 * Dr * Dr)/R4) - (D2r/R2))); | |
447 | D2Ndir.Multiply (offsetValue / R); | |
448 | ||
449 | if(IsDirectionChange) | |
450 | V2=-V2; | |
451 | ||
452 | V2.Add (Vec2d(D2Ndir)); | |
453 | ||
454 | // V1 = P' (U) : | |
455 | DNdir.Multiply(R); | |
456 | DNdir.Subtract (Ndir.Multiplied (Dr/R)); | |
457 | DNdir.Multiply (offsetValue/R2); | |
458 | theV1.Add (Vec2d(DNdir)); | |
459 | } | |
460 | else | |
461 | { | |
462 | // Same computation as IICURV in EUCLID-IS because the stability is | |
463 | // better. | |
7fd59977 | 464 | // V2 = P" (U) : |
7fd59977 | 465 | D2Ndir.Multiply (offsetValue/R); |
466 | D2Ndir.Subtract (DNdir.Multiplied (2.0 * offsetValue * Dr / R3)); | |
467 | D2Ndir.Add (Ndir.Multiplied | |
468 | (offsetValue * (((3.0 * Dr * Dr) / R5) - (D2r / R3)))); | |
32ca7a51 | 469 | |
470 | if(IsDirectionChange) | |
471 | V2=-V2; | |
472 | ||
7fd59977 | 473 | V2.Add (Vec2d(D2Ndir)); |
32ca7a51 | 474 | |
7fd59977 | 475 | // V1 = P' (U) |
32ca7a51 | 476 | DNdir.Multiply (offsetValue/R); |
477 | DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3)); | |
478 | theV1.Add (Vec2d(DNdir)); | |
479 | } | |
480 | ||
481 | //P (U) : | |
482 | D0(theU, P); | |
7fd59977 | 483 | } |
484 | ||
485 | ||
486 | //======================================================================= | |
487 | //function : D3 | |
488 | //purpose : | |
489 | //======================================================================= | |
490 | ||
32ca7a51 | 491 | void Geom2d_OffsetCurve::D3 (const Standard_Real theU, |
7fd59977 | 492 | Pnt2d& P, |
32ca7a51 | 493 | Vec2d& theV1, Vec2d& V2, Vec2d& V3) const { |
7fd59977 | 494 | |
495 | ||
496 | // P(u) = p(u) + Offset * Ndir / R | |
497 | // with R = || p' ^ Z|| and Ndir = P' ^ Z | |
498 | ||
499 | // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R - Ndir * (DR/R)) | |
500 | ||
501 | // P"(u) = p"(u) + (Offset / R) * (D2Ndir/DU - DNdir * (2.0 * Dr/ R**2) + | |
502 | // Ndir * ( (3.0 * Dr**2 / R**4) - (D2r / R**2))) | |
503 | ||
504 | //P"'(u) = p"'(u) + (Offset / R) * (D3Ndir - (3.0 * Dr/R**2 ) * D2Ndir - | |
505 | // (3.0 * D2r / R2) * DNdir) + (3.0 * Dr * Dr / R4) * DNdir - | |
506 | // (D3r/R2) * Ndir + (6.0 * Dr * Dr / R4) * Ndir + | |
507 | // (6.0 * Dr * D2r / R4) * Ndir - (15.0 * Dr* Dr* Dr /R6) * Ndir | |
508 | ||
32ca7a51 | 509 | const Standard_Real aTol = gp::Resolution(); |
7fd59977 | 510 | |
32ca7a51 | 511 | Standard_Boolean IsDirectionChange = Standard_False; |
7fd59977 | 512 | |
32ca7a51 | 513 | basisCurve->D3 (theU, P, theV1, V2, V3); |
514 | Vec2d V4 = basisCurve->DN (theU, 4); | |
515 | ||
516 | if(theV1.Magnitude() <= aTol) | |
517 | { | |
518 | const Standard_Real anUinfium = basisCurve->FirstParameter(); | |
519 | const Standard_Real anUsupremum = basisCurve->LastParameter(); | |
520 | ||
521 | const Standard_Real DivisionFactor = 1.e-3; | |
522 | Standard_Real du; | |
523 | if((anUsupremum >= RealLast()) || (anUinfium <= RealFirst())) | |
524 | du = 0.0; | |
525 | else | |
526 | du = anUsupremum-anUinfium; | |
527 | ||
528 | const Standard_Real aDelta = Max(du*DivisionFactor,MinStep); | |
529 | //Derivative is approximated by Taylor-series | |
530 | ||
531 | Standard_Integer anIndex = 1; //Derivative order | |
532 | Vec2d V; | |
533 | ||
534 | do | |
535 | { | |
536 | V = basisCurve->DN(theU,++anIndex); | |
537 | } | |
538 | while((V.Magnitude() <= aTol) && anIndex < maxDerivOrder); | |
539 | ||
540 | Standard_Real u; | |
541 | ||
542 | if(theU-anUinfium < aDelta) | |
543 | u = theU+aDelta; | |
544 | else | |
545 | u = theU-aDelta; | |
546 | ||
547 | Pnt2d P1, P2; | |
548 | basisCurve->D0(Min(theU, u),P1); | |
549 | basisCurve->D0(Max(theU, u),P2); | |
550 | ||
551 | Vec2d V1(P1,P2); | |
552 | Standard_Real aDirFactor = V.Dot(V1); | |
553 | ||
554 | if(aDirFactor < 0.0) | |
555 | { | |
556 | theV1 = -V; | |
557 | V2 = -basisCurve->DN (theU, anIndex + 1); | |
558 | V3 = -basisCurve->DN (theU, anIndex + 2); | |
559 | V4 = -basisCurve->DN (theU, anIndex + 3); | |
560 | ||
561 | IsDirectionChange = Standard_True; | |
562 | } | |
563 | else | |
564 | { | |
565 | theV1 = V; | |
566 | V2 = basisCurve->DN (theU, anIndex + 1); | |
567 | V3 = basisCurve->DN (theU, anIndex + 2); | |
568 | V4 = basisCurve->DN (theU, anIndex + 3); | |
569 | } | |
570 | }//if(V1.Magnitude() <= aTol) | |
571 | ||
572 | XY Ndir (theV1.Y(), -theV1.X()); | |
573 | XY DNdir (V2.Y(), -V2.X()); | |
574 | XY D2Ndir (V3.Y(), -V3.X()); | |
575 | XY D3Ndir (V4.Y(), -V4.X()); | |
576 | Standard_Real R2 = Ndir.SquareModulus(); | |
577 | Standard_Real R = Sqrt (R2); | |
578 | Standard_Real R3 = R2 * R; | |
579 | Standard_Real R4 = R2 * R2; | |
580 | Standard_Real R5 = R3 * R2; | |
581 | Standard_Real R6 = R3 * R3; | |
582 | Standard_Real R7 = R5 * R2; | |
583 | Standard_Real Dr = Ndir.Dot (DNdir); | |
584 | Standard_Real D2r = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir); | |
585 | Standard_Real D3r = Ndir.Dot (D3Ndir) + 3.0 * DNdir.Dot (D2Ndir); | |
586 | ||
587 | if (R7 <= gp::Resolution()) | |
588 | { | |
589 | //We try another computation but the stability is not very good | |
590 | //dixit ISG. | |
591 | ||
592 | if (R6 <= gp::Resolution()) | |
593 | Geom2d_UndefinedDerivative::Raise(); | |
594 | ||
595 | // V3 = P"' (U) : | |
596 | D3Ndir.Subtract (D2Ndir.Multiplied (3.0 * offsetValue * Dr / R2)); | |
597 | D3Ndir.Subtract ( | |
598 | (DNdir.Multiplied ((3.0 * offsetValue) * ((D2r/R2) + (Dr*Dr)/R4)))); | |
599 | D3Ndir.Add (Ndir.Multiplied ( | |
600 | (offsetValue * (6.0*Dr*Dr/R4 + 6.0*Dr*D2r/R4 - 15.0*Dr*Dr*Dr/R6 - D3r)))); | |
601 | D3Ndir.Multiply (offsetValue/R); | |
602 | ||
603 | if(IsDirectionChange) | |
604 | V3=-V3; | |
605 | ||
606 | V3.Add (Vec2d(D3Ndir)); | |
607 | ||
608 | ||
609 | // V2 = P" (U) : | |
610 | Standard_Real R4 = R2 * R2; | |
611 | D2Ndir.Subtract (DNdir.Multiplied (2.0 * Dr / R2)); | |
612 | D2Ndir.Subtract (Ndir.Multiplied (((3.0 * Dr * Dr)/R4) - (D2r/R2))); | |
613 | D2Ndir.Multiply (offsetValue / R); | |
614 | V2.Add (Vec2d(D2Ndir)); | |
615 | // V1 = P' (U) : | |
616 | DNdir.Multiply(R); | |
617 | DNdir.Subtract (Ndir.Multiplied (Dr/R)); | |
618 | DNdir.Multiply (offsetValue/R2); | |
619 | theV1.Add (Vec2d(DNdir)); | |
620 | } | |
621 | else | |
622 | { | |
623 | // Same computation as IICURV in EUCLID-IS because the stability is | |
624 | // better. | |
625 | // V3 = P"' (U) : | |
626 | D3Ndir.Multiply (offsetValue/R); | |
627 | D3Ndir.Subtract (D2Ndir.Multiplied (3.0 * offsetValue * Dr / R3)); | |
628 | D3Ndir.Subtract (DNdir.Multiplied ( | |
629 | ((3.0 * offsetValue) * ((D2r/R3) + (Dr*Dr)/R5))) ); | |
630 | D3Ndir.Add (Ndir.Multiplied ( | |
631 | (offsetValue * (6.0*Dr*Dr/R5 + 6.0*Dr*D2r/R5 - 15.0*Dr*Dr*Dr/R7 - D3r)))); | |
632 | ||
633 | if(IsDirectionChange) | |
634 | V3=-V3; | |
635 | ||
636 | V3.Add (Vec2d(D3Ndir)); | |
637 | ||
638 | // V2 = P" (U) : | |
639 | D2Ndir.Multiply (offsetValue/R); | |
640 | D2Ndir.Subtract (DNdir.Multiplied (2.0 * offsetValue * Dr / R3)); | |
641 | D2Ndir.Subtract (Ndir.Multiplied ( | |
642 | offsetValue * (((3.0 * Dr * Dr) / R5) - (D2r / R3)))); | |
643 | V2.Add (Vec2d(D2Ndir)); | |
644 | // V1 = P' (U) : | |
645 | DNdir.Multiply (offsetValue/R); | |
646 | DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3)); | |
647 | theV1.Add (Vec2d(DNdir)); | |
648 | } | |
649 | //P (U) : | |
650 | D0(theU, P); | |
651 | } | |
7fd59977 | 652 | |
653 | //======================================================================= | |
654 | //function : DN | |
655 | //purpose : | |
656 | //======================================================================= | |
657 | ||
658 | Vec2d Geom2d_OffsetCurve::DN (const Standard_Real U, | |
659 | const Standard_Integer N) const | |
660 | { | |
32ca7a51 | 661 | Standard_RangeError_Raise_if (N < 1, "Exception: Geom2d_OffsetCurve::DN(). N<1."); |
7fd59977 | 662 | |
663 | gp_Vec2d VN, VBidon; | |
664 | gp_Pnt2d PBidon; | |
665 | switch (N) { | |
666 | case 1: D1( U, PBidon, VN); break; | |
667 | case 2: D2( U, PBidon, VBidon, VN); break; | |
668 | case 3: D3( U, PBidon, VBidon, VBidon, VN); break; | |
669 | default: | |
32ca7a51 | 670 | Standard_NotImplemented::Raise("Exception: Derivative order is greater than 3. " |
671 | "Cannot compute of derivative."); | |
7fd59977 | 672 | } |
673 | ||
674 | return VN; | |
675 | } | |
676 | ||
677 | ||
678 | //======================================================================= | |
679 | //function : Value | |
680 | //purpose : | |
681 | //======================================================================= | |
682 | ||
32ca7a51 | 683 | void Geom2d_OffsetCurve::Value (const Standard_Real theU, |
684 | Pnt2d& theP, Pnt2d& thePbasis, | |
685 | Vec2d& theV1basis ) const | |
686 | { | |
687 | basisCurve->D1(theU, thePbasis, theV1basis); | |
688 | D0(theU,theP); | |
7fd59977 | 689 | } |
7fd59977 | 690 | |
691 | ||
692 | //======================================================================= | |
693 | //function : D1 | |
694 | //purpose : | |
695 | //======================================================================= | |
696 | ||
697 | void Geom2d_OffsetCurve::D1 (const Standard_Real U, | |
698 | Pnt2d& P, Pnt2d& Pbasis, | |
699 | Vec2d& V1, Vec2d& V1basis, | |
700 | Vec2d& V2basis ) const | |
701 | { | |
702 | // P(u) = p(u) + Offset * Ndir / R | |
703 | // with R = || p' ^ Z|| and Ndir = P' ^ Z | |
704 | ||
705 | // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R - Ndir * (DR/R)) | |
706 | ||
7fd59977 | 707 | basisCurve->D2 (U, Pbasis, V1basis, V2basis); |
708 | V1 = V1basis; | |
709 | Vec2d V2 = V2basis; | |
710 | Standard_Integer Index = 2; | |
32ca7a51 | 711 | while (V1.Magnitude() <= gp::Resolution() && Index <= maxDerivOrder) { |
7fd59977 | 712 | V1 = basisCurve->DN (U, Index); |
713 | Index++; | |
714 | } | |
715 | if (Index != 2) { | |
716 | V2 = basisCurve->DN (U, Index); | |
717 | } | |
718 | XY Ndir (V1.Y(), -V1.X()); | |
719 | XY DNdir (V2.Y(), -V2.X()); | |
720 | Standard_Real R2 = Ndir.SquareModulus(); | |
721 | Standard_Real R = Sqrt (R2); | |
722 | Standard_Real R3 = R * R2; | |
723 | Standard_Real Dr = Ndir.Dot (DNdir); | |
724 | if (R3 <= gp::Resolution()) { | |
725 | //We try another computation but the stability is not very good. | |
726 | if (R2 <= gp::Resolution()) { Geom2d_UndefinedDerivative::Raise(); } | |
727 | DNdir.Multiply(R); | |
728 | DNdir.Subtract (Ndir.Multiplied (Dr/R)); | |
729 | DNdir.Multiply (offsetValue / R2); | |
730 | V1.Add (Vec2d(DNdir)); | |
731 | } | |
732 | else { | |
733 | // Same computation as IICURV in EUCLID-IS because the stability is | |
734 | // better | |
735 | DNdir.Multiply (offsetValue/R); | |
736 | DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3)); | |
737 | V1.Add (Vec2d(DNdir)); | |
738 | } | |
739 | Ndir.Multiply (offsetValue/R); | |
740 | Ndir.Add (Pbasis.XY()); | |
741 | P.SetXY (Ndir); | |
742 | } | |
743 | ||
744 | ||
745 | //======================================================================= | |
746 | //function : D2 | |
747 | //purpose : | |
748 | //======================================================================= | |
749 | ||
750 | void Geom2d_OffsetCurve::D2 (const Standard_Real U, | |
751 | Pnt2d& P, Pnt2d& Pbasis, | |
752 | Vec2d& V1, Vec2d& V2, | |
753 | Vec2d& V1basis, Vec2d& V2basis, | |
754 | Vec2d& V3basis ) const | |
755 | { | |
756 | // P(u) = p(u) + Offset * Ndir / R | |
757 | // with R = || p' ^ Z|| and Ndir = P' ^ Z | |
758 | ||
759 | // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R - Ndir * (DR/R)) | |
760 | ||
761 | // P"(u) = p"(u) + (Offset / R) * (D2Ndir/DU - DNdir * (2.0 * Dr/ R**2) + | |
762 | // Ndir * ( (3.0 * Dr**2 / R**4) - (D2r / R**2))) | |
763 | ||
7fd59977 | 764 | basisCurve->D3 (U, Pbasis, V1basis, V2basis, V3basis); |
765 | Standard_Integer Index = 2; | |
766 | V1 = V1basis; | |
767 | V2 = V2basis; | |
768 | Vec2d V3 = V3basis; | |
32ca7a51 | 769 | while (V1.Magnitude() <= gp::Resolution() && Index <= maxDerivOrder) { |
7fd59977 | 770 | V1 = basisCurve->DN (U, Index); |
771 | Index++; | |
772 | } | |
773 | if (Index != 2) { | |
774 | V2 = basisCurve->DN (U, Index); | |
775 | V3 = basisCurve->DN (U, Index + 1); | |
776 | } | |
777 | XY Ndir (V1.Y(), -V1.X()); | |
778 | XY DNdir (V2.Y(), -V2.X()); | |
779 | XY D2Ndir (V3.Y(), -V3.X()); | |
780 | Standard_Real R2 = Ndir.SquareModulus(); | |
781 | Standard_Real R = Sqrt (R2); | |
782 | Standard_Real R3 = R2 * R; | |
783 | Standard_Real R4 = R2 * R2; | |
784 | Standard_Real R5 = R3 * R2; | |
785 | Standard_Real Dr = Ndir.Dot (DNdir); | |
786 | Standard_Real D2r = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir); | |
787 | if (R5 <= gp::Resolution()) { | |
788 | //We try another computation but the stability is not very good | |
789 | //dixit ISG. | |
790 | if (R4 <= gp::Resolution()) { Geom2d_UndefinedDerivative::Raise(); } | |
791 | // V2 = P" (U) : | |
792 | Standard_Real R4 = R2 * R2; | |
793 | D2Ndir.Subtract (DNdir.Multiplied (2.0 * Dr / R2)); | |
794 | D2Ndir.Subtract (Ndir.Multiplied (((3.0 * Dr * Dr)/R4) - (D2r/R2))); | |
795 | D2Ndir.Multiply (offsetValue / R); | |
796 | V2.Add (Vec2d(D2Ndir)); | |
797 | // V1 = P' (U) : | |
798 | DNdir.Multiply(R); | |
799 | DNdir.Subtract (Ndir.Multiplied (Dr/R)); | |
800 | DNdir.Multiply (offsetValue/R2); | |
801 | V1.Add (Vec2d(DNdir)); | |
802 | } | |
803 | else { | |
804 | // Same computation as IICURV in EUCLID-IS because the stability is | |
805 | // better. | |
806 | // V2 = P" (U) : | |
807 | D2Ndir.Multiply (offsetValue/R); | |
808 | D2Ndir.Subtract (DNdir.Multiplied (2.0 * offsetValue * Dr / R3)); | |
809 | D2Ndir.Subtract (Ndir.Multiplied ( | |
810 | offsetValue * (((3.0 * Dr * Dr) / R5) - (D2r / R3)) | |
811 | ) | |
812 | ); | |
813 | V2.Add (Vec2d(D2Ndir)); | |
814 | // V1 = P' (U) : | |
815 | DNdir.Multiply (offsetValue/R); | |
816 | DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3)); | |
817 | V1.Add (Vec2d(DNdir)); | |
818 | } | |
819 | //P (U) : | |
820 | Ndir.Multiply (offsetValue/R); | |
821 | Ndir.Add (Pbasis.XY()); | |
822 | P.SetXY (Ndir); | |
823 | } | |
824 | ||
825 | //======================================================================= | |
826 | //function : FirstParameter | |
827 | //purpose : | |
828 | //======================================================================= | |
829 | ||
830 | Standard_Real Geom2d_OffsetCurve::FirstParameter () const | |
831 | { | |
832 | return basisCurve->FirstParameter(); | |
833 | } | |
834 | ||
835 | //======================================================================= | |
836 | //function : LastParameter | |
837 | //purpose : | |
838 | //======================================================================= | |
839 | ||
840 | Standard_Real Geom2d_OffsetCurve::LastParameter () const | |
841 | { | |
842 | return basisCurve->LastParameter(); | |
843 | } | |
844 | ||
845 | ||
846 | //======================================================================= | |
847 | //function : Offset | |
848 | //purpose : | |
849 | //======================================================================= | |
850 | ||
851 | Standard_Real Geom2d_OffsetCurve::Offset () const { return offsetValue; } | |
852 | ||
853 | //======================================================================= | |
854 | //function : IsClosed | |
855 | //purpose : | |
856 | //======================================================================= | |
857 | ||
858 | Standard_Boolean Geom2d_OffsetCurve::IsClosed () const | |
859 | { | |
860 | gp_Pnt2d PF, PL; | |
861 | D0(FirstParameter(),PF); | |
862 | D0(LastParameter(),PL); | |
863 | return ( PF.Distance(PL) <= gp::Resolution()); | |
864 | } | |
865 | ||
866 | //======================================================================= | |
867 | //function : IsCN | |
868 | //purpose : | |
869 | //======================================================================= | |
870 | ||
871 | Standard_Boolean Geom2d_OffsetCurve::IsCN (const Standard_Integer N) const | |
872 | { | |
873 | Standard_RangeError_Raise_if (N < 0, " " ); | |
874 | return basisCurve->IsCN (N + 1); | |
875 | } | |
876 | ||
877 | //======================================================================= | |
878 | //function : IsPeriodic | |
879 | //purpose : | |
880 | //======================================================================= | |
881 | ||
882 | Standard_Boolean Geom2d_OffsetCurve::IsPeriodic () const | |
883 | { | |
884 | return basisCurve->IsPeriodic(); | |
885 | } | |
886 | ||
887 | //======================================================================= | |
888 | //function : Period | |
889 | //purpose : | |
890 | //======================================================================= | |
891 | ||
892 | Standard_Real Geom2d_OffsetCurve::Period() const | |
893 | { | |
894 | return basisCurve->Period(); | |
895 | } | |
896 | ||
897 | //======================================================================= | |
898 | //function : Transform | |
899 | //purpose : | |
900 | //======================================================================= | |
901 | ||
902 | void Geom2d_OffsetCurve::Transform (const Trsf2d& T) | |
903 | { | |
904 | basisCurve->Transform (T); | |
905 | offsetValue *= Abs(T.ScaleFactor()); | |
906 | } | |
907 | ||
908 | ||
909 | //======================================================================= | |
910 | //function : TransformedParameter | |
911 | //purpose : | |
912 | //======================================================================= | |
913 | ||
914 | Standard_Real Geom2d_OffsetCurve::TransformedParameter(const Standard_Real U, | |
915 | const gp_Trsf2d& T) const | |
916 | { | |
917 | return basisCurve->TransformedParameter(U,T); | |
918 | } | |
919 | ||
920 | //======================================================================= | |
921 | //function : ParametricTransformation | |
922 | //purpose : | |
923 | //======================================================================= | |
924 | ||
925 | Standard_Real Geom2d_OffsetCurve::ParametricTransformation(const gp_Trsf2d& T) const | |
926 | { | |
927 | return basisCurve->ParametricTransformation(T); | |
928 | } |