/*****************************************************************************/
/*                                                                           */
/*           F I S T :  Fast, Industrial-Strength Triangulation              */
/*                                                                           */
/*****************************************************************************/
/*                                                                           */
/* (C)                            Martin Held                                */
/* (C)            Universitaet Salzburg, Salzburg, Austria                   */
/*                                                                           */
/* This code is not in the public domain. All rights reserved! Please make   */
/* sure to read the full copyright statement contained in api_functions.cpp. */
/*                                                                           */
/*****************************************************************************/


#define BaseLength(u, v, base, delta) { \
   (base).x = (v).x - (u).x; \
   (base).y = (v).y - (u).y; \
   delta    = Abs(TO_MDOUBLE((base).x)) + Abs(TO_MDOUBLE((base).y)); }


#define SideLength(u, v, base, delta) { \
   (base).x = (v).x - (u).x; \
   (base).y = (v).y - (u).y; \
   delta    = (base).x * (base).x + (base).y * (base).y; }

#define PntPntDist(u, v, base, delta) { \
   (base).x = (v).x - (u).x; \
   (base).y = (v).y - (u).y; \
   delta    = sqrt((base).x * (base).x + (base).y * (base).y); }

#define PntPntDistSqd(u, v, base, delta) { \
   (base).x = TO_MDOUBLE((v).x) - TO_MDOUBLE((u).x); \
   (base).y = TO_MDOUBLE((v).y) - TO_MDOUBLE((u).y); \
   delta    = (base).x * (base).x + (base).y * (base).y; }


/*                                                                           */
/* this macro checks whether  i3,  which is collinear with  i1, i2,  is      */
/* between  i1, i2. note that we rely on the lexicographic sorting of the    */
/* points!                                                                   */
/*                                                                           */
#define InBetween(i1, i2, i3) \
   ((i1 <= i3)  &&  (i3 <= i2))


#define StrictlyInBetween(i1, i2, i3) \
   ((i1 < i3)  &&  (i3 < i2))


/*                                                                           */
/* this macro checks whether  i3,  which is collinear with  i1, i2,  is      */
/* between  i1, i2. note that we do _not_ rely on the lexicographic sorting  */
/* of the points!                                                            */
/*                                                                           */
#define InBetweenEps(data, i1, i2, i3) \
   (((Min(data->points[i1].x, data->points[i2].x) - EPS) <= data->points[i3].x)  && \
    ((Min(data->points[i1].y, data->points[i2].y) - EPS) <= data->points[i3].y)  && \
    ((Max(data->points[i1].x, data->points[i2].x) + EPS) >= data->points[i3].x)  && \
    ((Max(data->points[i1].y, data->points[i2].y) + EPS) >= data->points[i3].y))


/*                                                                           */
/* this macro computes the determinant  det(points[i],points[j],points[k])   */
/* in a consistent way.                                                      */
/*                                                                           */
#ifdef JRC_PREDICATE
double orient2dfast(point pa, point pb, point pc);
double orient2d(point pa, point pb, point pc);
#define StableDet2D(i, j, k, det) \
 { assert(InPointsList(i)); \
   assert(InPointsList(j)); \
   assert(InPointsList(k)); \
   \
   if ((i == j)  ||  (i == k)  ||   (j == k)) { \
      det = C_0_0; \
   } \
   else { \
      det = orient2d(points[i], points[j], points[k]); \
   } \
}
#else
#define StableDet2D(data, tmp, i, j, k, det) \
 { assert(InPointsList(data,i)); \
   assert(InPointsList(data,j)); \
   assert(InPointsList(data,k)); \
   \
   if ((i == j)  ||  (i == k)  ||   (j == k)) { \
      det = C_0_0; \
   } \
   else { \
	  tmp.numerics_h_p = data->points[i]; \
	  tmp.numerics_h_q = data->points[j]; \
	  tmp.numerics_h_r = data->points[k]; \
      \
      if (i < j) { \
         if (j < k)            /* i < j < k */ \
            det =  Det2D(tmp.numerics_h_p, tmp.numerics_h_q, tmp.numerics_h_r); \
         else if (i < k)       /* i < k < j */ \
            det = -Det2D(tmp.numerics_h_p, tmp.numerics_h_r, tmp.numerics_h_q); \
         else                  /* k < i < j */ \
            det =  Det2D(tmp.numerics_h_r, tmp.numerics_h_p, tmp.numerics_h_q); \
      } \
      else { \
         if (i < k)            /* j < i < k */ \
            det = -Det2D(tmp.numerics_h_q, tmp.numerics_h_p, tmp.numerics_h_r); \
         else if (j < k)       /* j < k < i */ \
            det =  Det2D(tmp.numerics_h_q, tmp.numerics_h_r, tmp.numerics_h_p); \
         else                  /* k < j < i */ \
            det = -Det2D(tmp.numerics_h_r, tmp.numerics_h_q, tmp.numerics_h_p); \
      } \
   } \
 }
#endif


/*                                                                           */
/* this macro sets ori to +1 if the points  i, j, k are given in CCW order,  */
/*                        -1 if the points  i, j, k are given in CW order,   */
/*                         0 if the points  i, j, k are collinear.           */
/*                                                                           */
#define Orientation(data, tmp, i, j, k, ori) \
 { StableDet2D(data, tmp, i, j, k, tmp.numerics_h_det); \
   \
   if (lt(tmp.numerics_h_det, EPS))        ori = -1; \
   else if (gt(tmp.numerics_h_det, EPS))   ori =  1; \
   else                                    ori =  0; \
}


/*                                                                           */
/* this macro sets ori to +1 if the points  i, j, k are given in CCW order,  */
/*                        -1 if the points  i, j, k are given in CW order,   */
/*                         0 if the points  i, j, k are collinear.           */
/*                                                                           */
#define OrientationThres(data, tmp, i, j, k, ori, thres)      \
 { StableDet2D(data, tmp, i, j, k, tmp.numerics_h_det); \
   \
   if (lt(tmp.numerics_h_det, thres))        ori = -1; \
   else if (gt(tmp.numerics_h_det, thres))   ori =  1; \
   else                                      ori =  0; \
}


/*                                                                           */
/* this macro checks whether  l  is in the cone defined by  i, j  and  j, k  */
/*                                                                           */
#define IsInConeStrict(data, tmp, i, j, k, l, convex, flag, thres)    \
 { assert(InPointsList(data, i)); \
   assert(InPointsList(data, j)); \
   assert(InPointsList(data, k)); \
   assert(InPointsList(data, l)); \
   \
   flag = true; \
   if (convex) { \
      if ((l != i)  &&  (i != j)) { \
         OrientationThres(data, tmp, i, j, l, tmp.numerics_h_ori1, thres);    \
         if (tmp.numerics_h_ori1 <= 0)              flag = false; \
      } \
      if ((l != k)  &&  (j != k)  &&  (flag == true)) { \
         OrientationThres(data, tmp, j, k, l, tmp.numerics_h_ori2, thres);    \
         if (tmp.numerics_h_ori2 <= 0)              flag = false; \
      } \
   } \
   else { \
      OrientationThres(data, tmp, i, j, l, tmp.numerics_h_ori1, thres);   \
      if (tmp.numerics_h_ori1 <= 0) { \
         OrientationThres(data, tmp, j, k, l, tmp.numerics_h_ori2, thres);    \
         if (tmp.numerics_h_ori2 <= 0)              flag = false; \
      } \
   } \
}


/*                                                                           */
/* this macro checks whether  l  is in the cone defined by  i, j  and  j, k  */
/*                                                                           */
#define IsInCone(data, tmp, i, j, k, l, convex, flag) \
 { assert(InPointsList(data, i)); \
   assert(InPointsList(data, j)); \
   assert(InPointsList(data, k)); \
   assert(InPointsList(data, l)); \
   \
   flag = true; \
   if (convex) { \
      if ((l != i)  &&  (i != j)) { \
         Orientation(data, tmp, i, j, l, tmp.numerics_h_ori1); \
         if (tmp.numerics_h_ori1 < 0)              flag = false; \
         else if (tmp.numerics_h_ori1 == 0) { \
            if (i < j) { \
               if (!InBetween(i, j, l))            flag = false; \
            } \
            else { \
               if (!InBetween(j, i, l))            flag = false; \
            } \
         } \
      } \
      if ((l != k)  &&  (j != k)  &&  (flag == true)) { \
         Orientation(data, tmp, j, k, l, tmp.numerics_h_ori2); \
         if (tmp.numerics_h_ori2 < 0)                flag = false; \
         else if (tmp.numerics_h_ori2 == 0) { \
            if (j < k) { \
               if (!InBetween(j, k, l))              flag = false; \
            } \
            else { \
               if (!InBetween(k, j, l))              flag = false; \
            } \
         } \
      } \
   } \
   else { \
      Orientation(data, tmp, i, j, l, tmp.numerics_h_ori1); \
      if (tmp.numerics_h_ori1 <= 0) { \
         Orientation(data, tmp, j, k, l, tmp.numerics_h_ori2); \
         if (tmp.numerics_h_ori2 <= 0) { \
            if (tmp.numerics_h_ori1 == 0) { \
               if (!InBetween(i, j, l))             flag = false; \
            } \
            else if (tmp.numerics_h_ori2 == 0) { \
               if (!InBetween(j, k, l))             flag = false; \
            } \
            else                                    flag = false; \
         } \
      } \
   } \
}



/*                                                                           */
/* this macro checks whether  l  is in the cone defined by  i, j  and  j, k  */
/*                                                                           */
#define IsInConeEps(data, tmp, i, j, k, l, convex, flag) \
 { assert(InPointsList(data, i)); \
   assert(InPointsList(data, j)); \
   assert(InPointsList(data, k)); \
   assert(InPointsList(data, l)); \
   \
   flag = true; \
   if (convex) { \
      if ((l != i)  &&  (i != j)) { \
         Orientation(data, tmp, i, j, l, tmp.numerics_h_ori1); \
         if (tmp.numerics_h_ori1 < 0)              flag = false; \
         else if (tmp.numerics_h_ori1 == 0) { \
            if (!InBetweenEps(data, i, j, l))      flag = false; \
         } \
      } \
      if ((l != k)  &&  (j != k)  &&  (flag == true)) { \
         Orientation(data, tmp, j, k, l, tmp.numerics_h_ori2); \
         if (tmp.numerics_h_ori2 < 0)              flag = false; \
         else if (tmp.numerics_h_ori2 == 0) { \
            if (!InBetweenEps(data, j, k, l))      flag = false; \
         } \
      } \
   } \
   else { \
      Orientation(data, tmp, i, j, l, tmp.numerics_h_ori1); \
      if (tmp.numerics_h_ori1 <= 0) { \
         Orientation(data, tmp, j, k, l, tmp.numerics_h_ori2); \
         if (tmp.numerics_h_ori2 <= 0) { \
            if (tmp.numerics_h_ori1 == 0) { \
               if (!InBetweenEps(data, i, j, l))    flag = false; \
            } \
            else if (tmp.numerics_h_ori2 == 0) { \
               if (!InBetweenEps(data, j, k, l))    flag = false; \
           } \
         else                                       flag = false; \
         } \
      } \
   } \
}




/*                                                                           */
/* this macro checks whether the diagonal  j, l  is in the cone defined by   */
/* i, j  and  j, k. if  left  then  j, k, l  forms the ear to be tested.     */
/* otherwise,  i, j, k  forms the ear to be tested.                          */
/*                                                                           */
#define DiagIsInConeEps(data, tmp, i, j, k, l, convex, left, flag) \
 { assert(InPointsList(data, i)); \
   assert(InPointsList(data, j)); \
   assert(InPointsList(data, k)); \
   assert(InPointsList(data, l)); \
   \
   flag = true; \
   if (convex) { \
      if (left) { \
         if ((l != i)  &&  (i != j)) { \
            Orientation(data, tmp, i, j, l,  tmp.numerics_h_ori1); \
            if (tmp.numerics_h_ori1 < 0)              flag = false; \
            else if (tmp.numerics_h_ori1 == 0) { \
               if (InBetweenEps(data, j, l, i))       flag = false; \
            } \
         } \
      } \
      else { \
         if ((l != k)  &&  (j != k)) { \
            Orientation(data, tmp, j, k, l, tmp.numerics_h_ori2); \
            if (tmp.numerics_h_ori2 < 0)              flag = false; \
            else if (tmp.numerics_h_ori2 == 0) { \
               if (InBetweenEps(data, j, l, k))       flag = false; \
            } \
         } \
      } \
   } \
}



#define DiagIsInCone(data, tmp, i, j, k, l, convex, left, flag) \
 { assert(InPointsList(data, i)); \
   assert(InPointsList(data, j)); \
   assert(InPointsList(data, k)); \
   assert(InPointsList(data, l)); \
   \
   flag = true; \
   if (convex) { \
      if (left) { \
         if ((l != i)  &&  (i != j)) { \
            Orientation(data, tmp, i, j, l, tmp.numerics_h_ori1); \
            if (tmp.numerics_h_ori1 < 0)              flag = false; \
            else if (tmp.numerics_h_ori1 == 0) { \
               if (InBetween(j, l, i))                flag = false; \
            } \
         } \
      } \
      else { \
         if ((l != k)  &&  (j != k)) { \
            Orientation(data, tmp, j, k, l, tmp.numerics_h_ori2); \
            if (tmp.numerics_h_ori2 < 0)              flag = false; \
            else if (tmp.numerics_h_ori2 == 0) { \
               if (InBetween(j, l, k))                flag = false; \
            } \
         } \
      } \
   } \
}


/*                                                                           */
/* returns                                                                   */
/*           -2 ... if angle is 360 degrees                                  */
/*           -1 ... if angle between 180 and 360 degrees                     */
/*            0 ... if angle is 180 degrees                                  */
/*            1 ... if angle between 0 and 180 degrees                       */
/*            2 ... if angle is 0 degrees                                    */
/*            3 ... if angle is about 180 degrees and triangle is very small */
/*                                                                           */
inline int IsConvexAngle(listdef *list, datadef *data, int i, int j, int k,
                         list_ind ind) {
   assert(InPointsList(data, i));
   assert(InPointsList(data, j));
   assert(InPointsList(data, k));
   tmp_data_def tmp;

   if (i == j) {
      if (j == k) {
         /*                                                                  */
         /* all three vertices are identical. we set the angle to -2.        */
         /* using -2 means to err on the safe side, as all the               */
         /* incarnations of this vertex will be clipped right at the         */
         /* start of the ear-clipping algorithm. thus, eventually there      */
         /* will be no other duplicates at this vertex position, and the     */
         /* regular classification of angles will yield the correct          */
         /* answer for j.                                                    */
         /*                                                                  */
         return -2;
      }
      else {
         /*                                                                  */
         /* two of the three vertices are identical; we set the angle to     */
         /* to 2 order to enable clipping of  j.                             */
         /*                                                                  */
         return 2;
      }
   }
   else if (j == k) {
      /*                                                                     */
      /* two vertices are identical. we could either determine the angle     */
      /* by means of yet another lengthy analysis, or simply set the         */
      /* angle to -1. using -1 means to err on the safe side, as all the     */
      /* incarnations of this vertex will be clipped right at the start      */
      /* of the ear-clipping algorithm. thus, eventually there will be no    */
      /* other duplicates at this vertex position, and the regular           */
      /* classification of angles will yield the correct answer for j.       */
      /*                                                                     */
      return -1;
   } else {
      StableDet2D(data, tmp, i, j, k, tmp.numerics_h_det);
      if (gt(tmp.numerics_h_det, EPS)) {
         return 1;
      }
      else if (lt(tmp.numerics_h_det, EPS)) {
         return -1;
      }
      else {
         /*                                                                  */
         /* 0, 180, or 360 degrees.                                          */
         /*                                                                  */
         VectorSub2D(data->points[i], data->points[j], tmp.numerics_h_p);
         VectorSub2D(data->points[k], data->points[j], tmp.numerics_h_q);
         tmp.numerics_h_dot = DotProduct2D(tmp.numerics_h_p, tmp.numerics_h_q);
         if (tmp.numerics_h_dot < C_0_0) {
            /*                                                               */
            /* 180 degrees.                                                  */
            /*                                                               */
            /*data->numerics_h_dot  = sqrt(Length2D(data->numerics_h_p));    */
            /*data->numerics_h_dot += sqrt(Length2D(data->numerics_h_q));    */
            if (tmp.numerics_h_det > C_0_0)         return 1;
            else                                    return 0;
         }
         else {
            /*                                                               */
            /* 0 or 360 degrees? this cannot be judged locally, and more     */
            /* work is needed.                                               */
            /*                                                               */
            return SpikeAngle(list, data, i, j, k, ind);
         }
      }
   }
}


/*                                                                           */
/* this macro checks whether point  i4  is inside of or on the boundary      */
/* of the triangle  i1, i2, i3                                               */
/*                                                                           */
#define PntInTriangle(data, tmp, i1, i2, i3, i4, inside) \
 { \
   inside = false; \
   Orientation(data, tmp, i2, i3, i4, tmp.numerics_h_ori1); \
   if (tmp.numerics_h_ori1 >= 0) { \
      Orientation(data, tmp, i1, i2, i4, tmp.numerics_h_ori1); \
      if (tmp.numerics_h_ori1 >= 0) { \
         Orientation(data, tmp, i3, i1, i4, tmp.numerics_h_ori1); \
         if (tmp.numerics_h_ori1 >= 0)  inside = true; \
      } \
   } \
}


/*                                                                           */
/* this macro checks whether point  i4  is inside of or on the boundary      */
/* of the triangle  i1, i2, i3. it also returns a classification if  i4  is  */
/* on the boundary of the triangle (except for the edge  i2, i3).            */
/*                                                                           */
#define VtxInTriangle(data, tmp, i1, i2, i3, i4, inside, classifier) \
 { \
   inside = false; \
   Orientation(data, tmp, i2, i3, i4, tmp.numerics_h_ori1); \
   if (tmp.numerics_h_ori1 >= 0) { \
      Orientation(data, tmp, i1, i2, i4, tmp.numerics_h_ori1); \
      if (tmp.numerics_h_ori1 > 0) { \
         Orientation(data, tmp, i3, i1, i4, tmp.numerics_h_ori1); \
         if (tmp.numerics_h_ori1 > 0) { \
            inside = true; \
            classifier = 0; \
         } \
         else if (tmp.numerics_h_ori1 == 0) { \
            inside = true; \
            classifier = 1; \
         } \
      } \
      else if (tmp.numerics_h_ori1 == 0) { \
         Orientation(data, tmp, i3, i1, i4, tmp.numerics_h_ori1); \
         if (tmp.numerics_h_ori1 > 0) { \
            inside = true; \
            classifier = 2; \
         } \
         else if (tmp.numerics_h_ori1 == 0) { \
            inside = true; \
            classifier = 3; \
         } \
      } \
   } \
}


#define ComputeRatio(A, B, C, RATIO)           \
{ \
   if (Min(B, C) <= ZERO) {                                     \
      RATIO = Min(B, C);                                        \
   }                                                            \
   else {                                                       \
      RATIO = CD_1_0 + (A - B - C) / (CD_2_0 * sqrt(B * C));    \
   }                                                            \
}

