From: jfa Date: Thu, 23 Dec 2021 10:02:29 +0000 (+0300) Subject: Call FrontTrack via its python interface X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=890b75c52ef73d2c9e4ecbd625e8fc6b59b8fbe1;p=modules%2Fsmesh.git Call FrontTrack via its python interface --- diff --git a/src/SMESH/SMESH_Homard.cxx b/src/SMESH/SMESH_Homard.cxx index 44dd23937..e02ea2ffe 100644 --- a/src/SMESH/SMESH_Homard.cxx +++ b/src/SMESH/SMESH_Homard.cxx @@ -35,58 +35,6 @@ #include #endif -//// - -#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 -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -namespace boofs = boost::filesystem; - -//// - namespace SMESHHOMARDImpl { @@ -2977,1737 +2925,4 @@ int HOMARD_Iteration::GetInfoCompute() const return _MessInfo; } - -/*! - * \brief Relocate nodes to lie on geometry - * \param [in] theInputMedFile - a MED file holding a mesh including nodes that will be - * moved onto the geometry - * \param [in] theOutputMedFile - a MED file to create, that will hold a modified mesh - * \param [in] theInputNodeFiles - an array of names of files describing groups of nodes that - * will be moved onto the geometry - * \param [in] theXaoFileName - a path to a file in XAO format holding the geometry and - * the geometrical groups. - * \param [in] theIsParallel - if \c true, all processors are used to treat boundary shapes - * in parallel. - */ -void FrontTrack::track( const std::string& theInputMedFile, - const std::string& theOutputMedFile, - const std::vector< std::string > & theInputNodeFiles, - const std::string& theXaoFileName, - bool theIsParallel ) -{ - // check arguments -#ifdef _DEBUG_ - std::cout << "FrontTrack::track" << std::endl; -#endif - - if ( theInputNodeFiles.empty() ) - return; - -#ifdef _DEBUG_ - std::cout << "Input MED file: " << theInputMedFile << std::endl; -#endif - if ( !FT_Utils::fileExists( theInputMedFile )) - throw std::invalid_argument( "Input MED file does not exist: " + theInputMedFile ); - -#ifdef _DEBUG_ - std::cout << "Output MED file: " << theOutputMedFile << std::endl; -#endif - if ( !FT_Utils::canWrite( theOutputMedFile )) - throw std::invalid_argument( "Can't create the output MED file: " + theOutputMedFile ); - - std::vector< std::string > theNodeFiles ; - for ( size_t i = 0; i < theInputNodeFiles.size(); ++i ) - { -#ifdef _DEBUG_ - std::cout << "Initial input node file #"<1 - // keep only files with more than 1 line: - std::ifstream fichier(theInputNodeFiles[i].c_str()); - std::string s; - unsigned int nb_lines = 0; - while(std::getline(fichier,s)) ++nb_lines; -// std::cout << ". nb_lines: " << nb_lines << std::endl; - if ( nb_lines >= 2 ) { theNodeFiles.push_back( theInputNodeFiles[i] ); } - } -#ifdef _DEBUG_ - for ( size_t i = 0; i < theNodeFiles.size(); ++i ) - { std::cout << "Valid input node file #"< - mfMesh( MEDCoupling::MEDFileUMesh::New( theInputMedFile )); - if ( mfMesh.isNull() ) - throw std::invalid_argument( "Failed to read the input MED file: " + theInputMedFile ); - - MEDCoupling::DataArrayDouble * nodeCoords = mfMesh->getCoords(); - if ( !nodeCoords || nodeCoords->empty() ) - throw std::invalid_argument( "No nodes in the input mesh" ); - - - // read a geometry - -#ifdef _DEBUG_ - std::cout << "Lecture de la geometrie" << std::endl; -#endif - XAO::Xao xao; - if ( !xao.importXAO( theXaoFileName ) || !xao.getGeometry() ) - throw std::invalid_argument( "Failed to read the XAO input file: " + theXaoFileName ); - -#ifdef _DEBUG_ - std::cout << "Conversion en BREP" << std::endl; -#endif - XAO::BrepGeometry* xaoGeom = dynamic_cast( xao.getGeometry() ); - if ( !xaoGeom || xaoGeom->getTopoDS_Shape().IsNull() ) - throw std::invalid_argument( "Failed to get a BREP shape from the XAO input file" ); - - - // read groups of nodes and associate them with boundary shapes using names (no projection so far) - -#ifdef _DEBUG_ - std::cout << "Lecture des groupes" << std::endl; -#endif - FT_NodeGroups nodeGroups; - nodeGroups.read( theNodeFiles, &xao, nodeCoords ); -#ifdef _DEBUG_ - std::cout << "Nombre de groupes : " << nodeGroups.nbOfGroups() << std::endl; -#endif - - // project nodes to the boundary shapes and change their coordinates - -#ifdef _DEBUG_ - std::cout << "Projection des noeuds, theIsParallel=" << theIsParallel << std::endl; -#endif - OSD_Parallel::For( 0, nodeGroups.nbOfGroups(), nodeGroups, !theIsParallel ); - - // save the modified mesh - -#ifdef _DEBUG_ - std::cout << "Ecriture du maillage" << std::endl; -#endif - const int erase = 2; - mfMesh->write( theOutputMedFile, /*mode=*/erase ); - - if ( !nodeGroups.isOK() ) - throw std::runtime_error("Unable to project some nodes"); -} - - //================================================================================ - /*! - * \brief Initialize FT_Projector's with all sub-shapes of given type - * \param [in] theMainShape - the shape to explore - * \param [in] theSubType - the type of sub-shapes - * \param [out] theProjectors - the projectors - */ - //================================================================================ - - void getProjectors( const TopoDS_Shape& theMainShape, - const TopAbs_ShapeEnum theSubType, - std::vector< FT_Projector > & theProjectors ) - { - TopTools_IndexedMapOfShape subShapes; - TopExp::MapShapes( theMainShape, theSubType, subShapes ); -#ifdef _DEBUG_ - std::cout << ". Nombre de subShapes : " << subShapes.Size() << std::endl; -#endif - - theProjectors.resize( subShapes.Size() ); - for ( int i = 1; i <= subShapes.Size(); ++i ) - theProjectors[ i-1 ].setBoundaryShape( subShapes( i )); - } - -//================================================================================ -/*! - * \brief Load node groups from files - * \param [in] theNodeFiles - an array of names of files describing groups of nodes that - * will be moved onto geometry - * \param [in] theXaoGeom - the whole geometry to project on - * \param [inout] theNodeCoords - array of node coordinates - */ -//================================================================================ - -void FT_NodeGroups::read( const std::vector< std::string >& theNodeFiles, - const XAO::Xao* theXao, - MEDCoupling::DataArrayDouble* theNodeCoords ) -{ - // get projectors for all boundary sub-shapes; - // index of a projector in the vector corresponds to a XAO index of a sub-shape - XAO::BrepGeometry* xaoGeom = dynamic_cast( theXao->getGeometry() ); - getProjectors( xaoGeom->getTopoDS_Shape(), TopAbs_EDGE, _projectors[0] ); - getProjectors( xaoGeom->getTopoDS_Shape(), TopAbs_FACE, _projectors[1] ); - - _nodesOnGeom.resize( theNodeFiles.size() ); - - // read node IDs and look for projectors to boundary sub-shapes by group name - FT_Utils::XaoGroups xaoGroups( theXao ); - for ( size_t i = 0; i < theNodeFiles.size(); ++i ) - { - _nodesOnGeom[i].read( theNodeFiles[i], xaoGroups, theNodeCoords, _projectors ); - } -} - -//================================================================================ -/*! - * \brief Project and move nodes of a given group of nodes - */ -//================================================================================ - -void FT_NodeGroups::projectAndMove( const int groupIndex ) -{ - _nodesOnGeom[ groupIndex ].projectAndMove(); -} - -//================================================================================ -/*! - * \brief Return true if all nodes were successfully relocated - */ -//================================================================================ - -bool FT_NodeGroups::isOK() const -{ - for ( size_t i = 0; i < _nodesOnGeom.size(); ++i ) - if ( ! _nodesOnGeom[ i ].isOK() ) - return false; - - return true; -} - -//================================================================================ -/*! - * \brief Print some statistics on node groups - */ -//================================================================================ - -void FT_NodeGroups::dumpStat() const -{ - for ( size_t i = 0; i < _nodesOnGeom.size(); ++i ) - { - std::cout << _nodesOnGeom[i].getShapeDim() << "D " - << _nodesOnGeom[i].nbNodes() << " nodes" << std::endl; - } -} - - /*! - * \brief Close a file at destruction - */ - struct FileCloser - { - FILE * _file; - - FileCloser( FILE * file ): _file( file ) {} - ~FileCloser() { if ( _file ) ::fclose( _file ); } - }; - -//================================================================================ -/*! - * \brief Read node ids from a file and find shapes for projection - * \param [in] theNodeFile - a name of file holding IDs of nodes that - * will be moved onto geometry - * \param [in] theXaoGroups - a tool returning FT_Projector's by XAO group name - * \param [inout] theNodeCoords - array of node coordinates - * \param [in] theAllProjectorsByDim - all projectors of 2 dimensions, ordered so that - * a vector index corresponds to a XAO sub-shape ID - */ -//================================================================================ - -void FT_NodesOnGeom::read( const std::string& theNodeFile, - const FT_Utils::XaoGroups& theXaoGroups, - MEDCoupling::DataArrayDouble* theNodeCoords, - std::vector< FT_Projector > * theAllProjectorsByDim ) -{ - _nodeCoords = theNodeCoords; - - FILE * file = ::fopen( theNodeFile.c_str(), "r" ); - if ( !file ) - throw std::invalid_argument( "Can't open an input node file: " + theNodeFile ); - - FileCloser fileCloser( file ); - - // ------------------------------------- - // get shape dimension by the file name - // ------------------------------------- - - // hope the file name is something like "frnD.**" with n in (1,2) - int dimPos = theNodeFile.size() - 5; - if ( theNodeFile[ dimPos ] == '2' ) - _shapeDim = 2; - else if ( theNodeFile[ dimPos ] == '1' ) - _shapeDim = 1; - else - throw std::invalid_argument( "Can't define dimension by node file name " + theNodeFile ); -#ifdef _DEBUG_ - std::cout << ". Dimension of the file " << theNodeFile << ": " << _shapeDim << std::endl; -#endif - - // ------------------------------------- - // read geom group names; several lines - // ------------------------------------- - - std::vector< std::string > geomNames; - - const int maxLineLen = 256; - char line[ maxLineLen ]; - - long int pos = ::ftell( file ); - while ( ::fgets( line, maxLineLen, file )) // read a line - { - if ( ::feof( file )) - { - return; // no nodes in the file - } - - // check if the line describes node ids in format 3I10 (e.g. " 120 1 43\n") - size_t lineLen = strlen( line ); - if ( lineLen >= 31 && - ::isdigit( line[9] ) && - line[10] == ' ' && - ::isdigit( line[19] ) && - line[20] == ' ' && - ::isdigit( line[29] ) && - ::isspace( line[30] )) - break; - - geomNames.push_back( line + 1 ); // skip the 1st white space - - pos = ::ftell( file ); // remember the position to return if the next line holds node ids - } - - ::fseek( file, pos, SEEK_SET ); // return to the 1st line holding nodes ids - - - // -------------- - // read node ids - // -------------- - - FT_NodeToMove nodeIds; - std::vector< int > ids; - - const int nbNodes = theNodeCoords->getNumberOfTuples(); // to check validity of node IDs - - while ( ::fgets( line, maxLineLen, file )) // read a line - { - // find node ids in the line - - char *beg = line, *end = 0; - long int id; - - ids.clear(); - while (( id = ::strtol( beg, &end, 10 )) && - ( beg != end )) - { - ids.push_back( id ); - if ( id > nbNodes ) - throw std::invalid_argument( "Too large node ID: " + FT_Utils::toStr( id )); - beg = end; - } - - if ( ids.size() >= 3 ) - { - std::vector< int >::iterator i = ids.begin(); - nodeIds._nodeToMove = *i; - nodeIds._neighborNodes.assign( ++i, ids.end() ); - - _nodes.push_back( nodeIds ); - } - - if ( ::feof( file )) - break; - } - - // ----------------------------------------------------------------- - // try to find FT_Projector's to boundary sub-shapes by group names - // ----------------------------------------------------------------- - - _allProjectors = & theAllProjectorsByDim[ _shapeDim - 1 ]; - - _projectors.reserve( geomNames.size() ); - std::vector< const FT_Projector* > projectors; - - for ( size_t i = 0; i < geomNames.size(); ++i ) - { - std::string & groupName = geomNames[i]; -#ifdef _DEBUG_ - std::cout << ". Group name: " << groupName << std::endl; -#endif - - // remove trailing white spaces - for ( int iC = groupName.size() - 1; iC >= 0; --iC ) - { - if ( ::isspace( groupName[iC] ) ) - groupName.resize( iC ); - else - break; - } - if ( groupName.empty() ) - continue; - - _groupNames.push_back( groupName ); // keep _groupNames for easier debug :) - - // get projectors by group name - theXaoGroups.getProjectors( groupName, _shapeDim, - theAllProjectorsByDim[ _shapeDim-1 ], projectors ); - } - - // ------------------------------ - // check the found FT_Projector's - // ------------------------------ - - if ( projectors.size() == 1 ) - { - _projectors.push_back( *projectors[ 0 ]); - } - else - { - Bnd_Box nodesBox; - for ( size_t i = 0; i < _nodes.size(); ++i ) - nodesBox.Add( getPoint( _nodes[i]._nodeToMove )); - - if ( projectors.size() > 1 ) - { - // more than one boundary shape; - // try to filter off unnecessary projectors using a bounding box of nodes - for ( size_t i = 0; i < projectors.size(); ++i ) - if ( !nodesBox.IsOut( projectors[ i ]->getBoundingBox() )) - _projectors.push_back( *projectors[ i ]); - } - - if ( _projectors.empty() ) - { - // select projectors using a bounding box of nodes - std::vector< FT_Projector > & allProjectors = *_allProjectors; - for ( size_t i = 0; i < allProjectors.size(); ++i ) - if ( !nodesBox.IsOut( allProjectors[ i ].getBoundingBox() )) - _projectors.push_back( allProjectors[ i ]); - - if ( _projectors.empty() && !_nodes.empty() ) - throw std::runtime_error("No boundary shape found for nodes in file " + theNodeFile ); - } - } - - // prepare for projection - create real projectors - for ( size_t i = 0; i < _projectors.size(); ++i ) - _projectors[ i ].prepareForProjection(); - -} - -//================================================================================ -/*! - * \brief Project nodes to the shapes and move them to new positions - */ -//================================================================================ - -void FT_NodesOnGeom::projectAndMove() -{ - _OK = true; -// -// 1. Préalables -// - // check if all the shapes are planar - bool isAllPlanar = true; - for ( size_t i = 0; i < _projectors.size() && isAllPlanar; ++i ) - isAllPlanar = _projectors[i].isPlanarBoundary(); - if ( isAllPlanar ) - return; - - // set nodes in the order suitable for optimal projection - putNodesInOrder(); - - // project and move nodes - - std::vector< FT_NodeToMove* > notProjectedNodes; - size_t iP, iProjector; - gp_Pnt newXyz; - -#ifdef _DEBUG_ - std::cout << ".. _projectors.size() = " << _projectors.size() << std::endl; - std::cout << ".. _nodesOrder.size() = " << _nodesOrder.size() << std::endl; -#endif -// -// 2. Calculs -// 2.1. Avec plusieurs shapes -// - if ( _projectors.size() > 1 ) - { - // the nodes are to be projected onto several boundary shapes; - // in addition to the projecting, classification on a shape is necessary - // in order to find out on which of the shapes a node is to be projected - - iProjector = 0; - for ( size_t i = 0; i < _nodesOrder.size(); ++i ) - { - FT_NodeToMove& nn = _nodes[ _nodesOrder[ i ]]; - gp_Pnt xyz = getPoint( nn._nodeToMove ); - gp_Pnt xyz1 = getPoint( nn._neighborNodes[0] ); - gp_Pnt xyz2 = getPoint( nn._neighborNodes[1] ); - double maxDist2 = xyz1.SquareDistance( xyz2 ) / 4.; - if ( _projectors[ iProjector ].projectAndClassify( xyz, maxDist2, newXyz, - nn._params, nn._nearParams )) - { - moveNode( nn._nodeToMove, newXyz ); - } - else // a node is not on iProjector-th shape, find the shape it is on - { - for ( iP = 1; iP < _projectors.size(); ++iP ) // check _projectors other than iProjector - { - iProjector = ( iProjector + 1 ) % _projectors.size(); - if ( _projectors[ iProjector ].projectAndClassify( xyz, maxDist2, newXyz, - nn._params, nn._nearParams )) - { - moveNode( nn._nodeToMove, newXyz ); - break; - } - } - if ( iP == _projectors.size() ) - { - notProjectedNodes.push_back( &nn ); - -#ifdef _DEBUG_ - std::cerr << "Warning: no shape found for node " << nn._nodeToMove << std::endl; - if ( !_groupNames.empty() ) - std::cerr << "Warning: group -- " << _groupNames[0] << std::endl; -#endif - } - } - } - } -// -// 2.2. Avec une seule shape -// - else // one shape - { - for ( size_t i = 0; i < _nodesOrder.size(); ++i ) - { - FT_NodeToMove& nn = _nodes[ _nodesOrder[ i ]]; - gp_Pnt xyz = getPoint( nn._nodeToMove ); - gp_Pnt xyz1 = getPoint( nn._neighborNodes[0] ); - gp_Pnt xyz2 = getPoint( nn._neighborNodes[1] ); - -// maxDist2 : le quart du carré de la distance entre les deux voisins du noeud à bouger - double maxDist2 = xyz1.SquareDistance( xyz2 ) / 4.; -#ifdef _DEBUG_ - std::cout << "\n.. maxDist2 = " << maxDist2 << " entre " << nn._neighborNodes[0] << " et " << nn._neighborNodes[1] << " - milieu " << nn._nodeToMove << " - d/2 = " << sqrt(maxDist2) << " - d = " << sqrt(xyz1.SquareDistance( xyz2 )) << std::endl; -#endif - if ( _projectors[ 0 ].project( xyz, maxDist2, newXyz, - nn._params, nn._nearParams )) - moveNode( nn._nodeToMove, newXyz ); - else - notProjectedNodes.push_back( &nn ); - } - } -// -// 3. Bilan -// - if ( !notProjectedNodes.empty() ) - { - // project nodes that are not projected by any of _projectors; - // a proper projector is selected by evaluation of a distance between neighbor nodes - // and a shape - - std::vector< FT_Projector > & projectors = *_allProjectors; - - iProjector = 0; - for ( size_t i = 0; i < notProjectedNodes.size(); ++i ) - { - FT_NodeToMove& nn = *notProjectedNodes[ i ]; - gp_Pnt xyz = getPoint( nn._nodeToMove ); - gp_Pnt xyz1 = getPoint( nn._neighborNodes[0] ); - gp_Pnt xyz2 = getPoint( nn._neighborNodes[1] ); - double maxDist2 = xyz1.SquareDistance( xyz2 ) / 4.; - double tol2 = 1e-6 * maxDist2; - - bool ok; - for ( iP = 0; iP < projectors.size(); ++iP ) - { - projectors[ iProjector ].prepareForProjection(); - projectors[ iProjector ].tryWithoutPrevSolution( true ); - - if (( ok = projectors[ iProjector ].isOnShape( xyz1, tol2, nn._params, nn._nearParams )) && - ( ok = projectors[ iProjector ].isOnShape( xyz2, tol2, nn._params, nn._params ))) - { - if ( nn._neighborNodes.size() == 4 ) - { - gp_Pnt xyz1 = getPoint( nn._neighborNodes[2] ); - gp_Pnt xyz2 = getPoint( nn._neighborNodes[3] ); - if (( ok = projectors[ iProjector ].isOnShape( xyz1, tol2, nn._params, nn._params ))) - ok = projectors[ iProjector ].isOnShape( xyz2, tol2, nn._params, nn._params ); - } - } - - if ( ok && projectors[iProjector].project( xyz, maxDist2, newXyz, nn._params, nn._params )) - { - moveNode( nn._nodeToMove, newXyz ); - break; - } - iProjector = ( iProjector + 1 ) % projectors.size(); - } - if ( iP == projectors.size() ) - { - _OK = false; - - std::cerr << "Error: not projected node " << nn._nodeToMove << std::endl; - } - } - } -} - -//================================================================================ -/*! - * \brief Put nodes in the order for optimal projection and set FT_NodeToMove::_nearParams - * to point to a FT_NodeToMove::_params of a node that will be projected earlier - */ -//================================================================================ - -void FT_NodesOnGeom::putNodesInOrder() -{ - if ( !_nodesOrder.empty() ) - return; - - // check if any of projectors can use parameters of a previously projected node on a shape - // to speed up projection - - bool isPrevSolutionUsed = false; - for ( size_t i = 0; i < _projectors.size() && !isPrevSolutionUsed; ++i ) - isPrevSolutionUsed = _projectors[i].canUsePrevSolution(); - - if ( !isPrevSolutionUsed ) - { - _nodesOrder.resize( _nodes.size() ); - for ( size_t i = 0; i < _nodesOrder.size(); ++i ) - _nodesOrder[ i ] = i; - return; - } - - // make a map to find a neighbor projected node - - // map of { FT_NodeToMove::_neighborNodes[i] } to { FT_NodeToMove* }; - // here we call FT_NodeToMove a 'link' as this data links a _neighborNodes[i] node to other nodes - typedef NCollection_DataMap< int, std::vector< FT_NodeToMove* > > TNodeIDToLinksMap; - TNodeIDToLinksMap neigborsMap; - - int mapSize = ( _shapeDim == 1 ) ? _nodes.size() + 1 : _nodes.size() * 3; - neigborsMap.Clear(); - neigborsMap.ReSize( mapSize ); - - std::vector< FT_NodeToMove* > linkVec, *linkVecPtr; - const int maxNbLinks = ( _shapeDim == 1 ) ? 2 : 6; // usual nb of links - - for ( size_t i = 0; i < _nodes.size(); ++i ) - { - FT_NodeToMove& nn = _nodes[i]; - for ( size_t iN = 0; iN < nn._neighborNodes.size(); ++iN ) - { - if ( !( linkVecPtr = neigborsMap.ChangeSeek( nn._neighborNodes[ iN ] ))) - { - linkVecPtr = neigborsMap.Bound( nn._neighborNodes[ iN ], linkVec ); - linkVecPtr->reserve( maxNbLinks ); - } - linkVecPtr->push_back( & nn ); - } - } - - // fill in _nodesOrder - - _nodesOrder.reserve( _nodes.size() ); - - std::list< FT_NodeToMove* > queue; - queue.push_back( &_nodes[0] ); - _nodes[0]._nearParams = _nodes[0]._params; // to avoid re-adding to the queue - - while ( !queue.empty() ) - { - FT_NodeToMove* nn = queue.front(); - queue.pop_front(); - - _nodesOrder.push_back( nn - & _nodes[0] ); - - // add neighbors to the queue and set their _nearParams = nn->_params - for ( size_t iN = 0; iN < nn->_neighborNodes.size(); ++iN ) - { - std::vector< FT_NodeToMove* >& linkVec = neigborsMap( nn->_neighborNodes[ iN ]); - for ( size_t iL = 0; iL < linkVec.size(); ++iL ) - { - FT_NodeToMove* nnn = linkVec[ iL ]; - if ( nnn != nn && nnn->_nearParams == 0 ) - { - nnn->_nearParams = nn->_params; - queue.push_back( nnn ); - } - } - } - } - _nodes[0]._nearParams = 0; // reset -} - -//================================================================================ -/*! - * \brief Get node coordinates. Node IDs count from a unit - */ -//================================================================================ - -gp_Pnt FT_NodesOnGeom::getPoint( const int nodeID ) -{ - const size_t dim = _nodeCoords->getNumberOfComponents(); - const double * xyz = _nodeCoords->getConstPointer() + ( dim * ( nodeID - 1 )); - return gp_Pnt( xyz[0], xyz[1], dim == 2 ? 0 : xyz[2] ); -} - -//================================================================================ -/*! - * \brief change node coordinates - */ -//================================================================================ - -void FT_NodesOnGeom::moveNode( const int nodeID, const gp_Pnt& newXyz ) -{ - const size_t dim = _nodeCoords->getNumberOfComponents(); - double z, *xyz = _nodeCoords->getPointer() + ( dim * ( nodeID - 1 )); - newXyz.Coord( xyz[0], xyz[1], dim == 2 ? z : xyz[2] ); -} - -//----------------------------------------------------------------------------- -/*! - * \brief Root class of a projector of a point to a boundary shape - */ -struct FT_RealProjector -{ - virtual ~FT_RealProjector() {} - - /*! - * \brief Project a point to a boundary shape - * \param [in] point - the point to project - * \param [out] newSolution - position on the shape (U or UV) found during the projection - * \param [in] prevSolution - position already found during the projection of a neighbor point - * \return gp_Pnt - the projection point - */ - virtual gp_Pnt project( const gp_Pnt& point, - double* newSolution, - const double* prevSolution = 0) = 0; - - /*! - * \brief Project a point to a boundary shape and check if the projection is within - * the shape boundary - * \param [in] point - the point to project - * \param [in] maxDist2 - the maximal allowed square distance between point and projection - * \param [out] projection - the projection point - * \param [out] newSolution - position on the shape (U or UV) found during the projection - * \param [in] prevSolution - position already found during the projection of a neighbor point - * \return bool - false if the projection point lies out of the shape boundary or - the distance the point and the projection is more than sqrt(maxDist2) - */ - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) = 0; - - // return true if a previously found solution can be used to speed up the projection - - virtual bool canUsePrevSolution() const { return false; } - - - double _dist; // distance between the point being projected and its projection -}; - -namespace // actual projection algorithms -{ - const double theEPS = 1e-12; - - //================================================================================ - /*! - * \brief Projector to any curve - */ - //================================================================================ - - struct CurveProjector : public FT_RealProjector - { - BRepAdaptor_Curve _curve; - double _tol; - ShapeAnalysis_Curve _projector; - double _uRange[2]; - - //----------------------------------------------------------------------------- - CurveProjector( const TopoDS_Edge& e, const double tol ): - _curve( e ), _tol( tol ) - { - BRep_Tool::Range( e, _uRange[0], _uRange[1] ); - } - - //----------------------------------------------------------------------------- - // project a point to the curve - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { -#ifdef _DEBUG_ - std::cout << ".. project a point to the curve prevSolution = " << prevSolution << std::endl; -#endif - gp_Pnt proj; - Standard_Real param; - - if ( prevSolution ) - { - _dist = _projector.NextProject( prevSolution[0], _curve, P, _tol, proj, param ); - } - else - { - _dist = _projector.Project( _curve, P, _tol, proj, param, false ); - } -#ifdef _DEBUG_ - std::cout << ".. _dist : " << _dist << std::endl; -#endif - proj = _curve.Value( param ); - - newSolution[0] = param; - - return proj; - } - - //----------------------------------------------------------------------------- - // project a point to a curve and check if the projection is within the curve boundary - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { -#ifdef _DEBUG_ - std::cout << ".. project a point to a curve and check " << std::endl; -#endif - projection = project( point, newSolution, prevSolution ); - return ( _uRange[0] < newSolution[0] && newSolution[0] < _uRange[1] && - _dist * _dist < maxDist2 ); - } - - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return true; } - }; - - //================================================================================ - /*! - * \brief Projector to a straight curve. Don't project, classify only - */ - //================================================================================ - - struct LineProjector : public FT_RealProjector - { - gp_Pnt _p0, _p1; - - //----------------------------------------------------------------------------- - LineProjector( TopoDS_Edge e ) - { - e.Orientation( TopAbs_FORWARD ); - _p0 = BRep_Tool::Pnt( TopExp::FirstVertex( e )); - _p1 = BRep_Tool::Pnt( TopExp::LastVertex ( e )); - } - - //----------------------------------------------------------------------------- - // does nothing - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { - return P; - } - //----------------------------------------------------------------------------- - // check if a point lies within the line segment - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - gp_Vec edge( _p0, _p1 ); - gp_Vec p0p ( _p0, point ); - double u = ( edge * p0p ) / edge.SquareMagnitude(); // param [0,1] on the edge - projection = ( 1. - u ) * _p0.XYZ() + u * _p1.XYZ(); // projection of the point on the edge - if ( u < 0 || 1 < u ) - return false; - - // check distance - return point.SquareDistance( projection ) < theEPS * theEPS; - } - }; - - //================================================================================ - /*! - * \brief Projector to a circular edge - */ - //================================================================================ - - struct CircleProjector : public FT_RealProjector - { - gp_Circ _circle; - double _uRange[2]; - - //----------------------------------------------------------------------------- - CircleProjector( const gp_Circ& c, const double f, const double l ): - _circle( c ) - { - _uRange[0] = f; - _uRange[1] = l; - } - - //----------------------------------------------------------------------------- - // project a point to the circle - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { - // assume that P is already on the the plane of circle, since - // it is in the middle of two points lying on the circle - - // move P to the circle - const gp_Pnt& O = _circle.Location(); - gp_Vec radiusVec( O, P ); - double radius = radiusVec.Magnitude(); - if ( radius < std::numeric_limits::min() ) - return P; // P in on the axe - - gp_Pnt proj = O.Translated( radiusVec.Multiplied( _circle.Radius() / radius )); - - _dist = _circle.Radius() - radius; - - return proj; - } - - //----------------------------------------------------------------------------- - // project and check if a projection lies within the circular edge - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - _dist = -1; - projection = project( point, newSolution ); - if ( _dist < 0 || // ? - _dist * _dist > maxDist2 ) - return false; - - newSolution[0] = ElCLib::Parameter( _circle, projection ); - return ( _uRange[0] < newSolution[0] && newSolution[0] < _uRange[1] ); - } - }; - - //================================================================================ - /*! - * \brief Projector to any surface - */ - //================================================================================ - - struct SurfaceProjector : public FT_RealProjector - { - ShapeAnalysis_Surface _projector; - double _tol; - BRepTopAdaptor_FClass2d* _classifier; - - //----------------------------------------------------------------------------- - SurfaceProjector( const TopoDS_Face& face, const double tol, BRepTopAdaptor_FClass2d* cls ): - _projector( BRep_Tool::Surface( face )), - _tol( tol ), - _classifier( cls ) - { - } - //----------------------------------------------------------------------------- - // delete _classifier - ~SurfaceProjector() - { - delete _classifier; - } - - //----------------------------------------------------------------------------- - // project a point to a surface - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { - gp_Pnt2d uv; - - if ( prevSolution ) - { - gp_Pnt2d prevUV( prevSolution[0], prevSolution[1] ); - uv = _projector.NextValueOfUV( prevUV, P, _tol ); - } - else - { - uv = _projector.ValueOfUV( P, _tol ); - } - - uv.Coord( newSolution[0], newSolution[1] ); - - gp_Pnt proj = _projector.Value( uv ); - - _dist = _projector.Gap(); - - return proj; - } - - //----------------------------------------------------------------------------- - // project a point to a surface and check if the projection is within the surface boundary - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - projection = project( point, newSolution, prevSolution ); - return ( _dist * _dist < maxDist2 ) && classify( newSolution ); - } - - //----------------------------------------------------------------------------- - // check if the projection is within the shape boundary - bool classify( const double* newSolution ) - { - TopAbs_State state = _classifier->Perform( gp_Pnt2d( newSolution[0], newSolution[1]) ); - return ( state != TopAbs_OUT ); - } - - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return true; } - }; - - //================================================================================ - /*! - * \brief Projector to a plane. Don't project, classify only - */ - //================================================================================ - - struct PlaneProjector : public SurfaceProjector - { - gp_Pln _plane; - bool _isRealPlane; // false means that a surface is planar but parametrization is different - - //----------------------------------------------------------------------------- - PlaneProjector( const gp_Pln& pln, - const TopoDS_Face& face, - BRepTopAdaptor_FClass2d* cls, - bool isRealPlane=true): - SurfaceProjector( face, 0, cls ), - _plane( pln ), - _isRealPlane( isRealPlane ) - {} - - //----------------------------------------------------------------------------- - // does nothing - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { - return P; - } - //----------------------------------------------------------------------------- - // check if a point lies within the boundry of the planar face - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - if ( _isRealPlane ) - { - ElSLib::PlaneParameters( _plane.Position(), point, newSolution[0], newSolution[1]); - projection = ElSLib::PlaneValue ( newSolution[0], newSolution[1], _plane.Position() ); - if ( projection.SquareDistance( point ) > theEPS * theEPS ) - return false; - - return SurfaceProjector::classify( newSolution ); - } - else - { - return SurfaceProjector::projectAndClassify( point, maxDist2, projection, - newSolution, prevSolution ); - } - } - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return false; } - }; - - //================================================================================ - /*! - * \brief Projector to a cylinder - */ - //================================================================================ - - struct CylinderProjector : public SurfaceProjector - { - gp_Cylinder _cylinder; - - //----------------------------------------------------------------------------- - CylinderProjector( const gp_Cylinder& c, - const TopoDS_Face& face, - BRepTopAdaptor_FClass2d* cls ): - SurfaceProjector( face, 0, cls ), - _cylinder( c ) - {} - - //----------------------------------------------------------------------------- - // project a point to the cylinder - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { - // project the point P to the cylinder axis -> Pp - const gp_Pnt& O = _cylinder.Position().Location(); - const gp_Dir& axe = _cylinder.Position().Direction(); - gp_Vec trsl = gp_Vec( axe ).Multiplied( gp_Vec( O, P ).Dot( axe )); - gp_Pnt Pp = O.Translated( trsl ); - - // move Pp to the cylinder - gp_Vec radiusVec( Pp, P ); - double radius = radiusVec.Magnitude(); - if ( radius < std::numeric_limits::min() ) - return P; // P in on the axe - - gp_Pnt proj = Pp.Translated( radiusVec.Multiplied( _cylinder.Radius() / radius )); - - _dist = _cylinder.Radius() - radius; - - return proj; - } - //----------------------------------------------------------------------------- - // project a point to the cylinder and check if the projection is within the surface boundary - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - ElSLib::CylinderParameters( _cylinder.Position(), _cylinder.Radius(), point, - newSolution[0], newSolution[1]); - projection = ElSLib::CylinderValue( newSolution[0], newSolution[1], - _cylinder.Position(), _cylinder.Radius() ); - _dist = point.Distance( projection ); - - return ( _dist * _dist < maxDist2 ) && SurfaceProjector::classify( newSolution ); - } - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return false; } - }; - - //================================================================================ - /*! - * \brief Projector to a cone - */ - //================================================================================ - - struct ConeProjector : public SurfaceProjector - { - gp_Cone _cone; - - //----------------------------------------------------------------------------- - ConeProjector( const gp_Cone& c, - const TopoDS_Face& face, - BRepTopAdaptor_FClass2d* cls ): - SurfaceProjector( face, 0, cls ), - _cone( c ) - {} - - //----------------------------------------------------------------------------- - // project a point to the cone - virtual gp_Pnt project( const gp_Pnt& point, - double* newSolution, - const double* prevSolution = 0) - { - ElSLib::ConeParameters( _cone.Position(), _cone.RefRadius(), _cone.SemiAngle(), - point, newSolution[0], newSolution[1]); - gp_Pnt proj = ElSLib::ConeValue( newSolution[0], newSolution[1], - _cone.Position(), _cone.RefRadius(), _cone.SemiAngle() ); - _dist = point.Distance( proj ); - - return proj; - } - - //----------------------------------------------------------------------------- - // project a point to the cone and check if the projection is within the surface boundary - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - projection = project( point, newSolution, prevSolution ); - - return ( _dist * _dist < maxDist2 ) && SurfaceProjector::classify( newSolution ); - } - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return false; } - }; - - //================================================================================ - /*! - * \brief Projector to a sphere - */ - //================================================================================ - - struct SphereProjector : public SurfaceProjector - { - gp_Sphere _sphere; - - //----------------------------------------------------------------------------- - SphereProjector( const gp_Sphere& s, - const TopoDS_Face& face, - BRepTopAdaptor_FClass2d* cls ): - SurfaceProjector( face, 0, cls ), - _sphere( s ) - {} - - //----------------------------------------------------------------------------- - // project a point to the sphere - virtual gp_Pnt project( const gp_Pnt& P, - double* newSolution, - const double* prevSolution = 0) - { - // move Pp to the Sphere - const gp_Pnt& O = _sphere.Location(); - gp_Vec radiusVec( O, P ); - double radius = radiusVec.Magnitude(); - if ( radius < std::numeric_limits::min() ) - return P; // P is on O - - gp_Pnt proj = O.Translated( radiusVec.Multiplied( _sphere.Radius() / radius )); - - _dist = _sphere.Radius() - radius; - - return proj; - } - - //----------------------------------------------------------------------------- - // project a point to the sphere and check if the projection is within the surface boundary - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - ElSLib::SphereParameters( _sphere.Position(), _sphere.Radius(), point, - newSolution[0], newSolution[1]); - projection = ElSLib::SphereValue( newSolution[0], newSolution[1], - _sphere.Position(), _sphere.Radius() ); - - return ( _dist * _dist < maxDist2 ) && SurfaceProjector::classify( newSolution ); - } - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return false; } - }; - - //================================================================================ - /*! - * \brief Projector to a torus - */ - //================================================================================ - - struct TorusProjector : public SurfaceProjector - { - gp_Torus _torus; - - //----------------------------------------------------------------------------- - TorusProjector( const gp_Torus& t, - const TopoDS_Face& face, - BRepTopAdaptor_FClass2d* cls ): - SurfaceProjector( face, 0, cls ), - _torus( t ) - {} - - //----------------------------------------------------------------------------- - // project a point to the torus - virtual gp_Pnt project( const gp_Pnt& point, - double* newSolution, - const double* prevSolution = 0) - { - ElSLib::TorusParameters( _torus.Position(), _torus.MajorRadius(), _torus.MinorRadius(), - point, newSolution[0], newSolution[1]); - gp_Pnt proj = ElSLib::TorusValue( newSolution[0], newSolution[1], - _torus.Position(), _torus.MajorRadius(), _torus.MinorRadius() ); - _dist = point.Distance( proj ); - - return proj; - } - - //----------------------------------------------------------------------------- - // project a point to the torus and check if the projection is within the surface boundary - virtual bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0) - { - projection = project( point, newSolution, prevSolution ); - - return ( _dist * _dist < maxDist2 ) && SurfaceProjector::classify( newSolution ); - } - //----------------------------------------------------------------------------- - // return true if a previously found solution can be used to speed up the projection - virtual bool canUsePrevSolution() const { return false; } - }; - - //================================================================================ - /*! - * \brief Check if a curve can be considered straight - */ - //================================================================================ - - bool isStraight( const GeomAdaptor_Curve& curve, const double tol ) - { - // rough check: evaluate how far from a straight line connecting the curve ends - // stand several internal points of the curve - - const double f = curve.FirstParameter(); - const double l = curve.LastParameter(); - const gp_Pnt pf = curve.Value( f ); - const gp_Pnt pl = curve.Value( l ); - const gp_Vec lineVec( pf, pl ); - const double lineLen2 = lineVec.SquareMagnitude(); - if ( lineLen2 < std::numeric_limits< double >::min() ) - return false; // E seems closed - - const double nbSamples = 7; - for ( int i = 0; i < nbSamples; ++i ) - { - const double r = ( i + 1 ) / nbSamples; - const gp_Pnt pi = curve.Value( f * r + l * ( 1 - r )); - const gp_Vec vi( pf, pi ); - const double h2 = lineVec.Crossed( vi ).SquareMagnitude() / lineLen2; - if ( h2 > tol * tol ) - return false; - } - - // thorough check - GCPnts_UniformDeflection divider( curve, tol ); - return ( divider.IsDone() && divider.NbPoints() < 3 ); - } -} - -//================================================================================ -/*! - * \brief Initialize with a boundary shape - */ -//================================================================================ - -FT_Projector::FT_Projector(const TopoDS_Shape& shape) -{ - _realProjector = 0; - setBoundaryShape( shape ); - _tryWOPrevSolution = false; -} - -//================================================================================ -/*! - * \brief Copy another projector - */ -//================================================================================ - -FT_Projector::FT_Projector(const FT_Projector& other) -{ - _realProjector = 0; - _shape = other._shape; - _bndBox = other._bndBox; - _tryWOPrevSolution = false; -} - -//================================================================================ -/*! - * \brief Destructor. Delete _realProjector - */ -//================================================================================ - -FT_Projector::~FT_Projector() -{ - delete _realProjector; -} - -//================================================================================ -/*! - * \brief Initialize with a boundary shape. Compute the bounding box - */ -//================================================================================ - -void FT_Projector::setBoundaryShape(const TopoDS_Shape& shape) -{ - delete _realProjector; _realProjector = 0; - _shape = shape; - if ( shape.IsNull() ) - return; - - BRepBndLib::Add( shape, _bndBox ); - _bndBox.Enlarge( 1e-5 * sqrt( _bndBox.SquareExtent() )); -} - -//================================================================================ -/*! - * \brief Create a real projector - */ -//================================================================================ - -void FT_Projector::prepareForProjection() -{ - if ( _shape.IsNull() || _realProjector ) - return; - - if ( _shape.ShapeType() == TopAbs_EDGE ) - { - const TopoDS_Edge& edge = TopoDS::Edge( _shape ); - - double tol = 1e-6 * sqrt( _bndBox.SquareExtent() ); - - double f,l; - Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, f,l ); - if ( curve.IsNull() ) - return; // degenerated edge - - GeomAdaptor_Curve acurve( curve, f, l ); - switch ( acurve.GetType() ) - { - case GeomAbs_Line: - _realProjector = new LineProjector( edge ); - break; - case GeomAbs_Circle: - _realProjector = new CircleProjector( acurve.Circle(), f, l ); - break; - case GeomAbs_BezierCurve: - case GeomAbs_BSplineCurve: - case GeomAbs_OffsetCurve: - case GeomAbs_OtherCurve: - if ( isStraight( acurve, tol )) - { - _realProjector = new LineProjector( edge ); - break; - } - case GeomAbs_Ellipse: - case GeomAbs_Hyperbola: - case GeomAbs_Parabola: - _realProjector = new CurveProjector( edge, tol ); - } - } - else if ( _shape.ShapeType() == TopAbs_FACE ) - { - TopoDS_Face face = TopoDS::Face( _shape ); - - Handle(Geom_Surface) surface = BRep_Tool::Surface( face ); - if ( surface.IsNull() ) - return; - - GeomAdaptor_Surface asurface( surface ); - Standard_Real tol = BRep_Tool::Tolerance( face ); - Standard_Real toluv = Min( asurface.UResolution( tol ), asurface.VResolution( tol )); - BRepTopAdaptor_FClass2d* classifier = new BRepTopAdaptor_FClass2d( face, toluv ); - - switch ( asurface.GetType() ) - { - case GeomAbs_Plane: - _realProjector = new PlaneProjector( asurface.Plane(), face, classifier ); - break; - case GeomAbs_Cylinder: - _realProjector = new CylinderProjector( asurface.Cylinder(), face, classifier ); - break; - case GeomAbs_Sphere: - _realProjector = new SphereProjector( asurface.Sphere(), face, classifier ); - break; - case GeomAbs_Cone: - _realProjector = new ConeProjector( asurface.Cone(), face, classifier ); - break; - case GeomAbs_Torus: - _realProjector = new TorusProjector( asurface.Torus(), face, classifier ); - break; - case GeomAbs_BezierSurface: - case GeomAbs_BSplineSurface: - case GeomAbs_SurfaceOfRevolution: - case GeomAbs_SurfaceOfExtrusion: - case GeomAbs_OffsetSurface: - case GeomAbs_OtherSurface: - GeomLib_IsPlanarSurface isPlaneCheck( surface, tol ); - if ( isPlaneCheck.IsPlanar() ) - { - _realProjector = new PlaneProjector( isPlaneCheck.Plan(), face, classifier, - /*isRealPlane=*/false); - } - else - { - _realProjector = new SurfaceProjector( face, tol, classifier ); - } - break; - } - - if ( !_realProjector ) - delete classifier; - } -} - -//================================================================================ -/*! - * \brief Return true if projection is not needed - */ -//================================================================================ - -bool FT_Projector::isPlanarBoundary() const -{ - return ( dynamic_cast< LineProjector* >( _realProjector ) || - dynamic_cast< PlaneProjector* >( _realProjector ) ); -} - -//================================================================================ -/*! - * \brief Check if a point lies on the boundary shape - * \param [in] point - the point to check - * \param [in] tol2 - a square tolerance allowing to decide whether a point is on the shape - * \param [in] newSolution - position on the shape (U or UV) of the point found - * during projecting - * \param [in] prevSolution - position on the shape (U or UV) of a neighbor point - * \return bool - \c true if the point lies on the boundary shape - * - * This method is used to select a shape by checking if all neighbor nodes of a node to move - * lie on a shape. - */ -//================================================================================ - -bool FT_Projector::isOnShape( const gp_Pnt& point, - const double tol2, - double* newSolution, - const double* prevSolution) -{ - if ( _bndBox.IsOut( point ) || !_realProjector ) - return false; - - gp_Pnt proj; - if ( isPlanarBoundary() ) - return projectAndClassify( point, tol2, proj, newSolution, prevSolution ); - - return project( point, tol2, proj, newSolution, prevSolution ); -} - -//================================================================================ -/*! - * \brief Project a point to the boundary shape - * \param [in] point - the point to project - * \param [in] maxDist2 - the maximal square distance between the point and the projection - * \param [out] projection - the projection - * \param [out] newSolution - position on the shape (U or UV) of the point found - * during projecting - * \param [in] prevSolution - already found position on the shape (U or UV) of a neighbor point - * \return bool - false if the distance between the point and the projection - * is more than sqrt(maxDist2) - * - * This method is used to project a node in the case where only one shape is found by name - */ -//================================================================================ - -bool FT_Projector::project( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution) -{ - if ( !_realProjector ) - return false; - - _realProjector->_dist = 1e100; - projection = _realProjector->project( point, newSolution, prevSolution ); - - bool ok = ( _realProjector->_dist * _realProjector->_dist < maxDist2 ); - if ( !ok && _tryWOPrevSolution && prevSolution ) - { - projection = _realProjector->project( point, newSolution ); - ok = ( _realProjector->_dist * _realProjector->_dist < maxDist2 ); - } - return ok; -} - -//================================================================================ -/*! - * \brief Project a point to the boundary shape and check if the projection lies within - * the shape boundary - * \param [in] point - the point to project - * \param [in] maxDist2 - the maximal square distance between the point and the projection - * \param [out] projection - the projection - * \param [out] newSolution - position on the shape (U or UV) of the point found - * during projecting - * \param [in] prevSolution - already found position on the shape (U or UV) of a neighbor point - * \return bool - false if the projection point lies out of the shape boundary or - * the distance between the point and the projection is more than sqrt(maxDist2) - * - * This method is used to project a node in the case where several shapes are selected for - * projection of a node group - */ -//================================================================================ - -bool FT_Projector::projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution) -{ - if ( _bndBox.IsOut( point ) || !_realProjector ) - return false; - - bool ok = _realProjector->projectAndClassify( point, maxDist2, projection, - newSolution, prevSolution ); - if ( !ok && _tryWOPrevSolution && prevSolution ) - ok = _realProjector->projectAndClassify( point, maxDist2, projection, newSolution ); - - return ok; -} - -//================================================================================ -/*! - * \brief Return true if a previously found solution can be used to speed up the projection - */ -//================================================================================ - -bool FT_Projector::canUsePrevSolution() const -{ - return ( _realProjector && _realProjector->canUsePrevSolution() ); -} - -//================================================================================ -/* - * \brief Check if a file exists - */ -//================================================================================ - -bool FT_Utils::fileExists( const std::string& path ) -{ - if ( path.empty() ) - return false; - - boost::system::error_code err; - bool res = boofs::exists( path, err ); - - return err ? false : res; -} - -//================================================================================ -/*! - * \brief Check if a file can be created/overwritten - */ -//================================================================================ - -bool FT_Utils::canWrite( const std::string& path ) -{ - if ( path.empty() ) - return false; - - bool can = false; -#ifdef WIN32 - - HANDLE file = CreateFile( path.c_str(), // name of the write - GENERIC_WRITE, // open for writing - 0, // do not share - NULL, // default security - OPEN_ALWAYS, // CREATE NEW or OPEN EXISTING - FILE_ATTRIBUTE_NORMAL, // normal file - NULL); // no attr. template - can = ( file != INVALID_HANDLE_VALUE ); - CloseHandle( file ); - -#else - - int file = ::open( path.c_str(), - O_WRONLY | O_CREAT, - S_IRUSR | S_IWUSR | S_IRGRP | S_IROTH ); // rw-r--r-- - can = ( file >= 0 ); - -#endif - - return can; -} - -//================================================================================ -/*! - * \brief Make a map of XAO groups - */ -//================================================================================ - -FT_Utils::XaoGroups::XaoGroups( const XAO::Xao* theXao ) -{ - XAO::Xao* xao = const_cast< XAO::Xao* >( theXao ); - - for ( int iG = 0; iG < theXao->countGroups(); ++iG ) - { - XAO::Group* group = xao->getGroup( iG ); - - if ( group->getDimension() == 1 ) - - _xaoGroups[ 0 ].insert( std::make_pair( group->getName(), group )); - - else if ( group->getDimension() == 2 ) - - _xaoGroups[ 1 ].insert( std::make_pair( group->getName(), group )); - } -} - -//================================================================================ -/*! - * \brief Return FT_Projector's by a group name - * \param [in] groupName - the group name - * \param [in] dim - the group dimension - * \param [in] allProjectors - the projector of all shapes of \a dim dimension - * \param [out] groupProjectors - projectors to shapes of the group - * \return int - number of found shapes - */ -//================================================================================ - -int FT_Utils::XaoGroups::getProjectors( const std::string& groupName, - const int dim, - const std::vector< FT_Projector > & allProjectors, - std::vector< const FT_Projector* > & groupProjectors) const -{ - // get namesake groups - - const TGroupByNameMap* groupMap = 0; - if ( dim == 1 ) - groupMap = &_xaoGroups[ 0 ]; - else if ( dim == 2 ) - groupMap = &_xaoGroups[ 1 ]; - else - return 0; - - TGroupByNameMap::const_iterator name2gr = groupMap->find( groupName ); - if ( name2gr == groupMap->end() ) - return 0; - - std::vector< XAO::Group* > groups; - groups.push_back( name2gr->second ); - - for ( ++name2gr; name2gr != groupMap->end(); ++name2gr ) - { - if ( name2gr->second->getName() == groupName ) - groups.push_back( name2gr->second ); - else - break; - } - - // get projectors - - int nbFound = 0; - for ( size_t i = 0; i < groups.size(); ++i ) - { - // IDs in XAO correspond to indices of allProjectors - std::set::iterator id = groups[i]->begin(), end = groups[i]->end(); - for ( ; id != end; ++id, ++nbFound ) - if ( *id < (int) allProjectors.size() ) - groupProjectors.push_back ( & allProjectors[ *id ]); - } - - return nbFound; -} - } // namespace SMESHHOMARDImpl /end/ diff --git a/src/SMESH/SMESH_Homard.hxx b/src/SMESH/SMESH_Homard.hxx index 416fe3937..45be48424 100644 --- a/src/SMESH/SMESH_Homard.hxx +++ b/src/SMESH/SMESH_Homard.hxx @@ -53,29 +53,14 @@ #include #include #include -#include #include #include -#include -#include -#include -#include - #if defined WIN32 #pragma warning ( disable: 4251 ) #endif -namespace MEDCoupling { - class DataArrayDouble; -} -namespace XAO { - class Xao; - class Group; - class BrepGeometry; -} - namespace SMESHHOMARDImpl { @@ -487,240 +472,6 @@ private: std::list _ListFieldInterpTSR; }; -// HOMARD/FrontTrack - -class FrontTrack -{ -public: - - /*! - * \brief Relocate nodes to lie on geometry - * \param [in] theInputMedFile - a MED file holding a mesh including nodes that will be - * moved onto the geometry - * \param [in] theOutputMedFile - a MED file to create, that will hold a modified mesh - * \param [in] theInputNodeFiles - an array of names of files describing groups of nodes that - * will be moved onto the geometry - * \param [in] theXaoFileName - a path to a file in XAO format holding the geometry and - * the geometrical groups. - * \param [in] theIsParallel - if \c true, all processors are used to treat boundary shapes - * in parallel. - */ - void track( const std::string& theInputMedFile, - const std::string& theOutputMedFile, - const std::vector< std::string > & theInputNodeFiles, - const std::string& theXaoFileName, - bool theIsParallel=true); - -}; - - -struct FT_RealProjector; - -/*! - * \brief Projector of a point to a boundary shape. Wrapper of a real projection algo - */ -class FT_Projector -{ -public: - - FT_Projector(const TopoDS_Shape& shape = TopoDS_Shape()); - FT_Projector(const FT_Projector& other); - ~FT_Projector(); - - // initialize with a boundary shape, compute the bounding box - void setBoundaryShape(const TopoDS_Shape& shape); - - // return the boundary shape - const TopoDS_Shape& getShape() const { return _shape; } - - // return the bounding box - const Bnd_Box getBoundingBox() const { return _bndBox; } - - - // create a real projector - void prepareForProjection(); - - // return true if a previously found solution can be used to speed up the projection - bool canUsePrevSolution() const; - - // return true if projection is not needed - bool isPlanarBoundary() const; - - - // switch a mode of usage of prevSolution. - // If projection fails, to try to project without usage of prevSolution. - // By default this mode is off - void tryWithoutPrevSolution( bool toTry ) { _tryWOPrevSolution = toTry; } - - // project a point to the boundary shape - bool project( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0); - - // project a point to the boundary shape and check if the projection is within the shape boundary - bool projectAndClassify( const gp_Pnt& point, - const double maxDist2, - gp_Pnt& projection, - double* newSolution, - const double* prevSolution = 0); - - // check if a point lies on the boundary shape - bool isOnShape( const gp_Pnt& point, - const double tol2, - double* newSolution, - const double* prevSolution = 0); - -private: - - FT_RealProjector* _realProjector; - Bnd_Box _bndBox; - TopoDS_Shape _shape; - bool _tryWOPrevSolution; -}; - -namespace FT_Utils -{ - // Check if a file exists - bool fileExists( const std::string& path ); - - // Check if a file can be created/overwritten - bool canWrite( const std::string& path ); - - // Transform anything printable to a string - template< typename T> std::string toStr( const T& t ) - { - std::ostringstream s; - s << t; - return s.str(); - } - - //-------------------------------------------------------------------------------------------- - /*! - * \brief Return projectors by group name - */ - struct XaoGroups - { - XaoGroups( const XAO::Xao* xao ); - - int getProjectors( const std::string& groupName, - const int dim, - const std::vector< FT_Projector > & allProjectors, - std::vector< const FT_Projector* > & groupProjectors ) const; - private: - - typedef std::multimap< std::string, XAO::Group* > TGroupByNameMap; - TGroupByNameMap _xaoGroups[ 2 ]; // by dim - }; -} // namespace FT_Utils - -/*! - * \brief Node group and geometry to project onto - */ -class FT_NodesOnGeom -{ -public: - - // read node IDs form a file and try to find a boundary sub-shape by name - void read( const std::string& nodesFile, - const FT_Utils::XaoGroups& xaoGroups, - MEDCoupling::DataArrayDouble* nodeCoords, - std::vector< FT_Projector > * allProjectorsByDim); - - // chose boundary shapes by evaluating distance between nodes and shapes - //void choseShape( const std::vector< FT_Utils::ShapeAndBndBox >& shapeAndBoxList ); - - // project nodes to the shapes and move them to new positions - void projectAndMove(); - - // return true if all nodes were successfully relocated - bool isOK() const { return _OK; } - - // return dimension of boundary shapes - int getShapeDim() const { return _shapeDim; } - - // return nb of nodes to move - int nbNodes() const { return _nodes.size(); } - - -private: - - // put nodes in the order for optimal projection - void putNodesInOrder(); - - // get node coordinates - gp_Pnt getPoint( const int nodeID ); - - // change node coordinates - void moveNode( const int nodeID, const gp_Pnt& xyz ); - - - // Ids of a node to move and its 2 or 4 neighbors - struct FT_NodeToMove - { - int _nodeToMove; - std::vector< int > _neighborNodes; - - double _params[2]; // parameters on shape (U or UV) found by projection - double *_nearParams; // _params of a neighbor already projected node - - FT_NodeToMove(): _nearParams(0) {} - }; - - std::vector< std::string > _groupNames; - int _shapeDim; // dimension of boundary shapes - std::vector< FT_NodeToMove > _nodes; // ids of nodes to move and their neighbors - std::vector< FT_Projector > _projectors; // FT_Projector's initialized with boundary shapes - std::vector< FT_Projector > * _allProjectors; // FT_Projector's for all shapes of _shapeDim - MEDCoupling::DataArrayDouble* _nodeCoords; - bool _OK; // projecting is successful - - // map of { FT_NodeToMove::_neighborNodes[i] } to { FT_NodeToMove* } - // this map is used to find neighbor nodes - typedef NCollection_DataMap< int, std::vector< FT_NodeToMove* > > TNodeIDToLinksMap; - TNodeIDToLinksMap _neigborsMap; - std::vector _nodesOrder; - -}; - -/*! - * \brief Container of node groups. - */ -class FT_NodeGroups -{ -public: - - // Load node groups from files - void read( const std::vector< std::string >& nodeFiles, - const XAO::Xao* xaoGeom, - MEDCoupling::DataArrayDouble* nodeCoords ); - - // return number of groups of nodes to move - int nbOfGroups() const { return _nodesOnGeom.size(); } - - // Move nodes of a group in parallel mode - void operator() ( const int groupIndex ) const - { - const_cast< FT_NodeGroups* >( this )->projectAndMove( groupIndex ); - } - - // Project and move nodes of a given group of nodes - void projectAndMove( const int groupIndex ); - - // return true if all nodes were successfully relocated - bool isOK() const; - - // print some statistics on node groups - void dumpStat() const; - -private: - - std::vector< FT_NodesOnGeom > _nodesOnGeom; - std::vector< FT_Projector > _projectors[2]; // curves and surfaces separately - -}; - }; // namespace SMESHHOMARDImpl #endif diff --git a/src/SMESHGUI/SMESHGUI_HomardAdaptDlg.cxx b/src/SMESHGUI/SMESHGUI_HomardAdaptDlg.cxx index c215ac0c8..dbd6ec225 100644 --- a/src/SMESHGUI/SMESHGUI_HomardAdaptDlg.cxx +++ b/src/SMESHGUI/SMESHGUI_HomardAdaptDlg.cxx @@ -423,7 +423,7 @@ bool SMESHGUI_HomardAdaptDlg::PushOnApply() if (anOutMed.isEmpty()) { // store in working directory and with default name QString aWorkingDir = myAdvOpt->workingDirectoryLineEdit->text().trimmed(); - QFileInfo aFileInfo (QDir(aWorkingDir), "Uniform_01_R.med"); + QFileInfo aFileInfo (QDir(aWorkingDir), "Uniform_R.med"); anOutMed = aFileInfo.absoluteFilePath(); // show it myArgs->mySelectOutMedFileLineEdit->setText(anOutMed); @@ -437,7 +437,7 @@ bool SMESHGUI_HomardAdaptDlg::PushOnApply() else { // Set file name without path for it to be created in current directory // (it will be iteration's dir, and it will be destroyed after) - aMeshFileOUT = "Uniform_01_R.med"; + aMeshFileOUT = "Uniform_R.med"; } myHomardGen->SetMeshFileOUT(aMeshFileOUT.c_str()); @@ -561,8 +561,7 @@ void SMESHGUI_HomardAdaptDlg::selectionChanged() myArgs->myInBrowserObject->setText(aMeshName); // Out mesh name default value - // TODO: add some suffix? "_R" or "_UnifRefin", or "_Uniform_01_R" - myArgs->myOutMeshNameLineEdit->setText(aMeshName); + myArgs->myOutMeshNameLineEdit->setText(aMeshName + "_Uniform_R"); // Output med file default value // Construct it from Input mesh name and working directory @@ -572,9 +571,9 @@ void SMESHGUI_HomardAdaptDlg::selectionChanged() } else { QString aWorkingDir = myAdvOpt->workingDirectoryLineEdit->text().trimmed(); - QFileInfo aFileInfo (QDir(aWorkingDir), aMeshName + QString("_Uniform_01_R.med")); + QFileInfo aFileInfo (QDir(aWorkingDir), aMeshName + QString("_Uniform_R.med")); for (int ii = 1; aFileInfo.exists(); ii++) { - QString anUniqueName = QString("%1_Uniform_01_R_%2.med").arg(aMeshName).arg(ii); + QString anUniqueName = QString("%1_Uniform_R_%2.med").arg(aMeshName).arg(ii); aFileInfo.setFile(QDir(aWorkingDir), anUniqueName); } myArgs->mySelectOutMedFileLineEdit->setText(aFileInfo.absoluteFilePath()); @@ -601,9 +600,8 @@ void SMESHGUI_HomardAdaptDlg::SetFileName() myArgs->mySelectInMedFileLineEdit->setText(fileName); // Out mesh name default value - // TODO: add some suffix? "_R" or "_UnifRefin", or "_Uniform_01_R" QString aMeshName = SMESH_HOMARD_QT_COMMUN::LireNomMaillage(fileName); - myArgs->myOutMeshNameLineEdit->setText(aMeshName); + myArgs->myOutMeshNameLineEdit->setText(aMeshName + "_Uniform_R"); // Output med file default value // Construct it from Input med file name and path @@ -613,9 +611,9 @@ void SMESHGUI_HomardAdaptDlg::SetFileName() if (lastdot != std::string::npos) fname = fname.substr(0, lastdot); QString fileNameOut = fname.c_str(); - QFileInfo aFileInfo (fileNameOut + QString("_Uniform_01_R.med")); + QFileInfo aFileInfo (fileNameOut + QString("_Uniform_R.med")); for (int ii = 1; aFileInfo.exists(); ii++) { - QString anUniqueName = QString("%1_Uniform_01_R_%2.med").arg(fileNameOut).arg(ii); + QString anUniqueName = QString("%1_Uniform_R_%2.med").arg(fileNameOut).arg(ii); aFileInfo.setFile(anUniqueName); } myArgs->mySelectOutMedFileLineEdit->setText(aFileInfo.absoluteFilePath()); diff --git a/src/SMESH_I/SMESH_Homard_i.cxx b/src/SMESH_I/SMESH_Homard_i.cxx index 19bfe7a68..130ca462d 100644 --- a/src/SMESH_I/SMESH_Homard_i.cxx +++ b/src/SMESH_I/SMESH_Homard_i.cxx @@ -34,6 +34,10 @@ #include "SALOME_LifeCycleCORBA.hxx" #include "SALOMEconfig.h" +// Have to be included before std headers +#include +#include + #include #include #include @@ -453,7 +457,7 @@ SMESHHOMARD::extrema* HOMARD_Cas_i::GetBoundingBox() std::vector LesExtremes = myHomardCas->GetBoundingBox(); ASSERT(LesExtremes.size() == 10); aResult->length(10); - for (int i = 0; i < (int)LesExtremes.size(); i++) { + for (unsigned int i = 0; i < LesExtremes.size(); i++) { aResult[i] = LesExtremes[i]; } return aResult._retn(); @@ -470,7 +474,7 @@ void HOMARD_Cas_i::SetGroups(const SMESHHOMARD::ListGroupType& ListGroup) { ASSERT(myHomardCas); std::list ListString; - for (int i = 0; i < ListGroup.length(); i++) + for (unsigned int i = 0; i < ListGroup.length(); i++) { ListString.push_back(std::string(ListGroup[i])); } @@ -957,6 +961,7 @@ HOMARD_Gen_i::HOMARD_Gen_i() : SALOME::GenericObj_i(SMESH_Gen_i::GetPOA()), //============================================================================= HOMARD_Gen_i::~HOMARD_Gen_i() { + MESSAGE ("HOMARD_Gen_i::~HOMARD_Gen_i()"); if (!myCase->_is_nil()) { CleanCase(); } @@ -1005,8 +1010,8 @@ CORBA::Long HOMARD_Gen_i::DeleteIteration(int numIter) MESSAGE ("DeleteIteration : numIter = " << numIter); if (numIter == 0) { - if (!CORBA::is_nil(myIteration1)) DeleteIteration(1); - myIteration0 = SMESHHOMARD::HOMARD_Iteration::_nil(); + if (CORBA::is_nil(myIteration1)) + myIteration0 = SMESHHOMARD::HOMARD_Iteration::_nil(); } else { if (!CORBA::is_nil(myIteration1)) { @@ -1088,18 +1093,28 @@ void HOMARD_Gen_i::AssociateCaseIter(int numIter, const char* labelIter) throw SALOME::SALOME_Exception(es); } - SMESHHOMARD::HOMARD_Iteration_var myIteration; - if (numIter == 0) myIteration = myIteration0; - else myIteration = myIteration1; - if (CORBA::is_nil(myIteration)) { - SALOME::ExceptionStruct es; - es.type = SALOME::BAD_PARAM; - es.text = "Invalid iteration"; - throw SALOME::SALOME_Exception(es); + if (numIter == 0) { + if (CORBA::is_nil(myIteration0)) { + SALOME::ExceptionStruct es; + es.type = SALOME::BAD_PARAM; + es.text = "Invalid iteration"; + throw SALOME::SALOME_Exception(es); + } + + myCase->AddIteration(myIteration0->GetName()); + myIteration0->SetCaseName("Case_1"); } + else { + if (CORBA::is_nil(myIteration1)) { + SALOME::ExceptionStruct es; + es.type = SALOME::BAD_PARAM; + es.text = "Invalid iteration"; + throw SALOME::SALOME_Exception(es); + } - myCase->AddIteration(myIteration->GetName()); - myIteration->SetCaseName("Case_1"); + myCase->AddIteration(myIteration1->GetName()); + myIteration1->SetCaseName("Case_1"); + } } //============================================================================= @@ -1261,7 +1276,7 @@ SMESHHOMARD::HOMARD_Cas_ptr HOMARD_Gen_i::CreateCaseOnMesh (const char* MeshName SMESHHOMARD::extrema_var aSeq = new SMESHHOMARD::extrema(); if (LesExtremes.size() != 10) { return 0; } aSeq->length(10); - for (int i = 0; i < LesExtremes.size(); i++) + for (unsigned int i = 0; i < LesExtremes.size(); i++) aSeq[i] = LesExtremes[i]; myCase->SetBoundingBox(aSeq); // Les groupes @@ -1343,7 +1358,7 @@ SMESHHOMARD::HOMARD_Cas_ptr HOMARD_Gen_i::CreateCase(const char* MeshName, SMESHHOMARD::extrema_var aSeq = new SMESHHOMARD::extrema(); if (LesExtremes.size() != 10) { return 0; } aSeq->length(10); - for (int i = 0; i < LesExtremes.size(); i++) + for (unsigned int i = 0; i < LesExtremes.size(); i++) aSeq[i] = LesExtremes[i]; myCase->SetBoundingBox(aSeq); // Les groupes @@ -1857,11 +1872,12 @@ void HOMARD_Gen_i::CleanCase() } // Delete all boundaries - std::map::const_iterator it_boundary; - for (it_boundary = _mesBoundarys.begin(); - it_boundary != _mesBoundarys.end(); ++it_boundary) { - DeleteBoundary((*it_boundary).first.c_str()); - } + //std::map::const_iterator it_boundary; + //for (it_boundary = _mesBoundarys.begin(); + // it_boundary != _mesBoundarys.end(); ++it_boundary) { + // DeleteBoundary((*it_boundary).first.c_str()); + //} + _mesBoundarys.clear(); // Delete iteration DeleteIteration(1); @@ -1877,6 +1893,7 @@ void HOMARD_Gen_i::CleanCase() if (!_CaseOnMedFile && !_TmpMeshFile.empty()) { SMESH_File aFile (_TmpMeshFile, false); if (aFile.exists()) aFile.remove(); + _TmpMeshFile = ""; } _SmeshMesh = SMESH::SMESH_Mesh::_nil(); } @@ -1904,14 +1921,7 @@ CORBA::Long HOMARD_Gen_i::ComputeAdap(SMESHHOMARD::HOMARD_Cas_var myCase, ASSERT(!CORBA::is_nil(myHypothesis)); // B. L'iteration parent - //const char* nomIterationParent = myIteration->GetIterParentName(); - SMESHHOMARD::HOMARD_Iteration_var myIterationParent = myIteration0; - ASSERT(!CORBA::is_nil(myIterationParent)); - // Si l'iteration parent n'est pas calculee, on le fait (recursivite amont) - //if (myIterationParent->GetState() == 1) { - // int codret = Compute(nomIterationParent); - // if (codret != 0) VERIFICATION("Pb au calcul de l'iteration precedente" == 0); - //} + ASSERT(!CORBA::is_nil(myIteration0)); // C. Le sous-répertoire de l'iteration precedente char* DirComputePa = ComputeDirPaManagement(myCase, myIteration); @@ -1923,9 +1933,9 @@ CORBA::Long HOMARD_Gen_i::ComputeAdap(SMESHHOMARD::HOMARD_Cas_var myCase, MESSAGE (". ConfType = " << ConfType); // D.3. Le maillage de depart - const char* NomMeshParent = myIterationParent->GetMeshName(); + const char* NomMeshParent = myIteration0->GetMeshName(); MESSAGE (". NomMeshParent = " << NomMeshParent); - const char* MeshFileParent = myIterationParent->GetMeshFile(); + const char* MeshFileParent = myIteration0->GetMeshFile(); MESSAGE (". MeshFileParent = " << MeshFileParent); // D.4. Le maillage associe a l'iteration @@ -2072,10 +2082,28 @@ CORBA::Long HOMARD_Gen_i::ComputeCAO(SMESHHOMARD::HOMARD_Cas_var myCase, // C. Lancement des projections MESSAGE (". Lancement des projections"); - // TODO? - SMESHHOMARDImpl::FrontTrack* myFrontTrack = new SMESHHOMARDImpl::FrontTrack(); - myFrontTrack->track(theInputMedFile, theOutputMedFile, - theInputNodeFiles, theXaoFileName, theIsParallel); + + assert(Py_IsInitialized()); + PyGILState_STATE gstate; + gstate = PyGILState_Ensure(); + PyRun_SimpleString("from FrontTrack import FrontTrack"); + // FrontTrack().track( fic_med_brut, fic_med_new, l_fr, xao_file ) + std::string pyCommand ("FrontTrack().track( \""); + pyCommand += theInputMedFile + "\", \"" + theOutputMedFile + "\", ["; + for (int i = 0; i < icpt; i++) { + if (i > 0) pyCommand += ", "; + pyCommand += "\""; + pyCommand += theInputNodeFiles[i]; + pyCommand += "\""; + } + pyCommand += "], \"" + theXaoFileName + "\", False )"; + MESSAGE (". Lancement des projections: pyCommand = " << pyCommand); + PyRun_SimpleString(pyCommand.c_str()); + PyGILState_Release(gstate); + + //SMESHHOMARDImpl::FrontTrack* myFrontTrack = new SMESHHOMARDImpl::FrontTrack(); + //myFrontTrack->track(theInputMedFile, theOutputMedFile, + // theInputNodeFiles, theXaoFileName, theIsParallel); // D. Transfert des coordonnées modifiées dans le fichier historique de HOMARD // On lance une exécution spéciale de HOMARD en attendant @@ -2367,7 +2395,8 @@ char* HOMARD_Gen_i::ComputeDirManagement(SMESHHOMARD::HOMARD_Cas_var myCase, //============================================================================= // Calcul d'une iteration : gestion du répertoire de calcul de l'iteration parent //============================================================================= -char* HOMARD_Gen_i::ComputeDirPaManagement(SMESHHOMARD::HOMARD_Cas_var myCase, SMESHHOMARD::HOMARD_Iteration_var myIteration) +char* HOMARD_Gen_i::ComputeDirPaManagement(SMESHHOMARD::HOMARD_Cas_var myCase, + SMESHHOMARD::HOMARD_Iteration_var myIteration) { MESSAGE ("ComputeDirPaManagement : répertoires pour le calcul"); // Le répertoire du cas @@ -2376,8 +2405,7 @@ char* HOMARD_Gen_i::ComputeDirPaManagement(SMESHHOMARD::HOMARD_Cas_var myCase, S // Le sous-répertoire de l'iteration precedente - SMESHHOMARD::HOMARD_Iteration_var myIterationParent = myIteration0; - const char* nomDirItPa = myIterationParent->GetDirNameLoc(); + const char* nomDirItPa = myIteration0->GetDirNameLoc(); std::stringstream DirComputePa; DirComputePa << nomDirCase << "/" << nomDirItPa; MESSAGE(". nomDirItPa = " << nomDirItPa); @@ -2671,6 +2699,8 @@ void HOMARD_Gen_i::PythonDump() for (it_boundary = _mesBoundarys.begin(); it_boundary != _mesBoundarys.end(); ++it_boundary) { SMESHHOMARD::HOMARD_Boundary_var maBoundary = (*it_boundary).second; + //MESSAGE ("PythonDump of boundary " << (*it_boundary).first << + // " : " << maBoundary->GetDumpPython()); pd << maBoundary->GetDumpPython(); }