OpenCascade/occt
1
0
Fork 0
mirror of https://git.dev.opencascade.org/repos/occt.git synced 2025-04-05 18:16:23 +03:00
Code Packages Releases Activity

0026339: [Regression in 6.9.0] Projecting a curve hangs

Browse Source 
Changed computation of point projection to more correct.
Calculation periodicity information added to cache.
Test case for issue CR26339

Small correction of test case for issue CR26339
This commit is contained in:
aml 2015-07-02 14:19:13 +03:00 committed by bugmaster
parent 4b65fc7750
commit ea3b6e72be
3 changed files with 254 additions and 134 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

346
src/ProjLib/ProjLib_ComputeApproxOnPolarSurface.cxx
View File

@ -81,204 +81,273 @@
//static Standard_Integer compteur = 0; //static Standard_Integer compteur = 0;
#endif #endif
struct aFuncStruct
{
Handle(Adaptor3d_HSurface) mySurf; // Surface where to project.
Handle(Adaptor3d_HCurve) myCurve; // Curve to project.
Handle(Adaptor2d_HCurve2d) myInitCurve2d; // Initial 2dcurve projection.
Standard_Real mySqProjOrtTol; // Used to filter non-orthogonal projected point.
Standard_Real myTolU;
Standard_Real myTolV;
Standard_Real myPeriod[2]; // U and V period correspondingly.
};
//=======================================================================
//function : aFuncValue
//purpose : compute functional value in (theU,theV) point
//=======================================================================
static Standard_Real anOrthogSqValue(const gp_Pnt& aBasePnt,
const Handle(Adaptor3d_HSurface)& Surf,
const Standard_Real theU,
const Standard_Real theV)
{
// Since find projection, formula is:
// F1 = Dot(S_U, Vec(aBasePnt, aProjPnt))
// F2 = Dot(S_V, Vec(aBasePnt, aProjPnt))
gp_Pnt aProjPnt;
gp_Vec aSu, aSv;
Surf->D1(theU, theV, aProjPnt, aSu, aSv);
gp_Vec aBaseVec(aBasePnt, aProjPnt);
if (aSu.SquareMagnitude() > Precision::SquareConfusion())
aSu.Normalize();
if (aSv.SquareMagnitude() > Precision::SquareConfusion())
aSv.Normalize();
Standard_Real aFirstPart = aSu.Dot(aBaseVec);
Standard_Real aSecondPart = aSv.Dot(aBaseVec);
return (aFirstPart * aFirstPart + aSecondPart * aSecondPart);
}
//======================================================================= //=======================================================================
//function : Value //function : Value
//purpose : (OCC217 - apo)- Compute Point2d that project on polar surface(<Surf>) 3D<Curve> //purpose : (OCC217 - apo)- Compute Point2d that project on polar surface(<Surf>) 3D<Curve>
// <InitCurve2d> use for calculate start 2D point. // <InitCurve2d> use for calculate start 2D point.
//======================================================================= //=======================================================================
static gp_Pnt2d Function_Value(const Standard_Real U, static gp_Pnt2d Function_Value(const Standard_Real theU,
const Handle(Adaptor3d_HSurface)& Surf, const aFuncStruct& theData)
const Handle(Adaptor3d_HCurve)& Curve,
const Handle(Adaptor2d_HCurve2d)& InitCurve2d,
//OCC217
const Standard_Real DistTol3d, const Standard_Real tolU, const Standard_Real tolV)
//const Standard_Real Tolerance)
{ {
//OCC217 gp_Pnt2d p2d = theData.myInitCurve2d->Value(theU) ;
//Standard_Real Tol3d = 100*Tolerance; gp_Pnt p = theData.myCurve->Value(theU);
gp_Pnt aSurfPnt = theData.mySurf->Value(p2d.X(), p2d.Y());
Standard_Real aSurfPntDist = aSurfPnt.SquareDistance(p);
gp_Pnt2d p2d = InitCurve2d->Value(U) ;
gp_Pnt p = Curve->Value(U);
// Curve->D0(U,p) ;
Standard_Real Uinf, Usup, Vinf, Vsup; Standard_Real Uinf, Usup, Vinf, Vsup;
Uinf = Surf->Surface().FirstUParameter(); Uinf = theData.mySurf->Surface().FirstUParameter();
Usup = Surf->Surface().LastUParameter(); Usup = theData.mySurf->Surface().LastUParameter();
Vinf = Surf->Surface().FirstVParameter(); Vinf = theData.mySurf->Surface().FirstVParameter();
Vsup = Surf->Surface().LastVParameter(); Vsup = theData.mySurf->Surface().LastVParameter();
// Check case when curve is close to co-parametrized isoline on surf.
if (Abs (p2d.X() - Uinf) < Precision::PConfusion() ||
Abs (p2d.X() - Usup) < Precision::PConfusion() )
{
// V isoline.
gp_Pnt aPnt;
theData.mySurf->D0(p2d.X(), theU, aPnt);
if (aPnt.SquareDistance(p) < aSurfPntDist)
p2d.SetY(theU);
}
if (Abs (p2d.Y() - Vinf) < Precision::PConfusion() ||
Abs (p2d.Y() - Vsup) < Precision::PConfusion() )
{
// U isoline.
gp_Pnt aPnt;
theData.mySurf->D0(theU, p2d.Y(), aPnt);
if (aPnt.SquareDistance(p) < aSurfPntDist)
p2d.SetX(theU);
}
Standard_Integer decalU = 0, decalV = 0; Standard_Integer decalU = 0, decalV = 0;
Standard_Real U0 = p2d.X(), V0 = p2d.Y(); Standard_Real U0 = p2d.X(), V0 = p2d.Y();
GeomAbs_SurfaceType Type = Surf->GetType(); GeomAbs_SurfaceType Type = theData.mySurf->GetType();
if((Type != GeomAbs_BSplineSurface) && if((Type != GeomAbs_BSplineSurface) &&
(Type != GeomAbs_BezierSurface) && (Type != GeomAbs_BezierSurface) &&
(Type != GeomAbs_OffsetSurface) ) { (Type != GeomAbs_OffsetSurface) )
{
// Analytical cases.
Standard_Real S = 0., T = 0.; Standard_Real S = 0., T = 0.;
switch (Type) { switch (Type)
// case GeomAbs_Plane: {
// {
// gp_Pln Plane = Surf->Plane();
// ElSLib::Parameters( Plane, p, S, T);
// break;
// }
case GeomAbs_Cylinder: case GeomAbs_Cylinder:
{ {
gp_Cylinder Cylinder = Surf->Cylinder(); gp_Cylinder Cylinder = theData.mySurf->Cylinder();
ElSLib::Parameters( Cylinder, p, S, T); ElSLib::Parameters( Cylinder, p, S, T);
if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1; if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1; if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
S += decalU*2*M_PI; S += decalU*2*M_PI;
break; break;
} }
case GeomAbs_Cone: case GeomAbs_Cone:
{ {
gp_Cone Cone = Surf->Cone(); gp_Cone Cone = theData.mySurf->Cone();
ElSLib::Parameters( Cone, p, S, T); ElSLib::Parameters( Cone, p, S, T);
if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1; if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1; if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
S += decalU*2*M_PI; S += decalU*2*M_PI;
break; break;
} }
case GeomAbs_Sphere: case GeomAbs_Sphere:
{ {
gp_Sphere Sphere = Surf->Sphere(); gp_Sphere Sphere = theData.mySurf->Sphere();
ElSLib::Parameters( Sphere, p, S, T); ElSLib::Parameters( Sphere, p, S, T);
if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1; if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1; if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
S += decalU*2*M_PI; S += decalU*2*M_PI;
if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1; if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1;
if(V0 > (Vsup+(Vsup-Vinf))) decalV = int((V0 - Vsup+(Vsup-Vinf))/(2*M_PI))+1; if(V0 > (Vsup+(Vsup-Vinf))) decalV = int((V0 - Vsup+(Vsup-Vinf))/(2*M_PI))+1;
T += decalV*2*M_PI; T += decalV*2*M_PI;
if(0.4*M_PI < Abs(U0 - S) && Abs(U0 - S) < 1.6*M_PI) { if(0.4*M_PI < Abs(U0 - S) && Abs(U0 - S) < 1.6*M_PI)
T = M_PI - T; {
if(U0 < S) T = M_PI - T;
S -= M_PI; if(U0 < S)
else S -= M_PI;
S += M_PI; else
} S += M_PI;
break; }
break;
} }
case GeomAbs_Torus: case GeomAbs_Torus:
{ {
gp_Torus Torus = Surf->Torus(); gp_Torus Torus = theData.mySurf->Torus();
ElSLib::Parameters( Torus, p, S, T); ElSLib::Parameters( Torus, p, S, T);
if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1; if(U0 < Uinf) decalU = -int((Uinf - U0)/(2*M_PI))-1;
if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1; if(U0 > Usup) decalU = int((U0 - Usup)/(2*M_PI))+1;
if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1; if(V0 < Vinf) decalV = -int((Vinf - V0)/(2*M_PI))-1;
if(V0 > Vsup) decalV = int((V0 - Vsup)/(2*M_PI))+1; if(V0 > Vsup) decalV = int((V0 - Vsup)/(2*M_PI))+1;
S += decalU*2*M_PI; T += decalV*2*M_PI; S += decalU*2*M_PI; T += decalV*2*M_PI;
break; break;
} }
default: default:
Standard_NoSuchObject::Raise("ProjLib_ComputeApproxOnPolarSurface::Value"); Standard_NoSuchObject::Raise("ProjLib_ComputeApproxOnPolarSurface::Value");
} }
return gp_Pnt2d(S, T); return gp_Pnt2d(S, T);
} }
////////////////// // Non-analytical case.
Standard_Real Dist2Min = RealLast(); Standard_Real Dist2Min = RealLast();
//OCC217 Standard_Real uperiod = theData.myPeriod[0],
//Standard_Real tolU,tolV ; vperiod = theData.myPeriod[1],
//tolU = Tolerance; u, v;
//tolV = Tolerance;
Standard_Real uperiod =0, vperiod = 0, u, v;
// U0 and V0 are the points within the initialized period // U0 and V0 are the points within the initialized period
// (periode with u and v), // (periode with u and v),
// U1 and V1 are the points for construction of tops // U1 and V1 are the points for construction of tops
if(Surf->IsUPeriodic() || Surf->IsUClosed()) { if(U0 < Uinf)
uperiod = Surf->LastUParameter() - Surf->FirstUParameter(); {
}
if(Surf->IsVPeriodic() || Surf->IsVClosed()) {
vperiod = Surf->LastVParameter() - Surf->FirstVParameter();
}
if(U0 < Uinf) {
if(!uperiod) if(!uperiod)
U0 = Uinf; U0 = Uinf;
else { else
{
decalU = int((Uinf - U0)/uperiod)+1; decalU = int((Uinf - U0)/uperiod)+1;
U0 += decalU*uperiod; U0 += decalU*uperiod;
} }
} }
if(U0 > Usup) { if(U0 > Usup)
{
if(!uperiod) if(!uperiod)
U0 = Usup; U0 = Usup;
else { else
{
decalU = -(int((U0 - Usup)/uperiod)+1); decalU = -(int((U0 - Usup)/uperiod)+1);
U0 += decalU*uperiod; U0 += decalU*uperiod;
} }
} }
if(V0 < Vinf) { if(V0 < Vinf)
{
if(!vperiod) if(!vperiod)
V0 = Vinf; V0 = Vinf;
else { else
{
decalV = int((Vinf - V0)/vperiod)+1; decalV = int((Vinf - V0)/vperiod)+1;
V0 += decalV*vperiod; V0 += decalV*vperiod;
} }
} }
if(V0 > Vsup) { if(V0 > Vsup)
{
if(!vperiod) if(!vperiod)
V0 = Vsup; V0 = Vsup;
else { else
{
decalV = -int((V0 - Vsup)/vperiod)-1; decalV = -int((V0 - Vsup)/vperiod)-1;
V0 += decalV*vperiod; V0 += decalV*vperiod;
} }
} }
// The surface around U0 is reduced // The surface around U0 is reduced.
Standard_Real uLittle = (Usup - Uinf)/10, vLittle = (Vsup - Vinf)/10; Standard_Real uLittle = (Usup - Uinf)/10, vLittle = (Vsup - Vinf)/10;
Standard_Real uInfLi = 0, vInfLi = 0,uSupLi = 0, vSupLi = 0; Standard_Real uInfLi = 0, vInfLi = 0,uSupLi = 0, vSupLi = 0;
if((U0 - Uinf) > uLittle) uInfLi = U0 - uLittle; else uInfLi = Uinf; if((U0 - Uinf) > uLittle) uInfLi = U0 - uLittle; else uInfLi = Uinf;
if((V0 - Vinf) > vLittle) vInfLi = V0 - vLittle; else vInfLi = Vinf; if((V0 - Vinf) > vLittle) vInfLi = V0 - vLittle; else vInfLi = Vinf;
if((Usup - U0) > uLittle) uSupLi = U0 + uLittle; else uSupLi = Usup; if((Usup - U0) > uLittle) uSupLi = U0 + uLittle; else uSupLi = Usup;
if((Vsup - V0) > vLittle) vSupLi = V0 + vLittle; else vSupLi = Vsup; if((Vsup - V0) > vLittle) vSupLi = V0 + vLittle; else vSupLi = Vsup;
// const Adaptor3d_Surface GAS = Surf->Surface();
GeomAdaptor_Surface SurfLittle; GeomAdaptor_Surface SurfLittle;
if (Type == GeomAbs_BSplineSurface) { if (Type == GeomAbs_BSplineSurface)
Handle(Geom_Surface) GBSS(Surf->Surface().BSpline()); {
Handle(Geom_Surface) GBSS(theData.mySurf->Surface().BSpline());
SurfLittle.Load(GBSS, uInfLi, uSupLi, vInfLi, vSupLi); SurfLittle.Load(GBSS, uInfLi, uSupLi, vInfLi, vSupLi);
} }
else if (Type == GeomAbs_BezierSurface) { else if (Type == GeomAbs_BezierSurface)
Handle(Geom_Surface) GS(Surf->Surface().Bezier()); {
Handle(Geom_Surface) GS(theData.mySurf->Surface().Bezier());
SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi); SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi);
} }
else if (Type == GeomAbs_OffsetSurface) { else if (Type == GeomAbs_OffsetSurface)
Handle(Geom_Surface) GS = GeomAdaptor::MakeSurface(Surf->Surface()); {
Handle(Geom_Surface) GS = GeomAdaptor::MakeSurface(theData.mySurf->Surface());
SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi); SurfLittle.Load(GS, uInfLi, uSupLi, vInfLi, vSupLi);
} }
else { else
{
Standard_NoSuchObject::Raise(""); Standard_NoSuchObject::Raise("");
} }
Extrema_GenLocateExtPS locext(p, SurfLittle, U0, V0, tolU, tolV); // Try to run simple search with initial point (U0, V0).
if (locext.IsDone()) { Extrema_GenLocateExtPS locext(p, SurfLittle, U0, V0, theData.myTolU, theData.myTolV);
Dist2Min = locext.SquareDistance(); if (locext.IsDone())
if (Dist2Min < DistTol3d * DistTol3d) { {
(locext.Point()).Parameter(u, v); locext.Point().Parameter(u, v);
gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod); Dist2Min = anOrthogSqValue(p, theData.mySurf, u, v);
return pnt; if (Dist2Min < theData.mySqProjOrtTol && // Point is projection.
} locext.SquareDistance() < aSurfPntDist + Precision::SquareConfusion()) // Point better than initial.
} {
Extrema_ExtPS ext(p, SurfLittle, tolU, tolV) ;
if (ext.IsDone() && ext.NbExt()>=1 ) {
Dist2Min = ext.SquareDistance(1);
Standard_Integer GoodValue = 1;
for ( Standard_Integer i = 2 ; i <= ext.NbExt() ; i++ )
if( Dist2Min > ext.SquareDistance(i)) {
Dist2Min = ext.SquareDistance(i);
GoodValue = i;
}
if (Dist2Min < DistTol3d * DistTol3d) {
(ext.Point(GoodValue)).Parameter(u,v);
gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod); gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod);
return pnt; return pnt;
} }
} }
// Perform whole param space search.
Extrema_ExtPS ext(p, SurfLittle, theData.myTolU, theData.myTolV);
if (ext.IsDone() && ext.NbExt() >= 1)
{
Dist2Min = ext.SquareDistance(1);
Standard_Integer GoodValue = 1;
for (Standard_Integer i = 2 ; i <= ext.NbExt() ; i++ )
{
if( Dist2Min > ext.SquareDistance(i))
{
Dist2Min = ext.SquareDistance(i);
GoodValue = i;
}
}
ext.Point(GoodValue).Parameter(u, v);
Dist2Min = anOrthogSqValue(p, theData.mySurf, u, v);
if (Dist2Min < theData.mySqProjOrtTol && // Point is projection.
ext.SquareDistance(GoodValue) < aSurfPntDist + Precision::SquareConfusion()) // Point better than initial.
{
gp_Pnt2d pnt(u - decalU*uperiod,v - decalV*vperiod);
return pnt;
}
}
// Both searches return bad values, use point from initial 2dcurve.
return p2d; return p2d;
} }
@ -290,11 +359,7 @@ static gp_Pnt2d Function_Value(const Standard_Real U,
class ProjLib_PolarFunction : public AppCont_Function class ProjLib_PolarFunction : public AppCont_Function
{ {
Handle(Adaptor3d_HCurve) myCurve; aFuncStruct myStruct;
Handle(Adaptor2d_HCurve2d) myInitialCurve2d;
Handle(Adaptor3d_HSurface) mySurface;
Standard_Real myTolU, myTolV;
Standard_Real myDistTol3d;
public : public :
@ -302,40 +367,53 @@ class ProjLib_PolarFunction : public AppCont_Function
const Handle(Adaptor3d_HSurface)& Surf, const Handle(Adaptor3d_HSurface)& Surf,
const Handle(Adaptor2d_HCurve2d)& InitialCurve2d, const Handle(Adaptor2d_HCurve2d)& InitialCurve2d,
const Standard_Real Tol3d) const Standard_Real Tol3d)
: myCurve(C),
myInitialCurve2d(InitialCurve2d),
mySurface(Surf),
myTolU(Surf->UResolution(Tol3d)),
myTolV(Surf->VResolution(Tol3d)),
myDistTol3d(100.0*Tol3d)
{ {
myNbPnt = 0; myNbPnt = 0;
myNbPnt2d = 1; myNbPnt2d = 1;
myStruct.myPeriod[0] = 0.0;
myStruct.myPeriod[1] = 0.0;
// Compute once information about periodicity.
if(Surf->IsUPeriodic() || Surf->IsUClosed())
{
myStruct.myPeriod[0] = Surf->LastUParameter() - Surf->FirstUParameter();
}
if(Surf->IsVPeriodic() || Surf->IsVClosed())
{
myStruct.myPeriod[1] = Surf->LastVParameter() - Surf->FirstVParameter();
}
myStruct.myCurve = C;
myStruct.myInitCurve2d = InitialCurve2d;
myStruct.mySurf = Surf;
myStruct.mySqProjOrtTol = 10000.0 * Tol3d * Tol3d;
myStruct.myTolU = Surf->UResolution(Tol3d);
myStruct.myTolV = Surf->VResolution(Tol3d);
} }
~ProjLib_PolarFunction() {} ~ProjLib_PolarFunction() {}
Standard_Real FirstParameter() const Standard_Real FirstParameter() const
{ {
return myCurve->FirstParameter(); return myStruct.myCurve->FirstParameter();
} }
Standard_Real LastParameter() const Standard_Real LastParameter() const
{ {
return myCurve->LastParameter(); return myStruct.myCurve->LastParameter();
} }
gp_Pnt2d Value(const Standard_Real t) const gp_Pnt2d Value(const Standard_Real t) const
{ {
return Function_Value return Function_Value(t, myStruct);
(t,mySurface,myCurve,myInitialCurve2d,myDistTol3d,myTolU,myTolV);
} }
Standard_Boolean Value(const Standard_Real theT, Standard_Boolean Value(const Standard_Real theT,
NCollection_Array1<gp_Pnt2d>& thePnt2d, NCollection_Array1<gp_Pnt2d>& thePnt2d,
NCollection_Array1<gp_Pnt>& /*thePnt*/) const NCollection_Array1<gp_Pnt>& /*thePnt*/) const
{ {
thePnt2d(1) = Function_Value(theT, mySurface, myCurve, myInitialCurve2d, myDistTol3d, myTolU, myTolV); thePnt2d(1) = Function_Value(theT, myStruct);
return Standard_True; return Standard_True;
} }
@ -1386,7 +1464,7 @@ Handle(Geom2d_BSplineCurve)
Standard_Real Tol3d = myTolerance; Standard_Real Tol3d = myTolerance;
Standard_Real DistTol3d = 1.0*Tol3d; Standard_Real DistTol3d = 1.0*Tol3d;
Standard_Real TolU = Surf->UResolution(Tol3d), TolV = Surf->VResolution(Tol3d); Standard_Real TolU = Surf->UResolution(Tol3d), TolV = Surf->VResolution(Tol3d);
Standard_Real Tol2d = Sqrt(TolU*TolU + TolV*TolV); Standard_Real Tol2d = Max(Sqrt(TolU*TolU + TolV*TolV), Precision::PConfusion());
Standard_Integer i; Standard_Integer i;
GeomAbs_SurfaceType TheTypeS = Surf->GetType(); GeomAbs_SurfaceType TheTypeS = Surf->GetType();

2
tests/bugs/moddata_1/bug22759
View File

@ -1,3 +1,5 @@
puts "TODO ?OCC26339 ALL: TEST INCOMPLETE"
puts "============" puts "============"
puts "OCC22759" puts "OCC22759"
puts "============" puts "============"

40
tests/bugs/moddata_3/bug26339 Normal file
View File

@ -0,0 +1,40 @@
puts "============"
puts "OCC26339"
puts "============"
puts ""
#######################################################################
# [Regression in 6.9.0] Projecting a curve hangs
#######################################################################
if { [regexp {Debug mode} [dversion]] } {
if { [regexp {Windows} [dversion]] } {
set max_time 10
} else {
set max_time 10
}
} else {
if { [regexp {Windows} [dversion]] } {
set max_time 3
} else {
set max_time 3
}
}
restore [locate_data_file bug26339_a_1886.brep] f
dchrono h reset
dchrono h start
fixshape r f 1e-5
dchrono h stop
set q [dchrono h show]
regexp {CPU user time: ([-0-9.+eE]+) seconds} $q full z
puts "$z"
if { $z > ${max_time} } {
puts "Elapsed time of projecting a curve is more than ${max_time} seconds - Faulty"
} else {
puts "Elapsed time of projecting a curve is less than ${max_time} seconds - OK"
}