// Created on: 1992-05-07
// Created by: Jacques GOUSSARD
// Copyright (c) 1992-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 <IntPatch_WLine.hxx>

static 
  Standard_Integer SetQuad(const Handle(Adaptor3d_HSurface)& theS,
                           GeomAbs_SurfaceType& theTS,
                           IntSurf_Quadric& theQuad);

//=======================================================================
//function : IntPatch_ImpImpIntersection
//purpose  : 
//=======================================================================
IntPatch_ImpImpIntersection::IntPatch_ImpImpIntersection ():
myDone(IntStatus_Fail)
{
}
//=======================================================================
//function : IntPatch_ImpImpIntersection
//purpose  : 
//=======================================================================
IntPatch_ImpImpIntersection::IntPatch_ImpImpIntersection
       (const Handle(Adaptor3d_HSurface)&  S1,
        const Handle(Adaptor3d_TopolTool)& D1,
        const Handle(Adaptor3d_HSurface)&  S2,
        const Handle(Adaptor3d_TopolTool)& D2,
        const Standard_Real TolArc,
        const Standard_Real TolTang,
        const Standard_Boolean theIsReqToKeepRLine)
{
  Perform(S1,D1,S2,D2,TolArc,TolTang, theIsReqToKeepRLine);
}
//=======================================================================
//function : Perform
//purpose  : 
//=======================================================================
void IntPatch_ImpImpIntersection::Perform(const Handle(Adaptor3d_HSurface)&  S1,
                                          const Handle(Adaptor3d_TopolTool)& D1,
                                          const Handle(Adaptor3d_HSurface)&  S2,
                                          const Handle(Adaptor3d_TopolTool)& D2,
                                          const Standard_Real TolArc,
                                          const Standard_Real TolTang,
                                          const Standard_Boolean theIsReqToKeepRLine)
{
  myDone = IntStatus_Fail;
  spnt.Clear();
  slin.Clear();

  Standard_Boolean isPostProcessingRequired = Standard_True;

  empt = Standard_True;
  tgte = Standard_False;
  oppo = Standard_False;

  Standard_Boolean all1 = Standard_False;
  Standard_Boolean all2 = Standard_False;
  Standard_Boolean SameSurf = Standard_False;
  Standard_Boolean multpoint = Standard_False;

  Standard_Boolean nosolonS1 = Standard_False;
  // indique s il y a des points sur restriction du carreau 1
  Standard_Boolean nosolonS2 = Standard_False;
  // indique s il y a des points sur restriction du carreau 2
  Standard_Integer i, nbpt, nbseg;
  IntPatch_SequenceOfSegmentOfTheSOnBounds edg1,edg2;
  IntPatch_SequenceOfPathPointOfTheSOnBounds pnt1,pnt2;
  //
  // On commence par intersecter les supports des surfaces
  IntSurf_Quadric quad1, quad2;
  IntPatch_ArcFunction AFunc;
  const Standard_Real Tolang = 1.e-8;
  GeomAbs_SurfaceType typs1, typs2;
  Standard_Boolean bEmpty = Standard_False;
  //
  const Standard_Integer iT1 = SetQuad(S1, typs1, quad1);
  const Standard_Integer iT2 = SetQuad(S2, typs2, quad2);
  //
  if (!iT1 || !iT2) {
    Standard_ConstructionError::Raise();
    return;
  }
  //
  const Standard_Boolean bReverse = iT1 > iT2;
  const Standard_Integer iTT = iT1*10 + iT2;
  //
  switch (iTT) {
    case 11: { // Plane/Plane
      if (!IntPP(quad1, quad2, Tolang, TolTang, SameSurf, slin)) {
        return;
      }
      break;
    }
    //
    case 12:
    case 21: { // Plane/Cylinder
      Standard_Real VMin, VMax, H;
      //
      const Handle(Adaptor3d_HSurface)& aSCyl = bReverse ? S2 : S1;
      VMin = aSCyl->FirstVParameter();
      VMax = aSCyl->LastVParameter();
      H = (Precision::IsNegativeInfinite(VMin) || 
           Precision::IsPositiveInfinite(VMax)) ? 0 : (VMax - VMin);
      //
      if (!IntPCy(quad1, quad2, Tolang, TolTang, bReverse, empt, slin, H)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 13:
    case 31: { // Plane/Cone
      if (!IntPCo(quad1, quad2, Tolang, TolTang, bReverse, empt, multpoint, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 14:
    case 41: { // Plane/Sphere
      if (!IntPSp(quad1, quad2, Tolang, TolTang, bReverse, empt, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 15:
    case 51: { // Plane/Torus
      if (!IntPTo(quad1, quad2, TolTang, bReverse, empt, slin)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 22:
      { // Cylinder/Cylinder
        Bnd_Box2d aBox1, aBox2;

        const Standard_Real aU1f = S1->FirstUParameter();
        Standard_Real aU1l = S1->LastUParameter();
        const Standard_Real aU2f = S2->FirstUParameter();
        Standard_Real aU2l = S2->LastUParameter();

        const Standard_Real anUperiod = 2.0*M_PI;

        if(aU1l - aU1f > anUperiod)
          aU1l = aU1f + anUperiod;

        if(aU2l - aU2f > anUperiod)
          aU2l = aU2f + anUperiod;

        aBox1.Add(gp_Pnt2d(aU1f, S1->FirstVParameter()));
        aBox1.Add(gp_Pnt2d(aU1l, S1->LastVParameter()));
        aBox2.Add(gp_Pnt2d(aU2f, S2->FirstVParameter()));
        aBox2.Add(gp_Pnt2d(aU2l, S2->LastVParameter()));

        // Resolution is too big if the cylinder radius is
        // too small. Therefore, we shall bounde its value above. 
        // Here, we use simple constant.
        const Standard_Real a2DTol = Min(1.0e-4, Min( S1->UResolution(TolTang),
                                            S2->UResolution(TolTang)));

        //The bigger range the bigger nunber of points in Walking-line (WLine)
        //we will be able to add and, consequently, we will obtain more
        //precise intersection line.
        //Every point of WLine is determined as function from U1-parameter,
        //where U1 is U-parameter on 1st quadric.
        //Therefore, we should use quadric with bigger range as 1st parameter
        //in IntCyCy() function.
        //On the other hand, there is no point in reversing in case of
        //analytical intersection (when result is line, ellipse, point...).
        //This result is independent of the arguments order.
        Standard_Boolean isReversed = ((aU2l - aU2f) < (aU1l - aU1f));

        if(isReversed)
        {
          myDone = IntCyCy(quad2, quad1, TolTang, a2DTol, aBox2, aBox1,
                                    Standard_True, empt, SameSurf, multpoint, slin, spnt);
        }
        else
        {
          myDone = IntCyCy(quad1, quad2, TolTang, a2DTol, aBox1, aBox2,
                                    Standard_False, empt, SameSurf, multpoint, slin, spnt);
        }

        if (myDone == IntPatch_ImpImpIntersection::IntStatus_Fail)
        {
          return;
        }

        bEmpty = empt;
        if(!slin.IsEmpty())
        {
          const Handle(IntPatch_WLine)& aWLine = 
                                    Handle(IntPatch_WLine)::DownCast(slin.Value(1));

          if(!aWLine.IsNull())
          {//No geometric solution
            isPostProcessingRequired = Standard_False;
          }
        }

        break;
      }
    //
    case 23:
    case 32: { // Cylinder/Cone
      if (!IntCyCo(quad1, quad2, TolTang, bReverse, empt, multpoint, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 24:
    case 42: { // Cylinder/Sphere
      if (!IntCySp(quad1, quad2, TolTang, bReverse, empt, multpoint, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 25:
    case 52: { // Cylinder/Torus
      if (!IntCyTo(quad1, quad2, TolTang, bReverse, empt, slin)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 33: { // Cone/Cone
      if (!IntCoCo(quad1, quad2, TolTang, empt, SameSurf, multpoint, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 34:
    case 43: { // Cone/Sphere
      if (!IntCoSp(quad1, quad2, TolTang, bReverse, empt, multpoint, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 35:
    case 53: { // Cone/Torus
      if (!IntCoTo(quad1, quad2, TolTang, bReverse, empt, slin)) {
        return;
      }
      break;
    }
    //
    case 44: { // Sphere/Sphere
      if (!IntSpSp(quad1, quad2, TolTang, empt, SameSurf, slin, spnt)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 45:
    case 54: { // Sphere/Torus
      if (!IntSpTo(quad1, quad2, TolTang, bReverse, empt, slin)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    case 55: { // Torus/Torus
      if (!IntToTo(quad1, quad2, TolTang, SameSurf, empt, slin)) {
        return;
      }
      bEmpty = empt;
      break;
    }
    //
    default: {
      Standard_ConstructionError::Raise();
      break;
    }
  }
  //
  if (bEmpty) {
    if (myDone == IntStatus_Fail)
      myDone = IntStatus_OK;

    return;
  }
  //

  if(isPostProcessingRequired)
  {
    if (!SameSurf) {
      AFunc.SetQuadric(quad2);
      AFunc.Set(S1);
    
      solrst.Perform(AFunc, D1, TolArc, TolTang);
      if (!solrst.IsDone()) {
        return;
      }

      if (solrst.AllArcSolution() && typs1 == typs2) {
        all1 = Standard_True;
      }
      nbpt = solrst.NbPoints();
      nbseg= solrst.NbSegments();
      for (i=1; i<= nbpt; i++) {
        pnt1.Append(solrst.Point(i));
      }
      for (i=1; i<= nbseg; i++) {
        edg1.Append(solrst.Segment(i));
      }
      nosolonS1 = (nbpt == 0) && (nbseg == 0);

      if (nosolonS1 && all1) {  // cas de face sans restrictions
        all1 = Standard_False;
      }
    }//if (!SameSurf) {
    else {
      nosolonS1 = Standard_True;
    }

    if (!SameSurf) {
      AFunc.SetQuadric(quad1);
      AFunc.Set(S2);

      solrst.Perform(AFunc, D2, TolArc, TolTang);
      if (!solrst.IsDone()) {
        return;
      }
    
      if (solrst.AllArcSolution() && typs1 == typs2) {
        all2 = Standard_True;
      }
      nbpt = solrst.NbPoints();
      nbseg= solrst.NbSegments();
      for (i=1; i<= nbpt; i++) {
        pnt2.Append(solrst.Point(i));
      }
    
      for (i=1; i<= nbseg; i++) {
        edg2.Append(solrst.Segment(i));
      }
      nosolonS2 = (nbpt == 0) && (nbseg == 0);

      if (nosolonS2 && all2) {  // cas de face sans restrictions
        all2 = Standard_False;
      }
    }// if (!SameSurf) {
    else {
      nosolonS2 = Standard_True;
    }
    //
    if (SameSurf || (all1 && all2)) {
      // faces "paralleles" parfaites
      empt = Standard_False;
      tgte = Standard_True;
      slin.Clear();
      spnt.Clear();

      gp_Pnt Ptreference;

      switch (typs1) {
      case GeomAbs_Plane:      {
        Ptreference = (S1->Plane()).Location();
      }
        break;
      case GeomAbs_Cylinder: {
        Ptreference = ElSLib::Value(0.,0.,S1->Cylinder());
      }
        break;
      case GeomAbs_Sphere:      {
        Ptreference = ElSLib::Value(M_PI/4.,M_PI/4.,S1->Sphere());
      }
        break;
      case GeomAbs_Cone:      {
        Ptreference = ElSLib::Value(0.,10.,S1->Cone());
      }
        break;
      case GeomAbs_Torus:      {
        Ptreference = ElSLib::Value(0.,0.,S1->Torus());
      }
        break;
      default: 
        break;
      }
      //
      oppo = quad1.Normale(Ptreference).Dot(quad2.Normale(Ptreference)) < 0.0;
      myDone = IntStatus_OK;
      return;
    }// if (SameSurf || (all1 && all2)) {

    if (!nosolonS1 || !nosolonS2) {
      empt = Standard_False;
      // C est la qu il faut commencer a bosser...
      PutPointsOnLine(S1,S2,pnt1, slin, Standard_True, D1, quad1,quad2,
                      multpoint,TolArc);
    
      PutPointsOnLine(S1,S2,pnt2, slin, Standard_False,D2, quad2,quad1,
                      multpoint,TolArc);

      if (edg1.Length() != 0) {
        ProcessSegments(edg1,slin,quad1,quad2,Standard_True,TolArc);
      }

      if (edg2.Length() != 0) {
        ProcessSegments(edg2,slin,quad1,quad2,Standard_False,TolArc);
      }

      if (edg1.Length() !=0 || edg2.Length() !=0) {
        //      ProcessRLine(slin,S1,S2,TolArc);
        ProcessRLine(slin,quad1,quad2,TolArc, theIsReqToKeepRLine);
      }
    }//if (!nosolonS1 || !nosolonS2) {
    else {
      empt = ((slin.Length()==0) && (spnt.Length()==0));
    }
  }
  
  Standard_Integer  nblin = slin.Length(),
                    aNbPnt = spnt.Length();
  //
  //modified by NIZNHY-PKV Tue Sep 06 10:03:35 2011f
  if (aNbPnt) {
    IntPatch_SequenceOfPoint aSIP;
    //
    for(i=1; i<=aNbPnt; ++i) {
      Standard_Real aU1, aV1, aU2, aV2;
      gp_Pnt2d aP2D;
      TopAbs_State aState1, aState2;
      //
      const IntPatch_Point& aIP=spnt(i);
      aIP.Parameters(aU1, aV1, aU2, aV2);
      //
      aP2D.SetCoord(aU1, aV1);
      aState1=D1->Classify(aP2D, TolArc);
      //
      aP2D.SetCoord(aU2, aV2);
      aState2=D2->Classify(aP2D, TolArc);
      //
      if(aState1!=TopAbs_OUT && aState2!=TopAbs_OUT) {
        aSIP.Append(aIP);
      }
    }
    //
    spnt.Clear();
    //
    aNbPnt=aSIP.Length();
    for(i=1; i<=aNbPnt; ++i) {
      const IntPatch_Point& aIP=aSIP(i);
      spnt.Append(aIP);
    }
    //
  }//  if (aNbPnt) {
  //modified by NIZNHY-PKV Tue Sep 06 10:18:20 2011t
  //
  for(i=1; i<=nblin; i++) {
    IntPatch_IType thetype = slin.Value(i)->ArcType();
    if(  (thetype ==  IntPatch_Ellipse)
       ||(thetype ==  IntPatch_Circle)
       ||(thetype ==  IntPatch_Lin)
       ||(thetype ==  IntPatch_Parabola)
       ||(thetype ==  IntPatch_Hyperbola)) { 
      Handle(IntPatch_GLine)& glin = *((Handle(IntPatch_GLine)*)&slin.Value(i));
    glin->ComputeVertexParameters(TolArc); 
    }
    else if(thetype == IntPatch_Analytic) { 
      Handle(IntPatch_ALine)& aligold = *((Handle(IntPatch_ALine)*)&slin.Value(i));
      aligold->ComputeVertexParameters(TolArc);
    }
    else if(thetype == IntPatch_Restriction) {
      Handle(IntPatch_RLine)& rlig = *((Handle(IntPatch_RLine)*)&slin.Value(i));
      rlig->ComputeVertexParameters(TolArc); 
    }
  }
  //
  //----------------------------------------------------------------
  //-- On place 2 vertex sur les courbes de GLine qui n en 
  //-- contiennent pas. 
  for(i=1; i<=nblin; i++) {
    gp_Pnt P;
    IntPatch_Point point;
    Standard_Real u1,v1,u2,v2; 
    if(slin.Value(i)->ArcType() == IntPatch_Circle) { 
      const Handle(IntPatch_GLine)& glin = *((Handle(IntPatch_GLine)*)&slin.Value(i));
      if(glin->NbVertex() == 0) { 
        gp_Circ Circ = glin->Circle();
        P=ElCLib::Value(0.0,Circ);
        quad1.Parameters(P,u1,v1);
        quad2.Parameters(P,u2,v2);
        point.SetValue(P,TolArc,Standard_False);
        point.SetParameters(u1,v1,u2,v2);
        point.SetParameter(0.0);
        glin->AddVertex(point);

        P=ElCLib::Value(0.0,Circ);
        quad1.Parameters(P,u1,v1);
        quad2.Parameters(P,u2,v2);
        point.SetValue(P,TolArc,Standard_False);
        point.SetParameters(u1,v1,u2,v2);
        point.SetParameter(M_PI+M_PI);
        glin->AddVertex(point);
      }
    }
    
    else if(slin.Value(i)->ArcType() == IntPatch_Ellipse) { 
      const Handle(IntPatch_GLine)& glin = *((Handle(IntPatch_GLine)*)&slin.Value(i));
      if(glin->NbVertex() == 0) { 
        gp_Elips Elips = glin->Ellipse();
        P=ElCLib::Value(0.0,Elips);
        quad1.Parameters(P,u1,v1);
        quad2.Parameters(P,u2,v2);
        point.SetValue(P,TolArc,Standard_False);
        point.SetParameters(u1,v1,u2,v2);
        point.SetParameter(0.0);
        glin->AddVertex(point);

        P=ElCLib::Value(0.0,Elips);
        quad1.Parameters(P,u1,v1);
        quad2.Parameters(P,u2,v2);
        point.SetValue(P,TolArc,Standard_False);
        point.SetParameters(u1,v1,u2,v2);
        point.SetParameter(M_PI+M_PI);
        glin->AddVertex(point);
      }
    }
  }
  myDone = IntStatus_OK;
}

//=======================================================================
//function : SetQuad
//purpose  : 
//=======================================================================
Standard_Integer SetQuad(const Handle(Adaptor3d_HSurface)& theS,
                         GeomAbs_SurfaceType& theTS,
                         IntSurf_Quadric& theQuad)
{
  theTS = theS->GetType();
  Standard_Integer iRet = 0;
  switch (theTS) {
  case GeomAbs_Plane:
    theQuad.SetValue(theS->Plane());
    iRet = 1;
    break;
  case GeomAbs_Cylinder:
    theQuad.SetValue(theS->Cylinder());
    iRet = 2;
    break;
  case GeomAbs_Cone:
    theQuad.SetValue(theS->Cone());
    iRet = 3;
    break;
  case GeomAbs_Sphere:
    theQuad.SetValue(theS->Sphere());
    iRet = 4;
    break;
  case GeomAbs_Torus:
    theQuad.SetValue(theS->Torus());
    iRet = 5;
    break;
  default:
    break;
  }
  //
  return iRet;
}