//============================================================================

// homoclip.cpp : polygon clipping in homogeneous coordinates

// Author: Maxim Garber
// Date:   02/09/2001

//============================================================================

/************************************************************
4D vector functions
************************************************************/
// returns p dot q
float DotProduct(float* p, float*q)
{return p[0]*q[0] + p[1]*q[1] + p[2]*q[2] + p[3]*q[3];}

// p = q
void VectorCopy(float* p, float* q)
{p[0] = q[0]; p[1] = q[1]; p[2] = q[2]; p[3] = q[3];}

// r = p - q 
void VectorSubtract(float* r, float* p, float* q)
{ r[0] = p[0] - q[0]; r[1] = p[1] - q[1]; r[2] = p[2] - q[2]; r[3] = p[3] - q[3];}

// r = p + q
void VectorAdd(float* r, float* p, float* q)
{ r[0] = p[0] + q[0]; r[1] = p[1] + q[1]; r[2] = p[2] + q[2]; r[3] = p[3] + q[3];}

// r = cp  where c is a scalar
void VectorMult(float* r, float* p, float c)
{r[0] = c*p[0]; r[1] = c*p[1]; r[2] = c*p[2]; r[3] = c*p[3];}


/************************************************************
			ClipPolyToPlane
Pre : InPts has NumInPts*4 points, plane is 4-vector which 
      represents a 4D plane equation of the form ax+by+cz+dw = 0.

Post: OutPts contains the points of the polygon after it has
      been clipped to the plane. Terurn value in number of 
      points in OutPts.
************************************************************/

int ClipPolyToPlane(float *InPts, const int NumInPts, float* plane, float* &OutPts)
{
  int NumOutPts = 0;

  // current and previous point in polygon traversal
  float *point, *prev  = InPts;     
  // the point's dot products with the clipping plane
  float  pointDotPlane, prevDotPlane = DotProduct(prev, plane);

  // traverse all points in polygon
  for(int i = 1; i < NumInPts+1; i++)
  {
    // get point and its dot product with the plane
    point = &InPts[4*(i % NumInPts)];
    pointDotPlane = DotProduct(point, plane);

    if( pointDotPlane <= 0)           // current point is inside
    {
      if(prevDotPlane > 0)            // previous point was outside
      {
        // intersect line segment with plane
        float slope[4], intersection[4];
        VectorSubtract(slope, point, prev);
        float t = -prevDotPlane / DotProduct(slope, plane);
        VectorMult(intersection, slope, t);
        VectorAdd(intersection, prev, intersection);
        // copy intersect point to output
        VectorCopy(&OutPts[4*NumOutPts], intersection);
        NumOutPts++;   
      }
      int outIndex = 4*NumOutPts;
      VectorCopy(&OutPts[outIndex], point); // copy current point to output 
      prev = &OutPts[outIndex];             // update previous point
      NumOutPts++;
    }
    else                           // current point is outside
    {
      if(prevDotPlane <= 0)        // previous point was inside
      {
        // intersect line segment with plane
        float slope[4], intersection[4];
        VectorSubtract(slope, point, prev);
        float t = -prevDotPlane / DotProduct(slope, plane);
        VectorMult(intersection, slope, t);
        VectorAdd(intersection, prev, intersection);
        // copy intersection to output
        VectorCopy(&OutPts[4*NumOutPts], intersection);
        NumOutPts++;
      }
      // update previous point
      prev = point;
    }
  prevDotPlane = DotProduct(prev, plane);
  }
  return(NumOutPts);
}


//----------------------------------------------------------------------------

// Clips a 4D homogeneous polygon defined by the packed array of float InPts.

// to the viewing frustum defined by w components of the verts. The clipped polygon is

// put in OutPts (which must be a different location than InPts) and the number

// of vertices in the clipped polygon is returned. InPts must have NumInPts*4

// floats (enough to contain poly). Regular orthographic NDC clipping can be 

// achieved by making the w coordinate of each point be 1. OutPts will be

// allocated and return filled with the clipped polygon (could contain as 

// many as (NumInPts+6)*4 floats)

//----------------------------------------------------------------------------
int ClipPolyToFrustum(float *InPts, const int NumInPts, float* &OutPts)

{

  OutPts  = new float[(NumInPts+6)*4];
  float *OutPts1 = new float[(NumInPts+6)*4];
  int NumOutPts;

  // CLIP InPts TO OutPts IN HOMOGENEOUS COORDS
  
  float  Plane[4] = {1,0,0,-1}; // plane: x <= w
  NumOutPts = ClipPolyToPlane(InPts, NumInPts, Plane, OutPts1);
  Plane[0] = 0; Plane[1] = 1;   // plane: y <= w
  NumOutPts = ClipPolyToPlane(OutPts1, NumOutPts, Plane, OutPts);
  Plane[1] = 0; Plane[2] = 1;   // plane: z <= w
  NumOutPts = ClipPolyToPlane(OutPts, NumOutPts, Plane, OutPts1);

  Plane[2] = 0; Plane[0] = -1;  // plane: x >= -w
  NumOutPts = ClipPolyToPlane(OutPts1, NumOutPts, Plane, OutPts);
  Plane[0] = 0; Plane[1] = -1;  // plane: x >= -w
  NumOutPts = ClipPolyToPlane(OutPts, NumOutPts, Plane, OutPts1);
  Plane[1] = 0; Plane[2] = -1;  // plane: x >= -w
  NumOutPts = ClipPolyToPlane(OutPts1, NumOutPts, Plane, OutPts);

  return(NumOutPts);

}