How do I efficiently check if a point is inside a rotated rectangle?

Question

Part for the sake of optimization, part for learning purposes, I will dare to ask: How can I most efficiently check whether a 2D point P is inside a 2D rotated rectangle XYZW, using C# or C++?

Currently, what I am doing is using a "point in triangle" algorithm found in the book Real Time Collision Detection, and running it twice (for the two triangles that make up the rectangle, say XYZ and XZW):

bool PointInTriangle(Vector2 A, Vector2 B, Vector2 C, Vector2 P)
{
 // Compute vectors        
 Vector2 v0 = C - A;
 Vector2 v1 = B - A;
 Vector2 v2 = P - A;

 // Compute dot products
 float dot00 = Vector2.Dot(v0, v0);
 float dot01 = Vector2.Dot(v0, v1);
 float dot02 = Vector2.Dot(v0, v2);
 float dot11 = Vector2.Dot(v1, v1);
 float dot12 = Vector2.Dot(v1, v2);

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

 // Check if point is in triangle
 if(u >= 0 && v >= 0 && (u + v) < 1)
    { return true; } else { return false; }
}


bool PointInRectangle(Vector2 X, Vector2 Y, Vector2 Z, Vector2 W, Vector2 P)
{
 if(PointInTriangle(X,Y,Z,P)) return true;
 if(PointInTriangle(X,Z,W,P)) return true;
 return false;
}

However, I have the feeling that there might be a cleaner and faster way. Specifically, to reduce the number of math operations.

Answer 1

An easy and straightforward optimization would be to change the final condition in PointInTriangle:

bool PointInRectangle(Vector2 A, Vector2 B, Vector2 C, Vector2 P) {
  ...
  if(u >= 0 && v >= 0 && u <= 1 && v <= 1)
      { return true; } else { return false; }
  }
}

the code was pretty much PointInRectangle already, the condition (u + v) < 1 was there to check if it is not in the "second" triangle of rectangle.

Alternatively, you could also do an isLeft (first code example on page, also greatly explained) test four times, and check that they all return results with the same sign (which one depends on whether the points were given in clockwise or counterclockwise order) for the point to be inside. This works for any other convex polygon too.

float isLeft( Point P0, Point P1, Point P2 )
{
    return ( (P1.x - P0.x) * (P2.y - P0.y) - (P2.x - P0.x) * (P1.y - P0.y) );
}
bool PointInRectangle(Vector2 X, Vector2 Y, Vector2 Z, Vector2 W, Vector2 P)
{
    return (isLeft(X, Y, P) > 0 && isLeft(Y, Z, P) > 0 && isLeft(Z, W, P) > 0 && isLeft(W, X, p) > 0);
}

Answer 2

Edit: The OPs comment has been skeptical about the efficiency of the suggested negative circular bound check to improve the algorithm in order do checks whether an arbitrary 2D point lies within a rotated and/or moving rectangle. Fiddling around a bit on my 2D game engine (OpenGL/C++), I supplement my answer by providing a performance benchmark of my algorithm against the OPs current point-in-rectangle-check algorithms (and variations).

I originally suggested to leave the algorithm in place (as it is nearly optimal), but simplify through mere game logic: (1) using a pre-processed circle around the original rectangle; (2) do a distance check and if the point lies within the given circle; (3) use the OPs or other any straightforward algorithm (I recommend the isLeft algorithm as provided in another answer). The logic behind my suggestion is that checking whether a point is within a circle is considerably more efficient than a boundary check of a rotated rectangle or any other polygon.

My initial scenario for a benchmark test is to run a large number of appearing and disappearing dots (whose position changes in every game-loop) in a constrained space that will be filled with around 20 rotating/moving squares. I have published a video (youtube link) for illustration purposes. Notice the parameters: number of randomly appearing dots, number or rectangles. I will benchmark with the following parameters:

OFF: Straightforward algorithm as provided by the OP without circle boundary negative checks

ON: Using per-processed (boundary) circles around the rectangles as a first exclusion check

ON + Stack: Creating circle boundaries at run-time within the loop on the stack

ON + Square Distance: Using square distances as a further optimization to avoid taking the more expensive square root algorithm (Pieter Geerkens).

Here is a summary of the various performances of different algorithms by showing the time it takes to iterate through the loop.

The x-axis shows an increased complexity by adding more dots (and thus slowing down the loop). (For example, at 1000 randomly appearing point checks in a confided space with 20 rectangles, the loop iterates and calls the algorithm 20000 times.) The y-axis shows the time it take (ms) to complete the entire loop using a high resolution performance timer. More than 20 ms would be problematic for a decent game as it would not take advantage of the high fps to interpolate a smooth animation and the game may appear thus 'rugged' at times.

Result 1: A pre-processed circular bound algorithm with a fast negative check within the loop improves the performance by 1900% compared to the regular algorithm (5% of the original loop time without a check). The result holds approximately proportional to the number of iterations within a loop, thus it does not matter if we check 10 or 10000 randomly appearing points. Thus, in this illustration one can increase the number of objects safely to 10k without feeling a performance loss.

Result 2: It has been suggested by a previous comment that the algorithm may be faster but memory intensive. However, note that storing a float for the pre-processed circle size takes merely 4 bytes. This should pose no real issue unless the O.P. plans to run simultaneously 100000+ objects. An alternative and memory efficient approach is to calculate the circle maximum size on the stack within the loop and letting it go out of scope with every iteration and thus having practically no memory usage for some unknown price of speed. Indeed, the result shows that this approach is indeed slower than using a pre-processed circle size, but it still shows a considerable performance improvement of around 1150% (i.e. 8% of the original processing time).

Result 3: I further improve the result 1 algorithm by using squared distances instead of actual distances and thus taking an computationally expensive square root operation. This only sligthtly boosts the performance (2400%). (Note: I also try hash tables for pre-processed arrays for square roots approximations with a similar but slightly worse result)

Result 4: I further check moving/colliding the rectangles around; however, this does not change the basic results (as expected) as the logical check remains essentially the same.

Result 5: I vary the number of rectangles and find that the algorithm becomes even more efficient the less crowdy the space is filled (not shown in demo). The result is also somewhat expected, as the probability decreases for a point to appear within tiny space between a circle and the object's boundaries. On the other extreme, I try to increase the number of rectangles too 100 within the same confined tiny space AND vary them dynamically in size at run time at within the loop (sin(iterator)). This still performs extremely well with increase in performance by 570% (or 15% of the original loop time).

Result 6: I test alternative algorithms suggested on here and find a very slight but not significant difference in performance (2%). The interesting and more simple IsLeft algorithm performs very well with a boost of performance by 17% (85% of the original calculation time) but nowhere near the efficiency of a quick negative check algorithm.

My point is to first consider lean design and game logic, especially when dealing with boundaries and collision events. The OPs current algorithm is already fairly efficient and a further optimization is not as critical as optimizing the underlying concept itself. Moreover, it is good to communicate the scope and purpose of the game, as the efficiency of an algorithm critically depends on them.

I suggest to always attempt to benchmark any complex algorithm during the game design stage as merely looking at the plain code may not reveal the truth about actual run-time performance. The suggested algorithm may not be here even necessary, if, for example, one wishes to merely test if the mouse cursor lies within a rectangle or not, or, when the majority of objects are already touching. If the majority of points checks are within the rectangle, the algorithm will be less efficient. (However, then it would be possible to establish an 'inner circle' boundary as a secondary negative check.) Circle/sphere boundary checks are very useful for any decent collision detection of a large number of objects that have naturally some space in between them.

Rec Points  Iter    OFF     ON     ON_Stack     ON_SqrDist  Ileft Algorithm (Wondra)
            (ms)    (ms)    (ms)    (ms)        (ms)        (ms)
20  10      200     0.29    0.02    0.04        0.02        0.17
20  100     2000    2.23    0.10    0.20        0.09        1.69
20  1000    20000   24.48   1.25    1.99        1.05        16.95
20  10000   200000  243.85  12.54   19.61       10.85       160.58

Answer 3

Defining a rectangle with 4 points makes it possible to make a trapezoid. If however, you would define it by x, y, width, height and a rotation around its middle, you could just rotate the point you're checking by the inverse rotation of your rectangle (around the same origin) and then check if it's in the original rectangle.

Answer 4

I didn't have the time to benchmark this, but my suggestion would be to store the transformation matrix that transforms the rectangle into the axis aligned square in the x- and y-range from 0 to 1. In other words store the matrix that transforms one corner of the rectangle into (0,0) and the opposite one into (1,1).

This would of course be more expensive if the rectangle is moved a lot and collision is checked rather seldomly, but if there are a lot more checks than updates to the rectangle it would at least be faster than the original approach of testing against two triangles, as the six dot products would be replaced with one matrix multiplication.

But as always the speed of this algorithm depends a lot on the kind of checks you expect to be performed. If most of the points aren't even close to the rectangle performing a simple distance check (e.g. (point.x - firstCorner.x) > aLargeDistance) might result in a large speedup, while it might even slow things down if almost all of the points are inside the rectangle.

EDIT: This is what my Rectangle-class would look like:

class Rectangle
{
public:
    Matrix3x3 _transform;

    Rectangle()
    {}

    void setCorners(Vector2 p_a, Vector2 p_b, Vector2 p_c)
    {
        // create a matrix from the two edges of the rectangle
        Vector2 edgeX = p_b - p_a;
        Vector2 edgeY = p_c - p_a;

        // and then create the inverse of that matrix because we want to 
        // transform points from world coordinates into "rectangle coordinates".
        float scaling = 1/(edgeX._x*edgeY._y - edgeY._x*edgeX._y);

        _transform._columns[0]._x = scaling * edgeY._y;
        _transform._columns[0]._y = - scaling * edgeX._y;
        _transform._columns[1]._x = - scaling * edgeY._x;
        _transform._columns[1]._y = scaling * edgeX._x;

        // the third column is the translation, which also has to be transformed into "rectangle space"
        _transform._columns[2]._x = -p_a._x * _transform._columns[0]._x - p_a._y * _transform._columns[1]._x;
        _transform._columns[2]._y = -p_a._x * _transform._columns[0]._y - p_a._y * _transform._columns[1]._y;
    }

    bool isInside(Vector2 p_point)
    {
        Vector2 test = _transform.transform(p_point);
        return  (test._x>=0)
                && (test._x<=1)
                && (test._y>=0)
                && (test._y<=1);
    }
};

This is the full listing of my benchmark:

#include <cstdlib>
#include <math.h>
#include <iostream>

#include <sys/time.h>

using namespace std;

class Vector2
{
public:
    float _x;
    float _y;

    Vector2()
    :_x(0)
    ,_y(0)
    {}

    Vector2(float p_x, float p_y)
        : _x (p_x)
        , _y (p_y)
        {}

    Vector2 operator-(const Vector2& p_other) const
    {
        return Vector2(_x-p_other._x, _y-p_other._y);
    }

    Vector2 operator+(const Vector2& p_other) const
    {
        return Vector2(_x+p_other._x, _y+p_other._y);
    }

    Vector2 operator*(float p_factor) const
    {
        return Vector2(_x*p_factor, _y*p_factor);
    }

    static float Dot(Vector2 p_a, Vector2 p_b)
    {
        return (p_a._x*p_b._x + p_a._y*p_b._y);
    }
};

bool PointInTriangle(Vector2 A, Vector2 B, Vector2 C, Vector2 P)
{
 // Compute vectors        
 Vector2 v0 = C - A;
 Vector2 v1 = B - A;
 Vector2 v2 = P - A;

 // Compute dot products
 float dot00 = Vector2::Dot(v0, v0);
 float dot01 = Vector2::Dot(v0, v1);
 float dot02 = Vector2::Dot(v0, v2);
 float dot11 = Vector2::Dot(v1, v1);
 float dot12 = Vector2::Dot(v1, v2);

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

 // Check if point is in triangle
 if(u >= 0 && v >= 0 && (u + v) < 1)
    { return true; } else { return false; }
}


bool PointInRectangle(Vector2 X, Vector2 Y, Vector2 Z, Vector2 W, Vector2 P)
{
 if(PointInTriangle(X,Y,Z,P)) return true;
 if(PointInTriangle(X,Z,W,P)) return true;
 return false;
}

class Matrix3x3
{
public:
    Vector2 _columns[3];

    Vector2 transform(Vector2 p_in)
    {
        return _columns[0] * p_in._x + _columns[1] * p_in._y + _columns[2];
    }
};

class Rectangle
{
public:
    Matrix3x3 _transform;

    Rectangle()
    {}

    void setCorners(Vector2 p_a, Vector2 p_b, Vector2 p_c)
    {
        // create a matrix from the two edges of the rectangle
        Vector2 edgeX = p_b - p_a;
        Vector2 edgeY = p_c - p_a;

        // and then create the inverse of that matrix because we want to 
        // transform points from world coordinates into "rectangle coordinates".
        float scaling = 1/(edgeX._x*edgeY._y - edgeY._x*edgeX._y);

        _transform._columns[0]._x = scaling * edgeY._y;
        _transform._columns[0]._y = - scaling * edgeX._y;
        _transform._columns[1]._x = - scaling * edgeY._x;
        _transform._columns[1]._y = scaling * edgeX._x;

        // the third column is the translation, which also has to be transformed into "rectangle space"
        _transform._columns[2]._x = -p_a._x * _transform._columns[0]._x - p_a._y * _transform._columns[1]._x;
        _transform._columns[2]._y = -p_a._x * _transform._columns[0]._y - p_a._y * _transform._columns[1]._y;
    }

    bool isInside(Vector2 p_point)
    {
        Vector2 test = _transform.transform(p_point);
        return  (test._x>=0)
                && (test._x<=1)
                && (test._y>=0)
                && (test._y<=1);
    }
};

void runTest(float& outA, float& outB)
{
    Rectangle r;
    r.setCorners(Vector2(0,0.5), Vector2(0.5,1), Vector2(0.5,0));

    int numTests = 10000;

    Vector2 points[numTests];

    Vector2 cornerA[numTests];
    Vector2 cornerB[numTests];
    Vector2 cornerC[numTests];
    Vector2 cornerD[numTests];

    bool results[numTests];
    bool resultsB[numTests];

    for (int i=0; i<numTests; ++i)
    {
        points[i]._x = rand() / ((float)RAND_MAX);
        points[i]._y = rand() / ((float)RAND_MAX);

        cornerA[i]._x = rand() / ((float)RAND_MAX);
        cornerA[i]._y = rand() / ((float)RAND_MAX);

        Vector2 edgeA;
        edgeA._x = rand() / ((float)RAND_MAX);
        edgeA._y = rand() / ((float)RAND_MAX);

        Vector2 edgeB;
        edgeB._x = rand() / ((float)RAND_MAX);
        edgeB._y = rand() / ((float)RAND_MAX);

        cornerB[i] = cornerA[i] + edgeA;
        cornerC[i] = cornerA[i] + edgeB;
        cornerD[i] = cornerA[i] + edgeA + edgeB;
    }

    struct timeval start, end;

    gettimeofday(&start, NULL);
    for (int i=0; i<numTests; ++i)
    {
        r.setCorners(cornerA[i], cornerB[i], cornerC[i]);
        results[i] = r.isInside(points[i]);
    }
    gettimeofday(&end, NULL);
    float elapsed = (end.tv_sec - start.tv_sec)*1000;
    elapsed += (end.tv_usec - start.tv_usec)*0.001;
    outA += elapsed;

    gettimeofday(&start, NULL);
    for (int i=0; i<numTests; ++i)
    {
        resultsB[i] = PointInRectangle(cornerA[i], cornerB[i], cornerC[i], cornerD[i], points[i]);
    }
    gettimeofday(&end, NULL);
    elapsed = (end.tv_sec - start.tv_sec)*1000;
    elapsed += (end.tv_usec - start.tv_usec)*0.001;
    outB += elapsed;
}

/*
 * 
 */
int main(int argc, char** argv) 
{
    float a = 0;
    float b = 0;

    for (int i=0; i<5000; i++)
    {
        runTest(a, b);
    }

    std::cout << "Result: " << a << " / " << b << std::endl;

    return 0;
}

The code certainly isn't beautiful, but I don't immediately see any major bugs. With that code I get results that indicate that my solution is about twice as fast if the rectangle is moved between each check. If it doesn't move then my code seems to be more than five times faster.

If you know how the code is going to be used you can even speed it up a little more by separating the transformation and the checks into the two dimensions. For example in a racing game it would probably be faster to check the coordinate first that points in the driving direction, because many obstacles will be in front or behind the car, but hardly any will be right or left of it.