+//================================================================================
+/*!
+ * \brief Return MaxTolerance( face ), probably cached
+ */
+//================================================================================
+
+double SMESH_MesherHelper::getFaceMaxTol( const TopoDS_Shape& face ) const
+{
+ int faceID = GetMeshDS()->ShapeToIndex( face );
+
+ SMESH_MesherHelper* me = const_cast< SMESH_MesherHelper* >( this );
+ double & tol = me->myFaceMaxTol.insert( make_pair( faceID, -1. )).first->second;
+ if ( tol < 0 )
+ tol = MaxTolerance( face );
+
+ return tol;
+}
+
+//================================================================================
+/*!
+ * \brief Return an angle between two EDGEs sharing a common VERTEX with reference
+ * of the FACE normal
+ * \return double - the angle (between -Pi and Pi), negative if the angle is concave,
+ * 1e100 in case of failure
+ * \waring Care about order of the EDGEs and their orientation to be as they are
+ * within the FACE! Don't pass degenerated EDGEs neither!
+ */
+//================================================================================
+
+double SMESH_MesherHelper::GetAngle( const TopoDS_Edge & theE1,
+ const TopoDS_Edge & theE2,
+ const TopoDS_Face & theFace,
+ const TopoDS_Vertex & theCommonV,
+ gp_Vec* theFaceNormal)
+{
+ double angle = 1e100;
+ try
+ {
+ double f,l;
+ Handle(Geom_Curve) c1 = BRep_Tool::Curve( theE1, f,l );
+ Handle(Geom_Curve) c2 = BRep_Tool::Curve( theE2, f,l );
+ Handle(Geom2d_Curve) c2d1 = BRep_Tool::CurveOnSurface( theE1, theFace, f,l );
+ Handle(Geom_Surface) surf = BRep_Tool::Surface( theFace );
+ double p1 = BRep_Tool::Parameter( theCommonV, theE1 );
+ double p2 = BRep_Tool::Parameter( theCommonV, theE2 );
+ if ( c1.IsNull() || c2.IsNull() )
+ return angle;
+ gp_XY uv = c2d1->Value( p1 ).XY();
+ gp_Vec du, dv; gp_Pnt p;
+ surf->D1( uv.X(), uv.Y(), p, du, dv );
+ gp_Vec vec1, vec2, vecRef = du ^ dv;
+ int nbLoops = 0;
+ double p1tmp = p1;
+ while ( vecRef.SquareMagnitude() < 1e-25 )
+ {
+ double dp = ( l - f ) / 1000.;
+ p1tmp += dp * (( Abs( p1 - f ) > Abs( p1 - l )) ? -1. : +1.);
+ uv = c2d1->Value( p1tmp ).XY();
+ surf->D1( uv.X(), uv.Y(), p, du, dv );
+ vecRef = du ^ dv;
+ if ( ++nbLoops > 10 )
+ {
+#ifdef _DEBUG_
+ cout << "SMESH_MesherHelper::GetAngle(): Captured in a sigularity" << endl;
+#endif
+ return angle;
+ }
+ }
+ if ( theFace.Orientation() == TopAbs_REVERSED )
+ vecRef.Reverse();
+ if ( theFaceNormal ) *theFaceNormal = vecRef;
+
+ c1->D1( p1, p, vec1 );
+ c2->D1( p2, p, vec2 );
+ // TopoDS_Face F = theFace;
+ // if ( F.Orientation() == TopAbs_INTERNAL )
+ // F.Orientation( TopAbs_FORWARD );
+ if ( theE1.Orientation() /*GetSubShapeOri( F, theE1 )*/ == TopAbs_REVERSED )
+ vec1.Reverse();
+ if ( theE2.Orientation() /*GetSubShapeOri( F, theE2 )*/ == TopAbs_REVERSED )
+ vec2.Reverse();
+ angle = vec1.AngleWithRef( vec2, vecRef );
+
+ if ( Abs ( angle ) >= 0.99 * M_PI )
+ {
+ BRep_Tool::Range( theE1, f, l );
+ p1 += 1e-7 * ( p1-f < l-p1 ? +1. : -1. );
+ c1->D1( p1, p, vec1 );
+ if ( theE1.Orientation() == TopAbs_REVERSED )
+ vec1.Reverse();
+ BRep_Tool::Range( theE2, f, l );
+ p2 += 1e-7 * ( p2-f < l-p2 ? +1. : -1. );
+ c2->D1( p2, p, vec2 );
+ if ( theE2.Orientation() == TopAbs_REVERSED )
+ vec2.Reverse();
+ angle = vec1.AngleWithRef( vec2, vecRef );
+ }
+ }
+ catch (...)
+ {
+ }
+ return angle;
+}
+