src/ProjLib/ProjLib_ComputeApproxOnPolarSurface.cxx

   1 // Created by: Bruno DUMORTIER
   2 // Copyright (c) 1995-1999 Matra Datavision
   3 // Copyright (c) 1999-2012 OPEN CASCADE SAS
   4 //
   5 // The content of this file is subject to the Open CASCADE Technology Public
   6 // License Version 6.5 (the "License"). You may not use the content of this file
   7 // except in compliance with the License. Please obtain a copy of the License
   8 // at http://www.opencascade.org and read it completely before using this file.
   9 //
  10 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
  11 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
  12 //
  13 // The Original Code and all software distributed under the License is
  14 // distributed on an "AS IS" basis, without warranty of any kind, and the
  15 // Initial Developer hereby disclaims all such warranties, including without
  16 // limitation, any warranties of merchantability, fitness for a particular
  17 // purpose or non-infringement. Please see the License for the specific terms
  18 // and conditions governing the rights and limitations under the License.
  19
  20 #include <ProjLib_ComputeApproxOnPolarSurface.hxx>
  21 #include <AppCont_Function2d.hxx>
  22 #include <ElSLib.hxx>
  23 #include <ElCLib.hxx>
  24 #include <BSplCLib.hxx>
  25 #include <PLib.hxx>
  26 #include <Standard_NoSuchObject.hxx>
  27 #include <Geom_UndefinedDerivative.hxx>
  28 #include <gp_Trsf.hxx>
  29 #include <Precision.hxx>
  30 #include <Approx_FitAndDivide2d.hxx>
  31 #include <math.hxx>
  32 #include <AppParCurves_MultiCurve.hxx>
  33 #include <Geom_Surface.hxx>
  34 #include <Geom2d_BSplineCurve.hxx>
  35 #include <Geom2d_BezierCurve.hxx>
  36 #include <Geom2d_Line.hxx>
  37 #include <Geom2d_Circle.hxx>
  38 #include <Geom2d_Ellipse.hxx>
  39 #include <Geom2d_Hyperbola.hxx>
  40 #include <Geom2d_Parabola.hxx>
  41 #include <Geom2d_TrimmedCurve.hxx>
  42 #include <Geom_BSplineSurface.hxx>
  43 #include <Geom_BezierSurface.hxx>
  44 #include <Geom_BSplineCurve.hxx>
  45 #include <Geom_BezierCurve.hxx>
  46 #include <Geom_TrimmedCurve.hxx>
  47
  48 #include <TColgp_Array1OfPnt2d.hxx>
  49 #include <TColgp_Array2OfPnt2d.hxx>
  50 #include <TColgp_Array1OfPnt.hxx>
  51 #include <TColgp_SequenceOfPnt2d.hxx>
  52 #include <TColStd_Array1OfReal.hxx>
  53 #include <TColStd_Array1OfInteger.hxx>
  54 #include <TColStd_SequenceOfReal.hxx>
  55 #include <TColStd_ListOfTransient.hxx>
  56
  57 #include <GeomAbs_SurfaceType.hxx>
  58 #include <GeomAbs_CurveType.hxx>
  59 #include <Handle_Adaptor3d_HCurve.hxx>
  60 #include <Handle_Adaptor3d_HSurface.hxx>
  61 #include <Adaptor3d_Surface.hxx>
  62 #include <Adaptor3d_Curve.hxx>
  63 #include <Adaptor3d_HSurface.hxx>
  64 #include <Adaptor3d_HCurve.hxx>
  65 #include <Adaptor2d_HCurve2d.hxx>
  66 #include <Geom2dAdaptor_Curve.hxx>
  67 #include <Geom2dAdaptor_HCurve.hxx>
  68 #include <GeomAdaptor_HCurve.hxx>
  69 #include <GeomAdaptor.hxx>
  70 #include <GeomAdaptor_Surface.hxx>
  71 #include <TColgp_SequenceOfPnt.hxx>
  72
  73 #include <gp_Pnt.hxx>
  74 #include <gp_Pnt2d.hxx>
  75 #include <gp_Vec2d.hxx>
  76 #include <Extrema_GenLocateExtPS.hxx>
  77 #include <Extrema_ExtPS.hxx>
  78 #include <GCPnts_QuasiUniformAbscissa.hxx>
  79 #include <Standard_DomainError.hxx>
  80 //#include <GeomLib_IsIso.hxx>
  81 //#include <GeomLib_CheckSameParameter.hxx>
  82
  83 #ifdef DEB
  84 #ifdef DRAW
  85 #include <DrawTrSurf.hxx>
  86 #endif
  87 //static Standard_Integer compteur = 0;
  88 #endif
  89
  90 //=======================================================================
  91 //function : Value
  92 //purpose  : (OCC217 - apo)- Compute Point2d that project on polar surface(<Surf>) 3D<Curve>
  93 //            <InitCurve2d> use for calculate start 2D point.
  94 //=======================================================================
  95
  96 static gp_Pnt2d Function_Value(const Standard_Real U,
  97                                const Handle(Adaptor3d_HSurface)& Surf,
  98                                const Handle(Adaptor3d_HCurve)& Curve,
  99                                const Handle(Adaptor2d_HCurve2d)& InitCurve2d,
 100                                //OCC217
 101                                const Standard_Real DistTol3d, const Standard_Real tolU, const Standard_Real tolV)
 102                                //const Standard_Real Tolerance)
 103 {
 104   //OCC217
 105   //Standard_Real Tol3d = 100*Tolerance;
 106
 107   gp_Pnt2d p2d = InitCurve2d->Value(U) ;
 108   gp_Pnt p = Curve->Value(U);
 109 //  Curve->D0(U,p) ;
 110   Standard_Real Uinf, Usup, Vinf, Vsup;
 111   Uinf = Surf->Surface().FirstUParameter();
 112   Usup = Surf->Surface().LastUParameter();
 113   Vinf = Surf->Surface().FirstVParameter();
 114   Vsup = Surf->Surface().LastVParameter();
 115   Standard_Integer decalU = 0, decalV = 0;
 116   Standard_Real U0 = p2d.X(), V0 = p2d.Y();
 117
 118   GeomAbs_SurfaceType Type = Surf->GetType();
 119   if((Type != GeomAbs_BSplineSurface) &&
 120      (Type != GeomAbs_BezierSurface)  &&
 121      (Type != GeomAbs_OffsetSurface)    ) {
 122     Standard_Real S, T;
 123     switch (Type) {
 124 //    case GeomAbs_Plane:
 125 //      {
 126 //      gp_Pln Plane = Surf->Plane();
 127 //      ElSLib::Parameters( Plane, p, S, T);
 128 //      break;
 129 //      }
 130     case GeomAbs_Cylinder:
 131       {
 132         gp_Cylinder Cylinder = Surf->Cylinder();
 133         ElSLib::Parameters( Cylinder, p, S, T);
 134         if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
 135         if(U0 > Usup) decalU =  int((U0 - Usup)/(2*M_PI))+1;
 136         S += decalU*2*M_PI;
 137         break;
 138       }
 139     case GeomAbs_Cone:
 140       {
 141         gp_Cone Cone = Surf->Cone();
 142         ElSLib::Parameters( Cone, p, S, T);
 143         if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
 144         if(U0 > Usup) decalU =  int((U0 - Usup)/(2*M_PI))+1;
 145         S += decalU*2*M_PI;
 146         break;
 147       }
 148     case GeomAbs_Sphere:
 149       {
 150         gp_Sphere Sphere = Surf->Sphere();
 151         ElSLib::Parameters( Sphere, p, S, T);
 152         if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
 153         if(U0 > Usup) decalU =  int((U0 - Usup)/(2*M_PI))+1;
 154         S += decalU*2*M_PI;
 155         if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1;
 156         if(V0 > (Vsup+(Vsup-Vinf))) decalV =  int((V0 - Vsup+(Vsup-Vinf))/(2*M_PI))+1;
 157         T += decalV*2*M_PI;
 158         if(0.4*M_PI < Abs(U0 - S) && Abs(U0 - S) < 1.6*M_PI) {
 159           T = M_PI - T;
 160           if(U0 < S)
 161             S -= M_PI;
 162           else
 163             S += M_PI;
 164         }
 165         break;
 166       }
 167     case GeomAbs_Torus:
 168       {
 169         gp_Torus Torus = Surf->Torus();
 170         ElSLib::Parameters( Torus, p, S, T);
 171         if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
 172         if(U0 > Usup) decalU =  int((U0 - Usup)/(2*M_PI))+1;
 173         if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1;
 174         if(V0 > Vsup) decalV =  int((V0 - Vsup)/(2*M_PI))+1;
 175         S += decalU*2*M_PI; T += decalV*2*M_PI;
 176         break;
 177       }
 178     default:
 179       Standard_NoSuchObject::Raise("ProjLib_ComputeApproxOnPolarSurface::Value");
 180     }
 181     return gp_Pnt2d(S, T);
 182   }
 183
 184   //////////////////
 185   Standard_Real Dist2Min = RealLast();
 186   //OCC217
 187   //Standard_Real tolU,tolV ;
 188   //tolU = Tolerance;
 189   //tolV = Tolerance;
 190
 191   Standard_Real uperiod =0, vperiod = 0, u, v;
 192   // U0 and V0 are the points within the initialized period
 193   // (periode with u and v),
 194   // U1 and V1 are the points for construction of tops
 195
 196   if(Surf->IsUPeriodic() || Surf->IsUClosed()) {
 197     uperiod =  Surf->LastUParameter() - Surf->FirstUParameter();
 198   }
 199   if(Surf->IsVPeriodic() || Surf->IsVClosed()) {
 200     vperiod = Surf->LastVParameter() - Surf->FirstVParameter();
 201   }
 202   if(U0 < Uinf)
 203     if(!uperiod)
 204       U0 = Uinf;
 205     else {
 206       decalU = int((Uinf - U0)/uperiod)+1;
 207       U0 += decalU*uperiod;
 208     }
 209   if(U0 > Usup)
 210     if(!uperiod)
 211       U0 = Usup;
 212     else {
 213       decalU = -(int((U0 - Usup)/uperiod)+1);
 214       U0 += decalU*uperiod;
 215     }
 216   if(V0 < Vinf)
 217     if(!vperiod)
 218       V0 = Vinf;
 219     else {
 220       decalV = int((Vinf - V0)/vperiod)+1;
 221       V0 += decalV*vperiod;
 222     }
 223   if(V0 > Vsup)
 224     if(!vperiod)
 225       V0 = Vsup;
 226     else {
 227       decalV = -int((V0 - Vsup)/vperiod)-1;
 228       V0 += decalV*vperiod;
 229     }
 230
 231   // The surface around U0 is reduced
 232   Standard_Real uLittle = (Usup - Uinf)/10, vLittle = (Vsup - Vinf)/10;
 233   Standard_Real uInfLi = 0, vInfLi = 0,uSupLi = 0, vSupLi = 0;
 234   if((U0 - Uinf) > uLittle) uInfLi = U0 - uLittle; else uInfLi = Uinf;
 235   if((V0 - Vinf) > vLittle) vInfLi = V0 - vLittle; else vInfLi = Vinf;
 236   if((Usup - U0) > uLittle) uSupLi = U0 + uLittle; else uSupLi = Usup;
 237   if((Vsup - V0) > vLittle) vSupLi = V0 + vLittle; else vSupLi = Vsup;
 238
 239   //  const Adaptor3d_Surface GAS = Surf->Surface();
 240
 241   GeomAdaptor_Surface SurfLittle;
 242   if (Type == GeomAbs_BSplineSurface) {
 243     Handle(Geom_Surface) GBSS(Surf->Surface().BSpline());
 244     SurfLittle.Load(GBSS, uInfLi, uSupLi, vInfLi, vSupLi);
 245   }
 246   else if (Type == GeomAbs_BezierSurface) {
 247     Handle(Geom_Surface) GS(Surf->Surface().Bezier());
 248     SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi);
 249   }
 250   else if (Type == GeomAbs_OffsetSurface) {
 251     Handle(Geom_Surface) GS = GeomAdaptor::MakeSurface(Surf->Surface());
 252     SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi);
 253   }
 254   else {
 255     Standard_NoSuchObject::Raise("");
 256   }
 257
 258   Extrema_GenLocateExtPS  locext(p, SurfLittle, U0, V0, tolU, tolV);
 259   if (locext.IsDone()) {
 260     Dist2Min = locext.SquareDistance();
 261     if (Dist2Min < DistTol3d * DistTol3d) {
 262       (locext.Point()).Parameter(u, v);
 263       gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod);
 264       return pnt;
 265     }
 266   }
 267
 268   Extrema_ExtPS  ext(p, SurfLittle, tolU, tolV) ;
 269   if (ext.IsDone() && ext.NbExt()>=1 ) {
 270     Dist2Min = ext.SquareDistance(1);
 271     Standard_Integer GoodValue = 1;
 272     for ( Standard_Integer i = 2 ; i <= ext.NbExt() ; i++ )
 273       if( Dist2Min > ext.SquareDistance(i)) {
 274         Dist2Min = ext.SquareDistance(i);
 275         GoodValue = i;
 276       }
 277     if (Dist2Min < DistTol3d * DistTol3d) {
 278       (ext.Point(GoodValue)).Parameter(u,v);
 279       gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod);
 280       return pnt;
 281     }
 282   }
 283
 284   return p2d;
 285 }
 286
 287
 288 //=======================================================================
 289 //function : ProjLib_PolarFunction
 290 //purpose  : (OCC217 - apo)- This class produce interface to call "gp_Pnt2d Function_Value(...)"
 291 //=======================================================================
 292
 293 class ProjLib_PolarFunction : public AppCont_Function2d
 294 {
 295   Handle(Adaptor3d_HCurve)           myCurve;
 296   Handle(Adaptor2d_HCurve2d)         myInitialCurve2d ;
 297   Handle(Adaptor3d_HSurface)         mySurface ;
 298   //OCC217
 299   Standard_Real                      myTolU, myTolV;
 300   Standard_Real                      myDistTol3d;
 301   //Standard_Real                    myTolerance ;
 302
 303   public :
 304
 305   ProjLib_PolarFunction(const Handle(Adaptor3d_HCurve) & C,
 306                         const Handle(Adaptor3d_HSurface)& Surf,
 307                         const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
 308                         //OCC217
 309                         const Standard_Real Tol3d) :
 310                         //const Standard_Real Tolerance) :
 311   myCurve(C),
 312   myInitialCurve2d(InitialCurve2d),
 313   mySurface(Surf),
 314   //OCC217
 315   myTolU(Surf->UResolution(Tol3d)),
 316   myTolV(Surf->VResolution(Tol3d)),
 317   myDistTol3d(100.0*Tol3d){} ;
 318   //myTolerance(Tolerance){} ;
 319
 320   ~ProjLib_PolarFunction() {}
 321
 322   Standard_Real FirstParameter() const
 323   {return (myCurve->FirstParameter()/*+1.e-9*/);}
 324
 325   Standard_Real LastParameter() const
 326   {return (myCurve->LastParameter()/*-1.e-9*/);}
 327
 328   gp_Pnt2d Value(const Standard_Real t) const {
 329     return Function_Value
 330       (t,mySurface,myCurve,myInitialCurve2d,myDistTol3d,myTolU,myTolV) ;  //OCC217
 331       //(t,mySurface,myCurve,myInitialCurve2d,myTolerance) ;
 332   }
 333
 334 //  Standard_Boolean D1(const Standard_Real t, gp_Pnt2d& P, gp_Vec2d& V) const
 335   Standard_Boolean D1(const Standard_Real , gp_Pnt2d& , gp_Vec2d& ) const
 336     {return Standard_False;}
 337 };
 338
 339 //=======================================================================
 340 //function : ProjLib_ComputeApproxOnPolarSurface
 341 //purpose  :
 342 //=======================================================================
 343
 344 ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface(){}
 345
 346
 347 //=======================================================================
 348 //function : ProjLib_ComputeApproxOnPolarSurface
 349 //purpose  :
 350 //=======================================================================
 351
 352 ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface
 353 (const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
 354  const Handle(Adaptor3d_HCurve)& Curve,
 355  const Handle(Adaptor3d_HSurface)& S,
 356  const Standard_Real tol3d):myProjIsDone(Standard_False)  //OCC217
 357  //const Standard_Real tolerance):myProjIsDone(Standard_False)
 358 {
 359   myTolerance = tol3d; //OCC217
 360   //myTolerance = Max(Tolerance,Precision::PApproximation());
 361   myBSpline = Perform(InitialCurve2d,Curve,S);
 362 }
 363 //=======================================================================
 364 //function : ProjLib_ComputeApproxOnPolarSurface
 365 //purpose  : Process the case of sewing
 366 //=======================================================================
 367
 368 ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface
 369 (const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
 370  const Handle(Adaptor2d_HCurve2d)& InitialCurve2dBis,
 371  const Handle(Adaptor3d_HCurve)& Curve,
 372  const Handle(Adaptor3d_HSurface)& S,
 373  const Standard_Real tol3d):myProjIsDone(Standard_False)  //OCC217
 374  //const Standard_Real tolerance):myProjIsDone(Standard_False)
 375 {// InitialCurve2d and InitialCurve2dBis are two pcurves of the sewing
 376   myTolerance = tol3d; //OCC217
 377   //myTolerance = Max(tolerance,Precision::PApproximation());
 378   Handle(Geom2d_BSplineCurve) bsc = Perform(InitialCurve2d,Curve,S);
 379   if(myProjIsDone) {
 380     gp_Pnt2d P2dproj, P2d, P2dBis;
 381     P2dproj = bsc->StartPoint();
 382     P2d    = InitialCurve2d->Value(InitialCurve2d->FirstParameter());
 383     P2dBis = InitialCurve2dBis->Value(InitialCurve2dBis->FirstParameter());
 384
 385     Standard_Real Dist, DistBis;
 386     Dist    = P2dproj.Distance(P2d);
 387     DistBis = P2dproj.Distance(P2dBis);
 388     if( Dist < DistBis)  {
 389       // myBSpline2d is the pcurve that is found. It is translated to obtain myCurve2d
 390       myBSpline = bsc;
 391       Handle(Geom2d_Geometry) GG = myBSpline->Translated(P2d, P2dBis);
 392       my2ndCurve = Handle(Geom2d_Curve)::DownCast(GG);
 393     }
 394     else {
 395       my2ndCurve = bsc;
 396       Handle(Geom2d_Geometry) GG = my2ndCurve->Translated(P2dBis, P2d);
 397       myBSpline = Handle(Geom2d_BSplineCurve)::DownCast(GG);
 398     }
 399   }
 400 }
 401
 402 //=======================================================================
 403 //function : ProjLib_ComputeApproxOnPolarSurface
 404 //purpose  : case without curve of initialization
 405 //=======================================================================
 406
 407 ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface
 408 (const Handle(Adaptor3d_HCurve)& Curve,
 409  const Handle(Adaptor3d_HSurface)& S,
 410  const Standard_Real tol3d):myProjIsDone(Standard_False)  //OCC217
 411  //const Standard_Real tolerance):myProjIsDone(Standard_False)
 412 {
 413   myTolerance = tol3d; //OCC217
 414   //myTolerance = Max(tolerance,Precision::PApproximation());
 415   const Handle_Adaptor2d_HCurve2d InitCurve2d ;
 416   myBSpline = Perform(InitCurve2d,Curve,S);
 417 }
 418
 419 static Handle(Geom2d_BSplineCurve) Concat(Handle(Geom2d_BSplineCurve) C1,
 420                                           Handle(Geom2d_BSplineCurve) C2)
 421 {
 422   Standard_Integer deg, deg1, deg2;
 423   deg1 = C1->Degree();
 424   deg2 = C2->Degree();
 425
 426   if ( deg1 < deg2) {
 427     C1->IncreaseDegree(deg2);
 428     deg = deg2;
 429   }
 430   else if ( deg2 < deg1) {
 431     C2->IncreaseDegree(deg1);
 432     deg = deg1;
 433   }
 434   else deg = deg1;
 435
 436   Standard_Integer np1,np2,nk1,nk2,np,nk;
 437   np1 = C1->NbPoles();
 438   nk1 = C1->NbKnots();
 439   np2 = C2->NbPoles();
 440   nk2 = C2->NbKnots();
 441   nk = nk1 + nk2 -1;
 442   np = np1 + np2 -1;
 443
 444   TColStd_Array1OfReal    K1(1,nk1); C1->Knots(K1);
 445   TColStd_Array1OfInteger M1(1,nk1); C1->Multiplicities(M1);
 446   TColgp_Array1OfPnt2d    P1(1,np1); C1->Poles(P1);
 447   TColStd_Array1OfReal    K2(1,nk2); C2->Knots(K2);
 448   TColStd_Array1OfInteger M2(1,nk2); C2->Multiplicities(M2);
 449   TColgp_Array1OfPnt2d    P2(1,np2); C2->Poles(P2);
 450
 451   // Compute the new BSplineCurve
 452   TColStd_Array1OfReal    K(1,nk);
 453   TColStd_Array1OfInteger M(1,nk);
 454   TColgp_Array1OfPnt2d    P(1,np);
 455
 456   Standard_Integer i, count = 0;
 457   // Set Knots and Mults
 458   for ( i = 1; i <= nk1; i++) {
 459     count++;
 460     K(count) = K1(i);
 461     M(count) = M1(i);
 462   }
 463   M(count) = deg;
 464   for ( i = 2; i <= nk2; i++) {
 465     count++;
 466     K(count) = K2(i);
 467     M(count) = M2(i);
 468   }
 469   // Set the Poles
 470   count = 0;
 471   for (i = 1; i <= np1; i++) {
 472     count++;
 473     P(count) = P1(i);
 474   }
 475   for (i = 2; i <= np2; i++) {
 476     count++;
 477     P(count) = P2(i);
 478   }
 479
 480   Handle(Geom2d_BSplineCurve) BS =
 481     new Geom2d_BSplineCurve(P,K,M,deg);
 482   return BS;
 483 }
 484
 485
 486 //=======================================================================
 487 //function : Perform
 488 //purpose  :
 489 //=======================================================================
 490 Handle(Geom2d_BSplineCurve) ProjLib_ComputeApproxOnPolarSurface::Perform
 491 (const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
 492  const Handle(Adaptor3d_HCurve)& Curve,
 493  const Handle(Adaptor3d_HSurface)& S)
 494 {
 495   //OCC217
 496   Standard_Real Tol3d = myTolerance;
 497   Standard_Real ParamTol = Precision::PApproximation();
 498
 499   Handle(Adaptor2d_HCurve2d) AHC2d = InitialCurve2d;
 500   Handle(Adaptor3d_HCurve) AHC = Curve;
 501
 502 // if the curve 3d is a BSpline with degree C0, it is cut into sections with degree C1
 503 // -> bug cts18237
 504   GeomAbs_CurveType typeCurve = Curve->GetType();
 505   if(typeCurve == GeomAbs_BSplineCurve) {
 506     TColStd_ListOfTransient LOfBSpline2d;
 507     Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
 508     Standard_Integer nbInter = Curve->NbIntervals(GeomAbs_C1);
 509     if(nbInter > 1) {
 510       Standard_Integer i, j;
 511       Handle(Geom_TrimmedCurve) GTC;
 512       Handle(Geom2d_TrimmedCurve) G2dTC;
 513       TColStd_Array1OfReal Inter(1,nbInter+1);
 514       Curve->Intervals(Inter,GeomAbs_C1);
 515       Standard_Real firstinter = Inter.Value(1), secondinter = Inter.Value(2);
 516       // initialization 3d
 517       GTC = new Geom_TrimmedCurve(BSC, firstinter, secondinter);
 518       AHC = new GeomAdaptor_HCurve(GTC);
 519
 520       // if there is an initialization curve:
 521       // - either this is a BSpline C0, with discontinuity at the same parameters of nodes
 522       // and the sections C1 are taken
 523       // - or this is a curve C1 and the sections of intrest are taken otherwise the curve is created.
 524
 525       // initialization 2d
 526       Standard_Integer nbInter2d;
 527       Standard_Boolean C2dIsToCompute;
 528       C2dIsToCompute = InitialCurve2d.IsNull();
 529       Handle(Geom2d_BSplineCurve) BSC2d;
 530       Handle(Geom2d_Curve) G2dC;
 531
 532       if(!C2dIsToCompute) {
 533         nbInter2d = InitialCurve2d->NbIntervals(GeomAbs_C1);
 534         TColStd_Array1OfReal Inter2d(1,nbInter2d+1);
 535         InitialCurve2d->Intervals(Inter2d,GeomAbs_C1);
 536         j = 1;
 537         for(i = 1,j = 1;i <= nbInter;i++)
 538           if(Abs(Inter.Value(i) - Inter2d.Value(j)) < ParamTol) { //OCC217
 539           //if(Abs(Inter.Value(i) - Inter2d.Value(j)) < myTolerance) {
 540             if (j > nbInter2d) break;
 541             j++;
 542           }
 543         if(j != (nbInter2d+1)) {
 544           C2dIsToCompute = Standard_True;
 545         }
 546       }
 547
 548       if(C2dIsToCompute) {
 549         AHC2d = BuildInitialCurve2d(AHC, S);
 550       }
 551       else {
 552         typeCurve = InitialCurve2d->GetType();
 553         switch (typeCurve) {
 554         case GeomAbs_Line: {
 555           G2dC = new Geom2d_Line(InitialCurve2d->Line());
 556           break;
 557         }
 558         case GeomAbs_Circle: {
 559           G2dC = new Geom2d_Circle(InitialCurve2d->Circle());
 560           break;
 561         }
 562         case GeomAbs_Ellipse: {
 563           G2dC = new Geom2d_Ellipse(InitialCurve2d->Ellipse());
 564           break;
 565         }
 566         case GeomAbs_Hyperbola: {
 567           G2dC = new Geom2d_Hyperbola(InitialCurve2d->Hyperbola());
 568           break;
 569         }
 570         case GeomAbs_Parabola: {
 571           G2dC = new Geom2d_Parabola(InitialCurve2d->Parabola());
 572           break;
 573         }
 574         case GeomAbs_BezierCurve: {
 575           G2dC = InitialCurve2d->Bezier();
 576           break;
 577         }
 578         case GeomAbs_BSplineCurve: {
 579           G2dC = InitialCurve2d->BSpline();
 580           break;
 581         }
 582         case GeomAbs_OtherCurve:
 583         default:
 584           break;
 585         }
 586         gp_Pnt2d fp2d = G2dC->Value(firstinter), lp2d = G2dC->Value(secondinter);
 587         gp_Pnt fps, lps, fpc, lpc;
 588         S->D0(fp2d.X(), fp2d.Y(), fps);
 589         S->D0(lp2d.X(), lp2d.Y(), lps);
 590         Curve->D0(firstinter, fpc);
 591         Curve->D0(secondinter, lpc);
 592         //OCC217
 593         if((fps.IsEqual(fpc, Tol3d)) &&
 594            (lps.IsEqual(lpc, Tol3d))) {
 595         //if((fps.IsEqual(fpc, myTolerance)) &&
 596         //   (lps.IsEqual(lpc, myTolerance))) {
 597           G2dTC = new Geom2d_TrimmedCurve(G2dC, firstinter, secondinter);
 598           Geom2dAdaptor_Curve G2dAC(G2dTC);
 599           AHC2d = new Geom2dAdaptor_HCurve(G2dAC);
 600           myProjIsDone = Standard_True;
 601         }
 602         else {
 603           AHC2d = BuildInitialCurve2d(AHC, S);
 604           C2dIsToCompute = Standard_True;
 605         }
 606       }
 607
 608       if(myProjIsDone) {
 609         BSC2d = ProjectUsingInitialCurve2d(AHC, S, AHC2d);
 610         if(BSC2d.IsNull()) return Handle(Geom2d_BSplineCurve)(); //IFV
 611         LOfBSpline2d.Append(BSC2d);
 612       }
 613       else {
 614         return Handle(Geom2d_BSplineCurve)();
 615       }
 616
 617
 618
 619       Standard_Real iinter, ip1inter;
 620       Standard_Integer nbK2d, deg;
 621       nbK2d = BSC2d->NbKnots(); deg = BSC2d->Degree();
 622
 623       for(i = 2;i <= nbInter;i++) {
 624         iinter = Inter.Value(i);
 625         ip1inter = Inter.Value(i+1);
 626         // general case 3d
 627         GTC->SetTrim(iinter, ip1inter);
 628         AHC = new GeomAdaptor_HCurve(GTC);
 629
 630         // general case 2d
 631         if(C2dIsToCompute) {
 632           AHC2d = BuildInitialCurve2d(AHC, S);
 633         }
 634         else {
 635           gp_Pnt2d fp2d = G2dC->Value(iinter), lp2d = G2dC->Value(ip1inter);
 636           gp_Pnt fps, lps, fpc, lpc;
 637           S->D0(fp2d.X(), fp2d.Y(), fps);
 638           S->D0(lp2d.X(), lp2d.Y(), lps);
 639           Curve->D0(iinter, fpc);
 640           Curve->D0(ip1inter, lpc);
 641           //OCC217
 642           if((fps.IsEqual(fpc, Tol3d)) &&
 643              (lps.IsEqual(lpc, Tol3d))) {
 644           //if((fps.IsEqual(fpc, myTolerance)) &&
 645           //   (lps.IsEqual(lpc, myTolerance))) {
 646             G2dTC->SetTrim(iinter, ip1inter);
 647             Geom2dAdaptor_Curve G2dAC(G2dTC);
 648             AHC2d = new Geom2dAdaptor_HCurve(G2dAC);
 649             myProjIsDone = Standard_True;
 650           }
 651           else {
 652             AHC2d = BuildInitialCurve2d(AHC, S);
 653           }
 654         }
 655         if(myProjIsDone) {
 656           BSC2d = ProjectUsingInitialCurve2d(AHC, S, AHC2d);
 657           if(BSC2d.IsNull()) {
 658             return Handle(Geom2d_BSplineCurve)();
 659           }
 660           LOfBSpline2d.Append(BSC2d);
 661           nbK2d += BSC2d->NbKnots() - 1;
 662           deg = Max(deg, BSC2d->Degree());
 663         }
 664         else {
 665           return Handle(Geom2d_BSplineCurve)();
 666         }
 667       }
 668
 669       Standard_Integer NbC = LOfBSpline2d.Extent();
 670       Handle(Geom2d_BSplineCurve) CurBS;
 671       CurBS = Handle(Geom2d_BSplineCurve)::DownCast(LOfBSpline2d.First());
 672       LOfBSpline2d.RemoveFirst();
 673       for (Standard_Integer ii = 2; ii <= NbC; ii++) {
 674         Handle(Geom2d_BSplineCurve) BS =
 675           Handle(Geom2d_BSplineCurve)::DownCast(LOfBSpline2d.First());
 676         CurBS = Concat(CurBS,BS);
 677         LOfBSpline2d.RemoveFirst();
 678       }
 679       return CurBS;
 680     }
 681   }
 682
 683   if(InitialCurve2d.IsNull()) {
 684     AHC2d = BuildInitialCurve2d(Curve, S);
 685     if(!myProjIsDone)
 686       return Handle(Geom2d_BSplineCurve)();
 687   }
 688   return ProjectUsingInitialCurve2d(AHC, S, AHC2d);
 689 }
 690
 691 //=======================================================================
 692 //function : ProjLib_BuildInitialCurve2d
 693 //purpose  :
 694 //=======================================================================
 695
 696 Handle(Adaptor2d_HCurve2d)
 697      ProjLib_ComputeApproxOnPolarSurface::
 698      BuildInitialCurve2d(const Handle(Adaptor3d_HCurve)&   Curve,
 699                          const Handle(Adaptor3d_HSurface)& Surf)
 700 {
 701   //  discretize the Curve with quasiuniform deflection
 702   //  density at least NbOfPnts points
 703   myProjIsDone = Standard_False;
 704
 705   //OCC217
 706   Standard_Real Tol3d = myTolerance;
 707   Standard_Real TolU = Surf->UResolution(Tol3d), TolV = Surf->VResolution(Tol3d);
 708   Standard_Real DistTol3d = 100.0*Tol3d;
 709
 710   // NO myTol is Tol2d !!!!
 711   //Standard_Real TolU = myTolerance, TolV = myTolerance;
 712   //Standard_Real Tol3d = 100*myTolerance; // At random Balthazar.
 713
 714   Standard_Integer NbOfPnts = 61;
 715   GCPnts_QuasiUniformAbscissa QUA(Curve->GetCurve(),NbOfPnts);
 716   TColgp_Array1OfPnt Pts(1,NbOfPnts);
 717   TColStd_Array1OfReal Param(1,NbOfPnts);
 718   Standard_Integer i, j;
 719   for( i = 1; i <= NbOfPnts ; i++ ) {
 720     Param(i) = QUA.Parameter(i);
 721     Pts(i) = Curve->Value(Param(i));
 722   }
 723
 724   TColgp_Array1OfPnt2d Pts2d(1,NbOfPnts);
 725   TColStd_Array1OfInteger Mult(1,NbOfPnts);
 726   Mult.Init(1);
 727   Mult(1) = Mult(NbOfPnts) = 2;
 728
 729   Standard_Real Vinf, Vsup;
 730   Vinf = Surf->Surface().FirstVParameter();
 731   Vsup = Surf->Surface().LastVParameter();
 732   GeomAbs_SurfaceType Type = Surf->GetType();
 733   if((Type != GeomAbs_BSplineSurface) && (Type != GeomAbs_BezierSurface) &&
 734      (Type != GeomAbs_OffsetSurface)) {
 735     Standard_Real S, T;
 736 //    Standard_Integer usens = 0, vsens = 0;
 737     // to know the position relatively to the period
 738     switch (Type) {
 739 //    case GeomAbs_Plane:
 740 //      {
 741 //      gp_Pln Plane = Surf->Plane();
 742 //      for ( i = 1 ; i <= NbOfPnts ; i++) {
 743 //        ElSLib::Parameters( Plane, Pts(i), S, T);
 744 //        Pts2d(i).SetCoord(S,T);
 745 //      }
 746 //      myProjIsDone = Standard_True;
 747 //      break;
 748 //      }
 749     case GeomAbs_Cylinder:
 750       {
 751 //      Standard_Real Sloc, Tloc;
 752         Standard_Real Sloc;
 753         Standard_Integer usens = 0;
 754         gp_Cylinder Cylinder = Surf->Cylinder();
 755         ElSLib::Parameters( Cylinder, Pts(1), S, T);
 756         Pts2d(1).SetCoord(S,T);
 757         for ( i = 2 ; i <= NbOfPnts ; i++) {
 758           Sloc = S;
 759           ElSLib::Parameters( Cylinder, Pts(i), S, T);
 760           if(Abs(Sloc - S) > M_PI)
 761             if(Sloc > S)
 762               usens++;
 763             else
 764               usens--;
 765           Pts2d(i).SetCoord(S+usens*2*M_PI,T);
 766         }
 767         myProjIsDone = Standard_True;
 768         break;
 769       }
 770     case GeomAbs_Cone:
 771       {
 772 //      Standard_Real Sloc, Tloc;
 773         Standard_Real Sloc;
 774         Standard_Integer usens = 0;
 775         gp_Cone Cone = Surf->Cone();
 776         ElSLib::Parameters( Cone, Pts(1), S, T);
 777         Pts2d(1).SetCoord(S,T);
 778         for ( i = 2 ; i <= NbOfPnts ; i++) {
 779           Sloc = S;
 780           ElSLib::Parameters( Cone, Pts(i), S, T);
 781           if(Abs(Sloc - S) > M_PI)
 782             if(Sloc > S)
 783               usens++;
 784             else
 785               usens--;
 786           Pts2d(i).SetCoord(S+usens*2*M_PI,T);
 787         }
 788         myProjIsDone = Standard_True;
 789         break;
 790       }
 791     case GeomAbs_Sphere:
 792       {
 793         Standard_Real Sloc, Tloc;
 794         Standard_Integer usens = 0, vsens = 0; //usens steps by half-period
 795         Standard_Boolean vparit = Standard_False;
 796         gp_Sphere Sphere = Surf->Sphere();
 797         ElSLib::Parameters( Sphere, Pts(1), S, T);
 798         Pts2d(1).SetCoord(S,T);
 799         for ( i = 2 ; i <= NbOfPnts ; i++) {
 800           Sloc = S;Tloc = T;
 801           ElSLib::Parameters( Sphere, Pts(i), S, T);
 802           if(1.6*M_PI < Abs(Sloc - S))
 803             if(Sloc > S)
 804               usens += 2;
 805             else
 806               usens -= 2;
 807           if(1.6*M_PI > Abs(Sloc - S) && Abs(Sloc - S) > 0.4*M_PI) {
 808             vparit = !vparit;
 809             if(Sloc > S)
 810               usens++;
 811             else
 812               usens--;
 813             if(Abs(Tloc - Vsup) < (Vsup - Vinf)/5)
 814               vsens++;
 815             else
 816               vsens--;
 817           }
 818           if(vparit) {
 819             Pts2d(i).SetCoord(S+usens*M_PI,(M_PI - T)*(vsens-1));
 820           }
 821           else {
 822             Pts2d(i).SetCoord(S+usens*M_PI,T+vsens*M_PI);
 823
 824           }
 825         }
 826         myProjIsDone = Standard_True;
 827         break;
 828       }
 829     case GeomAbs_Torus:
 830       {
 831         Standard_Real Sloc, Tloc;
 832         Standard_Integer usens = 0, vsens = 0;
 833         gp_Torus Torus = Surf->Torus();
 834         ElSLib::Parameters( Torus, Pts(1), S, T);
 835         Pts2d(1).SetCoord(S,T);
 836         for ( i = 2 ; i <= NbOfPnts ; i++) {
 837           Sloc = S; Tloc = T;
 838           ElSLib::Parameters( Torus, Pts(i), S, T);
 839           if(Abs(Sloc - S) > M_PI)
 840             if(Sloc > S)
 841               usens++;
 842             else
 843               usens--;
 844           if(Abs(Tloc - T) > M_PI)
 845             if(Tloc > T)
 846               vsens++;
 847             else
 848               vsens--;
 849           Pts2d(i).SetCoord(S+usens*2*M_PI,T+vsens*2*M_PI);
 850         }
 851         myProjIsDone = Standard_True;
 852         break;
 853       }
 854     default:
 855       Standard_NoSuchObject::Raise("ProjLib_ComputeApproxOnPolarSurface::BuildInitialCurve2d");
 856     }
 857   }
 858   else {
 859     Standard_Real Uinf = Surf->Surface().FirstUParameter();
 860     Standard_Real Usup = Surf->Surface().LastUParameter();
 861
 862     myProjIsDone = Standard_False;
 863     Standard_Real Dist2Min = 1.e+200, u = 0., v = 0.;
 864     gp_Pnt pntproj;
 865
 866     TColgp_SequenceOfPnt2d Sols;
 867     Standard_Boolean areManyZeros = Standard_False;
 868
 869     Curve->D0(Param.Value(1), pntproj) ;
 870     Extrema_ExtPS  aExtPS(pntproj, Surf->Surface(), TolU, TolV) ;
 871     Standard_Real aMinSqDist = RealLast();
 872     if (aExtPS.IsDone())
 873     {
 874       for (i = 1; i <= aExtPS.NbExt(); i++)
 875       {
 876         Standard_Real aSqDist = aExtPS.SquareDistance(i);
 877         if (aSqDist < aMinSqDist)
 878           aMinSqDist = aSqDist;
 879       }
 880     }
 881     if (aMinSqDist > DistTol3d * DistTol3d) //try to project with less tolerance
 882     {
 883       TolU = Min(TolU, Precision::PConfusion());
 884       TolV = Min(TolV, Precision::PConfusion());
 885       aExtPS.Initialize(Surf->Surface(),
 886                         Surf->Surface().FirstUParameter(), Surf->Surface().LastUParameter(),
 887                         Surf->Surface().FirstVParameter(), Surf->Surface().LastVParameter(),
 888                         TolU, TolV);
 889       aExtPS.Perform(pntproj);
 890     }
 891
 892     if( aExtPS.IsDone() && aExtPS.NbExt() >= 1 ) {
 893
 894       Standard_Integer GoodValue = 1;
 895
 896       for ( i = 1 ; i <= aExtPS.NbExt() ; i++ ) {
 897         if( aExtPS.SquareDistance(i) < DistTol3d * DistTol3d ) {
 898           if( aExtPS.SquareDistance(i) <= 1.e-18 ) {
 899             aExtPS.Point(i).Parameter(u,v);
 900             gp_Pnt2d p2d(u,v);
 901             Standard_Boolean isSame = Standard_False;
 902             for( j = 1; j <= Sols.Length(); j++ ) {
 903               if( p2d.SquareDistance( Sols.Value(j) ) <= 1.e-18 ) {
 904                 isSame = Standard_True;
 905                 break;
 906               }
 907             }
 908             if( !isSame ) Sols.Append( p2d );
 909           }
 910           if( Dist2Min > aExtPS.SquareDistance(i) ) {
 911             Dist2Min = aExtPS.SquareDistance(i);
 912             GoodValue = i;
 913           }
 914         }
 915       }
 916
 917       if( Sols.Length() > 1 ) areManyZeros = Standard_True;
 918
 919       if( Dist2Min <= DistTol3d * DistTol3d) {
 920         if( !areManyZeros ) {
 921           aExtPS.Point(GoodValue).Parameter(u,v);
 922           Pts2d(1).SetCoord(u,v);
 923           myProjIsDone = Standard_True;
 924         }
 925         else {
 926           Standard_Integer nbSols = Sols.Length();
 927           Standard_Real Dist2Max = -1.e+200;
 928           for( i = 1; i <= nbSols; i++ ) {
 929             const gp_Pnt2d& aP1 = Sols.Value(i);
 930             for( j = i+1; j <= nbSols; j++ ) {
 931               const gp_Pnt2d& aP2 = Sols.Value(j);
 932               Standard_Real aDist2 = aP1.SquareDistance(aP2);
 933               if( aDist2 > Dist2Max ) Dist2Max = aDist2;
 934             }
 935           }
 936       Standard_Real aMaxT2 = Max(TolU,TolV);
 937       aMaxT2 *= aMaxT2;
 938           if( Dist2Max > aMaxT2 ) {
 939             Standard_Integer tPp = 0;
 940             for( i = 1; i <= 5; i++ ) {
 941               Standard_Integer nbExtOk = 0;
 942               Standard_Integer indExt = 0;
 943               Standard_Integer iT = 1 + (NbOfPnts - 1)/5*i;
 944               Curve->D0( Param.Value(iT), pntproj );
 945               Extrema_ExtPS aTPS( pntproj, Surf->Surface(), TolU, TolV );
 946               Dist2Min = 1.e+200;
 947               if( aTPS.IsDone() && aTPS.NbExt() >= 1 ) {
 948                 for( j = 1 ; j <= aTPS.NbExt() ; j++ ) {
 949                   if( aTPS.SquareDistance(j) < DistTol3d * DistTol3d ) {
 950                     nbExtOk++;
 951                     if( aTPS.SquareDistance(j) < Dist2Min ) {
 952                       Dist2Min = aTPS.SquareDistance(j);
 953                       indExt = j;
 954                     }
 955                   }
 956                 }
 957               }
 958               if( nbExtOk == 1 ) {
 959                 tPp = iT;
 960                 aTPS.Point(indExt).Parameter(u,v);
 961                 break;
 962               }
 963             }
 964
 965             if( tPp != 0 ) {
 966               gp_Pnt2d aPp = gp_Pnt2d(u,v);
 967               gp_Pnt2d aPn;
 968               j = 1;
 969               Standard_Boolean isFound = Standard_False;
 970               while( !isFound ) {
 971                 Curve->D0( Param.Value(tPp+j), pntproj );
 972                 Extrema_ExtPS aTPS( pntproj, Surf->Surface(), TolU, TolV );
 973                 Dist2Min = 1.e+200;
 974                 Standard_Integer indExt = 0;
 975                 if( aTPS.IsDone() && aTPS.NbExt() >= 1 ) {
 976                   for( i = 1 ; i <= aTPS.NbExt() ; i++ ) {
 977                     if( aTPS.SquareDistance(i) < DistTol3d * DistTol3d && aTPS.SquareDistance(i) < Dist2Min ) {
 978                       Dist2Min = aTPS.SquareDistance(i);
 979                       indExt = i;
 980                       isFound = Standard_True;
 981                     }
 982                   }
 983                 }
 984                 if( isFound ) {
 985                   aTPS.Point(indExt).Parameter(u,v);
 986                   aPn = gp_Pnt2d(u,v);
 987                   break;
 988                 }
 989                 j++;
 990                 if( (tPp+j) > NbOfPnts ) break;
 991               }
 992
 993               if( isFound ) {
 994                 gp_Vec2d atV(aPp,aPn);
 995                 Standard_Boolean isChosen = Standard_False;
 996                 for( i = 1; i <= nbSols; i++ ) {
 997                   const gp_Pnt2d& aP1 = Sols.Value(i);
 998                   gp_Vec2d asV(aP1,aPp);
 999                   if( asV.Dot(atV) > 0. ) {
1000                     isChosen = Standard_True;
1001                     Pts2d(1).SetCoord(aP1.X(),aP1.Y());
1002                     myProjIsDone = Standard_True;
1003                     break;
1004                   }
1005                 }
1006                 if( !isChosen ) {
1007                   aExtPS.Point(GoodValue).Parameter(u,v);
1008                   Pts2d(1).SetCoord(u,v);
1009                   myProjIsDone = Standard_True;
1010                 }
1011               }
1012               else {
1013                 aExtPS.Point(GoodValue).Parameter(u,v);
1014                 Pts2d(1).SetCoord(u,v);
1015                 myProjIsDone = Standard_True;
1016               }
1017             }
1018             else {
1019               aExtPS.Point(GoodValue).Parameter(u,v);
1020               Pts2d(1).SetCoord(u,v);
1021               myProjIsDone = Standard_True;
1022             }
1023           }
1024           else {
1025             aExtPS.Point(GoodValue).Parameter(u,v);
1026             Pts2d(1).SetCoord(u,v);
1027             myProjIsDone = Standard_True;
1028           }
1029         }
1030       }
1031
1032       //  calculate the following points with GenLocate_ExtPS
1033       // (and store the result and each parameter in a sequence)
1034       Standard_Integer usens = 0, vsens = 0;
1035       // to know the position relatively to the period
1036       Standard_Real U0 = u, V0 = v, U1 = u, V1 = v, uperiod =0, vperiod = 0;
1037       // U0 and V0 are the points in the initialized period
1038       // (period with u and v),
1039       // U1 and V1 are the points for construction of poles
1040
1041
1042       if(Surf->IsUPeriodic() || Surf->IsUClosed()) {
1043         uperiod = Surf->LastUParameter() - Surf->FirstUParameter();
1044       }
1045
1046       if(Surf->IsVPeriodic() || Surf->IsVClosed()) {
1047         vperiod = Surf->LastVParameter() - Surf->FirstVParameter();
1048       }
1049
1050       for ( i = 2 ; i <= NbOfPnts ; i++)
1051         if(myProjIsDone) {
1052           myProjIsDone = Standard_False;
1053           Dist2Min = RealLast();
1054           Curve->D0(Param.Value(i), pntproj);
1055           Extrema_GenLocateExtPS  aLocateExtPS
1056             (pntproj, Surf->Surface(), U0, V0, TolU, TolV) ;
1057
1058           if (aLocateExtPS.IsDone())
1059             if (aLocateExtPS.SquareDistance() < DistTol3d * DistTol3d) {  //OCC217
1060             //if (aLocateExtPS.SquareDistance() < Tol3d * Tol3d) {
1061               (aLocateExtPS.Point()).Parameter(U0,V0);
1062               U1 = U0 + usens*uperiod;
1063               V1 = V0 + vsens*vperiod;
1064               Pts2d(i).SetCoord(U1,V1);
1065               myProjIsDone = Standard_True;
1066             }
1067           if(!myProjIsDone && uperiod) {
1068             Standard_Real Uinf, Usup, Uaux;
1069             Uinf = Surf->Surface().FirstUParameter();
1070             Usup = Surf->Surface().LastUParameter();
1071             if((Usup - U0) > (U0 - Uinf))
1072               Uaux = 2*Uinf - U0 + uperiod;
1073             else
1074               Uaux = 2*Usup - U0 - uperiod;
1075             Extrema_GenLocateExtPS  locext(pntproj,
1076                                            Surf->Surface(),
1077                                            Uaux, V0, TolU, TolV);
1078             if (locext.IsDone())
1079               if (locext.SquareDistance() < DistTol3d * DistTol3d) {  //OCC217
1080               //if (locext.SquareDistance() < Tol3d * Tol3d) {
1081                 (locext.Point()).Parameter(u,v);
1082                 if((Usup - U0) > (U0 - Uinf))
1083                   usens--;
1084                 else
1085                   usens++;
1086                 U0 = u; V0 = v;
1087                 U1 = U0 + usens*uperiod;
1088                 V1 = V0 + vsens*vperiod;
1089                 Pts2d(i).SetCoord(U1,V1);
1090                 myProjIsDone = Standard_True;
1091               }
1092           }
1093           if(!myProjIsDone && vperiod) {
1094             Standard_Real Vinf, Vsup, Vaux;
1095             Vinf = Surf->Surface().FirstVParameter();
1096             Vsup = Surf->Surface().LastVParameter();
1097             if((Vsup - V0) > (V0 - Vinf))
1098               Vaux = 2*Vinf - V0 + vperiod;
1099             else
1100               Vaux = 2*Vsup - V0 - vperiod;
1101             Extrema_GenLocateExtPS  locext(pntproj,
1102                                            Surf->Surface(),
1103                                            U0, Vaux, TolU, TolV) ;
1104             if (locext.IsDone())
1105               if (locext.SquareDistance() < DistTol3d * DistTol3d) {  //OCC217
1106               //if (locext.SquareDistance() < Tol3d * Tol3d) {
1107                 (locext.Point()).Parameter(u,v);
1108                 if((Vsup - V0) > (V0 - Vinf))
1109                   vsens--;
1110                 else
1111                   vsens++;
1112                 U0 = u; V0 = v;
1113                 U1 = U0 + usens*uperiod;
1114                 V1 = V0 + vsens*vperiod;
1115                 Pts2d(i).SetCoord(U1,V1);
1116                 myProjIsDone = Standard_True;
1117               }
1118           }
1119           if(!myProjIsDone && uperiod && vperiod) {
1120             Standard_Real Uaux, Vaux;
1121             if((Usup - U0) > (U0 - Uinf))
1122               Uaux = 2*Uinf - U0 + uperiod;
1123             else
1124               Uaux = 2*Usup - U0 - uperiod;
1125             if((Vsup - V0) > (V0 - Vinf))
1126               Vaux = 2*Vinf - V0 + vperiod;
1127             else
1128               Vaux = 2*Vsup - V0 - vperiod;
1129             Extrema_GenLocateExtPS  locext(pntproj,
1130                                            Surf->Surface(),
1131                                            Uaux, Vaux, TolU, TolV);
1132             if (locext.IsDone())
1133               if (locext.SquareDistance() < DistTol3d * DistTol3d) {
1134               //if (locext.SquareDistance() < Tol3d * Tol3d) {
1135                 (locext.Point()).Parameter(u,v);
1136                 if((Usup - U0) > (U0 - Uinf))
1137                   usens--;
1138                 else
1139                   usens++;
1140                 if((Vsup - V0) > (V0 - Vinf))
1141                   vsens--;
1142                 else
1143                   vsens++;
1144                 U0 = u; V0 = v;
1145                 U1 = U0 + usens*uperiod;
1146                 V1 = V0 + vsens*vperiod;
1147                 Pts2d(i).SetCoord(U1,V1);
1148                 myProjIsDone = Standard_True;
1149               }
1150           }
1151           if(!myProjIsDone) {
1152             Extrema_ExtPS ext(pntproj, Surf->Surface(), TolU, TolV) ;
1153             if (ext.IsDone()) {
1154               Dist2Min = ext.SquareDistance(1);
1155               Standard_Integer GoodValue = 1;
1156               for ( j = 2 ; j <= ext.NbExt() ; j++ )
1157                 if( Dist2Min > ext.SquareDistance(j)) {
1158                   Dist2Min = ext.SquareDistance(j);
1159                   GoodValue = j;
1160                 }
1161               if (Dist2Min < DistTol3d * DistTol3d) {
1162               //if (Dist2Min < Tol3d * Tol3d) {
1163                 (ext.Point(GoodValue)).Parameter(u,v);
1164                 if(uperiod)
1165                   if((U0 - u) > (2*uperiod/3)) {
1166                     usens++;
1167                   }
1168                   else
1169                     if((u - U0) > (2*uperiod/3)) {
1170                       usens--;
1171                     }
1172                 if(vperiod)
1173                   if((V0 - v) > (vperiod/2)) {
1174                     vsens++;
1175                   }
1176                   else
1177                     if((v - V0) > (vperiod/2)) {
1178                       vsens--;
1179                     }
1180                 U0 = u; V0 = v;
1181                 U1 = U0 + usens*uperiod;
1182                 V1 = V0 + vsens*vperiod;
1183                 Pts2d(i).SetCoord(U1,V1);
1184                 myProjIsDone = Standard_True;
1185               }
1186             }
1187           }
1188         }
1189         else break;
1190     }
1191   }
1192   // -- Pnts2d is transformed into Geom2d_BSplineCurve, with the help of Param and Mult
1193   if(myProjIsDone) {
1194     myBSpline = new Geom2d_BSplineCurve(Pts2d,Param,Mult,1);
1195     Geom2dAdaptor_Curve GAC(myBSpline);
1196     Handle(Adaptor2d_HCurve2d) IC2d = new Geom2dAdaptor_HCurve(GAC);
1197 #ifdef DEB
1198 //    char name [100];
1199 //    sprintf(name,"%s_%d","build",compteur++);
1200 //    DrawTrSurf::Set(name,myBSpline);
1201 #endif
1202     return IC2d;
1203   }
1204   else {
1205 //  Modified by Sergey KHROMOV - Thu Apr 18 10:57:50 2002 Begin
1206 //     Standard_NoSuchObject_Raise_if(1,"ProjLib_Compu: build echec");
1207 //  Modified by Sergey KHROMOV - Thu Apr 18 10:57:51 2002 End
1208     return Handle(Adaptor2d_HCurve2d)();
1209   }
1210   myProjIsDone = Standard_False;
1211 //  Modified by Sergey KHROMOV - Thu Apr 18 10:58:01 2002 Begin
1212 //   Standard_NoSuchObject_Raise_if(1,"ProjLib_ComputeOnPS: build echec");
1213 //  Modified by Sergey KHROMOV - Thu Apr 18 10:58:02 2002 End
1214   return Handle(Adaptor2d_HCurve2d)();
1215 }
1216
1217
1218
1219
1220 //=======================================================================
1221 //function : ProjLib_ProjectUsingInitialCurve2d
1222 //purpose  :
1223 //=======================================================================
1224 Handle(Geom2d_BSplineCurve)
1225      ProjLib_ComputeApproxOnPolarSurface::
1226      ProjectUsingInitialCurve2d(const Handle(Adaptor3d_HCurve)& Curve,
1227                                 const Handle(Adaptor3d_HSurface)& Surf,
1228                                 const Handle(Adaptor2d_HCurve2d)& InitCurve2d)
1229 {
1230   //OCC217
1231   Standard_Real Tol3d = myTolerance;
1232   Standard_Real DistTol3d = 1.0*Tol3d;
1233   Standard_Real TolU = Surf->UResolution(Tol3d), TolV = Surf->VResolution(Tol3d);
1234   Standard_Real Tol2d = Sqrt(TolU*TolU + TolV*TolV);
1235
1236   Standard_Integer i;
1237   GeomAbs_SurfaceType TheTypeS = Surf->GetType();
1238   GeomAbs_CurveType TheTypeC = Curve->GetType();
1239 //  Handle(Standard_Type) TheTypeS = Surf->DynamicType();
1240 //  Handle(Standard_Type) TheTypeC = Curve->DynamicType();   // si on a :
1241 //  if(TheTypeS == STANDARD_TYPE(Geom_BSplineSurface)) {
1242   if(TheTypeS == GeomAbs_Plane) {
1243     Standard_Real S, T;
1244     gp_Pln Plane = Surf->Plane();
1245     if(TheTypeC == GeomAbs_BSplineCurve) {
1246       Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
1247       TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1248       for(i = 1;i <= Curve->NbPoles();i++) {
1249         ElSLib::Parameters( Plane, BSC->Pole(i), S, T);
1250         Poles2d(i).SetCoord(S,T);
1251       }
1252       TColStd_Array1OfReal Knots(1, BSC->NbKnots());
1253       BSC->Knots(Knots);
1254       TColStd_Array1OfInteger Mults(1, BSC->NbKnots());
1255       BSC->Multiplicities(Mults);
1256       if(BSC->IsRational()) {
1257         TColStd_Array1OfReal Weights(1, BSC->NbPoles());
1258         BSC->Weights(Weights);
1259         return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1260                                        BSC->Degree(), BSC->IsPeriodic()) ;
1261       }
1262       return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1263                                      BSC->Degree(), BSC->IsPeriodic()) ;
1264
1265     }
1266     if(TheTypeC == GeomAbs_BezierCurve) {
1267       Handle(Geom_BezierCurve) BC = Curve->Bezier();
1268       TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1269       for(i = 1;i <= Curve->NbPoles();i++) {
1270         ElSLib::Parameters( Plane, BC->Pole(i), S, T);
1271         Poles2d(i).SetCoord(S,T);
1272       }
1273       TColStd_Array1OfReal Knots(1, 2);
1274       Knots.SetValue(1,0.0);
1275       Knots.SetValue(2,1.0);
1276       TColStd_Array1OfInteger Mults(1, 2);
1277       Mults.Init(BC->NbPoles());
1278       if(BC->IsRational()) {
1279         TColStd_Array1OfReal Weights(1, BC->NbPoles());
1280         BC->Weights(Weights);
1281         return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1282                                        BC->Degree(), BC->IsPeriodic()) ;
1283       }
1284       return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1285                                      BC->Degree(), BC->IsPeriodic()) ;
1286     }
1287   }
1288   if(TheTypeS == GeomAbs_BSplineSurface) {
1289     Handle(Geom_BSplineSurface) BSS = Surf->BSpline();
1290     if((BSS->MaxDegree() == 1) &&
1291        (BSS->NbUPoles() == 2) &&
1292        (BSS->NbVPoles() == 2)) {
1293       gp_Pnt p11 = BSS->Pole(1,1);
1294       gp_Pnt p12 = BSS->Pole(1,2);
1295       gp_Pnt p21 = BSS->Pole(2,1);
1296       gp_Pnt p22 = BSS->Pole(2,2);
1297       gp_Vec V1(p11,p12);
1298       gp_Vec V2(p21,p22);
1299       if(V1.IsEqual(V2,Tol3d,Tol3d/(p11.Distance(p12)*180/M_PI))){  //OCC217
1300       //if(V1.IsEqual(V2,myTolerance,myTolerance/(p11.Distance(p12)*180/M_PI))){
1301         //  so the polar surface is plane
1302         //  and if it is enough to projet the  poles of Curve
1303         Standard_Integer Dist2Min = IntegerLast();
1304         Standard_Real u,v;
1305         //OCC217
1306         //Standard_Real TolU = Surf->UResolution(myTolerance)
1307         //  , TolV = Surf->VResolution(myTolerance);
1308 //      gp_Pnt pntproj;
1309         if(TheTypeC == GeomAbs_BSplineCurve) {
1310           Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
1311           TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1312           for(i = 1;i <= Curve->NbPoles();i++) {
1313             myProjIsDone = Standard_False;
1314             Dist2Min = IntegerLast();
1315             Extrema_GenLocateExtPS  extrloc(BSC->Pole(i),Surf->Surface(),(p11.X()+p22.X())/2,
1316                                             (p11.Y()+p22.Y())/2,TolU,TolV) ;
1317             if (extrloc.IsDone()) {
1318               Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1319               if (Dist2Min < DistTol3d * DistTol3d) {  //OCC217
1320               //if (Dist2Min < myTolerance * myTolerance) {
1321                 (extrloc.Point()).Parameter(u,v);
1322                 Poles2d(i).SetCoord(u,v);
1323                 myProjIsDone = Standard_True;
1324               }
1325               else break;
1326             }
1327             else break;
1328             if(!myProjIsDone)
1329               break;
1330           }
1331           if(myProjIsDone) {
1332             TColStd_Array1OfReal Knots(1, BSC->NbKnots());
1333             BSC->Knots(Knots);
1334             TColStd_Array1OfInteger Mults(1, BSC->NbKnots());
1335             BSC->Multiplicities(Mults);
1336             if(BSC->IsRational()) {
1337               TColStd_Array1OfReal Weights(1, BSC->NbPoles());
1338               BSC->Weights(Weights);
1339               return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1340                                              BSC->Degree(), BSC->IsPeriodic()) ;
1341             }
1342             return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1343                                            BSC->Degree(), BSC->IsPeriodic()) ;
1344
1345
1346           }
1347         }
1348         if(TheTypeC == GeomAbs_BezierCurve) {
1349           Handle(Geom_BezierCurve) BC = Curve->Bezier();
1350           TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1351           for(i = 1;i <= Curve->NbPoles();i++) {
1352             Dist2Min = IntegerLast();
1353             Extrema_GenLocateExtPS  extrloc(BC->Pole(i),Surf->Surface(),0.5,
1354                                             0.5,TolU,TolV) ;
1355             if (extrloc.IsDone()) {
1356               Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1357               if (Dist2Min < DistTol3d * DistTol3d) {  //OCC217
1358               //if (Dist2Min < myTolerance * myTolerance) {
1359                 (extrloc.Point()).Parameter(u,v);
1360                 Poles2d(i).SetCoord(u,v);
1361                 myProjIsDone = Standard_True;
1362               }
1363               else break;
1364             }
1365             else break;
1366             if(myProjIsDone)
1367               myProjIsDone = Standard_False;
1368             else break;
1369           }
1370           if(myProjIsDone) {
1371             TColStd_Array1OfReal Knots(1, 2);
1372             Knots.SetValue(1,0.0);
1373             Knots.SetValue(2,1.0);
1374             TColStd_Array1OfInteger Mults(1, 2);
1375             Mults.Init(BC->NbPoles());
1376             if(BC->IsRational()) {
1377               TColStd_Array1OfReal Weights(1, BC->NbPoles());
1378               BC->Weights(Weights);
1379               return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1380                                                     BC->Degree(), BC->IsPeriodic()) ;
1381             }
1382             return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1383                                                   BC->Degree(), BC->IsPeriodic()) ;
1384           }
1385         }
1386       }
1387     }
1388   }
1389   else if(TheTypeS == GeomAbs_BezierSurface) {
1390     Handle(Geom_BezierSurface) BS = Surf->Bezier();
1391     if((BS->MaxDegree() == 1) &&
1392        (BS->NbUPoles() == 2) &&
1393        (BS->NbVPoles() == 2)) {
1394       gp_Pnt p11 = BS->Pole(1,1);
1395       gp_Pnt p12 = BS->Pole(1,2);
1396       gp_Pnt p21 = BS->Pole(2,1);
1397       gp_Pnt p22 = BS->Pole(2,2);
1398       gp_Vec V1(p11,p12);
1399       gp_Vec V2(p21,p22);
1400       if(V1.IsEqual(V2,Tol3d,Tol3d/(p11.Distance(p12)*180/M_PI))){ //OCC217
1401       //if (V1.IsEqual(V2,myTolerance,myTolerance/(p11.Distance(p12)*180/M_PI))){
1402         //    and if it is enough to project the poles of Curve
1403         Standard_Integer Dist2Min = IntegerLast();
1404         Standard_Real u,v;
1405         //OCC217
1406         //Standard_Real TolU = Surf->UResolution(myTolerance)
1407         //  , TolV = Surf->VResolution(myTolerance);
1408
1409 //      gp_Pnt pntproj;
1410         if(TheTypeC == GeomAbs_BSplineCurve) {
1411           Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
1412           TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1413           for(i = 1;i <= Curve->NbPoles();i++) {
1414             myProjIsDone = Standard_False;
1415             Dist2Min = IntegerLast();
1416             Extrema_GenLocateExtPS  extrloc(BSC->Pole(i),Surf->Surface(),(p11.X()+p22.X())/2,
1417                                             (p11.Y()+p22.Y())/2,TolU,TolV) ;
1418             if (extrloc.IsDone()) {
1419               Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1420               if (Dist2Min < DistTol3d * DistTol3d) {  //OCC217
1421               //if (Dist2Min < myTolerance * myTolerance) {
1422                 (extrloc.Point()).Parameter(u,v);
1423                 Poles2d(i).SetCoord(u,v);
1424                 myProjIsDone = Standard_True;
1425               }
1426               else break;
1427             }
1428             else break;
1429             if(!myProjIsDone)
1430               break;
1431           }
1432           if(myProjIsDone) {
1433             TColStd_Array1OfReal Knots(1, BSC->NbKnots());
1434             BSC->Knots(Knots);
1435             TColStd_Array1OfInteger Mults(1, BSC->NbKnots());
1436             BSC->Multiplicities(Mults);
1437             if(BSC->IsRational()) {
1438               TColStd_Array1OfReal Weights(1, BSC->NbPoles());
1439               BSC->Weights(Weights);
1440               return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1441                                                     BSC->Degree(), BSC->IsPeriodic()) ;
1442             }
1443             return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1444                                                   BSC->Degree(), BSC->IsPeriodic()) ;
1445
1446
1447           }
1448         }
1449         if(TheTypeC == GeomAbs_BezierCurve) {
1450           Handle(Geom_BezierCurve) BC = Curve->Bezier();
1451           TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1452           for(i = 1;i <= Curve->NbPoles();i++) {
1453             Dist2Min = IntegerLast();
1454             Extrema_GenLocateExtPS  extrloc(BC->Pole(i),Surf->Surface(),0.5,
1455                                             0.5,TolU,TolV) ;
1456             if (extrloc.IsDone()) {
1457               Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1458               if (Dist2Min < DistTol3d * DistTol3d) {  //OCC217
1459               //if (Dist2Min < myTolerance * myTolerance) {
1460                 (extrloc.Point()).Parameter(u,v);
1461                 Poles2d(i).SetCoord(u,v);
1462                 myProjIsDone = Standard_True;
1463               }
1464               else break;
1465             }
1466             else break;
1467             if(myProjIsDone)
1468               myProjIsDone = Standard_False;
1469             else break;
1470           }
1471           if(myProjIsDone) {
1472             TColStd_Array1OfReal Knots(1, 2);
1473             Knots.SetValue(1,0.0);
1474             Knots.SetValue(2,1.0);
1475             TColStd_Array1OfInteger Mults(1, 2);
1476             Mults.Init(BC->NbPoles());
1477             if(BC->IsRational()) {
1478               TColStd_Array1OfReal Weights(1, BC->NbPoles());
1479               BC->Weights(Weights);
1480               return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1481                                                     BC->Degree(), BC->IsPeriodic()) ;
1482             }
1483             return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1484                                                   BC->Degree(), BC->IsPeriodic()) ;
1485           }
1486         }
1487       }
1488     }
1489   }
1490
1491   ProjLib_PolarFunction F(Curve, Surf, InitCurve2d, Tol3d) ;  //OCC217
1492   //ProjLib_PolarFunction F(Curve, Surf, InitCurve2d, myTolerance) ;
1493
1494 #ifdef DEB
1495   Standard_Integer Nb = 50;
1496
1497   Standard_Real U, U1, U2;
1498   U1 = F.FirstParameter();
1499   U2 = F.LastParameter();
1500
1501   TColgp_Array1OfPnt2d    DummyPoles(1,Nb+1);
1502   TColStd_Array1OfReal    DummyKnots(1,Nb+1);
1503   TColStd_Array1OfInteger DummyMults(1,Nb+1);
1504   DummyMults.Init(1);
1505   DummyMults(1) = 2;
1506   DummyMults(Nb+1) = 2;
1507   for (Standard_Integer ij = 0; ij <= Nb; ij++) {
1508     U = (Nb-ij)*U1 + ij*U2;
1509     U /= Nb;
1510     DummyPoles(ij+1) = F.Value(U);
1511     DummyKnots(ij+1) = ij;
1512   }
1513   Handle(Geom2d_BSplineCurve) DummyC2d =
1514     new Geom2d_BSplineCurve(DummyPoles, DummyKnots, DummyMults, 1);
1515   Standard_CString Temp = "bs2d";
1516 #ifdef DRAW
1517   DrawTrSurf::Set(Temp,DummyC2d);
1518 #endif
1519 //  DrawTrSurf::Set((Standard_CString ) "bs2d",DummyC2d);
1520   Handle(Geom2dAdaptor_HCurve) DDD =
1521     Handle(Geom2dAdaptor_HCurve)::DownCast(InitCurve2d);
1522
1523   Temp = "initc2d";
1524 #ifdef DRAW
1525   DrawTrSurf::Set(Temp,DDD->ChangeCurve2d().Curve());
1526 #endif
1527 //  DrawTrSurf::Set((Standard_CString ) "initc2d",DDD->ChangeCurve2d().Curve());
1528 #endif
1529
1530   Standard_Integer Deg1,Deg2;
1531 //  Deg1 = 8;
1532 //  Deg2 = 8;
1533   Deg1 = 2; //IFV
1534   Deg2 = 8; //IFV
1535
1536   Approx_FitAndDivide2d Fit(F,Deg1,Deg2,Tol3d,Tol2d, //OCC217
1537   //Approx_FitAndDivide2d Fit(F,Deg1,Deg2,myTolerance,myTolerance,
1538                             Standard_True);
1539
1540   if(Fit.IsAllApproximated()) {
1541     Standard_Integer i;
1542     Standard_Integer NbCurves = Fit.NbMultiCurves();
1543     Standard_Integer MaxDeg = 0;
1544     // To transform the MultiCurve into BSpline, it is required that all
1545     // Bezier constituing it have the same degree -> Calculation of MaxDeg
1546     Standard_Integer NbPoles  = 1;
1547     for (i = 1; i <= NbCurves; i++) {
1548       Standard_Integer Deg = Fit.Value(i).Degree();
1549       MaxDeg = Max ( MaxDeg, Deg);
1550     }
1551
1552     NbPoles = MaxDeg * NbCurves + 1;               //Tops on the BSpline
1553     TColgp_Array1OfPnt2d  Poles( 1, NbPoles);
1554
1555     TColgp_Array1OfPnt2d TempPoles( 1, MaxDeg + 1);//to augment the degree
1556
1557     TColStd_Array1OfReal Knots( 1, NbCurves + 1);  //Nodes of the BSpline
1558
1559     Standard_Integer Compt = 1;
1560     for (i = 1; i <= NbCurves; i++) {
1561       Fit.Parameters(i, Knots(i), Knots(i+1));
1562       AppParCurves_MultiCurve MC = Fit.Value( i);       //Load the Ith Curve
1563       TColgp_Array1OfPnt2d Poles2d( 1, MC.Degree() + 1);//Retrieve the tops
1564       MC.Curve(1, Poles2d);
1565
1566       //Eventual augmentation of the degree
1567       Standard_Integer Inc = MaxDeg - MC.Degree();
1568       if ( Inc > 0) {
1569 //      BSplCLib::IncreaseDegree( Inc, Poles2d, PLib::NoWeights(),
1570         BSplCLib::IncreaseDegree( MaxDeg, Poles2d, PLib::NoWeights(),
1571                          TempPoles, PLib::NoWeights());
1572         //update of tops of the PCurve
1573         for (Standard_Integer j = 1 ; j <= MaxDeg + 1; j++) {
1574           Poles.SetValue( Compt, TempPoles( j));
1575           Compt++;
1576         }
1577       }
1578       else {
1579         //update of tops of the PCurve
1580         for (Standard_Integer j = 1 ; j <= MaxDeg + 1; j++) {
1581           Poles.SetValue( Compt, Poles2d( j));
1582           Compt++;
1583         }
1584       }
1585
1586       Compt--;
1587     }
1588
1589     //update of fields of  ProjLib_Approx
1590     Standard_Integer NbKnots = NbCurves + 1;
1591
1592     // The start and end nodes are not correct : Cf: opening of the interval
1593     //Knots( 1) -= 1.e-9;
1594     //Knots(NbKnots) += 1.e-9;
1595
1596
1597     TColStd_Array1OfInteger   Mults( 1, NbKnots);
1598     Mults.Init(MaxDeg);
1599     Mults.SetValue( 1, MaxDeg + 1);
1600     Mults.SetValue(NbKnots, MaxDeg + 1);
1601     myProjIsDone = Standard_True;
1602     Handle(Geom2d_BSplineCurve) Dummy =
1603       new Geom2d_BSplineCurve(Poles,Knots,Mults,MaxDeg);
1604
1605     // try to smoother the Curve GeomAbs_C1.
1606
1607     Standard_Boolean OK = Standard_True;
1608
1609     for (Standard_Integer ij = 2; ij < NbKnots; ij++) {
1610       OK = OK && Dummy->RemoveKnot(ij,MaxDeg-1,Tol3d);  //OCC217
1611       //OK = OK && Dummy->RemoveKnot(ij,MaxDeg-1,myTolerance);
1612     }
1613 #ifdef DEB
1614     if (!OK) {
1615       cout << "ProjLib_ComputeApproxOnPolarSurface : Smoothing echoue"<<endl;
1616     }
1617 #endif
1618     return Dummy;
1619   }
1620   return Handle(Geom2d_BSplineCurve)();
1621 }
1622
1623 //=======================================================================
1624 //function : BSpline
1625 //purpose  :
1626 //=======================================================================
1627
1628 Handle(Geom2d_BSplineCurve)
1629      ProjLib_ComputeApproxOnPolarSurface::BSpline() const
1630
1631 {
1632 //  Modified by Sergey KHROMOV - Thu Apr 18 11:16:46 2002 End
1633 //   Standard_NoSuchObject_Raise_if
1634 //     (!myProjIsDone,
1635 //      "ProjLib_ComputeApproxOnPolarSurface:BSpline");
1636 //  Modified by Sergey KHROMOV - Thu Apr 18 11:16:47 2002 End
1637   return myBSpline ;
1638 }
1639
1640 //=======================================================================
1641 //function : Curve2d
1642 //purpose  :
1643 //=======================================================================
1644
1645 Handle(Geom2d_Curve)
1646      ProjLib_ComputeApproxOnPolarSurface::Curve2d() const
1647
1648 {
1649   Standard_NoSuchObject_Raise_if
1650     (!myProjIsDone,
1651      "ProjLib_ComputeApproxOnPolarSurface:2ndCurve2d");
1652   return my2ndCurve ;
1653 }
1654
1655
1656 //=======================================================================
1657 //function : IsDone
1658 //purpose  :
1659 //=======================================================================
1660
1661 Standard_Boolean ProjLib_ComputeApproxOnPolarSurface::IsDone() const
1662
1663 {
1664   return myProjIsDone;
1665 }