Circle line-segment collision detection algorithm?
Solution 1
Taking
- E is the starting point of the ray,
- L is the end point of the ray,
- C is the center of sphere you're testing against
- r is the radius of that sphere
Compute:
d = L - E ( Direction vector of ray, from start to end )
f = E - C ( Vector from center sphere to ray start )
Then the intersection is found by..
Plugging:
P = E + t * d
This is a parametric equation:
Px = Ex + tdx
Py = Ey + tdy
into
(x - h)2 + (y - k)2 = r2
(h,k) = center of circle.
Note: We've simplified the problem to 2D here, the solution we get applies also in 3D
to get:
- Expand x2 - 2xh + h2 + y2 - 2yk + k2 - r2 = 0
- Plug
x = ex + tdx
y = ey + tdy
( ex + tdx )2 - 2( ex + tdx )h + h2 + ( ey + tdy )2 - 2( ey + tdy )k + k2 - r2 = 0 - Explode ex2 + 2extdx + t2dx2 - 2exh - 2tdxh + h2 + ey2 + 2eytdy + t2dy2 - 2eyk - 2tdyk + k2 - r2 = 0
- Group t2( dx2 + dy2 ) + 2t( exdx + eydy - dxh - dyk ) + ex2 + ey2 - 2exh - 2eyk + h2 + k2 - r2 = 0
- Finally,
t2( d · d ) + 2t( e · d - d · c ) + e · e - 2( e · c ) + c · c - r2 = 0
Where d is the vector d and · is the dot product. - And then, t2( d · d ) + 2t( d · ( e - c ) ) + ( e - c ) · ( e - c ) - r2 = 0
- Letting f = e - c t2( d · d ) + 2t( d · f ) + f · f - r2 = 0
So we get:
t2 * (d · d) + 2t*( f · d ) + ( f · f - r2 ) = 0
So solving the quadratic equation:
float a = d.Dot( d ) ;
float b = 2*f.Dot( d ) ;
float c = f.Dot( f ) - r*r ;
float discriminant = b*b-4*a*c;
if( discriminant < 0 )
{
// no intersection
}
else
{
// ray didn't totally miss sphere,
// so there is a solution to
// the equation.
discriminant = sqrt( discriminant );
// either solution may be on or off the ray so need to test both
// t1 is always the smaller value, because BOTH discriminant and
// a are nonnegative.
float t1 = (-b - discriminant)/(2*a);
float t2 = (-b + discriminant)/(2*a);
// 3x HIT cases:
// -o-> --|--> | | --|->
// Impale(t1 hit,t2 hit), Poke(t1 hit,t2>1), ExitWound(t1<0, t2 hit),
// 3x MISS cases:
// -> o o -> | -> |
// FallShort (t1>1,t2>1), Past (t1<0,t2<0), CompletelyInside(t1<0, t2>1)
if( t1 >= 0 && t1 <= 1 )
{
// t1 is the intersection, and it's closer than t2
// (since t1 uses -b - discriminant)
// Impale, Poke
return true ;
}
// here t1 didn't intersect so we are either started
// inside the sphere or completely past it
if( t2 >= 0 && t2 <= 1 )
{
// ExitWound
return true ;
}
// no intn: FallShort, Past, CompletelyInside
return false ;
}
Solution 2
No one seems to consider projection, am I completely off track here?
Project the vector AC
onto AB
. The projected vector, AD
, gives the new point D
.
If the distance between D
and C
is smaller than (or equal to) R
we have an intersection.
Like this:
Community Edit:
For anyone stumbling across this post later and wondering how such an algorithm can be implemented, here is a general implementation written in JavaScript using common vector manipulation functions.
/**
* Returns the distance from line segment AB to point C
*/
function distanceSegmentToPoint(A, B, C) {
// Compute vectors AC and AB
const AC = sub(C, A);
const AB = sub(B, A);
// Get point D by taking the projection of AC onto AB then adding the offset of A
const D = add(proj(AC, AB), A);
const AD = sub(D, A);
// D might not be on AB so calculate k of D down AB (aka solve AD = k * AB)
// We can use either component, but choose larger value to reduce the chance of dividing by zero
const k = Math.abs(AB.x) > Math.abs(AB.y) ? AD.x / AB.x : AD.y / AB.y;
// Check if D is off either end of the line segment
if (k <= 0.0) {
return Math.sqrt(hypot2(C, A));
} else if (k >= 1.0) {
return Math.sqrt(hypot2(C, B));
}
return Math.sqrt(hypot2(C, D));
}
For this implementation I used a couple common vector manipulation functions that you are likely to already have provided in whatever environment you might be working in. However if you do not already have these functions available, here is how they can be implemented.
// Define some common functions for working with vectors
const add = (a, b) => ({x: a.x + b.x, y: a.y + b.y});
const sub = (a, b) => ({x: a.x - b.x, y: a.y - b.y});
const dot = (a, b) => a.x * b.x + a.y * b.y;
const hypot2 = (a, b) => dot(sub(a, b), sub(a, b));
// Function for projecting some vector a onto b
function proj(a, b) {
const k = dot(a, b) / dot(b, b);
return {x: k * b.x, y: k * b.y};
}
Solution 3
I would use the algorithm to compute the distance between a point (circle center) and a line (line AB). This can then be used to determine the intersection points of the line with the circle.
Let say we have the points A, B, C. Ax and Ay are the x and y components of the A points. Same for B and C. The scalar R is the circle radius.
This algorithm requires that A, B and C are distinct points and that R is not 0.
Here is the algorithm
// compute the euclidean distance between A and B
LAB = sqrt( (Bx-Ax)²+(By-Ay)² )
// compute the direction vector D from A to B
Dx = (Bx-Ax)/LAB
Dy = (By-Ay)/LAB
// the equation of the line AB is x = Dx*t + Ax, y = Dy*t + Ay with 0 <= t <= LAB.
// compute the distance between the points A and E, where
// E is the point of AB closest the circle center (Cx, Cy)
t = Dx*(Cx-Ax) + Dy*(Cy-Ay)
// compute the coordinates of the point E
Ex = t*Dx+Ax
Ey = t*Dy+Ay
// compute the euclidean distance between E and C
LEC = sqrt((Ex-Cx)²+(Ey-Cy)²)
// test if the line intersects the circle
if( LEC < R )
{
// compute distance from t to circle intersection point
dt = sqrt( R² - LEC²)
// compute first intersection point
Fx = (t-dt)*Dx + Ax
Fy = (t-dt)*Dy + Ay
// compute second intersection point
Gx = (t+dt)*Dx + Ax
Gy = (t+dt)*Dy + Ay
}
// else test if the line is tangent to circle
else if( LEC == R )
// tangent point to circle is E
else
// line doesn't touch circle
Solution 4
Okay, I won't give you code, but since you have tagged this algorithm, I don't think that will matter to you. First, you have to get a vector perpendicular to the line.
You will have an unknown variable in y = ax + c
( c
will be unknown )
To solve for that, Calculate it's value when the line passes through the center of the circle.
That is,
Plug in the location of the circle center to the line equation and solve for c
.
Then calculate the intersection point of the original line and its normal.
This will give you the closest point on the line to the circle.
Calculate the distance between this point and the circle center (using the magnitude of the vector).
If this is less than the radius of the circle - voila, we have an intersection!
Solution 5
Another method uses the triangle ABC area formula. The intersection test is simpler and more efficient than the projection method, but finding the coordinates of the intersection point requires more work. At least it will be delayed to the point it is required.
The formula to compute the triangle area is : area = bh/2
where b is the base length and h is the height. We chose the segment AB to be the base so that h is the shortest distance from C, the circle center, to the line.
Since the triangle area can also be computed by a vector dot product we can determine h.
// compute the triangle area times 2 (area = area2/2)
area2 = abs( (Bx-Ax)*(Cy-Ay) - (Cx-Ax)(By-Ay) )
// compute the AB segment length
LAB = sqrt( (Bx-Ax)² + (By-Ay)² )
// compute the triangle height
h = area2/LAB
// if the line intersects the circle
if( h < R )
{
...
}
UPDATE 1 :
You could optimize the code by using the fast inverse square root computation described here to get a good approximation of 1/LAB.
Computing the intersection point is not that difficult. Here it goes
// compute the line AB direction vector components
Dx = (Bx-Ax)/LAB
Dy = (By-Ay)/LAB
// compute the distance from A toward B of closest point to C
t = Dx*(Cx-Ax) + Dy*(Cy-Ay)
// t should be equal to sqrt( (Cx-Ax)² + (Cy-Ay)² - h² )
// compute the intersection point distance from t
dt = sqrt( R² - h² )
// compute first intersection point coordinate
Ex = Ax + (t-dt)*Dx
Ey = Ay + (t-dt)*Dy
// compute second intersection point coordinate
Fx = Ax + (t+dt)*Dx
Fy = Ay + (t+dt)*Dy
If h = R then the line AB is tangent to the circle and the value dt = 0 and E = F. The point coordinates are those of E and F.
You should check that A is different of B and the segment length is not null if this may happen in your application.
Related videos on Youtube
Mizipzor
A crazy coder in a passionate hunt for greater wisdom. I take great interest in anything involving math and algorithms. Especially path finding, artificial life, cellular automata and emergent behavior.
Updated on March 23, 2022Comments
-
Mizipzor about 2 years
I have a line from A to B and a circle positioned at C with the radius R.
What is a good algorithm to use to check whether the line intersects the circle? And at what coordinate along the circles edge it occurred?
-
Jason S almost 15 yearsHmm. One question: are you talking about the infinite line through A and B, or the finite line segment from A to B?
-
Mizipzor almost 15 yearsIn this case, its the finite line segment. Is "line" called something else depending on if its finite or infinite?
-
chmike almost 15 yearsIs there a performance requirement ? Should it be a fast method ?
-
Mizipzor almost 15 yearsAt this point, no, all the algorithms here that Ive tried doesnt slow the application down noticeably.
-
MestreLion almost 10 years@Mizipzor yes, they are called something else: line segments. If you just say "line" it's implied an infinite one.
-
xakepp35 about 6 yearsIf you are doing collision detection in 2D game, then you will like to see my answer
-
-
Martijn almost 15 yearsIt's in 2D, not 3D; as you say, this doesn't really matter
-
Mizipzor almost 15 yearsThat was, in fact, what I wanted. I want the theory, a google search of line-circle collision algorithm turns up only code as far as I can see.
-
Martin almost 15 yearsI'm no mathematician so I thought I'd be better outlining a general approach and leaving it to others to figure out specific maths (although I does look rather trivial)
-
Jason S almost 15 yearsmathworld.wolfram.com/Point-LineDistance3-Dimensional.html and mathworld.wolfram.com/Point-LineDistance2-Dimensional.html are better & from a more reputable site
-
Martin almost 15 yearsI explained a little better about the closest point, and linked to mathworld instead of pbourke :)
-
Mizipzor almost 15 yearsOk, c is unknown in your equation, but what is "a"? The other answers seem to refer to that variable as "alpha" and "t". Although, this is just a linear function (y=kx+m), quite basic math, so I suddenly feel a bit rusty. Isnt k also unknown? Or do you mean we can assume m=0 and solve k? Wouldnt then m (that is, c) always be zero for our solved k?
-
a_m0d almost 15 yearsOh, sorry - I am using the simple equation of a line with a gradient and offset (the cartesian equation). I assumed that you were storing the line as such an equation - in which case you use the negative of the gradient for k. If you don't have the line stored like this, you could calculate k as (y2-y1)/(x2-x1)
-
a_m0d almost 15 yearsWe don't assume that m is zero; we calculate the gradient first (so that the equation of the line then looks like y=2x+m as an example), and then once we have the gradient we can solve for m by plugging in the center of the circle for y and x.
-
Mizipzor almost 15 yearsSeems to work if I do straight copy and paste, but Im looking to understand it to. In (x-h)^2+(y-k)^2=r^2 what is h and k? Is k to constant by which the line/ray increases on y over x? And what is t? Looking at the code it seems you have assumed its 1 (so its just "removed"). Do these formulas have a name or something? Maybe I can look them up in detail on Wolfram.
-
bobobobo almost 15 yearsh and k are the center of the circle that you're intersecting against. t is the parameter of the line equation. In the code, t1 and t2 are the solutions. t1 and t2 tell you "how far along the ray" the intersection happened.
-
Mizipzor almost 15 yearsI like the simplicity in this method. Maybe I could adapt some of the surrounding code to not need the actual collision point itself, Ill see what happens if I use A or B rather than the calculated point in between.
-
chmike almost 15 yearsI suggest another method below (with code) that uses the triangle area to compute the distance of the circle center to the line.
-
chmike almost 15 yearsVery elegant solution that requires only one square root computation and should in principle work in 3D. What would be the 3D equation and how do we find the roots ?
-
chmike almost 15 yearsOk, got it. The dot product is simply computed over the three elements (x,y,z) vectors. I will switch my code to this algorithm.
-
Admin over 14 yearsCan someone fully explain how we're going from (x - h)^2 + (y - k)^2 = r^2 to t^2 * (d DOT d) + (2t * (f DOT d)) + (f DOT f - r^2) = 0? I tried plugging in the parametric equivalent of x and y and got something I couldn't reduce. Thanks very much.
-
jeffreynolte over 13 yearsThere are many details to take into consideration: does D lie between AB? Is C perpendicular distance to the line larger than the radius? All of these involve magnitude of vector, i.e. square root.
-
jeffreynolte over 13 yearsYou missed one case since we are talking about a line segment: when the segment ends in the circle.
-
Juozas Kontvainis over 13 years@ADB actually my algorithm only works for infinite lines only, not line segments. There are many cases that it does not handle with line segments.
-
Ivaylo Strandjev over 12 yearsThe sign of b must be the opposite.
-
Ben about 12 yearsGood idea, but how do you then compute the two intersection points?
-
Admin almost 12 years+1 Awesome explanation! But I think this assumes a line, not a line segment. So, if the nearest point on this line to the circle's center wasn't between points A and B, it would still be counted.
-
Andrew Stacey almost 12 years(I just came across someone using this method and getting false negatives. Looking at the code, I saw that it is only testing one intersection point to see if it is on the line segment, but it needs to test both. I edited in the second test as well.)
-
bobobobo almost 12 years@AndrewStacey Good edit, I had thought that would be straight forward. I'd also like to add that the smaller of t1/t2 is going to be closer to the ray start position, but I left it out.
-
Andrew Stacey almost 12 years@bobobobo I figured one way of reading "t2 is the other intersection point" as "you also need to test t2" but as I'd just seen someone not realise that I thought it worth editing in (make the internet a better place, and all that). Your second sentence is misleading as the smaller of t1/t2 is not necessarily on the line (indeed, this was the problem I saw: the edge case of when one point is on the line segment and the other isn't; I know it's obvious but sometimes the obvious isn't obvious until it's been pointed out).
-
Derek 朕會功夫 almost 12 years
P = E + t * d
What ist
? -
Stonetip almost 12 yearst = Dx*(Cx-Ax) + Dy*(Cy-Ax) should read t = Dx*(Cx-Ax) + Dy*(Cy-Ay)
-
bobobobo over 11 yearsIn C#/XNA you can use
Ray.Intersects(BoundingSphere)
-
Nicolas Mommaerts about 11 yearsNot sure why, but the code doesn't seem to work for the Impale case. It does when I add if t1 <= 0 && t1 >= -1 && t2 <= 0 && t2 >= -1 as true condition, but then it also gives a false positive on one side of the finite line, when the circle is on the "infinite" part. I don't understand the math yet, but copy/pasters, beware.
-
ericsoco almost 11 yearsjust edited -- first line calculates triangle area using a cross product, not a dot product. verified with code here: stackoverflow.com/questions/2533011/…
-
ericsoco almost 11 yearsnote also, the first half of this answer tests for intersection with a line, not a line segment (as asked in the question).
-
Prabindh about 10 yearsIf possible, also post with some sample values so that we can quickly grasp the flow.
-
Rialgar over 9 yearsI just came up with a better way of understanding the final formula. (P-C)^2 equals r^2, since P is on the circle; C+f equals E by definition; E + td equals P by definition; so C + f + td equals P; so (f + td)^2 equals r^2; so t^2*d^2 + t*2*fd equals r^2; so we get our quadratic formula
-
Spider over 9 years@Mizipzor. bit late, but I dont agree with this point. Coz, it does not consider that AB segment line position at all. Your point is true when only AB is perpendicularly located against CD point! Consider AB segment line where CD is not perpendicular to AB at all. What is your opinion on that?
-
Admin over 9 years@Ben the two intersection points (let's call them
P1
&P2
,Px
in general) are points onAB
that are exactlyr
distance fromC
- simply calculate |CD| distance and from Pyth's you know that |CD|^2+|PxD|^2 = r^2 - making this trivial if you know how to calculate line's parametric equation (which I assume you should). -
Admin over 9 years@Spider it doesn't matter. In general, since this is a variant of sphere-line intersection problem, Mizipzor's strategy is perfectly valid.
CD
is a projection, it is perpendicular by definition. -
Prashant over 9 yearsif any line that is not intersecting the circle and its both point p1 and p2 are inside circle. in this case how your algorithm work??
-
chmike over 9 yearsYou have to test t-dt and t+dt. If t-dt < 0, than p1 is inside the circle. If t+dt > 1, than p2 is inside the circle. This is true if LEC < R of course.
-
Andrew over 9 yearsIt's an old question, but there's a good resource on this and related algorithms at this website: paulbourke.net/geometry/pointlineplane
-
punchcard about 9 yearsThanks. I liked this pgm comments as explanation because there was no use of the words "dot product" since my math is rusty. However t and dt are not between 0..1 so while changing this to python I changed t to be divided by LAB**2. My understanding is the first division by LAB is to project the center of the circle to the line AB, and the second division by LAB is to normalize it into the range 0..1. Also the dt needs to be divided by LAB so it is also normalized. Thus "if (t-dt >=0.0)" first intersection exists "if (t+dt <= 1.0)" second intersection exists. This worked with testing.
-
Mahmoud Ayman almost 9 years@chmike Can I please know why did you write in the part where LEC<R:
t+dt
once andt-dt
in another point? Why did you one time plus and the other subtract? -
chmike almost 9 yearsBecause the intersection point with the circle are at "distance"
t+dt
andt-dt
on the line.t
is the point on the line closest to the center of the circle. The intersection points with the circle are at a symmetric distance fromt
. The intersection points are at "distances"t-dt
andt+dt
. I quoted distance because it it's not the euclidian distance. To get the euclidian distance fromA
wheret=0
, you have to multiply the value byLAB
. -
Matt W over 8 yearsThis code calculates intersections on the line's vector. How can you determine if the intersection occurs on the line between the end points rather than outside the line?
-
Marconius over 8 years@Matt W You mean "How to determine if the intersection occurs outside the end points of the line section AB"? Just think about t as a measure of distance along the line. Point A is at
t=0
. Point B att=LAB
. When both intersection points (t1=t-td
andt2=t+td
) have negative values than intersections are outside the section (behind point A looking from the section side of the point). When t1 and t2 are bigger than LAB then they are outside too (this time behind B point). Intersection t1 (or t2) occurs between A and B only when t1 (or t2) it is between 0 and LAB. -
Gino over 8 yearssuggestions: First, have it handle the case where line segment is vertical (i.e. has infinite slope). Secondly, have it only return those points that actually fall within the range of the original line segment -- i believe it's happily returning all points that fall on the infinite line, even if those points lie outside of the line segment.
-
byJeevan about 8 yearsWith
underRadical
extra ')' -
arkon almost 8 years"E is the starting point of the ray, L is the end point of the ray" Rays only have one point, extending infinitely along a vector. What you're describing is a line segment.
-
sliders_alpha over 7 yearsWhat language is this? I need to read the documentation of Impale(), Past(), etc.
-
msumme over 7 yearsThe original question was about line-segments, not circle-line intersection, which is a much easier problem.
-
TommasoF about 7 yearsHave a look at this link math.stackexchange.com/questions/311921/… for theoretic explanation.
-
ickydime almost 7 yearsFor those asking what is t... You need to implement the code example where he finds t1 and t2. Then drop that into: P = E + t * d For example at // ExitWound you could then write (C#): var intersection = new Vector3(E.x + t1 * d.x, secondPoint.y, E.y + t1 * d.y); Worked great. Thanks!
-
Minestrone-Soup over 6 years@NicolasMommaerts really old question, but if people are having the same problem as Nicolas as I did, it's because you did f = C - E instead of f = E - C
-
ShawnFeatherly over 6 yearsA great explanation of this answer: scratchapixel.com/lessons/3d-basic-rendering/…
-
Anthony almost 6 yearsNote that this algorithm requires the line and the circle lie on the same plane, which is always the case in 2D.
-
unlut almost 6 yearsIsn't this answer in incomplete? It finds whether a circle and line intersect, not a line segment. In my opinion correct check should be something like:
(projectionPointInCircle and projectionPointOnLine) or endPointsInCircle
Last check is necessary since line segment may penetrate circle and finish before it goes past center. endPoints refer to starting and ending points of line segment. -
fishinear over 5 yearsThere is a small fault in the equations. t is divided by LAB, and therefore lies between 0 (at A) and 1 (at B). But dt is not divided by LAB. Either both t and dt should be divided by LAB, or neither should.
-
chmike over 5 years@fishinear you are right. I should have divided by LAB. I'll fix it.
-
fishinear over 5 years@chmike Actually, I think I was wrong in my remark. You divide D by LAB, not t, which means that t is actually between 0 and LAB (the comment in the code is incorrect in that respect). That means again that dt should not be divided by LAB either, and your original code was correct.
-
fishinear over 5 years@chmike Suppose Ax == Bx == Cx, then Dx is zero and Dy is equal to 1. That means that t must be equal to By - Ay at point B, not 1.
-
chmike over 5 years@fishinear you are right. If I want t to be 0 at A and 1 at B, then Dx = Bx - Ax and Dy = By - Ay. In this case tDx + Ax = x and tDy+Ay = y is correct. There is another error. The scalar product between AC and AB computed as Dx*(Cx-Ax)+Dy*(Cy-Ay) = LACLABcos θ. The value t for the point E is cos θ. I must thus divide the result of the scalar product by (LAC*LAB). This can be verified when C = B. In this case t = 1. We can then compute the coordinates of the point E. This solution to the question is not very efficient.
-
fishinear over 5 years@chmike Actually, in your original code, you divide Dx and Dy by LAB. That means that Dx*(Cx-Ax)+Dy*(Cy-Ay) is already divided by LAB, and therefore is equal to LAC*cos θ = LAE. And that is exactly what t needs to be if it needs to equal LAB at B. So your original code was actually correct, just not the comment that t is 1 at B. I still think it is an excellent solution.
-
chmike over 5 years@fishinear It should now be fixed. You're first comment confused me.
-
GMartigny over 5 years@unlut You're right. It even depend on the kind of calculation you use to project C on AB. I made a simulation if someone want to play around : desmos.com/calculator/gb7hk2dfwp
-
Mike about 5 yearsNote: This works well for lines, but does not work for line segments.
-
bobobobo over 4 years@Derek朕會功夫
t
is a parameter. InP = E + t * d
,P
is some point in space we will get to BY STARTING atE
, then movingt
units along the direction vectord
.t
is a variable that describes how far along the direction vectord
we are -
Chris almost 4 yearsIf the line segment intersects the circle but only slightly so it does not pass its center point, won't this function return false when it should return true? The t value might be outside the range 0..1.
-
WDUK almost 3 yearsAn image to go with this would be very helpful to understand it more
-
WDUK almost 3 years@GMartigny this does not work for line segments. It works for infinite lines, but the Q was about line segments wasn't it ? Image here showing it would detect intersection but the line segment doesn't overlap the circle. White line is the projection: imgur.com/IQh4lcJ
-
Denis Rizov about 2 yearsIn other words
P = E + d * t
is linear interpolation along the line (akat
is in range[0, 1]
). If the discriminant is>= 0
that means the infinite ray intersects with the circle. If we want to find if the non-finite line intersects with the circle, we need to find a solution of the equation in range[0, 1]