X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FSMESHUtils%2FSMESH_MAT2d.cxx;h=66ac2ff309a0487e65c931b1aa78ff2fe2299294;hp=f0f4d52a62966f469c77dd4e26312a0f40bd4b0d;hb=27c5228fcfe39a6405378da49ac700ee2adca6cc;hpb=fe7d1d57677486d8c546226dc2bf573fbfb6679d diff --git a/src/SMESHUtils/SMESH_MAT2d.cxx b/src/SMESHUtils/SMESH_MAT2d.cxx index f0f4d52a6..66ac2ff30 100644 --- a/src/SMESHUtils/SMESH_MAT2d.cxx +++ b/src/SMESHUtils/SMESH_MAT2d.cxx @@ -47,6 +47,7 @@ #include #ifdef _DEBUG_ +//#define _MYDEBUG_ #include "SMESH_File.hxx" #include "SMESH_Comment.hxx" #endif @@ -66,8 +67,8 @@ namespace struct InPoint { - int _a, _b; - double _param; + int _a, _b; // coordinates + double _param; // param on EDGE InPoint(int x, int y, double param) : _a(x), _b(y), _param(param) {} InPoint() : _a(0), _b(0), _param(0) {} @@ -76,6 +77,8 @@ namespace size_t index( const vector< InPoint >& inPoints ) const { return this - &inPoints[0]; } bool operator==( const InPoint& other ) const { return _a == other._a && _b == other._b; } + bool operator==( const TVDVertex* v ) const { return ( Abs( _a - v->x() ) < 1. && + Abs( _b - v->y() ) < 1. ); } }; // ------------------------------------------------------------------------------------- @@ -105,10 +108,12 @@ namespace // check if a TVDEdge begins at my end or ends at my start inline bool InSegment::isConnected( const TVDEdge* edge ) { - return ((Abs( edge->vertex0()->x() - _p1->_a ) < 1.&& - Abs( edge->vertex0()->y() - _p1->_b ) < 1. ) || - (Abs( edge->vertex1()->x() - _p0->_a ) < 1.&& - Abs( edge->vertex1()->y() - _p0->_b ) < 1. )); + return (( edge->vertex0() && edge->vertex1() ) + && + ((Abs( edge->vertex0()->x() - _p1->_a ) < 1.&& + Abs( edge->vertex0()->y() - _p1->_b ) < 1. ) || + (Abs( edge->vertex1()->x() - _p0->_a ) < 1.&& + Abs( edge->vertex1()->y() - _p0->_b ) < 1. ))); } // check if a MA TVDEdge is outside of a domain @@ -147,7 +152,7 @@ namespace // } // ------------------------------------------------------------------------------------- -#ifdef _DEBUG_ +#ifdef _MYDEBUG_ // writes segments into a txt file readable by voronoi_visualizer void inSegmentsToFile( vector< InSegment>& inSegments) { @@ -155,6 +160,7 @@ namespace return; const char* fileName = "/misc/dn25/salome/eap/salome/misc/Code/C++/MAdebug.txt"; SMESH_File file(fileName, false ); + file.remove(); file.openForWriting(); SMESH_Comment text; text << "0\n"; // nb points @@ -180,7 +186,7 @@ namespace if ( !edge->vertex1() ) cout << ") -> ( INF, INF"; else - cout << ") -> (" << edge->vertex1()->x() << ", " << edge->vertex1()->y(); + cout << ") -> ( " << edge->vertex1()->x() << ", " << edge->vertex1()->y(); cout << ")\t cell=" << edge->cell() << " iBnd=" << edge->color() << " twin=" << edge->twin() @@ -253,6 +259,7 @@ namespace boost { namespace { const int theNoBrachID = 0; // std::numeric_limits::max(); + double theScale[2]; // scale used in bndSegsToMesh() // ------------------------------------------------------------------------------------- /*! @@ -298,7 +305,7 @@ namespace if ( _edge ) // pass branch to an opposite BndSeg { size_t oppSegIndex = SMESH_MAT2d::Branch::getBndSegment( _edge->twin() ); - if ( oppSegIndex < bndSegs.size() /*&& bndSegs[ oppSegIndex ]._branchID == theNoBrachID*/ ) + if ( oppSegIndex < bndSegs.size() && bndSegs[ oppSegIndex ]._branchID == theNoBrachID ) bndSegs[ oppSegIndex ]._branchID = -branchID; } } @@ -349,7 +356,64 @@ namespace //================================================================================ /*! - * \brief Computes length of of TVDEdge + * \brief debug: to visually check found MA edges + */ + //================================================================================ + + void bndSegsToMesh( const vector< BndSeg >& bndSegs ) + { +#ifdef _MYDEBUG_ + if ( !getenv("bndSegsToMesh")) return; + map< const TVDVertex *, int > v2Node; + map< const TVDVertex *, int >::iterator v2n; + set< const TVDEdge* > addedEdges; + + const char* fileName = "/misc/dn25/salome/eap/salome/misc/Code/C++/MAedges.py"; + SMESH_File file(fileName, false ); + file.remove(); + file.openForWriting(); + SMESH_Comment text; + text << "import salome, SMESH\n"; + text << "salome.salome_init()\n"; + text << "from salome.smesh import smeshBuilder\n"; + text << "smesh = smeshBuilder.New(salome.myStudy)\n"; + text << "m=smesh.Mesh()\n"; + for ( size_t i = 0; i < bndSegs.size(); ++i ) + { + if ( !bndSegs[i]._edge ) + text << "# " << i << " NULL edge"; + else if ( !bndSegs[i]._edge->vertex0() || + !bndSegs[i]._edge->vertex1() ) + text << "# " << i << " INFINITE edge"; + else if ( addedEdges.insert( bndSegs[i]._edge ).second && + addedEdges.insert( bndSegs[i]._edge->twin() ).second ) + { + v2n = v2Node.insert( make_pair( bndSegs[i]._edge->vertex0(), v2Node.size() + 1 )).first; + int n0 = v2n->second; + if ( n0 == v2Node.size() ) + text << "n" << n0 << " = m.AddNode( " + << bndSegs[i]._edge->vertex0()->x() / theScale[0] << ", " + << bndSegs[i]._edge->vertex0()->y() / theScale[1] << ", 0 )\n"; + + v2n = v2Node.insert( make_pair( bndSegs[i]._edge->vertex1(), v2Node.size() + 1 )).first; + int n1 = v2n->second; + if ( n1 == v2Node.size() ) + text << "n" << n1 << " = m.AddNode( " + << bndSegs[i]._edge->vertex1()->x() / theScale[0] << ", " + << bndSegs[i]._edge->vertex1()->y() / theScale[1] << ", 0 )\n"; + + text << "e" << i << " = m.AddEdge([ n" << n0 << ", n" << n1 << " ])\n"; + } + } + text << "\n"; + file.write( text.c_str(), text.size() ); + cout << "Write " << fileName << endl; +#endif + } + + //================================================================================ + /*! + * \brief Computes length of a TVDEdge */ //================================================================================ @@ -593,6 +657,10 @@ namespace } } } + // debug + theScale[0] = scale[0]; + theScale[1] = scale[1]; + return true; } @@ -729,9 +797,12 @@ namespace continue; inPntChecked[ pInd ] = true; - const TVDEdge* edge = // a TVDEdge passing through an end of inSeg - is2nd ? maEdges.front()->prev() : maEdges.back()->next(); - while ( true ) + const TVDEdge* maE = is2nd ? maEdges.front() : maEdges.back(); + if ( inPnt == ( is2nd ? maE->vertex0() : maE->vertex1() )) + continue; + const TVDEdge* edge = // a secondary TVDEdge connecting inPnt and maE + is2nd ? maE->prev() : maE->next(); + while ( inSeg.isConnected( edge )) { if ( edge->is_primary() ) break; // this should not happen const TVDEdge* edge2 = edge->twin(); // we are in a neighbor cell, add MA edges to inPnt @@ -773,7 +844,7 @@ namespace inPoints[0]._edges.clear(); } - // Divide InSegment's into BndSeg's + // Divide InSegment's into BndSeg's (so that each BndSeg corresponds to one MA edge) vector< BndSeg > bndSegs; bndSegs.reserve( inSegments.size() * 3 ); @@ -791,25 +862,26 @@ namespace inPntChecked[ ip0 ] = false; // segments of InSegment's - size_t nbMaEdges = inSeg._edges.size(); + const size_t nbMaEdges = inSeg._edges.size(); switch ( nbMaEdges ) { case 0: // "around" circle center bndSegs.push_back( BndSeg( &inSeg, 0, inSeg._p1->_param )); break; case 1: bndSegs.push_back( BndSeg( &inSeg, inSeg._edges.back(), inSeg._p1->_param )); break; default: - vector< double > len; - len.push_back(0); - for ( e = inSeg._edges.rbegin(); e != inSeg._edges.rend(); ++e ) - len.push_back( len.back() + length( *e )); - + gp_XY inSegDir( inSeg._p1->_a - inSeg._p0->_a, + inSeg._p1->_b - inSeg._p0->_b ); + const double inSegLen2 = inSegDir.SquareModulus(); e = inSeg._edges.rbegin(); - for ( size_t l = 1; l < len.size(); ++e, ++l ) + for ( size_t iE = 1; iE < nbMaEdges; ++e, ++iE ) { - double dl = len[l] / len.back(); - double u = dl * inSeg._p1->_param + ( 1. - dl ) * inSeg._p0->_param; + gp_XY toMA( (*e)->vertex0()->x() - inSeg._p0->_a, + (*e)->vertex0()->y() - inSeg._p0->_b ); + double r = toMA * inSegDir / inSegLen2; + double u = r * inSeg._p1->_param + ( 1. - r ) * inSeg._p0->_param; bndSegs.push_back( BndSeg( &inSeg, *e, u )); } + bndSegs.push_back( BndSeg( &inSeg, *e, inSeg._p1->_param )); } // segments around 2nd concave point size_t ip1 = inSeg._p1->index( inPoints ); @@ -824,6 +896,8 @@ namespace for ( size_t i = 0; i < bndSegs.size(); ++i ) bndSegs[i].setIndexToEdge( i ); + bndSegsToMesh( bndSegs ); // debug: visually check found MA edges + // Find TVDEdge's of Branches and associate them with bndSegs @@ -838,7 +912,7 @@ namespace size_t i1st = 0; while ( i1st < bndSegs.size() && !bndSegs[i1st].hasOppositeEdge( noEdgeID )) ++i1st; - bndSegs[i1st].setBranch( branchID, bndSegs ); // set to the i-th and the opposite bndSeg + bndSegs[i1st].setBranch( branchID, bndSegs ); // set to the i-th and to the opposite bndSeg branchEdges[ branchID ].push_back( bndSegs[i1st]._edge ); for ( size_t i = i1st+1; i < bndSegs.size(); ++i ) @@ -865,7 +939,7 @@ namespace endType.insert( make_pair( bndSegs[i]._edge->vertex1(), SMESH_MAT2d::BE_BRANCH_POINT )); } - bndSegs[i].setBranch( branchID, bndSegs ); // set to i-th and the opposite bndSeg + bndSegs[i].setBranch( branchID, bndSegs ); // set to i-th and to the opposite bndSeg if ( bndSegs[i].hasOppositeEdge( noEdgeID )) branchEdges[ bndSegs[i].branchID() ].push_back( bndSegs[i]._edge ); } @@ -1053,15 +1127,10 @@ namespace { // no MA edge, bndSeg corresponds to an end point of a branch if ( bndPoints._maEdges.empty() ) - { - // should not get here according to algo design??? edgeInd = 0; - } else - { edgeInd = branchEdges[ brID ].size(); - dInd = bndSegs[ i ]._branchID > 0 ? +1 : -1; - } + dInd = bndSegs[ i ]._branchID > 0 ? +1 : -1; } bndPoints._maEdges.push_back( make_pair( br, ( 1 + edgeInd ) * dInd )); @@ -1070,7 +1139,7 @@ namespace iSeg = iSegEnd; - } // loop on all bndSegs + } // loop on all bndSegs to construct Boundary // Initialize branches @@ -1187,7 +1256,7 @@ void SMESH_MAT2d::MedialAxis::getPoints( const Branch* branch, //================================================================================ /*! * \brief Returns a BranchPoint corresponding to a given point on a geom EDGE - * \param [in] iGeomEdge - index of geom EDGE within a vector passed at MA construction + * \param [in] iEdge - index of geom EDGE within a vector passed at MA construction * \param [in] u - parameter of the point on EDGE curve * \param [out] p - the found BranchPoint * \return bool - is OK @@ -1242,9 +1311,9 @@ bool SMESH_MAT2d::Boundary::getBranchPoint( const std::size_t iEdge, const std::pair< const Branch*, int >& maE = points._maEdges[ i ]; bool maReverse = ( maE.second < 0 ); - p._branch = maE.first; - p._iEdge = ( maReverse ? -maE.second : maE.second ) - 1; // countered from 1 to store sign - p._edgeParam = maReverse ? ( 1. - edgeParam ) : edgeParam; + p._branch = maE.first; + p._iEdge = ( maReverse ? -maE.second : maE.second ) - 1; // countered from 1 to store sign + p._edgeParam = ( maE.first && maReverse ) ? ( 1. - edgeParam ) : edgeParam; return true; } @@ -1507,6 +1576,12 @@ bool SMESH_MAT2d::Branch::getBoundaryPoints(const BranchPoint& p, bool SMESH_MAT2d::Branch::getParameter(const BranchPoint & p, double & u ) const { + if ( this != p._branch && p._branch ) + return p._branch->getParameter( p, u ); + + if ( isRemoved() ) + return _proxyPoint._branch->getParameter( _proxyPoint, u ); + if ( p._iEdge > _params.size()-1 ) return false; if ( p._iEdge == _params.size()-1 )