Salome HOME
87af665b6af086270441b6ae5b2c3345bf2e1cdc
[modules/smesh.git] / src / StdMeshers / StdMeshers_FaceSide.cxx
1 // Copyright (C) 2007-2011  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
5 //
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License.
10 //
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14 // Lesser General Public License for more details.
15 //
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
19 //
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 //
22
23 // File      : StdMeshers_FaceSide.hxx
24 // Created   : Wed Jan 31 18:41:25 2007
25 // Author    : Edward AGAPOV (eap)
26 // Module    : SMESH
27 //
28 #include "StdMeshers_FaceSide.hxx"
29
30 #include "SMDS_EdgePosition.hxx"
31 #include "SMDS_MeshNode.hxx"
32 #include "SMESHDS_Mesh.hxx"
33 #include "SMESHDS_SubMesh.hxx"
34 #include "SMESH_Algo.hxx"
35 #include "SMESH_Mesh.hxx"
36 #include "SMESH_MesherHelper.hxx"
37 #include "SMESH_ComputeError.hxx"
38 #include "SMESH_Block.hxx"
39
40 #include <Adaptor2d_Curve2d.hxx>
41 #include <BRepAdaptor_CompCurve.hxx>
42 #include <BRep_Builder.hxx>
43 #include <BRep_Tool.hxx>
44 #include <GCPnts_AbscissaPoint.hxx>
45 #include <Geom2dAdaptor_Curve.hxx>
46 #include <TopExp.hxx>
47 #include <TopExp_Explorer.hxx>
48 #include <TopoDS.hxx>
49 #include <TopoDS_Face.hxx>
50 #include <TopoDS_Vertex.hxx>
51 #include <TopoDS_Wire.hxx>
52
53 #include <map>
54
55 #include "utilities.h"
56
57 //================================================================================
58 /*!
59  * \brief Constructor of a side of one edge
60   * \param theFace - the face
61   * \param theEdge - the edge
62  */
63 //================================================================================
64
65 StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
66                                          const TopoDS_Edge& theEdge,
67                                          SMESH_Mesh*        theMesh,
68                                          const bool         theIsForward,
69                                          const bool         theIgnoreMediumNodes)
70 {
71   list<TopoDS_Edge> edges(1,theEdge);
72   *this = StdMeshers_FaceSide( theFace, edges, theMesh, theIsForward, theIgnoreMediumNodes );
73 }
74
75 //================================================================================
76 /*!
77  * \brief Constructor of a side of several edges
78   * \param theFace - the face
79   * \param theEdge - the edge
80  */
81 //================================================================================
82
83 StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
84                                          list<TopoDS_Edge>& theEdges,
85                                          SMESH_Mesh*        theMesh,
86                                          const bool         theIsForward,
87                                          const bool         theIgnoreMediumNodes)
88 {
89   int nbEdges = theEdges.size();
90   myEdge.resize( nbEdges );
91   myEdgeID.resize( nbEdges );
92   myC2d.resize( nbEdges );
93   myC3dAdaptor.resize( nbEdges );
94   myFirst.resize( nbEdges );
95   myLast.resize( nbEdges );
96   myNormPar.resize( nbEdges );
97   myEdgeLength.resize( nbEdges );
98   myIsUniform.resize( nbEdges, true );
99   myLength = 0;
100   myNbPonits = myNbSegments = 0;
101   myMesh = theMesh;
102   myMissingVertexNodes = false;
103   myIgnoreMediumNodes = theIgnoreMediumNodes;
104   myDefaultPnt2d = gp_Pnt2d( 1e+100, 1e+100 );
105   if ( nbEdges == 0 ) return;
106
107   SMESHDS_Mesh* meshDS = theMesh->GetMeshDS();
108
109   int nbDegen = 0;
110   list<TopoDS_Edge>::iterator edge = theEdges.begin();
111   TopoDS_Iterator vExp;
112   for ( int index = 0; edge != theEdges.end(); ++index, ++edge )
113   {
114     int i = theIsForward ? index : nbEdges - index - 1;
115     myEdgeLength[i] = SMESH_Algo::EdgeLength( *edge );
116     if ( myEdgeLength[i] < DBL_MIN ) nbDegen++;
117     myLength += myEdgeLength[i];
118     myEdge[i] = *edge;
119     myEdgeID[i] = meshDS->ShapeToIndex( *edge );
120     if ( !theIsForward ) myEdge[i].Reverse();
121
122     if ( theFace.IsNull() )
123       BRep_Tool::Range( *edge, myFirst[i], myLast[i] );
124     else
125       myC2d[i] = BRep_Tool::CurveOnSurface( *edge, theFace, myFirst[i], myLast[i] );
126     if ( myEdge[i].Orientation() == TopAbs_REVERSED )
127       std::swap( myFirst[i], myLast[i] );
128
129     if ( SMESHDS_SubMesh* sm = meshDS->MeshElements( *edge )) {
130       int nbN = sm->NbNodes();
131       if ( theIgnoreMediumNodes ) {
132         SMDS_ElemIteratorPtr elemIt = sm->GetElements();
133         if ( elemIt->more() && elemIt->next()->IsQuadratic() )
134           nbN -= sm->NbElements();
135       }
136       myNbPonits += nbN;
137       myNbSegments += sm->NbElements();
138     }
139     // TopExp::FirstVertex() and TopExp::LastVertex() return NULL from INTERNAL edge
140     vExp.Initialize( *edge );
141     if ( vExp.Value().Orientation() == TopAbs_REVERSED ) vExp.Next();
142     if ( SMESH_Algo::VertexNode( TopoDS::Vertex( vExp.Value()), meshDS ))
143       myNbPonits += 1; // for the first end
144     else
145       myMissingVertexNodes = true;
146
147     // check if edge has non-uniform parametrization (issue 0020705)
148     if ( !myC2d[i].IsNull() && myEdgeLength[i] > DBL_MIN)
149     {
150       Geom2dAdaptor_Curve A2dC( myC2d[i],
151                                 std::min( myFirst[i], myLast[i] ),
152                                 std::max( myFirst[i], myLast[i] ));
153       double p2 = myFirst[i]+(myLast[i]-myFirst[i])/2., p4 = myFirst[i]+(myLast[i]-myFirst[i])/4.;
154       double d2 = GCPnts_AbscissaPoint::Length( A2dC, myFirst[i], p2 );
155       double d4 = GCPnts_AbscissaPoint::Length( A2dC, myFirst[i], p4 );
156       //cout<<"len = "<<len<<"  d2 = "<<d2<<"  fabs(2*d2/len-1.0) = "<<fabs(2*d2/len-1.0)<<endl;
157       myIsUniform[i] = !( fabs(2*d2/myEdgeLength[i]-1.0) > 0.01 || fabs(2*d4/d2-1.0) > 0.01 );
158       if ( !myIsUniform[i] )
159       {
160         double fp,lp;
161         TopLoc_Location L;
162         Handle(Geom_Curve) C3d = BRep_Tool::Curve(myEdge[i],L,fp,lp);
163         myC3dAdaptor[i].Load( C3d, fp,lp );
164       }
165     }
166   }
167   vExp.Initialize( theEdges.back() );
168   if ( vExp.Value().Orientation() != TopAbs_REVERSED ) vExp.Next();
169   if ( vExp.More() )
170   {
171     if ( SMESH_Algo::VertexNode( TopoDS::Vertex( vExp.Value()), meshDS ))
172       myNbPonits++; // for the last end
173     else
174       myMissingVertexNodes = true;
175   }
176   if ( nbEdges > 1 && myLength > DBL_MIN ) {
177     const double degenNormLen = 1.e-5;
178     double totLength = myLength;
179     if ( nbDegen )
180       totLength += myLength * degenNormLen * nbDegen;
181     double prevNormPar = 0;
182     for ( int i = 0; i < nbEdges; ++i ) {
183       if ( myEdgeLength[ i ] < DBL_MIN )
184         myEdgeLength[ i ] = myLength * degenNormLen;
185       myNormPar[ i ] = prevNormPar + myEdgeLength[i]/totLength;
186       prevNormPar = myNormPar[ i ];
187     }
188   }
189   myNormPar[nbEdges-1] = 1.;
190   //dump();
191 }
192
193 //================================================================================
194 /*!
195  * \brief Constructor of a side for vertex using data from other FaceSide
196   * \param theVertex - the vertex
197   * \param theSide - the side
198  */
199 //================================================================================
200
201 StdMeshers_FaceSide::StdMeshers_FaceSide(const SMDS_MeshNode* theNode,
202                                          const gp_Pnt2d thePnt2d,
203                                          const StdMeshers_FaceSide* theSide)
204 {
205   myC2d.resize(1);
206   myLength = 0;
207   myMesh = theSide->GetMesh();
208   myDefaultPnt2d = thePnt2d;
209
210   myPoints = theSide->GetUVPtStruct();
211   myNbPonits = myNbSegments = myPoints.size();
212   std::vector<uvPtStruct>::iterator it = myPoints.begin();
213   for(; it!=myPoints.end(); it++) {
214     (*it).u = thePnt2d.X();
215     (*it).v = thePnt2d.Y();
216     (*it).y = 0.0;
217     (*it).node = theNode;
218   }
219 }
220
221 //================================================================================
222 /*!
223  * \brief Return info on nodes on the side
224   * \retval UVPtStruct* - array of data structures
225  */
226 //================================================================================
227
228 const vector<UVPtStruct>& StdMeshers_FaceSide::GetUVPtStruct(bool   isXConst,
229                                                              double constValue) const
230 {
231   if ( myPoints.empty() ) {
232
233     if ( NbEdges() == 0 ) return myPoints;
234
235     SMESHDS_Mesh* meshDS = myMesh->GetMeshDS();
236     SMESH_MesherHelper helper(*myMesh);
237     bool paramOK;
238
239     // sort nodes of all edges putting them into a map
240
241     map< double, const SMDS_MeshNode*> u2node;
242     //int nbOnDegen = 0;
243     for ( int i = 0; i < myEdge.size(); ++i )
244     {
245       // Put 1st vertex node of a current edge
246       TopoDS_Vertex VV[2]; // TopExp::FirstVertex() returns NULL for INTERNAL edge
247       VV[0] = SMESH_MesherHelper::IthVertex( 0, myEdge[i]);
248       VV[1] = SMESH_MesherHelper::IthVertex( 1, myEdge[i]);
249       const SMDS_MeshNode* node = SMESH_Algo::VertexNode( VV[0], meshDS );
250       double prevNormPar = ( i == 0 ? 0 : myNormPar[ i-1 ]); // normalized param
251       if ( node ) { // internal nodes may be missing
252         u2node.insert( make_pair( prevNormPar, node ));
253       }
254       else if ( i == 0 ) {
255         MESSAGE(" NO NODE on VERTEX" );
256         return myPoints;
257       }
258
259       // Put internal nodes
260       SMESHDS_SubMesh* sm = meshDS->MeshElements( myEdge[i] );
261       if ( !sm ) continue;
262       vector< pair< double, const SMDS_MeshNode*> > u2nodeVec;
263       u2nodeVec.reserve( sm->NbNodes() );
264       SMDS_NodeIteratorPtr nItr = sm->GetNodes();
265       double paramSize = myLast[i] - myFirst[i];
266       double r = myNormPar[i] - prevNormPar;
267       if ( !myIsUniform[i] )
268         while ( nItr->more() )
269         {
270           const SMDS_MeshNode* node = nItr->next();
271           if ( myIgnoreMediumNodes && SMESH_MeshEditor::IsMedium( node, SMDSAbs_Edge ))
272             continue;
273           double u = helper.GetNodeU( myEdge[i], node, 0, &paramOK );
274           double aLenU = GCPnts_AbscissaPoint::Length
275             ( const_cast<GeomAdaptor_Curve&>( myC3dAdaptor[i]), myFirst[i], u );
276           if ( myEdgeLength[i] < aLenU ) // nonregression test "3D_mesh_NETGEN/G6"
277           {
278             u2nodeVec.clear();
279             break;
280           }
281           double normPar = prevNormPar + r*aLenU/myEdgeLength[i];
282           u2nodeVec.push_back( make_pair( normPar, node ));
283         }
284       nItr = sm->GetNodes();
285       if ( u2nodeVec.empty() )
286         while ( nItr->more() )
287         {
288           const SMDS_MeshNode* node = nItr->next();
289           if ( myIgnoreMediumNodes && SMESH_MeshEditor::IsMedium( node, SMDSAbs_Edge ))
290             continue;
291           double u = helper.GetNodeU( myEdge[i], node, 0, &paramOK );
292
293           // paramSize is signed so orientation is taken into account
294           double normPar = prevNormPar + r * ( u - myFirst[i] ) / paramSize;
295           u2nodeVec.push_back( make_pair( normPar, node ));
296         }
297       u2node.insert( u2nodeVec.begin(), u2nodeVec.end() );
298
299       // Put 2nd vertex node for a last edge
300       if ( i+1 == myEdge.size() ) {
301         node = SMESH_Algo::VertexNode( VV[1], meshDS );
302         if ( !node ) {
303           MESSAGE(" NO NODE on VERTEX" );
304           return myPoints;
305         }
306         u2node.insert( make_pair( 1., node ));
307       }
308     }
309     if ( u2node.size() != myNbPonits ) {
310       MESSAGE("Wrong node parameters on edges, u2node.size():"
311               <<u2node.size()<<" !=  myNbPonits:"<<myNbPonits);
312       return myPoints;
313     }
314
315     // fill array of UVPtStruct
316
317     vector<uvPtStruct>* points = const_cast<vector<uvPtStruct>*>( &myPoints );
318     points->resize( myNbPonits );
319
320     int EdgeIndex = 0;
321     double prevNormPar = 0, paramSize = myNormPar[ EdgeIndex ];
322     map< double, const SMDS_MeshNode*>::iterator u_node = u2node.begin();
323     for (int i = 0 ; u_node != u2node.end(); ++u_node, ++i ) {
324       UVPtStruct & uvPt = (*points)[i];
325       uvPt.node = u_node->second;
326       uvPt.x = uvPt.y = uvPt.normParam = u_node->first;
327       if ( isXConst ) uvPt.x = constValue;
328       else            uvPt.y = constValue;
329       const SMDS_EdgePosition* epos =
330         dynamic_cast<const SMDS_EdgePosition*>(uvPt.node->GetPosition());
331       if (( myNormPar[ EdgeIndex ] < uvPt.normParam ) ||
332           ( epos && uvPt.node->getshapeId() != myEdgeID[ EdgeIndex ])) // for myMissingVertexNodes
333       {
334         prevNormPar = myNormPar[ EdgeIndex ];
335         ++EdgeIndex;
336 #ifdef _DEBUG_
337         if ( EdgeIndex >= myEdge.size() ) {
338           dump("DEBUG");
339           MESSAGE ( "WRONg EdgeIndex " << 1+EdgeIndex
340                     << " myNormPar.size()="<<myNormPar.size()
341                     << " myNormPar["<< EdgeIndex<<"]="<< myNormPar[ EdgeIndex ]
342                     << " uvPt.normParam="<<uvPt.normParam );
343         }
344 #endif
345         paramSize = myNormPar[ EdgeIndex ] - prevNormPar;
346       }
347       if ( epos ) {
348         uvPt.param = epos->GetUParameter();
349       }
350       else {
351         double r = ( uvPt.normParam - prevNormPar )/ paramSize;
352 //         uvPt.param = myFirst[EdgeIndex] * ( 1 - r ) + myLast[EdgeIndex] * r;
353         uvPt.param = ( r > 0.5 ? myLast[EdgeIndex] : myFirst[EdgeIndex] );
354       }
355       if ( !myC2d[ EdgeIndex ].IsNull() ) {
356         gp_Pnt2d p = myC2d[ EdgeIndex ]->Value( uvPt.param );
357         uvPt.u = p.X();
358         uvPt.v = p.Y();
359       }
360       else {
361         uvPt.u = uvPt.v = 1e+100;
362       }
363     }
364   }
365   return myPoints;
366 }
367
368 //================================================================================
369 /*!
370  * \brief Falsificate info on nodes
371   * \param nbSeg - nb of segments on the side
372   * \retval UVPtStruct* - array of data structures
373  */
374 //================================================================================
375
376 const vector<UVPtStruct>& StdMeshers_FaceSide::SimulateUVPtStruct(int    nbSeg,
377                                                                   bool   isXConst,
378                                                                   double constValue) const
379 {
380   if ( myFalsePoints.empty() ) {
381
382     if ( NbEdges() == 0 ) return myFalsePoints;
383
384     vector<uvPtStruct>* points = const_cast<vector<uvPtStruct>*>( &myFalsePoints );
385     points->resize( nbSeg+1 );
386
387     int EdgeIndex = 0;
388     double prevNormPar = 0, paramSize = myNormPar[ EdgeIndex ];
389     for (int i = 0 ; i < myFalsePoints.size(); ++i ) {
390       double normPar = double(i) / double(nbSeg);
391       UVPtStruct & uvPt = (*points)[i];
392       uvPt.node = 0;
393       uvPt.x = uvPt.y = uvPt.param = uvPt.normParam = normPar;
394       if ( isXConst ) uvPt.x = constValue;
395       else            uvPt.y = constValue;
396       if ( myNormPar[ EdgeIndex ] < normPar ) {
397         prevNormPar = myNormPar[ EdgeIndex ];
398         ++EdgeIndex;
399         paramSize = myNormPar[ EdgeIndex ] - prevNormPar;
400       }
401       double r = ( normPar - prevNormPar )/ paramSize;
402       uvPt.param = myFirst[EdgeIndex] * ( 1 - r ) + myLast[EdgeIndex] * r;
403       if ( !myC2d[ EdgeIndex ].IsNull() ) {
404         gp_Pnt2d p = myC2d[ EdgeIndex ]->Value( uvPt.param );
405         uvPt.u = p.X();
406         uvPt.v = p.Y();
407       }
408       else {
409         uvPt.u = uvPt.v = 1e+100;
410       }
411     }
412   }
413   return myFalsePoints;
414 }
415 // gp_Pnt StdMeshers_FaceSide::Value(double U) const
416 // {
417 // }
418
419 //================================================================================
420 /*!
421  * \brief reverse order of vector elements
422   * \param vec - vector to reverse
423  */
424 //================================================================================
425
426 template <typename T > void reverse(vector<T> & vec)
427 {
428   for ( int f=0, r=vec.size()-1; f < r; ++f, --r )
429     std::swap( vec[f], vec[r] );
430 }
431
432 //================================================================================
433 /*!
434  * \brief Change orientation of side geometry
435  */
436 //================================================================================
437
438 void StdMeshers_FaceSide::Reverse()
439 {
440   int nbEdges = myEdge.size();
441   for ( int i = nbEdges-1; i >= 0; --i ) {
442     std::swap( myFirst[i], myLast[i] );
443     myEdge[i].Reverse();
444     if ( i > 0 ) // at the first loop 1. is overwritten
445       myNormPar[i] = 1 - myNormPar[i-1];
446   }
447   if ( nbEdges > 1 ) {
448     reverse( myEdge );
449     reverse( myEdgeID );
450     reverse( myC2d );
451     reverse( myC3dAdaptor );
452     reverse( myFirst );
453     reverse( myLast );
454     reverse( myNormPar );
455     reverse( myEdgeLength );
456     reverse( myIsUniform );
457   }
458   if ( nbEdges > 0 )
459   {
460     myNormPar[nbEdges-1]=1.;
461     myPoints.clear();
462     myFalsePoints.clear();
463   }
464 }
465
466 //================================================================================
467 /*!
468  * \brief Show side features
469  */
470 //================================================================================
471
472 void StdMeshers_FaceSide::dump(const char* msg) const
473 {
474 #ifdef _DEBUG_
475   if (msg) MESSAGE ( std::endl << msg );
476   MESSAGE_BEGIN ("NB EDGES: "<< myEdge.size() );
477   MESSAGE_ADD ( "nbPoints: "<<myNbPonits<<" vecSize: " << myPoints.size()<<" "<<myFalsePoints.size() );
478   for ( int i=0; i<myEdge.size(); ++i)
479   {
480     MESSAGE_ADD ( "\t"<<i+1 );
481     MESSAGE_ADD ( "\tEDGE: " );
482     if (myEdge[i].IsNull()) {
483       MESSAGE_ADD ( "NULL" );
484     }
485     else {
486       TopAbs::Print(myEdge[i].Orientation(),cout)<<" "<<myEdge[i].TShape().operator->()<<endl;
487       MESSAGE_ADD ( "\tV1: " << TopExp::FirstVertex( myEdge[i], 1).TShape().operator->()
488                  << "  V2: " << TopExp::LastVertex( myEdge[i], 1).TShape().operator->() );
489     }
490     MESSAGE_ADD ( "\tC2d: ");
491     
492     if (myC2d[i].IsNull()) {
493       MESSAGE_ADD ( "NULL" );
494     }
495     else {
496       MESSAGE_ADD ( myC2d[i].operator->() );
497     }
498       
499     MESSAGE_ADD ( "\tF: "<<myFirst[i]<< " L: "<< myLast[i] );
500     MESSAGE_END ( "\tnormPar: "<<myNormPar[i]<<endl );
501   }
502 #endif
503 }
504
505 //================================================================================
506 /*!
507  * \brief Creates a Adaptor2d_Curve2d to be used in SMESH_Block
508   * \retval Adaptor2d_Curve2d* - 
509  */
510 //================================================================================
511
512 struct Adaptor2dCurve2d : public Adaptor2d_Curve2d
513 {
514   const StdMeshers_FaceSide* mySide;
515   Adaptor2dCurve2d(const StdMeshers_FaceSide* faceSide):mySide(faceSide) {}
516   gp_Pnt2d Value(const Standard_Real U) const { return mySide->Value2d( U ); }
517   Standard_Real FirstParameter() const { return 0; }
518   Standard_Real LastParameter() const { return 1; }
519 };
520
521 Adaptor2d_Curve2d* StdMeshers_FaceSide::GetCurve2d() const
522 {
523   return new Adaptor2dCurve2d( this );
524 }
525
526 //================================================================================
527 /*!
528  * \brief Creates a fully functional Adaptor_Curve
529  */
530 //================================================================================
531
532 BRepAdaptor_CompCurve* StdMeshers_FaceSide::GetCurve3d() const
533 {
534   if ( myEdge.empty() )
535     return 0;
536
537   TopoDS_Wire aWire;
538   BRep_Builder aBuilder;
539   aBuilder.MakeWire(aWire);
540   for ( int i=0; i<myEdge.size(); ++i )
541     aBuilder.Add( aWire, myEdge[i] );
542
543   if ( myEdge.size() == 2 && FirstVertex().IsSame( LastVertex() ))
544     aWire.Closed(true); // issue 0021141
545
546   return new BRepAdaptor_CompCurve( aWire );
547 }
548
549
550 //================================================================================
551 /*!
552  * \brief Return 2D point by normalized parameter
553   * \param U - normalized parameter value
554   * \retval gp_Pnt2d - point
555  */
556 //================================================================================
557
558 gp_Pnt2d StdMeshers_FaceSide::Value2d(double U) const
559 {
560   if ( !myC2d[0].IsNull() ) {
561     int i = EdgeIndex( U );
562     double prevU = i ? myNormPar[ i-1 ] : 0;
563     double r = ( U - prevU )/ ( myNormPar[ i ] - prevU );
564
565     double par = myFirst[i] * ( 1 - r ) + myLast[i] * r;
566     
567     // check parametrization of curve
568     if( !myIsUniform[i] )
569     {
570       double aLen3dU = r * myEdgeLength[i] * ( myFirst[i]>myLast[i] ? -1. : 1.);
571       GCPnts_AbscissaPoint AbPnt
572         ( const_cast<GeomAdaptor_Curve&>( myC3dAdaptor[i]), aLen3dU, myFirst[i] );
573       if( AbPnt.IsDone() ) {
574         par = AbPnt.Parameter();
575       }
576     }
577     return myC2d[ i ]->Value(par);
578
579   }
580   return myDefaultPnt2d;
581 }
582
583 //================================================================================
584 /*!
585  * \brief Return wires of a face as StdMeshers_FaceSide's
586  */
587 //================================================================================
588
589 TSideVector StdMeshers_FaceSide::GetFaceWires(const TopoDS_Face& theFace,
590                                               SMESH_Mesh &       theMesh,
591                                               const bool         theIgnoreMediumNodes,
592                                               TError &           theError)
593 {
594   TopoDS_Vertex V1;
595   list< TopoDS_Edge > edges, internalEdges;
596   list< int > nbEdgesInWires;
597   int nbWires = SMESH_Block::GetOrderedEdges (theFace, V1, edges, nbEdgesInWires);
598
599   // split list of all edges into separate wires
600   TSideVector wires( nbWires );
601   list< int >::iterator nbE = nbEdgesInWires.begin();
602   list< TopoDS_Edge >::iterator from = edges.begin(), to = from;
603   for ( int iW = 0; iW < nbWires; ++iW, ++nbE )
604   {
605     std::advance( to, *nbE );
606     if ( *nbE == 0 ) // Issue 0020676
607     {
608       --nbWires;
609       --iW;
610       wires.resize( nbWires );
611       continue;
612     }
613     list< TopoDS_Edge > wireEdges( from, to );
614     // assure that there is a node on the first vertex
615     // as StdMeshers_FaceSide::GetUVPtStruct() requires
616     if ( wireEdges.front().Orientation() != TopAbs_INTERNAL ) // Issue 0020676
617     {
618       while ( !SMESH_Algo::VertexNode( TopExp::FirstVertex( wireEdges.front(), true),
619                                        theMesh.GetMeshDS()))
620       {
621         wireEdges.splice(wireEdges.end(), wireEdges,
622                          wireEdges.begin(), ++wireEdges.begin());
623         if ( from->IsSame( wireEdges.front() )) {
624           theError = TError
625             ( new SMESH_ComputeError(COMPERR_BAD_INPUT_MESH,"No nodes on vertices"));
626           return TSideVector(0);
627         }
628       }
629     }
630     else if ( *nbE > 1 ) // Issue 0020676 (Face_pb_netgen.brep) - several internal edges in a wire
631     {
632       internalEdges.splice( internalEdges.end(), wireEdges, ++wireEdges.begin(), wireEdges.end());
633     }
634
635     StdMeshers_FaceSide* wire = new StdMeshers_FaceSide( theFace, wireEdges, &theMesh,
636                                                          /*isForward=*/true, theIgnoreMediumNodes);
637     wires[ iW ] = StdMeshers_FaceSidePtr( wire );
638     from = to;
639   }
640   while ( !internalEdges.empty() )
641   {
642     StdMeshers_FaceSide* wire = new StdMeshers_FaceSide( theFace, internalEdges.back(), &theMesh,
643                                                          /*isForward=*/true, theIgnoreMediumNodes);
644     wires.push_back( StdMeshers_FaceSidePtr( wire ));
645     internalEdges.pop_back();
646   }
647   return wires;
648 }
649
650 //================================================================================
651 /*!
652  * \brief Return 1st vertex of the i-the edge
653  */
654 //================================================================================
655
656 TopoDS_Vertex StdMeshers_FaceSide::FirstVertex(int i) const
657 {
658   TopoDS_Vertex v;
659   if ( i < NbEdges() )
660   {
661     v = myEdge[i].Orientation() <= TopAbs_REVERSED ? // FORWARD || REVERSED
662         TopExp::FirstVertex( myEdge[i], 1 )        :
663         TopoDS::Vertex( TopoDS_Iterator( myEdge[i] ).Value() );
664   }
665   return v;
666 }
667
668 //================================================================================
669 /*!
670  * \brief Return last vertex of the i-the edge
671  */
672 //================================================================================
673
674 TopoDS_Vertex StdMeshers_FaceSide::LastVertex(int i) const
675 {
676   TopoDS_Vertex v;
677   if ( i < NbEdges() )
678   {
679     const TopoDS_Edge& edge = i<0 ? myEdge[ NbEdges() + i ] : myEdge[i];
680     if ( edge.Orientation() <= TopAbs_REVERSED ) // FORWARD || REVERSED
681       v = TopExp::LastVertex( edge, 1 );
682     else
683       for ( TopoDS_Iterator vIt( edge ); vIt.More(); vIt.Next() )
684         v = TopoDS::Vertex( vIt.Value() );
685   }
686   return v;
687 }
688