Adding tolerance to a point in polygon test

Question

I've been using this method which was taken from Game Coding Complete to detect whether a point is inside of a polygon. It works in almost every case, but is failing on a few edge cases, and I can't figure out the reason.

For example, given a polygon with vertices at (0,0) (0,100) and (100,100), the algorithm is returning:

• True for any point strictly inside the polygon
• False for any of the vertices
• False for (0, 50) which lies on one of the edges of the polygon
• True (?) for (50,50) which is also on one of the edges of the polygon

I'd actually like to relax the algorithm so that it returns true in all of these cases. In other words, it should return true for points that are strictly inside, for the vertices themselves, and for points on the edges of the polygon.

If possible I'd also like to give it enough tolerance so that it always tend towards "true" in face of floating point fluctuations. For example, I have another method, that given a line segment and a point, returns the closest location on the line segment to the given point.

Currently, given any point outside the polygon and one of its edges, there are cases where the result is categorized as being inside by the method above, while other points are considered outside. I'd like to give it enough tolerance so that it always returns true in this situation.

The way I've currently solved the problem is an hack, which consists of using an external library to inflate the polygon by a few pixels, and performing the tests on the inflated polygon, but I'd really like to replace this with a proper solution.

Answer 1

Since the algorithm already parses all edges of the polygon to see how many times a ray cast from test crosses them, I think it's reasonable to add a check to see whether test lies on the edge exactly (or within an epsilon).

In order to avoid too much additional complexity (the usual point - segment distance computation is really awful), I suggest approximating the edge to an extremely thin ellipse and see whether test is on that ellipse. This is the resulting code (only the inner loop shown):

oldPoint = polygon[points-1];
float oldSqDist = sqlength(oldPoint - test);

for (unsigned int i=0 ; i < points; i++)
{
    newPoint = polygon[i];
    float newSqDist = sqlength(newPoint - test);

    if (oldSqDist + newSqDist + 2.0f * std::sqrt(oldSqDist * newSqDist) - sqlength(newPoint - oldPoint) < EPSILON)
        return true;

    if (newPoint.x > oldPoint.x)
    {
        left = oldPoint;
        right = newPoint;
    }
    else
    {
        left = newPoint;
        right = oldPoint;
    }

    if ((newPoint.x < test.x) == (test.x <= oldPoint.x)
       && (test.y-left.y) * (right.x-left.x)
        < (right.y-left.y) * (test.x-left.x) )
    {
        inside=!inside;
    }

    oldPoint = newPoint;
    oldSqDist = newSqDist;
}

Where sqlength() is your favourite way of doing x * x + y * y. I tested the code with EPSILON = 1e-10f and got good results. In real life you should have EPSILON vary relatively to the size of your polygon.

Answer 2

A common approach to this sort of thing is to add epsilons to the floating-point comparisons. I'm not completely sure how to do this in the code you linked, since it seems that it's flipping the inside/outside bit as a result of some comparisons - so perhaps you'd need the epsilon to be positive when currently inside, and negative when outside, or vice versa. The idea would be to make the compares forgiving in the direction that tends to make the result 'inside'.

Answer 3

May I suggest another approach? As described in Real Time Collision Detection, to test if a point is inside a convex polygon oriented counter-clockwise, one fast way is to subdivide in half the polygon and check if the triangle formed by v0 (the first vertex of the poly), vk (k = n/2, with n the number of vertices) and p is oriented CCW. If so we can check again, incrementing k and see if the triangle is again CCW or not. If not the point lies inside the polygon, if not we must check again. The test is implemented as follow (just copying from the text):

int PointInConvexPolygon(Point p, int n, Point v[])
{
    int low = 0, high = n;
    do 
    {
        int mid = (low + high)/2;
        if (TriangleIsCCW(v[0], v[mid], p)) low = mid;
        else high = mid;
    } 
    while (low + 1 < high);
    if (low == 0 || high == n) return 0;
    return TriangleIsCCW(v[low], v[high], p);
}

http://pastebin.com/QkHi1Hrg

Here TriangleIsCCW checks if p lies on the left or on the right the segment formed by v[low] and v[high], and can be implemented as a method that checks if the determinant |ax ay 1| |bx by 1| |cx cy 1| with A (ax, ay), B (bx, by), C (cx, cy) the points checked. If the determinant is positive, the points are oriented CCW, if is 0 they are collinear. This can be used for particular cases as the points lying on the edge or on the vertex.