Salome HOME
Copyright update: 2016
[modules/smesh.git] / src / StdMeshers / StdMeshers_QuadFromMedialAxis_1D2D.cxx
index 680279e75c226a7341b0b1a2a854d059848eac09..46d661191d4799816baa7972924315dbde517156 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2015  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2016  CEA/DEN, EDF R&D, OPEN CASCADE
 //
 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
@@ -61,6 +61,8 @@
 #include <list>
 #include <vector>
 
+using namespace std;
+
 //================================================================================
 /*!
  * \brief 1D algo
@@ -551,8 +553,8 @@ namespace
                                allEdges, theShortEdges[ nbBranchPoints > 0 ] ))
         return false;
 
-      for ( size_t iS = 0; iS < theShortEdges[ nbBranchPoints ].size(); ++iS )
-        shortMap.Add( theShortEdges[ nbBranchPoints ][ iS ]);
+      for ( size_t iS = 0; iS < theShortEdges[ nbBranchPoints > 0 ].size(); ++iS )
+        shortMap.Add( theShortEdges[ nbBranchPoints > 0 ][ iS ]);
 
       ++nbBranchPoints;
     }
@@ -949,7 +951,7 @@ namespace
   {
     const SMDS_MeshNode* _node;
     double               _u;
-    int                  _edgeInd; // index in theSinuEdges vector
+    size_t               _edgeInd; // index in theSinuEdges vector
 
     NodePoint(): _node(0), _u(0), _edgeInd(-1) {}
     NodePoint(const SMDS_MeshNode* n, double u, size_t iEdge ): _node(n), _u(u), _edgeInd(iEdge) {}
@@ -1039,7 +1041,7 @@ namespace
     const vector< Handle(Geom_Curve) >& theCurves = theSinuFace._sinuCurves;
 
     double uMA;
-    SMESH_MAT2d::BoundaryPoint bp[2];
+    SMESH_MAT2d::BoundaryPoint bp[2]; // 2 sinuous sides
     const SMESH_MAT2d::Branch& branch = *theMA.getBranch(0);
     {
       // add to thePointsOnE NodePoint's of ends of theSinuEdges
@@ -1138,7 +1140,7 @@ namespace
       // projection is set to the BoundaryPoint of this projection
 
       // evaluate distance to neighbor projections
-      const double rShort = 0.2;
+      const double rShort = 0.33;
       bool isShortPrev[2], isShortNext[2], isPrevCloser[2];
       TMAPar2NPoints::iterator u2NPPrev = u2NP, u2NPNext = u2NP;
       --u2NPPrev; ++u2NPNext;
@@ -1251,8 +1253,8 @@ namespace
     const vector<TopoDS_Edge>&    theSinuEdges = theSinuFace._sinuEdges;
     const vector< Handle(Geom_Curve) >& curves = theSinuFace._sinuCurves;
 
-    SMESH_MAT2d::BoundaryPoint bp[2];
-    const SMESH_MAT2d::Branch& branch = *theMA.getBranch(0);
+    //SMESH_MAT2d::BoundaryPoint bp[2];
+    //const SMESH_MAT2d::Branch& branch = *theMA.getBranch(0);
 
     typedef TMAPar2NPoints::iterator TIterator;
 
@@ -1292,7 +1294,7 @@ namespace
         {
           // find an existing node on VERTEX among sameU2NP and get underlying EDGEs
           const SMDS_MeshNode* existingNode = 0;
-          set< int > edgeInds;
+          set< size_t > edgeInds;
           NodePoint* np;
           for ( size_t i = 0; i < sameU2NP.size(); ++i )
           {
@@ -1306,10 +1308,10 @@ namespace
 
           TIterator u2NPprev = sameU2NP.front();
           TIterator u2NPnext = sameU2NP.back() ;
-          if ( u2NPprev->first > 0. ) --u2NPprev;
-          if ( u2NPnext->first < 1. ) ++u2NPprev;
+          if ( u2NPprev->first < 0. ) ++u2NPprev;
+          if ( u2NPnext->first > 1. ) --u2NPnext;
 
-          set< int >::iterator edgeID = edgeInds.begin();
+          set< size_t >::iterator edgeID = edgeInds.begin();
           for ( ; edgeID != edgeInds.end(); ++edgeID )
           {
             // get U range on iEdge within which the equal points will be distributed
@@ -1322,8 +1324,8 @@ namespace
 
             if ( u0 == u1 )
             {
-              if ( np->_node ) --u2NPprev;
-              else             ++u2NPnext;
+              if ( u2NPprev != thePointsOnE.begin() ) --u2NPprev;
+              if ( u2NPnext != --thePointsOnE.end() ) ++u2NPnext;
               np = &get( u2NPprev->second, iSide );
               u0 = getUOnEdgeByPoint( *edgeID, np, theSinuFace );
               np = &get( u2NPnext->second, iSide );
@@ -1331,7 +1333,7 @@ namespace
             }
 
             // distribute points and create nodes
-            double du = ( u1 - u0 ) / ( sameU2NP.size() + !existingNode );
+            double du = ( u1 - u0 ) / ( sameU2NP.size() + 1 /*!existingNode*/ );
             double u  = u0 + du;
             for ( size_t i = 0; i < sameU2NP.size(); ++i )
             {
@@ -1503,8 +1505,8 @@ namespace
         theFace._quad->side[ 1 ] = StdMeshers_FaceSide::New( uvsNew );
       }
 
-      if ( theFace._quad->side[ 1 ].NbPoints() !=
-           theFace._quad->side[ 3 ].NbPoints())
+      if ( theFace._quad->side[ 1 ].GetUVPtStruct().empty() ||
+           theFace._quad->side[ 3 ].GetUVPtStruct().empty() )
         return false;
 
     } // if ( theFace.IsRing() )
@@ -1553,6 +1555,7 @@ namespace
     vector< int >                edgeIDs   ( theSinuEdges.size() ); // IDs in the main shape
     vector< bool >               isComputed( theSinuEdges.size() );
     curves.resize( theSinuEdges.size(), 0 );
+    bool                         allComputed = true;
     for ( size_t i = 0; i < theSinuEdges.size(); ++i )
     {
       curves[i] = BRep_Tool::Curve( theSinuEdges[i], f,l );
@@ -1561,6 +1564,8 @@ namespace
       SMESH_subMesh* sm = mesh->GetSubMesh( theSinuEdges[i] );
       edgeIDs   [i] = sm->GetId();
       isComputed[i] = ( !sm->IsEmpty() );
+      if ( !isComputed[i] )
+        allComputed = false;
     }
 
     const SMESH_MAT2d::Branch& branch = *theMA.getBranch(0);
@@ -1568,7 +1573,9 @@ namespace
 
     vector< std::size_t > edgeIDs1, edgeIDs2; // indices in theSinuEdges
     vector< SMESH_MAT2d::BranchPoint > divPoints;
-    branch.getOppositeGeomEdges( edgeIDs1, edgeIDs2, divPoints );
+    if ( !allComputed )
+      branch.getOppositeGeomEdges( edgeIDs1, edgeIDs2, divPoints );
+
     for ( size_t i = 0; i < edgeIDs1.size(); ++i )
       if ( isComputed[ edgeIDs1[i]] &&
            isComputed[ edgeIDs2[i]] )
@@ -1587,15 +1594,20 @@ namespace
           return false;
       }
 
-    // map param on MA to parameters of nodes on a pair of theSinuEdges
+    // map (param on MA) to (parameters of nodes on a pair of theSinuEdges)
     TMAPar2NPoints pointsOnE;
     vector<double> maParams;
+    set<int>       projectedEdges; // treated EDGEs which 'isComputed'
 
     // compute params of nodes on EDGEs by projecting division points from MA
 
     for ( size_t iEdgePair = 0; iEdgePair < edgeIDs1.size(); ++iEdgePair )
       // loop on pairs of opposite EDGEs
     {
+      if ( projectedEdges.count( edgeIDs1[ iEdgePair ]) ||
+           projectedEdges.count( edgeIDs2[ iEdgePair ]) )
+        continue;
+
       // --------------------------------------------------------------------------------
       if ( isComputed[ edgeIDs1[ iEdgePair ]] !=                    // one EDGE is meshed
            isComputed[ edgeIDs2[ iEdgePair ]])
@@ -1610,6 +1622,8 @@ namespace
         if ( !SMESH_Algo::GetSortedNodesOnEdge( meshDS, theSinuEdges[ iEdgeComputed ], /*skipMedium=*/true, nodeParams ))
           return false;
 
+        projectedEdges.insert( iEdgeComputed );
+
         SMESH_MAT2d::BoundaryPoint& bndPnt = bp[ 1-iSideComputed ];
         SMESH_MAT2d::BranchPoint brp;
         NodePoint npN, npB; // NodePoint's initialized by node and BoundaryPoint
@@ -1618,10 +1632,10 @@ namespace
 
         double maParam1st, maParamLast, maParam;
         if ( !theMA.getBoundary().getBranchPoint( iEdgeComputed, nodeParams.begin()->first, brp ))
-            return false;
+          return false;
         branch.getParameter( brp, maParam1st );
         if ( !theMA.getBoundary().getBranchPoint( iEdgeComputed, nodeParams.rbegin()->first, brp ))
-            return false;
+          return false;
         branch.getParameter( brp, maParamLast );
 
         map< double, const SMDS_MeshNode* >::iterator u2n = nodeParams.begin(), u2nEnd = nodeParams.end();
@@ -1641,28 +1655,6 @@ namespace
           npB = NodePoint( bndPnt );
           pos = pointsOnE.insert( hint, make_pair( maParam, make_pair( np0, np1 )));
         }
-
-        // move iEdgePair forward;
-        // find divPoints most close to max MA param
-        if ( edgeIDs1.size() > 1 )
-        {
-          maParamLast = pointsOnE.rbegin()->first;
-          int iClosest;
-          double minDist = 1.;
-          for ( ; iEdgePair < edgeIDs1.size()-1; ++iEdgePair )
-          {
-            branch.getParameter( divPoints[iEdgePair], maParam );
-            double d = Abs( maParamLast - maParam );
-            if ( d < minDist )
-              minDist = d, iClosest = iEdgePair;
-            else
-              break;
-          }
-          if ( Abs( maParamLast - 1. ) < minDist )
-            break; // the last pair treated
-          else
-            iEdgePair = iClosest;
-        }
       }
       // --------------------------------------------------------------------------------
       else if ( !isComputed[ edgeIDs1[ iEdgePair ]] &&         // none of EDGEs is meshed
@@ -1854,12 +1846,12 @@ namespace
     const double dksi = 0.5, deta = 0.5;
     const double  dksi2 = dksi*dksi, deta2 = deta*deta;
     double err = 0., g11, g22, g12;
-    int nbErr = 0;
+    //int nbErr = 0;
 
     FaceQuadStruct& q = *quad;
     UVPtStruct pNew;
 
-    double refArea = area( q.UVPt(0,0), q.UVPt(1,0), q.UVPt(1,1) );
+    //double refArea = area( q.UVPt(0,0), q.UVPt(1,0), q.UVPt(1,1) );
 
     for ( int iLoop = 0; iLoop < nbLoops; ++iLoop )
     {