From 7f18f75436555a166cd36bb245804705d2b531ed Mon Sep 17 00:00:00 2001
From: eap <eap@opencascade.com>
Date: Thu, 22 Apr 2010 13:32:59 +0000
Subject: [PATCH] Regress of 3D_mesh_NETGEN/G6 test    * avoid pb that for
 internal node GCPnts_AbscissaPoint::Length() return value larger than total
 edge length    * cash values used for work with non-uniformly paramtrized
 edges

---
 src/StdMeshers/StdMeshers_FaceSide.cxx | 102 +++++++++++--------------
 src/StdMeshers/StdMeshers_FaceSide.hxx |   8 +-
 2 files changed, 50 insertions(+), 60 deletions(-)

diff --git a/src/StdMeshers/StdMeshers_FaceSide.cxx b/src/StdMeshers/StdMeshers_FaceSide.cxx
index f95693298..32f798dee 100644
--- a/src/StdMeshers/StdMeshers_FaceSide.cxx
+++ b/src/StdMeshers/StdMeshers_FaceSide.cxx
@@ -92,9 +92,12 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
   int nbEdges = theEdges.size();
   myEdge.resize( nbEdges );
   myC2d.resize( nbEdges );
+  myC3dAdaptor.resize( nbEdges );
   myFirst.resize( nbEdges );
   myLast.resize( nbEdges );
   myNormPar.resize( nbEdges );
+  myEdgeLength.resize( nbEdges );
+  myIsUniform.resize( nbEdges );
   myLength = 0;
   myNbPonits = myNbSegments = 0;
   myMesh = theMesh;
@@ -104,7 +107,6 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
   if ( nbEdges == 0 ) return;
 
   SMESHDS_Mesh* meshDS = theMesh->GetMeshDS();
-  vector<double> len( nbEdges );
 
   int nbDegen = 0;
   list<TopoDS_Edge>::iterator edge = theEdges.begin();
@@ -112,9 +114,9 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
   for ( int index = 0; edge != theEdges.end(); ++index, ++edge )
   {
     int i = theIsForward ? index : nbEdges - index - 1;
-    len[i] = SMESH_Algo::EdgeLength( *edge );
-    if ( len[i] < DBL_MIN ) nbDegen++;
-    myLength += len[i];
+    myEdgeLength[i] = SMESH_Algo::EdgeLength( *edge );
+    if ( myEdgeLength[i] < DBL_MIN ) nbDegen++;
+    myLength += myEdgeLength[i];
     myEdge[i] = *edge;
     if ( !theIsForward ) myEdge[i].Reverse();
 
@@ -142,6 +144,24 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
       myNbPonits += 1; // for the first end
     else
       myMissingVertexNodes = true;
+
+    // check if edge has non-uniform parametrization (issue 0020705)
+    if ( !myC2d[i].IsNull() )
+    {
+      Geom2dAdaptor_Curve A2dC( myC2d[i] );
+      double p2 = myFirst[i]+(myLast[i]-myFirst[i])/2., p4 = myFirst[i]+(myLast[i]-myFirst[i])/4.;
+      double d2 = GCPnts_AbscissaPoint::Length( A2dC, myFirst[i], p2 );
+      double d4 = GCPnts_AbscissaPoint::Length( A2dC, myFirst[i], p4 );
+      //cout<<"len = "<<len<<"  d2 = "<<d2<<"  fabs(2*d2/len-1.0) = "<<fabs(2*d2/len-1.0)<<endl;
+      myIsUniform[i] = !( fabs(2*d2/myEdgeLength[i]-1.0) > 0.01 || fabs(2*d4/d2-1.0) > 0.01 );
+      if ( !myIsUniform[i] )
+      {
+        double fp,lp;
+        TopLoc_Location L;
+        Handle(Geom_Curve) C3d = BRep_Tool::Curve(myEdge[i],L,fp,lp);
+        myC3dAdaptor[i].Load( C3d, fp,lp );
+      }
+    }
   }
   vExp.Initialize( theEdges.back() );
   if ( vExp.Value().Orientation() != TopAbs_REVERSED ) vExp.Next();
@@ -159,9 +179,9 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace,
       totLength += myLength * degenNormLen * nbDegen;
     double prevNormPar = 0;
     for ( int i = 0; i < nbEdges; ++i ) {
-      if ( len[ i ] < DBL_MIN )
-        len[ i ] = myLength * degenNormLen;
-      myNormPar[ i ] = prevNormPar + len[i]/totLength;
+      if ( myEdgeLength[ i ] < DBL_MIN )
+        myEdgeLength[ i ] = myLength * degenNormLen;
+      myNormPar[ i ] = prevNormPar + myEdgeLength[i]/totLength;
       prevNormPar = myNormPar[ i ];
     }
   }
@@ -197,28 +217,6 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const SMDS_MeshNode* theNode,
   }
 }
 
-
-//=======================================================================
-//function : IsUniform
-//purpose  : auxilary function
-//=======================================================================
-bool IsUniform(const Handle(Geom2d_Curve)& C2d, double fp, double lp)
-{
-  //cout<<"IsUniform  fp = "<<fp<<"  lp = "<<lp<<endl;
-  if(C2d.IsNull())
-    return true;
-  Geom2dAdaptor_Curve A2dC(C2d);
-  double d1 = GCPnts_AbscissaPoint::Length( A2dC, fp, lp );
-  double d2 = GCPnts_AbscissaPoint::Length( A2dC, fp, fp+(lp-fp)/2. );
-  double d4 = GCPnts_AbscissaPoint::Length( A2dC, fp, fp+(lp-fp)/4. );
-  //cout<<"d1 = "<<d1<<"  d2 = "<<d2<<"  fabs(2*d2/d1-1.0) = "<<fabs(2*d2/d1-1.0)<<endl;
-  if( fabs(2*d2/d1-1.0) > 0.01 || fabs(2*d4/d2-1.0) > 0.01 )
-    return false;
-
-  return true;
-}
-
-
 //================================================================================
 /*!
  * \brief Return info on nodes on the side
@@ -234,6 +232,8 @@ const vector<UVPtStruct>& StdMeshers_FaceSide::GetUVPtStruct(bool   isXConst,
     if ( NbEdges() == 0 ) return myPoints;
 
     SMESHDS_Mesh* meshDS = myMesh->GetMeshDS();
+    SMESH_MesherHelper helper(*myMesh);
+    bool paramOK;
 
     // sort nodes of all edges putting them into a map
 
@@ -265,32 +265,27 @@ const vector<UVPtStruct>& StdMeshers_FaceSide::GetUVPtStruct(bool   isXConst,
         u2node.insert( make_pair( 1., node ));
       }
 
-      bool IsUni = IsUniform( myC2d[i], myFirst[i], myLast[i] );
-
       // put internal nodes
       SMESHDS_SubMesh* sm = meshDS->MeshElements( myEdge[i] );
       if ( !sm ) continue;
       SMDS_NodeIteratorPtr nItr = sm->GetNodes();
       double paramSize = myLast[i] - myFirst[i];
       double r = myNormPar[i] - prevNormPar;
-      while ( nItr->more() ) {
+      while ( nItr->more() )
+      {
         const SMDS_MeshNode* node = nItr->next();
         if ( myIgnoreMediumNodes && SMESH_MeshEditor::IsMedium( node, SMDSAbs_Edge ))
           continue;
-        const SMDS_EdgePosition* epos =
-          static_cast<const SMDS_EdgePosition*>(node->GetPosition().get());
-        double u = epos->GetUParameter();
-        // paramSize is signed so orientation is taken into account
+        double u = helper.GetNodeU( myEdge[i], node, &paramOK );
 
+        // paramSize is signed so orientation is taken into account
         double normPar = prevNormPar + r * ( u - myFirst[i] ) / paramSize;
-        if(!IsUni) {
-          double fp,lp;
-          TopLoc_Location L;
-          Handle(Geom_Curve) C3d = BRep_Tool::Curve(myEdge[i],L,fp,lp);
-          GeomAdaptor_Curve A3dC( C3d );
-          double aLen = GCPnts_AbscissaPoint::Length( A3dC, myFirst[i], myLast[i] );
-          double aLenU = GCPnts_AbscissaPoint::Length( A3dC, myFirst[i], u );
-          normPar = prevNormPar + r*aLenU/aLen;
+        if(!myIsUniform[i])
+        {
+          double aLenU = GCPnts_AbscissaPoint::Length
+            ( const_cast<GeomAdaptor_Curve&>( myC3dAdaptor[i]), myFirst[i], u );
+          if ( myEdgeLength[i] > aLenU ) // nonregression test "3D_mesh_NETGEN/G6"
+            normPar = prevNormPar + r*aLenU/myEdgeLength[i];
         }
 #ifdef _DEBUG_
         if ( normPar > 1 || normPar < 0) {
@@ -552,19 +547,11 @@ gp_Pnt2d StdMeshers_FaceSide::Value2d(double U) const
     double par = myFirst[i] * ( 1 - r ) + myLast[i] * r;
     
     // check parametrization of curve
-    if( !IsUniform( myC2d[i], myFirst[i], myLast[i] ) ) {
-      double fp,lp;
-      TopLoc_Location L;
-      Handle(Geom_Curve) C3d = BRep_Tool::Curve(myEdge[i],L,fp,lp);
-      fp = myFirst[i];
-      lp = myLast[i];
-      GeomAdaptor_Curve A3dC( C3d );
-      double aLen3d = GCPnts_AbscissaPoint::Length( A3dC, fp, lp );
-      double aLen3dU = aLen3d*r;
-      if(fp>lp) {
-        aLen3dU = -aLen3dU;
-      }
-      GCPnts_AbscissaPoint AbPnt( A3dC, aLen3dU, fp );
+    if( !myIsUniform[i] )
+    {
+      double aLen3dU = r * myEdgeLength[i] * ( myFirst[i]>myLast[i] ? -1. : 1.);
+      GCPnts_AbscissaPoint AbPnt
+        ( const_cast<GeomAdaptor_Curve&>( myC3dAdaptor[i]), aLen3dU, myFirst[i] );
       if( AbPnt.IsDone() ) {
         par = AbPnt.Parameter();
       }
@@ -572,7 +559,6 @@ gp_Pnt2d StdMeshers_FaceSide::Value2d(double U) const
     return myC2d[ i ]->Value(par);
 
   }
-  //return gp_Pnt2d( 1e+100, 1e+100 );
   return myDefaultPnt2d;
 }
 
diff --git a/src/StdMeshers/StdMeshers_FaceSide.hxx b/src/StdMeshers/StdMeshers_FaceSide.hxx
index 3fd8b81b3..9c1d05fbb 100644
--- a/src/StdMeshers/StdMeshers_FaceSide.hxx
+++ b/src/StdMeshers/StdMeshers_FaceSide.hxx
@@ -29,10 +29,11 @@
 #ifndef StdMeshers_FaceSide_HeaderFile
 #define StdMeshers_FaceSide_HeaderFile
 
-#include <gp_Pnt2d.hxx>
+#include <Geom2d_Curve.hxx>
+#include <GeomAdaptor_Curve.hxx>
 #include <TopoDS_Edge.hxx>
 #include <TopoDS_Vertex.hxx>
-#include <Geom2d_Curve.hxx>
+#include <gp_Pnt2d.hxx>
 
 #include "SMESH_StdMeshers.hxx"
 
@@ -200,8 +201,11 @@ protected:
   std::vector<uvPtStruct>           myPoints, myFalsePoints;
   std::vector<TopoDS_Edge>          myEdge;
   std::vector<Handle(Geom2d_Curve)> myC2d;
+  std::vector<GeomAdaptor_Curve>    myC3dAdaptor;
   std::vector<double>               myFirst, myLast;
   std::vector<double>               myNormPar;
+  std::vector<double>               myEdgeLength;
+  std::vector<double>               myIsUniform;
   double                            myLength;
   int                               myNbPonits, myNbSegments;
   SMESH_Mesh*                       myMesh;
-- 
2.39.2