[occt.git] / src / Extrema / Extrema_ExtCS.cxx

// Created on: 1995-07-18
// Created by: Modelistation
// Copyright (c) 1995-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.

//  Modified by skv - Thu Jul  7 12:29:34 2005 OCC9134

#include <Adaptor3d_Curve.hxx>
#include <Adaptor3d_Surface.hxx>
#include <Bnd_Box.hxx>
#include <BndLib_AddSurface.hxx>
#include <ElCLib.hxx>
#include <ElSLib.hxx>
#include <Extrema_ExtCS.hxx>
#include <Extrema_ExtPS.hxx>
#include <Extrema_GenExtCS.hxx>
#include <Extrema_POnCurv.hxx>
#include <Extrema_POnSurf.hxx>
#include <GeomAbs_CurveType.hxx>
#include <gp_Cone.hxx>
#include <gp_Cylinder.hxx>
#include <gp_Lin.hxx>
#include <gp_Pln.hxx>
#include <gp_Pnt.hxx>
#include <gp_Sphere.hxx>
#include <gp_Torus.hxx>
#include <Precision.hxx>
#include <Standard_OutOfRange.hxx>
#include <StdFail_NotDone.hxx>
#include <TColStd_Array1OfReal.hxx>
#include <Extrema_ExtPS.hxx>

Extrema_ExtCS::Extrema_ExtCS()
: myS(NULL),
  myDone(Standard_False),
  myIsPar(Standard_False),
  myuinf(0.0),
  myusup(0.0),
  myvinf(0.0),
  myvsup(0.0),
  mytolC(0.0),
  mytolS(0.0),
  myucinf(0.0),
  myucsup(0.0),
  myStype(GeomAbs_OtherSurface)
{
}

Extrema_ExtCS::Extrema_ExtCS(const Adaptor3d_Curve&   C,
  const Adaptor3d_Surface& S,
  const Standard_Real    TolC,
  const Standard_Real    TolS)

{
  Initialize(S, S.FirstUParameter(), S.LastUParameter(), 
    S.FirstVParameter(), S.LastVParameter(), 
    TolC, TolS);
  Perform(C, C.FirstParameter(), C.LastParameter());
}

Extrema_ExtCS::Extrema_ExtCS(const Adaptor3d_Curve&   C,
  const Adaptor3d_Surface& S,
  const Standard_Real    UCinf,
  const Standard_Real    UCsup,
  const Standard_Real    Uinf,	
  const Standard_Real    Usup,
  const Standard_Real    Vinf,	
  const Standard_Real    Vsup,
  const Standard_Real    TolC,
  const Standard_Real    TolS)

{
  Initialize(S, Uinf, Usup, Vinf, Vsup, TolC, TolS);
  Perform(C, UCinf, UCsup);
}


void Extrema_ExtCS::Initialize(const Adaptor3d_Surface& S,
  const Standard_Real    Uinf,	
  const Standard_Real    Usup,
  const Standard_Real    Vinf,	
  const Standard_Real    Vsup,
  const Standard_Real    TolC,
  const Standard_Real    TolS)
{
  myS = (Adaptor3d_SurfacePtr)&S;
  myIsPar = Standard_False;
  myuinf  = Uinf;
  myusup  = Usup;
  myvinf  = Vinf;
  myvsup  = Vsup;
  mytolC  = TolC;
  mytolS  = TolS;
  myStype  = myS->GetType();
}


void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
  const Standard_Real  Uinf,
  const Standard_Real  Usup)
{
  myucinf = Uinf;
  myucsup = Usup;
  myPOnS.Clear();
  myPOnC.Clear();
  mySqDist.Clear();
  Standard_Integer i, j;
  Standard_Integer NbT, NbU, NbV;
  NbT = 12; NbU = NbV = 10;
  GeomAbs_CurveType myCtype  = C.GetType();

  myDone = Standard_False;
  // Try analytic computation of extrema
  Standard_Boolean isComputeAnalytic = Standard_True;

  switch(myCtype) {

  case GeomAbs_Line: 
    {

      switch(myStype) {
      case GeomAbs_Sphere:
        myExtElCS.Perform(C.Line(), myS->Sphere());
        break;   
      case GeomAbs_Cylinder:
        myExtElCS.Perform(C.Line(), myS->Cylinder());
        break;
      case GeomAbs_Plane:
        myExtElCS.Perform(C.Line(), myS->Plane());
        if (myExtElCS.IsParallel())   break;
        Standard_FALLTHROUGH

      case GeomAbs_Torus:
      case GeomAbs_Cone:
      case GeomAbs_BezierSurface:
      case GeomAbs_BSplineSurface:
      case GeomAbs_SurfaceOfRevolution:
      case GeomAbs_SurfaceOfExtrusion:
      case GeomAbs_OffsetSurface:
      case GeomAbs_OtherSurface:
        {
          Standard_Real cfirst = myucinf, clast = myucsup;
          Standard_Real ufirst = myS->FirstUParameter(), ulast = myS->LastUParameter(), 
            vfirst = myS->FirstVParameter(), vlast = myS->LastVParameter();

          if (!(Precision::IsInfinite(ufirst) || Precision::IsInfinite(ulast) ||
            Precision::IsInfinite(vfirst) || Precision::IsInfinite(vlast)))
          {
            Standard_Real tmin = Precision::Infinite(), tmax = -tmin;
            Standard_Real xmin, ymin, zmin, xmax, ymax, zmax;
            Bnd_Box aSurfBox;
            BndLib_AddSurface::Add(*myS, ufirst, ulast, vfirst, vlast, Precision::Confusion(), aSurfBox);
            aSurfBox.Get(xmin, ymin, zmin, xmax, ymax, zmax);
            gp_Lin aLin = C.Line();
            Standard_Real aParOnLin;
            gp_Pnt aLimPntArray[8];

            aLimPntArray[0].SetCoord(xmin, ymin, zmin);
            aLimPntArray[1].SetCoord(xmax, ymin, zmin);
            aLimPntArray[2].SetCoord(xmin, ymax, zmin);
            aLimPntArray[3].SetCoord(xmax, ymax, zmin);
            aLimPntArray[4].SetCoord(xmin, ymin, zmax);
            aLimPntArray[5].SetCoord(xmax, ymin, zmax);
            aLimPntArray[6].SetCoord(xmin, ymax, zmax);
            aLimPntArray[7].SetCoord(xmax, ymax, zmax);

            for (i = 0; i <= 7; i++) {
              aParOnLin = ElCLib::Parameter(aLin, aLimPntArray[i]);
              tmin = Min(aParOnLin, tmin);
              tmax = Max(aParOnLin, tmax);
            }
            cfirst = Max(cfirst, tmin);
            clast = Min(clast, tmax);
          }

          if (myS->IsUPeriodic())
            NbU = 13;
          if (myS->IsVPeriodic())
            NbV = 13;

          if (clast - cfirst <= Precision::Confusion())
          {
            Standard_Real aCPar = (cfirst + clast) / 2.;
            gp_Pnt aPm = C.Value(aCPar);
            Extrema_ExtPS anExtPS(aPm, *myS, ufirst, ulast,
              vfirst, vlast, mytolS, mytolS, Extrema_ExtFlag_MIN);
            myDone = anExtPS.IsDone();
            if (myDone) {
              Standard_Integer NbExt = anExtPS.NbExt();
              Standard_Real T = aCPar, U, V;
              Extrema_POnCurv PC;
              Extrema_POnSurf PS;
              for (i = 1; i <= NbExt; i++) {
                PS = anExtPS.Point(i);
                PS.Parameter(U, V);
                AddSolution(C, T, U, V, PC.Value(), PS.Value(), anExtPS.SquareDistance(i));
              }
            }
            return;
          }

          Extrema_GenExtCS Ext(C, *myS, NbT, NbU, NbV, cfirst, clast, ufirst, ulast,
            vfirst, vlast, mytolC, mytolS);

          myDone = Ext.IsDone();
          if (myDone) {
            Standard_Integer NbExt = Ext.NbExt();
            Standard_Real T,U,V;
            Extrema_POnCurv PC;
            Extrema_POnSurf PS;
            for (i = 1; i <= NbExt; i++) {
              PC = Ext.PointOnCurve(i);
              PS = Ext.PointOnSurface(i);
              T = PC.Parameter();
              PS.Parameter(U, V);
              AddSolution(C, T, U, V, PC.Value(), PS.Value(), Ext.SquareDistance(i));
            }
          }
          return;
        }
      }
      break;
    }
    //  Modified by skv - Thu Jul  7 12:29:34 2005 OCC9134 Begin
  case GeomAbs_Circle:
    {
      if(myStype == GeomAbs_Cylinder) {
        myExtElCS.Perform(C.Circle(), myS->Cylinder());
        break;
      }
      else if(myStype == GeomAbs_Plane)
      {
        myExtElCS.Perform(C.Circle(), myS->Plane());
        break;
      }
      else if (myStype == GeomAbs_Sphere)
      {
        myExtElCS.Perform(C.Circle(), myS->Sphere());
        break;
      }
    }
    Standard_FALLTHROUGH
  case GeomAbs_Hyperbola: 
    {
      if(myCtype == GeomAbs_Hyperbola && myStype == GeomAbs_Plane) {
        //  Modified by skv - Thu Jul  7 12:29:34 2005 OCC9134 End
        myExtElCS.Perform(C.Hyperbola(), myS->Plane());
        break;
      }
    }
    Standard_FALLTHROUGH
  default:
    {
      isComputeAnalytic = Standard_False;
      break;
    }
  }

  if (isComputeAnalytic)
  {
    if (myExtElCS.IsDone())
    {
      myDone = Standard_True;
      myIsPar = myExtElCS.IsParallel();
      if (myIsPar)
      {
        mySqDist.Append(myExtElCS.SquareDistance(1));
      }
      else
      {
        Standard_Integer NbExt = myExtElCS.NbExt();
        for (i = 1; i <= NbExt; i++)
        {
          Extrema_POnCurv PC;
          Extrema_POnSurf PS;
          myExtElCS.Points(i, PC, PS);
          Standard_Real Ucurve = PC.Parameter();
          Standard_Real U, V;
          PS.Parameter(U, V);
          AddSolution(C, Ucurve, U, V, PC.Value(), PS.Value(), myExtElCS.SquareDistance(i));
        }

        if (mySqDist.Length() == 0 && NbExt > 0)
        {
          // Analytical extrema seem to be out of curve/surface boundaries.
          // Try extremity points of curve.
          gp_Pnt aPOnC[2], aPOnS[2];
          Standard_Real aT[2] = { myucinf, myucsup }, U[2], V[2];
          Standard_Real aDist[2] = { -1, -1 };
          for (i = 0; i < 2; ++i)
          {
            if (Precision::IsInfinite(aT[i]))
              continue;

            aPOnC[i] = C.Value(aT[i]);
            switch (myStype)
            {
              case GeomAbs_Plane:
              {
                ElSLib::Parameters(myS->Plane(), aPOnC[i], U[i], V[i]);
                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Plane());
                break;
              }
              case GeomAbs_Sphere:
              {
                ElSLib::Parameters(myS->Sphere(), aPOnC[i], U[i], V[i]);
                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Sphere());
                break;
              }
              case GeomAbs_Cylinder:
              {
                ElSLib::Parameters(myS->Cylinder(), aPOnC[i], U[i], V[i]);
                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Cylinder());
                break;
              }
              case GeomAbs_Torus:
              {
                ElSLib::Parameters(myS->Torus(), aPOnC[i], U[i], V[i]);
                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Torus());
                break;
              }
              case GeomAbs_Cone:
              {
                ElSLib::Parameters(myS->Cone(), aPOnC[i], U[i], V[i]);
                aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Cone());
                break;
              }
              default:
                continue;
            }

            aDist[i] = aPOnC[i].SquareDistance(aPOnS[i]);
          }

          Standard_Boolean bAdd[2] = {Standard_False, Standard_False};

          // Choose solution to add
          if (aDist[0] >= 0. && aDist[1] >= 0.)
          {
            Standard_Real aDiff = aDist[0] - aDist[1];
            // Both computed -> take only minimal
            if (Abs(aDiff) < Precision::Confusion())
              // Add both
              bAdd[0] = bAdd[1] = Standard_True;
            else if (aDiff < 0)
              // Add first
              bAdd[0] = Standard_True;
            else
              // Add second
              bAdd[1] = Standard_True;
          }
          else if (aDist[0] >= 0.)
            // Add first
            bAdd[0] = Standard_True;
          else if (aDist[1] >= 0.)
            // Add second
            bAdd[1] = Standard_True;

          for (i = 0; i < 2; ++i)
          {
            if (bAdd[i])
              AddSolution(C, aT[i], U[i], V[i], aPOnC[i], aPOnS[i], aDist[i]);
          }
        }
      }
      return;
    }
  }

  // Elementary extrema is not done, try generic solution
  Extrema_GenExtCS Ext;
  Ext.Initialize(*myS, NbU, NbV, mytolS);
  if (myCtype == GeomAbs_Hyperbola) {
    Standard_Real tmin = Max(-20., C.FirstParameter());
    Standard_Real tmax = Min(20., C.LastParameter());
    Ext.Perform(C, NbT, tmin, tmax, mytolC); // to avoid overflow
  }
  else {
    if ((myCtype == GeomAbs_Circle       && NbT < 13) ||
      (myCtype == GeomAbs_BSplineCurve && NbT < 13))
    {
      NbT = 13;
    }
    Ext.Perform(C, NbT, mytolC);
  }

  myDone = Ext.IsDone();
  if (myDone) {
    Standard_Integer NbExt = Ext.NbExt();
    Standard_Real T, U, V;
    Extrema_POnCurv PC;
    Extrema_POnSurf PS;
    for (i = 1; i <= NbExt; i++) {
      PC = Ext.PointOnCurve(i);
      PS = Ext.PointOnSurface(i);
      T = PC.Parameter();
      PS.Parameter(U, V);
      AddSolution(C, T, U, V, PC.Value(), PS.Value(), Ext.SquareDistance(i));
    }

    //Add sharp points
    Standard_Integer SolNumber = mySqDist.Length();
    Standard_Address CopyC = (Standard_Address)&C;
    Adaptor3d_Curve& aC = *(Adaptor3d_Curve*)CopyC;
    Standard_Integer NbIntervals = aC.NbIntervals(GeomAbs_C1);
    TColStd_Array1OfReal SharpPoints(1, NbIntervals + 1);
    aC.Intervals(SharpPoints, GeomAbs_C1);

    Extrema_ExtPS aProjPS;
    aProjPS.Initialize(*myS,
      myS->FirstUParameter(),
      myS->LastUParameter(),
      myS->FirstVParameter(),
      myS->LastVParameter(),
      mytolS,
      mytolS);

    for (i = 2; i < SharpPoints.Upper(); ++i)
    {
      T = SharpPoints(i);
      gp_Pnt aPnt = C.Value(T);
      aProjPS.Perform(aPnt);
      if (!aProjPS.IsDone())
        continue;
      Standard_Integer NbProj = aProjPS.NbExt(), jmin = 0;
      Standard_Real MinSqDist = RealLast();
      for (j = 1; j <= NbProj; j++)
      {
        Standard_Real aSqDist = aProjPS.SquareDistance(j);
        if (aSqDist < MinSqDist)
        {
          MinSqDist = aSqDist;
          jmin = j;
        }
      }
      if (jmin != 0)
      {
        aProjPS.Point(jmin).Parameter(U, V);
        AddSolution(C, T, U, V,
          aPnt, aProjPS.Point(jmin).Value(), MinSqDist);
      }
    }
    //Cut sharp solutions to keep only minimum and maximum
    Standard_Integer imin = SolNumber + 1, imax = mySqDist.Length();
    for (i = SolNumber + 1; i <= mySqDist.Length(); i++)
    {
      if (mySqDist(i) < mySqDist(imin))
        imin = i;
      if (mySqDist(i) > mySqDist(imax))
        imax = i;
    }
    if (mySqDist.Length() > SolNumber + 2)
    {
      Standard_Real MinSqDist = mySqDist(imin);
      Standard_Real MaxSqDist = mySqDist(imax);
      Extrema_POnCurv MinPC = myPOnC(imin);
      Extrema_POnCurv MaxPC = myPOnC(imax);
      Extrema_POnSurf MinPS = myPOnS(imin);
      Extrema_POnSurf MaxPS = myPOnS(imax);

      mySqDist.Remove(SolNumber + 1, mySqDist.Length());
      myPOnC.Remove(SolNumber + 1, myPOnC.Length());
      myPOnS.Remove(SolNumber + 1, myPOnS.Length());

      mySqDist.Append(MinSqDist);
      myPOnC.Append(MinPC);
      myPOnS.Append(MinPS);
      mySqDist.Append(MaxSqDist);
      myPOnC.Append(MaxPC);
      myPOnS.Append(MaxPS);
    }
  }
}


Standard_Boolean Extrema_ExtCS::IsDone() const
{
  return myDone;
}

Standard_Boolean Extrema_ExtCS::IsParallel() const
{
  if (!IsDone())
  {
    throw StdFail_NotDone();
  }

  return myIsPar;
}


Standard_Real Extrema_ExtCS::SquareDistance(const Standard_Integer N) const
{
  if (N < 1 || N > NbExt())
  {
    throw Standard_OutOfRange();
  }

  return mySqDist.Value(N);
}


Standard_Integer Extrema_ExtCS::NbExt() const
{
  if (!IsDone())
  {
    throw StdFail_NotDone();
  }

  return mySqDist.Length();
}


void Extrema_ExtCS::Points(const Standard_Integer N,
                           Extrema_POnCurv&       P1,
                           Extrema_POnSurf&       P2) const
{
  if (N < 1 || N > NbExt())
  {
    throw Standard_OutOfRange();
  }

  P1 = myPOnC.Value(N);
  P2 = myPOnS.Value(N);
}

Standard_Boolean Extrema_ExtCS::AddSolution(const Adaptor3d_Curve& theCurve,
  const Standard_Real aT,
  const Standard_Real aU,
  const Standard_Real aV,
  const gp_Pnt& PointOnCurve,
  const gp_Pnt& PointOnSurf,
  const Standard_Real SquareDist)
{
  Standard_Boolean Added = Standard_False;

  Standard_Real T = aT, U = aU, V = aV;

  if (theCurve.IsPeriodic())
    T = ElCLib::InPeriod(T, myucinf, myucinf + theCurve.Period());
  if (myS->IsUPeriodic())
    U = ElCLib::InPeriod(U, myuinf, myuinf + myS->UPeriod());
  if (myS->IsVPeriodic())
    V = ElCLib::InPeriod(V, myvinf, myvinf + myS->VPeriod());

  Extrema_POnCurv aPC;
  Extrema_POnSurf aPS;
  if ((myucinf-T) <= mytolC && (T-myucsup) <= mytolC &&
    (myuinf-U) <= mytolS && (U-myusup) <= mytolS &&
    (myvinf-V) <= mytolS && (V-myvsup) <= mytolS)
  {
    Standard_Boolean IsNewSolution = Standard_True;
    for (Standard_Integer j = 1; j <= mySqDist.Length(); j++)
    {
      aPC = myPOnC(j);
      aPS = myPOnS(j);
      Standard_Real Tj = aPC.Parameter();
      Standard_Real Uj, Vj;
      aPS.Parameter(Uj, Vj);
      if (Abs(T - Tj) <= mytolC &&
        Abs(U - Uj) <= mytolS &&
        Abs(V - Vj) <= mytolS)
      {
        IsNewSolution = Standard_False;
        break;
      }
    }
    if (IsNewSolution)
    {
      mySqDist.Append(SquareDist);
      aPC.SetValues(T, PointOnCurve);
      myPOnC.Append(aPC);
      myPOnS.Append(Extrema_POnSurf(U, V, PointOnSurf));
      Added = Standard_True;
    }
  }
  return Added;
}
Commit	Line	Data
b311480e	1	// Created on: 1995-07-18
	2	// Created by: Modelistation
	3	// Copyright (c) 1995-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
	17	// Modified by skv - Thu Jul 7 12:29:34 2005 OCC9134
	18
42cf5bc1	19	#include <Adaptor3d_Curve.hxx>
	20	#include <Adaptor3d_Surface.hxx>
	21	#include <Bnd_Box.hxx>
	22	#include <BndLib_AddSurface.hxx>
	23	#include <ElCLib.hxx>
	24	#include <ElSLib.hxx>
	25	#include <Extrema_ExtCS.hxx>
42cf5bc1	26	#include <Extrema_ExtPS.hxx>
7fd59977	27	#include <Extrema_GenExtCS.hxx>
42cf5bc1	28	#include <Extrema_POnCurv.hxx>
42cf5bc1	29	#include <Extrema_POnSurf.hxx>
7fd59977	30	#include <GeomAbs_CurveType.hxx>
7fd59977	31	#include <gp_Cone.hxx>
42cf5bc1	32	#include <gp_Cylinder.hxx>
	33	#include <gp_Lin.hxx>
	34	#include <gp_Pln.hxx>
	35	#include <gp_Pnt.hxx>
7fd59977	36	#include <gp_Sphere.hxx>
7fd59977	37	#include <gp_Torus.hxx>
42cf5bc1	38	#include <Precision.hxx>
42cf5bc1	39	#include <Standard_OutOfRange.hxx>
42cf5bc1	40	#include <StdFail_NotDone.hxx>
38308958	41	#include <TColStd_Array1OfReal.hxx>
e8e8b273	42	#include <Extrema_ExtPS.hxx>
7fd59977	43
d533dafb	44	Extrema_ExtCS::Extrema_ExtCS()
	45	: myS(NULL),
	46	myDone(Standard_False),
	47	myIsPar(Standard_False),
	48	myuinf(0.0),
	49	myusup(0.0),
	50	myvinf(0.0),
	51	myvsup(0.0),
	52	mytolC(0.0),
	53	mytolS(0.0),
	54	myucinf(0.0),
	55	myucsup(0.0),
	56	myStype(GeomAbs_OtherSurface)
7fd59977	57	{
7fd59977	58	}
	59
	60	Extrema_ExtCS::Extrema_ExtCS(const Adaptor3d_Curve& C,
29d778bf	61	const Adaptor3d_Surface& S,
	62	const Standard_Real TolC,
	63	const Standard_Real TolS)
7fd59977	64
	65	{
	66	Initialize(S, S.FirstUParameter(), S.LastUParameter(),
29d778bf	67	S.FirstVParameter(), S.LastVParameter(),
29d778bf	68	TolC, TolS);
7fd59977	69	Perform(C, C.FirstParameter(), C.LastParameter());
	70	}
	71
	72	Extrema_ExtCS::Extrema_ExtCS(const Adaptor3d_Curve& C,
29d778bf	73	const Adaptor3d_Surface& S,
	74	const Standard_Real UCinf,
	75	const Standard_Real UCsup,
	76	const Standard_Real Uinf,
	77	const Standard_Real Usup,
	78	const Standard_Real Vinf,
	79	const Standard_Real Vsup,
	80	const Standard_Real TolC,
	81	const Standard_Real TolS)
7fd59977	82
	83	{
	84	Initialize(S, Uinf, Usup, Vinf, Vsup, TolC, TolS);
	85	Perform(C, UCinf, UCsup);
	86	}
	87
	88
	89	void Extrema_ExtCS::Initialize(const Adaptor3d_Surface& S,
29d778bf	90	const Standard_Real Uinf,
	91	const Standard_Real Usup,
	92	const Standard_Real Vinf,
	93	const Standard_Real Vsup,
	94	const Standard_Real TolC,
	95	const Standard_Real TolS)
7fd59977	96	{
	97	myS = (Adaptor3d_SurfacePtr)&S;
	98	myIsPar = Standard_False;
	99	myuinf = Uinf;
	100	myusup = Usup;
	101	myvinf = Vinf;
	102	myvsup = Vsup;
	103	mytolC = TolC;
	104	mytolS = TolS;
	105	myStype = myS->GetType();
	106	}
	107
29d778bf	108
7fd59977	109	void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C,
29d778bf	110	const Standard_Real Uinf,
29d778bf	111	const Standard_Real Usup)
7fd59977	112	{
	113	myucinf = Uinf;
	114	myucsup = Usup;
	115	myPOnS.Clear();
	116	myPOnC.Clear();
	117	mySqDist.Clear();
38308958	118	Standard_Integer i, j;
7fd59977	119	Standard_Integer NbT, NbU, NbV;
79aa9b5c	120	NbT = 12; NbU = NbV = 10;
7fd59977	121	GeomAbs_CurveType myCtype = C.GetType();
7fd59977	122
03cca6f7	123	myDone = Standard_False;
	124	// Try analytic computation of extrema
	125	Standard_Boolean isComputeAnalytic = Standard_True;
7fd59977	126
	127	switch(myCtype) {
	128
	129	case GeomAbs_Line:
	130	{
29d778bf	131
7fd59977	132	switch(myStype) {
7fd59977	133	case GeomAbs_Sphere:
29d778bf	134	myExtElCS.Perform(C.Line(), myS->Sphere());
29d778bf	135	break;
7fd59977	136	case GeomAbs_Cylinder:
29d778bf	137	myExtElCS.Perform(C.Line(), myS->Cylinder());
29d778bf	138	break;
7fd59977	139	case GeomAbs_Plane:
29d778bf	140	myExtElCS.Perform(C.Line(), myS->Plane());
29d778bf	141	if (myExtElCS.IsParallel()) break;
b1811c1d	142	Standard_FALLTHROUGH
7fd59977	143
	144	case GeomAbs_Torus:
	145	case GeomAbs_Cone:
	146	case GeomAbs_BezierSurface:
	147	case GeomAbs_BSplineSurface:
	148	case GeomAbs_SurfaceOfRevolution:
	149	case GeomAbs_SurfaceOfExtrusion:
1896126e	150	case GeomAbs_OffsetSurface:
7fd59977	151	case GeomAbs_OtherSurface:
29d778bf	152	{
	153	Standard_Real cfirst = myucinf, clast = myucsup;
	154	Standard_Real ufirst = myS->FirstUParameter(), ulast = myS->LastUParameter(),
	155	vfirst = myS->FirstVParameter(), vlast = myS->LastVParameter();
7fd59977	156
4563cf3e	157	if (!(Precision::IsInfinite(ufirst) \|\| Precision::IsInfinite(ulast) \|\|
	158	Precision::IsInfinite(vfirst) \|\| Precision::IsInfinite(vlast)))
	159	{
	160	Standard_Real tmin = Precision::Infinite(), tmax = -tmin;
	161	Standard_Real xmin, ymin, zmin, xmax, ymax, zmax;
29d778bf	162	Bnd_Box aSurfBox;
1896126e	163	BndLib_AddSurface::Add(*myS, ufirst, ulast, vfirst, vlast, Precision::Confusion(), aSurfBox);
29d778bf	164	aSurfBox.Get(xmin, ymin, zmin, xmax, ymax, zmax);
29d778bf	165	gp_Lin aLin = C.Line();
4563cf3e	166	Standard_Real aParOnLin;
	167	gp_Pnt aLimPntArray[8];
	168
	169	aLimPntArray[0].SetCoord(xmin, ymin, zmin);
	170	aLimPntArray[1].SetCoord(xmax, ymin, zmin);
	171	aLimPntArray[2].SetCoord(xmin, ymax, zmin);
	172	aLimPntArray[3].SetCoord(xmax, ymax, zmin);
	173	aLimPntArray[4].SetCoord(xmin, ymin, zmax);
	174	aLimPntArray[5].SetCoord(xmax, ymin, zmax);
	175	aLimPntArray[6].SetCoord(xmin, ymax, zmax);
	176	aLimPntArray[7].SetCoord(xmax, ymax, zmax);
	177
	178	for (i = 0; i <= 7; i++) {
	179	aParOnLin = ElCLib::Parameter(aLin, aLimPntArray[i]);
	180	tmin = Min(aParOnLin, tmin);
	181	tmax = Max(aParOnLin, tmax);
29d778bf	182	}
29d778bf	183	cfirst = Max(cfirst, tmin);
4563cf3e	184	clast = Min(clast, tmax);
29d778bf	185	}
29d778bf	186
fa6d1712	187	if (myS->IsUPeriodic())
	188	NbU = 13;
	189	if (myS->IsVPeriodic())
	190	NbV = 13;
29d778bf	191
e8e8b273	192	if (clast - cfirst <= Precision::Confusion())
	193	{
	194	Standard_Real aCPar = (cfirst + clast) / 2.;
	195	gp_Pnt aPm = C.Value(aCPar);
	196	Extrema_ExtPS anExtPS(aPm, *myS, ufirst, ulast,
	197	vfirst, vlast, mytolS, mytolS, Extrema_ExtFlag_MIN);
	198	myDone = anExtPS.IsDone();
	199	if (myDone) {
	200	Standard_Integer NbExt = anExtPS.NbExt();
	201	Standard_Real T = aCPar, U, V;
	202	Extrema_POnCurv PC;
	203	Extrema_POnSurf PS;
	204	for (i = 1; i <= NbExt; i++) {
	205	PS = anExtPS.Point(i);
	206	PS.Parameter(U, V);
	207	AddSolution(C, T, U, V, PC.Value(), PS.Value(), anExtPS.SquareDistance(i));
	208	}
	209	}
	210	return;
	211	}
	212
29d778bf	213	Extrema_GenExtCS Ext(C, *myS, NbT, NbU, NbV, cfirst, clast, ufirst, ulast,
	214	vfirst, vlast, mytolC, mytolS);
	215
	216	myDone = Ext.IsDone();
	217	if (myDone) {
	218	Standard_Integer NbExt = Ext.NbExt();
	219	Standard_Real T,U,V;
	220	Extrema_POnCurv PC;
	221	Extrema_POnSurf PS;
	222	for (i = 1; i <= NbExt; i++) {
	223	PC = Ext.PointOnCurve(i);
	224	PS = Ext.PointOnSurface(i);
	225	T = PC.Parameter();
	226	PS.Parameter(U, V);
38308958	227	AddSolution(C, T, U, V, PC.Value(), PS.Value(), Ext.SquareDistance(i));
29d778bf	228	}
	229	}
	230	return;
29d778bf	231	}
7fd59977	232	}
	233	break;
	234	}
29d778bf	235	// Modified by skv - Thu Jul 7 12:29:34 2005 OCC9134 Begin
7fd59977	236	case GeomAbs_Circle:
	237	{
	238	if(myStype == GeomAbs_Cylinder) {
29d778bf	239	myExtElCS.Perform(C.Circle(), myS->Cylinder());
29d778bf	240	break;
7fd59977	241	}
9dfbbfe6	242	else if(myStype == GeomAbs_Plane)
	243	{
	244	myExtElCS.Perform(C.Circle(), myS->Plane());
	245	break;
	246	}
03cca6f7	247	else if (myStype == GeomAbs_Sphere)
	248	{
	249	myExtElCS.Perform(C.Circle(), myS->Sphere());
	250	break;
	251	}
7fd59977	252	}
b1811c1d	253	Standard_FALLTHROUGH
7fd59977	254	case GeomAbs_Hyperbola:
	255	{
	256	if(myCtype == GeomAbs_Hyperbola && myStype == GeomAbs_Plane) {
29d778bf	257	// Modified by skv - Thu Jul 7 12:29:34 2005 OCC9134 End
	258	myExtElCS.Perform(C.Hyperbola(), myS->Plane());
	259	break;
7fd59977	260	}
7fd59977	261	}
b1811c1d	262	Standard_FALLTHROUGH
7fd59977	263	default:
7fd59977	264	{
03cca6f7	265	isComputeAnalytic = Standard_False;
	266	break;
	267	}
	268	}
	269
	270	if (isComputeAnalytic)
	271	{
	272	if (myExtElCS.IsDone())
	273	{
	274	myDone = Standard_True;
	275	myIsPar = myExtElCS.IsParallel();
	276	if (myIsPar)
	277	{
	278	mySqDist.Append(myExtElCS.SquareDistance(1));
7fd59977	279	}
03cca6f7	280	else
	281	{
	282	Standard_Integer NbExt = myExtElCS.NbExt();
	283	for (i = 1; i <= NbExt; i++)
fa6d1712	284	{
03cca6f7	285	Extrema_POnCurv PC;
	286	Extrema_POnSurf PS;
	287	myExtElCS.Points(i, PC, PS);
	288	Standard_Real Ucurve = PC.Parameter();
	289	Standard_Real U, V;
29d778bf	290	PS.Parameter(U, V);
03cca6f7	291	AddSolution(C, Ucurve, U, V, PC.Value(), PS.Value(), myExtElCS.SquareDistance(i));
29d778bf	292	}
38308958	293
03cca6f7	294	if (mySqDist.Length() == 0 && NbExt > 0)
38308958	295	{
03cca6f7	296	// Analytical extrema seem to be out of curve/surface boundaries.
	297	// Try extremity points of curve.
	298	gp_Pnt aPOnC[2], aPOnS[2];
	299	Standard_Real aT[2] = { myucinf, myucsup }, U[2], V[2];
	300	Standard_Real aDist[2] = { -1, -1 };
	301	for (i = 0; i < 2; ++i)
38308958	302	{
03cca6f7	303	if (Precision::IsInfinite(aT[i]))
	304	continue;
	305
	306	aPOnC[i] = C.Value(aT[i]);
	307	switch (myStype)
38308958	308	{
03cca6f7	309	case GeomAbs_Plane:
	310	{
	311	ElSLib::Parameters(myS->Plane(), aPOnC[i], U[i], V[i]);
	312	aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Plane());
	313	break;
	314	}
	315	case GeomAbs_Sphere:
	316	{
	317	ElSLib::Parameters(myS->Sphere(), aPOnC[i], U[i], V[i]);
	318	aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Sphere());
	319	break;
	320	}
	321	case GeomAbs_Cylinder:
	322	{
	323	ElSLib::Parameters(myS->Cylinder(), aPOnC[i], U[i], V[i]);
	324	aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Cylinder());
	325	break;
	326	}
	327	case GeomAbs_Torus:
	328	{
	329	ElSLib::Parameters(myS->Torus(), aPOnC[i], U[i], V[i]);
	330	aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Torus());
	331	break;
	332	}
	333	case GeomAbs_Cone:
	334	{
	335	ElSLib::Parameters(myS->Cone(), aPOnC[i], U[i], V[i]);
	336	aPOnS[i] = ElSLib::Value(U[i], V[i], myS->Cone());
	337	break;
	338	}
	339	default:
	340	continue;
38308958	341	}
03cca6f7	342
03cca6f7	343	aDist[i] = aPOnC[i].SquareDistance(aPOnS[i]);
38308958	344	}
03cca6f7	345
	346	Standard_Boolean bAdd[2] = {Standard_False, Standard_False};
	347
	348	// Choose solution to add
	349	if (aDist[0] >= 0. && aDist[1] >= 0.)
38308958	350	{
03cca6f7	351	Standard_Real aDiff = aDist[0] - aDist[1];
	352	// Both computed -> take only minimal
	353	if (Abs(aDiff) < Precision::Confusion())
	354	// Add both
	355	bAdd[0] = bAdd[1] = Standard_True;
	356	else if (aDiff < 0)
	357	// Add first
	358	bAdd[0] = Standard_True;
	359	else
	360	// Add second
	361	bAdd[1] = Standard_True;
	362	}
	363	else if (aDist[0] >= 0.)
	364	// Add first
	365	bAdd[0] = Standard_True;
	366	else if (aDist[1] >= 0.)
	367	// Add second
	368	bAdd[1] = Standard_True;
	369
	370	for (i = 0; i < 2; ++i)
	371	{
	372	if (bAdd[i])
	373	AddSolution(C, aT[i], U[i], V[i], aPOnC[i], aPOnS[i], aDist[i]);
38308958	374	}
38308958	375	}
7fd59977	376	}
	377	return;
	378	}
7fd59977	379	}
29d778bf	380
03cca6f7	381	// Elementary extrema is not done, try generic solution
	382	Extrema_GenExtCS Ext;
	383	Ext.Initialize(*myS, NbU, NbV, mytolS);
	384	if (myCtype == GeomAbs_Hyperbola) {
	385	Standard_Real tmin = Max(-20., C.FirstParameter());
	386	Standard_Real tmax = Min(20., C.LastParameter());
	387	Ext.Perform(C, NbT, tmin, tmax, mytolC); // to avoid overflow
	388	}
	389	else {
	390	if ((myCtype == GeomAbs_Circle && NbT < 13) \|\|
	391	(myCtype == GeomAbs_BSplineCurve && NbT < 13))
	392	{
	393	NbT = 13;
	394	}
	395	Ext.Perform(C, NbT, mytolC);
	396	}
	397
	398	myDone = Ext.IsDone();
7fd59977	399	if (myDone) {
03cca6f7	400	Standard_Integer NbExt = Ext.NbExt();
	401	Standard_Real T, U, V;
	402	Extrema_POnCurv PC;
	403	Extrema_POnSurf PS;
	404	for (i = 1; i <= NbExt; i++) {
	405	PC = Ext.PointOnCurve(i);
	406	PS = Ext.PointOnSurface(i);
	407	T = PC.Parameter();
	408	PS.Parameter(U, V);
	409	AddSolution(C, T, U, V, PC.Value(), PS.Value(), Ext.SquareDistance(i));
7fd59977	410	}
03cca6f7	411
	412	//Add sharp points
	413	Standard_Integer SolNumber = mySqDist.Length();
	414	Standard_Address CopyC = (Standard_Address)&C;
	415	Adaptor3d_Curve& aC = (Adaptor3d_Curve)CopyC;
	416	Standard_Integer NbIntervals = aC.NbIntervals(GeomAbs_C1);
	417	TColStd_Array1OfReal SharpPoints(1, NbIntervals + 1);
	418	aC.Intervals(SharpPoints, GeomAbs_C1);
	419
	420	Extrema_ExtPS aProjPS;
	421	aProjPS.Initialize(*myS,
	422	myS->FirstUParameter(),
	423	myS->LastUParameter(),
	424	myS->FirstVParameter(),
	425	myS->LastVParameter(),
	426	mytolS,
	427	mytolS);
	428
	429	for (i = 2; i < SharpPoints.Upper(); ++i)
	430	{
	431	T = SharpPoints(i);
	432	gp_Pnt aPnt = C.Value(T);
	433	aProjPS.Perform(aPnt);
	434	if (!aProjPS.IsDone())
	435	continue;
	436	Standard_Integer NbProj = aProjPS.NbExt(), jmin = 0;
	437	Standard_Real MinSqDist = RealLast();
	438	for (j = 1; j <= NbProj; j++)
9dfbbfe6	439	{
03cca6f7	440	Standard_Real aSqDist = aProjPS.SquareDistance(j);
03cca6f7	441	if (aSqDist < MinSqDist)
9dfbbfe6	442	{
03cca6f7	443	MinSqDist = aSqDist;
03cca6f7	444	jmin = j;
9dfbbfe6	445	}
9dfbbfe6	446	}
03cca6f7	447	if (jmin != 0)
	448	{
	449	aProjPS.Point(jmin).Parameter(U, V);
	450	AddSolution(C, T, U, V,
	451	aPnt, aProjPS.Point(jmin).Value(), MinSqDist);
	452	}
	453	}
	454	//Cut sharp solutions to keep only minimum and maximum
	455	Standard_Integer imin = SolNumber + 1, imax = mySqDist.Length();
	456	for (i = SolNumber + 1; i <= mySqDist.Length(); i++)
	457	{
	458	if (mySqDist(i) < mySqDist(imin))
	459	imin = i;
	460	if (mySqDist(i) > mySqDist(imax))
	461	imax = i;
	462	}
	463	if (mySqDist.Length() > SolNumber + 2)
	464	{
	465	Standard_Real MinSqDist = mySqDist(imin);
	466	Standard_Real MaxSqDist = mySqDist(imax);
	467	Extrema_POnCurv MinPC = myPOnC(imin);
	468	Extrema_POnCurv MaxPC = myPOnC(imax);
	469	Extrema_POnSurf MinPS = myPOnS(imin);
	470	Extrema_POnSurf MaxPS = myPOnS(imax);
	471
	472	mySqDist.Remove(SolNumber + 1, mySqDist.Length());
	473	myPOnC.Remove(SolNumber + 1, myPOnC.Length());
	474	myPOnS.Remove(SolNumber + 1, myPOnS.Length());
	475
	476	mySqDist.Append(MinSqDist);
	477	myPOnC.Append(MinPC);
	478	myPOnS.Append(MinPS);
	479	mySqDist.Append(MaxSqDist);
	480	myPOnC.Append(MaxPC);
	481	myPOnS.Append(MaxPS);
7fd59977	482	}
7fd59977	483	}
7fd59977	484	}
	485
	486
	487	Standard_Boolean Extrema_ExtCS::IsDone() const
	488	{
	489	return myDone;
	490	}
	491
	492	Standard_Boolean Extrema_ExtCS::IsParallel() const
	493	{
638ad7f3	494	if (!IsDone())
	495	{
	496	throw StdFail_NotDone();
	497	}
	498
7fd59977	499	return myIsPar;
	500	}
	501
	502
	503	Standard_Real Extrema_ExtCS::SquareDistance(const Standard_Integer N) const
	504	{
638ad7f3	505	if (N < 1 \|\| N > NbExt())
	506	{
	507	throw Standard_OutOfRange();
	508	}
	509
7fd59977	510	return mySqDist.Value(N);
	511	}
	512
	513
	514	Standard_Integer Extrema_ExtCS::NbExt() const
	515	{
638ad7f3	516	if (!IsDone())
	517	{
	518	throw StdFail_NotDone();
	519	}
	520
	521	return mySqDist.Length();
7fd59977	522	}
	523
	524
	525
	526	void Extrema_ExtCS::Points(const Standard_Integer N,
638ad7f3	527	Extrema_POnCurv& P1,
638ad7f3	528	Extrema_POnSurf& P2) const
7fd59977	529	{
638ad7f3	530	if (N < 1 \|\| N > NbExt())
	531	{
	532	throw Standard_OutOfRange();
	533	}
	534
7fd59977	535	P1 = myPOnC.Value(N);
	536	P2 = myPOnS.Value(N);
	537	}
38308958	538
38308958	539	Standard_Boolean Extrema_ExtCS::AddSolution(const Adaptor3d_Curve& theCurve,
29d778bf	540	const Standard_Real aT,
	541	const Standard_Real aU,
	542	const Standard_Real aV,
	543	const gp_Pnt& PointOnCurve,
	544	const gp_Pnt& PointOnSurf,
	545	const Standard_Real SquareDist)
38308958	546	{
	547	Standard_Boolean Added = Standard_False;
	548
	549	Standard_Real T = aT, U = aU, V = aV;
29d778bf	550
38308958	551	if (theCurve.IsPeriodic())
	552	T = ElCLib::InPeriod(T, myucinf, myucinf + theCurve.Period());
	553	if (myS->IsUPeriodic())
	554	U = ElCLib::InPeriod(U, myuinf, myuinf + myS->UPeriod());
	555	if (myS->IsVPeriodic())
	556	V = ElCLib::InPeriod(V, myvinf, myvinf + myS->VPeriod());
	557
	558	Extrema_POnCurv aPC;
	559	Extrema_POnSurf aPS;
	560	if ((myucinf-T) <= mytolC && (T-myucsup) <= mytolC &&
29d778bf	561	(myuinf-U) <= mytolS && (U-myusup) <= mytolS &&
29d778bf	562	(myvinf-V) <= mytolS && (V-myvsup) <= mytolS)
38308958	563	{
	564	Standard_Boolean IsNewSolution = Standard_True;
	565	for (Standard_Integer j = 1; j <= mySqDist.Length(); j++)
	566	{
	567	aPC = myPOnC(j);
	568	aPS = myPOnS(j);
	569	Standard_Real Tj = aPC.Parameter();
	570	Standard_Real Uj, Vj;
	571	aPS.Parameter(Uj, Vj);
	572	if (Abs(T - Tj) <= mytolC &&
29d778bf	573	Abs(U - Uj) <= mytolS &&
29d778bf	574	Abs(V - Vj) <= mytolS)
38308958	575	{
	576	IsNewSolution = Standard_False;
	577	break;
	578	}
	579	}
	580	if (IsNewSolution)
	581	{
	582	mySqDist.Append(SquareDist);
	583	aPC.SetValues(T, PointOnCurve);
	584	myPOnC.Append(aPC);
	585	myPOnS.Append(Extrema_POnSurf(U, V, PointOnSurf));
	586	Added = Standard_True;
	587	}
	588	}
	589	return Added;
	590	}