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