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

0022883: Extrema can not find projection of 3D point on surface.

Browse Source 
Optimization of process of surface discretization: increase the number of parametric points in case of a complex surface geometry (BSpline and Bezier surfaces).
Small correction process of building subgrid in Extrema_GenExtPS::FindSolution().
Minor corrections (formatting, duplicate statements)
This commit is contained in:
ama 2012-05-05 17:04:19 +04:00
parent 8e0115e401
commit 5368adff54
6 changed files with 348 additions and 132 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

65
src/Extrema/Extrema_ExtPS.cxx
View File

@ -60,41 +60,43 @@ static Standard_Boolean IsoIsDeg (const Adaptor3d_Surface& S,
Standard_Real Step,D1NormMax; Standard_Real Step,D1NormMax;
if (IT == GeomAbs_IsoV) if (IT == GeomAbs_IsoV)
{ {
Step = (U2 - U1)/10; if( !Precision::IsInfinite(U1) && !Precision::IsInfinite(U2) )
if(Step < Precision::PConfusion()) {
return Standard_False;
}
if(Step < Precision::PConfusion()) {
return Standard_False;
}
D1NormMax=0.;
for (T=U1;T<=U2;T=T+Step)
{ {
S.D1(T,Param,P,D1U,D1V); Step = (U2 - U1)/10;
D1NormMax=Max(D1NormMax,D1U.Magnitude()); if(Step < Precision::PConfusion()) {
} return Standard_False;
}
D1NormMax=0.;
if (D1NormMax >TolMax || D1NormMax < TolMin ) for (T=U1;T<=U2;T=T+Step)
Along = Standard_False; {
S.D1(T,Param,P,D1U,D1V);
D1NormMax=Max(D1NormMax,D1U.Magnitude());
}
if (D1NormMax >TolMax || D1NormMax < TolMin )
Along = Standard_False;
}
} }
else else
{ {
Step = (V2 - V1)/10; if( !Precision::IsInfinite(V1) && !Precision::IsInfinite(V2) )
if(Step < Precision::PConfusion()) {
return Standard_False;
}
if(Step < Precision::PConfusion()) {
return Standard_False;
}
D1NormMax=0.;
for (T=V1;T<=V2;T=T+Step)
{ {
S.D1(Param,T,P,D1U,D1V); Step = (V2 - V1)/10;
D1NormMax=Max(D1NormMax,D1V.Magnitude()); if(Step < Precision::PConfusion()) {
return Standard_False;
}
D1NormMax=0.;
for (T=V1;T<=V2;T=T+Step)
{
S.D1(Param,T,P,D1U,D1V);
D1NormMax=Max(D1NormMax,D1V.Magnitude());
}
if (D1NormMax >TolMax || D1NormMax < TolMin )
Along = Standard_False;
} }
if (D1NormMax >TolMax || D1NormMax < TolMin )
Along = Standard_False;
} }
@ -240,14 +242,7 @@ void Extrema_ExtPS::Perform(const gp_Pnt& P)
myPoints.Clear(); myPoints.Clear();
mySqDist.Clear(); mySqDist.Clear();
Standard_Integer i; Standard_Integer i;
P11 = myS->Value(myuinf, myvinf);
P12 = myS->Value(myuinf, myvsup);
P21 = myS->Value(myusup, myvinf);
P22 = myS->Value(myusup, myvsup);
d11 = P.SquareDistance(P11);
d12 = P.SquareDistance(P12);
d21 = P.SquareDistance(P21);
d22 = P.SquareDistance(P22);
switch(mytype) { switch(mytype) {

13
src/Extrema/Extrema_FuncExtCS.cxx
View File

@ -23,6 +23,7 @@
#include <Extrema_FuncExtCS.ixx> #include <Extrema_FuncExtCS.ixx>
#include <gp_Vec.hxx> #include <gp_Vec.hxx>
#include <Standard_TypeMismatch.hxx> #include <Standard_TypeMismatch.hxx>
#include <Precision.hxx>
/*----------------------------------------------------------------------------- /*-----------------------------------------------------------------------------
Fonction permettant de rechercher une distance extremale entre une courbe C Fonction permettant de rechercher une distance extremale entre une courbe C
@ -209,7 +210,17 @@ Standard_Integer Extrema_FuncExtCS::GetStateNumber()
Value(UVSol, Sol); Value(UVSol, Sol);
cout <<"F(1)= "<<Sol(1)<<" F(2)= "<<Sol(2)<<" F(3)= "<<Sol(3)<<endl; cout <<"F(1)= "<<Sol(1)<<" F(2)= "<<Sol(2)<<" F(3)= "<<Sol(3)<<endl;
#endif #endif
//comparison of solution with previous solutions
Standard_Real tol2d = Precision::PConfusion() * Precision::PConfusion();
Standard_Integer i = 1, nbSol = mySqDist.Length();
for( ; i <= nbSol; i++)
{
Standard_Real aU = myPoint1(i).Parameter();
if( (myU - aU) * (myU - aU) <= tol2d )
break;
}
if (i <= nbSol)
return 0;
mySqDist.Append(myP1.SquareDistance(myP2)); mySqDist.Append(myP1.SquareDistance(myP2));
myPoint1.Append(Extrema_POnCurv(myt,myP1)); myPoint1.Append(Extrema_POnCurv(myt,myP1));
myPoint2.Append(Extrema_POnSurf(myU,myV,myP2)); myPoint2.Append(Extrema_POnSurf(myU,myV,myP2));

13
src/Extrema/Extrema_FuncExtPS.cxx
View File

@ -130,6 +130,19 @@ Standard_Boolean Extrema_FuncExtPS::Values (const math_Vector& UV,
Standard_Integer Extrema_FuncExtPS::GetStateNumber () Standard_Integer Extrema_FuncExtPS::GetStateNumber ()
{ {
if (!myPinit || !mySinit) Standard_TypeMismatch::Raise(); if (!myPinit || !mySinit) Standard_TypeMismatch::Raise();
//comparison of solution with previous solutions
Standard_Integer i = 1, nbSol = mySqDist.Length();
Standard_Real tol2d = Precision::PConfusion() * Precision::PConfusion();
for( ; i <= nbSol; i++)
{
Standard_Real aU, aV;
myPoint(i).Parameter(aU, aV);
if( ((myU - aU ) * (myU - aU ) + (myV - aV ) * (myV - aV )) <= tol2d )
break;
}
if( i <= nbSol)
return 0;
mySqDist.Append(myPs.SquareDistance(myP)); mySqDist.Append(myPs.SquareDistance(myP));
myPoint.Append(Extrema_POnSurf(myU,myV,myPs)); myPoint.Append(Extrema_POnSurf(myU,myV,myPs));
return 0; return 0;

14
src/Extrema/Extrema_GenExtPS.cdl
View File

@ -36,7 +36,8 @@ uses POnSurf from Extrema,
ExtAlgo from Extrema, ExtAlgo from Extrema,
HArray1OfSphere from Bnd, HArray1OfSphere from Bnd,
Vector from math, Vector from math,
HArray2OfPnt from TColgp HArray2OfPnt from TColgp,
HArray1OfReal from TColStd
raises NotDone from StdFail, raises NotDone from StdFail,
OutOfRange from Standard, OutOfRange from Standard,
@ -137,8 +138,15 @@ is
BuildTree(me : in out) BuildTree(me : in out)
is static private; is static private;
FindSolution(me: in out; P : Pnt from gp; UV : Vector from math; PasU, PasV : Real; f : ExtFlag from Extrema) FindSolution(me: in out; P : Pnt from gp; UV : Vector from math; theNoU, theNoV : Integer; f : ExtFlag from Extrema)
is static private; is static private;
GetGridPoints(me: in out; theSurf: Surface from Adaptor3d) is private;
---Purpose: Selection of points to build grid, depending on the type of surface
BuildGrid(me: in out) is private;
---Purpose: Creation of grid of parametric points
fields fields
myDone : Boolean; myDone : Boolean;
@ -158,5 +166,7 @@ fields
myS : SurfacePtr from Adaptor3d; myS : SurfacePtr from Adaptor3d;
myFlag : ExtFlag from Extrema; myFlag : ExtFlag from Extrema;
myAlgo : ExtAlgo from Extrema; myAlgo : ExtAlgo from Extrema;
myUParams : HArray1OfReal from TColStd;
myVParams : HArray1OfReal from TColStd;
end GenExtPS; end GenExtPS;

355
src/Extrema/Extrema_GenExtPS.cxx
View File

@ -37,7 +37,16 @@
#include <Extrema_ExtFlag.hxx> #include <Extrema_ExtFlag.hxx>
#include <Bnd_Array1OfSphere.hxx> #include <Bnd_Array1OfSphere.hxx>
#include <Bnd_HArray1OfSphere.hxx> #include <Bnd_HArray1OfSphere.hxx>
#include <Precision.hxx>
#include <Geom_OffsetSurface.hxx>
#include <Geom_RectangularTrimmedSurface.hxx>
#include <Geom_BSplineSurface.hxx>
#include <Geom_BezierSurface.hxx>
#include <Adaptor3d_HSurface.hxx>
#include <Adaptor3d_HCurve.hxx>
#include <Adaptor3d_Curve.hxx>
#include <Geom_BSplineCurve.hxx>
#include <Geom_BezierCurve.hxx>
//IMPLEMENT_HARRAY1(Extrema_HArray1OfSphere) //IMPLEMENT_HARRAY1(Extrema_HArray1OfSphere)
@ -322,7 +331,6 @@ void Extrema_GenExtPS::Initialize(const Adaptor3d_Surface& S,
const Standard_Real TolU, const Standard_Real TolU,
const Standard_Real TolV) const Standard_Real TolV)
{ {
myInit = Standard_True;
myS = (Adaptor3d_SurfacePtr)&S; myS = (Adaptor3d_SurfacePtr)&S;
myusample = NbU; myusample = NbU;
myvsample = NbV; myvsample = NbV;
@ -339,35 +347,191 @@ void Extrema_GenExtPS::Initialize(const Adaptor3d_Surface& S,
myF.Initialize(S); myF.Initialize(S);
mySphereUBTree.Nullify(); mySphereUBTree.Nullify();
myUParams.Nullify();
myVParams.Nullify();
myInit = Standard_False;
}
if(myAlgo == Extrema_ExtAlgo_Grad) inline static void fillParams(const TColStd_Array1OfReal& theKnots,
Standard_Integer theDegree,
Standard_Real theParMin,
Standard_Real theParMax,
Handle_TColStd_HArray1OfReal& theParams,
Standard_Integer theSample)
{
NCollection_Vector<Standard_Real> aParams;
Standard_Integer i = 1;
Standard_Real aPrevPar = theParMin;
aParams.Append(aPrevPar);
//calculation the array of parametric points depending on the knots array variation and degree of given surface
for ( ; i < theKnots.Length() && theKnots(i) < (theParMax - Precision::PConfusion()); i++ )
{ {
//If flag was changed and extrema not reinitialized Extrema would fail if ( theKnots(i+1) < theParMin + Precision::PConfusion())
mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1); continue;
Standard_Real PasU = myusup - myumin;
Standard_Real PasV = myvsup - myvmin;
Standard_Real U0 = PasU / myusample / 100.;
Standard_Real V0 = PasV / myvsample / 100.;
gp_Pnt P1;
PasU = (PasU - U0) / (myusample - 1);
PasV = (PasV - V0) / (myvsample - 1);
U0 = U0/2. + myumin;
V0 = V0/2. + myvmin;
// Calculation of distances Standard_Real aStep = (theKnots(i+1) - theKnots(i))/Max(theDegree,2);
Standard_Integer k =1;
for( ; k <= theDegree ; k++)
{
Standard_Real aPar = theKnots(i) + k * aStep;
if(aPar > theParMax - Precision::PConfusion())
break;
if(aPar > aPrevPar + Precision::PConfusion() )
{
aParams.Append(aPar);
aPrevPar = aPar;
}
}
}
aParams.Append(theParMax);
Standard_Integer nbPar = aParams.Length();
//in case of an insufficient number of points the grid will be built later
if (nbPar < theSample)
return;
theParams = new TColStd_HArray1OfReal(1, nbPar );
for( i = 0; i < nbPar; i++)
theParams->SetValue(i+1,aParams(i));
}
void Extrema_GenExtPS::GetGridPoints( const Adaptor3d_Surface& theSurf)
{
//creation parametric points for BSpline and Bezier surfaces
//with taking into account of Degree and NbKnots of BSpline or Bezier geometry
if(theSurf.GetType() == GeomAbs_OffsetSurface)
GetGridPoints(theSurf.BasisSurface()->Surface());
//parametric points for BSpline surfaces
else if( theSurf.GetType() == GeomAbs_BSplineSurface)
{
Handle(Geom_BSplineSurface) aBspl = theSurf.BSpline();
if(!aBspl.IsNull())
{
TColStd_Array1OfReal aUKnots(1, aBspl->NbUKnots());
aBspl->UKnots( aUKnots);
TColStd_Array1OfReal aVKnots(1, aBspl->NbVKnots());
aBspl->VKnots( aVKnots);
fillParams(aUKnots,aBspl->UDegree(),myumin, myusup, myUParams, myusample);
fillParams(aVKnots,aBspl->VDegree(),myvmin, myvsup, myVParams, myvsample);
}
}
//calculation parametric points for Bezier surfaces
else if(theSurf.GetType() == GeomAbs_BezierSurface)
{
Handle(Geom_BezierSurface) aBezier = theSurf.Bezier();
if(aBezier.IsNull())
return;
TColStd_Array1OfReal aUKnots(1,2);
TColStd_Array1OfReal aVKnots(1,2);
aBezier->Bounds(aUKnots(1), aUKnots(2), aVKnots(1), aVKnots(2));
fillParams(aUKnots, aBezier->UDegree(), myumin, myusup, myUParams, myusample);
fillParams(aVKnots, aBezier->VDegree(), myvmin, myvsup, myVParams, myvsample);
}
//creation points for surfaces based on BSpline or Bezier curves
else if(theSurf.GetType() == GeomAbs_SurfaceOfRevolution ||
theSurf.GetType() == GeomAbs_SurfaceOfExtrusion)
{
Handle(TColStd_HArray1OfReal) anArrKnots;
Standard_Integer aDegree = 0;
GeomAbs_CurveType aType = theSurf.BasisCurve()->Curve().GetType();
if(theSurf.BasisCurve()->Curve().GetType() == GeomAbs_BSplineCurve)
{
Handle(Geom_BSplineCurve) aBspl = theSurf.BasisCurve()->Curve().BSpline();
if(!aBspl.IsNull())
{
anArrKnots = new TColStd_HArray1OfReal(1,aBspl->NbKnots());
aBspl->Knots( anArrKnots->ChangeArray1() );
aDegree = aBspl->Degree();
}
}
if(theSurf.BasisCurve()->Curve().GetType() == GeomAbs_BezierCurve)
{
Handle(Geom_BezierCurve) aBez = theSurf.BasisCurve()->Curve().Bezier();
if(!aBez.IsNull())
{
anArrKnots = new TColStd_HArray1OfReal(1,2);
anArrKnots->SetValue(1, aBez->FirstParameter());
anArrKnots->SetValue(2, aBez->LastParameter());
aDegree = aBez->Degree();
}
}
if(anArrKnots.IsNull())
return;
if( theSurf.GetType() == GeomAbs_SurfaceOfRevolution )
fillParams( anArrKnots->Array1(), aDegree, myvmin, myvsup, myVParams, myvsample);
else
fillParams( anArrKnots->Array1(), aDegree, myumin, myusup, myUParams, myusample);
}
//update the number of points in sample
if(!myUParams.IsNull())
myusample = myUParams->Length();
if( !myVParams.IsNull())
myvsample = myVParams->Length();
}
void Extrema_GenExtPS::BuildGrid()
{
//if grid was already built do nothing
if(myInit)
return;
//build parametric grid in case of a complex surface geometry (BSpline and Bezier surfaces)
GetGridPoints(*myS);
//build grid in other cases
if( myUParams.IsNull() )
{
Standard_Real PasU = myusup - myumin;
Standard_Real U0 = PasU / myusample / 100.;
PasU = (PasU - U0) / (myusample - 1);
U0 = U0/2. + myumin;
myUParams = new TColStd_HArray1OfReal(1,myusample );
Standard_Integer NoU;
Standard_Real U = U0;
for ( NoU = 1 ; NoU <= myusample; NoU++, U += PasU)
myUParams->SetValue(NoU, U);
}
if( myVParams.IsNull())
{
Standard_Real PasV = myvsup - myvmin;
Standard_Real V0 = PasV / myvsample / 100.;
PasV = (PasV - V0) / (myvsample - 1);
V0 = V0/2. + myvmin;
myVParams = new TColStd_HArray1OfReal(1,myvsample );
Standard_Integer NoV;
Standard_Real V = V0;
for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV)
myVParams->SetValue(NoV, V);
}
//If flag was changed and extrema not reinitialized Extrema would fail
mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1);
// Calculation of distances
Standard_Integer NoU, NoV; Standard_Integer NoU, NoV;
Standard_Real U, V;
for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) { for ( NoU = 1 ; NoU <= myusample; NoU++ ) {
for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) { for ( NoV = 1 ; NoV <= myvsample; NoV++) {
P1 = myS->Value(U, V); gp_Pnt P1 = myS->Value(myUParams->Value(NoU), myVParams->Value(NoV));
mypoints->SetValue(NoU,NoV,P1); mypoints->SetValue(NoU,NoV,P1);
} }
} }
}
myInit = Standard_True;
}
//mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1);
/* /*
a- Constitution of the table of distances (TbDist(0,myusample+1,0,myvsample+1)): a- Constitution of the table of distances (TbDist(0,myusample+1,0,myvsample+1)):
--------------------------------------------------------------- ---------------------------------------------------------------
@ -375,9 +539,6 @@ a- Constitution of the table of distances (TbDist(0,myusample+1,0,myvsample+1)):
// Parameterisation of the sample // Parameterisation of the sample
}
void Extrema_GenExtPS::BuildTree() void Extrema_GenExtPS::BuildTree()
{ {
// if tree already exists, assume it is already correctly filled // if tree already exists, assume it is already correctly filled
@ -394,16 +555,26 @@ void Extrema_GenExtPS::BuildTree()
U0 = U0/2. + myumin; U0 = U0/2. + myumin;
V0 = V0/2. + myvmin; V0 = V0/2. + myvmin;
//build grid of parametric points
myUParams = new TColStd_HArray1OfReal(1,myusample );
myVParams = new TColStd_HArray1OfReal(1,myvsample );
Standard_Integer NoU, NoV;
Standard_Real U = U0, V = V0;
for ( NoU = 1 ; NoU <= myusample; NoU++, U += PasU)
myUParams->SetValue(NoU, U);
for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV)
myVParams->SetValue(NoV, V);
// Calculation of distances // Calculation of distances
mySphereUBTree = new Extrema_UBTreeOfSphere; mySphereUBTree = new Extrema_UBTreeOfSphere;
Extrema_UBTreeFillerOfSphere aFiller(*mySphereUBTree); Extrema_UBTreeFillerOfSphere aFiller(*mySphereUBTree);
Standard_Integer i = 0; Standard_Integer i = 0;
Standard_Real U, V;
mySphereArray = new Bnd_HArray1OfSphere(0, myusample * myvsample); mySphereArray = new Bnd_HArray1OfSphere(0, myusample * myvsample);
Standard_Integer NoU, NoV;
for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) { for ( NoU = 1; NoU <= myusample; NoU++ ) {
for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) { for ( NoV = 1; NoV <= myvsample; NoV++) {
P1 = myS->Value(U, V); P1 = myS->Value(myUParams->Value(NoU), myVParams->Value(NoV));
Bnd_Sphere aSph(P1.XYZ(), 0/*mytolu < mytolv ? mytolu : mytolv*/, NoU, NoV); Bnd_Sphere aSph(P1.XYZ(), 0/*mytolu < mytolv ? mytolu : mytolv*/, NoU, NoV);
aFiller.Add(i, aSph); aFiller.Add(i, aSph);
mySphereArray->SetValue( i, aSph ); mySphereArray->SetValue( i, aSph );
@ -413,7 +584,10 @@ void Extrema_GenExtPS::BuildTree()
aFiller.Fill(); aFiller.Fill();
} }
void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, const Standard_Real PasU, const Standard_Real PasV, const Extrema_ExtFlag f) void Extrema_GenExtPS::FindSolution(const gp_Pnt& P,
const math_Vector& UV,
const Standard_Integer theNoU,
const Standard_Integer theNoV, const Extrema_ExtFlag f)
{ {
math_Vector Tol(1,2); math_Vector Tol(1,2);
Tol(1) = mytolu; Tol(1) = mytolu;
@ -438,6 +612,7 @@ void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, cons
math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter); math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
Standard_Boolean ToResolveOnSubgrid = Standard_False; Standard_Boolean ToResolveOnSubgrid = Standard_False;
Standard_Boolean NewSolution = Standard_False;
if (f == Extrema_ExtFlag_MIN) if (f == Extrema_ExtFlag_MIN)
{ {
if(S.IsDone()) if(S.IsDone())
@ -447,46 +622,46 @@ void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, cons
gp_Pnt PSol = myS->Value(root(1), root(2)); gp_Pnt PSol = myS->Value(root(1), root(2));
DistSol = P.SquareDistance(PSol); DistSol = P.SquareDistance(PSol);
if(Abs(errors(1)) > eps || Abs(errors(2)) > eps || DistStart < DistSol) if(Abs(errors(1)) > eps || Abs(errors(2)) > eps || DistStart < DistSol)
{
//try to improve solution on subgrid of sample points //try to improve solution on subgrid of sample points
ToResolveOnSubgrid = Standard_True; ToResolveOnSubgrid = Standard_True;
}
} }
else else
{
//keep found roots and try to find solution
Standard_Integer nbExt = myF.NbExt();
Standard_Integer k = 1;
for( ; k <= nbExt; k++)
{
Standard_Real aD = myF.SquareDistance(k);
if(aD < DistSol)
{
DistSol = aD;
myF.Point(k).Parameter(UV(1),UV(2));
NewSolution = Standard_True;
}
}
ToResolveOnSubgrid = Standard_True; ToResolveOnSubgrid = Standard_True;
}
if (ToResolveOnSubgrid) if (ToResolveOnSubgrid)
{ {
Standard_Real u1 = Max(UV(1) - PasU, myumin), u2 = Min(UV(1) + PasU, myusup); //extension of interval to find new solution
Standard_Real v1 = Max(UV(2) - PasV, myvmin), v2 = Min(UV(2) + PasV, myvsup); Standard_Real u1 = theNoU == 1 ? myumin : myUParams->Value(theNoU-1),
u2 = theNoU == myusample ? myusup : myUParams->Value(theNoU + 1);
if(u2 - u1 < 2.*PasU) { Standard_Real v1 = theNoV == 1 ? myvmin : myVParams->Value(theNoV-1),
if(Abs(u1 - myumin) < 1.e-9) { v2 = theNoV == myvsample ? myvsup : myVParams->Value(theNoV + 1);
u2 = u1 + 2.*PasU;
u2 = Min(u2, myusup);
}
if(Abs(u2 - myusup) < 1.e-9) {
u1 = u2 - 2.*PasU;
u1 = Max(u1, myumin);
}
}
if(v2 - v1 < 2.*PasV) {
if(Abs(v1 - myvmin) < 1.e-9) {
v2 = v1 + 2.*PasV;
v2 = Min(v2, myvsup);
}
if(Abs(v2 - myvsup) < 1.e-9) {
v1 = v2 - 2.*PasV;
v1 = Max(v1, myvmin);
}
}
Standard_Real du = (u2 - u1)/(nbsubsample-1); Standard_Real du = (u2 - u1)/(nbsubsample-1);
Standard_Real dv = (v2 - v1)/(nbsubsample-1); Standard_Real dv = (v2 - v1)/(nbsubsample-1);
Standard_Real u, v; Standard_Real u, v;
Standard_Real dist; Standard_Real dist;
Standard_Boolean NewSolution = Standard_False; //try to find solution on subgrid
Standard_Integer Nu, Nv; Standard_Integer Nu, Nv;
Standard_Integer minU = 0;
Standard_Integer minV = 0;
for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) { for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) {
for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) { for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) {
gp_Pnt Puv = myS->Value(u, v); gp_Pnt Puv = myS->Value(u, v);
@ -497,13 +672,20 @@ void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, cons
UV(2) = v; UV(2) = v;
NewSolution = Standard_True; NewSolution = Standard_True;
DistSol = dist; DistSol = dist;
minU = Nu;
minV = Nv;
} }
} }
} }
if(NewSolution) { if(NewSolution) {
//try to precise //try to precise
UVinf(1) = u1 + (minU == 1 ? 0 : minU - 2) * du;
UVinf(2) = v1 + (minV == 1 ? 0 : minV - 2) * dv;
UVsup(1) = u1 + (minU == nbsubsample ? nbsubsample : minU +1) * du;;
UVsup(2) = v1 + (minV == nbsubsample ? nbsubsample : minV +1) * dv;
math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter); math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
} }
} //end of if (ToResolveOnSubgrid) } //end of if (ToResolveOnSubgrid)
} //end of if (f == Extrema_ExtFlag_MIN) } //end of if (f == Extrema_ExtFlag_MIN)
@ -518,40 +700,33 @@ void Extrema_GenExtPS::SetFlag(const Extrema_ExtFlag F)
void Extrema_GenExtPS::SetAlgo(const Extrema_ExtAlgo A) void Extrema_GenExtPS::SetAlgo(const Extrema_ExtAlgo A)
{ {
if(myAlgo != A)
myInit = Standard_False;
myAlgo = A; myAlgo = A;
} }
void Extrema_GenExtPS::Perform(const gp_Pnt& P) void Extrema_GenExtPS::Perform(const gp_Pnt& P)
{ {
//BuildTree();
myDone = Standard_False; myDone = Standard_False;
myF.SetPoint(P); myF.SetPoint(P);
Standard_Real PasU = myusup - myumin;
Standard_Real PasV = myvsup - myvmin;
Standard_Real U, U0 = PasU / myusample / 100.;
Standard_Real V, V0 = PasV / myvsample / 100.;
PasU = (PasU - U0) / (myusample - 1);
PasV = (PasV - V0) / (myvsample - 1);
U0 = U0/2.+myumin;
V0 = V0/2.+myvmin;
//math_Vector Tol(1, 2);
math_Vector UV(1,2); math_Vector UV(1,2);
if(myAlgo == Extrema_ExtAlgo_Grad) if(myAlgo == Extrema_ExtAlgo_Grad)
{ {
BuildGrid();
Standard_Integer NoU,NoV; Standard_Integer NoU,NoV;
TColStd_Array2OfReal TheDist(0, myusample+1, 0, myvsample+1); TColStd_Array2OfReal TheDist(0, myusample+1, 0, myvsample+1);
for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) { for ( NoU = 1 ; NoU <= myusample; NoU++) {
for ( NoV = 1; NoV <= myvsample; NoV++) {
TheDist(NoU, NoV) = P.SquareDistance(mypoints->Value(NoU, NoV)); TheDist(NoU, NoV) = P.SquareDistance(mypoints->Value(NoU, NoV));
} }
} }
TColStd_Array2OfInteger TbSel(0,myusample+1,0,myvsample+1);
TbSel.Init(0);
Standard_Real Dist; Standard_Real Dist;
@ -567,7 +742,6 @@ void Extrema_GenExtPS::Perform(const gp_Pnt& P)
} }
for (NoU = 1; NoU <= myusample; NoU++) { for (NoU = 1; NoU <= myusample; NoU++) {
for (NoV = 1; NoV <= myvsample; NoV++) { for (NoV = 1; NoV <= myvsample; NoV++) {
if (TbSel(NoU,NoV) == 0) {
Dist = TheDist(NoU,NoV); Dist = TheDist(NoU,NoV);
if ((TheDist(NoU-1,NoV-1) >= Dist) && if ((TheDist(NoU-1,NoV-1) >= Dist) &&
(TheDist(NoU-1,NoV ) >= Dist) && (TheDist(NoU-1,NoV ) >= Dist) &&
@ -578,11 +752,12 @@ void Extrema_GenExtPS::Perform(const gp_Pnt& P)
(TheDist(NoU+1,NoV ) >= Dist) && (TheDist(NoU+1,NoV ) >= Dist) &&
(TheDist(NoU+1,NoV+1) >= Dist)) { (TheDist(NoU+1,NoV+1) >= Dist)) {
//Create array of UV vectors to calculate min //Create array of UV vectors to calculate min
UV(1) = U0 + (NoU - 1) * PasU;
UV(2) = V0 + (NoV - 1) * PasV; UV(1) = myUParams->Value(NoU);
FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN); UV(2) = myVParams->Value(NoV);
FindSolution(P, UV, NoU, NoV, Extrema_ExtFlag_MIN);
} }
}
} }
} }
} }
@ -599,7 +774,6 @@ void Extrema_GenExtPS::Perform(const gp_Pnt& P)
} }
for (NoU = 1; NoU <= myusample; NoU++) { for (NoU = 1; NoU <= myusample; NoU++) {
for (NoV = 1; NoV <= myvsample; NoV++) { for (NoV = 1; NoV <= myvsample; NoV++) {
if (TbSel(NoU,NoV) == 0) {
Dist = TheDist(NoU,NoV); Dist = TheDist(NoU,NoV);
if ((TheDist(NoU-1,NoV-1) <= Dist) && if ((TheDist(NoU-1,NoV-1) <= Dist) &&
(TheDist(NoU-1,NoV ) <= Dist) && (TheDist(NoU-1,NoV ) <= Dist) &&
@ -610,12 +784,11 @@ void Extrema_GenExtPS::Perform(const gp_Pnt& P)
(TheDist(NoU+1,NoV ) <= Dist) && (TheDist(NoU+1,NoV ) <= Dist) &&
(TheDist(NoU+1,NoV+1) <= Dist)) { (TheDist(NoU+1,NoV+1) <= Dist)) {
//Create array of UV vectors to calculate max //Create array of UV vectors to calculate max
UV(1) = U0 + (NoU - 1) * PasU; UV(1) = myUParams->Value(NoU);
UV(2) = V0 + (NoV - 1) * PasV; UV(2) = myVParams->Value(NoV);
FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX); FindSolution(P, UV, NoU, NoV, Extrema_ExtFlag_MAX);
} }
} }
}
} }
} }
} }
@ -632,10 +805,11 @@ void Extrema_GenExtPS::Perform(const gp_Pnt& P)
//TODO: check if no solution in binary tree //TODO: check if no solution in binary tree
Bnd_Sphere& aSph = aSelector.Sphere(); Bnd_Sphere& aSph = aSelector.Sphere();
UV(1) = U0 + (aSph.U() - 1) * PasU; UV(1) = myUParams->Value(aSph.U());//U0 + (aSph.U() - 1) * PasU;
UV(2) = V0 + (aSph.V() - 1) * PasV; UV(2) = myVParams->Value(aSph.V());//V0 + (aSph.V() - 1) * PasV;
FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN); //FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
FindSolution(P, UV, aSph.U(), aSph.V(), Extrema_ExtFlag_MIN);
} }
if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX) if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
{ {
@ -646,11 +820,12 @@ void Extrema_GenExtPS::Perform(const gp_Pnt& P)
Standard_Integer aNbSel = mySphereUBTree->Select( aSelector ); Standard_Integer aNbSel = mySphereUBTree->Select( aSelector );
//TODO: check if no solution in binary tree //TODO: check if no solution in binary tree
Bnd_Sphere& aSph = aSelector.Sphere(); Bnd_Sphere& aSph = aSelector.Sphere();
UV(1) = myUParams->Value(aSph.U());
UV(2) = myVParams->Value(aSph.V());
//UV(1) = U0 + (aSph.U() - 1) * PasU;
//UV(2) = V0 + (aSph.V() - 1) * PasV;
UV(1) = U0 + (aSph.U() - 1) * PasU; FindSolution(P, UV, aSph.U(),aSph.V(), Extrema_ExtFlag_MAX);
UV(2) = V0 + (aSph.V() - 1) * PasV;
FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
} }
} }
} }

20
src/math/math_FunctionSetRoot.cxx
View File

@ -747,7 +747,7 @@ void math_FunctionSetRoot::Perform(math_FunctionSetWithDerivatives& F,
// s'il on est dejas sur la solution, il faut leurer ce test pour eviter // s'il on est dejas sur la solution, il faut leurer ce test pour eviter
// de faire une seconde iteration... // de faire une seconde iteration...
Save(0) = Max (F2, EpsSqrt); Save(0) = Max (F2, EpsSqrt);
Standard_Real aTol_Func = Epsilon(F2);
FSR_DEBUG("=== Mode Debug de Function Set Root" << endl) FSR_DEBUG("=== Mode Debug de Function Set Root" << endl)
FSR_DEBUG(" F2 Initial = " << F2) FSR_DEBUG(" F2 Initial = " << F2)
@ -1038,12 +1038,18 @@ void math_FunctionSetRoot::Perform(math_FunctionSetWithDerivatives& F,
//fin du test solution //fin du test solution
// Analyse de la progression... // Analyse de la progression...
if (F2 < PreviousMinimum) { //comparison of current minimum and previous minimum
if ((F2 - PreviousMinimum) <= aTol_Func){
if (Kount > 5) { if (Kount > 5) {
// L'historique est il bon ? // L'historique est il bon ?
if (F2 >= 0.95*Save(Kount - 5)) { if (F2 >= 0.95*Save(Kount - 5)) {
if (!ChangeDirection) ChangeDirection = Standard_True; if (!ChangeDirection) ChangeDirection = Standard_True;
else return; // si un gain inf a 5% on sort else
{
Done = Standard_True;
State = F.GetStateNumber();
return; // si un gain inf a 5% on sort
}
} }
else ChangeDirection = Standard_False; //Si oui on recommence else ChangeDirection = Standard_False; //Si oui on recommence
} }
@ -1084,9 +1090,15 @@ void math_FunctionSetRoot::Perform(math_FunctionSetWithDerivatives& F,
return; return;
} }
} }
else return; // y a plus d'issues else
{
State = F.GetStateNumber();
return; // y a plus d'issues
}
} }
} }
State = F.GetStateNumber();
} }