0024023: Revamp the OCCT Handle -- downcast (automatic)
[occt.git] / src / TopOpeBRepTool / TopOpeBRepTool_CurveTool.cxx
CommitLineData
b311480e 1// Created on: 1993-06-24
2// Created by: Jean Yves LEBEY
3// Copyright (c) 1993-1999 Matra Datavision
973c2be1 4// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 5//
973c2be1 6// This file is part of Open CASCADE Technology software library.
b311480e 7//
d5f74e42 8// This library is free software; you can redistribute it and/or modify it under
9// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 10// by the Free Software Foundation, with special exception defined in the file
11// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12// distribution for complete text of the license and disclaimer of any warranty.
b311480e 13//
973c2be1 14// Alternatively, this file may be used under the terms of Open CASCADE
15// commercial license or contractual agreement.
7fd59977 16
17#include <TopOpeBRepTool_CurveTool.ixx>
18
19#include <ProjLib_ProjectedCurve.hxx>
20#include <BRepApprox_Approx.hxx>
21#include <BRepApprox_ApproxLine.hxx>
22#include <TopoDS.hxx>
23#include <TopoDS_Face.hxx>
24#include <Geom_BSplineCurve.hxx>
25#include <Geom2d_Line.hxx>
26#include <Geom2d_Circle.hxx>
27#include <Geom2d_Ellipse.hxx>
28#include <Geom2d_Parabola.hxx>
29#include <Geom2d_Hyperbola.hxx>
30#include <Geom2d_BSplineCurve.hxx>
31#include <TColStd_Array1OfReal.hxx>
32#include <GeomLib_CheckBSplineCurve.hxx>
33#include <GeomLib_Check2dBSplineCurve.hxx>
34#include <TColStd_Array1OfInteger.hxx>
35#include <BRepAdaptor_Surface.hxx>
36#include <GeomAdaptor_Curve.hxx>
37#include <BRepAdaptor_HSurface.hxx>
38#include <GeomAdaptor_HCurve.hxx>
39#include <GeomAbs_SurfaceType.hxx>
40#include <gp_Pln.hxx>
41#include <gp_Pnt2d.hxx>
42#include <ElSLib.hxx>
43#include <Standard_NotImplemented.hxx>
44#include <Standard_ProgramError.hxx>
45#include <OSD_Chronometer.hxx>
46
47#include <GeomTools_Curve2dSet.hxx>
48#include <TColStd_Array1OfBoolean.hxx>
49#include <gp_Lin.hxx>
50#include <gce_MakeLin.hxx>
51#include <gp_Lin2d.hxx>
52#include <gce_MakeLin2d.hxx>
53#include <BRepTools.hxx>
54#include <Precision.hxx>
55#include <gce_MakeCirc.hxx>
56#include <gp_Vec.hxx>
57//#include <Approx.hxx>
58#include <Geom_Surface.hxx>
59#include <TopLoc_Location.hxx>
60#include <Geom_RectangularTrimmedSurface.hxx>
61#include <BRep_Tool.hxx>
62
0797d9d3 63#ifdef OCCT_DEBUG
7fd59977 64#include <TopOpeBRepTool_KRO.hxx>
65TOPKRO KRO_CURVETOOL_APPRO("approximation");
1d0a9d4d 66extern Standard_Boolean TopOpeBRepTool_GettraceKRO();
67extern Standard_Boolean TopOpeBRepTool_GettracePCURV();
68extern Standard_Boolean TopOpeBRepTool_GettraceCHKBSPL();
7fd59977 69#endif
70//#define DRAW
71//#define IFV
72#define CurveImprovement
73#ifdef DRAW
74#include <DrawTrSurf.hxx>
ec357c5c 75#include <Geom2d_Curve.hxx>
7fd59977 76static Standard_Integer NbCalls = 0;
77#endif
78//=======================================================================
79//function : CurveTool
80//purpose :
81//=======================================================================
82
83TopOpeBRepTool_CurveTool::TopOpeBRepTool_CurveTool()
84{}
85
86//=======================================================================
87//function : CurveTool
88//purpose :
89//=======================================================================
90
91TopOpeBRepTool_CurveTool::TopOpeBRepTool_CurveTool
92(const TopOpeBRepTool_OutCurveType O)
93{
94 TopOpeBRepTool_GeomTool GT(O);
95 SetGeomTool(GT);
96}
97
98//=======================================================================
99//function : CurveTool
100//purpose :
101//=======================================================================
102
103TopOpeBRepTool_CurveTool::TopOpeBRepTool_CurveTool
104(const TopOpeBRepTool_GeomTool& GT)
105{
106 SetGeomTool(GT);
107}
108
109//=======================================================================
110//function : ChangeGeomTool
111//purpose :
112//=======================================================================
113
114TopOpeBRepTool_GeomTool& TopOpeBRepTool_CurveTool::ChangeGeomTool()
115{
116 return myGeomTool;
117}
118
119//=======================================================================
120//function : GetGeomTool
121//purpose :
122//=======================================================================
123
124const TopOpeBRepTool_GeomTool& TopOpeBRepTool_CurveTool::GetGeomTool()const
125{
126 return myGeomTool;
127}
128
129//=======================================================================
130//function : SetGeomTool
131//purpose :
132//=======================================================================
133
134void TopOpeBRepTool_CurveTool::SetGeomTool
135(const TopOpeBRepTool_GeomTool& GT)
136{
137 myGeomTool.Define(GT);
138}
139
140//-----------------------------------------------------------------------
141//function : MakePCurve
142//purpose :
143//-----------------------------------------------------------------------
144Standard_EXPORT Handle(Geom2d_Curve) MakePCurve(const ProjLib_ProjectedCurve& PC)
145{
146 Handle(Geom2d_Curve) C2D;
147 switch (PC.GetType()) {
148 case GeomAbs_Line : C2D = new Geom2d_Line(PC.Line()); break;
149 case GeomAbs_Circle : C2D = new Geom2d_Circle(PC.Circle()); break;
150 case GeomAbs_Ellipse : C2D = new Geom2d_Ellipse(PC.Ellipse()); break;
151 case GeomAbs_Parabola : C2D = new Geom2d_Parabola(PC.Parabola()); break;
152 case GeomAbs_Hyperbola : C2D = new Geom2d_Hyperbola(PC.Hyperbola()); break;
153 case GeomAbs_BSplineCurve : C2D = PC.BSpline(); break;
154 case GeomAbs_BezierCurve : case GeomAbs_OtherCurve : default :
155 Standard_NotImplemented::Raise("CurveTool::MakePCurve");
156 break;
157 }
158 return C2D;
159}
160
161//------------------------------------------------------------------
162static Standard_Boolean CheckApproxResults
163 (const BRepApprox_Approx& Approx)
164//------------------------------------------------------------------
165{
166 const AppParCurves_MultiBSpCurve& amc = Approx.Value(1);
167 Standard_Integer np = amc.NbPoles();
168 Standard_Integer nc = amc.NbCurves();
169 if (np < 2 || nc < 1) return Standard_False;
170
171 // check the knots for coincidence
172 const TColStd_Array1OfReal& knots = amc.Knots();
173 for (Standard_Integer i = knots.Lower(); i < knots.Upper(); i++) {
174 if (knots(i+1) - knots(i) <= Epsilon(knots(i))) {
175 return Standard_False;
176 }
177 }
178 return Standard_True;
179}
180
181//------------------------------------------------------------------
182static Standard_Boolean CheckPCurve
183 (const Handle(Geom2d_Curve)& aPC, const TopoDS_Face& aFace)
184//------------------------------------------------------------------
185// check if points of the pcurve are out of the face bounds
186{
187 const Standard_Integer NPoints = 23;
188 Standard_Real umin,umax,vmin,vmax;
189
190 BRepTools::UVBounds(aFace, umin, umax, vmin, vmax);
191 Standard_Real tolU = Max ((umax-umin)*0.01, Precision::Confusion());
192 Standard_Real tolV = Max ((vmax-vmin)*0.01, Precision::Confusion());
193 Standard_Real fp = aPC->FirstParameter();
194 Standard_Real lp = aPC->LastParameter();
195 Standard_Real step = (lp-fp)/(NPoints+1);
196
197 // adjust domain for periodic surfaces
198 TopLoc_Location aLoc;
199 Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace, aLoc);
200 if (aSurf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)))
201 aSurf = (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurf))->BasisSurface();
202
203 gp_Pnt2d pnt = aPC->Value((fp+lp)/2);
204 Standard_Real u,v;
205 pnt.Coord(u,v);
206
207 if (aSurf->IsUPeriodic()) {
208 Standard_Real aPer = aSurf->UPeriod();
209 Standard_Integer nshift = (Standard_Integer) ((u-umin)/aPer);
210 if (u < umin+aPer*nshift) nshift--;
211 umin += aPer*nshift;
212 umax += aPer*nshift;
213 }
214 if (aSurf->IsVPeriodic()) {
215 Standard_Real aPer = aSurf->VPeriod();
216 Standard_Integer nshift = (Standard_Integer) ((v-vmin)/aPer);
217 if (v < vmin+aPer*nshift) nshift--;
218 vmin += aPer*nshift;
219 vmax += aPer*nshift;
220 }
221
222 Standard_Integer i;
223 for (i=1; i <= NPoints; i++) {
224 Standard_Real p = fp + i * step;
225 pnt = aPC->Value(p);
226 pnt.Coord(u,v);
227 if (umin-u > tolU || u-umax > tolU ||
228 vmin-v > tolV || v-vmax > tolV)
229 return Standard_False;
230 }
231 return Standard_True;
232}
233
234//------------------------------------------------------------------
235static Handle(Geom_Curve) MakeCurve3DfromWLineApprox
236(const BRepApprox_Approx& Approx,const Standard_Integer )
237//------------------------------------------------------------------
238{
239 const AppParCurves_MultiBSpCurve& amc = Approx.Value(1);
240 Standard_Integer np = amc.NbPoles();
241 //Standard_Integer nc = amc.NbCurves();
242 TColgp_Array1OfPnt poles3d(1,np);
243 Standard_Integer ic = 1;
244 //for (ic=1; ic<=nc; ic++) {
245 //if (ic == CI) {
246 amc.Curve(ic,poles3d);
247 //}
248 //}
249
250 const TColStd_Array1OfReal& knots = amc.Knots();
251 const TColStd_Array1OfInteger& mults = amc.Multiplicities();
252 Standard_Integer degree = amc.Degree();
253 Handle(Geom_Curve) C3D = new Geom_BSplineCurve(poles3d,knots,mults,degree);
254 return C3D;
255}
256
257//------------------------------------------------------------------
258static Handle(Geom2d_Curve) MakeCurve2DfromWLineApproxAndPlane
259(const BRepApprox_Approx& Approx,const gp_Pln& Pl)
260//------------------------------------------------------------------
261{
262 const AppParCurves_MultiBSpCurve& amc = Approx.Value(1);
263 Standard_Integer np = amc.NbPoles();
264 TColgp_Array1OfPnt2d poles2d(1,np);
265 TColgp_Array1OfPnt poles3d(1,np);
266 amc.Curve(1,poles3d);
267 for(Standard_Integer i=1; i<=np; i++) {
268 Standard_Real U,V; ElSLib::Parameters(Pl,poles3d.Value(i),U,V);
269 poles2d.SetValue(i,gp_Pnt2d(U,V));
270 }
271 const TColStd_Array1OfReal& knots = amc.Knots();
272 const TColStd_Array1OfInteger& mults = amc.Multiplicities();
273 Standard_Integer degree = amc.Degree();
274 Handle(Geom2d_Curve) C2D = new Geom2d_BSplineCurve(poles2d,knots,mults,degree);
275 return C2D;
276}
277
278//------------------------------------------------------------------
279static Handle(Geom2d_Curve) MakeCurve2DfromWLineApprox
280(const BRepApprox_Approx& Approx,const Standard_Integer CI)
281//------------------------------------------------------------------
282{
283 const AppParCurves_MultiBSpCurve& amc = Approx.Value(1);
284 Standard_Integer np = amc.NbPoles();
285 TColgp_Array1OfPnt2d poles2d(1,np);
286 Standard_Integer nc = amc.NbCurves();
287 for (Standard_Integer ic=1; ic<=nc; ic++) if (ic == CI) amc.Curve(ic,poles2d);
288 const TColStd_Array1OfReal& knots = amc.Knots();
289 const TColStd_Array1OfInteger& mults = amc.Multiplicities();
290 Standard_Integer degree = amc.Degree();
291 Handle(Geom2d_Curve) C2D = new Geom2d_BSplineCurve(poles2d,knots,mults,degree);
292 return C2D;
293}
294
295
296//=======================================================================
297//function : MakeCurves
298//purpose :
299//=======================================================================
300
301Standard_Boolean TopOpeBRepTool_CurveTool::MakeCurves
302(const Standard_Real parmin, const Standard_Real parmax,
303 const Handle(Geom_Curve)& C3D,
304 const Handle(Geom2d_Curve)& PC1, const Handle(Geom2d_Curve)& PC2,
305 const TopoDS_Shape& S1, const TopoDS_Shape& S2,
306 Handle(Geom_Curve)& C3Dnew,
307 Handle(Geom2d_Curve)& PC1new, Handle(Geom2d_Curve)& PC2new,
308 Standard_Real& TolReached3d, Standard_Real& TolReached2d) const
309{
310 Standard_Boolean CompC3D = myGeomTool.CompC3D();
311
312 //cout << "MakeCurves begin" << endl;
313 if (!CompC3D) return Standard_False;
314
315 Standard_Boolean CompPC1 = myGeomTool.CompPC1();
316 Standard_Boolean CompPC2 = myGeomTool.CompPC2();
317 Standard_Real tol3d,tol2d;
318 Standard_Boolean RelativeTol;
319 myGeomTool.GetTolerances(tol3d,tol2d,RelativeTol);
320 Standard_Integer NbPntMax = myGeomTool.NbPntMax();
321
0797d9d3 322#ifdef OCCT_DEBUG
7fd59977 323 if (TopOpeBRepTool_GettraceKRO()) KRO_CURVETOOL_APPRO.Start();
324#endif
325
326//----------------------------------
327///*
328#ifdef IFV
329 char name[16];
330 char *nm = &name[0];
331 sprintf(name,"C3D_%d",++NbCalls);
332 DrawTrSurf::Set(nm, C3D);
333 sprintf(name,"PC1_%d",NbCalls);
334 DrawTrSurf::Set(nm, PC1);
335 sprintf(name,"PC2_%d",NbCalls);
336 DrawTrSurf::Set(nm, PC2);
337#endif
338//*/
339//---------------------------------------------
340
341 Standard_Integer iparmin = (Standard_Integer)parmin;
342 Standard_Integer iparmax = (Standard_Integer)parmax;
343
c5f3a425 344 Handle(Geom_BSplineCurve)& HC3D = Handle(Geom_BSplineCurve)::DownCast (C3D);
345 Handle(Geom2d_BSplineCurve)& HPC1 = Handle(Geom2d_BSplineCurve)::DownCast (PC1);
346 Handle(Geom2d_BSplineCurve)& HPC2 = Handle(Geom2d_BSplineCurve)::DownCast (PC2);
7fd59977 347
348//--------------------- IFV - "improving" initial curves
349#ifdef CurveImprovement
350 Standard_Integer nbpol = HC3D->NbPoles();
351 //cout <<"nbpol = " << nbpol << endl;
352 if(nbpol > 100) {
353 TColgp_Array1OfPnt PolC3D(1, nbpol);
354 TColgp_Array1OfPnt2d PolPC1(1, nbpol);
355 TColgp_Array1OfPnt2d PolPC2(1, nbpol);
356 TColStd_Array1OfBoolean IsValid(1, nbpol);
357 IsValid.Init(Standard_True);
358 Standard_Real tol = Max(1.e-10, 100.*tol3d*tol3d); //tol *= tol; - square distance
359 Standard_Real tl2d = tol*(tol2d*tol2d)/(tol3d*tol3d);
360 Standard_Real def = tol;
361 Standard_Real def2d = tol2d;
362 HC3D->Poles(PolC3D);
363 if(CompPC1) HPC1->Poles(PolPC1);
364 if(CompPC2) HPC2->Poles(PolPC2);
365
366 Standard_Integer ip = 1, NbPol = 1;
367 Standard_Real d, d1, d2;
368 gp_Pnt P = PolC3D(ip);
369 gp_Pnt2d P1, P2;
370 if(CompPC1) P1 = PolPC1(ip);
371 if(CompPC2) P2 = PolPC2(ip);
372
373 for(ip = 2; ip <= nbpol; ip++) {
374 if( IsValid(ip) ) {
375 d = P.SquareDistance(PolC3D(ip));
376 if(CompPC1) {d1 = P1.SquareDistance(PolPC1(ip));} else {d1 = 0.;}
377 if(CompPC2) {d2 = P2.SquareDistance(PolPC2(ip));} else {d2 = 0.;}
378 if(d > tol || d1 > tl2d || d2 > tl2d ) {
379 Standard_Real dd=0.;
380 if(ip < nbpol ) dd = P.Distance(PolC3D(ip+1));
381 if(ip < nbpol && dd < 10.*tol) {
382 gce_MakeLin mkL(P, PolC3D(ip+1));
383 if(mkL.IsDone()) {
384 gp_Lin L = mkL.Value();
385 d = L.SquareDistance(PolC3D(ip));
386 if(CompPC1) {
387 gp_Lin2d L1 = gce_MakeLin2d(P1, PolPC1(ip+1));
388 d1 = L1.SquareDistance(PolPC1(ip));
389 }
390 else d1 = 0.;
391 if(CompPC2) {
392 gp_Lin2d L2 = gce_MakeLin2d(P2, PolPC2(ip+1));
393 d2 = L2.SquareDistance(PolPC2(ip));
394 }
395 else d2 = 0.;
396
397 if(d > def || d1 > def2d || d2 > def2d ) {
398 NbPol++;
399 P = PolC3D(ip);
400 if(CompPC1) P1 = PolPC1(ip);
401 if(CompPC2) P2 = PolPC2(ip);
402 }
403 else {
404 IsValid(ip) = Standard_False;
405 }
406 }
407 else {
408 IsValid(ip+1) = Standard_False;
409 NbPol++;
410 P = PolC3D(ip);
411 if(CompPC1) P1 = PolPC1(ip);
412 if(CompPC2) P2 = PolPC2(ip);
413 }
414 }
415 else {
416 NbPol++;
417 P = PolC3D(ip);
418 if(CompPC1) P1 = PolPC1(ip);
419 if(CompPC2) P2 = PolPC2(ip);
420 }
421 }
422 else {
423 IsValid(ip) = Standard_False;
424 }
425 }
426 }
427
428 if(NbPol < 2) {IsValid(nbpol) = Standard_True; NbPol++;}
429
430 if(NbPol < nbpol) {
431#ifdef IFV
432 cout << "NbPol < nbpol : " << NbPol << " " << nbpol << endl;
433#endif
434 TColgp_Array1OfPnt Polc3d(1, NbPol);
435 TColgp_Array1OfPnt2d Polpc1(1, NbPol);
436 TColgp_Array1OfPnt2d Polpc2(1, NbPol);
437 TColStd_Array1OfReal knots(1,NbPol);
438 TColStd_Array1OfInteger mults(1,NbPol);
439 mults.Init(1); mults(1) = 2; mults(NbPol) = 2;
440 Standard_Integer count = 0;
441 for(ip = 1; ip <= nbpol; ip++) {
442 if(IsValid(ip)) {
443 count++;
444 Polc3d(count) = PolC3D(ip);
445 if(CompPC1) Polpc1(count) = PolPC1(ip);
446 if(CompPC2) Polpc2(count) = PolPC2(ip);
447 knots(count) = count;
448 }
449 }
450
451 Polc3d(NbPol) = PolC3D(nbpol);
452 if(CompPC1) Polpc1(NbPol) = PolPC1(nbpol);
453 if(CompPC2) Polpc2(NbPol) = PolPC2(nbpol);
454
455 HC3D = new Geom_BSplineCurve(Polc3d, knots, mults, 1);
456 if(CompPC1) HPC1 = new Geom2d_BSplineCurve(Polpc1, knots, mults, 1);
457 if(CompPC2) HPC2 = new Geom2d_BSplineCurve(Polpc2, knots, mults, 1);
458 iparmax = NbPol;
459
460#ifdef IFV
461 sprintf(name,"C3Dmod_%d",NbCalls);
462 nm = &name[0];
463 DrawTrSurf::Set(nm, HC3D);
464 sprintf(name,"PC1mod_%d",NbCalls);
465 nm = &name[0];
466 DrawTrSurf::Set(nm, HPC1);
467 sprintf(name,"PC2mod_%d",NbCalls);
468 nm = &name[0];
469 DrawTrSurf::Set(nm, HPC2);
470#endif
471
472 }
473 }
474//--------------- IFV - end "improving"
475#endif
476
477
478 BRepApprox_Approx Approx;
479
480 Standard_Integer degmin = 4;
481 Standard_Integer degmax = 8;
482 Approx_ParametrizationType parametrization = Approx_ChordLength;
483
484 Standard_Integer npol = HC3D->NbPoles();
485 TColgp_Array1OfPnt Polc3d(1, npol);
486 TColStd_Array1OfReal par(1, npol);
487 HC3D->Poles(Polc3d);
488 gp_Pnt P = Polc3d(1);
489
490 Standard_Boolean IsBad = Standard_False;
491 Standard_Integer ip;
492 Standard_Real ksi = 0., kc, kcprev = 0.;
493 Standard_Real dist;
494 par(1) = 0.;
495 for(ip = 2; ip <= npol; ip++) {
496 dist = P.Distance(Polc3d(ip));
497
498 if(dist < Precision::Confusion()) {
499 IsBad = Standard_True;
500 break;
501 }
502
503 par(ip) = par(ip-1) + dist;
504 P = Polc3d(ip);
505 }
506
507 if(!IsBad) {
508
509 TColStd_Array1OfReal Curvature(1, npol);
510
511 if(npol > 3) {
512 Standard_Integer ic = 1;
513 for(ip = 2; ip <= npol-1; ip += npol-3) {
514 gp_Vec v1(Polc3d(ip-1),Polc3d(ip));
515 gp_Vec v2(Polc3d(ip),Polc3d(ip+1));
516 if(!v1.IsParallel(v2, 1.e-4)) {
517 gce_MakeCirc mc(Polc3d(ip-1),Polc3d(ip),Polc3d(ip+1));
518 gp_Circ cir = mc.Value();
519 kc = 1./cir.Radius();
520 }
521 else kc = 0.;
522
523 Curvature(ic) = kc;
524 ic = npol;
525 }
526 }
527 else if(npol == 3) {
528 ip = 2;
529 gp_Vec v1(Polc3d(ip-1),Polc3d(ip));
530 gp_Vec v2(Polc3d(ip),Polc3d(ip+1));
531 if(!v1.IsParallel(v2, 1.e-4)) {
532 gce_MakeCirc mc(Polc3d(ip-1),Polc3d(ip),Polc3d(ip+1));
533 gp_Circ cir = mc.Value();
534 kc = 1./cir.Radius();
535 }
536 else kc = 0.;
537 Curvature(1) = Curvature(npol) = kc;
538 }
539 else {
540 Curvature(1) = Curvature(npol) = 0.;
541 }
542
543 ip = 1;
544 Standard_Real dt = par(ip+1) - par(ip);
545 Standard_Real dx = (Polc3d(ip+1).X() - Polc3d(ip).X())/dt,
546 dy = (Polc3d(ip+1).Y() - Polc3d(ip).Y())/dt,
547 dz = (Polc3d(ip+1).Z() - Polc3d(ip).Z())/dt;
548 Standard_Real dx1, dy1, dz1, d2x, d2y, d2z, d2t;
549
550 for(ip = 2; ip <= npol-1; ip++) {
551 dt = par(ip+1) - par(ip);
552 dx1 = (Polc3d(ip+1).X() - Polc3d(ip).X())/dt;
553 dy1 = (Polc3d(ip+1).Y() - Polc3d(ip).Y())/dt,
554 dz1 = (Polc3d(ip+1).Z() - Polc3d(ip).Z())/dt;
555 d2t = 2./(par(ip+1) - par(ip-1));
556 d2x = d2t*(dx1 - dx);
557 d2y = d2t*(dy1 - dy);
558 d2z = d2t*(dz1 - dz);
559 Curvature(ip) = Sqrt(d2x*d2x + d2y*d2y + d2z*d2z);
560 dx = dx1; dy = dy1; dz = dz1;
561 }
562
563 Standard_Real crit = 1000.; // empirical criterion !!!
564
565 dt = par(2) - par(1);
566 kcprev = (Curvature(2) - Curvature(1))/dt;
567 for(ip = 2; ip <= npol-1; ip++) {
568 dt = par(ip+1) - par(ip);
569 kc = (Curvature(ip+1) - Curvature(ip))/dt;
570 ksi = ksi + Abs(kc - kcprev);
571 if(ksi > crit) {IsBad = Standard_True;break;}
572 kc = kcprev;
573 }
574
575 }
576 //cout << NbCalls << " ksi = " << ksi << endl;
577 //cout << "IsBad = " << IsBad << endl;
578
579 if(IsBad){
580 Standard_Real tt = Min(10.*tol3d,Precision::Approximation());
581 tol2d = tt * tol2d / tol3d;
582 tol3d = tt;
583 NbPntMax = 40;
584 degmax = 4;
585 parametrization = Approx_Centripetal;
586 }
587
588 Standard_Integer nitmax = 0; // use projection only
589 Standard_Boolean withtangency = Standard_True;
590
591 Standard_Boolean compminmaxUV = Standard_True;
592 BRepAdaptor_Surface BAS1(TopoDS::Face(S1),compminmaxUV);
593 BRepAdaptor_Surface BAS2(TopoDS::Face(S2),compminmaxUV);
594
595
596 Handle(BRepApprox_ApproxLine) AL;
597 AL = new BRepApprox_ApproxLine(HC3D,HPC1,HPC2);
598
599 Approx.SetParameters(tol3d,tol2d,RelativeTol,degmin,degmax,nitmax,NbPntMax,withtangency,
600 parametrization);
601
602 if (CompC3D && CompPC1 && BAS1.GetType() == GeomAbs_Plane) {
603 //-- The curve X,Y,Z and U2,V2 is approximated
604 Approx.Perform(BAS1,BAS2,AL,CompC3D,Standard_False,CompPC2,iparmin,iparmax);
605 }
606
607 else if(CompC3D && CompPC2 && BAS2.GetType() == GeomAbs_Plane) {
608 //-- The curve X,Y,Z and U1,V1 is approximated
609 Approx.Perform(BAS1,BAS2,AL,CompC3D,CompPC1,Standard_False,iparmin,iparmax);
610 }
611
612 else {
613 Approx.Perform(BAS1,BAS2,AL,CompC3D,CompPC1,CompPC2,iparmin,iparmax);
614 }
615
616 // MSV Nov 9, 2001: do not raise exception in the case of failure,
617 // but return status
618
619 Standard_Boolean done = Approx.IsDone();
620 done = done && CheckApproxResults(Approx);
621
622 if (done) {
623 if (CompC3D && CompPC1 && BAS1.GetType() == GeomAbs_Plane) {
624 C3Dnew = ::MakeCurve3DfromWLineApprox(Approx,1);
625 PC1new = ::MakeCurve2DfromWLineApproxAndPlane(Approx,BAS1.Plane());
626 if (CompPC2) PC2new = ::MakeCurve2DfromWLineApprox(Approx,2);
627 }
628 else if(CompC3D && CompPC2 && BAS2.GetType() == GeomAbs_Plane) {
629 C3Dnew = ::MakeCurve3DfromWLineApprox(Approx,1);
630 if (CompPC1) PC1new = ::MakeCurve2DfromWLineApprox(Approx,2);
631 PC2new = ::MakeCurve2DfromWLineApproxAndPlane(Approx,BAS2.Plane());
632 }
633 else {
634 if (CompC3D) C3Dnew = ::MakeCurve3DfromWLineApprox(Approx,1);
635 if (CompPC1) PC1new = ::MakeCurve2DfromWLineApprox(Approx,2);
636 if (CompPC2) {
637 Standard_Integer i32 = (CompPC1) ? 3 : 2;
638 PC2new = ::MakeCurve2DfromWLineApprox(Approx,i32);
639 }
640 }
641
642 // check the pcurves relatively the faces bounds
643 if (CompPC1)
644 done = done && CheckPCurve(PC1new,TopoDS::Face(S1));
645 if (CompPC2)
646 done = done && CheckPCurve(PC2new,TopoDS::Face(S2));
647 }
648
649 if (!done) {
650 if (CompC3D) C3Dnew.Nullify();
651 if (CompPC1) PC1new.Nullify();
652 if (CompPC2) PC2new.Nullify();
653 return Standard_False;
654 }
655
656#ifdef IFV
657 sprintf(name,"C3Dnew_%d", NbCalls);
658 nm = &name[0];
659 DrawTrSurf::Set(nm, C3Dnew);
660 if (CompPC1) {
661 sprintf(name,"PC1new_%d", NbCalls);
662 nm = &name[0];
663 DrawTrSurf::Set(nm, PC1new);
664 }
665 if (CompPC2) {
666 sprintf(name,"PC2new_%d", NbCalls);
667 nm = &name[0];
668 DrawTrSurf::Set(nm, PC2new);
669 }
670
671#endif
672
673 TolReached3d = Approx.TolReached3d();
674 TolReached2d = Approx.TolReached2d();
675#ifdef IFV
676 cout << NbCalls << " : Tol = " << TolReached3d << " " << TolReached2d << endl;
677#endif
678
679 Standard_Boolean bf, bl;
680
c5f3a425 681 Handle(Geom_BSplineCurve) Curve (Handle(Geom_BSplineCurve)::DownCast(C3Dnew));
7fd59977 682 if(!Curve.IsNull()) {
683 GeomLib_CheckBSplineCurve cbsc(Curve, 1.e-7, 0.1);
684 cbsc.NeedTangentFix(bf, bl);
685
0797d9d3 686#ifdef OCCT_DEBUG
7fd59977 687 if (TopOpeBRepTool_GettraceCHKBSPL()) {
688 if(bf || bl) {
689 cout<<"Problem orientation GeomLib_CheckBSplineCurve : First = "<<bf;
690 cout<<" Last = "<<bl<<endl;
691 }
692 }
693#endif
694 cbsc.FixTangent(bf, bl);
695 }
696
c5f3a425 697 Handle(Geom2d_BSplineCurve) Curve2df (Handle(Geom2d_BSplineCurve)::DownCast(PC1new));
7fd59977 698 if(!Curve2df.IsNull()) {
699 GeomLib_Check2dBSplineCurve cbsc2df(Curve2df, 1.e-7, 0.1);
700 cbsc2df.NeedTangentFix(bf, bl);
0797d9d3 701#ifdef OCCT_DEBUG
7fd59977 702 if (TopOpeBRepTool_GettraceCHKBSPL()) {
703 if(bf || bl) {
704 cout<<"Problem orientation GeomLib_CheckBSplineCurve : First = "<<bf;
705 cout<<" Last = "<<bl<<endl;
706 }
707 }
708#endif
709 cbsc2df.FixTangent(bf, bl);
710 }
711
c5f3a425 712 Handle(Geom2d_BSplineCurve) Curve2ds (Handle(Geom2d_BSplineCurve)::DownCast(PC2new));
7fd59977 713 if(!Curve2ds.IsNull()) {
714 GeomLib_Check2dBSplineCurve cbsc2ds(Curve2ds, 1.e-7, 0.1);
715 cbsc2ds.NeedTangentFix(bf, bl);
0797d9d3 716#ifdef OCCT_DEBUG
7fd59977 717 if (TopOpeBRepTool_GettraceCHKBSPL()) {
718 if(bf || bl) {
719 cout<<"Problem orientation GeomLib_CheckBSplineCurve : First = "<<bf;
720 cout<<" Last = "<<bl<<endl;
721 }
722 }
723#endif
724 cbsc2ds.FixTangent(bf, bl);
725 }
726
0797d9d3 727#ifdef OCCT_DEBUG
7fd59977 728 if (TopOpeBRepTool_GettraceKRO()) KRO_CURVETOOL_APPRO.Stop();
729#endif
730// cout << "MakeCurves end" << endl;
731
732 return Standard_True;
733}
734
735
736//=======================================================================
737//function : MakeBSpline1fromPnt
738//purpose :
739//=======================================================================
740
741Handle(Geom_Curve) TopOpeBRepTool_CurveTool::MakeBSpline1fromPnt
742(const TColgp_Array1OfPnt& Points)
743{
744 Standard_Integer Degree = 1;
745
746 Standard_Integer i,nbpoints = Points.Length();
747 Standard_Integer nbknots = nbpoints - Degree +1;
748
749 // First compute the parameters
750 // Standard_Real length = 0.;
751 // TColStd_Array1OfReal parameters(1,nbpoints);
752 // for (i = 1; i < nbpoints; i++) {
753 // parameters(i) = length;
754 // Standard_Real dist = Points(i).Distance(Points(i+1));
755 // length += dist;
756 // }
757 // parameters(nbpoints) = length;
758
759 // knots and multiplicities
760 TColStd_Array1OfReal knots(1,nbknots);
761 TColStd_Array1OfInteger mults(1,nbknots);
762 mults.Init(1);
763 mults(1) = mults(nbknots) = Degree + 1;
764
765 // knots(1) = 0;
766 // for (i=2;i<nbknots;i++) knots(i) = (parameters(i) + parameters(i+1)) /2.;
767 // knots(nbknots) = length;
768
769 // take point index as parameter : JYL 01/AUG/94
770 for (i = 1; i <= nbknots; i++) knots(i) = (Standard_Real)i;
771 Handle(Geom_Curve) C = new Geom_BSplineCurve(Points,knots,mults,Degree);
772 return C;
773}
774
775
776//=======================================================================
777//function : MakeBSpline1fromPnt2d
778//purpose :
779//=======================================================================
780
781Handle(Geom2d_Curve) TopOpeBRepTool_CurveTool::MakeBSpline1fromPnt2d
782(const TColgp_Array1OfPnt2d& Points)
783{
784 Standard_Integer Degree = 1;
785
786 Standard_Integer i,nbpoints = Points.Length();
787 Standard_Integer nbknots = nbpoints - Degree +1;
788
789 // First compute the parameters
790 // Standard_Real length = 0;
791 // TColStd_Array1OfReal parameters(1,nbpoints);
792 // for (i = 1; i < nbpoints; i++) {
793 // parameters(i) = length;
794 // Standard_Real dist = Points(i).Distance(Points(i+1));
795 // length += dist;
796 // }
797 // parameters(nbpoints) = length;
798
799 // knots and multiplicities
800 TColStd_Array1OfReal knots(1,nbknots);
801 TColStd_Array1OfInteger mults(1,nbknots);
802 mults.Init(1);
803 mults(1) = mults(nbknots) = Degree + 1;
804
805 // knots(1) = 0;
806 // for (i=2;i<nbknots;i++) knots(i) = (parameters(i) + parameters(i+1)) /2.;
807 // knots(nbknots) = length;
808
809 // take point index as parameter : JYL 01/AUG/94
810 for (i = 1; i <= nbknots; i++) knots(i) = (Standard_Real)i;
811 Handle(Geom2d_Curve) C = new Geom2d_BSplineCurve(Points,knots,mults,Degree);
812 return C;
813}
814
815//=======================================================================
816//function : IsProjectable
817//purpose :
818//=======================================================================
819
820Standard_Boolean TopOpeBRepTool_CurveTool::IsProjectable
821(const TopoDS_Shape& S, const Handle(Geom_Curve)& C3D)
822{
823 const TopoDS_Face& F = TopoDS::Face(S);
824 Standard_Boolean compminmaxUV = Standard_False;
825 BRepAdaptor_Surface BAS(F,compminmaxUV);
826 GeomAbs_SurfaceType suty = BAS.GetType();
827 GeomAdaptor_Curve GAC(C3D);
828 GeomAbs_CurveType cuty = GAC.GetType();
829
830 // --------
831 // avoid projection of 3d curve on surface in case
832 // of a quadric (ellipse,hyperbola,parabola) on a cone.
833 // Projection fails when the curve in not fully inside the UV domain
834 // of the cone : only part of 2d curve is built.
835 // NYI : projection of quadric on cone (crossing cone domain)
836 // --------
837
838 Standard_Boolean projectable = Standard_True;
839 if ( suty == GeomAbs_Cone ) {
840 if( (cuty == GeomAbs_Ellipse) ||
841 ( cuty == GeomAbs_Hyperbola) ||
842 ( cuty == GeomAbs_Parabola) ) {
843 projectable = Standard_False;
844 }
845 }
846 else if ( suty == GeomAbs_Cylinder ) {
847 if (cuty == GeomAbs_Ellipse) {
848 projectable = Standard_False;
849 }
850 }
851 else if ( suty == GeomAbs_Sphere ) {
852 if (cuty == GeomAbs_Circle) {
853 projectable = Standard_False;
854 }
855 }
856 else if ( suty == GeomAbs_Torus ) {
857 if (cuty == GeomAbs_Circle) {
858 projectable = Standard_False;
859 }
860 }
861
0797d9d3 862#ifdef OCCT_DEBUG
7fd59977 863 if (TopOpeBRepTool_GettracePCURV()) {
864 cout<<"--- IsProjectable : ";
865 if (projectable) cout<<"projectable"<<endl;
866 else cout<<"NOT projectable"<<endl;
867 }
868#endif
869
870 return projectable;
871}
872
873
874//=======================================================================
875//function : MakePCurveOnFace
876//purpose :
877//=======================================================================
878
879Handle(Geom2d_Curve) TopOpeBRepTool_CurveTool::MakePCurveOnFace
880(const TopoDS_Shape& S,
881 const Handle(Geom_Curve)& C3D,
882 Standard_Real& TolReached2d,
883 const Standard_Real first, const Standard_Real last)
884
885{
886 Standard_Boolean trim = Standard_False;
887 if (first < last)
888 trim = Standard_True;
889
890 const TopoDS_Face& F = TopoDS::Face(S);
891 Standard_Boolean compminmaxUV = Standard_False;
892 BRepAdaptor_Surface BAS(F,compminmaxUV);
893 GeomAdaptor_Curve GAC;
894 if (trim) GAC.Load(C3D,first,last);
895 else GAC.Load(C3D);
896 Handle(BRepAdaptor_HSurface) BAHS = new BRepAdaptor_HSurface(BAS);
897 Handle(GeomAdaptor_HCurve) BAHC = new GeomAdaptor_HCurve(GAC);
898 ProjLib_ProjectedCurve projcurv(BAHS,BAHC);
899 Handle(Geom2d_Curve) C2D = ::MakePCurve(projcurv);
900 TolReached2d = projcurv.GetTolerance();
901
902 Standard_Real UMin, UMax, VMin, VMax;
903 BRepTools::UVBounds(F,UMin, UMax, VMin, VMax);
904
905 Standard_Real f = GAC.FirstParameter();
906 Standard_Real l = GAC.LastParameter();
907 Standard_Real t =(f+l)*.5;
908 gp_Pnt2d pC2D; C2D->D0(t,pC2D);
909 Standard_Real u2 = pC2D.X();
910 Standard_Real v2 = pC2D.Y();
911
912 if (BAS.GetType() == GeomAbs_Sphere) {
913 // MSV: consider quasiperiodic shift of pcurve
914 Standard_Real VFirst = BAS.FirstVParameter();
915 Standard_Real VLast = BAS.LastVParameter();
916 Standard_Boolean mincond = v2 < VFirst;
917 Standard_Boolean maxcond = v2 > VLast;
918 if (mincond || maxcond) {
919 Handle(Geom2d_Curve) PCT = Handle(Geom2d_Curve)::DownCast(C2D->Copy());
920 // make mirror relative to the isoline of apex -PI/2 or PI/2
921 gp_Trsf2d aTrsf;
c6541a0c
D
922 gp_Pnt2d po(0,-M_PI/2);
923 if (maxcond) po.SetY(M_PI/2);
7fd59977 924 aTrsf.SetMirror(gp_Ax2d(po, gp_Dir2d(1,0)));
925 PCT->Transform(aTrsf);
926 // add translation along U direction on PI
c6541a0c 927 gp_Vec2d vec(M_PI,0);
7fd59977 928 Standard_Real UFirst = BAS.FirstUParameter();
c6541a0c 929 if (u2-UFirst-M_PI > -1e-7) vec.Reverse();
7fd59977 930 PCT->Translate(vec);
931 C2D = PCT;
932 // recompute the test point
933 C2D->D0(t,pC2D);
934 u2 = pC2D.X();
935 v2 = pC2D.Y();
936 }
937 }
938
939 Standard_Real du = 0.;
940 if (BAHS->IsUPeriodic()) {
941 //modified by NIZHNY-MZV Thu Mar 30 10:03:15 2000
942 Standard_Boolean mincond = (UMin - u2 > 1e-7) ? Standard_True : Standard_False;
943 Standard_Boolean maxcond = (u2 - UMax > 1e-7) ? Standard_True : Standard_False;
944 Standard_Boolean decalu = mincond || maxcond;
945 if (decalu) du = ( mincond ) ? BAHS->UPeriod() : -BAHS->UPeriod();
946 //Standard_Boolean decalu = ( u2 < UMin || u2 > UMax);
947 //if (decalu) du = ( u2 < UMin ) ? BAHS->UPeriod() : -BAHS->UPeriod();
948 }
949 Standard_Real dv = 0.;
950 if (BAHS->IsVPeriodic()) {
951 //modified by NIZHNY-MZV Thu Mar 30 10:06:24 2000
952 Standard_Boolean mincond = (VMin - v2 > 1e-7) ? Standard_True : Standard_False;
953 Standard_Boolean maxcond = (v2 - VMax > 1e-7) ? Standard_True : Standard_False;
954 Standard_Boolean decalv = mincond || maxcond;
955 if (decalv) dv = ( mincond ) ? BAHS->VPeriod() : -BAHS->VPeriod();
956 //Standard_Boolean decalv = ( v2 < VMin || v2 > VMax);
957 //if (decalv) dv = ( v2 < VMin ) ? BAHS->VPeriod() : -BAHS->VPeriod();
958 }
959
960 if ( du != 0. || dv != 0.) {
961 Handle(Geom2d_Curve) PCT = Handle(Geom2d_Curve)::DownCast(C2D->Copy());
962 PCT->Translate(gp_Vec2d(du,dv));
963 C2D = PCT;
964 }
965
966 return C2D;
967}