src/BRepBlend/BRepBlend_RstRstConstRad.cxx

   1 // Created on: 1997-02-10
   2 // Created by: Jacques GOUSSARD
   3 // Copyright (c) 1997-1999 Matra Datavision
   4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   5 //
   6 // This file is part of Open CASCADE Technology software library.
   7 //
   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
  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.
  13 //
  14 // Alternatively, this file may be used under the terms of Open CASCADE
  15 // commercial license or contractual agreement.
  16
  17
  18 #include <Adaptor2d_HCurve2d.hxx>
  19 #include <Adaptor3d_HCurve.hxx>
  20 #include <Adaptor3d_HSurface.hxx>
  21 #include <Blend_Point.hxx>
  22 #include <BlendFunc.hxx>
  23 #include <BRepBlend_RstRstConstRad.hxx>
  24 #include <ElCLib.hxx>
  25 #include <GeomFill.hxx>
  26 #include <gp.hxx>
  27 #include <gp_Circ.hxx>
  28 #include <gp_Pnt.hxx>
  29 #include <gp_Pnt2d.hxx>
  30 #include <gp_Vec.hxx>
  31 #include <gp_Vec2d.hxx>
  32 #include <math_Gauss.hxx>
  33 #include <math_Matrix.hxx>
  34 #include <math_SVD.hxx>
  35 #include <Precision.hxx>
  36 #include <Standard_DomainError.hxx>
  37 #include <Standard_NotImplemented.hxx>
  38
  39 #define Eps 1.e-15
  40
  41 static void t3dto2d(Standard_Real& a,
  42                     Standard_Real& b,
  43                     const gp_Vec& A,
  44                     const gp_Vec& B,
  45                     const gp_Vec& C)
  46 {
  47   Standard_Real AB = A.Dot(B);
  48   Standard_Real AC = A.Dot(C);
  49   Standard_Real BC = B.Dot(C);
  50   Standard_Real BB = B.Dot(B);
  51   Standard_Real CC = C.Dot(C);
  52   Standard_Real deno = (BB*CC-BC*BC);
  53   a = (AB*CC-AC*BC)/deno;
  54   b = (AC*BB-AB*BC)/deno;
  55 }
  56
  57 //=======================================================================
  58 //function : BRepBlend_RstRstConstRad
  59 //purpose  :
  60 //=======================================================================
  61
  62 BRepBlend_RstRstConstRad::BRepBlend_RstRstConstRad
  63 (const Handle(Adaptor3d_HSurface)& Surf1,
  64  const Handle(Adaptor2d_HCurve2d)& Rst1,
  65  const Handle(Adaptor3d_HSurface)& Surf2,
  66  const Handle(Adaptor2d_HCurve2d)& Rst2,
  67  const Handle(Adaptor3d_HCurve)&   CGuide):
  68  surf1(Surf1), surf2(Surf2), rst1(Rst1), rst2(Rst2),
  69  cons1(Rst1, Surf1), cons2(Rst2, Surf2),
  70  guide(CGuide), tguide(CGuide),
  71  istangent(Standard_True), maxang(RealFirst()), minang(RealLast()),
  72  distmin(RealLast()),
  73  mySShape(BlendFunc_Rational)
  74 {}
  75
  76 //=======================================================================
  77 //function : NbVariables
  78 //purpose  :
  79 //=======================================================================
  80
  81 Standard_Integer BRepBlend_RstRstConstRad::NbVariables() const
  82 {
  83   return 2;
  84 }
  85
  86 //=======================================================================
  87 //function : NbEquations
  88 //purpose  :
  89 //=======================================================================
  90
  91 Standard_Integer BRepBlend_RstRstConstRad::NbEquations() const
  92 {
  93   return 2;
  94 }
  95
  96 //=======================================================================
  97 //function : Value
  98 //purpose  :
  99 //=======================================================================
 100
 101 Standard_Boolean BRepBlend_RstRstConstRad::Value(const math_Vector& X,
 102                                                  math_Vector&       F)
 103 {
 104   ptrst1 = cons1.Value(X(1));
 105   ptrst2 = cons2.Value(X(2));
 106
 107   F(1)   = nplan.XYZ().Dot(ptrst1.XYZ()) + theD;
 108   F(2)   = nplan.XYZ().Dot(ptrst2.XYZ()) + theD;
 109
 110   return Standard_True;
 111 }
 112
 113 //=======================================================================
 114 //function : Derivatives
 115 //purpose  :
 116 //=======================================================================
 117
 118 Standard_Boolean BRepBlend_RstRstConstRad::Derivatives(const math_Vector& X,
 119                                                        math_Matrix&       D)
 120 {
 121   gp_Vec d11, d21;
 122
 123   cons1.D1(X(1), ptrst1, d11);
 124   cons2.D1(X(2), ptrst2, d21);
 125
 126   D(1,1) = nplan.Dot(d11);
 127   D(1,2) = 0.;
 128
 129   D(2,1) = 0.;
 130   D(2,2) = nplan.Dot(d21);
 131
 132   return Standard_True;
 133 }
 134
 135 //=======================================================================
 136 //function : Values
 137 //purpose  :
 138 //=======================================================================
 139
 140 Standard_Boolean BRepBlend_RstRstConstRad::Values(const math_Vector& X,
 141                                                   math_Vector&       F,
 142                                                   math_Matrix&       D)
 143 {
 144   Value(X, F);
 145   Derivatives(X, D);
 146
 147   return Standard_True;
 148 }
 149
 150 //=======================================================================
 151 //function : Set
 152 //purpose  :
 153 //=======================================================================
 154
 155 void BRepBlend_RstRstConstRad::Set(const Handle(Adaptor3d_HSurface)& SurfRef1,
 156                                    const Handle(Adaptor2d_HCurve2d)& RstRef1,
 157                                    const Handle(Adaptor3d_HSurface)& SurfRef2,
 158                                    const Handle(Adaptor2d_HCurve2d)& RstRef2)
 159 {
 160   surfref1 = SurfRef1;
 161   surfref2 = SurfRef2;
 162   rstref1  = RstRef1;
 163   rstref2  = RstRef2;
 164 }
 165
 166 //=======================================================================
 167 //function : Set
 168 //purpose  :
 169 //=======================================================================
 170
 171 void BRepBlend_RstRstConstRad::Set(const Standard_Real Param)
 172 {
 173   tguide->D2(Param, ptgui, d1gui, d2gui);
 174   normtg = d1gui.Magnitude();
 175   nplan  = d1gui.Normalized();
 176   theD   = - (nplan.XYZ().Dot(ptgui.XYZ()));
 177 }
 178
 179 //=======================================================================
 180 //function : Set
 181 //purpose  :
 182 //=======================================================================
 183
 184 void BRepBlend_RstRstConstRad::Set(const Standard_Real First,
 185                                    const Standard_Real Last)
 186 {
 187  tguide = guide->Trim(First, Last, 1.e-12);
 188 }
 189
 190 //=======================================================================
 191 //function : GetTolerance
 192 //purpose  :
 193 //=======================================================================
 194
 195 void BRepBlend_RstRstConstRad::GetTolerance(math_Vector&        Tolerance,
 196                                             const Standard_Real Tol) const
 197 {
 198   Tolerance(1) = cons1.Resolution(Tol);
 199   Tolerance(2) = cons2.Resolution(Tol);
 200 }
 201
 202 //=======================================================================
 203 //function : GetBounds
 204 //purpose  :
 205 //=======================================================================
 206
 207 void BRepBlend_RstRstConstRad::GetBounds(math_Vector& InfBound,
 208                                          math_Vector& SupBound) const
 209 {
 210   InfBound(1) = cons1.FirstParameter();
 211   InfBound(2) = cons2.FirstParameter();
 212   SupBound(1) = cons1.LastParameter();
 213   SupBound(2) = cons2.LastParameter();
 214
 215 }
 216
 217 //=======================================================================
 218 //function : IsSolution
 219 //purpose  :
 220 //=======================================================================
 221
 222 Standard_Boolean BRepBlend_RstRstConstRad::IsSolution(const math_Vector&  Sol,
 223                                                       const Standard_Real Tol)
 224
 225
 226 {
 227   math_Vector valsol(1, 2), secmember(1, 2);
 228   math_Matrix gradsol(1, 2, 1, 2);
 229
 230   gp_Vec dnplan, d1urst1, d1vrst1, d1urst2, d1vrst2, d11, d21, temp;
 231   gp_Pnt bid;
 232
 233   Standard_Real Cosa, Sina, Angle;
 234
 235   Values(Sol, valsol, gradsol);
 236
 237   if (Abs(valsol(1)) <= Tol &&
 238       Abs(valsol(2)) <= Tol ) {
 239
 240     // Calculation of tangents
 241     prmrst1  = Sol(1);
 242     pt2drst1 = rst1->Value(prmrst1);
 243     prmrst2  = Sol(2);
 244     pt2drst2 = rst2->Value(prmrst2);
 245
 246     cons1.D1(Sol(1), ptrst1, d11);
 247     cons2.D1(Sol(2), ptrst2, d21);
 248
 249     dnplan.SetLinearForm(1./normtg, d2gui,
 250                          -1./normtg * (nplan.Dot(d2gui)), nplan);
 251
 252     temp.SetXYZ(ptrst1.XYZ() - ptgui.XYZ());
 253     secmember(1) = normtg - dnplan.Dot(temp);
 254
 255     temp.SetXYZ(ptrst2.XYZ() - ptgui.XYZ());
 256     secmember(2) = normtg - dnplan.Dot(temp);
 257
 258     math_Gauss Resol(gradsol);
 259
 260     if (Resol.IsDone()) {
 261       Resol.Solve(secmember);
 262       istangent = Standard_False;
 263     }
 264     else {
 265       math_SVD SingRS (gradsol);
 266       if (SingRS.IsDone()) {
 267         math_Vector DEDT(1,3);
 268         DEDT = secmember;
 269         SingRS.Solve(DEDT, secmember, 1.e-6);
 270         istangent = Standard_False;
 271       }
 272       else istangent = Standard_True;
 273     }
 274
 275
 276     if (!istangent) {
 277       tgrst1 = secmember(1) * d11;
 278       tgrst2 = secmember(2) * d21;
 279
 280       Standard_Real a, b;
 281       surf1->D1(pt2drst1.X(), pt2drst1.Y(), bid, d1urst1, d1vrst1);
 282       t3dto2d(a, b, tgrst1, d1urst1, d1vrst1);
 283       tg2drst1.SetCoord(a, b);
 284       surf2->D1(pt2drst2.X(), pt2drst2.Y(), bid, d1urst2, d1vrst2);
 285       t3dto2d(a, b, tgrst1, d1urst2, d1vrst2);
 286       tg2drst2.SetCoord(a, b);
 287     }
 288
 289     gp_Pnt Center;
 290     gp_Vec NotUsed;
 291     Standard_Boolean IsCenter;
 292
 293     IsCenter = CenterCircleRst1Rst2(ptrst1, ptrst2, nplan, Center, NotUsed);
 294
 295     if (!IsCenter) return Standard_False;
 296
 297     gp_Vec n1(Center, ptrst1) , n2(Center, ptrst2);
 298
 299     n1.Normalize();
 300     n2.Normalize();
 301
 302     Cosa = n1.Dot(n2);
 303     Sina = nplan.Dot(n1.Crossed(n2));
 304
 305     if (choix%2 != 0) {
 306       Sina = -Sina;  //nplan is changed into -nplan
 307     }
 308
 309     Angle = ACos(Cosa);
 310     if (Sina < 0.) {
 311       Angle = 2.*M_PI - Angle;
 312     }
 313
 314     if (Angle > maxang) {maxang = Angle;}
 315     if (Angle < minang) {minang = Angle;}
 316     distmin = Min( distmin, ptrst1.Distance(ptrst2));
 317
 318     return Standard_True;
 319   }
 320   istangent = Standard_True;
 321   return Standard_False;
 322 }
 323
 324 //=======================================================================
 325 //function : GetMinimalDistance
 326 //purpose  :
 327 //=======================================================================
 328
 329 Standard_Real BRepBlend_RstRstConstRad::GetMinimalDistance() const
 330 {
 331   return distmin;
 332 }
 333
 334 //=======================================================================
 335 //function : PointOnRst1
 336 //purpose  :
 337 //=======================================================================
 338
 339 const gp_Pnt& BRepBlend_RstRstConstRad::PointOnRst1() const
 340 {
 341   return ptrst1;
 342 }
 343
 344 //=======================================================================
 345 //function : PointOnRst2
 346 //purpose  :
 347 //=======================================================================
 348
 349 const gp_Pnt& BRepBlend_RstRstConstRad::PointOnRst2() const
 350 {
 351   return ptrst2;
 352 }
 353
 354 //=======================================================================
 355 //function : Pnt2dOnRst1
 356 //purpose  :
 357 //=======================================================================
 358
 359 const gp_Pnt2d& BRepBlend_RstRstConstRad::Pnt2dOnRst1() const
 360 {
 361   return pt2drst1;
 362 }
 363
 364 //=======================================================================
 365 //function : Pnt2dOnRst2
 366 //purpose  :
 367 //=======================================================================
 368
 369 const gp_Pnt2d& BRepBlend_RstRstConstRad::Pnt2dOnRst2() const
 370 {
 371   return pt2drst2;
 372 }
 373
 374 //=======================================================================
 375 //function : ParameterOnRst1
 376 //purpose  :
 377 //=======================================================================
 378
 379 Standard_Real BRepBlend_RstRstConstRad::ParameterOnRst1() const
 380 {
 381   return prmrst1;
 382 }
 383
 384 //=======================================================================
 385 //function : ParameterOnRst2
 386 //purpose  :
 387 //=======================================================================
 388
 389 Standard_Real BRepBlend_RstRstConstRad::ParameterOnRst2() const
 390 {
 391   return prmrst2;
 392 }
 393 //=======================================================================
 394 //function : IsTangencyPoint
 395 //purpose  :
 396 //=======================================================================
 397
 398 Standard_Boolean BRepBlend_RstRstConstRad::IsTangencyPoint() const
 399 {
 400   return istangent;
 401 }
 402
 403 //=======================================================================
 404 //function : TangentOnRst1
 405 //purpose  :
 406 //=======================================================================
 407
 408 const gp_Vec& BRepBlend_RstRstConstRad::TangentOnRst1() const
 409 {
 410   if (istangent) {throw Standard_DomainError();}
 411   return tgrst1;
 412 }
 413
 414 //=======================================================================
 415 //function : Tangent2dOnRst1
 416 //purpose  :
 417 //=======================================================================
 418
 419 const gp_Vec2d& BRepBlend_RstRstConstRad::Tangent2dOnRst1() const
 420 {
 421   if (istangent) {throw Standard_DomainError();}
 422   return tg2drst1;
 423 }
 424
 425 //=======================================================================
 426 //function : TangentOnRst2
 427 //purpose  :
 428 //=======================================================================
 429
 430 const gp_Vec& BRepBlend_RstRstConstRad::TangentOnRst2() const
 431 {
 432   if (istangent) {throw Standard_DomainError();}
 433   return tgrst2;
 434 }
 435
 436 //=======================================================================
 437 //function : Tangent2dOnRst2
 438 //purpose  :
 439 //=======================================================================
 440
 441 const gp_Vec2d& BRepBlend_RstRstConstRad::Tangent2dOnRst2() const
 442 {
 443   if (istangent) {throw Standard_DomainError();}
 444   return tg2drst2;
 445 }
 446
 447 //=======================================================================
 448 //function : Decroch
 449 //purpose  :
 450 //=======================================================================
 451
 452 Blend_DecrochStatus BRepBlend_RstRstConstRad::Decroch(const math_Vector& Sol,
 453                                                       gp_Vec&            NRst1,
 454                                                       gp_Vec&            TgRst1,
 455                                                       gp_Vec&            NRst2,
 456                                                       gp_Vec&            TgRst2)const
 457 {
 458   gp_Vec NRst1InPlane, NRst2InPlane;
 459   gp_Pnt PtTmp1, PtTmp2, Center;
 460   gp_Vec d1u, d1v, centptrst, NotUsed;
 461   Standard_Real norm, unsurnorm;
 462   Standard_Real u,v;
 463
 464   rstref1->Value(Sol(1)).Coord(u, v);
 465   surfref1->D1(u, v,PtTmp1,d1u,d1v);
 466   // Normal to the reference surface 1
 467   NRst1     = d1u.Crossed(d1v);
 468   rstref2->Value(Sol(2)).Coord(u, v);
 469   surfref2->D1(u, v, PtTmp2, d1u, d1v);
 470   // Normal to the reference surface 2
 471   NRst2     = d1u.Crossed(d1v);
 472
 473   CenterCircleRst1Rst2(PtTmp1, PtTmp2, nplan, Center, NotUsed);
 474
 475   norm      = nplan.Crossed(NRst1).Magnitude();
 476   unsurnorm = 1. / norm;
 477
 478   NRst1InPlane.SetLinearForm(nplan.Dot(NRst1) * unsurnorm, nplan, -unsurnorm, NRst1);
 479
 480   centptrst.SetXYZ(PtTmp1.XYZ() - Center.XYZ());
 481
 482   if (centptrst.Dot(NRst1InPlane) < 0.) NRst1InPlane.Reverse();
 483
 484   TgRst1    = nplan.Crossed(centptrst);
 485
 486   norm      = nplan.Crossed(NRst2).Magnitude();
 487   unsurnorm = 1./ norm;
 488   NRst2InPlane.SetLinearForm(nplan.Dot(NRst2) * unsurnorm, nplan, -unsurnorm, NRst2);
 489   centptrst.SetXYZ(PtTmp2.XYZ() - Center.XYZ());
 490
 491
 492   if (centptrst.Dot(NRst2InPlane) < 0.) NRst2InPlane.Reverse();
 493
 494   TgRst2 = nplan.Crossed(centptrst);
 495
 496   if (choix %2 != 0) {
 497     TgRst1.Reverse();
 498     TgRst2.Reverse();
 499   }
 500
 501   // The vectors are returned
 502   if (NRst1InPlane.Dot(TgRst1) > -1.e-10) {
 503     if (NRst2InPlane.Dot(TgRst2) < 1.e-10) {
 504       return Blend_DecrochBoth;
 505     }
 506     else {
 507       return Blend_DecrochRst1;
 508     }
 509   }
 510   else {
 511     if (NRst2InPlane.Dot(TgRst2) < 1.e-10) {
 512       return Blend_DecrochRst2;
 513     }
 514     else {
 515       return Blend_NoDecroch;
 516     }
 517   }
 518
 519 }
 520
 521 //=======================================================================
 522 //function : Set
 523 //purpose  :
 524 //=======================================================================
 525
 526 void BRepBlend_RstRstConstRad::Set(const Standard_Real    Radius,
 527                                    const Standard_Integer Choix)
 528 {
 529   choix = Choix;
 530   ray   = Abs(Radius);
 531
 532 }
 533
 534 //=======================================================================
 535 //function : Set
 536 //purpose  :
 537 //=======================================================================
 538
 539 void BRepBlend_RstRstConstRad::Set(const BlendFunc_SectionShape TypeSection)
 540 {
 541   mySShape = TypeSection;
 542 }
 543
 544
 545
 546 //=======================================================================
 547 //function : CenterCircleRst1Rst2
 548 //purpose  : Calculate the center of the circle passing by two points of restrictions
 549 //=======================================================================
 550 Standard_Boolean  BRepBlend_RstRstConstRad::CenterCircleRst1Rst2(const gp_Pnt&       PtRst1,
 551                                                                  const gp_Pnt&       PtRst2,
 552                                                                  const gp_Vec&       np,
 553                                                                  gp_Pnt&             Center,
 554                                                                  gp_Vec&             VdMed) const
 555 {
 556
 557   gp_Vec rst1rst2(PtRst1, PtRst2);
 558   gp_Vec   vdmedNor; //,NRst1;    vdmedNor directong vector of the perpendicular bisector
 559   Standard_Real norm2;
 560   Standard_Real Dist;// distance between the middle of PtRst1,PtRst2 and Center
 561
 562   // Calculate the center of the circle
 563   VdMed = rst1rst2.Crossed(np);
 564   norm2 = rst1rst2.SquareMagnitude();
 565   Dist  = ray * ray - 0.25 * norm2;
 566
 567   if (choix > 2) {
 568     VdMed.Reverse();
 569   }
 570
 571   if (Dist < - 1.E-07) return Standard_False;
 572
 573   if (Dist > 1.E-07) {
 574     Dist = sqrt(Dist);
 575     vdmedNor = VdMed.Normalized();
 576     Center.SetXYZ(0.5 * rst1rst2.XYZ() + PtRst1.XYZ() + Dist * vdmedNor.XYZ());
 577   }
 578   else
 579   {
 580     Center.SetXYZ(0.5 * rst1rst2.XYZ() + PtRst1.XYZ());
 581   }
 582
 583   return Standard_True;
 584
 585 }
 586
 587
 588
 589
 590
 591
 592 //=======================================================================
 593 //function : Section
 594 //purpose  :
 595 //=======================================================================
 596
 597 void BRepBlend_RstRstConstRad::Section(const Standard_Real Param,
 598                                        const Standard_Real U,
 599                                        const Standard_Real V,
 600                                        Standard_Real&      Pdeb,
 601                                        Standard_Real&      Pfin,
 602                                        gp_Circ&               C)
 603 {
 604   gp_Vec ns, np, NotUsed;
 605   gp_Pnt Center;
 606
 607   tguide->D1(Param, ptgui, d1gui);
 608   np       = d1gui.Normalized();
 609   ptrst1   = cons1.Value(U);
 610   ptrst2   = cons2.Value(V);
 611
 612   CenterCircleRst1Rst2(ptrst1, ptrst2, np, Center, NotUsed);
 613
 614   C.SetRadius(Abs(ray));
 615   ns = gp_Vec(Center, ptrst1).Normalized();
 616
 617   if (choix%2 != 0) {
 618     np.Reverse();
 619   }
 620
 621   C.SetPosition(gp_Ax2(Center, np, ns));
 622   Pdeb = 0; //ElCLib::Parameter(C, pts);
 623   Pfin = ElCLib::Parameter(C, ptrst2);
 624
 625   // Test of angles negative and almost null : Special Case
 626   if (Pfin > 1.5 * M_PI) {
 627     np.Reverse();
 628     C.SetPosition(gp_Ax2(Center, np, ns));
 629     Pfin = ElCLib::Parameter(C, ptrst2);
 630   }
 631   if (Pfin < Precision::PConfusion()) Pfin += Precision::PConfusion();
 632 }
 633
 634 //=======================================================================
 635 //function : IsRational
 636 //purpose  :
 637 //=======================================================================
 638
 639 Standard_Boolean BRepBlend_RstRstConstRad::IsRational () const
 640 {
 641   return  (mySShape==BlendFunc_Rational || mySShape==BlendFunc_QuasiAngular);
 642 }
 643
 644 //=======================================================================
 645 //function : GetSectionSize
 646 //purpose  :
 647 //=======================================================================
 648
 649 Standard_Real BRepBlend_RstRstConstRad::GetSectionSize() const
 650 {
 651   return maxang * Abs(ray);
 652 }
 653
 654 //=======================================================================
 655 //function : GetMinimalWeight
 656 //purpose  :
 657 //=======================================================================
 658
 659 void BRepBlend_RstRstConstRad::GetMinimalWeight(TColStd_Array1OfReal& Weights) const
 660 {
 661   BlendFunc::GetMinimalWeights(mySShape, myTConv, minang, maxang, Weights );
 662   // It is supposed that it does not depend on the Radius!
 663 }
 664
 665 //=======================================================================
 666 //function : NbIntervals
 667 //purpose  :
 668 //=======================================================================
 669
 670 Standard_Integer BRepBlend_RstRstConstRad::NbIntervals (const GeomAbs_Shape S) const
 671 {
 672   return guide->NbIntervals(BlendFunc::NextShape(S));
 673 }
 674
 675 //=======================================================================
 676 //function : Intervals
 677 //purpose  :
 678 //=======================================================================
 679
 680 void BRepBlend_RstRstConstRad::Intervals (TColStd_Array1OfReal& T,
 681                                           const GeomAbs_Shape S) const
 682 {
 683   guide->Intervals(T, BlendFunc::NextShape(S));
 684 }
 685
 686 //=======================================================================
 687 //function : GetShape
 688 //purpose  :
 689 //=======================================================================
 690
 691 void BRepBlend_RstRstConstRad::GetShape (Standard_Integer& NbPoles,
 692                                          Standard_Integer& NbKnots,
 693                                          Standard_Integer& Degree,
 694                                          Standard_Integer& NbPoles2d)
 695 {
 696   NbPoles2d = 2;
 697   BlendFunc::GetShape(mySShape, maxang, NbPoles, NbKnots, Degree, myTConv);
 698 }
 699
 700 //=======================================================================
 701 //function : GetTolerance
 702 //purpose  : Determine Tolerances to be used in approximations.
 703 //=======================================================================
 704
 705 void BRepBlend_RstRstConstRad::GetTolerance(const Standard_Real BoundTol,
 706                                             const Standard_Real SurfTol,
 707                                             const Standard_Real AngleTol,
 708                                             math_Vector& Tol3d,
 709                                             math_Vector& Tol1d) const
 710 {
 711   Standard_Integer low = Tol3d.Lower(), up = Tol3d.Upper();
 712   Standard_Real Tol;
 713   Tol= GeomFill::GetTolerance(myTConv, minang, Abs(ray),
 714                                AngleTol, SurfTol);
 715   Tol1d.Init(SurfTol);
 716   Tol3d.Init(SurfTol);
 717   Tol3d(low+1) = Tol3d(up-1) = Min(Tol, SurfTol);
 718   Tol3d(low)   = Tol3d(up)   = Min(Tol, BoundTol);
 719 }
 720
 721 //=======================================================================
 722 //function : Knots
 723 //purpose  :
 724 //=======================================================================
 725
 726 void BRepBlend_RstRstConstRad::Knots(TColStd_Array1OfReal& TKnots)
 727 {
 728   GeomFill::Knots(myTConv, TKnots);
 729 }
 730
 731 //=======================================================================
 732 //function : Mults
 733 //purpose  :
 734 //=======================================================================
 735
 736 void BRepBlend_RstRstConstRad::Mults(TColStd_Array1OfInteger& TMults)
 737 {
 738   GeomFill::Mults(myTConv, TMults);
 739 }
 740
 741 //=======================================================================
 742 //function : Section
 743 //purpose  :
 744 //=======================================================================
 745
 746 void BRepBlend_RstRstConstRad::Section(const Blend_Point& P,
 747                                        TColgp_Array1OfPnt& Poles,
 748                                        TColgp_Array1OfPnt2d& Poles2d,
 749                                        TColStd_Array1OfReal& Weights)
 750 {
 751   gp_Vec ns, ns2, NotUsed;
 752   gp_Pnt Center;
 753   Standard_Real u, v;
 754
 755   Standard_Real prm    = P.Parameter();
 756   Standard_Integer low = Poles.Lower();
 757   Standard_Integer upp = Poles.Upper();
 758
 759   tguide->D1(prm,ptgui, d1gui);
 760   nplan  = d1gui.Normalized();
 761
 762   u = P.ParameterOnC1();
 763   v = P.ParameterOnC2();
 764
 765   gp_Pnt2d  pt2d1 = rst1->Value(u);
 766   gp_Pnt2d  pt2d2 = rst2->Value(v);
 767
 768   ptrst1 = cons1.Value(u);
 769   ptrst2 = cons2.Value(v);
 770   distmin = Min (distmin, ptrst1.Distance(ptrst2));
 771
 772   Poles2d(Poles2d.Lower()).SetCoord(pt2d1.X(),pt2d1.Y());
 773   Poles2d(Poles2d.Upper()).SetCoord(pt2d2.X(),pt2d2.Y());
 774
 775   // Linear case
 776   if (mySShape == BlendFunc_Linear) {
 777     Poles(low)   = ptrst1;
 778     Poles(upp)   = ptrst2;
 779     Weights(low) = 1.0;
 780     Weights(upp) = 1.0;
 781     return;
 782   }
 783
 784   // Calculate the center of the circle
 785   CenterCircleRst1Rst2(ptrst1, ptrst2, nplan, Center, NotUsed);
 786
 787   // normals to the section with points
 788   ns  = gp_Vec(Center, ptrst1).Normalized();
 789   ns2 = gp_Vec(Center, ptrst2).Normalized();
 790
 791   if (choix%2 != 0) {
 792     nplan.Reverse();
 793   }
 794
 795   GeomFill::GetCircle(myTConv,
 796                       ns, ns2,
 797                       nplan, ptrst1, ptrst2,
 798                       Abs(ray), Center,
 799                       Poles, Weights);
 800 }
 801
 802 //=======================================================================
 803 //function : Section
 804 //purpose  :
 805 //=======================================================================
 806
 807 Standard_Boolean BRepBlend_RstRstConstRad::Section(const Blend_Point& P,
 808                                                    TColgp_Array1OfPnt& Poles,
 809                                                    TColgp_Array1OfVec& DPoles,
 810                                                    TColgp_Array1OfPnt2d& Poles2d,
 811                                                    TColgp_Array1OfVec2d& DPoles2d,
 812                                                    TColStd_Array1OfReal& Weights,
 813                                                    TColStd_Array1OfReal& DWeights)
 814 {
 815
 816   gp_Vec d11, d21;
 817   gp_Vec  dnplan, d1n1, d1n2;//,np2, dnp2;
 818   gp_Vec temp, tgct;
 819   gp_Vec d1urst, d1vrst;
 820   gp_Pnt Center, NotUsed;
 821
 822   Standard_Real norm2, normmed, Dist;
 823
 824   math_Vector sol(1, 2), valsol(1, 2), secmember(1, 2);
 825   math_Matrix gradsol(1, 2, 1, 2);
 826
 827   Standard_Real prm       = P.Parameter();
 828   Standard_Integer low    = Poles.Lower();
 829   Standard_Integer upp    = Poles.Upper();
 830   Standard_Boolean istgt;
 831
 832   tguide->D2(prm, ptgui, d1gui, d2gui);
 833   normtg = d1gui.Magnitude();
 834   nplan  = d1gui.Normalized();
 835   dnplan.SetLinearForm(1./normtg, d2gui,
 836                        -1./normtg * (nplan.Dot(d2gui)), nplan);
 837
 838   sol(1)   = prmrst1 = P.ParameterOnC1();
 839   sol(2)   = prmrst2 = P.ParameterOnC2();
 840   pt2drst1 = rst1->Value(prmrst1);
 841   pt2drst2 = rst2->Value(prmrst2);
 842
 843   Values(sol, valsol, gradsol);
 844
 845   cons1.D1(sol(1), ptrst1, d11);
 846   cons2.D1(sol(2), ptrst2, d21);
 847
 848   temp.SetXYZ(ptrst1.XYZ() - ptgui.XYZ());
 849   secmember(1) = normtg - dnplan.Dot(temp);
 850
 851   temp.SetXYZ(ptrst2.XYZ() - ptgui.XYZ());
 852   secmember(2) = normtg - dnplan.Dot(temp);
 853
 854   math_Gauss Resol(gradsol, 1.e-9);
 855
 856   if (Resol.IsDone()) {
 857     istgt = Standard_False;
 858     Resol.Solve(secmember);
 859   }
 860   else {
 861     math_SVD SingRS (gradsol);
 862     if (SingRS.IsDone()) {
 863       math_Vector DEDT(1,2);
 864       DEDT = secmember;
 865       SingRS.Solve(DEDT, secmember, 1.e-6);
 866       istgt = Standard_False;
 867     }
 868     else istgt = Standard_True;
 869   }
 870
 871   gp_Vec med;
 872   gp_Vec rst1rst2(ptrst1, ptrst2);
 873   Standard_Boolean IsCenter;
 874
 875   IsCenter = CenterCircleRst1Rst2(ptrst1, ptrst2, nplan, Center, med);
 876   if (!IsCenter) return Standard_False;
 877
 878   normmed = med.Magnitude();
 879   med.Normalize();
 880   gp_Vec n1(Center, ptrst1), n2(Center, ptrst2);
 881   n1.Normalize();
 882   n2.Normalize();
 883
 884   if (!istgt) {
 885     // secmember contains derivatives of parameters on curves
 886     // compared to t
 887     tgrst1 = secmember(1) * d11;
 888     tgrst2 = secmember(2) * d21;
 889
 890     gp_Vec d1rst1rst2;
 891
 892     norm2      = rst1rst2.SquareMagnitude();
 893     d1rst1rst2 = tgrst2 - tgrst1;
 894     Dist       = ray * ray - 0.25 * norm2;
 895
 896     if (Dist >  1.E-07) {
 897       gp_Vec d1P1P2CrosNp, dmed;
 898       d1P1P2CrosNp = d1rst1rst2.Crossed(nplan) + rst1rst2.Crossed(dnplan);
 899       // derivative of the perpendicular bisector
 900       dmed = d1P1P2CrosNp - med.Dot(d1P1P2CrosNp) * med;
 901       dmed /= normmed;
 902       Dist = sqrt(Dist);
 903       Standard_Real d1Dist = - (0.25 / Dist) * rst1rst2.Dot(d1rst1rst2);
 904
 905       if  (choix > 2) {
 906         dmed.Reverse();
 907       }
 908
 909       // the derivative of coefficient Dist is located in dmed
 910       dmed.SetLinearForm(Dist, dmed, d1Dist, med);
 911       d1rst1rst2 *= 0.5;
 912       // derivative of the Normal to the curve in P1
 913       d1n1 = - (dmed + d1rst1rst2) / ray;
 914
 915       // derivative of the Normal to the curve in P2
 916       d1n2 = (d1rst1rst2 - dmed) / ray;
 917     }
 918     else {
 919       d1rst1rst2 *= 0.5;
 920       // Normal to the curve in P1
 921       d1n1 = - d1rst1rst2 / ray;
 922
 923       // Normal to the curve in P2
 924       d1n2 = d1rst1rst2 / ray;
 925     }
 926   }
 927
 928   // Tops 2d
 929
 930   Poles2d(Poles2d.Lower()).SetCoord(pt2drst1.X(), pt2drst1.Y());
 931   Poles2d(Poles2d.Upper()).SetCoord(pt2drst2.X(), pt2drst2.Y());
 932   if (!istgt) {
 933     Standard_Real a, b;
 934     surf1->D1(pt2drst1.X(), pt2drst1.Y(), NotUsed, d1urst, d1vrst);
 935     t3dto2d(a,b,tgrst1, d1urst, d1vrst);
 936     DPoles2d(Poles2d.Lower()).SetCoord(a, b);
 937
 938     surf2->D1(pt2drst2.X(), pt2drst2.Y(), NotUsed, d1urst, d1vrst);
 939     t3dto2d(a, b, tgrst2, d1urst, d1vrst);
 940     DPoles2d(Poles2d.Upper()).SetCoord(a, b);
 941   }
 942
 943   // Linear case
 944   if (mySShape == BlendFunc_Linear) {
 945     Poles(low)   = ptrst1;
 946     Poles(upp)   = ptrst2;
 947     Weights(low) = 1.0;
 948     Weights(upp) = 1.0;
 949     if (!istgt) {
 950       DPoles(low)   = tgrst1;
 951       DPoles(upp)   = tgrst2;
 952       DWeights(low) = 0.0;
 953       DWeights(upp) = 0.0;
 954     }
 955     return (!istgt);
 956   }
 957
 958   // Case of the circle
 959   // tangent to the center of the circle
 960   if (!istgt) {
 961     tgct.SetLinearForm(-ray, d1n1, tgrst1);
 962   }
 963
 964
 965   if (choix%2 != 0) {
 966     nplan.Reverse();
 967     dnplan.Reverse();
 968   }
 969
 970   if (!istgt) {
 971     return GeomFill::GetCircle(myTConv,
 972                               n1, n2,
 973                               d1n1, d1n2,
 974                               nplan, dnplan,
 975                               ptrst1, ptrst2,
 976                               tgrst1, tgrst2,
 977                               Abs(ray), 0,
 978                               Center, tgct,
 979                               Poles,
 980                               DPoles,
 981                               Weights,
 982                               DWeights);
 983   }
 984   else {
 985     GeomFill::GetCircle(myTConv,
 986                        n1, n2,
 987                        nplan, ptrst1, ptrst2,
 988                        Abs(ray), Center,
 989                        Poles, Weights);
 990     return Standard_False;
 991   }
 992 }
 993
 994 //=======================================================================
 995 //function : Section
 996 //purpose  :
 997 //=======================================================================
 998
 999 Standard_Boolean BRepBlend_RstRstConstRad::Section
1000 (const Blend_Point&,
1001  TColgp_Array1OfPnt&,
1002  TColgp_Array1OfVec&,
1003  TColgp_Array1OfVec&,
1004  TColgp_Array1OfPnt2d&,
1005  TColgp_Array1OfVec2d&,
1006  TColgp_Array1OfVec2d&,
1007  TColStd_Array1OfReal&,
1008  TColStd_Array1OfReal&,
1009  TColStd_Array1OfReal&)
1010 {
1011   return Standard_False;
1012 }
1013
1014
1015 void BRepBlend_RstRstConstRad::Resolution(const Standard_Integer IC2d,
1016                                           const Standard_Real Tol,
1017                                           Standard_Real& TolU,
1018                                           Standard_Real& TolV) const
1019 {
1020   if(IC2d == 1){
1021     TolU = surf1->UResolution(Tol);
1022     TolV = surf1->VResolution(Tol);
1023   }
1024   else {
1025     TolU = surf2->UResolution(Tol);
1026     TolV = surf2->VResolution(Tol);
1027   }
1028 }
1029