src/IntPatch/IntPatch_Intersection.cxx

   1 // Created by: Modelization
   2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   3 //
   4 // This file is part of Open CASCADE Technology software library.
   5 //
   6 // This library is free software; you can redistribute it and/or modify it under
   7 // the terms of the GNU Lesser General Public License version 2.1 as published
   8 // by the Free Software Foundation, with special exception defined in the file
   9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
  10 // distribution for complete text of the license and disclaimer of any warranty.
  11 //
  12 // Alternatively, this file may be used under the terms of Open CASCADE
  13 // commercial license or contractual agreement.
  14
  15
  16 #include <Adaptor3d_HSurface.hxx>
  17 #include <Adaptor3d_TopolTool.hxx>
  18 #include <IntPatch_ALine.hxx>
  19 #include <IntPatch_ALineToWLine.hxx>
  20 #include <IntPatch_GLine.hxx>
  21 #include <IntPatch_ImpImpIntersection.hxx>
  22 #include <IntPatch_ImpPrmIntersection.hxx>
  23 #include <IntPatch_Intersection.hxx>
  24 #include <IntPatch_Line.hxx>
  25 #include <IntPatch_Point.hxx>
  26 #include <IntPatch_PrmPrmIntersection.hxx>
  27 #include <IntPatch_RLine.hxx>
  28 #include <IntPatch_WLine.hxx>
  29 #include <IntSurf_Quadric.hxx>
  30 #include <Standard_ConstructionError.hxx>
  31 #include <Standard_DomainError.hxx>
  32 #include <Standard_OutOfRange.hxx>
  33 #include <StdFail_NotDone.hxx>
  34
  35 #include <stdio.h>
  36 #define DEBUG 0
  37 static const Standard_Integer aNbPointsInALine = 200;
  38
  39 //=======================================================================
  40 //function : MinMax
  41 //purpose  : Replaces theParMIN = MIN(theParMIN, theParMAX),
  42 //                    theParMAX = MAX(theParMIN, theParMAX).
  43 //=======================================================================
  44 static inline void MinMax(Standard_Real& theParMIN, Standard_Real& theParMAX)
  45 {
  46   if(theParMIN > theParMAX)
  47   {
  48     const Standard_Real aTmp = theParMAX;
  49     theParMAX = theParMIN;
  50     theParMIN = aTmp;
  51   }
  52 }
  53
  54 //=======================================================================
  55 //function : IsSeam
  56 //purpose  : Returns:
  57 //            0 - if interval [theU1, theU2] does not intersect the "seam-edge"
  58 //                or if "seam-edge" do not exist;
  59 //            1 - if interval (theU1, theU2) intersect the "seam-edge".
  60 //            2 - if theU1 or/and theU2 lie ON the "seam-edge"
  61 //
  62 //ATTENTION!!!
  63 //  If (theU1 == theU2) then this function will return only both 0 or 2.
  64 //=======================================================================
  65 static Standard_Integer IsSeam( const Standard_Real theU1,
  66                                 const Standard_Real theU2,
  67                                 const Standard_Real thePeriod)
  68 {
  69   if(IsEqual(thePeriod, 0.0))
  70     return 0;
  71
  72   //If interval [theU1, theU2] intersect seam-edge then there exists an integer
  73   //number N such as
  74   //    (theU1 <= T*N <= theU2) <=> (theU1/T <= N <= theU2/T),
  75   //where T is the period.
  76   //I.e. the inerval [theU1/T, theU2/T] must contain at least one
  77   //integer number. In this case, Floor(theU1/T) and Floor(theU2/T)
  78   //return different values or theU1/T is strictly integer number.
  79   //Examples:
  80   //  1. theU1/T==2.8, theU2/T==3.5 => Floor(theU1/T) == 2, Floor(theU2/T) == 3.
  81   //  2. theU1/T==2.0, theU2/T==2.6 => Floor(theU1/T) == Floor(theU2/T) == 2.
  82
  83   const Standard_Real aVal1 = theU1/thePeriod,
  84                       aVal2 = theU2/thePeriod;
  85   const Standard_Integer aPar1 = static_cast<Standard_Integer>(Floor(aVal1));
  86   const Standard_Integer aPar2 = static_cast<Standard_Integer>(Floor(aVal2));
  87   if(aPar1 != aPar2)
  88   {//Interval (theU1, theU2] intersects seam-edge
  89     if(IsEqual(aVal2, static_cast<Standard_Real>(aPar2)))
  90     {//aVal2 is an integer number => theU2 lies ON the "seam-edge"
  91       return 2;
  92     }
  93
  94     return 1;
  95   }
  96
  97   //Here, aPar1 == aPar2.
  98
  99   if(IsEqual(aVal1, static_cast<Standard_Real>(aPar1)))
 100   {//aVal1 is an integer number => theU1 lies ON the "seam-edge"
 101     return 2;
 102   }
 103
 104   //If aVal2 is a true integer number then always (aPar1 != aPar2).
 105
 106   return 0;
 107 }
 108
 109 //=======================================================================
 110 //function : IsSeamOrBound
 111 //purpose  : Returns TRUE if segment [thePtf, thePtl] intersects "seam-edge"
 112 //            (if it exist) or surface boundaries and both thePtf and thePtl do
 113 //            not match "seam-edge" or boundaries.
 114 //           Point thePtmid lies in this segment. If thePtmid match
 115 //            "seam-edge" or boundaries strictly (without any tolerance) then
 116 //            the function will return TRUE.
 117 //            See comments in function body for detail information.
 118 //=======================================================================
 119 static Standard_Boolean IsSeamOrBound(const IntSurf_PntOn2S& thePtf,
 120                                       const IntSurf_PntOn2S& thePtl,
 121                                       const IntSurf_PntOn2S& thePtmid,
 122                                       const Standard_Real theU1Period,
 123                                       const Standard_Real theU2Period,
 124                                       const Standard_Real theV1Period,
 125                                       const Standard_Real theV2Period,
 126                                       const Standard_Real theUfSurf1,
 127                                       const Standard_Real theUlSurf1,
 128                                       const Standard_Real theVfSurf1,
 129                                       const Standard_Real theVlSurf1,
 130                                       const Standard_Real theUfSurf2,
 131                                       const Standard_Real theUlSurf2,
 132                                       const Standard_Real theVfSurf2,
 133                                       const Standard_Real theVlSurf2)
 134 {
 135   Standard_Real aU11 = 0.0, aU12 = 0.0, aV11 = 0.0, aV12 = 0.0;
 136   Standard_Real aU21 = 0.0, aU22 = 0.0, aV21 = 0.0, aV22 = 0.0;
 137   thePtf.Parameters(aU11, aV11, aU12, aV12);
 138   thePtl.Parameters(aU21, aV21, aU22, aV22);
 139
 140   MinMax(aU11, aU21);
 141   MinMax(aV11, aV21);
 142   MinMax(aU12, aU22);
 143   MinMax(aV12, aV22);
 144
 145   if((aU11 - theUfSurf1)*(aU21 - theUfSurf1) < 0.0)
 146   {//Interval [aU11, aU21] intersects theUfSurf1
 147     return Standard_True;
 148   }
 149
 150   if((aU11 - theUlSurf1)*(aU21 - theUlSurf1) < 0.0)
 151   {//Interval [aU11, aU21] intersects theUlSurf1
 152     return Standard_True;
 153   }
 154
 155   if((aV11 - theVfSurf1)*(aV21 - theVfSurf1) < 0.0)
 156   {//Interval [aV11, aV21] intersects theVfSurf1
 157     return Standard_True;
 158   }
 159
 160   if((aV11 - theVlSurf1)*(aV21 - theVlSurf1) < 0.0)
 161   {//Interval [aV11, aV21] intersects theVlSurf1
 162     return Standard_True;
 163   }
 164
 165   if((aU12 - theUfSurf2)*(aU22 - theUfSurf2) < 0.0)
 166   {//Interval [aU12, aU22] intersects theUfSurf2
 167     return Standard_True;
 168   }
 169
 170   if((aU12 - theUlSurf2)*(aU22 - theUlSurf2) < 0.0)
 171   {//Interval [aU12, aU22] intersects theUlSurf2
 172     return Standard_True;
 173   }
 174
 175   if((aV12 - theVfSurf2)*(aV22 - theVfSurf2) < 0.0)
 176   {//Interval [aV12, aV22] intersects theVfSurf2
 177     return Standard_True;
 178   }
 179
 180   if((aV12 - theVlSurf2)*(aV22 - theVlSurf2) < 0.0)
 181   {//Interval [aV12, aV22] intersects theVlSurf2
 182     return Standard_True;
 183   }
 184
 185   if(IsSeam(aU11, aU21, theU1Period))
 186     return Standard_True;
 187
 188   if(IsSeam(aV11, aV21, theV1Period))
 189     return Standard_True;
 190
 191   if(IsSeam(aU12, aU22, theU2Period))
 192     return Standard_True;
 193
 194   if(IsSeam(aV12, aV22, theV2Period))
 195     return Standard_True;
 196
 197   /*
 198     The segment [thePtf, thePtl] does not intersect the boundaries and
 199     the seam-edge of the surfaces.
 200     Nevertheless, following situation is possible:
 201
 202                   seam or
 203                    bound
 204                      |
 205         thePtf  *    |
 206                      |
 207                      * thePtmid
 208           thePtl  *  |
 209                      |
 210
 211     This case must be processed, too.
 212   */
 213
 214   Standard_Real aU1 = 0.0, aU2 = 0.0, aV1 = 0.0, aV2 = 0.0;
 215   thePtmid.Parameters(aU1, aV1, aU2, aV2);
 216
 217   if(IsEqual(aU1, theUfSurf1) || IsEqual(aU1, theUlSurf1))
 218     return Standard_True;
 219
 220   if(IsEqual(aU2, theUfSurf2) || IsEqual(aU2, theUlSurf2))
 221     return Standard_True;
 222
 223   if(IsEqual(aV1, theVfSurf1) || IsEqual(aV1, theVlSurf1))
 224     return Standard_True;
 225
 226   if(IsEqual(aV2, theVfSurf2) || IsEqual(aV2, theVlSurf2))
 227     return Standard_True;
 228
 229   if(IsSeam(aU1, aU1, theU1Period))
 230     return Standard_True;
 231
 232   if(IsSeam(aU2, aU2, theU2Period))
 233     return Standard_True;
 234
 235   if(IsSeam(aV1, aV1, theV1Period))
 236     return Standard_True;
 237
 238   if(IsSeam(aV2, aV2, theV2Period))
 239     return Standard_True;
 240
 241   return Standard_False;
 242 }
 243
 244 //=======================================================================
 245 //function : JoinWLines
 246 //purpose  : joins all WLines from theSlin to one if it is possible and
 247 //            records the result into theSlin again.
 248 //            Lines will be kept to be splitted if:
 249 //              a) they are separated (has no common points);
 250 //              b) resulted line (after joining) go through
 251 //                 seam-edges or surface boundaries.
 252 //
 253 //          In addition, if points in theSPnt lies at least in one of
 254 //          the line in theSlin, this point will be deleted.
 255 //=======================================================================
 256 static void JoinWLines(IntPatch_SequenceOfLine& theSlin,
 257                 IntPatch_SequenceOfPoint& theSPnt,
 258                 const Standard_Real theTol3D,
 259                 const Standard_Real theU1Period,
 260                 const Standard_Real theU2Period,
 261                 const Standard_Real theV1Period,
 262                 const Standard_Real theV2Period,
 263                 const Standard_Real theUfSurf1,
 264                 const Standard_Real theUlSurf1,
 265                 const Standard_Real theVfSurf1,
 266                 const Standard_Real theVlSurf1,
 267                 const Standard_Real theUfSurf2,
 268                 const Standard_Real theUlSurf2,
 269                 const Standard_Real theVfSurf2,
 270                 const Standard_Real theVlSurf2)
 271 {
 272   if(theSlin.Length() == 0)
 273     return;
 274
 275   for(Standard_Integer aNumOfLine1 = 1; aNumOfLine1 <= theSlin.Length(); aNumOfLine1++)
 276   {
 277     Handle(IntPatch_WLine) aWLine1 (Handle(IntPatch_WLine)::DownCast(theSlin.Value(aNumOfLine1)));
 278
 279     if(aWLine1.IsNull())
 280     {//We must have failed to join not-point-lines
 281       return;
 282     }
 283
 284     const Standard_Integer aNbPntsWL1 = aWLine1->NbPnts();
 285     const IntSurf_PntOn2S& aPntFWL1 = aWLine1->Point(1);
 286     const IntSurf_PntOn2S& aPntLWL1 = aWLine1->Point(aNbPntsWL1);
 287
 288     for(Standard_Integer aNPt = 1; aNPt <= theSPnt.Length(); aNPt++)
 289     {
 290       const IntSurf_PntOn2S aPntCur = theSPnt.Value(aNPt).PntOn2S();
 291
 292       if( aPntCur.IsSame(aPntFWL1, Precision::Confusion()) ||
 293         aPntCur.IsSame(aPntLWL1, Precision::Confusion()))
 294       {
 295         theSPnt.Remove(aNPt);
 296         aNPt--;
 297       }
 298     }
 299
 300     Standard_Boolean hasBeenRemoved = Standard_False;
 301     for(Standard_Integer aNumOfLine2 = aNumOfLine1 + 1; aNumOfLine2 <= theSlin.Length(); aNumOfLine2++)
 302     {
 303       Handle(IntPatch_WLine) aWLine2 (Handle(IntPatch_WLine)::DownCast(theSlin.Value(aNumOfLine2)));
 304
 305       const Standard_Integer aNbPntsWL2 = aWLine2->NbPnts();
 306
 307       const IntSurf_PntOn2S& aPntFWL1 = aWLine1->Point(1);
 308       const IntSurf_PntOn2S& aPntLWL1 = aWLine1->Point(aNbPntsWL1);
 309
 310       const IntSurf_PntOn2S& aPntFWL2 = aWLine2->Point(1);
 311       const IntSurf_PntOn2S& aPntLWL2 = aWLine2->Point(aNbPntsWL2);
 312
 313       if(aPntFWL1.IsSame(aPntFWL2, Precision::Confusion()))
 314       {
 315         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(2);
 316         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(2);
 317         if(!IsSeamOrBound(aPt1, aPt2, aPntFWL1, theU1Period, theU2Period,
 318                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
 319                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
 320                           theVfSurf2, theVlSurf2))
 321         {
 322           aWLine1->ClearVertexes();
 323           for(Standard_Integer aNPt = 1; aNPt <= aNbPntsWL2; aNPt++)
 324           {
 325             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
 326             aWLine1->Curve()->InsertBefore(1, aPt);
 327           }
 328
 329           aWLine1->ComputeVertexParameters(theTol3D);
 330
 331           theSlin.Remove(aNumOfLine2);
 332           aNumOfLine2--;
 333           hasBeenRemoved = Standard_True;
 334
 335           continue;
 336         }
 337       }
 338
 339       if(aPntFWL1.IsSame(aPntLWL2, Precision::Confusion()))
 340       {
 341         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(2);
 342         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(aNbPntsWL2-1);
 343         if(!IsSeamOrBound(aPt1, aPt2, aPntFWL1, theU1Period, theU2Period,
 344                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
 345                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
 346                           theVfSurf2, theVlSurf2))
 347         {
 348           aWLine1->ClearVertexes();
 349           for(Standard_Integer aNPt = aNbPntsWL2; aNPt >= 1; aNPt--)
 350           {
 351             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
 352             aWLine1->Curve()->InsertBefore(1, aPt);
 353           }
 354
 355           aWLine1->ComputeVertexParameters(theTol3D);
 356
 357           theSlin.Remove(aNumOfLine2);
 358           aNumOfLine2--;
 359           hasBeenRemoved = Standard_True;
 360
 361           continue;
 362         }
 363       }
 364
 365       if(aPntLWL1.IsSame(aPntFWL2, Precision::Confusion()))
 366       {
 367         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(aNbPntsWL1-1);
 368         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(2);
 369         if(!IsSeamOrBound(aPt1, aPt2, aPntLWL1, theU1Period, theU2Period,
 370                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
 371                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
 372                           theVfSurf2, theVlSurf2))
 373         {
 374           aWLine1->ClearVertexes();
 375           for(Standard_Integer aNPt = 1; aNPt <= aNbPntsWL2; aNPt++)
 376           {
 377             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
 378             aWLine1->Curve()->Add(aPt);
 379           }
 380
 381           aWLine1->ComputeVertexParameters(theTol3D);
 382
 383           theSlin.Remove(aNumOfLine2);
 384           aNumOfLine2--;
 385           hasBeenRemoved = Standard_True;
 386
 387           continue;
 388         }
 389       }
 390
 391       if(aPntLWL1.IsSame(aPntLWL2, Precision::Confusion()))
 392       {
 393         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(aNbPntsWL1-1);
 394         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(aNbPntsWL2-1);
 395         if(!IsSeamOrBound(aPt1, aPt2, aPntLWL1, theU1Period, theU2Period,
 396                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
 397                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
 398                           theVfSurf2, theVlSurf2))
 399         {
 400           aWLine1->ClearVertexes();
 401           for(Standard_Integer aNPt = aNbPntsWL2; aNPt >= 1; aNPt--)
 402           {
 403             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
 404             aWLine1->Curve()->Add(aPt);
 405           }
 406
 407           aWLine1->ComputeVertexParameters(theTol3D);
 408
 409           theSlin.Remove(aNumOfLine2);
 410           aNumOfLine2--;
 411           hasBeenRemoved = Standard_True;
 412
 413           continue;
 414         }
 415       }
 416     }
 417
 418     if(hasBeenRemoved)
 419       aNumOfLine1--;
 420   }
 421 }
 422
 423 //======================================================================
 424 // function: SequenceOfLine
 425 //======================================================================
 426 const IntPatch_SequenceOfLine& IntPatch_Intersection::SequenceOfLine() const { return(slin); }
 427
 428 //======================================================================
 429 // function: IntPatch_Intersection
 430 //======================================================================
 431 IntPatch_Intersection::IntPatch_Intersection ()
 432  : done(Standard_False),
 433    //empt, tgte, oppo,
 434    myTolArc(0.0), myTolTang(0.0),
 435    myUVMaxStep(0.0), myFleche(0.0),
 436    myIsStartPnt(Standard_False)
 437    //myU1Start, myV1Start, myU2Start, myV2Start
 438 {
 439 }
 440
 441 //======================================================================
 442 // function: IntPatch_Intersection
 443 //======================================================================
 444 IntPatch_Intersection::IntPatch_Intersection(const Handle(Adaptor3d_HSurface)&  S1,
 445                                              const Handle(Adaptor3d_TopolTool)& D1,
 446                                              const Handle(Adaptor3d_HSurface)&  S2,
 447                                              const Handle(Adaptor3d_TopolTool)& D2,
 448                                              const Standard_Real TolArc,
 449                                              const Standard_Real TolTang)
 450  : done(Standard_False),
 451    //empt, tgte, oppo,
 452    myTolArc(TolArc), myTolTang(TolTang),
 453    myUVMaxStep(0.0), myFleche(0.0),
 454    myIsStartPnt(Standard_False)
 455    //myU1Start, myV1Start, myU2Start, myV2Start
 456 {
 457   if(myTolArc<1e-8) myTolArc=1e-8;
 458   if(myTolTang<1e-8) myTolTang=1e-8;
 459   if(myTolArc>0.5) myTolArc=0.5;
 460   if(myTolTang>0.5) myTolTang=0.5;
 461   Perform(S1,D1,S2,D2,TolArc,TolTang);
 462 }
 463
 464 //======================================================================
 465 // function: IntPatch_Intersection
 466 //======================================================================
 467 IntPatch_Intersection::IntPatch_Intersection(const Handle(Adaptor3d_HSurface)&  S1,
 468                                              const Handle(Adaptor3d_TopolTool)& D1,
 469                                              const Standard_Real TolArc,
 470                                              const Standard_Real TolTang)
 471  : done(Standard_False),
 472    //empt, tgte, oppo,
 473    myTolArc(TolArc), myTolTang(TolTang),
 474    myUVMaxStep(0.0), myFleche(0.0),
 475    myIsStartPnt(Standard_False)
 476    //myU1Start, myV1Start, myU2Start, myV2Start
 477 {
 478   Perform(S1,D1,TolArc,TolTang);
 479 }
 480
 481 //======================================================================
 482 // function: SetTolerances
 483 //======================================================================
 484 void IntPatch_Intersection::SetTolerances(const Standard_Real TolArc,
 485                                           const Standard_Real TolTang,
 486                                           const Standard_Real UVMaxStep,
 487                                           const Standard_Real Fleche)
 488 {
 489   myTolArc     = TolArc;
 490   myTolTang    = TolTang;
 491   myUVMaxStep  = UVMaxStep;
 492   myFleche     = Fleche;
 493   if(myTolArc<1e-8) myTolArc=1e-8;
 494   if(myTolTang<1e-8) myTolTang=1e-8;
 495   if(myTolArc>0.5) myTolArc=0.5;
 496   if(myTolTang>0.5) myTolTang=0.5;
 497   if(myFleche<1.0e-3) myFleche=1e-3;
 498   if(myUVMaxStep<1.0e-3) myUVMaxStep=1e-3;
 499   if(myFleche>10) myFleche=10;
 500   if(myUVMaxStep>0.5) myUVMaxStep=0.5;
 501 }
 502
 503 //======================================================================
 504 // function: Perform
 505 //======================================================================
 506 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  S1,
 507                                     const Handle(Adaptor3d_TopolTool)& D1,
 508                                     const Standard_Real TolArc,
 509                                     const Standard_Real TolTang)
 510 {
 511   myTolArc = TolArc;
 512   myTolTang = TolTang;
 513   if(myFleche == 0.0)  myFleche = 0.01;
 514   if(myUVMaxStep==0.0) myUVMaxStep = 0.01;
 515
 516   done = Standard_True;
 517   spnt.Clear();
 518   slin.Clear();
 519
 520   empt = Standard_True;
 521   tgte = Standard_False;
 522   oppo = Standard_False;
 523
 524   switch (S1->GetType())
 525   {
 526     case GeomAbs_Plane:
 527     case GeomAbs_Cylinder:
 528     case GeomAbs_Sphere:
 529     case GeomAbs_Cone:
 530     case GeomAbs_Torus: break;
 531     default:
 532     {
 533       IntPatch_PrmPrmIntersection interpp;
 534       interpp.Perform(S1,D1,TolArc,TolTang,myFleche,myUVMaxStep);
 535       if (interpp.IsDone())
 536       {
 537         done = Standard_True;
 538         tgte = Standard_False;
 539         empt = interpp.IsEmpty();
 540         const Standard_Integer nblm = interpp.NbLines();
 541         for (Standard_Integer i=1; i<=nblm; i++) slin.Append(interpp.Line(i));
 542       }
 543     }
 544     break;
 545   }
 546 }
 547
 548 /////////////////////////////////////////////////////////////////////////////
 549 //  These several support functions provide methods which can help basic   //
 550 //  algorithm to intersect infinite surfaces of the following types:       //
 551 //                                                                         //
 552 //  a.) SurfaceOfExtrusion;                                                //
 553 //  b.) SurfaceOfRevolution;                                               //
 554 //  c.) OffsetSurface.                                                     //
 555 //                                                                         //
 556 /////////////////////////////////////////////////////////////////////////////
 557 #include <TColgp_Array1OfXYZ.hxx>
 558 #include <TColgp_Array1OfPnt.hxx>
 559 #include <TColgp_SequenceOfPnt.hxx>
 560 #include <Extrema_ExtPS.hxx>
 561 #include <Extrema_POnSurf.hxx>
 562 #include <Geom2d_Curve.hxx>
 563 #include <Geom2dAPI_InterCurveCurve.hxx>
 564 #include <GeomAdaptor.hxx>
 565 #include <GeomAdaptor_HCurve.hxx>
 566 #include <GeomAdaptor_Curve.hxx>
 567 #include <GeomAdaptor_Surface.hxx>
 568 #include <GeomAdaptor_HSurface.hxx>
 569 #include <Geom_Plane.hxx>
 570 #include <ProjLib_ProjectOnPlane.hxx>
 571 #include <GeomProjLib.hxx>
 572 #include <ElCLib.hxx>
 573 #include <Geom_TrimmedCurve.hxx>
 574 #include <Geom_Surface.hxx>
 575 #include <Geom_SurfaceOfLinearExtrusion.hxx>
 576 #include <Geom_OffsetSurface.hxx>
 577 #include <Geom_SurfaceOfRevolution.hxx>
 578 #include <Geom_RectangularTrimmedSurface.hxx>
 579
 580 //===============================================================
 581 //function: FUN_GetMinMaxXYZPnt
 582 //===============================================================
 583 static void FUN_GetMinMaxXYZPnt( const Handle(Adaptor3d_HSurface)& S,
 584                                  gp_Pnt& pMin, gp_Pnt& pMax )
 585 {
 586   const Standard_Real DU = 0.25 * Abs(S->LastUParameter() - S->FirstUParameter());
 587   const Standard_Real DV = 0.25 * Abs(S->LastVParameter() - S->FirstVParameter());
 588   Standard_Real tMinXYZ = RealLast();
 589   Standard_Real tMaxXYZ = -tMinXYZ;
 590   gp_Pnt PUV, ptMax, ptMin;
 591   for(Standard_Real U = S->FirstUParameter(); U <= S->LastUParameter(); U += DU)
 592   {
 593     for(Standard_Real V = S->FirstVParameter(); V <= S->LastVParameter(); V += DV)
 594     {
 595       S->D0(U,V,PUV);
 596       const Standard_Real cXYZ = PUV.XYZ().Modulus();
 597       if(cXYZ > tMaxXYZ) { tMaxXYZ = cXYZ; ptMax = PUV; }
 598       if(cXYZ < tMinXYZ) { tMinXYZ = cXYZ; ptMin = PUV; }
 599     }
 600   }
 601   pMin = ptMin;
 602   pMax = ptMax;
 603 }
 604 //==========================================================================
 605 //function: FUN_TrimInfSurf
 606 //==========================================================================
 607 static void FUN_TrimInfSurf(const gp_Pnt& Pmin,
 608                             const gp_Pnt& Pmax,
 609                             const Handle(Adaptor3d_HSurface)& InfSurf,
 610                             const Standard_Real& AlternativeTrimPrm,
 611                             Handle(Adaptor3d_HSurface)& TrimS)
 612 {
 613   Standard_Real TP = AlternativeTrimPrm;
 614   Extrema_ExtPS ext1(Pmin, InfSurf->Surface(), 1.e-7, 1.e-7);
 615   Extrema_ExtPS ext2(Pmax, InfSurf->Surface(), 1.e-7, 1.e-7);
 616   if(ext1.IsDone() || ext2.IsDone())
 617   {
 618     Standard_Real Umax = -1.e+100, Umin = 1.e+100, Vmax = -1.e+100, Vmin = 1.e+100, cU, cV;
 619     if(ext1.IsDone())
 620     {
 621       for(Standard_Integer i = 1; i <= ext1.NbExt(); i++)
 622       {
 623         const Extrema_POnSurf & pons = ext1.Point(i);
 624         pons.Parameter(cU,cV);
 625         if(cU > Umax) Umax = cU;
 626         if(cU < Umin) Umin = cU;
 627         if(cV > Vmax) Vmax = cV;
 628         if(cV < Vmin) Vmin = cV;
 629       }
 630     }
 631     if(ext2.IsDone())
 632     {
 633       for(Standard_Integer i = 1; i <= ext2.NbExt(); i++)
 634       {
 635         const Extrema_POnSurf & pons = ext2.Point(i);
 636         pons.Parameter(cU,cV);
 637         if(cU > Umax) Umax = cU;
 638         if(cU < Umin) Umin = cU;
 639         if(cV > Vmax) Vmax = cV;
 640         if(cV < Vmin) Vmin = cV;
 641       }
 642     }
 643     TP = Max(Abs(Umin),Max(Abs(Umax),Max(Abs(Vmin),Abs(Vmax))));
 644   }
 645   if(TP == 0.) { TrimS = InfSurf; return; }
 646   else
 647   {
 648     const Standard_Boolean Uinf = Precision::IsNegativeInfinite(InfSurf->FirstUParameter());
 649     const Standard_Boolean Usup = Precision::IsPositiveInfinite(InfSurf->LastUParameter());
 650     const Standard_Boolean Vinf = Precision::IsNegativeInfinite(InfSurf->FirstVParameter());
 651     const Standard_Boolean Vsup = Precision::IsPositiveInfinite(InfSurf->LastVParameter());
 652     Handle(Adaptor3d_HSurface) TmpSS;
 653     Standard_Integer IsTrimed = 0;
 654     const Standard_Real tp = 1000.0 * TP;
 655     if(Vinf && Vsup) { TrimS = InfSurf->VTrim(-tp, tp, 1.0e-7); IsTrimed = 1; }
 656     if(Vinf && !Vsup){ TrimS = InfSurf->VTrim(-tp, InfSurf->LastVParameter(), 1.0e-7); IsTrimed = 1; }
 657     if(!Vinf && Vsup){ TrimS = InfSurf->VTrim(InfSurf->FirstVParameter(), tp, 1.0e-7); IsTrimed = 1; }
 658     if(IsTrimed)
 659     {
 660       TmpSS = TrimS;
 661       if(Uinf && Usup)  TrimS = TmpSS->UTrim(-tp, tp, 1.0e-7);
 662       if(Uinf && !Usup) TrimS = TmpSS->UTrim(-tp, InfSurf->LastUParameter(), 1.0e-7);
 663       if(!Uinf && Usup) TrimS = TmpSS->UTrim(InfSurf->FirstUParameter(), tp, 1.0e-7);
 664     }
 665     else
 666     {
 667       if(Uinf && Usup)  TrimS = InfSurf->UTrim(-tp, tp, 1.0e-7);
 668       if(Uinf && !Usup) TrimS = InfSurf->UTrim(-tp, InfSurf->LastUParameter(), 1.0e-7);
 669       if(!Uinf && Usup) TrimS = InfSurf->UTrim(InfSurf->FirstUParameter(), tp, 1.0e-7);
 670     }
 671   }
 672 }
 673 //================================================================================
 674 //function: FUN_GetUiso
 675 //================================================================================
 676 static void FUN_GetUiso(const Handle(Geom_Surface)& GS,
 677                         const GeomAbs_SurfaceType&  T,
 678                         const Standard_Real&        FirstV,
 679                         const Standard_Real&        LastV,
 680                         const Standard_Boolean&     IsVC,
 681                         const Standard_Boolean&     IsVP,
 682                         const Standard_Real&        U,
 683                         Handle(Geom_Curve)&         I)
 684 {
 685   if(T !=  GeomAbs_OffsetSurface)
 686   {
 687     Handle(Geom_Curve) gc = GS->UIso(U);
 688     if(IsVP && (FirstV == 0.0 && LastV == (2.*M_PI))) I = gc;
 689     else
 690     {
 691       Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstV,LastV);
 692       //szv:I = Handle(Geom_Curve)::DownCast(gtc);
 693       I = gtc;
 694     }
 695   }
 696   else//OffsetSurface
 697   {
 698     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (GS);
 699     const Handle(Geom_Surface) bs = gos->BasisSurface();
 700     Handle(Geom_Curve) gcbs = bs->UIso(U);
 701     GeomAdaptor_Curve gac(gcbs);
 702     const GeomAbs_CurveType GACT = gac.GetType();
 703     if(IsVP || IsVC || GACT == GeomAbs_BSplineCurve || GACT == GeomAbs_BezierCurve || Abs(LastV - FirstV) < 1.e+5)
 704     {
 705       Handle(Geom_Curve) gc = gos->UIso(U);
 706       if(IsVP && (FirstV == 0.0 && LastV == (2*M_PI))) I = gc;
 707       else
 708       {
 709         Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstV,LastV);
 710         //szv:I = Handle(Geom_Curve)::DownCast(gtc);
 711         I = gtc;
 712       }
 713     }
 714     else//Offset Line, Parab, Hyperb
 715     {
 716       Standard_Real VmTr, VMTr;
 717       if(GACT != GeomAbs_Hyperbola)
 718       {
 719         if(FirstV >= 0. && LastV >= 0.){ VmTr = FirstV; VMTr = ((LastV - FirstV) > 1.e+4) ? (FirstV + 1.e+4) : LastV; }
 720         else if(FirstV < 0. && LastV < 0.){ VMTr = LastV; VmTr = ((FirstV - LastV) < -1.e+4) ? (LastV - 1.e+4) : FirstV; }
 721         else { VmTr = (FirstV < -1.e+4) ? -1.e+4 : FirstV; VMTr = (LastV > 1.e+4) ? 1.e+4 : LastV; }
 722       }
 723       else//Hyperbola
 724       {
 725         if(FirstV >= 0. && LastV >= 0.)
 726         {
 727           if(FirstV > 4.) return;
 728           VmTr = FirstV; VMTr = (LastV > 4.) ? 4. : LastV;
 729         }
 730         else if(FirstV < 0. && LastV < 0.)
 731         {
 732           if(LastV < -4.) return;
 733           VMTr = LastV; VmTr = (FirstV < -4.) ? -4. : FirstV;
 734         }
 735         else { VmTr = (FirstV < -4.) ? -4. : FirstV; VMTr = (LastV > 4.) ? 4. : LastV; }
 736       }
 737       //Make trimmed surface
 738       Handle(Geom_RectangularTrimmedSurface) rts = new Geom_RectangularTrimmedSurface(gos,VmTr,VMTr,Standard_True);
 739       I = rts->UIso(U);
 740     }
 741   }
 742 }
 743 //================================================================================
 744 //function: FUN_GetViso
 745 //================================================================================
 746 static void FUN_GetViso(const Handle(Geom_Surface)& GS,
 747                         const GeomAbs_SurfaceType&  T,
 748                         const Standard_Real&        FirstU,
 749                         const Standard_Real&        LastU,
 750                         const Standard_Boolean&     IsUC,
 751                         const Standard_Boolean&     IsUP,
 752                         const Standard_Real&        V,
 753                         Handle(Geom_Curve)&         I)
 754 {
 755   if(T !=  GeomAbs_OffsetSurface)
 756   {
 757     Handle(Geom_Curve) gc = GS->VIso(V);
 758     if(IsUP && (FirstU == 0.0 && LastU == (2*M_PI))) I = gc;
 759     else
 760     {
 761       Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstU,LastU);
 762       //szv:I = Handle(Geom_Curve)::DownCast(gtc);
 763       I = gtc;
 764     }
 765   }
 766   else//OffsetSurface
 767   {
 768     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (GS);
 769     const Handle(Geom_Surface) bs = gos->BasisSurface();
 770     Handle(Geom_Curve) gcbs = bs->VIso(V);
 771     GeomAdaptor_Curve gac(gcbs);
 772     const GeomAbs_CurveType GACT = gac.GetType();
 773     if(IsUP || IsUC || GACT == GeomAbs_BSplineCurve || GACT == GeomAbs_BezierCurve || Abs(LastU - FirstU) < 1.e+5)
 774     {
 775       Handle(Geom_Curve) gc = gos->VIso(V);
 776       if(IsUP && (FirstU == 0.0 && LastU == (2*M_PI))) I = gc;
 777       else
 778       {
 779         Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstU,LastU);
 780         //szv:I = Handle(Geom_Curve)::DownCast(gtc);
 781         I = gtc;
 782       }
 783     }
 784     else//Offset Line, Parab, Hyperb
 785     {
 786       Standard_Real UmTr, UMTr;
 787       if(GACT != GeomAbs_Hyperbola)
 788       {
 789         if(FirstU >= 0. && LastU >= 0.){ UmTr = FirstU; UMTr = ((LastU - FirstU) > 1.e+4) ? (FirstU + 1.e+4) : LastU; }
 790         else if(FirstU < 0. && LastU < 0.){ UMTr = LastU; UmTr = ((FirstU - LastU) < -1.e+4) ? (LastU - 1.e+4) : FirstU; }
 791         else { UmTr = (FirstU < -1.e+4) ? -1.e+4 : FirstU; UMTr = (LastU > 1.e+4) ? 1.e+4 : LastU; }
 792       }
 793       else//Hyperbola
 794       {
 795         if(FirstU >= 0. && LastU >= 0.)
 796         {
 797           if(FirstU > 4.) return;
 798           UmTr = FirstU; UMTr = (LastU > 4.) ? 4. : LastU;
 799         }
 800         else if(FirstU < 0. && LastU < 0.)
 801         {
 802           if(LastU < -4.) return;
 803           UMTr = LastU; UmTr = (FirstU < -4.) ? -4. : FirstU;
 804         }
 805         else { UmTr = (FirstU < -4.) ? -4. : FirstU; UMTr = (LastU > 4.) ? 4. : LastU; }
 806       }
 807       //Make trimmed surface
 808       Handle(Geom_RectangularTrimmedSurface) rts = new Geom_RectangularTrimmedSurface(gos,UmTr,UMTr,Standard_True);
 809       I = rts->VIso(V);
 810     }
 811   }
 812 }
 813 //================================================================================
 814 //function: FUN_PL_Intersection
 815 //================================================================================
 816 static void FUN_PL_Intersection(const Handle(Adaptor3d_HSurface)& S1,
 817                                 const GeomAbs_SurfaceType&        T1,
 818                                 const Handle(Adaptor3d_HSurface)& S2,
 819                                 const GeomAbs_SurfaceType&        T2,
 820                                 Standard_Boolean&                 IsOk,
 821                                 TColgp_SequenceOfPnt&             SP,
 822                                 gp_Vec&                           DV)
 823 {
 824   IsOk = Standard_False;
 825   // 1. Check: both surfaces have U(V)isos - lines.
 826   DV = gp_Vec(0.,0.,1.);
 827   Standard_Boolean isoS1isLine[2] = {0, 0};
 828   Standard_Boolean isoS2isLine[2] = {0, 0};
 829   Handle(Geom_Curve) C1, C2;
 830   const GeomAdaptor_Surface & gas1 = *(GeomAdaptor_Surface*)(&(S1->Surface()));
 831   const GeomAdaptor_Surface & gas2 = *(GeomAdaptor_Surface*)(&(S2->Surface()));
 832   const Handle(Geom_Surface) gs1 = gas1.Surface();
 833   const Handle(Geom_Surface) gs2 = gas2.Surface();
 834   Standard_Real MS1[2], MS2[2];
 835   MS1[0] = 0.5 * (S1->LastUParameter() + S1->FirstUParameter());
 836   MS1[1] = 0.5 * (S1->LastVParameter() + S1->FirstVParameter());
 837   MS2[0] = 0.5 * (S2->LastUParameter() + S2->FirstUParameter());
 838   MS2[1] = 0.5 * (S2->LastVParameter() + S2->FirstVParameter());
 839   if(T1 == GeomAbs_SurfaceOfExtrusion) isoS1isLine[0] = Standard_True;
 840   else if(!S1->IsVPeriodic() && !S1->IsVClosed()) {
 841     if(T1 != GeomAbs_OffsetSurface) C1 = gs1->UIso(MS1[0]);
 842     else {
 843       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs1);
 844       const Handle(Geom_Surface) bs = gos->BasisSurface();
 845       C1 = bs->UIso(MS1[0]);
 846     }
 847     GeomAdaptor_Curve gac(C1);
 848     if(gac.GetType() == GeomAbs_Line) isoS1isLine[0] = Standard_True;
 849   }
 850   if(!S1->IsUPeriodic() && !S1->IsUClosed()) {
 851     if(T1 != GeomAbs_OffsetSurface) C1 = gs1->VIso(MS1[1]);
 852     else {
 853       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs1);
 854       const Handle(Geom_Surface) bs = gos->BasisSurface();
 855       C1 = bs->VIso(MS1[1]);
 856     }
 857     GeomAdaptor_Curve gac(C1);
 858     if(gac.GetType() == GeomAbs_Line) isoS1isLine[1] = Standard_True;
 859   }
 860   if(T2 == GeomAbs_SurfaceOfExtrusion) isoS2isLine[0] = Standard_True;
 861   else if(!S2->IsVPeriodic() && !S2->IsVClosed()) {
 862     if(T2 != GeomAbs_OffsetSurface) C2 = gs2->UIso(MS2[0]);
 863     else {
 864       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs2);
 865       const Handle(Geom_Surface) bs = gos->BasisSurface();
 866       C2 = bs->UIso(MS2[0]);
 867     }
 868     GeomAdaptor_Curve gac(C2);
 869     if(gac.GetType() == GeomAbs_Line) isoS2isLine[0] = Standard_True;
 870   }
 871   if(!S2->IsUPeriodic() && !S2->IsUClosed()) {
 872     if(T2 != GeomAbs_OffsetSurface) C2 = gs2->VIso(MS2[1]);
 873     else {
 874       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs2);
 875       const Handle(Geom_Surface) bs = gos->BasisSurface();
 876       C2 = bs->VIso(MS2[1]);
 877     }
 878     GeomAdaptor_Curve gac(C2);
 879     if(gac.GetType() == GeomAbs_Line) isoS2isLine[1] = Standard_True;
 880   }
 881   Standard_Boolean IsBothLines = ((isoS1isLine[0] || isoS1isLine[1]) &&
 882                                   (isoS2isLine[0] || isoS2isLine[1]));
 883   if(!IsBothLines){
 884     return;
 885   }
 886   // 2. Check: Uiso lines of both surfaces are collinear.
 887   gp_Pnt puvS1, puvS2;
 888   gp_Vec derS1[2], derS2[2];
 889   S1->D1(MS1[0], MS1[1], puvS1, derS1[0], derS1[1]);
 890   S2->D1(MS2[0], MS2[1], puvS2, derS2[0], derS2[1]);
 891   C1.Nullify(); C2.Nullify();
 892   Standard_Integer iso = 0;
 893   if(isoS1isLine[0] && isoS2isLine[0] &&
 894      derS1[1].IsParallel(derS2[1],Precision::Angular())) {
 895     iso = 1;
 896     FUN_GetViso(gs1,T1,S1->FirstUParameter(),S1->LastUParameter(),
 897                 S1->IsUClosed(),S1->IsUPeriodic(),MS1[1],C1);
 898     FUN_GetViso(gs2,T2,S2->FirstUParameter(),S2->LastUParameter(),
 899                 S2->IsUClosed(),S2->IsUPeriodic(),MS2[1],C2);
 900   }
 901   else if(isoS1isLine[0] && isoS2isLine[1] &&
 902           derS1[1].IsParallel(derS2[0],Precision::Angular())) {
 903     iso = 1;
 904     FUN_GetViso(gs1,T1,S1->FirstUParameter(),S1->LastUParameter(),
 905                 S1->IsUClosed(),S1->IsUPeriodic(),MS1[1],C1);
 906     FUN_GetUiso(gs2,T2,S2->FirstVParameter(),S2->LastVParameter(),
 907                 S2->IsVClosed(),S2->IsVPeriodic(),MS2[0],C2);
 908   }
 909   else if(isoS1isLine[1] && isoS2isLine[0] &&
 910           derS1[0].IsParallel(derS2[1],Precision::Angular())) {
 911     iso = 0;
 912     FUN_GetUiso(gs1,T1,S1->FirstVParameter(),S1->LastVParameter(),
 913                 S1->IsVClosed(),S1->IsVPeriodic(),MS1[0],C1);
 914     FUN_GetViso(gs2,T2,S2->FirstUParameter(),S2->LastUParameter(),
 915                 S2->IsUClosed(),S2->IsUPeriodic(),MS2[1],C2);
 916   }
 917   else if(isoS1isLine[1] && isoS2isLine[1] &&
 918           derS1[0].IsParallel(derS2[0],Precision::Angular())) {
 919     iso = 0;
 920     FUN_GetUiso(gs1,T1,S1->FirstVParameter(),S1->LastVParameter(),
 921                 S1->IsVClosed(),S1->IsVPeriodic(),MS1[0],C1);
 922     FUN_GetUiso(gs2,T2,S2->FirstVParameter(),S2->LastVParameter(),
 923                 S2->IsVClosed(),S2->IsVPeriodic(),MS2[0],C2);
 924   }
 925   else {
 926     IsOk = Standard_False;
 927     return;
 928   }
 929   IsOk = Standard_True;
 930   // 3. Make intersections of V(U)isos
 931   if(C1.IsNull() || C2.IsNull()) return;
 932   DV = derS1[iso];
 933   Handle(Geom_Plane) GPln = new Geom_Plane(gp_Pln(puvS1,gp_Dir(DV)));
 934   Handle(Geom_Curve) C1Prj =
 935     GeomProjLib::ProjectOnPlane(C1,GPln,gp_Dir(DV),Standard_True);
 936   Handle(Geom_Curve) C2Prj =
 937     GeomProjLib::ProjectOnPlane(C2,GPln,gp_Dir(DV),Standard_True);
 938   if(C1Prj.IsNull() || C2Prj.IsNull()) return;
 939   Handle(Geom2d_Curve) C1Prj2d =
 940     GeomProjLib::Curve2d(C1Prj,Handle(Geom_Surface)::DownCast (GPln));
 941   Handle(Geom2d_Curve) C2Prj2d =
 942     GeomProjLib::Curve2d(C2Prj,Handle(Geom_Surface)::DownCast (GPln));
 943   Geom2dAPI_InterCurveCurve ICC(C1Prj2d,C2Prj2d,1.0e-7);
 944   if(ICC.NbPoints() > 0 )
 945   {
 946     for(Standard_Integer ip = 1; ip <= ICC.NbPoints(); ip++)
 947     {
 948       gp_Pnt2d P = ICC.Point(ip);
 949       gp_Pnt P3d = ElCLib::To3d(gp_Ax2(puvS1,gp_Dir(DV)),P);
 950       SP.Append(P3d);
 951     }
 952   }
 953 }
 954 //================================================================================
 955 //function: FUN_NewFirstLast
 956 //================================================================================
 957 static void FUN_NewFirstLast(const GeomAbs_CurveType& ga_ct,
 958                              const Standard_Real&     Fst,
 959                              const Standard_Real&     Lst,
 960                              const Standard_Real&     TrVal,
 961                              Standard_Real&           NewFst,
 962                              Standard_Real&           NewLst,
 963                              Standard_Boolean&        NeedTr)
 964 {
 965   NewFst = Fst; NewLst = Lst; NeedTr = Standard_False;
 966   switch (ga_ct)
 967   {
 968     case GeomAbs_Line:
 969     case GeomAbs_Parabola:
 970     {
 971       if(Abs(Lst - Fst) > TrVal)
 972       {
 973         if(Fst >= 0. && Lst >= 0.)
 974         {
 975           NewFst = Fst;
 976           NewLst = ((Fst + TrVal) < Lst) ? (Fst + TrVal) : Lst;
 977         }
 978         if(Fst < 0. && Lst < 0.)
 979         {
 980           NewLst = Lst;
 981           NewFst = ((Lst - TrVal) > Fst) ? (Lst - TrVal) : Fst;
 982         }
 983         else
 984         {
 985           NewFst = (Fst < -TrVal) ? -TrVal : Fst;
 986           NewLst = (Lst > TrVal) ? TrVal : Lst;
 987         }
 988         NeedTr = Standard_True;
 989       }
 990       break;
 991     }
 992     case GeomAbs_Hyperbola:
 993     {
 994       if(Abs(Lst - Fst) > 10.)
 995       {
 996         if(Fst >= 0. && Lst >= 0.)
 997         {
 998           if(Fst > 4.) return;
 999           NewFst = Fst;
1000           NewLst = (Lst > 4.) ? 4. : Lst;
1001         }
1002         if(Fst < 0. && Lst < 0.)
1003         {
1004           if(Lst < -4.) return;
1005           NewLst = Lst;
1006           NewFst = (Fst < -4.) ? -4. : Fst;
1007         }
1008         else
1009         {
1010           NewFst = (Fst < -4.) ? -4. : Fst;
1011           NewLst = (Lst > 4.) ? 4. : Lst;
1012         }
1013         NeedTr = Standard_True;
1014       }
1015       break;
1016     }
1017     default:
1018     break;
1019   }
1020 }
1021 //================================================================================
1022 //function: FUN_TrimBothSurf
1023 //================================================================================
1024 static void FUN_TrimBothSurf(const Handle(Adaptor3d_HSurface)& S1,
1025                              const GeomAbs_SurfaceType&        T1,
1026                              const Handle(Adaptor3d_HSurface)& S2,
1027                              const GeomAbs_SurfaceType&        T2,
1028                              const Standard_Real&              TV,
1029                              Handle(Adaptor3d_HSurface)&       NS1,
1030                              Handle(Adaptor3d_HSurface)&       NS2)
1031 {
1032   const GeomAdaptor_Surface & gas1 = *(GeomAdaptor_Surface*)(&(S1->Surface()));
1033   const GeomAdaptor_Surface & gas2 = *(GeomAdaptor_Surface*)(&(S2->Surface()));
1034   const Handle(Geom_Surface) gs1 = gas1.Surface();
1035   const Handle(Geom_Surface) gs2 = gas2.Surface();
1036   const Standard_Real UM1 = 0.5 * (S1->LastUParameter() + S1->FirstUParameter());
1037   const Standard_Real UM2 = 0.5 * (S2->LastUParameter() + S2->FirstUParameter());
1038   const Standard_Real VM1 = 0.5 * (S1->LastVParameter() + S1->FirstVParameter());
1039   const Standard_Real VM2 = 0.5 * (S2->LastVParameter() + S2->FirstVParameter());
1040   Handle(Geom_Curve) visoS1, visoS2, uisoS1, uisoS2;
1041   if(T1 != GeomAbs_OffsetSurface){ visoS1 = gs1->VIso(VM1); uisoS1 = gs1->UIso(UM1); }
1042   else
1043   {
1044     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs1);
1045     const Handle(Geom_Surface) bs = gos->BasisSurface();
1046     visoS1 = bs->VIso(VM1); uisoS1 = bs->UIso(UM1);
1047   }
1048   if(T2 != GeomAbs_OffsetSurface){ visoS2 = gs2->VIso(VM2); uisoS2 = gs2->UIso(UM2); }
1049   else
1050   {
1051     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs2);
1052     const Handle(Geom_Surface) bs = gos->BasisSurface();
1053     visoS2 = bs->VIso(VM2); uisoS2 = bs->UIso(UM2);
1054   }
1055   if(uisoS1.IsNull() || uisoS2.IsNull() || visoS1.IsNull() || visoS2.IsNull()){ NS1 = S1; NS2 = S2; return; }
1056   GeomAdaptor_Curve gau1(uisoS1);
1057   GeomAdaptor_Curve gav1(visoS1);
1058   GeomAdaptor_Curve gau2(uisoS2);
1059   GeomAdaptor_Curve gav2(visoS2);
1060   GeomAbs_CurveType GA_U1 = gau1.GetType();
1061   GeomAbs_CurveType GA_V1 = gav1.GetType();
1062   GeomAbs_CurveType GA_U2 = gau2.GetType();
1063   GeomAbs_CurveType GA_V2 = gav2.GetType();
1064   Standard_Boolean TrmU1 = Standard_False;
1065   Standard_Boolean TrmV1 = Standard_False;
1066   Standard_Boolean TrmU2 = Standard_False;
1067   Standard_Boolean TrmV2 = Standard_False;
1068   Standard_Real V1S1,V2S1,U1S1,U2S1, V1S2,V2S2,U1S2,U2S2;
1069   FUN_NewFirstLast(GA_U1,S1->FirstVParameter(),S1->LastVParameter(),TV,V1S1,V2S1,TrmV1);
1070   FUN_NewFirstLast(GA_V1,S1->FirstUParameter(),S1->LastUParameter(),TV,U1S1,U2S1,TrmU1);
1071   FUN_NewFirstLast(GA_U2,S2->FirstVParameter(),S2->LastVParameter(),TV,V1S2,V2S2,TrmV2);
1072   FUN_NewFirstLast(GA_V2,S2->FirstUParameter(),S2->LastUParameter(),TV,U1S2,U2S2,TrmU2);
1073   if(TrmV1) NS1 = S1->VTrim(V1S1, V2S1, 1.0e-7);
1074   if(TrmV2) NS2 = S2->VTrim(V1S2, V2S2, 1.0e-7);
1075   if(TrmU1)
1076   {
1077     if(TrmV1)
1078     {
1079       Handle(Adaptor3d_HSurface) TS = NS1;
1080       NS1 = TS->UTrim(U1S1, U2S1, 1.0e-7);
1081     }
1082     else NS1 = S1->UTrim(U1S1, U2S1, 1.0e-7);
1083   }
1084   if(TrmU2)
1085   {
1086     if(TrmV2)
1087     {
1088       Handle(Adaptor3d_HSurface) TS = NS2;
1089       NS2 = TS->UTrim(U1S2, U2S2, 1.0e-7);
1090     }
1091     else NS2 = S2->UTrim(U1S2, U2S2, 1.0e-7);
1092   }
1093 }
1094
1095 //=======================================================================
1096 //function : Perform
1097 //purpose  :
1098 //=======================================================================
1099 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  theS1,
1100                                     const Handle(Adaptor3d_TopolTool)& theD1,
1101                                     const Handle(Adaptor3d_HSurface)&  theS2,
1102                                     const Handle(Adaptor3d_TopolTool)& theD2,
1103                                     const Standard_Real TolArc,
1104                                     const Standard_Real TolTang,
1105                                     const Standard_Boolean isGeomInt,
1106                                     const Standard_Boolean theIsReqToKeepRLine)
1107 {
1108   myTolArc = TolArc;
1109   myTolTang = TolTang;
1110   if(myFleche <= Precision::PConfusion())
1111     myFleche = 0.01;
1112   if(myUVMaxStep <= Precision::PConfusion())
1113     myUVMaxStep = 0.01;
1114
1115   done = Standard_False;
1116   spnt.Clear();
1117   slin.Clear();
1118   empt = Standard_True;
1119   tgte = Standard_False;
1120   oppo = Standard_False;
1121
1122   GeomAbs_SurfaceType typs1 = theS1->GetType();
1123   GeomAbs_SurfaceType typs2 = theS2->GetType();
1124
1125   //treatment of the cases with cone or torus
1126   Standard_Boolean TreatAsBiParametric = Standard_False;
1127   Standard_Integer bGeomGeom = 0;
1128   //
1129   if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone ||
1130       typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1131     gp_Ax1 aCTAx, aGeomAx;
1132     GeomAbs_SurfaceType aCTType;
1133     Standard_Boolean bToCheck;
1134     //
1135     const Handle(Adaptor3d_HSurface)& aCTSurf =
1136       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS1 : theS2;
1137     const Handle(Adaptor3d_HSurface)& aGeomSurf =
1138       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS2 : theS1;
1139     //
1140     aCTType = aCTSurf->GetType();
1141     bToCheck = Standard_False;
1142     //
1143     if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone) {
1144       const gp_Cone aCon1 = (aCTType == GeomAbs_Cone) ?
1145         aCTSurf->Cone() : aGeomSurf->Cone();
1146       Standard_Real a1 = Abs(aCon1.SemiAngle());
1147       bToCheck = (a1 < 0.02) || (a1 > 1.55);
1148       //
1149       if (typs1 == typs2) {
1150         const gp_Cone aCon2 = aGeomSurf->Cone();
1151         Standard_Real a2 = Abs(aCon2.SemiAngle());
1152         bToCheck = bToCheck || (a2 < 0.02) || (a2 > 1.55);
1153         //
1154         if (a1 > 1.55 && a2 > 1.55) {//quasi-planes: if same domain, treat as canonic
1155           const gp_Ax1 A1 = aCon1.Axis(), A2 = aCon2.Axis();
1156           if (A1.IsParallel(A2,Precision::Angular())) {
1157             const gp_Pnt Apex1 = aCon1.Apex(), Apex2 = aCon2.Apex();
1158             const gp_Pln Plan1( Apex1, A1.Direction() );
1159             if (Plan1.Distance( Apex2 ) <= Precision::Confusion()) {
1160               bToCheck = Standard_False;
1161             }
1162           }
1163         }
1164       }
1165       //
1166       TreatAsBiParametric = bToCheck;
1167       if (aCTType == GeomAbs_Cone) {
1168         aCTAx = aCon1.Axis();
1169       }
1170     }
1171     //
1172     if (typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1173       const gp_Torus aTor1 = (aCTType == GeomAbs_Torus) ?
1174         aCTSurf->Torus() : aGeomSurf->Torus();
1175       bToCheck = aTor1.MajorRadius() > aTor1.MinorRadius();
1176       if (typs1 == typs2) {
1177         const gp_Torus aTor2 = aGeomSurf->Torus();
1178         bToCheck = aTor2.MajorRadius() > aTor2.MinorRadius();
1179       }
1180       //
1181       if (aCTType == GeomAbs_Torus) {
1182         aCTAx = aTor1.Axis();
1183       }
1184     }
1185     //
1186     if (bToCheck) {
1187       const gp_Lin aL1(aCTAx);
1188       //
1189       switch (aGeomSurf->GetType()) {
1190       case GeomAbs_Plane: {
1191         aGeomAx = aGeomSurf->Plane().Axis();
1192         if (aCTType == GeomAbs_Cone) {
1193           bGeomGeom = 1;
1194           if (Abs(aCTSurf->Cone().SemiAngle()) < 0.02) {
1195             Standard_Real ps = Abs(aCTAx.Direction().Dot(aGeomAx.Direction()));
1196             if(ps < 0.015) {
1197               bGeomGeom = 0;
1198             }
1199           }
1200         }
1201         else {
1202           if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) ||
1203               (aCTAx.IsNormal(aGeomAx, Precision::Angular()) &&
1204                (aGeomSurf->Plane().Distance(aCTAx.Location()) < Precision::Confusion()))) {
1205             bGeomGeom = 1;
1206           }
1207         }
1208         bToCheck = Standard_False;
1209         break;
1210       }
1211       case GeomAbs_Sphere: {
1212         if (aL1.Distance(aGeomSurf->Sphere().Location()) < Precision::Confusion()) {
1213           bGeomGeom = 1;
1214         }
1215         bToCheck = Standard_False;
1216         break;
1217       }
1218       case GeomAbs_Cylinder:
1219         aGeomAx = aGeomSurf->Cylinder().Axis();
1220         break;
1221       case GeomAbs_Cone:
1222         aGeomAx = aGeomSurf->Cone().Axis();
1223         break;
1224       case GeomAbs_Torus:
1225         aGeomAx = aGeomSurf->Torus().Axis();
1226         break;
1227       default:
1228         bToCheck = Standard_False;
1229         break;
1230       }
1231       //
1232       if (bToCheck) {
1233         if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) &&
1234             (aL1.Distance(aGeomAx.Location()) <= Precision::Confusion())) {
1235           bGeomGeom = 1;
1236         }
1237       }
1238       //
1239       if (bGeomGeom == 1) {
1240         TreatAsBiParametric = Standard_False;
1241       }
1242     }
1243   }
1244   //
1245
1246   if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite()) {
1247     TreatAsBiParametric= Standard_False;
1248   }
1249
1250 //  Modified by skv - Mon Sep 26 14:58:30 2005 Begin
1251 //   if(TreatAsBiParametric) { typs1 = typs2 = GeomAbs_BezierSurface; }
1252   if(TreatAsBiParametric)
1253   {
1254     if (typs1 == GeomAbs_Cone && typs2 == GeomAbs_Plane)
1255       typs1 = GeomAbs_BezierSurface; // Using Imp-Prm Intersector
1256     else if (typs1 == GeomAbs_Plane && typs2 == GeomAbs_Cone)
1257       typs2 = GeomAbs_BezierSurface; // Using Imp-Prm Intersector
1258     else {
1259       // Using Prm-Prm Intersector
1260       typs1 = GeomAbs_BezierSurface;
1261       typs2 = GeomAbs_BezierSurface;
1262     }
1263   }
1264 //  Modified by skv - Mon Sep 26 14:58:30 2005 End
1265
1266   // Surface type definition
1267   Standard_Integer ts1 = 0;
1268   switch (typs1)
1269   {
1270     case GeomAbs_Plane:
1271     case GeomAbs_Cylinder:
1272     case GeomAbs_Sphere:
1273     case GeomAbs_Cone: ts1 = 1; break;
1274     case GeomAbs_Torus: ts1 = bGeomGeom; break;
1275     default: break;
1276   }
1277
1278   Standard_Integer ts2 = 0;
1279   switch (typs2)
1280   {
1281     case GeomAbs_Plane:
1282     case GeomAbs_Cylinder:
1283     case GeomAbs_Sphere:
1284     case GeomAbs_Cone: ts2 = 1; break;
1285     case GeomAbs_Torus: ts2 = bGeomGeom; break;
1286     default: break;
1287   }
1288   //
1289   // treatment of the cases with torus and any other geom surface
1290   //
1291   // Possible intersection types: 1. ts1 == ts2 == 1 <Geom-Geom>
1292   //                              2. ts1 != ts2      <Geom-Param>
1293   //                              3. ts1 == ts2 == 0 <Param-Param>
1294
1295   // Geom - Geom
1296   if(ts1 == ts2 && ts1 == 1)
1297   {
1298     const Standard_Boolean RestrictLine = Standard_True;
1299     IntSurf_ListOfPntOn2S ListOfPnts;
1300     ListOfPnts.Clear();
1301     if(isGeomInt)
1302     {
1303       if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite())
1304       {
1305         GeomGeomPerfom( theS1, theD1, theS2, theD2, TolArc,
1306                         TolTang, ListOfPnts, RestrictLine,
1307                         typs1, typs2, theIsReqToKeepRLine);
1308       }
1309       else
1310       {
1311         GeomGeomPerfomTrimSurf( theS1, theD1, theS2, theD2,
1312                                 TolArc, TolTang, ListOfPnts, RestrictLine,
1313                                 typs1, typs2, theIsReqToKeepRLine);
1314       }
1315     }
1316     else
1317     {
1318       ParamParamPerfom(theS1, theD1, theS2, theD2,
1319               TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1320     }
1321   }
1322
1323   // Geom - Param
1324   if(ts1 != ts2)
1325   {
1326     GeomParamPerfom(theS1, theD1, theS2, theD2, ts1 == 0, typs1, typs2);
1327   }
1328
1329   // Param - Param
1330   if(ts1 == ts2 && ts1 == 0)
1331   {
1332     const Standard_Boolean RestrictLine = Standard_True;
1333     IntSurf_ListOfPntOn2S ListOfPnts;
1334     ListOfPnts.Clear();
1335
1336     ParamParamPerfom(theS1, theD1, theS2, theD2, TolArc,
1337                         TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1338   }
1339 }
1340
1341 //=======================================================================
1342 //function : Perform
1343 //purpose  :
1344 //=======================================================================
1345 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  theS1,
1346                                     const Handle(Adaptor3d_TopolTool)& theD1,
1347                                     const Handle(Adaptor3d_HSurface)&  theS2,
1348                                     const Handle(Adaptor3d_TopolTool)& theD2,
1349                                     const Standard_Real TolArc,
1350                                     const Standard_Real TolTang,
1351                                     IntSurf_ListOfPntOn2S& ListOfPnts,
1352                                     const Standard_Boolean RestrictLine,
1353                                     const Standard_Boolean isGeomInt)
1354 {
1355   myTolArc = TolArc;
1356   myTolTang = TolTang;
1357   if(myFleche <= Precision::PConfusion())
1358     myFleche = 0.01;
1359   if(myUVMaxStep <= Precision::PConfusion())
1360     myUVMaxStep = 0.01;
1361
1362   done = Standard_False;
1363   spnt.Clear();
1364   slin.Clear();
1365   empt = Standard_True;
1366   tgte = Standard_False;
1367   oppo = Standard_False;
1368
1369   GeomAbs_SurfaceType typs1 = theS1->GetType();
1370   GeomAbs_SurfaceType typs2 = theS2->GetType();
1371   //
1372   //treatment of the cases with cone or torus
1373   Standard_Boolean TreatAsBiParametric = Standard_False;
1374   Standard_Integer bGeomGeom = 0;
1375   //
1376   if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone ||
1377       typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1378     gp_Ax1 aCTAx, aGeomAx;
1379     GeomAbs_SurfaceType aCTType;
1380     Standard_Boolean bToCheck;
1381     //
1382     const Handle(Adaptor3d_HSurface)& aCTSurf =
1383       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS1 : theS2;
1384     const Handle(Adaptor3d_HSurface)& aGeomSurf =
1385       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS2 : theS1;
1386     //
1387     aCTType = aCTSurf->GetType();
1388     bToCheck = Standard_False;
1389     //
1390     if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone) {
1391       const gp_Cone aCon1 = (aCTType == GeomAbs_Cone) ?
1392         aCTSurf->Cone() : aGeomSurf->Cone();
1393       Standard_Real a1 = Abs(aCon1.SemiAngle());
1394       bToCheck = (a1 < 0.02) || (a1 > 1.55);
1395       //
1396       if (typs1 == typs2) {
1397         const gp_Cone aCon2 = aGeomSurf->Cone();
1398         Standard_Real a2 = Abs(aCon2.SemiAngle());
1399         bToCheck = bToCheck || (a2 < 0.02) || (a2 > 1.55);
1400         //
1401         if (a1 > 1.55 && a2 > 1.55) {//quasi-planes: if same domain, treat as canonic
1402           const gp_Ax1 A1 = aCon1.Axis(), A2 = aCon2.Axis();
1403           if (A1.IsParallel(A2,Precision::Angular())) {
1404             const gp_Pnt Apex1 = aCon1.Apex(), Apex2 = aCon2.Apex();
1405             const gp_Pln Plan1( Apex1, A1.Direction() );
1406             if (Plan1.Distance( Apex2 ) <= Precision::Confusion()) {
1407               bToCheck = Standard_False;
1408             }
1409           }
1410         }
1411       }
1412       //
1413       TreatAsBiParametric = bToCheck;
1414       if (aCTType == GeomAbs_Cone) {
1415         aCTAx = aCon1.Axis();
1416       }
1417     }
1418     //
1419     if (typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1420       const gp_Torus aTor1 = (aCTType == GeomAbs_Torus) ?
1421         aCTSurf->Torus() : aGeomSurf->Torus();
1422       bToCheck = aTor1.MajorRadius() > aTor1.MinorRadius();
1423       if (typs1 == typs2) {
1424         const gp_Torus aTor2 = aGeomSurf->Torus();
1425         bToCheck = aTor2.MajorRadius() > aTor2.MinorRadius();
1426       }
1427       //
1428       if (aCTType == GeomAbs_Torus) {
1429         aCTAx = aTor1.Axis();
1430       }
1431     }
1432     //
1433     if (bToCheck) {
1434       const gp_Lin aL1(aCTAx);
1435       //
1436       switch (aGeomSurf->GetType()) {
1437       case GeomAbs_Plane: {
1438         aGeomAx = aGeomSurf->Plane().Axis();
1439         if (aCTType == GeomAbs_Cone) {
1440           bGeomGeom = 1;
1441           if (Abs(aCTSurf->Cone().SemiAngle()) < 0.02) {
1442             Standard_Real ps = Abs(aCTAx.Direction().Dot(aGeomAx.Direction()));
1443             if(ps < 0.015) {
1444               bGeomGeom = 0;
1445             }
1446           }
1447         }
1448         else {
1449           if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) ||
1450               (aCTAx.IsNormal(aGeomAx, Precision::Angular()) &&
1451                (aGeomSurf->Plane().Distance(aCTAx.Location()) < Precision::Confusion()))) {
1452             bGeomGeom = 1;
1453           }
1454         }
1455         bToCheck = Standard_False;
1456         break;
1457       }
1458       case GeomAbs_Sphere: {
1459         if (aL1.Distance(aGeomSurf->Sphere().Location()) < Precision::Confusion()) {
1460           bGeomGeom = 1;
1461         }
1462         bToCheck = Standard_False;
1463         break;
1464       }
1465       case GeomAbs_Cylinder:
1466         aGeomAx = aGeomSurf->Cylinder().Axis();
1467         break;
1468       case GeomAbs_Cone:
1469         aGeomAx = aGeomSurf->Cone().Axis();
1470         break;
1471       case GeomAbs_Torus:
1472         aGeomAx = aGeomSurf->Torus().Axis();
1473         break;
1474       default:
1475         bToCheck = Standard_False;
1476         break;
1477       }
1478       //
1479       if (bToCheck) {
1480         if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) &&
1481             (aL1.Distance(aGeomAx.Location()) <= Precision::Confusion())) {
1482           bGeomGeom = 1;
1483         }
1484       }
1485       //
1486       if (bGeomGeom == 1) {
1487         TreatAsBiParametric = Standard_False;
1488       }
1489     }
1490   }
1491   //
1492
1493   if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite()) {
1494     TreatAsBiParametric= Standard_False;
1495   }
1496
1497   if(TreatAsBiParametric)
1498   {
1499     // Using Prm-Prm Intersector
1500     typs1 = GeomAbs_BezierSurface;
1501     typs2 = GeomAbs_BezierSurface;
1502   }
1503
1504   // Surface type definition
1505   Standard_Integer ts1 = 0;
1506   switch (typs1)
1507   {
1508     case GeomAbs_Plane:
1509     case GeomAbs_Cylinder:
1510     case GeomAbs_Sphere:
1511     case GeomAbs_Cone: ts1 = 1; break;
1512     case GeomAbs_Torus: ts1 = bGeomGeom; break;
1513     default: break;
1514   }
1515
1516   Standard_Integer ts2 = 0;
1517   switch (typs2)
1518   {
1519     case GeomAbs_Plane:
1520     case GeomAbs_Cylinder:
1521     case GeomAbs_Sphere:
1522     case GeomAbs_Cone: ts2 = 1; break;
1523     case GeomAbs_Torus: ts2 = bGeomGeom; break;
1524     default: break;
1525   }
1526   //
1527   // Possible intersection types: 1. ts1 == ts2 == 1 <Geom-Geom>
1528   //                              2. ts1 != ts2      <Geom-Param>
1529   //                              3. ts1 == ts2 == 0 <Param-Param>
1530
1531   if(!isGeomInt)
1532   {
1533     ParamParamPerfom(theS1, theD1, theS2, theD2,
1534                 TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1535   }
1536   else if(ts1 != ts2)
1537   {
1538     GeomParamPerfom(theS1, theD1, theS2, theD2, ts1 == 0, typs1, typs2);
1539   }
1540   else if (ts1 == 0)
1541   {
1542     ParamParamPerfom(theS1, theD1, theS2, theD2,
1543                 TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1544   }
1545   else if(ts1 == 1)
1546   {
1547     if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite())
1548     {
1549       GeomGeomPerfom(theS1, theD1, theS2, theD2, TolArc,
1550                       TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1551     }
1552     else
1553     {
1554       GeomGeomPerfomTrimSurf(theS1, theD1, theS2, theD2,
1555               TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1556     }
1557   }
1558 }
1559
1560 //=======================================================================
1561 //function : ParamParamPerfom
1562 //purpose  :
1563 //=======================================================================
1564 void IntPatch_Intersection::ParamParamPerfom(const Handle(Adaptor3d_HSurface)&  theS1,
1565                                              const Handle(Adaptor3d_TopolTool)& theD1,
1566                                              const Handle(Adaptor3d_HSurface)&  theS2,
1567                                              const Handle(Adaptor3d_TopolTool)& theD2,
1568                                              const Standard_Real TolArc,
1569                                              const Standard_Real TolTang,
1570                                              IntSurf_ListOfPntOn2S& ListOfPnts,
1571                                              const Standard_Boolean RestrictLine,
1572                                              const GeomAbs_SurfaceType typs1,
1573                                              const GeomAbs_SurfaceType typs2)
1574 {
1575   IntPatch_PrmPrmIntersection interpp;
1576   if(!theD1->DomainIsInfinite() && !theD2->DomainIsInfinite())
1577   {
1578     Standard_Boolean ClearFlag = Standard_True;
1579     if(!ListOfPnts.IsEmpty())
1580     {
1581       interpp.Perform(theS1,theD1,theS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep, ListOfPnts, RestrictLine);
1582       ClearFlag = Standard_False;
1583     }
1584     interpp.Perform(theS1,theD1,theS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep,ClearFlag);   //double call!!!!!!!
1585   }
1586   else if((theD1->DomainIsInfinite()) ^ (theD2->DomainIsInfinite()))
1587   {
1588     gp_Pnt pMaxXYZ, pMinXYZ;
1589     if(theD1->DomainIsInfinite())
1590     {
1591       FUN_GetMinMaxXYZPnt( theS2, pMinXYZ, pMaxXYZ );
1592       const Standard_Real MU = Max(Abs(theS2->FirstUParameter()),Abs(theS2->LastUParameter()));
1593       const Standard_Real MV = Max(Abs(theS2->FirstVParameter()),Abs(theS2->LastVParameter()));
1594       const Standard_Real AP = Max(MU, MV);
1595       Handle(Adaptor3d_HSurface) SS;
1596       FUN_TrimInfSurf(pMinXYZ, pMaxXYZ, theS1, AP, SS);
1597       interpp.Perform(SS,theD1,theS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep);
1598     }
1599     else
1600     {
1601       FUN_GetMinMaxXYZPnt( theS1, pMinXYZ, pMaxXYZ );
1602       const Standard_Real MU = Max(Abs(theS1->FirstUParameter()),Abs(theS1->LastUParameter()));
1603       const Standard_Real MV = Max(Abs(theS1->FirstVParameter()),Abs(theS1->LastVParameter()));
1604       const Standard_Real AP = Max(MU, MV);
1605       Handle(Adaptor3d_HSurface) SS;
1606       FUN_TrimInfSurf(pMinXYZ, pMaxXYZ, theS2, AP, SS);
1607       interpp.Perform(theS1, theD1, SS, theD2,TolArc,TolTang,myFleche,myUVMaxStep);
1608     }
1609   }//(theD1->DomainIsInfinite()) ^ (theD2->DomainIsInfinite())
1610   else
1611   {
1612     if(typs1 == GeomAbs_OtherSurface || typs2 == GeomAbs_OtherSurface)
1613     {
1614       done = Standard_False;
1615       return;
1616     }
1617
1618     Standard_Boolean IsPLInt = Standard_False;
1619     TColgp_SequenceOfPnt sop;
1620     gp_Vec v;
1621     FUN_PL_Intersection(theS1,typs1,theS2,typs2,IsPLInt,sop,v);
1622
1623     if(IsPLInt)
1624     {
1625       if(sop.Length() > 0)
1626       {
1627         for(Standard_Integer ip = 1; ip <= sop.Length(); ip++)
1628         {
1629           gp_Lin lin(sop.Value(ip),gp_Dir(v));
1630           Handle(IntPatch_GLine) gl = new IntPatch_GLine(lin,Standard_False);
1631           slin.Append(Handle(IntPatch_Line)::DownCast (gl));
1632         }
1633
1634         done = Standard_True;
1635       }
1636       else
1637         done = Standard_False;
1638
1639       return;
1640     }// 'COLLINEAR LINES'
1641     else
1642     {
1643       Handle(Adaptor3d_HSurface) nS1 = theS1;
1644       Handle(Adaptor3d_HSurface) nS2 = theS2;
1645       FUN_TrimBothSurf(theS1,typs1,theS2,typs2,1.e+8,nS1,nS2);
1646       interpp.Perform(nS1,theD1,nS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep);
1647     }// 'NON - COLLINEAR LINES'
1648   }// both domains are infinite
1649
1650   if (interpp.IsDone())
1651   {
1652     done = Standard_True;
1653     tgte = Standard_False;
1654     empt = interpp.IsEmpty();
1655
1656     for(Standard_Integer i = 1; i <= interpp.NbLines(); i++)
1657     {
1658       if(interpp.Line(i)->ArcType() != IntPatch_Walking)
1659         slin.Append(interpp.Line(i));
1660     }
1661
1662     for (Standard_Integer i = 1; i <= interpp.NbLines(); i++)
1663     {
1664       if(interpp.Line(i)->ArcType() == IntPatch_Walking)
1665         slin.Append(interpp.Line(i));
1666     }
1667   }
1668 }
1669
1670 //=======================================================================
1671 ////function : GeomGeomPerfom
1672 //purpose  :
1673 //=======================================================================
1674 void IntPatch_Intersection::GeomGeomPerfom(const Handle(Adaptor3d_HSurface)& theS1,
1675                                            const Handle(Adaptor3d_TopolTool)& theD1,
1676                                            const Handle(Adaptor3d_HSurface)& theS2,
1677                                            const Handle(Adaptor3d_TopolTool)& theD2,
1678                                            const Standard_Real TolArc,
1679                                            const Standard_Real TolTang,
1680                                            IntSurf_ListOfPntOn2S& ListOfPnts,
1681                                            const Standard_Boolean RestrictLine,
1682                                            const GeomAbs_SurfaceType theTyps1,
1683                                            const GeomAbs_SurfaceType theTyps2,
1684                                            const Standard_Boolean theIsReqToKeepRLine)
1685 {
1686   IntPatch_ImpImpIntersection interii(theS1,theD1,theS2,theD2,
1687                                       myTolArc,myTolTang, theIsReqToKeepRLine);
1688   const Standard_Boolean anIS = interii.IsDone();
1689   if (anIS)
1690   {
1691     done = anIS;
1692     empt = interii.IsEmpty();
1693     if (!empt)
1694     {
1695       tgte = interii.TangentFaces();
1696       if (tgte)
1697         oppo = interii.OppositeFaces();
1698
1699       for (Standard_Integer i = 1; i <= interii.NbLines(); i++)
1700       {
1701         const Handle(IntPatch_Line)& line = interii.Line(i);
1702         if (line->ArcType() == IntPatch_Analytic)
1703         {
1704           const GeomAbs_SurfaceType aTyps1 = theS1->GetType();
1705           const GeomAbs_SurfaceType aTyps2 = theS2->GetType();
1706           IntSurf_Quadric Quad1,Quad2;
1707
1708           switch(aTyps1)
1709           {
1710           case GeomAbs_Plane:
1711             Quad1.SetValue(theS1->Plane());
1712             break;
1713
1714           case GeomAbs_Cylinder:
1715             Quad1.SetValue(theS1->Cylinder());
1716             break;
1717
1718           case GeomAbs_Sphere:
1719             Quad1.SetValue(theS1->Sphere());
1720             break;
1721
1722           case GeomAbs_Cone:
1723             Quad1.SetValue(theS1->Cone());
1724             break;
1725
1726           case GeomAbs_Torus:
1727             Quad1.SetValue(theS1->Torus());
1728             break;
1729
1730           default:
1731             break;
1732           }
1733
1734           switch(aTyps2)
1735           {
1736           case GeomAbs_Plane:
1737             Quad2.SetValue(theS2->Plane());
1738             break;
1739           case GeomAbs_Cylinder:
1740             Quad2.SetValue(theS2->Cylinder());
1741             break;
1742
1743           case GeomAbs_Sphere:
1744             Quad2.SetValue(theS2->Sphere());
1745             break;
1746
1747           case GeomAbs_Cone:
1748             Quad2.SetValue(theS2->Cone());
1749             break;
1750
1751           case GeomAbs_Torus:
1752             Quad2.SetValue(theS2->Torus());
1753             break;
1754
1755           default:
1756             break;
1757           }
1758
1759           IntPatch_ALineToWLine AToW(Quad1,Quad2,0.01,0.05,aNbPointsInALine);
1760           Handle(IntPatch_Line) wlin=AToW.MakeWLine((*((Handle(IntPatch_ALine) *)(&line))));
1761           slin.Append(wlin);
1762         }
1763         else
1764           slin.Append(interii.Line(i));
1765       }
1766
1767       for (Standard_Integer i = 1; i <= interii.NbPnts(); i++)
1768       {
1769         spnt.Append(interii.Point(i));
1770       }
1771     }
1772   }
1773   else
1774     ParamParamPerfom(theS1, theD1, theS2, theD2,
1775                 TolArc, TolTang, ListOfPnts, RestrictLine, theTyps1, theTyps2);
1776 }
1777
1778 //=======================================================================
1779 //function : GeomParamPerfom
1780 //purpose  :
1781 //=======================================================================
1782 void IntPatch_Intersection::
1783   GeomParamPerfom(const Handle(Adaptor3d_HSurface)&  theS1,
1784                   const Handle(Adaptor3d_TopolTool)& theD1,
1785                   const Handle(Adaptor3d_HSurface)&  theS2,
1786                   const Handle(Adaptor3d_TopolTool)& theD2,
1787                   const Standard_Boolean isNotAnalitical,
1788                   const GeomAbs_SurfaceType typs1,
1789                   const GeomAbs_SurfaceType typs2)
1790 {
1791   IntPatch_ImpPrmIntersection interip;
1792   if (myIsStartPnt)
1793   {
1794     if (isNotAnalitical/*ts1 == 0*/)
1795       interip.SetStartPoint(myU1Start,myV1Start);
1796     else
1797       interip.SetStartPoint(myU2Start,myV2Start);
1798   }
1799
1800   if(theD1->DomainIsInfinite() && theD2->DomainIsInfinite())
1801   {
1802     Standard_Boolean IsPLInt = Standard_False;
1803     TColgp_SequenceOfPnt sop;
1804     gp_Vec v;
1805     FUN_PL_Intersection(theS1,typs1,theS2,typs2,IsPLInt,sop,v);
1806
1807     if(IsPLInt)
1808     {
1809       if(sop.Length() > 0)
1810       {
1811         for(Standard_Integer ip = 1; ip <= sop.Length(); ip++)
1812         {
1813           gp_Lin lin(sop.Value(ip),gp_Dir(v));
1814           Handle(IntPatch_GLine) gl = new IntPatch_GLine(lin,Standard_False);
1815           slin.Append(Handle(IntPatch_Line)::DownCast (gl));
1816         }
1817
1818         done = Standard_True;
1819       }
1820       else
1821         done = Standard_False;
1822
1823       return;
1824     }
1825     else
1826     {
1827       Handle(Adaptor3d_HSurface) nS1 = theS1;
1828       Handle(Adaptor3d_HSurface) nS2 = theS2;
1829       FUN_TrimBothSurf(theS1,typs1,theS2,typs2,1.e+5,nS1,nS2);
1830       interip.Perform(nS1,theD1,nS2,theD2,myTolArc,myTolTang,myFleche,myUVMaxStep);
1831     }
1832   }
1833   else
1834     interip.Perform(theS1,theD1,theS2,theD2,myTolArc,myTolTang,myFleche,myUVMaxStep);
1835
1836   if (interip.IsDone())
1837   {
1838     done = Standard_True;
1839     empt = interip.IsEmpty();
1840
1841     if (!empt)
1842     {
1843       const Standard_Integer aNbLines = interip.NbLines();
1844       for(Standard_Integer i = 1; i <= aNbLines; i++)
1845       {
1846         if(interip.Line(i)->ArcType() != IntPatch_Walking)
1847           slin.Append(interip.Line(i));
1848       }
1849
1850       for(Standard_Integer i = 1; i <= aNbLines; i++)
1851       {
1852         if(interip.Line(i)->ArcType() == IntPatch_Walking)
1853           slin.Append(interip.Line(i));
1854       }
1855
1856       for (Standard_Integer i = 1; i <= interip.NbPnts(); i++)
1857         spnt.Append(interip.Point(i));
1858     }
1859   }
1860 }
1861
1862 //=======================================================================
1863 //function : GeomGeomPerfomTrimSurf
1864 //purpose  : This function returns ready walking-line (which is not need
1865 //            in convertation) as an intersection line between two
1866 //            trimmed surfaces.
1867 //=======================================================================
1868 void IntPatch_Intersection::
1869   GeomGeomPerfomTrimSurf( const Handle(Adaptor3d_HSurface)& theS1,
1870                           const Handle(Adaptor3d_TopolTool)& theD1,
1871                           const Handle(Adaptor3d_HSurface)& theS2,
1872                           const Handle(Adaptor3d_TopolTool)& theD2,
1873                           const Standard_Real theTolArc,
1874                           const Standard_Real theTolTang,
1875                           IntSurf_ListOfPntOn2S& theListOfPnts,
1876                           const Standard_Boolean RestrictLine,
1877                           const GeomAbs_SurfaceType theTyps1,
1878                           const GeomAbs_SurfaceType theTyps2,
1879                           const Standard_Boolean theIsReqToKeepRLine)
1880 {
1881   IntSurf_Quadric Quad1,Quad2;
1882
1883   if((theTyps1 == GeomAbs_Cylinder) && (theTyps2 == GeomAbs_Cylinder))
1884   {
1885     IntPatch_ImpImpIntersection anInt;
1886     anInt.Perform(theS1, theD1, theS2, theD2, myTolArc,
1887                   myTolTang, Standard_True, theIsReqToKeepRLine);
1888
1889     done = anInt.IsDone();
1890
1891     if(done)
1892     {
1893       empt = anInt.IsEmpty();
1894       if (!empt)
1895       {
1896         tgte = anInt.TangentFaces();
1897         if (tgte)
1898           oppo = anInt.OppositeFaces();
1899
1900         const Standard_Integer aNbLin = anInt.NbLines();
1901         const Standard_Integer aNbPts = anInt.NbPnts();
1902
1903         for(Standard_Integer aLID = 1; aLID <= aNbLin; aLID++)
1904         {
1905           const Handle(IntPatch_Line)& aLine = anInt.Line(aLID);
1906           slin.Append(aLine);
1907         }
1908
1909         for(Standard_Integer aPID = 1; aPID <= aNbPts; aPID++)
1910         {
1911           const IntPatch_Point& aPoint = anInt.Point(aPID);
1912           spnt.Append(aPoint);
1913         }
1914
1915         JoinWLines( slin, spnt, theTolTang,
1916                     theS1->IsUPeriodic()? theS1->UPeriod() : 0.0,
1917                     theS2->IsUPeriodic()? theS2->UPeriod() : 0.0,
1918                     theS1->IsVPeriodic()? theS1->VPeriod() : 0.0,
1919                     theS2->IsVPeriodic()? theS2->VPeriod() : 0.0,
1920                     theS1->FirstUParameter(), theS1->LastUParameter(),
1921                     theS1->FirstVParameter(), theS1->LastVParameter(),
1922                     theS2->FirstUParameter(), theS2->LastUParameter(),
1923                     theS2->FirstVParameter(), theS2->LastVParameter());
1924       }
1925     }
1926   }
1927   else
1928   {
1929     GeomGeomPerfom(theS1, theD1, theS2, theD2,
1930             theTolArc, theTolTang, theListOfPnts,
1931             RestrictLine, theTyps1, theTyps2, theIsReqToKeepRLine);
1932   }
1933 }
1934
1935
1936 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  S1,
1937                                     const Handle(Adaptor3d_TopolTool)& D1,
1938                                     const Handle(Adaptor3d_HSurface)&  S2,
1939                                     const Handle(Adaptor3d_TopolTool)& D2,
1940                                     const Standard_Real U1,
1941                                     const Standard_Real V1,
1942                                     const Standard_Real U2,
1943                                     const Standard_Real V2,
1944                                     const Standard_Real TolArc,
1945                                     const Standard_Real TolTang)
1946 {
1947   myTolArc = TolArc;
1948   myTolTang = TolTang;
1949   if(myFleche == 0.0) {
1950 #if DEBUG
1951     //cout<<" -- IntPatch_Intersection::myFleche fixe par defaut a 0.01 --"<<endl;
1952     //cout<<" -- Utiliser la Methode SetTolerances( ... ) "<<endl;
1953 #endif
1954     myFleche = 0.01;
1955   }
1956   if(myUVMaxStep==0.0) {
1957 #if DEBUG
1958     //cout<<" -- IntPatch_Intersection::myUVMaxStep fixe par defaut a 0.01 --"<<endl;
1959     //cout<<" -- Utiliser la Methode SetTolerances( ... ) "<<endl;
1960 #endif
1961     myUVMaxStep = 0.01;
1962   }
1963
1964   done = Standard_False;
1965   spnt.Clear();
1966   slin.Clear();
1967
1968   empt = Standard_True;
1969   tgte = Standard_False;
1970   oppo = Standard_False;
1971
1972   const GeomAbs_SurfaceType typs1 = S1->GetType();
1973   const GeomAbs_SurfaceType typs2 = S2->GetType();
1974
1975   if(   typs1==GeomAbs_Plane
1976      || typs1==GeomAbs_Cylinder
1977      || typs1==GeomAbs_Sphere
1978      || typs1==GeomAbs_Cone
1979      || typs2==GeomAbs_Plane
1980      || typs2==GeomAbs_Cylinder
1981      || typs2==GeomAbs_Sphere
1982      || typs2==GeomAbs_Cone)
1983   {
1984     myIsStartPnt = Standard_True;
1985     myU1Start = U1; myV1Start = V1; myU2Start = U2; myV2Start = V2;
1986     Perform(S1,D1,S2,D2,TolArc,TolTang);
1987     myIsStartPnt = Standard_False;
1988   }
1989   else
1990   {
1991     IntPatch_PrmPrmIntersection interpp;
1992     interpp.Perform(S1,D1,S2,D2,U1,V1,U2,V2,TolArc,TolTang,myFleche,myUVMaxStep);
1993     if (interpp.IsDone())
1994     {
1995       done = Standard_True;
1996       tgte = Standard_False;
1997       empt = interpp.IsEmpty();
1998       const Standard_Integer nblm = interpp.NbLines();
1999       Standard_Integer i = 1;
2000       for (; i<=nblm; i++) slin.Append(interpp.Line(i));
2001     }
2002   }
2003 }
2004 //======================================================================
2005 #include <IntPatch_IType.hxx>
2006 #include <IntPatch_LineConstructor.hxx>
2007 #include <Adaptor2d_HCurve2d.hxx>
2008 #include <Geom_Curve.hxx>
2009 #define MAXR 200
2010
2011
2012 //void IntPatch_Intersection__MAJ_R(Handle(Adaptor2d_HCurve2d) *R1,
2013 //                                     Handle(Adaptor2d_HCurve2d) *R2,
2014 //                                     int *NR1,
2015 //                                     int *NR2,
2016 //                                     Standard_Integer nbR1,
2017 //                                     Standard_Integer nbR2,
2018 //                                     const IntPatch_Point& VTX)
2019 void IntPatch_Intersection__MAJ_R(Handle(Adaptor2d_HCurve2d) *,
2020                                   Handle(Adaptor2d_HCurve2d) *,
2021                                   int *,
2022                                   int *,
2023                                   Standard_Integer ,
2024                                   Standard_Integer ,
2025                                   const IntPatch_Point& )
2026 {
2027   /*
2028   if(VTX.IsOnDomS1()) {
2029
2030     //-- long unsigned ptr= *((long unsigned *)(((Handle(Standard_Transient) *)(&(VTX.ArcOnS1())))));
2031     for(Standard_Integer i=0; i<nbR1;i++) {
2032       if(VTX.ArcOnS1()==R1[i]) {
2033         NR1[i]++;
2034         printf("\n ******************************");
2035         return;
2036       }
2037     }
2038     printf("\n R Pas trouvee  (IntPatch)\n");
2039
2040   }
2041   */
2042 }
2043
2044
2045 //void IntPatch_Intersection::Dump(const Standard_Integer Mode,
2046 void IntPatch_Intersection::Dump(const Standard_Integer ,
2047                                  const Handle(Adaptor3d_HSurface)&  S1,
2048                                  const Handle(Adaptor3d_TopolTool)& D1,
2049                                  const Handle(Adaptor3d_HSurface)&  S2,
2050                                  const Handle(Adaptor3d_TopolTool)& D2) const
2051 {
2052
2053   //-- ----------------------------------------------------------------------
2054   //--  construction de la liste des restrictions & vertex
2055   //--
2056   int NR1[MAXR],NR2[MAXR];
2057   Handle(Adaptor2d_HCurve2d) R1[MAXR],R2[MAXR];
2058   Standard_Integer nbR1=0,nbR2=0;
2059   for(D1->Init();D1->More() && nbR1<MAXR; D1->Next()) {
2060     R1[nbR1]=D1->Value();
2061     NR1[nbR1]=0;
2062     nbR1++;
2063   }
2064   for(D2->Init();D2->More() && nbR2<MAXR; D2->Next()) {
2065     R2[nbR2]=D2->Value();
2066     NR2[nbR2]=0;
2067     nbR2++;
2068   }
2069
2070   printf("\nDUMP_INT:  ----empt:%2ud  tgte:%2ud  oppo:%2ud ---------------------------------",empt,tgte,empt);
2071   Standard_Integer i,j,nbr1,nbr2,nbgl,nbgc,nbge,nbgp,nbgh,nbl,nbr,nbg,nbw,nba;
2072   nbl=nbr=nbg=nbw=nba=nbgl=nbge=nbr1=nbr2=nbgc=nbgp=nbgh=0;
2073   nbl=NbLines();
2074   for(i=1;i<=nbl;i++) {
2075     const Handle(IntPatch_Line)& line=Line(i);
2076     const IntPatch_IType IType=line->ArcType();
2077     if(IType == IntPatch_Walking) nbw++;
2078     else     if(IType == IntPatch_Restriction) {
2079       nbr++;
2080       Handle(IntPatch_RLine) rlin (Handle(IntPatch_RLine)::DownCast (line));
2081       if(rlin->IsArcOnS1()) nbr1++;
2082       if(rlin->IsArcOnS2()) nbr2++;
2083     }
2084     else     if(IType == IntPatch_Analytic) nba++;
2085     else     {
2086       nbg++;
2087       if(IType == IntPatch_Lin) nbgl++;
2088       else if(IType == IntPatch_Circle) nbgc++;
2089       else if(IType == IntPatch_Parabola) nbgp++;
2090       else if(IType == IntPatch_Hyperbola) nbgh++;
2091       else if(IType == IntPatch_Ellipse) nbge++;
2092     }
2093   }
2094
2095
2096   printf("\nDUMP_INT:Lines:%2d Wlin:%2d Restr:%2d(On1:%2d On2:%2d) Ana:%2d Geom:%2d(L:%2d C:%2d E:%2d H:%2d P:%2d)",
2097          nbl,nbw,nbr,nbr1,nbr2,nba,nbg,nbgl,nbgc,nbge,nbgh,nbgp);
2098
2099   IntPatch_LineConstructor LineConstructor(2);
2100
2101   Standard_Integer nbllc=0;
2102   nbw=nbr=nbg=nba=0;
2103   Standard_Integer nbva,nbvw,nbvr,nbvg;
2104   nbva=nbvr=nbvw=nbvg=0;
2105   for (j=1; j<=nbl; j++) {
2106     Standard_Integer v,nbvtx;
2107     const Handle(IntPatch_Line)& intersLinej = Line(j);
2108     Standard_Integer NbLines;
2109     LineConstructor.Perform(SequenceOfLine(),intersLinej,S1,D1,S2,D2,1e-7);
2110     NbLines = LineConstructor.NbLines();
2111
2112     for(Standard_Integer k=1;k<=NbLines;k++) {
2113       nbllc++;
2114       const Handle(IntPatch_Line)& LineK = LineConstructor.Line(k);
2115       if (LineK->ArcType() == IntPatch_Analytic) {
2116         Handle(IntPatch_ALine) alin (Handle(IntPatch_ALine)::DownCast (LineK));
2117         nbvtx=alin->NbVertex();
2118         nbva+=nbvtx;        nba++;
2119         for(v=1;v<=nbvtx;v++) {
2120           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,alin->Vertex(v));
2121         }
2122       }
2123       else if (LineK->ArcType() == IntPatch_Restriction) {
2124         Handle(IntPatch_RLine) rlin (Handle(IntPatch_RLine)::DownCast (LineK));
2125         nbvtx=rlin->NbVertex();
2126         nbvr+=nbvtx;        nbr++;
2127         for(v=1;v<=nbvtx;v++) {
2128           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,rlin->Vertex(v));
2129         }
2130       }
2131       else if (LineK->ArcType() == IntPatch_Walking) {
2132         Handle(IntPatch_WLine) wlin (Handle(IntPatch_WLine)::DownCast (LineK));
2133         nbvtx=wlin->NbVertex();
2134         nbvw+=nbvtx;        nbw++;
2135         for(v=1;v<=nbvtx;v++) {
2136           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,wlin->Vertex(v));
2137         }
2138       }
2139       else {
2140         Handle(IntPatch_GLine) glin (Handle(IntPatch_GLine)::DownCast (LineK));
2141         nbvtx=glin->NbVertex();
2142         nbvg+=nbvtx;        nbg++;
2143         for(v=1;v<=nbvtx;v++) {
2144           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,glin->Vertex(v));
2145         }
2146       }
2147     }
2148   }
2149   printf("\nDUMP_LC :Lines:%2d WLin:%2d Restr:%2d Ana:%2d Geom:%2d",
2150          nbllc,nbw,nbr,nba,nbg);
2151   printf("\nDUMP_LC :vtx          :%2d     r:%2d    :%2d     :%2d",
2152          nbvw,nbvr,nbva,nbvg);
2153
2154
2155
2156    printf("\n");
2157 }