src/Geom/Geom_OffsetCurve.cxx

   1 // Created on: 1991-06-25
   2 // Created by: JCV
   3 // Copyright (c) 1991-1999 Matra Datavision
   4 // Copyright (c) 1999-2012 OPEN CASCADE SAS
   5 //
   6 // The content of this file is subject to the Open CASCADE Technology Public
   7 // License Version 6.5 (the "License"). You may not use the content of this file
   8 // except in compliance with the License. Please obtain a copy of the License
   9 // at http://www.opencascade.org and read it completely before using this file.
  10 //
  11 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
  12 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
  13 //
  14 // The Original Code and all software distributed under the License is
  15 // distributed on an "AS IS" basis, without warranty of any kind, and the
  16 // Initial Developer hereby disclaims all such warranties, including without
  17 // limitation, any warranties of merchantability, fitness for a particular
  18 // purpose or non-infringement. Please see the License for the specific terms
  19 // and conditions governing the rights and limitations under the License.
  20
  21
  22 // 24-Aug-95 : xab removed C1 and C2 test : appeller  D1 et D2
  23 //             avec discernement !
  24 // 19-09-97  : JPI correction derivee seconde
  25
  26
  27 #include <Geom_OffsetCurve.ixx>
  28
  29
  30 #include <gp.hxx>
  31 #include <gp_XYZ.hxx>
  32 #include <Geom_Line.hxx>
  33 #include <Geom_Circle.hxx>
  34 #include <Geom_Ellipse.hxx>
  35 #include <Geom_Hyperbola.hxx>
  36 #include <Geom_Parabola.hxx>
  37 #include <Geom_BezierCurve.hxx>
  38 #include <Geom_TrimmedCurve.hxx>
  39 #include <Geom_BSplineCurve.hxx>
  40 #include <Geom_Geometry.hxx>
  41
  42 #include <Geom_UndefinedDerivative.hxx>
  43 #include <Geom_UndefinedValue.hxx>
  44 #include <Standard_ConstructionError.hxx>
  45 #include <Standard_RangeError.hxx>
  46 #include <Standard_NotImplemented.hxx>
  47
  48 typedef Geom_OffsetCurve         OffsetCurve;
  49 typedef Handle(Geom_OffsetCurve) Handle(OffsetCurve);
  50 typedef Geom_Curve               Curve;
  51 typedef Handle(Geom_Curve)       Handle(Curve);
  52 typedef Handle(Geom_Geometry)    Handle(Geometry);
  53
  54 typedef gp_Dir  Dir;
  55 typedef gp_Pnt  Pnt;
  56 typedef gp_Trsf Trsf;
  57 typedef gp_Vec  Vec;
  58 typedef gp_XYZ  XYZ;
  59
  60
  61
  62
  63
  64 static const int MaxDegree = 9;
  65    //ordre de derivation maximum pour la recherche de la premiere
  66    //derivee non nulle
  67
  68
  69
  70 //=======================================================================
  71 //function : Copy
  72 //purpose  :
  73 //=======================================================================
  74
  75 Handle(Geom_Geometry) Geom_OffsetCurve::Copy () const {
  76
  77  Handle(OffsetCurve) C;
  78  C = new OffsetCurve (basisCurve, offsetValue, direction);
  79  return C;
  80 }
  81
  82
  83
  84 //=======================================================================
  85 //function : Geom_OffsetCurve
  86 //purpose  :
  87 //=======================================================================
  88
  89 Geom_OffsetCurve::Geom_OffsetCurve (const Handle(Curve)& C,
  90                                     const Standard_Real           Offset,
  91                                     const Dir&           V      )
  92  : direction(V), offsetValue(Offset)
  93 {
  94   if (C->DynamicType() == STANDARD_TYPE(Geom_OffsetCurve)) {
  95     Handle(OffsetCurve) OC = Handle(OffsetCurve)::DownCast(C);
  96     SetBasisCurve (OC->BasisCurve());
  97
  98     Standard_Real PrevOff = OC->Offset();
  99     gp_Vec V1(OC->Direction());
 100     gp_Vec V2(direction);
 101     gp_Vec Vdir(PrevOff*V1 + offsetValue*V2);
 102
 103     if (Offset >= 0.) {
 104       offsetValue = Vdir.Magnitude();
 105       direction.SetXYZ(Vdir.XYZ());
 106     } else {
 107       offsetValue = -Vdir.Magnitude();
 108       direction.SetXYZ((-Vdir).XYZ());
 109     }
 110   }
 111   else {
 112     SetBasisCurve(C);
 113   }
 114 }
 115
 116
 117 //=======================================================================
 118 //function : Reverse
 119 //purpose  :
 120 //=======================================================================
 121
 122 void Geom_OffsetCurve::Reverse ()
 123 {
 124   basisCurve->Reverse();
 125   offsetValue = -offsetValue;
 126 }
 127
 128
 129 //=======================================================================
 130 //function : ReversedParameter
 131 //purpose  :
 132 //=======================================================================
 133
 134 Standard_Real Geom_OffsetCurve::ReversedParameter( const Standard_Real U) const
 135 {
 136   return basisCurve->ReversedParameter( U);
 137 }
 138
 139 //=======================================================================
 140 //function : Direction
 141 //purpose  :
 142 //=======================================================================
 143
 144 const gp_Dir& Geom_OffsetCurve::Direction () const
 145   { return direction; }
 146
 147 //=======================================================================
 148 //function : SetDirection
 149 //purpose  :
 150 //=======================================================================
 151
 152 void Geom_OffsetCurve::SetDirection (const Dir& V)
 153   { direction = V; }
 154
 155 //=======================================================================
 156 //function : SetOffsetValue
 157 //purpose  :
 158 //=======================================================================
 159
 160 void Geom_OffsetCurve::SetOffsetValue (const Standard_Real D)
 161   { offsetValue = D; }
 162
 163
 164 //=======================================================================
 165 //function : IsPeriodic
 166 //purpose  :
 167 //=======================================================================
 168
 169 Standard_Boolean Geom_OffsetCurve::IsPeriodic () const
 170 {
 171   return basisCurve->IsPeriodic();
 172 }
 173
 174 //=======================================================================
 175 //function : Period
 176 //purpose  :
 177 //=======================================================================
 178
 179 Standard_Real Geom_OffsetCurve::Period () const
 180 {
 181   return basisCurve->Period();
 182 }
 183
 184 //=======================================================================
 185 //function : SetBasisCurve
 186 //purpose  :
 187 //=======================================================================
 188
 189 void Geom_OffsetCurve::SetBasisCurve (const Handle(Curve)& C)
 190 {
 191   Handle(Curve) aBasisCurve = Handle(Curve)::DownCast(C->Copy());
 192
 193   // Basis curve must be at least C1
 194   if (aBasisCurve->Continuity() == GeomAbs_C0)
 195   {
 196     // For B-splines it is sometimes possible to increase continuity by removing
 197     // unnecessarily duplicated knots
 198     if (aBasisCurve->IsKind(STANDARD_TYPE(Geom_BSplineCurve)))
 199     {
 200       Handle(Geom_BSplineCurve) aBCurve = Handle(Geom_BSplineCurve)::DownCast(aBasisCurve);
 201       Standard_Integer degree = aBCurve->Degree();
 202       Standard_Real Toler = Precision::Confusion();
 203       Standard_Integer start = aBCurve->IsPeriodic() ? 1 :  aBCurve->FirstUKnotIndex(),
 204                        finish = aBCurve->IsPeriodic() ? aBCurve->NbKnots() :  aBCurve->LastUKnotIndex();
 205       for (Standard_Integer i = start; i <= finish; i++)
 206       {
 207         Standard_Integer mult = aBCurve->Multiplicity(i);
 208         if ( mult == degree )
 209           aBCurve->RemoveKnot(i,degree - 1, Toler);
 210       }
 211     }
 212
 213     // Raise exception if still C0
 214     if (aBasisCurve->Continuity() == GeomAbs_C0)
 215       Standard_ConstructionError::Raise("Offset on C0 curve");
 216   }
 217
 218   basisCurve = aBasisCurve;
 219 }
 220
 221
 222
 223 //=======================================================================
 224 //function : BasisCurve
 225 //purpose  :
 226 //=======================================================================
 227
 228 Handle(Curve) Geom_OffsetCurve::BasisCurve () const
 229 {
 230   return basisCurve;
 231 }
 232
 233
 234 //=======================================================================
 235 //function : Continuity
 236 //purpose  :
 237 //=======================================================================
 238
 239 GeomAbs_Shape Geom_OffsetCurve::Continuity () const {
 240
 241   GeomAbs_Shape OffsetShape=GeomAbs_C0;
 242   switch (basisCurve->Continuity()) {
 243     case GeomAbs_C0 : OffsetShape = GeomAbs_C0;       break;
 244     case GeomAbs_C1 : OffsetShape = GeomAbs_C0;       break;
 245     case GeomAbs_C2 : OffsetShape = GeomAbs_C1;       break;
 246     case GeomAbs_C3 : OffsetShape = GeomAbs_C2;       break;
 247     case GeomAbs_CN : OffsetShape = GeomAbs_CN;       break;
 248     case GeomAbs_G1 : OffsetShape = GeomAbs_G1;       break;
 249     case GeomAbs_G2 : OffsetShape = GeomAbs_G2;       break;
 250   }
 251   return OffsetShape;
 252 }
 253
 254
 255 //=======================================================================
 256 //function : D0
 257 //purpose  :
 258 //=======================================================================
 259
 260 void Geom_OffsetCurve::D0 (const Standard_Real U, Pnt& P) const
 261 {
 262   gp_Pnt PBasis;
 263   gp_Vec VBasis;
 264   D0(U,P,PBasis,VBasis);
 265 }
 266
 267
 268 //=======================================================================
 269 //function : D1
 270 //purpose  :
 271 //=======================================================================
 272
 273 void Geom_OffsetCurve::D1 (const Standard_Real U, Pnt& P, Vec& V1) const
 274 {
 275   gp_Pnt PBasis;
 276   gp_Vec V1Basis,V2Basis;
 277   D1(U,P,PBasis,V1,V1Basis,V2Basis);
 278 }
 279
 280
 281
 282 //=======================================================================
 283 //function : D2
 284 //purpose  :
 285 //=======================================================================
 286
 287 void Geom_OffsetCurve::D2 (const Standard_Real U, Pnt& P, Vec& V1, Vec& V2) const
 288 {
 289   gp_Pnt PBasis;
 290   gp_Vec V1Basis,V2Basis,V3Basis;
 291   D2(U,P,PBasis,V1,V2,V1Basis,V2Basis,V3Basis);
 292 }
 293
 294
 295 //=======================================================================
 296 //function : D3
 297 //purpose  :
 298 //=======================================================================
 299
 300 void Geom_OffsetCurve::D3 (const Standard_Real U, Pnt& P, Vec& V1, Vec& V2, Vec& V3)
 301 const {
 302
 303
 304    // P(u) = p(u) + Offset * Ndir / R
 305    // with R = || p' ^ V|| and Ndir = P' ^ direction (local normal direction)
 306
 307    // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R -  Ndir * (DR/R))
 308
 309    // P"(u) = p"(u) + (Offset / R) * (D2Ndir/DU - DNdir * (2.0 * Dr/ R**2) +
 310    //         Ndir * ( (3.0 * Dr**2 / R**4) - (D2r / R**2)))
 311
 312    //P"'(u) = p"'(u) + (Offset / R) * (D3Ndir - (3.0 * Dr/R**2) * D2Ndir -
 313    //         (3.0 * D2r / R2) * DNdir + (3.0 * Dr * Dr / R4) * DNdir -
 314    //         (D3r/R2) * Ndir + (6.0 * Dr * Dr / R4) * Ndir +
 315    //         (6.0 * Dr * D2r / R4) * Ndir - (15.0 * Dr* Dr* Dr /R6) * Ndir
 316
 317
 318
 319   basisCurve->D3 (U, P, V1, V2, V3);
 320   Vec V4 = basisCurve->DN (U, 4);
 321   Standard_Integer Index = 2;
 322   while (V1.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 323     V1 = basisCurve->DN (U, Index);
 324     Index++;
 325   }
 326   if (Index != 2) {
 327     V2 = basisCurve->DN (U, Index);
 328     V3 = basisCurve->DN (U, Index + 1);
 329     V4 = basisCurve->DN (U, Index + 2);
 330   }
 331   XYZ OffsetDir = direction.XYZ();
 332   XYZ Ndir      = (V1.XYZ()).Crossed (OffsetDir);
 333   XYZ DNdir     = (V2.XYZ()).Crossed (OffsetDir);
 334   XYZ D2Ndir    = (V3.XYZ()).Crossed (OffsetDir);
 335   XYZ D3Ndir    = (V4.XYZ()).Crossed (OffsetDir);
 336   Standard_Real R2  = Ndir.SquareModulus();
 337   Standard_Real R   = Sqrt (R2);
 338   Standard_Real R3  = R2 * R;
 339   Standard_Real R4  = R2 * R2;
 340   Standard_Real R5  = R3 * R2;
 341   Standard_Real R6  = R3 * R3;
 342   Standard_Real R7  = R5 * R2;
 343   Standard_Real Dr  = Ndir.Dot (DNdir);
 344   Standard_Real D2r = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir);
 345   Standard_Real D3r = Ndir.Dot (D3Ndir) + 3.0 * DNdir.Dot (D2Ndir);
 346   if (R7 <= gp::Resolution()) {
 347     if (R6 <= gp::Resolution())  Geom_UndefinedDerivative::Raise();
 348     // V3 = P"' (U) :
 349     D3Ndir.Subtract (D2Ndir.Multiplied (3.0 * Dr / R2));
 350     D3Ndir.Subtract (DNdir.Multiplied (3.0 * ((D2r/R2) + (Dr*Dr/R4))));
 351     D3Ndir.Add (Ndir.Multiplied (6.0*Dr*Dr/R4 + 6.0*Dr*D2r/R4 -
 352                                  15.0*Dr*Dr*Dr/R6 - D3r));
 353     D3Ndir.Multiply (offsetValue/R);
 354     V3.Add (Vec(D3Ndir));
 355     // V2 = P" (U) :
 356     Standard_Real R4 = R2 * R2;
 357     D2Ndir.Subtract (DNdir.Multiplied (2.0 * Dr / R2));
 358     D2Ndir.Subtract (Ndir.Multiplied ((3.0 * Dr * Dr / R4) - (D2r / R2)));
 359     D2Ndir.Multiply (offsetValue / R);
 360     V2.Add (Vec(D2Ndir));
 361     // V1 = P' (U) :
 362     DNdir.Multiply(R);
 363     DNdir.Subtract (Ndir.Multiplied (Dr/R));
 364     DNdir.Multiply (offsetValue/R2);
 365     V1.Add (Vec(DNdir));
 366   }
 367   else {
 368     // V3 = P"' (U) :
 369     D3Ndir.Divide (R);
 370     D3Ndir.Subtract (D2Ndir.Multiplied (3.0 * Dr / R3));
 371     D3Ndir.Subtract (DNdir.Multiplied ((3.0 * ((D2r/R3) + (Dr*Dr)/R5))));
 372     D3Ndir.Add (Ndir.Multiplied (6.0*Dr*Dr/R5 + 6.0*Dr*D2r/R5 -
 373                 15.0*Dr*Dr*Dr/R7 - D3r));
 374     D3Ndir.Multiply (offsetValue);
 375     V3.Add (Vec(D3Ndir));
 376     // V2 = P" (U) :
 377     D2Ndir.Divide (R);
 378     D2Ndir.Subtract (DNdir.Multiplied (2.0 * Dr / R3));
 379     D2Ndir.Subtract (Ndir.Multiplied ((3.0 * Dr * Dr / R5) - (D2r / R3)));
 380     D2Ndir.Multiply (offsetValue);
 381     V2.Add (Vec(D2Ndir));
 382     // V1 = P' (U) :
 383     DNdir.Multiply (offsetValue/R);
 384     DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3));
 385     V1.Add (Vec(DNdir));
 386   }
 387   //P (U) :
 388   Ndir.Multiply (offsetValue/R);
 389   Ndir.Add (P.XYZ());
 390   P.SetXYZ (Ndir);
 391 }
 392
 393
 394 //=======================================================================
 395 //function : DN
 396 //purpose  :
 397 //=======================================================================
 398
 399 Vec Geom_OffsetCurve::DN (const Standard_Real U, const Standard_Integer N) const {
 400
 401
 402   if (N < 1)                  Standard_RangeError::Raise();
 403   XYZ OffsetDir = direction.XYZ();
 404   Pnt P;
 405   Vec V1, V2, dummy;
 406   if (N == 1) {
 407     basisCurve->D2 (U, P, V1, V2);
 408     Standard_Integer Index = 2;
 409     while (V1.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 410       V1 = basisCurve->DN (U, Index);
 411       Index++;
 412     }
 413     if (Index != 2)  V2 = basisCurve->DN (U, Index);
 414     XYZ Ndir  = (V1.XYZ()).Crossed (OffsetDir);
 415     XYZ DNdir = (V2.XYZ()).Crossed (OffsetDir);
 416     Standard_Real R2 = Ndir.SquareModulus();
 417     Standard_Real R  = Sqrt (R2);
 418     Standard_Real R3 = R * R2;
 419     Standard_Real Dr = Ndir.Dot (DNdir);
 420     if (R3 <= gp::Resolution()) {
 421       if (R2 <= gp::Resolution())  Geom_UndefinedDerivative::Raise();
 422       Ndir.Multiply (Dr/R);
 423       DNdir.Multiply(R);
 424       DNdir.Subtract (Ndir);
 425       DNdir.Multiply (offsetValue/R2);
 426       V1.Add (Vec(DNdir));
 427     }
 428     else {
 429       Ndir.Multiply (offsetValue * Dr / R3);
 430       DNdir.Multiply (offsetValue/R);
 431       DNdir.Subtract (Ndir);
 432       V1.Add (Vec(DNdir));
 433     }
 434     dummy = V1;
 435   }
 436   else if (N == 2) {
 437     Vec V3;
 438     basisCurve->D3 (U, P, V1, V2, V3);
 439     Standard_Integer Index = 2;
 440     while (V1.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 441       V1 = basisCurve->DN (U, Index);
 442       Index++;
 443     }
 444     if (Index != 2) {
 445       V2 = basisCurve->DN (U, Index);
 446       V3 = basisCurve->DN (U, Index + 1);
 447     }
 448     XYZ Ndir   = (V1.XYZ()).Crossed (OffsetDir);
 449     XYZ DNdir  = (V2.XYZ()).Crossed (OffsetDir);
 450     XYZ D2Ndir = (V3.XYZ()).Crossed (OffsetDir);
 451     Standard_Real R2  = Ndir.SquareModulus();
 452     Standard_Real R   = Sqrt (R2);
 453     Standard_Real R3  = R2 * R;   Standard_Real R4 = R2 * R2;  Standard_Real R5 = R3 * R2;
 454     Standard_Real Dr  = Ndir.Dot (DNdir);
 455     Standard_Real D2r = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir);
 456     if (R5 <= gp::Resolution()) {
 457       if (R4 <= gp::Resolution())  Geom_UndefinedDerivative::Raise();
 458       Ndir.Multiply ((3.0 * Dr * Dr / R4) - (D2r/R2));
 459       DNdir.Multiply (2.0 * Dr / R2);
 460       D2Ndir.Subtract (DNdir);
 461       D2Ndir.Subtract (Ndir);
 462       D2Ndir.Multiply (offsetValue / R);
 463       V2.Add (Vec(D2Ndir));
 464     }
 465     else {
 466       Ndir.Multiply ((3.0 * Dr * Dr / R4) - (D2r / R2));
 467       DNdir.Multiply (2.0 * Dr / R2);
 468       D2Ndir.Divide (R);
 469       D2Ndir.Subtract (DNdir);
 470       D2Ndir.Subtract (Ndir);
 471       D2Ndir.Multiply (offsetValue);
 472       V2.Add (Vec(D2Ndir));
 473     }
 474     dummy = V2;
 475   }
 476   else if (N == 3) {
 477     Vec V3;
 478     basisCurve->D3 (U, P, V1, V2, V3);
 479     Vec V4 = basisCurve->DN (U, 4);
 480     Standard_Integer Index = 2;
 481     while (V1.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 482       V1 = basisCurve->DN (U, Index);
 483       Index++;
 484     }
 485     if (Index != 2) {
 486       V2 = basisCurve->DN (U, Index);
 487       V3 = basisCurve->DN (U, Index + 1);
 488       V4 = basisCurve->DN (U, Index + 2);
 489     }
 490     XYZ Ndir = (V1.XYZ()).Crossed (OffsetDir);
 491     XYZ DNdir = (V2.XYZ()).Crossed (OffsetDir);
 492     XYZ D2Ndir = (V3.XYZ()).Crossed (OffsetDir);
 493     XYZ D3Ndir = (V4.XYZ()).Crossed (OffsetDir);
 494     Standard_Real R2  = Ndir.SquareModulus();
 495     Standard_Real R   = Sqrt (R2);  Standard_Real R3 = R2 * R;  Standard_Real R4 = R2 * R2;
 496     Standard_Real R5  = R3 * R2;  Standard_Real R6 = R3 * R3;  Standard_Real R7 = R5 * R2;
 497     Standard_Real Dr  = Ndir.Dot (DNdir);
 498     Standard_Real D2r = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir);
 499     Standard_Real D3r = Ndir.Dot (D3Ndir) + 3.0 * DNdir.Dot (D2Ndir);
 500     if (R7 <= gp::Resolution()) {
 501       if (R6 <= gp::Resolution()) Geom_UndefinedDerivative::Raise();
 502       D2Ndir.Multiply (3.0 * Dr / R2);
 503       DNdir.Multiply (3.0 * ((D2r/R2) + (Dr*Dr)/R4));
 504       Ndir.Multiply (6.0*Dr*Dr/R4 + 6.0*Dr*D2r/R4 - 15.0*Dr*Dr*Dr/R6 - D3r);
 505       D3Ndir.Subtract (D2Ndir);
 506       D3Ndir.Subtract (DNdir);
 507       D3Ndir.Add (Ndir);
 508       D3Ndir.Multiply (offsetValue/R);
 509       V3.Add (Vec(D3Ndir));
 510     }
 511     else {
 512       D2Ndir.Multiply (3.0 * Dr / R3);
 513       DNdir.Multiplied (3.0 * ((D2r/R3) + (Dr*Dr/R5)));
 514       Ndir.Multiply (6.0*Dr*Dr/R5 + 6.0*Dr*D2r/R5 - 15.0*Dr*Dr*Dr/R7 - D3r);
 515       D3Ndir.Divide (R);
 516       D3Ndir.Subtract (D2Ndir);
 517       D3Ndir.Subtract (DNdir);
 518       D3Ndir.Add (Ndir);
 519       D3Ndir.Multiply (offsetValue);
 520       V3.Add (Vec(D3Ndir));
 521     }
 522     dummy = V3;
 523   }
 524   else { Standard_NotImplemented::Raise();  }
 525
 526   return dummy;
 527 }
 528
 529 //=======================================================================
 530 //function : D0
 531 //purpose  :
 532 //=======================================================================
 533
 534 void  Geom_OffsetCurve::D0(const Standard_Real U,
 535                            gp_Pnt& P,
 536                            gp_Pnt& Pbasis,
 537                            gp_Vec& V1basis)const
 538 {
 539
 540   basisCurve->D1 (U, Pbasis, V1basis);
 541   Standard_Integer Index = 2;
 542   while (V1basis.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 543     V1basis = basisCurve->DN (U, Index);
 544     Index++;
 545   }
 546   XYZ Ndir = (V1basis.XYZ()).Crossed (direction.XYZ());
 547   Standard_Real R = Ndir.Modulus();
 548   if (R <= gp::Resolution())  Geom_UndefinedValue::Raise();
 549   Ndir.Multiply (offsetValue/R);
 550   Ndir.Add (Pbasis.XYZ());
 551   P.SetXYZ(Ndir);
 552 }
 553
 554 //=======================================================================
 555 //function : D1
 556 //purpose  :
 557 //=======================================================================
 558
 559 void Geom_OffsetCurve::D1 ( const Standard_Real U,
 560                             Pnt& P , Pnt& PBasis ,
 561                             Vec& V1, Vec& V1basis, Vec& V2basis) const {
 562
 563    // P(u) = p(u) + Offset * Ndir / R
 564    // with R = || p' ^ V|| and Ndir = P' ^ direction (local normal direction)
 565
 566    // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R -  Ndir * (DR/R))
 567
 568   basisCurve->D2 (U, PBasis, V1basis, V2basis);
 569   V1 = V1basis;
 570   Vec V2 = V2basis;
 571   Standard_Integer Index = 2;
 572   while (V1.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 573     V1 = basisCurve->DN (U, Index);
 574     Index++;
 575   }
 576   if (Index != 2)  V2 = basisCurve->DN (U, Index);
 577   XYZ OffsetDir = direction.XYZ();
 578   XYZ Ndir  = (V1.XYZ()).Crossed (OffsetDir);
 579   XYZ DNdir = (V2.XYZ()).Crossed (OffsetDir);
 580   Standard_Real R2 = Ndir.SquareModulus();
 581   Standard_Real R  = Sqrt (R2);
 582   Standard_Real R3 = R * R2;
 583   Standard_Real Dr = Ndir.Dot (DNdir);
 584   if (R3 <= gp::Resolution()) {
 585     //We try another computation but the stability is not very good.
 586     if (R2 <= gp::Resolution()) Geom_UndefinedDerivative::Raise();
 587     DNdir.Multiply(R);
 588     DNdir.Subtract (Ndir.Multiplied (Dr/R));
 589     DNdir.Multiply (offsetValue/R2);
 590     V1.Add (Vec(DNdir));
 591   }
 592   else {
 593     // Same computation as IICURV in EUCLID-IS because the stability is
 594     // better
 595     DNdir.Multiply (offsetValue/R);
 596     DNdir.Subtract (Ndir.Multiplied (offsetValue * Dr/R3));
 597     V1.Add (Vec(DNdir));
 598   }
 599   Ndir.Multiply (offsetValue/R);
 600   Ndir.Add (PBasis.XYZ());
 601   P.SetXYZ (Ndir);
 602 }
 603
 604
 605 //=======================================================================
 606 //function : D2
 607 //purpose  :
 608 //=======================================================================
 609
 610 void Geom_OffsetCurve::D2 (const Standard_Real U,
 611                            Pnt& P      , Pnt& PBasis ,
 612                            Vec& V1     , Vec& V2     ,
 613                            Vec& V1basis, Vec& V2basis, Vec& V3basis) const {
 614
 615    // P(u) = p(u) + Offset * Ndir / R
 616    // with R = || p' ^ V|| and Ndir = P' ^ direction (local normal direction)
 617
 618    // P'(u) = p'(u) + (Offset / R**2) * (DNdir/DU * R -  Ndir * (DR/R))
 619
 620    // P"(u) = p"(u) + (Offset / R) * (D2Ndir/DU - DNdir * (2.0 * Dr/ R**2) +
 621    //         Ndir * ( (3.0 * Dr**2 / R**4) - (D2r / R**2)))
 622
 623   basisCurve->D3 (U, PBasis, V1basis, V2basis, V3basis);
 624   Standard_Integer Index = 2;
 625   V1     = V1basis;
 626   V2     = V2basis;
 627   Vec V3 = V3basis;
 628   while (V1.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 629     V1 = basisCurve->DN (U, Index);
 630     Index++;
 631   }
 632   if (Index != 2) {
 633     V2 = basisCurve->DN (U, Index);
 634     V3 = basisCurve->DN (U, Index + 1);
 635   }
 636   XYZ OffsetDir = direction.XYZ();
 637   XYZ Ndir   = (V1.XYZ()).Crossed (OffsetDir);
 638   XYZ DNdir  = (V2.XYZ()).Crossed (OffsetDir);
 639   XYZ D2Ndir = (V3.XYZ()).Crossed (OffsetDir);
 640   Standard_Real R2    = Ndir.SquareModulus();
 641   Standard_Real R     = Sqrt (R2);
 642   Standard_Real R3    = R2 * R;
 643   Standard_Real R4    = R2 * R2;
 644   Standard_Real R5    = R3 * R2;
 645   Standard_Real Dr    = Ndir.Dot (DNdir);
 646   Standard_Real D2r   = Ndir.Dot (D2Ndir) + DNdir.Dot (DNdir);
 647   if (R5 <= gp::Resolution()) {
 648     //We try another computation but the stability is not very good
 649     //dixit ISG.
 650     if (R4 <= gp::Resolution())  Geom_UndefinedDerivative::Raise();
 651     // V2 = P" (U) :
 652     Standard_Real R4 = R2 * R2;
 653     D2Ndir.Subtract (DNdir.Multiplied (2.0 * Dr / R2));
 654     D2Ndir.Add (Ndir.Multiplied (((3.0 * Dr * Dr)/R4) - (D2r/R2)));
 655     D2Ndir.Multiply (offsetValue / R);
 656     V2.Add (Vec(D2Ndir));
 657     // V1 = P' (U) :
 658     DNdir.Multiply(R);
 659     DNdir.Subtract (Ndir.Multiplied (Dr/R));
 660     DNdir.Multiply (offsetValue/R2);
 661     V1.Add (Vec(DNdir));
 662   }
 663   else {
 664     // Same computation as IICURV in EUCLID-IS because the stability is
 665     // better.
 666     // V2 = P" (U) :
 667     D2Ndir.Multiply (offsetValue/R);
 668     D2Ndir.Subtract (DNdir.Multiplied (2.0 * offsetValue * Dr / R3));
 669     D2Ndir.Add (Ndir.Multiplied (
 670                      offsetValue * (((3.0 * Dr * Dr) / R5) - (D2r / R3))
 671                                      )
 672                     );
 673     V2.Add (Vec(D2Ndir));
 674     // V1 = P' (U) :
 675     DNdir.Multiply (offsetValue/R);
 676     DNdir.Subtract (Ndir.Multiplied (offsetValue*Dr/R3));
 677     V1.Add (Vec(DNdir));
 678   }
 679   //P (U) :
 680   Ndir.Multiply (offsetValue/R);
 681   Ndir.Add (PBasis.XYZ());
 682   P.SetXYZ (Ndir);
 683 }
 684
 685
 686 //=======================================================================
 687 //function : FirstParameter
 688 //purpose  :
 689 //=======================================================================
 690
 691 Standard_Real Geom_OffsetCurve::FirstParameter () const {
 692
 693    return basisCurve->FirstParameter();
 694 }
 695
 696
 697 //=======================================================================
 698 //function : LastParameter
 699 //purpose  :
 700 //=======================================================================
 701
 702 Standard_Real Geom_OffsetCurve::LastParameter () const {
 703
 704    return basisCurve->LastParameter();
 705 }
 706
 707
 708 //=======================================================================
 709 //function : Offset
 710 //purpose  :
 711 //=======================================================================
 712
 713 Standard_Real Geom_OffsetCurve::Offset () const { return offsetValue; }
 714
 715 //=======================================================================
 716 //function : Value
 717 //purpose  :
 718 //=======================================================================
 719
 720 void Geom_OffsetCurve::Value (
 721 const Standard_Real U, Pnt& P, Pnt& PBasis,  Vec& V1basis) const {
 722
 723   if (basisCurve->Continuity() == GeomAbs_C0)  Geom_UndefinedValue::Raise();
 724   basisCurve->D1 (U, PBasis, V1basis);
 725   Standard_Integer Index = 2;
 726   while (V1basis.Magnitude() <= gp::Resolution() && Index <= MaxDegree) {
 727     V1basis = basisCurve->DN (U, Index);
 728     Index++;
 729   }
 730   XYZ Ndir = (V1basis.XYZ()).Crossed (direction.XYZ());
 731   Standard_Real R = Ndir.Modulus();
 732   if (R <= gp::Resolution())  Geom_UndefinedValue::Raise();
 733   Ndir.Multiply (offsetValue/R);
 734   Ndir.Add (PBasis.XYZ());
 735   P.SetXYZ (Ndir);
 736 }
 737
 738
 739 //=======================================================================
 740 //function : IsClosed
 741 //purpose  :
 742 //=======================================================================
 743
 744 Standard_Boolean Geom_OffsetCurve::IsClosed () const
 745 {
 746   gp_Pnt PF,PL;
 747   D0(FirstParameter(),PF);
 748   D0(LastParameter(),PL);
 749   return ( PF.Distance(PL) <= gp::Resolution());
 750 }
 751
 752
 753
 754 //=======================================================================
 755 //function : IsCN
 756 //purpose  :
 757 //=======================================================================
 758
 759 Standard_Boolean Geom_OffsetCurve::IsCN (const Standard_Integer N) const {
 760
 761    Standard_RangeError_Raise_if (N < 0, " ");
 762    return basisCurve->IsCN (N + 1);
 763 }
 764
 765
 766 //=======================================================================
 767 //function : Transform
 768 //purpose  :
 769 //=======================================================================
 770
 771 void Geom_OffsetCurve::Transform (const Trsf& T) {
 772
 773   basisCurve->Transform (T);
 774   direction.Transform(T);
 775   offsetValue *= T.ScaleFactor();
 776 }
 777
 778 //=======================================================================
 779 //function : TransformedParameter
 780 //purpose  :
 781 //=======================================================================
 782
 783 Standard_Real Geom_OffsetCurve::TransformedParameter(const Standard_Real U,
 784                                                      const gp_Trsf& T) const
 785 {
 786   return basisCurve->TransformedParameter(U,T);
 787 }
 788
 789 //=======================================================================
 790 //function : ParametricTransformation
 791 //purpose  :
 792 //=======================================================================
 793
 794 Standard_Real Geom_OffsetCurve::ParametricTransformation(const gp_Trsf& T)
 795 const
 796 {
 797   return basisCurve->ParametricTransformation(T);
 798 }