Method InterpolateFunctionsUsingMVC

Common/DataModel/vtkPolygon.cxx:524–622

------------------------------------------------------------------------------ Compute interpolation weights using mean value coordinate.

//------------------------------------------------------------------------------
// Compute interpolation weights using mean value coordinate.
void vtkPolygon::InterpolateFunctionsUsingMVC(const double x[3], double* weights)
{
  int numPts = this->Points->GetNumberOfPoints();

  // Begin by initializing weights.
  for (int i = 0; i < numPts; i++)
  {
    weights[i] = 0.0;
  }

  // create local array for storing point-to-vertex vectors and distances
  std::vector<double> dist(numPts);
  std::vector<double> uVec(3 * numPts);
  static const double eps = 0.00000001;
  for (int i = 0; i < numPts; i++)
  {
    double pt[3];
    this->Points->GetPoint(i, pt);

    // point-to-vertex vector
    uVec[3 * i] = pt[0] - x[0];
    uVec[3 * i + 1] = pt[1] - x[1];
    uVec[3 * i + 2] = pt[2] - x[2];

    // distance
    dist[i] = vtkMath::Norm(uVec.data() + 3 * i);

    // handle special case when the point is really close to a vertex
    if (dist[i] < eps)
    {
      weights[i] = 1.0;
      return;
    }

    uVec[3 * i] /= dist[i];
    uVec[3 * i + 1] /= dist[i];
    uVec[3 * i + 2] /= dist[i];
  }

  // Now loop over all vertices to compute weight
  // w_i = ( tan(theta_i/2) + tan(theta_(i+1)/2) ) / dist_i
  // To do consider the simplification of
  // tan(alpha/2) = (1-cos(alpha))/sin(alpha)
  //              = (d0*d1 - cross(u0, u1))/(2*dot(u0,u1))
  std::vector<double> tanHalfTheta(numPts);
  for (int i = 0; i < numPts; i++)
  {
    int i1 = i + 1;
    if (i1 == numPts)
    {
      i1 = 0;
    }

    double* u0 = uVec.data() + 3 * i;
    double* u1 = uVec.data() + 3 * i1;

    double l = sqrt(vtkMath::Distance2BetweenPoints(u0, u1));
    double theta = 2.0 * asin(l / 2.0);