1 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 // This file is part of Open CASCADE Technology software library.
5 // This library is free software; you can redistribute it and/or modify it under
6 // the terms of the GNU Lesser General Public License version 2.1 as published
7 // by the Free Software Foundation, with special exception defined in the file
8 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
9 // distribution for complete text of the license and disclaimer of any warranty.
11 // Alternatively, this file may be used under the terms of Open CASCADE
12 // commercial license or contractual agreement.
14 // pdn 04.12.98 Add method using Adaptor_Curve
15 //:j8 abv 10.12.98 TR10 r0501_db.stp #9423
16 //pdn 25.12.98 private method ProjectAct
18 //:s5 abv 22.04.99 Adding debug printouts in catch {} blocks
19 // abv 14.05.99 S4174: Adding method for exact computing of the boundary box
20 // gka 21.06.99 S4208: adding method NextProject(Adaptor_Curve)
21 // msv 30.05.00 correct IsPlanar for a conic curve
23 #include <Adaptor3d_Curve.hxx>
24 #include <Bnd_Box2d.hxx>
26 #include <Extrema_ExtPC.hxx>
27 #include <Extrema_LocateExtPC.hxx>
28 #include <Geom2d_BezierCurve.hxx>
29 #include <Geom2d_BoundedCurve.hxx>
30 #include <Geom2d_BSplineCurve.hxx>
31 #include <Geom2d_Conic.hxx>
32 #include <Geom2d_Curve.hxx>
33 #include <Geom2d_Line.hxx>
34 #include <Geom2d_OffsetCurve.hxx>
35 #include <Geom2d_TrimmedCurve.hxx>
36 #include <Geom2dAdaptor_Curve.hxx>
37 #include <Geom2dInt_Geom2dCurveTool.hxx>
38 #include <Geom_BezierCurve.hxx>
39 #include <Geom_BoundedCurve.hxx>
40 #include <Geom_BSplineCurve.hxx>
41 #include <Geom_Circle.hxx>
42 #include <Geom_Conic.hxx>
43 #include <Geom_Curve.hxx>
44 #include <Geom_Line.hxx>
45 #include <Geom_OffsetCurve.hxx>
46 #include <Geom_TrimmedCurve.hxx>
47 #include <GeomAdaptor_Curve.hxx>
50 #include <Precision.hxx>
51 #include <ShapeAnalysis.hxx>
52 #include <ShapeAnalysis_Curve.hxx>
53 #include <ShapeExtend_ComplexCurve.hxx>
54 #include <Standard_ErrorHandler.hxx>
55 #include <Standard_Failure.hxx>
56 #include <TColgp_SequenceOfPnt.hxx>
58 //=======================================================================
59 //function : ProjectOnSegments
61 //=======================================================================
62 static void ProjectOnSegments (const Adaptor3d_Curve& AC, const gp_Pnt& P3D,
63 const Standard_Integer nbseg,
64 Standard_Real& uMin, Standard_Real& uMax,
65 Standard_Real& distmin, gp_Pnt& proj, Standard_Real& param)
67 // On considere <nbseg> points sur [uMin,uMax]
68 // Quel est le plus proche. Et quel est le nouvel intervalle
69 // (il ne peut pas deborder l ancien)
70 Standard_Real u, dist2, delta = (nbseg == 0)? 0 : (uMax-uMin)/nbseg; //szv#4:S4163:12Mar99 anti-exception
71 Standard_Real distmin2 = distmin * distmin;
72 Standard_Boolean aHasChanged = Standard_False;
73 for (Standard_Integer i = 0; i <= nbseg; i ++) {
74 u = uMin + (delta * i);
75 gp_Pnt PU = AC.Value (u);
76 dist2 = PU.SquareDistance (P3D);
77 if (dist2 < distmin2) { distmin2 = dist2; proj = PU; param = u; aHasChanged = Standard_True; }
80 distmin = Sqrt (distmin2);
82 cout<<"ShapeAnalysis_Geom:Project, param="<<param<<" -> distmin="<<distmin<<endl;
85 uMax = Min (uMax, param+delta);
86 uMin = Max (uMin, param-delta);
90 //=======================================================================
93 //=======================================================================
95 Standard_Real ShapeAnalysis_Curve::Project(const Handle(Geom_Curve)& C3D,
97 const Standard_Real preci,
100 const Standard_Boolean AdjustToEnds) const
102 Standard_Real uMin = C3D->FirstParameter();
103 Standard_Real uMax = C3D->LastParameter();
104 if (uMin < uMax) return Project (C3D,P3D,preci,proj,param,uMin,uMax,AdjustToEnds);
105 else return Project (C3D,P3D,preci,proj,param,uMax,uMin,AdjustToEnds);
108 //=======================================================================
111 //=======================================================================
113 Standard_Real ShapeAnalysis_Curve::Project(const Handle(Geom_Curve)& C3D,
115 const Standard_Real preci,
117 Standard_Real& param,
118 const Standard_Real cf,
119 const Standard_Real cl,
120 const Standard_Boolean AdjustToEnds) const
122 Standard_Real distmin;
123 Standard_Real uMin = (cf < cl ? cf : cl);
124 Standard_Real uMax = (cf < cl ? cl : cf);
126 GeomAdaptor_Curve GAC(C3D, uMin, uMax);
127 if (C3D->IsKind(STANDARD_TYPE(Geom_BoundedCurve))) {
128 Standard_Real prec = ( AdjustToEnds ? preci : Precision::Confusion() ); //:j8 abv 10 Dec 98: tr10_r0501_db.stp #9423: protection against densing of points near one end
129 gp_Pnt LowBound = GAC.Value(uMin);
130 gp_Pnt HigBound = GAC.Value(uMax);
131 distmin = LowBound.Distance(P3D);
132 if (distmin <= prec) {
137 distmin = HigBound.Distance(P3D);
138 if (distmin <= prec) {
145 if (!C3D->IsClosed()) {
146 //modified by rln on 16/12/97 after CSR# PRO11641 entity 20767
147 //the VIso was not closed (according to C3D->IsClosed()) while it "almost"
148 //was (the distance between ends of the curve was a little bit more than
149 //Precision::Confusion())
150 //in that case value 0.1 was too much and this method returned not correct parameter
153 // modified by pdn on 01.07.98 after BUC60195 entity 1952 (Min() added)
154 Standard_Real delta = Min (GAC.Resolution (preci), (uMax - uMin) * 0.1);
157 GAC.Load(C3D, uMin, uMax);
160 return ProjectAct(GAC, P3D, preci, proj, param);
163 //=======================================================================
166 //=======================================================================
168 Standard_Real ShapeAnalysis_Curve::Project(const Adaptor3d_Curve& C3D,
170 const Standard_Real preci,
172 Standard_Real& param,
173 const Standard_Boolean AdjustToEnds) const
177 Standard_Real uMin = C3D.FirstParameter();
178 Standard_Real uMax = C3D.LastParameter();
180 if (Precision::IsInfinite(uMin) && Precision::IsInfinite(uMax))
181 return ProjectAct(C3D, P3D, preci, proj, param);
183 Standard_Real distmin_L = Precision::Infinite(), distmin_H = Precision::Infinite();
184 Standard_Real prec = ( AdjustToEnds ? preci : Precision::Confusion() ); //:j8 abv 10 Dec 98: tr10_r0501_db.stp #9423: protection against densing of points near one end
185 gp_Pnt LowBound = C3D.Value(uMin);
186 gp_Pnt HigBound = C3D.Value(uMax);
187 distmin_L = LowBound.Distance(P3D);
188 distmin_H = HigBound.Distance(P3D);
190 if (distmin_L <= prec) {
196 if (distmin_H <= prec) {
202 Standard_Real distProj = ProjectAct(C3D, P3D, preci, proj, param);
203 if( distProj < distmin_L + Precision::Confusion() && distProj < distmin_H + Precision::Confusion())
206 if( distmin_L < distmin_H)
218 //=======================================================================
219 //function : ProjectAct
221 //=======================================================================
223 Standard_Real ShapeAnalysis_Curve::ProjectAct(const Adaptor3d_Curve& C3D,
225 const Standard_Real preci,
227 Standard_Real& param) const
230 Standard_Boolean OK = Standard_False;
234 Extrema_ExtPC myExtPC(P3D,C3D);
235 Standard_Real dist2Min = RealLast() , dist2;
236 Standard_Integer index = 0;
237 if ( myExtPC.IsDone() && ( myExtPC.NbExt() > 0) )
239 for ( Standard_Integer i = 1; i <= myExtPC.NbExt(); i++)
241 if (!myExtPC.IsMin(i))
244 dist2 = myExtPC.SquareDistance(i);
245 if ( dist2 < dist2Min)
247 dist2Min = dist2; index = i;
253 param = (myExtPC.Point(index)).Parameter();
254 proj = (myExtPC.Point(index)).Value();
259 catch(Standard_Failure const& anException) {
262 cout << "\nWarning: ShapeAnalysis_Curve::ProjectAct(): Exception in Extrema_ExtPC: ";
263 anException.Print(cout); cout << endl;
269 //szv#4:S4163:12Mar99 moved
270 Standard_Real uMin = C3D.FirstParameter(), uMax = C3D.LastParameter();
271 Standard_Boolean closed = Standard_False; // si on franchit les bornes ...
272 Standard_Real distmin = Precision::Infinite(), valclosed = 0.;
273 Standard_Real aModParam = param;
274 Standard_Real aModMin = distmin;
276 // PTV 29.05.2002 remember the old solution, cause it could be better
277 Standard_Real anOldParam =0.;
278 Standard_Boolean IsHaveOldSol = Standard_False;
281 IsHaveOldSol = Standard_True;
284 distmin = proj.Distance (P3D);
286 if (distmin > preci) OK = Standard_False;
287 // Cas TrimmedCurve a cheval. Voir la courbe de base.
288 // Si fermee, passer une periode
289 if (C3D.IsClosed()) {
290 closed = Standard_True;
291 valclosed = uMax - uMin; //szv#4:S4163:12Mar99 optimized
296 // BUG NICOLAS - Si le point est sur la courbe 0 Solutions
297 // Cela fonctionne avec ElCLib
299 // D une facon generale, on essaie de TOUJOURS retourner un resultat
300 // MEME PAS BIEN BON. L appelant pourra decider alors quoi faire
303 switch(C3D.GetType()) {
306 const gp_Circ& aCirc = C3D.Circle();
307 proj = aCirc.Position().Location();
308 if(aCirc.Radius() <= gp::Resolution() ||
309 P3D.SquareDistance(proj) <= gp::Resolution() ) {
310 param = C3D.FirstParameter();
311 proj = proj.XYZ() + aCirc.XAxis().Direction().XYZ() * aCirc.Radius();
314 param = ElCLib::Parameter(aCirc, P3D);
315 proj = ElCLib::Value(param, aCirc);
317 closed = Standard_True;
321 case GeomAbs_Hyperbola:
323 param = ElCLib::Parameter(C3D.Hyperbola(), P3D);
324 proj = ElCLib::Value(param, C3D.Hyperbola());
327 case GeomAbs_Parabola:
329 param = ElCLib::Parameter(C3D.Parabola(), P3D);
330 proj = ElCLib::Value(param, C3D.Parabola());
335 param = ElCLib::Parameter(C3D.Line(), P3D);
336 proj = ElCLib::Value(param, C3D.Line());
339 case GeomAbs_Ellipse:
341 param = ElCLib::Parameter(C3D.Ellipse(), P3D);
342 proj = ElCLib::Value(param, C3D.Ellipse());
343 closed = Standard_True;
350 // on ne va quand meme pas se laisser abattre ... ???
351 // on tente ceci : 21 points sur la courbe, quel est le plus proche
352 distmin = Precision::Infinite();
353 ProjectOnSegments (C3D,P3D,25,uMin,uMax,distmin,proj,param);
354 if (distmin <= preci)
356 Extrema_LocateExtPC aProjector (P3D, C3D, param/*U0*/, uMin, uMax, preci/*TolU*/);
357 if (aProjector.IsDone())
359 param = aProjector.Point().Parameter();
360 proj = aProjector.Point().Value();
361 Standard_Real aDistNewton = P3D.Distance(proj);
362 if (aDistNewton < aModMin)
365 // cout <<"newton failed"<<endl;
366 ProjectOnSegments (C3D,P3D,40,uMin,uMax,distmin,proj,param);
367 if (distmin <= preci) return distmin;
368 ProjectOnSegments (C3D,P3D,20,uMin,uMax,distmin,proj,param);
369 if (distmin <= preci) return distmin;
370 ProjectOnSegments (C3D,P3D,25,uMin,uMax,distmin,proj,param);
371 if (distmin <= preci) return distmin;
372 ProjectOnSegments (C3D,P3D,40,uMin,uMax,distmin,proj,param);
373 if (distmin <= preci) return distmin;
374 // soyons raisonnable ...
375 if(distmin > aModMin) {
385 // p = PPOC.LowerDistanceParameter(); cf try
386 if ( closed && ( param < uMin || param > uMax ) )
387 param += ShapeAnalysis::AdjustByPeriod ( param, 0.5 * ( uMin + uMax ), valclosed );
390 // PTV 29.05.2002 Compare old solution and new;
391 Standard_Real adist1, adist2;
392 adist1 = anOldProj.SquareDistance(P3D);
393 adist2 = proj.SquareDistance (P3D);
394 if (adist1 < adist2) {
399 return proj.Distance (P3D);
403 //=======================================================================
404 //function : NextProject
405 //purpose : Newton algo for projecting point on curve (S4030)
406 //=======================================================================
408 Standard_Real ShapeAnalysis_Curve::NextProject(const Standard_Real paramPrev,
409 const Handle(Geom_Curve)& C3D,
411 const Standard_Real preci,
413 Standard_Real& param,
414 const Standard_Real cf,
415 const Standard_Real cl,
416 const Standard_Boolean AdjustToEnds) const
418 Standard_Real uMin = (cf < cl ? cf : cl);
419 Standard_Real uMax = (cf < cl ? cl : cf);
420 Standard_Real distmin = Precision::Infinite();
421 GeomAdaptor_Curve GAC(C3D, uMin, uMax);
422 if (C3D->IsKind(STANDARD_TYPE(Geom_BoundedCurve))) {
423 Standard_Real prec = ( AdjustToEnds ? preci : Precision::Confusion() ); //:j8 abv 10 Dec 98: tr10_r0501_db.stp #9423: protection against densing of points near one end
424 gp_Pnt LowBound = GAC.Value(uMin);
425 gp_Pnt HigBound = GAC.Value(uMax);
426 distmin = LowBound.Distance(P3D);
427 if (distmin <= prec) {
432 distmin = HigBound.Distance(P3D);
433 if (distmin <= prec) {
440 if (!C3D->IsClosed()) {
441 //modified by rln on 16/12/97 after CSR# PRO11641 entity 20767
442 //the VIso was not closed (according to C3D->IsClosed()) while it "almost"
443 //was (the distance between ends of the curve was a little bit more than
444 //Precision::Confusion())
445 //in that case value 0.1 was too much and this method returned not correct parameter
448 // modified by pdn on 01.07.98 after BUC60195 entity 1952 (Min() added)
449 Standard_Real delta = Min (GAC.Resolution (preci), (uMax - uMin) * 0.1);
452 GAC.Load(C3D, uMin, uMax);
454 return NextProject ( paramPrev, GAC, P3D, preci, proj, param );
457 //=======================================================================
458 //function : NextProject
460 //=======================================================================
462 Standard_Real ShapeAnalysis_Curve::NextProject(const Standard_Real paramPrev,
463 const Adaptor3d_Curve& C3D,
465 const Standard_Real preci,
467 Standard_Real& param) const
469 Standard_Real uMin = C3D.FirstParameter();
470 Standard_Real uMax = C3D.LastParameter();
472 Extrema_LocateExtPC aProjector (P3D, C3D, paramPrev/*U0*/, uMin, uMax, preci/*TolU*/);
473 if (aProjector.IsDone()){
474 param = aProjector.Point().Parameter();
475 proj = aProjector.Point().Value();
476 return P3D.Distance(proj);
478 return Project(C3D, P3D, preci, proj, param, Standard_False);
481 //=======================================================================
482 //function : AdjustParameters
483 //purpose : Copied from StepToTopoDS_GeometricTuul::UpdateParam3d (Aug 2001)
484 //=======================================================================
486 Standard_Boolean ShapeAnalysis_Curve::ValidateRange (const Handle(Geom_Curve)& theCurve,
487 Standard_Real& First, Standard_Real& Last,
488 const Standard_Real preci) const
490 // First et/ou Last peuvent etre en dehors des bornes naturelles de la courbe.
491 // On donnera alors la valeur en bout a First et/ou Last
493 Standard_Real cf = theCurve->FirstParameter();
494 Standard_Real cl = theCurve->LastParameter();
495 // Standard_Real preci = BRepAPI::Precision();
497 if (theCurve->IsKind(STANDARD_TYPE(Geom_BoundedCurve)) && !theCurve->IsClosed()) {
500 cout << "Update Edge First Parameter to Curve First Parameter" << endl;
504 else if (First > cl) {
506 cout << "Update Edge First Parameter to Curve Last Parameter" << endl;
512 cout << "Update Edge Last Parameter to Curve First Parameter" << endl;
516 else if (Last > cl) {
518 cout << "Update Edge Last Parameter to Curve Last Parameter" << endl;
524 // 15.11.2002 PTV OCC966
525 if (ShapeAnalysis_Curve::IsPeriodic(theCurve)) {
526 ElCLib::AdjustPeriodic(cf,cl,Precision::PConfusion(),First,Last); //:a7 abv 11 Feb 98: preci -> PConfusion()
528 else if (First < Last) {
531 else if (theCurve->IsClosed()) {
532 // l'un des points projecte se trouve sur l'origine du parametrage
533 // de la courbe 3D. L algo a donne cl +- preci au lieu de cf ou vice-versa
534 // DANGER precision 3d applique a une espace 1d
536 // Last = cf au lieu de Last = cl
537 if (Abs(Last - cf) < Precision::PConfusion() /*preci*/) Last = cl ;
538 // First = cl au lieu de First = cf
539 else if (Abs(First - cl) < Precision::PConfusion() /*preci*/) First = cf;
541 // on se trouve dans un cas ou l origine est traversee
542 // illegal sur une courbe fermee non periodique
543 // on inverse quand meme les parametres !!!!!!
545 //:S4136 abv 20 Apr 99: r0701_ug.stp #6230: add check in 3d
546 if ( theCurve->Value(First).Distance(theCurve->Value(cf)) < preci ) First = cf;
547 if ( theCurve->Value(Last).Distance(theCurve->Value(cl)) < preci ) Last = cl;
548 if ( First > Last ) {
550 cout << "Warning : parameter range of edge crossing non periodic curve origin" << endl;
552 Standard_Real tmp = First;
558 // The curve is closed within the 3D tolerance
559 else if (theCurve->IsKind(STANDARD_TYPE(Geom_BSplineCurve))) {
560 Handle(Geom_BSplineCurve) aBSpline =
561 Handle(Geom_BSplineCurve)::DownCast(theCurve);
562 if (aBSpline->StartPoint().Distance(aBSpline->EndPoint()) <= preci ) {
563 //:S4136 <= BRepAPI::Precision()) {
564 // l'un des points projecte se trouve sur l'origine du parametrage
565 // de la courbe 3D. L algo a donne cl +- preci au lieu de cf ou vice-versa
566 // DANGER precision 3d applique a une espace 1d
568 // Last = cf au lieu de Last = cl
569 if (Abs(Last - cf) < Precision::PConfusion() /*preci*/) Last = cl ;
570 // First = cl au lieu de First = cf
571 else if (Abs(First - cl) < Precision::PConfusion() /*preci*/) First = cf;
573 // on se trouve dans un cas ou l origine est traversee
574 // illegal sur une courbe fermee non periodique
575 // on inverse quand meme les parametres !!!!!!
578 cout << "Warning : parameter range of edge crossing non periodic curve origin" << endl;
580 Standard_Real tmp = First;
585 //abv 15.03.00 #72 bm1_pe_t4 protection of exceptions in draw
586 else if ( First > Last ) {
588 cout << "Warning: parameter range is bad; curve reversed" << endl;
590 First = theCurve->ReversedParameter ( First );
591 Last = theCurve->ReversedParameter ( Last );
594 //:j9 abv 11 Dec 98: PRO7747 #4875, after :j8: else
595 if (First == Last) { //gka 10.07.1998 file PRO7656 entity 33334
596 First = cf; Last = cl;
597 return Standard_False;
602 cout << "UpdateParam3d Failed" << endl;
603 cout << " - Curve Type : " << theCurve->DynamicType() << endl;
604 cout << " - Param 1 : " << First << endl;
605 cout << " - Param 2 : " << Last << endl;
607 //abv 15.03.00 #72 bm1_pe_t4 protection of exceptions in draw
608 if ( First > Last ) {
610 cout << "Warning: parameter range is bad; curve reversed" << endl;
612 First = theCurve->ReversedParameter ( First );
613 Last = theCurve->ReversedParameter ( Last );
616 //pdn 11.01.99 #144 bm1_pe_t4 protection of exceptions in draw
618 First -= Precision::PConfusion();
619 Last += Precision::PConfusion();
621 return Standard_False;
623 return Standard_True;
626 //=======================================================================
627 //function : FillBndBox
628 //purpose : WORK-AROUND for methods brom BndLib which give not exact bounds
629 //=======================================================================
631 // search for extremum using Newton
632 static Standard_Integer SearchForExtremum (const Handle(Geom2d_Curve)& C2d,
633 const Standard_Real First,
634 const Standard_Real Last,
639 Standard_Real prevpar;
640 Standard_Integer nbOut = 0;
641 for (Standard_Integer i = 0; i <10; i++) {
645 C2d->D2 ( par, res, D1, D2 );
646 Standard_Real Det = ( D2 * dir );
647 if ( Abs ( Det ) < 1e-10 ) return Standard_True;
649 par -= ( D1 * dir ) / Det;
650 if ( Abs ( par - prevpar ) < Precision::PConfusion() ) return Standard_True;
654 if (nbOut++ > 2 || prevpar == First)
655 return Standard_False;
660 if (nbOut++ > 2 || prevpar == Last)
661 return Standard_False;
665 return Standard_True;
668 void ShapeAnalysis_Curve::FillBndBox (const Handle(Geom2d_Curve)& C2d,
669 const Standard_Real First,
670 const Standard_Real Last,
671 const Standard_Integer NPoints,
672 const Standard_Boolean Exact,
673 Bnd_Box2d &Box) const
676 Standard_Integer nseg = (NPoints < 2 ? 1 : NPoints - 1);
677 Standard_Real step = (Last - First) / nseg;
678 for (Standard_Integer i = 0; i <= nseg; i++) {
679 Standard_Real par = First + i * step;
680 gp_Pnt2d pnt = C2d->Value(par);
686 // We should solve the task on intervals of C2 continuity.
687 Geom2dAdaptor_Curve anAC(C2d, First, Last);
688 Standard_Integer nbInt = anAC.NbIntervals(GeomAbs_C2);
689 // If we have only 1 interval then use input NPoints parameter to get samples.
690 Standard_Integer nbSamples = (nbInt < 2 ? NPoints - 1 : nbInt);
691 TColStd_Array1OfReal aParams(1, nbSamples + 1);
692 if (nbSamples == nbInt)
693 anAC.Intervals(aParams, GeomAbs_C2);
695 Standard_Real step = (Last - First) / nbSamples;
696 for (Standard_Integer i = 0; i <= nbSamples; i++)
697 aParams(i+1) = First + i * step;
699 for (Standard_Integer i = 1; i <= nbSamples + 1; i++) {
700 Standard_Real aPar1 = aParams(i);
701 gp_Pnt2d aPnt = C2d->Value(aPar1);
703 if (i <= nbSamples) {
704 Standard_Real aPar2 = aParams(i + 1);
705 Standard_Real par = (aPar1 + aPar2) * 0.5;
707 Standard_Real parextr = par;
708 if (SearchForExtremum(C2d, aPar1, aPar2, gp_Vec2d(1, 0), parextr, pextr)) {
712 if (SearchForExtremum(C2d, aPar1, aPar2, gp_Vec2d(0, 1), parextr, pextr)) {
719 //=======================================================================
720 //function : SelectForwardSeam
722 //=======================================================================
724 Standard_Integer ShapeAnalysis_Curve::SelectForwardSeam(const Handle(Geom2d_Curve)& C1,
725 const Handle(Geom2d_Curve)& C2) const
727 // SelectForward est destine a devenir un outil distinct
728 // Il est sans doute optimisable !
730 Standard_Integer theCurveIndice = 0;
732 Handle(Geom2d_Line) L1 = Handle(Geom2d_Line)::DownCast(C1);
734 // if we have BoundedCurve, create a line from C1
735 Handle(Geom2d_BoundedCurve) BC1 = Handle(Geom2d_BoundedCurve)::DownCast(C1);
736 if (BC1.IsNull()) return theCurveIndice;
737 gp_Pnt2d StartBC1 = BC1->StartPoint();
738 gp_Pnt2d EndBC1 = BC1->EndPoint();
739 gp_Vec2d VecBC1(StartBC1, EndBC1);
740 if (VecBC1.SquareMagnitude() < gp::Resolution()) return theCurveIndice;
741 L1 = new Geom2d_Line(StartBC1, VecBC1);
744 Handle(Geom2d_Line) L2 = Handle(Geom2d_Line)::DownCast(C2);
746 // if we have BoundedCurve, creates a line from C2
747 Handle(Geom2d_BoundedCurve) BC2 = Handle(Geom2d_BoundedCurve)::DownCast(C2);
748 if (BC2.IsNull()) return theCurveIndice;
749 gp_Pnt2d StartBC2 = BC2->StartPoint();
750 gp_Pnt2d EndBC2 = BC2->EndPoint();
751 gp_Vec2d VecBC2(StartBC2, EndBC2);
752 if (VecBC2.SquareMagnitude() < gp::Resolution()) return theCurveIndice;
753 L2 = new Geom2d_Line(StartBC2, VecBC2);
756 Standard_Boolean UdirPos, UdirNeg, VdirPos, VdirNeg;
757 UdirPos = UdirNeg = VdirPos = VdirNeg = Standard_False;
759 gp_Dir2d theDir = L1->Direction();
760 gp_Pnt2d theLoc1 = L1->Location();
761 gp_Pnt2d theLoc2 = L2->Location();
763 if (theDir.X() > 0.) {
764 UdirPos = Standard_True; //szv#4:S4163:12Mar99 Udir unused
765 } else if (theDir.X() < 0.) {
766 UdirNeg = Standard_True; //szv#4:S4163:12Mar99 Udir unused
767 } else if (theDir.Y() > 0.) {
768 VdirPos = Standard_True; //szv#4:S4163:12Mar99 Vdir unused
769 } else if (theDir.Y() < 0.) {
770 VdirNeg = Standard_True; //szv#4:S4163:12Mar99 Vdir unused
774 // max of Loc1.X() Loc2.X()
775 if (theLoc1.X() > theLoc2.X()) theCurveIndice = 1;
776 else theCurveIndice = 2;
777 } else if (VdirNeg) {
778 if (theLoc1.X() > theLoc2.X()) theCurveIndice = 2;
779 else theCurveIndice = 1;
780 } else if (UdirPos) {
781 // min of Loc1.X() Loc2.X()
782 if (theLoc1.Y() < theLoc2.Y()) theCurveIndice = 1;
783 else theCurveIndice = 2;
784 } else if (UdirNeg) {
785 if (theLoc1.Y() < theLoc2.Y()) theCurveIndice = 2;
786 else theCurveIndice = 1;
789 return theCurveIndice;
793 //=============================================================================
794 // Static methods for IsPlanar
796 //=============================================================================
798 static gp_XYZ GetAnyNormal ( gp_XYZ orig )
801 if ( Abs ( orig.Z() ) < Precision::Confusion() )
802 Norm.SetCoord ( 0, 0, 1 );
804 Norm.SetCoord ( orig.Z(), 0, -orig.X() );
805 Standard_Real nrm = Norm.Modulus();
806 if ( nrm < Precision::Confusion() ) Norm.SetCoord ( 0, 0, 1 );
807 else Norm = Norm / nrm;
812 //=======================================================================
813 //function : GetSamplePoints
815 //=======================================================================
816 static void AppendControlPoles (TColgp_SequenceOfPnt& seq,
817 const Handle(Geom_Curve)& curve)
819 if ( curve->IsKind(STANDARD_TYPE(Geom_Line))) {
820 seq.Append(curve->Value(0));
821 seq.Append(curve->Value(1));
822 } else if ( curve->IsKind(STANDARD_TYPE(Geom_Conic))) {
823 seq.Append(curve->Value(0));
824 seq.Append(curve->Value(M_PI/2));
825 seq.Append(curve->Value(M_PI));
826 } else if ( curve->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
827 //DeclareAndCast(Geom_TrimmedCurve, Trimmed, curve);
828 Handle(Geom_TrimmedCurve) Trimmed = Handle(Geom_TrimmedCurve)::DownCast (curve);
829 // AppendControlPoles(seq,Trimmed->BasisCurve());
830 Handle(Geom_Curve) aBaseCrv = Trimmed->BasisCurve();
831 Standard_Boolean done = Standard_False;
832 if ( aBaseCrv->IsKind(STANDARD_TYPE(Geom_BSplineCurve))) {
835 Handle(Geom_Geometry) Ctmp = aBaseCrv->Copy();
836 Handle(Geom_BSplineCurve) bslp = Handle(Geom_BSplineCurve)::DownCast(Ctmp);
837 bslp->Segment(curve->FirstParameter(), curve->LastParameter());
838 AppendControlPoles(seq,bslp);
839 done = Standard_True;
841 catch (Standard_Failure) {
844 else if ( aBaseCrv->IsKind(STANDARD_TYPE(Geom_BezierCurve))) {
847 Handle(Geom_Geometry) Ctmp = aBaseCrv->Copy();
848 Handle(Geom_BezierCurve) bz = Handle(Geom_BezierCurve)::DownCast(Ctmp);
849 bz->Segment(curve->FirstParameter(), curve->LastParameter());
850 AppendControlPoles(seq,bz);
851 done = Standard_True;
853 catch (Standard_Failure) {
857 seq.Append(curve->Value(curve->FirstParameter()));
858 seq.Append(curve->Value((curve->FirstParameter() + curve->LastParameter())/2.));
859 seq.Append(curve->Value(curve->LastParameter()));
861 } else if ( curve->IsKind(STANDARD_TYPE(Geom_OffsetCurve))) {
862 //DeclareAndCast(Geom_OffsetCurve, OffsetC, curve);
863 Handle(Geom_OffsetCurve) OffsetC = Handle(Geom_OffsetCurve)::DownCast (curve);
864 // AppendControlPoles(seq,OffsetC->BasisCurve());
865 seq.Append(curve->Value(curve->FirstParameter()));
866 seq.Append(curve->Value((curve->FirstParameter() + curve->LastParameter())/2.));
867 seq.Append(curve->Value(curve->LastParameter()));
868 } else if ( curve->IsKind(STANDARD_TYPE(Geom_BSplineCurve))) {
869 //DeclareAndCast(Geom_BSplineCurve, BSpline, curve);
870 Handle(Geom_BSplineCurve) BSpline = Handle(Geom_BSplineCurve)::DownCast (curve);
871 TColgp_Array1OfPnt Poles(1,BSpline->NbPoles());
872 BSpline->Poles(Poles);
873 for(Standard_Integer i = 1; i <= BSpline->NbPoles(); i++)
874 seq.Append(Poles(i));
875 } else if ( curve->IsKind(STANDARD_TYPE(Geom_BezierCurve))) {
876 //DeclareAndCast(Geom_BezierCurve, Bezier, curve);
877 //Handle(Geom_BezierCurve) Bezier = Handle(Geom_BezierCurve)::DownCast(curve);
878 Handle(Geom_BezierCurve) Bezier = Handle(Geom_BezierCurve)::DownCast (curve);
879 TColgp_Array1OfPnt Poles(1,Bezier->NbPoles());
880 Bezier->Poles(Poles);
881 for(Standard_Integer i = 1; i <= Bezier->NbPoles(); i++)
882 seq.Append(Poles(i));
888 //=======================================================================
889 //function : IsPlanar
890 //purpose : Detects if points lie in some plane and returns normal
891 //=======================================================================
893 Standard_Boolean ShapeAnalysis_Curve::IsPlanar (const TColgp_Array1OfPnt& pnts,
895 const Standard_Real preci)
897 Standard_Real precision = (preci > 0.0)? preci : Precision::Confusion();
898 Standard_Boolean noNorm = (Normal.SquareModulus() == 0);
900 if (pnts.Length() < 3) {
901 gp_XYZ N1 = pnts(1).XYZ()-pnts(2).XYZ();
903 Normal = GetAnyNormal(N1);
904 return Standard_True;
906 return Abs(N1*Normal) < Precision::Confusion();
911 //define a center point
912 gp_XYZ aCenter(0,0,0);
913 Standard_Integer i = 1;
914 for (; i <= pnts.Length(); i++)
915 aCenter += pnts(i).XYZ();
916 aCenter/=pnts.Length();
919 aMaxDir = pnts(1).XYZ() - aCenter;
920 Normal = (pnts(pnts.Length()).XYZ() - aCenter) ^ aMaxDir;
922 for ( i = 1; i < pnts.Length(); i++) {
923 gp_XYZ aTmpDir = pnts(i+1).XYZ() - aCenter;
924 if(aTmpDir.SquareModulus() > aMaxDir.SquareModulus())
927 gp_XYZ aDelta = (pnts(i).XYZ() - aCenter) ^ (pnts(i+1).XYZ() - aCenter);
928 if(Normal*aDelta < 0)
934 // check if points are linear
935 Standard_Real nrm = Normal.Modulus();
936 if ( nrm < Precision::Confusion() ) {
937 Normal = GetAnyNormal(aMaxDir);
938 return Standard_True;
940 Normal = Normal / nrm;
942 Standard_Real mind = RealLast(), maxd = -RealLast(), dev;
943 for (Standard_Integer i = 1; i <= pnts.Length(); i++) {
944 dev = pnts(i).XYZ() * Normal;
945 if (dev < mind) mind = dev;
946 if (dev > maxd) maxd = dev;
949 return ((maxd - mind) <= precision);
953 //=======================================================================
954 //function : IsPlanar
956 //=======================================================================
958 Standard_Boolean ShapeAnalysis_Curve::IsPlanar (const Handle(Geom_Curve)& curve,
960 const Standard_Real preci)
962 Standard_Real precision = (preci > 0.0)? preci : Precision::Confusion();
963 Standard_Boolean noNorm = (Normal.SquareModulus() == 0);
965 if (curve->IsKind(STANDARD_TYPE(Geom_Line))) {
966 //DeclareAndCast(Geom_Line, Line, curve);
967 Handle(Geom_Line) Line = Handle(Geom_Line)::DownCast (curve);
968 gp_XYZ N1 = Line->Position().Direction().XYZ();
970 Normal = GetAnyNormal(N1);
971 return Standard_True;
973 return Abs(N1*Normal) < Precision::Confusion();
976 if (curve->IsKind(STANDARD_TYPE(Geom_Conic))) {
977 //DeclareAndCast(Geom_Conic, Conic, curve);
978 Handle(Geom_Conic) Conic = Handle(Geom_Conic)::DownCast (curve);
979 gp_XYZ N1 = Conic->Axis().Direction().XYZ();
982 return Standard_True;
984 gp_XYZ aVecMul = N1^Normal;
985 return aVecMul.SquareModulus() < Precision::SquareConfusion();
988 if (curve->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
989 //DeclareAndCast(Geom_TrimmedCurve, Trimmed, curve);
990 Handle(Geom_TrimmedCurve) Trimmed = Handle(Geom_TrimmedCurve)::DownCast (curve);
991 return IsPlanar(Trimmed->BasisCurve(),Normal,precision);
994 if (curve->IsKind(STANDARD_TYPE(Geom_OffsetCurve))) {
995 //DeclareAndCast(Geom_OffsetCurve, OffsetC, curve);
996 Handle(Geom_OffsetCurve) OffsetC = Handle(Geom_OffsetCurve)::DownCast (curve);
997 return IsPlanar(OffsetC->BasisCurve(),Normal,precision);
1000 if (curve->IsKind(STANDARD_TYPE(Geom_BSplineCurve))) {
1001 //DeclareAndCast(Geom_BSplineCurve, BSpline, curve);
1002 Handle(Geom_BSplineCurve) BSpline = Handle(Geom_BSplineCurve)::DownCast (curve);
1003 TColgp_Array1OfPnt Poles(1,BSpline->NbPoles());
1004 BSpline->Poles(Poles);
1005 return IsPlanar(Poles,Normal,precision);
1008 if (curve->IsKind(STANDARD_TYPE(Geom_BezierCurve))) {
1009 //DeclareAndCast(Geom_BezierCurve, Bezier, curve);
1010 Handle(Geom_BezierCurve) Bezier = Handle(Geom_BezierCurve)::DownCast (curve);
1011 TColgp_Array1OfPnt Poles(1,Bezier->NbPoles());
1012 Bezier->Poles(Poles);
1013 return IsPlanar(Poles,Normal,precision);
1016 if (curve->IsKind(STANDARD_TYPE(ShapeExtend_ComplexCurve))) {
1017 //DeclareAndCast(ShapeExtend_ComplexCurve, Complex, curve);
1018 Handle(ShapeExtend_ComplexCurve) Complex = Handle(ShapeExtend_ComplexCurve)::DownCast (curve);
1019 TColgp_SequenceOfPnt sequence;
1020 Standard_Integer i; // svv Jan11 2000 : porting on DEC
1021 for (i = 1; i <= Complex->NbCurves(); i++)
1022 AppendControlPoles(sequence,Complex->Curve(i));
1023 TColgp_Array1OfPnt Poles(1,sequence.Length());
1024 for (i=1; i <= sequence.Length(); i++) Poles(i) = sequence(i);
1025 return IsPlanar(Poles,Normal,precision);
1028 return Standard_False;
1031 //=======================================================================
1032 //function : GetSamplePoints
1034 //=======================================================================
1036 Standard_Boolean ShapeAnalysis_Curve::GetSamplePoints (const Handle(Geom_Curve)& curve,
1037 const Standard_Real first,
1038 const Standard_Real last,
1039 TColgp_SequenceOfPnt& seq)
1041 Standard_Real adelta = curve->LastParameter() - curve->FirstParameter();
1043 return Standard_False;
1045 Standard_Integer aK = (Standard_Integer)ceil ((last - first) / adelta);
1046 Standard_Integer nbp =100*aK;
1047 if(curve->IsKind(STANDARD_TYPE(Geom_Line)))
1049 else if(curve->IsKind(STANDARD_TYPE(Geom_Circle)))
1052 else if (curve->IsKind(STANDARD_TYPE(Geom_BSplineCurve))) {
1053 Handle(Geom_BSplineCurve) aBspl = Handle(Geom_BSplineCurve)::DownCast(curve);
1055 nbp = aBspl->NbKnots() * aBspl->Degree()*aK;
1056 if(nbp < 2.0) nbp=2;
1058 else if (curve->IsKind(STANDARD_TYPE(Geom_BezierCurve))) {
1059 Handle(Geom_BezierCurve) aB = Handle(Geom_BezierCurve)::DownCast(curve);
1060 nbp = 3 + aB->NbPoles();
1062 else if(curve->IsKind(STANDARD_TYPE(Geom_OffsetCurve))) {
1063 Handle(Geom_OffsetCurve) aC = Handle(Geom_OffsetCurve)::DownCast(curve);
1064 return GetSamplePoints(aC->BasisCurve(),first,last,seq);
1066 else if(curve->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
1067 Handle(Geom_TrimmedCurve) aC = Handle(Geom_TrimmedCurve)::DownCast(curve);
1068 return GetSamplePoints(aC->BasisCurve(),first,last,seq);
1071 GeomAdaptor_Curve GAC(curve);
1072 Standard_Real step = ( last - first ) / (Standard_Real)( nbp - 1 );
1073 for (Standard_Integer i = 0; i < nbp - 1; ++i)
1074 seq.Append(GAC.Value(first + step * i));
1075 seq.Append(GAC.Value(last));
1076 return Standard_True;
1078 //=======================================================================
1079 //function : GetSamplePoints
1081 //=======================================================================
1083 Standard_Boolean ShapeAnalysis_Curve::GetSamplePoints (const Handle(Geom2d_Curve)& curve,
1084 const Standard_Real first,
1085 const Standard_Real last,
1086 TColgp_SequenceOfPnt2d& seq)
1088 //:abv 05.06.02: TUBE.stp
1089 // Use the same distribution of points as BRepTopAdaptor_FClass2d for consistency
1090 Geom2dAdaptor_Curve C ( curve, first, last );
1091 Standard_Integer nbs = Geom2dInt_Geom2dCurveTool::NbSamples(C);
1092 //-- Attention aux bsplines rationnelles de degree 3. (bouts de cercles entre autres)
1093 if (nbs > 2) nbs*=4;
1094 Standard_Real step = ( last - first ) / (Standard_Real)( nbs - 1 );
1095 for (Standard_Integer i = 0; i < nbs - 1; ++i)
1096 seq.Append(C.Value(first + step * i));
1097 seq.Append(C.Value(last));
1098 return Standard_True;
1103 if ( curve->IsKind(STANDARD_TYPE(Geom2d_Line))) {
1104 seq.Append(curve->Value(first));
1105 seq.Append(curve->Value(last));
1106 return Standard_True;
1108 else if(curve->IsKind(STANDARD_TYPE(Geom2d_Conic))) {
1109 step = Min ( M_PI, last-first ) / 19; //:abv 05.06.02 TUBE.stp #19209...: M_PI/16
1110 // if( step>(last-first) ) {
1111 // seq.Append(curve->Value(first));
1112 // seq.Append(curve->Value((last+first)/2));
1113 // seq.Append(curve->Value(last));
1114 // return Standard_True;
1117 Standard_Real par=first;
1118 for(i=0; par<last; i++) {
1119 seq.Append(curve->Value(par));
1122 seq.Append(curve->Value(last));
1123 return Standard_True;
1126 else if ( curve->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve))) {
1127 DeclareAndCast(Geom2d_TrimmedCurve, Trimmed, curve);
1128 return GetControlPoints(Trimmed->BasisCurve(), seq, first, last);
1130 else if ( curve->IsKind(STANDARD_TYPE(Geom2d_BSplineCurve))) {
1131 DeclareAndCast(Geom2d_BSplineCurve, aBs, curve);
1132 TColStd_SequenceOfReal aSeqParam;
1134 aSeqParam.Append(first);
1135 for(i=1; i<=aBs->NbKnots(); i++) {
1136 if( aBs->Knot(i)>first && aBs->Knot(i)<last )
1137 aSeqParam.Append(aBs->Knot(i));
1139 aSeqParam.Append(last);
1140 Standard_Integer NbPoints=aBs->Degree();
1141 if( (aSeqParam.Length()-1)*NbPoints>10 ) {
1142 for(i=1; i<aSeqParam.Length(); i++) {
1143 Standard_Real FirstPar = aSeqParam.Value(i);
1144 Standard_Real LastPar = aSeqParam.Value(i+1);
1145 step = (LastPar-FirstPar)/NbPoints;
1146 for(Standard_Integer k=0; k<NbPoints; k++ ) {
1147 aBs->D0(FirstPar+step*k, Ptmp);
1151 aBs->D0(last, Ptmp);
1153 return Standard_True;
1156 step = (last-first)/10;
1157 for(Standard_Integer k=0; k<=10; k++ ) {
1158 aBs->D0(first+step*k, Ptmp);
1161 return Standard_True;
1165 else if ( curve->IsKind(STANDARD_TYPE(Geom2d_BezierCurve))) {
1166 DeclareAndCast(Geom2d_BezierCurve, aBz, curve);
1168 Standard_Integer NbPoints=aBz->Degree();
1169 step = (last-first)/NbPoints;
1170 for(Standard_Integer k=0; k<NbPoints; k++ ) {
1171 aBz->D0(first+step*k, Ptmp);
1174 aBz->D0(last, Ptmp);
1176 return Standard_True;
1179 return Standard_False;
1183 //=======================================================================
1184 //function : IsClosed
1186 //=======================================================================
1188 Standard_Boolean ShapeAnalysis_Curve::IsClosed(const Handle(Geom_Curve)& theCurve,
1189 const Standard_Real preci)
1191 if (theCurve->IsClosed())
1192 return Standard_True;
1194 Standard_Real prec = Max (preci, Precision::Confusion());
1197 f = theCurve->FirstParameter();
1198 l = theCurve->LastParameter();
1200 if (Precision::IsInfinite (f) || Precision::IsInfinite (l))
1201 return Standard_False;
1203 Standard_Real aClosedVal = theCurve->Value(f).SquareDistance(theCurve->Value(l));
1204 Standard_Real preci2 = prec*prec;
1206 return (aClosedVal <= preci2);
1208 //=======================================================================
1209 //function : IsPeriodic
1211 //=======================================================================
1213 Standard_Boolean ShapeAnalysis_Curve::IsPeriodic(const Handle(Geom_Curve)& theCurve)
1215 // 15.11.2002 PTV OCC966
1216 // remove regressions in DE tests (diva, divb, divc, toe3) in KAS:dev
1217 // ask IsPeriodic on BasisCurve
1218 Handle(Geom_Curve) aTmpCurve = theCurve;
1219 while ( (aTmpCurve->IsKind(STANDARD_TYPE(Geom_OffsetCurve))) ||
1220 (aTmpCurve->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) ) {
1221 if (aTmpCurve->IsKind(STANDARD_TYPE(Geom_OffsetCurve)))
1222 aTmpCurve = Handle(Geom_OffsetCurve)::DownCast(aTmpCurve)->BasisCurve();
1223 if (aTmpCurve->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)))
1224 aTmpCurve = Handle(Geom_TrimmedCurve)::DownCast(aTmpCurve)->BasisCurve();
1226 return aTmpCurve->IsPeriodic();
1229 Standard_Boolean ShapeAnalysis_Curve::IsPeriodic(const Handle(Geom2d_Curve)& theCurve)
1231 // 15.11.2002 PTV OCC966
1232 // remove regressions in DE tests (diva, divb, divc, toe3) in KAS:dev
1233 // ask IsPeriodic on BasisCurve
1234 Handle(Geom2d_Curve) aTmpCurve = theCurve;
1235 while ( (aTmpCurve->IsKind(STANDARD_TYPE(Geom2d_OffsetCurve))) ||
1236 (aTmpCurve->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve))) ) {
1237 if (aTmpCurve->IsKind(STANDARD_TYPE(Geom2d_OffsetCurve)))
1238 aTmpCurve = Handle(Geom2d_OffsetCurve)::DownCast(aTmpCurve)->BasisCurve();
1239 if (aTmpCurve->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)))
1240 aTmpCurve = Handle(Geom2d_TrimmedCurve)::DownCast(aTmpCurve)->BasisCurve();
1242 return aTmpCurve->IsPeriodic();