-// File: Extrema_GenExtPS.cxx
-// Created: Tue Jul 18 08:21:34 1995
-// Author: Modelistation
-// <model@metrox>
+// Created on: 1995-07-18
+// Created by: Modelistation
+// Copyright (c) 1995-1999 Matra Datavision
+// Copyright (c) 1999-2012 OPEN CASCADE SAS
+//
+// The content of this file is subject to the Open CASCADE Technology Public
+// License Version 6.5 (the "License"). You may not use the content of this file
+// except in compliance with the License. Please obtain a copy of the License
+// at http://www.opencascade.org and read it completely before using this file.
+//
+// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
+// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
+//
+// The Original Code and all software distributed under the License is
+// distributed on an "AS IS" basis, without warranty of any kind, and the
+// Initial Developer hereby disclaims all such warranties, including without
+// limitation, any warranties of merchantability, fitness for a particular
+// purpose or non-infringement. Please see the License for the specific terms
+// and conditions governing the rights and limitations under the License.
+
// Modified by skv - Thu Sep 30 15:21:07 2004 OCC593
#include <Extrema_ExtFlag.hxx>
#include <Bnd_Array1OfSphere.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)
return Standard_False;
}
+/*
+ * This function computes the point on surface parameters on edge.
+ * if it coincides with theParam0 or theParam1, it is returned.
+ */
+static Extrema_POnSurfParams ComputeEdgeParameters
+ (const Standard_Boolean IsUEdge,
+ const Extrema_POnSurfParams &theParam0,
+ const Extrema_POnSurfParams &theParam1,
+ const Adaptor3d_SurfacePtr &theSurf,
+ const gp_Pnt &thePoint,
+ const Standard_Real theDiffTol)
+{
+ const Standard_Real aSqrDist01 =
+ theParam0.Value().SquareDistance(theParam1.Value());
+
+ if (aSqrDist01 <= theDiffTol) {
+ // The points are confused. Get the first point and change its type.
+ return theParam0;
+ } else {
+ const Standard_Real aDiffDist =
+ Abs(theParam0.GetSqrDistance() - theParam1.GetSqrDistance());
+
+ if (aDiffDist >= aSqrDist01 - theDiffTol) {
+ // The shortest distance is one of the nodes.
+ if (theParam0.GetSqrDistance() > theParam1.GetSqrDistance()) {
+ // The shortest distance is the point 1.
+ return theParam1;
+ } else {
+ // The shortest distance is the point 0.
+ return theParam0;
+ }
+ } else {
+ // The shortest distance is inside the edge.
+ gp_XYZ aPoP(thePoint.XYZ().Subtracted(theParam0.Value().XYZ()));
+ gp_XYZ aPoP1(theParam1.Value().XYZ().Subtracted(theParam0.Value().XYZ()));
+ Standard_Real aRatio = aPoP.Dot(aPoP1)/aSqrDist01;
+ Standard_Real aU[2];
+ Standard_Real aV[2];
+
+ theParam0.Parameter(aU[0], aV[0]);
+ theParam1.Parameter(aU[1], aV[1]);
+
+ Standard_Real aUPar = aU[0];
+ Standard_Real aVPar = aV[0];
+ if (IsUEdge) {
+ aUPar += aRatio*(aU[1] - aU[0]);
+ } else {
+ aVPar += aRatio*(aV[1] - aV[0]);
+ }
+
+ Extrema_POnSurfParams aParam(aUPar, aVPar, theSurf->Value(aUPar, aVPar));
+ Standard_Integer anIndices[2];
+
+ theParam0.GetIndices(anIndices[0], anIndices[1]);
+ aParam.SetElementType(IsUEdge ? Extrema_UIsoEdge : Extrema_VIsoEdge);
+ aParam.SetSqrDistance(thePoint.SquareDistance(aParam.Value()));
+ aParam.SetIndices(anIndices[0], anIndices[1]);
+
+ return aParam;
+ }
+ }
+}
//=============================================================================
S using sampling (NbU,NbV).
Method:
- The algorithm bases on the hypothesis that sampling is precise enough
- pour que, s'il existe N distances extremales entre le point et la surface,
- alors il existe aussi N extrema entre le point et la grille.
+ The algorithm bases on the hypothesis that sampling is precise enough,
+ if there exist N extreme distances between the point and the surface,
+ so there also exist N extrema between the point and the grid.
So, the algorithm consists in starting from extrema of the grid to find the
extrema of the surface.
The extrema are calculated by the algorithm math_FunctionSetRoot with the
- update table TbSel.
-----------------------------------------------------------------------------*/
-static Standard_Boolean IsoIsDeg (const Adaptor3d_Surface& S,
- const Standard_Real Param,
- const GeomAbs_IsoType IT,
- const Standard_Real TolMin,
- const Standard_Real TolMax)
-{
- Standard_Real U1=0.,U2=0.,V1=0.,V2=0.,T;
- Standard_Boolean Along = Standard_True;
- U1 = S.FirstUParameter();
- U2 = S.LastUParameter();
- V1 = S.FirstVParameter();
- V2 = S.LastVParameter();
- gp_Vec D1U,D1V;
- gp_Pnt P;
- Standard_Real Step,D1NormMax;
- if (IT == GeomAbs_IsoV)
- {
- Step = (U2 - U1)/10;
- D1NormMax=0.;
- for (T=U1;T<=U2;T=T+Step)
- {
- S.D1(T,Param,P,D1U,D1V);
- D1NormMax=Max(D1NormMax,D1U.Magnitude());
- }
-
- if (D1NormMax >TolMax || D1NormMax < TolMin )
- Along = Standard_False;
- }
- else
- {
- Step = (V2 - V1)/10;
- 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;
-
-
- }
- return Along;
-}
-//----------------------------------------------------------
Extrema_GenExtPS::Extrema_GenExtPS()
{
myDone = Standard_False;
const Standard_Real TolU,
const Standard_Real TolV)
{
- myInit = Standard_True;
myS = (Adaptor3d_SurfacePtr)&S;
myusample = NbU;
myvsample = NbV;
myF.Initialize(S);
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
- mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1);
- 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;
+ if ( theKnots(i+1) < theParMin + Precision::PConfusion())
+ continue;
-// Calcul des 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;
+ 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(const gp_Pnt &thePoint)
+{
Standard_Integer NoU, NoV;
- Standard_Real U, V;
- for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
- for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
- P1 = myS->Value(U, V);
- mypoints->SetValue(NoU,NoV,P1);
+
+ //if grid was already built skip its creation
+ if (!myInit) {
+ //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 Extrema_HArray2OfPOnSurfParams
+ (0, myusample + 1, 0, myvsample + 1);
+ // Calculation of distances
+
+ for ( NoU = 1 ; NoU <= myusample; NoU++ ) {
+ for ( NoV = 1 ; NoV <= myvsample; NoV++) {
+ gp_Pnt aP1 = myS->Value(myUParams->Value(NoU), myVParams->Value(NoV));
+ Extrema_POnSurfParams aParam
+ (myUParams->Value(NoU), myVParams->Value(NoV), aP1);
+
+ aParam.SetElementType(Extrema_Node);
+ aParam.SetIndices(NoU, NoV);
+ myPoints->SetValue(NoU, NoV, aParam);
+ }
+ }
+
+ // Fill boundary with negative square distance.
+ // It is used for computation of Maximum.
+ for (NoV = 0; NoV <= myvsample + 1; NoV++) {
+ myPoints->ChangeValue(0, NoV).SetSqrDistance(-1.);
+ myPoints->ChangeValue(myusample + 1, NoV).SetSqrDistance(-1.);
+ }
+
+ for (NoU = 1; NoU <= myusample; NoU++) {
+ myPoints->ChangeValue(NoU, 0).SetSqrDistance(-1.);
+ myPoints->ChangeValue(NoU, myvsample + 1).SetSqrDistance(-1.);
+ }
+
+ myInit = Standard_True;
}
+
+ // Compute distances to mesh.
+ // Step 1. Compute distances to nodes.
+ for ( NoU = 1 ; NoU <= myusample; NoU++ ) {
+ for ( NoV = 1 ; NoV <= myvsample; NoV++) {
+ Extrema_POnSurfParams &aParam = myPoints->ChangeValue(NoU, NoV);
+
+ aParam.SetSqrDistance(thePoint.SquareDistance(aParam.Value()));
+ }
}
- //mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1);
+ // For search of minimum compute distances to mesh.
+ if(myFlag == Extrema_ExtFlag_MIN || myFlag == Extrema_ExtFlag_MINMAX)
+ {
+ // This is the tolerance of difference of squared values.
+ // No need to set it too small.
+ const Standard_Real aDiffTol = mytolu + mytolv;
+
+ // Step 2. Compute distances to edges.
+ // Assume UEdge(i, j) = { Point(i, j); Point(i + 1, j ) }
+ // Assume VEdge(i, j) = { Point(i, j); Point(i, j + 1) }
+ Handle(Extrema_HArray2OfPOnSurfParams) aUEdgePntParams =
+ new Extrema_HArray2OfPOnSurfParams(1, myusample - 1, 1, myvsample);
+ Handle(Extrema_HArray2OfPOnSurfParams) aVEdgePntParams =
+ new Extrema_HArray2OfPOnSurfParams(1, myusample, 1, myvsample - 1);
+
+ for ( NoU = 1 ; NoU <= myusample; NoU++ ) {
+ for ( NoV = 1 ; NoV <= myvsample; NoV++) {
+ const Extrema_POnSurfParams &aParam0 = myPoints->Value(NoU, NoV);
+
+ if (NoU < myusample) {
+ // Compute parameters to UEdge.
+ const Extrema_POnSurfParams &aParam1 = myPoints->Value(NoU + 1, NoV);
+ Extrema_POnSurfParams aUEdgeParam = ComputeEdgeParameters
+ (Standard_True, aParam0, aParam1, myS, thePoint, aDiffTol);
+
+ aUEdgePntParams->SetValue(NoU, NoV, aUEdgeParam);
+ }
-/*
-a- Constitution du tableau des distances (TbDist(0,myusample+1,0,myvsample+1)):
- ---------------------------------------------------------------
-*/
+ if (NoV < myvsample) {
+ // Compute parameters to VEdge.
+ const Extrema_POnSurfParams &aParam1 = myPoints->Value(NoU, NoV + 1);
+ Extrema_POnSurfParams aVEdgeParam = ComputeEdgeParameters
+ (Standard_False, aParam0, aParam1, myS, thePoint, aDiffTol);
-// Parametrage de l echantillon
+ aVEdgePntParams->SetValue(NoU, NoV, aVEdgeParam);
+ }
+ }
+ }
+
+ // Step 3. Compute distances to faces.
+ // Assume myFacePntParams(i, j) =
+ // { Point(i, j); Point(i + 1, j); Point(i + 1, j + 1); Point(i, j + 1) }
+ // Or
+ // { UEdge(i, j); VEdge(i + 1, j); UEdge(i, j + 1); VEdge(i, j) }
+ myFacePntParams =
+ new Extrema_HArray2OfPOnSurfParams(0, myusample, 0, myvsample);
+ Standard_Real aSqrDist01;
+ Standard_Real aDiffDist;
+ Standard_Boolean isOut;
+
+ for ( NoU = 1 ; NoU < myusample; NoU++ ) {
+ for ( NoV = 1 ; NoV < myvsample; NoV++) {
+ const Extrema_POnSurfParams &aUE0 = aUEdgePntParams->Value(NoU, NoV);
+ const Extrema_POnSurfParams &aUE1 = aUEdgePntParams->Value(NoU, NoV+1);
+ const Extrema_POnSurfParams &aVE0 = aVEdgePntParams->Value(NoU, NoV);
+ const Extrema_POnSurfParams &aVE1 = aVEdgePntParams->Value(NoU+1, NoV);
+
+ aSqrDist01 = aUE0.Value().SquareDistance(aUE1.Value());
+ aDiffDist = Abs(aUE0.GetSqrDistance() - aUE1.GetSqrDistance());
+ isOut = Standard_False;
+
+ if (aDiffDist >= aSqrDist01 - aDiffTol) {
+ // The projection is outside the face.
+ isOut = Standard_True;
+ } else {
+ aSqrDist01 = aVE0.Value().SquareDistance(aVE1.Value());
+ aDiffDist = Abs(aVE0.GetSqrDistance() - aVE1.GetSqrDistance());
+
+ if (aDiffDist >= aSqrDist01 - aDiffTol) {
+ // The projection is outside the face.
+ isOut = Standard_True;
+ }
+ }
+ if (isOut) {
+ // Get the closest point on an edge.
+ const Extrema_POnSurfParams &aUEMin =
+ aUE0.GetSqrDistance() < aUE1.GetSqrDistance() ? aUE0 : aUE1;
+ const Extrema_POnSurfParams &aVEMin =
+ aVE0.GetSqrDistance() < aVE1.GetSqrDistance() ? aVE0 : aVE1;
+ const Extrema_POnSurfParams &aEMin =
+ aUEMin.GetSqrDistance() < aVEMin.GetSqrDistance() ? aUEMin : aVEMin;
+
+ myFacePntParams->SetValue(NoU, NoV, aEMin);
+ } else {
+ // Find closest point inside the face.
+ Standard_Real aU[2];
+ Standard_Real aV[2];
+ Standard_Real aUPar;
+ Standard_Real aVPar;
+
+ // Compute U parameter.
+ aUE0.Parameter(aU[0], aV[0]);
+ aUE1.Parameter(aU[1], aV[1]);
+ aUPar = 0.5*(aU[0] + aU[1]);
+ // Compute V parameter.
+ aVE0.Parameter(aU[0], aV[0]);
+ aVE1.Parameter(aU[1], aV[1]);
+ aVPar = 0.5*(aV[0] + aV[1]);
+
+ Extrema_POnSurfParams aParam(aUPar, aVPar, myS->Value(aUPar, aVPar));
+
+ aParam.SetElementType(Extrema_Face);
+ aParam.SetSqrDistance(thePoint.SquareDistance(aParam.Value()));
+ aParam.SetIndices(NoU, NoV);
+ myFacePntParams->SetValue(NoU, NoV, aParam);
+ }
+ }
+ }
+
+ // Fill boundary with RealLast square distance.
+ for (NoV = 0; NoV <= myvsample; NoV++) {
+ myFacePntParams->ChangeValue(0, NoV).SetSqrDistance(RealLast());
+ myFacePntParams->ChangeValue(myusample, NoV).SetSqrDistance(RealLast());
+ }
+ for (NoU = 1; NoU < myusample; NoU++) {
+ myFacePntParams->ChangeValue(NoU, 0).SetSqrDistance(RealLast());
+ myFacePntParams->ChangeValue(NoU, myvsample).SetSqrDistance(RealLast());
+ }
+ }
}
+/*
+a- Constitution of the table of distances (TbDist(0,myusample+1,0,myvsample+1)):
+ ---------------------------------------------------------------
+*/
+
+// Parameterisation of the sample
+
void Extrema_GenExtPS::BuildTree()
{
// if tree already exists, assume it is already correctly filled
U0 = U0/2. + myumin;
V0 = V0/2. + myvmin;
- // Calcul des distances
+ //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
mySphereUBTree = new Extrema_UBTreeOfSphere;
Extrema_UBTreeFillerOfSphere aFiller(*mySphereUBTree);
Standard_Integer i = 0;
- Standard_Real U, V;
+
mySphereArray = new Bnd_HArray1OfSphere(0, myusample * myvsample);
- Standard_Integer NoU, NoV;
- for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
- for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
- P1 = myS->Value(U, V);
+
+ for ( NoU = 1; NoU <= myusample; NoU++ ) {
+ for ( NoV = 1; NoV <= myvsample; NoV++) {
+ P1 = myS->Value(myUParams->Value(NoU), myVParams->Value(NoV));
Bnd_Sphere aSph(P1.XYZ(), 0/*mytolu < mytolv ? mytolu : mytolv*/, NoU, NoV);
aFiller.Add(i, aSph);
mySphereArray->SetValue( i, aSph );
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 Extrema_POnSurfParams &theParams)
{
math_Vector Tol(1,2);
Tol(1) = mytolu;
Tol(2) = mytolv;
+ math_Vector UV(1, 2);
+
+ theParams.Parameter(UV(1), UV(2));
+
math_Vector UVinf(1,2), UVsup(1,2);
UVinf(1) = myumin;
UVinf(2) = myvmin;
math_Vector errors(1,2);
math_Vector root(1, 2);
- Standard_Real eps = 1.e-9;
- Standard_Integer nbsubsample = 11;
Standard_Integer aNbMaxIter = 100;
- if (myF.HasDegIso())
- aNbMaxIter = 150;
-
+ gp_Pnt PStart = theParams.Value();
+
math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
- if(S.IsDone()) {
- root = S.Root();
- myF.Value(root, errors);
- if(f == Extrema_ExtFlag_MIN)
- {
- if(Abs(errors(1)) > eps || Abs(errors(2)) > eps) {
- //try to improve solution on subgrid of sample points
- gp_Pnt PSol = myS->Value(root(1), root(2));
- Standard_Real DistSol = P.SquareDistance(PSol);
-
- Standard_Real u1 = Max(UV(1) - PasU, myumin), u2 = Min(UV(1) + PasU, myusup);
- Standard_Real v1 = Max(UV(2) - PasV, myvmin), v2 = Min(UV(2) + PasV, myvsup);
-
- if(u2 - u1 < 2.*PasU) {
- if(Abs(u1 - myumin) < 1.e-9) {
- 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 dv = (v2 - v1)/(nbsubsample-1);
- Standard_Real u, v;
- Standard_Real dist;
-
- Standard_Boolean NewSolution = Standard_False;
- Standard_Integer Nu, Nv;
- for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) {
- for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) {
- gp_Pnt Puv = myS->Value(u, v);
- dist = P.SquareDistance(Puv);
-
- if(dist < DistSol) {
- UV(1) = u;
- UV(2) = v;
- NewSolution = Standard_True;
- DistSol = dist;
- }
- }
- }
-
- if(NewSolution) {
- //try to precise
- math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
- }
-
- }
- }
- }
+
myDone = Standard_True;
}
void Extrema_GenExtPS::SetAlgo(const Extrema_ExtAlgo A)
{
+ if(myAlgo != A)
+ myInit = Standard_False;
myAlgo = A;
+
}
void Extrema_GenExtPS::Perform(const gp_Pnt& P)
{
- //BuildTree();
myDone = Standard_False;
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);
-
+
if(myAlgo == Extrema_ExtAlgo_Grad)
{
+ BuildGrid(P);
Standard_Integer NoU,NoV;
- 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) {
- TheDist(NoU, NoV) = P.SquareDistance(mypoints->Value(NoU, NoV));
- }
- }
-
- TColStd_Array2OfInteger TbSel(0,myusample+1,0,myvsample+1);
- TbSel.Init(0);
-
- Standard_Real Dist;
-
if(myFlag == Extrema_ExtFlag_MIN || myFlag == Extrema_ExtFlag_MINMAX)
{
- for (NoV = 0; NoV <= myvsample+1; NoV++) {
- TheDist(0,NoV) = RealLast();
- TheDist(myusample+1,NoV) = RealLast();
- }
- for (NoU = 1; NoU <= myusample; NoU++) {
- TheDist(NoU,0) = RealLast();
- TheDist(NoU,myvsample+1) = RealLast();
- }
- for (NoU = 1; NoU <= myusample; NoU++) {
- for (NoV = 1; NoV <= myvsample; NoV++) {
- if (TbSel(NoU,NoV) == 0) {
- Dist = TheDist(NoU,NoV);
- if ((TheDist(NoU-1,NoV-1) >= Dist) &&
- (TheDist(NoU-1,NoV ) >= Dist) &&
- (TheDist(NoU-1,NoV+1) >= Dist) &&
- (TheDist(NoU ,NoV-1) >= Dist) &&
- (TheDist(NoU ,NoV+1) >= Dist) &&
- (TheDist(NoU+1,NoV-1) >= Dist) &&
- (TheDist(NoU+1,NoV ) >= Dist) &&
- (TheDist(NoU+1,NoV+1) >= Dist)) {
- //Create array of UV vectors to calculate min
- UV(1) = U0 + (NoU - 1) * PasU;
- UV(2) = V0 + (NoV - 1) * PasV;
- FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
+ Extrema_ElementType anElemType;
+ Standard_Integer iU;
+ Standard_Integer iV;
+ Standard_Integer iU2;
+ Standard_Integer iV2;
+ Standard_Boolean isMin;
+ Standard_Integer i;
+
+ for (NoU = 1; NoU < myusample; NoU++) {
+ for (NoV = 1; NoV < myvsample; NoV++) {
+ const Extrema_POnSurfParams &aParam =
+ myFacePntParams->Value(NoU, NoV);
+
+ isMin = Standard_False;
+ anElemType = aParam.GetElementType();
+
+ if (anElemType == Extrema_Face) {
+ isMin = Standard_True;
+ } else {
+ // Check if it is a boundary edge or corner vertex.
+ aParam.GetIndices(iU, iV);
+
+ if (anElemType == Extrema_UIsoEdge) {
+ isMin = (iV == 1 || iV == myvsample);
+ } else if (anElemType == Extrema_VIsoEdge) {
+ isMin = (iU == 1 || iU == myusample);
+ } else if (anElemType == Extrema_Node) {
+ isMin = (iU == 1 || iU == myusample) &&
+ (iV == 1 || iV == myvsample);
}
+
+ if (!isMin) {
+ // This is a middle element.
+ if (anElemType == Extrema_UIsoEdge ||
+ (anElemType == Extrema_Node && (iU == 1 || iU == myusample))) {
+ // Check the down face.
+ const Extrema_POnSurfParams &aDownParam =
+ myFacePntParams->Value(NoU, NoV - 1);
+
+ if (aDownParam.GetElementType() == anElemType) {
+ aDownParam.GetIndices(iU2, iV2);
+ isMin = (iU == iU2 && iV == iV2);
+ }
+ } else if (anElemType == Extrema_VIsoEdge ||
+ (anElemType == Extrema_Node && (iV == 1 || iV == myvsample))) {
+ // Check the right face.
+ const Extrema_POnSurfParams &aRightParam =
+ myFacePntParams->Value(NoU - 1, NoV);
+
+ if (aRightParam.GetElementType() == anElemType) {
+ aRightParam.GetIndices(iU2, iV2);
+ isMin = (iU == iU2 && iV == iV2);
+ }
+ } else if (iU == NoU && iV == NoV) {
+ // Check the lower-left node. For this purpose it is necessary
+ // to check lower-left, lower and left faces.
+ isMin = Standard_True;
+
+ const Extrema_POnSurfParams *anOtherParam[3] =
+ { &myFacePntParams->Value(NoU, NoV - 1), // Down
+ &myFacePntParams->Value(NoU - 1, NoV - 1), // Lower-left
+ &myFacePntParams->Value(NoU - 1, NoV) }; // Left
+
+ for (i = 0; i < 3 && isMin; i++) {
+ if (anOtherParam[i]->GetElementType() == Extrema_Node) {
+ anOtherParam[i]->GetIndices(iU2, iV2);
+ isMin = (iU == iU2 && iV == iV2);
+ } else {
+ isMin = Standard_False;
+ }
+ }
+ }
+ }
+ }
+
+ if (isMin) {
+ FindSolution(P, aParam);
}
}
}
if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
{
- for (NoV = 0; NoV <= myvsample+1; NoV++) {
- TheDist(0,NoV) = RealFirst();
- TheDist(myusample+1,NoV) = RealFirst();
- }
- for (NoU = 1; NoU <= myusample; NoU++) {
- TheDist(NoU,0) = RealFirst();
- TheDist(NoU,myvsample+1) = RealFirst();
- }
+ Standard_Real Dist;
+
for (NoU = 1; NoU <= myusample; NoU++) {
for (NoV = 1; NoV <= myvsample; NoV++) {
- if (TbSel(NoU,NoV) == 0) {
- Dist = TheDist(NoU,NoV);
- if ((TheDist(NoU-1,NoV-1) <= Dist) &&
- (TheDist(NoU-1,NoV ) <= Dist) &&
- (TheDist(NoU-1,NoV+1) <= Dist) &&
- (TheDist(NoU ,NoV-1) <= Dist) &&
- (TheDist(NoU ,NoV+1) <= Dist) &&
- (TheDist(NoU+1,NoV-1) <= Dist) &&
- (TheDist(NoU+1,NoV ) <= Dist) &&
- (TheDist(NoU+1,NoV+1) <= Dist)) {
- //Create array of UV vectors to calculate max
- UV(1) = U0 + (NoU - 1) * PasU;
- UV(2) = V0 + (NoV - 1) * PasV;
- FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
- }
+ Dist = myPoints->Value(NoU, NoV).GetSqrDistance();
+ if ((myPoints->Value(NoU-1,NoV-1).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU-1,NoV ).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU-1,NoV+1).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU ,NoV-1).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU ,NoV+1).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU+1,NoV-1).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU+1,NoV ).GetSqrDistance() <= Dist) &&
+ (myPoints->Value(NoU+1,NoV+1).GetSqrDistance() <= Dist)) {
+ // Find maximum.
+ FindSolution(P, myPoints->Value(NoU, NoV));
}
}
}
Bnd_SphereUBTreeSelectorMin aSelector(mySphereArray, aSol);
//aSelector.SetMaxDist( RealLast() );
aSelector.DefineCheckPoint( P );
- Standard_Integer aNbSel = mySphereUBTree->Select( aSelector );
+ mySphereUBTree->Select( aSelector );
//TODO: check if no solution in binary tree
Bnd_Sphere& aSph = aSelector.Sphere();
+ Standard_Real aU = myUParams->Value(aSph.U());
+ Standard_Real aV = myVParams->Value(aSph.V());
+ Extrema_POnSurfParams aParams(aU, aV, myS->Value(aU, aV));
- UV(1) = U0 + (aSph.U() - 1) * PasU;
- UV(2) = V0 + (aSph.V() - 1) * PasV;
-
- FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
+ aParams.SetSqrDistance(P.SquareDistance(aParams.Value()));
+ aParams.SetIndices(aSph.U(), aSph.V());
+ FindSolution(P, aParams);
}
if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
{
Bnd_SphereUBTreeSelectorMax aSelector(mySphereArray, aSol);
//aSelector.SetMaxDist( RealLast() );
aSelector.DefineCheckPoint( P );
- Standard_Integer aNbSel = mySphereUBTree->Select( aSelector );
+ mySphereUBTree->Select( aSelector );
//TODO: check if no solution in binary tree
Bnd_Sphere& aSph = aSelector.Sphere();
+ Standard_Real aU = myUParams->Value(aSph.U());
+ Standard_Real aV = myVParams->Value(aSph.V());
+ Extrema_POnSurfParams aParams(aU, aV, myS->Value(aU, aV));
+ aParams.SetSqrDistance(P.SquareDistance(aParams.Value()));
+ aParams.SetIndices(aSph.U(), aSph.V());
- UV(1) = U0 + (aSph.U() - 1) * PasU;
- UV(2) = V0 + (aSph.V() - 1) * PasV;
-
- FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
+ FindSolution(P, aParams);
}
}
}
projects / occt.git / blobdiff