From: eap Date: Wed, 10 Aug 2011 10:18:40 +0000 (+0000) Subject: 0020511: EDF 1101 SMESH : Add CGNS to Mesh Format Supported X-Git-Tag: V6_4_0a1~124 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=418c46e9629912d2663714e7bd740000628a0c32;p=modules%2Fsmesh.git 0020511: EDF 1101 SMESH : Add CGNS to Mesh Format Supported Move general utils independent of SMESH data structures to SMESHUtils --- diff --git a/src/SMESH/Makefile.am b/src/SMESH/Makefile.am index e6e345333..040c9849f 100644 --- a/src/SMESH/Makefile.am +++ b/src/SMESH/Makefile.am @@ -40,17 +40,10 @@ salomeinclude_HEADERS = \ SMESH_3D_Algo.hxx \ SMESH_Group.hxx \ SMESH_MeshEditor.hxx \ - SMESH_Block.hxx \ SMESH_Pattern.hxx \ - SMESH_TypeDefs.hxx \ SMESH_MesherHelper.hxx \ - SMESH_Octree.hxx \ - SMESH_OctreeNode.hxx \ - SMESH_Comment.hxx \ - SMESH_ComputeError.hxx \ - SMESH_File.hxx \ - SMESH_SMESH.hxx \ - SMESH_ProxyMesh.hxx + SMESH_ProxyMesh.hxx \ + SMESH_SMESH.hxx # Libraries targets @@ -69,14 +62,10 @@ dist_libSMESHimpl_la_SOURCES = \ SMESH_3D_Algo.cxx \ SMESH_Group.cxx \ SMESH_MeshEditor.cxx \ - SMESH_Block.cxx \ SMESH_Pattern.cxx \ SMESH_HypoFilter.cxx \ - SMESH_MesherHelper.cxx \ - SMESH_Octree.cxx \ - SMESH_OctreeNode.cxx \ - SMESH_File.cxx \ - SMESH_ProxyMesh.cxx + SMESH_ProxyMesh.cxx \ + SMESH_MesherHelper.cxx # additionnal information to compile and link file libSMESHimpl_la_CPPFLAGS = \ @@ -93,8 +82,10 @@ libSMESHimpl_la_CPPFLAGS = \ -I$(srcdir)/../DriverMED \ -I$(srcdir)/../DriverUNV \ -I$(srcdir)/../DriverSTL \ + -I$(srcdir)/../DriverCGNS \ -I$(srcdir)/../SMDS \ - -I$(srcdir)/../SMESHDS + -I$(srcdir)/../SMESHDS \ + -I$(srcdir)/../SMESHUtils libSMESHimpl_la_LDFLAGS = \ ../SMESHDS/libSMESHDS.la \ @@ -103,5 +94,7 @@ libSMESHimpl_la_LDFLAGS = \ ../DriverSTL/libMeshDriverSTL.la \ ../DriverMED/libMeshDriverMED.la \ ../DriverUNV/libMeshDriverUNV.la \ + ../DriverCGNS/libMeshDriverCGNS.la \ + ../SMESHUtils/libSMESHUtils.la \ $(GEOM_LDFLAGS) -lNMTTools \ $(CAS_LDPATH) -lTKShHealing -lTKPrim -lTKG2d diff --git a/src/SMESH/SMESH_Block.cxx b/src/SMESH/SMESH_Block.cxx deleted file mode 100644 index b036993a5..000000000 --- a/src/SMESH/SMESH_Block.cxx +++ /dev/null @@ -1,1737 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// File : SMESH_Pattern.hxx -// Created : Mon Aug 2 10:30:00 2004 -// Author : Edward AGAPOV (eap) -// -#include "SMESH_Block.hxx" - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include "SMDS_MeshNode.hxx" -#include "SMDS_MeshVolume.hxx" -#include "SMDS_VolumeTool.hxx" -#include "utilities.h" - -#include - -using namespace std; - -//#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 -//purpose : -//======================================================================= - -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; -} - -//======================================================================= -//function : SMESH_Block::TEdge::Point -//purpose : -//======================================================================= - -gp_XYZ SMESH_Block::TEdge::Point( const gp_XYZ& theParams ) const -{ - 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); - } -} - -//======================================================================= -//function : SMESH_Block::TFace::GetUV -//purpose : -//======================================================================= - -gp_XY SMESH_Block::TFace::GetUV( const gp_XYZ& theParams ) const -{ - gp_XY uv(0.,0.); - 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 ] ); - u = ( 1 - u ) * myFirst[ iE ] + u * myLast[ iE ]; - uv += Ecoef * myC2d[ iE ]->Value( u ).XY(); - // corner addition - uv -= Vcoef * myCorner[ iE ]; - } - return uv; -} - -//======================================================================= -//function : SMESH_Block::TFace::Point -//purpose : -//======================================================================= - -gp_XYZ SMESH_Block::TFace::Point( const gp_XYZ& theParams ) const -{ - 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; -} - -//======================================================================= -//function : GetShapeCoef -//purpose : -//======================================================================= - -double* SMESH_Block::GetShapeCoef (const int theShapeID) -{ - static double shapeCoef[][3] = { - // V000, V100, V010, V110 - { -1,-1,-1 }, { 1,-1,-1 }, { -1, 1,-1 }, { 1, 1,-1 }, - // V001, V101, V011, V111, - { -1,-1, 1 }, { 1,-1, 1 }, { -1, 1, 1 }, { 1, 1, 1 }, - // Ex00, Ex10, Ex01, Ex11, - { 0,-1,-1 }, { 0, 1,-1 }, { 0,-1, 1 }, { 0, 1, 1 }, - // E0y0, E1y0, E0y1, E1y1, - { -1, 0,-1 }, { 1, 0,-1 }, { -1, 0, 1 }, { 1, 0, 1 }, - // E00z, E10z, E01z, E11z, - { -1,-1, 0 }, { 1,-1, 0 }, { -1, 1, 0 }, { 1, 1, 0 }, - // Fxy0, Fxy1, Fx0z, Fx1z, F0yz, F1yz, - { 0, 0,-1 }, { 0, 0, 1 }, { 0,-1, 0 }, { 0, 1, 0 }, { -1, 0, 0 }, { 1, 0, 0 }, - // ID_Shell - { 0, 0, 0 } - }; - if ( theShapeID < ID_V000 || theShapeID > ID_F1yz ) - return shapeCoef[ ID_Shell - 1 ]; - - return shapeCoef[ theShapeID - 1 ]; -} - -//======================================================================= -//function : ShellPoint -//purpose : return coordinates of a point in shell -//======================================================================= - -bool SMESH_Block::ShellPoint( const gp_XYZ& theParams, gp_XYZ& thePoint ) const -{ - thePoint.SetCoord( 0., 0., 0. ); - for ( int shapeID = ID_V000; shapeID < ID_Shell; shapeID++ ) - { - // coef - double* coefs = GetShapeCoef( shapeID ); - double k = 1; - for ( int iCoef = 0; iCoef < 3; iCoef++ ) { - if ( coefs[ iCoef ] != 0 ) { - if ( coefs[ iCoef ] < 0 ) - k *= ( 1. - theParams.Coord( iCoef + 1 )); - else - k *= theParams.Coord( iCoef + 1 ); - } - } - // 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& thePointOnShape, - gp_XYZ& thePoint ) -{ - if ( thePointOnShape.size() < ID_F1yz ) - return false; - - const double x = theParams.X(), y = theParams.Y(), z = theParams.Z(); - const double x1 = 1. - x, y1 = 1. - y, z1 = 1. - z; - const vector& 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 : -//======================================================================= - -Standard_Integer SMESH_Block::NbVariables() const -{ - return 3; -} - -//======================================================================= -//function : NbEquations -//purpose : -//======================================================================= - -Standard_Integer SMESH_Block::NbEquations() const -{ - return 1; -} - -//======================================================================= -//function : Value -//purpose : -//======================================================================= - -Standard_Boolean SMESH_Block::Value(const math_Vector& theXYZ, math_Vector& theFxyz) -{ - gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) ); - if ( params.IsEqual( myParam, DBL_MIN )) { // same param - theFxyz( 1 ) = funcValue( myValues[ SQUARE_DIST ]); - } - else { - ShellPoint( params, P ); - gp_Vec dP( P - myPoint ); - theFxyz(1) = funcValue( dP.SquareMagnitude() ); - } - return true; -} - -//======================================================================= -//function : Derivatives -//purpose : -//======================================================================= - -Standard_Boolean SMESH_Block::Derivatives(const math_Vector& XYZ,math_Matrix& Df) -{ - 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 : -//======================================================================= - -Standard_Boolean SMESH_Block::Values(const math_Vector& theXYZ, - math_Vector& theFxyz, - math_Matrix& theDf) -{ - gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) ); - if ( params.IsEqual( myParam, DBL_MIN )) { // same param - 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 - MESSAGE ( "PARAM GUESS: " << params.X() << " "<< params.Y() << " "<< params.X() ); - myNbIterations++; // how many times call ShellPoint() -#endif - ShellPoint( params, P ); - - 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 ( sqDist < myValues[ SQUARE_DIST ] ) // a better guess - { - // 3 partial derivatives - 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; - 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 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 { - double Dis = plnNorm * P - plnNorm * myPoint; - theDf( 1, iP + 1 ) = Dis/Direc; - } -#endif - } -#ifdef DEBUG_PARAM_COMPUTE - MESSAGE ( "F = " << theFxyz(1) << " DRV: " << theDf(1,1) << " " << theDf(1,2) << " " << theDf(1,3) ); - myNbIterations +=3; // how many times call ShellPoint() -#endif - - // store better values - 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 - - double loopTol = 10 * myTolerance; - int nbLoops = 0; - while ( distance() > loopTol && nbLoops <= 3 ) - { - paramSearch.Perform ( *static_cast(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(); - MESSAGE ( " ------ SOLUTION: ( "<< myParam.X() <<" "<< myParam.Y() <<" "<< myParam.Z() <<" )"< 0 ) - theParams.SetCoord( myFaceIndex, myFaceParam ); - - return true; -} - -//======================================================================= -//function : ComputeParameters -//purpose : compute point parameters in the block -//======================================================================= - -bool SMESH_Block::ComputeParameters(const gp_Pnt& thePoint, - gp_XYZ& theParams, - const int theShapeID, - const gp_XYZ& theParamsHint) -{ - if ( VertexParameters( theShapeID, theParams )) - return true; - - 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 ); - - // 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 - bool needGrid = false; - gp_XYZ par000( 0, 0, 0 ), par111( 1, 1, 1 ); - double zero = DBL_MIN * DBL_MIN; - for ( int iEdge = 0, iParam = 1; iParam <= 3 && !needGrid; iParam++ ) - { - if ( isOnFace && coef[ iParam - 1 ] != 0 ) { - 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 ); - gp_Vec v01( p0, p1 ), v0P( p0, thePoint ); - double len2 = v01.SquareMagnitude(); - double par = 0; - if ( len2 > zero ) { - par = v0P.Dot( v01 ) / len2; - if ( par < 0 || par > 1 ) { // projection falls out of line ends => needGrid - needGrid = true; - break; - } - } - sumParam += par; - } - 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 ( hasHint ) - { - start = theParamsHint; - } - else if ( myGridComputed ) - { - double minDist = DBL_MAX; - gp_XYZ* bestParam = 0; - for ( int iNode = 0; iNode < 27; iNode++ ) { - TxyzPair & prmPtn = my3x3x3GridNodes[ iNode ]; - double dist = ( thePoint.XYZ() - prmPtn.second ).SquareModulus(); - if ( dist < minDist ) { - minDist = dist; - bestParam = & prmPtn.first; - } - } - start = *bestParam; - } - - myFaceIndex = -1; - myFaceParam = 0.; - if ( isOnFace ) { - // put a point on the face - for ( int iCoord = 0; iCoord < 3; iCoord++ ) - if ( coef[ iCoord ] ) { - myFaceIndex = iCoord + 1; - myFaceParam = ( coef[ iCoord ] < 0.5 ) ? 0.0 : 1.0; - start.SetCoord( myFaceIndex, myFaceParam ); - } - } - -#ifdef DEBUG_PARAM_COMPUTE - MESSAGE ( " #### POINT " < sqDistance ) { // solution get worse - if ( ++nbGetWorst > 2 ) - return computeParameters( thePoint, theParams, solution ); - } -#ifdef DEBUG_PARAM_COMPUTE - MESSAGE ( "PARAMS: ( " << params.X() <<" "<< params.Y() <<" "<< params.Z() <<" )" ); - MESSAGE ( "DIST: " << sqrt( sqDist ) ); -#endif - - if ( sqDist < sqDistance ) { // get better - sqDistance = sqDist; - solution = params; - nbGetWorst = 0; - if ( sqDistance < sqTolerance ) // a solution found - break; - } - - // 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++; - } -#ifdef DEBUG_PARAM_COMPUTE - myNbIterations += nbLoops*4; // how many times ShellPoint called - mySumDist += sqrt( sqDistance ); - MESSAGE ( " ------ SOLUTION: ( "< 0 ) - theParams.SetCoord( myFaceIndex, myFaceParam ); - - return true; -} - -//======================================================================= -//function : VertexParameters -//purpose : return parameters of a vertex given by TShapeID -//======================================================================= - -bool SMESH_Block::VertexParameters(const int theVertexID, gp_XYZ& theParams) -{ - 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; -} - -//======================================================================= -//function : DumpShapeID -//purpose : debug an id of a block sub-shape -//======================================================================= - -#define CASEDUMP(id,strm) case id: strm << #id; break; - -ostream& SMESH_Block::DumpShapeID (const int id, ostream& stream) -{ - switch ( id ) { - CASEDUMP( ID_V000, stream ); - CASEDUMP( ID_V100, stream ); - CASEDUMP( ID_V010, stream ); - CASEDUMP( ID_V110, stream ); - CASEDUMP( ID_V001, stream ); - CASEDUMP( ID_V101, stream ); - CASEDUMP( ID_V011, stream ); - CASEDUMP( ID_V111, stream ); - CASEDUMP( ID_Ex00, stream ); - CASEDUMP( ID_Ex10, stream ); - CASEDUMP( ID_Ex01, stream ); - CASEDUMP( ID_Ex11, stream ); - CASEDUMP( ID_E0y0, stream ); - CASEDUMP( ID_E1y0, stream ); - CASEDUMP( ID_E0y1, stream ); - CASEDUMP( ID_E1y1, stream ); - CASEDUMP( ID_E00z, stream ); - CASEDUMP( ID_E10z, stream ); - CASEDUMP( ID_E01z, stream ); - CASEDUMP( ID_E11z, stream ); - CASEDUMP( ID_Fxy0, stream ); - CASEDUMP( ID_Fxy1, stream ); - CASEDUMP( ID_Fx0z, stream ); - CASEDUMP( ID_Fx1z, stream ); - CASEDUMP( ID_F0yz, stream ); - CASEDUMP( ID_F1yz, stream ); - CASEDUMP( ID_Shell, stream ); - default: stream << "ID_INVALID"; - } - return stream; -} - -//======================================================================= -//function : GetShapeIDByParams -//purpose : define an id of the block sub-shape by normlized point coord -//======================================================================= - -int SMESH_Block::GetShapeIDByParams ( const gp_XYZ& theCoord ) -{ - // id ( 0 - 26 ) computation: - - // vertex ( 0 - 7 ) : id = 1*x + 2*y + 4*z - - // edge || X ( 8 - 11 ) : id = 8 + 1*y + 2*z - // edge || Y ( 12 - 15 ): id = 1*x + 12 + 2*z - // edge || Z ( 16 - 19 ): id = 1*x + 2*y + 16 - - // face || XY ( 20 - 21 ): id = 8 + 12 + 1*z - 0 - // face || XZ ( 22 - 23 ): id = 8 + 1*y + 16 - 2 - // face || YZ ( 24 - 25 ): id = 1*x + 12 + 16 - 4 - - static int iAddBnd[] = { 1, 2, 4 }; - static int iAddNotBnd[] = { 8, 12, 16 }; - static int iFaceSubst[] = { 0, 2, 4 }; - - int id = 0; - int iOnBoundary = 0; - for ( int iCoord = 0; iCoord < 3; iCoord++ ) - { - double val = theCoord.Coord( iCoord + 1 ); - if ( val == 0.0 ) - iOnBoundary++; - else if ( val == 1.0 ) - id += iAddBnd[ iOnBoundary++ ]; - else - id += iAddNotBnd[ iCoord ]; - } - if ( iOnBoundary == 1 ) // face - id -= iFaceSubst[ (id - 20) / 4 ]; - else if ( iOnBoundary == 0 ) // shell - id = 26; - - if ( id > 26 || id < 0 ) { - MESSAGE( "GetShapeIDByParams() = " << id - <<" "<< theCoord.X() <<" "<< theCoord.Y() <<" "<< theCoord.Z() ); - } - - return id + 1; // shape ids start at 1 -} - -//================================================================================ -/*! - * \brief Return number of wires and a list of oredered edges. - * \param theFace - the face to process - * \param theFirstVertex - the vertex of the outer wire to set first in the returned - * list ( theFirstVertex may be NULL ) - * \param theEdges - all ordered edges of theFace (outer edges goes first). - * \param theNbEdgesInWires - nb of edges (== nb of vertices in closed wire) in each wire - * \param theShapeAnalysisAlgo - if true, ShapeAnalysis::OuterWire() is used to find - * the outer wire else BRepTools::OuterWire() is used. - * \retval int - nb of wires - * - * Always try to set a seam edge first. - * BRepTools::OuterWire() fails e.g. in the case of issue 0020184, - * ShapeAnalysis::OuterWire() fails in the case of issue 0020452 - */ -//================================================================================ - -int SMESH_Block::GetOrderedEdges (const TopoDS_Face& theFace, - TopoDS_Vertex theFirstVertex, - list< TopoDS_Edge >& theEdges, - list< int > & theNbEdgesInWires, - const bool theShapeAnalysisAlgo) -{ - // put wires in a list, so that an outer wire comes first - list aWireList; - TopoDS_Wire anOuterWire = - theShapeAnalysisAlgo ? ShapeAnalysis::OuterWire( theFace ) : BRepTools::OuterWire( theFace ); - for ( TopoDS_Iterator wIt (theFace); wIt.More(); wIt.Next() ) - if ( wIt.Value().ShapeType() == TopAbs_WIRE ) // it can be internal vertex! - { - if ( !anOuterWire.IsSame( wIt.Value() )) - aWireList.push_back( TopoDS::Wire( wIt.Value() )); - else - aWireList.push_front( TopoDS::Wire( wIt.Value() )); - } - - // loop on edges of wires - theNbEdgesInWires.clear(); - list::iterator wlIt = aWireList.begin(); - for ( ; wlIt != aWireList.end(); wlIt++ ) - { - int iE; - BRepTools_WireExplorer wExp( *wlIt, theFace ); - for ( iE = 0; wExp.More(); wExp.Next(), iE++ ) - { - TopoDS_Edge edge = wExp.Current(); - // commented for issue 0020557, other related ones: 0020526, PAL19080 - // edge = TopoDS::Edge( edge.Oriented( wExp.Orientation() )); - theEdges.push_back( edge ); - } - if ( iE == 0 ) // wExp returns nothing if e.g. the wire contains one internal edge - { // Issue 0020676 - for ( TopoDS_Iterator e( *wlIt ); e.More(); e.Next(), ++iE ) - theEdges.push_back( TopoDS::Edge( e.Value() )); - } - theNbEdgesInWires.push_back( iE ); - iE = 0; - if ( wlIt == aWireList.begin() && theEdges.size() > 1 ) { // the outer wire - // orient closed edges - list< TopoDS_Edge >::iterator eIt, eIt2; - for ( eIt = theEdges.begin(); eIt != theEdges.end(); eIt++ ) - { - TopoDS_Edge& edge = *eIt; - if ( TopExp::FirstVertex( edge ).IsSame( TopExp::LastVertex( edge ) )) - { - eIt2 = eIt; - bool isNext = ( eIt2 == theEdges.begin() ); - TopoDS_Edge edge2 = isNext ? *(++eIt2) : *(--eIt2); - double f1,l1,f2,l2; - Handle(Geom2d_Curve) c1 = BRep_Tool::CurveOnSurface( edge, theFace, f1,l1 ); - Handle(Geom2d_Curve) c2 = BRep_Tool::CurveOnSurface( edge2, theFace, f2,l2 ); - gp_Pnt2d pf = c1->Value( edge.Orientation() == TopAbs_FORWARD ? f1 : l1 ); - gp_Pnt2d pl = c1->Value( edge.Orientation() == TopAbs_FORWARD ? l1 : f1 ); - bool isFirst = ( edge2.Orientation() == TopAbs_FORWARD ? isNext : !isNext ); - gp_Pnt2d p2 = c2->Value( isFirst ? f2 : l2 ); - 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() ) { - 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()); - TopExp::Vertices( theEdges.front(), vv[0], vv[1], true ); - if ( iE++ > theNbEdgesInWires.back() ) { -#ifdef _DEBUG_ - gp_Pnt p = BRep_Tool::Pnt( theFirstVertex ); - MESSAGE ( " : Warning : vertex "<< theFirstVertex.TShape().operator->() - << " ( " << p.X() << " " << p.Y() << " " << p.Z() << " )" - << " not found in outer wire of face "<< theFace.TShape().operator->() - << " with vertices: " ); - 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 ); - MESSAGE_ADD ( i << " " << v.TShape().operator->() << " " - << p.X() << " " << p.Y() << " " << p.Z() << " " << std::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& 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 vFxy0, vFxy1; - - V000 = theNode000Index; - V001 = theNode001Index; - - // get faces sharing V000 and V001 - list 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::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 : 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_IndexedMapOfOrientedShape& theShapeIDMap ) -{ - MESSAGE(" ::LoadBlockShapes()"); - return ( FindBlockShapes( theShell, theVertex000, theVertex001, theShapeIDMap ) && - LoadBlockShapes( theShapeIDMap )); -} - -//======================================================================= -//function : LoadBlockShapes -//purpose : add sub-shapes of theBlock to theShapeIDMap so that they get -// IDs acoording to enum TShapeID -//======================================================================= - -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; - TopoDS_Shape E0y0, E1y0, E0y1, E1y1; - TopoDS_Shape E00z, E10z, E01z, E11z; - // 6 faces - TopoDS_Shape Fxy0, Fx0z, F0yz, Fxy1, Fx1z, F1yz; - - // nb of faces bound to a vertex in TopTools_IndexedDataMapOfShapeListOfShape - // filled by TopExp::MapShapesAndAncestors() - const int NB_FACES_BY_VERTEX = 6; - - TopTools_IndexedDataMapOfShapeListOfShape 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; - } - - V000 = theVertex000; - V001 = theVertex001; - - if ( V000.IsNull() ) { - // find vertex 000 - the one with smallest coordinates - double minVal = DBL_MAX, minX, val; - for ( int i = 1; i <= 8; i++ ) { - const TopoDS_Vertex& v = TopoDS::Vertex( vfMap.FindKey( i )); - gp_Pnt P = BRep_Tool::Pnt( v ); - val = P.X() + P.Y() + P.Z(); - if ( val < minVal || ( val == minVal && P.X() < minX )) { - V000 = v; - minVal = val; - minX = P.X(); - } - } - // find vertex 001 - the one on the most vertical edge passing through V000 - TopTools_IndexedDataMapOfShapeListOfShape 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; - TopTools_ListIteratorOfListOfShape eIt ( veMap.FindFromKey( V000 )); - for ( ; eIt.More(); eIt.Next() ) { - const TopoDS_Edge& e = TopoDS::Edge( eIt.Value() ); - TopoDS_Vertex v = TopExp::FirstVertex( e ); - if ( v.IsSame( V000 )) - v = TopExp::LastVertex( e ); - val = dir001 * gp_Vec( p000, BRep_Tool::Pnt( v )).Normalized(); - if ( val > maxVal ) { - V001 = v; - maxVal = val; - } - } - } - - // find the bottom (Fxy0), Fx0z and F0yz faces - - const TopTools_ListOfShape& f000List = vfMap.FindFromKey( V000 ); - const TopTools_ListOfShape& f001List = vfMap.FindFromKey( V001 ); - if (f000List.Extent() != NB_FACES_BY_VERTEX || - f001List.Extent() != NB_FACES_BY_VERTEX ) { - MESSAGE(" LoadBlockShapes() " << f000List.Extent() << " " << f001List.Extent()); - return false; - } - TopTools_ListIteratorOfListOfShape f001It, f000It ( f000List ); - int i, j, iFound1, iFound2; - for ( j = 0; f000It.More(); f000It.Next(), j++ ) - { - if ( NB_FACES_BY_VERTEX == 6 && j % 2 ) continue; // each face encounters twice - const TopoDS_Shape& F = f000It.Value(); - for ( i = 0, f001It.Initialize( f001List ); f001It.More(); f001It.Next(), i++ ) { - if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice - if ( F.IsSame( f001It.Value() )) - break; - } - if ( f001It.More() ) // Fx0z or F0yz found - if ( Fx0z.IsNull() ) { - Fx0z = F; - iFound1 = i; - } else { - F0yz = F; - iFound2 = i; - } - else // F is the bottom face - Fxy0 = F; - } - if ( Fxy0.IsNull() || Fx0z.IsNull() || F0yz.IsNull() ) { - MESSAGE( Fxy0.IsNull() <<" "<< Fx0z.IsNull() <<" "<< F0yz.IsNull() ); - return false; - } - - // choose the top face (Fxy1) - for ( i = 0, f001It.Initialize( f001List ); f001It.More(); f001It.Next(), i++ ) { - if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice - if ( i != iFound1 && i != iFound2 ) - break; - } - Fxy1 = f001It.Value(); - if ( Fxy1.IsNull() ) { - MESSAGE(" LoadBlockShapes() error "); - return false; - } - - // find bottom edges and veritices - list< TopoDS_Edge > eList; - list< int > nbVertexInWires; - GetOrderedEdges( TopoDS::Face( Fxy0 ), TopoDS::Vertex( V000 ), eList, nbVertexInWires ); - if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) { - MESSAGE(" LoadBlockShapes() error "); - return false; - } - list< TopoDS_Edge >::iterator elIt = eList.begin(); - for ( i = 0; elIt != eList.end(); elIt++, i++ ) - switch ( i ) { - case 0: E0y0 = *elIt; V010 = TopExp::LastVertex( *elIt, true ); break; - case 1: Ex10 = *elIt; V110 = TopExp::LastVertex( *elIt, true ); break; - case 2: E1y0 = *elIt; V100 = TopExp::LastVertex( *elIt, true ); break; - case 3: Ex00 = *elIt; break; - default:; - } - if ( i != 4 || E0y0.IsNull() || Ex10.IsNull() || E1y0.IsNull() || Ex00.IsNull() ) { - MESSAGE(" LoadBlockShapes() error, eList.size()=" << eList.size()); - return false; - } - - - // find top edges and veritices - eList.clear(); - GetOrderedEdges( TopoDS::Face( Fxy1 ), TopoDS::Vertex( V001 ), eList, nbVertexInWires ); - if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) { - MESSAGE(" LoadBlockShapes() error "); - return false; - } - for ( i = 0, elIt = eList.begin(); elIt != eList.end(); elIt++, i++ ) - switch ( i ) { - case 0: Ex01 = *elIt; V101 = TopExp::LastVertex( *elIt, true ); break; - case 1: E1y1 = *elIt; V111 = TopExp::LastVertex( *elIt, true ); break; - case 2: Ex11 = *elIt; V011 = TopExp::LastVertex( *elIt, true ); break; - case 3: E0y1 = *elIt; break; - default:; - } - if ( i != 4 || Ex01.IsNull() || E1y1.IsNull() || Ex11.IsNull() || E0y1.IsNull() ) { - MESSAGE(" LoadBlockShapes() error, eList.size()=" << eList.size()); - return false; - } - - // swap Fx0z and F0yz if necessary - TopExp_Explorer exp( Fx0z, TopAbs_VERTEX ); - for ( ; exp.More(); exp.Next() ) // Fx0z shares V101 and V100 - if ( V101.IsSame( exp.Current() ) || V100.IsSame( exp.Current() )) - break; // V101 or V100 found - if ( !exp.More() ) { // not found - std::swap( Fx0z, F0yz); - } - - // find Fx1z and F1yz faces - const TopTools_ListOfShape& f111List = vfMap.FindFromKey( V111 ); - const TopTools_ListOfShape& f110List = vfMap.FindFromKey( V110 ); - if (f111List.Extent() != NB_FACES_BY_VERTEX || - f110List.Extent() != NB_FACES_BY_VERTEX ) { - MESSAGE(" LoadBlockShapes() " << f111List.Extent() << " " << f110List.Extent()); - return false; - } - TopTools_ListIteratorOfListOfShape f111It, f110It ( f110List); - for ( j = 0 ; f110It.More(); f110It.Next(), j++ ) { - if ( NB_FACES_BY_VERTEX == 6 && j % 2 ) continue; // each face encounters twice - const TopoDS_Shape& F = f110It.Value(); - for ( i = 0, f111It.Initialize( f111List ); f111It.More(); f111It.Next(), i++ ) { - if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice - if ( F.IsSame( f111It.Value() )) { // Fx1z or F1yz found - if ( Fx1z.IsNull() ) - Fx1z = F; - else - F1yz = F; - } - } - } - if ( Fx1z.IsNull() || F1yz.IsNull() ) { - MESSAGE(" LoadBlockShapes() error "); - return false; - } - - // swap Fx1z and F1yz if necessary - for ( exp.Init( Fx1z, TopAbs_VERTEX ); exp.More(); exp.Next() ) - if ( V010.IsSame( exp.Current() ) || V011.IsSame( exp.Current() )) - break; - if ( !exp.More() ) { - std::swap( Fx1z, F1yz); - } - - // find vertical edges - for ( exp.Init( Fx0z, TopAbs_EDGE ); exp.More(); exp.Next() ) { - const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() ); - const TopoDS_Shape& vFirst = TopExp::FirstVertex( edge, true ); - if ( vFirst.IsSame( V001 )) - E00z = edge; - else if ( vFirst.IsSame( V100 )) - E10z = edge; - } - if ( E00z.IsNull() || E10z.IsNull() ) { - MESSAGE(" LoadBlockShapes() error "); - return false; - } - for ( exp.Init( Fx1z, TopAbs_EDGE ); exp.More(); exp.Next() ) { - const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() ); - const TopoDS_Shape& vFirst = TopExp::FirstVertex( edge, true ); - if ( vFirst.IsSame( V111 )) - E11z = edge; - else if ( vFirst.IsSame( V010 )) - E01z = edge; - } - if ( E01z.IsNull() || E11z.IsNull() ) { - MESSAGE(" LoadBlockShapes() error "); - return false; - } - - // load shapes in theShapeIDMap - - theShapeIDMap.Clear(); - - theShapeIDMap.Add(V000.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V100.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V010.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V110.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V001.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V101.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V011.Oriented( TopAbs_FORWARD )); - theShapeIDMap.Add(V111.Oriented( TopAbs_FORWARD )); - - theShapeIDMap.Add(Ex00); - theShapeIDMap.Add(Ex10); - theShapeIDMap.Add(Ex01); - theShapeIDMap.Add(Ex11); - - theShapeIDMap.Add(E0y0); - theShapeIDMap.Add(E1y0); - theShapeIDMap.Add(E0y1); - theShapeIDMap.Add(E1y1); - - theShapeIDMap.Add(E00z); - theShapeIDMap.Add(E10z); - theShapeIDMap.Add(E01z); - theShapeIDMap.Add(E11z); - - theShapeIDMap.Add(Fxy0); - theShapeIDMap.Add(Fxy1); - theShapeIDMap.Add(Fx0z); - theShapeIDMap.Add(Fx1z); - theShapeIDMap.Add(F0yz); - theShapeIDMap.Add(F1yz); - - theShapeIDMap.Add(theShell); - - 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++ ) - { - const TopoDS_Shape& S = theShapeIDMap( shapeID ); - switch ( S.ShapeType() ) - { - case TopAbs_VERTEX: { - - 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 ); - TEdge& tEdge = myEdge[ shapeID - ID_FirstE ]; - tEdge.Set( shapeID, - new BRepAdaptor_Curve( edge ), - IsForwardEdge( edge, theShapeIDMap )); - break; - } - case TopAbs_FACE: { - - if ( !LoadFace( TopoDS::Face( S ), shapeID, theShapeIDMap )) - return false; - break; - } - default: break; - } - } // loop on shapes in theShapeIDMap - - 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 -// u0 means "|| u, v == 0" -//======================================================================= - -void SMESH_Block::GetFaceEdgesIDs (const int faceID, vector< int >& edgeVec ) -{ - edgeVec.resize( 4 ); - switch ( faceID ) { - case ID_Fxy0: - edgeVec[ 0 ] = ID_Ex00; - edgeVec[ 1 ] = ID_Ex10; - edgeVec[ 2 ] = ID_E0y0; - edgeVec[ 3 ] = ID_E1y0; - break; - case ID_Fxy1: - edgeVec[ 0 ] = ID_Ex01; - edgeVec[ 1 ] = ID_Ex11; - edgeVec[ 2 ] = ID_E0y1; - edgeVec[ 3 ] = ID_E1y1; - break; - case ID_Fx0z: - edgeVec[ 0 ] = ID_Ex00; - edgeVec[ 1 ] = ID_Ex01; - edgeVec[ 2 ] = ID_E00z; - edgeVec[ 3 ] = ID_E10z; - break; - case ID_Fx1z: - edgeVec[ 0 ] = ID_Ex10; - edgeVec[ 1 ] = ID_Ex11; - edgeVec[ 2 ] = ID_E01z; - edgeVec[ 3 ] = ID_E11z; - break; - case ID_F0yz: - edgeVec[ 0 ] = ID_E0y0; - edgeVec[ 1 ] = ID_E0y1; - edgeVec[ 2 ] = ID_E00z; - edgeVec[ 3 ] = ID_E01z; - break; - case ID_F1yz: - edgeVec[ 0 ] = ID_E1y0; - edgeVec[ 1 ] = ID_E1y1; - edgeVec[ 2 ] = ID_E10z; - edgeVec[ 3 ] = ID_E11z; - break; - default: - 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 ); - } -} diff --git a/src/SMESH/SMESH_Block.hxx b/src/SMESH/SMESH_Block.hxx deleted file mode 100644 index 3771907ec..000000000 --- a/src/SMESH/SMESH_Block.hxx +++ /dev/null @@ -1,391 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// File : SMESH_Block.hxx -// Created : Tue Nov 30 12:42:18 2004 -// Author : Edward AGAPOV (eap) -// -#ifndef SMESH_Block_HeaderFile -#define SMESH_Block_HeaderFile - -#include "SMESH_SMESH.hxx" - -//#include -//#include -//#include - -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include -#include - -class SMDS_MeshVolume; -class SMDS_MeshNode; -class Adaptor3d_Surface; -class Adaptor2d_Curve2d; -class Adaptor3d_Curve; -class gp_Pnt; - -// ========================================================= -// class calculating coordinates of 3D points by normalized -// parameters inside the block and vice versa -// ========================================================= - -class SMESH_EXPORT SMESH_Block: public math_FunctionSetWithDerivatives -{ - public: - enum TShapeID { - // ---------------------------- - // Ids of the block sub-shapes - // ---------------------------- - ID_NONE = 0, - - ID_V000 = 1, ID_V100, ID_V010, ID_V110, ID_V001, ID_V101, ID_V011, ID_V111, - - ID_Ex00, ID_Ex10, ID_Ex01, ID_Ex11, - ID_E0y0, ID_E1y0, ID_E0y1, ID_E1y1, - ID_E00z, ID_E10z, ID_E01z, ID_E11z, - - ID_Fxy0, ID_Fxy1, ID_Fx0z, ID_Fx1z, ID_F0yz, ID_F1yz, - - ID_Shell - }; - enum { // to use TShapeID for indexing certain type subshapes - - ID_FirstV = ID_V000, ID_FirstE = ID_Ex00, ID_FirstF = ID_Fxy0 - - }; - - - public: - // ------------------------------------------------- - // Block topology in terms of block sub-shapes' ids - // ------------------------------------------------- - - static int NbVertices() { return 8; } - static int NbEdges() { return 12; } - static int NbFaces() { return 6; } - static int NbSubShapes() { return ID_Shell; } - // to avoid magic numbers when allocating memory for subshapes - - static inline bool IsVertexID( int theShapeID ) - { return ( theShapeID >= ID_V000 && theShapeID <= ID_V111 ); } - - static inline bool IsEdgeID( int theShapeID ) - { return ( theShapeID >= ID_Ex00 && theShapeID <= ID_E11z ); } - - static inline bool IsFaceID( int theShapeID ) - { return ( theShapeID >= ID_Fxy0 && theShapeID <= ID_F1yz ); } - - static int ShapeIndex( int theShapeID ) - { - if ( IsVertexID( theShapeID )) return theShapeID - ID_V000; - if ( IsEdgeID( theShapeID )) return theShapeID - ID_Ex00; - if ( IsFaceID( theShapeID )) return theShapeID - ID_Fxy0; - return 0; - } - // return index [0-...] for each type of sub-shapes, - // for example : - // ShapeIndex( ID_Ex00 ) == 0 - // ShapeIndex( ID_Ex10 ) == 1 - - static void GetFaceEdgesIDs (const int faceID, std::vector< int >& edgeVec ); - // return edges IDs of a face in the order u0, u1, 0v, 1v - - static void GetEdgeVertexIDs (const int edgeID, std::vector< int >& vertexVec ); - // return vertex IDs of an edge - - static int GetCoordIndOnEdge (const int theEdgeID) - { return (theEdgeID < ID_E0y0) ? 1 : (theEdgeID < ID_E00z) ? 2 : 3; } - // return an index of a coordinate which varies along the edge - - static double* GetShapeCoef (const int theShapeID); - // for theShapeID( TShapeID ), returns 3 coefficients used - // to compute an addition of an on-theShape point to coordinates - // of an in-shell point. If an in-shell point has parameters (Px,Py,Pz), - // then the addition of a point P is computed as P*kx*ky*kz and ki is - // defined by the returned coef like this: - // ki = (coef[i] == 0) ? 1 : (coef[i] < 0) ? 1 - Pi : Pi - - static int GetShapeIDByParams ( const gp_XYZ& theParams ); - // define an id of the block sub-shape by point parameters - - static std::ostream& DumpShapeID (const int theBlockShapeID, std::ostream& stream); - // DEBUG: dump an id of a block sub-shape - - - public: - // --------------- - // Initialization - // --------------- - - SMESH_Block(); - - bool LoadBlockShapes(const TopoDS_Shell& theShell, - const TopoDS_Vertex& theVertex000, - const TopoDS_Vertex& theVertex001, - TopTools_IndexedMapOfOrientedShape& theShapeIDMap ); - // Initialize block geometry with theShell, - // add sub-shapes of theBlock to theShapeIDMap so that they get - // IDs acoording to enum TShapeID - - bool LoadBlockShapes(const TopTools_IndexedMapOfOrientedShape& theShapeIDMap); - // Initialize block geometry with shapes from theShapeIDMap - - bool LoadMeshBlock(const SMDS_MeshVolume* theVolume, - const int theNode000Index, - const int theNode001Index, - std::vector& theOrderedNodes); - // prepare to work with theVolume and - // return nodes in theVolume corners in the order of TShapeID enum - - bool LoadFace(const TopoDS_Face& theFace, - const int theFaceID, - const TopTools_IndexedMapOfOrientedShape& theShapeIDMap); - // Load face geometry. - // It is enough to compute params or coordinates on the face. - // Face subshapes must be loaded into theShapeIDMap before - - static bool Insert(const TopoDS_Shape& theShape, - const int theShapeID, - TopTools_IndexedMapOfOrientedShape& theShapeIDMap); - // Insert theShape into theShapeIDMap with theShapeID, - // Not yet set shapes preceding theShapeID are filled with compounds - // Return true if theShape was successfully bound to theShapeID - - static bool FindBlockShapes(const TopoDS_Shell& theShell, - const TopoDS_Vertex& theVertex000, - const TopoDS_Vertex& theVertex001, - TopTools_IndexedMapOfOrientedShape& theShapeIDMap ); - // add sub-shapes of theBlock to theShapeIDMap so that they get - // IDs acoording to enum TShapeID - -public: - // --------------------------------- - // Define coordinates by parameters - // --------------------------------- - - bool VertexPoint( const int theVertexID, gp_XYZ& thePoint ) const { - if ( !IsVertexID( theVertexID )) return false; - thePoint = myPnt[ theVertexID - ID_FirstV ]; return true; - } - // return vertex coordinates, parameters are defined by theVertexID - - bool EdgePoint( const int theEdgeID, const gp_XYZ& theParams, gp_XYZ& thePoint ) const { - if ( !IsEdgeID( theEdgeID )) return false; - thePoint = myEdge[ theEdgeID - ID_FirstE ].Point( theParams ); return true; - } - // return coordinates of a point on edge - - bool EdgeU( const int theEdgeID, const gp_XYZ& theParams, double& theU ) const { - if ( !IsEdgeID( theEdgeID )) return false; - theU = myEdge[ theEdgeID - ID_FirstE ].GetU( theParams ); return true; - } - // return parameter on edge by in-block parameters - - bool FacePoint( const int theFaceID, const gp_XYZ& theParams, gp_XYZ& thePoint ) const { - if ( !IsFaceID ( theFaceID )) return false; - thePoint = myFace[ theFaceID - ID_FirstF ].Point( theParams ); return true; - } - // return coordinates of a point on face - - bool FaceUV( const int theFaceID, const gp_XYZ& theParams, gp_XY& theUV ) const { - if ( !IsFaceID ( theFaceID )) return false; - theUV = myFace[ theFaceID - ID_FirstF ].GetUV( theParams ); return true; - } - // return UV coordinates on a face by in-block parameters - - bool ShellPoint( const gp_XYZ& theParams, gp_XYZ& thePoint ) const; - // return coordinates of a point in shell - - static bool ShellPoint(const gp_XYZ& theParams, - const std::vector& thePointOnShape, - gp_XYZ& thePoint ); - // computes coordinates of a point in shell by points on sub-shapes - // and point parameters. - // thePointOnShape[ subShapeID ] must be a point on a subShape; - // thePointOnShape.size() == ID_Shell, thePointOnShape[0] not used - - - public: - // --------------------------------- - // Define parameters by coordinates - // --------------------------------- - - bool ComputeParameters (const gp_Pnt& thePoint, - gp_XYZ& theParams, - const int theShapeID = ID_Shell, - const gp_XYZ& theParamsHint = gp_XYZ(-1,-1,-1)); - // compute point parameters in the block. - // Note: for edges, it is better to use EdgeParameters() - - bool VertexParameters(const int theVertexID, gp_XYZ& theParams); - // return parameters of a vertex given by TShapeID - - bool EdgeParameters(const int theEdgeID, const double theU, gp_XYZ& theParams); - // return parameters of a point given by theU on edge - - - public: - // --------------- - // Block geomerty - // --------------- - - - - public: - // --------- - // Services - // --------- - - static bool IsForwardEdge (const TopoDS_Edge & theEdge, - const TopTools_IndexedMapOfOrientedShape& theShapeIDMap) { - int v1ID = theShapeIDMap.FindIndex( TopExp::FirstVertex( theEdge ).Oriented( TopAbs_FORWARD )); - int v2ID = theShapeIDMap.FindIndex( TopExp::LastVertex( theEdge ).Oriented( TopAbs_FORWARD )); - return ( v1ID < v2ID ); - } - // Return true if an in-block parameter increases along theEdge curve - - static int GetOrderedEdges (const TopoDS_Face& theFace, - TopoDS_Vertex theFirstVertex, - std::list< TopoDS_Edge >& theEdges, - std::list< int > & theNbEdgesInWires, - const bool theShapeAnalysisAlgo=false); - // Return nb wires and a list of oredered edges. - // It is used to assign indices to subshapes. - // theFirstVertex may be NULL. - // Always try to set a seam edge first - // if (theShapeAnalysisAlgo) then ShapeAnalysis::OuterWire() is used to find the outer - // wire else BRepTools::OuterWire() is used - - public: - // ----------------------------------------------------------- - // Methods of math_FunctionSetWithDerivatives used internally - // to define parameters by coordinates - // ----------------------------------------------------------- - Standard_Integer NbVariables() const; - Standard_Integer NbEquations() const; - Standard_Boolean Value(const math_Vector& X,math_Vector& F) ; - Standard_Boolean Derivatives(const math_Vector& X,math_Matrix& D) ; - Standard_Boolean Values(const math_Vector& X,math_Vector& F,math_Matrix& D) ; - Standard_Integer GetStateNumber (); - - protected: - - /*! - * \brief Call it after geometry initialisation - */ - void init(); - - // Note: to compute params of a point on a face, it is enough to set - // TFace, TEdge's and points for that face only - - // Note 2: curve adaptors need to have only Value(double), FirstParameter() and - // LastParameter() defined to be used by Block algoritms - - class SMESH_EXPORT TEdge { - int myCoordInd; - double myFirst; - double myLast; - Adaptor3d_Curve* myC3d; - // if mesh volume - gp_XYZ myNodes[2]; - public: - void Set( const int edgeID, Adaptor3d_Curve* curve, const bool isForward ); - void Set( const int edgeID, const gp_XYZ& node1, const gp_XYZ& node2 ); - Adaptor3d_Curve* GetCurve() const { return myC3d; } - double EndParam(int i) const { return i ? myLast : myFirst; } - int CoordInd() const { return myCoordInd; } - const gp_XYZ& NodeXYZ(int i) const { return i ? myNodes[1] : myNodes[0]; } - gp_XYZ Point( const gp_XYZ& theParams ) const; // Return coord by params - double GetU( const gp_XYZ& theParams ) const; // Return U by params - TEdge(): myC3d(0) {} - ~TEdge(); - }; - - class SMESH_EXPORT TFace { - // 4 edges in the order u0, u1, 0v, 1v - int myCoordInd[ 4 ]; - double myFirst [ 4 ]; - double myLast [ 4 ]; - Adaptor2d_Curve2d* myC2d [ 4 ]; - // 4 corner points in the order 00, 10, 11, 01 - gp_XY myCorner [ 4 ]; - // surface - Adaptor3d_Surface* myS; - // if mesh volume - gp_XYZ myNodes[4]; - public: - void Set( const int faceID, Adaptor3d_Surface* S, // must be in GetFaceEdgesIDs() order: - Adaptor2d_Curve2d* c2d[4], const bool isForward[4] ); - void Set( const int faceID, const TEdge& edgeU0, const TEdge& edgeU1 ); - gp_XY GetUV( const gp_XYZ& theParams ) const; - gp_XYZ Point( const gp_XYZ& theParams ) const; - int GetUInd() const { return myCoordInd[ 0 ]; } - int GetVInd() const { return myCoordInd[ 2 ]; } - void GetCoefs( int i, const gp_XYZ& theParams, double& eCoef, double& vCoef ) const; - TFace(): myS(0) { myC2d[0]=myC2d[1]=myC2d[2]=myC2d[3]=0; } - ~TFace(); - }; - - // geometry in the order as in TShapeID: - // 8 vertices - gp_XYZ myPnt[ 8 ]; - // 12 edges - TEdge myEdge[ 12 ]; - // 6 faces - TFace myFace[ 6 ]; - - // for param computation - - enum { SQUARE_DIST = 0, DRV_1, DRV_2, DRV_3 }; - double distance () const { return sqrt( myValues[ SQUARE_DIST ]); } - double funcValue(double sqDist) const { return mySquareFunc ? sqDist : sqrt(sqDist); } - bool computeParameters(const gp_Pnt& thePoint, gp_XYZ& theParams, const gp_XYZ& theParamsHint); - - int myFaceIndex; - double myFaceParam; - int myNbIterations; - double mySumDist; - double myTolerance; - bool mySquareFunc; - - gp_XYZ myPoint; // the given point - gp_XYZ myParam; // the best parameters guess - double myValues[ 4 ]; // values computed at myParam: square distance and 3 derivatives - - typedef std::pair TxyzPair; - TxyzPair my3x3x3GridNodes[ 27 ]; // to compute the first param guess - bool myGridComputed; -}; - - -#endif diff --git a/src/SMESH/SMESH_Comment.hxx b/src/SMESH/SMESH_Comment.hxx deleted file mode 100644 index b817d8ee6..000000000 --- a/src/SMESH/SMESH_Comment.hxx +++ /dev/null @@ -1,76 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// SMESH SMESH : implementaion of SMESH idl descriptions -// File : SMESH_Comment.hxx -// Created : Wed Mar 14 18:28:45 2007 -// Author : Edward AGAPOV (eap) -// Module : SMESH -// $Header: -// -#ifndef SMESH_Comment_HeaderFile -#define SMESH_Comment_HeaderFile - -# include -# include - -using namespace std; - -/*! - * \brief Class to generate string from any type - */ -class SMESH_Comment : public string -{ - ostringstream _s ; - -public : - - SMESH_Comment():string("") {} - - SMESH_Comment(const SMESH_Comment& c):string() { - _s << c.c_str() ; - this->string::operator=( _s.str() ); - } - - SMESH_Comment & operator=(const SMESH_Comment& c) { - _s << c.c_str() ; - this->string::operator=( _s.str() ); - return *this; - } - - template - SMESH_Comment( const T &anything ) { - _s << anything ; - this->string::operator=( _s.str() ); - } - - template - SMESH_Comment & operator<<( const T &anything ) { - _s << anything ; - this->string::operator=( _s.str() ); - return *this ; - } - - operator char*() const { - return (char*)c_str(); - } -}; - - -#endif diff --git a/src/SMESH/SMESH_ComputeError.hxx b/src/SMESH/SMESH_ComputeError.hxx deleted file mode 100644 index e5072c155..000000000 --- a/src/SMESH/SMESH_ComputeError.hxx +++ /dev/null @@ -1,106 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// File : SMESH_Hypothesis.hxx -// Author : Edward AGAPOV (eap) -// Module : SMESH -// -#ifndef SMESH_ComputeError_HeaderFile -#define SMESH_ComputeError_HeaderFile - -#include -#include -#include - -class SMESH_Algo; -class SMDS_MeshElement; -struct SMESH_ComputeError; - -typedef boost::shared_ptr SMESH_ComputeErrorPtr; - -// ============================================================= - -enum SMESH_ComputeErrorName -{ - // If you modify it, pls update SMESH_ComputeError::CommonName() below. - // Positive values are for algo specific errors - COMPERR_OK = -1, - COMPERR_BAD_INPUT_MESH = -2, //!< wrong mesh on lower submesh - COMPERR_STD_EXCEPTION = -3, //!< some std exception raised - COMPERR_OCC_EXCEPTION = -4, //!< OCC exception raised - COMPERR_SLM_EXCEPTION = -5, //!< SALOME exception raised - COMPERR_EXCEPTION = -6, //!< other exception raised - COMPERR_MEMORY_PB = -7, //!< std::bad_alloc exception - COMPERR_ALGO_FAILED = -8, //!< algo failed for some reason - COMPERR_BAD_SHAPE = -9, //!< bad geometry - COMPERR_WARNING = -10 //!< algo reports error but sub-mesh is computed anyway -}; - -// ============================================================= -/*! - * \brief Contains an algorithm and description of an occured error - */ -// ============================================================= - -struct SMESH_ComputeError -{ - int myName; //!< SMESH_ComputeErrorName or anything algo specific - std::string myComment; - const SMESH_Algo* myAlgo; - - std::list myBadElements; //!< to explain COMPERR_BAD_INPUT_MESH - - static SMESH_ComputeErrorPtr New( int error = COMPERR_OK, - std::string comment = "", - const SMESH_Algo* algo = 0) - { return SMESH_ComputeErrorPtr( new SMESH_ComputeError( error, comment, algo )); } - - SMESH_ComputeError(int error = COMPERR_OK, - std::string comment = "", - const SMESH_Algo* algo = 0) - :myName(error),myComment(comment),myAlgo(algo) {} - - bool IsOK() { return myName == COMPERR_OK; } - bool IsKO() { return myName != COMPERR_OK && myName != COMPERR_WARNING; } - bool IsCommon() { return myName < 0; } - inline std::string CommonName() const; - -}; - -#define _case2char(err) case err: return #err; - -std::string SMESH_ComputeError::CommonName() const -{ - switch( myName ) { - _case2char(COMPERR_OK ); - _case2char(COMPERR_BAD_INPUT_MESH); - _case2char(COMPERR_STD_EXCEPTION ); - _case2char(COMPERR_OCC_EXCEPTION ); - _case2char(COMPERR_SLM_EXCEPTION ); - _case2char(COMPERR_EXCEPTION ); - _case2char(COMPERR_MEMORY_PB ); - _case2char(COMPERR_ALGO_FAILED ); - _case2char(COMPERR_BAD_SHAPE ); - _case2char(COMPERR_WARNING ); - default:; - } - return ""; -} - -#endif diff --git a/src/SMESH/SMESH_File.cxx b/src/SMESH/SMESH_File.cxx deleted file mode 100644 index d630cdfde..000000000 --- a/src/SMESH/SMESH_File.cxx +++ /dev/null @@ -1,239 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// File : SMESH_File.cxx -// Created : Wed Mar 10 11:23:25 2010 -// Author : Edward AGAPOV (eap) -// -#include "SMESH_File.hxx" -#include "utilities.h" - -#include -#include -#include -#include -#include - -#include -#include - -#ifdef WIN32 -#include -#else -#include -#endif - -//================================================================================ -/*! - * \brief Creator opening the file for reading by default - */ -//================================================================================ - -SMESH_File::SMESH_File(const std::string& name, bool open) - :_name(name), _size(-1), _file(0), _map(0), _pos(0), _end(0) -{ - if ( open ) this->open(); -} - -//================================================================================ -/*! - * \brief Destructor closing the file - */ -//================================================================================ - -SMESH_File::~SMESH_File() -{ - close(); -} - -//================================================================================ -/*! - * \brief Open file for reading. Return true if there is something to read - */ -//================================================================================ - -bool SMESH_File::open() -{ - int length = size(); - if ( !_map && length > 0 ) - { -#ifdef WNT - _file = CreateFile(_name.data(), GENERIC_READ, FILE_SHARE_READ, - NULL, OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL); - bool ok = ( _file != INVALID_HANDLE_VALUE ); -#else - _file = ::open(_name.data(), O_RDONLY ); - bool ok = ( _file > 0 ); -#endif - if ( ok ) - { -#ifdef WNT - _mapObj = CreateFileMapping(_file, NULL, PAGE_READONLY, 0, (DWORD)length, NULL); - _map = (void*) MapViewOfFile( _mapObj, FILE_MAP_READ, 0, 0, 0 ); -#else - _map = ::mmap(0,length,PROT_READ,MAP_PRIVATE,_file,0); - if ( _map == MAP_FAILED ) _map = NULL; -#endif - if ( _map != NULL ) - { - _size = length; - _pos = (char*) _map; - _end = _pos + _size; - } - else - { -#ifdef WNT - CloseHandle(_mapObj); - CloseHandle(_file); -#else - ::close(_file); -#endif - } - } - } - return _pos; -} - -//================================================================================ -/*! - * \brief Close the file - */ -//================================================================================ - -void SMESH_File::close() -{ - if ( _map != NULL ) - { -#ifdef WNT - UnmapViewOfFile(_map); - CloseHandle(_mapObj); - CloseHandle(_file); -#else - ::munmap(_map, _size); - ::close(_file); -#endif - _map = NULL; - _pos = _end = 0; - _size = -1; - } -} - -//================================================================================ -/*! - * \brief Remove the file - */ -//================================================================================ - -bool SMESH_File::remove() -{ - close(); - try { - OSD_Path filePath(TCollection_AsciiString((char*)_name.data())); - OSD_File(filePath).Remove(); - } - catch ( Standard_ProgramError ) { - MESSAGE("Can't remove file: " << _name << " ; file does not exist or permission denied"); - return false; - } - return true; -} - -//================================================================================ -/*! - * \brief Return file size - */ -//================================================================================ - -int SMESH_File::size() const -{ - if ( _size >= 0 ) return _size; // size of open file - - int size = -1; - int file = ::open( _name.data(), O_RDONLY ); - if ( file > 0 ) - { - struct stat status; - int err = fstat( file, &status); - if ( !err ) - size = status.st_size; - ::close( file ); - } - return size; -} - -//================================================================================ -/*! - * \brief Set cursor to the given position - */ -//================================================================================ - -void SMESH_File::setPos(const char* pos) -{ - if ( pos > (const char*)_map && pos < _end ) - _pos = (char*) pos; -} - -//================================================================================ -/*! - * \brief Skip till current line end and return the skipped string - */ -//================================================================================ - -std::string SMESH_File::getLine() -{ - std::string line; - const char* p = _pos; - while ( !eof() ) - if ( *(++_pos) == '\n' ) - break; - line.append( p, _pos ); - if ( !eof() ) _pos++; - return line; -} - -//================================================================================ -/*! - * \brief Move cursor to the file beginning - */ -//================================================================================ - -void SMESH_File::rewind() -{ - _pos = (char*) _map; -} - -//================================================================================ -/*! - * \brief Fill vector by reading out integers from file. Vector size gives number - * of integers to read - */ -//================================================================================ - -bool SMESH_File::getInts(std::vector& ints) -{ - int i = 0; - while ( i < ints.size() ) - { - while ( !isdigit( *_pos ) && !eof()) ++_pos; - if ( eof() ) break; - if ( _pos[-1] == '-' ) --_pos; - ints[ i++ ] = strtol( _pos, (char**)&_pos, 10 ); - } - return ( i == ints.size() ); -} diff --git a/src/SMESH/SMESH_File.hxx b/src/SMESH/SMESH_File.hxx deleted file mode 100644 index 28400df80..000000000 --- a/src/SMESH/SMESH_File.hxx +++ /dev/null @@ -1,95 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// File : SMESH_File.hxx -// Created : Wed Mar 10 10:33:04 2010 -// Author : Edward AGAPOV (eap) -// -#ifndef __SMESH_File_HXX__ -#define __SMESH_File_HXX__ - -#include "SMESH_SMESH.hxx" - -#include -#include - -#ifdef WNT -#include -#else -#include -#endif - -/*! - * \brief High level util for effective file reading and other file operations - */ -class SMESH_EXPORT SMESH_File -{ -public: - - SMESH_File(const std::string& name, bool open=true); - - ~SMESH_File(); - - std::string getName() const { return _name; } - - bool open(); - - void close(); - - bool remove(); - - int size() const; - - // ------------------------ - // Access to file contents - // ------------------------ - - operator const char*() const { return _pos; } - - bool operator++() { return ++_pos < _end; } - - void operator +=(int posDelta) { _pos+=posDelta; } - - bool eof() const { return _pos >= _end; } - - const char* getPos() const { return _pos; } - - void setPos(const char* pos); - - std::string getLine(); - - void rewind(); - - bool getInts(std::vector& ids); - -private: - - std::string _name; //!< file name - int _size; //!< file size -#ifdef WNT - HANDLE _file, _mapObj; -#else - int _file; -#endif - void* _map; - const char* _pos; //!< current position - const char* _end; //!< position after file end -}; - -#endif diff --git a/src/SMESH/SMESH_Octree.cxx b/src/SMESH/SMESH_Octree.cxx deleted file mode 100644 index 11e345ea1..000000000 --- a/src/SMESH/SMESH_Octree.cxx +++ /dev/null @@ -1,187 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// SMESH SMESH_Octree : global Octree implementation -// File : SMESH_Octree.cxx -// Created : Tue Jan 16 16:00:00 2007 -// Author : Nicolas Geimer & Aurélien Motteux(OCC) -// Module : SMESH -// -#include "SMESH_Octree.hxx" - -//=========================================================================== -/*! - * Constructor. limit must be provided at tree root construction. - * limit will be deleted by SMESH_Octree. - */ -//=========================================================================== - -SMESH_Octree::SMESH_Octree (SMESH_Octree::Limit* limit): - myChildren(NULL), - myFather(NULL), - myIsLeaf( false ), - myLimit( limit ), - myLevel(0), - myBox(NULL) -{ -} - -//================================================================================ -/*! - * \brief Compute the Octree - */ -//================================================================================ - -void SMESH_Octree::compute() -{ - if ( myLevel==0 ) - { - myBox = buildRootBox(); - if ( myLimit->myMinBoxSize > 0. && maxSize() <= myLimit->myMinBoxSize ) - myIsLeaf = true; - else - buildChildren(); - } -} - -//====================================== -/*! - * \brief SMESH_Octree Destructor - */ -//====================================== - -SMESH_Octree::~SMESH_Octree () -{ - if(myChildren != NULL) - { - if(!isLeaf()) - { - for(int i = 0; i<8; i++) - delete myChildren[i]; - delete[] myChildren; - myChildren = 0; - } - } - if ( myBox ) - delete myBox; - myBox = 0; - if ( level() == 0 ) - delete myLimit; - myLimit = 0; -} - -//================================================================= -/*! - * \brief Build the 8 children boxes and call buildChildrenData() - */ -//================================================================= - -void SMESH_Octree::buildChildren() -{ - if ( isLeaf() ) return; - - myChildren = new SMESH_Octree*[8]; - - gp_XYZ min = myBox->CornerMin(); - gp_XYZ max = myBox->CornerMax(); - gp_XYZ HSize = (max - min)/2.; - gp_XYZ mid = min + HSize; - gp_XYZ childHsize = HSize/2.; - - // get the whole model size - double rootSize = 0; - { - SMESH_Octree* root = this; - while ( root->myLevel > 0 ) - root = root->myFather; - rootSize = root->maxSize(); - } - Standard_Real XminChild, YminChild, ZminChild; - gp_XYZ minChild; - for (int i = 0; i < 8; i++) - { - // We build the eight boxes, we need 2 points to do that: - // Min and Mid - // In binary, we can write i from 0 to 7 - // For instance : - // 5 is 101, it corresponds here in coordinates to ZYX - // If coordinate is 0 in Y-> box from Ymin to Ymid - // If coordinate is 1 in Y-> box from Ymid to Ymax - // Same scheme for X and Z - // I need the minChild to build the Bnd_B3d box. - - XminChild = (i%2==0)?min.X():mid.X(); - YminChild = ((i%4)/2==0)?min.Y():mid.Y(); - ZminChild = (i<4)?min.Z():mid.Z(); - minChild.SetCoord(XminChild, YminChild, ZminChild); - - // The child is of the same type than its father (For instance, a SMESH_OctreeNode) - // We allocate the memory we need for the child - myChildren[i] = allocateOctreeChild(); - // and we assign to him its box. - myChildren[i]->myFather = this; - myChildren[i]->myLimit = myLimit; - myChildren[i]->myLevel = myLevel + 1; - myChildren[i]->myBox = new Bnd_B3d(minChild+childHsize,childHsize); - myChildren[i]->myBox->Enlarge( rootSize * 1e-10 ); - if ( myLimit->myMinBoxSize > 0. && myChildren[i]->maxSize() <= myLimit->myMinBoxSize ) - myChildren[i]->myIsLeaf = true; - } - - // After building the 8 boxes, we put the data into the children. - buildChildrenData(); - - //After we pass to the next level of the Octree - for (int i = 0; i<8; i++) - myChildren[i]->buildChildren(); -} - -//================================================================================ -/*! - * \brief Tell if Octree is a leaf or not - * An inheriting class can influence it via myIsLeaf protected field - */ -//================================================================================ - -bool SMESH_Octree::isLeaf() const -{ - return myIsLeaf || ((myLimit->myMaxLevel > 0) ? (level() >= myLimit->myMaxLevel) : false ); -} - -//=========================================================================== -/*! - * \brief Compute the bigger dimension of my box - */ -//=========================================================================== - -double SMESH_Octree::maxSize() const -{ - if ( myBox ) - { - gp_XYZ min = myBox->CornerMin(); - gp_XYZ max = myBox->CornerMax(); - gp_XYZ Size = (max - min); - double returnVal = (Size.X()>Size.Y())?Size.X():Size.Y(); - return (returnVal>Size.Z())?returnVal:Size.Z(); - } - return 0.; -} diff --git a/src/SMESH/SMESH_Octree.hxx b/src/SMESH/SMESH_Octree.hxx deleted file mode 100644 index afcb1ecaf..000000000 --- a/src/SMESH/SMESH_Octree.hxx +++ /dev/null @@ -1,123 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// SMESH SMESH_Octree : global Octree implementation -// File : SMESH_Octree.hxx -// Created : Tue Jan 16 16:00:00 2007 -// Author : Nicolas Geimer & Aurélien Motteux (OCC) -// Module : SMESH -// -#ifndef _SMESH_OCTREE_HXX_ -#define _SMESH_OCTREE_HXX_ - -#include - -class SMESH_Octree { - -public: - - // Data limiting the tree height - struct Limit { - // MaxLevel of the Octree - int myMaxLevel; - // Minimal size of the Box - double myMinBoxSize; - - // Default: - // maxLevel-> 8^8 = 16777216 terminal trees - // minSize -> box size not checked - Limit(int maxLevel=8, double minSize=0.):myMaxLevel(maxLevel),myMinBoxSize(minSize) {} - virtual ~Limit() {} // it can be inherited - }; - - // Constructor. limit must be provided at tree root construction. - // limit will be deleted by SMESH_Octree - SMESH_Octree (Limit* limit=0); - - // Destructor - virtual ~SMESH_Octree (); - - // Compute the Octree. Must be called by constructor of inheriting class - void compute(); - - // Tell if Octree is a leaf or not. - // An inheriting class can influence it via myIsLeaf protected field - bool isLeaf() const; - - // Return its level - int level() const { return myLevel; } - - // Get box to the 3d Bounding Box of the Octree - const Bnd_B3d& getBox() const { return *myBox; } - - // Compute the bigger dimension of my box - double maxSize() const; - - // Return index of a child the given point is in - inline int getChildIndex(double x, double y, double z, const gp_XYZ& boxMiddle)const; - -protected: - // Return box of the whole tree - virtual Bnd_B3d* buildRootBox() = 0; - - // Constructor for children - virtual SMESH_Octree* allocateOctreeChild() const = 0; - - // Build the data in the 8 children - virtual void buildChildrenData() = 0; - - // members - - // Array of 8 Octree children - SMESH_Octree** myChildren; - - // Point the father, set to NULL for the level 0 - SMESH_Octree* myFather; - - // Tell us if the Octree is a leaf or not - bool myIsLeaf; - - // Tree limit - const Limit* myLimit; - -private: - // Build the 8 children boxes recursively - void buildChildren(); - - // Level of the Octree - int myLevel; - - Bnd_B3d* myBox; -}; - -//================================================================================ -/*! - * \brief Return index of a child the given point is in - */ -//================================================================================ - -inline int SMESH_Octree::getChildIndex(double x, double y, double z, const gp_XYZ& mid) const -{ - return (x > mid.X()) + ( y > mid.Y())*2 + (z > mid.Z())*4; -} - -#endif diff --git a/src/SMESH/SMESH_OctreeNode.cxx b/src/SMESH/SMESH_OctreeNode.cxx deleted file mode 100644 index f735011e6..000000000 --- a/src/SMESH/SMESH_OctreeNode.cxx +++ /dev/null @@ -1,435 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// SMESH SMESH_OctreeNode : Octree with Nodes set -// inherites global class SMESH_Octree -// File : SMESH_OctreeNode.cxx -// Created : Tue Jan 16 16:00:00 2007 -// Author : Nicolas Geimer & Aurelien Motteux (OCC) -// Module : SMESH -// -#include "SMESH_OctreeNode.hxx" - -#include "SMDS_SetIterator.hxx" -#include - -using namespace std; - -//=============================================================== -/*! - * \brief Constructor : Build all the Octree using Compute() - * \param theNodes - Set of nodes, the Octree is built from this nodes - * \param maxLevel - Maximum level for the leaves - * \param maxNbNodes - Maximum number of nodes, a leaf can contain - * \param minBoxSize - Minimal size of the Octree Box - */ -//================================================================ -SMESH_OctreeNode::SMESH_OctreeNode (const TIDSortedNodeSet & theNodes, const int maxLevel, - const int maxNbNodes , const double minBoxSize ) - :SMESH_Octree( new SMESH_Octree::Limit( maxLevel,minBoxSize)), - myMaxNbNodes(maxNbNodes), - myNodes(theNodes) -{ - compute(); -} - -//================================================================================ -/*! - * \brief Constructor used to allocate a child - */ -//================================================================================ - -SMESH_OctreeNode::SMESH_OctreeNode (int maxNbNodes): - SMESH_Octree(), myMaxNbNodes(maxNbNodes) -{ -} - -//================================================================================== -/*! - * \brief Construct an empty SMESH_OctreeNode used by SMESH_Octree::buildChildren() - */ -//================================================================================== - -SMESH_Octree* SMESH_OctreeNode::allocateOctreeChild() const -{ - return new SMESH_OctreeNode(myMaxNbNodes); -} - -//====================================== -/*! - * \brief Compute the first bounding box - * - * We take the max/min coord of the nodes - */ -//====================================== - -Bnd_B3d* SMESH_OctreeNode::buildRootBox() -{ - Bnd_B3d* box = new Bnd_B3d; - TIDSortedNodeSet::iterator it = myNodes.begin(); - for (; it != myNodes.end(); it++) { - const SMDS_MeshNode* n1 = *it; - gp_XYZ p1( n1->X(), n1->Y(), n1->Z() ); - box->Add(p1); - } - if ( myNodes.size() <= myMaxNbNodes ) - myIsLeaf = true; - - return box; -} - -//==================================================================================== -/*! - * \brief Tells us if Node is inside the current box with the precision "precision" - * \param Node - Node - * \param precision - The box is enlarged with this precision - * \retval bool - True if Node is in the box within precision - */ -//==================================================================================== - -const bool SMESH_OctreeNode::isInside (const gp_XYZ& p, const double precision) -{ - if (precision <= 0.) - return !(getBox().IsOut(p)); - Bnd_B3d BoxWithPrecision = getBox(); - BoxWithPrecision.Enlarge(precision); - return ! BoxWithPrecision.IsOut(p); -} - -//================================================ -/*! - * \brief Set the data of the children - * Shares the father's data with each of his child - */ -//================================================ -void SMESH_OctreeNode::buildChildrenData() -{ - gp_XYZ min = getBox().CornerMin(); - gp_XYZ max = getBox().CornerMax(); - gp_XYZ mid = (min + max)/2.; - - TIDSortedNodeSet::iterator it = myNodes.begin(); - while (it != myNodes.end()) - { - const SMDS_MeshNode* n1 = *it; - int ChildBoxNum = getChildIndex( n1->X(), n1->Y(), n1->Z(), mid ); - SMESH_OctreeNode* myChild = dynamic_cast (myChildren[ChildBoxNum]); - myChild->myNodes.insert(myChild->myNodes.end(),n1); - myNodes.erase( it ); - it = myNodes.begin(); - } - for (int i = 0; i < 8; i++) - { - SMESH_OctreeNode* myChild = dynamic_cast (myChildren[i]); - if ( myChild->myNodes.size() <= myMaxNbNodes ) - myChild->myIsLeaf = true; - } -} - -//=================================================================== -/*! - * \brief Return in Result a list of Nodes potentials to be near Node - * \param Node - Node - * \param precision - precision used - * \param Result - list of Nodes potentials to be near Node - */ -//==================================================================== -void SMESH_OctreeNode::NodesAround (const SMDS_MeshNode * Node, - list* Result, - const double precision) -{ - gp_XYZ p(Node->X(), Node->Y(), Node->Z()); - if (isInside(p, precision)) - { - if (isLeaf()) - { - Result->insert(Result->end(), myNodes.begin(), myNodes.end()); - } - else - { - for (int i = 0; i < 8; i++) - { - SMESH_OctreeNode* myChild = dynamic_cast (myChildren[i]); - myChild->NodesAround(Node, Result, precision); - } - } - } -} - -//================================================================================ -/*! - * \brief Return in dist2Nodes nodes mapped to their square distance from Node - * \param node - node to find nodes closest to - * \param dist2Nodes - map of found nodes and their distances - * \param precision - radius of a sphere to check nodes inside - * \retval bool - true if an exact overlapping found - */ -//================================================================================ - -bool SMESH_OctreeNode::NodesAround(const gp_XYZ &node, - map& dist2Nodes, - double precision) -{ - if ( !dist2Nodes.empty() ) - precision = min ( precision, sqrt( dist2Nodes.begin()->first )); - else if ( precision == 0. ) - precision = maxSize() / 2; - - //gp_XYZ p(node->X(), node->Y(), node->Z()); - if (isInside(node, precision)) - { - if (!isLeaf()) - { - // first check a child containing node - gp_XYZ mid = (getBox().CornerMin() + getBox().CornerMax()) / 2.; - int nodeChild = getChildIndex( node.X(), node.Y(), node.Z(), mid ); - if ( ((SMESH_OctreeNode*) myChildren[nodeChild])->NodesAround(node, dist2Nodes, precision)) - return true; - - for (int i = 0; i < 8; i++) - if ( i != nodeChild ) - if (((SMESH_OctreeNode*) myChildren[i])->NodesAround(node, dist2Nodes, precision)) - return true; - } - else if ( NbNodes() > 0 ) - { - double minDist = precision * precision; - gp_Pnt p1 ( node.X(), node.Y(), node.Z() ); - TIDSortedNodeSet::iterator nIt = myNodes.begin(); - for ( ; nIt != myNodes.end(); ++nIt ) - { - gp_Pnt p2 ( (*nIt)->X(), (*nIt)->Y(), (*nIt)->Z() ); - double dist2 = p1.SquareDistance( p2 ); - if ( dist2 < minDist ) - dist2Nodes.insert( make_pair( minDist = dist2, *nIt )); - } -// if ( dist2Nodes.size() > 1 ) // leave only closest node in dist2Nodes -// dist2Nodes.erase( ++dist2Nodes.begin(), dist2Nodes.end()); - - return ( sqrt( minDist) <= precision * 1e-12 ); - } - } - return false; -} - -//============================= -/*! - * \brief Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance - * Search for all the nodes in theSetOfNodes - * Static Method : no need to create an SMESH_OctreeNode - * \param theSetOfNodes - set of nodes we look at, modified during research - * \param theGroupsOfNodes - list of nodes closed to each other returned - * \param theTolerance - Precision used, default value is 0.00001 - * \param maxLevel - Maximum level for SMESH_OctreeNode constructed, default value is -1 (Infinite) - * \param maxNbNodes - maximum Nodes in a Leaf of the SMESH_OctreeNode constructed, default value is 5 - */ -//============================= -void SMESH_OctreeNode::FindCoincidentNodes (TIDSortedNodeSet& theSetOfNodes, - list< list< const SMDS_MeshNode*> >* theGroupsOfNodes, - const double theTolerance, - const int maxLevel, - const int maxNbNodes) -{ - SMESH_OctreeNode theOctreeNode(theSetOfNodes, maxLevel, maxNbNodes, theTolerance); - theOctreeNode.FindCoincidentNodes (&theSetOfNodes, theTolerance, theGroupsOfNodes); -} - -//============================= -/*! - * \brief Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance - * Search for all the nodes in theSetOfNodes - * \note The Octree itself is also modified by this method - * \param theSetOfNodes - set of nodes we look at, modified during research - * \param theTolerance - Precision used - * \param theGroupsOfNodes - list of nodes closed to each other returned - */ -//============================= -void SMESH_OctreeNode::FindCoincidentNodes ( TIDSortedNodeSet* theSetOfNodes, - const double theTolerance, - list< list< const SMDS_MeshNode*> >* theGroupsOfNodes) -{ - TIDSortedNodeSet::iterator it1 = theSetOfNodes->begin(); - list::iterator it2; - - while (it1 != theSetOfNodes->end()) - { - const SMDS_MeshNode * n1 = *it1; - - list ListOfCoincidentNodes;// Initialize the lists via a declaration, it's enough - - list * groupPtr = 0; - - // Searching for Nodes around n1 and put them in ListofCoincidentNodes. - // Found nodes are also erased from theSetOfNodes - FindCoincidentNodes(n1, theSetOfNodes, &ListOfCoincidentNodes, theTolerance); - - // We build a list {n1 + his neigbours} and add this list in theGroupsOfNodes - for (it2 = ListOfCoincidentNodes.begin(); it2 != ListOfCoincidentNodes.end(); it2++) - { - const SMDS_MeshNode* n2 = *it2; - if ( !groupPtr ) - { - theGroupsOfNodes->push_back( list() ); - groupPtr = & theGroupsOfNodes->back(); - groupPtr->push_back( n1 ); - } - if (groupPtr->front() > n2) - groupPtr->push_front( n2 ); - else - groupPtr->push_back( n2 ); - } - if (groupPtr != 0) - groupPtr->sort(); - - theSetOfNodes->erase(it1); - it1 = theSetOfNodes->begin(); - } -} - -//====================================================================================== -/*! - * \brief Return a list of nodes closed to Node and remove it from SetOfNodes - * \note The Octree itself is also modified by this method - * \param Node - We're searching the nodes next to him. - * \param SetOfNodes - set of nodes in which we erase the found nodes - * \param Result - list of nodes closed to Node - * \param precision - Precision used - */ -//====================================================================================== -void SMESH_OctreeNode::FindCoincidentNodes (const SMDS_MeshNode * Node, - TIDSortedNodeSet* SetOfNodes, - list* Result, - const double precision) -{ - gp_XYZ p(Node->X(), Node->Y(), Node->Z()); - bool isInsideBool = isInside(p, precision); - - if (isInsideBool) - { - // I'm only looking in the leaves, since all the nodes are stored there. - if (isLeaf()) - { - gp_Pnt p1 (Node->X(), Node->Y(), Node->Z()); - - TIDSortedNodeSet myNodesCopy = myNodes; - TIDSortedNodeSet::iterator it = myNodesCopy.begin(); - double tol2 = precision * precision; - bool squareBool; - - while (it != myNodesCopy.end()) - { - const SMDS_MeshNode* n2 = *it; - // We're only looking at nodes with a superior Id. - // JFA: Why? - //if (Node->GetID() < n2->GetID()) - if (Node->GetID() != n2->GetID()) // JFA: for bug 0020185 - { - gp_Pnt p2 (n2->X(), n2->Y(), n2->Z()); - // Distance optimized computation - squareBool = (p1.SquareDistance( p2 ) <= tol2); - - // If n2 inside the SquareDistance, we add it in Result and remove it from SetOfNodes and myNodes - if (squareBool) - { - Result->insert(Result->begin(), n2); - SetOfNodes->erase( n2 ); - myNodes.erase( n2 ); - } - } - //myNodesCopy.erase( it ); - //it = myNodesCopy.begin(); - it++; - } - if (Result->size() > 0) - myNodes.erase(Node); // JFA: for bug 0020185 - } - else - { - // If I'm not a leaf, I'm going to see my children ! - for (int i = 0; i < 8; i++) - { - SMESH_OctreeNode* myChild = dynamic_cast (myChildren[i]); - myChild->FindCoincidentNodes(Node, SetOfNodes, Result, precision); - } - } - } -} - -//================================================================================ -/*! - * \brief Update data according to node movement - */ -//================================================================================ - -void SMESH_OctreeNode::UpdateByMoveNode( const SMDS_MeshNode* node, const gp_Pnt& toPnt ) -{ - if ( isLeaf() ) - { - TIDSortedNodeSet::iterator pNode = myNodes.find( node ); - bool nodeInMe = ( pNode != myNodes.end() ); - - bool pointInMe = isInside( toPnt.Coord(), 1e-10 ); - - if ( pointInMe != nodeInMe ) - { - if ( pointInMe ) - myNodes.insert( node ); - else - myNodes.erase( node ); - } - } - else if ( myChildren ) - { - gp_XYZ mid = (getBox().CornerMin() + getBox().CornerMax()) / 2.; - int nodeChild = getChildIndex( node->X(), node->Y(), node->Z(), mid ); - int pointChild = getChildIndex( toPnt.X(), toPnt.Y(), toPnt.Z(), mid ); - if ( nodeChild != pointChild ) - { - ((SMESH_OctreeNode*) myChildren[ nodeChild ])->UpdateByMoveNode( node, toPnt ); - ((SMESH_OctreeNode*) myChildren[ pointChild ])->UpdateByMoveNode( node, toPnt ); - } - } -} - -//================================================================================ -/*! - * \brief Return iterator over children - */ -//================================================================================ -SMESH_OctreeNodeIteratorPtr SMESH_OctreeNode::GetChildrenIterator() -{ - return SMESH_OctreeNodeIteratorPtr - ( new SMDS_SetIterator< SMESH_OctreeNode*, SMESH_Octree** > - ( myChildren, (( isLeaf() || !myChildren ) ? myChildren : &myChildren[ 8 ] ))); -} - -//================================================================================ -/*! - * \brief Return nodes iterator - */ -//================================================================================ -SMDS_NodeIteratorPtr SMESH_OctreeNode::GetNodeIterator() -{ - return SMDS_NodeIteratorPtr - ( new SMDS_SetIterator< SMDS_pNode, TIDSortedNodeSet::const_iterator > - ( myNodes.begin(), myNodes.size() ? myNodes.end() : myNodes.begin())); -} diff --git a/src/SMESH/SMESH_OctreeNode.hxx b/src/SMESH/SMESH_OctreeNode.hxx deleted file mode 100644 index e8f5beac2..000000000 --- a/src/SMESH/SMESH_OctreeNode.hxx +++ /dev/null @@ -1,136 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// - -// SMESH SMESH_OctreeNode : Octree with Nodes set -// inherites global class SMESH_Octree -// File : SMESH_OctreeNode.hxx -// Created : Tue Jan 16 16:00:00 2007 -// Author : Nicolas Geimer & Aurelien Motteux (OCC) -// Module : SMESH -// -#ifndef _SMESH_OCTREENODE_HXX_ -#define _SMESH_OCTREENODE_HXX_ - -#include "SMESH_Octree.hxx" -#include -#include "SMDS_MeshNode.hxx" - -#include -#include -#include - -#include "SMDS_ElemIterator.hxx" - -//forward declaration -class SMDS_MeshNode; -class SMESH_OctreeNode; - -typedef SMDS_Iterator SMESH_OctreeNodeIterator; -typedef boost::shared_ptr SMESH_OctreeNodeIteratorPtr; -typedef std::set< const SMDS_MeshNode*, TIDCompare > TIDSortedNodeSet; - -class SMESH_OctreeNode : public SMESH_Octree { - -public: - - // Constructor - SMESH_OctreeNode (const TIDSortedNodeSet& theNodes, const int maxLevel = 8, - const int maxNbNodes = 5, const double minBoxSize = 0.); - -//============================= -/*! - * \brief Empty destructor - */ -//============================= - virtual ~SMESH_OctreeNode () {}; - - // Tells us if Node is inside the current box with the precision "precision" - virtual const bool isInside(const gp_XYZ& p, const double precision = 0.); - - // Return in Result a list of Nodes potentials to be near Node - void NodesAround(const SMDS_MeshNode * Node, - std::list* Result, - const double precision = 0.); - - // Return in dist2Nodes nodes mapped to their square distance from Node - bool NodesAround(const gp_XYZ& node, - std::map& dist2Nodes, - double precision); - - // Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance - // Search for all the nodes in nodes - void FindCoincidentNodes ( TIDSortedNodeSet* nodes, - const double theTolerance, - std::list< std::list< const SMDS_MeshNode*> >* theGroupsOfNodes); - - // Static method that return in theGroupsOfNodes a list of group of nodes close to each other within - // theTolerance search for all the nodes in nodes - static void FindCoincidentNodes ( TIDSortedNodeSet& nodes, - std::list< std::list< const SMDS_MeshNode*> >* theGroupsOfNodes, - const double theTolerance = 0.00001, - const int maxLevel = -1, - const int maxNbNodes = 5); - /*! - * \brief Update data according to node movement - */ - void UpdateByMoveNode( const SMDS_MeshNode* node, const gp_Pnt& toPnt ); - /*! - * \brief Return iterator over children - */ - SMESH_OctreeNodeIteratorPtr GetChildrenIterator(); - /*! - * \brief Return nodes iterator - */ - SMDS_NodeIteratorPtr GetNodeIterator(); - /*! - * \brief Return nb nodes in a tree - */ - int NbNodes() const { return myNodes.size(); } - -protected: - - SMESH_OctreeNode (int maxNbNodes ); - - // Compute the bounding box of the whole set of nodes myNodes - virtual Bnd_B3d* buildRootBox(); - - // Shares the father's data with each of his child - virtual void buildChildrenData(); - - // Construct an empty SMESH_OctreeNode used by SMESH_Octree::buildChildren() - virtual SMESH_Octree* allocateOctreeChild() const; - - // Return in result a list of nodes closed to Node and remove it from SetOfNodes - void FindCoincidentNodes( const SMDS_MeshNode * Node, - TIDSortedNodeSet* SetOfNodes, - std::list* Result, - const double precision); - - // The max number of nodes a leaf box can contain - int myMaxNbNodes; - - // The set of nodes inside the box of the Octree (Empty if Octree is not a leaf) - TIDSortedNodeSet myNodes; - -}; - -#endif diff --git a/src/SMESH/SMESH_TypeDefs.hxx b/src/SMESH/SMESH_TypeDefs.hxx deleted file mode 100644 index e337d56ff..000000000 --- a/src/SMESH/SMESH_TypeDefs.hxx +++ /dev/null @@ -1,162 +0,0 @@ -// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE -// -// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com -// -// File : SMESH_TypeDefs.hxx -// Created : Thu Jan 27 18:38:33 2011 -// Author : Edward AGAPOV (eap) - - -#ifndef __SMESH_TypeDefs_HXX__ -#define __SMESH_TypeDefs_HXX__ - -#include "SMESH_SMESH.hxx" - -#include - -#include - -#include -#include -#include - -typedef std::map > TElemOfElemListMap; -typedef std::map TNodeNodeMap; - -//!< Set of elements sorted by ID, to be used to assure predictability of edition -typedef std::set< const SMDS_MeshElement*, TIDCompare > TIDSortedElemSet; -typedef std::set< const SMDS_MeshNode*, TIDCompare > TIDSortedNodeSet; - -typedef pair< const SMDS_MeshNode*, const SMDS_MeshNode* > NLink; - - -//======================================================================= -/*! - * \brief A sorted pair of nodes - */ -//======================================================================= - -struct SMESH_TLink: public NLink -{ - SMESH_TLink(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2 ):NLink( n1, n2 ) - { if ( n1->GetID() < n2->GetID() ) std::swap( first, second ); } - SMESH_TLink(const NLink& link ):NLink( link ) - { if ( first->GetID() < second->GetID() ) std::swap( first, second ); } - const SMDS_MeshNode* node1() const { return first; } - const SMDS_MeshNode* node2() const { return second; } -}; - -//======================================================================= -/*! - * \brief SMESH_TLink knowing its orientation - */ -//======================================================================= - -struct SMESH_OrientedLink: public SMESH_TLink -{ - bool _reversed; - SMESH_OrientedLink(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2 ) - : SMESH_TLink( n1, n2 ), _reversed( n1 != node1() ) {} -}; - -//------------------------------------------ -/*! - * \brief SMDS_MeshNode -> gp_XYZ convertor - */ -//------------------------------------------ -struct SMESH_TNodeXYZ : public gp_XYZ -{ - const SMDS_MeshNode* _node; - SMESH_TNodeXYZ( const SMDS_MeshElement* e):gp_XYZ(0,0,0),_node(0) { - if (e) { - assert( e->GetType() == SMDSAbs_Node ); - _node = static_cast(e); - SetCoord( _node->X(), _node->Y(), _node->Z() ); - } - } - double Distance(const SMDS_MeshNode* n) const { return (SMESH_TNodeXYZ( n )-*this).Modulus(); } - double SquareDistance(const SMDS_MeshNode* n) const { return (SMESH_TNodeXYZ( n )-*this).SquareModulus(); } - bool operator==(const SMESH_TNodeXYZ& other) const { return _node == other._node; } -}; - -// -------------------------------------------------------------------------------- -// class SMESH_SequenceOfElemPtr -#include - -class SMDS_MeshElement; - -typedef const SMDS_MeshElement* SMDS_MeshElementPtr; - -DEFINE_BASECOLLECTION (SMESH_BaseCollectionElemPtr, SMDS_MeshElementPtr) -DEFINE_SEQUENCE (SMESH_SequenceOfElemPtr, SMESH_BaseCollectionElemPtr, SMDS_MeshElementPtr) - - -// -------------------------------------------------------------------------------- -// class SMESH_SequenceOfNode -typedef const SMDS_MeshNode* SMDS_MeshNodePtr; - -DEFINE_BASECOLLECTION (SMESH_BaseCollectionNodePtr, SMDS_MeshNodePtr) -DEFINE_SEQUENCE(SMESH_SequenceOfNode, - SMESH_BaseCollectionNodePtr, SMDS_MeshNodePtr) - -// -------------------------------------------------------------------------------- -// #include "SMESHDS_DataMapOfShape.hxx" - -// #include - -// #include - -/// Class SMESH_IndexedMapOfShape - -// DEFINE_BASECOLLECTION (SMESH_BaseCollectionShape, TopoDS_Shape) -// DEFINE_INDEXEDMAP (SMESH_IndexedMapOfShape, SMESH_BaseCollectionShape, TopoDS_Shape) - -/// Class SMESH_IndexedDataMapOfShapeIndexedMapOfShape - -// DEFINE_BASECOLLECTION (SMESH_BaseCollectionIndexedMapOfShape, SMESH_IndexedMapOfShape) -// DEFINE_INDEXEDDATAMAP (SMESH_IndexedDataMapOfShapeIndexedMapOfShape, -// SMESH_BaseCollectionIndexedMapOfShape, TopoDS_Shape, -// SMESH_IndexedMapOfShape) - -// -------------------------------------------------------------------------------- -// class SMESH_DataMapOfElemPtrSequenceOfElemPtr - -// SMESH_EXPORT -// inline Standard_Integer HashCode(SMDS_MeshElementPtr theElem, -// const Standard_Integer theUpper) -// { -// void* anElem = (void*) theElem; -// return HashCode(anElem,theUpper); -// } - -// SMESH_EXPORT -// inline Standard_Boolean IsEqual(SMDS_MeshElementPtr theOne, -// SMDS_MeshElementPtr theTwo) -// { -// return theOne == theTwo; -// } - -// DEFINE_BASECOLLECTION (SMESH_BaseCollectionSequenceOfElemPtr, SMESH_SequenceOfElemPtr) -// DEFINE_DATAMAP (SMESH_DataMapOfElemPtrSequenceOfElemPtr, -// SMESH_BaseCollectionSequenceOfElemPtr, -// SMDS_MeshElementPtr, SMESH_SequenceOfElemPtr) - -#endif diff --git a/src/SMESHUtils/Makefile.am b/src/SMESHUtils/Makefile.am new file mode 100644 index 000000000..8125be45c --- /dev/null +++ b/src/SMESHUtils/Makefile.am @@ -0,0 +1,57 @@ +# Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +# +# This library is free software; you can redistribute it and/or +# modify it under the terms of the GNU Lesser General Public +# License as published by the Free Software Foundation; either +# version 2.1 of the License. +# +# This library is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +# Lesser General Public License for more details. +# +# You should have received a copy of the GNU Lesser General Public +# 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +# + +# File : Makefile.in +# Module : SMESH +# +include $(top_srcdir)/adm_local/unix/make_common_starter.am + +# header files +salomeinclude_HEADERS = \ + SMESH_Block.hxx \ + SMESH_TypeDefs.hxx \ + SMESH_Octree.hxx \ + SMESH_OctreeNode.hxx \ + SMESH_Comment.hxx \ + SMESH_ComputeError.hxx \ + SMESH_File.hxx \ + SMESH_Utils.hxx + +# Libraries targets + +lib_LTLIBRARIES = libSMESHUtils.la + +dist_libSMESHUtils_la_SOURCES = \ + SMESH_Block.cxx \ + SMESH_Octree.cxx \ + SMESH_OctreeNode.cxx \ + SMESH_File.cxx + +# additionnal information to compile and link file +libSMESHUtils_la_CPPFLAGS = \ + $(KERNEL_CXXFLAGS) \ + $(CAS_CPPFLAGS) \ + $(VTK_INCLUDES) \ + $(BOOST_CPPFLAGS) \ + -I$(srcdir)/../SMDS \ + -I$(srcdir)/../SMESHDS + +libSMESHUtils_la_LDFLAGS = \ + ../SMESHDS/libSMESHDS.la \ + $(CAS_LDPATH) -lTKShHealing -lTKPrim -lTKG2d diff --git a/src/SMESHUtils/SMESH_Block.cxx b/src/SMESHUtils/SMESH_Block.cxx new file mode 100644 index 000000000..b036993a5 --- /dev/null +++ b/src/SMESHUtils/SMESH_Block.cxx @@ -0,0 +1,1737 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : SMESH_Pattern.hxx +// Created : Mon Aug 2 10:30:00 2004 +// Author : Edward AGAPOV (eap) +// +#include "SMESH_Block.hxx" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include "SMDS_MeshNode.hxx" +#include "SMDS_MeshVolume.hxx" +#include "SMDS_VolumeTool.hxx" +#include "utilities.h" + +#include + +using namespace std; + +//#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 +//purpose : +//======================================================================= + +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; +} + +//======================================================================= +//function : SMESH_Block::TEdge::Point +//purpose : +//======================================================================= + +gp_XYZ SMESH_Block::TEdge::Point( const gp_XYZ& theParams ) const +{ + 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); + } +} + +//======================================================================= +//function : SMESH_Block::TFace::GetUV +//purpose : +//======================================================================= + +gp_XY SMESH_Block::TFace::GetUV( const gp_XYZ& theParams ) const +{ + gp_XY uv(0.,0.); + 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 ] ); + u = ( 1 - u ) * myFirst[ iE ] + u * myLast[ iE ]; + uv += Ecoef * myC2d[ iE ]->Value( u ).XY(); + // corner addition + uv -= Vcoef * myCorner[ iE ]; + } + return uv; +} + +//======================================================================= +//function : SMESH_Block::TFace::Point +//purpose : +//======================================================================= + +gp_XYZ SMESH_Block::TFace::Point( const gp_XYZ& theParams ) const +{ + 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; +} + +//======================================================================= +//function : GetShapeCoef +//purpose : +//======================================================================= + +double* SMESH_Block::GetShapeCoef (const int theShapeID) +{ + static double shapeCoef[][3] = { + // V000, V100, V010, V110 + { -1,-1,-1 }, { 1,-1,-1 }, { -1, 1,-1 }, { 1, 1,-1 }, + // V001, V101, V011, V111, + { -1,-1, 1 }, { 1,-1, 1 }, { -1, 1, 1 }, { 1, 1, 1 }, + // Ex00, Ex10, Ex01, Ex11, + { 0,-1,-1 }, { 0, 1,-1 }, { 0,-1, 1 }, { 0, 1, 1 }, + // E0y0, E1y0, E0y1, E1y1, + { -1, 0,-1 }, { 1, 0,-1 }, { -1, 0, 1 }, { 1, 0, 1 }, + // E00z, E10z, E01z, E11z, + { -1,-1, 0 }, { 1,-1, 0 }, { -1, 1, 0 }, { 1, 1, 0 }, + // Fxy0, Fxy1, Fx0z, Fx1z, F0yz, F1yz, + { 0, 0,-1 }, { 0, 0, 1 }, { 0,-1, 0 }, { 0, 1, 0 }, { -1, 0, 0 }, { 1, 0, 0 }, + // ID_Shell + { 0, 0, 0 } + }; + if ( theShapeID < ID_V000 || theShapeID > ID_F1yz ) + return shapeCoef[ ID_Shell - 1 ]; + + return shapeCoef[ theShapeID - 1 ]; +} + +//======================================================================= +//function : ShellPoint +//purpose : return coordinates of a point in shell +//======================================================================= + +bool SMESH_Block::ShellPoint( const gp_XYZ& theParams, gp_XYZ& thePoint ) const +{ + thePoint.SetCoord( 0., 0., 0. ); + for ( int shapeID = ID_V000; shapeID < ID_Shell; shapeID++ ) + { + // coef + double* coefs = GetShapeCoef( shapeID ); + double k = 1; + for ( int iCoef = 0; iCoef < 3; iCoef++ ) { + if ( coefs[ iCoef ] != 0 ) { + if ( coefs[ iCoef ] < 0 ) + k *= ( 1. - theParams.Coord( iCoef + 1 )); + else + k *= theParams.Coord( iCoef + 1 ); + } + } + // 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& thePointOnShape, + gp_XYZ& thePoint ) +{ + if ( thePointOnShape.size() < ID_F1yz ) + return false; + + const double x = theParams.X(), y = theParams.Y(), z = theParams.Z(); + const double x1 = 1. - x, y1 = 1. - y, z1 = 1. - z; + const vector& 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 : +//======================================================================= + +Standard_Integer SMESH_Block::NbVariables() const +{ + return 3; +} + +//======================================================================= +//function : NbEquations +//purpose : +//======================================================================= + +Standard_Integer SMESH_Block::NbEquations() const +{ + return 1; +} + +//======================================================================= +//function : Value +//purpose : +//======================================================================= + +Standard_Boolean SMESH_Block::Value(const math_Vector& theXYZ, math_Vector& theFxyz) +{ + gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) ); + if ( params.IsEqual( myParam, DBL_MIN )) { // same param + theFxyz( 1 ) = funcValue( myValues[ SQUARE_DIST ]); + } + else { + ShellPoint( params, P ); + gp_Vec dP( P - myPoint ); + theFxyz(1) = funcValue( dP.SquareMagnitude() ); + } + return true; +} + +//======================================================================= +//function : Derivatives +//purpose : +//======================================================================= + +Standard_Boolean SMESH_Block::Derivatives(const math_Vector& XYZ,math_Matrix& Df) +{ + 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 : +//======================================================================= + +Standard_Boolean SMESH_Block::Values(const math_Vector& theXYZ, + math_Vector& theFxyz, + math_Matrix& theDf) +{ + gp_XYZ P, params( theXYZ(1), theXYZ(2), theXYZ(3) ); + if ( params.IsEqual( myParam, DBL_MIN )) { // same param + 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 + MESSAGE ( "PARAM GUESS: " << params.X() << " "<< params.Y() << " "<< params.X() ); + myNbIterations++; // how many times call ShellPoint() +#endif + ShellPoint( params, P ); + + 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 ( sqDist < myValues[ SQUARE_DIST ] ) // a better guess + { + // 3 partial derivatives + 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; + 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 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 { + double Dis = plnNorm * P - plnNorm * myPoint; + theDf( 1, iP + 1 ) = Dis/Direc; + } +#endif + } +#ifdef DEBUG_PARAM_COMPUTE + MESSAGE ( "F = " << theFxyz(1) << " DRV: " << theDf(1,1) << " " << theDf(1,2) << " " << theDf(1,3) ); + myNbIterations +=3; // how many times call ShellPoint() +#endif + + // store better values + 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 + + double loopTol = 10 * myTolerance; + int nbLoops = 0; + while ( distance() > loopTol && nbLoops <= 3 ) + { + paramSearch.Perform ( *static_cast(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(); + MESSAGE ( " ------ SOLUTION: ( "<< myParam.X() <<" "<< myParam.Y() <<" "<< myParam.Z() <<" )"< 0 ) + theParams.SetCoord( myFaceIndex, myFaceParam ); + + return true; +} + +//======================================================================= +//function : ComputeParameters +//purpose : compute point parameters in the block +//======================================================================= + +bool SMESH_Block::ComputeParameters(const gp_Pnt& thePoint, + gp_XYZ& theParams, + const int theShapeID, + const gp_XYZ& theParamsHint) +{ + if ( VertexParameters( theShapeID, theParams )) + return true; + + 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 ); + + // 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 + bool needGrid = false; + gp_XYZ par000( 0, 0, 0 ), par111( 1, 1, 1 ); + double zero = DBL_MIN * DBL_MIN; + for ( int iEdge = 0, iParam = 1; iParam <= 3 && !needGrid; iParam++ ) + { + if ( isOnFace && coef[ iParam - 1 ] != 0 ) { + 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 ); + gp_Vec v01( p0, p1 ), v0P( p0, thePoint ); + double len2 = v01.SquareMagnitude(); + double par = 0; + if ( len2 > zero ) { + par = v0P.Dot( v01 ) / len2; + if ( par < 0 || par > 1 ) { // projection falls out of line ends => needGrid + needGrid = true; + break; + } + } + sumParam += par; + } + 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 ( hasHint ) + { + start = theParamsHint; + } + else if ( myGridComputed ) + { + double minDist = DBL_MAX; + gp_XYZ* bestParam = 0; + for ( int iNode = 0; iNode < 27; iNode++ ) { + TxyzPair & prmPtn = my3x3x3GridNodes[ iNode ]; + double dist = ( thePoint.XYZ() - prmPtn.second ).SquareModulus(); + if ( dist < minDist ) { + minDist = dist; + bestParam = & prmPtn.first; + } + } + start = *bestParam; + } + + myFaceIndex = -1; + myFaceParam = 0.; + if ( isOnFace ) { + // put a point on the face + for ( int iCoord = 0; iCoord < 3; iCoord++ ) + if ( coef[ iCoord ] ) { + myFaceIndex = iCoord + 1; + myFaceParam = ( coef[ iCoord ] < 0.5 ) ? 0.0 : 1.0; + start.SetCoord( myFaceIndex, myFaceParam ); + } + } + +#ifdef DEBUG_PARAM_COMPUTE + MESSAGE ( " #### POINT " < sqDistance ) { // solution get worse + if ( ++nbGetWorst > 2 ) + return computeParameters( thePoint, theParams, solution ); + } +#ifdef DEBUG_PARAM_COMPUTE + MESSAGE ( "PARAMS: ( " << params.X() <<" "<< params.Y() <<" "<< params.Z() <<" )" ); + MESSAGE ( "DIST: " << sqrt( sqDist ) ); +#endif + + if ( sqDist < sqDistance ) { // get better + sqDistance = sqDist; + solution = params; + nbGetWorst = 0; + if ( sqDistance < sqTolerance ) // a solution found + break; + } + + // 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++; + } +#ifdef DEBUG_PARAM_COMPUTE + myNbIterations += nbLoops*4; // how many times ShellPoint called + mySumDist += sqrt( sqDistance ); + MESSAGE ( " ------ SOLUTION: ( "< 0 ) + theParams.SetCoord( myFaceIndex, myFaceParam ); + + return true; +} + +//======================================================================= +//function : VertexParameters +//purpose : return parameters of a vertex given by TShapeID +//======================================================================= + +bool SMESH_Block::VertexParameters(const int theVertexID, gp_XYZ& theParams) +{ + 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; +} + +//======================================================================= +//function : DumpShapeID +//purpose : debug an id of a block sub-shape +//======================================================================= + +#define CASEDUMP(id,strm) case id: strm << #id; break; + +ostream& SMESH_Block::DumpShapeID (const int id, ostream& stream) +{ + switch ( id ) { + CASEDUMP( ID_V000, stream ); + CASEDUMP( ID_V100, stream ); + CASEDUMP( ID_V010, stream ); + CASEDUMP( ID_V110, stream ); + CASEDUMP( ID_V001, stream ); + CASEDUMP( ID_V101, stream ); + CASEDUMP( ID_V011, stream ); + CASEDUMP( ID_V111, stream ); + CASEDUMP( ID_Ex00, stream ); + CASEDUMP( ID_Ex10, stream ); + CASEDUMP( ID_Ex01, stream ); + CASEDUMP( ID_Ex11, stream ); + CASEDUMP( ID_E0y0, stream ); + CASEDUMP( ID_E1y0, stream ); + CASEDUMP( ID_E0y1, stream ); + CASEDUMP( ID_E1y1, stream ); + CASEDUMP( ID_E00z, stream ); + CASEDUMP( ID_E10z, stream ); + CASEDUMP( ID_E01z, stream ); + CASEDUMP( ID_E11z, stream ); + CASEDUMP( ID_Fxy0, stream ); + CASEDUMP( ID_Fxy1, stream ); + CASEDUMP( ID_Fx0z, stream ); + CASEDUMP( ID_Fx1z, stream ); + CASEDUMP( ID_F0yz, stream ); + CASEDUMP( ID_F1yz, stream ); + CASEDUMP( ID_Shell, stream ); + default: stream << "ID_INVALID"; + } + return stream; +} + +//======================================================================= +//function : GetShapeIDByParams +//purpose : define an id of the block sub-shape by normlized point coord +//======================================================================= + +int SMESH_Block::GetShapeIDByParams ( const gp_XYZ& theCoord ) +{ + // id ( 0 - 26 ) computation: + + // vertex ( 0 - 7 ) : id = 1*x + 2*y + 4*z + + // edge || X ( 8 - 11 ) : id = 8 + 1*y + 2*z + // edge || Y ( 12 - 15 ): id = 1*x + 12 + 2*z + // edge || Z ( 16 - 19 ): id = 1*x + 2*y + 16 + + // face || XY ( 20 - 21 ): id = 8 + 12 + 1*z - 0 + // face || XZ ( 22 - 23 ): id = 8 + 1*y + 16 - 2 + // face || YZ ( 24 - 25 ): id = 1*x + 12 + 16 - 4 + + static int iAddBnd[] = { 1, 2, 4 }; + static int iAddNotBnd[] = { 8, 12, 16 }; + static int iFaceSubst[] = { 0, 2, 4 }; + + int id = 0; + int iOnBoundary = 0; + for ( int iCoord = 0; iCoord < 3; iCoord++ ) + { + double val = theCoord.Coord( iCoord + 1 ); + if ( val == 0.0 ) + iOnBoundary++; + else if ( val == 1.0 ) + id += iAddBnd[ iOnBoundary++ ]; + else + id += iAddNotBnd[ iCoord ]; + } + if ( iOnBoundary == 1 ) // face + id -= iFaceSubst[ (id - 20) / 4 ]; + else if ( iOnBoundary == 0 ) // shell + id = 26; + + if ( id > 26 || id < 0 ) { + MESSAGE( "GetShapeIDByParams() = " << id + <<" "<< theCoord.X() <<" "<< theCoord.Y() <<" "<< theCoord.Z() ); + } + + return id + 1; // shape ids start at 1 +} + +//================================================================================ +/*! + * \brief Return number of wires and a list of oredered edges. + * \param theFace - the face to process + * \param theFirstVertex - the vertex of the outer wire to set first in the returned + * list ( theFirstVertex may be NULL ) + * \param theEdges - all ordered edges of theFace (outer edges goes first). + * \param theNbEdgesInWires - nb of edges (== nb of vertices in closed wire) in each wire + * \param theShapeAnalysisAlgo - if true, ShapeAnalysis::OuterWire() is used to find + * the outer wire else BRepTools::OuterWire() is used. + * \retval int - nb of wires + * + * Always try to set a seam edge first. + * BRepTools::OuterWire() fails e.g. in the case of issue 0020184, + * ShapeAnalysis::OuterWire() fails in the case of issue 0020452 + */ +//================================================================================ + +int SMESH_Block::GetOrderedEdges (const TopoDS_Face& theFace, + TopoDS_Vertex theFirstVertex, + list< TopoDS_Edge >& theEdges, + list< int > & theNbEdgesInWires, + const bool theShapeAnalysisAlgo) +{ + // put wires in a list, so that an outer wire comes first + list aWireList; + TopoDS_Wire anOuterWire = + theShapeAnalysisAlgo ? ShapeAnalysis::OuterWire( theFace ) : BRepTools::OuterWire( theFace ); + for ( TopoDS_Iterator wIt (theFace); wIt.More(); wIt.Next() ) + if ( wIt.Value().ShapeType() == TopAbs_WIRE ) // it can be internal vertex! + { + if ( !anOuterWire.IsSame( wIt.Value() )) + aWireList.push_back( TopoDS::Wire( wIt.Value() )); + else + aWireList.push_front( TopoDS::Wire( wIt.Value() )); + } + + // loop on edges of wires + theNbEdgesInWires.clear(); + list::iterator wlIt = aWireList.begin(); + for ( ; wlIt != aWireList.end(); wlIt++ ) + { + int iE; + BRepTools_WireExplorer wExp( *wlIt, theFace ); + for ( iE = 0; wExp.More(); wExp.Next(), iE++ ) + { + TopoDS_Edge edge = wExp.Current(); + // commented for issue 0020557, other related ones: 0020526, PAL19080 + // edge = TopoDS::Edge( edge.Oriented( wExp.Orientation() )); + theEdges.push_back( edge ); + } + if ( iE == 0 ) // wExp returns nothing if e.g. the wire contains one internal edge + { // Issue 0020676 + for ( TopoDS_Iterator e( *wlIt ); e.More(); e.Next(), ++iE ) + theEdges.push_back( TopoDS::Edge( e.Value() )); + } + theNbEdgesInWires.push_back( iE ); + iE = 0; + if ( wlIt == aWireList.begin() && theEdges.size() > 1 ) { // the outer wire + // orient closed edges + list< TopoDS_Edge >::iterator eIt, eIt2; + for ( eIt = theEdges.begin(); eIt != theEdges.end(); eIt++ ) + { + TopoDS_Edge& edge = *eIt; + if ( TopExp::FirstVertex( edge ).IsSame( TopExp::LastVertex( edge ) )) + { + eIt2 = eIt; + bool isNext = ( eIt2 == theEdges.begin() ); + TopoDS_Edge edge2 = isNext ? *(++eIt2) : *(--eIt2); + double f1,l1,f2,l2; + Handle(Geom2d_Curve) c1 = BRep_Tool::CurveOnSurface( edge, theFace, f1,l1 ); + Handle(Geom2d_Curve) c2 = BRep_Tool::CurveOnSurface( edge2, theFace, f2,l2 ); + gp_Pnt2d pf = c1->Value( edge.Orientation() == TopAbs_FORWARD ? f1 : l1 ); + gp_Pnt2d pl = c1->Value( edge.Orientation() == TopAbs_FORWARD ? l1 : f1 ); + bool isFirst = ( edge2.Orientation() == TopAbs_FORWARD ? isNext : !isNext ); + gp_Pnt2d p2 = c2->Value( isFirst ? f2 : l2 ); + 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() ) { + 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()); + TopExp::Vertices( theEdges.front(), vv[0], vv[1], true ); + if ( iE++ > theNbEdgesInWires.back() ) { +#ifdef _DEBUG_ + gp_Pnt p = BRep_Tool::Pnt( theFirstVertex ); + MESSAGE ( " : Warning : vertex "<< theFirstVertex.TShape().operator->() + << " ( " << p.X() << " " << p.Y() << " " << p.Z() << " )" + << " not found in outer wire of face "<< theFace.TShape().operator->() + << " with vertices: " ); + 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 ); + MESSAGE_ADD ( i << " " << v.TShape().operator->() << " " + << p.X() << " " << p.Y() << " " << p.Z() << " " << std::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& 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 vFxy0, vFxy1; + + V000 = theNode000Index; + V001 = theNode001Index; + + // get faces sharing V000 and V001 + list 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::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 : 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_IndexedMapOfOrientedShape& theShapeIDMap ) +{ + MESSAGE(" ::LoadBlockShapes()"); + return ( FindBlockShapes( theShell, theVertex000, theVertex001, theShapeIDMap ) && + LoadBlockShapes( theShapeIDMap )); +} + +//======================================================================= +//function : LoadBlockShapes +//purpose : add sub-shapes of theBlock to theShapeIDMap so that they get +// IDs acoording to enum TShapeID +//======================================================================= + +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; + TopoDS_Shape E0y0, E1y0, E0y1, E1y1; + TopoDS_Shape E00z, E10z, E01z, E11z; + // 6 faces + TopoDS_Shape Fxy0, Fx0z, F0yz, Fxy1, Fx1z, F1yz; + + // nb of faces bound to a vertex in TopTools_IndexedDataMapOfShapeListOfShape + // filled by TopExp::MapShapesAndAncestors() + const int NB_FACES_BY_VERTEX = 6; + + TopTools_IndexedDataMapOfShapeListOfShape 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; + } + + V000 = theVertex000; + V001 = theVertex001; + + if ( V000.IsNull() ) { + // find vertex 000 - the one with smallest coordinates + double minVal = DBL_MAX, minX, val; + for ( int i = 1; i <= 8; i++ ) { + const TopoDS_Vertex& v = TopoDS::Vertex( vfMap.FindKey( i )); + gp_Pnt P = BRep_Tool::Pnt( v ); + val = P.X() + P.Y() + P.Z(); + if ( val < minVal || ( val == minVal && P.X() < minX )) { + V000 = v; + minVal = val; + minX = P.X(); + } + } + // find vertex 001 - the one on the most vertical edge passing through V000 + TopTools_IndexedDataMapOfShapeListOfShape 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; + TopTools_ListIteratorOfListOfShape eIt ( veMap.FindFromKey( V000 )); + for ( ; eIt.More(); eIt.Next() ) { + const TopoDS_Edge& e = TopoDS::Edge( eIt.Value() ); + TopoDS_Vertex v = TopExp::FirstVertex( e ); + if ( v.IsSame( V000 )) + v = TopExp::LastVertex( e ); + val = dir001 * gp_Vec( p000, BRep_Tool::Pnt( v )).Normalized(); + if ( val > maxVal ) { + V001 = v; + maxVal = val; + } + } + } + + // find the bottom (Fxy0), Fx0z and F0yz faces + + const TopTools_ListOfShape& f000List = vfMap.FindFromKey( V000 ); + const TopTools_ListOfShape& f001List = vfMap.FindFromKey( V001 ); + if (f000List.Extent() != NB_FACES_BY_VERTEX || + f001List.Extent() != NB_FACES_BY_VERTEX ) { + MESSAGE(" LoadBlockShapes() " << f000List.Extent() << " " << f001List.Extent()); + return false; + } + TopTools_ListIteratorOfListOfShape f001It, f000It ( f000List ); + int i, j, iFound1, iFound2; + for ( j = 0; f000It.More(); f000It.Next(), j++ ) + { + if ( NB_FACES_BY_VERTEX == 6 && j % 2 ) continue; // each face encounters twice + const TopoDS_Shape& F = f000It.Value(); + for ( i = 0, f001It.Initialize( f001List ); f001It.More(); f001It.Next(), i++ ) { + if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice + if ( F.IsSame( f001It.Value() )) + break; + } + if ( f001It.More() ) // Fx0z or F0yz found + if ( Fx0z.IsNull() ) { + Fx0z = F; + iFound1 = i; + } else { + F0yz = F; + iFound2 = i; + } + else // F is the bottom face + Fxy0 = F; + } + if ( Fxy0.IsNull() || Fx0z.IsNull() || F0yz.IsNull() ) { + MESSAGE( Fxy0.IsNull() <<" "<< Fx0z.IsNull() <<" "<< F0yz.IsNull() ); + return false; + } + + // choose the top face (Fxy1) + for ( i = 0, f001It.Initialize( f001List ); f001It.More(); f001It.Next(), i++ ) { + if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice + if ( i != iFound1 && i != iFound2 ) + break; + } + Fxy1 = f001It.Value(); + if ( Fxy1.IsNull() ) { + MESSAGE(" LoadBlockShapes() error "); + return false; + } + + // find bottom edges and veritices + list< TopoDS_Edge > eList; + list< int > nbVertexInWires; + GetOrderedEdges( TopoDS::Face( Fxy0 ), TopoDS::Vertex( V000 ), eList, nbVertexInWires ); + if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) { + MESSAGE(" LoadBlockShapes() error "); + return false; + } + list< TopoDS_Edge >::iterator elIt = eList.begin(); + for ( i = 0; elIt != eList.end(); elIt++, i++ ) + switch ( i ) { + case 0: E0y0 = *elIt; V010 = TopExp::LastVertex( *elIt, true ); break; + case 1: Ex10 = *elIt; V110 = TopExp::LastVertex( *elIt, true ); break; + case 2: E1y0 = *elIt; V100 = TopExp::LastVertex( *elIt, true ); break; + case 3: Ex00 = *elIt; break; + default:; + } + if ( i != 4 || E0y0.IsNull() || Ex10.IsNull() || E1y0.IsNull() || Ex00.IsNull() ) { + MESSAGE(" LoadBlockShapes() error, eList.size()=" << eList.size()); + return false; + } + + + // find top edges and veritices + eList.clear(); + GetOrderedEdges( TopoDS::Face( Fxy1 ), TopoDS::Vertex( V001 ), eList, nbVertexInWires ); + if ( nbVertexInWires.size() != 1 || nbVertexInWires.front() != 4 ) { + MESSAGE(" LoadBlockShapes() error "); + return false; + } + for ( i = 0, elIt = eList.begin(); elIt != eList.end(); elIt++, i++ ) + switch ( i ) { + case 0: Ex01 = *elIt; V101 = TopExp::LastVertex( *elIt, true ); break; + case 1: E1y1 = *elIt; V111 = TopExp::LastVertex( *elIt, true ); break; + case 2: Ex11 = *elIt; V011 = TopExp::LastVertex( *elIt, true ); break; + case 3: E0y1 = *elIt; break; + default:; + } + if ( i != 4 || Ex01.IsNull() || E1y1.IsNull() || Ex11.IsNull() || E0y1.IsNull() ) { + MESSAGE(" LoadBlockShapes() error, eList.size()=" << eList.size()); + return false; + } + + // swap Fx0z and F0yz if necessary + TopExp_Explorer exp( Fx0z, TopAbs_VERTEX ); + for ( ; exp.More(); exp.Next() ) // Fx0z shares V101 and V100 + if ( V101.IsSame( exp.Current() ) || V100.IsSame( exp.Current() )) + break; // V101 or V100 found + if ( !exp.More() ) { // not found + std::swap( Fx0z, F0yz); + } + + // find Fx1z and F1yz faces + const TopTools_ListOfShape& f111List = vfMap.FindFromKey( V111 ); + const TopTools_ListOfShape& f110List = vfMap.FindFromKey( V110 ); + if (f111List.Extent() != NB_FACES_BY_VERTEX || + f110List.Extent() != NB_FACES_BY_VERTEX ) { + MESSAGE(" LoadBlockShapes() " << f111List.Extent() << " " << f110List.Extent()); + return false; + } + TopTools_ListIteratorOfListOfShape f111It, f110It ( f110List); + for ( j = 0 ; f110It.More(); f110It.Next(), j++ ) { + if ( NB_FACES_BY_VERTEX == 6 && j % 2 ) continue; // each face encounters twice + const TopoDS_Shape& F = f110It.Value(); + for ( i = 0, f111It.Initialize( f111List ); f111It.More(); f111It.Next(), i++ ) { + if ( NB_FACES_BY_VERTEX == 6 && i % 2 ) continue; // each face encounters twice + if ( F.IsSame( f111It.Value() )) { // Fx1z or F1yz found + if ( Fx1z.IsNull() ) + Fx1z = F; + else + F1yz = F; + } + } + } + if ( Fx1z.IsNull() || F1yz.IsNull() ) { + MESSAGE(" LoadBlockShapes() error "); + return false; + } + + // swap Fx1z and F1yz if necessary + for ( exp.Init( Fx1z, TopAbs_VERTEX ); exp.More(); exp.Next() ) + if ( V010.IsSame( exp.Current() ) || V011.IsSame( exp.Current() )) + break; + if ( !exp.More() ) { + std::swap( Fx1z, F1yz); + } + + // find vertical edges + for ( exp.Init( Fx0z, TopAbs_EDGE ); exp.More(); exp.Next() ) { + const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() ); + const TopoDS_Shape& vFirst = TopExp::FirstVertex( edge, true ); + if ( vFirst.IsSame( V001 )) + E00z = edge; + else if ( vFirst.IsSame( V100 )) + E10z = edge; + } + if ( E00z.IsNull() || E10z.IsNull() ) { + MESSAGE(" LoadBlockShapes() error "); + return false; + } + for ( exp.Init( Fx1z, TopAbs_EDGE ); exp.More(); exp.Next() ) { + const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() ); + const TopoDS_Shape& vFirst = TopExp::FirstVertex( edge, true ); + if ( vFirst.IsSame( V111 )) + E11z = edge; + else if ( vFirst.IsSame( V010 )) + E01z = edge; + } + if ( E01z.IsNull() || E11z.IsNull() ) { + MESSAGE(" LoadBlockShapes() error "); + return false; + } + + // load shapes in theShapeIDMap + + theShapeIDMap.Clear(); + + theShapeIDMap.Add(V000.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V100.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V010.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V110.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V001.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V101.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V011.Oriented( TopAbs_FORWARD )); + theShapeIDMap.Add(V111.Oriented( TopAbs_FORWARD )); + + theShapeIDMap.Add(Ex00); + theShapeIDMap.Add(Ex10); + theShapeIDMap.Add(Ex01); + theShapeIDMap.Add(Ex11); + + theShapeIDMap.Add(E0y0); + theShapeIDMap.Add(E1y0); + theShapeIDMap.Add(E0y1); + theShapeIDMap.Add(E1y1); + + theShapeIDMap.Add(E00z); + theShapeIDMap.Add(E10z); + theShapeIDMap.Add(E01z); + theShapeIDMap.Add(E11z); + + theShapeIDMap.Add(Fxy0); + theShapeIDMap.Add(Fxy1); + theShapeIDMap.Add(Fx0z); + theShapeIDMap.Add(Fx1z); + theShapeIDMap.Add(F0yz); + theShapeIDMap.Add(F1yz); + + theShapeIDMap.Add(theShell); + + 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++ ) + { + const TopoDS_Shape& S = theShapeIDMap( shapeID ); + switch ( S.ShapeType() ) + { + case TopAbs_VERTEX: { + + 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 ); + TEdge& tEdge = myEdge[ shapeID - ID_FirstE ]; + tEdge.Set( shapeID, + new BRepAdaptor_Curve( edge ), + IsForwardEdge( edge, theShapeIDMap )); + break; + } + case TopAbs_FACE: { + + if ( !LoadFace( TopoDS::Face( S ), shapeID, theShapeIDMap )) + return false; + break; + } + default: break; + } + } // loop on shapes in theShapeIDMap + + 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 +// u0 means "|| u, v == 0" +//======================================================================= + +void SMESH_Block::GetFaceEdgesIDs (const int faceID, vector< int >& edgeVec ) +{ + edgeVec.resize( 4 ); + switch ( faceID ) { + case ID_Fxy0: + edgeVec[ 0 ] = ID_Ex00; + edgeVec[ 1 ] = ID_Ex10; + edgeVec[ 2 ] = ID_E0y0; + edgeVec[ 3 ] = ID_E1y0; + break; + case ID_Fxy1: + edgeVec[ 0 ] = ID_Ex01; + edgeVec[ 1 ] = ID_Ex11; + edgeVec[ 2 ] = ID_E0y1; + edgeVec[ 3 ] = ID_E1y1; + break; + case ID_Fx0z: + edgeVec[ 0 ] = ID_Ex00; + edgeVec[ 1 ] = ID_Ex01; + edgeVec[ 2 ] = ID_E00z; + edgeVec[ 3 ] = ID_E10z; + break; + case ID_Fx1z: + edgeVec[ 0 ] = ID_Ex10; + edgeVec[ 1 ] = ID_Ex11; + edgeVec[ 2 ] = ID_E01z; + edgeVec[ 3 ] = ID_E11z; + break; + case ID_F0yz: + edgeVec[ 0 ] = ID_E0y0; + edgeVec[ 1 ] = ID_E0y1; + edgeVec[ 2 ] = ID_E00z; + edgeVec[ 3 ] = ID_E01z; + break; + case ID_F1yz: + edgeVec[ 0 ] = ID_E1y0; + edgeVec[ 1 ] = ID_E1y1; + edgeVec[ 2 ] = ID_E10z; + edgeVec[ 3 ] = ID_E11z; + break; + default: + 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 ); + } +} diff --git a/src/SMESHUtils/SMESH_Block.hxx b/src/SMESHUtils/SMESH_Block.hxx new file mode 100644 index 000000000..a6b1e72d1 --- /dev/null +++ b/src/SMESHUtils/SMESH_Block.hxx @@ -0,0 +1,391 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : SMESH_Block.hxx +// Created : Tue Nov 30 12:42:18 2004 +// Author : Edward AGAPOV (eap) +// +#ifndef SMESH_Block_HeaderFile +#define SMESH_Block_HeaderFile + +#include "SMESH_Utils.hxx" + +//#include +//#include +//#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +class SMDS_MeshVolume; +class SMDS_MeshNode; +class Adaptor3d_Surface; +class Adaptor2d_Curve2d; +class Adaptor3d_Curve; +class gp_Pnt; + +// ========================================================= +// class calculating coordinates of 3D points by normalized +// parameters inside the block and vice versa +// ========================================================= + +class SMESHUtils_EXPORT SMESH_Block: public math_FunctionSetWithDerivatives +{ + public: + enum TShapeID { + // ---------------------------- + // Ids of the block sub-shapes + // ---------------------------- + ID_NONE = 0, + + ID_V000 = 1, ID_V100, ID_V010, ID_V110, ID_V001, ID_V101, ID_V011, ID_V111, + + ID_Ex00, ID_Ex10, ID_Ex01, ID_Ex11, + ID_E0y0, ID_E1y0, ID_E0y1, ID_E1y1, + ID_E00z, ID_E10z, ID_E01z, ID_E11z, + + ID_Fxy0, ID_Fxy1, ID_Fx0z, ID_Fx1z, ID_F0yz, ID_F1yz, + + ID_Shell + }; + enum { // to use TShapeID for indexing certain type subshapes + + ID_FirstV = ID_V000, ID_FirstE = ID_Ex00, ID_FirstF = ID_Fxy0 + + }; + + + public: + // ------------------------------------------------- + // Block topology in terms of block sub-shapes' ids + // ------------------------------------------------- + + static int NbVertices() { return 8; } + static int NbEdges() { return 12; } + static int NbFaces() { return 6; } + static int NbSubShapes() { return ID_Shell; } + // to avoid magic numbers when allocating memory for subshapes + + static inline bool IsVertexID( int theShapeID ) + { return ( theShapeID >= ID_V000 && theShapeID <= ID_V111 ); } + + static inline bool IsEdgeID( int theShapeID ) + { return ( theShapeID >= ID_Ex00 && theShapeID <= ID_E11z ); } + + static inline bool IsFaceID( int theShapeID ) + { return ( theShapeID >= ID_Fxy0 && theShapeID <= ID_F1yz ); } + + static int ShapeIndex( int theShapeID ) + { + if ( IsVertexID( theShapeID )) return theShapeID - ID_V000; + if ( IsEdgeID( theShapeID )) return theShapeID - ID_Ex00; + if ( IsFaceID( theShapeID )) return theShapeID - ID_Fxy0; + return 0; + } + // return index [0-...] for each type of sub-shapes, + // for example : + // ShapeIndex( ID_Ex00 ) == 0 + // ShapeIndex( ID_Ex10 ) == 1 + + static void GetFaceEdgesIDs (const int faceID, std::vector< int >& edgeVec ); + // return edges IDs of a face in the order u0, u1, 0v, 1v + + static void GetEdgeVertexIDs (const int edgeID, std::vector< int >& vertexVec ); + // return vertex IDs of an edge + + static int GetCoordIndOnEdge (const int theEdgeID) + { return (theEdgeID < ID_E0y0) ? 1 : (theEdgeID < ID_E00z) ? 2 : 3; } + // return an index of a coordinate which varies along the edge + + static double* GetShapeCoef (const int theShapeID); + // for theShapeID( TShapeID ), returns 3 coefficients used + // to compute an addition of an on-theShape point to coordinates + // of an in-shell point. If an in-shell point has parameters (Px,Py,Pz), + // then the addition of a point P is computed as P*kx*ky*kz and ki is + // defined by the returned coef like this: + // ki = (coef[i] == 0) ? 1 : (coef[i] < 0) ? 1 - Pi : Pi + + static int GetShapeIDByParams ( const gp_XYZ& theParams ); + // define an id of the block sub-shape by point parameters + + static std::ostream& DumpShapeID (const int theBlockShapeID, std::ostream& stream); + // DEBUG: dump an id of a block sub-shape + + + public: + // --------------- + // Initialization + // --------------- + + SMESH_Block(); + + bool LoadBlockShapes(const TopoDS_Shell& theShell, + const TopoDS_Vertex& theVertex000, + const TopoDS_Vertex& theVertex001, + TopTools_IndexedMapOfOrientedShape& theShapeIDMap ); + // Initialize block geometry with theShell, + // add sub-shapes of theBlock to theShapeIDMap so that they get + // IDs acoording to enum TShapeID + + bool LoadBlockShapes(const TopTools_IndexedMapOfOrientedShape& theShapeIDMap); + // Initialize block geometry with shapes from theShapeIDMap + + bool LoadMeshBlock(const SMDS_MeshVolume* theVolume, + const int theNode000Index, + const int theNode001Index, + std::vector& theOrderedNodes); + // prepare to work with theVolume and + // return nodes in theVolume corners in the order of TShapeID enum + + bool LoadFace(const TopoDS_Face& theFace, + const int theFaceID, + const TopTools_IndexedMapOfOrientedShape& theShapeIDMap); + // Load face geometry. + // It is enough to compute params or coordinates on the face. + // Face subshapes must be loaded into theShapeIDMap before + + static bool Insert(const TopoDS_Shape& theShape, + const int theShapeID, + TopTools_IndexedMapOfOrientedShape& theShapeIDMap); + // Insert theShape into theShapeIDMap with theShapeID, + // Not yet set shapes preceding theShapeID are filled with compounds + // Return true if theShape was successfully bound to theShapeID + + static bool FindBlockShapes(const TopoDS_Shell& theShell, + const TopoDS_Vertex& theVertex000, + const TopoDS_Vertex& theVertex001, + TopTools_IndexedMapOfOrientedShape& theShapeIDMap ); + // add sub-shapes of theBlock to theShapeIDMap so that they get + // IDs acoording to enum TShapeID + +public: + // --------------------------------- + // Define coordinates by parameters + // --------------------------------- + + bool VertexPoint( const int theVertexID, gp_XYZ& thePoint ) const { + if ( !IsVertexID( theVertexID )) return false; + thePoint = myPnt[ theVertexID - ID_FirstV ]; return true; + } + // return vertex coordinates, parameters are defined by theVertexID + + bool EdgePoint( const int theEdgeID, const gp_XYZ& theParams, gp_XYZ& thePoint ) const { + if ( !IsEdgeID( theEdgeID )) return false; + thePoint = myEdge[ theEdgeID - ID_FirstE ].Point( theParams ); return true; + } + // return coordinates of a point on edge + + bool EdgeU( const int theEdgeID, const gp_XYZ& theParams, double& theU ) const { + if ( !IsEdgeID( theEdgeID )) return false; + theU = myEdge[ theEdgeID - ID_FirstE ].GetU( theParams ); return true; + } + // return parameter on edge by in-block parameters + + bool FacePoint( const int theFaceID, const gp_XYZ& theParams, gp_XYZ& thePoint ) const { + if ( !IsFaceID ( theFaceID )) return false; + thePoint = myFace[ theFaceID - ID_FirstF ].Point( theParams ); return true; + } + // return coordinates of a point on face + + bool FaceUV( const int theFaceID, const gp_XYZ& theParams, gp_XY& theUV ) const { + if ( !IsFaceID ( theFaceID )) return false; + theUV = myFace[ theFaceID - ID_FirstF ].GetUV( theParams ); return true; + } + // return UV coordinates on a face by in-block parameters + + bool ShellPoint( const gp_XYZ& theParams, gp_XYZ& thePoint ) const; + // return coordinates of a point in shell + + static bool ShellPoint(const gp_XYZ& theParams, + const std::vector& thePointOnShape, + gp_XYZ& thePoint ); + // computes coordinates of a point in shell by points on sub-shapes + // and point parameters. + // thePointOnShape[ subShapeID ] must be a point on a subShape; + // thePointOnShape.size() == ID_Shell, thePointOnShape[0] not used + + + public: + // --------------------------------- + // Define parameters by coordinates + // --------------------------------- + + bool ComputeParameters (const gp_Pnt& thePoint, + gp_XYZ& theParams, + const int theShapeID = ID_Shell, + const gp_XYZ& theParamsHint = gp_XYZ(-1,-1,-1)); + // compute point parameters in the block. + // Note: for edges, it is better to use EdgeParameters() + + bool VertexParameters(const int theVertexID, gp_XYZ& theParams); + // return parameters of a vertex given by TShapeID + + bool EdgeParameters(const int theEdgeID, const double theU, gp_XYZ& theParams); + // return parameters of a point given by theU on edge + + + public: + // --------------- + // Block geomerty + // --------------- + + + + public: + // --------- + // Services + // --------- + + static bool IsForwardEdge (const TopoDS_Edge & theEdge, + const TopTools_IndexedMapOfOrientedShape& theShapeIDMap) { + int v1ID = theShapeIDMap.FindIndex( TopExp::FirstVertex( theEdge ).Oriented( TopAbs_FORWARD )); + int v2ID = theShapeIDMap.FindIndex( TopExp::LastVertex( theEdge ).Oriented( TopAbs_FORWARD )); + return ( v1ID < v2ID ); + } + // Return true if an in-block parameter increases along theEdge curve + + static int GetOrderedEdges (const TopoDS_Face& theFace, + TopoDS_Vertex theFirstVertex, + std::list< TopoDS_Edge >& theEdges, + std::list< int > & theNbEdgesInWires, + const bool theShapeAnalysisAlgo=false); + // Return nb wires and a list of oredered edges. + // It is used to assign indices to subshapes. + // theFirstVertex may be NULL. + // Always try to set a seam edge first + // if (theShapeAnalysisAlgo) then ShapeAnalysis::OuterWire() is used to find the outer + // wire else BRepTools::OuterWire() is used + + public: + // ----------------------------------------------------------- + // Methods of math_FunctionSetWithDerivatives used internally + // to define parameters by coordinates + // ----------------------------------------------------------- + Standard_Integer NbVariables() const; + Standard_Integer NbEquations() const; + Standard_Boolean Value(const math_Vector& X,math_Vector& F) ; + Standard_Boolean Derivatives(const math_Vector& X,math_Matrix& D) ; + Standard_Boolean Values(const math_Vector& X,math_Vector& F,math_Matrix& D) ; + Standard_Integer GetStateNumber (); + + protected: + + /*! + * \brief Call it after geometry initialisation + */ + void init(); + + // Note: to compute params of a point on a face, it is enough to set + // TFace, TEdge's and points for that face only + + // Note 2: curve adaptors need to have only Value(double), FirstParameter() and + // LastParameter() defined to be used by Block algoritms + + class SMESHUtils_EXPORT TEdge { + int myCoordInd; + double myFirst; + double myLast; + Adaptor3d_Curve* myC3d; + // if mesh volume + gp_XYZ myNodes[2]; + public: + void Set( const int edgeID, Adaptor3d_Curve* curve, const bool isForward ); + void Set( const int edgeID, const gp_XYZ& node1, const gp_XYZ& node2 ); + Adaptor3d_Curve* GetCurve() const { return myC3d; } + double EndParam(int i) const { return i ? myLast : myFirst; } + int CoordInd() const { return myCoordInd; } + const gp_XYZ& NodeXYZ(int i) const { return i ? myNodes[1] : myNodes[0]; } + gp_XYZ Point( const gp_XYZ& theParams ) const; // Return coord by params + double GetU( const gp_XYZ& theParams ) const; // Return U by params + TEdge(): myC3d(0) {} + ~TEdge(); + }; + + class SMESHUtils_EXPORT TFace { + // 4 edges in the order u0, u1, 0v, 1v + int myCoordInd[ 4 ]; + double myFirst [ 4 ]; + double myLast [ 4 ]; + Adaptor2d_Curve2d* myC2d [ 4 ]; + // 4 corner points in the order 00, 10, 11, 01 + gp_XY myCorner [ 4 ]; + // surface + Adaptor3d_Surface* myS; + // if mesh volume + gp_XYZ myNodes[4]; + public: + void Set( const int faceID, Adaptor3d_Surface* S, // must be in GetFaceEdgesIDs() order: + Adaptor2d_Curve2d* c2d[4], const bool isForward[4] ); + void Set( const int faceID, const TEdge& edgeU0, const TEdge& edgeU1 ); + gp_XY GetUV( const gp_XYZ& theParams ) const; + gp_XYZ Point( const gp_XYZ& theParams ) const; + int GetUInd() const { return myCoordInd[ 0 ]; } + int GetVInd() const { return myCoordInd[ 2 ]; } + void GetCoefs( int i, const gp_XYZ& theParams, double& eCoef, double& vCoef ) const; + TFace(): myS(0) { myC2d[0]=myC2d[1]=myC2d[2]=myC2d[3]=0; } + ~TFace(); + }; + + // geometry in the order as in TShapeID: + // 8 vertices + gp_XYZ myPnt[ 8 ]; + // 12 edges + TEdge myEdge[ 12 ]; + // 6 faces + TFace myFace[ 6 ]; + + // for param computation + + enum { SQUARE_DIST = 0, DRV_1, DRV_2, DRV_3 }; + double distance () const { return sqrt( myValues[ SQUARE_DIST ]); } + double funcValue(double sqDist) const { return mySquareFunc ? sqDist : sqrt(sqDist); } + bool computeParameters(const gp_Pnt& thePoint, gp_XYZ& theParams, const gp_XYZ& theParamsHint); + + int myFaceIndex; + double myFaceParam; + int myNbIterations; + double mySumDist; + double myTolerance; + bool mySquareFunc; + + gp_XYZ myPoint; // the given point + gp_XYZ myParam; // the best parameters guess + double myValues[ 4 ]; // values computed at myParam: square distance and 3 derivatives + + typedef std::pair TxyzPair; + TxyzPair my3x3x3GridNodes[ 27 ]; // to compute the first param guess + bool myGridComputed; +}; + + +#endif diff --git a/src/SMESHUtils/SMESH_Comment.hxx b/src/SMESHUtils/SMESH_Comment.hxx new file mode 100644 index 000000000..b817d8ee6 --- /dev/null +++ b/src/SMESHUtils/SMESH_Comment.hxx @@ -0,0 +1,76 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// SMESH SMESH : implementaion of SMESH idl descriptions +// File : SMESH_Comment.hxx +// Created : Wed Mar 14 18:28:45 2007 +// Author : Edward AGAPOV (eap) +// Module : SMESH +// $Header: +// +#ifndef SMESH_Comment_HeaderFile +#define SMESH_Comment_HeaderFile + +# include +# include + +using namespace std; + +/*! + * \brief Class to generate string from any type + */ +class SMESH_Comment : public string +{ + ostringstream _s ; + +public : + + SMESH_Comment():string("") {} + + SMESH_Comment(const SMESH_Comment& c):string() { + _s << c.c_str() ; + this->string::operator=( _s.str() ); + } + + SMESH_Comment & operator=(const SMESH_Comment& c) { + _s << c.c_str() ; + this->string::operator=( _s.str() ); + return *this; + } + + template + SMESH_Comment( const T &anything ) { + _s << anything ; + this->string::operator=( _s.str() ); + } + + template + SMESH_Comment & operator<<( const T &anything ) { + _s << anything ; + this->string::operator=( _s.str() ); + return *this ; + } + + operator char*() const { + return (char*)c_str(); + } +}; + + +#endif diff --git a/src/SMESHUtils/SMESH_ComputeError.hxx b/src/SMESHUtils/SMESH_ComputeError.hxx new file mode 100644 index 000000000..e5072c155 --- /dev/null +++ b/src/SMESHUtils/SMESH_ComputeError.hxx @@ -0,0 +1,106 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : SMESH_Hypothesis.hxx +// Author : Edward AGAPOV (eap) +// Module : SMESH +// +#ifndef SMESH_ComputeError_HeaderFile +#define SMESH_ComputeError_HeaderFile + +#include +#include +#include + +class SMESH_Algo; +class SMDS_MeshElement; +struct SMESH_ComputeError; + +typedef boost::shared_ptr SMESH_ComputeErrorPtr; + +// ============================================================= + +enum SMESH_ComputeErrorName +{ + // If you modify it, pls update SMESH_ComputeError::CommonName() below. + // Positive values are for algo specific errors + COMPERR_OK = -1, + COMPERR_BAD_INPUT_MESH = -2, //!< wrong mesh on lower submesh + COMPERR_STD_EXCEPTION = -3, //!< some std exception raised + COMPERR_OCC_EXCEPTION = -4, //!< OCC exception raised + COMPERR_SLM_EXCEPTION = -5, //!< SALOME exception raised + COMPERR_EXCEPTION = -6, //!< other exception raised + COMPERR_MEMORY_PB = -7, //!< std::bad_alloc exception + COMPERR_ALGO_FAILED = -8, //!< algo failed for some reason + COMPERR_BAD_SHAPE = -9, //!< bad geometry + COMPERR_WARNING = -10 //!< algo reports error but sub-mesh is computed anyway +}; + +// ============================================================= +/*! + * \brief Contains an algorithm and description of an occured error + */ +// ============================================================= + +struct SMESH_ComputeError +{ + int myName; //!< SMESH_ComputeErrorName or anything algo specific + std::string myComment; + const SMESH_Algo* myAlgo; + + std::list myBadElements; //!< to explain COMPERR_BAD_INPUT_MESH + + static SMESH_ComputeErrorPtr New( int error = COMPERR_OK, + std::string comment = "", + const SMESH_Algo* algo = 0) + { return SMESH_ComputeErrorPtr( new SMESH_ComputeError( error, comment, algo )); } + + SMESH_ComputeError(int error = COMPERR_OK, + std::string comment = "", + const SMESH_Algo* algo = 0) + :myName(error),myComment(comment),myAlgo(algo) {} + + bool IsOK() { return myName == COMPERR_OK; } + bool IsKO() { return myName != COMPERR_OK && myName != COMPERR_WARNING; } + bool IsCommon() { return myName < 0; } + inline std::string CommonName() const; + +}; + +#define _case2char(err) case err: return #err; + +std::string SMESH_ComputeError::CommonName() const +{ + switch( myName ) { + _case2char(COMPERR_OK ); + _case2char(COMPERR_BAD_INPUT_MESH); + _case2char(COMPERR_STD_EXCEPTION ); + _case2char(COMPERR_OCC_EXCEPTION ); + _case2char(COMPERR_SLM_EXCEPTION ); + _case2char(COMPERR_EXCEPTION ); + _case2char(COMPERR_MEMORY_PB ); + _case2char(COMPERR_ALGO_FAILED ); + _case2char(COMPERR_BAD_SHAPE ); + _case2char(COMPERR_WARNING ); + default:; + } + return ""; +} + +#endif diff --git a/src/SMESHUtils/SMESH_File.cxx b/src/SMESHUtils/SMESH_File.cxx new file mode 100644 index 000000000..d630cdfde --- /dev/null +++ b/src/SMESHUtils/SMESH_File.cxx @@ -0,0 +1,239 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : SMESH_File.cxx +// Created : Wed Mar 10 11:23:25 2010 +// Author : Edward AGAPOV (eap) +// +#include "SMESH_File.hxx" +#include "utilities.h" + +#include +#include +#include +#include +#include + +#include +#include + +#ifdef WIN32 +#include +#else +#include +#endif + +//================================================================================ +/*! + * \brief Creator opening the file for reading by default + */ +//================================================================================ + +SMESH_File::SMESH_File(const std::string& name, bool open) + :_name(name), _size(-1), _file(0), _map(0), _pos(0), _end(0) +{ + if ( open ) this->open(); +} + +//================================================================================ +/*! + * \brief Destructor closing the file + */ +//================================================================================ + +SMESH_File::~SMESH_File() +{ + close(); +} + +//================================================================================ +/*! + * \brief Open file for reading. Return true if there is something to read + */ +//================================================================================ + +bool SMESH_File::open() +{ + int length = size(); + if ( !_map && length > 0 ) + { +#ifdef WNT + _file = CreateFile(_name.data(), GENERIC_READ, FILE_SHARE_READ, + NULL, OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL); + bool ok = ( _file != INVALID_HANDLE_VALUE ); +#else + _file = ::open(_name.data(), O_RDONLY ); + bool ok = ( _file > 0 ); +#endif + if ( ok ) + { +#ifdef WNT + _mapObj = CreateFileMapping(_file, NULL, PAGE_READONLY, 0, (DWORD)length, NULL); + _map = (void*) MapViewOfFile( _mapObj, FILE_MAP_READ, 0, 0, 0 ); +#else + _map = ::mmap(0,length,PROT_READ,MAP_PRIVATE,_file,0); + if ( _map == MAP_FAILED ) _map = NULL; +#endif + if ( _map != NULL ) + { + _size = length; + _pos = (char*) _map; + _end = _pos + _size; + } + else + { +#ifdef WNT + CloseHandle(_mapObj); + CloseHandle(_file); +#else + ::close(_file); +#endif + } + } + } + return _pos; +} + +//================================================================================ +/*! + * \brief Close the file + */ +//================================================================================ + +void SMESH_File::close() +{ + if ( _map != NULL ) + { +#ifdef WNT + UnmapViewOfFile(_map); + CloseHandle(_mapObj); + CloseHandle(_file); +#else + ::munmap(_map, _size); + ::close(_file); +#endif + _map = NULL; + _pos = _end = 0; + _size = -1; + } +} + +//================================================================================ +/*! + * \brief Remove the file + */ +//================================================================================ + +bool SMESH_File::remove() +{ + close(); + try { + OSD_Path filePath(TCollection_AsciiString((char*)_name.data())); + OSD_File(filePath).Remove(); + } + catch ( Standard_ProgramError ) { + MESSAGE("Can't remove file: " << _name << " ; file does not exist or permission denied"); + return false; + } + return true; +} + +//================================================================================ +/*! + * \brief Return file size + */ +//================================================================================ + +int SMESH_File::size() const +{ + if ( _size >= 0 ) return _size; // size of open file + + int size = -1; + int file = ::open( _name.data(), O_RDONLY ); + if ( file > 0 ) + { + struct stat status; + int err = fstat( file, &status); + if ( !err ) + size = status.st_size; + ::close( file ); + } + return size; +} + +//================================================================================ +/*! + * \brief Set cursor to the given position + */ +//================================================================================ + +void SMESH_File::setPos(const char* pos) +{ + if ( pos > (const char*)_map && pos < _end ) + _pos = (char*) pos; +} + +//================================================================================ +/*! + * \brief Skip till current line end and return the skipped string + */ +//================================================================================ + +std::string SMESH_File::getLine() +{ + std::string line; + const char* p = _pos; + while ( !eof() ) + if ( *(++_pos) == '\n' ) + break; + line.append( p, _pos ); + if ( !eof() ) _pos++; + return line; +} + +//================================================================================ +/*! + * \brief Move cursor to the file beginning + */ +//================================================================================ + +void SMESH_File::rewind() +{ + _pos = (char*) _map; +} + +//================================================================================ +/*! + * \brief Fill vector by reading out integers from file. Vector size gives number + * of integers to read + */ +//================================================================================ + +bool SMESH_File::getInts(std::vector& ints) +{ + int i = 0; + while ( i < ints.size() ) + { + while ( !isdigit( *_pos ) && !eof()) ++_pos; + if ( eof() ) break; + if ( _pos[-1] == '-' ) --_pos; + ints[ i++ ] = strtol( _pos, (char**)&_pos, 10 ); + } + return ( i == ints.size() ); +} diff --git a/src/SMESHUtils/SMESH_File.hxx b/src/SMESHUtils/SMESH_File.hxx new file mode 100644 index 000000000..22794e9e6 --- /dev/null +++ b/src/SMESHUtils/SMESH_File.hxx @@ -0,0 +1,95 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : SMESH_File.hxx +// Created : Wed Mar 10 10:33:04 2010 +// Author : Edward AGAPOV (eap) +// +#ifndef __SMESH_File_HXX__ +#define __SMESH_File_HXX__ + +#include "SMESH_Utils.hxx" + +#include +#include + +#ifdef WNT +#include +#else +#include +#endif + +/*! + * \brief High level util for effective file reading and other file operations + */ +class SMESHUtils_EXPORT SMESH_File +{ +public: + + SMESH_File(const std::string& name, bool open=true); + + ~SMESH_File(); + + std::string getName() const { return _name; } + + bool open(); + + void close(); + + bool remove(); + + int size() const; + + // ------------------------ + // Access to file contents + // ------------------------ + + operator const char*() const { return _pos; } + + bool operator++() { return ++_pos < _end; } + + void operator +=(int posDelta) { _pos+=posDelta; } + + bool eof() const { return _pos >= _end; } + + const char* getPos() const { return _pos; } + + void setPos(const char* pos); + + std::string getLine(); + + void rewind(); + + bool getInts(std::vector& ids); + +private: + + std::string _name; //!< file name + int _size; //!< file size +#ifdef WNT + HANDLE _file, _mapObj; +#else + int _file; +#endif + void* _map; + const char* _pos; //!< current position + const char* _end; //!< position after file end +}; + +#endif diff --git a/src/SMESHUtils/SMESH_Octree.cxx b/src/SMESHUtils/SMESH_Octree.cxx new file mode 100644 index 000000000..11e345ea1 --- /dev/null +++ b/src/SMESHUtils/SMESH_Octree.cxx @@ -0,0 +1,187 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// SMESH SMESH_Octree : global Octree implementation +// File : SMESH_Octree.cxx +// Created : Tue Jan 16 16:00:00 2007 +// Author : Nicolas Geimer & Aurélien Motteux(OCC) +// Module : SMESH +// +#include "SMESH_Octree.hxx" + +//=========================================================================== +/*! + * Constructor. limit must be provided at tree root construction. + * limit will be deleted by SMESH_Octree. + */ +//=========================================================================== + +SMESH_Octree::SMESH_Octree (SMESH_Octree::Limit* limit): + myChildren(NULL), + myFather(NULL), + myIsLeaf( false ), + myLimit( limit ), + myLevel(0), + myBox(NULL) +{ +} + +//================================================================================ +/*! + * \brief Compute the Octree + */ +//================================================================================ + +void SMESH_Octree::compute() +{ + if ( myLevel==0 ) + { + myBox = buildRootBox(); + if ( myLimit->myMinBoxSize > 0. && maxSize() <= myLimit->myMinBoxSize ) + myIsLeaf = true; + else + buildChildren(); + } +} + +//====================================== +/*! + * \brief SMESH_Octree Destructor + */ +//====================================== + +SMESH_Octree::~SMESH_Octree () +{ + if(myChildren != NULL) + { + if(!isLeaf()) + { + for(int i = 0; i<8; i++) + delete myChildren[i]; + delete[] myChildren; + myChildren = 0; + } + } + if ( myBox ) + delete myBox; + myBox = 0; + if ( level() == 0 ) + delete myLimit; + myLimit = 0; +} + +//================================================================= +/*! + * \brief Build the 8 children boxes and call buildChildrenData() + */ +//================================================================= + +void SMESH_Octree::buildChildren() +{ + if ( isLeaf() ) return; + + myChildren = new SMESH_Octree*[8]; + + gp_XYZ min = myBox->CornerMin(); + gp_XYZ max = myBox->CornerMax(); + gp_XYZ HSize = (max - min)/2.; + gp_XYZ mid = min + HSize; + gp_XYZ childHsize = HSize/2.; + + // get the whole model size + double rootSize = 0; + { + SMESH_Octree* root = this; + while ( root->myLevel > 0 ) + root = root->myFather; + rootSize = root->maxSize(); + } + Standard_Real XminChild, YminChild, ZminChild; + gp_XYZ minChild; + for (int i = 0; i < 8; i++) + { + // We build the eight boxes, we need 2 points to do that: + // Min and Mid + // In binary, we can write i from 0 to 7 + // For instance : + // 5 is 101, it corresponds here in coordinates to ZYX + // If coordinate is 0 in Y-> box from Ymin to Ymid + // If coordinate is 1 in Y-> box from Ymid to Ymax + // Same scheme for X and Z + // I need the minChild to build the Bnd_B3d box. + + XminChild = (i%2==0)?min.X():mid.X(); + YminChild = ((i%4)/2==0)?min.Y():mid.Y(); + ZminChild = (i<4)?min.Z():mid.Z(); + minChild.SetCoord(XminChild, YminChild, ZminChild); + + // The child is of the same type than its father (For instance, a SMESH_OctreeNode) + // We allocate the memory we need for the child + myChildren[i] = allocateOctreeChild(); + // and we assign to him its box. + myChildren[i]->myFather = this; + myChildren[i]->myLimit = myLimit; + myChildren[i]->myLevel = myLevel + 1; + myChildren[i]->myBox = new Bnd_B3d(minChild+childHsize,childHsize); + myChildren[i]->myBox->Enlarge( rootSize * 1e-10 ); + if ( myLimit->myMinBoxSize > 0. && myChildren[i]->maxSize() <= myLimit->myMinBoxSize ) + myChildren[i]->myIsLeaf = true; + } + + // After building the 8 boxes, we put the data into the children. + buildChildrenData(); + + //After we pass to the next level of the Octree + for (int i = 0; i<8; i++) + myChildren[i]->buildChildren(); +} + +//================================================================================ +/*! + * \brief Tell if Octree is a leaf or not + * An inheriting class can influence it via myIsLeaf protected field + */ +//================================================================================ + +bool SMESH_Octree::isLeaf() const +{ + return myIsLeaf || ((myLimit->myMaxLevel > 0) ? (level() >= myLimit->myMaxLevel) : false ); +} + +//=========================================================================== +/*! + * \brief Compute the bigger dimension of my box + */ +//=========================================================================== + +double SMESH_Octree::maxSize() const +{ + if ( myBox ) + { + gp_XYZ min = myBox->CornerMin(); + gp_XYZ max = myBox->CornerMax(); + gp_XYZ Size = (max - min); + double returnVal = (Size.X()>Size.Y())?Size.X():Size.Y(); + return (returnVal>Size.Z())?returnVal:Size.Z(); + } + return 0.; +} diff --git a/src/SMESHUtils/SMESH_Octree.hxx b/src/SMESHUtils/SMESH_Octree.hxx new file mode 100644 index 000000000..afcb1ecaf --- /dev/null +++ b/src/SMESHUtils/SMESH_Octree.hxx @@ -0,0 +1,123 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// SMESH SMESH_Octree : global Octree implementation +// File : SMESH_Octree.hxx +// Created : Tue Jan 16 16:00:00 2007 +// Author : Nicolas Geimer & Aurélien Motteux (OCC) +// Module : SMESH +// +#ifndef _SMESH_OCTREE_HXX_ +#define _SMESH_OCTREE_HXX_ + +#include + +class SMESH_Octree { + +public: + + // Data limiting the tree height + struct Limit { + // MaxLevel of the Octree + int myMaxLevel; + // Minimal size of the Box + double myMinBoxSize; + + // Default: + // maxLevel-> 8^8 = 16777216 terminal trees + // minSize -> box size not checked + Limit(int maxLevel=8, double minSize=0.):myMaxLevel(maxLevel),myMinBoxSize(minSize) {} + virtual ~Limit() {} // it can be inherited + }; + + // Constructor. limit must be provided at tree root construction. + // limit will be deleted by SMESH_Octree + SMESH_Octree (Limit* limit=0); + + // Destructor + virtual ~SMESH_Octree (); + + // Compute the Octree. Must be called by constructor of inheriting class + void compute(); + + // Tell if Octree is a leaf or not. + // An inheriting class can influence it via myIsLeaf protected field + bool isLeaf() const; + + // Return its level + int level() const { return myLevel; } + + // Get box to the 3d Bounding Box of the Octree + const Bnd_B3d& getBox() const { return *myBox; } + + // Compute the bigger dimension of my box + double maxSize() const; + + // Return index of a child the given point is in + inline int getChildIndex(double x, double y, double z, const gp_XYZ& boxMiddle)const; + +protected: + // Return box of the whole tree + virtual Bnd_B3d* buildRootBox() = 0; + + // Constructor for children + virtual SMESH_Octree* allocateOctreeChild() const = 0; + + // Build the data in the 8 children + virtual void buildChildrenData() = 0; + + // members + + // Array of 8 Octree children + SMESH_Octree** myChildren; + + // Point the father, set to NULL for the level 0 + SMESH_Octree* myFather; + + // Tell us if the Octree is a leaf or not + bool myIsLeaf; + + // Tree limit + const Limit* myLimit; + +private: + // Build the 8 children boxes recursively + void buildChildren(); + + // Level of the Octree + int myLevel; + + Bnd_B3d* myBox; +}; + +//================================================================================ +/*! + * \brief Return index of a child the given point is in + */ +//================================================================================ + +inline int SMESH_Octree::getChildIndex(double x, double y, double z, const gp_XYZ& mid) const +{ + return (x > mid.X()) + ( y > mid.Y())*2 + (z > mid.Z())*4; +} + +#endif diff --git a/src/SMESHUtils/SMESH_OctreeNode.cxx b/src/SMESHUtils/SMESH_OctreeNode.cxx new file mode 100644 index 000000000..f735011e6 --- /dev/null +++ b/src/SMESHUtils/SMESH_OctreeNode.cxx @@ -0,0 +1,435 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// SMESH SMESH_OctreeNode : Octree with Nodes set +// inherites global class SMESH_Octree +// File : SMESH_OctreeNode.cxx +// Created : Tue Jan 16 16:00:00 2007 +// Author : Nicolas Geimer & Aurelien Motteux (OCC) +// Module : SMESH +// +#include "SMESH_OctreeNode.hxx" + +#include "SMDS_SetIterator.hxx" +#include + +using namespace std; + +//=============================================================== +/*! + * \brief Constructor : Build all the Octree using Compute() + * \param theNodes - Set of nodes, the Octree is built from this nodes + * \param maxLevel - Maximum level for the leaves + * \param maxNbNodes - Maximum number of nodes, a leaf can contain + * \param minBoxSize - Minimal size of the Octree Box + */ +//================================================================ +SMESH_OctreeNode::SMESH_OctreeNode (const TIDSortedNodeSet & theNodes, const int maxLevel, + const int maxNbNodes , const double minBoxSize ) + :SMESH_Octree( new SMESH_Octree::Limit( maxLevel,minBoxSize)), + myMaxNbNodes(maxNbNodes), + myNodes(theNodes) +{ + compute(); +} + +//================================================================================ +/*! + * \brief Constructor used to allocate a child + */ +//================================================================================ + +SMESH_OctreeNode::SMESH_OctreeNode (int maxNbNodes): + SMESH_Octree(), myMaxNbNodes(maxNbNodes) +{ +} + +//================================================================================== +/*! + * \brief Construct an empty SMESH_OctreeNode used by SMESH_Octree::buildChildren() + */ +//================================================================================== + +SMESH_Octree* SMESH_OctreeNode::allocateOctreeChild() const +{ + return new SMESH_OctreeNode(myMaxNbNodes); +} + +//====================================== +/*! + * \brief Compute the first bounding box + * + * We take the max/min coord of the nodes + */ +//====================================== + +Bnd_B3d* SMESH_OctreeNode::buildRootBox() +{ + Bnd_B3d* box = new Bnd_B3d; + TIDSortedNodeSet::iterator it = myNodes.begin(); + for (; it != myNodes.end(); it++) { + const SMDS_MeshNode* n1 = *it; + gp_XYZ p1( n1->X(), n1->Y(), n1->Z() ); + box->Add(p1); + } + if ( myNodes.size() <= myMaxNbNodes ) + myIsLeaf = true; + + return box; +} + +//==================================================================================== +/*! + * \brief Tells us if Node is inside the current box with the precision "precision" + * \param Node - Node + * \param precision - The box is enlarged with this precision + * \retval bool - True if Node is in the box within precision + */ +//==================================================================================== + +const bool SMESH_OctreeNode::isInside (const gp_XYZ& p, const double precision) +{ + if (precision <= 0.) + return !(getBox().IsOut(p)); + Bnd_B3d BoxWithPrecision = getBox(); + BoxWithPrecision.Enlarge(precision); + return ! BoxWithPrecision.IsOut(p); +} + +//================================================ +/*! + * \brief Set the data of the children + * Shares the father's data with each of his child + */ +//================================================ +void SMESH_OctreeNode::buildChildrenData() +{ + gp_XYZ min = getBox().CornerMin(); + gp_XYZ max = getBox().CornerMax(); + gp_XYZ mid = (min + max)/2.; + + TIDSortedNodeSet::iterator it = myNodes.begin(); + while (it != myNodes.end()) + { + const SMDS_MeshNode* n1 = *it; + int ChildBoxNum = getChildIndex( n1->X(), n1->Y(), n1->Z(), mid ); + SMESH_OctreeNode* myChild = dynamic_cast (myChildren[ChildBoxNum]); + myChild->myNodes.insert(myChild->myNodes.end(),n1); + myNodes.erase( it ); + it = myNodes.begin(); + } + for (int i = 0; i < 8; i++) + { + SMESH_OctreeNode* myChild = dynamic_cast (myChildren[i]); + if ( myChild->myNodes.size() <= myMaxNbNodes ) + myChild->myIsLeaf = true; + } +} + +//=================================================================== +/*! + * \brief Return in Result a list of Nodes potentials to be near Node + * \param Node - Node + * \param precision - precision used + * \param Result - list of Nodes potentials to be near Node + */ +//==================================================================== +void SMESH_OctreeNode::NodesAround (const SMDS_MeshNode * Node, + list* Result, + const double precision) +{ + gp_XYZ p(Node->X(), Node->Y(), Node->Z()); + if (isInside(p, precision)) + { + if (isLeaf()) + { + Result->insert(Result->end(), myNodes.begin(), myNodes.end()); + } + else + { + for (int i = 0; i < 8; i++) + { + SMESH_OctreeNode* myChild = dynamic_cast (myChildren[i]); + myChild->NodesAround(Node, Result, precision); + } + } + } +} + +//================================================================================ +/*! + * \brief Return in dist2Nodes nodes mapped to their square distance from Node + * \param node - node to find nodes closest to + * \param dist2Nodes - map of found nodes and their distances + * \param precision - radius of a sphere to check nodes inside + * \retval bool - true if an exact overlapping found + */ +//================================================================================ + +bool SMESH_OctreeNode::NodesAround(const gp_XYZ &node, + map& dist2Nodes, + double precision) +{ + if ( !dist2Nodes.empty() ) + precision = min ( precision, sqrt( dist2Nodes.begin()->first )); + else if ( precision == 0. ) + precision = maxSize() / 2; + + //gp_XYZ p(node->X(), node->Y(), node->Z()); + if (isInside(node, precision)) + { + if (!isLeaf()) + { + // first check a child containing node + gp_XYZ mid = (getBox().CornerMin() + getBox().CornerMax()) / 2.; + int nodeChild = getChildIndex( node.X(), node.Y(), node.Z(), mid ); + if ( ((SMESH_OctreeNode*) myChildren[nodeChild])->NodesAround(node, dist2Nodes, precision)) + return true; + + for (int i = 0; i < 8; i++) + if ( i != nodeChild ) + if (((SMESH_OctreeNode*) myChildren[i])->NodesAround(node, dist2Nodes, precision)) + return true; + } + else if ( NbNodes() > 0 ) + { + double minDist = precision * precision; + gp_Pnt p1 ( node.X(), node.Y(), node.Z() ); + TIDSortedNodeSet::iterator nIt = myNodes.begin(); + for ( ; nIt != myNodes.end(); ++nIt ) + { + gp_Pnt p2 ( (*nIt)->X(), (*nIt)->Y(), (*nIt)->Z() ); + double dist2 = p1.SquareDistance( p2 ); + if ( dist2 < minDist ) + dist2Nodes.insert( make_pair( minDist = dist2, *nIt )); + } +// if ( dist2Nodes.size() > 1 ) // leave only closest node in dist2Nodes +// dist2Nodes.erase( ++dist2Nodes.begin(), dist2Nodes.end()); + + return ( sqrt( minDist) <= precision * 1e-12 ); + } + } + return false; +} + +//============================= +/*! + * \brief Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance + * Search for all the nodes in theSetOfNodes + * Static Method : no need to create an SMESH_OctreeNode + * \param theSetOfNodes - set of nodes we look at, modified during research + * \param theGroupsOfNodes - list of nodes closed to each other returned + * \param theTolerance - Precision used, default value is 0.00001 + * \param maxLevel - Maximum level for SMESH_OctreeNode constructed, default value is -1 (Infinite) + * \param maxNbNodes - maximum Nodes in a Leaf of the SMESH_OctreeNode constructed, default value is 5 + */ +//============================= +void SMESH_OctreeNode::FindCoincidentNodes (TIDSortedNodeSet& theSetOfNodes, + list< list< const SMDS_MeshNode*> >* theGroupsOfNodes, + const double theTolerance, + const int maxLevel, + const int maxNbNodes) +{ + SMESH_OctreeNode theOctreeNode(theSetOfNodes, maxLevel, maxNbNodes, theTolerance); + theOctreeNode.FindCoincidentNodes (&theSetOfNodes, theTolerance, theGroupsOfNodes); +} + +//============================= +/*! + * \brief Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance + * Search for all the nodes in theSetOfNodes + * \note The Octree itself is also modified by this method + * \param theSetOfNodes - set of nodes we look at, modified during research + * \param theTolerance - Precision used + * \param theGroupsOfNodes - list of nodes closed to each other returned + */ +//============================= +void SMESH_OctreeNode::FindCoincidentNodes ( TIDSortedNodeSet* theSetOfNodes, + const double theTolerance, + list< list< const SMDS_MeshNode*> >* theGroupsOfNodes) +{ + TIDSortedNodeSet::iterator it1 = theSetOfNodes->begin(); + list::iterator it2; + + while (it1 != theSetOfNodes->end()) + { + const SMDS_MeshNode * n1 = *it1; + + list ListOfCoincidentNodes;// Initialize the lists via a declaration, it's enough + + list * groupPtr = 0; + + // Searching for Nodes around n1 and put them in ListofCoincidentNodes. + // Found nodes are also erased from theSetOfNodes + FindCoincidentNodes(n1, theSetOfNodes, &ListOfCoincidentNodes, theTolerance); + + // We build a list {n1 + his neigbours} and add this list in theGroupsOfNodes + for (it2 = ListOfCoincidentNodes.begin(); it2 != ListOfCoincidentNodes.end(); it2++) + { + const SMDS_MeshNode* n2 = *it2; + if ( !groupPtr ) + { + theGroupsOfNodes->push_back( list() ); + groupPtr = & theGroupsOfNodes->back(); + groupPtr->push_back( n1 ); + } + if (groupPtr->front() > n2) + groupPtr->push_front( n2 ); + else + groupPtr->push_back( n2 ); + } + if (groupPtr != 0) + groupPtr->sort(); + + theSetOfNodes->erase(it1); + it1 = theSetOfNodes->begin(); + } +} + +//====================================================================================== +/*! + * \brief Return a list of nodes closed to Node and remove it from SetOfNodes + * \note The Octree itself is also modified by this method + * \param Node - We're searching the nodes next to him. + * \param SetOfNodes - set of nodes in which we erase the found nodes + * \param Result - list of nodes closed to Node + * \param precision - Precision used + */ +//====================================================================================== +void SMESH_OctreeNode::FindCoincidentNodes (const SMDS_MeshNode * Node, + TIDSortedNodeSet* SetOfNodes, + list* Result, + const double precision) +{ + gp_XYZ p(Node->X(), Node->Y(), Node->Z()); + bool isInsideBool = isInside(p, precision); + + if (isInsideBool) + { + // I'm only looking in the leaves, since all the nodes are stored there. + if (isLeaf()) + { + gp_Pnt p1 (Node->X(), Node->Y(), Node->Z()); + + TIDSortedNodeSet myNodesCopy = myNodes; + TIDSortedNodeSet::iterator it = myNodesCopy.begin(); + double tol2 = precision * precision; + bool squareBool; + + while (it != myNodesCopy.end()) + { + const SMDS_MeshNode* n2 = *it; + // We're only looking at nodes with a superior Id. + // JFA: Why? + //if (Node->GetID() < n2->GetID()) + if (Node->GetID() != n2->GetID()) // JFA: for bug 0020185 + { + gp_Pnt p2 (n2->X(), n2->Y(), n2->Z()); + // Distance optimized computation + squareBool = (p1.SquareDistance( p2 ) <= tol2); + + // If n2 inside the SquareDistance, we add it in Result and remove it from SetOfNodes and myNodes + if (squareBool) + { + Result->insert(Result->begin(), n2); + SetOfNodes->erase( n2 ); + myNodes.erase( n2 ); + } + } + //myNodesCopy.erase( it ); + //it = myNodesCopy.begin(); + it++; + } + if (Result->size() > 0) + myNodes.erase(Node); // JFA: for bug 0020185 + } + else + { + // If I'm not a leaf, I'm going to see my children ! + for (int i = 0; i < 8; i++) + { + SMESH_OctreeNode* myChild = dynamic_cast (myChildren[i]); + myChild->FindCoincidentNodes(Node, SetOfNodes, Result, precision); + } + } + } +} + +//================================================================================ +/*! + * \brief Update data according to node movement + */ +//================================================================================ + +void SMESH_OctreeNode::UpdateByMoveNode( const SMDS_MeshNode* node, const gp_Pnt& toPnt ) +{ + if ( isLeaf() ) + { + TIDSortedNodeSet::iterator pNode = myNodes.find( node ); + bool nodeInMe = ( pNode != myNodes.end() ); + + bool pointInMe = isInside( toPnt.Coord(), 1e-10 ); + + if ( pointInMe != nodeInMe ) + { + if ( pointInMe ) + myNodes.insert( node ); + else + myNodes.erase( node ); + } + } + else if ( myChildren ) + { + gp_XYZ mid = (getBox().CornerMin() + getBox().CornerMax()) / 2.; + int nodeChild = getChildIndex( node->X(), node->Y(), node->Z(), mid ); + int pointChild = getChildIndex( toPnt.X(), toPnt.Y(), toPnt.Z(), mid ); + if ( nodeChild != pointChild ) + { + ((SMESH_OctreeNode*) myChildren[ nodeChild ])->UpdateByMoveNode( node, toPnt ); + ((SMESH_OctreeNode*) myChildren[ pointChild ])->UpdateByMoveNode( node, toPnt ); + } + } +} + +//================================================================================ +/*! + * \brief Return iterator over children + */ +//================================================================================ +SMESH_OctreeNodeIteratorPtr SMESH_OctreeNode::GetChildrenIterator() +{ + return SMESH_OctreeNodeIteratorPtr + ( new SMDS_SetIterator< SMESH_OctreeNode*, SMESH_Octree** > + ( myChildren, (( isLeaf() || !myChildren ) ? myChildren : &myChildren[ 8 ] ))); +} + +//================================================================================ +/*! + * \brief Return nodes iterator + */ +//================================================================================ +SMDS_NodeIteratorPtr SMESH_OctreeNode::GetNodeIterator() +{ + return SMDS_NodeIteratorPtr + ( new SMDS_SetIterator< SMDS_pNode, TIDSortedNodeSet::const_iterator > + ( myNodes.begin(), myNodes.size() ? myNodes.end() : myNodes.begin())); +} diff --git a/src/SMESHUtils/SMESH_OctreeNode.hxx b/src/SMESHUtils/SMESH_OctreeNode.hxx new file mode 100644 index 000000000..e8f5beac2 --- /dev/null +++ b/src/SMESHUtils/SMESH_OctreeNode.hxx @@ -0,0 +1,136 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// SMESH SMESH_OctreeNode : Octree with Nodes set +// inherites global class SMESH_Octree +// File : SMESH_OctreeNode.hxx +// Created : Tue Jan 16 16:00:00 2007 +// Author : Nicolas Geimer & Aurelien Motteux (OCC) +// Module : SMESH +// +#ifndef _SMESH_OCTREENODE_HXX_ +#define _SMESH_OCTREENODE_HXX_ + +#include "SMESH_Octree.hxx" +#include +#include "SMDS_MeshNode.hxx" + +#include +#include +#include + +#include "SMDS_ElemIterator.hxx" + +//forward declaration +class SMDS_MeshNode; +class SMESH_OctreeNode; + +typedef SMDS_Iterator SMESH_OctreeNodeIterator; +typedef boost::shared_ptr SMESH_OctreeNodeIteratorPtr; +typedef std::set< const SMDS_MeshNode*, TIDCompare > TIDSortedNodeSet; + +class SMESH_OctreeNode : public SMESH_Octree { + +public: + + // Constructor + SMESH_OctreeNode (const TIDSortedNodeSet& theNodes, const int maxLevel = 8, + const int maxNbNodes = 5, const double minBoxSize = 0.); + +//============================= +/*! + * \brief Empty destructor + */ +//============================= + virtual ~SMESH_OctreeNode () {}; + + // Tells us if Node is inside the current box with the precision "precision" + virtual const bool isInside(const gp_XYZ& p, const double precision = 0.); + + // Return in Result a list of Nodes potentials to be near Node + void NodesAround(const SMDS_MeshNode * Node, + std::list* Result, + const double precision = 0.); + + // Return in dist2Nodes nodes mapped to their square distance from Node + bool NodesAround(const gp_XYZ& node, + std::map& dist2Nodes, + double precision); + + // Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance + // Search for all the nodes in nodes + void FindCoincidentNodes ( TIDSortedNodeSet* nodes, + const double theTolerance, + std::list< std::list< const SMDS_MeshNode*> >* theGroupsOfNodes); + + // Static method that return in theGroupsOfNodes a list of group of nodes close to each other within + // theTolerance search for all the nodes in nodes + static void FindCoincidentNodes ( TIDSortedNodeSet& nodes, + std::list< std::list< const SMDS_MeshNode*> >* theGroupsOfNodes, + const double theTolerance = 0.00001, + const int maxLevel = -1, + const int maxNbNodes = 5); + /*! + * \brief Update data according to node movement + */ + void UpdateByMoveNode( const SMDS_MeshNode* node, const gp_Pnt& toPnt ); + /*! + * \brief Return iterator over children + */ + SMESH_OctreeNodeIteratorPtr GetChildrenIterator(); + /*! + * \brief Return nodes iterator + */ + SMDS_NodeIteratorPtr GetNodeIterator(); + /*! + * \brief Return nb nodes in a tree + */ + int NbNodes() const { return myNodes.size(); } + +protected: + + SMESH_OctreeNode (int maxNbNodes ); + + // Compute the bounding box of the whole set of nodes myNodes + virtual Bnd_B3d* buildRootBox(); + + // Shares the father's data with each of his child + virtual void buildChildrenData(); + + // Construct an empty SMESH_OctreeNode used by SMESH_Octree::buildChildren() + virtual SMESH_Octree* allocateOctreeChild() const; + + // Return in result a list of nodes closed to Node and remove it from SetOfNodes + void FindCoincidentNodes( const SMDS_MeshNode * Node, + TIDSortedNodeSet* SetOfNodes, + std::list* Result, + const double precision); + + // The max number of nodes a leaf box can contain + int myMaxNbNodes; + + // The set of nodes inside the box of the Octree (Empty if Octree is not a leaf) + TIDSortedNodeSet myNodes; + +}; + +#endif diff --git a/src/SMESHUtils/SMESH_TypeDefs.hxx b/src/SMESHUtils/SMESH_TypeDefs.hxx new file mode 100644 index 000000000..404565812 --- /dev/null +++ b/src/SMESHUtils/SMESH_TypeDefs.hxx @@ -0,0 +1,175 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// File : SMESH_TypeDefs.hxx +// Created : Thu Jan 27 18:38:33 2011 +// Author : Edward AGAPOV (eap) + + +#ifndef __SMESH_TypeDefs_HXX__ +#define __SMESH_TypeDefs_HXX__ + +#include "SMESH_Utils.hxx" + +#include + +#include + +#include +#include +#include + +typedef std::map > TElemOfElemListMap; +typedef std::map TNodeNodeMap; + +//!< Set of elements sorted by ID, to be used to assure predictability of edition +typedef std::set< const SMDS_MeshElement*, TIDCompare > TIDSortedElemSet; +typedef std::set< const SMDS_MeshNode*, TIDCompare > TIDSortedNodeSet; + +typedef pair< const SMDS_MeshNode*, const SMDS_MeshNode* > NLink; + + +namespace SMESHUtils +{ + /*! + * \brief Enforce freeing memory allocated by std::vector + */ + template + void FreeVector(TVECTOR& vec) + { + TVECTOR v2; + vec.swap( v2 ); + } +} + +//======================================================================= +/*! + * \brief A sorted pair of nodes + */ +//======================================================================= + +struct SMESH_TLink: public NLink +{ + SMESH_TLink(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2 ):NLink( n1, n2 ) + { if ( n1->GetID() < n2->GetID() ) std::swap( first, second ); } + SMESH_TLink(const NLink& link ):NLink( link ) + { if ( first->GetID() < second->GetID() ) std::swap( first, second ); } + const SMDS_MeshNode* node1() const { return first; } + const SMDS_MeshNode* node2() const { return second; } +}; + +//======================================================================= +/*! + * \brief SMESH_TLink knowing its orientation + */ +//======================================================================= + +struct SMESH_OrientedLink: public SMESH_TLink +{ + bool _reversed; + SMESH_OrientedLink(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2 ) + : SMESH_TLink( n1, n2 ), _reversed( n1 != node1() ) {} +}; + +//------------------------------------------ +/*! + * \brief SMDS_MeshNode -> gp_XYZ convertor + */ +//------------------------------------------ +struct SMESH_TNodeXYZ : public gp_XYZ +{ + const SMDS_MeshNode* _node; + SMESH_TNodeXYZ( const SMDS_MeshElement* e):gp_XYZ(0,0,0),_node(0) { + if (e) { + assert( e->GetType() == SMDSAbs_Node ); + _node = static_cast(e); + SetCoord( _node->X(), _node->Y(), _node->Z() ); + } + } + double Distance(const SMDS_MeshNode* n) const { return (SMESH_TNodeXYZ( n )-*this).Modulus(); } + double SquareDistance(const SMDS_MeshNode* n) const { return (SMESH_TNodeXYZ( n )-*this).SquareModulus(); } + bool operator==(const SMESH_TNodeXYZ& other) const { return _node == other._node; } +}; + +// -------------------------------------------------------------------------------- +// class SMESH_SequenceOfElemPtr +#include + +class SMDS_MeshElement; + +typedef const SMDS_MeshElement* SMDS_MeshElementPtr; + +DEFINE_BASECOLLECTION (SMESH_BaseCollectionElemPtr, SMDS_MeshElementPtr) +DEFINE_SEQUENCE (SMESH_SequenceOfElemPtr, SMESH_BaseCollectionElemPtr, SMDS_MeshElementPtr) + + +// -------------------------------------------------------------------------------- +// class SMESH_SequenceOfNode +typedef const SMDS_MeshNode* SMDS_MeshNodePtr; + +DEFINE_BASECOLLECTION (SMESH_BaseCollectionNodePtr, SMDS_MeshNodePtr) +DEFINE_SEQUENCE(SMESH_SequenceOfNode, + SMESH_BaseCollectionNodePtr, SMDS_MeshNodePtr) + +// -------------------------------------------------------------------------------- +// #include "SMESHDS_DataMapOfShape.hxx" + +// #include + +// #include + +/// Class SMESH_IndexedMapOfShape + +// DEFINE_BASECOLLECTION (SMESH_BaseCollectionShape, TopoDS_Shape) +// DEFINE_INDEXEDMAP (SMESH_IndexedMapOfShape, SMESH_BaseCollectionShape, TopoDS_Shape) + +/// Class SMESH_IndexedDataMapOfShapeIndexedMapOfShape + +// DEFINE_BASECOLLECTION (SMESH_BaseCollectionIndexedMapOfShape, SMESH_IndexedMapOfShape) +// DEFINE_INDEXEDDATAMAP (SMESH_IndexedDataMapOfShapeIndexedMapOfShape, +// SMESH_BaseCollectionIndexedMapOfShape, TopoDS_Shape, +// SMESH_IndexedMapOfShape) + +// -------------------------------------------------------------------------------- +// class SMESH_DataMapOfElemPtrSequenceOfElemPtr + +// SMESHUtils_EXPORT +// inline Standard_Integer HashCode(SMDS_MeshElementPtr theElem, +// const Standard_Integer theUpper) +// { +// void* anElem = (void*) theElem; +// return HashCode(anElem,theUpper); +// } + +// SMESHUtils_EXPORT +// inline Standard_Boolean IsEqual(SMDS_MeshElementPtr theOne, +// SMDS_MeshElementPtr theTwo) +// { +// return theOne == theTwo; +// } + +// DEFINE_BASECOLLECTION (SMESH_BaseCollectionSequenceOfElemPtr, SMESH_SequenceOfElemPtr) +// DEFINE_DATAMAP (SMESH_DataMapOfElemPtrSequenceOfElemPtr, +// SMESH_BaseCollectionSequenceOfElemPtr, +// SMDS_MeshElementPtr, SMESH_SequenceOfElemPtr) + +#endif diff --git a/src/SMESHUtils/SMESH_Utils.hxx b/src/SMESHUtils/SMESH_Utils.hxx new file mode 100755 index 000000000..de8f2a0f3 --- /dev/null +++ b/src/SMESHUtils/SMESH_Utils.hxx @@ -0,0 +1,40 @@ +// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// 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.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : SMESH_Utils.hxx +// Author : Alexander A. BORODIN +// Module : SMESH +// +#ifndef _SMESH_Utils_hxx_ +#define _SMESH_Utils_hxx_ + +#ifdef WNT + #if defined SMESHUtils_EXPORTS + #define SMESHUtils_EXPORT __declspec( dllexport ) + #else + #define SMESHUtils_EXPORT __declspec( dllimport ) + #endif +#else + #define SMESHUtils_EXPORT +#endif + +#endif