src/IntPatch/IntPatch_ImpPrmIntersection.cxx

   1 // Created on: 1992-05-07
   2 // Created by: Jacques GOUSSARD
   3 // Copyright (c) 1992-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 #include <IntPatch_ImpPrmIntersection.hxx>
  18
  19 #include <Adaptor3d_Surface.hxx>
  20 #include <Adaptor3d_TopolTool.hxx>
  21 #include <ElCLib.hxx>
  22 #include <ElSLib.hxx>
  23 #include <IntPatch_ArcFunction.hxx>
  24 #include <IntPatch_PointLine.hxx>
  25 #include <IntPatch_RLine.hxx>
  26 #include <IntPatch_RstInt.hxx>
  27 #include <IntPatch_SpecialPoints.hxx>
  28 #include <IntPatch_TheIWLineOfTheIWalking.hxx>
  29 #include <IntPatch_TheIWalking.hxx>
  30 #include <IntPatch_TheSurfFunction.hxx>
  31 #include <IntPatch_WLine.hxx>
  32 #include <IntSurf.hxx>
  33 #include <IntSurf_Quadric.hxx>
  34 #include <IntSurf_QuadricTool.hxx>
  35 #include <IntSurf_SequenceOfPathPoint.hxx>
  36 #include <TColStd_Array1OfInteger.hxx>
  37 #include <TopAbs_Orientation.hxx>
  38 #include <TopTrans_CurveTransition.hxx>
  39 #include <math_Matrix.hxx>
  40 #include <math_Vector.hxx>
  41
  42 #ifndef OCCT_DEBUG
  43 #define No_Standard_RangeError
  44 #define No_Standard_OutOfRange
  45 #endif
  46
  47 static Standard_Boolean DecomposeResult(const Handle(IntPatch_PointLine)& theLine,
  48                                         const Standard_Boolean       IsReversed,
  49                                         const IntSurf_Quadric&       theQuad,
  50                                         const Handle(Adaptor3d_TopolTool)& thePDomain,
  51                                         const Handle(Adaptor3d_Surface)&  theQSurf,
  52                                         const Handle(Adaptor3d_Surface)&  theOtherSurf,
  53                                         const Standard_Real                theArcTol,
  54                                         const Standard_Real                theTolTang,
  55                                         IntPatch_SequenceOfLine&           theLines);
  56 static
  57   void ComputeTangency (const IntPatch_TheSOnBounds& solrst,
  58   IntSurf_SequenceOfPathPoint& seqpdep,
  59   const Handle(Adaptor3d_TopolTool)& Domain,
  60   IntPatch_TheSurfFunction& Func,
  61   const Handle(Adaptor3d_Surface)& PSurf,
  62   TColStd_Array1OfInteger& Destination);
  63 static
  64   void Recadre(const Standard_Boolean ,
  65   GeomAbs_SurfaceType typeS1,
  66   GeomAbs_SurfaceType typeS2,
  67   IntPatch_Point&  pt,
  68   const Handle(IntPatch_TheIWLineOfTheIWalking)& iwline,
  69   Standard_Integer Param,
  70   Standard_Real U1,
  71   Standard_Real V1,
  72   Standard_Real U2,
  73   Standard_Real V2);
  74
  75 static
  76   Standard_Boolean IsCoincide(IntPatch_TheSurfFunction& theFunc,
  77                               const Handle(IntPatch_PointLine)& theLine,
  78                               const Handle(Adaptor2d_Curve2d)& theArc,
  79                               const Standard_Boolean isTheSurface1Using,
  80                               const Standard_Real theToler3D,
  81                               const Standard_Real theToler2D,
  82                               const Standard_Real thePeriod);
  83
  84 //=======================================================================
  85 //function : IsSeamOrPole
  86 //purpose  :
  87 //=======================================================================
  88 static IntPatch_SpecPntType IsSeamOrPole(const Handle(Adaptor3d_Surface)& theQSurf,
  89                                          const Handle(IntSurf_LineOn2S)& theLine,
  90                                          const Standard_Boolean IsReversed,
  91                                          const Standard_Integer theRefIndex,
  92                                          const Standard_Real theTol3D,
  93                                          const Standard_Real theDeltaMax)
  94 {
  95   if((theRefIndex < 1) || (theRefIndex >= theLine->NbPoints()))
  96     return IntPatch_SPntNone;
  97
  98   //Parameters on Quadric and on parametric for reference point
  99   Standard_Real aUQRef, aVQRef, aUPRef, aVPRef;
 100   Standard_Real aUQNext, aVQNext, aUPNext, aVPNext;
 101
 102   const gp_Pnt &aP3d = theLine->Value(theRefIndex + 1).Value();
 103
 104   if(IsReversed)
 105   {
 106     theLine->Value(theRefIndex).Parameters  (aUPRef, aVPRef, aUQRef, aVQRef);
 107     theLine->Value(theRefIndex+1).Parameters(aUPNext, aVPNext, aUQNext, aVQNext);
 108   }
 109   else
 110   {
 111     theLine->Value(theRefIndex).Parameters  (aUQRef, aVQRef, aUPRef, aVPRef);
 112     theLine->Value(theRefIndex+1).Parameters(aUQNext, aVQNext, aUPNext, aVPNext);
 113   }
 114
 115   const GeomAbs_SurfaceType aType = theQSurf->GetType();
 116
 117   if ((aType == GeomAbs_Cone) &&
 118       (theQSurf->Cone().Apex().SquareDistance(aP3d) < theTol3D*theTol3D))
 119   {
 120     return IntPatch_SPntPoleSeamU;
 121   }
 122   else if (aType == GeomAbs_Sphere)
 123   {
 124     const Standard_Real aSqTol = theTol3D*theTol3D;
 125     gp_Pnt aP(ElSLib::Value(0.0, M_PI_2, theQSurf->Sphere()));
 126     if (aP.SquareDistance(aP3d) < aSqTol)
 127     {
 128       return IntPatch_SPntPoleSeamU;
 129     }
 130
 131     aP = ElSLib::Value(0.0, -M_PI_2, theQSurf->Sphere());
 132     if (aP.SquareDistance(aP3d) < aSqTol)
 133     {
 134       return IntPatch_SPntPoleSeamU;
 135     }
 136   }
 137
 138
 139   const Standard_Real aDeltaU = Abs(aUQRef - aUQNext);
 140
 141   if((aType != GeomAbs_Torus) && (aDeltaU < theDeltaMax))
 142     return IntPatch_SPntNone;
 143
 144   switch(aType)
 145   {
 146   case GeomAbs_Cylinder:
 147     return IntPatch_SPntSeamU;
 148
 149   case GeomAbs_Torus:
 150     {
 151       const Standard_Real aDeltaV = Abs(aVQRef - aVQNext);
 152
 153       if((aDeltaU >= theDeltaMax) && (aDeltaV >= theDeltaMax))
 154         return IntPatch_SPntSeamUV;
 155
 156       if(aDeltaU >= theDeltaMax)
 157         return IntPatch_SPntSeamU;
 158
 159       if(aDeltaV >= theDeltaMax)
 160         return IntPatch_SPntSeamV;
 161     }
 162
 163     break;
 164   case GeomAbs_Sphere:
 165   case GeomAbs_Cone:
 166     return IntPatch_SPntPoleSeamU;
 167   default:
 168     break;
 169   }
 170
 171   return IntPatch_SPntNone;
 172 }
 173
 174 //=======================================================================
 175 //function : IntPatch_ImpPrmIntersection
 176 //purpose  :
 177 //=======================================================================
 178 IntPatch_ImpPrmIntersection::IntPatch_ImpPrmIntersection ()
 179   : done(Standard_False),
 180   empt(Standard_False),
 181   myIsStartPnt(Standard_False),
 182   myUStart(0.0),
 183   myVStart(0.0)
 184 { }
 185
 186
 187 //=======================================================================
 188 //function : IntPatch_ImpPrmIntersection
 189 //purpose  :
 190 //=======================================================================
 191
 192 IntPatch_ImpPrmIntersection::IntPatch_ImpPrmIntersection
 193   (const Handle(Adaptor3d_Surface)&    Surf1,
 194   const Handle(Adaptor3d_TopolTool)&   D1,
 195   const Handle(Adaptor3d_Surface)&    Surf2,
 196   const Handle(Adaptor3d_TopolTool)&   D2,
 197   const Standard_Real    TolArc,
 198   const Standard_Real    TolTang,
 199   const Standard_Real    Fleche,
 200   const Standard_Real    Pas)
 201   : done(Standard_False),
 202   empt(Standard_False),
 203   myIsStartPnt(Standard_False),
 204   myUStart(0.0),
 205   myVStart(0.0)
 206 {
 207   Perform(Surf1,D1,Surf2,D2,TolArc,TolTang,Fleche,Pas);
 208 }
 209
 210
 211 //=======================================================================
 212 //function : SetStartPoint
 213 //purpose  :
 214 //=======================================================================
 215
 216 void IntPatch_ImpPrmIntersection::SetStartPoint(const Standard_Real U,
 217   const Standard_Real V)
 218 {
 219   myIsStartPnt = Standard_True;
 220   myUStart = U; myVStart = V;
 221 }
 222
 223 //=======================================================================
 224 //function : ComputeTangency
 225 //purpose  :
 226 //=======================================================================
 227 void ComputeTangency (const IntPatch_TheSOnBounds& solrst,
 228   IntSurf_SequenceOfPathPoint& seqpdep,
 229   const Handle(Adaptor3d_TopolTool)& Domain,
 230   IntPatch_TheSurfFunction& Func,
 231   const Handle(Adaptor3d_Surface)& PSurf,
 232   TColStd_Array1OfInteger& Destination)
 233 {
 234   Standard_Integer i,k, NbPoints, seqlength;
 235   Standard_Real theparam,test;
 236   Standard_Boolean fairpt, ispassing;
 237   TopAbs_Orientation arcorien,vtxorien;
 238   Handle(Adaptor2d_Curve2d) thearc;
 239   Handle(Adaptor3d_HVertex) vtx,vtxbis;
 240   //Standard_Boolean ispassing;
 241   IntPatch_ThePathPointOfTheSOnBounds PStart;
 242   IntSurf_PathPoint PPoint;
 243   gp_Vec vectg;
 244   gp_Dir2d dirtg;
 245   gp_Pnt ptbid;
 246   gp_Vec d1u,d1v,v1,v2;
 247   gp_Pnt2d p2d;
 248   gp_Vec2d d2d;
 249   //
 250   double aX[2], aF[1], aD[1][2];
 251   math_Vector X(aX, 1, 2);
 252   math_Vector F(aF, 1, 1);
 253   math_Matrix D(aD, 1, 1, 1, 2);
 254   //
 255   seqlength = 0;
 256   NbPoints = solrst.NbPoints();
 257   for (i=1; i<= NbPoints; i++) {
 258     if (Destination(i) == 0) {
 259       PStart = solrst.Point(i);
 260       thearc   = PStart.Arc();
 261       theparam = PStart.Parameter();
 262       arcorien = Domain->Orientation(thearc);
 263       ispassing = (arcorien == TopAbs_INTERNAL ||
 264         arcorien == TopAbs_EXTERNAL);
 265
 266       thearc->D0(theparam,p2d);
 267       X(1) = p2d.X();
 268       X(2) = p2d.Y();
 269       PPoint.SetValue(PStart.Value(),X(1),X(2));
 270
 271       Func.Values(X,F,D);
 272       if (Func.IsTangent()) {
 273         PPoint.SetTangency(Standard_True);
 274         Destination(i) = seqlength+1;
 275         if (!PStart.IsNew()) {
 276           vtx = PStart.Vertex();
 277           for (k=i+1; k<=NbPoints; k++) {
 278             if (Destination(k) ==0) {
 279               PStart = solrst.Point(k);
 280               if (!PStart.IsNew()) {
 281                 vtxbis = PStart.Vertex();
 282                 if (Domain->Identical(vtx,vtxbis)) {
 283                   thearc   = PStart.Arc();
 284                   theparam = PStart.Parameter();
 285                   arcorien = Domain->Orientation(thearc);
 286                   ispassing = ispassing && (arcorien == TopAbs_INTERNAL ||
 287                     arcorien == TopAbs_EXTERNAL);
 288
 289                   thearc->D0(theparam,p2d);
 290                   PPoint.AddUV(p2d.X(),p2d.Y());
 291                   Destination(k) = seqlength+1;
 292                 }
 293               }
 294             }
 295           }
 296         }
 297         PPoint.SetPassing(ispassing);
 298         seqpdep.Append(PPoint);
 299         seqlength++;
 300       }
 301       else { // on a un point de depart potentiel
 302
 303         vectg = Func.Direction3d();
 304         dirtg = Func.Direction2d();
 305
 306         PSurf->D1(X(1),X(2),ptbid,d1u,d1v);
 307         thearc->D1(theparam,p2d,d2d);
 308         v2.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 309         v1 = d1u.Crossed(d1v);
 310
 311         test = vectg.Dot(v1.Crossed(v2));
 312         if (PStart.IsNew()) {
 313           if ((test < 0. && arcorien == TopAbs_FORWARD) ||
 314             (test > 0. && arcorien == TopAbs_REVERSED)) {
 315               vectg.Reverse();
 316               dirtg.Reverse();
 317           }
 318           PPoint.SetDirections(vectg,dirtg);
 319           PPoint.SetPassing(ispassing);
 320           Destination(i) = seqlength+1;
 321           seqpdep.Append(PPoint);
 322           seqlength++;
 323         }
 324         else { // traiter la transition complexe
 325           gp_Dir bidnorm(1.,1.,1.);
 326           Standard_Real tole = 1.e-8;
 327           TopAbs_Orientation LocTrans;
 328           TopTrans_CurveTransition comptrans;
 329           comptrans.Reset(vectg,bidnorm,0.);
 330           if (arcorien == TopAbs_FORWARD ||
 331             arcorien == TopAbs_REVERSED) {
 332               // pour essai
 333
 334               vtx = PStart.Vertex();
 335               vtxorien = Domain->Orientation(vtx);
 336               if (Abs(test) <= tole) {
 337                 LocTrans = TopAbs_EXTERNAL; // et pourquoi pas INTERNAL
 338               }
 339               else {
 340                 if (((test > 0.)&& arcorien == TopAbs_FORWARD) ||
 341                   ((test < 0.)&& arcorien == TopAbs_REVERSED)){
 342                     LocTrans = TopAbs_FORWARD;
 343                 }
 344                 else {
 345                   LocTrans = TopAbs_REVERSED;
 346                 }
 347                 if (arcorien == TopAbs_REVERSED) {v2.Reverse();}
 348               }
 349
 350               comptrans.Compare(tole,v2,bidnorm,0.,LocTrans,vtxorien);
 351           }
 352           Destination(i) = seqlength+1;
 353           for (k= i+1; k<=NbPoints; k++) {
 354             if (Destination(k) == 0) {
 355               PStart = solrst.Point(k);
 356               if (!PStart.IsNew()) {
 357                 vtxbis = PStart.Vertex();
 358                 if (Domain->Identical(vtx,vtxbis)) {
 359                   thearc   = PStart.Arc();
 360                   theparam = PStart.Parameter();
 361                   arcorien = Domain->Orientation(thearc);
 362
 363                   PPoint.AddUV(X(1),X(2));
 364
 365                   thearc->D1(theparam,p2d,d2d);
 366                   PPoint.AddUV(p2d.X(),p2d.Y());
 367
 368                   if (arcorien == TopAbs_FORWARD ||
 369                     arcorien == TopAbs_REVERSED) {
 370                       ispassing = Standard_False;
 371                       v2.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 372
 373                       test = vectg.Dot(v1.Crossed(v2));
 374                       vtxorien = Domain->Orientation(PStart.Vertex());
 375                       if (Abs(test) <= tole) {
 376                         LocTrans = TopAbs_EXTERNAL; // et pourquoi pas INTERNAL
 377                       }
 378                       else {
 379                         if (((test > 0.)&& arcorien == TopAbs_FORWARD) ||
 380                           ((test < 0.)&& arcorien == TopAbs_REVERSED)){
 381                             LocTrans = TopAbs_FORWARD;
 382                         }
 383                         else {
 384                           LocTrans = TopAbs_REVERSED;
 385                         }
 386                         if (arcorien == TopAbs_REVERSED) {v2.Reverse();}
 387                       }
 388
 389                       comptrans.Compare(tole,v2,bidnorm,0.,LocTrans,vtxorien);
 390                   }
 391                   Destination(k) = seqlength+1;
 392                 }
 393               }
 394             }
 395           }
 396           fairpt = Standard_True;
 397           if (!ispassing) {
 398             TopAbs_State Before = comptrans.StateBefore();
 399             TopAbs_State After  = comptrans.StateAfter();
 400             if ((Before == TopAbs_UNKNOWN)||(After == TopAbs_UNKNOWN)) {
 401               fairpt = Standard_False;
 402             }
 403             else if (Before == TopAbs_IN) {
 404               if (After == TopAbs_IN) {
 405                 ispassing = Standard_True;
 406               }
 407               else {
 408                 vectg.Reverse();
 409                 dirtg.Reverse();
 410               }
 411             }
 412             else {
 413               if (After !=TopAbs_IN) {
 414                 fairpt = Standard_False;
 415               }
 416             }
 417           }
 418           if (fairpt) {
 419             PPoint.SetDirections(vectg,dirtg);
 420             PPoint.SetPassing(ispassing);
 421             seqpdep.Append(PPoint);
 422             seqlength++;
 423           }
 424           else { // il faut remettre en "ordre" si on ne garde pas le point.
 425             for (k=i; k <=NbPoints ; k++) {
 426               if (Destination(k)==seqlength + 1) {
 427                 Destination(k) = -Destination(k);
 428               }
 429             }
 430           }
 431         }
 432       }
 433     }
 434   }
 435 }
 436 //=======================================================================
 437 //function : Recadre
 438 //purpose  :
 439 //=======================================================================
 440 void Recadre(const Standard_Boolean ,
 441   GeomAbs_SurfaceType typeS1,
 442   GeomAbs_SurfaceType typeS2,
 443   IntPatch_Point&  pt,
 444   const Handle(IntPatch_TheIWLineOfTheIWalking)& iwline,
 445   Standard_Integer Param,
 446   Standard_Real U1,
 447   Standard_Real V1,
 448   Standard_Real U2,
 449   Standard_Real V2)
 450 {
 451   Standard_Real U1p,V1p,U2p,V2p;
 452   iwline->Line()->Value(Param).Parameters(U1p,V1p,U2p,V2p);
 453   switch(typeS1)
 454   {
 455   case GeomAbs_Torus:
 456     while(V1<(V1p-1.5*M_PI)) V1+=M_PI+M_PI;
 457     while(V1>(V1p+1.5*M_PI)) V1-=M_PI+M_PI;
 458     Standard_FALLTHROUGH
 459   case GeomAbs_Cylinder:
 460   case GeomAbs_Cone:
 461   case GeomAbs_Sphere:
 462     while(U1<(U1p-1.5*M_PI)) U1+=M_PI+M_PI;
 463     while(U1>(U1p+1.5*M_PI)) U1-=M_PI+M_PI;
 464   default:
 465     break;
 466   }
 467   switch(typeS2)
 468   {
 469   case GeomAbs_Torus:
 470     while(V2<(V2p-1.5*M_PI)) V2+=M_PI+M_PI;
 471     while(V2>(V2p+1.5*M_PI)) V2-=M_PI+M_PI;
 472     Standard_FALLTHROUGH
 473   case GeomAbs_Cylinder:
 474   case GeomAbs_Cone:
 475   case GeomAbs_Sphere:
 476     while(U2<(U2p-1.5*M_PI)) U2+=M_PI+M_PI;
 477     while(U2>(U2p+1.5*M_PI)) U2-=M_PI+M_PI;
 478   default:
 479     break;
 480   }
 481   pt.SetParameters(U1,V1,U2,V2);
 482 }
 483
 484 //=======================================================================
 485 //function : Perform
 486 //purpose  :
 487 //=======================================================================
 488 void IntPatch_ImpPrmIntersection::Perform (const Handle(Adaptor3d_Surface)& Surf1,
 489   const Handle(Adaptor3d_TopolTool)& D1,
 490   const Handle(Adaptor3d_Surface)& Surf2,
 491   const Handle(Adaptor3d_TopolTool)& D2,
 492   const Standard_Real TolArc,
 493   const Standard_Real TolTang,
 494   const Standard_Real Fleche,
 495   const Standard_Real Pas)
 496 {
 497   Standard_Boolean reversed, procf, procl, dofirst, dolast;
 498   Standard_Integer indfirst = 0, indlast = 0, ind2, NbSegm;
 499   Standard_Integer NbPointIns, NbPointRst, Nblines, Nbpts, NbPointDep;
 500   Standard_Real U1,V1,U2,V2,paramf,paraml,currentparam;
 501
 502   IntPatch_TheSegmentOfTheSOnBounds thesegm;
 503   IntSurf_PathPoint PPoint;
 504
 505   Handle(IntPatch_RLine) rline;
 506   Handle(IntPatch_WLine) wline;
 507   IntPatch_ThePathPointOfTheSOnBounds PStart,PStartf,PStartl;
 508   IntPatch_Point ptdeb,ptfin,ptbis;
 509
 510   IntPatch_IType typ;
 511   IntSurf_Transition TLine,TArc;
 512   IntSurf_TypeTrans trans1,trans2;
 513   gp_Pnt valpt,ptbid;
 514   gp_Vec tgline,tgrst,norm1,norm2,d1u,d1v;
 515   gp_Dir DirNormale;
 516   gp_Vec VecNormale;
 517
 518   gp_Pnt2d p2d;
 519   gp_Vec2d d2d;
 520
 521   Handle(Adaptor2d_Curve2d) currentarc;
 522   GeomAbs_SurfaceType typeS1, typeS2;
 523   IntSurf_Quadric Quad;
 524   IntPatch_TheSurfFunction Func;
 525   IntPatch_ArcFunction AFunc;
 526   //
 527   typeS1 = Surf1->GetType();
 528   typeS2 = Surf2->GetType();
 529
 530   paramf =0.;
 531   paraml =0.;
 532   trans1 = IntSurf_Undecided;
 533   trans2 = IntSurf_Undecided;
 534   //
 535   done = Standard_False;
 536   empt = Standard_True;
 537   slin.Clear();
 538   spnt.Clear();
 539   //
 540   reversed = Standard_False;
 541   switch (typeS1)
 542   {
 543   case GeomAbs_Plane:
 544     Quad.SetValue(Surf1->Plane());
 545     break;
 546
 547   case GeomAbs_Cylinder:
 548     Quad.SetValue(Surf1->Cylinder());
 549     break;
 550
 551   case GeomAbs_Sphere:
 552     Quad.SetValue(Surf1->Sphere());
 553     break;
 554
 555   case GeomAbs_Cone:
 556     Quad.SetValue(Surf1->Cone());
 557     break;
 558
 559   default:
 560     {
 561       reversed = Standard_True;
 562       switch (typeS2)
 563       {
 564       case GeomAbs_Plane:
 565         Quad.SetValue(Surf2->Plane());
 566         break;
 567
 568       case GeomAbs_Cylinder:
 569         Quad.SetValue(Surf2->Cylinder());
 570         break;
 571
 572       case GeomAbs_Sphere:
 573         Quad.SetValue(Surf2->Sphere());
 574         break;
 575
 576       case GeomAbs_Cone:
 577         Quad.SetValue(Surf2->Cone());
 578         break;
 579       default:
 580         {
 581           throw Standard_ConstructionError();
 582           break;
 583         }
 584       }
 585     }
 586     break;
 587   }
 588   //
 589   Standard_Real aLocalPas = Pas;
 590    GeomAbs_SurfaceType aSType = reversed ? Surf1->GetType() : Surf2->GetType();
 591
 592   if (aSType == GeomAbs_BezierSurface || aSType == GeomAbs_BSplineSurface)
 593   {
 594     Standard_Real aMinRes = Precision::Infinite();
 595     GeomAbs_Shape aCont = GeomAbs_C0;
 596     Standard_Integer aMaxDeg = 0;
 597     const Standard_Real aLimRes = 1.e-10;
 598
 599     if (reversed)
 600     {
 601       aMinRes = Min(Surf1->UResolution(Precision::Confusion()),
 602         Surf1->VResolution(Precision::Confusion()));
 603       aCont = (GeomAbs_Shape)Min(Surf1->UContinuity(), Surf1->VContinuity());
 604       aMaxDeg = Max(Surf1->UDegree(), Surf1->VDegree());
 605     }
 606     else
 607     {
 608       aMinRes = Min(Surf2->UResolution(Precision::Confusion()),
 609         Surf2->VResolution(Precision::Confusion()));
 610       aCont = (GeomAbs_Shape)Min(Surf2->UContinuity(), Surf2->VContinuity());
 611       aMaxDeg = Max(Surf2->UDegree(), Surf2->VDegree());
 612     }
 613
 614     if (aMinRes < aLimRes && aCont > GeomAbs_C0 && aMaxDeg > 3)
 615     {
 616       aLocalPas = Min(Pas, 0.0001);
 617     }
 618   }
 619
 620   Func.SetImplicitSurface(Quad);
 621   Func.Set(IntSurf_QuadricTool::Tolerance(Quad));
 622   AFunc.SetQuadric(Quad);
 623   //
 624   if (!reversed) {
 625     Func.Set(Surf2);
 626     AFunc.Set(Surf2);
 627   }
 628   else {
 629     Func.Set(Surf1);
 630     AFunc.Set(Surf1);
 631   }
 632   //
 633   if (!reversed) {
 634     solrst.Perform(AFunc,D2,TolArc,TolTang);
 635   }
 636   else {
 637     solrst.Perform(AFunc,D1,TolArc,TolTang);
 638   }
 639   if (!solrst.IsDone()) {
 640     return;
 641   }
 642   //
 643   IntSurf_SequenceOfPathPoint seqpdep;
 644   IntSurf_SequenceOfInteriorPoint seqpins;
 645   //
 646   NbPointRst = solrst.NbPoints();
 647   TColStd_Array1OfInteger Destination(1,NbPointRst+1); Destination.Init(0);
 648   if (NbPointRst) {
 649     if (!reversed) {
 650       ComputeTangency(solrst,seqpdep,D2,Func,Surf2,Destination);
 651     }
 652     else {
 653       ComputeTangency(solrst,seqpdep,D1,Func,Surf1,Destination);
 654     }
 655   }
 656   //
 657   Standard_Boolean SearchIns = Standard_True;
 658   if(Quad.TypeQuadric() == GeomAbs_Plane && solrst.NbSegments() > 0)
 659   {
 660     //For such kind of cases it is possible that whole surface is on one side of plane,
 661     //plane only touches surface and does not cross it,
 662     //so no inner points exist.
 663     SearchIns = Standard_False;
 664     Handle(Adaptor3d_TopolTool) T;
 665     if(reversed)
 666     {
 667       T = D1;
 668     }
 669     else
 670     {
 671       T = D2;
 672     }
 673     Standard_Integer aNbSamples = 0;
 674     aNbSamples = T->NbSamples();
 675     gp_Pnt2d s2d;
 676     gp_Pnt s3d;
 677     Standard_Real aValf[1], aUVap[2];
 678     math_Vector Valf(aValf,1,1), UVap(aUVap,1,2);
 679     T->SamplePoint(1,s2d, s3d);
 680     UVap(1)=s2d.X();
 681     UVap(2)=s2d.Y();
 682     Func.Value(UVap,Valf);
 683     Standard_Real rvalf = Sign(1.,Valf(1));
 684     for(Standard_Integer i = 2; i <= aNbSamples; ++i)
 685     {
 686       T->SamplePoint(i,s2d, s3d);
 687       UVap(1)=s2d.X();
 688       UVap(2)=s2d.Y();
 689       Func.Value(UVap,Valf);
 690       if(rvalf * Valf(1) < 0.)
 691       {
 692         SearchIns = Standard_True;
 693         break;
 694       }
 695     }
 696   }
 697   // Recherche des points interieurs
 698   NbPointIns = 0;
 699   if(SearchIns) {
 700     if (!reversed) {
 701       if (myIsStartPnt)
 702         solins.Perform(Func,Surf2,myUStart,myVStart);
 703       else
 704         solins.Perform(Func,Surf2,D2,TolTang);
 705     }
 706     else {
 707       if (myIsStartPnt)
 708         solins.Perform(Func,Surf1,myUStart,myVStart);
 709       else
 710         solins.Perform(Func,Surf1,D1,TolTang);
 711     }
 712     NbPointIns = solins.NbPoints();
 713     for (Standard_Integer i=1; i <= NbPointIns; i++) {
 714       seqpins.Append(solins.Value(i));
 715     }
 716   }
 717   //
 718   NbPointDep=seqpdep.Length();
 719   //
 720   if (NbPointDep || NbPointIns) {
 721     IntPatch_TheIWalking iwalk(TolTang, Fleche, aLocalPas);
 722     iwalk.Perform(seqpdep, seqpins, Func, reversed ? Surf1 : Surf2, reversed);
 723
 724     if(!iwalk.IsDone()) {
 725       return;
 726     }
 727
 728     Standard_Real Vmin, Vmax, TolV = 1.e-14;
 729     if (!reversed) { //Surf1 is quadric
 730       Vmin = Surf1->FirstVParameter();
 731       Vmax = Surf1->LastVParameter();
 732     }
 733     else { //Surf2 is quadric
 734       Vmin = Surf2->FirstVParameter();
 735       Vmax = Surf2->LastVParameter();
 736     }
 737     //
 738     Nblines = iwalk.NbLines();
 739     for (Standard_Integer j=1; j<=Nblines; j++) {
 740       const Handle(IntPatch_TheIWLineOfTheIWalking)&  iwline  = iwalk.Value(j);
 741       const Handle(IntSurf_LineOn2S)&                 thelin  = iwline->Line();
 742
 743       Nbpts = thelin->NbPoints();
 744       if(Nbpts>=2) {
 745         Standard_Integer k = 0;
 746         tgline = iwline->TangentVector(k);
 747         if(k>=1 && k<=Nbpts) { } else { k=Nbpts>>1; }
 748         valpt = thelin->Value(k).Value();
 749
 750         if (!reversed) {
 751           thelin->Value(k).ParametersOnS2(U2,V2);
 752           norm1 = Quad.Normale(valpt);
 753           Surf2->D1(U2,V2,ptbid,d1u,d1v);
 754           norm2 = d1u.Crossed(d1v);
 755         }
 756         else {
 757           thelin->Value(k).ParametersOnS1(U2,V2);
 758           norm2 = Quad.Normale(valpt);
 759           Surf1->D1(U2,V2,ptbid,d1u,d1v);
 760           norm1 = d1u.Crossed(d1v);
 761         }
 762         if (tgline.DotCross(norm2,norm1) > 0.) {
 763           trans1 = IntSurf_Out;
 764           trans2 = IntSurf_In;
 765         }
 766         else {
 767           trans1 = IntSurf_In;
 768           trans2 = IntSurf_Out;
 769         }
 770
 771         //
 772         Standard_Real AnU1,AnU2,AnV2;
 773
 774         GeomAbs_SurfaceType typQuad = Quad.TypeQuadric();
 775         Standard_Boolean arecadr=Standard_False;
 776         valpt = thelin->Value(1).Value();
 777         Quad.Parameters(valpt,AnU1,V1);
 778
 779         if((V1 < Vmin) && (Vmin-V1 < TolV)) V1 = Vmin;
 780         if((V1 > Vmax) && (V1-Vmax < TolV)) V1 = Vmax;
 781
 782         if(reversed) {
 783           thelin->SetUV(1,Standard_False,AnU1,V1); //-- on va lire u2,v2
 784           thelin->Value(1).ParametersOnS1(AnU2,AnV2);
 785         }
 786         else {
 787           thelin->SetUV(1,Standard_True,AnU1,V1);  //-- on va lire u1,v1
 788           thelin->Value(1).ParametersOnS2(AnU2,AnV2);
 789         }
 790
 791         if(typQuad==GeomAbs_Cylinder ||
 792           typQuad==GeomAbs_Cone ||
 793           typQuad==GeomAbs_Sphere) {
 794             arecadr=Standard_True;
 795         }
 796         //
 797         for (k=2; k<=Nbpts; ++k) {
 798           valpt = thelin->Value(k).Value();
 799           Quad.Parameters(valpt,U1,V1);
 800           //
 801           if((V1 < Vmin) && (Vmin-V1 < TolV)) {
 802             V1 = Vmin;
 803           }
 804           if((V1 > Vmax) && (V1-Vmax < TolV)) {
 805             V1 = Vmax;
 806           }
 807           //
 808           if(arecadr) {
 809             //modified by NIZNHY-PKV Fri Mar 28 15:06:01 2008f
 810             Standard_Real aCf, aTwoPI;
 811             //
 812             aCf=0.;
 813             aTwoPI=M_PI+M_PI;
 814             if ((U1-AnU1) >  1.5*M_PI) {
 815               while ((U1-AnU1) > (1.5*M_PI+aCf*aTwoPI)) {
 816                 aCf=aCf+1.;
 817               }
 818               U1=U1-aCf*aTwoPI;
 819             }
 820             //
 821             else {
 822               while ((U1-AnU1) < (-1.5*M_PI-aCf*aTwoPI)) {
 823                 aCf=aCf+1.;
 824               }
 825               U1=U1+aCf*aTwoPI;
 826             }
 827             // was:
 828             //if ((U1-AnU1) >  1.5*M_PI) {
 829             //  U1-=M_PI+M_PI;
 830             //}
 831             //else if ((U1-AnU1) < -1.5*M_PI) {
 832             //  U1+=M_PI+M_PI;
 833             //}
 834             //modified by NIZNHY-PKV Fri Mar 28 15:06:11 2008t
 835           }
 836           //
 837           if(reversed) {
 838             thelin->SetUV(k,Standard_False,U1,V1);
 839
 840             thelin->Value(k).ParametersOnS1(U2,V2);
 841             switch(typeS1) {
 842             case GeomAbs_Cylinder:
 843             case GeomAbs_Cone:
 844             case GeomAbs_Sphere:
 845             case GeomAbs_Torus:
 846               while(U2<(AnU2-1.5*M_PI)) U2+=M_PI+M_PI;
 847               while(U2>(AnU2+1.5*M_PI)) U2-=M_PI+M_PI;
 848               break;
 849             default:
 850               break;
 851             }
 852             if(typeS2==GeomAbs_Torus) {
 853               while(V2<(AnV2-1.5*M_PI)) V2+=M_PI+M_PI;
 854               while(V2>(AnV2+1.5*M_PI)) V2-=M_PI+M_PI;
 855             }
 856             thelin->SetUV(k,Standard_True,U2,V2);
 857           }
 858           else {
 859             thelin->SetUV(k,Standard_True,U1,V1);
 860
 861             thelin->Value(k).ParametersOnS2(U2,V2);
 862             switch(typeS2) {
 863             case GeomAbs_Cylinder:
 864             case GeomAbs_Cone:
 865             case GeomAbs_Sphere:
 866             case GeomAbs_Torus:
 867               while(U2<(AnU2-1.5*M_PI)) U2+=M_PI+M_PI;
 868               while(U2>(AnU2+1.5*M_PI)) U2-=M_PI+M_PI;
 869               break;
 870             default:
 871               break;
 872             }
 873             if(typeS2==GeomAbs_Torus) {
 874               while(V2<(AnV2-1.5*M_PI)) V2+=M_PI+M_PI;
 875               while(V2>(AnV2+1.5*M_PI)) V2-=M_PI+M_PI;
 876             }
 877             thelin->SetUV(k,Standard_False,U2,V2);
 878
 879           }
 880
 881           AnU1=U1;
 882           AnU2=U2;
 883           AnV2=V2;
 884         }
 885         // <-A
 886         wline = new IntPatch_WLine(thelin,Standard_False,trans1,trans2);
 887         wline->SetCreatingWayInfo(IntPatch_WLine::IntPatch_WLImpPrm);
 888
 889 #ifdef INTPATCH_IMPPRMINTERSECTION_DEBUG
 890         wline->Dump(0);
 891 #endif
 892
 893         if (   iwline->HasFirstPoint()
 894           && iwline->IsTangentAtBegining() == Standard_False)
 895         {
 896           indfirst = iwline->FirstPointIndex();
 897           PPoint = seqpdep(indfirst);
 898           tgline = PPoint.Direction3d();
 899           Standard_Integer themult = PPoint.Multiplicity();
 900           for (Standard_Integer i=NbPointRst; i>=1; i--) {
 901             if (Destination(i) == indfirst) {
 902               if (!reversed) { //-- typeS1 = Pln || Cyl || Sph || Cone
 903                 Quad.Parameters(PPoint.Value(),U1,V1);
 904
 905                 if((V1 < Vmin) && (Vmin-V1 < TolV)) V1 = Vmin;
 906                 if((V1 > Vmax) && (V1-Vmax < TolV)) V1 = Vmax;
 907
 908                 PPoint.Parameters(themult,U2,V2);
 909                 Surf2->D1(U2,V2,ptbid,d1u,d1v); //-- @@@@
 910               }
 911               else {  //-- typeS1 != Pln && Cyl && Sph && Cone
 912                 Quad.Parameters(PPoint.Value(),U2,V2);
 913
 914                 if((V2 < Vmin) && (Vmin-V2 < TolV)) V2 = Vmin;
 915                 if((V2 > Vmax) && (V2-Vmax < TolV)) V2 = Vmax;
 916
 917                 PPoint.Parameters(themult,U1,V1);
 918                 Surf1->D1(U1,V1,ptbid,d1u,d1v); //-- @@@@
 919               }
 920
 921               VecNormale = d1u.Crossed(d1v);
 922               //-- Modif du 27 Septembre 94 (Recadrage des pts U,V)
 923               ptdeb.SetValue(PPoint.Value(),TolArc,Standard_False);
 924               ptdeb.SetParameters(U1,V1,U2,V2);
 925               ptdeb.SetParameter(1.);
 926
 927               Recadre(reversed,typeS1,typeS2,ptdeb,iwline,1,U1,V1,U2,V2);
 928
 929               currentarc = solrst.Point(i).Arc();
 930               currentparam = solrst.Point(i).Parameter();
 931               currentarc->D1(currentparam,p2d,d2d);
 932               tgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 933
 934               Standard_Real squaremagnitudeVecNormale = VecNormale.SquareMagnitude();
 935               if(squaremagnitudeVecNormale > 1e-13) {
 936                 DirNormale=VecNormale;
 937                 IntSurf::MakeTransition(tgline,tgrst,DirNormale,TLine,TArc);
 938               }
 939               else {
 940                 TLine.SetValue(Standard_True,IntSurf_Undecided);
 941                 TArc.SetValue(Standard_True,IntSurf_Undecided);
 942               }
 943
 944               ptdeb.SetArc(reversed,currentarc,currentparam,TLine,TArc);
 945               if (!solrst.Point(i).IsNew()) {
 946                 ptdeb.SetVertex(reversed,solrst.Point(i).Vertex());
 947               }
 948               wline->AddVertex(ptdeb);
 949               if (themult == 0) {
 950                 wline->SetFirstPoint(wline->NbVertex());
 951               }
 952
 953               themult--;
 954             }
 955           }
 956         }
 957         else if (iwline->IsTangentAtBegining())
 958         {
 959           gp_Pnt psol = thelin->Value(1).Value();
 960           thelin->Value(1).ParametersOnS1(U1,V1);
 961           thelin->Value(1).ParametersOnS2(U2,V2);
 962           ptdeb.SetValue(psol,TolArc,Standard_True);
 963           ptdeb.SetParameters(U1,V1,U2,V2);
 964           ptdeb.SetParameter(1.);
 965           wline->AddVertex(ptdeb);
 966           wline->SetFirstPoint(wline->NbVertex());
 967         }
 968         else
 969         {
 970           gp_Pnt psol = thelin->Value(1).Value();
 971           thelin->Value(1).ParametersOnS1(U1,V1);
 972           thelin->Value(1).ParametersOnS2(U2,V2);
 973           ptdeb.SetValue(psol,TolArc,Standard_False);
 974           ptdeb.SetParameters(U1,V1,U2,V2);
 975           ptdeb.SetParameter(1.);
 976           wline->AddVertex(ptdeb);
 977           wline->SetFirstPoint(wline->NbVertex());
 978         }
 979
 980
 981         if (   iwline->HasLastPoint()
 982           && iwline->IsTangentAtEnd() == Standard_False)
 983         {
 984           indlast = iwline->LastPointIndex();
 985           PPoint = seqpdep(indlast);
 986           tgline = PPoint.Direction3d().Reversed();
 987           Standard_Integer themult = PPoint.Multiplicity();
 988           for (Standard_Integer i=NbPointRst; i >=1; i--) {
 989             if (Destination(i) == indlast) {
 990               if (!reversed) {
 991                 Quad.Parameters(PPoint.Value(),U1,V1);
 992
 993                 if((V1 < Vmin) && (Vmin-V1 < TolV)) V1 = Vmin;
 994                 if((V1 > Vmax) && (V1-Vmax < TolV)) V1 = Vmax;
 995
 996                 PPoint.Parameters(themult,U2,V2);
 997                 Surf2->D1(U2,V2,ptbid,d1u,d1v); //-- @@@@
 998                 VecNormale = d1u.Crossed(d1v);                    //-- @@@@
 999               }
1000               else {
1001                 Quad.Parameters(PPoint.Value(),U2,V2);
1002
1003                 if((V2 < Vmin) && (Vmin-V2 < TolV)) V2 = Vmin;
1004                 if((V2 > Vmax) && (V2-Vmax < TolV)) V2 = Vmax;
1005
1006                 PPoint.Parameters(themult,U1,V1);
1007                 Surf1->D1(U1,V1,ptbid,d1u,d1v); //-- @@@@
1008                 VecNormale = d1u.Crossed(d1v);                    //-- @@@@
1009               }
1010
1011               ptfin.SetValue(PPoint.Value(),TolArc,Standard_False);
1012               ptfin.SetParameters(U1,V1,U2,V2);
1013               ptfin.SetParameter(Nbpts);
1014
1015               Recadre(reversed,typeS1,typeS2,ptfin,iwline,Nbpts-1,U1,V1,U2,V2);
1016
1017               currentarc = solrst.Point(i).Arc();
1018               currentparam = solrst.Point(i).Parameter();
1019               currentarc->D1(currentparam,p2d,d2d);
1020               tgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1021
1022
1023               Standard_Real squaremagnitudeVecNormale = VecNormale.SquareMagnitude();
1024               if(squaremagnitudeVecNormale > 1e-13) {
1025                 DirNormale=VecNormale;
1026                 IntSurf::MakeTransition(tgline,tgrst,DirNormale,TLine,TArc);
1027               }
1028               else {
1029                 TLine.SetValue(Standard_True,IntSurf_Undecided);
1030                 TArc.SetValue(Standard_True,IntSurf_Undecided);
1031               }
1032
1033
1034               ptfin.SetArc(reversed,currentarc,currentparam,TLine,TArc);
1035               if (!solrst.Point(i).IsNew()) {
1036                 ptfin.SetVertex(reversed,solrst.Point(i).Vertex());
1037               }
1038               wline->AddVertex(ptfin);
1039               if (themult == 0) {
1040                 wline->SetLastPoint(wline->NbVertex());
1041               }
1042
1043               themult--;
1044             }
1045           }
1046         }
1047         else if (iwline->IsTangentAtEnd())
1048         {
1049           gp_Pnt psol = thelin->Value(Nbpts).Value();
1050           thelin->Value(Nbpts).ParametersOnS1(U1,V1);
1051           thelin->Value(Nbpts).ParametersOnS2(U2,V2);
1052           ptfin.SetValue(psol,TolArc,Standard_True);
1053           ptfin.SetParameters(U1,V1,U2,V2);
1054           ptfin.SetParameter(Nbpts);
1055           wline->AddVertex(ptfin);
1056           wline->SetLastPoint(wline->NbVertex());
1057         }
1058         else
1059         {
1060           gp_Pnt psol = thelin->Value(Nbpts).Value();
1061           thelin->Value(Nbpts).ParametersOnS1(U1,V1);
1062           thelin->Value(Nbpts).ParametersOnS2(U2,V2);
1063           ptfin.SetValue(psol,TolArc,Standard_False);
1064           ptfin.SetParameters(U1,V1,U2,V2);
1065           ptfin.SetParameter(Nbpts);
1066           wline->AddVertex(ptfin);
1067           wline->SetLastPoint(wline->NbVertex());
1068         }
1069         //
1070         // Il faut traiter les points de passage.
1071         slin.Append(wline);
1072       }// if(Nbpts>=2) {
1073     }// for (j=1; j<=Nblines; j++) {
1074
1075     // ON GERE LES RACCORDS ENTRE LIGNES. ELLE NE PEUVENT SE RACCORDER
1076     // QUE SUR DES POINTS DE TANGENCE
1077
1078
1079     Nblines = slin.Length();
1080     for (Standard_Integer j=1; j<=Nblines-1; j++) {
1081       dofirst = dolast = Standard_False;
1082       const  Handle(IntPatch_Line)& slinj = slin(j);
1083       Handle(IntPatch_WLine) wlin1 (Handle(IntPatch_WLine)::DownCast (slinj));
1084       if (wlin1->HasFirstPoint()) {
1085         ptdeb = wlin1->FirstPoint(indfirst);
1086         if (ptdeb.IsTangencyPoint()) {
1087           dofirst = Standard_True;
1088         }
1089       }
1090       if (wlin1->HasLastPoint()) {
1091         ptfin = wlin1->LastPoint(indlast);
1092         if (ptfin.IsTangencyPoint()) {
1093           dolast = Standard_True;
1094         }
1095       }
1096
1097       if (dofirst || dolast) {
1098         for (Standard_Integer k=j+1; k<=Nblines;k++) {
1099           const  Handle(IntPatch_Line)& slink = slin(k);
1100           Handle(IntPatch_WLine) wlin2 (Handle(IntPatch_WLine)::DownCast (slink));
1101           if (wlin2->HasFirstPoint()) {
1102             ptbis = wlin2->FirstPoint(ind2);
1103             if (ptbis.IsTangencyPoint()) {
1104               if (dofirst ) {
1105                 if (ptdeb.Value().Distance(ptbis.Value()) <= TolArc) {
1106                   ptdeb.SetMultiple(Standard_True);
1107                   if (!ptbis.IsMultiple()) {
1108                     ptbis.SetMultiple(Standard_True);
1109                     wlin2->Replace(ind2,ptbis);
1110                   }
1111                 }
1112               }
1113               if (dolast ) {
1114                 if (ptfin.Value().Distance(ptbis.Value()) <= TolArc) {
1115                   ptfin.SetMultiple(Standard_True);
1116                   if (!ptbis.IsMultiple()) {
1117                     ptbis.SetMultiple(Standard_True);
1118                     wlin2->Replace(ind2,ptbis);
1119                   }
1120                 }
1121               }
1122             }
1123           }
1124           if (wlin2->HasLastPoint()) {
1125             ptbis = wlin2->LastPoint(ind2);
1126             if (ptbis.IsTangencyPoint()) {
1127               if (dofirst ) {
1128                 if (ptdeb.Value().Distance(ptbis.Value()) <= TolArc) {
1129                   ptdeb.SetMultiple(Standard_True);
1130                   if (!ptbis.IsMultiple()) {
1131                     ptbis.SetMultiple(Standard_True);
1132                     wlin2->Replace(ind2,ptbis);
1133                   }
1134                 }
1135               }
1136               if (dolast ) {
1137                 if (ptfin.Value().Distance(ptbis.Value()) <= TolArc) {
1138                   ptfin.SetMultiple(Standard_True);
1139                   if (!ptbis.IsMultiple()) {
1140                     ptbis.SetMultiple(Standard_True);
1141                     wlin2->Replace(ind2,ptbis);
1142                   }
1143                 }
1144               }
1145             }
1146           }
1147         }
1148         if(dofirst)
1149           wlin1->Replace(indfirst,ptdeb);
1150         if(dolast)
1151           wlin1->Replace(indlast,ptfin);
1152       }
1153     }
1154   }// if (seqpdep.Length() != 0 || seqpins.Length() != 0) {
1155   //
1156   // Treatment the segments
1157   NbSegm = solrst.NbSegments();
1158   if (NbSegm) {
1159     for(Standard_Integer i=1; i<=NbSegm; i++) {
1160       thesegm = solrst.Segment(i);
1161       //Check if segment is degenerated
1162       if(thesegm.HasFirstPoint() && thesegm.HasLastPoint())
1163       {
1164         Standard_Real tol2 = Precision::Confusion();
1165         tol2 *= tol2;
1166         const gp_Pnt& aPf = thesegm.FirstPoint().Value();
1167         const gp_Pnt& aPl = thesegm.LastPoint().Value();
1168         if(aPf.SquareDistance(aPl) <= tol2)
1169         {
1170           //segment can be degenerated - check inner point
1171           paramf = thesegm.FirstPoint().Parameter();
1172           paraml = thesegm.LastPoint().Parameter();
1173           gp_Pnt2d _p2d =
1174             thesegm.Curve()->Value(.57735 * paramf + 0.42265 * paraml);
1175           gp_Pnt aPm;
1176           if(reversed)
1177           {
1178             Surf1->D0(_p2d.X(), _p2d.Y(), aPm);
1179           }
1180           else
1181           {
1182             Surf2->D0(_p2d.X(), _p2d.Y(), aPm);
1183           }
1184           if(aPm.SquareDistance(aPf) <= tol2)
1185           {
1186             //Degenerated
1187             continue;
1188           }
1189         }
1190       }
1191
1192
1193       //----------------------------------------------------------------------
1194       // on cree une ligne d intersection contenant uniquement le segment.
1195       // VOIR POUR LA TRANSITION DE LA LIGNE
1196       // On ajoute aussi un polygone pour le traitement des intersections
1197       // entre ligne et restrictions de la surface implicite (PutVertexOnLine)
1198       //----------------------------------------------------------------------
1199       //-- Calcul de la transition sur la rline (12 fev 97)
1200       //-- reversed a le sens de OnFirst
1201       //--
1202       dofirst = dolast  = Standard_False;
1203       procf = Standard_False;
1204       procl = Standard_False;
1205       IntSurf_Transition TLineUnk,TArcUnk;
1206
1207       IntPatch_Point _thepointAtBeg;
1208       IntPatch_Point _thepointAtEnd;
1209
1210       Standard_Boolean TransitionOK=Standard_False;
1211
1212       if(thesegm.HasFirstPoint()) {
1213         Standard_Real _u1,_v1,_u2,_v2;
1214
1215         dofirst = Standard_True;
1216         PStartf = thesegm.FirstPoint();
1217         paramf = PStartf.Parameter();
1218
1219         gp_Pnt2d _p2d     = thesegm.Curve()->Value(paramf);
1220         Handle(Adaptor3d_HVertex) _vtx;
1221         if(PStartf.IsNew()==Standard_False)
1222           _vtx= PStartf.Vertex();
1223         const gp_Pnt& _Pp = PStartf.Value();
1224         _thepointAtBeg.SetValue(_Pp,PStartf.Tolerance(),Standard_False);
1225         if (!reversed) { //-- typeS1 = Pln || Cyl || Sph || Cone
1226           Quad.Parameters(_Pp,_u1,_v1);
1227           _u2=_p2d.X(); _v2=_p2d.Y();
1228         }
1229         else {  //-- typeS1 != Pln && Cyl && Sph && Cone
1230           Quad.Parameters(_Pp,_u2,_v2);
1231           _u1=_p2d.X(); _v1=_p2d.Y();
1232         }
1233         _thepointAtBeg.SetParameters(_u1,_v1,_u2,_v2);
1234         _thepointAtBeg.SetParameter(paramf);
1235         if(PStartf.IsNew()==Standard_False)
1236           _thepointAtBeg.SetVertex(reversed,_vtx);
1237         _thepointAtBeg.SetArc(reversed,thesegm.Curve(),paramf,TLineUnk,TArcUnk);
1238
1239
1240         gp_Vec d1u1,d1v1,d1u2,d1v2; gp_Vec2d _d2d;
1241         Surf1->D1(_u1,_v1,ptbid,d1u1,d1v1);
1242         norm1 = d1u1.Crossed(d1v1);
1243         Surf2->D1(_u2,_v2,ptbid,d1u2,d1v2);
1244         norm2 = d1u2.Crossed(d1v2);
1245
1246         thesegm.Curve()->D1(paramf,_p2d,_d2d);
1247         if(reversed) {
1248           tgline.SetLinearForm(_d2d.X(),d1u1,_d2d.Y(),d1v1);
1249         }
1250         else {
1251           tgline.SetLinearForm(_d2d.X(),d1u2,_d2d.Y(),d1v2);
1252         }
1253         _u1=tgline.DotCross(norm2,norm1);
1254         TransitionOK=Standard_True;
1255         if (_u1 > 0.00000001) {
1256           trans1 = IntSurf_Out;
1257           trans2 = IntSurf_In;
1258         }
1259         else if(_u1 < -0.00000001) {
1260           trans1 = IntSurf_In;
1261           trans2 = IntSurf_Out;
1262         }
1263         else {
1264           TransitionOK=Standard_False;
1265         }
1266       }
1267       if(thesegm.HasLastPoint()) {
1268         Standard_Real _u1,_v1,_u2,_v2;
1269
1270         dolast = Standard_True;
1271         PStartl = thesegm.LastPoint();
1272         paraml = PStartl.Parameter();
1273
1274         gp_Pnt2d _p2d = thesegm.Curve()->Value(paraml);
1275         Handle(Adaptor3d_HVertex) _vtx;
1276         if(PStartl.IsNew()==Standard_False)
1277           _vtx = PStartl.Vertex();
1278         const gp_Pnt& _Pp = PStartl.Value();
1279         IntPatch_Point _thepoint;
1280         _thepointAtEnd.SetValue(_Pp,PStartl.Tolerance(),Standard_False);
1281         if (!reversed) { //-- typeS1 = Pln || Cyl || Sph || Cone
1282           Quad.Parameters(_Pp,_u1,_v1);
1283           _u2=_p2d.X(); _v2=_p2d.Y();
1284         }
1285         else {  //-- typeS1 != Pln && Cyl && Sph && Cone
1286           Quad.Parameters(_Pp,_u2,_v2);
1287           _u1=_p2d.X(); _v1=_p2d.Y();
1288         }
1289         _thepointAtEnd.SetParameters(_u1,_v1,_u2,_v2);
1290         _thepointAtEnd.SetParameter(paraml);
1291         if(PStartl.IsNew()==Standard_False)
1292           _thepointAtEnd.SetVertex(reversed,_vtx);
1293         _thepointAtEnd.SetArc(reversed,thesegm.Curve(),paraml,TLineUnk,TArcUnk);
1294
1295
1296
1297         gp_Vec d1u1,d1v1,d1u2,d1v2; gp_Vec2d _d2d;
1298         Surf1->D1(_u1,_v1,ptbid,d1u1,d1v1);
1299         norm1 = d1u1.Crossed(d1v1);
1300         Surf2->D1(_u2,_v2,ptbid,d1u2,d1v2);
1301         norm2 = d1u2.Crossed(d1v2);
1302
1303         thesegm.Curve()->D1(paraml,_p2d,_d2d);
1304         if(reversed) {
1305           tgline.SetLinearForm(_d2d.X(),d1u1,_d2d.Y(),d1v1);
1306         }
1307         else {
1308           tgline.SetLinearForm(_d2d.X(),d1u2,_d2d.Y(),d1v2);
1309         }
1310         _u1=tgline.DotCross(norm2,norm1);
1311         TransitionOK=Standard_True;
1312         if (_u1 > 0.00000001) {
1313           trans1 = IntSurf_Out;
1314           trans2 = IntSurf_In;
1315         }
1316         else if(_u1 < -0.00000001) {
1317           trans1 = IntSurf_In;
1318           trans2 = IntSurf_Out;
1319         }
1320         else {
1321           TransitionOK=Standard_False;
1322         }
1323       }
1324       if(TransitionOK==Standard_False) {
1325         //-- rline = new IntPatch_RLine (thesegm.Curve(),reversed,Standard_False);
1326         rline =  new IntPatch_RLine (Standard_False);
1327         if(reversed) {
1328           rline->SetArcOnS1(thesegm.Curve());
1329         }
1330         else {
1331           rline->SetArcOnS2(thesegm.Curve());
1332         }
1333       }
1334       else {
1335         //-- rline = new IntPatch_RLine (thesegm.Curve(),reversed,Standard_False,trans1,trans2);
1336         rline =  new IntPatch_RLine (Standard_False,trans1,trans2);
1337         if(reversed) {
1338           rline->SetArcOnS1(thesegm.Curve());
1339         }
1340         else {
1341           rline->SetArcOnS2(thesegm.Curve());
1342         }
1343       }
1344
1345       //------------------------------
1346       //-- Ajout des points
1347       //--
1348       if (thesegm.HasFirstPoint()) {
1349         rline->AddVertex(_thepointAtBeg);
1350         rline->SetFirstPoint(rline->NbVertex());
1351       }
1352
1353       if (thesegm.HasLastPoint()) {
1354         rline->AddVertex(_thepointAtEnd);
1355         rline->SetLastPoint(rline->NbVertex());
1356       }
1357
1358       // Polygone sur restriction solution
1359       if (dofirst && dolast) {
1360         Standard_Real prm;
1361         gp_Pnt ptpoly;
1362         IntSurf_PntOn2S p2s;
1363         Handle(IntSurf_LineOn2S) Thelin = new IntSurf_LineOn2S ();
1364         Handle(Adaptor2d_Curve2d) arcsegm = thesegm.Curve();
1365         Standard_Integer nbsample = 100;
1366
1367         if (!reversed) {
1368           for (Standard_Integer j=1; j<=nbsample; j++) {
1369             prm = paramf + (j-1)*(paraml-paramf)/(nbsample-1);
1370             arcsegm->D0(prm,p2d);
1371             Surf2->D0(p2d.X(),p2d.Y(),ptpoly);
1372
1373             Quad.Parameters(ptpoly,U1,V1);
1374             p2s.SetValue(ptpoly,U1,V1,p2d.X(),p2d.Y());
1375             Thelin->Add(p2s);
1376           }
1377         }
1378         else {
1379           for (Standard_Integer j=1; j<=nbsample; j++) {
1380             prm = paramf + (j-1)*(paraml-paramf)/(nbsample-1);
1381             arcsegm->D0(prm,p2d);
1382             Surf1->D0(p2d.X(),p2d.Y(),ptpoly);
1383
1384             Quad.Parameters(ptpoly,U2,V2);
1385             p2s.SetValue(ptpoly,p2d.X(),p2d.Y(),U2,V2);
1386             Thelin->Add(p2s);
1387           }
1388         }
1389         rline->Add(Thelin);
1390       }
1391
1392       if (dofirst || dolast) {
1393         Nblines = slin.Length();
1394         for (Standard_Integer j=1; j<=Nblines; j++) {
1395           const Handle(IntPatch_Line)& slinj = slin(j);
1396           typ = slinj->ArcType();
1397           if (typ == IntPatch_Walking) {
1398             Nbpts = Handle(IntPatch_WLine)::DownCast (slinj)->NbVertex();
1399           }
1400           else {
1401             Nbpts = Handle(IntPatch_RLine)::DownCast (slinj)->NbVertex();
1402           }
1403           for (Standard_Integer k=1; k<=Nbpts;k++) {
1404             if (typ == IntPatch_Walking) {
1405               ptdeb = Handle(IntPatch_WLine)::DownCast (slinj)->Vertex(k);
1406             }
1407             else {
1408               ptdeb = Handle(IntPatch_RLine)::DownCast (slinj)->Vertex(k);
1409             }
1410             if (dofirst) {
1411
1412               if (ptdeb.Value().Distance(PStartf.Value()) <=TolArc) {
1413                 ptdeb.SetMultiple(Standard_True);
1414                 if (typ == IntPatch_Walking) {
1415                   Handle(IntPatch_WLine)::DownCast (slinj)->Replace(k,ptdeb);
1416                 }
1417                 else {
1418                   Handle(IntPatch_RLine)::DownCast (slinj)->Replace(k,ptdeb);
1419                 }
1420                 ptdeb.SetParameter(paramf);
1421                 rline->AddVertex(ptdeb);
1422                 if (!procf){
1423                   procf=Standard_True;
1424                   rline->SetFirstPoint(rline->NbVertex());
1425                 }
1426               }
1427             }
1428             if (dolast) {
1429               if(dofirst) { //-- on recharge le ptdeb
1430                 if (typ == IntPatch_Walking) {
1431                   ptdeb = Handle(IntPatch_WLine)::DownCast (slinj)->Vertex(k);
1432                 }
1433                 else {
1434                   ptdeb = Handle(IntPatch_RLine)::DownCast (slinj)->Vertex(k);
1435                 }
1436               }
1437               if (ptdeb.Value().Distance(PStartl.Value()) <=TolArc) {
1438                 ptdeb.SetMultiple(Standard_True);
1439                 if (typ == IntPatch_Walking) {
1440                   Handle(IntPatch_WLine)::DownCast (slinj)->Replace(k,ptdeb);
1441                 }
1442                 else {
1443                   Handle(IntPatch_RLine)::DownCast (slinj)->Replace(k,ptdeb);
1444                 }
1445                 ptdeb.SetParameter(paraml);
1446                 rline->AddVertex(ptdeb);
1447                 if (!procl){
1448                   procl=Standard_True;
1449                   rline->SetLastPoint(rline->NbVertex());
1450                 }
1451               }
1452             }
1453           }
1454         }
1455       }
1456       slin.Append(rline);
1457     }
1458   }// if (NbSegm)
1459   //
1460   // on traite les restrictions de la surface implicite
1461
1462   for (Standard_Integer i=1, aNbLin = slin.Length(); i<=aNbLin; i++)
1463   {
1464     Handle(IntPatch_PointLine) aL = Handle(IntPatch_PointLine)::DownCast(slin(i));
1465
1466     if (!reversed)
1467       IntPatch_RstInt::PutVertexOnLine(aL,Surf1,D1,Surf2,Standard_True,TolTang);
1468     else
1469       IntPatch_RstInt::PutVertexOnLine(aL,Surf2,D2,Surf1,Standard_False,TolTang);
1470
1471     if (aL->NbPnts() <= 2)
1472     {
1473       Standard_Boolean aCond = aL->NbPnts() < 2;
1474       if (!aCond)
1475         aCond = (aL->Point(1).IsSame(aL->Point(2), Precision::Confusion()));
1476
1477       if (aCond)
1478       {
1479         slin.Remove(i);
1480         i--;
1481         aNbLin--;
1482         continue;
1483       }
1484     }
1485
1486     if(aL->ArcType() == IntPatch_Walking)
1487     {
1488       const Handle(IntPatch_WLine) aWL = Handle(IntPatch_WLine)::DownCast(aL);
1489       slin.Append(aWL);
1490       slin.Remove(i);
1491       i--;
1492       aNbLin--;
1493     }
1494   }
1495
1496   // Now slin is filled as follows: lower indices correspond to Restriction line,
1497   // after (higher indices) - only Walking-line.
1498
1499   const Standard_Real aTol3d = Max(Func.Tolerance(), TolTang);
1500   const Handle(Adaptor3d_Surface)& aQSurf = (reversed) ? Surf2 : Surf1;
1501   const Handle(Adaptor3d_Surface)& anOtherSurf = (reversed) ? Surf1 : Surf2;
1502
1503   for (Standard_Integer i = 1; i <= slin.Length(); i++)
1504   {
1505     const Handle(IntPatch_PointLine)& aL1 = Handle(IntPatch_PointLine)::DownCast(slin(i));
1506     const Handle(IntPatch_RLine)& aRL1 = Handle(IntPatch_RLine)::DownCast(aL1);
1507
1508     if(aRL1.IsNull())
1509     {
1510       //Walking-Walking cases are not supported
1511       break;
1512     }
1513
1514     const Handle(Adaptor2d_Curve2d)& anArc = aRL1->IsArcOnS1() ?
1515                                               aRL1->ArcOnS1() :
1516                                               aRL1->ArcOnS2();
1517     if(anArc->GetType() != GeomAbs_Line)
1518     {
1519       //Restriction line must be isoline.
1520       //Other cases are not supported by
1521       //existing algorithms.
1522
1523       break;
1524     }
1525
1526     Standard_Boolean isFirstDeleted = Standard_False;
1527
1528     for(Standard_Integer j = i + 1; j <= slin.Length(); j++)
1529     {
1530       Handle(IntPatch_PointLine) aL2 = Handle(IntPatch_PointLine)::DownCast(slin(j));
1531       Handle(IntPatch_RLine) aRL2 = Handle(IntPatch_RLine)::DownCast(aL2);
1532
1533       //Here aL1 (i-th line) is Restriction-line and aL2 (j-th line) is
1534       //Restriction or Walking
1535
1536       if(!aRL2.IsNull())
1537       {
1538         const Handle(Adaptor2d_Curve2d)& anArc2 = aRL2->IsArcOnS1() ?
1539                                                    aRL2->ArcOnS1() :
1540                                                    aRL2->ArcOnS2();
1541         if(anArc2->GetType() != GeomAbs_Line)
1542         {
1543           //Restriction line must be isoline.
1544           //Other cases are not supported by
1545           //existing algorithms.
1546
1547           continue;
1548         }
1549       }
1550
1551       //aDir can be equal to one of following four values only
1552       //(because Reastriction line is boundary of rectangular surface):
1553       //either {0, 1} or {0, -1} or {1, 0} or {-1, 0}.
1554       const gp_Dir2d aDir = anArc->Line().Direction();
1555
1556       Standard_Real aTol2d = anOtherSurf->UResolution(aTol3d),
1557                     aPeriod = anOtherSurf->IsVPeriodic() ? anOtherSurf->VPeriod() : 0.0;
1558
1559       if(Abs(aDir.X()) < 0.5)
1560       {//Restriction directs along V-direction
1561         aTol2d = anOtherSurf->VResolution(aTol3d);
1562         aPeriod = anOtherSurf->IsUPeriodic() ? anOtherSurf->UPeriod() : 0.0;
1563       }
1564
1565       const Standard_Boolean isCoincide = IsCoincide(Func, aL2, anArc, aRL1->IsArcOnS1(),
1566                                                       aTol3d, aTol2d, aPeriod);
1567
1568       if(isCoincide)
1569       {
1570         if(aRL2.IsNull())
1571         {//Delete Walking-line
1572           slin.Remove(j);
1573           j--;
1574         }
1575         else
1576         {//Restriction-Restriction
1577           const Handle(Adaptor2d_Curve2d)& anArc2 = aRL2->IsArcOnS1() ?
1578                                                      aRL2->ArcOnS1() :
1579                                                      aRL2->ArcOnS2();
1580
1581           const Standard_Real aRange2 = anArc2->LastParameter() -
1582                                         anArc2->FirstParameter();
1583           const Standard_Real aRange1 = anArc->LastParameter() -
1584                                         anArc->FirstParameter();
1585
1586           if(aRange2 > aRange1)
1587           {
1588             isFirstDeleted = Standard_True;
1589             break;
1590           }
1591           else
1592           {//Delete j-th line
1593             slin.Remove(j);
1594             j--;
1595           }
1596         }
1597       }
1598     } //for(Standard_Integer j = i + 1; j <= slin.Length(); j++)
1599
1600     if(isFirstDeleted)
1601     {//Delete i-th line
1602       slin.Remove(i--);
1603     }
1604   }//for (Standard_Integer i = 1; i <= slin.Length(); i++)
1605
1606   empt = (slin.Length() == 0 && spnt.Length() == 0);
1607   done = Standard_True;
1608
1609
1610   if(slin.Length() == 0)
1611     return;
1612
1613   Standard_Boolean isDecomposeRequired =  (Quad.TypeQuadric() == GeomAbs_Cone) ||
1614                                           (Quad.TypeQuadric() == GeomAbs_Sphere) ||
1615                                           (Quad.TypeQuadric() == GeomAbs_Cylinder) ||
1616                                           (Quad.TypeQuadric() == GeomAbs_Torus);
1617
1618   if(!isDecomposeRequired)
1619     return;
1620
1621   // post processing for cones and spheres
1622
1623   const Handle(Adaptor3d_TopolTool)& PDomain = (reversed) ? D1 : D2;
1624
1625   IntPatch_SequenceOfLine dslin;
1626   Standard_Boolean isDecompose = Standard_False;
1627   for(Standard_Integer i = 1; i <= slin.Length(); i++ )
1628   {
1629     if(DecomposeResult( Handle(IntPatch_PointLine)::DownCast(slin(i)),
1630                                         reversed, Quad, PDomain, aQSurf,
1631                                         anOtherSurf, TolArc, aTol3d, dslin))
1632     {
1633       isDecompose = Standard_True;
1634     }
1635   }
1636
1637   if(!isDecompose)
1638     return;
1639
1640   slin.Clear();
1641   for(Standard_Integer i = 1; i <= dslin.Length(); i++ )
1642     slin.Append(dslin(i));
1643 }
1644
1645 // correct U parameter of the start point of line on Quadric
1646 // (change 0->2PI or vs, if necessary)
1647 static Standard_Real AdjustUFirst(Standard_Real U1,Standard_Real U2)
1648 {
1649   Standard_Real u = U1;
1650
1651   // case: no adjustment
1652   if( U1 > 0. && U1 < (2.*M_PI) )
1653     return u;
1654
1655   // case: near '0'
1656   if( U1 == 0. || fabs(U1) <= 1.e-9 ) {
1657     if( U2 > 0. && U2 < (2.*M_PI) )
1658       u = ( U2 < ((2.*M_PI)-U2) ) ? 0. : (2.*M_PI);
1659     else {
1660       Standard_Real uu = U2;
1661       if( U2 > (2.*M_PI) )
1662         while( uu > (2.*M_PI) )
1663           uu -= (2.*M_PI);
1664       else
1665         while( uu < 0.)
1666           uu += (2.*M_PI);
1667
1668       u = ( uu < ((2.*M_PI)-uu) ) ? 0. : (2.*M_PI);
1669     }
1670   }
1671   // case: near '2PI'
1672   else if( U1 == (2.*M_PI) || fabs((2.*M_PI)-fabs(U1)) <= 1.e-9 ) {
1673     if( U2 > 0. && U2 < (2.*M_PI) )
1674       u = ( U2 < ((2.*M_PI)-U2) ) ? 0. : (2.*M_PI);
1675     else {
1676       Standard_Real uu = U2;
1677       if( U2 > (2.*M_PI) )
1678         while( uu > (2.*M_PI) )
1679           uu -= (2.*M_PI);
1680       else
1681         while( uu < 0.)
1682           uu += (2.*M_PI);
1683
1684       u = ( uu < ((2.*M_PI)-uu) ) ? 0. : (2.*M_PI);
1685     }
1686   }
1687   // case: '<0. || >2PI'
1688   else {
1689     if(U1 < 0.)
1690       while(u < 0.)
1691         u += 2.*M_PI;
1692     if(U1 > (2.*M_PI))
1693       while(u > (2.*M_PI))
1694         u -= (2.*M_PI);
1695   }
1696
1697   return u;
1698 }
1699
1700 // collect vertices, reject equals
1701 static Handle(IntSurf_LineOn2S) GetVertices(const Handle(IntPatch_PointLine)& thePLine,
1702                                             const Standard_Real           TOL3D,
1703                                             const Standard_Real           TOL2D)
1704 {
1705   //  Standard_Real TOL3D = 1.e-12, TOL2D = 1.e-8;
1706
1707   Handle(IntSurf_LineOn2S) vertices = new IntSurf_LineOn2S();
1708
1709   Standard_Real U1 = 0., U2 = 0., V1 = 0., V2 = 0.;
1710   Standard_Integer i = 0, k = 0;
1711   Standard_Integer NbVrt = thePLine->NbVertex();
1712
1713   TColStd_Array1OfInteger anVrts(1,NbVrt);
1714   anVrts.Init(0);
1715
1716   // check equal vertices
1717   for(i = 1; i <= NbVrt; i++) {
1718
1719     if( anVrts(i) == -1 ) continue;
1720
1721     const IntPatch_Point& Pi = thePLine->Vertex(i);
1722
1723     for(k = (i+1); k <= NbVrt; k++) {
1724
1725       if( anVrts(k) == -1 ) continue;
1726
1727       const IntPatch_Point& Pk = thePLine->Vertex(k);
1728
1729       if(Pi.Value().Distance(Pk.Value()) <= TOL3D) {
1730         // suggest the points are equal;
1731         // test 2d parameters on surface
1732         Standard_Boolean sameU1 = Standard_False;
1733         Standard_Boolean sameV1 = Standard_False;
1734         Standard_Boolean sameU2 = Standard_False;
1735         Standard_Boolean sameV2 = Standard_False;
1736
1737         Pi.ParametersOnS1(U1,V1);
1738         Pk.ParametersOnS1(U2,V2);
1739         if(fabs(U1-U2) <= TOL2D) sameU1 = Standard_True;
1740         if(fabs(V1-V2) <= TOL2D) sameV1 = Standard_True;
1741
1742         Pi.ParametersOnS2(U1,V1);
1743         Pk.ParametersOnS2(U2,V2);
1744         if(fabs(U1-U2) <= TOL2D) sameU2 = Standard_True;
1745         if(fabs(V1-V2) <= TOL2D) sameV2 = Standard_True;
1746
1747         if((sameU1 && sameV1) && (sameU2 && sameV2))
1748           anVrts(k) = -1;
1749       }
1750     }
1751   }
1752
1753   // copy further processed vertices
1754   for(i = 1; i <= NbVrt; i++) {
1755     if( anVrts(i) == -1 ) continue;
1756     vertices->Add(thePLine->Vertex(i).PntOn2S());
1757   }
1758   return vertices;
1759 }
1760
1761 static void SearchVertices(const Handle(IntSurf_LineOn2S)& Line,
1762   const Handle(IntSurf_LineOn2S)& Vertices,
1763   TColStd_Array1OfInteger&        PTypes)
1764 {
1765   Standard_Integer nbp = Line->NbPoints(), nbv = Vertices->NbPoints();
1766   Standard_Integer ip = 0, iv = 0;
1767   for(ip = 1; ip <= nbp; ip++) {
1768     const IntSurf_PntOn2S& aP = Line->Value(ip);
1769     Standard_Integer type = 0;
1770     for(iv = 1; iv <= nbv; iv++) {
1771       const IntSurf_PntOn2S& aV = Vertices->Value(iv);
1772       if(aP.IsSame(aV, Precision::Confusion(), Precision::PConfusion())) {
1773         type = iv;
1774         break;
1775       }
1776     }
1777     PTypes(ip) = type;
1778   }
1779 }
1780
1781 static inline Standard_Boolean IsSeamParameter(const Standard_Real U,
1782   const Standard_Real TOL2D)
1783 {
1784   return (fabs(U) <= TOL2D || fabs(2.*M_PI - U) <= TOL2D);
1785 }
1786
1787 static inline Standard_Real AdjustU(const Standard_Real U)
1788 {
1789   Standard_Real u = U, DBLPI = 2.*M_PI;
1790   if(u < 0. || u > DBLPI) {
1791     if(u < 0.)
1792       while(u < 0.)
1793         u += DBLPI;
1794     else
1795       while(u > DBLPI)
1796         u -= DBLPI;
1797   }
1798   return u;
1799 }
1800
1801 static inline void Correct2DBounds(const Standard_Real UF,
1802   const Standard_Real UL,
1803   const Standard_Real VF,
1804   const Standard_Real VL,
1805   const Standard_Real TOL2D,
1806   Standard_Real&      U,
1807   Standard_Real&      V)
1808 {
1809   Standard_Real Eps = 1.e-16;
1810   Standard_Real dUF = fabs(U - UF);
1811   Standard_Real dUL = fabs(U - UL);
1812   Standard_Real dVF = fabs(V - VF);
1813   Standard_Real dVL = fabs(V - VL);
1814   if(dUF <= TOL2D && dUF > Eps) U = UF;
1815   if(dUL <= TOL2D && dUL > Eps) U = UL;
1816   if(dVF <= TOL2D && dVF > Eps) V = VF;
1817   if(dVL <= TOL2D && dVL > Eps) V = VL;
1818 }
1819
1820 static void AdjustLine(Handle(IntSurf_LineOn2S)& Line,
1821   const Standard_Boolean    IsReversed,
1822   const Handle(Adaptor3d_Surface)&         QSurf,
1823   const Standard_Real       TOL2D)
1824 {
1825   Standard_Real VF = QSurf->FirstVParameter();
1826   Standard_Real VL = QSurf->LastVParameter();
1827   Standard_Real UF = QSurf->FirstUParameter();
1828   Standard_Real UL = QSurf->LastUParameter();
1829
1830   Standard_Integer nbp = Line->NbPoints(), ip = 0;
1831   Standard_Real U = 0., V = 0.;
1832   for(ip = 1; ip <= nbp; ip++) {
1833     if(IsReversed) {
1834       Line->Value(ip).ParametersOnS2(U,V);
1835       U = AdjustU(U);
1836       Correct2DBounds(UF,UL,VF,VL,TOL2D,U,V);
1837       Line->SetUV(ip,Standard_False,U,V);
1838     }
1839     else {
1840       Line->Value(ip).ParametersOnS1(U,V);
1841       U = AdjustU(U);
1842       Correct2DBounds(UF,UL,VF,VL,TOL2D,U,V);
1843       Line->SetUV(ip,Standard_True,U,V);
1844     }
1845   }
1846 }
1847
1848 static Standard_Boolean InsertSeamVertices(Handle(IntSurf_LineOn2S)&       Line,
1849   const Standard_Boolean          IsReversed,
1850   Handle(IntSurf_LineOn2S)&       Vertices,
1851   const TColStd_Array1OfInteger&  PTypes,
1852   const Standard_Real             TOL2D)
1853 {
1854   Standard_Boolean result = Standard_False;
1855   Standard_Integer ip = 0, nbp = Line->NbPoints();
1856   Standard_Real U = 0., V = 0.;
1857   for(ip = 1; ip <= nbp; ip++) {
1858     Standard_Integer ipt = PTypes(ip);
1859     if(ipt != 0) {
1860       const IntSurf_PntOn2S& aP = Line->Value(ip);
1861       if(IsReversed)
1862         aP.ParametersOnS2(U,V); // S2 - quadric
1863       else
1864         aP.ParametersOnS1(U,V); // S1 - quadric
1865       U = AdjustU(U);
1866       if(IsSeamParameter(U,TOL2D)) {
1867         if(ip == 1 || ip == nbp) {
1868           Standard_Real U1 = 0., V1 = 0.;
1869           Standard_Integer ipp = (ip == 1) ? (ip+1) : (ip-1);
1870           if(IsReversed)
1871             Line->Value(ipp).ParametersOnS2(U1,V1); // S2 - quadric
1872           else
1873             Line->Value(ipp).ParametersOnS1(U1,V1); // S1 - quadric
1874           Standard_Real u = AdjustUFirst(U,U1);
1875           if(fabs(u-U) >= 1.5*M_PI) {
1876             Standard_Real U2 = 0., V2 = 0.;
1877             if(IsReversed) {
1878               Line->Value(ip).ParametersOnS1(U2,V2); // prm
1879               Line->SetUV(ip,Standard_False,u,V);
1880               Line->SetUV(ip,Standard_True,U2,V2);
1881             }
1882             else {
1883               Line->Value(ip).ParametersOnS2(U2,V2); // prm
1884               Line->SetUV(ip,Standard_True,u,V);
1885               Line->SetUV(ip,Standard_False,U2,V2);
1886             }
1887           }
1888         }
1889         else {
1890           Standard_Integer ipp = ip - 1;
1891           Standard_Integer ipn = ip + 1;
1892           Standard_Real U1 = 0., V1 = 0., U2 = 0., V2 = 0.;
1893           if(IsReversed) {
1894             Line->Value(ipp).ParametersOnS2(U1,V1); // quad
1895             Line->Value(ipn).ParametersOnS2(U2,V2); // quad
1896           }
1897           else {
1898             Line->Value(ipp).ParametersOnS1(U1,V1); // quad
1899             Line->Value(ipn).ParametersOnS1(U2,V2); // quad
1900           }
1901           U1 = AdjustU(U1);
1902           U2 = AdjustU(U2);
1903           Standard_Boolean pnearZero = (fabs(U1) < fabs(2.*M_PI-U1)) ? Standard_True : Standard_False;
1904           Standard_Boolean cnearZero = (fabs(U) < fabs(2.*M_PI-U)) ? Standard_True : Standard_False;
1905           if(pnearZero == cnearZero) {
1906             if(!IsSeamParameter(U2,TOL2D) && !IsSeamParameter(U1,TOL2D)) {
1907               Standard_Real nU = (cnearZero) ? (2.*M_PI) : 0.;
1908               IntSurf_PntOn2S nP;
1909               nP.SetValue(aP.Value());
1910               Standard_Real U3 = 0., V3 = 0.;
1911               if(IsReversed) {
1912                 Line->Value(ip).ParametersOnS1(U3,V3); // prm
1913                 nP.SetValue(Standard_False,nU,V);
1914                 nP.SetValue(Standard_True,U3,V3);
1915               }
1916               else {
1917                 Line->Value(ip).ParametersOnS2(U3,V3); // prm
1918                 nP.SetValue(Standard_True,nU,V);
1919                 nP.SetValue(Standard_False,U3,V3);
1920               }
1921               Line->InsertBefore(ipn,nP);
1922               Vertices->Add(nP);
1923               result = Standard_True;
1924               break;
1925             }
1926           }
1927           else {
1928             if(!IsSeamParameter(U2,TOL2D) && !IsSeamParameter(U1,TOL2D)) {
1929               Standard_Real nU = (cnearZero) ? (2.*M_PI) : 0.;
1930               IntSurf_PntOn2S nP;
1931               nP.SetValue(aP.Value());
1932               Standard_Real U3 = 0., V3 = 0.;
1933               if(IsReversed) {
1934                 Line->Value(ip).ParametersOnS1(U3,V3); // prm
1935                 nP.SetValue(Standard_False,nU,V);
1936                 nP.SetValue(Standard_True,U3,V3);
1937               }
1938               else {
1939                 Line->Value(ip).ParametersOnS2(U3,V3); // prm
1940                 nP.SetValue(Standard_True,nU,V);
1941                 nP.SetValue(Standard_False,U3,V3);
1942               }
1943               Line->InsertBefore(ip,nP);
1944               Vertices->Add(nP);
1945               result = Standard_True;
1946               break;
1947             }
1948             else {
1949               // Line->InsertBefore(ip,Line->Value(ipn));
1950               // Line->RemovePoint(ip+2);
1951               // result = Standard_True;
1952               // std::cout << "swap vertex " << std::endl;
1953               // break;
1954             }
1955           }
1956         }
1957       }
1958     }
1959   }
1960   return result;
1961 }
1962
1963 static void ToSmooth( const Handle(IntSurf_LineOn2S)& Line,
1964                       const Standard_Boolean          IsReversed,
1965                       const IntSurf_Quadric&          Quad,
1966                       const Standard_Boolean          IsFirst,
1967                       Standard_Real&                  D3D)
1968 {
1969   if(Line->NbPoints() <= 10)
1970     return;
1971
1972   D3D = 0.;
1973   Standard_Integer NbTestPnts = Line->NbPoints() / 5;
1974   if(NbTestPnts < 5) NbTestPnts = 5;
1975
1976   Standard_Integer startp = (IsFirst) ? 2 : (Line->NbPoints() - NbTestPnts - 2);
1977   Standard_Integer ip = 0;
1978   Standard_Real Uc = 0., Vc = 0., Un = 0., Vn = 0., DDU = 0., DDV = 0.;
1979
1980   for(ip = startp; ip <= NbTestPnts; ip++) {
1981     if(IsReversed) {
1982       Line->Value(ip).ParametersOnS2(Uc,Vc); // S2 - quadric
1983       Line->Value(ip+1).ParametersOnS2(Un,Vn);
1984     }
1985     else {
1986       Line->Value(ip).ParametersOnS1(Uc,Vc); // S1 - quadric
1987       Line->Value(ip+1).ParametersOnS1(Un,Vn);
1988     }
1989     DDU += fabs(fabs(Uc)-fabs(Un));
1990     DDV += fabs(fabs(Vc)-fabs(Vn));
1991
1992     if(ip > startp) {
1993       Standard_Real DP = Line->Value(ip).Value().Distance(Line->Value(ip-1).Value());
1994       D3D += DP;
1995     }
1996   }
1997
1998   DDU /= (Standard_Real) NbTestPnts + 1;
1999   DDV /= (Standard_Real) NbTestPnts + 1;
2000
2001   D3D /= (Standard_Real) NbTestPnts + 1;
2002
2003
2004   Standard_Integer Index1 = (IsFirst) ? 1 : (Line->NbPoints());
2005   Standard_Integer Index2 = (IsFirst) ? 2 : (Line->NbPoints()-1);
2006   Standard_Integer Index3 = (IsFirst) ? 3 : (Line->NbPoints()-2);
2007
2008   Standard_Boolean doU = Standard_False;
2009
2010   Standard_Real U1 = 0., U2 = 0., V1 = 0., V2 = 0., U3 = 0., V3 = 0.;
2011
2012   if(IsReversed) {
2013     Line->Value(Index1).ParametersOnS2(U1,V1); // S2 - quadric
2014     Line->Value(Index2).ParametersOnS2(U2,V2);
2015     Line->Value(Index3).ParametersOnS2(U3,V3);
2016   }
2017   else {
2018     Line->Value(Index1).ParametersOnS1(U1,V1); // S1 - quadric
2019     Line->Value(Index2).ParametersOnS1(U2,V2);
2020     Line->Value(Index3).ParametersOnS1(U3,V3);
2021   }
2022
2023   if(!doU && Quad.TypeQuadric() == GeomAbs_Sphere) {
2024     if(fabs(fabs(U1)-fabs(U2)) > (M_PI/16.)) doU = Standard_True;
2025
2026     if(doU && (fabs(U1) <= 1.e-9 || fabs(U1-2.*M_PI) <= 1.e-9)) {
2027       if(fabs(V1-M_PI/2.) <= 1.e-9 || fabs(V1+M_PI/2.) <= 1.e-9) {}
2028       else {
2029         doU = Standard_False;
2030       }
2031     }
2032   }
2033
2034   if(Quad.TypeQuadric() == GeomAbs_Cone) {
2035     Standard_Real Uapx = 0., Vapx = 0.;
2036     Quad.Parameters(Quad.Cone().Apex(),Uapx,Vapx);
2037
2038     if(fabs(fabs(U1)-fabs(U2)) > M_PI/32.) doU = Standard_True;
2039
2040     if(doU && (fabs(U1) <= 1.e-9 || fabs(U1-2.*M_PI) <= 1.e-9)) {
2041       if(fabs(V1-Vapx) <= 1.e-9) {}
2042       else {
2043         doU = Standard_False;
2044       }
2045     }
2046   }
2047
2048   if(doU) {
2049     Standard_Real dU = Min((DDU/10.),5.e-8);
2050     Standard_Real U = (U2 > U3) ? (U2 + dU) : (U2 - dU);
2051     if(IsReversed)
2052       Line->SetUV(Index1,Standard_False,U,V1);
2053     else
2054       Line->SetUV(Index1,Standard_True,U,V1);
2055     U1 = U;
2056   }
2057 }
2058
2059 static Standard_Boolean TestMiddleOnPrm(const IntSurf_PntOn2S& aP,
2060                                         const IntSurf_PntOn2S& aV,
2061                                         const Standard_Boolean IsReversed,
2062                                         const Standard_Real    ArcTol,
2063                                         const Handle(Adaptor3d_TopolTool)&  PDomain)
2064
2065 {
2066   Standard_Boolean result = Standard_False;
2067   Standard_Real Up = 0., Vp = 0., Uv = 0., Vv = 0.;
2068   if(IsReversed) {
2069     aP.ParametersOnS1(Up,Vp); //S1 - parametric
2070     aV.ParametersOnS1(Uv,Vv);
2071   }
2072   else {
2073     aP.ParametersOnS2(Up,Vp); // S2 - parametric
2074     aV.ParametersOnS2(Uv,Vv);
2075   }
2076   Standard_Real Um = (Up + Uv)*0.5, Vm = (Vp + Vv)*0.5;
2077   gp_Pnt2d a2DPntM(Um,Vm);
2078   TopAbs_State PosM = PDomain->Classify(a2DPntM,ArcTol);
2079   if(PosM == TopAbs_ON || PosM == TopAbs_IN )
2080     result = Standard_True;
2081   return result;
2082 }
2083
2084 static void VerifyVertices( const Handle(IntSurf_LineOn2S)&    Line,
2085                             const Standard_Boolean             IsReversed,
2086                             const Handle(IntSurf_LineOn2S)&    Vertices,
2087                             const Standard_Real                TOL2D,
2088                             const Standard_Real       ArcTol,
2089                             const Handle(Adaptor3d_TopolTool)& PDomain,
2090                             IntSurf_PntOn2S&          VrtF,
2091                             Standard_Boolean&         AddFirst,
2092                             IntSurf_PntOn2S&          VrtL,
2093                             Standard_Boolean&         AddLast)
2094 {
2095   Standard_Integer nbp = Line->NbPoints(), nbv = Vertices->NbPoints();
2096   Standard_Integer FIndexSame = 0, FIndexNear = 0, LIndexSame = 0, LIndexNear = 0;
2097   const IntSurf_PntOn2S& aPF = Line->Value(1);
2098   const IntSurf_PntOn2S& aPL = Line->Value(nbp);
2099   Standard_Real UF = 0., VF = 0., UL = 0., VL = 0.;
2100   if(IsReversed) {
2101     aPF.ParametersOnS2(UF,VF);
2102     aPL.ParametersOnS2(UL,VL);
2103   }
2104   else {
2105     aPF.ParametersOnS1(UF,VF);
2106     aPL.ParametersOnS1(UL,VL);
2107   }
2108   gp_Pnt2d a2DPF(UF,VF);
2109   gp_Pnt2d a2DPL(UL,VL);
2110   Standard_Real DistMinF = 1.e+100, DistMinL = 1.e+100;
2111   Standard_Integer FConjugated = 0, LConjugated = 0;
2112
2113   Standard_Integer iv = 0;
2114
2115   for(iv = 1; iv <= nbv; iv++) {
2116     Standard_Real Uv = 0., Vv = 0.;
2117     if(IsReversed) {
2118       Vertices->Value(iv).ParametersOnS2(Uv,Vv);
2119       Uv = AdjustU(Uv);
2120       Vertices->SetUV(iv,Standard_False,Uv,Vv);
2121     }
2122     else {
2123       Vertices->Value(iv).ParametersOnS1(Uv,Vv);
2124       Uv = AdjustU(Uv);
2125       Vertices->SetUV(iv,Standard_True,Uv,Vv);
2126     }
2127   }
2128
2129   for(iv = 1; iv <= nbv; iv++) {
2130     const IntSurf_PntOn2S& aV = Vertices->Value(iv);
2131     if(aPF.IsSame(aV, Precision::Confusion(), Precision::PConfusion())) {
2132       FIndexSame = iv;
2133       break;
2134     }
2135     else {
2136       Standard_Real Uv = 0., Vv = 0.;
2137       if(IsReversed)
2138         aV.ParametersOnS2(Uv,Vv);
2139       else
2140         aV.ParametersOnS1(Uv,Vv);
2141       gp_Pnt2d a2DV(Uv,Vv);
2142       Standard_Real Dist = a2DV.Distance(a2DPF);
2143       if(Dist < DistMinF) {
2144         DistMinF = Dist;
2145         FIndexNear = iv;
2146         if(FConjugated != 0)
2147           FConjugated = 0;
2148       }
2149       if(IsSeamParameter(Uv,TOL2D)) {
2150         Standard_Real Ucv = (fabs(Uv) < fabs(2.*M_PI-Uv)) ? (2.*M_PI) : 0.;
2151         gp_Pnt2d a2DCV(Ucv,Vv);
2152         Standard_Real CDist = a2DCV.Distance(a2DPF);
2153         if(CDist < DistMinF) {
2154           DistMinF = CDist;
2155           FConjugated = iv;
2156           FIndexNear = iv;
2157         }
2158       }
2159     }
2160   }
2161
2162   for(iv = 1; iv <= nbv; iv++) {
2163     const IntSurf_PntOn2S& aV = Vertices->Value(iv);
2164     if(aPL.IsSame(aV, Precision::Confusion(), Precision::PConfusion())) {
2165       LIndexSame = iv;
2166       break;
2167     }
2168     else {
2169       Standard_Real Uv = 0., Vv = 0.;
2170       if(IsReversed)
2171         aV.ParametersOnS2(Uv,Vv);
2172       else
2173         aV.ParametersOnS1(Uv,Vv);
2174       gp_Pnt2d a2DV(Uv,Vv);
2175       Standard_Real Dist = a2DV.Distance(a2DPL);
2176       if(Dist < DistMinL) {
2177         DistMinL = Dist;
2178         LIndexNear = iv;
2179         if(LConjugated != 0)
2180           LConjugated = 0;
2181       }
2182       if(IsSeamParameter(Uv,TOL2D)) {
2183         Standard_Real Ucv = (fabs(Uv) < fabs(2.*M_PI-Uv)) ? (2.*M_PI) : 0.;
2184         gp_Pnt2d a2DCV(Ucv,Vv);
2185         Standard_Real CDist = a2DCV.Distance(a2DPL);
2186         if(CDist < DistMinL) {
2187           DistMinL = CDist;
2188           LConjugated = iv;
2189           LIndexNear = iv;
2190         }
2191       }
2192     }
2193   }
2194
2195   AddFirst = Standard_False;
2196   AddLast  = Standard_False;
2197
2198   if(FIndexSame == 0) {
2199     if(FIndexNear != 0) {
2200       const IntSurf_PntOn2S& aV = Vertices->Value(FIndexNear);
2201       Standard_Real Uv = 0., Vv = 0.;
2202       if(IsReversed)
2203         aV.ParametersOnS2(Uv,Vv);
2204       else
2205         aV.ParametersOnS1(Uv,Vv);
2206       if(IsSeamParameter(Uv,TOL2D)) {
2207         Standard_Real Ucv = (fabs(Uv) < fabs(2.*M_PI-Uv)) ? (2.*M_PI) : 0.;
2208         Standard_Boolean test = TestMiddleOnPrm(aPF,aV,IsReversed,ArcTol,PDomain);
2209         if(test) {
2210           VrtF.SetValue(aV.Value());
2211           if(IsReversed) {
2212             Standard_Real U2 = 0., V2 = 0.;
2213             aV.ParametersOnS1(U2,V2); // S1 - prm
2214             VrtF.SetValue(Standard_True,U2,V2);
2215             if(FConjugated == 0)
2216               VrtF.SetValue(Standard_False,Uv,Vv);
2217             else
2218               VrtF.SetValue(Standard_False,Ucv,Vv);
2219           }
2220           else {
2221             Standard_Real U2 = 0., V2 = 0.;
2222             aV.ParametersOnS2(U2,V2); // S2 - prm
2223             VrtF.SetValue(Standard_False,U2,V2);
2224             if(FConjugated == 0)
2225               VrtF.SetValue(Standard_True,Uv,Vv);
2226             else
2227               VrtF.SetValue(Standard_True,Ucv,Vv);
2228           }
2229           Standard_Real Dist3D = VrtF.Value().Distance(aPF.Value());
2230           if(Dist3D > 1.5e-7 && DistMinF > TOL2D) {
2231             AddFirst = Standard_True;
2232           }
2233         }
2234       }
2235       else {
2236         // to do: analyze internal vertex
2237       }
2238     }
2239   }
2240
2241   if(LIndexSame == 0) {
2242     if(LIndexNear != 0) {
2243       const IntSurf_PntOn2S& aV = Vertices->Value(LIndexNear);
2244       Standard_Real Uv = 0., Vv = 0.;
2245       if(IsReversed)
2246         aV.ParametersOnS2(Uv,Vv);
2247       else
2248         aV.ParametersOnS1(Uv,Vv);
2249       if(IsSeamParameter(Uv,TOL2D)) {
2250         Standard_Real Ucv = (fabs(Uv) < fabs(2.*M_PI-Uv)) ? (2.*M_PI) : 0.;
2251         Standard_Boolean test = TestMiddleOnPrm(aPL,aV,IsReversed,ArcTol,PDomain);
2252         if(test) {
2253           VrtL.SetValue(aV.Value());
2254           if(IsReversed) {
2255             Standard_Real U2 = 0., V2 = 0.;
2256             aV.ParametersOnS1(U2,V2); // S1 - prm
2257             VrtL.SetValue(Standard_True,U2,V2);
2258             if(LConjugated == 0)
2259               VrtL.SetValue(Standard_False,Uv,Vv);
2260             else
2261               VrtL.SetValue(Standard_False,Ucv,Vv);
2262           }
2263           else {
2264             Standard_Real U2 = 0., V2 = 0.;
2265             aV.ParametersOnS2(U2,V2); // S2 - prm
2266             VrtL.SetValue(Standard_False,U2,V2);
2267             if(LConjugated == 0)
2268               VrtL.SetValue(Standard_True,Uv,Vv);
2269             else
2270               VrtL.SetValue(Standard_True,Ucv,Vv);
2271           }
2272           Standard_Real Dist3D = VrtL.Value().Distance(aPL.Value());
2273           if(Dist3D > 1.5e-7 && DistMinL > TOL2D) {
2274             AddLast = Standard_True;
2275           }
2276         }
2277       }
2278       else {
2279         // to do: analyze internal vertex
2280       }
2281     }
2282   }
2283 }
2284
2285 static Standard_Boolean AddVertices(Handle(IntSurf_LineOn2S)& Line,
2286   const IntSurf_PntOn2S&    VrtF,
2287   const Standard_Boolean    AddFirst,
2288   const IntSurf_PntOn2S&    VrtL,
2289   const Standard_Boolean    AddLast,
2290   const Standard_Real       D3DF,
2291   const Standard_Real       D3DL)
2292 {
2293   Standard_Boolean result = Standard_False;
2294   if(AddFirst) {
2295     Standard_Real DF = Line->Value(1).Value().Distance(VrtF.Value());
2296     if((D3DF*2.) > DF && DF > 1.5e-7) {
2297       Line->InsertBefore(1,VrtF);
2298       result = Standard_True;
2299     }
2300   }
2301   if(AddLast) {
2302     Standard_Real DL = Line->Value(Line->NbPoints()).Value().Distance(VrtL.Value());
2303     if((D3DL*2.) > DL && DL > 1.5e-7) {
2304       Line->Add(VrtL);
2305       result = Standard_True;
2306     }
2307   }
2308   return result;
2309 }
2310
2311
2312 static void PutIntVertices(const Handle(IntPatch_PointLine)&    Line,
2313   Handle(IntSurf_LineOn2S)& Result,
2314   Standard_Boolean          theIsReversed,
2315   Handle(IntSurf_LineOn2S)& Vertices,
2316   const Standard_Real       ArcTol)
2317 {
2318   Standard_Integer nbp = Result->NbPoints(), nbv = Vertices->NbPoints();
2319
2320   if(nbp < 3)
2321     return;
2322
2323   const Handle(IntPatch_RLine) aRLine = Handle(IntPatch_RLine)::DownCast(Line);
2324
2325   Standard_Integer ip = 0, iv = 0;
2326   gp_Pnt aPnt;
2327   IntPatch_Point thePnt;
2328   Standard_Real U1 = 0., V1 = 0., U2 = 0., V2 = 0.;
2329
2330   for(ip = 2; ip <= (nbp-1); ip++) {
2331     const IntSurf_PntOn2S& aP = Result->Value(ip);
2332     for(iv = 1; iv <= nbv; iv++) {
2333       const IntSurf_PntOn2S& aV = Vertices->Value(iv);
2334       if(aP.IsSame(aV, Precision::Confusion(), Precision::PConfusion())) {
2335         aPnt = Result->Value(ip).Value();
2336         Result->Value(ip).ParametersOnS1(U1,V1);
2337         Result->Value(ip).ParametersOnS2(U2,V2);
2338         thePnt.SetValue(aPnt,ArcTol,Standard_False);
2339         thePnt.SetParameters(U1,V1,U2,V2);
2340
2341         Standard_Real aParam = (Standard_Real)ip;
2342
2343         if(!aRLine.IsNull())
2344         {
2345           //In fact, aRLine is always on the parametric surface.
2346           //If (theIsReversed == TRUE) then (U1, V1) - point on
2347           //parametric surface, otherwise - point on quadric.
2348           const Handle(Adaptor2d_Curve2d)& anArc = aRLine->IsArcOnS1() ?
2349                                                     aRLine->ArcOnS1() :
2350                                                     aRLine->ArcOnS2();
2351
2352           const gp_Lin2d aLin(anArc->Line());
2353           gp_Pnt2d aPSurf;
2354
2355           if(theIsReversed)
2356           {
2357             aPSurf.SetCoord(U1, V1);
2358           }
2359           else
2360           {
2361             aPSurf.SetCoord(U2, V2);
2362           }
2363
2364           aParam = ElCLib::Parameter(aLin, aPSurf);
2365         }
2366
2367         thePnt.SetParameter(aParam);
2368         Line->AddVertex(thePnt);
2369       }
2370     }
2371   }
2372 }
2373
2374 static Standard_Boolean HasInternals(Handle(IntSurf_LineOn2S)& Line,
2375   Handle(IntSurf_LineOn2S)& Vertices)
2376 {
2377   Standard_Integer nbp = Line->NbPoints(), nbv = Vertices->NbPoints();
2378   Standard_Integer ip = 0, iv = 0;
2379   Standard_Boolean result = Standard_False;
2380
2381   if(nbp < 3)
2382     return result;
2383
2384   for(ip = 2; ip <= (nbp-1); ip++) {
2385     const IntSurf_PntOn2S& aP = Line->Value(ip);
2386     for(iv = 1; iv <= nbv; iv++) {
2387       const IntSurf_PntOn2S& aV = Vertices->Value(iv);
2388       if(aP.IsSame(aV, Precision::Confusion(), Precision::PConfusion())) {
2389         result = Standard_True;
2390         break;
2391       }
2392     }
2393     if(result)
2394       break;
2395   }
2396
2397   return result;
2398 }
2399 static Handle(IntPatch_WLine) MakeSplitWLine (Handle(IntPatch_WLine)&        WLine,
2400   Standard_Boolean         Tang,
2401   IntSurf_TypeTrans        Trans1,
2402   IntSurf_TypeTrans        Trans2,
2403   Standard_Real            ArcTol,
2404   Standard_Integer         ParFirst,
2405   Standard_Integer         ParLast)
2406 {
2407   Handle(IntSurf_LineOn2S) SLine = WLine->Curve();
2408   Handle(IntSurf_LineOn2S) sline = new IntSurf_LineOn2S();
2409
2410   Standard_Integer ip = 0;
2411   for(ip = ParFirst; ip <= ParLast; ip++)
2412     sline->Add(SLine->Value(ip));
2413
2414   Handle(IntPatch_WLine) wline = new IntPatch_WLine(sline,Tang,Trans1,Trans2);
2415   wline->SetCreatingWayInfo(IntPatch_WLine::IntPatch_WLImpPrm);
2416
2417   gp_Pnt aSPnt;
2418   IntPatch_Point TPntF,TPntL;
2419   Standard_Real uu1 = 0., vv1 = 0., uu2 = 0., vv2 = 0.;
2420
2421   aSPnt = sline->Value(1).Value();
2422   sline->Value(1).ParametersOnS1(uu1,vv1);
2423   sline->Value(1).ParametersOnS2(uu2,vv2);
2424   TPntF.SetValue(aSPnt,ArcTol,Standard_False);
2425   TPntF.SetParameters(uu1,vv1,uu2,vv2);
2426   TPntF.SetParameter(1.);
2427   wline->AddVertex(TPntF);
2428   wline->SetFirstPoint(1);
2429
2430   aSPnt =  sline->Value(sline->NbPoints()).Value();
2431   sline->Value(sline->NbPoints()).ParametersOnS1(uu1,vv1);
2432   sline->Value(sline->NbPoints()).ParametersOnS2(uu2,vv2);
2433   TPntL.SetValue(aSPnt,ArcTol,Standard_False);
2434   TPntL.SetParameters(uu1,vv1,uu2,vv2);
2435   TPntL.SetParameter((Standard_Real)sline->NbPoints());
2436   wline->AddVertex(TPntL);
2437   wline->SetLastPoint(wline->NbVertex());
2438
2439   return wline;
2440 }
2441
2442 static Standard_Boolean SplitOnSegments(Handle(IntPatch_WLine)&        WLine,
2443   Standard_Boolean         Tang,
2444   IntSurf_TypeTrans        Trans1,
2445   IntSurf_TypeTrans        Trans2,
2446   Standard_Real            ArcTol,
2447   IntPatch_SequenceOfLine& Segments)
2448 {
2449   Standard_Boolean result = Standard_False;
2450   Segments.Clear();
2451
2452   Standard_Integer nbv = WLine->NbVertex();
2453   if(nbv > 3) {
2454     Standard_Integer iv = 0;
2455     for(iv = 1; iv < nbv; iv++) {
2456       Standard_Integer firstPar =
2457                     (Standard_Integer) WLine->Vertex(iv).ParameterOnLine();
2458       Standard_Integer lastPar  =
2459                     (Standard_Integer) WLine->Vertex(iv+1).ParameterOnLine();
2460       if((lastPar - firstPar) <= 1)
2461         continue;
2462       else {
2463         Handle(IntPatch_WLine) splitwline = MakeSplitWLine(WLine,Tang,Trans1,Trans2,
2464                                                                 ArcTol,firstPar,lastPar);
2465         Segments.Append(splitwline);
2466         if(!result)
2467           result = Standard_True;
2468       }
2469     }
2470   }
2471   return result;
2472 }
2473
2474 //=======================================================================
2475 //function : IsPointOnBoundary
2476 //purpose  : Returns TRUE if point <theParam> matches <theBoundary +/- thePeriod>
2477 //            with given tolerance criterion.
2478 //            For not-periodic case, thePeriod must be equal to 0.0.
2479 //=======================================================================
2480 static Standard_Boolean IsPointOnBoundary(const Standard_Real theToler2D,
2481                                           const Standard_Real theBoundary,
2482                                           const Standard_Real thePeriod,
2483                                           const Standard_Real theParam)
2484 {
2485   Standard_Real aDelta = Abs(theParam - theBoundary);
2486   if (thePeriod != 0.0)
2487   {
2488     aDelta = fmod(aDelta, thePeriod);
2489
2490     // 0 <= aDelta < thePeriod
2491     return ((aDelta < theToler2D) || ((thePeriod - aDelta) < theToler2D));
2492   }
2493
2494   // Here, thePeriod == 0.0, aDelta > 0.0
2495
2496   return (aDelta < theToler2D);
2497 }
2498
2499 //=======================================================================
2500 //function : DetectOfBoundaryAchievement
2501 //purpose  : Can change values of theNewLine (by adding the computed point on boundary,
2502 //            which parameter will be adjusted) and theIsOnBoundary variables.
2503 //=======================================================================
2504 static void DetectOfBoundaryAchievement(const Handle(Adaptor3d_Surface)& theQSurf, // quadric
2505                                         const Standard_Boolean theIsReversed,
2506                                         const Handle(IntSurf_LineOn2S)& theSourceLine,
2507                                         const Standard_Integer thePointIndex,
2508                                         const Standard_Real theToler2D,
2509                                         Handle(IntSurf_LineOn2S)& theNewLine,
2510                                         Standard_Boolean& theIsOnBoundary)
2511 {
2512   const Standard_Real aUPeriod = theQSurf->IsUPeriodic() ? theQSurf->UPeriod() : 0.0,
2513                       aVPeriod = theQSurf->IsVPeriodic() ? theQSurf->VPeriod() : 0.0;
2514   const Standard_Real aUf = theQSurf->FirstUParameter(),
2515                       aUl = theQSurf->LastUParameter(),
2516                       aVf = theQSurf->FirstVParameter(),
2517                       aVl = theQSurf->LastVParameter();
2518
2519   const IntSurf_PntOn2S &aPPrev = theSourceLine->Value(thePointIndex - 1),
2520                         &aPCurr = theSourceLine->Value(thePointIndex);
2521   Standard_Real aUPrev, aVPrev, aUCurr, aVCurr;
2522   if (theIsReversed)
2523   {
2524     aPPrev.ParametersOnS2(aUPrev, aVPrev);   // S2 - quadric, set U,V by Pnt3D
2525     aPCurr.ParametersOnS2(aUCurr, aVCurr);   // S2 - quadric, set U,V by Pnt3D
2526   }
2527   else
2528   {
2529     aPPrev.ParametersOnS1(aUPrev, aVPrev);    // S1 - quadric, set U,V by Pnt3D
2530     aPCurr.ParametersOnS1(aUCurr, aVCurr);    // S1 - quadric, set U,V by Pnt3D
2531   }
2532
2533   // Ignore cases when the WLine goes along the surface boundary completely.
2534
2535   if (IsPointOnBoundary(theToler2D, aUf, aUPeriod, aUCurr) &&
2536       !IsPointOnBoundary(theToler2D, aUf, aUPeriod, aUPrev))
2537   {
2538     theIsOnBoundary = Standard_True;
2539   }
2540   else if (IsPointOnBoundary(theToler2D, aUl, aUPeriod, aUCurr) &&
2541            !IsPointOnBoundary(theToler2D, aUl, aUPeriod, aUPrev))
2542   {
2543     theIsOnBoundary = Standard_True;
2544   }
2545   else if (IsPointOnBoundary(theToler2D, aVf, aVPeriod, aVCurr) &&
2546            !IsPointOnBoundary(theToler2D, aVf, aVPeriod, aVPrev))
2547   {
2548     theIsOnBoundary = Standard_True;
2549   }
2550   else if (IsPointOnBoundary(theToler2D, aVl, aVPeriod, aVCurr) &&
2551            !IsPointOnBoundary(theToler2D, aVl, aVPeriod, aVPrev))
2552   {
2553     theIsOnBoundary = Standard_True;
2554   }
2555
2556   if (theIsOnBoundary)
2557   {
2558     // Adjust, to avoid bad jumping of the WLine.
2559
2560     const Standard_Real aDu = (aUPrev - aUCurr);
2561     const Standard_Real aDv = (aVPrev - aVCurr);
2562     if (aUPeriod > 0.0 && (2.0*Abs(aDu) > aUPeriod))
2563     {
2564       aUCurr += Sign(aUPeriod, aDu);
2565     }
2566
2567     if (aVPeriod > 0.0 && (2.0*Abs(aDv) > aVPeriod))
2568     {
2569       aVCurr += Sign(aVPeriod, aDv);
2570     }
2571
2572     IntSurf_PntOn2S aPoint = aPCurr;
2573     aPoint.SetValue(!theIsReversed, aUCurr, aVCurr);
2574     theNewLine->Add(aPoint);
2575   }
2576 }
2577 //=======================================================================
2578 //function : DecomposeResult
2579 //purpose  : Split <theLine> in the places where it passes through seam edge
2580 //            or singularity (apex of cone or pole of sphere).
2581 //            This passage is detected by jump of U-parameter
2582 //            from point to point.
2583 //=======================================================================
2584 static Standard_Boolean DecomposeResult(const Handle(IntPatch_PointLine)& theLine,
2585                                         const Standard_Boolean       IsReversed,
2586                                         const IntSurf_Quadric&       theQuad,
2587                                         const Handle(Adaptor3d_TopolTool)& thePDomain,
2588                                         const Handle(Adaptor3d_Surface)&  theQSurf, //quadric
2589                                         const Handle(Adaptor3d_Surface)&  thePSurf, //parametric
2590                                         const Standard_Real                theArcTol,
2591                                         const Standard_Real                theTolTang,
2592                                         IntPatch_SequenceOfLine&           theLines)
2593 {
2594   if(theLine->ArcType() == IntPatch_Restriction)
2595   {
2596     const Handle(IntPatch_RLine)& aRL = Handle(IntPatch_RLine)::DownCast(theLine);
2597     if(!aRL.IsNull())
2598     {
2599       const Handle(Adaptor2d_Curve2d)& anArc = aRL->IsArcOnS1() ?
2600                                         aRL->ArcOnS1() :
2601                                         aRL->ArcOnS2();
2602       if(anArc->GetType() != GeomAbs_Line)
2603       {
2604         //Restriction line must be isoline.
2605         //Other cases are not supported by
2606         //existing algorithms.
2607
2608         return Standard_False;
2609       }
2610     }
2611   }
2612
2613   const Standard_Real aDeltaUmax = M_PI_2;
2614   const Standard_Real aTOL3D = 1.e-10,
2615                       aTOL2D = Precision::PConfusion(),
2616                       aTOL2DS = Precision::PConfusion();
2617
2618   const Handle(IntSurf_LineOn2S)& aSLine = theLine->Curve();
2619
2620   if(aSLine->NbPoints() <= 2)
2621   {
2622     return Standard_False;
2623   }
2624
2625   //Deletes repeated vertices
2626   Handle(IntSurf_LineOn2S) aVLine = GetVertices(theLine,aTOL3D,aTOL2D);
2627
2628   Handle(IntSurf_LineOn2S) aSSLine(aSLine);
2629
2630   if(aSSLine->NbPoints() <= 1)
2631     return Standard_False;
2632
2633   AdjustLine(aSSLine,IsReversed,theQSurf,aTOL2D);
2634
2635   if(theLine->ArcType() == IntPatch_Walking)
2636   {
2637     Standard_Boolean isInserted = Standard_True;
2638     while(isInserted)
2639     {
2640       const Standard_Integer aNbPnts = aSSLine->NbPoints();
2641       TColStd_Array1OfInteger aPTypes(1,aNbPnts);
2642       SearchVertices(aSSLine,aVLine,aPTypes);
2643       isInserted = InsertSeamVertices(aSSLine,IsReversed,aVLine,aPTypes,aTOL2D);
2644     }
2645   }
2646
2647   const Standard_Integer aLindex = aSSLine->NbPoints();
2648   Standard_Integer aFindex = 1, aBindex = 0;
2649
2650   // build WLine parts (if any)
2651   Standard_Boolean flNextLine = Standard_True;
2652   Standard_Boolean hasBeenDecomposed = Standard_False;
2653   IntPatch_SpecPntType aPrePointExist = IntPatch_SPntNone;
2654
2655   IntSurf_PntOn2S PrePoint;
2656   while(flNextLine)
2657   {
2658     // reset variables
2659     flNextLine = Standard_False;
2660     Standard_Boolean isDecomposited = Standard_False;
2661
2662     Handle(IntSurf_LineOn2S) sline = new IntSurf_LineOn2S();
2663
2664     //if((Lindex-Findex+1) <= 2 )
2665     if((aLindex <= aFindex) && !aPrePointExist)
2666     {
2667       //break of "while(flNextLine)" cycle
2668       break;
2669     }
2670
2671     if(aPrePointExist)
2672     {
2673       const IntSurf_PntOn2S& aRefPt = aSSLine->Value(aFindex);
2674
2675       const Standard_Real aURes = theQSurf->UResolution(theArcTol),
2676                           aVRes = theQSurf->VResolution(theArcTol);
2677
2678       const Standard_Real aTol2d = (aPrePointExist == IntPatch_SPntPole) ? -1.0 :
2679               (aPrePointExist == IntPatch_SPntSeamV)? aVRes :
2680               (aPrePointExist == IntPatch_SPntSeamUV)? Max(aURes, aVRes) : aURes;
2681
2682       if(IntPatch_SpecialPoints::ContinueAfterSpecialPoint(theQSurf, thePSurf, aRefPt,
2683                                                               aPrePointExist, aTol2d,
2684                                                               PrePoint, IsReversed))
2685       {
2686         sline->Add(PrePoint);
2687
2688         //Avoid adding duplicate points.
2689         for (;aFindex <= aLindex; aFindex++)
2690         {
2691           if (!PrePoint.IsSame(aSSLine->Value(aFindex), theTolTang))
2692           {
2693             break;
2694           }
2695         }
2696       }
2697       else
2698       {
2699         //break of "while(flNextLine)" cycle
2700         break;
2701       }
2702     }
2703
2704     aPrePointExist = IntPatch_SPntNone;
2705
2706     // analyze other points
2707     for(Standard_Integer k = aFindex; k <= aLindex; k++)
2708     {
2709       if( k == aFindex )
2710       {
2711         PrePoint = aSSLine->Value(k);
2712         sline->Add(PrePoint);
2713         continue;
2714       }
2715
2716       //Check whether the current point is on the boundary of theQSurf.
2717       //If that is TRUE then the Walking-line will be decomposed in this point.
2718       //However, this boundary is not singular-point (like seam or pole of sphere).
2719       //Therefore, its processing will be simplified.
2720       Standard_Boolean isOnBoundary = Standard_False;
2721
2722       // Values of sline and isOnBoundary can be changed by this function
2723       DetectOfBoundaryAchievement(theQSurf, IsReversed, aSSLine,
2724                                   k, aTOL2D, sline, isOnBoundary);
2725
2726       aPrePointExist = IsSeamOrPole(theQSurf, aSSLine, IsReversed,
2727                                     k - 1, theTolTang, aDeltaUmax);
2728
2729       if (isOnBoundary && (aPrePointExist != IntPatch_SPntPoleSeamU))
2730       {
2731         // If the considered point is on seam then its UV-parameters
2732         // are defined to within the surface period. Therefore, we can
2733         // trust already computed parameters of this point.
2734         // However, if this point (which is on the surface boundary) is
2735         // a sphere pole or cone apex then its (point's) parameters
2736         // have to be recomputed in the code below
2737         // (see IntPatch_SpecialPoints::AddSingularPole() method).
2738         // E.g. see "bugs modalg_6 bug26684_2" test case.
2739
2740         aPrePointExist = IntPatch_SPntNone;
2741       }
2742
2743       if (aPrePointExist != IntPatch_SPntNone)
2744       {
2745         aBindex = k;
2746         isDecomposited = Standard_True;
2747         ////
2748         const IntSurf_PntOn2S& aRefPt = aSSLine->Value(aBindex-1);
2749
2750         Standard_Real aCompareTol3D = Precision::Confusion();
2751         Standard_Real aCompareTol2D = Precision::PConfusion();
2752
2753         IntSurf_PntOn2S aNewPoint = aRefPt;
2754         IntPatch_SpecPntType aLastType = IntPatch_SPntNone;
2755
2756         if(aPrePointExist == IntPatch_SPntSeamUV)
2757         {
2758           aPrePointExist = IntPatch_SPntNone;
2759           aLastType = IntPatch_SPntSeamUV;
2760           IntPatch_SpecialPoints::AddCrossUVIsoPoint(theQSurf, thePSurf,
2761                                                         aRefPt, theTolTang,
2762                                                         aNewPoint, IsReversed);
2763         }
2764         else if(aPrePointExist == IntPatch_SPntSeamV)
2765         {//WLine goes through seam
2766           aPrePointExist = IntPatch_SPntNone;
2767           aLastType = IntPatch_SPntSeamV;
2768
2769           //Not quadric point
2770           Standard_Real aU0 = 0.0, aV0 = 0.0;
2771           //Quadric point
2772           Standard_Real aUQuadRef = 0.0, aVQuadRef = 0.0;
2773
2774           if(IsReversed)
2775           {
2776             aRefPt.Parameters(aU0, aV0, aUQuadRef, aVQuadRef);
2777           }
2778           else
2779           {
2780             aRefPt.Parameters(aUQuadRef, aVQuadRef, aU0, aV0);
2781           }
2782
2783           math_Vector aTol(1, 3), aStartPoint(1,3),
2784             anInfBound(1, 3), aSupBound(1, 3);
2785
2786           //Parameters on parametric surface
2787           Standard_Real aUp = 0.0, aVp = 0.0, aUq = 0.0, aVq = 0.0;
2788           if(IsReversed)
2789           {
2790             aSSLine->Value(k).Parameters(aUp, aVp, aUq, aVq);
2791           }
2792           else
2793           {
2794             aSSLine->Value(k).Parameters(aUq, aVq, aUp, aVp);
2795           }
2796
2797           aTol(1) = thePSurf->UResolution(theArcTol);
2798           aTol(2) = thePSurf->VResolution(theArcTol);
2799           aTol(3) = theQSurf->UResolution(theArcTol);
2800           aStartPoint(1) = 0.5*(aU0 + aUp);
2801           aStartPoint(2) = 0.5*(aV0 + aVp);
2802           aStartPoint(3) = 0.5*(aUQuadRef + aUq);
2803           anInfBound(1) = thePSurf->FirstUParameter();
2804           anInfBound(2) = thePSurf->FirstVParameter();
2805           anInfBound(3) = theQSurf->FirstUParameter();
2806           aSupBound(1) = thePSurf->LastUParameter();
2807           aSupBound(2) = thePSurf->LastVParameter();
2808           aSupBound(3) = theQSurf->LastUParameter();
2809
2810           IntPatch_SpecialPoints::
2811                       AddPointOnUorVIso(theQSurf, thePSurf, aRefPt, Standard_False, 0.0,
2812                                         aTol, aStartPoint, anInfBound, aSupBound,
2813                                         aNewPoint, IsReversed);
2814         }
2815         else if(aPrePointExist == IntPatch_SPntPoleSeamU)
2816         {
2817           aPrePointExist = IntPatch_SPntNone;
2818
2819           IntPatch_Point aVert;
2820           aVert.SetValue(aRefPt);
2821           aVert.SetTolerance(theTolTang);
2822
2823           if(IntPatch_SpecialPoints::
2824                       AddSingularPole(theQSurf, thePSurf, aRefPt,
2825                                       aVert, aNewPoint, IsReversed))
2826           {
2827             aPrePointExist = IntPatch_SPntPole;
2828             aLastType = IntPatch_SPntPole;
2829             if (isOnBoundary)
2830             {
2831               // It is necessary to replace earlier added point on
2832               // the surface boundary with the pole. For that,
2833               // here we delete excess point. New point will be added later.
2834               isOnBoundary = Standard_False;
2835               sline->RemovePoint(sline->NbPoints());
2836             }
2837
2838             aCompareTol2D = -1.0;
2839           } //if(IntPatch_AddSpecialPoints::AddSingularPole(...))
2840           else
2841           {//Pole is not an intersection point
2842             aPrePointExist = IntPatch_SPntSeamU;
2843           }
2844         }
2845
2846         if(aPrePointExist == IntPatch_SPntSeamU)
2847         {//WLine goes through seam
2848           aPrePointExist = IntPatch_SPntNone;
2849           aLastType = IntPatch_SPntSeamU;
2850
2851           //Not quadric point
2852           Standard_Real aU0 = 0.0, aV0 = 0.0;
2853           //Quadric point
2854           Standard_Real aUQuadRef = 0.0, aVQuadRef = 0.0;
2855
2856           if(IsReversed)
2857           {
2858             aRefPt.Parameters(aU0, aV0, aUQuadRef, aVQuadRef);
2859           }
2860           else
2861           {
2862             aRefPt.Parameters(aUQuadRef, aVQuadRef, aU0, aV0);
2863           }
2864
2865           math_Vector aTol(1, 3), aStartPoint(1,3),
2866                       anInfBound(1, 3), aSupBound(1, 3);
2867
2868           //Parameters on parametric surface
2869           Standard_Real aUp = 0.0, aVp = 0.0, aUq = 0.0, aVq = 0.0;
2870           if (IsReversed)
2871           {
2872             aSSLine->Value(k).Parameters(aUp, aVp, aUq, aVq);
2873           }
2874           else
2875           {
2876             aSSLine->Value(k).Parameters(aUq, aVq, aUp, aVp);
2877           }
2878
2879           aTol(1) = thePSurf->UResolution(theArcTol);
2880           aTol(2) = thePSurf->VResolution(theArcTol);
2881           aTol(3) = theQSurf->VResolution(theArcTol);
2882           aStartPoint(1) = 0.5*(aU0 + aUp);
2883           aStartPoint(2) = 0.5*(aV0 + aVp);
2884           aStartPoint(3) = 0.5*(aVQuadRef + aVq);
2885           anInfBound(1) = thePSurf->FirstUParameter();
2886           anInfBound(2) = thePSurf->FirstVParameter();
2887           anInfBound(3) = theQSurf->FirstVParameter();
2888           aSupBound(1) = thePSurf->LastUParameter();
2889           aSupBound(2) = thePSurf->LastVParameter();
2890           aSupBound(3) = theQSurf->LastVParameter();
2891
2892           IntPatch_SpecialPoints::
2893                 AddPointOnUorVIso(theQSurf, thePSurf, aRefPt, Standard_True, 0.0, aTol,
2894                                   aStartPoint, anInfBound, aSupBound, aNewPoint,
2895                                   IsReversed);
2896         }
2897
2898         if(!aNewPoint.IsSame(aRefPt, aCompareTol3D, aCompareTol2D))
2899         {
2900           if (isOnBoundary)
2901             break;
2902
2903           sline->Add(aNewPoint);
2904           aPrePointExist = aLastType;
2905           PrePoint = aNewPoint;
2906         }
2907         else
2908         {
2909           if (isOnBoundary || (sline->NbPoints() == 1))
2910           {
2911             //FIRST point of the sline is the pole of the quadric.
2912             //Therefore, there is no point in decomposition.
2913
2914             // If the considered point is on surface boundary then
2915             // it is already marked as vertex. So, decomposition is
2916             // not required, too.
2917
2918             PrePoint = aRefPt;
2919             aPrePointExist = isOnBoundary ? IntPatch_SPntNone : aLastType;
2920           }
2921         }
2922
2923         ////
2924         break;
2925       } //if (aPrePointExist != IntPatch_SPntNone) cond.
2926
2927       PrePoint = aSSLine->Value(k);
2928
2929       if (isOnBoundary)
2930       {
2931         aBindex = k;
2932         isDecomposited = Standard_True;
2933         aPrePointExist = IntPatch_SPntNone;
2934         break;
2935       }
2936       else
2937       {
2938         sline->Add(aSSLine->Value(k));
2939       }
2940     } //for(Standard_Integer k = aFindex; k <= aLindex; k++)
2941
2942     //Creation of new line as part of existing theLine.
2943     //This part is defined by sline.
2944
2945     if(sline->NbPoints() == 1)
2946     {
2947       flNextLine = Standard_True;
2948
2949       if (aFindex < aBindex)
2950         aFindex = aBindex;
2951
2952       //Go to the next part of aSSLine
2953       //because we cannot create the line
2954       //with single point.
2955
2956       continue;
2957     }
2958
2959     IntSurf_PntOn2S aVF, aVL;
2960     Standard_Boolean addVF = Standard_False, addVL = Standard_False;
2961     VerifyVertices(sline,IsReversed,aVLine,aTOL2DS,theArcTol,
2962                                     thePDomain,aVF,addVF,aVL,addVL);
2963
2964     Standard_Boolean hasInternals = HasInternals(sline,aVLine);
2965
2966     Standard_Real D3F = 0., D3L = 0.;
2967     ToSmooth(sline,IsReversed,theQuad,Standard_True,D3F);
2968     ToSmooth(sline,IsReversed,theQuad,Standard_False,D3L);
2969
2970     //if(D3F <= 1.5e-7 && sline->NbPoints() >=3) {
2971     //  D3F = sline->Value(2).Value().Distance(sline->Value(3).Value());
2972     //}
2973     //if(D3L <= 1.5e-7 && sline->NbPoints() >=3) {
2974     //  D3L = sline->Value(sline->NbPoints()-1).Value().Distance(sline->
2975     //                                  Value(sline->NbPoints()-2).Value());
2976     //}
2977
2978     if(addVF || addVL)
2979     {
2980       Standard_Boolean isAdded = AddVertices(sline,aVF,addVF,aVL,addVL,D3F,D3L);
2981       if(isAdded)
2982       {
2983         ToSmooth(sline,IsReversed,theQuad,Standard_True,D3F);
2984         ToSmooth(sline,IsReversed,theQuad,Standard_False,D3L);
2985       }
2986     }
2987
2988     if(theLine->ArcType() == IntPatch_Walking)
2989     {
2990       IntPatch_Point aTPntF, aTPntL;
2991
2992       Handle(IntPatch_WLine) wline =
2993                           new IntPatch_WLine(sline,Standard_False,
2994                           theLine->TransitionOnS1(),theLine->TransitionOnS2());
2995       wline->SetCreatingWayInfo(IntPatch_WLine::IntPatch_WLImpPrm);
2996
2997       Standard_Real aU1 = 0.0, aV1 = 0.0, aU2 = 0.0, aV2 = 0.0;
2998       gp_Pnt aSPnt(sline->Value(1).Value());
2999       sline->Value(1).Parameters(aU1, aV1, aU2, aV2);
3000       aTPntF.SetValue(aSPnt,theArcTol,Standard_False);
3001       aTPntF.SetParameters(aU1, aV1, aU2, aV2);
3002       aTPntF.SetParameter(1.0);
3003       wline->AddVertex(aTPntF);
3004       wline->SetFirstPoint(1);
3005
3006       if(hasInternals)
3007       {
3008         PutIntVertices(wline,sline,IsReversed,aVLine,theArcTol);
3009       }
3010
3011       aSPnt =  sline->Value(sline->NbPoints()).Value();
3012       sline->Value(sline->NbPoints()).Parameters(aU1, aV1, aU2, aV2);
3013       aTPntL.SetValue(aSPnt,theArcTol,Standard_False);
3014       aTPntL.SetParameters(aU1, aV1, aU2, aV2);
3015       aTPntL.SetParameter(sline->NbPoints());
3016       wline->AddVertex(aTPntL);
3017       wline->SetLastPoint(wline->NbVertex());
3018
3019       IntPatch_SequenceOfLine segm;
3020       Standard_Boolean isSplited = SplitOnSegments(wline,Standard_False,
3021                       theLine->TransitionOnS1(),theLine->TransitionOnS2(),theArcTol,segm);
3022
3023       if(!isSplited)
3024       {
3025         theLines.Append(wline);
3026       }
3027       else
3028       {
3029         Standard_Integer nbsegms = segm.Length();
3030         Standard_Integer iseg = 0;
3031         for(iseg = 1; iseg <= nbsegms; iseg++)
3032           theLines.Append(segm(iseg));
3033       }
3034     }
3035     else
3036     {//theLine->ArcType() == IntPatch_Restriction
3037       if(!isDecomposited && !hasBeenDecomposed)
3038       {
3039         //The line has not been changed
3040         theLines.Append(Handle(IntPatch_RLine)::DownCast(theLine));
3041         return hasBeenDecomposed;
3042       }
3043
3044       IntPatch_Point aTPnt;
3045       gp_Pnt2d aPSurf;
3046       gp_Pnt aSPnt;
3047
3048       Handle(IntPatch_RLine) aRLine = new IntPatch_RLine(*Handle(IntPatch_RLine)::DownCast(theLine));
3049
3050       aRLine->ClearVertexes();
3051       aRLine->SetCurve(sline);
3052
3053       if(hasInternals)
3054       {
3055         PutIntVertices(aRLine,sline,IsReversed,aVLine,theArcTol);
3056       }
3057
3058       const Handle(Adaptor2d_Curve2d)& anArc = aRLine->IsArcOnS1() ?
3059                                                 aRLine->ArcOnS1() :
3060                                                 aRLine->ArcOnS2();
3061
3062       Standard_Real aFPar = anArc->FirstParameter(),
3063                     aLPar = anArc->LastParameter();
3064
3065       const IntSurf_PntOn2S &aRFirst = sline->Value(1),
3066                             &aRLast = sline->Value(sline->NbPoints());
3067
3068       const gp_Lin2d aLin(anArc->Line());
3069
3070       for(Standard_Integer aFLIndex = 0; aFLIndex < 2; aFLIndex++)
3071       {
3072         Standard_Real aU1 = 0.0, aV1 = 0.0, aU2 = 0.0, aV2 = 0.0;
3073         if(aFLIndex == 0)
3074         {
3075           aRFirst.Parameters(aU1, aV1, aU2, aV2);
3076           aSPnt.SetXYZ(aRFirst.Value().XYZ());
3077         }
3078         else
3079         {
3080           aRLast.Parameters(aU1, aV1, aU2, aV2);
3081           aSPnt.SetXYZ(aRLast.Value().XYZ());
3082         }
3083
3084         if(IsReversed)
3085         {
3086           aPSurf.SetCoord(aU1, aV1);
3087         }
3088         else
3089         {
3090           aPSurf.SetCoord(aU2, aV2);
3091         }
3092
3093         Standard_Real aPar = ElCLib::Parameter(aLin, aPSurf);
3094
3095         if(aFLIndex == 0)
3096         {
3097           aFPar = Max(aFPar, aPar);
3098           aPar = aFPar;
3099         }
3100         else
3101         {
3102           aLPar = Min(aLPar, aPar);
3103           aPar = aLPar;
3104         }
3105
3106         aTPnt.SetParameter(aPar);
3107         aTPnt.SetValue(aSPnt,theArcTol,Standard_False);
3108         aTPnt.SetParameters(aU1, aV1, aU2, aV2);
3109
3110         aRLine->AddVertex(aTPnt);
3111       }
3112
3113       if(aLPar - aFPar > Precision::PConfusion())
3114       {
3115         aRLine->SetFirstPoint(1);
3116         aRLine->SetLastPoint(aRLine->NbVertex());
3117
3118         anArc->Trim(aFPar, aLPar, theArcTol);
3119
3120         theLines.Append(aRLine);
3121       }
3122     }
3123
3124     if(isDecomposited)
3125     {
3126       aFindex = aBindex;
3127       flNextLine = hasBeenDecomposed = Standard_True;
3128     }
3129   }
3130
3131   return hasBeenDecomposed;
3132 }
3133
3134 //=======================================================================
3135 //function : CheckSegmSegm
3136 //purpose  : Returns TRUE if the segment [theParF, theParL] is included
3137 //            in the segment [theRefParF, theRefParL] segment.
3138 //=======================================================================
3139 static Standard_Boolean CheckSegmSegm(const Standard_Real theRefParF,
3140                                       const Standard_Real theRefParL,
3141                                       const Standard_Real theParF,
3142                                       const Standard_Real theParL)
3143 {
3144   if((theParF < theRefParF) || (theParF > theRefParL))
3145   {
3146     return Standard_False;
3147   }
3148
3149   if((theParL < theRefParF) || (theParL > theRefParL))
3150   {
3151     return Standard_False;
3152   }
3153
3154   return Standard_True;
3155 }
3156
3157 //=======================================================================
3158 //function : IsCoincide
3159 //purpose  : Check, if theLine is coincided with theArc (in 2d-space).
3160 //
3161 // Attention!!!
3162 //            Cases when theArc is not 2d-line adaptor are suppored by
3163 //          TopOpeBRep classes only (i.e. are archaic).
3164 //=======================================================================
3165 Standard_Boolean IsCoincide(IntPatch_TheSurfFunction& theFunc,
3166                             const Handle(IntPatch_PointLine)& theLine,
3167                             const Handle(Adaptor2d_Curve2d)& theArc,
3168                             const Standard_Boolean isTheSurface1Using, //Surf1 is parametric?
3169                             const Standard_Real theToler3D,
3170                             const Standard_Real theToler2D,
3171                             const Standard_Real thePeriod) // Period of parametric surface in direction which is perpendicular to theArc direction.
3172 {
3173   const Standard_Real aCoeffs[] = { 0.02447174185,  0.09549150281, 0.20610737385, 0.34549150281, /*Sin(x)*Sin(x)*/
3174                                     0.5, 0.65450849719, 0.79389262615 };
3175   if(theLine->ArcType() == IntPatch_Restriction)
3176   {//Restriction-restriction processing
3177     const Handle(IntPatch_RLine)& aRL2 = Handle(IntPatch_RLine)::DownCast(theLine);
3178     const Handle(Adaptor2d_Curve2d)& anArc = aRL2->IsArcOnS1() ? aRL2->ArcOnS1() : aRL2->ArcOnS2();
3179
3180     if(anArc->GetType() != GeomAbs_Line)
3181     {
3182       //Restriction line must be isoline.
3183       //Other cases are not supported by
3184       //existing algorithms.
3185
3186       return Standard_False;
3187     }
3188
3189     const gp_Lin2d aLin1(theArc->Line()),
3190                    aLin2(anArc->Line());
3191
3192     if(!aLin1.Direction().IsParallel(aLin2.Direction(), Precision::Angular()))
3193     {
3194       return Standard_False;
3195     }
3196
3197     const Standard_Real aDist =
3198             theArc->Line().Distance(anArc->Line());
3199     if((aDist < theToler2D) || (Abs(aDist - thePeriod) < theToler2D))
3200     {
3201       const Standard_Real aRf = theArc->FirstParameter(),
3202                           aRl = theArc->LastParameter();
3203       const Standard_Real aParf = anArc->FirstParameter(),
3204                           aParl = anArc->LastParameter();
3205       const gp_Pnt2d aP1(ElCLib::Value(aParf, aLin2)),
3206                      aP2(ElCLib::Value(aParl, aLin2));
3207
3208       Standard_Real aParam1 = ElCLib::Parameter(aLin1, aP1),
3209                     aParam2 = ElCLib::Parameter(aLin1, aP2);
3210
3211       if(CheckSegmSegm(aRf, aRl, aParam1, aParam2))
3212         return Standard_True;
3213
3214       //Lines are parallel. Therefore, there is no point in
3215       //projecting points to another line in order to check
3216       //if segment second line is included in segment of first one.
3217
3218       return CheckSegmSegm(aParam1, aParam2, aRf, aRl);
3219     }
3220
3221     return Standard_False;
3222   }
3223
3224   const Standard_Integer aNbPnts = theLine->NbPnts();
3225   const Standard_Real aUAf = theArc->FirstParameter(),
3226                       aUAl = theArc->LastParameter();
3227   const gp_Lin2d anArcLin(theArc->Line());
3228
3229   math_Vector aX(1, 2), aVal(1, 1);
3230
3231   for(Standard_Integer aPtID = 1; aPtID <= aNbPnts; aPtID++)
3232   {
3233     Standard_Real aUf = 0.0, aVf = 0.0;
3234     if(isTheSurface1Using)
3235       theLine->Point(aPtID).ParametersOnS1(aUf, aVf);
3236     else
3237       theLine->Point(aPtID).ParametersOnS2(aUf, aVf);
3238
3239     //Take 2d-point in parametric surface (because theArc is
3240     //2d-line in parametric surface).
3241     const gp_Pnt2d aPloc(aUf, aVf);
3242
3243     const Standard_Real aRParam = ElCLib::Parameter(anArcLin, aPloc);
3244
3245     if((aRParam < aUAf) || (aRParam > aUAl))
3246       return Standard_False;
3247
3248     const gp_Pnt2d aPmin(ElCLib::Value(aRParam, anArcLin));
3249
3250     const Standard_Real aDist = aPloc.Distance(aPmin);
3251     if((aDist < theToler2D) || (Abs(aDist - thePeriod) < theToler2D))
3252     {//Considered point is in Restriction line.
3253      //Go to the next point.
3254       continue;
3255     }
3256
3257     //Check if intermediate points between aPloc and theArc are
3258     //intersection point (i.e. if aPloc is in tangent zone between
3259     //two intersected surfaces).
3260
3261     const Standard_Real aUl = aPmin.X(), aVl = aPmin.Y();
3262
3263     Standard_Real aU, aV;
3264     Standard_Real dU = aUl - aUf, dV = aVl - aVf;
3265     for(Standard_Integer i = 0; i < 7; i++)
3266     {
3267       aU = aUf + aCoeffs[i] * dU;
3268       aV = aVf + aCoeffs[i] * dV;
3269
3270       aX.Value(1) = aU;
3271       aX.Value(2) = aV;
3272
3273       if(!theFunc.Value(aX, aVal))
3274       {
3275         return Standard_False;
3276       }
3277
3278       if(Abs(theFunc.Root()) > theToler3D)
3279       {
3280         return Standard_False;
3281       }
3282     }
3283   }
3284
3285   return Standard_True;
3286 }