
#include <d3tess.h>
#include <d3mincentroid.h>



typedef point2<double> pt2;
typedef point3<double> pt3;

#define SHOW(x) #x << '=' << (x)



d3mincentroid::d3mincentroid(d3tess & tess_)
  : d3minoperator(tess_)
{
}

boolc d3mincentroid::eval(uintc a, uintc b)
{
  if (tess.isconvex(a,b)==false)
    return false;

  simplexD2linked const & A(tess.vi[a]);
  simplexD2linked const & B(tess.vi[b]);

  uint ai = A.niInverse(b);
  uint bi = B.niInverse(a);

  uint p[4];
  p[0] = A.pi[ai];
  p[3] = B.pi[bi];
  p[1] = A.pi[(ai+1)%3];
  p[2] = A.pi[(ai+2)%3];

  vector<pt3> const & pt(tess.pt);

  pt2 P0(pt[p[0]]);
  pt2 P1(pt[p[1]]);
  pt2 P2(pt[p[2]]);
  pt2 P3(pt[p[3]]);

  // Calculate the centroid.
  pt2 C0(P0);
  C0 += P1;
  C0 += P2;
  C0 *= (1.0/3.0);

  // Calculate the average radius of points to the centroid.
  double radius(0.0);
  radius += (C0-P0).distance();
  radius += (C0-P1).distance();
  radius += (C0-P2).distance();
  radius /= 3.0;

  if ((C0-P3).dot() > radius*radius )
    return false;

  tess.flip(a,b);

  return true;
}



