src/IntPatch/IntPatch_ImpImpIntersection_0.gxx

   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 //  Modified by skv - Thu Jan 15 15:57:15 2004 OCC4455
  18
  19 #include <IntPatch_ThePathPointOfTheSOnBounds.hxx>
  20 #include <IntPatch_TheSegmentOfTheSOnBounds.hxx>
  21 #include <IntPatch_RLine.hxx>
  22 #include <IntSurf.hxx>
  23 #include <TColStd_Array1OfInteger.hxx>
  24 #include <TColStd_SequenceOfReal.hxx>
  25 #include <IntPatch_GLine.hxx>
  26 #include <Extrema_ExtPC.hxx>
  27 #include <GeomAdaptor_Curve.hxx>
  28 #include <Geom_Ellipse.hxx>
  29 #include <Geom_Parabola.hxx>
  30 #include <Geom_Hyperbola.hxx>
  31
  32
  33 static void PutPointsOnLine(const Handle(Adaptor3d_HSurface)& S1,
  34                             const Handle(Adaptor3d_HSurface)& S2,
  35                             const IntPatch_SequenceOfPathPointOfTheSOnBounds&,
  36                             const IntPatch_SequenceOfLine&,
  37                             const Standard_Boolean,
  38                             const Handle(Adaptor3d_TopolTool)&,
  39                             const IntSurf_Quadric&,
  40                             const IntSurf_Quadric&,
  41                             const Standard_Boolean,
  42                             const Standard_Real);
  43
  44 static Standard_Boolean MultiplePoint (const IntPatch_SequenceOfPathPointOfTheSOnBounds&,
  45                                        const Handle(Adaptor3d_TopolTool)&,
  46                                        const IntSurf_Quadric&,
  47                                        const gp_Vec&,
  48                                        const IntPatch_SequenceOfLine&,
  49                                        TColStd_Array1OfInteger&,
  50                                        TColStd_Array1OfInteger&,
  51                                        const Standard_Integer,
  52                                        const Standard_Boolean);
  53
  54 static Standard_Boolean PointOnSecondDom (const IntPatch_SequenceOfPathPointOfTheSOnBounds&,
  55                                           const Handle(Adaptor3d_TopolTool)&,
  56                                           const IntSurf_Quadric&,
  57                                           const gp_Vec&,
  58                                           const gp_Vec&,
  59                                           const Handle(IntPatch_Line)&,
  60                                           TColStd_Array1OfInteger&,
  61                                           const Standard_Integer);
  62
  63 static Standard_Boolean SingleLine (const gp_Pnt&,
  64                                     const Handle(IntPatch_Line)&,
  65                                     const Standard_Real,
  66                                     Standard_Real&,
  67                                     gp_Vec&);
  68
  69
  70 static Standard_Boolean FindLine (gp_Pnt&,
  71                                   const IntPatch_SequenceOfLine&,
  72                                   const Standard_Real,
  73                                   Standard_Real&,
  74                                   gp_Vec&,
  75                                   Standard_Integer&,
  76                                   Standard_Integer,
  77                                   const Handle(Adaptor2d_HCurve2d)&,
  78                                   Standard_Real&,
  79                                   gp_Pnt& pointonarc,
  80                                   const IntSurf_Quadric&);
  81
  82 static void ProcessSegments (const IntPatch_SequenceOfSegmentOfTheSOnBounds&,
  83                              IntPatch_SequenceOfLine&,
  84                              const IntSurf_Quadric&,
  85                              const IntSurf_Quadric&,
  86                              const Standard_Boolean,
  87                              const Standard_Real);
  88
  89 static void ProcessRLine (IntPatch_SequenceOfLine&,
  90                           const IntSurf_Quadric&,
  91                           const IntSurf_Quadric&,
  92                           const Standard_Real,
  93                           const Standard_Boolean theIsReqToKeepRLine);
  94
  95 //-- le calcul de dist est completement faux ds la routine ci dessous a revoir (lbr le 18 nov 97)
  96 Standard_Boolean IntersectionWithAnArc(gp_Pnt& PSurf,
  97                                        const Handle(IntPatch_ALine)& alin,
  98                                        Standard_Real& para,
  99                                        const Handle(Adaptor2d_HCurve2d)& thearc,
 100                                        Standard_Real& _theparameteronarc,
 101                                        gp_Pnt& thepointonarc,
 102                                        const IntSurf_Quadric& QuadSurf,
 103                                        const Standard_Real u0alin,
 104                                        const Standard_Real u1alin,
 105                                        Standard_Real& actualdist) {
 106   Standard_Real dtheta,theta;
 107 #ifdef OCCT_DEBUG
 108   //Standard_Real u,v,A,B,C,cost,sint,sign;
 109 #endif
 110   //-- recherche bete du point le plus proche de thearc->Value(...)
 111   dtheta = (u1alin-u0alin)*0.01;
 112   Standard_Real du=0.000000001;
 113   Standard_Real distmin = RealLast();
 114
 115   Standard_Real thetamin = 0.;
 116
 117   Standard_Real theparameteronarc = _theparameteronarc;
 118   for(Standard_Real _theta=u0alin+dtheta; _theta<=u1alin-dtheta; _theta+=dtheta) {
 119     gp_Pnt P=alin->Value(_theta);
 120     Standard_Real d=P.Distance(PSurf);
 121     if(d<distmin) {
 122       thetamin=_theta;
 123       distmin=d;
 124     }
 125   }
 126
 127   Standard_Real bestpara =0., besttheta =0., bestdist =0., distinit =0. ;
 128
 129   //-- Distance initiale
 130   {
 131     gp_Pnt pp0 = alin->Value(thetamin);
 132     Standard_Real ua0,va0;
 133     QuadSurf.Parameters(pp0,ua0,va0);
 134     gp_Pnt2d p2d;
 135     gp_Vec2d d2d;
 136     thearc->D1(theparameteronarc,p2d,d2d);
 137     gp_Vec2d PaPr(gp_Pnt2d(ua0,va0),p2d);
 138     distinit=PaPr.Magnitude();
 139   }
 140   theta = thetamin;
 141   //-- recherche a partir de theta et theparameteronarc
 142   Standard_Boolean cpasok=Standard_True;
 143   Standard_Integer nbiter=0;
 144   Standard_Real drmax = (thearc->LastParameter() - thearc->FirstParameter())*0.05;
 145   Standard_Real damax = (u1alin-u0alin)*0.05;
 146
 147
 148
 149   bestdist = RealLast();
 150
 151   do {
 152     Standard_Real ua0,va0,ua1,va1;
 153     //-- alin->Curve().InternalUVValue(theta,ua0,va0,A,B,C,cost,sint,sign);
 154     //-- alin->Curve().InternalUVValue(theta+du,ua1,va1,A,B,C,cost,sint,sign);
 155     gp_Pnt pp0 = alin->Value(theta);
 156     gp_Pnt pp1 = alin->Value(theta+du);
 157     QuadSurf.Parameters(pp0,ua0,va0);
 158     QuadSurf.Parameters(pp1,ua1,va1);
 159
 160
 161     gp_Vec2d D1a((ua1-ua0)/du,(va1-va0)/du);
 162     gp_Pnt2d p2d;
 163     gp_Vec2d d2d;
 164     thearc->D1(theparameteronarc,p2d,d2d);
 165     gp_Vec2d PaPr(gp_Pnt2d(ua0,va0),p2d);
 166
 167     Standard_Real pbd=PaPr.Magnitude();
 168     if(bestdist>pbd) {
 169       bestdist = pbd;
 170       bestpara = theparameteronarc;
 171       besttheta = theta;
 172     }
 173
 174     D1a.SetCoord(-D1a.X(),-D1a.Y());
 175
 176     Standard_Real d  = D1a.X()    *   d2d.Y()   -   D1a.Y()     * d2d.X();
 177
 178     Standard_Real da = (-PaPr.X())*   d2d.Y()   -   (-PaPr.Y()) * d2d.X();
 179     Standard_Real dr =  D1a.X()   * (-PaPr.Y()) -   D1a.Y()     * (-PaPr.X());
 180     if(Abs(d)>1e-15) {
 181       da/=d;
 182       dr/=d;
 183     }
 184     else {
 185       if(Abs(PaPr.X())>Abs(PaPr.Y())) {
 186         Standard_Real xx=PaPr.X();
 187         xx*=0.5;
 188         if(D1a.X()) {
 189           da = -xx/D1a.X();
 190         }
 191         if(d2d.X()) {
 192           dr = -xx/d2d.X();
 193         }
 194       }
 195       else {
 196         Standard_Real yy=PaPr.Y();
 197         yy*=0.5;
 198         if(D1a.Y()) {
 199           da = -yy/D1a.Y();
 200         }
 201         if(d2d.Y()) {
 202           dr = -yy/d2d.Y();
 203         }
 204       }
 205     }
 206 //--     Standard_Real da = -PaPr.Dot(D1a);
 207 //--     Standard_Real dr = -PaPr.Dot(d2d);
 208
 209     if(da<-damax) da=-damax;
 210     else if(da>damax) da=damax;
 211     if(dr<-drmax) dr=-drmax;
 212     else if(dr>drmax) dr=drmax;
 213
 214     if(Abs(da)<1e-10 && Abs(dr)<1e-10) {
 215       para = theta;
 216       PSurf = alin->Value(para);
 217       _theparameteronarc=theparameteronarc;
 218       thepointonarc = alin->Value(para);
 219       cpasok=Standard_False;
 220 //--      printf("\nt:%d",nbiter);
 221       actualdist = bestdist;
 222       return(Standard_True);
 223     }
 224     else {
 225       theta+=da;
 226       theparameteronarc+=dr;
 227       if(   theparameteronarc>thearc->LastParameter() ) {
 228         theparameteronarc = thearc->LastParameter();
 229       }
 230       if(  theparameteronarc<thearc->FirstParameter() ) {
 231         theparameteronarc = thearc->FirstParameter();
 232       }
 233       if( theta < u0alin) {
 234         theta = u0alin;
 235       }
 236       if( theta > u1alin-du) {
 237         theta = u1alin-du-du;
 238       }
 239     }
 240     nbiter++;
 241   }
 242   while(cpasok && nbiter<20);
 243   if(bestdist < distinit) {
 244     para = besttheta;
 245     PSurf = alin->Value(para);
 246     _theparameteronarc=bestpara;
 247     thepointonarc = alin->Value(para);
 248 //--     printf("\nT:%d",nbiter);
 249     actualdist=bestdist;
 250     return(Standard_True);
 251   }
 252 //--   printf("\nF:%d",nbiter);
 253   return(Standard_False);
 254 }
 255
 256
 257
 258
 259
 260 //-- ======================================================================
 261 static void Recadre(const Handle(Adaptor3d_HSurface)& myHS1,
 262                     const Handle(Adaptor3d_HSurface)& myHS2,
 263                     Standard_Real& u1,
 264                     Standard_Real& v1,
 265                     Standard_Real& u2,
 266                     Standard_Real& v2) {
 267   Standard_Real f,l,lmf,fpls2;
 268   GeomAbs_SurfaceType typs1 = myHS1->GetType();
 269   GeomAbs_SurfaceType typs2 = myHS2->GetType();
 270
 271   Standard_Boolean myHS1IsUPeriodic,myHS1IsVPeriodic;
 272   switch (typs1) {
 273   case GeomAbs_Cylinder:
 274   case GeomAbs_Cone:
 275   case GeomAbs_Sphere:
 276     {
 277       myHS1IsUPeriodic = Standard_True;
 278       myHS1IsVPeriodic = Standard_False;
 279       break;
 280     }
 281   case GeomAbs_Torus:
 282     {
 283       myHS1IsUPeriodic = myHS1IsVPeriodic = Standard_True;
 284       break;
 285     }
 286   default:
 287      {
 288        //-- Le cas de biparametrees periodiques est gere en amont
 289        myHS1IsUPeriodic = myHS1IsVPeriodic = Standard_False;
 290        break;
 291      }
 292   }
 293
 294   Standard_Boolean myHS2IsUPeriodic,myHS2IsVPeriodic;
 295   switch (typs2) {
 296   case GeomAbs_Cylinder:
 297   case GeomAbs_Cone:
 298   case GeomAbs_Sphere:
 299     {
 300       myHS2IsUPeriodic = Standard_True;
 301       myHS2IsVPeriodic = Standard_False;
 302       break;
 303     }
 304   case GeomAbs_Torus:
 305     {
 306       myHS2IsUPeriodic = myHS2IsVPeriodic = Standard_True;
 307       break;
 308     }
 309   default:
 310      {
 311        //-- Le cas de biparametrees periodiques est gere en amont
 312        myHS2IsUPeriodic = myHS2IsVPeriodic = Standard_False;
 313        break;
 314      }
 315   }
 316   if(myHS1IsUPeriodic) {
 317     lmf = M_PI+M_PI; //-- myHS1->UPeriod();
 318     f = myHS1->FirstUParameter();
 319     l = myHS1->LastUParameter();
 320     fpls2=0.5*(f+l);
 321     while((u1 < f)&&((fpls2-u1) > (u1+lmf-fpls2)   )) { u1+=lmf; }
 322     while((u1 > l)&&((u1-fpls2) > (fpls2-(u1-lmf)) )) { u1-=lmf; }
 323
 324   }
 325   if(myHS1IsVPeriodic) {
 326     lmf = M_PI+M_PI; //-- myHS1->VPeriod();
 327     f = myHS1->FirstVParameter();
 328     l = myHS1->LastVParameter();
 329     fpls2=0.5*(f+l);
 330     while((v1 < f)&&((fpls2-v1) > (v1+lmf-fpls2)   )) { v1+=lmf; }
 331     while((v1 > l)&&((v1-fpls2) > (fpls2-(v1-lmf)) )) { v1-=lmf; }
 332     //--    while(v1 < f) { v1+=lmf; }
 333     //--    while(v1 > l) { v1-=lmf; }
 334   }
 335   if(myHS2IsUPeriodic) {
 336     lmf = M_PI+M_PI; //-- myHS2->UPeriod();
 337     f = myHS2->FirstUParameter();
 338     l = myHS2->LastUParameter();
 339     fpls2=0.5*(f+l);
 340     while((u2 < f)&&((fpls2-u2) > (u2+lmf-fpls2)   )) { u2+=lmf; }
 341     while((u2 > l)&&((u2-fpls2) > (fpls2-(u2-lmf)) )) { u2-=lmf; }
 342     //-- while(u2 < f) { u2+=lmf; }
 343     //-- while(u2 > l) { u2-=lmf; }
 344   }
 345   if(myHS2IsVPeriodic) {
 346     lmf = M_PI+M_PI; //-- myHS2->VPeriod();
 347     f = myHS2->FirstVParameter();
 348     l = myHS2->LastVParameter();
 349     fpls2=0.5*(f+l);
 350     while((v2 < f)&&((fpls2-v2) > (v2+lmf-fpls2)   )) { v2+=lmf; }
 351     while((v2 > l)&&((v2-fpls2) > (fpls2-(v2-lmf)) )) { v2-=lmf; }
 352     //-- while(v2 < f) { v2+=lmf; }
 353     //-- while(v2 > l) { v2-=lmf; }
 354   }
 355 }
 356 //=======================================================================
 357 //function : PutPointsOnLine
 358 //purpose  :
 359 //=======================================================================
 360 void PutPointsOnLine(const Handle(Adaptor3d_HSurface)& S1,
 361                      const Handle(Adaptor3d_HSurface)& S2,
 362                      const IntPatch_SequenceOfPathPointOfTheSOnBounds& listpnt,
 363                      const IntPatch_SequenceOfLine& slin,
 364                      const Standard_Boolean OnFirst,
 365                      const Handle(Adaptor3d_TopolTool)& Domain,
 366                      const IntSurf_Quadric& QuadSurf,
 367                      const IntSurf_Quadric& OtherQuad,
 368                      const Standard_Boolean multpoint,
 369                      const Standard_Real Tolarc) {
 370
 371   // Traitement des point (de listpnt) de depart. On les replace sur
 372   // la ligne d intersection, en leur affectant la transition correcte sur
 373   // cette ligne.
 374   Standard_Integer nbpnt = listpnt.Length();
 375   Standard_Integer nblin=slin.Length();
 376
 377   if (!slin.Length() || !nbpnt) {
 378     return;
 379   }
 380   //
 381   Standard_Integer i,k;
 382   Standard_Integer linenumber;
 383   Standard_Real paraint = 0.,currentparameter,tolerance;
 384   Standard_Real U1,V1,U2,V2;
 385   Standard_Boolean goon;
 386
 387
 388   gp_Pnt Psurf, ptbid;
 389   gp_Vec Normale, Vtgint, Vtgrst;
 390
 391   gp_Vec d1u,d1v;
 392   gp_Pnt2d p2d;
 393   gp_Vec2d d2d;
 394
 395   IntSurf_Transition Transline,Transarc;
 396
 397   Handle(Adaptor2d_HCurve2d) currentarc;
 398   Handle(Adaptor3d_HVertex) vtx,vtxbis;
 399
 400   IntPatch_Point solpnt;
 401   IntPatch_ThePathPointOfTheSOnBounds currentpointonrst;
 402   IntPatch_IType TheType;
 403
 404   TColStd_Array1OfInteger UsedLine(1,nblin);
 405   TColStd_Array1OfInteger Done(1,nbpnt);
 406   for(i=1;i<=nbpnt;i++) Done(i) = 0; //-- Initialisation a la main
 407
 408   for (i=1; i<=nbpnt; i++) {
 409
 410     if (Done(i) != 1) {
 411       currentpointonrst = listpnt.Value(i);
 412       Psurf   = currentpointonrst.Value();   // Point dans l espace
 413       tolerance = currentpointonrst.Tolerance();
 414
 415       // On recherche d abord si on a correspondance avec un "point multiple"
 416       UsedLine.Init(0);
 417
 418       goon = Standard_True;
 419       if (multpoint) {
 420 #if 1
 421         Normale = QuadSurf.Normale(Psurf);     // Normale a la surface au point
 422         currentarc = currentpointonrst.Arc();
 423         currentparameter = currentpointonrst.Parameter();
 424         currentarc->D1(currentparameter,p2d,d2d);
 425         QuadSurf.D1(p2d.X(),p2d.Y(),ptbid,d1u,d1v);
 426         Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 427 #endif
 428         goon = MultiplePoint(listpnt,Domain,QuadSurf,Normale,slin,Done, UsedLine,
 429                              i,OnFirst);
 430       }
 431       if (goon) {
 432         Standard_Boolean linefound;
 433
 434         for(Standard_Integer indiceline = 1; indiceline <=slin.Length(); indiceline++) {
 435           if( UsedLine(indiceline) != 0 )
 436             continue;
 437           linenumber = indiceline;
 438
 439           //-- Attention , les points peuvent etre deplaces
 440           //-- il faut reprendre le point original
 441           currentpointonrst = listpnt.Value(i);
 442           currentarc = currentpointonrst.Arc();
 443           currentparameter = currentpointonrst.Parameter();
 444           Psurf   = currentpointonrst.Value();   // Point dans l espace
 445           tolerance = currentpointonrst.Tolerance();
 446           //--
 447
 448
 449           //  Modified by skv - Thu Jan 15 15:57:15 2004 OCC4455 Begin
 450           if (! currentpointonrst.IsNew()) {
 451             Handle(Adaptor3d_HVertex) aVtx    = currentpointonrst.Vertex();
 452             Standard_Real aVtxTol = aVtx->Resolution(currentarc);
 453             Standard_Real aTolAng = 0.01*tolerance;
 454
 455             tolerance = Max(tolerance, aVtxTol);
 456
 457             gp_Vec aNorm1  = QuadSurf.Normale(Psurf);
 458             gp_Vec aNorm2  = OtherQuad.Normale(Psurf);
 459             //
 460             if (aNorm1.Magnitude()>gp::Resolution() &&
 461                 aNorm2.Magnitude()>gp::Resolution()) {
 462
 463               if (aNorm1.IsParallel(aNorm2, aTolAng))
 464                 tolerance = Sqrt(tolerance);
 465             }//
 466           }
 467           //  Modified by skv - Thu Jan 15 15:57:15 2004 OCC4455 End
 468           gp_Pnt pointonarc;
 469           Vtgint.SetCoord(0,0,0);
 470           linefound = FindLine(Psurf,slin,tolerance,paraint,Vtgint,linenumber,indiceline,
 471                                currentarc,currentparameter,pointonarc,QuadSurf);
 472           if (linefound) {
 473
 474 #if 1
 475             Normale = QuadSurf.Normale(Psurf);     // Normale a la surface au point
 476             currentarc = currentpointonrst.Arc();
 477             //-- currentparameter = currentpointonrst.Parameter();
 478             currentarc->D1(currentparameter,p2d,d2d);
 479             QuadSurf.D1(p2d.X(),p2d.Y(),ptbid,d1u,d1v);
 480             Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 481 #endif
 482
 483
 484
 485             const Handle(IntPatch_Line)& lin = slin.Value(linenumber);
 486             TheType = lin->ArcType();
 487
 488             if (!OnFirst) { // on cherche la correspondance entre point sur domaine
 489               // de la premiere surface et point sur domaine de la
 490               // deuxieme surface
 491
 492               goon = PointOnSecondDom (listpnt, Domain, QuadSurf, Normale,
 493                                        Vtgint, lin, Done, i);
 494             }
 495
 496             if (goon) {
 497               //-- Modification du 4 avril 97    tolerance->Tolarc
 498               //-- on replace sur le vertex la tolerance d entree et
 499               //-- non la tolerance qui a servi au FindLine
 500               solpnt.SetValue(Psurf,Tolarc,Standard_False);
 501
 502               U1 = p2d.X(); V1 = p2d.Y();
 503               OtherQuad.Parameters(Psurf,U2,V2);
 504
 505               if (OnFirst) {
 506                 Recadre(S1,S2,U1,V1,U2,V2);
 507                 solpnt.SetParameters(U1,V1,U2,V2);
 508               }
 509               else {
 510                 Recadre(S1,S2,U2,V2,U1,V1);
 511                 solpnt.SetParameters(U2,V2,U1,V1);
 512               }
 513               solpnt.SetParameter(paraint);
 514
 515               if (! currentpointonrst.IsNew()) {
 516                 vtx = currentpointonrst.Vertex();
 517                 solpnt.SetVertex(OnFirst,vtx);
 518               }
 519               else {
 520                 //-- goon = Standard_False; ????
 521               }
 522
 523               if(Normale.SquareMagnitude()<1e-16) {
 524                 Transline.SetValue(Standard_True,IntSurf_Undecided);
 525                 Transarc.SetValue(Standard_True,IntSurf_Undecided);
 526               }
 527               else {
 528                 IntSurf::MakeTransition(Vtgint,Vtgrst,Normale,Transline,Transarc);
 529               }
 530               solpnt.SetArc(OnFirst,currentarc, currentparameter,
 531                             Transline,Transarc);
 532               if (TheType == IntPatch_Analytic) {
 533                 (*((Handle(IntPatch_ALine)*)&lin))->AddVertex(solpnt);
 534               }
 535               else {
 536                 (*((Handle(IntPatch_GLine)*)&lin))->AddVertex(solpnt);
 537               }
 538               Done(i) = 1;
 539
 540               if (goon) {
 541                 for (k=i+1; k<= nbpnt; k++) {
 542                   if (Done(k) != 1) {
 543                     currentpointonrst = listpnt.Value(k);
 544                     if (!currentpointonrst.IsNew()) {
 545                       vtxbis = currentpointonrst.Vertex();
 546                       if(vtx.IsNull()) {
 547                       }
 548                       else if (Domain->Identical(vtx, vtxbis)) {
 549                         solpnt.SetVertex(OnFirst,vtxbis);
 550                         currentarc = currentpointonrst.Arc();
 551                         currentparameter = currentpointonrst.Parameter();
 552
 553                         //                    currentarc->D1(currentparameter,ptbid,Vtgrst);
 554                         currentarc->D1(currentparameter,p2d,d2d);
 555                         Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 556                         if(Normale.SquareMagnitude()<1e-16) {
 557                           Transline.SetValue(Standard_True,IntSurf_Undecided);
 558                           Transarc.SetValue(Standard_True,IntSurf_Undecided);
 559                         }
 560                         else {
 561                           IntSurf::MakeTransition(Vtgint,Vtgrst,Normale,
 562                                                   Transline,Transarc);
 563                         }
 564                         solpnt.SetArc(OnFirst,currentarc,currentparameter,
 565                                       Transline,Transarc);
 566                         if (TheType == IntPatch_Analytic) {
 567                           (*((Handle(IntPatch_ALine)*)&lin))->AddVertex(solpnt);
 568                         }
 569                         else {
 570                           (*((Handle(IntPatch_GLine)*)&lin))->AddVertex(solpnt);
 571                         }
 572                         Done(k) = 1;
 573                       }
 574                     }
 575                   }
 576                 }
 577               }
 578             }
 579           }
 580           else {
 581             Done(i) = 1; // il faudra tester si IsNew ou pas
 582             // et traiter en consequence
 583           }
 584         }
 585       }
 586     }
 587   }
 588 }
 589
 590
 591 Standard_Boolean  MultiplePoint (const IntPatch_SequenceOfPathPointOfTheSOnBounds& listpnt,
 592                                  const Handle(Adaptor3d_TopolTool)& Domain,
 593                                  const IntSurf_Quadric& QuadSurf,
 594                                  const gp_Vec&    Normale,
 595                                  const IntPatch_SequenceOfLine& slin,
 596                                  TColStd_Array1OfInteger& Done,
 597                                  TColStd_Array1OfInteger& UsedLine,
 598                                  const Standard_Integer Index,
 599                                  const Standard_Boolean OnFirst) {
 600
 601 // Traitement des points "multiples".
 602
 603
 604   Standard_Integer k,ii,jj,nbvtx;
 605   Standard_Integer nblin = slin.Length();
 606   IntPatch_IType TheType;
 607
 608
 609   IntSurf_Transition Transline,Transarc;
 610
 611
 612   IntPatch_Point intpt;
 613   Handle(Adaptor2d_HCurve2d) currentarc;
 614   Handle(Adaptor3d_HVertex) vtx,vtxbis;
 615
 616   Standard_Integer nbpnt = listpnt.Length();
 617   IntPatch_ThePathPointOfTheSOnBounds currentpointonrst = listpnt.Value(Index);
 618   IntPatch_ThePathPointOfTheSOnBounds otherpt;
 619   gp_Pnt Point = currentpointonrst.Value();
 620   TColStd_Array1OfInteger localdone(1,nbpnt); localdone.Init(0);
 621   for (ii=1; ii<=nbpnt; ii++) {
 622     localdone(ii)=Done(ii);
 623   }
 624
 625   Standard_Real currentparameter;
 626   Standard_Real Paraint;
 627   gp_Vec Vtgint,Vtgrst;
 628   gp_Pnt ptbid;
 629
 630   gp_Vec d1u,d1v;
 631   gp_Pnt2d p2d;
 632   gp_Vec2d d2d;
 633
 634   Standard_Boolean goon;
 635
 636   Standard_Boolean Retvalue = Standard_True;
 637
 638   for (ii = 1; ii <= nblin; ii++) {
 639     const Handle(IntPatch_Line)& slinValueii = slin.Value(ii);
 640     TheType = slinValueii->ArcType();
 641     if (TheType == IntPatch_Analytic) {
 642       nbvtx = (*((Handle(IntPatch_ALine)*)&slinValueii))->NbVertex();
 643     }
 644     else {
 645       nbvtx = (*((Handle(IntPatch_GLine)*)&slinValueii))->NbVertex();
 646     }
 647     jj = 1;
 648     while (jj <= nbvtx) {
 649       if (TheType == IntPatch_Analytic) {
 650         intpt = (*((Handle(IntPatch_ALine)*)&slinValueii))->Vertex(jj);
 651       }
 652       else {
 653         intpt = (*((Handle(IntPatch_GLine)*)&slinValueii))->Vertex(jj);
 654       }
 655       if (intpt.IsMultiple() &&
 656           (( OnFirst && !intpt.IsOnDomS1()) ||
 657            (!OnFirst && !intpt.IsOnDomS2()))) {
 658         if (Point.Distance(intpt.Value()) <= intpt.Tolerance()) {
 659           Retvalue = Standard_False;
 660           Standard_Boolean foo = SingleLine(Point,slinValueii,
 661                                             intpt.Tolerance(),Paraint,Vtgint);
 662           if (!foo) {
 663             return Standard_False;   // ne doit pas se produire
 664           }
 665
 666           if (!currentpointonrst.IsNew()) {
 667             goon = Standard_True;
 668             vtx = currentpointonrst.Vertex();
 669             intpt.SetVertex(OnFirst,vtx);
 670           }
 671           else {
 672             goon = Standard_False;
 673           }
 674           currentarc = currentpointonrst.Arc();
 675           currentparameter = currentpointonrst.Parameter();
 676           currentarc->D1(currentparameter,p2d,d2d);
 677           QuadSurf.D1(p2d.X(),p2d.Y(),ptbid,d1u,d1v);
 678           Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 679
 680           //-- Si la normale est nulle (apex d un cone) On simule une transition UNKNOWN
 681           if(Normale.SquareMagnitude()<1e-16) {
 682             Transline.SetValue(Standard_True,IntSurf_Undecided);
 683             Transarc.SetValue(Standard_True,IntSurf_Undecided);
 684           }
 685           else {
 686             IntSurf::MakeTransition(Vtgint,Vtgrst,Normale,Transline,Transarc);
 687           }
 688
 689           //-- Avant, on ne mettait pas ce point (17 nov 97)
 690           //--printf("\n ImpImp_0  : Point(%g,%g,%g)  intpt(%g,%g,%g) \n",
 691           //--  Point.X(),Point.Y(),Point.Z(),intpt.Value().X(),intpt.Value().Y(),intpt.Value().Z());
 692           intpt.SetValue(Point);
 693
 694           intpt.SetArc(OnFirst,currentarc,currentparameter,
 695                        Transline,Transarc);
 696
 697
 698           if (TheType == IntPatch_Analytic) {
 699             (*((Handle(IntPatch_ALine)*)&slinValueii))->Replace(jj,intpt);
 700           }
 701           else {
 702             (*((Handle(IntPatch_GLine)*)&slinValueii))->Replace(jj,intpt);
 703           }
 704           localdone(Index) = 1;
 705           if (goon) {
 706             for (k=Index+1; k<= nbpnt; k++) {
 707               if (Done(k) != 1) {
 708                 otherpt= listpnt.Value(k);
 709                 if (!otherpt.IsNew()) {
 710                   vtxbis = otherpt.Vertex();
 711                   if (Domain->Identical(vtx, vtxbis)) {
 712                     intpt.SetVertex(OnFirst,vtxbis);
 713                     currentarc = otherpt.Arc();
 714                     currentparameter = otherpt.Parameter();
 715
 716                     currentarc->D1(currentparameter,p2d,d2d);
 717                     Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 718                     if(Normale.SquareMagnitude()<1e-16) {
 719                       Transline.SetValue(Standard_True,IntSurf_Undecided);
 720                       Transarc.SetValue(Standard_True,IntSurf_Undecided);
 721                     }
 722                     else {
 723                       IntSurf::MakeTransition(Vtgint,Vtgrst,Normale,
 724                                               Transline,Transarc);
 725                     }
 726                     intpt.SetArc(OnFirst,currentarc,currentparameter,
 727                                  Transline,Transarc);
 728                     if (TheType == IntPatch_Analytic) {
 729                       (*((Handle(IntPatch_ALine)*)&slinValueii))->AddVertex(intpt);
 730                     }
 731                     else {
 732                       (*((Handle(IntPatch_GLine)*)&slinValueii))->AddVertex(intpt);
 733                     }
 734                     UsedLine(ii) = 1;
 735                     Retvalue = Standard_True;
 736                     localdone(k) = 1;
 737                   }
 738                 }
 739               }
 740             }
 741           }
 742 //--           jj = nbvtx +1;
 743         }
 744 //--         else {
 745           jj = jj+1;
 746 //--        }
 747       }
 748       else {
 749         jj = jj+1;
 750       }
 751     }
 752   }
 753
 754   for (ii=1; ii<=nbpnt;ii++) {
 755     Done(ii) = localdone(ii);
 756   }
 757
 758   return Retvalue;
 759 }
 760
 761
 762
 763 Standard_Boolean PointOnSecondDom (const IntPatch_SequenceOfPathPointOfTheSOnBounds& listpnt,
 764                                    const Handle(Adaptor3d_TopolTool)& Domain,
 765                                    const IntSurf_Quadric& QuadSurf,
 766                                    const gp_Vec&    Normale,
 767                                    const gp_Vec&    Vtgint,
 768                                    const Handle(IntPatch_Line)& lin,
 769                                    TColStd_Array1OfInteger& Done,
 770                                    const Standard_Integer Index)
 771
 772
 773 // Duplication des points sur domaine de l autre surface.
 774 // On sait que le vertex sous-jacent est PntRef
 775
 776
 777 {
 778
 779   Standard_Integer k,jj,nbvtx;
 780   IntPatch_IType TheType;
 781
 782   IntSurf_Transition Transline,Transarc;
 783   IntPatch_Point intpt;
 784   Handle(Adaptor2d_HCurve2d) currentarc;
 785   Handle(Adaptor3d_HVertex) vtx,vtxbis;
 786   gp_Pnt ptbid;
 787   gp_Vec Vtgrst;
 788
 789   gp_Vec d1u,d1v;
 790   gp_Pnt2d p2d;
 791   gp_Vec2d d2d;
 792
 793   Standard_Integer nbpnt = listpnt.Length();
 794   IntPatch_ThePathPointOfTheSOnBounds currentpointonrst = listpnt.Value(Index);
 795   Standard_Real currentparameter;
 796
 797   Standard_Boolean goon;
 798   Standard_Boolean Retvalue = Standard_True;
 799
 800   TheType = lin->ArcType();
 801   if (TheType == IntPatch_Analytic) {
 802     nbvtx = (*((Handle(IntPatch_ALine)*)&lin))->NbVertex();
 803   }
 804   else {
 805     nbvtx = (*((Handle(IntPatch_GLine)*)&lin))->NbVertex();
 806   }
 807   jj = 1;
 808   while (jj <= nbvtx) {
 809     if (TheType == IntPatch_Analytic) {
 810       intpt = (*((Handle(IntPatch_ALine)*)&lin))->Vertex(jj);
 811     }
 812     else {
 813       intpt = (*((Handle(IntPatch_GLine)*)&lin))->Vertex(jj);
 814     }
 815     if (!intpt.IsOnDomS2()) {
 816       if (currentpointonrst.Value().Distance(intpt.Value()) <=
 817           intpt.Tolerance()) {
 818         Retvalue = Standard_False;
 819         if (!currentpointonrst.IsNew()) {
 820           goon = Standard_True;
 821           vtx = currentpointonrst.Vertex();
 822           intpt.SetVertex(Standard_False,vtx);
 823         }
 824         else {
 825           goon = Standard_False;
 826         }
 827         currentarc = currentpointonrst.Arc();
 828         currentparameter = currentpointonrst.Parameter();
 829         currentarc->D1(currentparameter,p2d,d2d);
 830         QuadSurf.D1(p2d.X(),p2d.Y(),ptbid,d1u,d1v);
 831         Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 832         if(Normale.SquareMagnitude()<1e-16) {
 833           Transline.SetValue(Standard_True,IntSurf_Undecided);
 834           Transarc.SetValue(Standard_True,IntSurf_Undecided);
 835         }
 836         else {
 837           IntSurf::MakeTransition(Vtgint,Vtgrst,Normale,Transline,Transarc);
 838         }
 839         intpt.SetArc(Standard_False,currentarc,currentparameter,
 840                      Transline,Transarc);
 841         if (TheType == IntPatch_Analytic) {
 842           (*((Handle(IntPatch_ALine)*)&lin))->Replace(jj,intpt);
 843         }
 844         else {
 845           (*((Handle(IntPatch_GLine)*)&lin))->Replace(jj,intpt);
 846         }
 847         Done(Index) = 1;
 848
 849         if (goon) {
 850           for (k=Index+1; k<= nbpnt; k++) {
 851             if (Done(k) != 1) {
 852               currentpointonrst = listpnt.Value(k);
 853               if (!currentpointonrst.IsNew()) {
 854                 vtxbis = currentpointonrst.Vertex();
 855                   if (Domain->Identical(vtx, vtxbis)) {
 856                   intpt.SetVertex(Standard_False,vtxbis);
 857                   currentarc = currentpointonrst.Arc();
 858                   currentparameter = currentpointonrst.Parameter();
 859                   currentarc->D1(currentparameter,p2d,d2d);
 860                   Vtgrst.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
 861                   if(Normale.SquareMagnitude()<1e-16) {
 862                     Transline.SetValue(Standard_True,IntSurf_Undecided);
 863                     Transarc.SetValue(Standard_True,IntSurf_Undecided);
 864                   }
 865                   else {
 866                     IntSurf::MakeTransition(Vtgint,Vtgrst,Normale,
 867                                             Transline,Transarc);
 868                   }
 869                   intpt.SetArc(Standard_False,currentarc,currentparameter,
 870                                Transline,Transarc);
 871                   if (TheType == IntPatch_Analytic) {
 872                     (*((Handle(IntPatch_ALine)*)&lin))->AddVertex(intpt);
 873                   }
 874                   else {
 875                     (*((Handle(IntPatch_GLine)*)&lin))->AddVertex(intpt);
 876                   }
 877                   Done(k) = 1;
 878                 }
 879               }
 880             }
 881           }
 882         }
 883         //-- jj = nbvtx + 1;
 884         jj++;
 885       }
 886       else {
 887         jj = jj+1;
 888       }
 889     }
 890     else {
 891       jj = jj+1;
 892     }
 893     if (TheType == IntPatch_Analytic) {
 894       nbvtx = (*((Handle(IntPatch_ALine)*)&lin))->NbVertex();
 895     }
 896     else {
 897       nbvtx = (*((Handle(IntPatch_GLine)*)&lin))->NbVertex();
 898     }
 899   }
 900   return Retvalue;
 901 }
 902
 903
 904
 905 Standard_Boolean FindLine (gp_Pnt& Psurf,
 906                            const IntPatch_SequenceOfLine& slin,
 907                            const Standard_Real Tol,
 908                            Standard_Real& Paraint,
 909                            gp_Vec& Vtgtint,
 910                            Standard_Integer& Range,
 911                            Standard_Integer OnlyThisLine,
 912                            const Handle(Adaptor2d_HCurve2d)& thearc,
 913                            Standard_Real& theparameteronarc,
 914                            gp_Pnt& thepointonarc,
 915                            const IntSurf_Quadric& QuadSurf)
 916 {
 917
 918 // Traitement du point de depart ayant pour representation Psurf
 919 // dans l espace. On recherche la ligne d intersection contenant ce point.
 920 // On a en sortie la ligne, et le parametre et sa tangente du point sur
 921 // la ligne d intersection.
 922
 923   Standard_Real distmin = RealLast();
 924   Standard_Real dist,para;
 925   Standard_Real lower,upper;
 926   gp_Pnt pt;
 927   Standard_Integer i;
 928   IntPatch_IType typarc;
 929
 930   Standard_Integer nblin = slin.Length();
 931   for (i=1; i<=nblin; i++) {
 932     if(OnlyThisLine) { i=OnlyThisLine; nblin=0; }
 933     const Handle(IntPatch_Line)& lin = slin.Value(i);
 934     typarc = lin->ArcType();
 935     if (typarc == IntPatch_Analytic) {
 936       Standard_Boolean foo;
 937       lower = (*((Handle(IntPatch_ALine)*)&lin))->FirstParameter(foo);
 938       upper = (*((Handle(IntPatch_ALine)*)&lin))->LastParameter(foo);
 939     }
 940     else {
 941       if ((*((Handle(IntPatch_GLine)*)&lin))->HasFirstPoint()) {
 942         lower = (*((Handle(IntPatch_GLine)*)&lin))->FirstPoint().ParameterOnLine();
 943       }
 944       else {
 945         lower = RealFirst();
 946       }
 947       if ((*((Handle(IntPatch_GLine)*)&lin))->HasLastPoint()) {
 948         upper = (*((Handle(IntPatch_GLine)*)&lin))->LastPoint().ParameterOnLine();
 949       }
 950       else {
 951         upper = RealLast();
 952       }
 953     }
 954
 955     switch (typarc) {
 956     case IntPatch_Lin :
 957       {
 958         para = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Line(),Psurf);
 959         if (para <= upper && para >= lower) {
 960           pt = ElCLib::Value(para,(*((Handle(IntPatch_GLine)*)&lin))->Line());
 961           dist = Psurf.Distance(pt);
 962           if (dist< distmin) {
 963             distmin = dist;
 964             Paraint = para;
 965             Range = i;
 966           }
 967         }
 968       }
 969       break;
 970     case IntPatch_Circle :
 971       {
 972         para = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Circle(),Psurf);
 973         if ((para <= upper && para >= lower) ||
 974             (para + 2.*M_PI <=upper && para + 2.*M_PI >= lower) ||
 975             (para - 2.*M_PI <=upper && para - 2.*M_PI >= lower)) {
 976           pt = ElCLib::Value(para,(*((Handle(IntPatch_GLine)*)&lin))->Circle());
 977           dist = Psurf.Distance(pt);
 978           if (dist< distmin) {
 979             distmin = dist;
 980             Paraint = para;
 981             Range = i;
 982           }
 983         }
 984       }
 985       break;
 986     case IntPatch_Ellipse :
 987       {
 988         para = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Ellipse(),Psurf);
 989         if ((para <= upper && para >= lower) ||
 990             (para + 2.*M_PI <=upper && para + 2.*M_PI >= lower) ||
 991             (para - 2.*M_PI <=upper && para - 2.*M_PI >= lower)) {
 992           pt = ElCLib::Value(para,(*((Handle(IntPatch_GLine)*)&lin))->Ellipse());
 993           dist = Psurf.Distance(pt);
 994           if (dist< distmin) {
 995             distmin = dist;
 996             Paraint = para;
 997             Range = i;
 998           }
 999         }
1000       }
1001       break;
1002     case IntPatch_Parabola :
1003       {
1004
1005 #if 0
1006         para = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Parabola(),Psurf);
1007         if (para <= upper && para >= lower) {
1008           pt = ElCLib::Value(para,(*((Handle(IntPatch_GLine)*)&lin))->Parabola());
1009           dist = Psurf.Distance(pt);
1010           if (dist< distmin) {
1011             distmin = dist;
1012             Paraint = para;
1013             Range = i;
1014           }
1015         }
1016 #else
1017         //-- Le calcul du parametre sur une parabole est mal fait ds ElCLib. Il ne tient pas compte
1018         //-- de la meilleure facon de calculer (axe X ou axe Y). Bilan : Si la parabole est tres
1019         //-- pointue (focal de l'ordre de 1e-2 et si le point est a un parametre grand, ca foire. )
1020         //-- On ne peut pas modifier faciolement ds ElCLib car on ne passe pas la focale. ...
1021         const gp_Parab& Parab=(*((Handle(IntPatch_GLine)*)&lin))->Parabola();
1022         para = ElCLib::Parameter(Parab,Psurf);
1023         if (para <= upper && para >= lower) {
1024           Standard_Integer amelioration=0;
1025           //-- cout<<"\n ****** \n";
1026           do {
1027             Standard_Real parabis = para+0.0000001;
1028
1029             pt = ElCLib::Value(para,Parab);
1030             dist = Psurf.Distance(pt);
1031
1032             gp_Pnt ptbis = ElCLib::Value(parabis,Parab);
1033             Standard_Real distbis = Psurf.Distance(ptbis);
1034
1035             Standard_Real ddist = distbis-dist;
1036
1037             //--cout<<" para: "<<para<<"    dist:"<<dist<<"   ddist:"<<ddist<<endl;
1038
1039             if (dist< distmin) {
1040               distmin = dist;
1041               Paraint = para;
1042               Range = i;
1043             }
1044             if(dist<1.0e-9 && dist>-1.0e-9) { amelioration=100; }
1045
1046             if(ddist>1.0e-9 || ddist<-1.0e-9 ) {
1047               para=para-dist*(parabis-para)/ddist;
1048             }
1049             else {
1050               amelioration=100;
1051             }
1052           }
1053           while(++amelioration < 5);
1054         }
1055
1056
1057 #endif
1058       }
1059       break;
1060     case IntPatch_Hyperbola :
1061       {
1062         para = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Hyperbola(),Psurf);
1063         if (para <= upper && para >= lower) {
1064           pt = ElCLib::Value(para,(*((Handle(IntPatch_GLine)*)&lin))->Hyperbola());
1065           dist = Psurf.Distance(pt);
1066           if (dist< distmin) {
1067             distmin = dist;
1068             Paraint = para;
1069             Range = i;
1070           }
1071         }
1072       }
1073       break;
1074
1075     case IntPatch_Analytic :
1076       {
1077         const Handle(IntPatch_ALine)& alin = (*((Handle(IntPatch_ALine)*)&lin));
1078         Standard_Boolean ok = alin->FindParameter(Psurf,para);
1079         if (ok) {
1080           pt = alin->Value(para);
1081           dist = Psurf.Distance(pt);
1082           if (dist< distmin) {
1083             distmin = dist;
1084             Paraint = para;
1085             Range = i;
1086           }
1087         }
1088         else {
1089           //-- le point n a pas ete trouve par bete projection.
1090           //-- on essaie l intersection avec la restriction en 2d
1091           Standard_Real theparamonarc = theparameteronarc;
1092 //#ifdef OCCT_DEBUG
1093 //        Standard_Real anpara=para;
1094 //#endif
1095           gp_Pnt CopiePsurf=Psurf;
1096           Standard_Boolean ok=IntersectionWithAnArc(CopiePsurf,alin,para,thearc,theparamonarc,thepointonarc,QuadSurf,lower,upper,dist);
1097
1098           //--printf("\nIntersectionWithAnArc   %d \n Psurf(%g,%g,%g)->(%g,%g,%g)  dist=%g\n para(%g)->(%g)\n paraonarc(%g)->(%g)",
1099           //--   ok,Psurf.X(),Psurf.Y(),Psurf.Z(),thepointonarc.X(),thepointonarc.Y(),thepointonarc.Z(),dist,
1100           //--    anpara,para,theparameteronarc,theparamonarc);
1101           dist = CopiePsurf.Distance(Psurf);
1102           if(ok) {
1103             if(dist<Tol) {
1104               theparameteronarc = theparamonarc;
1105               Psurf = thepointonarc;
1106               distmin = dist;
1107               Paraint = para;
1108               Range = i;
1109             }
1110           }
1111         }
1112       }
1113       break;
1114
1115     case IntPatch_Walking: // impossible . c est pour eviter les warnings
1116       {
1117       }
1118     case IntPatch_Restriction: // impossible . c est pour eviter les warnings
1119       {
1120       }
1121     }
1122   }
1123
1124   if (distmin > Tol) {
1125     return Standard_False;
1126   }
1127
1128   typarc = slin.Value(Range)->ArcType();
1129
1130   // Calcul de la tangente.
1131   switch (typarc) {
1132   case IntPatch_Lin :
1133     Vtgtint = (*((Handle(IntPatch_GLine)*)&slin(Range)))->Line().Direction();
1134     break;
1135   case IntPatch_Circle :
1136     Vtgtint = ElCLib::DN(Paraint,(*((Handle(IntPatch_GLine)*)&slin(Range)))->Circle(),1);
1137     break;
1138   case IntPatch_Ellipse :
1139     Vtgtint = ElCLib::DN(Paraint,(*((Handle(IntPatch_GLine)*)&slin(Range)))->Ellipse(),1);
1140     break;
1141   case IntPatch_Parabola :
1142     Vtgtint = ElCLib::DN(Paraint,(*((Handle(IntPatch_GLine)*)&slin(Range)))->Parabola(),1);
1143     break;
1144   case IntPatch_Hyperbola :
1145     Vtgtint = ElCLib::DN(Paraint,(*((Handle(IntPatch_GLine)*)&slin(Range)))->Hyperbola(),1);
1146     break;
1147
1148   case IntPatch_Analytic:
1149     {
1150       const Handle(IntPatch_ALine)& alin = (*((Handle(IntPatch_ALine)*)&slin(Range)));
1151       Standard_Boolean abid = alin->D1(Paraint,pt,Vtgtint);
1152       if (!abid) {
1153         Standard_Real domaininf,domainsup,paramproche;
1154         Standard_Boolean boolbid;
1155         domaininf = alin->FirstParameter(boolbid);
1156         domainsup = alin->LastParameter(boolbid);
1157         if(Paraint>=domaininf && Paraint<=domainsup) {
1158           Standard_Real DeltaParam = 0.001 * (domainsup-domaininf);
1159           if(Paraint-domaininf >= domainsup-Paraint) {
1160             //-- On decale le point vers le parametre le plus eloigne.
1161             DeltaParam = -DeltaParam;
1162           }
1163           Standard_Integer kountbid = 0;
1164           Standard_Boolean bornok = Standard_True;
1165           paramproche = Paraint;
1166           do {
1167             paramproche+=DeltaParam;
1168             kountbid++;
1169             gp_Pnt ptbid;
1170             if(paramproche>=domaininf && paramproche<=domainsup) {
1171               abid = alin->D1(paramproche,ptbid,Vtgtint);
1172             }
1173             else {
1174               bornok = Standard_False;
1175             }
1176           }
1177           while(abid==Standard_False && kountbid<5 && bornok);
1178           //-- Attention aux points de tangence (croisement de 4 lignes )
1179           bornok = Standard_True;
1180           kountbid = 0;
1181           gp_Vec OVtgtint(0.0,0.0,0.0);
1182           paramproche = Paraint;
1183           do {
1184             paramproche-=DeltaParam;
1185             kountbid++;
1186             gp_Pnt ptbid;
1187             if(paramproche>=domaininf && paramproche<=domainsup) {
1188               abid = alin->D1(paramproche,ptbid,OVtgtint);
1189             }
1190             else {
1191               bornok = Standard_False;
1192             }
1193           }
1194           while(abid==Standard_False && kountbid<5 && bornok);
1195           if(bornok) {
1196             paramproche = Vtgtint.Dot(OVtgtint);
1197             if(paramproche<=0.0) abid = Standard_False;
1198           }
1199         }
1200         if(!abid) {
1201           //-- cout << "Pb sur Calcul de derivee 111 " << endl;
1202           Vtgtint.SetCoord(0.,0.,0.);
1203         }
1204       }
1205     }
1206     break;
1207   case IntPatch_Walking: // impossible . c est pour eviter les warnings
1208     {
1209     }
1210   case IntPatch_Restriction: // impossible . c est pour eviter les warnings
1211     {
1212     }
1213
1214   }
1215   return Standard_True;
1216 }
1217
1218
1219 Standard_Boolean  SingleLine (const gp_Pnt& Psurf,
1220                               const Handle(IntPatch_Line)& lin,
1221                               const Standard_Real Tol,
1222                               Standard_Real& Paraint,
1223                               gp_Vec& Vtgtint) {
1224
1225 // Traitement du point de depart ayant pour representation Psurf
1226 // dans l espace. On le replace sur la ligne d intersection; On a en sortie
1227 // son parametre et sa tangente sur la ligne d intersection.
1228 // La fonction renvoie False si le point projete est a une distance
1229 // superieure a Tol du point a projeter.
1230
1231   IntPatch_IType typarc = lin->ArcType();
1232
1233   Standard_Real parproj = 0.;
1234   gp_Vec tgint;
1235   gp_Pnt ptproj;
1236   Standard_Boolean retvalue;
1237
1238
1239   switch (typarc) {
1240   case IntPatch_Lin :
1241     parproj = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Line(),Psurf);
1242     ElCLib::D1(parproj,(*((Handle(IntPatch_GLine)*)&lin))->Line(),ptproj,tgint);
1243     break;
1244   case IntPatch_Circle :
1245     parproj = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Circle(),Psurf);
1246     ElCLib::D1(parproj,(*((Handle(IntPatch_GLine)*)&lin))->Circle(),ptproj,tgint);
1247     break;
1248   case IntPatch_Ellipse :
1249     parproj = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Ellipse(),Psurf);
1250     ElCLib::D1(parproj,(*((Handle(IntPatch_GLine)*)&lin))->Ellipse(),ptproj,tgint);
1251     break;
1252   case IntPatch_Parabola :
1253     parproj = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Parabola(),Psurf);
1254     ElCLib::D1(parproj,(*((Handle(IntPatch_GLine)*)&lin))->Parabola(),ptproj,tgint);
1255     break;
1256   case IntPatch_Hyperbola :
1257     parproj = ElCLib::Parameter((*((Handle(IntPatch_GLine)*)&lin))->Hyperbola(),Psurf);
1258     ElCLib::D1(parproj,(*((Handle(IntPatch_GLine)*)&lin))->Hyperbola(),ptproj,tgint);
1259     break;
1260   case IntPatch_Analytic :
1261     {
1262       const Handle(IntPatch_ALine)& alin = (*((Handle(IntPatch_ALine)*)&lin));
1263       Standard_Boolean ok = alin->FindParameter(Psurf,parproj);
1264       if (ok) {
1265         gp_Pnt ptbid;
1266         Standard_Boolean bid = alin->D1(parproj,ptbid,tgint);
1267         if (!bid) {
1268           Standard_Real domaininf,domainsup,paramproche;
1269           Standard_Boolean boolbid;
1270           domaininf = alin->FirstParameter(boolbid);
1271           domainsup = alin->LastParameter(boolbid);
1272           if(parproj>=domaininf && parproj<=domainsup) {
1273             Standard_Real DeltaParam = 0.001 * (domainsup-domaininf);
1274             if(parproj-domaininf >= domainsup-parproj) {
1275               //-- On decale le point vers le parametre le plus eloigne.
1276               DeltaParam = -DeltaParam;
1277             }
1278             Standard_Integer kountbid = 0;
1279             paramproche = parproj;
1280             do {
1281               paramproche+=DeltaParam;
1282               kountbid++;
1283               bid = alin->D1(paramproche,ptbid,tgint);
1284             }
1285             while(bid==Standard_False && kountbid<5);
1286             ptproj = Psurf;
1287           }
1288           if(!bid) {
1289             //-- cout << "Pb sur Calcul de derivee ALine " << endl;
1290             tgint.SetCoord(0.,0.,0.);
1291             return(Standard_False);
1292           }
1293         }
1294         else {
1295           ptproj = Psurf;
1296         }
1297       }
1298       else {
1299         //-- cout << "---- Pb sur ligne analytique dans SingleLine" << endl;
1300         //-- cout << "     Find Parameter"<<endl;
1301         return Standard_False;
1302       }
1303     }
1304     break;
1305   case IntPatch_Walking: // impossible . c est pour eviter les warnings
1306     {
1307     }
1308   case IntPatch_Restriction: // impossible . c est pour eviter les warnings
1309     {
1310     }
1311   }
1312
1313   if (Psurf.Distance(ptproj) <= Tol) {
1314     Paraint = parproj;
1315     Vtgtint = tgint;
1316     retvalue = Standard_True;
1317   }
1318   else {
1319     retvalue = Standard_False;
1320   }
1321   return retvalue;
1322 }
1323
1324
1325 void ProcessSegments (const IntPatch_SequenceOfSegmentOfTheSOnBounds& listedg,
1326                       IntPatch_SequenceOfLine& slin,
1327                       const IntSurf_Quadric& Quad1,
1328                       const IntSurf_Quadric& Quad2,
1329                       const Standard_Boolean OnFirst,
1330                       const Standard_Real TolArc) {
1331
1332   Standard_Integer i,j,k;
1333   Standard_Integer nbedg = listedg.Length();
1334   Standard_Integer Nblines,Nbpts;
1335
1336   Handle(Adaptor2d_HCurve2d) arcRef;
1337   IntPatch_Point ptvtx, newptvtx;
1338
1339   Handle(IntPatch_RLine)  rline; //-- On fait rline = new ... par la suite
1340
1341   IntPatch_TheSegmentOfTheSOnBounds thesegsol;
1342   IntPatch_ThePathPointOfTheSOnBounds PStartf,PStartl;
1343   Standard_Boolean dofirst,dolast,procf,procl;
1344
1345   Standard_Real paramf =0.,paraml =0.,U1 =0.,V1 =0.,U2 =0.,V2 =0.;
1346
1347   IntPatch_IType typ;
1348   IntSurf_TypeTrans trans1,trans2;
1349   IntSurf_Transition TRest,TArc;
1350   gp_Vec tgline,norm1,norm2,tgarc;
1351   gp_Pnt valpt;
1352
1353   gp_Vec d1u,d1v;
1354   gp_Pnt2d p2d;
1355   gp_Vec2d d2d;
1356
1357
1358   for (i = 1; i <= nbedg; i++) {
1359     Standard_Boolean EdgeDegenere=Standard_False;
1360     thesegsol = listedg.Value(i);
1361     arcRef = thesegsol.Curve();
1362
1363
1364     rline = new IntPatch_RLine(Standard_False);
1365     if(OnFirst) {
1366       rline->SetArcOnS1(arcRef);
1367     }
1368     else {
1369       rline->SetArcOnS2(arcRef);
1370     }
1371
1372 // Traitement des points debut/fin du segment solution.
1373
1374     dofirst = Standard_False;
1375     dolast  = Standard_False;
1376     procf = Standard_False;
1377     procl = Standard_False;
1378
1379     if (thesegsol.HasFirstPoint()) {
1380       dofirst = Standard_True;
1381       PStartf = thesegsol.FirstPoint();
1382       paramf = PStartf.Parameter();
1383     }
1384     if (thesegsol.HasLastPoint()) {
1385       dolast = Standard_True;
1386       PStartl = thesegsol.LastPoint();
1387       paraml = PStartl.Parameter();
1388     }
1389
1390     if (dofirst && dolast) { // determination de la transition de la ligne
1391       arcRef->D1(0.5*(paramf+paraml),p2d,d2d);
1392       if (OnFirst) {
1393         Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1394       }
1395       else {
1396         Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1397       }
1398       tgline.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1399
1400       if(d1u.Magnitude()<1e-7) { //-- edge degenere ?
1401          EdgeDegenere=Standard_True;
1402         for(Standard_Integer edg=0;edg<=10;edg++) {
1403           arcRef->D1(paramf+(paraml-paramf)*edg*0.1,p2d,d2d);
1404           if (OnFirst) {
1405             Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1406           }
1407           else {
1408             Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1409           }
1410
1411           if(d1u.Magnitude()>1e-7) {
1412             EdgeDegenere=Standard_False;
1413           }
1414         }
1415         rline = new IntPatch_RLine(Standard_False);
1416         if(OnFirst) {
1417           rline->SetArcOnS1(arcRef);
1418         }
1419         else {
1420           rline->SetArcOnS2(arcRef);
1421         }
1422       }
1423       else {
1424         norm2 = Quad2.Normale(valpt);
1425         norm1 = Quad1.Normale(valpt);
1426
1427         if (tgline.DotCross(norm2,norm1) > 0.000000001) {
1428           trans1 = IntSurf_Out;
1429           trans2 = IntSurf_In;
1430         }
1431         else if (tgline.DotCross(norm2,norm1) < -0.000000001){
1432           trans1 = IntSurf_In;
1433           trans2 = IntSurf_Out;
1434         }
1435         else {
1436           trans1 = trans2 = IntSurf_Undecided;
1437         }
1438         rline = new IntPatch_RLine(Standard_False,trans1,trans2);
1439         if(OnFirst) {
1440           rline->SetArcOnS1(arcRef);
1441         }
1442         else {
1443           rline->SetArcOnS2(arcRef);
1444         }
1445       }
1446     }
1447     else {
1448       rline = new IntPatch_RLine(Standard_False);
1449       if(OnFirst) {
1450         rline->SetArcOnS1(arcRef);
1451       }
1452       else {
1453         rline->SetArcOnS2(arcRef);
1454       }
1455     }
1456
1457     if (dofirst || dolast) {
1458       Nblines = slin.Length();
1459       for (j=1; j<=Nblines; j++) {
1460         const Handle(IntPatch_Line)& slinj = slin(j);
1461         typ = slinj->ArcType();
1462         if (typ == IntPatch_Analytic) {
1463           Nbpts = (*((Handle(IntPatch_ALine)*)&slinj))->NbVertex();
1464         }
1465         else if (typ == IntPatch_Restriction) {
1466           Nbpts = (*((Handle(IntPatch_RLine)*)&slinj))->NbVertex();
1467         }
1468         else {
1469           Nbpts = (*((Handle(IntPatch_GLine)*)&slinj))->NbVertex();
1470         }
1471         for (k=1; k<=Nbpts;k++) {
1472           if (typ == IntPatch_Analytic) {
1473             ptvtx = (*((Handle(IntPatch_ALine)*)&slinj))->Vertex(k);
1474           }
1475           else if (typ == IntPatch_Restriction) {
1476             ptvtx = (*((Handle(IntPatch_RLine)*)&slinj))->Vertex(k);
1477           }
1478           else {
1479             ptvtx = (*((Handle(IntPatch_GLine)*)&slinj))->Vertex(k);
1480           }
1481
1482           if (EdgeDegenere==Standard_False && dofirst) {
1483             if (ptvtx.Value().Distance(PStartf.Value()) <=TolArc) {
1484               ptvtx.SetMultiple(Standard_True);
1485               if (typ == IntPatch_Analytic) {
1486                 (*((Handle(IntPatch_ALine)*)&slinj))->Replace(k,ptvtx);
1487               }
1488               else if (typ == IntPatch_Restriction) {
1489                 (*((Handle(IntPatch_RLine)*)&slinj))->Replace(k,ptvtx);
1490               }
1491               else {
1492                 (*((Handle(IntPatch_GLine)*)&slinj))->Replace(k,ptvtx);
1493               }
1494               newptvtx = ptvtx;
1495               newptvtx.SetParameter(paramf);
1496               //Recalcul des  transitions si point sur restriction
1497
1498               arcRef->D1(paramf,p2d,d2d);
1499               if (OnFirst) {
1500                 Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1501               }
1502               else {
1503                 Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1504               }
1505               tgline.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1506               if (ptvtx.IsOnDomS1()) {
1507                 const Handle(Adaptor2d_HCurve2d)& thearc = ptvtx.ArcOnS1();
1508                 thearc->D1(ptvtx.ParameterOnArc1(),p2d,d2d);
1509                 Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1510                 tgarc.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1511                 norm1 = d1u.Crossed(d1v);
1512                 if(norm1.SquareMagnitude()<1e-16) {
1513                   TRest.SetValue(Standard_True,IntSurf_Undecided);
1514                   TArc.SetValue(Standard_True,IntSurf_Undecided);
1515                 }
1516                 else {
1517                   IntSurf::MakeTransition(tgline,tgarc,norm1,TRest,TArc);
1518                 }
1519                 newptvtx.SetArc(Standard_True,thearc,ptvtx.ParameterOnArc1(),
1520                                 TRest,TArc);
1521               }
1522               if (ptvtx.IsOnDomS2()) {
1523                 const Handle(Adaptor2d_HCurve2d)& thearc = ptvtx.ArcOnS2();
1524                 thearc->D1(ptvtx.ParameterOnArc2(),p2d,d2d);
1525                 Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1526                 tgarc.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1527                 norm2 = d1u.Crossed(d1v);
1528                 if(norm2.SquareMagnitude()<1e-16) {
1529                   TRest.SetValue(Standard_True,IntSurf_Undecided);
1530                   TArc.SetValue(Standard_True,IntSurf_Undecided);
1531                 }
1532                 else {
1533                   IntSurf::MakeTransition(tgline,tgarc,norm2,TRest,TArc);
1534                 }
1535                 newptvtx.SetArc(Standard_False,thearc,ptvtx.ParameterOnArc2(),
1536                                 TRest,TArc);
1537               }
1538
1539               rline->AddVertex(newptvtx);
1540               if (!procf){
1541                 procf=Standard_True;
1542                 rline->SetFirstPoint(rline->NbVertex());
1543               }
1544             }
1545           }
1546           if (EdgeDegenere==Standard_False && dolast) {
1547             if (ptvtx.Value().Distance(PStartl.Value()) <=TolArc) {
1548               ptvtx.SetMultiple(Standard_True);
1549               if (typ == IntPatch_Analytic) {
1550                 (*((Handle(IntPatch_ALine)*)&slinj))->Replace(k,ptvtx);
1551               }
1552               else if (typ == IntPatch_Restriction) {
1553                 (*((Handle(IntPatch_RLine)*)&slinj))->Replace(k,ptvtx);
1554               }
1555               else {
1556                 (*((Handle(IntPatch_GLine)*)&slinj))->Replace(k,ptvtx);
1557               }
1558
1559               newptvtx = ptvtx;
1560               newptvtx.SetParameter(paraml);
1561               //Recalcul des  transitions si point sur restriction
1562
1563               arcRef->D1(paraml,p2d,d2d);
1564               if (OnFirst) {
1565                 Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1566               }
1567               else {
1568                 Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1569               }
1570               tgline.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1571               if (ptvtx.IsOnDomS1()) {
1572                 const Handle(Adaptor2d_HCurve2d)& thearc = ptvtx.ArcOnS1();
1573                 thearc->D1(ptvtx.ParameterOnArc1(),p2d,d2d);
1574                 Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1575                 tgarc.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1576                 norm1 = d1u.Crossed(d1v);
1577                 if(norm1.SquareMagnitude()<1e-16) {
1578                   TRest.SetValue(Standard_True,IntSurf_Undecided);
1579                   TArc.SetValue(Standard_True,IntSurf_Undecided);
1580                 }
1581                 else {
1582                   IntSurf::MakeTransition(tgline,tgarc,norm1,TRest,TArc);
1583                 }
1584                 newptvtx.SetArc(Standard_True,thearc,ptvtx.ParameterOnArc1(),
1585                                 TRest,TArc);
1586               }
1587               if (ptvtx.IsOnDomS2()) {
1588                 const Handle(Adaptor2d_HCurve2d)& thearc = ptvtx.ArcOnS2();
1589                 thearc->D1(ptvtx.ParameterOnArc2(),p2d,d2d);
1590                 Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1591                 tgarc.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
1592                 norm2 = d1u.Crossed(d1v);
1593                 if(norm2.SquareMagnitude()<1e-16) {
1594                   TRest.SetValue(Standard_True,IntSurf_Undecided);
1595                   TArc.SetValue(Standard_True,IntSurf_Undecided);
1596                 }
1597                 else {
1598                   IntSurf::MakeTransition(tgline,tgarc,norm2,TRest,TArc);
1599                 }
1600                 newptvtx.SetArc(Standard_False,thearc,ptvtx.ParameterOnArc2(),
1601                                 TRest,TArc);
1602               }
1603
1604               rline->AddVertex(newptvtx);
1605               if (!procl){
1606                 procl=Standard_True;
1607                 rline->SetLastPoint(rline->NbVertex());
1608               }
1609             }
1610           }
1611         }
1612 // Si on a traite le pt debut et/ou fin, on ne doit pas recommencer si
1613 // il (ils) correspond(ent) a un point multiple.
1614
1615         if (procf) {
1616           dofirst = Standard_False;
1617         }
1618         if (procl) {
1619           dolast  = Standard_False;
1620         }
1621       }
1622     }
1623 // Si on n a pas trouve le point debut et./ou fin sur une des lignes
1624 // d intersection, il faut quand-meme le placer sur la restriction solution
1625
1626     if (dofirst) {
1627       ptvtx.SetValue(PStartf.Value(),PStartf.Tolerance(),Standard_False);
1628       Quad1.Parameters(PStartf.Value(),U1,V1);
1629       Quad2.Parameters(PStartf.Value(),U2,V2);
1630       ptvtx.SetParameters(U1,V1,U2,V2);
1631       ptvtx.SetParameter(paramf);
1632       if (! PStartf.IsNew()) {
1633         IntSurf_Transition Transline;
1634         IntSurf_Transition Transarc;
1635         ptvtx.SetVertex(OnFirst,PStartf.Vertex());
1636         ptvtx.SetArc(OnFirst,PStartf.Arc(),PStartf.Parameter(),
1637                      Transline,Transarc);
1638       }
1639
1640       rline->AddVertex(ptvtx);
1641       rline->SetFirstPoint(rline->NbVertex());
1642     }
1643     if (dolast) {
1644       ptvtx.SetValue(PStartl.Value(),PStartl.Tolerance(),Standard_False);
1645       Quad1.Parameters(PStartl.Value(),U1,V1);
1646       Quad2.Parameters(PStartl.Value(),U2,V2);
1647       ptvtx.SetParameters(U1,V1,U2,V2);
1648       ptvtx.SetParameter(paraml);
1649       if (! PStartl.IsNew()) {
1650         IntSurf_Transition Transline;
1651         IntSurf_Transition Transarc;
1652
1653         ptvtx.SetVertex(OnFirst,PStartl.Vertex());
1654         ptvtx.SetArc(OnFirst,PStartl.Arc(),PStartl.Parameter(),
1655                      Transline,Transarc);
1656       }
1657
1658       rline->AddVertex(ptvtx);
1659       rline->SetLastPoint(rline->NbVertex());
1660     }
1661     slin.Append(rline);
1662   }
1663 }
1664
1665 inline const gp_Pnt& PointValue(const Handle(IntPatch_RLine) theRLine,
1666                                 const Standard_Integer theIndex)
1667 {
1668   return theRLine->Point(theIndex).Value();
1669 }
1670
1671 inline const gp_Pnt& VertexValue( const Handle(IntPatch_RLine) theRLine,
1672                                   const Standard_Integer theIndex)
1673 {
1674   return theRLine->Vertex(theIndex).Value();
1675 }
1676
1677 static Standard_Real SquareDistance(const Handle(IntPatch_GLine)& theGLine,
1678                                     const gp_Pnt& theP,
1679                                     Extrema_ExtPC& theExtr)
1680 {
1681   Standard_Real aSQDist = RealLast();
1682   switch(theGLine->ArcType())
1683   {
1684   case IntPatch_Lin:
1685     aSQDist = theGLine->Line().SquareDistance(theP);
1686     break;
1687   case IntPatch_Circle:
1688     aSQDist = theGLine->Circle().SquareDistance(theP);
1689     break;
1690   default:
1691     theExtr.Perform(theP);
1692     if(!theExtr.IsDone() || !theExtr.NbExt())
1693     {
1694       //Lines are not overlapped
1695       return aSQDist;
1696     }
1697
1698     aSQDist = theExtr.SquareDistance(1);
1699     const Standard_Integer aNbExtr = theExtr.NbExt();
1700     for ( Standard_Integer i = 2; i <= aNbExtr; i++)
1701     {
1702       const Standard_Real aSQD = theExtr.SquareDistance(i);
1703       if (aSQD < aSQDist)
1704       {
1705         aSQDist = aSQD;
1706       }
1707     }
1708   }
1709
1710   return aSQDist;
1711 }
1712
1713 static Standard_Boolean IsRLineGood(const IntSurf_Quadric& Quad1,
1714                                     const IntSurf_Quadric& Quad2,
1715                                     const Handle(IntPatch_GLine) theGLine,
1716                                     const Handle(IntPatch_RLine) theRLine,
1717                                     const Standard_Real theTol)
1718 {
1719   const Standard_Real aSQTol = theTol*theTol;
1720   const IntPatch_IType aGType = theGLine->ArcType();
1721   Standard_Integer aNbPntsM1 = 0;
1722
1723   const gp_Pnt& (*Value) (const Handle(IntPatch_RLine), const Standard_Integer);
1724
1725   if(theRLine->HasPolygon())
1726   {
1727     aNbPntsM1 = theRLine->NbPnts()-1;
1728     Value = PointValue;
1729   }
1730   else
1731   {
1732     aNbPntsM1 = theRLine->NbVertex()-1;
1733     Value = VertexValue;
1734   }
1735
1736   if(aNbPntsM1 < 1)
1737     return Standard_False;
1738
1739   Extrema_ExtPC anExtr;
1740   GeomAdaptor_Curve anAC;
1741   Handle(Geom_Curve) aCurv;
1742
1743   if(aGType == IntPatch_Ellipse)
1744     aCurv = new Geom_Ellipse(theGLine->Ellipse());
1745   else if(aGType == IntPatch_Parabola)
1746     aCurv = new Geom_Parabola(theGLine->Parabola());
1747   else if(aGType == IntPatch_Hyperbola)
1748     aCurv = new Geom_Hyperbola(theGLine->Hyperbola());
1749
1750   if(!aCurv.IsNull())
1751   {
1752     const Standard_Real anUinf = aCurv->FirstParameter(),
1753                         anUsup = aCurv->LastParameter();
1754     anAC.Load(aCurv, anUinf, anUsup);
1755     anExtr.Initialize(anAC, anUinf, anUsup);
1756   }
1757
1758   if(aNbPntsM1 == 1)
1759   {
1760     gp_Pnt aP1(Value(theRLine, 1)), aP2(Value(theRLine, 2));
1761
1762     if(aP1.SquareDistance(aP2) < aSQTol)
1763     {
1764       //RLine is degenerated
1765       return Standard_False;
1766     }
1767
1768     gp_Pnt aPMid;
1769     if(theRLine->IsArcOnS1())
1770     {
1771       const Handle(Adaptor2d_HCurve2d)& anAC2d = theRLine->ArcOnS1();
1772       const Standard_Real aParF = anAC2d->FirstParameter(),
1773                           aParL = anAC2d->LastParameter();
1774       gp_Pnt2d aP2d(anAC2d->Value(0.5*(aParF+aParL)));
1775       aPMid = Quad1.Value(aP2d.X(), aP2d.Y());
1776     }
1777     else
1778     {
1779       const Handle(Adaptor2d_HCurve2d)& anAC2d = theRLine->ArcOnS2();
1780       const Standard_Real aParF = anAC2d->FirstParameter(),
1781                           aParL = anAC2d->LastParameter();
1782       gp_Pnt2d aP2d(anAC2d->Value(0.5*(aParF+aParL)));
1783       aPMid = Quad2.Value(aP2d.X(), aP2d.Y());
1784     }
1785
1786     const Standard_Real aSQDist = SquareDistance(theGLine, aPMid, anExtr);
1787     return (aSQDist > aSQTol);
1788   }
1789
1790   for(Standard_Integer i = 2; i <= aNbPntsM1; i++)
1791   {
1792     const gp_Pnt aP(Value(theRLine, i));
1793     const Standard_Real aSQDist = SquareDistance(theGLine, aP, anExtr);
1794
1795     if(aSQDist > aSQTol)
1796       return Standard_True;
1797   }
1798
1799   return Standard_False;
1800 }
1801
1802
1803 void ProcessRLine (IntPatch_SequenceOfLine& slin,
1804 //                 const Handle(Adaptor3d_HSurface)& Surf1,
1805 //                 const Handle(Adaptor3d_HSurface)& Surf2,
1806                    const IntSurf_Quadric& Quad1,
1807                    const IntSurf_Quadric& Quad2,
1808                    const Standard_Real _TolArc,
1809                    const Standard_Boolean theIsReqToKeepRLine) {
1810
1811 // On cherche a placer sur les restrictions solutions les points "multiples"
1812 // des autres lignes d intersection
1813 // Pas forcemment le plus efficace : on rique de projeter plusieurs fois
1814 // le meme point sur la meme restriction...
1815
1816   Standard_Real TolArc=100.0*_TolArc;
1817   if(TolArc>0.1) TolArc=0.1;
1818
1819   Standard_Integer i,j,k;
1820   Standard_Integer Nblin,Nbvtx, Nbpt;
1821
1822   Standard_Boolean OnFirst = Standard_False,project = Standard_False,keeppoint = Standard_False;
1823
1824   Handle(Adaptor2d_HCurve2d) arcref;
1825   Standard_Real paramproj,paramf,paraml;
1826
1827   TColgp_SequenceOfPnt seq_Pnt3d;
1828   TColStd_SequenceOfReal seq_Real;
1829
1830   gp_Pnt ptproj,toproj,valpt;
1831
1832   gp_Pnt2d p2d;
1833   gp_Vec2d d2d;
1834   gp_Vec d1u,d1v,tgrest,tgarc,norm;
1835   IntSurf_Transition TRest,TArc;
1836 #ifndef OCCT_DEBUG
1837   Standard_Real U =0.,V =0.;
1838 #else
1839   Standard_Real U,V;
1840 #endif
1841   IntPatch_Point Ptvtx,newptvtx;
1842
1843   IntPatch_IType typ1,typ2;
1844
1845
1846   Nblin = slin.Length();
1847   for (i=1; i<=Nblin; i++) {
1848     const Handle(IntPatch_Line)& slini = slin(i);
1849     typ1 = slini->ArcType();
1850
1851     Standard_Boolean HasToDeleteRLine = Standard_False;
1852     if (typ1 == IntPatch_Restriction) {
1853       seq_Pnt3d.Clear();
1854       seq_Real.Clear();
1855
1856       for (j=1; j<=Nblin; j++) {
1857         const  Handle(IntPatch_Line)& slinj = slin(j);
1858         Nbpt = seq_Pnt3d.Length();      // important que ce soit ici
1859         typ2 = slinj->ArcType();
1860         if (typ2 != IntPatch_Restriction) {
1861           //-- arcref = (*((Handle(IntPatch_RLine)*)&slini))->Arc();
1862           //-- OnFirst = (*((Handle(IntPatch_RLine)*)&slini))->IsOnFirstSurface();
1863
1864           //-- DES CHOSES A FAIRE ICI
1865           if((*((Handle(IntPatch_RLine)*)&slini))->IsArcOnS1()) {
1866             OnFirst=Standard_True;
1867             arcref= (*((Handle(IntPatch_RLine)*)&slini))->ArcOnS1();
1868           }
1869           else if((*((Handle(IntPatch_RLine)*)&slini))->IsArcOnS2()) {
1870             arcref= (*((Handle(IntPatch_RLine)*)&slini))->ArcOnS2();
1871             OnFirst=Standard_False;
1872           }
1873           if ((*((Handle(IntPatch_RLine)*)&slini))->HasFirstPoint()) {
1874             paramf = (*((Handle(IntPatch_RLine)*)&slini))->FirstPoint().ParameterOnLine();
1875           }
1876           else {
1877             // cout << "Pas de param debut sur rst solution" << endl;
1878             paramf = RealFirst();
1879           }
1880           if ((*((Handle(IntPatch_RLine)*)&slini))->HasLastPoint()) {
1881             paraml = (*((Handle(IntPatch_RLine)*)&slini))->LastPoint().ParameterOnLine();
1882           }
1883           else {
1884             // cout << "Pas de param debut sur rst solution" << endl;
1885             paraml = RealLast();
1886           }
1887
1888           if (typ2 == IntPatch_Analytic) {
1889             Nbvtx = (*((Handle(IntPatch_ALine)*)&slinj))->NbVertex();
1890           }
1891           else {
1892             Nbvtx = (*((Handle(IntPatch_GLine)*)&slinj))->NbVertex();
1893           }
1894
1895
1896           Standard_Boolean EdgeDegenere=Standard_True;
1897           for(Standard_Integer edg=0;EdgeDegenere && edg<=10;edg++) {
1898             arcref->D1(paramf+(paraml-paramf)*edg*0.1,p2d,d2d);
1899             if (OnFirst) {
1900               Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1901             }
1902             else {
1903               Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
1904             }
1905             if(d1u.Magnitude()>1e-7) {
1906               EdgeDegenere=Standard_False;
1907             }
1908           }
1909
1910           for (k=1; EdgeDegenere==Standard_False && k<=Nbvtx; k++) {
1911             if (typ2 == IntPatch_Analytic) {
1912               Ptvtx = (*((Handle(IntPatch_ALine)*)&slinj))->Vertex(k);
1913             }
1914             else {
1915               Ptvtx = (*((Handle(IntPatch_GLine)*)&slinj))->Vertex(k);
1916             }
1917             if ((OnFirst && !Ptvtx.IsOnDomS1()) ||
1918                 (!OnFirst && !Ptvtx.IsOnDomS2())) {
1919               // Si OnFirst && OnDomS1, c est qu on est a une extremite
1920               // ca doit etre traite par Process Segment...
1921               project = Standard_True;
1922               keeppoint = Standard_False;
1923               toproj = Ptvtx.Value();
1924
1925               Standard_Integer jj;
1926               for (jj = 1; jj <= Nbpt; jj++) {
1927               //for (Standard_Integer jj = 1; jj <= Nbpt; jj++) {
1928                 if (toproj.Distance(seq_Pnt3d(jj)) < _TolArc) {
1929                   project = Standard_False;
1930                   break;
1931                 }
1932               }
1933               if (project) { //-- il faut projeter pour trouver le point sur la rline.
1934                 if (OnFirst) {
1935                   Ptvtx.ParametersOnS1(U,V);
1936                 }
1937                 else {
1938                   Ptvtx.ParametersOnS2(U,V);
1939                 }
1940
1941                 project = IntPatch_HInterTool::Project(arcref,gp_Pnt2d(U,V),
1942                                                 paramproj,p2d);
1943
1944                 if (project) {
1945                   if (OnFirst) {
1946                     ptproj = Quad1.Value(p2d.X(),p2d.Y());
1947                   }
1948                   else {
1949                     ptproj = Quad2.Value(p2d.X(),p2d.Y());
1950                   }
1951                   if ((toproj.Distance(ptproj) <=100*TolArc) &&
1952                       (paramproj >= paramf) && (paramproj <= paraml)){
1953                     newptvtx = Ptvtx;
1954                     newptvtx.SetParameter(paramproj);
1955                     keeppoint = Standard_True;
1956                     seq_Pnt3d.Append(toproj);
1957                     seq_Real.Append(paramproj);
1958
1959                     //-- verifier que si la restriction arcref est trouvee, elle porte ce vertex
1960                     for (int ri=1; ri<=Nblin; ri++) {
1961                       const Handle(IntPatch_Line)& slinri = slin(ri);
1962                       if (slinri->ArcType() == IntPatch_Restriction) {
1963                         if(OnFirst && (*((Handle(IntPatch_RLine)*)&slinri))->IsArcOnS1()) {
1964                           if(arcref == (*((Handle(IntPatch_RLine)*)&slinri))->ArcOnS1()) {
1965                             (*((Handle(IntPatch_RLine)*)&slinri))->AddVertex(newptvtx);
1966                             //printf("\n ImpImpIntersection_0.gxx CAS1 \n");
1967                           }
1968                         }
1969                         else if(OnFirst==Standard_False && (*((Handle(IntPatch_RLine)*)&slinri))->IsArcOnS2()) {
1970                           if(arcref == (*((Handle(IntPatch_RLine)*)&slinri))->ArcOnS2()) {
1971                             (*((Handle(IntPatch_RLine)*)&slinri))->AddVertex(newptvtx);
1972                             //printf("\n ImpImpIntersection_0.gxx CAS2 \n");
1973                           }
1974                         }
1975                       }
1976                     }
1977                     // -- --------------------------------------------------
1978                   }
1979                 }
1980               }
1981               else {
1982                 keeppoint = Standard_True;
1983                 newptvtx = Ptvtx;
1984                 newptvtx.SetParameter(seq_Real(jj));
1985               }
1986               if (keeppoint) {
1987                 Ptvtx.SetMultiple(Standard_True);
1988                 newptvtx.SetMultiple(Standard_True);
1989
1990                 if (typ2 == IntPatch_Analytic) {
1991                   (*((Handle(IntPatch_ALine)*)&slinj))->Replace(k,Ptvtx);
1992                 }
1993                 else {
1994                   (*((Handle(IntPatch_GLine)*)&slinj))->Replace(k,Ptvtx);
1995                 }
1996
1997                 if (Ptvtx.IsOnDomS1() || Ptvtx.IsOnDomS2()) {
1998
1999                   arcref->D1(newptvtx.ParameterOnLine(),p2d,d2d);
2000
2001                   if (OnFirst) { // donc OnDomS2
2002                     Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
2003                     tgrest.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
2004
2005                     const Handle(Adaptor2d_HCurve2d)& thearc = Ptvtx.ArcOnS2();
2006                     thearc->D1(Ptvtx.ParameterOnArc2(),p2d,d2d);
2007                     Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
2008                     tgarc.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
2009                     norm = d1u.Crossed(d1v); //Quad2.Normale(valpt);
2010                     if(norm.SquareMagnitude()<1e-16) {
2011                       TRest.SetValue(Standard_True,IntSurf_Undecided);
2012                       TArc.SetValue(Standard_True,IntSurf_Undecided);
2013                     }
2014                     else {
2015                       IntSurf::MakeTransition(tgrest,tgarc,norm,TRest,TArc);
2016                     }
2017                     newptvtx.SetArc(Standard_False,thearc,
2018                                     Ptvtx.ParameterOnArc2(),TRest,TArc);
2019
2020                   }
2021                   else { // donc OnDomS1
2022                     Quad2.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
2023                     tgrest.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
2024
2025                     const Handle(Adaptor2d_HCurve2d)& thearc = Ptvtx.ArcOnS1();
2026                     thearc->D1(Ptvtx.ParameterOnArc1(),p2d,d2d);
2027                     Quad1.D1(p2d.X(),p2d.Y(),valpt,d1u,d1v);
2028                     tgarc.SetLinearForm(d2d.X(),d1u,d2d.Y(),d1v);
2029                     norm = d1u.Crossed(d1v); //Quad1.Normale(valpt);
2030                     if(norm.SquareMagnitude()<1e-16) {
2031                       TRest.SetValue(Standard_True,IntSurf_Undecided);
2032                       TArc.SetValue(Standard_True,IntSurf_Undecided);
2033                     }
2034                     else {
2035                       IntSurf::MakeTransition(tgrest,tgarc,norm,TRest,TArc);
2036                     }
2037                     newptvtx.SetArc(Standard_True,thearc,
2038                                     Ptvtx.ParameterOnArc1(),TRest,TArc);
2039                   }
2040                 } //-- if (Ptvtx.IsOnDomS1() || Ptvtx.IsOnDomS2())
2041
2042                 (*((Handle(IntPatch_RLine)*)&slini))->AddVertex(newptvtx);
2043
2044               } //-- if (keeppoint)
2045             } //-- if ((OnFirst && !Ptvtx.IsOnDomS1())||(!OnFirst && !Ptvtx.IsOnDomS2()))
2046           } //-- boucle sur les vertex
2047
2048           if(!theIsReqToKeepRLine)
2049           {
2050             Handle(IntPatch_GLine) aGL = Handle(IntPatch_GLine)::DownCast(slinj);
2051
2052             if(!aGL.IsNull())
2053             {
2054               HasToDeleteRLine = !IsRLineGood(Quad1, Quad2, aGL,
2055                                         Handle(IntPatch_RLine)::DownCast(slini), TolArc);
2056             }
2057
2058             if(HasToDeleteRLine)
2059             {
2060               break;
2061             }
2062           }
2063         } //-- if (typ2 != IntPatch_Restriction)
2064       } //-- for (j=1; j<=Nblin; j++)
2065     } //-- if (typ1 == IntPatch_Restriction)
2066
2067     if(HasToDeleteRLine)
2068     {
2069       slin.Remove(i);
2070       i--;
2071       Nblin = slin.Length();
2072       continue;
2073     }
2074   } //-- for (i=1; i<=Nblin; i++)
2075 }