[occt.git] / src / BRepMesh / BRepMesh_GeomTool.cxx

// Created on: 1993-09-29
// Created by: Isabelle GRIGNON
// Copyright (c) 1993-1999 Matra Datavision
// Copyright (c) 1999-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.

#include <BRepMesh_GeomTool.hxx>

#include <BRepMesh_DefaultRangeSplitter.hxx>

#include <TopAbs_Orientation.hxx>
#include <CSLib.hxx>
#include <Precision.hxx>
#include <Adaptor3d_IsoCurve.hxx>
#include <Adaptor3d_HCurve.hxx>
#include <BRepAdaptor_Curve.hxx>
#include <BRepAdaptor_HSurface.hxx>
#include <Geom2d_Curve.hxx>
#include <BRep_Tool.hxx>

namespace
{
  void ComputeErrFactors (const Standard_Real               theDeflection,
                          const Handle(Adaptor3d_HSurface)& theFace,
                          Standard_Real&                    theErrFactorU,
                          Standard_Real&                    theErrFactorV)
  {
    theErrFactorU = theDeflection * 10.;
    theErrFactorV = theDeflection * 10.;

    switch (theFace->GetType ())
    {
    case GeomAbs_Cylinder:
    case GeomAbs_Cone:
    case GeomAbs_Sphere:
    case GeomAbs_Torus:
      break;

    case GeomAbs_SurfaceOfExtrusion:
    case GeomAbs_SurfaceOfRevolution:
      {
        Handle(Adaptor3d_HCurve) aCurve = theFace->BasisCurve ();
        if (aCurve->GetType () == GeomAbs_BSplineCurve && aCurve->Degree () > 2)
        {
          theErrFactorV /= (aCurve->Degree () * aCurve->NbKnots ());
        }
        break;
      }
    case GeomAbs_BezierSurface:
      {
        if (theFace->UDegree () > 2)
        {
          theErrFactorU /= (theFace->UDegree ());
        }
        if (theFace->VDegree () > 2)
        {
          theErrFactorV /= (theFace->VDegree ());
        }
        break;
      }
    case GeomAbs_BSplineSurface:
      {
        if (theFace->UDegree () > 2)
        {
          theErrFactorU /= (theFace->UDegree () * theFace->NbUKnots ());
        }
        if (theFace->VDegree () > 2)
        {
          theErrFactorV /= (theFace->VDegree () *  theFace->NbVKnots ());
        }
        break;
      }

    case GeomAbs_Plane:
    default:
      theErrFactorU = theErrFactorV = 1.;
    }
  }

  void AdjustCellsCounts (const Handle(Adaptor3d_HSurface)& theFace,
                          const Standard_Integer            theNbVertices,
                          Standard_Integer&                 theCellsCountU,
                          Standard_Integer&                 theCellsCountV)
  {
    const GeomAbs_SurfaceType aType = theFace->GetType ();
    if (aType == GeomAbs_OtherSurface)
    {
      // fallback to the default behavior
      theCellsCountU = theCellsCountV = -1;
      return;
    }

    Standard_Real aSqNbVert = theNbVertices;
    if (aType == GeomAbs_Plane)
    {
      theCellsCountU = theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
    }
    else if (aType == GeomAbs_Cylinder || aType == GeomAbs_Cone)
    {
      theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
    }
    else if (aType == GeomAbs_SurfaceOfExtrusion || aType == GeomAbs_SurfaceOfRevolution)
    {
      Handle (Adaptor3d_HCurve) aCurve = theFace->BasisCurve ();
      if (aCurve->GetType () == GeomAbs_Line ||
         (aCurve->GetType () == GeomAbs_BSplineCurve && aCurve->Degree () < 2))
      {
        // planar, cylindrical, conical cases
        if (aType == GeomAbs_SurfaceOfExtrusion)
          theCellsCountU = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
        else
          theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
      }
      if (aType == GeomAbs_SurfaceOfExtrusion)
      {
        // V is always a line
        theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
      }
    }
    else if (aType == GeomAbs_BezierSurface || aType == GeomAbs_BSplineSurface)
    {
      if (theFace->UDegree () < 2)
      {
        theCellsCountU = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
      }
      if (theFace->VDegree () < 2)
      {
        theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
      }
    }

    theCellsCountU = Max (theCellsCountU, 2);
    theCellsCountV = Max (theCellsCountV, 2);
  }
}

//=======================================================================
//function : Constructor
//purpose  :
//=======================================================================
BRepMesh_GeomTool::BRepMesh_GeomTool(
  const BRepAdaptor_Curve& theCurve,
  const Standard_Real      theFirstParam,
  const Standard_Real      theLastParam,
  const Standard_Real      theLinDeflection,
  const Standard_Real      theAngDeflection,
  const Standard_Integer   theMinPointsNb,
  const Standard_Real      theMinSize)
  : myEdge(&theCurve.Edge()),
    myIsoType(GeomAbs_NoneIso)
{
  myDiscretTool.Initialize(theCurve, theFirstParam, theLastParam,
    theAngDeflection, theLinDeflection, theMinPointsNb, 
    Precision::PConfusion(), theMinSize);
}

//=======================================================================
//function : Constructor
//purpose  :
//=======================================================================
BRepMesh_GeomTool::BRepMesh_GeomTool(
  const Handle(BRepAdaptor_HSurface)& theSurface,
  const GeomAbs_IsoType               theIsoType,
  const Standard_Real                 theParamIso,
  const Standard_Real                 theFirstParam,
  const Standard_Real                 theLastParam,
  const Standard_Real                 theLinDeflection,
  const Standard_Real                 theAngDeflection,
  const Standard_Integer              theMinPointsNb,
  const Standard_Real                 theMinSize)
  : myEdge(NULL),
    myIsoType(theIsoType)
{
  Adaptor3d_IsoCurve aIso(theSurface, theIsoType, theParamIso,
    theFirstParam, theLastParam);

  myDiscretTool.Initialize(aIso, theFirstParam, theLastParam,
    theAngDeflection, theLinDeflection, theMinPointsNb,
    Precision::PConfusion(), theMinSize);
}

//=======================================================================
//function : Value
//purpose  :
//=======================================================================
Standard_Boolean BRepMesh_GeomTool::Value(
  const Standard_Integer              theIndex,
  const Handle(BRepAdaptor_HSurface)& theSurface,
  Standard_Real&                      theParam,
  gp_Pnt&                             thePoint,
  gp_Pnt2d&                           theUV) const
{
  if (theIndex < 1 || theIndex > NbPoints())
    return Standard_False;

  if (myEdge == NULL)
    return Standard_False;

  thePoint = myDiscretTool.Value(theIndex);
  theParam = myDiscretTool.Parameter(theIndex);

  const TopoDS_Face& aFace = ((BRepAdaptor_Surface*)&(theSurface->Surface()))->Face();

  Standard_Real aFirst, aLast;
  Handle(Geom2d_Curve) aCurve = 
    BRep_Tool::CurveOnSurface(*myEdge, aFace, aFirst, aLast);

  aCurve->D0(theParam, theUV);

  return Standard_True;
}

//=======================================================================
//function : Value
//purpose  :
//=======================================================================
Standard_Boolean BRepMesh_GeomTool::Value(
  const Standard_Integer theIndex,
  const Standard_Real    theIsoParam,
  Standard_Real&         theParam,
  gp_Pnt&                thePoint,
  gp_Pnt2d&              theUV) const
{
  if (theIndex < 1 || theIndex > NbPoints())
    return Standard_False;

  thePoint = myDiscretTool.Value(theIndex);
  theParam = myDiscretTool.Parameter(theIndex);

  if (myIsoType == GeomAbs_IsoU)
    theUV.SetCoord(theIsoParam, theParam);
  else
    theUV.SetCoord(theParam, theIsoParam);

  return Standard_True;
}

//=======================================================================
//function : Normal
//purpose  :
//=======================================================================
Standard_Boolean BRepMesh_GeomTool::Normal(
  const Handle(BRepAdaptor_HSurface)& theSurface,
  const Standard_Real                 theParamU,
  const Standard_Real                 theParamV,
  gp_Pnt&                             thePoint,
  gp_Dir&                             theNormal)
{
  Standard_Boolean isOK = Standard_True;
  gp_Vec aD1U, aD1V;

  theSurface->D1(theParamU, theParamV, thePoint, aD1U, aD1V);

  CSLib_DerivativeStatus aStatus;
  CSLib::Normal(aD1U, aD1V, Precision::Angular(), aStatus, theNormal);
  if (aStatus != CSLib_Done)
  {
    gp_Vec aD2U,aD2V,aD2UV;
    theSurface->D2(theParamU, theParamV, thePoint, aD1U, aD1V, aD2U, aD2V, aD2UV);
    CSLib_NormalStatus aNormalStatus;
    CSLib::Normal(aD1U, aD1V, aD2U, aD2V, aD2UV, Precision::Angular(), 
      isOK, aNormalStatus, theNormal);
  }

  if (!isOK)
    return Standard_False;

  const TopoDS_Face& aFace = ((BRepAdaptor_Surface*)&(theSurface->Surface()))->Face();
  TopAbs_Orientation aOri = aFace.Orientation();
  if (aOri == TopAbs_REVERSED)
    theNormal.Reverse();

  return Standard_True;
}

//=============================================================================
//function : IntLinLin
//purpose  : 
//=============================================================================
BRepMesh_GeomTool::IntFlag BRepMesh_GeomTool::IntLinLin(
  const gp_XY&  theStartPnt1,
  const gp_XY&  theEndPnt1,
  const gp_XY&  theStartPnt2,
  const gp_XY&  theEndPnt2,
  gp_XY&        theIntPnt,
  Standard_Real (&theParamOnSegment)[2])
{
  gp_XY aVec1    = theEndPnt1   - theStartPnt1;
  gp_XY aVec2    = theEndPnt2   - theStartPnt2;
  gp_XY aVecO1O2 = theStartPnt2 - theStartPnt1;
    
  Standard_Real aCrossD1D2 = aVec1    ^ aVec2;
  Standard_Real aCrossD1D3 = aVecO1O2 ^ aVec2;

  const Standard_Real aPrec = gp::Resolution();
  // Are edgegs codirectional
  if ( Abs( aCrossD1D2 ) < aPrec )
  {
    // Just a parallel case?
    if( Abs( aCrossD1D3 ) < aPrec )
      return BRepMesh_GeomTool::Same;
    else
      return BRepMesh_GeomTool::NoIntersection;
  }

  theParamOnSegment[0] = aCrossD1D3 / aCrossD1D2;
  theIntPnt = theStartPnt1 + theParamOnSegment[0] * aVec1;

  Standard_Real aCrossD2D3 = aVecO1O2 ^ aVec1;
  theParamOnSegment[1] = aCrossD2D3 / aCrossD1D2;

  return BRepMesh_GeomTool::Cross;
}

//=============================================================================
//function : IntSegSeg
//purpose  : 
//=============================================================================
BRepMesh_GeomTool::IntFlag BRepMesh_GeomTool::IntSegSeg(
  const gp_XY&           theStartPnt1,
  const gp_XY&           theEndPnt1,
  const gp_XY&           theStartPnt2,
  const gp_XY&           theEndPnt2,
  const Standard_Boolean isConsiderEndPointTouch,
  const Standard_Boolean isConsiderPointOnSegment,
  gp_Pnt2d&              theIntPnt)
{
  Standard_Integer aPointHash[] = {
    classifyPoint(theStartPnt1, theEndPnt1, theStartPnt2),
    classifyPoint(theStartPnt1, theEndPnt1, theEndPnt2  ),
    classifyPoint(theStartPnt2, theEndPnt2, theStartPnt1),
    classifyPoint(theStartPnt2, theEndPnt2, theEndPnt1  )
  };

  // Consider case when edges have shared vertex
  if ( aPointHash[0] < 0 || aPointHash[1] < 0 )
  {
    if ( isConsiderEndPointTouch )
      return BRepMesh_GeomTool::EndPointTouch;

    return BRepMesh_GeomTool::NoIntersection;
  }

  Standard_Integer aPosHash = 
    aPointHash[0] + aPointHash[1] + aPointHash[2] + aPointHash[3];

  /*=========================================*/
  /*  1) hash code == 1:

                    0+
                    /
           0      1/         0
           +======+==========+
  
      2) hash code == 2:

           0    1        1   0
        a) +----+========+---+

           0       1   1     0
        b) +-------+===+=====+

                                             */
  /*=========================================*/
  if ( aPosHash == 1 )
  {
    if (isConsiderPointOnSegment)
    {
      if (aPointHash[0] == 1)
        theIntPnt = theStartPnt1;
      else if (aPointHash[1] == 1)
        theIntPnt = theEndPnt1;
      else if (aPointHash[2] == 1)
        theIntPnt = theStartPnt2;
      else
        theIntPnt = theEndPnt2;

      return BRepMesh_GeomTool::PointOnSegment;
    }

    return BRepMesh_GeomTool::NoIntersection;
  }
  else if ( aPosHash == 2 )
    return BRepMesh_GeomTool::Glued;

  Standard_Real aParam[2];
  IntFlag aIntFlag = IntLinLin(theStartPnt1, theEndPnt1, 
    theStartPnt2, theEndPnt2, theIntPnt.ChangeCoord(), aParam);

  if (aIntFlag == BRepMesh_GeomTool::NoIntersection)
    return BRepMesh_GeomTool::NoIntersection;

  if (aIntFlag == BRepMesh_GeomTool::Same)
  {
    if ( aPosHash < -2 )
      return BRepMesh_GeomTool::Same;
    else if ( aPosHash == -1 )
      return BRepMesh_GeomTool::Glued;

    return BRepMesh_GeomTool::NoIntersection;
  }

  // Cross
  // Intersection is out of segments ranges
  const Standard_Real aPrec    = Precision::PConfusion();
  const Standard_Real aEndPrec = 1 - aPrec;
  for (Standard_Integer i = 0; i < 2; ++i)
  {
    if(aParam[i] < aPrec || aParam[i] > aEndPrec )
      return BRepMesh_GeomTool::NoIntersection;
  }
 
  return BRepMesh_GeomTool::Cross;
}

//=============================================================================
//function : CellsCount
//purpose  :
//=============================================================================
std::pair<Standard_Integer, Standard_Integer> BRepMesh_GeomTool::CellsCount (
  const Handle (Adaptor3d_HSurface)&   theSurface,
  const Standard_Integer               theVerticesNb,
  const Standard_Real                  theDeflection,
  const BRepMesh_DefaultRangeSplitter* theRangeSplitter)
{
  if (theRangeSplitter == NULL)
    return std::pair<Standard_Integer, Standard_Integer>(-1, -1);

  const GeomAbs_SurfaceType aType = theSurface->GetType ();

  Standard_Real anErrFactorU, anErrFactorV;
  ComputeErrFactors(theDeflection, theSurface, anErrFactorU, anErrFactorV);

  const std::pair<Standard_Real, Standard_Real>& aRangeU = theRangeSplitter->GetRangeU();
  const std::pair<Standard_Real, Standard_Real>& aRangeV = theRangeSplitter->GetRangeV();
  const std::pair<Standard_Real, Standard_Real>& aDelta  = theRangeSplitter->GetDelta ();

  Standard_Integer aCellsCountU, aCellsCountV;
  if (aType == GeomAbs_Torus)
  {
    aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
      (aRangeU.second - aRangeU.first) / aDelta.first)));
    aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
      (aRangeV.second - aRangeV.first) / aDelta.second)));
  }
  else if (aType == GeomAbs_Cylinder)
  {
    aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
      (aRangeU.second - aRangeU.first) / aDelta.first /
      (aRangeV.second - aRangeV.first))));
    aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
      (aRangeV.second - aRangeV.first) / anErrFactorV)));
  }
  else
  {
    aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
      (aRangeU.second - aRangeU.first) / aDelta.first  / anErrFactorU)));
    aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
      (aRangeV.second - aRangeV.first) / aDelta.second / anErrFactorV)));
  }

  AdjustCellsCounts(theSurface, theVerticesNb, aCellsCountU, aCellsCountV);
  return std::pair<Standard_Integer, Standard_Integer>(aCellsCountU, aCellsCountV);
}

//=============================================================================
//function : classifyPoint
//purpose  : 
//=============================================================================
Standard_Integer BRepMesh_GeomTool::classifyPoint(
  const gp_XY& thePoint1,
  const gp_XY& thePoint2,
  const gp_XY& thePointToCheck)
{
  gp_XY aP1 = thePoint2       - thePoint1;
  gp_XY aP2 = thePointToCheck - thePoint1;
  
  const Standard_Real aPrec   = Precision::PConfusion();
  const Standard_Real aSqPrec = aPrec * aPrec;
  Standard_Real aDist = Abs(aP1 ^ aP2);
  if (aDist > aPrec)
  {
    aDist = (aDist * aDist) / aP1.SquareModulus();
    if (aDist > aSqPrec)
      return 0; //out
  }
    
  gp_XY aMult = aP1.Multiplied(aP2);
  if ( aMult.X() < 0.0 || aMult.Y() < 0.0 )
    return 0; //out
    
  if (aP1.SquareModulus() < aP2.SquareModulus())
    return 0; //out
    
  if (thePointToCheck.IsEqual(thePoint1, aPrec) || 
      thePointToCheck.IsEqual(thePoint2, aPrec))
  {
    return -1; //coinsides with an end point
  }
    
  return 1;
}
Commit	Line	Data
b311480e	1	// Created on: 1993-09-29
	2	// Created by: Isabelle GRIGNON
	3	// Copyright (c) 1993-1999 Matra Datavision
973c2be1	4	// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e	5	//
973c2be1	6	// This file is part of Open CASCADE Technology software library.
b311480e	7	//
d5f74e42	8	// This library is free software; you can redistribute it and/or modify it under
d5f74e42	9	// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1	10	// by the Free Software Foundation, with special exception defined in the file
	11	// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
	12	// distribution for complete text of the license and disclaimer of any warranty.
b311480e	13	//
973c2be1	14	// Alternatively, this file may be used under the terms of Open CASCADE
973c2be1	15	// commercial license or contractual agreement.
7fd59977	16
fc9b36d6	17	#include <BRepMesh_GeomTool.hxx>
7fd59977	18
7bd071ed	19	#include <BRepMesh_DefaultRangeSplitter.hxx>
7bd071ed	20
7fd59977	21	#include <TopAbs_Orientation.hxx>
7fd59977	22	#include <CSLib.hxx>
	23	#include <Precision.hxx>
	24	#include <Adaptor3d_IsoCurve.hxx>
7bd071ed	25	#include <Adaptor3d_HCurve.hxx>
7fd59977	26	#include <BRepAdaptor_Curve.hxx>
fc9b36d6	27	#include <BRepAdaptor_HSurface.hxx>
	28	#include <Geom2d_Curve.hxx>
	29	#include <BRep_Tool.hxx>
7fd59977	30
7bd071ed	31	namespace
	32	{
	33	void ComputeErrFactors (const Standard_Real theDeflection,
	34	const Handle(Adaptor3d_HSurface)& theFace,
	35	Standard_Real& theErrFactorU,
	36	Standard_Real& theErrFactorV)
	37	{
	38	theErrFactorU = theDeflection * 10.;
	39	theErrFactorV = theDeflection * 10.;
	40
	41	switch (theFace->GetType ())
	42	{
	43	case GeomAbs_Cylinder:
	44	case GeomAbs_Cone:
	45	case GeomAbs_Sphere:
	46	case GeomAbs_Torus:
	47	break;
	48
	49	case GeomAbs_SurfaceOfExtrusion:
	50	case GeomAbs_SurfaceOfRevolution:
	51	{
	52	Handle(Adaptor3d_HCurve) aCurve = theFace->BasisCurve ();
	53	if (aCurve->GetType () == GeomAbs_BSplineCurve && aCurve->Degree () > 2)
	54	{
	55	theErrFactorV /= (aCurve->Degree () * aCurve->NbKnots ());
	56	}
	57	break;
	58	}
	59	case GeomAbs_BezierSurface:
	60	{
	61	if (theFace->UDegree () > 2)
	62	{
	63	theErrFactorU /= (theFace->UDegree ());
	64	}
	65	if (theFace->VDegree () > 2)
	66	{
	67	theErrFactorV /= (theFace->VDegree ());
	68	}
	69	break;
	70	}
	71	case GeomAbs_BSplineSurface:
	72	{
	73	if (theFace->UDegree () > 2)
	74	{
	75	theErrFactorU /= (theFace->UDegree () * theFace->NbUKnots ());
	76	}
	77	if (theFace->VDegree () > 2)
	78	{
	79	theErrFactorV /= (theFace->VDegree () * theFace->NbVKnots ());
	80	}
	81	break;
	82	}
	83
	84	case GeomAbs_Plane:
	85	default:
	86	theErrFactorU = theErrFactorV = 1.;
	87	}
	88	}
	89
	90	void AdjustCellsCounts (const Handle(Adaptor3d_HSurface)& theFace,
	91	const Standard_Integer theNbVertices,
	92	Standard_Integer& theCellsCountU,
	93	Standard_Integer& theCellsCountV)
	94	{
95	const GeomAbs_SurfaceType aType = theFace->GetType ();
96	if (aType == GeomAbs_OtherSurface)
97	{
98	// fallback to the default behavior
99	theCellsCountU = theCellsCountV = -1;
100	return;
101	}
102
103	Standard_Real aSqNbVert = theNbVertices;
104	if (aType == GeomAbs_Plane)
105	{
106	theCellsCountU = theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
107	}
108	else if (aType == GeomAbs_Cylinder \|\| aType == GeomAbs_Cone)
109	{
110	theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
111	}
112	else if (aType == GeomAbs_SurfaceOfExtrusion \|\| aType == GeomAbs_SurfaceOfRevolution)
113	{
114	Handle (Adaptor3d_HCurve) aCurve = theFace->BasisCurve ();
115	if (aCurve->GetType () == GeomAbs_Line \|\|
116	(aCurve->GetType () == GeomAbs_BSplineCurve && aCurve->Degree () < 2))
117	{
118	// planar, cylindrical, conical cases
119	if (aType == GeomAbs_SurfaceOfExtrusion)
120	theCellsCountU = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
121	else
122	theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
123	}
124	if (aType == GeomAbs_SurfaceOfExtrusion)
125	{
126	// V is always a line
127	theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
128	}
129	}
130	else if (aType == GeomAbs_BezierSurface \|\| aType == GeomAbs_BSplineSurface)
131	{
132	if (theFace->UDegree () < 2)
133	{
134	theCellsCountU = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
135	}
136	if (theFace->VDegree () < 2)
137	{
138	theCellsCountV = (Standard_Integer)Ceiling (Pow (2, Log10 (aSqNbVert)));
139	}
140	}
141
142	theCellsCountU = Max (theCellsCountU, 2);
143	theCellsCountV = Max (theCellsCountV, 2);
144	}
145	}
146
fc9b36d6	147	//=======================================================================
	148	//function : Constructor
	149	//purpose :
	150	//=======================================================================
	151	BRepMesh_GeomTool::BRepMesh_GeomTool(
	152	const BRepAdaptor_Curve& theCurve,
	153	const Standard_Real theFirstParam,
	154	const Standard_Real theLastParam,
	155	const Standard_Real theLinDeflection,
	156	const Standard_Real theAngDeflection,
74da0216	157	const Standard_Integer theMinPointsNb,
74da0216	158	const Standard_Real theMinSize)
fc9b36d6	159	: myEdge(&theCurve.Edge()),
fc9b36d6	160	myIsoType(GeomAbs_NoneIso)
7fd59977	161	{
fc9b36d6	162	myDiscretTool.Initialize(theCurve, theFirstParam, theLastParam,
74da0216	163	theAngDeflection, theLinDeflection, theMinPointsNb,
74da0216	164	Precision::PConfusion(), theMinSize);
7fd59977	165	}
fc9b36d6	166
	167	//=======================================================================
	168	//function : Constructor
	169	//purpose :
	170	//=======================================================================
	171	BRepMesh_GeomTool::BRepMesh_GeomTool(
	172	const Handle(BRepAdaptor_HSurface)& theSurface,
	173	const GeomAbs_IsoType theIsoType,
	174	const Standard_Real theParamIso,
	175	const Standard_Real theFirstParam,
	176	const Standard_Real theLastParam,
	177	const Standard_Real theLinDeflection,
	178	const Standard_Real theAngDeflection,
74da0216	179	const Standard_Integer theMinPointsNb,
74da0216	180	const Standard_Real theMinSize)
fc9b36d6	181	: myEdge(NULL),
fc9b36d6	182	myIsoType(theIsoType)
7fd59977	183	{
fc9b36d6	184	Adaptor3d_IsoCurve aIso(theSurface, theIsoType, theParamIso,
	185	theFirstParam, theLastParam);
	186
	187	myDiscretTool.Initialize(aIso, theFirstParam, theLastParam,
74da0216	188	theAngDeflection, theLinDeflection, theMinPointsNb,
74da0216	189	Precision::PConfusion(), theMinSize);
7fd59977	190	}
7fd59977	191
fc9b36d6	192	//=======================================================================
	193	//function : Value
	194	//purpose :
	195	//=======================================================================
	196	Standard_Boolean BRepMesh_GeomTool::Value(
	197	const Standard_Integer theIndex,
7bd071ed	198	const Handle(BRepAdaptor_HSurface)& theSurface,
fc9b36d6	199	Standard_Real& theParam,
7bd071ed	200	gp_Pnt& thePoint,
7bd071ed	201	gp_Pnt2d& theUV) const
7fd59977	202	{
fc9b36d6	203	if (theIndex < 1 \|\| theIndex > NbPoints())
	204	return Standard_False;
	205
	206	if (myEdge == NULL)
	207	return Standard_False;
	208
	209	thePoint = myDiscretTool.Value(theIndex);
	210	theParam = myDiscretTool.Parameter(theIndex);
	211
7bd071ed	212	const TopoDS_Face& aFace = ((BRepAdaptor_Surface*)&(theSurface->Surface()))->Face();
	213
	214	Standard_Real aFirst, aLast;
	215	Handle(Geom2d_Curve) aCurve =
	216	BRep_Tool::CurveOnSurface(*myEdge, aFace, aFirst, aLast);
	217
	218	aCurve->D0(theParam, theUV);
	219
fc9b36d6	220	return Standard_True;
7fd59977	221	}
7fd59977	222
fc9b36d6	223	//=======================================================================
	224	//function : Value
	225	//purpose :
	226	//=======================================================================
7bd071ed	227	Standard_Boolean BRepMesh_GeomTool::Value(
	228	const Standard_Integer theIndex,
	229	const Standard_Real theIsoParam,
	230	Standard_Real& theParam,
	231	gp_Pnt& thePoint,
	232	gp_Pnt2d& theUV) const
fc9b36d6	233	{
	234	if (theIndex < 1 \|\| theIndex > NbPoints())
	235	return Standard_False;
	236
	237	thePoint = myDiscretTool.Value(theIndex);
	238	theParam = myDiscretTool.Parameter(theIndex);
	239
	240	if (myIsoType == GeomAbs_IsoU)
	241	theUV.SetCoord(theIsoParam, theParam);
	242	else
	243	theUV.SetCoord(theParam, theIsoParam);
	244
	245	return Standard_True;
7fd59977	246	}
7fd59977	247
fc9b36d6	248	//=======================================================================
	249	//function : Normal
	250	//purpose :
	251	//=======================================================================
	252	Standard_Boolean BRepMesh_GeomTool::Normal(
	253	const Handle(BRepAdaptor_HSurface)& theSurface,
	254	const Standard_Real theParamU,
	255	const Standard_Real theParamV,
	256	gp_Pnt& thePoint,
	257	gp_Dir& theNormal)
7fd59977	258	{
fc9b36d6	259	Standard_Boolean isOK = Standard_True;
	260	gp_Vec aD1U, aD1V;
	261
	262	theSurface->D1(theParamU, theParamV, thePoint, aD1U, aD1V);
	263
	264	CSLib_DerivativeStatus aStatus;
	265	CSLib::Normal(aD1U, aD1V, Precision::Angular(), aStatus, theNormal);
	266	if (aStatus != CSLib_Done)
	267	{
	268	gp_Vec aD2U,aD2V,aD2UV;
	269	theSurface->D2(theParamU, theParamV, thePoint, aD1U, aD1V, aD2U, aD2V, aD2UV);
	270	CSLib_NormalStatus aNormalStatus;
	271	CSLib::Normal(aD1U, aD1V, aD2U, aD2V, aD2UV, Precision::Angular(),
	272	isOK, aNormalStatus, theNormal);
7fd59977	273	}
fc9b36d6	274
	275	if (!isOK)
	276	return Standard_False;
	277
	278	const TopoDS_Face& aFace = ((BRepAdaptor_Surface*)&(theSurface->Surface()))->Face();
	279	TopAbs_Orientation aOri = aFace.Orientation();
	280	if (aOri == TopAbs_REVERSED)
	281	theNormal.Reverse();
	282
	283	return Standard_True;
7fd59977	284	}
7fd59977	285
fc9b36d6	286	//=============================================================================
	287	//function : IntLinLin
	288	//purpose :
	289	//=============================================================================
	290	BRepMesh_GeomTool::IntFlag BRepMesh_GeomTool::IntLinLin(
	291	const gp_XY& theStartPnt1,
	292	const gp_XY& theEndPnt1,
	293	const gp_XY& theStartPnt2,
	294	const gp_XY& theEndPnt2,
	295	gp_XY& theIntPnt,
	296	Standard_Real (&theParamOnSegment)[2])
7fd59977	297	{
fc9b36d6	298	gp_XY aVec1 = theEndPnt1 - theStartPnt1;
	299	gp_XY aVec2 = theEndPnt2 - theStartPnt2;
	300	gp_XY aVecO1O2 = theStartPnt2 - theStartPnt1;
	301
	302	Standard_Real aCrossD1D2 = aVec1 ^ aVec2;
	303	Standard_Real aCrossD1D3 = aVecO1O2 ^ aVec2;
	304
	305	const Standard_Real aPrec = gp::Resolution();
	306	// Are edgegs codirectional
	307	if ( Abs( aCrossD1D2 ) < aPrec )
	308	{
	309	// Just a parallel case?
	310	if( Abs( aCrossD1D3 ) < aPrec )
	311	return BRepMesh_GeomTool::Same;
	312	else
	313	return BRepMesh_GeomTool::NoIntersection;
	314	}
	315
	316	theParamOnSegment[0] = aCrossD1D3 / aCrossD1D2;
	317	theIntPnt = theStartPnt1 + theParamOnSegment[0] * aVec1;
	318
	319	Standard_Real aCrossD2D3 = aVecO1O2 ^ aVec1;
	320	theParamOnSegment[1] = aCrossD2D3 / aCrossD1D2;
	321
	322	return BRepMesh_GeomTool::Cross;
7fd59977	323	}
7fd59977	324
fc9b36d6	325	//=============================================================================
	326	//function : IntSegSeg
	327	//purpose :
	328	//=============================================================================
	329	BRepMesh_GeomTool::IntFlag BRepMesh_GeomTool::IntSegSeg(
	330	const gp_XY& theStartPnt1,
	331	const gp_XY& theEndPnt1,
	332	const gp_XY& theStartPnt2,
	333	const gp_XY& theEndPnt2,
	334	const Standard_Boolean isConsiderEndPointTouch,
	335	const Standard_Boolean isConsiderPointOnSegment,
	336	gp_Pnt2d& theIntPnt)
7fd59977	337	{
fc9b36d6	338	Standard_Integer aPointHash[] = {
	339	classifyPoint(theStartPnt1, theEndPnt1, theStartPnt2),
	340	classifyPoint(theStartPnt1, theEndPnt1, theEndPnt2 ),
	341	classifyPoint(theStartPnt2, theEndPnt2, theStartPnt1),
	342	classifyPoint(theStartPnt2, theEndPnt2, theEndPnt1 )
	343	};
	344
	345	// Consider case when edges have shared vertex
71316196	346	if ( aPointHash[0] < 0 \|\| aPointHash[1] < 0 )
fc9b36d6	347	{
71316196	348	if ( isConsiderEndPointTouch )
fc9b36d6	349	return BRepMesh_GeomTool::EndPointTouch;
71316196	350
71316196	351	return BRepMesh_GeomTool::NoIntersection;
fc9b36d6	352	}
	353
	354	Standard_Integer aPosHash =
	355	aPointHash[0] + aPointHash[1] + aPointHash[2] + aPointHash[3];
	356
	357	/=========================================/
	358	/* 1) hash code == 1:
	359
	360	0+
	361	/
	362	0 1/ 0
	363	+======+==========+
	364
	365	2) hash code == 2:
	366
	367	0 1 1 0
	368	a) +----+========+---+
	369
	370	0 1 1 0
	371	b) +-------+===+=====+
	372
	373	*/
	374	/=========================================/
	375	if ( aPosHash == 1 )
	376	{
	377	if (isConsiderPointOnSegment)
	378	{
	379	if (aPointHash[0] == 1)
	380	theIntPnt = theStartPnt1;
	381	else if (aPointHash[1] == 1)
	382	theIntPnt = theEndPnt1;
	383	else if (aPointHash[2] == 1)
	384	theIntPnt = theStartPnt2;
	385	else
	386	theIntPnt = theEndPnt2;
	387
	388	return BRepMesh_GeomTool::PointOnSegment;
	389	}
	390
	391	return BRepMesh_GeomTool::NoIntersection;
	392	}
	393	else if ( aPosHash == 2 )
	394	return BRepMesh_GeomTool::Glued;
	395
	396	Standard_Real aParam[2];
	397	IntFlag aIntFlag = IntLinLin(theStartPnt1, theEndPnt1,
	398	theStartPnt2, theEndPnt2, theIntPnt.ChangeCoord(), aParam);
	399
	400	if (aIntFlag == BRepMesh_GeomTool::NoIntersection)
	401	return BRepMesh_GeomTool::NoIntersection;
	402
	403	if (aIntFlag == BRepMesh_GeomTool::Same)
	404	{
	405	if ( aPosHash < -2 )
	406	return BRepMesh_GeomTool::Same;
	407	else if ( aPosHash == -1 )
	408	return BRepMesh_GeomTool::Glued;
	409
	410	return BRepMesh_GeomTool::NoIntersection;
	411	}
	412
	413	// Cross
	414	// Intersection is out of segments ranges
	415	const Standard_Real aPrec = Precision::PConfusion();
416	const Standard_Real aEndPrec = 1 - aPrec;
417	for (Standard_Integer i = 0; i < 2; ++i)
418	{
7bd071ed	419	if(aParam[i] < aPrec \|\| aParam[i] > aEndPrec )
fc9b36d6	420	return BRepMesh_GeomTool::NoIntersection;
	421	}
	422
	423	return BRepMesh_GeomTool::Cross;
	424	}
	425
	426	//=============================================================================
7bd071ed	427	//function : CellsCount
	428	//purpose :
	429	//=============================================================================
	430	std::pair<Standard_Integer, Standard_Integer> BRepMesh_GeomTool::CellsCount (
	431	const Handle (Adaptor3d_HSurface)& theSurface,
	432	const Standard_Integer theVerticesNb,
	433	const Standard_Real theDeflection,
	434	const BRepMesh_DefaultRangeSplitter* theRangeSplitter)
	435	{
	436	if (theRangeSplitter == NULL)
	437	return std::pair<Standard_Integer, Standard_Integer>(-1, -1);
	438
	439	const GeomAbs_SurfaceType aType = theSurface->GetType ();
	440
	441	Standard_Real anErrFactorU, anErrFactorV;
	442	ComputeErrFactors(theDeflection, theSurface, anErrFactorU, anErrFactorV);
	443
	444	const std::pair<Standard_Real, Standard_Real>& aRangeU = theRangeSplitter->GetRangeU();
	445	const std::pair<Standard_Real, Standard_Real>& aRangeV = theRangeSplitter->GetRangeV();
	446	const std::pair<Standard_Real, Standard_Real>& aDelta = theRangeSplitter->GetDelta ();
	447
	448	Standard_Integer aCellsCountU, aCellsCountV;
	449	if (aType == GeomAbs_Torus)
	450	{
	451	aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
	452	(aRangeU.second - aRangeU.first) / aDelta.first)));
	453	aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
	454	(aRangeV.second - aRangeV.first) / aDelta.second)));
	455	}
	456	else if (aType == GeomAbs_Cylinder)
	457	{
	458	aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
	459	(aRangeU.second - aRangeU.first) / aDelta.first /
	460	(aRangeV.second - aRangeV.first))));
	461	aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
	462	(aRangeV.second - aRangeV.first) / anErrFactorV)));
	463	}
	464	else
	465	{
	466	aCellsCountU = (Standard_Integer)Ceiling(Pow(2, Log10(
	467	(aRangeU.second - aRangeU.first) / aDelta.first / anErrFactorU)));
	468	aCellsCountV = (Standard_Integer)Ceiling(Pow(2, Log10(
	469	(aRangeV.second - aRangeV.first) / aDelta.second / anErrFactorV)));
	470	}
	471
	472	AdjustCellsCounts(theSurface, theVerticesNb, aCellsCountU, aCellsCountV);
	473	return std::pair<Standard_Integer, Standard_Integer>(aCellsCountU, aCellsCountV);
	474	}
	475
	476	//=============================================================================
fc9b36d6	477	//function : classifyPoint
	478	//purpose :
	479	//=============================================================================
	480	Standard_Integer BRepMesh_GeomTool::classifyPoint(
	481	const gp_XY& thePoint1,
	482	const gp_XY& thePoint2,
	483	const gp_XY& thePointToCheck)
7fd59977	484	{
fc9b36d6	485	gp_XY aP1 = thePoint2 - thePoint1;
	486	gp_XY aP2 = thePointToCheck - thePoint1;
	487
	488	const Standard_Real aPrec = Precision::PConfusion();
	489	const Standard_Real aSqPrec = aPrec * aPrec;
	490	Standard_Real aDist = Abs(aP1 ^ aP2);
	491	if (aDist > aPrec)
	492	{
	493	aDist = (aDist * aDist) / aP1.SquareModulus();
	494	if (aDist > aSqPrec)
	495	return 0; //out
7fd59977	496	}
fc9b36d6	497
	498	gp_XY aMult = aP1.Multiplied(aP2);
	499	if ( aMult.X() < 0.0 \|\| aMult.Y() < 0.0 )
	500	return 0; //out
	501
	502	if (aP1.SquareModulus() < aP2.SquareModulus())
	503	return 0; //out
	504
	505	if (thePointToCheck.IsEqual(thePoint1, aPrec) \|\|
	506	thePointToCheck.IsEqual(thePoint2, aPrec))
	507	{
	508	return -1; //coinsides with an end point
7fd59977	509	}
fc9b36d6	510
fc9b36d6	511	return 1;
7fd59977	512	}