1 // Created on: 1993-09-07
2 // Created by: Bruno DUMORTIER
3 // Copyright (c) 1993-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
6 // This file is part of Open CASCADE Technology software library.
8 // This library is free software; you can redistribute it and / or modify it
9 // under the terms of the GNU Lesser General Public version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
17 // modified by NIZHNY-OFV Thu Jan 20 11:04:19 2005
19 #include <ProjLib_ComputeApprox.hxx>
21 #include <GeomAbs_SurfaceType.hxx>
22 #include <GeomAbs_CurveType.hxx>
23 #include <AppCont_Function2d.hxx>
24 #include <Convert_CompBezierCurves2dToBSplineCurve2d.hxx>
27 #include <BSplCLib.hxx>
28 #include <Standard_NoSuchObject.hxx>
29 #include <Geom_UndefinedDerivative.hxx>
31 #include <gp_Trsf.hxx>
32 #include <Precision.hxx>
33 #include <Approx_FitAndDivide2d.hxx>
34 #include <AppParCurves_MultiCurve.hxx>
35 #include <Handle_Adaptor3d_HCurve.hxx>
36 #include <Adaptor3d_HCurve.hxx>
37 #include <Handle_Adaptor3d_HSurface.hxx>
38 #include <Adaptor3d_HSurface.hxx>
39 #include <TColgp_Array1OfPnt2d.hxx>
40 #include <TColgp_Array1OfPnt.hxx>
41 #include <TColStd_Array1OfReal.hxx>
42 #include <TColStd_Array1OfInteger.hxx>
43 #include <Geom_BSplineCurve.hxx>
44 #include <Geom_BezierCurve.hxx>
45 #include <Geom2d_BSplineCurve.hxx>
46 #include <Geom2d_BezierCurve.hxx>
49 #include <DrawTrSurf.hxx>
52 static Standard_Boolean AffichValue = Standard_False;
56 void Parameters(const Handle(Adaptor3d_HCurve)& myCurve,
57 const Handle(Adaptor3d_HSurface)& mySurface,
59 const Standard_Integer iFirst,
60 const Standard_Real aTolU,
64 //=======================================================================
67 //=======================================================================
69 static inline Standard_Boolean IsEqual(Standard_Real Check,Standard_Real With,Standard_Real Toler)
71 return ((Abs(Check - With) < Toler) ? Standard_True : Standard_False);
75 //=======================================================================
78 //=======================================================================
80 static gp_Pnt2d Function_Value(const Standard_Real U,
81 const Handle(Adaptor3d_HCurve)& myCurve,
82 const Handle(Adaptor3d_HSurface)& mySurface,
83 const Standard_Real U1,
84 const Standard_Real U2,
85 const Standard_Real V1,
86 const Standard_Real V2,
87 const Standard_Boolean UCouture,
88 const Standard_Boolean VCouture )
90 Standard_Real S = 0., T = 0.;
92 gp_Pnt P3d = myCurve->Value(U);
93 GeomAbs_SurfaceType SType = mySurface->GetType();
99 gp_Pln Plane = mySurface->Plane();
100 ElSLib::Parameters( Plane, P3d, S, T);
103 case GeomAbs_Cylinder:
105 gp_Cylinder Cylinder = mySurface->Cylinder();
106 ElSLib::Parameters( Cylinder, P3d, S, T);
111 gp_Cone Cone = mySurface->Cone();
112 ElSLib::Parameters( Cone, P3d, S, T);
117 gp_Sphere Sphere = mySurface->Sphere();
118 ElSLib::Parameters(Sphere, P3d, S, T);
123 gp_Torus Torus = mySurface->Torus();
124 ElSLib::Parameters( Torus, P3d, S, T);
128 Standard_NoSuchObject::Raise("ProjLib_ComputeApprox::Value");
133 S = ElCLib::InPeriod(S, U1, U2);
137 if(SType == GeomAbs_Sphere) {
138 if ( Abs( S - U1 ) > M_PI ) {
143 S = ElCLib::InPeriod(S, U1, U2);
146 T = ElCLib::InPeriod(T, V1, V2);
149 return gp_Pnt2d(S, T);
151 //=======================================================================
154 //=======================================================================
155 static Standard_Boolean Function_D1( const Standard_Real U,
158 const Handle(Adaptor3d_HCurve)& myCurve,
159 const Handle(Adaptor3d_HSurface)& mySurface,
160 const Standard_Real U1,
161 const Standard_Real U2,
162 const Standard_Real V1,
163 const Standard_Real V2,
164 const Standard_Boolean UCouture,
165 const Standard_Boolean VCouture )
168 Standard_Real dU, dV;
170 P = Function_Value(U,myCurve,mySurface,U1,U2,V1,V2,UCouture,VCouture);
172 GeomAbs_SurfaceType Type = mySurface->GetType();
177 case GeomAbs_Cylinder:
183 myCurve->D1(U,P3d,T);
184 mySurface->D1(P.X(),P.Y(),P3d,D1U,D1V);
188 Standard_Real Nu = D1U.SquareMagnitude();
189 Standard_Real Nv = D1V.SquareMagnitude();
191 if ( Nu < Epsilon(1.) || Nv < Epsilon(1.))
192 return Standard_False;
196 D = gp_Vec2d( dU, dV);
201 return Standard_False;
204 return Standard_True;
207 //=======================================================================
208 //function : Function_SetUVBounds
210 //=======================================================================
211 static void Function_SetUVBounds(Standard_Real& myU1,
215 Standard_Boolean& UCouture,
216 Standard_Boolean& VCouture,
217 const Handle(Adaptor3d_HCurve)& myCurve,
218 const Handle(Adaptor3d_HSurface)& mySurface)
220 Standard_Real W1, W2, W;
223 W1 = myCurve->FirstParameter();
224 W2 = myCurve->LastParameter ();
226 // on ouvre l`intervalle
229 P1 = myCurve->Value(W1);
230 P2 = myCurve->Value(W2);
231 P = myCurve->Value(W);
233 switch ( mySurface->GetType()) {
236 gp_Cone Cone = mySurface->Cone();
237 VCouture = Standard_False;
239 switch( myCurve->GetType() ){
240 case GeomAbs_Parabola:
241 case GeomAbs_Hyperbola:
242 case GeomAbs_Ellipse:{
243 Standard_Real U1, U2, V1, V2, U , V;
244 ElSLib::Parameters( Cone, P1, U1, V1);
245 ElSLib::Parameters( Cone, P2, U2, V2);
246 ElSLib::Parameters( Cone, P , U , V );
249 if ( ( U1 < U && U < U2 ) && !myCurve->IsClosed() ) {
250 UCouture = Standard_False;
253 UCouture = Standard_True;
254 myU2 = myU1 + 2*M_PI;
260 Standard_Real U1, V1, U , V, Delta = 0., d = 0., pmin = W1, pmax = W1, dmax = 0., Uf, Ul;
261 ElSLib::Parameters( Cone, P1, U1, V1);
262 ElSLib::Parameters( Cone, P2, Ul, V1);
263 myU1 = U1; myU2 = U1; Uf = U1;
264 Standard_Real Step = .1;
265 Standard_Integer nbp = (Standard_Integer)((W2 - W1) / Step + 1);
267 Step = (W2 - W1) / (nbp - 1);
268 Standard_Boolean isclandper = (!(myCurve->IsClosed()) && !(myCurve->IsPeriodic()));
269 for(Standard_Real par = W1 + Step; par <= W2; par += Step) {
270 if(!isclandper) par += Step;
271 P = myCurve->Value(par);
272 ElSLib::Parameters( Cone, P, U, V);
276 if( ( (IsEqual(U,(2*M_PI),1.e-10) && (U1 >= 0. && U1 <= M_PI)) &&
277 (IsEqual(U,Ul,1.e-10) && !IsEqual(Uf,0.,1.e-10)) ) && isclandper ) U = 0.;
278 else Delta -= 2*M_PI;
283 if( ( (IsEqual(U,0.,1.e-10) && (U1 >= M_PI && U1 <= (2*M_PI))) &&
284 (IsEqual(U,Ul,1.e-10) && !IsEqual(Uf,(2*M_PI),1.e-10)) ) && isclandper ) U = 2*M_PI;
285 else Delta += 2*M_PI;
289 dmax = Max(dmax, Abs(d));
290 if(U < myU1) {myU1 = U; pmin = par;}
291 if(U > myU2) {myU2 = U; pmax = par;}
295 if(!(Abs(pmin - W1) <= Precision::PConfusion() || Abs(pmin - W2) <= Precision::PConfusion()) ) myU1 -= dmax*.5;
296 if(!(Abs(pmax - W1) <= Precision::PConfusion() || Abs(pmax - W2) <= Precision::PConfusion()) ) myU2 += dmax*.5;
298 if((myU1 >=0. && myU1 <= 2*M_PI) && (myU2 >=0. && myU2 <= 2*M_PI) ) UCouture = Standard_False;
300 U = ( myU1 + myU2 ) /2.;
303 UCouture = Standard_True;
307 }// switch curve type
311 case GeomAbs_Cylinder: {
312 gp_Cylinder Cylinder = mySurface->Cylinder();
313 VCouture = Standard_False;
315 if (myCurve->GetType() == GeomAbs_Ellipse) {
317 Standard_Real U1, U2, V1, V2, U , V;
318 ElSLib::Parameters( Cylinder, P1, U1, V1);
319 ElSLib::Parameters( Cylinder, P2, U2, V2);
320 ElSLib::Parameters( Cylinder, P , U , V );
324 if ( !myCurve->IsClosed()) {
325 if ( myU1 < U && U < myU2) {
326 U = ( myU1 + myU2 ) /2.;
331 U = ( myU1 + myU2 ) /2.;
341 UCouture = Standard_True;
347 myCurve->D1(W1,P3d,T);
348 mySurface->D1(U1,U2,P3d,D1U,D1V);
349 Standard_Real dU = T.Dot(D1U);
351 UCouture = Standard_True;
353 myU2 = myU1 + 2*M_PI;
362 Standard_Real U1, V1, U , V;
363 ElSLib::Parameters( Cylinder, P1, U1, V1);
364 Standard_Real Step = .1, Delta = 0.;
365 Standard_Real eps = M_PI, dmax = 0., d = 0.;
366 Standard_Integer nbp = (Standard_Integer)((W2 - W1) / Step + 1);
368 Step = (W2 - W1) / (nbp - 1);
369 myU1 = U1; myU2 = U1;
370 Standard_Real pmin = W1, pmax = W1, plim = W2+.1*Step;
371 for(Standard_Real par = W1 + Step; par <= plim; par += Step) {
372 P = myCurve->Value(par);
373 ElSLib::Parameters( Cylinder, P, U, V);
388 dmax = Max(dmax, Abs(d));
389 if(U < myU1) {myU1 = U; pmin = par;}
390 if(U > myU2) {myU2 = U; pmax = par;}
394 if(!(Abs(pmin - W1) <= Precision::PConfusion() ||
395 Abs(pmin - W2) <= Precision::PConfusion()) ) myU1 -= dmax*.5;
396 if(!(Abs(pmax - W1) <= Precision::PConfusion() ||
397 Abs(pmax - W2) <= Precision::PConfusion()) ) myU2 += dmax*.5;
399 if((myU1 >=0. && myU1 <= 2*M_PI) &&
400 (myU2 >=0. && myU2 <= 2*M_PI) ) {
401 UCouture = Standard_False;
404 U = ( myU1 + myU2 ) /2.;
407 UCouture = Standard_True;
413 case GeomAbs_Sphere:{
414 VCouture = Standard_False;
415 gp_Sphere SP = mySurface->Sphere();
416 if ( myCurve->GetType() == GeomAbs_Circle) {
417 UCouture = Standard_True;
419 // on cherche a savoir le nombre de fois que la couture est
421 // si 0 ou 2 fois : la PCurve est fermee et dans l`intervalle
422 // [Uc-PI, Uc+PI] (Uc: U du centre du cercle)
423 // si 1 fois : la PCurve est ouverte et dans l`intervalle
426 // pour determiner le nombre de solution, on resoud le systeme
427 // x^2 + y^2 + z^2 = R^2 (1)
428 // A x + B y + C z + D = 0 (2)
431 // REM : (1) (2) : equation du cercle
432 // (1) (3) (4) : equation de la couture.
433 Standard_Integer NbSolutions = 0;
434 Standard_Real A, B, C, D, R, Tol = 1.e-10;
435 Standard_Real U1, U2, V1, V2, aTPC;
438 aTPC=Precision::PConfusion();
440 gp_Circ Circle = myCurve->Circle();
441 Trsf.SetTransformation(SP.Position());
442 Circle.Transform(Trsf);
445 gp_Pln Plane( gp_Ax3(Circle.Position()));
446 Plane.Coefficients(A,B,C,D);
451 if ( ( R - Abs(D/A)) > Tol) NbSolutions = 2;
452 else if ( Abs(R - Abs(D/A))< Tol) NbSolutions = 1;
453 else NbSolutions = 0;
458 Standard_Real delta = R*R*(A*A+C*C) - D*D;
460 if ( Abs(delta) < Tol*Tol) {
461 if ( A*D > 0.) NbSolutions = 1;
463 else if ( delta > 0) {
468 if ( xx > Tol) NbSolutions++;
471 if ( xx > Tol) NbSolutions++;
477 Standard_Real UU = 0.;
478 ElSLib::Parameters(SP, P1, U1, V1);
479 Standard_Real eps = 2.*Epsilon(1.);
480 Standard_Real dt = Max(Precision::PConfusion(), 0.01*(W2-W1));
483 //May be U1 must be equal 2*PI?
484 gp_Pnt Pd = myCurve->Value(W1+dt);
485 Standard_Real ud, vd;
486 ElSLib::Parameters(SP, Pd, ud, vd);
487 if(Abs(U1 - ud) > M_PI)
492 else if(Abs(2.*M_PI - U1) < eps)
495 gp_Pnt Pd = myCurve->Value(W1+dt);
496 Standard_Real ud, vd;
497 ElSLib::Parameters(SP, Pd, ud, vd);
498 if(Abs(U1 - ud) > M_PI)
504 ElSLib::Parameters(SP, P2, U2, V1);
507 //May be U2 must be equal 2*PI?
508 gp_Pnt Pd = myCurve->Value(W2-dt);
509 Standard_Real ud, vd;
510 ElSLib::Parameters(SP, Pd, ud, vd);
511 if(Abs(U2 - ud) > M_PI)
516 else if(Abs(2.*M_PI - U2) < eps)
519 gp_Pnt Pd = myCurve->Value(W2-dt);
520 Standard_Real ud, vd;
521 ElSLib::Parameters(SP, Pd, ud, vd);
522 if(Abs(U2 - ud) > M_PI)
528 ElSLib::Parameters(SP, P, UU, V1);
529 Standard_Real UUmi = Min(Min(U1,UU),Min(UU,U2));
530 Standard_Real UUma = Max(Max(U1,UU),Max(UU,U2));
531 Standard_Boolean reCalc = ((UUmi >= 0. && UUmi <= M_PI) && (UUma >= 0. && UUma <= M_PI));
533 P2 = myCurve->Value(W1+M_PI/8);
534 ElSLib::Parameters(SP,P2,U2,V2);
536 if ( NbSolutions == 1) {
537 if ( Abs(U1-U2) > M_PI) { // on traverse la couture
547 else { // on ne traverse pas la couture
558 else { // 0 ou 2 solutions
559 gp_Pnt Center = Circle.Location();
561 ElSLib::SphereParameters(gp_Ax3(gp::XOY()),1,Center, U, V);
566 // eval the VCouture.
567 if ( (C==0) || Abs(Abs(D/C)-R) > 1.e-10) {
568 VCouture = Standard_False;
571 VCouture = Standard_True;
572 UCouture = Standard_True;
576 myV2 = 3 * M_PI / 2.;
579 myV1 = -3 * M_PI / 2.;
583 // si P1.Z() vaut +/- R on est sur le sommet : pas significatif.
584 gp_Pnt pp = P1.Transformed(Trsf);
586 if ( Abs( Abs(pp.Z()) - R) < Tol) {
587 gp_Pnt Center = Circle.Location();
589 ElSLib::SphereParameters(gp_Ax3(gp::XOY()),1,Center, U, V);
592 VCouture = Standard_False;
597 myV1 = -1.e+100; myV2 = 1.e+100;
598 Standard_Real UU1 = myU1, UU2 = myU2;
599 if((Abs(UU1) <= (2.*M_PI) && Abs(UU2) <= (2.*M_PI)) && NbSolutions == 1 && reCalc) {
600 gp_Pnt Center = Circle.Location();
602 ElSLib::SphereParameters(gp_Ax3(gp::XOY()),1,Center, U, V);
604 myU1 = Min(UU1,myU1);
605 myU2 = myU1 + 2.*M_PI;
609 }//if ( myCurve->GetType() == GeomAbs_Circle)
612 Standard_Real U1, V1, U , V;
613 ElSLib::Parameters( SP, P1, U1, V1);
614 Standard_Real Step = .1, Delta = 0.;
615 Standard_Real eps = M_PI, dmax = 0., d = 0.;
616 Standard_Integer nbp = (Standard_Integer)((W2 - W1) / Step + 1);
618 Step = (W2 - W1) / (nbp - 1);
619 myU1 = U1; myU2 = U1;
620 Standard_Real pmin = W1, pmax = W1, plim = W2+.1*Step;
621 for(Standard_Real par = W1 + Step; par <= plim; par += Step) {
622 P = myCurve->Value(par);
623 ElSLib::Parameters( SP, P, U, V);
638 dmax = Max(dmax, Abs(d));
639 if(U < myU1) {myU1 = U; pmin = par;}
640 if(U > myU2) {myU2 = U; pmax = par;}
644 if(!(Abs(pmin - W1) <= Precision::PConfusion() ||
645 Abs(pmin - W2) <= Precision::PConfusion()) ) myU1 -= dmax*.5;
646 if(!(Abs(pmax - W1) <= Precision::PConfusion() ||
647 Abs(pmax - W2) <= Precision::PConfusion()) ) myU2 += dmax*.5;
649 if((myU1 >=0. && myU1 <= 2*M_PI) &&
650 (myU2 >=0. && myU2 <= 2*M_PI) ) {
653 UCouture = Standard_False;
656 U = ( myU1 + myU2 ) /2.;
659 UCouture = Standard_True;
662 VCouture = Standard_False;
668 gp_Torus TR = mySurface->Torus();
669 Standard_Real U1, V1, U , V;
670 ElSLib::Parameters( TR, P1, U1, V1);
671 Standard_Real Step = .1, DeltaU = 0., DeltaV = 0.;
672 Standard_Real eps = M_PI, dmaxU = 0., dU = 0., dmaxV = 0., dV = 0.;
673 Standard_Integer nbp = (Standard_Integer)((W2 - W1) / Step + 1);
675 Step = (W2 - W1) / (nbp - 1);
676 myU1 = U1; myU2 = U1;
677 myV1 = V1; myV2 = V1;
678 Standard_Real pminU = W1, pmaxU = W1, pminV = W1, pmaxV = W1,
680 for(Standard_Real par = W1 + Step; par <= plim; par += Step) {
681 P = myCurve->Value(par);
682 ElSLib::Parameters( TR, P, U, V);
711 dmaxU = Max(dmaxU, Abs(dU));
712 dmaxV = Max(dmaxV, Abs(dV));
713 if(U < myU1) {myU1 = U; pminU = par;}
714 if(U > myU2) {myU2 = U; pmaxU = par;}
715 if(V < myV1) {myV1 = V; pminV = par;}
716 if(V > myV2) {myV2 = V; pmaxV = par;}
721 if(!(Abs(pminU - W1) <= Precision::PConfusion() ||
722 Abs(pminU - W2) <= Precision::PConfusion()) ) myU1 -= dmaxU*.5;
723 if(!(Abs(pmaxU - W1) <= Precision::PConfusion() ||
724 Abs(pmaxU - W2) <= Precision::PConfusion()) ) myU2 += dmaxU*.5;
725 if(!(Abs(pminV - W1) <= Precision::PConfusion() ||
726 Abs(pminV - W2) <= Precision::PConfusion()) ) myV1 -= dmaxV*.5;
727 if(!(Abs(pmaxV - W1) <= Precision::PConfusion() ||
728 Abs(pmaxV - W2) <= Precision::PConfusion()) ) myV2 += dmaxV*.5;
730 if((myU1 >=0. && myU1 <= 2*M_PI) &&
731 (myU2 >=0. && myU2 <= 2*M_PI) ) {
734 UCouture = Standard_False;
737 U = ( myU1 + myU2 ) /2.;
740 UCouture = Standard_True;
742 if((myV1 >=0. && myV1 <= 2*M_PI) &&
743 (myV2 >=0. && myV2 <= 2*M_PI) ) {
744 VCouture = Standard_False;
747 V = ( myV1 + myV2 ) /2.;
750 VCouture = Standard_True;
758 UCouture = Standard_False;
759 VCouture = Standard_False;
766 //=======================================================================
767 //classn : ProjLib_Function
769 //=======================================================================
770 class ProjLib_Function : public AppCont_Function2d
772 Handle(Adaptor3d_HCurve) myCurve;
773 Handle(Adaptor3d_HSurface) mySurface;
777 Standard_Real myU1,myU2,myV1,myV2;
778 Standard_Boolean UCouture,VCouture;
780 ProjLib_Function(const Handle(Adaptor3d_HCurve)& C,
781 const Handle(Adaptor3d_HSurface)& S) :
782 myCurve(C), mySurface(S),
787 UCouture(Standard_False),
788 VCouture(Standard_False)
789 {Function_SetUVBounds(myU1,myU2,myV1,myV2,UCouture,VCouture,myCurve,mySurface);}
791 Standard_Real FirstParameter() const
792 {return (myCurve->FirstParameter() + 1.e-9);}
794 Standard_Real LastParameter() const
795 {return (myCurve->LastParameter() -1.e-9);}
798 gp_Pnt2d Value(const Standard_Real t) const
799 {return Function_Value(t,myCurve,mySurface,myU1,myU2,myV1,myV2,UCouture,VCouture);}
801 Standard_Boolean D1(const Standard_Real t, gp_Pnt2d& P, gp_Vec2d& V) const
802 {return Function_D1(t,P,V,myCurve,mySurface,myU1,myU2,myV1,myV2,UCouture,VCouture);}
805 //=======================================================================
806 //function : ProjLib_ComputeApprox
808 //=======================================================================
810 ProjLib_ComputeApprox::ProjLib_ComputeApprox
811 (const Handle(Adaptor3d_HCurve) & C,
812 const Handle(Adaptor3d_HSurface) & S,
813 const Standard_Real Tol )
815 // if the surface is a plane and the curve a BSpline or a BezierCurve,
816 // don`t make an Approx but only the projection of the poles.
818 myTolerance = Max(Precision::PApproximation(),Tol);
819 Standard_Integer NbKnots, NbPoles ;
820 GeomAbs_CurveType CType = C->GetType();
821 GeomAbs_SurfaceType SType = S->GetType();
823 Standard_Boolean SurfIsAnal = (SType != GeomAbs_BSplineSurface) &&
824 (SType != GeomAbs_BezierSurface) &&
825 (SType != GeomAbs_OtherSurface) ;
827 Standard_Boolean CurvIsAnal = (CType != GeomAbs_BSplineCurve) &&
828 (CType != GeomAbs_BezierCurve) &&
829 (CType != GeomAbs_OtherCurve) ;
831 Standard_Boolean simplecase = SurfIsAnal && CurvIsAnal;
833 if (CType == GeomAbs_BSplineCurve &&
834 SType == GeomAbs_Plane ) {
836 // get the poles and eventually the weights
837 Handle(Geom_BSplineCurve) BS = C->BSpline();
838 NbPoles = BS->NbPoles();
839 TColgp_Array1OfPnt P3d( 1, NbPoles);
840 TColgp_Array1OfPnt2d Poles( 1, NbPoles);
841 TColStd_Array1OfReal Weights( 1, NbPoles);
842 if ( BS->IsRational()) BS->Weights(Weights);
844 gp_Pln Plane = S->Plane();
846 for ( Standard_Integer i = 1; i <= NbPoles; i++) {
847 ElSLib::Parameters( Plane, P3d(i), U, V);
848 Poles.SetValue(i,gp_Pnt2d(U,V));
850 NbKnots = BS->NbKnots();
851 TColStd_Array1OfReal Knots(1,NbKnots);
852 TColStd_Array1OfInteger Mults(1,NbKnots);
854 BS->Multiplicities(Mults) ;
855 // get the knots and mults if BSplineCurve
856 if ( BS->IsRational()) {
857 myBSpline = new Geom2d_BSplineCurve(Poles,
865 myBSpline = new Geom2d_BSplineCurve(Poles,
872 else if (CType == GeomAbs_BezierCurve &&
873 SType == GeomAbs_Plane ) {
875 // get the poles and eventually the weights
876 Handle(Geom_BezierCurve) BezierCurvePtr = C->Bezier() ;
877 NbPoles = BezierCurvePtr->NbPoles();
878 TColgp_Array1OfPnt P3d( 1, NbPoles);
879 TColgp_Array1OfPnt2d Poles( 1, NbPoles);
880 TColStd_Array1OfReal Weights( 1, NbPoles);
881 if ( BezierCurvePtr->IsRational()) {
882 BezierCurvePtr->Weights(Weights);
884 BezierCurvePtr->Poles( P3d);
886 // project the 3D-Poles on the plane
888 gp_Pln Plane = S->Plane();
890 for ( Standard_Integer i = 1; i <= NbPoles; i++) {
891 ElSLib::Parameters( Plane, P3d(i), U, V);
892 Poles.SetValue(i,gp_Pnt2d(U,V));
894 if ( BezierCurvePtr->IsRational()) {
895 myBezier = new Geom2d_BezierCurve(Poles, Weights);
898 myBezier = new Geom2d_BezierCurve(Poles);
902 ProjLib_Function F( C, S);
906 Standard_Integer Nb = 20;
907 Standard_Real U1, U2, dU, U;
908 U1 = F.FirstParameter();
909 U2 = F.LastParameter();
910 dU = ( U2 - U1) / Nb;
911 TColStd_Array1OfInteger Mults(1,Nb+1);
912 TColStd_Array1OfReal Knots(1,Nb+1);
913 TColgp_Array1OfPnt2d Poles(1,Nb+1);
914 for ( Standard_Integer i = 1; i <= Nb+1; i++) {
916 Poles(i) = F.Value(U);
924 char* ResultName = "Result";
925 DrawTrSurf::Set(ResultName,new Geom2d_BSplineCurve(Poles,Knots,Mults,1));
926 // DrawTrSurf::Set("Result",new Geom2d_BSplineCurve(Poles,Knots,Mults,1));
932 Standard_Integer Deg1, Deg2;
942 Approx_FitAndDivide2d Fit(F,Deg1,Deg2,myTolerance,myTolerance,
944 if(Fit.IsAllApproximated()) {
946 Standard_Integer NbCurves = Fit.NbMultiCurves();
948 // on essaie de rendre la courbe au moins C1
949 Convert_CompBezierCurves2dToBSplineCurve2d Conv;
952 Standard_Real Tol3d,Tol2d;
953 for (i = 1; i <= NbCurves; i++) {
954 Fit.Error(i,Tol3d, Tol2d);
955 myTolerance = Max(myTolerance, Tol2d);
956 AppParCurves_MultiCurve MC = Fit.Value( i); //Charge la Ieme Curve
957 TColgp_Array1OfPnt2d Poles2d( 1, MC.Degree() + 1);//Recupere les poles
958 MC.Curve(1, Poles2d);
960 Conv.AddCurve(Poles2d);
963 //mise a jour des fields de ProjLib_Approx
966 NbPoles = Conv.NbPoles();
967 NbKnots = Conv.NbKnots();
970 if(NbPoles <= 0 || NbPoles > 100000)
972 if(NbKnots <= 0 || NbKnots > 100000)
975 TColgp_Array1OfPnt2d NewPoles(1,NbPoles);
976 TColStd_Array1OfReal NewKnots(1,NbKnots);
977 TColStd_Array1OfInteger NewMults(1,NbKnots);
979 Conv.KnotsAndMults(NewKnots,NewMults);
980 Conv.Poles(NewPoles);
982 BSplCLib::Reparametrize(C->FirstParameter(),
986 // il faut recadrer les poles de debut et de fin:
987 // ( Car pour les problemes de couture, on a du ouvrir l`intervalle
988 // de definition de la courbe.)
989 // On choisit de calculer ces poles par prolongement de la courbe
995 U = C->FirstParameter() - 1.e-9;
1001 BSplCLib::NoWeights(),
1005 NewPoles.SetValue(1,P);
1006 U = C->LastParameter() + 1.e-9;
1012 BSplCLib::NoWeights(),
1016 NewPoles.SetValue(NbPoles,P);
1017 myBSpline = new Geom2d_BSplineCurve (NewPoles,
1023 Standard_Integer NbCurves = Fit.NbMultiCurves();
1025 Standard_Real Tol3d,Tol2d;
1026 Fit.Error(NbCurves,Tol3d, Tol2d);
1027 myTolerance = Tol2d;
1032 Standard_Real UFirst = F.FirstParameter();
1033 gp_Pnt P3d = C->Value( UFirst );
1034 Standard_Real u = 0., v = 0.;
1039 gp_Pln Plane = S->Plane();
1040 ElSLib::Parameters( Plane, P3d, u, v );
1043 case GeomAbs_Cylinder:
1045 gp_Cylinder Cylinder = S->Cylinder();
1046 ElSLib::Parameters( Cylinder, P3d, u, v );
1051 gp_Cone Cone = S->Cone();
1052 ElSLib::Parameters( Cone, P3d, u, v );
1055 case GeomAbs_Sphere:
1057 gp_Sphere Sphere = S->Sphere();
1058 ElSLib::Parameters( Sphere, P3d, u, v );
1063 gp_Torus Torus = S->Torus();
1064 ElSLib::Parameters( Torus, P3d, u, v );
1068 Standard_NoSuchObject::Raise("ProjLib_ComputeApprox::Value");
1070 Standard_Boolean ToMirror = Standard_False;
1071 Standard_Real du = 0., dv = 0.;
1072 Standard_Integer number;
1075 if (SType == GeomAbs_Sphere && Abs(u-F.myU1) > M_PI)
1077 ToMirror = Standard_True;
1081 Standard_Real newV = ElCLib::InPeriod( v, F.myV1, F.myV2 );
1082 number = (Standard_Integer) (Floor((newV-v)/(F.myV2-F.myV1)));
1083 dv -= number*(F.myV2-F.myV1);
1085 if (F.UCouture || (F.VCouture && SType == GeomAbs_Sphere))
1087 gp_Pnt2d P2d = F.Value( UFirst );
1088 number = (Standard_Integer) (Floor((P2d.X()-u)/M_PI + Epsilon(M_PI)));
1092 if (!myBSpline.IsNull())
1094 if (du != 0. || dv != 0.)
1095 myBSpline->Translate( gp_Vec2d(du,dv) );
1098 gp_Ax2d Axe( gp_Pnt2d(0.,0.), gp_Dir2d(1.,0.) );
1099 myBSpline->Mirror( Axe );
1105 //=======================================================================
1106 //function : BSpline
1108 //=======================================================================
1110 Handle(Geom2d_BSplineCurve) ProjLib_ComputeApprox::BSpline() const
1116 //=======================================================================
1119 //=======================================================================
1121 Handle(Geom2d_BezierCurve) ProjLib_ComputeApprox::Bezier() const
1128 //=======================================================================
1129 //function : Tolerance
1131 //=======================================================================
1133 Standard_Real ProjLib_ComputeApprox::Tolerance() const