src/IntCurveSurface/IntCurveSurface_Inter.gxx

   1 // File:        IntCurveSurface_Inter.gxx
   2 // Created:     Fri  Apr 9 15:05:43 1993
   3 // Author:      Laurent BUCHARD
   4 //              <lbr@topsn3>
   5 #ifndef DEB
   6 #define No_Standard_RangeError
   7 #define No_Standard_OutOfRange
   8 #endif
   9
  10
  11 #define  TOLTANGENCY         0.00000001
  12 #define  TOLERANCE_ANGULAIRE 0.00000001
  13 #define  TOLERANCE           0.00000001
  14
  15 #define NBSAMPLESONCIRCLE  32
  16 #define NBSAMPLESONELLIPSE 32
  17 #define NBSAMPLESONPARAB   16
  18 #define NBSAMPLESONHYPR    32
  19
  20 #include <ElSLib.hxx>
  21 #include <Intf_SectionPoint.hxx>
  22 #include <Intf_TangentZone.hxx>
  23 #include <Intf_Tool.hxx>
  24 #include <math_FunctionSetRoot.hxx>
  25 #include <IntCurveSurface_IntersectionPoint.hxx>
  26 #include <IntCurveSurface_IntersectionSegment.hxx>
  27 #include <IntAna_IntConicQuad.hxx>
  28 #include <IntAna_Quadric.hxx>
  29 #include <gp_Vec.hxx>
  30 #include <gp_Dir.hxx>
  31 #include <Precision.hxx>
  32 #include <IntAna_IntLinTorus.hxx>
  33
  34 #include <GeomAbs_Shape.hxx>
  35 #include <GeomAbs_CurveType.hxx>
  36 #include <TColStd_Array1OfReal.hxx>
  37 #include <TColgp_Array1OfPnt.hxx>
  38 #include <gp_Pnt.hxx>
  39 #include <gp_Pln.hxx>
  40 #include <gp_Lin.hxx>
  41 #include <GProp_PEquation.hxx>
  42 #include <GProp_PGProps.hxx>
  43 #include <GProp_PrincipalProps.hxx>
  44 #include <Extrema_ExtElC.hxx>
  45 #include <Extrema_POnCurv.hxx>
  46
  47 #if 0
  48 //-- jgv patch (from)
  49 #include <Extrema_ExtPS.hxx>
  50 //-- jgv patch (to)
  51 #endif
  52
  53 #include <ProjLib_Plane.hxx>
  54 #include <IntAna2d_AnaIntersection.hxx>
  55 #include <gp_Lin2d.hxx>
  56 #include <IntAna2d_IntPoint.hxx>
  57 #include <gp_Parab2d.hxx>
  58 #include <IntAna2d_Conic.hxx>
  59 #include <IntAna2d_AnaIntersection.hxx>
  60 #include <gp_Hypr2d.hxx>
  61 #include <Adaptor3d_HCurve.hxx>
  62 #include <Adaptor3d_HSurface.hxx>
  63
  64 #if 0
  65 //-- jgv patch (from)
  66 #include <BndLib_Add3dCurve.hxx>
  67 #include <GeomAdaptor_Curve.hxx>
  68 #include <Geom_Line.hxx>
  69 #include <Geom_Circle.hxx>
  70 #include <Geom_Ellipse.hxx>
  71 #include <Geom_Hyperbola.hxx>
  72 #include <Geom_Parabola.hxx>
  73 #include <Geom_TrimmedCurve.hxx>
  74 #include <GeomAdaptor_Surface.hxx>
  75 #include <Geom_Plane.hxx>
  76 #include <Geom_CylindricalSurface.hxx>
  77 #include <Geom_ConicalSurface.hxx>
  78 #include <Geom_ToroidalSurface.hxx>
  79 #include <Geom_SphericalSurface.hxx>
  80 #include <Geom_SurfaceOfRevolution.hxx>
  81 #include <Geom_SurfaceOfLinearExtrusion.hxx>
  82 #include <Geom_OffsetSurface.hxx>
  83 #include <Geom_RectangularTrimmedSurface.hxx>
  84 //-- jgv patch (to)
  85 #endif
  86
  87 #include <TColgp_Array2OfPnt.hxx>
  88 #include <TColStd_HArray1OfReal.hxx>
  89 #include <Adaptor3d_TopolTool.hxx>
  90 #include <ElCLib.hxx>
  91
  92 //================================================================================
  93 static void EstLimForInfExtr(const gp_Lin&   Line,
  94                              const TheSurface& surface,
  95                              const Standard_Boolean IsOffSurf,
  96                              const Standard_Integer nbsu,
  97                              const Standard_Boolean U1inf,
  98                              const Standard_Boolean U2inf,
  99                              const Standard_Boolean V1inf,
 100                              const Standard_Boolean V2inf,
 101                              Standard_Real& U1new,
 102                              Standard_Real& U2new,
 103                              Standard_Real& V1new,
 104                              Standard_Real& V2new,
 105                              Standard_Boolean& NoIntersection);
 106
 107 //================================================================================
 108 static void EstLimForInfRevl(const gp_Lin&   Line,
 109                              const TheSurface& surface,
 110                              const Standard_Boolean U1inf,
 111                              const Standard_Boolean U2inf,
 112                              const Standard_Boolean V1inf,
 113                              const Standard_Boolean V2inf,
 114                              Standard_Real& U1new,
 115                              Standard_Real& U2new,
 116                              Standard_Real& V1new,
 117                              Standard_Real& V2new,
 118                              Standard_Boolean& NoIntersection);
 119
 120 //================================================================================
 121 static void EstLimForInfOffs(const gp_Lin&   Line,
 122                              const TheSurface& surface,
 123                              const Standard_Integer nbsu,
 124                              const Standard_Boolean U1inf,
 125                              const Standard_Boolean U2inf,
 126                              const Standard_Boolean V1inf,
 127                              const Standard_Boolean V2inf,
 128                              Standard_Real& U1new,
 129                              Standard_Real& U2new,
 130                              Standard_Real& V1new,
 131                              Standard_Real& V2new,
 132                              Standard_Boolean& NoIntersection);
 133
 134 //================================================================================
 135 static void EstLimForInfSurf(Standard_Real& U1new,
 136                              Standard_Real& U2new,
 137                              Standard_Real& V1new,
 138                              Standard_Real& V2new);
 139
 140 //================================================================================
 141 static void SectionPointToParameters(const Intf_SectionPoint& Sp,
 142                                      const IntCurveSurface_ThePolyhedron& Surf,
 143                                      const IntCurveSurface_ThePolygon&    Curv,
 144                                      Standard_Real& u,
 145                                      Standard_Real& v,
 146                                      Standard_Real& w);
 147 //================================================================================
 148 static void IntCurveSurface_ComputeTransitions(const TheCurve& curve,
 149                                                const Standard_Real w,
 150                                                IntCurveSurface_TransitionOnCurve&   TransOnCurve,
 151                                                const TheSurface& surface,
 152                                                const Standard_Real u,
 153                                                const Standard_Real v);
 154 //================================================================================
 155 static void IntCurveSurface_ComputeParamsOnQuadric(const TheSurface& surface,
 156                                                    const gp_Pnt& P,
 157                                                    Standard_Real& u,
 158                                                    Standard_Real& v);
 159 //================================================================================
 160 #if 0
 161 static Handle(Geom_Curve) GetCurve(const Handle(Adaptor3d_HCurve) AdCurve);
 162 //================================================================================
 163 static Handle(Geom_Surface) GetSurface(const Handle(Adaptor3d_HSurface) AdSurface);
 164 #endif
 165 //================================================================================
 166 //==
 167 IntCurveSurface_Inter::IntCurveSurface_Inter() {
 168 }
 169
 170 static Standard_Boolean DoTrim(const TheCurve& curve,
 171                                const TheSurface& surface)
 172 {
 173   Standard_Boolean isAnaCurve = Standard_False, isAnaSurface = Standard_False;
 174   GeomAbs_CurveType CurveType = TheCurveTool::GetType(curve);
 175   switch (CurveType) {
 176   case GeomAbs_Line:
 177   case GeomAbs_Circle:
 178   case GeomAbs_Ellipse:
 179   case GeomAbs_Hyperbola:
 180   case GeomAbs_Parabola: isAnaCurve = Standard_True; break;
 181   default: break;
 182   }
 183   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
 184   switch (SurfaceType) {
 185   case GeomAbs_Plane:
 186   case GeomAbs_Cylinder:
 187   case GeomAbs_Cone:
 188   case GeomAbs_Sphere:
 189   case GeomAbs_Torus: isAnaSurface = Standard_True; break;
 190   default: break;
 191   }
 192   Standard_Boolean result = (isAnaCurve && isAnaSurface) ? Standard_False : Standard_True;
 193   if(result) {
 194     Standard_Boolean isUClosed = (TheSurfaceTool::IsUClosed(surface) || TheSurfaceTool::IsUPeriodic(surface));
 195     Standard_Boolean isVClosed = (TheSurfaceTool::IsVClosed(surface) || TheSurfaceTool::IsVPeriodic(surface));
 196     if(isUClosed && isVClosed)
 197       result = Standard_False;
 198   }
 199   return result;
 200 }
 201
 202 // modified by NIZHNY-MKK  Tue Jul 26 14:41:59 2005
 203 // static void DoSurface(const TheSurface&   surface,
 204 void IntCurveSurface_Inter::DoSurface(const TheSurface&   surface,
 205                                       const Standard_Real u0,
 206                                       const Standard_Real u1,
 207                                       const Standard_Real v0,
 208                                       const Standard_Real v1,
 209                                       TColgp_Array2OfPnt& pntsOnSurface,
 210                                       Bnd_Box&            boxSurface,
 211                                       Standard_Real&      gap)
 212 {
 213   Standard_Integer iU = 0, iV = 0;
 214   Standard_Real U = 0., V = 0;
 215 // modified by NIZHNY-MKK  Mon Oct  3 17:38:45 2005
 216 //   Standard_Real dU = fabs(u1-u0)/50., dV = fabs(v1-v0)/50.;
 217   Standard_Real dU = (u1-u0)/50., dV = (v1-v0)/50.;
 218   gp_Pnt aPnt;
 219
 220   for(iU = 0; iU < 50; iU++) {
 221
 222     if(iU == 0)
 223       U = u0;
 224     else if(iU == 49)
 225       U = u1;
 226     else
 227       U = u0 + dU * ((Standard_Real)iU);
 228
 229     for(iV = 0; iV < 50; iV++) {
 230
 231       if(iV == 0)
 232         V = v0;
 233       else if(iV == 49)
 234         V = v1;
 235       else
 236         V = v0 + dV * ((Standard_Real)iV);
 237
 238       TheSurfaceTool::D0(surface,U,V,aPnt);
 239       boxSurface.Add(aPnt);
 240       pntsOnSurface.SetValue(iU+1,iV+1,aPnt);
 241     }
 242   }
 243   Standard_Real Ures = TheSurfaceTool::UResolution(surface,dU);
 244   Standard_Real Vres = TheSurfaceTool::VResolution(surface,dV);
 245   gap = Max(Ures,Vres);
 246 }
 247
 248 static void DoCurve(const TheCurve& curve,
 249                     Bnd_Box&        boxCurve)
 250 {
 251   Standard_Real CF = TheCurveTool::FirstParameter(curve);
 252   Standard_Real CL = TheCurveTool::LastParameter(curve);
 253   Standard_Real C = 0., dC = fabs(CL-CF)/50.;
 254   Standard_Integer iC = 0;
 255   gp_Pnt aPnt;
 256
 257   for(iC = 0; iC < 50; iC++) {
 258     if(iC == 0)
 259       C = CF;
 260     else if(iC == 49)
 261       C = CL;
 262     else
 263       C = CF + dC *((Standard_Real)iC);
 264
 265     TheCurveTool::D0(curve,C,aPnt);
 266     boxCurve.Add(aPnt);
 267   }
 268 }
 269 static void DoCommon(TColStd_Array1OfReal& Coords,
 270                      Standard_Real&        Cmin,
 271                      Standard_Real&        Cmax)
 272 {
 273   Standard_Integer i = 0, j = 0;
 274   for(i = 1; i <= 4; i++) {
 275     for(j = i+1; j <= 4; j++) {
 276       if(Coords(j) > Coords(i)) {
 277         Standard_Real c = Coords(i);
 278         Coords.SetValue(i,Coords(j));
 279         Coords.SetValue(j,c);
 280       }
 281     }
 282   }
 283   Cmax = Coords(2);
 284   Cmin = Coords(3);
 285 }
 286
 287 static void DoCommonBox(const Bnd_Box&        boxSurface,
 288                         const Bnd_Box&        boxCurve,
 289                         TColStd_Array1OfReal& X,
 290                         TColStd_Array1OfReal& Y,
 291                         TColStd_Array1OfReal& Z)
 292 {
 293   Standard_Real SBXmin = 0., SBYmin = 0., SBZmin = 0.;
 294   Standard_Real SBXmax = 0., SBYmax = 0., SBZmax = 0.;
 295   boxSurface.Get(SBXmin,SBYmin,SBZmin,SBXmax,SBYmax,SBZmax);
 296
 297   Standard_Real CBXmin = 0., CBYmin = 0., CBZmin = 0.;
 298   Standard_Real CBXmax = 0., CBYmax = 0., CBZmax = 0.;
 299   boxCurve.Get(CBXmin,CBYmin,CBZmin,CBXmax,CBYmax,CBZmax);
 300
 301   TColStd_Array1OfReal Coord(1,4);
 302
 303   Coord(1) = SBXmin; Coord(2) = SBXmax; Coord(3) = CBXmin; Coord(4) = CBXmax;
 304   Standard_Real CXmin = SBXmin, CXmax = SBXmax;
 305   DoCommon(Coord,CXmin,CXmax);
 306
 307   Coord(1) = SBYmin; Coord(2) = SBYmax; Coord(3) = CBYmin; Coord(4) = CBYmax;
 308   Standard_Real CYmin = SBYmin, CYmax = SBYmax;
 309   DoCommon(Coord,CYmin,CYmax);
 310
 311   Coord(1) = SBZmin; Coord(2) = SBZmax; Coord(3) = CBZmin; Coord(4) = CBZmax;
 312   Standard_Real CZmin = SBZmin, CZmax = SBZmax;
 313   DoCommon(Coord,CZmin,CZmax);
 314
 315   X.SetValue(1,CXmin); X.SetValue(2,CXmax);
 316   Y.SetValue(1,CYmin); Y.SetValue(2,CYmax);
 317   Z.SetValue(1,CZmin); Z.SetValue(2,CZmax);
 318 }
 319
 320 // modified by NIZHNY-MKK  Tue Jul 26 14:41:42 2005
 321 void IntCurveSurface_Inter::DoNewBounds(
 322 // static void DoNewBounds(
 323                                         const TheSurface&           surface,
 324                                         const Standard_Real         u0,
 325                                         const Standard_Real         u1,
 326                                         const Standard_Real         v0,
 327                                         const Standard_Real         v1,
 328                                         const TColgp_Array2OfPnt&   pntsOnSurface,
 329                                         const TColStd_Array1OfReal& X,
 330                                         const TColStd_Array1OfReal& Y,
 331                                         const TColStd_Array1OfReal& Z,
 332                                         TColStd_Array1OfReal&       Bounds)
 333 {
 334   Bounds.SetValue(1,u0);
 335   Bounds.SetValue(2,u1);
 336   Bounds.SetValue(3,v0);
 337   Bounds.SetValue(4,v1);
 338
 339   Standard_Boolean isUClosed = (TheSurfaceTool::IsUClosed(surface) || TheSurfaceTool::IsUPeriodic(surface));
 340   Standard_Boolean isVClosed = (TheSurfaceTool::IsVClosed(surface) || TheSurfaceTool::IsVPeriodic(surface));
 341   Standard_Boolean checkU = (isUClosed) ? Standard_False : Standard_True;
 342   Standard_Boolean checkV = (isVClosed) ? Standard_False : Standard_True;
 343
 344   Standard_Integer i = 0, j = 0, k = 0, iU = 0, iV = 0;
 345   Standard_Integer iUmin = 50, iVmin = 50, iUmax = 1, iVmax = 1;
 346
 347   for(i = 1; i <= 2; i++) {
 348     for(j = 1; j <= 2; j++) {
 349       for(k = 1; k <= 2; k++) {
 350         gp_Pnt aPoint(X(i),Y(j),Z(k));
 351         Standard_Real DistMin = 1.e+100;
 352         Standard_Integer diU = 0, diV = 0;
 353         for(iU = 1; iU <= 50; iU++) {
 354           for(iV = 1; iV <= 50; iV++) {
 355             const gp_Pnt aP = pntsOnSurface.Value(iU,iV);
 356             Standard_Real dist = aP.SquareDistance(aPoint);
 357             if(dist < DistMin) {
 358               DistMin = dist;
 359               diU = iU;
 360               diV = iV;
 361             }
 362           }
 363         }
 364         if(diU > 0 && diU < iUmin) iUmin = diU;
 365         if(diU > 0 && diU > iUmax) iUmax = diU;
 366         if(diV > 0 && diV < iVmin) iVmin = diV;
 367         if(diV > 0 && diV > iVmax) iVmax = diV;
 368       }
 369     }
 370   }
 371
 372 // modified by NIZHNY-MKK  Mon Oct  3 17:38:43 2005
 373 //   Standard_Real dU = fabs(u1-u0)/50., dV = fabs(v1-v0)/50.;
 374   Standard_Real dU = (u1-u0)/50., dV = (v1-v0)/50.;
 375
 376   Standard_Real USmin = u0 + dU * ((Standard_Real)(iUmin - 1));
 377   Standard_Real USmax = u0 + dU * ((Standard_Real)(iUmax - 1));
 378   Standard_Real VSmin = v0 + dV * ((Standard_Real)(iVmin - 1));
 379   Standard_Real VSmax = v0 + dV * ((Standard_Real)(iVmax - 1));
 380
 381   if(USmin > USmax) {
 382     Standard_Real tmp = USmax;
 383     USmax = USmin;
 384     USmin = tmp;
 385   }
 386   if(VSmin > VSmax) {
 387     Standard_Real tmp = VSmax;
 388     VSmax = VSmin;
 389     VSmin = tmp;
 390   }
 391
 392   USmin -= 1.5*dU;
 393   if(USmin < u0) USmin = u0;
 394   USmax += 1.5*dU;
 395   if(USmax > u1) USmax = u1;
 396   VSmin -= 1.5*dV;
 397   if(VSmin < v0) VSmin = v0;
 398   VSmax += 1.5*dV;
 399   if(VSmax > v1) VSmax = v1;
 400
 401   if(checkU) {
 402     Bounds.SetValue(1,USmin);
 403     Bounds.SetValue(2,USmax);
 404   }
 405   if(checkV) {
 406     Bounds.SetValue(3,VSmin);
 407     Bounds.SetValue(4,VSmax);
 408   }
 409 }
 410
 411 //================================================================================
 412 //==          P e r f o r m     g e n e r a l
 413 //==
 414 //==          Decompose la surface si besoin est
 415 //================================================================================
 416 void IntCurveSurface_Inter::Perform(const TheCurve&   curve,
 417                                     const TheSurface& surface) {
 418   ResetFields();
 419   done = Standard_True;
 420   Standard_Integer NbUOnS = TheSurfaceTool::NbUIntervals(surface,GeomAbs_C2);
 421   Standard_Integer NbVOnS = TheSurfaceTool::NbVIntervals(surface,GeomAbs_C2);
 422 #ifdef DEB
 423   Standard_Integer NbOnC  = TheCurveTool::NbIntervals(curve,GeomAbs_C2);
 424 #else
 425   TheCurveTool::NbIntervals(curve,GeomAbs_C2);
 426 #endif
 427   Standard_Real U0,U1,V0,V1;
 428
 429   if(NbUOnS > 1) {
 430     TColStd_Array1OfReal TabU(1,NbUOnS+1);
 431     TheSurfaceTool::UIntervals(surface,TabU,GeomAbs_C2);
 432     for(Standard_Integer iu = 1;iu <= NbUOnS; iu++) {
 433       U0 = TabU.Value(iu);
 434       U1 = TabU.Value(iu+1);
 435       if(NbVOnS > 1) {
 436         TColStd_Array1OfReal TabV(1,NbVOnS+1);
 437         TheSurfaceTool::VIntervals(surface,TabV,GeomAbs_C2);
 438         for(Standard_Integer iv = 1;iv <= NbVOnS; iv++) {
 439           V0 = TabV.Value(iv);
 440           V1 = TabV.Value(iv+1);
 441           Perform(curve,surface,U0,V0,U1,V1);
 442         }
 443       }
 444       else { //------ composite en U     non composite en V
 445         V0 = TheSurfaceTool::FirstVParameter(surface);
 446         V1 = TheSurfaceTool::LastVParameter(surface);
 447         Perform(curve,surface,U0,V0,U1,V1);
 448       }
 449     }
 450   }
 451   else if(NbVOnS > 1) { //---------- non composite en U  composite en V
 452     U0 = TheSurfaceTool::FirstUParameter(surface);
 453     U1 = TheSurfaceTool::LastUParameter(surface);
 454     TColStd_Array1OfReal TabV(1,NbVOnS+1);
 455     TheSurfaceTool::VIntervals(surface,TabV,GeomAbs_C2);
 456     for(Standard_Integer iv = 1;iv <= NbVOnS; iv++) {
 457       V0 = TabV.Value(iv);
 458       V1 = TabV.Value(iv+1);
 459       Perform(curve,surface,U0,V0,U1,V1);
 460     }
 461   }
 462   else {
 463     V0 = TheSurfaceTool::FirstVParameter(surface);
 464     V1 = TheSurfaceTool::LastVParameter(surface);
 465     U0 = TheSurfaceTool::FirstUParameter(surface);
 466     U1 = TheSurfaceTool::LastUParameter(surface);
 467
 468     //-- ofv: begin
 469     Standard_Boolean doTrim = DoTrim(curve,surface);
 470     if(doTrim) {
 471       TColgp_Array2OfPnt aPS(1,50,1,50);
 472       Bnd_Box SB;
 473       Standard_Real g = 1.e-7;
 474       DoSurface(surface,U0,U1,V0,V1,aPS,SB,g);
 475       Bnd_Box CB;
 476       DoCurve(curve,CB);
 477       CB.Enlarge(g);
 478       TColStd_Array1OfReal X(1,2), Y(1,2), Z(1,2);
 479       DoCommonBox(SB,CB,X,Y,Z);
 480       TColStd_Array1OfReal B(1,4);
 481       DoNewBounds(surface,U0,U1,V0,V1,aPS,X,Y,Z,B);
 482       U0 = B(1); U1 = B(2); V0 = B(3); V1 = B(4);
 483     }
 484     //-- ofv: end
 485
 486 #if 0
 487     //-- jgv patch (from)
 488     //Computing of local bounds
 489     Standard_Real LocalU0 = U0, LocalU1 = U1, LocalV0 = V0, LocalV1 = V1;
 490     Standard_Real Umin = RealLast(), Vmin = RealLast(), Umax = RealFirst(), Vmax = RealFirst();
 491     Bnd_Box CurveBox;
 492     Standard_Integer i, j, k;
 493     //Making GeomAdaptor_Curve
 494     Standard_Real f = TheCurveTool::FirstParameter(curve);
 495     Standard_Real l = TheCurveTool::LastParameter(curve);
 496     GeomAbs_CurveType CurveType = TheCurveTool::GetType(curve);
 497     Handle(Geom_Curve) theCurve;
 498     switch (CurveType)
 499       {
 500       case GeomAbs_Line:
 501         theCurve = new Geom_Line( TheCurveTool::Line(curve) );
 502         break;
 503       case GeomAbs_Circle:
 504         theCurve = new Geom_Circle( TheCurveTool::Circle(curve) );
 505         break;
 506       case GeomAbs_Ellipse:
 507         theCurve = new Geom_Ellipse( TheCurveTool::Ellipse(curve) );
 508         break;
 509       case GeomAbs_Hyperbola:
 510         theCurve = new Geom_Hyperbola( TheCurveTool::Hyperbola(curve) );
 511         break;
 512       case GeomAbs_Parabola:
 513         theCurve = new Geom_Parabola( TheCurveTool::Parabola(curve) );
 514         break;
 515       case GeomAbs_BezierCurve:
 516         theCurve = TheCurveTool::Bezier(curve);
 517         break;
 518       case GeomAbs_BSplineCurve:
 519         theCurve = TheCurveTool::BSpline(curve);
 520         break;
 521       }
 522     if (!theCurve.IsNull())
 523       {
 524         GeomAdaptor_Curve GACurve( theCurve, f, l );
 525         BndLib_Add3dCurve::Add( GACurve, Precision::Confusion(), CurveBox );
 526       }
 527     else
 528       {
 529         Standard_Integer nbp = TheCurveTool::NbSamples(curve,f,l);
 530         Standard_Real delta = (f-l)/nbp;
 531         for (i = 0; i <= nbp; i++)
 532           {
 533             gp_Pnt aPnt = TheCurveTool::Value(curve, f + i*delta);
 534             CurveBox.Add(aPnt);
 535           }
 536       }
 537     Standard_Real X[2], Y[2], Z[2];
 538     CurveBox.Get( X[0], Y[0], Z[0], X[1], Y[1], Z[1] );
 539     Standard_Real Ures = TheSurfaceTool::UResolution(surface, Precision::Confusion());
 540     Standard_Real Vres = TheSurfaceTool::VResolution(surface, Precision::Confusion());
 541     //Making GeomAdaptor_Surface
 542     GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
 543     Handle(Geom_Surface) theSurface;
 544     switch (SurfaceType)
 545       {
 546       case GeomAbs_Plane:
 547         theSurface = new Geom_Plane( TheSurfaceTool::Plane(surface) );
 548         break;
 549       case GeomAbs_Cylinder:
 550         theSurface = new Geom_CylindricalSurface( TheSurfaceTool::Cylinder(surface) );
 551         break;
 552       case GeomAbs_Cone:
 553         theSurface = new Geom_ConicalSurface( TheSurfaceTool::Cone(surface) );
 554         break;
 555       case GeomAbs_Torus:
 556         theSurface = new Geom_ToroidalSurface( TheSurfaceTool::Torus(surface) );
 557         break;
 558       case GeomAbs_Sphere:
 559         theSurface = new Geom_SphericalSurface( TheSurfaceTool::Sphere(surface) );
 560         break;
 561       case GeomAbs_BezierSurface:
 562         theSurface = TheSurfaceTool::Bezier(surface);
 563         break;
 564       case GeomAbs_BSplineSurface:
 565         theSurface = TheSurfaceTool::BSpline(surface);
 566         break;
 567       case GeomAbs_SurfaceOfRevolution:
 568         {
 569           gp_Ax1 Axis = TheSurfaceTool::AxeOfRevolution(surface);
 570           Handle(Adaptor3d_HCurve) AdBC = TheSurfaceTool::BasisCurve(surface);
 571           Handle(Geom_Curve) BC = GetCurve(AdBC);
 572           if (!BC.IsNull())
 573             theSurface = new Geom_SurfaceOfRevolution( BC, Axis );
 574           break;
 575         }
 576       case GeomAbs_SurfaceOfExtrusion:
 577         {
 578           gp_Dir Direction = TheSurfaceTool::Direction(surface);
 579           Handle(Adaptor3d_HCurve) AdBC = TheSurfaceTool::BasisCurve(surface);
 580           Handle(Geom_Curve) BC = GetCurve(AdBC);
 581           if (!BC.IsNull())
 582             theSurface = new Geom_SurfaceOfLinearExtrusion( BC, Direction );
 583           break;
 584         }
 585       case GeomAbs_OffsetSurface:
 586         {
 587           Standard_Real OffsetValue = TheSurfaceTool::OffsetValue(surface);
 588           Handle(Adaptor3d_HSurface) AdBS = TheSurfaceTool::BasisSurface(surface);
 589           Handle(Geom_Surface) BS = GetSurface(AdBS);
 590           if (!BS.IsNull())
 591             theSurface = new Geom_OffsetSurface( BS, OffsetValue );
 592           break;
 593         }
 594       }
 595     if (!theSurface.IsNull())
 596       {
 597         GeomAdaptor_Surface GASurface( theSurface );
 598         Extrema_ExtPS Projector;
 599         Projector.Initialize( GASurface, U0, U1, V0, V1, Ures, Vres );
 600         for (i = 0; i <= 1; i++)
 601           for (j = 0; j <= 1; j++)
 602             for (k = 0; k <= 1; k++)
 603               {
 604                 gp_Pnt aPoint( X[i], Y[j], Z[k] );
 605                 Projector.Perform( aPoint );
 606                 if (Projector.IsDone() && Projector.NbExt() > 0)
 607                   {
 608                     Standard_Real mindist = RealLast();
 609                     Standard_Integer minind, ind;
 610                     for (ind = 1; ind <= Projector.NbExt(); ind++)
 611                       {
 612                         Standard_Real dist = Projector.Value(ind);
 613                         if (dist < mindist)
 614                           {
 615                             mindist = dist;
 616                             minind  = ind;
 617                           }
 618                       }
 619                     Extrema_POnSurf pons = Projector.Point(minind);
 620                     Standard_Real U, V;
 621                     pons.Parameter( U, V );
 622                     if (U < Umin) Umin = U;
 623                     if (U > Umax) Umax = U;
 624                     if (V < Vmin) Vmin = V;
 625                     if (V > Vmax) Vmax = V;
 626                   }
 627               }
 628         Umin -= Ures; Umax += Ures; Vmin -= Vres; Vmax += Vres;
 629         if (Umin > U0 && Umin <= U1) LocalU0 = Umin;
 630         if (Umax < U1 && Umax >= U0) LocalU1 = Umax;
 631         if (Vmin > V0 && Vmin <= V1) LocalV0 = Vmin;
 632         if (Vmax < V1 && Vmax >= V0) LocalV1 = Vmax;
 633         U0 = LocalU0; U1 = LocalU1; V0 = LocalV0; V1 = LocalV1;
 634       }
 635     //-- jgv patch (to)
 636 #endif
 637     Perform(curve,surface,U0,V0,U1,V1);
 638   }
 639 }
 640 //================================================================================
 641 //================================================================================
 642 void IntCurveSurface_Inter::Perform(const TheCurve&   curve,
 643                                     const TheSurface& surface,
 644                                     const Standard_Real U1,const Standard_Real V1,
 645                                     const Standard_Real U2,const Standard_Real V2) {
 646
 647   GeomAbs_CurveType CurveType = TheCurveTool::GetType(curve);
 648
 649   switch(CurveType) {
 650   case GeomAbs_Line:
 651     {
 652       PerformConicSurf(TheCurveTool::Line(curve),curve,surface,U1,V1,U2,V2);
 653       break;
 654     }
 655   case GeomAbs_Circle:
 656     {
 657       PerformConicSurf(TheCurveTool::Circle(curve),curve,surface,U1,V1,U2,V2);
 658       break;
 659     }
 660   case GeomAbs_Ellipse:
 661     {
 662       PerformConicSurf(TheCurveTool::Ellipse(curve),curve,surface,U1,V1,U2,V2);
 663       break;
 664     }
 665   case GeomAbs_Parabola:
 666     {
 667       PerformConicSurf(TheCurveTool::Parabola(curve),curve,surface,U1,V1,U2,V2);
 668       break;
 669     }
 670   case GeomAbs_Hyperbola:
 671     {
 672       PerformConicSurf(TheCurveTool::Hyperbola(curve),curve,surface,U1,V1,U2,V2);
 673       break;
 674     }
 675   default:
 676     {
 677       Standard_Integer nbIntervalsOnCurve = TheCurveTool::NbIntervals(curve,GeomAbs_C2);
 678       GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
 679       if(   (SurfaceType != GeomAbs_Plane)
 680          && (SurfaceType != GeomAbs_Cylinder)
 681          && (SurfaceType != GeomAbs_Cone)
 682          && (SurfaceType != GeomAbs_Sphere)) {
 683
 684         if(nbIntervalsOnCurve > 1) {
 685           TColStd_Array1OfReal TabW(1,nbIntervalsOnCurve+1);
 686           TheCurveTool::Intervals(curve,TabW,GeomAbs_C2);
 687           for(Standard_Integer i = 1; i<=nbIntervalsOnCurve; i++) {
 688             Standard_Real u1,u2;
 689             u1 = TabW.Value(i);
 690             u2 = TabW.Value(i+1);
 691
 692             Handle(TColStd_HArray1OfReal) aPars;
 693             Standard_Real defl = 0.1;
 694             Standard_Integer NbMin = 10;
 695             TheCurveTool::SamplePars(curve, u1, u2, defl, NbMin, aPars);
 696
 697 //          IntCurveSurface_ThePolygon  polygon(curve,u1,u2,TheCurveTool::NbSamples(curve,u1,u2));
 698             IntCurveSurface_ThePolygon  polygon(curve, aPars->Array1());
 699             InternalPerform(curve,polygon,surface,U1,V1,U2,V2);
 700           }
 701         }
 702         else {
 703           Standard_Real u1,u2;
 704           u1 = TheCurveTool::FirstParameter(curve);
 705           u2 = TheCurveTool::LastParameter(curve);
 706
 707           Handle(TColStd_HArray1OfReal) aPars;
 708           Standard_Real defl = 0.1;
 709           Standard_Integer NbMin = 10;
 710           TheCurveTool::SamplePars(curve, u1, u2, defl, NbMin, aPars);
 711
 712 //        IntCurveSurface_ThePolygon       polygon(curve,TheCurveTool::NbSamples(curve,u1,u2));
 713           IntCurveSurface_ThePolygon  polygon(curve, aPars->Array1());
 714           InternalPerform(curve,polygon,surface,U1,V1,U2,V2);
 715         }
 716       }
 717       else {   //-- la surface est une quadrique
 718         InternalPerformCurveQuadric(curve,surface);
 719       }
 720     }
 721   }
 722 }
 723 //================================================================================
 724 void IntCurveSurface_Inter::Perform(const TheCurve&   curve,
 725                                     const IntCurveSurface_ThePolygon& polygon,
 726                                     const TheSurface& surface) {
 727   ResetFields();
 728   done = Standard_True;
 729   Standard_Real u1,v1,u2,v2;
 730   u1 = TheSurfaceTool::FirstUParameter(surface);
 731   v1 = TheSurfaceTool::FirstVParameter(surface);
 732   u2 = TheSurfaceTool::LastUParameter(surface);
 733   v2 = TheSurfaceTool::LastVParameter(surface);
 734   Standard_Integer nbsu,nbsv;
 735   nbsu = TheSurfaceTool::NbSamplesU(surface,u1,u2);
 736   nbsv = TheSurfaceTool::NbSamplesV(surface,v1,v2);
 737   if(nbsu>40) nbsu=40;
 738   if(nbsv>40) nbsv=40;
 739   IntCurveSurface_ThePolyhedron    polyhedron(surface,nbsu,nbsv,u1,v1,u2,v2);
 740   Perform(curve,polygon,surface,polyhedron);
 741 }
 742 //================================================================================
 743 void IntCurveSurface_Inter::Perform(const TheCurve&   curve,
 744                                     const TheSurface& surface,
 745                                     const IntCurveSurface_ThePolyhedron& polyhedron) {
 746   ResetFields();
 747   done = Standard_True;
 748   Standard_Real u1 = TheCurveTool::FirstParameter(curve);
 749   Standard_Real u2 = TheCurveTool::LastParameter(curve);
 750   IntCurveSurface_ThePolygon polygon(curve,TheCurveTool::NbSamples(curve,u1,u2));
 751   Perform(curve,polygon,surface,polyhedron);
 752 }
 753 //================================================================================
 754 void IntCurveSurface_Inter::Perform(const TheCurve&       curve,
 755                                     const IntCurveSurface_ThePolygon&     polygon,
 756                                     const TheSurface&     surface,
 757                                     const IntCurveSurface_ThePolyhedron&  polyhedron) {
 758   ResetFields();
 759   done = Standard_True;
 760   Standard_Real u1,v1,u2,v2;
 761   u1 = TheSurfaceTool::FirstUParameter(surface);
 762   v1 = TheSurfaceTool::FirstVParameter(surface);
 763   u2 = TheSurfaceTool::LastUParameter(surface);
 764   v2 = TheSurfaceTool::LastVParameter(surface);
 765   InternalPerform(curve,polygon,surface,polyhedron,u1,v1,u2,v2);
 766 }
 767
 768 void IntCurveSurface_Inter::Perform(const TheCurve&       curve,
 769                                     const IntCurveSurface_ThePolygon&     polygon,
 770                                     const TheSurface&     surface,
 771                                     const IntCurveSurface_ThePolyhedron&  polyhedron,
 772                                     Bnd_BoundSortBox& BndBSB) {
 773   ResetFields();
 774   done = Standard_True;
 775   Standard_Real u1,v1,u2,v2;
 776   u1 = TheSurfaceTool::FirstUParameter(surface);
 777   v1 = TheSurfaceTool::FirstVParameter(surface);
 778   u2 = TheSurfaceTool::LastUParameter(surface);
 779   v2 = TheSurfaceTool::LastVParameter(surface);
 780   InternalPerform(curve,polygon,surface,polyhedron,u1,v1,u2,v2,BndBSB);
 781 }
 782 //================================================================================
 783 //==          C a l c u l   d u   p o i n t   a p p r o c h e                   ==
 784 //==              p u i s   d u   p o i n t   E x a c t                         ==
 785 //================================================================================
 786 void IntCurveSurface_Inter::InternalPerform(const TheCurve&       curve,
 787                                             const IntCurveSurface_ThePolygon&     polygon,
 788                                             const TheSurface&     surface,
 789                                             const IntCurveSurface_ThePolyhedron&  polyhedron,
 790                                             const Standard_Real u0,
 791                                             const Standard_Real v0,
 792                                             const Standard_Real u1,
 793                                             const Standard_Real v1,
 794                                             Bnd_BoundSortBox& BSB) {
 795   IntCurveSurface_TheInterference  interference(polygon,polyhedron,BSB);
 796   IntCurveSurface_TheCSFunction    theicsfunction(surface,curve);
 797   IntCurveSurface_TheExactInter    intersectionExacte(theicsfunction,TOLTANGENCY);
 798   math_FunctionSetRoot rsnld(intersectionExacte.Function());
 799
 800 //  Standard_Real    u,v,w,winit;
 801   Standard_Real    u,v,w;
 802   gp_Pnt P;
 803   Standard_Real winf = polygon.InfParameter();
 804   Standard_Real wsup = polygon.SupParameter();
 805   Standard_Integer NbSectionPoints = interference.NbSectionPoints();
 806   Standard_Integer NbTangentZones  = interference.NbTangentZones();
 807
 808
 809   //-- Les interferences renvoient parfois de nombreuses fois (>20) les memes points
 810
 811   Standard_Integer i,NbStartPoints=NbSectionPoints;
 812   for(i=1; i<= NbTangentZones; i++) {
 813     const Intf_TangentZone& TZ = interference.ZoneValue(i);
 814     Standard_Integer nbpnts = TZ.NumberOfPoints();
 815     NbStartPoints+=nbpnts;
 816   }
 817
 818   if(NbStartPoints) {
 819     Standard_Real *TabU = new Standard_Real [NbStartPoints+1];
 820     Standard_Real *TabV = new Standard_Real [NbStartPoints+1];
 821     Standard_Real *TabW = new Standard_Real [NbStartPoints+1];
 822     Standard_Integer IndexPoint=0;
 823
 824     for(i=1; i<= NbSectionPoints; i++) {
 825       const Intf_SectionPoint& SP = interference.PntValue(i);
 826       SectionPointToParameters(SP,polyhedron,polygon,u,v,w);
 827       TabU[IndexPoint]=u;
 828       TabV[IndexPoint]=v;
 829       TabW[IndexPoint]=w;
 830       IndexPoint++;
 831     }
 832     for(i=1; i<= NbTangentZones; i++) {
 833       const Intf_TangentZone& TZ = interference.ZoneValue(i);
 834       Standard_Integer nbpnts = TZ.NumberOfPoints();
 835       for(Standard_Integer j=1; j<=nbpnts; j++) {
 836         const Intf_SectionPoint& SP = TZ.GetPoint(j);
 837         SectionPointToParameters(SP,polyhedron,polygon,u,v,w);
 838         TabU[IndexPoint]=u;
 839         TabV[IndexPoint]=v;
 840         TabW[IndexPoint]=w;
 841         IndexPoint++;
 842       }
 843     }
 844
 845     //-- Tri
 846     Standard_Real su=0,sv=0,sw=0,ptol;
 847     ptol = 10*Precision::PConfusion();
 848
 849     //-- Tri suivant la variable W
 850     Standard_Boolean Triok;
 851     do {
 852       Triok=Standard_True;
 853       Standard_Integer im1;
 854       for(i=1,im1=0;i<NbStartPoints;im1++,i++) {
 855         if(TabW[i] < TabW[im1]) {
 856           Standard_Real t=TabW[i]; TabW[i]=TabW[im1]; TabW[im1]=t;
 857           t=TabU[i]; TabU[i]=TabU[im1]; TabU[im1]=t;
 858           t=TabV[i]; TabV[i]=TabV[im1]; TabV[im1]=t;
 859           Triok=Standard_False;
 860         }
 861       }
 862     }
 863     while(Triok==Standard_False);
 864
 865     //-- On trie pour des meme W suivant U
 866     do {
 867       Triok=Standard_True;
 868       Standard_Integer im1;
 869       for(i=1,im1=0;i<NbStartPoints;im1++,i++) {
 870 // modified by NIZHNY-MKK  Mon Oct  3 17:38:49 2005
 871 //      if(Abs(TabW[i]-TabW[im1])<ptol) {
 872         if((TabW[i]-TabW[im1])<ptol) {
 873           TabW[i]=TabW[im1];
 874           if(TabU[i] < TabU[im1]) {
 875             Standard_Real t=TabU[i]; TabU[i]=TabU[im1]; TabU[im1]=t;
 876             t=TabV[i]; TabV[i]=TabV[im1]; TabV[im1]=t;
 877             Triok=Standard_False;
 878           }
 879         }
 880       }
 881     }
 882     while(Triok==Standard_False);
 883
 884     //-- On trie pour des meme U et W suivant V
 885     do {
 886       Triok=Standard_True;
 887       Standard_Integer im1;
 888       for(i=1,im1=0;i<NbStartPoints;im1++,i++) {
 889 // modified by NIZHNY-MKK  Mon Oct  3 17:38:52 2005
 890 //      if((Abs(TabW[i]-TabW[im1])<ptol) && (Abs(TabU[i]-TabU[im1])<ptol)) {
 891         if(((TabW[i]-TabW[im1])<ptol) && ((TabU[i]-TabU[im1])<ptol)) {
 892           TabU[i]=TabU[im1];
 893           if(TabV[i] < TabV[im1]) {
 894             Standard_Real t=TabV[i]; TabV[i]=TabV[im1]; TabV[im1]=t;
 895             Triok=Standard_False;
 896           }
 897         }
 898       }
 899     }
 900     while(Triok==Standard_False);
 901
 902
 903     for(i=0;i<NbStartPoints; i++) {
 904       u=TabU[i]; v=TabV[i]; w=TabW[i];
 905       if(i==0) {
 906         su=u-1;
 907       }
 908       if(Abs(u-su)>ptol || Abs(v-sv)>ptol || Abs(w-sw)>ptol) {
 909         intersectionExacte.Perform(u,v,w,rsnld,u0,u1,v0,v1,winf,wsup);
 910         if(intersectionExacte.IsDone()) {
 911           if(!intersectionExacte.IsEmpty()) {
 912             P=intersectionExacte.Point();
 913             w=intersectionExacte.ParameterOnCurve();
 914             intersectionExacte.ParameterOnSurface(u,v);
 915             AppendPoint(curve,w,surface,u,v);
 916           }
 917         }
 918       }
 919       su=TabU[i]; sv=TabV[i]; sw=TabW[i];
 920     }
 921     delete [] TabW;
 922     delete [] TabV;
 923     delete [] TabU;
 924   }
 925 }
 926
 927 void IntCurveSurface_Inter::InternalPerform(const TheCurve&       curve,
 928                                             const IntCurveSurface_ThePolygon&     polygon,
 929                                             const TheSurface&     surface,
 930                                             const IntCurveSurface_ThePolyhedron&  polyhedron,
 931                                             const Standard_Real u0,
 932                                             const Standard_Real v0,
 933                                             const Standard_Real u1,
 934                                             const Standard_Real v1) {
 935   IntCurveSurface_TheInterference  interference(polygon,polyhedron);
 936   IntCurveSurface_TheCSFunction    theicsfunction(surface,curve);
 937   IntCurveSurface_TheExactInter    intersectionExacte(theicsfunction,TOLTANGENCY);
 938   math_FunctionSetRoot rsnld(intersectionExacte.Function());
 939
 940 //  Standard_Real    u,v,w,winit;
 941   Standard_Real    u,v,w;
 942   gp_Pnt P;
 943   Standard_Real winf = polygon.InfParameter();
 944   Standard_Real wsup = polygon.SupParameter();
 945   Standard_Integer NbSectionPoints = interference.NbSectionPoints();
 946   Standard_Integer NbTangentZones  = interference.NbTangentZones();
 947
 948
 949   //-- Les interferences renvoient parfois de nombreuses fois (>20) les memes points
 950
 951   Standard_Integer i,NbStartPoints=NbSectionPoints;
 952   for(i=1; i<= NbTangentZones; i++) {
 953     const Intf_TangentZone& TZ = interference.ZoneValue(i);
 954     Standard_Integer nbpnts = TZ.NumberOfPoints();
 955     NbStartPoints+=nbpnts;
 956   }
 957
 958   if(NbStartPoints) {
 959     Standard_Real *TabU = new Standard_Real [NbStartPoints+1];
 960     Standard_Real *TabV = new Standard_Real [NbStartPoints+1];
 961     Standard_Real *TabW = new Standard_Real [NbStartPoints+1];
 962     Standard_Integer IndexPoint=0;
 963
 964     for(i=1; i<= NbSectionPoints; i++) {
 965       const Intf_SectionPoint& SP = interference.PntValue(i);
 966       SectionPointToParameters(SP,polyhedron,polygon,u,v,w);
 967       TabU[IndexPoint]=u;
 968       TabV[IndexPoint]=v;
 969       TabW[IndexPoint]=w;
 970       IndexPoint++;
 971     }
 972     for(i=1; i<= NbTangentZones; i++) {
 973       const Intf_TangentZone& TZ = interference.ZoneValue(i);
 974       Standard_Integer nbpnts = TZ.NumberOfPoints();
 975       for(Standard_Integer j=1; j<=nbpnts; j++) {
 976         const Intf_SectionPoint& SP = TZ.GetPoint(j);
 977         SectionPointToParameters(SP,polyhedron,polygon,u,v,w);
 978         TabU[IndexPoint]=u;
 979         TabV[IndexPoint]=v;
 980         TabW[IndexPoint]=w;
 981         IndexPoint++;
 982       }
 983     }
 984
 985     //-- Tri
 986     Standard_Real su=0,sv=0,sw=0,ptol;
 987     ptol = 10*Precision::PConfusion();
 988
 989     //-- Tri suivant la variable W
 990     Standard_Boolean Triok;
 991     do {
 992       Triok=Standard_True;
 993       Standard_Integer im1;
 994       for(i=1,im1=0;i<NbStartPoints;im1++,i++) {
 995         if(TabW[i] < TabW[im1]) {
 996           Standard_Real t=TabW[i]; TabW[i]=TabW[im1]; TabW[im1]=t;
 997           t=TabU[i]; TabU[i]=TabU[im1]; TabU[im1]=t;
 998           t=TabV[i]; TabV[i]=TabV[im1]; TabV[im1]=t;
 999           Triok=Standard_False;
1000         }
1001       }
1002     }
1003     while(Triok==Standard_False);
1004
1005     //-- On trie pour des meme W suivant U
1006     do {
1007       Triok=Standard_True;
1008       Standard_Integer im1;
1009       for(i=1,im1=0;i<NbStartPoints;im1++,i++) {
1010 // modified by NIZHNY-MKK  Mon Oct  3 17:38:56 2005
1011 //      if(Abs(TabW[i]-TabW[im1])<ptol) {
1012         if((TabW[i]-TabW[im1])<ptol) {
1013           TabW[i]=TabW[im1];
1014           if(TabU[i] < TabU[im1]) {
1015             Standard_Real t=TabU[i]; TabU[i]=TabU[im1]; TabU[im1]=t;
1016             t=TabV[i]; TabV[i]=TabV[im1]; TabV[im1]=t;
1017             Triok=Standard_False;
1018           }
1019         }
1020       }
1021     }
1022     while(Triok==Standard_False);
1023
1024     //-- On trie pour des meme U et W suivant V
1025     do {
1026       Triok=Standard_True;
1027       Standard_Integer im1;
1028       for(i=1,im1=0;i<NbStartPoints;im1++,i++) {
1029 // modified by NIZHNY-MKK  Mon Oct  3 17:38:58 2005
1030 //      if((Abs(TabW[i]-TabW[im1])<ptol) && (Abs(TabU[i]-TabU[im1])<ptol)) {
1031         if(((TabW[i]-TabW[im1])<ptol) && ((TabU[i]-TabU[im1])<ptol)) {
1032           TabU[i]=TabU[im1];
1033           if(TabV[i] < TabV[im1]) {
1034             Standard_Real t=TabV[i]; TabV[i]=TabV[im1]; TabV[im1]=t;
1035             Triok=Standard_False;
1036           }
1037         }
1038       }
1039     }
1040     while(Triok==Standard_False);
1041
1042
1043     for(i=0;i<NbStartPoints; i++) {
1044       u=TabU[i]; v=TabV[i]; w=TabW[i];
1045       if(i==0) {
1046         su=u-1;
1047       }
1048       if(Abs(u-su)>ptol || Abs(v-sv)>ptol || Abs(w-sw)>ptol) {
1049         intersectionExacte.Perform(u,v,w,rsnld,u0,u1,v0,v1,winf,wsup);
1050         if(intersectionExacte.IsDone()) {
1051           if(!intersectionExacte.IsEmpty()) {
1052             P=intersectionExacte.Point();
1053             w=intersectionExacte.ParameterOnCurve();
1054             intersectionExacte.ParameterOnSurface(u,v);
1055             AppendPoint(curve,w,surface,u,v);
1056           }
1057         }
1058       }
1059       su=TabU[i]; sv=TabV[i]; sw=TabW[i];
1060     }
1061     delete [] TabW;
1062     delete [] TabV;
1063     delete [] TabU;
1064   }
1065 }
1066 //================================================================================
1067 void IntCurveSurface_Inter::InternalPerformCurveQuadric(const TheCurve&                 curve,
1068                                                         const TheSurface&               surface) {
1069   IntCurveSurface_TheQuadCurvExactInter QuadCurv(surface,curve);
1070   if(QuadCurv.IsDone()) {
1071     Standard_Integer NbRoots = QuadCurv.NbRoots();
1072 #ifdef DEB
1073     Standard_Integer NbInter = QuadCurv.NbIntervals();
1074 #else
1075     QuadCurv.NbIntervals();
1076 #endif
1077     Standard_Real u,v,w;
1078     for(Standard_Integer i = 1; i<= NbRoots; i++) {
1079       w = QuadCurv.Root(i);
1080       IntCurveSurface_ComputeParamsOnQuadric(surface,TheCurveTool::Value(curve,w),u,v);
1081       AppendPoint(curve,w,surface,u,v);
1082     }
1083     //-- Intervals non traites .............................................
1084   }
1085 }
1086
1087 //================================================================================
1088 void IntCurveSurface_Inter::InternalPerform(const TheCurve&       curve,
1089                                             const IntCurveSurface_ThePolygon&  polygon,
1090                                             const TheSurface&     surface,
1091                                             const Standard_Real U1,
1092                                             const Standard_Real V1,
1093                                             const Standard_Real U2,
1094                                             const Standard_Real V2) {
1095   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1096   if(   (SurfaceType != GeomAbs_Plane)
1097      && (SurfaceType != GeomAbs_Cylinder)
1098      && (SurfaceType != GeomAbs_Cone)
1099      && (SurfaceType != GeomAbs_Sphere) ) {
1100     if(SurfaceType != GeomAbs_BSplineSurface) {
1101       Standard_Integer nbsu,nbsv;
1102       nbsu = TheSurfaceTool::NbSamplesU(surface,U1,U2);
1103       nbsv = TheSurfaceTool::NbSamplesV(surface,V1,V2);
1104       if(nbsu>40) nbsu=40;
1105       if(nbsv>40) nbsv=40;
1106       IntCurveSurface_ThePolyhedron polyhedron(surface,nbsu,nbsv,U1,V1,U2,V2);
1107       InternalPerform(curve,polygon,surface,polyhedron,U1,V1,U2,V2);
1108     }
1109     else {
1110       Handle(Adaptor3d_HSurface) aS = TheSurfaceTool::UTrim(surface, U1, U2, 1.e-9);
1111       aS = aS->VTrim(V1, V2, 1.e-9);
1112       Handle(Adaptor3d_TopolTool) aTopTool = new Adaptor3d_TopolTool(aS);
1113       Standard_Real defl = 0.1;
1114       aTopTool->SamplePnts(defl, 10, 10);
1115
1116       Standard_Integer nbpu = aTopTool->NbSamplesU();
1117       Standard_Integer nbpv = aTopTool->NbSamplesV();
1118       TColStd_Array1OfReal Upars(1, nbpu), Vpars(1, nbpv);
1119       aTopTool->UParameters(Upars);
1120       aTopTool->VParameters(Vpars);
1121
1122       IntCurveSurface_ThePolyhedron polyhedron(surface,Upars, Vpars);
1123       InternalPerform(curve,polygon,surface,polyhedron,U1,V1,U2,V2);
1124     }
1125   }
1126   else {
1127     IntCurveSurface_TheQuadCurvExactInter QuadCurv(surface,curve);
1128     if(QuadCurv.IsDone()) {
1129       Standard_Integer NbRoots = QuadCurv.NbRoots();
1130 #ifdef DEB
1131       Standard_Integer NbInter = QuadCurv.NbIntervals();
1132 #else
1133       QuadCurv.NbIntervals();
1134 #endif
1135       Standard_Real u,v,w;
1136       for(Standard_Integer i = 1; i<= NbRoots; i++) {
1137         w = QuadCurv.Root(i);
1138         IntCurveSurface_ComputeParamsOnQuadric(surface,TheCurveTool::Value(curve,w),u,v);
1139         AppendPoint(curve,w,surface,u,v);
1140       }
1141       //-- Intervalles non traites .............................................
1142     }
1143   } //-- Fin : la Surface  est une quadrique
1144 }
1145 //================================================================================
1146 void IntCurveSurface_Inter::PerformConicSurf(const gp_Lin&   Line,
1147                                              const TheCurve&       curve,
1148                                              const TheSurface&     surface,
1149                                              const Standard_Real U1,
1150                                              const Standard_Real V1,
1151                                              const Standard_Real U2,
1152                                              const Standard_Real V2) {
1153
1154
1155   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1156   switch(SurfaceType) {
1157   case GeomAbs_Plane:
1158     {
1159       IntAna_IntConicQuad LinPlane(Line,TheSurfaceTool::Plane(surface),TOLERANCE_ANGULAIRE);
1160       AppendIntAna(curve,surface,LinPlane);
1161       break;
1162     }
1163   case GeomAbs_Cylinder:
1164     {
1165       IntAna_IntConicQuad LinCylinder(Line,TheSurfaceTool::Cylinder(surface));
1166       AppendIntAna(curve,surface,LinCylinder);
1167       break;
1168     }
1169   case GeomAbs_Sphere:
1170     {
1171       IntAna_IntConicQuad LinSphere(Line,TheSurfaceTool::Sphere(surface));
1172       AppendIntAna(curve,surface,LinSphere);
1173       break;
1174     }
1175   case GeomAbs_Torus:
1176     {
1177       IntAna_IntLinTorus intlintorus(Line,TheSurfaceTool::Torus(surface));
1178       if(intlintorus.IsDone()) {
1179         Standard_Integer nbp = intlintorus.NbPoints();
1180         Standard_Real fi,theta,w;
1181         for(Standard_Integer i = 1; i<= nbp; i++) {
1182 #ifndef DEB
1183           gp_Pnt P;
1184           P = intlintorus.Value(i);
1185 #else
1186           gp_Pnt P(intlintorus.Value(i));
1187 #endif
1188           w = intlintorus.ParamOnLine(i);
1189           intlintorus.ParamOnTorus(i,fi,theta);
1190           AppendPoint(curve,w,surface,fi,theta);
1191         }
1192         break;
1193       }
1194     }  //-- Si Done retourne False, On passe dans Default !!
1195   case GeomAbs_Cone:
1196     {
1197       //OCC516(apo)->
1198       const Standard_Real correction = 1.E+5*Precision::Angular();
1199       gp_Cone cn = TheSurfaceTool::Cone(surface);
1200       if(Abs(cn.SemiAngle()) < M_PI/2.0 - correction){
1201         IntAna_IntConicQuad LinCone(Line,cn);
1202         AppendIntAna(curve,surface,LinCone);
1203         break;
1204       }//<-OCC516(apo)
1205     }
1206   default:
1207     {
1208       Standard_Integer nbsu,nbsv;
1209       nbsu = TheSurfaceTool::NbSamplesU(surface,U1,U2);
1210       nbsv = TheSurfaceTool::NbSamplesV(surface,V1,V2);
1211
1212       Standard_Boolean U1inf = Precision::IsInfinite(U1);
1213       Standard_Boolean U2inf = Precision::IsInfinite(U2);
1214       Standard_Boolean V1inf = Precision::IsInfinite(V1);
1215       Standard_Boolean V2inf = Precision::IsInfinite(V2);
1216
1217       Standard_Real U1new=U1, U2new=U2, V1new=V1, V2new=V2;
1218
1219       Standard_Boolean NoIntersection = Standard_False;
1220
1221       if(U1inf || U2inf || V1inf || V2inf ) {
1222
1223         if(SurfaceType == GeomAbs_SurfaceOfExtrusion) {
1224
1225           EstLimForInfExtr(Line, surface, Standard_False, nbsu,
1226                            U1inf, U2inf, V1inf, V2inf,
1227                            U1new, U2new, V1new, V2new, NoIntersection);
1228
1229         }
1230         else if(SurfaceType == GeomAbs_SurfaceOfRevolution) {
1231
1232           EstLimForInfRevl(Line, surface,
1233                            U1inf, U2inf, V1inf, V2inf,
1234                            U1new, U2new, V1new, V2new, NoIntersection);
1235
1236         }
1237         else if(SurfaceType == GeomAbs_OffsetSurface) {
1238
1239           EstLimForInfOffs(Line, surface, nbsu,
1240                            U1inf, U2inf, V1inf, V2inf,
1241                            U1new, U2new, V1new, V2new, NoIntersection);
1242         }
1243         else {
1244
1245           EstLimForInfSurf(U1new, U2new, V1new, V2new);
1246         }
1247
1248       }
1249
1250
1251       if(NoIntersection) return;
1252
1253
1254       // modified by NIZHNY-OFV  Mon Aug 20 14:56:47 2001 (60963 begin)
1255       if(nbsu<20) nbsu=20;
1256       if(nbsv<20) nbsv=20;
1257       // modified by NIZHNY-OFV  Mon Aug 20 14:57:06 2001 (60963 end)
1258       IntCurveSurface_ThePolyhedron polyhedron(surface,nbsu,nbsv,U1new,V1new,U2new,V2new);
1259       Intf_Tool                         bndTool;
1260       Bnd_Box                          boxLine;
1261       bndTool.LinBox(Line,polyhedron.Bounding(),boxLine);
1262       for(Standard_Integer nbseg=1; nbseg<= bndTool.NbSegments(); nbseg++) {
1263         Standard_Real pinf = bndTool.BeginParam(nbseg);
1264         Standard_Real psup = bndTool.EndParam(nbseg);
1265         if((psup - pinf)<1e-10) { pinf-=1e-10; psup+=1e-10; }
1266         IntCurveSurface_ThePolygon polygon(curve, pinf,psup,2);
1267         InternalPerform(curve,polygon,surface,polyhedron,U1new,V1new,U2new,V2new);
1268       }
1269     }
1270   }
1271 }
1272 //================================================================================
1273 void IntCurveSurface_Inter::PerformConicSurf(const gp_Circ&        Circle,
1274                                              const TheCurve&       curve,
1275                                              const TheSurface&     surface,
1276                                              const Standard_Real U1,
1277                                              const Standard_Real V1,
1278                                              const Standard_Real U2,
1279                                              const Standard_Real V2) {
1280
1281   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1282   switch(SurfaceType) {
1283   case GeomAbs_Plane:
1284     {
1285       IntAna_IntConicQuad CircPlane(Circle,TheSurfaceTool::Plane(surface),TOLERANCE_ANGULAIRE,TOLERANCE);
1286       AppendIntAna(curve,surface,CircPlane);
1287       break;
1288     }
1289   case GeomAbs_Cylinder:
1290     {
1291       IntAna_IntConicQuad CircCylinder(Circle,TheSurfaceTool::Cylinder(surface));
1292       AppendIntAna(curve,surface,CircCylinder);
1293       break;
1294     }
1295   case GeomAbs_Cone:
1296     {
1297       IntAna_IntConicQuad CircCone(Circle,TheSurfaceTool::Cone(surface));
1298       AppendIntAna(curve,surface,CircCone);
1299       break;
1300     }
1301   case GeomAbs_Sphere:
1302     {
1303       IntAna_IntConicQuad CircSphere(Circle,TheSurfaceTool::Sphere(surface));
1304       AppendIntAna(curve,surface,CircSphere);
1305       break;
1306     }
1307   default:
1308     {
1309       IntCurveSurface_ThePolygon polygon(curve,NBSAMPLESONCIRCLE);
1310       InternalPerform(curve,polygon,surface,U1,V1,U2,V2);
1311     }
1312   }
1313 }
1314 //================================================================================
1315 void IntCurveSurface_Inter::PerformConicSurf(const gp_Elips&       Ellipse,
1316                                              const TheCurve&       curve,
1317                                              const TheSurface&     surface,
1318                                              const Standard_Real U1,
1319                                              const Standard_Real V1,
1320                                              const Standard_Real U2,
1321                                              const Standard_Real V2) {
1322
1323   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1324   switch(SurfaceType) {
1325   case GeomAbs_Plane:
1326     {
1327       IntAna_IntConicQuad EllipsePlane(Ellipse,TheSurfaceTool::Plane(surface),TOLERANCE_ANGULAIRE,TOLERANCE);
1328       AppendIntAna(curve,surface,EllipsePlane);
1329       break;
1330     }
1331   case GeomAbs_Cylinder:
1332     {
1333       IntAna_IntConicQuad EllipseCylinder(Ellipse,TheSurfaceTool::Cylinder(surface));
1334       AppendIntAna(curve,surface,EllipseCylinder);
1335       break;
1336     }
1337   case GeomAbs_Cone:
1338      {
1339       IntAna_IntConicQuad EllipseCone(Ellipse,TheSurfaceTool::Cone(surface));
1340       AppendIntAna(curve,surface,EllipseCone);
1341       break;
1342     }
1343   case GeomAbs_Sphere:
1344     {
1345       IntAna_IntConicQuad EllipseSphere(Ellipse,TheSurfaceTool::Sphere(surface));
1346       AppendIntAna(curve,surface,EllipseSphere);
1347       break;
1348     }
1349   default:
1350     {
1351       IntCurveSurface_ThePolygon polygon(curve,NBSAMPLESONELLIPSE);
1352       InternalPerform(curve,polygon,surface,U1,V1,U2,V2);
1353     }
1354   }
1355 }
1356 //================================================================================
1357 void IntCurveSurface_Inter::PerformConicSurf(const gp_Parab&       Parab,
1358                                              const TheCurve&       curve,
1359                                              const TheSurface&     surface,
1360                                              const Standard_Real U1,
1361                                              const Standard_Real V1,
1362                                              const Standard_Real U2,
1363                                              const Standard_Real V2) {
1364
1365   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1366   switch(SurfaceType) {
1367   case GeomAbs_Plane:
1368     {
1369       IntAna_IntConicQuad ParabPlane(Parab,TheSurfaceTool::Plane(surface),TOLERANCE_ANGULAIRE);
1370       AppendIntAna(curve,surface,ParabPlane);
1371       break;
1372     }
1373   case GeomAbs_Cylinder:
1374     {
1375       IntAna_IntConicQuad ParabCylinder(Parab,TheSurfaceTool::Cylinder(surface));
1376       AppendIntAna(curve,surface,ParabCylinder);
1377       break;
1378     }
1379   case GeomAbs_Cone:
1380      {
1381       IntAna_IntConicQuad ParabCone(Parab,TheSurfaceTool::Cone(surface));
1382       AppendIntAna(curve,surface,ParabCone);
1383       break;
1384     }
1385   case GeomAbs_Sphere:
1386     {
1387       IntAna_IntConicQuad ParabSphere(Parab,TheSurfaceTool::Sphere(surface));
1388       AppendIntAna(curve,surface,ParabSphere);
1389       break;
1390     }
1391   default:
1392     {
1393       Standard_Integer nbsu,nbsv;
1394       nbsu = TheSurfaceTool::NbSamplesU(surface,U1,U2);
1395       nbsv = TheSurfaceTool::NbSamplesV(surface,V1,V2);
1396       if(nbsu>40) nbsu=40;
1397       if(nbsv>40) nbsv=40;
1398       IntCurveSurface_ThePolyhedron polyhedron(surface,nbsu,nbsv,U1,V1,U2,V2);
1399       Intf_Tool                         bndTool;
1400       Bnd_Box                          boxParab;
1401       bndTool.ParabBox(Parab,polyhedron.Bounding(),boxParab);
1402       for(Standard_Integer nbseg=1; nbseg<= bndTool.NbSegments(); nbseg++) {
1403         IntCurveSurface_ThePolygon polygon(curve,
1404                                            bndTool.BeginParam(nbseg),
1405                                            bndTool.EndParam(nbseg),
1406                                            NBSAMPLESONPARAB);
1407         InternalPerform(curve,polygon,surface,polyhedron,U1,V1,U2,V2);
1408       }
1409     }
1410   }
1411 }
1412 //================================================================================
1413 void IntCurveSurface_Inter::PerformConicSurf(const gp_Hypr&        Hypr,
1414                                              const TheCurve&       curve,
1415                                              const TheSurface&     surface,
1416                                              const Standard_Real U1,
1417                                              const Standard_Real V1,
1418                                              const Standard_Real U2,
1419                                              const Standard_Real V2) {
1420
1421   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1422   switch(SurfaceType) {
1423   case GeomAbs_Plane:
1424     {
1425       IntAna_IntConicQuad HyprPlane(Hypr,TheSurfaceTool::Plane(surface),TOLERANCE_ANGULAIRE);
1426       AppendIntAna(curve,surface,HyprPlane);
1427       break;
1428     }
1429   case GeomAbs_Cylinder:
1430     {
1431       IntAna_IntConicQuad HyprCylinder(Hypr,TheSurfaceTool::Cylinder(surface));
1432       AppendIntAna(curve,surface,HyprCylinder);
1433       break;
1434     }
1435   case GeomAbs_Cone:
1436      {
1437       IntAna_IntConicQuad HyprCone(Hypr,TheSurfaceTool::Cone(surface));
1438       AppendIntAna(curve,surface,HyprCone);
1439       break;
1440     }
1441   case GeomAbs_Sphere:
1442     {
1443       IntAna_IntConicQuad HyprSphere(Hypr,TheSurfaceTool::Sphere(surface));
1444       AppendIntAna(curve,surface,HyprSphere);
1445       break;
1446     }
1447   default:
1448     {
1449       Standard_Integer nbsu,nbsv;
1450       nbsu = TheSurfaceTool::NbSamplesU(surface,U1,U2);
1451       nbsv = TheSurfaceTool::NbSamplesV(surface,V1,V2);
1452       if(nbsu>40) nbsu=40;
1453       if(nbsv>40) nbsv=40;
1454       IntCurveSurface_ThePolyhedron polyhedron(surface,nbsu,nbsv,U1,V1,U2,V2);
1455       Intf_Tool                         bndTool;
1456       Bnd_Box                          boxHypr;
1457       bndTool.HyprBox(Hypr,polyhedron.Bounding(),boxHypr);
1458       for(Standard_Integer nbseg=1; nbseg<= bndTool.NbSegments(); nbseg++) {
1459         IntCurveSurface_ThePolygon polygon(curve,
1460                                            bndTool.BeginParam(nbseg),
1461                                            bndTool.EndParam(nbseg),
1462                                            NBSAMPLESONHYPR);
1463         InternalPerform(curve,polygon,surface,polyhedron,U1,V1,U2,V2);
1464       }
1465     }
1466   }
1467 }
1468 //================================================================================
1469 void IntCurveSurface_Inter::AppendIntAna(const TheCurve& curve,
1470                                          const TheSurface& surface,
1471                                          const IntAna_IntConicQuad& intana_ConicQuad) {
1472   if(intana_ConicQuad.IsDone()) {
1473     if(intana_ConicQuad.IsInQuadric()) {
1474       //-- cout<<" Courbe Dans la Quadrique !!! Non Traite !!!"<<endl;
1475     }
1476     else if(intana_ConicQuad.IsParallel()) {
1477       //-- cout<<" Courbe // a la Quadrique !!! Non Traite !!!"<<endl;
1478     }
1479     else {
1480       Standard_Integer nbp = intana_ConicQuad.NbPoints();
1481       Standard_Real u,v,w;
1482       for(Standard_Integer i = 1; i<= nbp; i++) {
1483         gp_Pnt P(intana_ConicQuad.Point(i));
1484         w = intana_ConicQuad.ParamOnConic(i);
1485         IntCurveSurface_ComputeParamsOnQuadric(surface,P,u,v);
1486         AppendPoint(curve,w,surface,u,v);
1487       }
1488     }
1489   }
1490   else {
1491     //-- cout<<" IntAna Conic Quad Not Done  "<<endl;
1492   }
1493 }
1494 //================================================================================
1495 void IntCurveSurface_Inter::AppendPoint(const TheCurve& curve,
1496                                         const Standard_Real lw,
1497                                         const TheSurface& surface,
1498                                         const Standard_Real su,
1499                                         const Standard_Real sv) {
1500
1501   Standard_Real W0 = TheCurveTool::FirstParameter(curve);
1502   Standard_Real W1 = TheCurveTool::LastParameter(curve);
1503   Standard_Real U0 = TheSurfaceTool::FirstUParameter(surface);
1504   Standard_Real U1 = TheSurfaceTool::LastUParameter(surface);
1505   Standard_Real V0 = TheSurfaceTool::FirstVParameter(surface);
1506   Standard_Real V1 = TheSurfaceTool::LastVParameter(surface);
1507   //-- Test si la courbe est periodique
1508   Standard_Real w = lw, u = su, v = sv;
1509
1510   GeomAbs_CurveType aCType = TheCurveTool::GetType(curve);
1511
1512   if(TheCurveTool::IsPeriodic(curve)
1513      || aCType == GeomAbs_Circle
1514      || aCType == GeomAbs_Ellipse) {
1515     w = ElCLib::InPeriod(w, W0, W0 + TheCurveTool::Period(curve));
1516   }
1517
1518   if((W0 - w) >= TOLTANGENCY || (w - W1) >= TOLTANGENCY) return;
1519
1520   GeomAbs_SurfaceType aSType = TheSurfaceTool::GetType(surface);
1521   if (TheSurfaceTool::IsUPeriodic(surface)
1522       || aSType == GeomAbs_Cylinder
1523       || aSType == GeomAbs_Cone
1524       || aSType == GeomAbs_Sphere) {
1525     u = ElCLib::InPeriod(u, U0, U0 + TheSurfaceTool::UPeriod(surface));
1526   }
1527
1528   if (TheSurfaceTool::IsVPeriodic(surface)) {
1529     v = ElCLib::InPeriod(v, V0, V0 + TheSurfaceTool::VPeriod(surface));
1530   }
1531
1532   if((U0 - u) >= TOLTANGENCY || (u - U1) >= TOLTANGENCY) return;
1533   if((V0 - v) >= TOLTANGENCY || (v - V1) >= TOLTANGENCY) return;
1534
1535   IntCurveSurface_TransitionOnCurve   TransOnCurve;
1536   IntCurveSurface_ComputeTransitions(curve,w,TransOnCurve,
1537                                      surface,u,v);
1538   gp_Pnt P(TheCurveTool::Value(curve,w));
1539   IntCurveSurface_IntersectionPoint IP(P,u,v,w,TransOnCurve);
1540   Append(IP); //-- invoque la methode de IntCurveSurface_Intersection.
1541
1542 }
1543
1544 //================================================================================
1545 void IntCurveSurface_Inter::AppendSegment(const TheCurve& ,
1546                                           const Standard_Real ,
1547                                           const Standard_Real ,
1548                                           const TheSurface& ) {
1549   cout<<" !!! Not Yet Implemented    IntCurveSurface_Inter::Append(const IntCurveSurf ...)"<<endl;
1550 }
1551
1552 //================================================================================
1553 //== P o i n t   d   i n t e r f e r e n c e   - - >     U , V       e t   W    ==
1554 //================================================================================
1555 void SectionPointToParameters(const Intf_SectionPoint&              Sp,
1556                               const IntCurveSurface_ThePolyhedron&  Polyhedron,
1557                               const IntCurveSurface_ThePolygon&     Polygon,
1558                               Standard_Real& U,
1559                               Standard_Real& V,
1560                               Standard_Real& W) {
1561
1562   Intf_PIType       typ;
1563   Standard_Integer  Adr1,Adr2;
1564   Standard_Real     Param,u,v;
1565   gp_Pnt P(Sp.Pnt());
1566
1567   Standard_Integer Pt1,Pt2,Pt3;
1568   Standard_Real u1,v1,param;
1569   //----------------------------------------------------------------------
1570   //--          Calcul des parametres approches sur la surface          --
1571   //----------------------------------------------------------------------
1572
1573   Sp.InfoSecond(typ,Adr1,Adr2,Param);
1574   switch(typ) {
1575   case Intf_VERTEX:   //-- Adr1 est le numero du vertex
1576     {
1577       Polyhedron.Parameters(Adr1,u1,v1);
1578       break;
1579     }
1580   case Intf_EDGE:
1581     {
1582       Polyhedron.Parameters(Adr1,u1,v1);
1583       Polyhedron.Parameters(Adr2,u,v);
1584       u1+= Param * (u-u1);
1585       v1+= Param * (v-v1);
1586       break;
1587     }
1588   case Intf_FACE:
1589     {
1590       Standard_Real ua,va,ub,vb,uc,vc,ca,cb,cc,cabc;
1591       Polyhedron.Triangle(Adr1,Pt1,Pt2,Pt3);
1592       gp_Pnt PA(Polyhedron.Point(Pt1));
1593       gp_Pnt PB(Polyhedron.Point(Pt2));
1594       gp_Pnt PC(Polyhedron.Point(Pt3));
1595       Polyhedron.Parameters(Pt1,ua,va);
1596       Polyhedron.Parameters(Pt2,ub,vb);
1597       Polyhedron.Parameters(Pt3,uc,vc);
1598       gp_Vec Normale(gp_Vec(PA,PB).Crossed(gp_Vec(PA,PC)));
1599       cc = (gp_Vec(PA,PB).Crossed(gp_Vec(PA,P))).Dot(Normale);
1600       ca = (gp_Vec(PB,PC).Crossed(gp_Vec(PB,P))).Dot(Normale);
1601       cb = (gp_Vec(PC,PA).Crossed(gp_Vec(PC,P))).Dot(Normale);
1602       cabc = ca + cb + cc;
1603
1604       ca/=cabc;     cb/=cabc;       cc/=cabc;
1605
1606       u1 = ca * ua + cb * ub + cc * uc;
1607       v1 = ca * va + cb * vb + cc * vc;
1608       break;
1609     }
1610   default:
1611     {
1612       cout<<" Default dans SectionPointToParameters "<<endl;
1613       break;
1614     }
1615   }
1616   //----------------------------------------------------------------------
1617   //--                Calcul du point approche sur la Curve             --
1618   //----------------------------------------------------------------------
1619   Standard_Integer SegIndex;
1620
1621   Sp.InfoFirst(typ,SegIndex,param);
1622   W = Polygon.ApproxParamOnCurve(SegIndex,param);
1623   if(param>1.0 || param<0.0) {
1624     //-- IntCurveSurface_ThePolyhedronTool::Dump(Polyhedron);
1625     //-- IntCurveSurface_ThePolygonTool::Dump(Polygon);
1626     }
1627   U = u1;
1628   V = v1;
1629 }
1630 //================================================================================
1631 void IntCurveSurface_ComputeTransitions(const TheCurve& curve,
1632                                         const Standard_Real w,
1633                                         IntCurveSurface_TransitionOnCurve&   TransOnCurve,
1634                                         const TheSurface& surface,
1635                                         const Standard_Real u,
1636                                         const Standard_Real v) {
1637
1638   gp_Vec NSurf,D1U,D1V;//TgCurv;
1639   gp_Pnt Psurf;
1640   Standard_Real CosDir;
1641
1642
1643   TheSurfaceTool::D1(surface,u,v,Psurf,D1U,D1V);
1644   NSurf = D1U.Crossed(D1V);
1645   TheCurveTool::D1(curve,w,Psurf,D1U);
1646   Standard_Real Norm = NSurf.Magnitude();
1647   if(Norm>TOLERANCE_ANGULAIRE) {
1648     D1U.Normalize();
1649     CosDir = NSurf.Dot(D1U);
1650     CosDir/=Norm;
1651     if( -CosDir > TOLERANCE_ANGULAIRE) {
1652       //--  --Curve--->    <----Surface----
1653       TransOnCurve   = IntCurveSurface_In;
1654     }
1655     else if(CosDir > TOLERANCE_ANGULAIRE) {
1656       //--  --Curve--->  ----Surface-->
1657       TransOnCurve   = IntCurveSurface_Out;
1658     }
1659     else {
1660       TransOnCurve   = IntCurveSurface_Tangent;
1661     }
1662   }
1663   else {
1664     TransOnCurve   = IntCurveSurface_Tangent;
1665   }
1666 }
1667 //================================================================================
1668 void IntCurveSurface_ComputeParamsOnQuadric(const TheSurface& surface,
1669                                             const gp_Pnt& P,
1670                                             Standard_Real& u,
1671                                             Standard_Real& v) {
1672   GeomAbs_SurfaceType SurfaceType = TheSurfaceTool::GetType(surface);
1673   switch(SurfaceType) {
1674   case GeomAbs_Plane:
1675     {
1676       ElSLib::Parameters(TheSurfaceTool::Plane(surface),P,u,v);
1677         break;
1678       }
1679   case GeomAbs_Cylinder:
1680     {
1681       ElSLib::Parameters(TheSurfaceTool::Cylinder(surface),P,u,v);
1682         break;
1683       }
1684   case GeomAbs_Cone:
1685     {
1686       ElSLib::Parameters(TheSurfaceTool::Cone(surface),P,u,v);
1687         break;
1688       }
1689   case GeomAbs_Sphere:
1690     {
1691       ElSLib::Parameters(TheSurfaceTool::Sphere(surface),P,u,v);
1692         break;
1693       }
1694 #ifndef DEB
1695   default:      break;
1696 #endif
1697   }
1698 }
1699 //================================================================================
1700
1701
1702
1703
1704 /*
1705
1706   Standard_Real u1,v1,u2,v2;
1707   u1 = TheSurfaceTool::FirstUParameter(surface);
1708   v1 = TheSurfaceTool::FirstVParameter(surface);
1709   u2 = TheSurfaceTool::LastUParameter(surface);
1710   v2 = TheSurfaceTool::LastVParameter(surface);
1711   Standard_Integer nbsu,nbsv;
1712   nbsu = TheSurfaceTool::NbSamplesU(surface,u1,u2);
1713   nbsv = TheSurfaceTool::NbSamplesV(surface,v1,v2);
1714   IntCurveSurface_ThePolyhedron    polyhedron(surface,nbsu,nbsv,u1,v1,u2,v2);
1715
1716 */
1717 //====================================================================================
1718 // Estimation of limits for infinite surfaces
1719 //====================================================================================
1720
1721 //================================================================================
1722 static void EstLimForInfExtr(const gp_Lin&   Line,
1723                              const TheSurface& surface,
1724                              const Standard_Boolean IsOffSurf,
1725                              const Standard_Integer nbsu,
1726                              const Standard_Boolean U1inf,
1727                              const Standard_Boolean U2inf,
1728                              const Standard_Boolean V1inf,
1729                              const Standard_Boolean V2inf,
1730                              Standard_Real& U1new,
1731                              Standard_Real& U2new,
1732                              Standard_Real& V1new,
1733                              Standard_Real& V2new,
1734                              Standard_Boolean& NoIntersection)
1735
1736 {
1737
1738   NoIntersection = Standard_False;
1739
1740   Handle(Adaptor3d_HSurface) aBasSurf;
1741
1742   if(IsOffSurf) aBasSurf = TheSurfaceTool::BasisSurface(surface);
1743
1744   gp_Dir aDirOfExt;
1745
1746   if(IsOffSurf) aDirOfExt = aBasSurf->Direction();
1747   else aDirOfExt = TheSurfaceTool::Direction(surface);
1748
1749   Standard_Real tolang = TOLERANCE_ANGULAIRE;
1750
1751   if(aDirOfExt.IsParallel(Line.Direction(), tolang)) {
1752     NoIntersection = Standard_True;
1753     return;
1754   }
1755
1756   if((V1inf || V2inf) && !(U1inf || U2inf)) {
1757
1758     Standard_Real vmin = RealLast(), vmax = -vmin;
1759     gp_Lin aL;
1760     Standard_Real step = (U2new - U1new) / nbsu;
1761     Standard_Real u = U1new, v;
1762     gp_Pnt aP;
1763     Extrema_POnCurv aP1, aP2;
1764     Standard_Integer i;
1765
1766     for(i = 0; i <= nbsu; i++) {
1767
1768       TheSurfaceTool::D0(surface, u, 0., aP);
1769       aL.SetLocation(aP);
1770       aL.SetDirection(aDirOfExt);
1771
1772       Extrema_ExtElC aExtr(aL, Line, tolang);
1773
1774       if(!aExtr.IsDone()) return;
1775
1776       if(aExtr.IsParallel()) {
1777         NoIntersection = Standard_True;
1778         return;
1779       }
1780
1781       aExtr.Points(1, aP1, aP2);
1782       v = aP1.Parameter();
1783       vmin = Min(vmin, v);
1784       vmax = Max(vmax, v);
1785
1786       u += step;
1787
1788     }
1789
1790     vmin = vmin - Abs(vmin) - 10.;
1791     vmax = vmax + Abs(vmax) + 10.;
1792
1793     V1new = Max(V1new, vmin);
1794     V2new = Min(V2new, vmax);
1795
1796   }
1797   else if(U1inf || U2inf) {
1798
1799     Standard_Real umin = RealLast(), umax = -umin;
1800     Standard_Real u0 = Min(Max(0., U1new), U2new);
1801     Standard_Real v0 = Min(Max(0., V1new), V2new);
1802     gp_Pnt aP;
1803     TheSurfaceTool::D0(surface, u0, v0, aP);
1804     gp_Pln aRefPln(aP, aDirOfExt);
1805
1806     Handle(Adaptor3d_HCurve) aBasCurv;
1807
1808     if(IsOffSurf) aBasCurv = aBasSurf->BasisCurve();
1809     else aBasCurv = TheSurfaceTool::BasisCurve(surface);
1810
1811     ProjLib_Plane Projector(aRefPln);
1812
1813     Projector.Project(Line);
1814
1815     if(!Projector.IsDone()) return;
1816
1817     gp_Lin2d Line2d = Projector.Line();
1818
1819     GeomAbs_CurveType aCurvTyp = aBasCurv->GetType();
1820
1821     if(aCurvTyp == GeomAbs_Line) {
1822
1823       Projector.Project(aBasCurv->Line());
1824
1825       if(!Projector.IsDone()) return;
1826
1827       gp_Lin2d aL2d = Projector.Line();
1828
1829       IntAna2d_AnaIntersection anInter(Line2d, aL2d);
1830
1831       if(!anInter.IsDone()) return;
1832
1833       if(anInter.IsEmpty() || anInter.IdenticalElements() ||
1834                               anInter.ParallelElements()     ) {
1835         NoIntersection = Standard_True;
1836         return;
1837       }
1838
1839       const IntAna2d_IntPoint& anIntPnt = anInter.Point(1);
1840       umin = umax = anIntPnt.ParamOnSecond();
1841     }
1842     else if(aCurvTyp == GeomAbs_Parabola || aCurvTyp == GeomAbs_Hyperbola) {
1843
1844       IntAna2d_Conic aCon(Line2d);
1845       IntAna2d_AnaIntersection anInter;
1846
1847       if(aCurvTyp == GeomAbs_Parabola) {
1848         Projector.Project(aBasCurv->Parabola());
1849         if(!Projector.IsDone()) return;
1850
1851         const gp_Parab2d& aP2d = Projector.Parabola();
1852
1853         anInter.Perform(aP2d, aCon);
1854       }
1855       else {
1856         Projector.Project(aBasCurv->Hyperbola());
1857         if(!Projector.IsDone()) return;
1858
1859         const gp_Hypr2d& aH2d = Projector.Hyperbola();
1860         anInter.Perform(aH2d, aCon);
1861       }
1862
1863       if(!anInter.IsDone()) return;
1864
1865       if(anInter.IsEmpty()) {
1866         NoIntersection = Standard_True;
1867         return;
1868       }
1869
1870       Standard_Integer i, nbint = anInter.NbPoints();
1871       for(i = 1; i <= nbint; i++) {
1872
1873         const IntAna2d_IntPoint& anIntPnt = anInter.Point(i);
1874
1875         umin = Min(anIntPnt.ParamOnFirst(), umin);
1876         umax = Max(anIntPnt.ParamOnFirst(), umax);
1877
1878       }
1879
1880     }
1881     else {
1882       return;
1883     }
1884
1885     umin = umin - Abs(umin) - 10;
1886     umax = umax + Abs(umax) + 10;
1887
1888     U1new = Max(U1new, umin);
1889     U2new = Min(U2new, umax);
1890
1891     if(V1inf || V2inf) {
1892       EstLimForInfExtr(Line, surface, IsOffSurf, nbsu,
1893                        Standard_False, Standard_False, V1inf, V2inf,
1894                        U1new, U2new, V1new, V2new, NoIntersection);
1895     }
1896   }
1897
1898   return;
1899
1900 }
1901
1902 //================================================================================
1903 //=======================================================================
1904 //function : ProjectIntersectAndEstLim
1905 //purpose  : project <theLine> and it's X-axe symmetric line to <thePln> and
1906 //           intersect resulting curve with <theBasCurvProj>.
1907 //           Then estimate max and min parameters of intersection on
1908 //           <theBasCurvProj>.
1909 //           Is called from EstLimForInfRevl()
1910 //=======================================================================
1911
1912 static void ProjectIntersectAndEstLim(const gp_Lin&        theLine,
1913                                       const gp_Pln&        thePln,
1914                                       const ProjLib_Plane& theBasCurvProj,
1915                                       Standard_Real&       theVmin,
1916                                       Standard_Real&       theVmax,
1917                                       Standard_Boolean&    theNoIntersection)
1918 {
1919   ProjLib_Plane aLineProj( thePln, theLine );
1920   if (!aLineProj.IsDone()) {
1921 #ifdef DEB
1922     cout
1923       << "Info: IntCurveSurface_Inter::ProjectIntersectAndEstLim(), !aLineProj.IsDone()"
1924         << endl;
1925 #endif
1926     return;
1927   }
1928   gp_Lin2d aLin2d = aLineProj.Line();
1929
1930   // make a second line X-axe symmetric to the first one
1931   gp_Pnt2d aP1 = aLin2d.Location();
1932   gp_Pnt2d aP2 (aP1.XY() + aLin2d.Direction().XY());
1933   gp_Pnt2d aP1sym (aP1.X(), -aP1.Y());
1934   gp_Pnt2d aP2sym (aP2.X(), -aP2.Y());
1935   gp_Lin2d aLin2dsym (aP1sym, gp_Vec2d(aP1sym,aP2sym));
1936
1937   // intersect projections
1938   IntAna2d_Conic aCon    (aLin2d);
1939   IntAna2d_Conic aConSym (aLin2dsym);
1940   IntAna2d_AnaIntersection anIntersect;
1941   IntAna2d_AnaIntersection anIntersectSym;
1942
1943   switch (theBasCurvProj.GetType()) {
1944   case GeomAbs_Line:
1945     anIntersectSym.Perform(theBasCurvProj.Line(), aConSym);
1946     anIntersect.Perform(theBasCurvProj.Line(), aCon); break;
1947   case GeomAbs_Hyperbola:
1948     anIntersectSym.Perform(theBasCurvProj.Hyperbola(), aConSym);
1949     anIntersect.Perform(theBasCurvProj.Hyperbola(), aCon); break;
1950   case GeomAbs_Parabola:
1951     anIntersectSym.Perform(theBasCurvProj.Parabola(), aConSym);
1952     anIntersect.Perform(theBasCurvProj.Parabola(), aCon); break;
1953   default:
1954     return; // not infinite curve
1955   }
1956
1957   // retrieve params of intersections
1958   Standard_Integer aNbIntPnt = anIntersect.IsDone() ? anIntersect.NbPoints() : 0;
1959   Standard_Integer aNbIntPntSym = anIntersectSym.IsDone() ? anIntersectSym.NbPoints() : 0;
1960   Standard_Integer iPnt, aNbPnt = Max (aNbIntPnt, aNbIntPntSym);
1961
1962   if (aNbPnt == 0) {
1963     theNoIntersection = Standard_True;
1964     return;
1965   }
1966   Standard_Real aParam;
1967   for (iPnt = 1; iPnt <= aNbPnt; iPnt++)
1968   {
1969     if (iPnt <= aNbIntPnt) {
1970       const IntAna2d_IntPoint& aIntPnt = anIntersect.Point(iPnt);
1971       aParam = aIntPnt.ParamOnFirst();
1972       theVmin = Min (theVmin, aParam );
1973       theVmax = Max (theVmax, aParam );
1974     }
1975     if (iPnt <= aNbIntPntSym) {
1976       const IntAna2d_IntPoint& aIntPnt = anIntersectSym.Point(iPnt);
1977       aParam = aIntPnt.ParamOnFirst();
1978       theVmin = Min (theVmin, aParam );
1979       theVmax = Max (theVmax, aParam );
1980     }
1981   }
1982
1983   return;
1984 }
1985 //=======================================================================
1986 //function : EstLimForInfRevl
1987 //purpose  : Estimate V1 and V2 to pass to InternalPerform() if they are
1988 //           infinite for Surface of Revolution
1989 //           Algo: intersect projections of Line and basis curve on the
1990 //           plane passing through revolution axe
1991 //=======================================================================
1992
1993 static void EstLimForInfRevl(const gp_Lin&   Line,
1994                              const TheSurface& surface,
1995                              const Standard_Boolean U1inf,
1996                              const Standard_Boolean U2inf,
1997                              const Standard_Boolean V1inf,
1998                              const Standard_Boolean V2inf,
1999                              Standard_Real& U1new,
2000                              Standard_Real& U2new,
2001                              Standard_Real& V1new,
2002                              Standard_Real& V2new,
2003                              Standard_Boolean& NoIntersection)
2004 {
2005
2006   NoIntersection = Standard_False;
2007
2008   if (U1inf || U2inf) {
2009     if (U1inf)
2010       U1new = Max (0., U1new);
2011     else
2012       U2new = Min (2 * M_PI, U2new);
2013     if (! V1inf && !V2inf) return;
2014   }
2015
2016   Handle(Adaptor3d_HCurve) aBasisCurve = TheSurfaceTool::BasisCurve(surface);
2017   gp_Ax1 aRevAx = TheSurfaceTool::AxeOfRevolution(surface);
2018   gp_Vec aXVec = aRevAx.Direction();
2019   Standard_Real aTolAng = Precision::Angular();
2020
2021   // make plane to project a basis curve
2022   gp_Pnt O = aRevAx.Location();
2023   Standard_Real aU = 0.;
2024   gp_Pnt P = aBasisCurve->Value(aU);
2025   while (O.SquareDistance(P) <= Precision::PConfusion() ||
2026          aXVec.IsParallel( gp_Vec(O,P), aTolAng)) {
2027     aU += 1.;
2028     P = aBasisCurve->Value(aU);
2029     if (aU > 3)
2030       // basis curve is a line coinciding with aXVec, P is any not on aXVec
2031       P = gp_Pnt(aU, aU+1, aU+2);
2032   }
2033   gp_Vec aNVec = aXVec ^ gp_Vec(O,P);
2034   gp_Pln aPln (gp_Ax3 (O, aNVec ,aXVec));
2035
2036   // project basic curve
2037   ProjLib_Plane aBasCurvProj(aPln);
2038   switch (aBasisCurve->GetType()) {
2039   case GeomAbs_Line:
2040     aBasCurvProj.Project(aBasisCurve->Line     ()); break;
2041   case GeomAbs_Hyperbola:
2042     aBasCurvProj.Project(aBasisCurve->Hyperbola()); break;
2043   case GeomAbs_Parabola:
2044     aBasCurvProj.Project(aBasisCurve->Parabola ()); break;
2045   default:
2046     return; // not infinite curve
2047   }
2048   if (!aBasCurvProj.IsDone()) {
2049 #ifdef DEB
2050     cout << "Info: IntCurveSurface_Inter::EstLimForInfRevl(), !aBasCurvProj.IsDone()" << endl;
2051 #endif
2052     return;
2053   }
2054   // make plane to project Line
2055   if (aXVec.IsParallel( Line.Direction(), aTolAng)) {
2056     P = Line.Location();
2057     while (O.SquareDistance(P) <= Precision::PConfusion()) {
2058       aU += 1.;
2059       P = gp_Pnt(aU, aU+1, aU+2); // any not on aXVec
2060     }
2061     aNVec = aXVec ^ gp_Vec(O,P);
2062   } else
2063     aNVec = aXVec.Crossed( Line.Direction() );
2064
2065   aPln  = gp_Pln (gp_Ax3 (O, aNVec ,aXVec));
2066
2067   // make a second plane perpendicular to the first one, rotated around aXVec
2068   gp_Pln aPlnPrp = aPln.Rotated (gp_Ax1 (O,aXVec), M_PI/2.);
2069
2070   // project Line and it's X-axe symmetric one to plane and intersect
2071   // resulting curve with projection of Basic Curev
2072   Standard_Real aVmin = RealLast(), aVmax = -aVmin;
2073   Standard_Boolean aNoInt1 = Standard_False, aNoInt2 = Standard_False;
2074   ProjectIntersectAndEstLim (Line, aPln,    aBasCurvProj, aVmin, aVmax, aNoInt1);
2075   ProjectIntersectAndEstLim (Line, aPlnPrp, aBasCurvProj, aVmin, aVmax, aNoInt2);
2076
2077   if (aNoInt1 && aNoInt2) {
2078     NoIntersection = Standard_True;
2079     return;
2080   }
2081
2082   aVmin = aVmin - Abs(aVmin) - 10;
2083   aVmax = aVmax + Abs(aVmax) + 10;
2084
2085   if (V1inf) V1new = aVmin;
2086   if (V2inf) V2new = aVmax;
2087
2088   //cout << "EstLimForInfRevl: Vmin " << V1new << " Vmax " << V2new << endl;
2089
2090   return;
2091 }
2092
2093 //================================================================================
2094 static void EstLimForInfOffs(const gp_Lin&   Line,
2095                              const TheSurface& surface,
2096                              const Standard_Integer nbsu,
2097                              const Standard_Boolean U1inf,
2098                              const Standard_Boolean U2inf,
2099                              const Standard_Boolean V1inf,
2100                              const Standard_Boolean V2inf,
2101                              Standard_Real& U1new,
2102                              Standard_Real& U2new,
2103                              Standard_Real& V1new,
2104                              Standard_Real& V2new,
2105                              Standard_Boolean& NoIntersection)
2106 {
2107
2108   NoIntersection = Standard_False;
2109
2110   const Handle(Adaptor3d_HSurface)& aBasSurf = TheSurfaceTool::BasisSurface(surface);
2111   Standard_Real anOffVal = TheSurfaceTool::OffsetValue(surface);
2112
2113   GeomAbs_SurfaceType aTypeOfBasSurf = aBasSurf->GetType();
2114
2115   //  case for plane, cylinder and cone - make equivalent surface;
2116   if(aTypeOfBasSurf == GeomAbs_Plane) {
2117
2118     gp_Pln aPln = aBasSurf->Plane();
2119     gp_Vec aT = aPln.Position().XDirection()^aPln.Position().YDirection();
2120     aT *= anOffVal;
2121     aPln.Translate(aT);
2122     IntAna_IntConicQuad LinPlane(Line,aPln,TOLERANCE_ANGULAIRE);
2123
2124     if(!LinPlane.IsDone()) return;
2125
2126     if(LinPlane.IsParallel() || LinPlane.IsInQuadric()) {
2127
2128       NoIntersection = Standard_True;
2129       return;
2130
2131     }
2132
2133     Standard_Real u, v;
2134     ElSLib::Parameters(aPln, LinPlane.Point(1), u, v);
2135     U1new = Max(U1new, u - 10.);
2136     U2new = Min(U2new, u + 10.);
2137     V1new = Max(V1new, v - 10.);
2138     V2new = Min(V2new, v + 10.);
2139
2140   }
2141   else if(aTypeOfBasSurf == GeomAbs_Cylinder) {
2142
2143     gp_Cylinder aCyl = aBasSurf->Cylinder();
2144
2145     Standard_Real aR = aCyl.Radius();
2146     gp_Ax3 anA = aCyl.Position();
2147
2148     if (anA.Direct())
2149       aR += anOffVal;
2150     else
2151       aR -= anOffVal;
2152
2153     if ( aR >= TOLTANGENCY ) {
2154       aCyl.SetRadius(aR);
2155     }
2156     else if ( aR <= -TOLTANGENCY ){
2157       anA.Rotate(gp_Ax1(anA.Location(), anA.Direction()), M_PI);
2158       aCyl.SetPosition(anA);
2159 // modified by NIZHNY-MKK  Mon Oct  3 17:37:54 2005
2160 //       aCyl.SetRadius(Abs(aR));
2161       aCyl.SetRadius(-aR);
2162     }
2163     else {
2164
2165       NoIntersection = Standard_True;
2166       return;
2167
2168     }
2169
2170     IntAna_IntConicQuad LinCylinder(Line, aCyl);
2171
2172     if(!LinCylinder.IsDone()) return;
2173
2174     if(LinCylinder.IsParallel() || LinCylinder.IsInQuadric()) {
2175
2176       NoIntersection = Standard_True;
2177       return;
2178
2179     }
2180
2181     Standard_Integer i, nbp = LinCylinder.NbPoints();
2182     Standard_Real vmin = RealLast(), vmax = -vmin, u, v;
2183
2184     for(i = 1; i <= nbp; i++) {
2185
2186       ElSLib::Parameters(aCyl, LinCylinder.Point(i), u, v);
2187       vmin = Min(vmin, v);
2188       vmax = Max(vmax, v);
2189
2190     }
2191
2192     V1new = Max(V1new, vmin - Abs(vmin) - 10.);
2193     V2new = Min(V2new, vmax + Abs(vmax) + 10.);
2194
2195   }
2196   else if(aTypeOfBasSurf == GeomAbs_Cone) {
2197
2198     gp_Cone aCon = aBasSurf->Cone();
2199     Standard_Real anAng = aCon.SemiAngle();
2200     Standard_Real aR = aCon.RefRadius() + anOffVal * Cos(anAng);
2201     gp_Ax3 anA = aCon.Position();
2202     if ( aR >= 0.) {
2203
2204       gp_Vec aZ( anA.Direction());
2205       aZ *= - anOffVal * Sin(anAng);
2206       anA.Translate(aZ);
2207       aCon.SetPosition(anA);
2208       aCon.SetRadius(aR);
2209       aCon.SetSemiAngle(anAng);
2210
2211     }
2212     else {
2213
2214       return;
2215
2216     }
2217
2218     IntAna_IntConicQuad LinCone(Line, aCon);
2219
2220     if(!LinCone.IsDone()) return;
2221
2222     if(LinCone.IsParallel() || LinCone.IsInQuadric()) {
2223
2224       NoIntersection = Standard_True;
2225       return;
2226
2227     }
2228
2229     Standard_Integer i, nbp = LinCone.NbPoints();
2230     Standard_Real vmin = RealLast(), vmax = -vmin, u, v;
2231
2232     for(i = 1; i <= nbp; i++) {
2233
2234       ElSLib::Parameters(aCon, LinCone.Point(i), u, v);
2235       vmin = Min(vmin, v);
2236       vmax = Max(vmax, v);
2237
2238     }
2239
2240     V1new = Max(V1new, vmin - Abs(vmin) - 10.);
2241     V2new = Min(V2new, vmax + Abs(vmax) + 10.);
2242
2243   }
2244   else if(aTypeOfBasSurf == GeomAbs_SurfaceOfExtrusion) {
2245
2246     Standard_Real anU1 = U1new, anU2 = U2new;
2247
2248     EstLimForInfExtr(Line, surface, Standard_True, nbsu,
2249                      U1inf, U2inf, V1inf, V2inf,
2250                      U1new, U2new, V1new, V2new, NoIntersection);
2251
2252     if(NoIntersection) return;
2253
2254     if(U1inf || U2inf) {
2255
2256       GeomAbs_CurveType aBasCurvType = aBasSurf->BasisCurve()->GetType();
2257       if(aBasCurvType == GeomAbs_Line) {
2258         U1new = Max(anU1, -1.e10);
2259         U2new = Min(anU2,  1.e10);
2260       }
2261       else if(aBasCurvType == GeomAbs_Parabola) {
2262         gp_Parab aPrb  = aBasSurf->BasisCurve()->Parabola();
2263         Standard_Real aF = aPrb.Focal();
2264         Standard_Real dU = 2.e5 * Sqrt(aF);
2265         U1new = Max(anU1, -dU);
2266         U2new = Min(anU2,  dU);
2267       }
2268       else if(aBasCurvType == GeomAbs_Hyperbola) {
2269         U1new = Max(anU1, -30.);
2270         U2new = Min(anU2,  30.);
2271       }
2272       else {
2273         U1new = Max(anU1, -1.e10);
2274         U2new = Min(anU2,  1.e10);
2275       }
2276
2277     }
2278
2279   }
2280   else if(aTypeOfBasSurf == GeomAbs_SurfaceOfRevolution) {
2281
2282     GeomAbs_CurveType aBasCurvType = aBasSurf->BasisCurve()->GetType();
2283     if(aBasCurvType == GeomAbs_Line) {
2284       V1new = Max(V1new, -1.e10);
2285       V2new = Min(V2new,  1.e10);
2286     }
2287     else if(aBasCurvType == GeomAbs_Parabola) {
2288       gp_Parab aPrb  = aBasSurf->BasisCurve()->Parabola();
2289       Standard_Real aF = aPrb.Focal();
2290       Standard_Real dV = 2.e5 * Sqrt(aF);
2291       V1new = Max(V1new, -dV);
2292       V2new = Min(V2new,  dV);
2293     }
2294     else if(aBasCurvType == GeomAbs_Hyperbola) {
2295      V1new  = Max(V1new, -30.);
2296      V2new  = Min(V2new,  30.);
2297     }
2298     else {
2299      V1new  = Max(V1new, -1.e10);
2300      V2new  = Min(V2new,  1.e10);
2301     }
2302
2303   }
2304   else {
2305
2306     V1new = Max(V1new, -1.e10);
2307     V2new = Min(V2new,  1.e10);
2308
2309   }
2310
2311 }
2312
2313 //================================================================================
2314 static void EstLimForInfSurf(Standard_Real& U1new,
2315                              Standard_Real& U2new,
2316                              Standard_Real& V1new,
2317                              Standard_Real& V2new)
2318 {
2319     U1new = Max(U1new, -1.e10);
2320     U2new = Min(U2new,  1.e10);
2321     V1new = Max(V1new, -1.e10);
2322     V2new = Min(V2new,  1.e10);
2323 }
2324
2325 #if 0
2326 //-- jgv patch (from)
2327 static Handle(Geom_Curve) GetCurve(const Handle(Adaptor3d_HCurve) AdCurve)
2328 {
2329   Handle(Geom_Curve) theCurve;
2330   GeomAbs_CurveType CurveType = AdCurve->GetType();
2331   switch (CurveType)
2332     {
2333     case GeomAbs_Line:
2334       theCurve = new Geom_Line( AdCurve->Line() );
2335       break;
2336     case GeomAbs_Circle:
2337       theCurve = new Geom_Circle( AdCurve->Circle() );
2338       break;
2339     case GeomAbs_Ellipse:
2340       theCurve = new Geom_Ellipse( AdCurve->Ellipse() );
2341       break;
2342     case GeomAbs_Hyperbola:
2343       theCurve = new Geom_Hyperbola( AdCurve->Hyperbola() );
2344       break;
2345     case GeomAbs_Parabola:
2346       theCurve = new Geom_Parabola( AdCurve->Parabola() );
2347       break;
2348     case GeomAbs_BezierCurve:
2349       theCurve = AdCurve->Bezier();
2350       break;
2351     case GeomAbs_BSplineCurve:
2352       theCurve = AdCurve->BSpline();
2353       break;
2354     }
2355   if (!theCurve.IsNull())
2356     {
2357       Standard_Real f = AdCurve->FirstParameter();
2358       Standard_Real l = AdCurve->LastParameter();
2359       theCurve = new Geom_TrimmedCurve( theCurve, f, l );
2360     }
2361   return theCurve;
2362 }
2363
2364 static Handle(Geom_Surface) GetSurface(const Handle(Adaptor3d_HSurface) AdSurface)
2365 {
2366   Handle(Geom_Surface) theSurface;
2367   GeomAbs_SurfaceType SurfaceType = AdSurface->GetType();
2368   switch (SurfaceType)
2369     {
2370     case GeomAbs_Plane:
2371       theSurface = new Geom_Plane( AdSurface->Plane() );
2372       break;
2373     case GeomAbs_Cylinder:
2374       theSurface = new Geom_CylindricalSurface( AdSurface->Cylinder() );
2375       break;
2376     case GeomAbs_Cone:
2377       theSurface = new Geom_ConicalSurface( AdSurface->Cone() );
2378       break;
2379     case GeomAbs_Torus:
2380       theSurface = new Geom_ToroidalSurface( AdSurface->Torus() );
2381       break;
2382     case GeomAbs_Sphere:
2383       theSurface = new Geom_SphericalSurface( AdSurface->Sphere() );
2384       break;
2385     case GeomAbs_BezierSurface:
2386       theSurface = AdSurface->Bezier();
2387       break;
2388     case GeomAbs_BSplineSurface:
2389       theSurface = AdSurface->BSpline();
2390       break;
2391     case GeomAbs_SurfaceOfRevolution:
2392       {
2393         gp_Ax1 Axis = AdSurface->AxeOfRevolution();
2394         Handle(Adaptor3d_HCurve) AdBC = AdSurface->BasisCurve();
2395         Handle(Geom_Curve) BC = GetCurve(AdBC);
2396         if (!BC.IsNull())
2397           theSurface = new Geom_SurfaceOfRevolution( BC, Axis );
2398         break;
2399       }
2400     case GeomAbs_SurfaceOfExtrusion:
2401       {
2402         gp_Dir Direction = AdSurface->Direction();
2403         Handle(Adaptor3d_HCurve) AdBC = AdSurface->BasisCurve();
2404         Handle(Geom_Curve) BC = GetCurve(AdBC);
2405         if (!BC.IsNull())
2406           theSurface = new Geom_SurfaceOfLinearExtrusion( BC, Direction );
2407         break;
2408       }
2409     }
2410   if (!theSurface.IsNull())
2411     {
2412       Standard_Real uf = AdSurface->FirstUParameter();
2413       Standard_Real ul = AdSurface->LastUParameter();
2414       Standard_Real vf = AdSurface->FirstVParameter();
2415       Standard_Real vl = AdSurface->LastVParameter();
2416       theSurface = new Geom_RectangularTrimmedSurface( theSurface, uf, ul, vf, vl );
2417     }
2418   return theSurface;
2419 }
2420 //-- jgv patch (to)
2421 #endif