c++: Find the closest point on a mesh to a query point

Question

Here is what I need: Given a point(x,y,z) in 3d space, and a mesh compose of some vertices(x,y,z), to calculate and return the close point coordinate on that mesh. The function probably like this:

bool closePointOnMesh(const Point& queryPoint, const Mesh& myMesh, float maxDistance);

I have done some searching, and probably I will choose octree to reduce the calculation.

But there are still many details that I can't understand:

1: How the octree node been subdivided, so each node contains may contains 0~some triangles? It is easier to subdivided the cell further based on vertices and just store vertices directly.

2: How the octree structure helps to reduce the calculation, I know if the cell is empty I will just disregard it. But do I need to get all the closest point within each triangle face in a octree cell to the queryPoint, so I finally get the most closest point of all? that sound still heavy. Beside it will be more easier if I just iter through all the triangles, get the closest point from them, which means no need for the octree???

3: Is there a fast way to get the closest point to a point within a triangle face?

4: how the maxDistance limit helps to reduce the calculation?

Answer 1

For #3, here's some code on how to get the closest point of a triangle. It projects the point onto the triangle's plane, and then clamps the barycentric coordinates to [0,1], and uses those values computes the closest point.

Copied below:

vector3 closesPointOnTriangle( const vector3 *triangle, const vector3 &sourcePosition )
{
vector3 edge0 = triangle[1] - triangle[0];
vector3 edge1 = triangle[2] - triangle[0];
vector3 v0 = triangle[0] - sourcePosition;

float a = edge0.dot( edge0 );
float b = edge0.dot( edge1 );
float c = edge1.dot( edge1 );
float d = edge0.dot( v0 );
float e = edge1.dot( v0 );

float det = a*c - b*b;
float s = b*e - c*d;
float t = b*d - a*e;

if ( s + t < det )
{
    if ( s < 0.f )
    {
        if ( t < 0.f )
        {
            if ( d < 0.f )
            {
                s = clamp( -d/a, 0.f, 1.f );
                t = 0.f;
            }
            else
            {
                s = 0.f;
                t = clamp( -e/c, 0.f, 1.f );
            }
        }
        else
        {
            s = 0.f;
            t = clamp( -e/c, 0.f, 1.f );
        }
    }
    else if ( t < 0.f )
    {
        s = clamp( -d/a, 0.f, 1.f );
        t = 0.f;
    }
    else
    {
        float invDet = 1.f / det;
        s *= invDet;
        t *= invDet;
    }
}
else
{
    if ( s < 0.f )
    {
        float tmp0 = b+d;
        float tmp1 = c+e;
        if ( tmp1 > tmp0 )
        {
            float numer = tmp1 - tmp0;
            float denom = a-2*b+c;
            s = clamp( numer/denom, 0.f, 1.f );
            t = 1-s;
        }
        else
        {
            t = clamp( -e/c, 0.f, 1.f );
            s = 0.f;
        }
    }
    else if ( t < 0.f )
    {
        if ( a+d > b+e )
        {
            float numer = c+e-b-d;
            float denom = a-2*b+c;
            s = clamp( numer/denom, 0.f, 1.f );
            t = 1-s;
        }
        else
        {
            s = clamp( -e/c, 0.f, 1.f );
            t = 0.f;
        }
    }
    else
    {
        float numer = c+e-b-d;
        float denom = a-2*b+c;
        s = clamp( numer/denom, 0.f, 1.f );
        t = 1.f - s;
    }
}

return triangle[0] + s * edge0 + t * edge1;
}