0022048: Visualization, AIS_InteractiveContext - single object selection should alway...
[occt.git] / src / ProjLib / ProjLib_ComputeApproxOnPolarSurface.cxx
CommitLineData
b311480e 1// Created by: Bruno DUMORTIER
2// Copyright (c) 1995-1999 Matra Datavision
973c2be1 3// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 4//
973c2be1 5// This file is part of Open CASCADE Technology software library.
b311480e 6//
d5f74e42 7// This library is free software; you can redistribute it and/or modify it under
8// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 9// by the Free Software Foundation, with special exception defined in the file
10// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11// distribution for complete text of the license and disclaimer of any warranty.
b311480e 12//
973c2be1 13// Alternatively, this file may be used under the terms of Open CASCADE
14// commercial license or contractual agreement.
b311480e 15
7fd59977 16#include <ProjLib_ComputeApproxOnPolarSurface.hxx>
368cdde6 17#include <AppCont_Function.hxx>
7fd59977 18#include <ElSLib.hxx>
19#include <ElCLib.hxx>
20#include <BSplCLib.hxx>
21#include <PLib.hxx>
22#include <Standard_NoSuchObject.hxx>
23#include <Geom_UndefinedDerivative.hxx>
24#include <gp_Trsf.hxx>
25#include <Precision.hxx>
26#include <Approx_FitAndDivide2d.hxx>
27#include <math.hxx>
28#include <AppParCurves_MultiCurve.hxx>
29#include <Geom_Surface.hxx>
30#include <Geom2d_BSplineCurve.hxx>
31#include <Geom2d_BezierCurve.hxx>
32#include <Geom2d_Line.hxx>
33#include <Geom2d_Circle.hxx>
34#include <Geom2d_Ellipse.hxx>
35#include <Geom2d_Hyperbola.hxx>
36#include <Geom2d_Parabola.hxx>
37#include <Geom2d_TrimmedCurve.hxx>
38#include <Geom_BSplineSurface.hxx>
39#include <Geom_BezierSurface.hxx>
40#include <Geom_BSplineCurve.hxx>
41#include <Geom_BezierCurve.hxx>
42#include <Geom_TrimmedCurve.hxx>
43
44#include <TColgp_Array1OfPnt2d.hxx>
45#include <TColgp_Array2OfPnt2d.hxx>
46#include <TColgp_Array1OfPnt.hxx>
47#include <TColgp_SequenceOfPnt2d.hxx>
48#include <TColStd_Array1OfReal.hxx>
49#include <TColStd_Array1OfInteger.hxx>
50#include <TColStd_SequenceOfReal.hxx>
51#include <TColStd_ListOfTransient.hxx>
52
53#include <GeomAbs_SurfaceType.hxx>
54#include <GeomAbs_CurveType.hxx>
7fd59977 55#include <Adaptor3d_Surface.hxx>
56#include <Adaptor3d_Curve.hxx>
57#include <Adaptor3d_HSurface.hxx>
58#include <Adaptor3d_HCurve.hxx>
59#include <Adaptor2d_HCurve2d.hxx>
60#include <Geom2dAdaptor_Curve.hxx>
61#include <Geom2dAdaptor_HCurve.hxx>
62#include <GeomAdaptor_HCurve.hxx>
63#include <GeomAdaptor.hxx>
64#include <GeomAdaptor_Surface.hxx>
65#include <TColgp_SequenceOfPnt.hxx>
66
67#include <gp_Pnt.hxx>
68#include <gp_Pnt2d.hxx>
69#include <gp_Vec2d.hxx>
70#include <Extrema_GenLocateExtPS.hxx>
71#include <Extrema_ExtPS.hxx>
72#include <GCPnts_QuasiUniformAbscissa.hxx>
73#include <Standard_DomainError.hxx>
74//#include <GeomLib_IsIso.hxx>
75//#include <GeomLib_CheckSameParameter.hxx>
76
0797d9d3 77#ifdef OCCT_DEBUG
7fd59977 78#ifdef DRAW
79#include <DrawTrSurf.hxx>
ec357c5c 80#include <Geom2d_Curve.hxx>
7fd59977 81#endif
82//static Standard_Integer compteur = 0;
83#endif
84
ea3b6e72 85struct aFuncStruct
86{
87 Handle(Adaptor3d_HSurface) mySurf; // Surface where to project.
88 Handle(Adaptor3d_HCurve) myCurve; // Curve to project.
89 Handle(Adaptor2d_HCurve2d) myInitCurve2d; // Initial 2dcurve projection.
90 Standard_Real mySqProjOrtTol; // Used to filter non-orthogonal projected point.
91 Standard_Real myTolU;
92 Standard_Real myTolV;
93 Standard_Real myPeriod[2]; // U and V period correspondingly.
94};
95
7afe616f 96//=======================================================================
97//function : computePeriodicity
98//purpose : Compute period information on adaptor.
99//=======================================================================
100static void computePeriodicity(const Handle(Adaptor3d_HSurface)& theSurf,
101 Standard_Real &theUPeriod,
102 Standard_Real &theVPeriod)
103{
104 theUPeriod = 0.0;
105 theVPeriod = 0.0;
106
107 // Compute once information about periodicity.
108 // Param space may be reduced in case of rectangular trimmed surface,
109 // in this case really trimmed bounds should be set as unperiodic.
110 Standard_Real aTrimF, aTrimL, aBaseF, aBaseL, aDummyF, aDummyL;
111 Handle(Geom_Surface) aS = GeomAdaptor::MakeSurface(theSurf->Surface(), Standard_False); // Not trim.
112 // U param space.
113 if (theSurf->IsUPeriodic())
114 {
115 theUPeriod = theSurf->UPeriod();
116 }
117 else if(theSurf->IsUClosed())
118 {
119 theUPeriod = theSurf->LastUParameter() - theSurf->FirstUParameter();
120 }
121 if (theUPeriod != 0.0)
122 {
123 aTrimF = theSurf->FirstUParameter(); // Trimmed first
124 aTrimL = theSurf->LastUParameter(); // Trimmed last
125 aS->Bounds(aBaseF, aBaseL, aDummyF, aDummyL); // Non-trimmed values.
126 if (Abs (aBaseF - aTrimF) + Abs (aBaseL - aTrimL) > Precision::PConfusion())
127 {
128 // Param space reduced.
129 theUPeriod = 0.0;
130 }
131 }
132
133 // V param space.
134 if (theSurf->IsVPeriodic())
135 {
136 theVPeriod = theSurf->VPeriod();
137 }
138 else if(theSurf->IsVClosed())
139 {
140 theVPeriod = theSurf->LastVParameter() - theSurf->FirstVParameter();
141 }
142 if (theVPeriod != 0.0)
143 {
144 aTrimF = theSurf->FirstVParameter(); // Trimmed first
145 aTrimL = theSurf->LastVParameter(); // Trimmed last
146 aS->Bounds(aDummyF, aDummyL, aBaseF, aBaseL); // Non-trimmed values.
147 if (Abs (aBaseF - aTrimF) + Abs (aBaseL - aTrimL) > Precision::PConfusion())
148 {
149 // Param space reduced.
150 theVPeriod = 0.0;
151 }
152 }
153}
154
ea3b6e72 155//=======================================================================
156//function : aFuncValue
157//purpose : compute functional value in (theU,theV) point
158//=======================================================================
159static Standard_Real anOrthogSqValue(const gp_Pnt& aBasePnt,
160 const Handle(Adaptor3d_HSurface)& Surf,
161 const Standard_Real theU,
162 const Standard_Real theV)
163{
164 // Since find projection, formula is:
165 // F1 = Dot(S_U, Vec(aBasePnt, aProjPnt))
166 // F2 = Dot(S_V, Vec(aBasePnt, aProjPnt))
167
168 gp_Pnt aProjPnt;
169 gp_Vec aSu, aSv;
170
171 Surf->D1(theU, theV, aProjPnt, aSu, aSv);
172 gp_Vec aBaseVec(aBasePnt, aProjPnt);
173
174 if (aSu.SquareMagnitude() > Precision::SquareConfusion())
175 aSu.Normalize();
176
177 if (aSv.SquareMagnitude() > Precision::SquareConfusion())
178 aSv.Normalize();
179
180 Standard_Real aFirstPart = aSu.Dot(aBaseVec);
181 Standard_Real aSecondPart = aSv.Dot(aBaseVec);
182 return (aFirstPart * aFirstPart + aSecondPart * aSecondPart);
183}
184
7fd59977 185//=======================================================================
186//function : Value
187//purpose : (OCC217 - apo)- Compute Point2d that project on polar surface(<Surf>) 3D<Curve>
188// <InitCurve2d> use for calculate start 2D point.
189//=======================================================================
ea3b6e72 190static gp_Pnt2d Function_Value(const Standard_Real theU,
191 const aFuncStruct& theData)
7fd59977 192{
ea3b6e72 193 gp_Pnt2d p2d = theData.myInitCurve2d->Value(theU) ;
194 gp_Pnt p = theData.myCurve->Value(theU);
195 gp_Pnt aSurfPnt = theData.mySurf->Value(p2d.X(), p2d.Y());
196 Standard_Real aSurfPntDist = aSurfPnt.SquareDistance(p);
7fd59977 197
7fd59977 198 Standard_Real Uinf, Usup, Vinf, Vsup;
ea3b6e72 199 Uinf = theData.mySurf->Surface().FirstUParameter();
200 Usup = theData.mySurf->Surface().LastUParameter();
201 Vinf = theData.mySurf->Surface().FirstVParameter();
202 Vsup = theData.mySurf->Surface().LastVParameter();
203
204 // Check case when curve is close to co-parametrized isoline on surf.
205 if (Abs (p2d.X() - Uinf) < Precision::PConfusion() ||
206 Abs (p2d.X() - Usup) < Precision::PConfusion() )
207 {
208 // V isoline.
209 gp_Pnt aPnt;
210 theData.mySurf->D0(p2d.X(), theU, aPnt);
211 if (aPnt.SquareDistance(p) < aSurfPntDist)
212 p2d.SetY(theU);
213 }
214
215 if (Abs (p2d.Y() - Vinf) < Precision::PConfusion() ||
216 Abs (p2d.Y() - Vsup) < Precision::PConfusion() )
217 {
218 // U isoline.
219 gp_Pnt aPnt;
220 theData.mySurf->D0(theU, p2d.Y(), aPnt);
221 if (aPnt.SquareDistance(p) < aSurfPntDist)
222 p2d.SetX(theU);
223 }
224
7fd59977 225 Standard_Integer decalU = 0, decalV = 0;
226 Standard_Real U0 = p2d.X(), V0 = p2d.Y();
ea3b6e72 227
228 GeomAbs_SurfaceType Type = theData.mySurf->GetType();
7fd59977 229 if((Type != GeomAbs_BSplineSurface) &&
230 (Type != GeomAbs_BezierSurface) &&
ea3b6e72 231 (Type != GeomAbs_OffsetSurface) )
232 {
233 // Analytical cases.
1d47d8d0 234 Standard_Real S = 0., T = 0.;
ea3b6e72 235 switch (Type)
236 {
7fd59977 237 case GeomAbs_Cylinder:
238 {
ea3b6e72 239 gp_Cylinder Cylinder = theData.mySurf->Cylinder();
240 ElSLib::Parameters( Cylinder, p, S, T);
241 if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
242 if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
243 S += decalU*2*M_PI;
244 break;
7fd59977 245 }
246 case GeomAbs_Cone:
247 {
ea3b6e72 248 gp_Cone Cone = theData.mySurf->Cone();
249 ElSLib::Parameters( Cone, p, S, T);
250 if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
251 if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
252 S += decalU*2*M_PI;
253 break;
7fd59977 254 }
255 case GeomAbs_Sphere:
256 {
ea3b6e72 257 gp_Sphere Sphere = theData.mySurf->Sphere();
258 ElSLib::Parameters( Sphere, p, S, T);
259 if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
260 if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
261 S += decalU*2*M_PI;
262 if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1;
263 if(V0 > (Vsup+(Vsup-Vinf))) decalV = int((V0 - Vsup+(Vsup-Vinf))/(2*M_PI))+1;
264 T += decalV*2*M_PI;
265 if(0.4*M_PI < Abs(U0 - S) && Abs(U0 - S) < 1.6*M_PI)
266 {
267 T = M_PI - T;
268 if(U0 < S)
269 S -= M_PI;
270 else
271 S += M_PI;
272 }
273 break;
7fd59977 274 }
275 case GeomAbs_Torus:
276 {
ea3b6e72 277 gp_Torus Torus = theData.mySurf->Torus();
278 ElSLib::Parameters( Torus, p, S, T);
279 if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
280 if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
281 if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1;
282 if(V0 > Vsup) decalV = int((V0 - Vsup)/(2*M_PI))+1;
283 S += decalU*2*M_PI; T += decalV*2*M_PI;
284 break;
7fd59977 285 }
286 default:
9775fa61 287 throw Standard_NoSuchObject("ProjLib_ComputeApproxOnPolarSurface::Value");
7fd59977 288 }
289 return gp_Pnt2d(S, T);
290 }
ea3b6e72 291
292 // Non-analytical case.
7fd59977 293 Standard_Real Dist2Min = RealLast();
ea3b6e72 294 Standard_Real uperiod = theData.myPeriod[0],
295 vperiod = theData.myPeriod[1],
296 u, v;
7afe616f 297
298 // U0 and V0 are the points within the initialized period.
ea3b6e72 299 if(U0 < Uinf)
300 {
7fd59977 301 if(!uperiod)
302 U0 = Uinf;
ea3b6e72 303 else
304 {
7fd59977 305 decalU = int((Uinf - U0)/uperiod)+1;
306 U0 += decalU*uperiod;
307 }
eafb234b 308 }
ea3b6e72 309 if(U0 > Usup)
310 {
7fd59977 311 if(!uperiod)
312 U0 = Usup;
ea3b6e72 313 else
314 {
7fd59977 315 decalU = -(int((U0 - Usup)/uperiod)+1);
316 U0 += decalU*uperiod;
317 }
eafb234b 318 }
ea3b6e72 319 if(V0 < Vinf)
320 {
7fd59977 321 if(!vperiod)
322 V0 = Vinf;
ea3b6e72 323 else
324 {
7fd59977 325 decalV = int((Vinf - V0)/vperiod)+1;
326 V0 += decalV*vperiod;
327 }
eafb234b 328 }
ea3b6e72 329 if(V0 > Vsup)
330 {
7fd59977 331 if(!vperiod)
332 V0 = Vsup;
ea3b6e72 333 else
334 {
7fd59977 335 decalV = -int((V0 - Vsup)/vperiod)-1;
336 V0 += decalV*vperiod;
337 }
eafb234b 338 }
7afe616f 339
340 // The surface around (U0,V0) is reduced.
7fd59977 341 Standard_Real uLittle = (Usup - Uinf)/10, vLittle = (Vsup - Vinf)/10;
342 Standard_Real uInfLi = 0, vInfLi = 0,uSupLi = 0, vSupLi = 0;
343 if((U0 - Uinf) > uLittle) uInfLi = U0 - uLittle; else uInfLi = Uinf;
344 if((V0 - Vinf) > vLittle) vInfLi = V0 - vLittle; else vInfLi = Vinf;
345 if((Usup - U0) > uLittle) uSupLi = U0 + uLittle; else uSupLi = Usup;
346 if((Vsup - V0) > vLittle) vSupLi = V0 + vLittle; else vSupLi = Vsup;
7fd59977 347
348 GeomAdaptor_Surface SurfLittle;
ea3b6e72 349 if (Type == GeomAbs_BSplineSurface)
350 {
351 Handle(Geom_Surface) GBSS(theData.mySurf->Surface().BSpline());
7fd59977 352 SurfLittle.Load(GBSS, uInfLi, uSupLi, vInfLi, vSupLi);
353 }
ea3b6e72 354 else if (Type == GeomAbs_BezierSurface)
355 {
356 Handle(Geom_Surface) GS(theData.mySurf->Surface().Bezier());
7fd59977 357 SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi);
358 }
ea3b6e72 359 else if (Type == GeomAbs_OffsetSurface)
360 {
361 Handle(Geom_Surface) GS = GeomAdaptor::MakeSurface(theData.mySurf->Surface());
7fd59977 362 SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi);
363 }
ea3b6e72 364 else
365 {
2d2b3d53 366 throw Standard_NoSuchObject ("ProjLib_ComputeApproxOnPolarSurface::ProjectUsingInitialCurve2d() - unknown surface type");
7fd59977 367 }
368
ea3b6e72 369 // Try to run simple search with initial point (U0, V0).
e5260e1d 370 Extrema_GenLocateExtPS locext(SurfLittle, theData.myTolU, theData.myTolV);
371 locext.Perform(p, U0, V0);
ea3b6e72 372 if (locext.IsDone())
373 {
374 locext.Point().Parameter(u, v);
375 Dist2Min = anOrthogSqValue(p, theData.mySurf, u, v);
376 if (Dist2Min < theData.mySqProjOrtTol && // Point is projection.
377 locext.SquareDistance() < aSurfPntDist + Precision::SquareConfusion()) // Point better than initial.
378 {
7fd59977 379 gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod);
380 return pnt;
381 }
ea3b6e72 382 }
383
384 // Perform whole param space search.
385 Extrema_ExtPS ext(p, SurfLittle, theData.myTolU, theData.myTolV);
386 if (ext.IsDone() && ext.NbExt() >= 1)
387 {
7fd59977 388 Dist2Min = ext.SquareDistance(1);
389 Standard_Integer GoodValue = 1;
ea3b6e72 390 for (Standard_Integer i = 2 ; i <= ext.NbExt() ; i++ )
391 {
392 if( Dist2Min > ext.SquareDistance(i))
393 {
394 Dist2Min = ext.SquareDistance(i);
395 GoodValue = i;
7fd59977 396 }
ea3b6e72 397 }
398 ext.Point(GoodValue).Parameter(u, v);
399 Dist2Min = anOrthogSqValue(p, theData.mySurf, u, v);
400 if (Dist2Min < theData.mySqProjOrtTol && // Point is projection.
401 ext.SquareDistance(GoodValue) < aSurfPntDist + Precision::SquareConfusion()) // Point better than initial.
402 {
7fd59977 403 gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod);
404 return pnt;
405 }
406 }
ea3b6e72 407
408 // Both searches return bad values, use point from initial 2dcurve.
7fd59977 409 return p2d;
410}
411
412
413//=======================================================================
414//function : ProjLib_PolarFunction
415//purpose : (OCC217 - apo)- This class produce interface to call "gp_Pnt2d Function_Value(...)"
416//=======================================================================
417
368cdde6 418class ProjLib_PolarFunction : public AppCont_Function
7fd59977 419{
ea3b6e72 420 aFuncStruct myStruct;
368cdde6 421
7fd59977 422 public :
368cdde6 423
7fd59977 424 ProjLib_PolarFunction(const Handle(Adaptor3d_HCurve) & C,
368cdde6 425 const Handle(Adaptor3d_HSurface)& Surf,
426 const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
427 const Standard_Real Tol3d)
368cdde6 428 {
429 myNbPnt = 0;
430 myNbPnt2d = 1;
ea3b6e72 431
7afe616f 432 computePeriodicity(Surf, myStruct.myPeriod[0], myStruct.myPeriod[1]);
ea3b6e72 433
434 myStruct.myCurve = C;
435 myStruct.myInitCurve2d = InitialCurve2d;
436 myStruct.mySurf = Surf;
437 myStruct.mySqProjOrtTol = 10000.0 * Tol3d * Tol3d;
438 myStruct.myTolU = Surf->UResolution(Tol3d);
439 myStruct.myTolV = Surf->VResolution(Tol3d);
368cdde6 440 }
441
7fd59977 442 ~ProjLib_PolarFunction() {}
368cdde6 443
7fd59977 444 Standard_Real FirstParameter() const
368cdde6 445 {
ea3b6e72 446 return myStruct.myCurve->FirstParameter();
368cdde6 447 }
448
7fd59977 449 Standard_Real LastParameter() const
368cdde6 450 {
ea3b6e72 451 return myStruct.myCurve->LastParameter();
7fd59977 452 }
368cdde6 453
ea3b6e72 454 gp_Pnt2d Value(const Standard_Real t) const
fa0f5a55 455 {
ea3b6e72 456 return Function_Value(t, myStruct);
fa0f5a55 457 }
458
368cdde6 459 Standard_Boolean Value(const Standard_Real theT,
460 NCollection_Array1<gp_Pnt2d>& thePnt2d,
461 NCollection_Array1<gp_Pnt>& /*thePnt*/) const
462 {
ea3b6e72 463 thePnt2d(1) = Function_Value(theT, myStruct);
368cdde6 464 return Standard_True;
465 }
466
467 Standard_Boolean D1(const Standard_Real /*theT*/,
468 NCollection_Array1<gp_Vec2d>& /*theVec2d*/,
469 NCollection_Array1<gp_Vec>& /*theVec*/) const
7fd59977 470 {return Standard_False;}
471};
472
473//=======================================================================
474//function : ProjLib_ComputeApproxOnPolarSurface
475//purpose :
476//=======================================================================
477
cbff1e55 478ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface()
479: myProjIsDone(Standard_False),
480 myTolerance (-1.0)
481{
482}
7fd59977 483
484
485//=======================================================================
486//function : ProjLib_ComputeApproxOnPolarSurface
487//purpose :
488//=======================================================================
489
490ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface
cbff1e55 491 (const Handle(Adaptor2d_HCurve2d)& theInitialCurve2d,
492 const Handle(Adaptor3d_HCurve)& theCurve,
493 const Handle(Adaptor3d_HSurface)& theSurface,
494 const Standard_Real theTolerance3D)
495: myProjIsDone(Standard_False),
496 myTolerance (theTolerance3D)
7fd59977 497{
cbff1e55 498 myBSpline = Perform(theInitialCurve2d, theCurve, theSurface);
7fd59977 499}
500//=======================================================================
501//function : ProjLib_ComputeApproxOnPolarSurface
502//purpose : Process the case of sewing
503//=======================================================================
504
505ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface
cbff1e55 506 (const Handle(Adaptor2d_HCurve2d)& theInitialCurve2d,
507 const Handle(Adaptor2d_HCurve2d)& theInitialCurve2dBis,
508 const Handle(Adaptor3d_HCurve)& theCurve,
509 const Handle(Adaptor3d_HSurface)& theSurface,
510 const Standard_Real theTolerance3D)
511: myProjIsDone(Standard_False),
512 myTolerance (theTolerance3D)
513{
514 // InitialCurve2d and InitialCurve2dBis are two pcurves of the sewing
515 Handle(Geom2d_BSplineCurve) bsc =
516 Perform(theInitialCurve2d, theCurve, theSurface);
517
7fd59977 518 if(myProjIsDone) {
519 gp_Pnt2d P2dproj, P2d, P2dBis;
520 P2dproj = bsc->StartPoint();
cbff1e55 521 P2d = theInitialCurve2d->Value(theInitialCurve2d->FirstParameter());
522 P2dBis = theInitialCurve2dBis->Value(theInitialCurve2dBis->FirstParameter());
7fd59977 523
524 Standard_Real Dist, DistBis;
525 Dist = P2dproj.Distance(P2d);
526 DistBis = P2dproj.Distance(P2dBis);
527 if( Dist < DistBis) {
528 // myBSpline2d is the pcurve that is found. It is translated to obtain myCurve2d
529 myBSpline = bsc;
530 Handle(Geom2d_Geometry) GG = myBSpline->Translated(P2d, P2dBis);
531 my2ndCurve = Handle(Geom2d_Curve)::DownCast(GG);
532 }
533 else {
534 my2ndCurve = bsc;
535 Handle(Geom2d_Geometry) GG = my2ndCurve->Translated(P2dBis, P2d);
536 myBSpline = Handle(Geom2d_BSplineCurve)::DownCast(GG);
537 }
538 }
539}
540
541//=======================================================================
542//function : ProjLib_ComputeApproxOnPolarSurface
543//purpose : case without curve of initialization
544//=======================================================================
545
546ProjLib_ComputeApproxOnPolarSurface::ProjLib_ComputeApproxOnPolarSurface
cbff1e55 547 (const Handle(Adaptor3d_HCurve)& theCurve,
548 const Handle(Adaptor3d_HSurface)& theSurface,
549 const Standard_Real theTolerance3D)
550: myProjIsDone(Standard_False),
551 myTolerance (theTolerance3D)
7fd59977 552{
cbff1e55 553 const Handle(Adaptor2d_HCurve2d) anInitCurve2d;
554 myBSpline = Perform(anInitCurve2d, theCurve, theSurface);
7fd59977 555}
556
cbff1e55 557//=======================================================================
558//function : Concat
559//purpose :
560//=======================================================================
561
7fd59977 562static Handle(Geom2d_BSplineCurve) Concat(Handle(Geom2d_BSplineCurve) C1,
460f4f69 563 Handle(Geom2d_BSplineCurve) C2,
564 Standard_Real theUJump,
565 Standard_Real theVJump)
7fd59977 566{
567 Standard_Integer deg, deg1, deg2;
568 deg1 = C1->Degree();
569 deg2 = C2->Degree();
570
571 if ( deg1 < deg2) {
572 C1->IncreaseDegree(deg2);
573 deg = deg2;
574 }
575 else if ( deg2 < deg1) {
576 C2->IncreaseDegree(deg1);
577 deg = deg1;
578 }
579 else deg = deg1;
580
581 Standard_Integer np1,np2,nk1,nk2,np,nk;
582 np1 = C1->NbPoles();
583 nk1 = C1->NbKnots();
584 np2 = C2->NbPoles();
585 nk2 = C2->NbKnots();
586 nk = nk1 + nk2 -1;
587 np = np1 + np2 -1;
588
589 TColStd_Array1OfReal K1(1,nk1); C1->Knots(K1);
590 TColStd_Array1OfInteger M1(1,nk1); C1->Multiplicities(M1);
591 TColgp_Array1OfPnt2d P1(1,np1); C1->Poles(P1);
592 TColStd_Array1OfReal K2(1,nk2); C2->Knots(K2);
593 TColStd_Array1OfInteger M2(1,nk2); C2->Multiplicities(M2);
594 TColgp_Array1OfPnt2d P2(1,np2); C2->Poles(P2);
595
596 // Compute the new BSplineCurve
597 TColStd_Array1OfReal K(1,nk);
598 TColStd_Array1OfInteger M(1,nk);
599 TColgp_Array1OfPnt2d P(1,np);
600
601 Standard_Integer i, count = 0;
602 // Set Knots and Mults
603 for ( i = 1; i <= nk1; i++) {
604 count++;
605 K(count) = K1(i);
606 M(count) = M1(i);
607 }
608 M(count) = deg;
609 for ( i = 2; i <= nk2; i++) {
610 count++;
611 K(count) = K2(i);
612 M(count) = M2(i);
613 }
614 // Set the Poles
615 count = 0;
616 for (i = 1; i <= np1; i++) {
617 count++;
618 P(count) = P1(i);
619 }
620 for (i = 2; i <= np2; i++) {
621 count++;
460f4f69 622 P(count).SetX(P2(i).X() + theUJump);
623 P(count).SetY(P2(i).Y() + theVJump);
7fd59977 624 }
625
626 Handle(Geom2d_BSplineCurve) BS =
627 new Geom2d_BSplineCurve(P,K,M,deg);
628 return BS;
629}
630
631
632//=======================================================================
633//function : Perform
634//purpose :
635//=======================================================================
636Handle(Geom2d_BSplineCurve) ProjLib_ComputeApproxOnPolarSurface::Perform
637(const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
638 const Handle(Adaptor3d_HCurve)& Curve,
639 const Handle(Adaptor3d_HSurface)& S)
640{
641 //OCC217
642 Standard_Real Tol3d = myTolerance;
643 Standard_Real ParamTol = Precision::PApproximation();
644
645 Handle(Adaptor2d_HCurve2d) AHC2d = InitialCurve2d;
646 Handle(Adaptor3d_HCurve) AHC = Curve;
647
648// if the curve 3d is a BSpline with degree C0, it is cut into sections with degree C1
649// -> bug cts18237
650 GeomAbs_CurveType typeCurve = Curve->GetType();
651 if(typeCurve == GeomAbs_BSplineCurve) {
652 TColStd_ListOfTransient LOfBSpline2d;
653 Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
654 Standard_Integer nbInter = Curve->NbIntervals(GeomAbs_C1);
655 if(nbInter > 1) {
656 Standard_Integer i, j;
657 Handle(Geom_TrimmedCurve) GTC;
658 Handle(Geom2d_TrimmedCurve) G2dTC;
659 TColStd_Array1OfReal Inter(1,nbInter+1);
660 Curve->Intervals(Inter,GeomAbs_C1);
661 Standard_Real firstinter = Inter.Value(1), secondinter = Inter.Value(2);
662 // initialization 3d
663 GTC = new Geom_TrimmedCurve(BSC, firstinter, secondinter);
664 AHC = new GeomAdaptor_HCurve(GTC);
665
666 // if there is an initialization curve:
667 // - either this is a BSpline C0, with discontinuity at the same parameters of nodes
668 // and the sections C1 are taken
669 // - or this is a curve C1 and the sections of intrest are taken otherwise the curve is created.
670
671 // initialization 2d
672 Standard_Integer nbInter2d;
673 Standard_Boolean C2dIsToCompute;
674 C2dIsToCompute = InitialCurve2d.IsNull();
675 Handle(Geom2d_BSplineCurve) BSC2d;
676 Handle(Geom2d_Curve) G2dC;
677
678 if(!C2dIsToCompute) {
679 nbInter2d = InitialCurve2d->NbIntervals(GeomAbs_C1);
680 TColStd_Array1OfReal Inter2d(1,nbInter2d+1);
681 InitialCurve2d->Intervals(Inter2d,GeomAbs_C1);
682 j = 1;
683 for(i = 1,j = 1;i <= nbInter;i++)
684 if(Abs(Inter.Value(i) - Inter2d.Value(j)) < ParamTol) { //OCC217
685 //if(Abs(Inter.Value(i) - Inter2d.Value(j)) < myTolerance) {
686 if (j > nbInter2d) break;
687 j++;
688 }
689 if(j != (nbInter2d+1)) {
690 C2dIsToCompute = Standard_True;
691 }
692 }
693
694 if(C2dIsToCompute) {
695 AHC2d = BuildInitialCurve2d(AHC, S);
696 }
697 else {
698 typeCurve = InitialCurve2d->GetType();
699 switch (typeCurve) {
700 case GeomAbs_Line: {
701 G2dC = new Geom2d_Line(InitialCurve2d->Line());
702 break;
703 }
704 case GeomAbs_Circle: {
705 G2dC = new Geom2d_Circle(InitialCurve2d->Circle());
706 break;
707 }
708 case GeomAbs_Ellipse: {
709 G2dC = new Geom2d_Ellipse(InitialCurve2d->Ellipse());
710 break;
711 }
712 case GeomAbs_Hyperbola: {
713 G2dC = new Geom2d_Hyperbola(InitialCurve2d->Hyperbola());
714 break;
715 }
716 case GeomAbs_Parabola: {
717 G2dC = new Geom2d_Parabola(InitialCurve2d->Parabola());
718 break;
719 }
720 case GeomAbs_BezierCurve: {
721 G2dC = InitialCurve2d->Bezier();
722 break;
723 }
724 case GeomAbs_BSplineCurve: {
725 G2dC = InitialCurve2d->BSpline();
726 break;
727 }
728 case GeomAbs_OtherCurve:
729 default:
730 break;
731 }
732 gp_Pnt2d fp2d = G2dC->Value(firstinter), lp2d = G2dC->Value(secondinter);
733 gp_Pnt fps, lps, fpc, lpc;
734 S->D0(fp2d.X(), fp2d.Y(), fps);
735 S->D0(lp2d.X(), lp2d.Y(), lps);
736 Curve->D0(firstinter, fpc);
737 Curve->D0(secondinter, lpc);
738 //OCC217
739 if((fps.IsEqual(fpc, Tol3d)) &&
740 (lps.IsEqual(lpc, Tol3d))) {
741 //if((fps.IsEqual(fpc, myTolerance)) &&
742 // (lps.IsEqual(lpc, myTolerance))) {
743 G2dTC = new Geom2d_TrimmedCurve(G2dC, firstinter, secondinter);
744 Geom2dAdaptor_Curve G2dAC(G2dTC);
745 AHC2d = new Geom2dAdaptor_HCurve(G2dAC);
746 myProjIsDone = Standard_True;
747 }
748 else {
749 AHC2d = BuildInitialCurve2d(AHC, S);
750 C2dIsToCompute = Standard_True;
751 }
752 }
753
754 if(myProjIsDone) {
755 BSC2d = ProjectUsingInitialCurve2d(AHC, S, AHC2d);
756 if(BSC2d.IsNull()) return Handle(Geom2d_BSplineCurve)(); //IFV
757 LOfBSpline2d.Append(BSC2d);
758 }
759 else {
760 return Handle(Geom2d_BSplineCurve)();
761 }
762
763
764
765 Standard_Real iinter, ip1inter;
766 Standard_Integer nbK2d, deg;
767 nbK2d = BSC2d->NbKnots(); deg = BSC2d->Degree();
768
769 for(i = 2;i <= nbInter;i++) {
770 iinter = Inter.Value(i);
771 ip1inter = Inter.Value(i+1);
772 // general case 3d
773 GTC->SetTrim(iinter, ip1inter);
774 AHC = new GeomAdaptor_HCurve(GTC);
775
776 // general case 2d
777 if(C2dIsToCompute) {
778 AHC2d = BuildInitialCurve2d(AHC, S);
779 }
780 else {
781 gp_Pnt2d fp2d = G2dC->Value(iinter), lp2d = G2dC->Value(ip1inter);
782 gp_Pnt fps, lps, fpc, lpc;
783 S->D0(fp2d.X(), fp2d.Y(), fps);
784 S->D0(lp2d.X(), lp2d.Y(), lps);
785 Curve->D0(iinter, fpc);
786 Curve->D0(ip1inter, lpc);
787 //OCC217
788 if((fps.IsEqual(fpc, Tol3d)) &&
789 (lps.IsEqual(lpc, Tol3d))) {
790 //if((fps.IsEqual(fpc, myTolerance)) &&
791 // (lps.IsEqual(lpc, myTolerance))) {
792 G2dTC->SetTrim(iinter, ip1inter);
793 Geom2dAdaptor_Curve G2dAC(G2dTC);
794 AHC2d = new Geom2dAdaptor_HCurve(G2dAC);
795 myProjIsDone = Standard_True;
796 }
797 else {
798 AHC2d = BuildInitialCurve2d(AHC, S);
799 }
800 }
801 if(myProjIsDone) {
802 BSC2d = ProjectUsingInitialCurve2d(AHC, S, AHC2d);
803 if(BSC2d.IsNull()) {
804 return Handle(Geom2d_BSplineCurve)();
805 }
806 LOfBSpline2d.Append(BSC2d);
807 nbK2d += BSC2d->NbKnots() - 1;
808 deg = Max(deg, BSC2d->Degree());
809 }
810 else {
811 return Handle(Geom2d_BSplineCurve)();
812 }
813 }
814
7afe616f 815 Standard_Real anUPeriod, anVPeriod;
816 computePeriodicity(S, anUPeriod, anVPeriod);
7fd59977 817 Standard_Integer NbC = LOfBSpline2d.Extent();
818 Handle(Geom2d_BSplineCurve) CurBS;
819 CurBS = Handle(Geom2d_BSplineCurve)::DownCast(LOfBSpline2d.First());
820 LOfBSpline2d.RemoveFirst();
460f4f69 821 for (Standard_Integer ii = 2; ii <= NbC; ii++)
822 {
823 Handle(Geom2d_BSplineCurve) BS =
824 Handle(Geom2d_BSplineCurve)::DownCast(LOfBSpline2d.First());
825
826 //Check for period jump in point of contact.
827 gp_Pnt2d aC1End = CurBS->Pole(CurBS->NbPoles()); // End of C1.
828 gp_Pnt2d aC2Beg = BS->Pole(1); // Beginning of C2.
829 Standard_Real anUJump = 0.0, anVJump = 0.0;
830
7afe616f 831 if (anUPeriod > 0.0 &&
832 Abs (aC1End.X() - aC2Beg.X()) > (anUPeriod ) / 2.01)
460f4f69 833 {
7afe616f 834 Standard_Real aMultCoeff = aC2Beg.X() < aC1End.X() ? 1.0 : -1.0;
835 anUJump = (anUPeriod) * aMultCoeff;
460f4f69 836 }
837
7afe616f 838 if (anVPeriod &&
839 Abs (aC1End.Y() - aC2Beg.Y()) > (anVPeriod) / 2.01)
460f4f69 840 {
7afe616f 841 Standard_Real aMultCoeff = aC2Beg.Y() < aC1End.Y() ? 1.0 : -1.0;
842 anVJump = (anVPeriod) * aMultCoeff;
460f4f69 843 }
844
845 CurBS = Concat(CurBS,BS, anUJump, anVJump);
846 LOfBSpline2d.RemoveFirst();
7fd59977 847 }
848 return CurBS;
849 }
850 }
851
852 if(InitialCurve2d.IsNull()) {
853 AHC2d = BuildInitialCurve2d(Curve, S);
854 if(!myProjIsDone)
855 return Handle(Geom2d_BSplineCurve)();
856 }
857 return ProjectUsingInitialCurve2d(AHC, S, AHC2d);
858}
859
860//=======================================================================
861//function : ProjLib_BuildInitialCurve2d
862//purpose :
863//=======================================================================
864
865Handle(Adaptor2d_HCurve2d)
866 ProjLib_ComputeApproxOnPolarSurface::
867 BuildInitialCurve2d(const Handle(Adaptor3d_HCurve)& Curve,
868 const Handle(Adaptor3d_HSurface)& Surf)
869{
870 // discretize the Curve with quasiuniform deflection
871 // density at least NbOfPnts points
872 myProjIsDone = Standard_False;
873
874 //OCC217
875 Standard_Real Tol3d = myTolerance;
876 Standard_Real TolU = Surf->UResolution(Tol3d), TolV = Surf->VResolution(Tol3d);
877 Standard_Real DistTol3d = 100.0*Tol3d;
878
7afe616f 879 Standard_Real uperiod = 0.0, vperiod = 0.0;
880 computePeriodicity(Surf, uperiod, vperiod);
ee9451ab 881
7fd59977 882 // NO myTol is Tol2d !!!!
883 //Standard_Real TolU = myTolerance, TolV = myTolerance;
884 //Standard_Real Tol3d = 100*myTolerance; // At random Balthazar.
885
886 Standard_Integer NbOfPnts = 61;
887 GCPnts_QuasiUniformAbscissa QUA(Curve->GetCurve(),NbOfPnts);
888 TColgp_Array1OfPnt Pts(1,NbOfPnts);
889 TColStd_Array1OfReal Param(1,NbOfPnts);
890 Standard_Integer i, j;
891 for( i = 1; i <= NbOfPnts ; i++ ) {
892 Param(i) = QUA.Parameter(i);
893 Pts(i) = Curve->Value(Param(i));
894 }
895
896 TColgp_Array1OfPnt2d Pts2d(1,NbOfPnts);
897 TColStd_Array1OfInteger Mult(1,NbOfPnts);
898 Mult.Init(1);
899 Mult(1) = Mult(NbOfPnts) = 2;
900
ee9451ab 901 Standard_Real Uinf, Usup, Vinf, Vsup;
902 Uinf = Surf->Surface().FirstUParameter();
903 Usup = Surf->Surface().LastUParameter();
7fd59977 904 Vinf = Surf->Surface().FirstVParameter();
905 Vsup = Surf->Surface().LastVParameter();
906 GeomAbs_SurfaceType Type = Surf->GetType();
907 if((Type != GeomAbs_BSplineSurface) && (Type != GeomAbs_BezierSurface) &&
908 (Type != GeomAbs_OffsetSurface)) {
909 Standard_Real S, T;
910// Standard_Integer usens = 0, vsens = 0;
911 // to know the position relatively to the period
912 switch (Type) {
913// case GeomAbs_Plane:
914// {
915// gp_Pln Plane = Surf->Plane();
916// for ( i = 1 ; i <= NbOfPnts ; i++) {
917// ElSLib::Parameters( Plane, Pts(i), S, T);
918// Pts2d(i).SetCoord(S,T);
919// }
920// myProjIsDone = Standard_True;
921// break;
922// }
923 case GeomAbs_Cylinder:
924 {
925// Standard_Real Sloc, Tloc;
eafb234b 926 Standard_Real Sloc;
927 Standard_Integer usens = 0;
928 gp_Cylinder Cylinder = Surf->Cylinder();
929 ElSLib::Parameters( Cylinder, Pts(1), S, T);
930 Pts2d(1).SetCoord(S,T);
931 for ( i = 2 ; i <= NbOfPnts ; i++) {
932 Sloc = S;
933 ElSLib::Parameters( Cylinder, Pts(i), S, T);
934 if(Abs(Sloc - S) > M_PI) {
935 if(Sloc > S)
936 usens++;
937 else
938 usens--;
939 }
940 Pts2d(i).SetCoord(S+usens*2*M_PI,T);
941 }
942 myProjIsDone = Standard_True;
943 break;
7fd59977 944 }
945 case GeomAbs_Cone:
946 {
947// Standard_Real Sloc, Tloc;
eafb234b 948 Standard_Real Sloc;
949 Standard_Integer usens = 0;
950 gp_Cone Cone = Surf->Cone();
951 ElSLib::Parameters( Cone, Pts(1), S, T);
952 Pts2d(1).SetCoord(S,T);
953 for ( i = 2 ; i <= NbOfPnts ; i++) {
954 Sloc = S;
955 ElSLib::Parameters( Cone, Pts(i), S, T);
956 if(Abs(Sloc - S) > M_PI) {
957 if(Sloc > S)
958 usens++;
959 else
960 usens--;
961 }
962 Pts2d(i).SetCoord(S+usens*2*M_PI,T);
963 }
964 myProjIsDone = Standard_True;
965 break;
7fd59977 966 }
967 case GeomAbs_Sphere:
968 {
969 Standard_Real Sloc, Tloc;
970 Standard_Integer usens = 0, vsens = 0; //usens steps by half-period
971 Standard_Boolean vparit = Standard_False;
972 gp_Sphere Sphere = Surf->Sphere();
973 ElSLib::Parameters( Sphere, Pts(1), S, T);
974 Pts2d(1).SetCoord(S,T);
975 for ( i = 2 ; i <= NbOfPnts ; i++) {
976 Sloc = S;Tloc = T;
977 ElSLib::Parameters( Sphere, Pts(i), S, T);
eafb234b 978 if(1.6*M_PI < Abs(Sloc - S)) {
7fd59977 979 if(Sloc > S)
980 usens += 2;
981 else
982 usens -= 2;
eafb234b 983 }
c6541a0c 984 if(1.6*M_PI > Abs(Sloc - S) && Abs(Sloc - S) > 0.4*M_PI) {
7fd59977 985 vparit = !vparit;
986 if(Sloc > S)
987 usens++;
988 else
989 usens--;
990 if(Abs(Tloc - Vsup) < (Vsup - Vinf)/5)
991 vsens++;
992 else
993 vsens--;
994 }
995 if(vparit) {
c6541a0c 996 Pts2d(i).SetCoord(S+usens*M_PI,(M_PI - T)*(vsens-1));
7fd59977 997 }
998 else {
c6541a0c 999 Pts2d(i).SetCoord(S+usens*M_PI,T+vsens*M_PI);
7fd59977 1000
1001 }
1002 }
1003 myProjIsDone = Standard_True;
1004 break;
1005 }
1006 case GeomAbs_Torus:
1007 {
1008 Standard_Real Sloc, Tloc;
1009 Standard_Integer usens = 0, vsens = 0;
1010 gp_Torus Torus = Surf->Torus();
1011 ElSLib::Parameters( Torus, Pts(1), S, T);
1012 Pts2d(1).SetCoord(S,T);
1013 for ( i = 2 ; i <= NbOfPnts ; i++) {
1014 Sloc = S; Tloc = T;
1015 ElSLib::Parameters( Torus, Pts(i), S, T);
eafb234b 1016 if(Abs(Sloc - S) > M_PI) {
7fd59977 1017 if(Sloc > S)
1018 usens++;
1019 else
1020 usens--;
eafb234b 1021 }
1022 if(Abs(Tloc - T) > M_PI) {
7fd59977 1023 if(Tloc > T)
1024 vsens++;
1025 else
1026 vsens--;
eafb234b 1027 }
c6541a0c 1028 Pts2d(i).SetCoord(S+usens*2*M_PI,T+vsens*2*M_PI);
7fd59977 1029 }
1030 myProjIsDone = Standard_True;
1031 break;
1032 }
1033 default:
9775fa61 1034 throw Standard_NoSuchObject("ProjLib_ComputeApproxOnPolarSurface::BuildInitialCurve2d");
7fd59977 1035 }
1036 }
1037 else {
7fd59977 1038 myProjIsDone = Standard_False;
1039 Standard_Real Dist2Min = 1.e+200, u = 0., v = 0.;
1040 gp_Pnt pntproj;
1041
1042 TColgp_SequenceOfPnt2d Sols;
1043 Standard_Boolean areManyZeros = Standard_False;
1044
1045 Curve->D0(Param.Value(1), pntproj) ;
1046 Extrema_ExtPS aExtPS(pntproj, Surf->Surface(), TolU, TolV) ;
f7e3c52f 1047 Standard_Real aMinSqDist = RealLast();
1048 if (aExtPS.IsDone())
1049 {
1050 for (i = 1; i <= aExtPS.NbExt(); i++)
1051 {
1052 Standard_Real aSqDist = aExtPS.SquareDistance(i);
1053 if (aSqDist < aMinSqDist)
1054 aMinSqDist = aSqDist;
1055 }
1056 }
1057 if (aMinSqDist > DistTol3d * DistTol3d) //try to project with less tolerance
1058 {
1059 TolU = Min(TolU, Precision::PConfusion());
1060 TolV = Min(TolV, Precision::PConfusion());
1061 aExtPS.Initialize(Surf->Surface(),
1062 Surf->Surface().FirstUParameter(), Surf->Surface().LastUParameter(),
1063 Surf->Surface().FirstVParameter(), Surf->Surface().LastVParameter(),
1064 TolU, TolV);
1065 aExtPS.Perform(pntproj);
1066 }
7fd59977 1067
1068 if( aExtPS.IsDone() && aExtPS.NbExt() >= 1 ) {
1069
1070 Standard_Integer GoodValue = 1;
1071
1072 for ( i = 1 ; i <= aExtPS.NbExt() ; i++ ) {
1073 if( aExtPS.SquareDistance(i) < DistTol3d * DistTol3d ) {
1074 if( aExtPS.SquareDistance(i) <= 1.e-18 ) {
1075 aExtPS.Point(i).Parameter(u,v);
1076 gp_Pnt2d p2d(u,v);
1077 Standard_Boolean isSame = Standard_False;
1078 for( j = 1; j <= Sols.Length(); j++ ) {
1079 if( p2d.SquareDistance( Sols.Value(j) ) <= 1.e-18 ) {
1080 isSame = Standard_True;
1081 break;
1082 }
1083 }
1084 if( !isSame ) Sols.Append( p2d );
1085 }
1086 if( Dist2Min > aExtPS.SquareDistance(i) ) {
1087 Dist2Min = aExtPS.SquareDistance(i);
1088 GoodValue = i;
1089 }
1090 }
1091 }
1092
1093 if( Sols.Length() > 1 ) areManyZeros = Standard_True;
1094
1095 if( Dist2Min <= DistTol3d * DistTol3d) {
1096 if( !areManyZeros ) {
1097 aExtPS.Point(GoodValue).Parameter(u,v);
1098 Pts2d(1).SetCoord(u,v);
1099 myProjIsDone = Standard_True;
1100 }
1101 else {
1102 Standard_Integer nbSols = Sols.Length();
1103 Standard_Real Dist2Max = -1.e+200;
1104 for( i = 1; i <= nbSols; i++ ) {
1105 const gp_Pnt2d& aP1 = Sols.Value(i);
1106 for( j = i+1; j <= nbSols; j++ ) {
1107 const gp_Pnt2d& aP2 = Sols.Value(j);
1108 Standard_Real aDist2 = aP1.SquareDistance(aP2);
1109 if( aDist2 > Dist2Max ) Dist2Max = aDist2;
1110 }
1111 }
1112 Standard_Real aMaxT2 = Max(TolU,TolV);
1113 aMaxT2 *= aMaxT2;
1114 if( Dist2Max > aMaxT2 ) {
1115 Standard_Integer tPp = 0;
1116 for( i = 1; i <= 5; i++ ) {
1117 Standard_Integer nbExtOk = 0;
1118 Standard_Integer indExt = 0;
1119 Standard_Integer iT = 1 + (NbOfPnts - 1)/5*i;
1120 Curve->D0( Param.Value(iT), pntproj );
1121 Extrema_ExtPS aTPS( pntproj, Surf->Surface(), TolU, TolV );
1122 Dist2Min = 1.e+200;
1123 if( aTPS.IsDone() && aTPS.NbExt() >= 1 ) {
1124 for( j = 1 ; j <= aTPS.NbExt() ; j++ ) {
1125 if( aTPS.SquareDistance(j) < DistTol3d * DistTol3d ) {
1126 nbExtOk++;
1127 if( aTPS.SquareDistance(j) < Dist2Min ) {
1128 Dist2Min = aTPS.SquareDistance(j);
1129 indExt = j;
1130 }
1131 }
1132 }
1133 }
1134 if( nbExtOk == 1 ) {
1135 tPp = iT;
1136 aTPS.Point(indExt).Parameter(u,v);
1137 break;
1138 }
1139 }
1140
1141 if( tPp != 0 ) {
1142 gp_Pnt2d aPp = gp_Pnt2d(u,v);
1143 gp_Pnt2d aPn;
1144 j = 1;
1145 Standard_Boolean isFound = Standard_False;
1146 while( !isFound ) {
1147 Curve->D0( Param.Value(tPp+j), pntproj );
1148 Extrema_ExtPS aTPS( pntproj, Surf->Surface(), TolU, TolV );
1149 Dist2Min = 1.e+200;
1150 Standard_Integer indExt = 0;
1151 if( aTPS.IsDone() && aTPS.NbExt() >= 1 ) {
1152 for( i = 1 ; i <= aTPS.NbExt() ; i++ ) {
1153 if( aTPS.SquareDistance(i) < DistTol3d * DistTol3d && aTPS.SquareDistance(i) < Dist2Min ) {
1154 Dist2Min = aTPS.SquareDistance(i);
1155 indExt = i;
1156 isFound = Standard_True;
1157 }
1158 }
1159 }
1160 if( isFound ) {
1161 aTPS.Point(indExt).Parameter(u,v);
1162 aPn = gp_Pnt2d(u,v);
1163 break;
1164 }
1165 j++;
1166 if( (tPp+j) > NbOfPnts ) break;
1167 }
1168
1169 if( isFound ) {
1170 gp_Vec2d atV(aPp,aPn);
1171 Standard_Boolean isChosen = Standard_False;
1172 for( i = 1; i <= nbSols; i++ ) {
1173 const gp_Pnt2d& aP1 = Sols.Value(i);
1174 gp_Vec2d asV(aP1,aPp);
1175 if( asV.Dot(atV) > 0. ) {
1176 isChosen = Standard_True;
1177 Pts2d(1).SetCoord(aP1.X(),aP1.Y());
1178 myProjIsDone = Standard_True;
1179 break;
1180 }
1181 }
1182 if( !isChosen ) {
1183 aExtPS.Point(GoodValue).Parameter(u,v);
1184 Pts2d(1).SetCoord(u,v);
1185 myProjIsDone = Standard_True;
1186 }
1187 }
1188 else {
1189 aExtPS.Point(GoodValue).Parameter(u,v);
1190 Pts2d(1).SetCoord(u,v);
1191 myProjIsDone = Standard_True;
1192 }
1193 }
1194 else {
1195 aExtPS.Point(GoodValue).Parameter(u,v);
1196 Pts2d(1).SetCoord(u,v);
1197 myProjIsDone = Standard_True;
1198 }
1199 }
1200 else {
1201 aExtPS.Point(GoodValue).Parameter(u,v);
1202 Pts2d(1).SetCoord(u,v);
1203 myProjIsDone = Standard_True;
1204 }
1205 }
1206 }
1207
1208 // calculate the following points with GenLocate_ExtPS
1209 // (and store the result and each parameter in a sequence)
1210 Standard_Integer usens = 0, vsens = 0;
1211 // to know the position relatively to the period
ee9451ab 1212 Standard_Real U0 = u, V0 = v, U1 = u, V1 = v;
7fd59977 1213 // U0 and V0 are the points in the initialized period
1214 // (period with u and v),
1215 // U1 and V1 are the points for construction of poles
1216
7fd59977 1217 for ( i = 2 ; i <= NbOfPnts ; i++)
1218 if(myProjIsDone) {
1219 myProjIsDone = Standard_False;
1220 Dist2Min = RealLast();
1221 Curve->D0(Param.Value(i), pntproj);
e5260e1d 1222 Extrema_GenLocateExtPS aLocateExtPS(Surf->Surface(), TolU, TolV);
1223 aLocateExtPS.Perform(pntproj, U0, V0);
1224
7fd59977 1225 if (aLocateExtPS.IsDone())
15173a08 1226 {
1227 if (aLocateExtPS.SquareDistance() < DistTol3d * DistTol3d)
1228 { //OCC217
1229 //if (aLocateExtPS.SquareDistance() < Tol3d * Tol3d) {
7fd59977 1230 (aLocateExtPS.Point()).Parameter(U0,V0);
1231 U1 = U0 + usens*uperiod;
1232 V1 = V0 + vsens*vperiod;
1233 Pts2d(i).SetCoord(U1,V1);
1234 myProjIsDone = Standard_True;
1235 }
15173a08 1236 else
1237 {
1238 Extrema_ExtPS aGlobalExtr(pntproj, Surf->Surface(), TolU, TolV);
1239 if (aGlobalExtr.IsDone())
1240 {
1241 Standard_Real LocalMinSqDist = RealLast();
1242 Standard_Integer imin = 0;
1243 for (Standard_Integer isol = 1; isol <= aGlobalExtr.NbExt(); isol++)
1244 {
1245 Standard_Real aSqDist = aGlobalExtr.SquareDistance(isol);
1246 if (aSqDist < LocalMinSqDist)
1247 {
1248 LocalMinSqDist = aSqDist;
1249 imin = isol;
1250 }
1251 }
1252 if (LocalMinSqDist < DistTol3d * DistTol3d)
1253 {
1254 Standard_Real LocalU, LocalV;
1255 aGlobalExtr.Point(imin).Parameter(LocalU, LocalV);
1256 if (uperiod > 0. && Abs(U0 - LocalU) >= uperiod/2.)
1257 {
1258 if (LocalU > U0)
1259 usens = -1;
1260 else
1261 usens = 1;
1262 }
1263 if (vperiod > 0. && Abs(V0 - LocalV) >= vperiod/2.)
1264 {
1265 if (LocalV > V0)
1266 vsens = -1;
1267 else
1268 vsens = 1;
1269 }
1270 U0 = LocalU; V0 = LocalV;
1271 U1 = U0 + usens*uperiod;
1272 V1 = V0 + vsens*vperiod;
1273 Pts2d(i).SetCoord(U1,V1);
1274 myProjIsDone = Standard_True;
7a8c6a36 1275
1276 if((i == 2) && (!IsEqual(uperiod, 0.0) || !IsEqual(vperiod, 0.0)))
1277 {//Make 1st point more precise for periodic surfaces
1278 const Standard_Integer aSize = 3;
1279 const gp_Pnt2d aP(Pts2d(2));
1280 Standard_Real aUpar[aSize], aVpar[aSize];
1281 Pts2d(1).Coord(aUpar[1], aVpar[1]);
1282 aUpar[0] = aUpar[1] - uperiod;
1283 aUpar[2] = aUpar[1] + uperiod;
1284 aVpar[0] = aVpar[1] - vperiod;
1285 aVpar[2] = aVpar[1] + vperiod;
1286
1287 Standard_Real aSQdistMin = RealLast();
1288 Standard_Integer aBestUInd = 1, aBestVInd = 1;
1289 const Standard_Integer aSizeU = IsEqual(uperiod, 0.0) ? 1 : aSize,
1290 aSizeV = IsEqual(vperiod, 0.0) ? 1 : aSize;
1291 for(Standard_Integer uInd = 0; uInd < aSizeU; uInd++)
1292 {
1293 for(Standard_Integer vInd = 0; vInd < aSizeV; vInd++)
1294 {
1295 Standard_Real aSQdist = aP.SquareDistance(gp_Pnt2d(aUpar[uInd], aVpar[vInd]));
1296 if(aSQdist < aSQdistMin)
1297 {
1298 aSQdistMin = aSQdist;
1299 aBestUInd = uInd;
1300 aBestVInd = vInd;
1301 }
1302 }
1303 }
1304
1305 Pts2d(1).SetCoord(aUpar[aBestUInd], aVpar[aBestVInd]);
1306 }//if(i == 2) condition
15173a08 1307 }
1308 }
1309 }
1310 }
7fd59977 1311 if(!myProjIsDone && uperiod) {
51740958 1312 Standard_Real aUinf, aUsup, Uaux;
1313 aUinf = Surf->Surface().FirstUParameter();
1314 aUsup = Surf->Surface().LastUParameter();
1315 if((aUsup - U0) > (U0 - aUinf))
1316 Uaux = 2*aUinf - U0 + uperiod;
7fd59977 1317 else
51740958 1318 Uaux = 2*aUsup - U0 - uperiod;
e5260e1d 1319
1320 Extrema_GenLocateExtPS locext(Surf->Surface(), TolU, TolV);
1321 locext.Perform(pntproj, Uaux, V0);
1322
7fd59977 1323 if (locext.IsDone())
1324 if (locext.SquareDistance() < DistTol3d * DistTol3d) { //OCC217
1325 //if (locext.SquareDistance() < Tol3d * Tol3d) {
1326 (locext.Point()).Parameter(u,v);
51740958 1327 if((aUsup - U0) > (U0 - aUinf))
7fd59977 1328 usens--;
1329 else
1330 usens++;
1331 U0 = u; V0 = v;
1332 U1 = U0 + usens*uperiod;
1333 V1 = V0 + vsens*vperiod;
1334 Pts2d(i).SetCoord(U1,V1);
1335 myProjIsDone = Standard_True;
1336 }
1337 }
1338 if(!myProjIsDone && vperiod) {
51740958 1339 Standard_Real aVinf, aVsup, Vaux;
1340 aVinf = Surf->Surface().FirstVParameter();
1341 aVsup = Surf->Surface().LastVParameter();
1342 if((aVsup - V0) > (V0 - aVinf))
1343 Vaux = 2*aVinf - V0 + vperiod;
7fd59977 1344 else
51740958 1345 Vaux = 2*aVsup - V0 - vperiod;
e5260e1d 1346
1347 Extrema_GenLocateExtPS locext(Surf->Surface(), TolU, TolV);
1348 locext.Perform(pntproj, U0, Vaux);
1349
7fd59977 1350 if (locext.IsDone())
1351 if (locext.SquareDistance() < DistTol3d * DistTol3d) { //OCC217
1352 //if (locext.SquareDistance() < Tol3d * Tol3d) {
1353 (locext.Point()).Parameter(u,v);
51740958 1354 if((aVsup - V0) > (V0 - aVinf))
7fd59977 1355 vsens--;
1356 else
1357 vsens++;
1358 U0 = u; V0 = v;
1359 U1 = U0 + usens*uperiod;
1360 V1 = V0 + vsens*vperiod;
1361 Pts2d(i).SetCoord(U1,V1);
1362 myProjIsDone = Standard_True;
1363 }
1364 }
1365 if(!myProjIsDone && uperiod && vperiod) {
1366 Standard_Real Uaux, Vaux;
1367 if((Usup - U0) > (U0 - Uinf))
1368 Uaux = 2*Uinf - U0 + uperiod;
1369 else
1370 Uaux = 2*Usup - U0 - uperiod;
1371 if((Vsup - V0) > (V0 - Vinf))
1372 Vaux = 2*Vinf - V0 + vperiod;
1373 else
1374 Vaux = 2*Vsup - V0 - vperiod;
e5260e1d 1375
1376 Extrema_GenLocateExtPS locext(Surf->Surface(), TolU, TolV);
1377 locext.Perform(pntproj, Uaux, Vaux);
1378
7fd59977 1379 if (locext.IsDone())
1380 if (locext.SquareDistance() < DistTol3d * DistTol3d) {
1381 //if (locext.SquareDistance() < Tol3d * Tol3d) {
1382 (locext.Point()).Parameter(u,v);
1383 if((Usup - U0) > (U0 - Uinf))
1384 usens--;
1385 else
1386 usens++;
1387 if((Vsup - V0) > (V0 - Vinf))
1388 vsens--;
1389 else
1390 vsens++;
1391 U0 = u; V0 = v;
1392 U1 = U0 + usens*uperiod;
1393 V1 = V0 + vsens*vperiod;
1394 Pts2d(i).SetCoord(U1,V1);
1395 myProjIsDone = Standard_True;
1396 }
1397 }
1398 if(!myProjIsDone) {
1399 Extrema_ExtPS ext(pntproj, Surf->Surface(), TolU, TolV) ;
1400 if (ext.IsDone()) {
1401 Dist2Min = ext.SquareDistance(1);
51740958 1402 Standard_Integer aGoodValue = 1;
7fd59977 1403 for ( j = 2 ; j <= ext.NbExt() ; j++ )
1404 if( Dist2Min > ext.SquareDistance(j)) {
1405 Dist2Min = ext.SquareDistance(j);
51740958 1406 aGoodValue = j;
7fd59977 1407 }
1408 if (Dist2Min < DistTol3d * DistTol3d) {
1409 //if (Dist2Min < Tol3d * Tol3d) {
51740958 1410 (ext.Point(aGoodValue)).Parameter(u,v);
eafb234b 1411 if(uperiod) {
7fd59977 1412 if((U0 - u) > (2*uperiod/3)) {
1413 usens++;
1414 }
1415 else
1416 if((u - U0) > (2*uperiod/3)) {
1417 usens--;
1418 }
eafb234b 1419 }
1420 if(vperiod) {
7fd59977 1421 if((V0 - v) > (vperiod/2)) {
1422 vsens++;
1423 }
1424 else
1425 if((v - V0) > (vperiod/2)) {
1426 vsens--;
1427 }
eafb234b 1428 }
7fd59977 1429 U0 = u; V0 = v;
1430 U1 = U0 + usens*uperiod;
1431 V1 = V0 + vsens*vperiod;
1432 Pts2d(i).SetCoord(U1,V1);
1433 myProjIsDone = Standard_True;
1434 }
1435 }
1436 }
1437 }
1438 else break;
1439 }
1440 }
1441 // -- Pnts2d is transformed into Geom2d_BSplineCurve, with the help of Param and Mult
1442 if(myProjIsDone) {
1443 myBSpline = new Geom2d_BSplineCurve(Pts2d,Param,Mult,1);
ee9451ab 1444 //jgv: put the curve into parametric range
1445 gp_Pnt2d MidPoint = myBSpline->Value(0.5*(myBSpline->FirstParameter() + myBSpline->LastParameter()));
1446 Standard_Real TestU = MidPoint.X(), TestV = MidPoint.Y();
1447 Standard_Real sense = 0.;
1448 if (uperiod)
1449 {
1450 if (TestU < Uinf - TolU)
1451 sense = 1.;
1452 else if (TestU > Usup + TolU)
1453 sense = -1;
1454 while (TestU < Uinf - TolU || TestU > Usup + TolU)
1455 TestU += sense * uperiod;
1456 }
1457 if (vperiod)
1458 {
1459 sense = 0.;
1460 if (TestV < Vinf - TolV)
1461 sense = 1.;
1462 else if (TestV > Vsup + TolV)
1463 sense = -1.;
1464 while (TestV < Vinf - TolV || TestV > Vsup + TolV)
1465 TestV += sense * vperiod;
1466 }
1467 gp_Vec2d Offset(TestU - MidPoint.X(), TestV - MidPoint.Y());
1468 if (Abs(Offset.X()) > gp::Resolution() ||
1469 Abs(Offset.Y()) > gp::Resolution())
1470 myBSpline->Translate(Offset);
1471 //////////////////////////////////////////
7fd59977 1472 Geom2dAdaptor_Curve GAC(myBSpline);
1473 Handle(Adaptor2d_HCurve2d) IC2d = new Geom2dAdaptor_HCurve(GAC);
0797d9d3 1474#ifdef OCCT_DEBUG
7fd59977 1475// char name [100];
1476// sprintf(name,"%s_%d","build",compteur++);
1477// DrawTrSurf::Set(name,myBSpline);
1478#endif
1479 return IC2d;
1480 }
1481 else {
1482// Modified by Sergey KHROMOV - Thu Apr 18 10:57:50 2002 Begin
1483// Standard_NoSuchObject_Raise_if(1,"ProjLib_Compu: build echec");
1484// Modified by Sergey KHROMOV - Thu Apr 18 10:57:51 2002 End
1485 return Handle(Adaptor2d_HCurve2d)();
1486 }
d3f26155 1487// myProjIsDone = Standard_False;
7fd59977 1488// Modified by Sergey KHROMOV - Thu Apr 18 10:58:01 2002 Begin
1489// Standard_NoSuchObject_Raise_if(1,"ProjLib_ComputeOnPS: build echec");
1490// Modified by Sergey KHROMOV - Thu Apr 18 10:58:02 2002 End
7fd59977 1491}
1492
1493
1494
1495
1496//=======================================================================
1497//function : ProjLib_ProjectUsingInitialCurve2d
1498//purpose :
1499//=======================================================================
1500Handle(Geom2d_BSplineCurve)
1501 ProjLib_ComputeApproxOnPolarSurface::
1502 ProjectUsingInitialCurve2d(const Handle(Adaptor3d_HCurve)& Curve,
1503 const Handle(Adaptor3d_HSurface)& Surf,
1504 const Handle(Adaptor2d_HCurve2d)& InitCurve2d)
1505{
1506 //OCC217
1507 Standard_Real Tol3d = myTolerance;
1508 Standard_Real DistTol3d = 1.0*Tol3d;
1509 Standard_Real TolU = Surf->UResolution(Tol3d), TolV = Surf->VResolution(Tol3d);
ea3b6e72 1510 Standard_Real Tol2d = Max(Sqrt(TolU*TolU + TolV*TolV), Precision::PConfusion());
7fd59977 1511
1512 Standard_Integer i;
1513 GeomAbs_SurfaceType TheTypeS = Surf->GetType();
1514 GeomAbs_CurveType TheTypeC = Curve->GetType();
1515// Handle(Standard_Type) TheTypeS = Surf->DynamicType();
1516// Handle(Standard_Type) TheTypeC = Curve->DynamicType(); // si on a :
1517// if(TheTypeS == STANDARD_TYPE(Geom_BSplineSurface)) {
1518 if(TheTypeS == GeomAbs_Plane) {
1519 Standard_Real S, T;
1520 gp_Pln Plane = Surf->Plane();
1521 if(TheTypeC == GeomAbs_BSplineCurve) {
1522 Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
1523 TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1524 for(i = 1;i <= Curve->NbPoles();i++) {
1525 ElSLib::Parameters( Plane, BSC->Pole(i), S, T);
1526 Poles2d(i).SetCoord(S,T);
1527 }
1528 TColStd_Array1OfReal Knots(1, BSC->NbKnots());
1529 BSC->Knots(Knots);
1530 TColStd_Array1OfInteger Mults(1, BSC->NbKnots());
1531 BSC->Multiplicities(Mults);
1532 if(BSC->IsRational()) {
1533 TColStd_Array1OfReal Weights(1, BSC->NbPoles());
1534 BSC->Weights(Weights);
1535 return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1536 BSC->Degree(), BSC->IsPeriodic()) ;
1537 }
1538 return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1539 BSC->Degree(), BSC->IsPeriodic()) ;
1540
1541 }
1542 if(TheTypeC == GeomAbs_BezierCurve) {
1543 Handle(Geom_BezierCurve) BC = Curve->Bezier();
1544 TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1545 for(i = 1;i <= Curve->NbPoles();i++) {
1546 ElSLib::Parameters( Plane, BC->Pole(i), S, T);
1547 Poles2d(i).SetCoord(S,T);
1548 }
1549 TColStd_Array1OfReal Knots(1, 2);
1550 Knots.SetValue(1,0.0);
1551 Knots.SetValue(2,1.0);
1552 TColStd_Array1OfInteger Mults(1, 2);
1553 Mults.Init(BC->NbPoles());
1554 if(BC->IsRational()) {
1555 TColStd_Array1OfReal Weights(1, BC->NbPoles());
1556 BC->Weights(Weights);
1557 return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1558 BC->Degree(), BC->IsPeriodic()) ;
1559 }
1560 return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1561 BC->Degree(), BC->IsPeriodic()) ;
1562 }
1563 }
1564 if(TheTypeS == GeomAbs_BSplineSurface) {
1565 Handle(Geom_BSplineSurface) BSS = Surf->BSpline();
1566 if((BSS->MaxDegree() == 1) &&
1567 (BSS->NbUPoles() == 2) &&
1568 (BSS->NbVPoles() == 2)) {
1569 gp_Pnt p11 = BSS->Pole(1,1);
1570 gp_Pnt p12 = BSS->Pole(1,2);
1571 gp_Pnt p21 = BSS->Pole(2,1);
1572 gp_Pnt p22 = BSS->Pole(2,2);
1573 gp_Vec V1(p11,p12);
1574 gp_Vec V2(p21,p22);
c6541a0c
D
1575 if(V1.IsEqual(V2,Tol3d,Tol3d/(p11.Distance(p12)*180/M_PI))){ //OCC217
1576 //if(V1.IsEqual(V2,myTolerance,myTolerance/(p11.Distance(p12)*180/M_PI))){
7fd59977 1577 // so the polar surface is plane
1578 // and if it is enough to projet the poles of Curve
1579 Standard_Integer Dist2Min = IntegerLast();
1580 Standard_Real u,v;
1581 //OCC217
1582 //Standard_Real TolU = Surf->UResolution(myTolerance)
1583 // , TolV = Surf->VResolution(myTolerance);
1584// gp_Pnt pntproj;
1585 if(TheTypeC == GeomAbs_BSplineCurve) {
1586 Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
1587 TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1588 for(i = 1;i <= Curve->NbPoles();i++) {
1589 myProjIsDone = Standard_False;
1590 Dist2Min = IntegerLast();
e5260e1d 1591
1592 Extrema_GenLocateExtPS extrloc(Surf->Surface(), TolU, TolV);
1593 extrloc.Perform(BSC->Pole(i), (p11.X()+p22.X())/2, (p11.Y()+p22.Y())/2);
1594
7fd59977 1595 if (extrloc.IsDone()) {
1596 Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1597 if (Dist2Min < DistTol3d * DistTol3d) { //OCC217
1598 //if (Dist2Min < myTolerance * myTolerance) {
1599 (extrloc.Point()).Parameter(u,v);
1600 Poles2d(i).SetCoord(u,v);
1601 myProjIsDone = Standard_True;
1602 }
1603 else break;
1604 }
1605 else break;
1606 if(!myProjIsDone)
1607 break;
1608 }
1609 if(myProjIsDone) {
1610 TColStd_Array1OfReal Knots(1, BSC->NbKnots());
1611 BSC->Knots(Knots);
1612 TColStd_Array1OfInteger Mults(1, BSC->NbKnots());
1613 BSC->Multiplicities(Mults);
1614 if(BSC->IsRational()) {
1615 TColStd_Array1OfReal Weights(1, BSC->NbPoles());
1616 BSC->Weights(Weights);
1617 return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1618 BSC->Degree(), BSC->IsPeriodic()) ;
1619 }
1620 return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1621 BSC->Degree(), BSC->IsPeriodic()) ;
1622
1623
1624 }
1625 }
1626 if(TheTypeC == GeomAbs_BezierCurve) {
1627 Handle(Geom_BezierCurve) BC = Curve->Bezier();
1628 TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1629 for(i = 1;i <= Curve->NbPoles();i++) {
1630 Dist2Min = IntegerLast();
e5260e1d 1631
1632 Extrema_GenLocateExtPS extrloc(Surf->Surface(), TolU, TolV);
1633 extrloc.Perform(BC->Pole(i), 0.5, 0.5);
1634
7fd59977 1635 if (extrloc.IsDone()) {
1636 Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1637 if (Dist2Min < DistTol3d * DistTol3d) { //OCC217
1638 //if (Dist2Min < myTolerance * myTolerance) {
1639 (extrloc.Point()).Parameter(u,v);
1640 Poles2d(i).SetCoord(u,v);
1641 myProjIsDone = Standard_True;
1642 }
1643 else break;
1644 }
1645 else break;
1646 if(myProjIsDone)
1647 myProjIsDone = Standard_False;
1648 else break;
1649 }
1650 if(myProjIsDone) {
1651 TColStd_Array1OfReal Knots(1, 2);
1652 Knots.SetValue(1,0.0);
1653 Knots.SetValue(2,1.0);
1654 TColStd_Array1OfInteger Mults(1, 2);
1655 Mults.Init(BC->NbPoles());
1656 if(BC->IsRational()) {
1657 TColStd_Array1OfReal Weights(1, BC->NbPoles());
1658 BC->Weights(Weights);
1659 return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1660 BC->Degree(), BC->IsPeriodic()) ;
1661 }
1662 return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1663 BC->Degree(), BC->IsPeriodic()) ;
1664 }
1665 }
1666 }
1667 }
1668 }
1669 else if(TheTypeS == GeomAbs_BezierSurface) {
1670 Handle(Geom_BezierSurface) BS = Surf->Bezier();
1671 if((BS->MaxDegree() == 1) &&
1672 (BS->NbUPoles() == 2) &&
1673 (BS->NbVPoles() == 2)) {
1674 gp_Pnt p11 = BS->Pole(1,1);
1675 gp_Pnt p12 = BS->Pole(1,2);
1676 gp_Pnt p21 = BS->Pole(2,1);
1677 gp_Pnt p22 = BS->Pole(2,2);
1678 gp_Vec V1(p11,p12);
1679 gp_Vec V2(p21,p22);
c6541a0c
D
1680 if(V1.IsEqual(V2,Tol3d,Tol3d/(p11.Distance(p12)*180/M_PI))){ //OCC217
1681 //if (V1.IsEqual(V2,myTolerance,myTolerance/(p11.Distance(p12)*180/M_PI))){
7fd59977 1682 // and if it is enough to project the poles of Curve
1683 Standard_Integer Dist2Min = IntegerLast();
1684 Standard_Real u,v;
1685 //OCC217
1686 //Standard_Real TolU = Surf->UResolution(myTolerance)
1687 // , TolV = Surf->VResolution(myTolerance);
1688
1689// gp_Pnt pntproj;
1690 if(TheTypeC == GeomAbs_BSplineCurve) {
1691 Handle(Geom_BSplineCurve) BSC = Curve->BSpline();
1692 TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1693 for(i = 1;i <= Curve->NbPoles();i++) {
1694 myProjIsDone = Standard_False;
1695 Dist2Min = IntegerLast();
e5260e1d 1696
1697 Extrema_GenLocateExtPS extrloc(Surf->Surface(), TolU, TolV);
1698 extrloc.Perform(BSC->Pole(i), (p11.X()+p22.X())/2, (p11.Y()+p22.Y())/2);
1699
7fd59977 1700 if (extrloc.IsDone()) {
1701 Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1702 if (Dist2Min < DistTol3d * DistTol3d) { //OCC217
1703 //if (Dist2Min < myTolerance * myTolerance) {
1704 (extrloc.Point()).Parameter(u,v);
1705 Poles2d(i).SetCoord(u,v);
1706 myProjIsDone = Standard_True;
1707 }
1708 else break;
1709 }
1710 else break;
1711 if(!myProjIsDone)
1712 break;
1713 }
1714 if(myProjIsDone) {
1715 TColStd_Array1OfReal Knots(1, BSC->NbKnots());
1716 BSC->Knots(Knots);
1717 TColStd_Array1OfInteger Mults(1, BSC->NbKnots());
1718 BSC->Multiplicities(Mults);
1719 if(BSC->IsRational()) {
1720 TColStd_Array1OfReal Weights(1, BSC->NbPoles());
1721 BSC->Weights(Weights);
1722 return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1723 BSC->Degree(), BSC->IsPeriodic()) ;
1724 }
1725 return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1726 BSC->Degree(), BSC->IsPeriodic()) ;
1727
1728
1729 }
1730 }
1731 if(TheTypeC == GeomAbs_BezierCurve) {
1732 Handle(Geom_BezierCurve) BC = Curve->Bezier();
1733 TColgp_Array1OfPnt2d Poles2d(1,Curve->NbPoles());
1734 for(i = 1;i <= Curve->NbPoles();i++) {
1735 Dist2Min = IntegerLast();
e5260e1d 1736
1737 Extrema_GenLocateExtPS extrloc(Surf->Surface(), TolU, TolV);
1738 extrloc.Perform(BC->Pole(i), 0.5, 0.5);
1739
7fd59977 1740 if (extrloc.IsDone()) {
1741 Dist2Min = (Standard_Integer ) extrloc.SquareDistance();
1742 if (Dist2Min < DistTol3d * DistTol3d) { //OCC217
1743 //if (Dist2Min < myTolerance * myTolerance) {
1744 (extrloc.Point()).Parameter(u,v);
1745 Poles2d(i).SetCoord(u,v);
1746 myProjIsDone = Standard_True;
1747 }
1748 else break;
1749 }
1750 else break;
1751 if(myProjIsDone)
1752 myProjIsDone = Standard_False;
1753 else break;
1754 }
1755 if(myProjIsDone) {
1756 TColStd_Array1OfReal Knots(1, 2);
1757 Knots.SetValue(1,0.0);
1758 Knots.SetValue(2,1.0);
1759 TColStd_Array1OfInteger Mults(1, 2);
1760 Mults.Init(BC->NbPoles());
1761 if(BC->IsRational()) {
1762 TColStd_Array1OfReal Weights(1, BC->NbPoles());
1763 BC->Weights(Weights);
1764 return new Geom2d_BSplineCurve(Poles2d, Weights, Knots, Mults,
1765 BC->Degree(), BC->IsPeriodic()) ;
1766 }
1767 return new Geom2d_BSplineCurve(Poles2d, Knots, Mults,
1768 BC->Degree(), BC->IsPeriodic()) ;
1769 }
1770 }
1771 }
1772 }
1773 }
1774
1775 ProjLib_PolarFunction F(Curve, Surf, InitCurve2d, Tol3d) ; //OCC217
1776 //ProjLib_PolarFunction F(Curve, Surf, InitCurve2d, myTolerance) ;
1777
0797d9d3 1778#ifdef OCCT_DEBUG
7fd59977 1779 Standard_Integer Nb = 50;
1780
1781 Standard_Real U, U1, U2;
1782 U1 = F.FirstParameter();
1783 U2 = F.LastParameter();
1784
1785 TColgp_Array1OfPnt2d DummyPoles(1,Nb+1);
1786 TColStd_Array1OfReal DummyKnots(1,Nb+1);
1787 TColStd_Array1OfInteger DummyMults(1,Nb+1);
1788 DummyMults.Init(1);
1789 DummyMults(1) = 2;
1790 DummyMults(Nb+1) = 2;
1791 for (Standard_Integer ij = 0; ij <= Nb; ij++) {
1792 U = (Nb-ij)*U1 + ij*U2;
1793 U /= Nb;
1794 DummyPoles(ij+1) = F.Value(U);
1795 DummyKnots(ij+1) = ij;
1796 }
1797 Handle(Geom2d_BSplineCurve) DummyC2d =
1798 new Geom2d_BSplineCurve(DummyPoles, DummyKnots, DummyMults, 1);
7fd59977 1799#ifdef DRAW
96a95605 1800 Standard_CString Temp = "bs2d";
7fd59977 1801 DrawTrSurf::Set(Temp,DummyC2d);
1802#endif
1803// DrawTrSurf::Set((Standard_CString ) "bs2d",DummyC2d);
1804 Handle(Geom2dAdaptor_HCurve) DDD =
1805 Handle(Geom2dAdaptor_HCurve)::DownCast(InitCurve2d);
1806
7fd59977 1807#ifdef DRAW
96a95605 1808 Temp = "initc2d";
7fd59977 1809 DrawTrSurf::Set(Temp,DDD->ChangeCurve2d().Curve());
1810#endif
1811// DrawTrSurf::Set((Standard_CString ) "initc2d",DDD->ChangeCurve2d().Curve());
1812#endif
1813
1814 Standard_Integer Deg1,Deg2;
1815// Deg1 = 8;
1816// Deg2 = 8;
1817 Deg1 = 2; //IFV
1818 Deg2 = 8; //IFV
1819
1820 Approx_FitAndDivide2d Fit(F,Deg1,Deg2,Tol3d,Tol2d, //OCC217
1821 //Approx_FitAndDivide2d Fit(F,Deg1,Deg2,myTolerance,myTolerance,
1822 Standard_True);
1823
1824 if(Fit.IsAllApproximated()) {
51740958 1825 Standard_Integer j;
7fd59977 1826 Standard_Integer NbCurves = Fit.NbMultiCurves();
1827 Standard_Integer MaxDeg = 0;
1828 // To transform the MultiCurve into BSpline, it is required that all
1829 // Bezier constituing it have the same degree -> Calculation of MaxDeg
1830 Standard_Integer NbPoles = 1;
51740958 1831 for (j = 1; j <= NbCurves; j++) {
1832 Standard_Integer Deg = Fit.Value(j).Degree();
7fd59977 1833 MaxDeg = Max ( MaxDeg, Deg);
1834 }
1835
1836 NbPoles = MaxDeg * NbCurves + 1; //Tops on the BSpline
1837 TColgp_Array1OfPnt2d Poles( 1, NbPoles);
1838
1839 TColgp_Array1OfPnt2d TempPoles( 1, MaxDeg + 1);//to augment the degree
1840
1841 TColStd_Array1OfReal Knots( 1, NbCurves + 1); //Nodes of the BSpline
1842
1843 Standard_Integer Compt = 1;
1844 for (i = 1; i <= NbCurves; i++) {
1845 Fit.Parameters(i, Knots(i), Knots(i+1));
1846 AppParCurves_MultiCurve MC = Fit.Value( i); //Load the Ith Curve
1847 TColgp_Array1OfPnt2d Poles2d( 1, MC.Degree() + 1);//Retrieve the tops
1848 MC.Curve(1, Poles2d);
1849
1850 //Eventual augmentation of the degree
1851 Standard_Integer Inc = MaxDeg - MC.Degree();
1852 if ( Inc > 0) {
1853// BSplCLib::IncreaseDegree( Inc, Poles2d, PLib::NoWeights(),
0e14656b 1854 BSplCLib::IncreaseDegree( MaxDeg, Poles2d, BSplCLib::NoWeights(),
1855 TempPoles, BSplCLib::NoWeights());
7fd59977 1856 //update of tops of the PCurve
51740958 1857 for (Standard_Integer k = 1 ; k <= MaxDeg + 1; k++) {
1858 Poles.SetValue( Compt, TempPoles( k));
7fd59977 1859 Compt++;
1860 }
1861 }
1862 else {
1863 //update of tops of the PCurve
51740958 1864 for (Standard_Integer k = 1 ; k <= MaxDeg + 1; k++) {
1865 Poles.SetValue( Compt, Poles2d( k));
7fd59977 1866 Compt++;
1867 }
1868 }
1869
1870 Compt--;
1871 }
1872
1873 //update of fields of ProjLib_Approx
1874 Standard_Integer NbKnots = NbCurves + 1;
368cdde6 1875
7fd59977 1876 TColStd_Array1OfInteger Mults( 1, NbKnots);
1877 Mults.Init(MaxDeg);
1878 Mults.SetValue( 1, MaxDeg + 1);
1879 Mults.SetValue(NbKnots, MaxDeg + 1);
1880 myProjIsDone = Standard_True;
1881 Handle(Geom2d_BSplineCurve) Dummy =
1882 new Geom2d_BSplineCurve(Poles,Knots,Mults,MaxDeg);
1883
1884 // try to smoother the Curve GeomAbs_C1.
1885
1886 Standard_Boolean OK = Standard_True;
1887
1888 for (Standard_Integer ij = 2; ij < NbKnots; ij++) {
1889 OK = OK && Dummy->RemoveKnot(ij,MaxDeg-1,Tol3d); //OCC217
1890 //OK = OK && Dummy->RemoveKnot(ij,MaxDeg-1,myTolerance);
1891 }
0797d9d3 1892#ifdef OCCT_DEBUG
7fd59977 1893 if (!OK) {
1894 cout << "ProjLib_ComputeApproxOnPolarSurface : Smoothing echoue"<<endl;
1895 }
1896#endif
1897 return Dummy;
1898 }
1899 return Handle(Geom2d_BSplineCurve)();
1900}
1901
1902//=======================================================================
1903//function : BSpline
1904//purpose :
1905//=======================================================================
1906
1907Handle(Geom2d_BSplineCurve)
1908 ProjLib_ComputeApproxOnPolarSurface::BSpline() const
1909
1910{
1911// Modified by Sergey KHROMOV - Thu Apr 18 11:16:46 2002 End
1912// Standard_NoSuchObject_Raise_if
1913// (!myProjIsDone,
1914// "ProjLib_ComputeApproxOnPolarSurface:BSpline");
1915// Modified by Sergey KHROMOV - Thu Apr 18 11:16:47 2002 End
1916 return myBSpline ;
1917}
1918
1919//=======================================================================
1920//function : Curve2d
1921//purpose :
1922//=======================================================================
1923
1924Handle(Geom2d_Curve)
1925 ProjLib_ComputeApproxOnPolarSurface::Curve2d() const
1926
1927{
1928 Standard_NoSuchObject_Raise_if
1929 (!myProjIsDone,
1930 "ProjLib_ComputeApproxOnPolarSurface:2ndCurve2d");
1931 return my2ndCurve ;
1932}
1933
1934
1935//=======================================================================
1936//function : IsDone
1937//purpose :
1938//=======================================================================
1939
1940Standard_Boolean ProjLib_ComputeApproxOnPolarSurface::IsDone() const
1941
1942{
1943 return myProjIsDone;
1944}