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