src/IntCurveSurface/IntCurveSurface_Inter.gxx

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