I know that there are already many questions about how to check whether a line segment contains a point or not, but my question is specifically about which method will most probably be most efficient for a performance-critical application.
Looking around at various SO answers, there seem to be three general approaches:
Using distance:
float distance(Vector2 a, Vector2 b)
{
return std::sqrtf(std::powf((a.x - b.x), 2.f) + std::powf((a.y - b.y), 2.f));
}
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
return distance(lineP, point) + distance(point,lineQ) == distance(lineP, lineQ);
}
This one I write off quickly since three square root operations doesn't exactly scream fast.
Calculating the slope of the line segment and then checking if the point is within the bounds:
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
float slope = (lineP.y - lineQ.y) / (lineP.x - lineQ.x);
float c = lineP.y - lineP.x * slope;
float minX = std::min(lineP.x, lineQ.x);
float maxX = (minX == lineP.x) ? lineQ.x : lineP.x;
return point.x * slope + c == point.y && point.x >= minX && point.x <= maxX;
}
This is the first approach that I came up with when writing such a function but I am unsure if it is the most efficient in terms of performance.
Using cross product and dot product:
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
float crossproduct = (point.y - lineP.y) * (lineQ.x - lineP.x) - (point.x - lineP.x) * (lineQ.y - lineP.y)
if (fabs(crossproduct) > std::numeric_limits<float>::epsilon())
{
return false;
}
float dotproduct = (point.x - lineP.x) * (lineQ.x - lineP.x) + (point.y - lineP.y)*(lineQ.y - lineP.y)
if (dotproduct < 0.f)
{
return false;
}
float squaredlengthba = (lineQ.x - lineP.x)*(lineQ.x - lineP.x) + (lineQ.y - lineP.y)*(lineQ.y - lineP.y)
if (dotproduct > squaredlengthba)
{
return false;
}
return true;
}
This one seems longer than the line slope solution but is different in that it doesn't use division, and I've been told that division is 10 times slower than multiplication. So maybe this is the superior solution?
Or is there another, faster way?
You should absolutely profile, ideally in real conditions, instead of trying to analyze the code to guess which version is fastest. However, there are several other things you should consider here.
Your distance function calls std::powf to square a value X instead of just doing X * X. I'm unsure whether a compiler will optimize this, but it just looks wrong to me anyway.
Your containsPoint function calls distance three times, which means three square roots as you mention. However, you can easily simplify this to just one square root by finding a solution to the equation sqrt(a) + sqrt(b) = sqrt(c). Solving for c, you can find c = a + b + 2*sqrt(a*b)
With this, the code becomes:
float distanceSquared(Vector2 a, Vector2 b)
{
auto dx = (a.x - b.x);
auto dy = (a.y - b.y);
return dx*dx + dy*dy;
}
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
auto d2PPoint = distanceSquared(lineP, point);
auto d2QPoint = distanceSquared(lineQ, point);
auto d2PQ = distanceSquared(lineP, lineQ);
return d2PQ == d2PPoint + d2QPoint + 2*std::sqrt(d2PPoint*d2QPoint);
}
Also, you probably want to use an approximate equality instead of a straight comparison with ==.
This method simply does not work for a vertical line, though you can of course add a special case for it. It introduces a branch, but it's an unlikely one. Also, to check if x is between a and b without knowing whether a < b or a > b, you can multiply the differences: check if (x-a)*(x-b) <= 0. Additionally, you can do without the division by just rescaling things.
The code becomes:
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
float dxPQ = lineP.x - lineQ.x;
if (dxPQ == 0.f) {
return (point.y - lineP.y) * (point.y - lineQ.y) <= 0.f;
}
float dyPQ = lineP.y - lineQ.y;
float c = lineP.y * dxPQ - lineP.x * dyPQ;
bool xBetween = ((point.x - lineP.x) * (point.x - lineQ.x) <= 0.f);
return xBetween && point.x * dyPQ + c * dxPQ == point.y * dxPQ;
}
Again, you probably want to use an approximate equality for the last part.
Using the dot product is a rather strange way of checking whether the point is within bounds. Instead, you can check the x and y coordinates separately:
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
bool xBetween = ((point.x - lineP.x) * (point.x - lineQ.x) <= 0.f);
bool yBetween = ((point.y - lineP.y) * (point.y - lineQ.y) <= 0.f);
float crossproduct = (point.y - lineP.y) * (lineQ.x - lineP.x) - (point.x - lineP.x) * (lineQ.y - lineP.y);
return xBetween && yBetween && std::abs(crossproduct) <= std::numeric_limits<float>::epsilon();
}
For some reason, you used an approximate equality here, but not in the other versions. In this version, you should use an epsilon proportional to the length of the line, since the cross product will be as well. It does introduce a square root, unless you compare the square of the cross product (in which case you can also drop the std::abs). Something like this:
bool containsPoint(Vector2 lineP, Vector2 lineQ, Vector2 point)
{
bool xBetween = ((point.x - lineP.x) * (point.x - lineQ.x) <= 0.f);
bool yBetween = ((point.y - lineP.y) * (point.y - lineQ.y) <= 0.f);
if (!xBetween || !yBetween) { // early return or can be moved to the end
return false;
}
float dxPQ = lineQ.x - lineP.x;
float dyPQ = lineQ.y - lineP.y;
float lineLengthSquared = dxPQ*dxPQ + dyPQ*dyPQ;
float crossproduct = (point.y - lineP.y) * dxPQ - (point.x - lineP.x) * dyPQ;
float tolerance = std::numeric_limits<float>::epsilon();
return crossproduct * crossproduct <= lineLengthSquared * tolerance * tolerance;
}
This version looks to me like the most correct and the most understandable (maybe put a few comments here and there). I would not be surprised if it's also the most performant, with or without the branch.