判断一个点是否在三角形内部 以及重心计算

dot product 和 cross product 的区别

• A dot product of two vectors is also called the scalar product. It is the product of the magnitude of the two vectors and the cosine of the angle that they form with each other. (It is also called the inner product or the projection product.)

• A cross product of two vectors is also called the vector product. It is the product of the magnitude of the two vectors and the sine of the angle that they form with each other.
  
  The difference between the dot product and the cross product of two vectors is that the result of the dot product is a scalar quantity, whereas the result of the cross product is a vector quantity.

判断一个点是否在三角形内部 —— cross product法

原文

For a point to be inside the traingle A B C it must be below AB and left of BC and right of AC. If any one of these tests fails we can return early.

how do we tell if a point is on the correct side of a line?

• If you take the cross product of [B-A] and [p-A], you’ll get a vector pointing out of the screen.
• if you take the cross product of [B-A] and [p’-A] you’ll get a vector pointing into the screen.
• In fact， if you cross [B-A] with the vector from A to any point above the line AB, the resulting vector points out of the screen
• while using any point below AB yields a vector pointing into the screen.
  The only question remaining is: how do we know what direction the cross product should point in?
  Answer is : Using point c, any point p where [B-A] cross [p-A] does not point in the same direction as [B-A] cross [C-A] isn’t inside the triangle.

实现伪代码：

function SameSide(p1,p2, a,b)
    cp1 = CrossProduct(b-a, p1-a)
    cp2 = CrossProduct(b-a, p2-a)
    if DotProduct(cp1, cp2) >= 0 then return true
    else return false

function PointInTriangle(p, a,b,c)
    if SameSide(p,a, b,c) and SameSide(p,b, a,c)
        and SameSide(p,c, a,b) then return true
    else return false

判断一个点是否在三角形内部 —— Barycentric法

Now we can get to any point on the plane just by starting at A and walking some distance along (C - A) and then from there walking some more in the direction (B - A).

we can now describe any point on the plane as
P = A + u * (C - A) + v * (B - A)

Notice now that if u or v < 0 then we’ve walked in the wrong direction and must be outside the triangle. Also if u or v > 1 then we’ve walked too far in a direction and are outside the triangle. Finally if u + v > 1 then we’ve crossed the edge BC again leaving the triangle.

Given u and v we can easily calculate the point P with the above equation, but how can we go in the reverse direction and calculate u and v from a given point P? Time for some math!

实现伪代码：

// Compute vectors        
v0 = C - A
v1 = B - A
v2 = P - A

// Compute dot products
dot00 = dot(v0, v0)
dot01 = dot(v0, v1)
dot02 = dot(v0, v2)
dot11 = dot(v1, v1)
dot12 = dot(v1, v2)

// Compute barycentric coordinates
invDenom = 1 / (dot00 * dot11 - dot01 * dot01)
u = (dot11 * dot02 - dot01 * dot12) * invDenom
v = (dot00 * dot12 - dot01 * dot02) * invDenom

// Check if point is in triangle
return (u >= 0) && (v >= 0) && (u + v < 1)

重心计算

原文

It follows the method in shirley, where you compute the area of the triangles formed by embedding the point P inside the triangle.

code

Vector Triangle::getBarycentricCoordinatesAt( const Vector & P ) const
{
  Vector bary ;

  // The area of a triangle is 
  real areaABC = DOT( normal, CROSS( (b - a), (c - a) )  ) ;
  real areaPBC = DOT( normal, CROSS( (b - P), (c - P) )  ) ;
  real areaPCA = DOT( normal, CROSS( (c - P), (a - P) )  ) ;

  bary.x = areaPBC / areaABC ; // alpha
  bary.y = areaPCA / areaABC ; // beta
  bary.z = 1.0f - bary.x - bary.y ; // gamma

  return bary ;
}

另一种方法：
Transcribed from Christer Ericson’s Real-Time Collision Detection (which, incidentally, is an excellent book)

// Compute barycentric coordinates (u, v, w) for
// point p with respect to triangle (a, b, c)
void Barycentric(Point p, Point a, Point b, Point c, float &u, float &v, float &w)
{
    Vector v0 = b - a, v1 = c - a, v2 = p - a;
    float d00 = Dot(v0, v0);
    float d01 = Dot(v0, v1);
    float d11 = Dot(v1, v1);
    float d20 = Dot(v2, v0);
    float d21 = Dot(v2, v1);
    float denom = d00 * d11 - d01 * d01;
    v = (d11 * d20 - d01 * d21) / denom;
    w = (d00 * d21 - d01 * d20) / denom;
    u = 1.0f - v - w;
}