// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
-// See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
// File : SMESH_Pattern.hxx
// Created : Mon Aug 2 10:30:00 2004
#include "SMESH_Block.hxx"
+#include <BRepAdaptor_Curve.hxx>
+#include <BRepAdaptor_Curve2d.hxx>
+#include <BRepAdaptor_Surface.hxx>
#include <BRepTools.hxx>
#include <BRepTools_WireExplorer.hxx>
+#include <BRep_Builder.hxx>
#include <BRep_Tool.hxx>
-#include <TopAbs_ShapeEnum.hxx>
+#include <Bnd_Box.hxx>
+#include <Extrema_ExtPC.hxx>
+#include <Extrema_POnCurv.hxx>
+#include <Geom2d_Curve.hxx>
#include <TopExp_Explorer.hxx>
-#include <TopLoc_Location.hxx>
#include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
#include <TopTools_ListIteratorOfListOfShape.hxx>
#include <TopTools_ListOfShape.hxx>
#include <TopoDS.hxx>
+#include <TopoDS_Compound.hxx>
#include <TopoDS_Face.hxx>
#include <TopoDS_Iterator.hxx>
#include <TopoDS_Wire.hxx>
#include <gp_Trsf.hxx>
#include <gp_Vec.hxx>
#include <math_FunctionSetRoot.hxx>
+#include <math_Matrix.hxx>
+#include <math_Vector.hxx>
+#include "SMDS_MeshNode.hxx"
+#include "SMDS_MeshVolume.hxx"
+#include "SMDS_VolumeTool.hxx"
#include "utilities.h"
#include <list>
using namespace std;
-#define SQRT_FUNC 1
+//#define DEBUG_PARAM_COMPUTE
+
+//================================================================================
+/*!
+ * \brief Set edge data
+ * \param edgeID - block subshape ID
+ * \param curve - edge geometry
+ * \param isForward - is curve orientation coincides with edge orientation in the block
+ */
+//================================================================================
+
+void SMESH_Block::TEdge::Set( const int edgeID, Adaptor3d_Curve* curve, const bool isForward )
+{
+ myCoordInd = SMESH_Block::GetCoordIndOnEdge( edgeID );
+ if ( myC3d ) delete myC3d;
+ myC3d = curve;
+ myFirst = curve->FirstParameter();
+ myLast = curve->LastParameter();
+ if ( !isForward )
+ std::swap( myFirst, myLast );
+}
+
+//================================================================================
+/*!
+ * \brief Set coordinates of nodes at edge ends to work with mesh block
+ * \param edgeID - block subshape ID
+ * \param node1 - coordinates of node with lower ID
+ * \param node2 - coordinates of node with upper ID
+ */
+//================================================================================
+
+void SMESH_Block::TEdge::Set( const int edgeID, const gp_XYZ& node1, const gp_XYZ& node2 )
+{
+ myCoordInd = SMESH_Block::GetCoordIndOnEdge( edgeID );
+ myNodes[ 0 ] = node1;
+ myNodes[ 1 ] = node2;
+
+ if ( myC3d ) delete myC3d;
+ myC3d = 0;
+}
//=======================================================================
//function : SMESH_Block::TEdge::GetU
double SMESH_Block::TEdge::GetU( const gp_XYZ& theParams ) const
{
double u = theParams.Coord( myCoordInd );
+ if ( !myC3d ) // if mesh block
+ return u;
return ( 1 - u ) * myFirst + u * myLast;
}
gp_XYZ SMESH_Block::TEdge::Point( const gp_XYZ& theParams ) const
{
- gp_XYZ p = myC3d->Value( GetU( theParams )).XYZ();
- if ( myTrsf.Form() != gp_Identity )
- myTrsf.Transforms( p );
- return p;
+ double u = GetU( theParams );
+ if ( myC3d ) return myC3d->Value( u ).XYZ();
+ // mesh block
+ return myNodes[0] * ( 1 - u ) + myNodes[1] * u;
+}
+
+//================================================================================
+/*!
+ * \brief Destructor
+ */
+//================================================================================
+
+SMESH_Block::TEdge::~TEdge()
+{
+ if ( myC3d ) delete myC3d;
+}
+
+//================================================================================
+/*!
+ * \brief Set face data
+ * \param faceID - block subshape ID
+ * \param S - face surface geometry
+ * \param c2d - 4 pcurves in the order as returned by GetFaceEdgesIDs(faceID)
+ * \param isForward - orientation of pcurves comparing with block edge direction
+ */
+//================================================================================
+
+void SMESH_Block::TFace::Set( const int faceID,
+ Adaptor3d_Surface* S,
+ Adaptor2d_Curve2d* c2D[4],
+ const bool isForward[4] )
+{
+ if ( myS ) delete myS;
+ myS = S;
+ // pcurves
+ vector< int > edgeIdVec;
+ GetFaceEdgesIDs( faceID, edgeIdVec );
+ for ( int iE = 0; iE < edgeIdVec.size(); iE++ ) // loop on 4 edges
+ {
+ myCoordInd[ iE ] = GetCoordIndOnEdge( edgeIdVec[ iE ] );
+ if ( myC2d[ iE ]) delete myC2d[ iE ];
+ myC2d[ iE ] = c2D[ iE ];
+ myFirst[ iE ] = myC2d[ iE ]->FirstParameter();
+ myLast [ iE ] = myC2d[ iE ]->LastParameter();
+ if ( !isForward[ iE ])
+ std::swap( myFirst[ iE ], myLast[ iE ] );
+ }
+ // 2d corners
+ myCorner[ 0 ] = myC2d[ 0 ]->Value( myFirst[0] ).XY();
+ myCorner[ 1 ] = myC2d[ 0 ]->Value( myLast[0] ).XY();
+ myCorner[ 2 ] = myC2d[ 1 ]->Value( myLast[1] ).XY();
+ myCorner[ 3 ] = myC2d[ 1 ]->Value( myFirst[1] ).XY();
+}
+
+//================================================================================
+/*!
+ * \brief Set face data to work with mesh block
+ * \param faceID - block subshape ID
+ * \param edgeU0 - filled data of edge u0 = GetFaceEdgesIDs(faceID)[ 0 ]
+ * \param edgeU1 - filled data of edge u1 = GetFaceEdgesIDs(faceID)[ 1 ]
+ */
+//================================================================================
+
+void SMESH_Block::TFace::Set( const int faceID, const TEdge& edgeU0, const TEdge& edgeU1 )
+{
+ vector< int > edgeIdVec;
+ GetFaceEdgesIDs( faceID, edgeIdVec );
+ myNodes[ 0 ] = edgeU0.NodeXYZ( 1 );
+ myNodes[ 1 ] = edgeU0.NodeXYZ( 0 );
+ myNodes[ 2 ] = edgeU1.NodeXYZ( 0 );
+ myNodes[ 3 ] = edgeU1.NodeXYZ( 1 );
+ myCoordInd[ 0 ] = GetCoordIndOnEdge( edgeIdVec[ 0 ] );
+ myCoordInd[ 1 ] = GetCoordIndOnEdge( edgeIdVec[ 1 ] );
+ myCoordInd[ 2 ] = GetCoordIndOnEdge( edgeIdVec[ 2 ] );
+ myCoordInd[ 3 ] = GetCoordIndOnEdge( edgeIdVec[ 3 ] );
+ if ( myS ) delete myS;
+ myS = 0;
+}
+
+//================================================================================
+/*!
+ * \brief Destructor
+ */
+//================================================================================
+
+SMESH_Block::TFace::~TFace()
+{
+ if ( myS ) delete myS;
+ for ( int i = 0 ; i < 4; ++i )
+ if ( myC2d[ i ]) delete myC2d[ i ];
+}
+
+//=======================================================================
+//function : SMESH_Block::TFace::GetCoefs
+//purpose : return coefficients for addition of [0-3]-th edge and vertex
+//=======================================================================
+
+void SMESH_Block::TFace::GetCoefs(int iE,
+ const gp_XYZ& theParams,
+ double& Ecoef,
+ double& Vcoef ) const
+{
+ double dU = theParams.Coord( GetUInd() );
+ double dV = theParams.Coord( GetVInd() );
+ switch ( iE ) {
+ case 0:
+ Ecoef = ( 1 - dV ); // u0
+ Vcoef = ( 1 - dU ) * ( 1 - dV ); break; // 00
+ case 1:
+ Ecoef = dV; // u1
+ Vcoef = dU * ( 1 - dV ); break; // 10
+ case 2:
+ Ecoef = ( 1 - dU ); // 0v
+ Vcoef = dU * dV ; break; // 11
+ case 3:
+ Ecoef = dU ; // 1v
+ Vcoef = ( 1 - dU ) * dV ; break; // 01
+ default: ASSERT(0);
+ }
}
//=======================================================================
gp_XY SMESH_Block::TFace::GetUV( const gp_XYZ& theParams ) const
{
gp_XY uv(0.,0.);
- double dU = theParams.Coord( GetUInd() );
- double dV = theParams.Coord( GetVInd() );
for ( int iE = 0; iE < 4; iE++ ) // loop on 4 edges
{
double Ecoef = 0, Vcoef = 0;
- switch ( iE ) {
- case 0:
- Ecoef = ( 1 - dV ); // u0
- Vcoef = ( 1 - dU ) * ( 1 - dV ); break; // 00
- case 1:
- Ecoef = dV; // u1
- Vcoef = dU * ( 1 - dV ); break; // 10
- case 2:
- Ecoef = ( 1 - dU ); // 0v
- Vcoef = dU * dV ; break; // 11
- case 3:
- Ecoef = dU ; // 1v
- Vcoef = ( 1 - dU ) * dV ; break; // 01
- default:;
- }
+ GetCoefs( iE, theParams, Ecoef, Vcoef );
// edge addition
double u = theParams.Coord( myCoordInd[ iE ] );
u = ( 1 - u ) * myFirst[ iE ] + u * myLast[ iE ];
gp_XYZ SMESH_Block::TFace::Point( const gp_XYZ& theParams ) const
{
- gp_XY uv = GetUV( theParams );
- gp_XYZ p = myS->Value( uv.X(), uv.Y() ).XYZ();
- if ( myTrsf.Form() != gp_Identity )
- myTrsf.Transforms( p );
+ gp_XYZ p(0.,0.,0.);
+ if ( !myS ) // if mesh block
+ {
+ for ( int iE = 0; iE < 4; iE++ ) // loop on 4 edges
+ {
+ double Ecoef = 0, Vcoef = 0;
+ GetCoefs( iE, theParams, Ecoef, Vcoef );
+ // edge addition
+ double u = theParams.Coord( myCoordInd[ iE ] );
+ int i1 = 0, i2 = 1;
+ switch ( iE ) {
+ case 1: i1 = 3; i2 = 2; break;
+ case 2: i1 = 1; i2 = 2; break;
+ case 3: i1 = 0; i2 = 3; break;
+ }
+ p += Ecoef * ( myNodes[ i1 ] * ( 1 - u ) + myNodes[ i2 ] * u );
+ // corner addition
+ p -= Vcoef * myNodes[ iE ];
+ }
+
+ }
+ else // shape block
+ {
+ gp_XY uv = GetUV( theParams );
+ p = myS->Value( uv.X(), uv.Y() ).XYZ();
+ }
return p;
}
k *= theParams.Coord( iCoef + 1 );
}
}
- // point on a shape
- gp_XYZ Ps;
- if ( shapeID < ID_Ex00 ) // vertex
- VertexPoint( shapeID, Ps );
- else if ( shapeID < ID_Fxy0 ) { // edge
- EdgePoint( shapeID, theParams, Ps );
- k = -k;
- } else // face
- FacePoint( shapeID, theParams, Ps );
-
- thePoint += k * Ps;
+ // add point on a shape
+ if ( fabs( k ) > DBL_MIN )
+ {
+ gp_XYZ Ps;
+ if ( shapeID < ID_Ex00 ) // vertex
+ VertexPoint( shapeID, Ps );
+ else if ( shapeID < ID_Fxy0 ) { // edge
+ EdgePoint( shapeID, theParams, Ps );
+ k = -k;
+ } else // face
+ FacePoint( shapeID, theParams, Ps );
+
+ thePoint += k * Ps;
+ }
}
return true;
}
+//=======================================================================
+//function : ShellPoint
+//purpose : computes coordinates of a point in shell by points on sub-shapes;
+// thePointOnShape[ subShapeID ] must be a point on a subShape
+//=======================================================================
+
+bool SMESH_Block::ShellPoint(const gp_XYZ& theParams,
+ const vector<gp_XYZ>& thePointOnShape,
+ gp_XYZ& thePoint )
+{
+ if ( thePointOnShape.size() < ID_F1yz )
+ return false;
+
+ double x = theParams.X(), y = theParams.Y(), z = theParams.Z();
+ double x1 = 1. - x, y1 = 1. - y, z1 = 1. - z;
+ const vector<gp_XYZ>& p = thePointOnShape;
+
+ thePoint =
+ x1 * p[ID_F0yz] + x * p[ID_F1yz] +
+ y1 * p[ID_Fx0z] + y * p[ID_Fx1z] +
+ z1 * p[ID_Fxy0] + z * p[ID_Fxy1] +
+ x1 * (y1 * (z1 * p[ID_V000] + z * p[ID_V001]) +
+ y * (z1 * p[ID_V010] + z * p[ID_V011])) +
+ x * (y1 * (z1 * p[ID_V100] + z * p[ID_V101]) +
+ y * (z1 * p[ID_V110] + z * p[ID_V111]));
+ thePoint -=
+ x1 * (y1 * p[ID_E00z] + y * p[ID_E01z]) +
+ x * (y1 * p[ID_E10z] + y * p[ID_E11z]) +
+ y1 * (z1 * p[ID_Ex00] + z * p[ID_Ex01]) +
+ y * (z1 * p[ID_Ex10] + z * p[ID_Ex11]) +
+ z1 * (x1 * p[ID_E0y0] + x * p[ID_E1y0]) +
+ z * (x1 * p[ID_E0y1] + x * p[ID_E1y1]);
+
+ return true;
+}
+
+//=======================================================================
+//function : Constructor
+//purpose :
+//=======================================================================
+
+SMESH_Block::SMESH_Block():
+ myNbIterations(0),
+ mySumDist(0.),
+ myTolerance(-1.) // to be re-initialized
+{
+}
+
+
//=======================================================================
//function : NbVariables
//purpose :
{
gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) );
if ( params.IsEqual( myParam, DBL_MIN )) { // same param
- theFxyz( 1 ) = myValues[ 0 ];
+ theFxyz( 1 ) = funcValue( myValues[ SQUARE_DIST ]);
}
else {
ShellPoint( params, P );
gp_Vec dP( P - myPoint );
- theFxyz(1) = SQRT_FUNC ? dP.SquareMagnitude() : dP.Magnitude();
+ theFxyz(1) = funcValue( dP.SquareMagnitude() );
}
return true;
}
Standard_Boolean SMESH_Block::Derivatives(const math_Vector& XYZ,math_Matrix& Df)
{
- MESSAGE( "SMESH_Block::Derivatives()");
math_Vector F(1,3);
return Values(XYZ,F,Df);
}
+//=======================================================================
+//function : GetStateNumber
+//purpose :
+//=======================================================================
+
+Standard_Integer SMESH_Block::GetStateNumber ()
+{
+ return 0; //myValues[0] < 1e-1;
+}
+
//=======================================================================
//function : Values
//purpose :
math_Vector& theFxyz,
math_Matrix& theDf)
{
-// MESSAGE( endl<<"SMESH_Block::Values( "<<theXYZ(1)<<", "<<theXYZ(2)<<", "<<theXYZ(3)<<")");
-
gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) );
if ( params.IsEqual( myParam, DBL_MIN )) { // same param
- theFxyz( 1 ) = myValues[ 0 ];
- theDf( 1,1 ) = myValues[ 1 ];
- theDf( 1,2 ) = myValues[ 2 ];
- theDf( 1,3 ) = myValues[ 3 ];
+ theFxyz( 1 ) = funcValue( myValues[ SQUARE_DIST ] );
+ theDf( 1, DRV_1 ) = myValues[ DRV_1 ];
+ theDf( 1, DRV_2 ) = myValues[ DRV_2 ];
+ theDf( 1, DRV_3 ) = myValues[ DRV_3 ];
return true;
}
-
+#ifdef DEBUG_PARAM_COMPUTE
+ cout << "PARAM GUESS: " << params.X() << " "<< params.Y() << " "<< params.X() << endl;
+ myNbIterations++; // how many times call ShellPoint()
+#endif
ShellPoint( params, P );
- //myNbIterations++; // how many time call ShellPoint()
-
- gp_Vec dP( P - myPoint );
- theFxyz(1) = SQRT_FUNC ? dP.SquareMagnitude() : dP.Magnitude();
- if ( theFxyz(1) < 1e-6 ) {
- myParam = params;
- myValues[ 0 ]= 0;
- theDf( 1,1 ) = 0;
- theDf( 1,2 ) = 0;
- theDf( 1,3 ) = 0;
+
+ gp_Vec dP( myPoint, P );
+ double sqDist = dP.SquareMagnitude();
+ theFxyz(1) = funcValue( sqDist );
+
+ if ( sqDist < myTolerance * myTolerance ) { // a solution found
+ myParam = params;
+ myValues[ SQUARE_DIST ] = sqDist;
+ theFxyz(1) = theDf( 1,1 ) = theDf( 1,2 ) = theDf( 1,3 ) = 0;
return true;
}
- if ( theFxyz(1) < myValues[0] )
+ if ( sqDist < myValues[ SQUARE_DIST ] ) // a better guess
{
// 3 partial derivatives
- gp_Vec drv[ 3 ];
+ gp_Vec drv[ 3 ]; // where we move with a small step in each direction
for ( int iP = 1; iP <= 3; iP++ ) {
+ if ( iP == myFaceIndex ) {
+ drv[ iP - 1 ] = gp_Vec(0,0,0);
+ continue;
+ }
gp_XYZ Pi;
- params.SetCoord( iP, theXYZ( iP ) + 0.001 );
+ bool onEdge = ( theXYZ( iP ) + 0.001 > 1. );
+ if ( onEdge )
+ params.SetCoord( iP, theXYZ( iP ) - 0.001 );
+ else
+ params.SetCoord( iP, theXYZ( iP ) + 0.001 );
ShellPoint( params, Pi );
params.SetCoord( iP, theXYZ( iP ) ); // restore params
gp_Vec dPi ( P, Pi );
+ if ( onEdge ) dPi *= -1.;
double mag = dPi.Magnitude();
if ( mag > DBL_MIN )
dPi /= mag;
drv[ iP - 1 ] = dPi;
}
for ( int iP = 0; iP < 3; iP++ ) {
- if ( iP == myFaceIndex )
- theDf( 1, iP + 1 ) = myFaceParam;
+#if 1
+ theDf( 1, iP + 1 ) = dP * drv[iP];
+#else
+ // Distance from P to plane passing through myPoint and defined
+ // by the 2 other derivative directions:
+ // like IntAna_IntConicQuad::Perform (const gp_Lin& L, const gp_Pln& P)
+ // where L is (P -> myPoint), P is defined by the 2 other derivative direction
+ int iPrev = ( iP ? iP - 1 : 2 );
+ int iNext = ( iP == 2 ? 0 : iP + 1 );
+ gp_Vec plnNorm = drv[ iPrev ].Crossed( drv [ iNext ] );
+ double Direc = plnNorm * drv[ iP ];
+ if ( Abs(Direc) <= DBL_MIN )
+ theDf( 1, iP + 1 ) = dP * drv[ iP ];
else {
- // like IntAna_IntConicQuad::Perform (const gp_Lin& L, const gp_Pln& P)
- // where L is (P -> myPoint), P is defined by the 2 other derivative direction
- int iPrev = ( iP ? iP - 1 : 2 );
- int iNext = ( iP == 2 ? 0 : iP + 1 );
- gp_Vec plnNorm = drv[ iPrev ].Crossed( drv [ iNext ] );
- double Direc = plnNorm * drv[ iP ];
- if ( Abs(Direc) <= DBL_MIN )
- theDf( 1, iP + 1 ) = dP * drv[ iP ];
- else {
- double Dis = plnNorm * P - plnNorm * myPoint;
- theDf( 1, iP + 1 ) = Dis/Direc;
- }
+ double Dis = plnNorm * P - plnNorm * myPoint;
+ theDf( 1, iP + 1 ) = Dis/Direc;
}
+#endif
}
- //myNbIterations +=3; // how many time call ShellPoint()
+#ifdef DEBUG_PARAM_COMPUTE
+ cout << "F = " << theFxyz(1) <<
+ " DRV: " << theDf(1,1) << " " << theDf(1,2) << " " << theDf(1,3) << endl;
+ myNbIterations +=3; // how many times call ShellPoint()
+#endif
// store better values
- myParam = params;
- myValues[0]= theFxyz(1);
- myValues[1]= theDf(1,1);
- myValues[2]= theDf(1,2);
- myValues[3]= theDf(1,3);
+ myParam = params;
+ myValues[SQUARE_DIST]= sqDist;
+ myValues[DRV_1] = theDf(1,DRV_1);
+ myValues[DRV_2] = theDf(1,DRV_2);
+ myValues[DRV_3] = theDf(1,DRV_3);
+ }
+
+ return true;
+}
+
+//============================================================================
+//function : computeParameters
+//purpose : compute point parameters in the block using math_FunctionSetRoot
+//============================================================================
+
+bool SMESH_Block::computeParameters(const gp_Pnt& thePoint,
+ gp_XYZ& theParams,
+ const gp_XYZ& theParamsHint)
+{
+ myPoint = thePoint.XYZ();
+
+ myParam.SetCoord( -1,-1,-1 );
+ myValues[ SQUARE_DIST ] = 1e100;
+
+ math_Vector low ( 1, 3, 0.0 );
+ math_Vector up ( 1, 3, 1.0 );
+ math_Vector tol ( 1, 3, 1e-4 );
+ math_Vector start( 1, 3, 0.0 );
+ start( 1 ) = theParamsHint.X();
+ start( 2 ) = theParamsHint.Y();
+ start( 3 ) = theParamsHint.Z();
+
+ math_FunctionSetRoot paramSearch( *this, tol );
+
+ mySquareFunc = 0; // large approaching steps
+ //if ( hasHint ) mySquareFunc = 1; // small approaching steps
-// SCRUTE( theFxyz(1) );
-// SCRUTE( theDf( 1,1 ));
-// SCRUTE( theDf( 1,2 ));
-// SCRUTE( theDf( 1,3 ));
+ double loopTol = 10 * myTolerance;
+ int nbLoops = 0;
+ while ( distance() > loopTol && nbLoops <= 3 )
+ {
+ paramSearch.Perform ( *static_cast<math_FunctionSetWithDerivatives*>(this),
+ start, low, up );
+ start( 1 ) = myParam.X();
+ start( 2 ) = myParam.Y();
+ start( 3 ) = myParam.Z();
+ mySquareFunc = !mySquareFunc;
+ nbLoops++;
}
+#ifdef DEBUG_PARAM_COMPUTE
+ mySumDist += distance();
+ cout << " ------ SOLUTION: ( "<< myParam.X() <<" "<< myParam.Y() <<" "<< myParam.Z() <<" )"<<endl
+ << " ------ DIST : " << distance() << "\t Tol=" << myTolerance << "\t Nb LOOPS=" << nbLoops << endl
+ << " ------ NB IT: " << myNbIterations << ", SUM DIST: " << mySumDist << endl;
+#endif
+
+ theParams = myParam;
+
+ if ( myFaceIndex > 0 )
+ theParams.SetCoord( myFaceIndex, myFaceParam );
return true;
}
bool SMESH_Block::ComputeParameters(const gp_Pnt& thePoint,
gp_XYZ& theParams,
- const int theShapeID)
+ const int theShapeID,
+ const gp_XYZ& theParamsHint)
{
-// MESSAGE( endl<<"SMESH_Block::ComputeParameters( "
-// <<thePoint.X()<<", "<<thePoint.Y()<<", "<<thePoint.Z()<<")");
-
- myPoint = thePoint.XYZ();
+ if ( VertexParameters( theShapeID, theParams ))
+ return true;
- myParam.SetCoord( -1,-1,-1 );
- myValues[0] = 1e100;
+ if ( IsEdgeID( theShapeID )) {
+ TEdge& e = myEdge[ theShapeID - ID_FirstE ];
+ Adaptor3d_Curve* curve = e.GetCurve();
+ Extrema_ExtPC anExtPC( thePoint, *curve, curve->FirstParameter(), curve->LastParameter() );
+ int i, nb = anExtPC.IsDone() ? anExtPC.NbExt() : 0;
+ for ( i = 1; i <= nb; i++ ) {
+ if ( anExtPC.IsMin( i ))
+ return EdgeParameters( theShapeID, anExtPC.Point( i ).Parameter(), theParams );
+ }
+ return false;
+ }
const bool isOnFace = IsFaceID( theShapeID );
double * coef = GetShapeCoef( theShapeID );
- // the first guess
- math_Vector start( 1, 3, 0.0 );
- if ( !myGridComputed )
+ // Find the first guess paremeters
+
+ gp_XYZ start(0, 0, 0);
+
+ bool hasHint = ( 0 <= theParamsHint.X() && theParamsHint.X() <= 1 &&
+ 0 <= theParamsHint.Y() && theParamsHint.Y() <= 1 &&
+ 0 <= theParamsHint.Y() && theParamsHint.Y() <= 1 );
+ if ( !hasHint && !myGridComputed )
{
// define the first guess by thePoint projection on lines
// connecting vertices
iEdge += 4;
continue;
}
+ double sumParam = 0;
for ( int iE = 0; iE < 4; iE++, iEdge++ ) { // loop on 4 parallel edges
gp_Pnt p0 = myEdge[ iEdge ].Point( par000 );
gp_Pnt p1 = myEdge[ iEdge ].Point( par111 );
double par = 0;
if ( len2 > zero ) {
par = v0P.Dot( v01 ) / len2;
- if ( par < 0 || par > 1 ) {
+ if ( par < 0 || par > 1 ) { // projection falls out of line ends => needGrid
needGrid = true;
break;
}
}
- start( iParam ) += par;
+ sumParam += par;
}
- start( iParam ) /= 4.;
+ start.SetCoord( iParam, sumParam / 4.);
}
if ( needGrid ) {
// compute nodes of 3 x 3 x 3 grid
int iNode = 0;
+ Bnd_Box box;
for ( double x = 0.25; x < 0.9; x += 0.25 )
for ( double y = 0.25; y < 0.9; y += 0.25 )
for ( double z = 0.25; z < 0.9; z += 0.25 ) {
TxyzPair & prmPtn = my3x3x3GridNodes[ iNode++ ];
prmPtn.first.SetCoord( x, y, z );
ShellPoint( prmPtn.first, prmPtn.second );
+ box.Add( gp_Pnt( prmPtn.second ));
}
myGridComputed = true;
+ myTolerance = sqrt( box.SquareExtent() ) * 1e-5;
}
}
- if ( myGridComputed ) {
+
+ if ( hasHint )
+ {
+ start = theParamsHint;
+ }
+ else if ( myGridComputed )
+ {
double minDist = DBL_MAX;
gp_XYZ* bestParam = 0;
for ( int iNode = 0; iNode < 27; iNode++ ) {
bestParam = & prmPtn.first;
}
}
- start( 1 ) = bestParam->X();
- start( 2 ) = bestParam->Y();
- start( 3 ) = bestParam->Z();
+ start = *bestParam;
}
- int myFaceIndex = -1;
+ int myFaceIndex = -1;
+ double myFaceParam = 0.;
if ( isOnFace ) {
// put a point on the face
for ( int iCoord = 0; iCoord < 3; iCoord++ )
if ( coef[ iCoord ] ) {
- myFaceIndex = iCoord;
- myFaceParam = ( coef[ myFaceIndex ] < 0.5 ) ? 0.0 : 1.0;
- start( iCoord + 1 ) = myFaceParam;
+ myFaceIndex = iCoord + 1;
+ myFaceParam = ( coef[ iCoord ] < 0.5 ) ? 0.0 : 1.0;
+ start.SetCoord( myFaceIndex, myFaceParam );
}
}
- math_Vector low ( 1, 3, 0.0 );
- math_Vector up ( 1, 3, 1.0 );
- math_Vector tol ( 1, 3, 1e-4 );
- math_FunctionSetRoot paramSearch( *this, tol );
- int nbLoops = 0;
- while ( myValues[0] > 1e-1 && nbLoops++ < 10 ) {
- paramSearch.Perform ( *this, start, low, up );
- if ( !paramSearch.IsDone() ) {
- //MESSAGE( " !paramSearch.IsDone() " );
+#ifdef DEBUG_PARAM_COMPUTE
+ cout << " #### POINT " <<thePoint.X()<<" "<<thePoint.Y()<<" "<<thePoint.Z()<<" ####"<< endl;
+#endif
+
+ if ( myTolerance < 0 ) myTolerance = 1e-6;
+
+ const double parDelta = 1e-4;
+ const double sqTolerance = myTolerance * myTolerance;
+
+ gp_XYZ solution = start, params = start;
+ double sqDistance = 1e100;
+ int nbLoops = 0, nbGetWorst = 0;
+
+ while ( nbLoops <= 100 )
+ {
+ gp_XYZ P, Pi;
+ ShellPoint( params, P );
+
+ gp_Vec dP( thePoint, P );
+ double sqDist = dP.SquareMagnitude();
+
+ if ( sqDist > sqDistance ) { // solution get worse
+ if ( ++nbGetWorst > 2 )
+ return computeParameters( thePoint, theParams, solution );
}
- else {
- //MESSAGE( " NB ITERATIONS: " << paramSearch.NbIterations() );
+#ifdef DEBUG_PARAM_COMPUTE
+ cout << "PARAMS: ( " << params.X() <<" "<< params.Y() <<" "<< params.Z() <<" )"<< endl;
+ cout << "DIST: " << sqrt( sqDist ) << endl;
+#endif
+
+ if ( sqDist < sqDistance ) { // get better
+ sqDistance = sqDist;
+ solution = params;
+ nbGetWorst = 0;
+ if ( sqDistance < sqTolerance ) // a solution found
+ break;
}
- start( 1 ) = myParam.X();
- start( 2 ) = myParam.Y();
- start( 3 ) = myParam.Z();
- //MESSAGE( "Distance: " << ( SQRT_FUNC ? sqrt(myValues[0]) : myValues[0] ));
+
+ // look for a next better solution
+ for ( int iP = 1; iP <= 3; iP++ ) {
+ if ( iP == myFaceIndex )
+ continue;
+ // see where we move with a small (=parDelta) step in this direction
+ gp_XYZ nearParams = params;
+ bool onEdge = ( params.Coord( iP ) + parDelta > 1. );
+ if ( onEdge )
+ nearParams.SetCoord( iP, params.Coord( iP ) - parDelta );
+ else
+ nearParams.SetCoord( iP, params.Coord( iP ) + parDelta );
+ ShellPoint( nearParams, Pi );
+ gp_Vec dPi ( P, Pi );
+ if ( onEdge ) dPi *= -1.;
+ // modify a parameter
+ double mag = dPi.Magnitude();
+ if ( mag < DBL_MIN )
+ continue;
+ gp_Vec dir = dPi / mag; // dir we move modifying the parameter
+ double dist = dir * dP; // where we should get to
+ double dPar = dist / mag * parDelta; // predict parameter change
+ double curPar = params.Coord( iP );
+ double par = curPar - dPar; // new parameter value
+ while ( par > 1 || par < 0 ) {
+ dPar /= 2.;
+ par = curPar - dPar;
+ }
+ params.SetCoord( iP, par );
+ }
+
+ nbLoops++;
}
-// MESSAGE( endl << myParam.X() << " " << myParam.Y() << " " << myParam.Z() << endl);
-// mySumDist += myValues[0];
-// MESSAGE( " TOTAL NB ITERATIONS: " << myNbIterations <<
-// " DIST: " << ( SQRT_FUNC ? sqrt(mySumDist) : mySumDist ));
+#ifdef DEBUG_PARAM_COMPUTE
+ myNbIterations += nbLoops*4; // how many times ShellPoint called
+ mySumDist += sqrt( sqDistance );
+ cout << " ------ SOLUTION: ( "<<solution.X()<<" "<<solution.Y()<<" "<<solution.Z()<<" )"<<endl
+ << " ------ DIST : " << sqrt( sqDistance ) << "\t Tol=" << myTolerance << "\t Nb LOOPS=" << nbLoops << endl
+ << " ------ NB IT: " << myNbIterations << ", SUM DIST: " << mySumDist << endl;
+#endif
+ theParams = solution;
- theParams = myParam;
+ if ( myFaceIndex > 0 )
+ theParams.SetCoord( myFaceIndex, myFaceParam );
return true;
}
//=======================================================================
-//function : GetStateNumber
-//purpose :
+//function : VertexParameters
+//purpose : return parameters of a vertex given by TShapeID
//=======================================================================
-Standard_Integer SMESH_Block::GetStateNumber ()
+bool SMESH_Block::VertexParameters(const int theVertexID, gp_XYZ& theParams)
{
-// MESSAGE( endl<<"SMESH_Block::GetStateNumber( "<<myParam.X()<<", "<<
-// myParam.Y()<<", "<<myParam.Z()<<") DISTANCE: " << myValues[0]);
- return myValues[0] < 1e-1;
+ switch ( theVertexID ) {
+ case ID_V000: theParams.SetCoord(0., 0., 0.); return true;
+ case ID_V100: theParams.SetCoord(1., 0., 0.); return true;
+ case ID_V110: theParams.SetCoord(1., 1., 0.); return true;
+ case ID_V010: theParams.SetCoord(0., 1., 0.); return true;
+ default:;
+ }
+ return false;
+}
+
+//=======================================================================
+//function : EdgeParameters
+//purpose : return parameters of a point given by theU on edge
+//=======================================================================
+
+bool SMESH_Block::EdgeParameters(const int theEdgeID, const double theU, gp_XYZ& theParams)
+{
+ if ( IsEdgeID( theEdgeID )) {
+ vector< int > vertexVec;
+ GetEdgeVertexIDs( theEdgeID, vertexVec );
+ VertexParameters( vertexVec[0], theParams );
+ TEdge& e = myEdge[ theEdgeID - ID_Ex00 ];
+ double param = ( theU - e.EndParam(0) ) / ( e.EndParam(1) - e.EndParam(0) );
+ theParams.SetCoord( e.CoordInd(), param );
+ return true;
+ }
+ return false;
}
//=======================================================================
return id + 1; // shape ids start at 1
}
-
//=======================================================================
-//function : getOrderedEdges
+//function : GetOrderedEdges
//purpose : return nb wires and a list of oredered edges
//=======================================================================
-static int getOrderedEdges (const TopoDS_Face& theFace,
- const TopoDS_Vertex& theFirstVertex,
- list< TopoDS_Edge >& theEdges,
- list< int > & theNbVertexInWires)
+int SMESH_Block::GetOrderedEdges (const TopoDS_Face& theFace,
+ TopoDS_Vertex theFirstVertex,
+ list< TopoDS_Edge >& theEdges,
+ list< int > & theNbVertexInWires)
{
// put wires in a list, so that an outer wire comes first
list<TopoDS_Wire> aWireList;
isFirst = ( p2.SquareDistance( pf ) < p2.SquareDistance( pl ));
if ( isNext ? isFirst : !isFirst )
edge.Reverse();
+ // to make a seam go first
+ if ( theFirstVertex.IsNull() )
+ theFirstVertex = TopExp::FirstVertex( edge, true );
}
}
// rotate theEdges until it begins from theFirstVertex
- if ( ! theFirstVertex.IsNull() )
- while ( !theFirstVertex.IsSame( TopExp::FirstVertex( theEdges.front(), true )))
+ if ( ! theFirstVertex.IsNull() ) {
+ TopoDS_Vertex vv[2];
+ TopExp::Vertices( theEdges.front(), vv[0], vv[1], true );
+ // on closed face, make seam edge the first in the list
+ while ( !vv[0].IsSame( theFirstVertex ) || vv[0].IsSame( vv[1] ))
{
theEdges.splice(theEdges.end(), theEdges,
- theEdges.begin(), ++ theEdges.begin());
- if ( iE++ > theNbVertexInWires.back() )
+ theEdges.begin(), ++theEdges.begin());
+ TopExp::Vertices( theEdges.front(), vv[0], vv[1], true );
+ if ( iE++ > theNbVertexInWires.back() ) {
+#ifdef _DEBUG_
+ gp_Pnt p = BRep_Tool::Pnt( theFirstVertex );
+ cout << " : Warning : vertex "<< theFirstVertex.TShape().operator->()
+ << " ( " << p.X() << " " << p.Y() << " " << p.Z() << " )"
+ << " not found in outer wire of face "<< theFace.TShape().operator->()
+ << " with vertices: " << endl;
+ wExp.Init( *wlIt, theFace );
+ for ( int i = 0; wExp.More(); wExp.Next(), i++ )
+ {
+ TopoDS_Edge edge = wExp.Current();
+ edge = TopoDS::Edge( edge.Oriented( wExp.Orientation() ));
+ TopoDS_Vertex v = TopExp::FirstVertex( edge, true );
+ gp_Pnt p = BRep_Tool::Pnt( v );
+ cout << i << " " << v.TShape().operator->() << " "
+ << p.X() << " " << p.Y() << " " << p.Z() << " " << endl;
+ }
+#endif
break; // break infinite loop
+ }
}
- }
+ }
+ } // end outer wire
}
return aWireList.size();
}
+//================================================================================
+/*!
+ * \brief Call it after geometry initialisation
+ */
+//================================================================================
+
+void SMESH_Block::init()
+{
+ myNbIterations = 0;
+ mySumDist = 0;
+ myGridComputed = false;
+}
+
+//=======================================================================
+//function : LoadMeshBlock
+//purpose : prepare to work with theVolume
+//=======================================================================
+
+#define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
+
+bool SMESH_Block::LoadMeshBlock(const SMDS_MeshVolume* theVolume,
+ const int theNode000Index,
+ const int theNode001Index,
+ vector<const SMDS_MeshNode*>& theOrderedNodes)
+{
+ MESSAGE(" ::LoadMeshBlock()");
+ init();
+
+ SMDS_VolumeTool vTool;
+ if (!vTool.Set( theVolume ) || vTool.NbNodes() != 8 ||
+ !vTool.IsLinked( theNode000Index, theNode001Index )) {
+ MESSAGE(" Bad arguments ");
+ return false;
+ }
+ vTool.SetExternalNormal();
+ // In terms of indices used for access to nodes and faces in SMDS_VolumeTool:
+ int V000, V100, V010, V110, V001, V101, V011, V111; // 8 vertices
+ int Fxy0, Fxy1; // bottom and top faces
+ // vertices of faces
+ vector<int> vFxy0, vFxy1;
+
+ V000 = theNode000Index;
+ V001 = theNode001Index;
+
+ // get faces sharing V000 and V001
+ list<int> fV000, fV001;
+ int i, iF, iE, iN;
+ for ( iF = 0; iF < vTool.NbFaces(); ++iF ) {
+ const int* nid = vTool.GetFaceNodesIndices( iF );
+ for ( iN = 0; iN < 4; ++iN )
+ if ( nid[ iN ] == V000 ) {
+ fV000.push_back( iF );
+ } else if ( nid[ iN ] == V001 ) {
+ fV001.push_back( iF );
+ }
+ }
+
+ // find the bottom (Fxy0), the top (Fxy1) faces
+ list<int>::iterator fIt1, fIt2, Fxy0Pos;
+ for ( fIt1 = fV000.begin(); fIt1 != fV000.end(); fIt1++) {
+ fIt2 = std::find( fV001.begin(), fV001.end(), *fIt1 );
+ if ( fIt2 != fV001.end() ) { // *fIt1 is in the both lists
+ fV001.erase( fIt2 ); // erase Fx0z or F0yz from fV001
+ } else { // *fIt1 is in fV000 only
+ Fxy0Pos = fIt1; // points to Fxy0
+ }
+ }
+ Fxy0 = *Fxy0Pos;
+ Fxy1 = fV001.front();
+ const SMDS_MeshNode** nn = vTool.GetNodes();
+
+ // find bottom veritices, their order is that a face normal is external
+ vFxy0.resize(4);
+ const int* nid = vTool.GetFaceNodesIndices( Fxy0 );
+ for ( i = 0; i < 4; ++i )
+ if ( nid[ i ] == V000 )
+ break;
+ for ( iN = 0; iN < 4; ++iN, ++i ) {
+ if ( i == 4 ) i = 0;
+ vFxy0[ iN ] = nid[ i ];
+ }
+ // find top veritices, their order is that a face normal is external
+ vFxy1.resize(4);
+ nid = vTool.GetFaceNodesIndices( Fxy1 );
+ for ( i = 0; i < 4; ++i )
+ if ( nid[ i ] == V001 )
+ break;
+ for ( iN = 0; iN < 4; ++iN, ++i ) {
+ if ( i == 4 ) i = 0;
+ vFxy1[ iN ] = nid[ i ];
+ }
+ // find indices of the rest veritices
+ V100 = vFxy0[3];
+ V010 = vFxy0[1];
+ V110 = vFxy0[2];
+ V101 = vFxy1[1];
+ V011 = vFxy1[3];
+ V111 = vFxy1[2];
+
+ // set points coordinates
+ myPnt[ ID_V000 - 1 ] = gpXYZ( nn[ V000 ] );
+ myPnt[ ID_V100 - 1 ] = gpXYZ( nn[ V100 ] );
+ myPnt[ ID_V010 - 1 ] = gpXYZ( nn[ V010 ] );
+ myPnt[ ID_V110 - 1 ] = gpXYZ( nn[ V110 ] );
+ myPnt[ ID_V001 - 1 ] = gpXYZ( nn[ V001 ] );
+ myPnt[ ID_V101 - 1 ] = gpXYZ( nn[ V101 ] );
+ myPnt[ ID_V011 - 1 ] = gpXYZ( nn[ V011 ] );
+ myPnt[ ID_V111 - 1 ] = gpXYZ( nn[ V111 ] );
+
+ // fill theOrderedNodes
+ theOrderedNodes.resize( 8 );
+ theOrderedNodes[ 0 ] = nn[ V000 ];
+ theOrderedNodes[ 1 ] = nn[ V100 ];
+ theOrderedNodes[ 2 ] = nn[ V010 ];
+ theOrderedNodes[ 3 ] = nn[ V110 ];
+ theOrderedNodes[ 4 ] = nn[ V001 ];
+ theOrderedNodes[ 5 ] = nn[ V101 ];
+ theOrderedNodes[ 6 ] = nn[ V011 ];
+ theOrderedNodes[ 7 ] = nn[ V111 ];
+
+ // fill edges
+ vector< int > vertexVec;
+ for ( iE = 0; iE < NbEdges(); ++iE ) {
+ GetEdgeVertexIDs(( iE + ID_FirstE ), vertexVec );
+ myEdge[ iE ].Set(( iE + ID_FirstE ),
+ myPnt[ vertexVec[0] - 1 ],
+ myPnt[ vertexVec[1] - 1 ]);
+ }
+
+ // fill faces' corners
+ for ( iF = ID_Fxy0; iF < ID_Shell; ++iF )
+ {
+ TFace& tFace = myFace[ iF - ID_FirstF ];
+ vector< int > edgeIdVec(4, -1);
+ GetFaceEdgesIDs( iF, edgeIdVec );
+ tFace.Set( iF, myEdge[ edgeIdVec [ 0 ] - ID_Ex00], myEdge[ edgeIdVec [ 1 ] - ID_Ex00]);
+ }
+
+ return true;
+}
//=======================================================================
//function : LoadBlockShapes
-//purpose : add sub-shapes of theBlock to theShapeIDMap so that they get
+//purpose : Initialize block geometry with theShell,
+// add sub-shapes of theBlock to theShapeIDMap so that they get
// IDs acoording to enum TShapeID
//=======================================================================
bool SMESH_Block::LoadBlockShapes(const TopoDS_Shell& theShell,
const TopoDS_Vertex& theVertex000,
const TopoDS_Vertex& theVertex001,
-// TopTools_IndexedMapOfShape& theShapeIDMap
TopTools_IndexedMapOfOrientedShape& theShapeIDMap )
{
MESSAGE(" ::LoadBlockShapes()");
+ return ( FindBlockShapes( theShell, theVertex000, theVertex001, theShapeIDMap ) &&
+ LoadBlockShapes( theShapeIDMap ));
+}
- myShell = theShell;
- myNbIterations = 0;
- mySumDist = 0;
- myGridComputed = false;
+//=======================================================================
+//function : LoadBlockShapes
+//purpose : add sub-shapes of theBlock to theShapeIDMap so that they get
+// IDs acoording to enum TShapeID
+//=======================================================================
- // 6 vertices
+bool SMESH_Block::FindBlockShapes(const TopoDS_Shell& theShell,
+ const TopoDS_Vertex& theVertex000,
+ const TopoDS_Vertex& theVertex001,
+ TopTools_IndexedMapOfOrientedShape& theShapeIDMap )
+{
+ MESSAGE(" ::FindBlockShapes()");
+
+ // 8 vertices
TopoDS_Shape V000, V100, V010, V110, V001, V101, V011, V111;
// 12 edges
TopoDS_Shape Ex00, Ex10, Ex01, Ex11;
const int NB_FACES_BY_VERTEX = 6;
TopTools_IndexedDataMapOfShapeListOfShape vfMap;
- TopExp::MapShapesAndAncestors( myShell, TopAbs_VERTEX, TopAbs_FACE, vfMap );
+ TopExp::MapShapesAndAncestors( theShell, TopAbs_VERTEX, TopAbs_FACE, vfMap );
if ( vfMap.Extent() != 8 ) {
MESSAGE(" Wrong nb of vertices in the block: " << vfMap.Extent() );
return false;
}
// find vertex 001 - the one on the most vertical edge passing through V000
TopTools_IndexedDataMapOfShapeListOfShape veMap;
- TopExp::MapShapesAndAncestors( myShell, TopAbs_VERTEX, TopAbs_EDGE, veMap );
+ TopExp::MapShapesAndAncestors( theShell, TopAbs_VERTEX, TopAbs_EDGE, veMap );
gp_Vec dir001 = gp::DZ();
gp_Pnt p000 = BRep_Tool::Pnt( TopoDS::Vertex( V000 ));
double maxVal = -DBL_MAX;
// find bottom edges and veritices
list< TopoDS_Edge > eList;
list< int > nbVertexInWires;
- getOrderedEdges( TopoDS::Face( Fxy0 ), TopoDS::Vertex( V000 ), eList, nbVertexInWires );
+ GetOrderedEdges( TopoDS::Face( Fxy0 ), TopoDS::Vertex( V000 ), eList, nbVertexInWires );
if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) {
MESSAGE(" LoadBlockShapes() error ");
return false;
// find top edges and veritices
eList.clear();
- getOrderedEdges( TopoDS::Face( Fxy1 ), TopoDS::Vertex( V001 ), eList, nbVertexInWires );
+ GetOrderedEdges( TopoDS::Face( Fxy1 ), TopoDS::Vertex( V001 ), eList, nbVertexInWires );
if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) {
MESSAGE(" LoadBlockShapes() error ");
return false;
if ( V101.IsSame( exp.Current() ) || V100.IsSame( exp.Current() ))
break; // V101 or V100 found
if ( !exp.More() ) { // not found
- TopoDS_Shape f = Fx0z; Fx0z = F0yz; F0yz = f;
+ std::swap( Fx0z, F0yz);
}
// find Fx1z and F1yz faces
if ( V010.IsSame( exp.Current() ) || V011.IsSame( exp.Current() ))
break;
if ( !exp.More() ) {
- TopoDS_Shape f = Fx1z; Fx1z = F1yz; F1yz = f;
+ std::swap( Fx1z, F1yz);
}
// find vertical edges
theShapeIDMap.Add(F0yz);
theShapeIDMap.Add(F1yz);
- theShapeIDMap.Add(myShell);
+ theShapeIDMap.Add(theShell);
- if ( theShapeIDMap.Extent() != 27 ) {
- MESSAGE("LoadBlockShapes() " << theShapeIDMap.Extent() );
- return false;
- }
+ return true;
+}
+
+//================================================================================
+/*!
+ * \brief Initialize block geometry with shapes from theShapeIDMap
+ * \param theShapeIDMap - map of block subshapes
+ * \retval bool - is a success
+ */
+//================================================================================
+
+bool SMESH_Block::LoadBlockShapes(const TopTools_IndexedMapOfOrientedShape& theShapeIDMap)
+{
+ init();
// store shapes geometry
for ( int shapeID = 1; shapeID < theShapeIDMap.Extent(); shapeID++ )
{
case TopAbs_VERTEX: {
- if ( shapeID > ID_V111 ) {
- MESSAGE(" shapeID =" << shapeID );
- return false;
- }
- myPnt[ shapeID - ID_V000 ] =
- BRep_Tool::Pnt( TopoDS::Vertex( S )).XYZ();
+ if ( !IsVertexID( ID_V111 )) return false;
+ myPnt[ shapeID - ID_V000 ] = BRep_Tool::Pnt( TopoDS::Vertex( S )).XYZ();
break;
}
case TopAbs_EDGE: {
+ if ( !IsEdgeID( shapeID )) return false;
const TopoDS_Edge& edge = TopoDS::Edge( S );
- if ( shapeID < ID_Ex00 || shapeID > ID_E11z || edge.IsNull() ) {
- MESSAGE(" shapeID =" << shapeID );
- return false;
- }
- TEdge& tEdge = myEdge[ shapeID - ID_Ex00 ];
- tEdge.myCoordInd = GetCoordIndOnEdge( shapeID );
- TopLoc_Location loc;
- tEdge.myC3d = BRep_Tool::Curve( edge, loc, tEdge.myFirst, tEdge.myLast );
- if ( !IsForwardEdge( edge, theShapeIDMap ))
- Swap( tEdge.myFirst, tEdge.myLast );
- tEdge.myTrsf = loc;
+ TEdge& tEdge = myEdge[ shapeID - ID_FirstE ];
+ tEdge.Set( shapeID,
+ new BRepAdaptor_Curve( edge ),
+ IsForwardEdge( edge, theShapeIDMap ));
break;
}
case TopAbs_FACE: {
- const TopoDS_Face& face = TopoDS::Face( S );
- if ( shapeID < ID_Fxy0 || shapeID > ID_F1yz || face.IsNull() ) {
- MESSAGE(" shapeID =" << shapeID );
+ if ( !LoadFace( TopoDS::Face( S ), shapeID, theShapeIDMap ))
return false;
- }
- TFace& tFace = myFace[ shapeID - ID_Fxy0 ];
- // pcurves
- vector< int > edgeIdVec(4, -1);
- GetFaceEdgesIDs( shapeID, edgeIdVec );
- for ( int iE = 0; iE < 4; iE++ ) // loop on 4 edges
- {
- const TopoDS_Edge& edge = TopoDS::Edge( theShapeIDMap( edgeIdVec[ iE ]));
- tFace.myCoordInd[ iE ] = GetCoordIndOnEdge( edgeIdVec[ iE ] );
- tFace.myC2d[ iE ] =
- BRep_Tool::CurveOnSurface( edge, face, tFace.myFirst[iE], tFace.myLast[iE] );
- if ( !IsForwardEdge( edge, theShapeIDMap ))
- Swap( tFace.myFirst[ iE ], tFace.myLast[ iE ] );
- }
- // 2d corners
- tFace.myCorner[ 0 ] = tFace.myC2d[ 0 ]->Value( tFace.myFirst[0] ).XY();
- tFace.myCorner[ 1 ] = tFace.myC2d[ 0 ]->Value( tFace.myLast[0] ).XY();
- tFace.myCorner[ 2 ] = tFace.myC2d[ 1 ]->Value( tFace.myLast[1] ).XY();
- tFace.myCorner[ 3 ] = tFace.myC2d[ 1 ]->Value( tFace.myFirst[1] ).XY();
- // sufrace
- TopLoc_Location loc;
- tFace.myS = BRep_Tool::Surface( face, loc );
- tFace.myTrsf = loc;
break;
}
default: break;
return true;
}
+//================================================================================
+/*!
+ * \brief Load face geometry
+ * \param theFace - face
+ * \param theFaceID - face in-block ID
+ * \param theShapeIDMap - map of block subshapes
+ * \retval bool - is a success
+ *
+ * It is enough to compute params or coordinates on the face.
+ * Face subshapes must be loaded into theShapeIDMap before
+ */
+//================================================================================
+
+bool SMESH_Block::LoadFace(const TopoDS_Face& theFace,
+ const int theFaceID,
+ const TopTools_IndexedMapOfOrientedShape& theShapeIDMap)
+{
+ if ( !IsFaceID( theFaceID ) ) return false;
+ // pcurves
+ Adaptor2d_Curve2d* c2d[4];
+ bool isForward[4];
+ vector< int > edgeIdVec;
+ GetFaceEdgesIDs( theFaceID, edgeIdVec );
+ for ( int iE = 0; iE < edgeIdVec.size(); iE++ ) // loop on 4 edges
+ {
+ if ( edgeIdVec[ iE ] > theShapeIDMap.Extent() )
+ return false;
+ const TopoDS_Edge& edge = TopoDS::Edge( theShapeIDMap( edgeIdVec[ iE ]));
+ c2d[ iE ] = new BRepAdaptor_Curve2d( edge, theFace );
+ isForward[ iE ] = IsForwardEdge( edge, theShapeIDMap );
+ }
+ TFace& tFace = myFace[ theFaceID - ID_FirstF ];
+ tFace.Set( theFaceID, new BRepAdaptor_Surface( theFace ), c2d, isForward );
+ return true;
+}
+
+//================================================================================
+/*!
+ * \brief/ Insert theShape into theShapeIDMap with theShapeID
+ * \param theShape - shape to insert
+ * \param theShapeID - shape in-block ID
+ * \param theShapeIDMap - map of block subshapes
+ */
+//================================================================================
+
+bool SMESH_Block::Insert(const TopoDS_Shape& theShape,
+ const int theShapeID,
+ TopTools_IndexedMapOfOrientedShape& theShapeIDMap)
+{
+ if ( !theShape.IsNull() && theShapeID > 0 )
+ {
+ if ( theShapeIDMap.Contains( theShape ))
+ return ( theShapeIDMap.FindIndex( theShape ) == theShapeID );
+
+ if ( theShapeID <= theShapeIDMap.Extent() ) {
+ theShapeIDMap.Substitute( theShapeID, theShape );
+ }
+ else {
+ while ( theShapeIDMap.Extent() < theShapeID - 1 ) {
+ TopoDS_Compound comp;
+ BRep_Builder().MakeCompound( comp );
+ theShapeIDMap.Add( comp );
+ }
+ theShapeIDMap.Add( theShape );
+ }
+ return true;
+ }
+ return false;
+}
+
//=======================================================================
//function : GetFaceEdgesIDs
//purpose : return edges IDs in the order u0, u1, 0v, 1v
void SMESH_Block::GetFaceEdgesIDs (const int faceID, vector< int >& edgeVec )
{
+ edgeVec.resize( 4 );
switch ( faceID ) {
case ID_Fxy0:
edgeVec[ 0 ] = ID_Ex00;
MESSAGE(" GetFaceEdgesIDs(), wrong face ID: " << faceID );
}
}
+
+//=======================================================================
+//function : GetEdgeVertexIDs
+//purpose : return vertex IDs of an edge
+//=======================================================================
+
+void SMESH_Block::GetEdgeVertexIDs (const int edgeID, vector< int >& vertexVec )
+{
+ vertexVec.resize( 2 );
+ switch ( edgeID ) {
+
+ case ID_Ex00:
+ vertexVec[ 0 ] = ID_V000;
+ vertexVec[ 1 ] = ID_V100;
+ break;
+ case ID_Ex10:
+ vertexVec[ 0 ] = ID_V010;
+ vertexVec[ 1 ] = ID_V110;
+ break;
+ case ID_Ex01:
+ vertexVec[ 0 ] = ID_V001;
+ vertexVec[ 1 ] = ID_V101;
+ break;
+ case ID_Ex11:
+ vertexVec[ 0 ] = ID_V011;
+ vertexVec[ 1 ] = ID_V111;
+ break;
+
+ case ID_E0y0:
+ vertexVec[ 0 ] = ID_V000;
+ vertexVec[ 1 ] = ID_V010;
+ break;
+ case ID_E1y0:
+ vertexVec[ 0 ] = ID_V100;
+ vertexVec[ 1 ] = ID_V110;
+ break;
+ case ID_E0y1:
+ vertexVec[ 0 ] = ID_V001;
+ vertexVec[ 1 ] = ID_V011;
+ break;
+ case ID_E1y1:
+ vertexVec[ 0 ] = ID_V101;
+ vertexVec[ 1 ] = ID_V111;
+ break;
+
+ case ID_E00z:
+ vertexVec[ 0 ] = ID_V000;
+ vertexVec[ 1 ] = ID_V001;
+ break;
+ case ID_E10z:
+ vertexVec[ 0 ] = ID_V100;
+ vertexVec[ 1 ] = ID_V101;
+ break;
+ case ID_E01z:
+ vertexVec[ 0 ] = ID_V010;
+ vertexVec[ 1 ] = ID_V011;
+ break;
+ case ID_E11z:
+ vertexVec[ 0 ] = ID_V110;
+ vertexVec[ 1 ] = ID_V111;
+ break;
+ default:
+ vertexVec.resize(0);
+ MESSAGE(" GetEdgeVertexIDs(), wrong edge ID: " << edgeID );
+ }
+}