How do you detect where two line segments intersect? [closed]

Question

How do I determine whether or not two lines intersect, and if they do, at what x,y point?

Scour

It might help to think of the edges of the rectangle as separate lines instead of the complete polygon. – Chrestomathy

Answer 1

There’s a nice approach to this problem that uses vector cross products. Define the 2-dimensional vector cross product v × w to be v_x w_y − v_y w_x.

Suppose the two line segments run from p to p + r and from q to q + s. Then any point on the first line is representable as p + t r (for a scalar parameter t) and any point on the second line as q + u s (for a scalar parameter u).

The two lines intersect if we can find t and u such that:

p + t r = q + u s

Cross both sides with s, getting

(p + t r) × s = (q + u s) × s

And since s × s = 0, this means

t (r × s) = (q − p) × s

And therefore, solving for t:

t = (q − p) × s / (r × s)

In the same way, we can solve for u:

(p + t r) × r = (q + u s) × r

u (s × r) = (p − q) × r

u = (p − q) × r / (s × r)

To reduce the number of computation steps, it's convenient to rewrite this as follows (remembering that s × r = − r × s):

u = (q − p) × r / (r × s)

Now there are four cases:

1. If r × s = 0 and (q − p) × r = 0, then the two lines are collinear.

   In this case, express the endpoints of the second segment (q and q + s) in terms of the equation of the first line segment (p + t r):

   t₀ = (q − p) · r / (r · r)

   t₁ = (q + s − p) · r / (r · r) = t₀ + s · r / (r · r)

   If the interval between t₀ and t₁ intersects the interval [0, 1] then the line segments are collinear and overlapping; otherwise they are collinear and disjoint.

   Note that if s and r point in opposite directions, then s · r < 0 and so the interval to be checked is [t₁, t₀] rather than [t₀, t₁].

2. If r × s = 0 and (q − p) × r ≠ 0, then the two lines are parallel and non-intersecting.

3. If r × s ≠ 0 and 0 ≤ t ≤ 1 and 0 ≤ u ≤ 1, the two line segments meet at the point p + t r = q + u s.

4. Otherwise, the two line segments are not parallel but do not intersect.

Credit: this method is the 2-dimensional specialization of the 3D line intersection algorithm from the article "Intersection of two lines in three-space" by Ronald Goldman, published in Graphics Gems, page 304. In three dimensions, the usual case is that the lines are skew (neither parallel nor intersecting) in which case the method gives the points of closest approach of the two lines.

Shellishellie

This is how it's done in every serious 2D library I've seen. It relies on point-line duality in homogenous space which I still have only a tenuous grasp on. – Cosh

This is essentially the same technique as mine, but I use the dot product instead of cross product. In this case, I believe the efficiency is approximentally identical. – Zebrass

Unfortunately, after more tests it seems that this solution also give incorrect results. If you take a spacial case in which one of the considered segments is a point that lies on the other segment the r × s = 0 test says that segments are parallel. – Impasto

As I wrote above, "There are two cases: if (q − p) × r = 0 too, then the lines are collinear, otherwise they never intersect" You've found the collinear case. – Shellishellie

Thanks Garets, for your comment. However I think that your fucntion still doesn't give correct results for some special cases. First case, as I've said before(but was also not correct) if you will take two segments such as S1 = ((0,0),(0,0)) and S2 = ((5,5),(10,7)) your function will result the example as collinear and parallel. S1 is a special case when segment is a point. I've also found an example when false results are given when one of the segments has length = 1 – Impasto

I might be wrong with the length = 1 case (do not have this data to check), but the case which I've given above works as I said. Or I implemented it in a wrong way... Also in this case S1 = ((11,11),(-1,-1)) S2 = ((0,0), (0,10)) function says that those two segments are collinear. But their cross point lies only in (0,0). – Impasto

Converting your example into my notation, I get p=(11,11), r=(-12,-12), q=(0,0), s=(0,10), r×s=-120, t=11/12, u=0. Since r×s is non-zero, the segments are not parallel. – Shellishellie

Hi, shouldn't it be noted that t and u can be negative? So -1 ≤ t ≤ 1 and -1 ≤ u ≤ 1 ? Just for future reference... – Androgynous

@myrkos: No. The first line segment runs "from p to p + r" so when it's represented in parametric terms as "p + tr" then the segment corresponds to 0 ≤ t ≤ 1. Similarly for the other segment. – Shellishellie

@GarethRees How can I know if two line segments are near collinear? – Immunogenetics

@Josh: the angle θ between the vectors r and s is given by the formula θ = arccos(r·s / |r||s|) – Shellishellie

Gareth, I feel like I must be missing something, but how do you divide (a vector) by a vector? Your solutions for t and u end with / (r × s), but (r × s) is a vector, right? A vector (0, 0, rx * sy - ry * sx). And the left-hand side is similarly a vector parallel to the z axis. So... do I just divide the z component by the other z component? Is the formula for t actually |(q − p) × s| / |(r × s)|? – Footed

@LarsH: see the first paragraph. – Shellishellie

@GarethRees: ok thanks, I skipped over that. (Sure enough, I was missing something.) – Footed

@GarethRees Just so I understand this right, if either t or u are greater than 1, then they lines intersect outside of the region we are interested in? – Cahier

@Joe: yes, that's right. – Shellishellie

@GarethRees Brilliant! One last thing, where you say "if (q − p) × r = 0" which means they are collinear, does the same apply for (q - p) x s? I ask because you say two cases, but only give one equation :) (If I'm understanding it correctly, anyway!) – Cahier

@Joe: If r is parallel to s, and also r is parallel to q − p, then s must also be parallel to q − p (being parallel is an equivalence relation). The two cases I referred to are the cases (q − p) × r = 0 and (q − p) × r ≠ 0. – Shellishellie

@GarethRees Thanks for clearing that up! (Even if it was 3 years after your OP!) Much appreciated – Cahier

@GarethRees I hate to contradict what you've said, but I'm running into to situation where I'm having negative values for t or u, see this picture for where it's occuring: i.imgur.com/0owwP.png In this image, the t value for P + tR is going to be negative, since the intersection point is the opposite direction of PR ! – Cahier

The nodes should be P + R etc. – Cahier

@Joe: "the intersection point is on the original pair of line segments if 0 ≤ t ≤ 1 and 0 ≤ u ≤ 1." – Shellishellie

I feel like the reason that line intersection tests are expensive is that it uses at least two divisions. It's pretty clear we'll need to do at least six (or thereabouts) multiplies but the necessity of computing two reciprocals is... well, troubling isn't the right word, but I really want to know if it's possible to avoid division. – Hummingbird

Oh right, r cross s only needs to have its reciprocal computed once. This appears to be about as good as we'll get. – Hummingbird

Want to add, line segment from (1,3) to (5,1) can be represented by the vector (4,-2) – Bouldin

For those interested, here is a simple C# implementation, taking PointF start and end coordinates for lines, that seems to be working: ideone.com/PnPJgb – Harbird

What would be a good epsilon value to test whether r × s = 0 and (q − p) × r = 0 ? – Susquehanna

@peoro: There's no need to use an epsilon here. We're just trying to avoid division by zero, so exact equality is the appropriate test. – Shellishellie

@Mat: I think there's 2 issues w/ your C#. One, for collinearity you need to test BOTH rxs and CmPxr=0, so the latter should be nested in the former. Two, the test for collinearity overlap looks wrong (there is no mention of D, and there are examples where the D location is deciding factor) – Asphyxiate

I put together a JavaScript implementation following @Matt. I made corrections for the errors pointed out by Tekito. – Broaddus

@Mat you need an Epsilon too - if the endpoints are same or one of them lies on the other segment, the algorithm decides by pure chance (= division mistake). Also collinearity test does not work as mentioned above. Try for example (1, 0) (-1, 0)(2, 0)(0, 0) what results false. – Signal

@Broaddus I believe there is a bug in your example; at line 38 it should be var u = uNumerator / -denominator; This is because the denominator for u is (s × r) = -(r × s). – Decompress

This seems to fail for the special case when both lines are collinear, overlap or touch and ||q-p|| is greater than ||r|| and ||s||. E.g. p,p+r = (1,0), (5, 0) and q, q+s = (15,0) and (5, 0). But this can be fixed easily by swapping p and p+r when r.s < 0. – Explicable

Re "this is the magnitude of the 3-dimensional cross product" - since magnitudes are always positive, I think it would be more correct to say that it's the z component of the 3D cross product. – Former

Also, to be totally clear, I would add in [0, 1] when you say The two lines intersect if we can find t and u such that. You do specify this later, but the initial statement could be misleading. – Former

@Broaddus There is an inconsistency in your implementation altough its the best I found here: When two segments share an endpoint you correctly return true (sometimes), but when two collinear segments meet in the shared midle point then you return false: The inconsitency is in the <= vs < usage. Example (0,0) (0,2) (0,2) (0,4): the segments meet at (0,2) but return false: where it returns true in case (0,0) (0,1) (0,1) (1,1). It returns also true in a variation of example 1: (0,0) (0,2) (0,2) (0,3) – Marchand

@Marchand Good catch. I don't think it was a case of <= vs <, because that would break other examples. I added in a special case for equal points as I think it makes the code more readable than trying to include those cases in the exclusive or. Check it out! – Broaddus

@Broaddus A code must be reliable. And yes there is no readable code for the lines segement intersection that works. At least none that works. Finally I took the code from geomalgorithms.com which was ported from the book "Computational Geometry in C". This works. Further I learned that touching segements should be treated as intersetion, and who need should check and remove that afterwords. For the post removal a boolean return is not sufficient, so the intersection (-interval) is needed. – Marchand

@GarethRees Thanks for your nice post! I've made a "terse" C# implementation at codeproject.com, where I try to stick to your algorithm and naming. However, I think that @Explicable has a point with his edge case. Given r=p2-p;s=q2-q, if I were to add something like if (r*s < 0) {Vector.Swap(p, p2); r=p2-p)} I seem to detect the "non-prober" intersection, which I won't if I don't swap. – Dogmatist

I appreciate this post although I think a slight improvement would be to make it so the line segments were represented by (a, a+b), and the other as (m, m+n). That way in the equations it's a little easier to keep straight which vector is which. Or even (p, p+q), (r, r+s). But now I've got the mnemonic "public relations, quicksilver" to keep the pairs straight. – Workout

How do you calculate t0 and t1 in the colinear case? They look like scalars, but on the right-hand side of their equations you get a vector, correct? – Forfend

@mucaho: The dot product of two vectors returns a scalar. – Shellishellie

@GarethRees Aha, there is a difference between · (dot product) and * (multiplication). Seems like t0 / t1 are the scalars obtained by projecting the vector q - p / q + s - p onto r and are expressed in parametric units along r. source – Forfend

Why is this better than solving algebraically? i.e. Set the line equations equal to each other, solve for the intersection, and see whether that point lies on both of the segments. y = m1 * x + b1; y = m2 * x + b2 Set them equal: m1 * x + b1 = m2 * x + b2. Solve for x: x = (b2 - b1) / (m1 - m2). Plug back into one of the originals to find y... etc. – Gaullism

@cp.engr: That's not a complete solution. How are you representing the line segments, and how are you determining whether the intersection lies on both of the segments? And how are you handling the collinear and parallel cases (when m1 = m2)? – Shellishellie

I didn't mean to suggest that my comment contained those details. True, you would need some extra cases. Lines are collinear when m1 = m2 and b1 = b2. Lines are parallel when m1 = m2 and b1 != b2. You'd also need a separate case for vertical lines. Lines segments are defined by their two endpoints. In the normal case, having found the intersection of the lines, test whether that point is on both line segments by checking for xMin <= xIntersect <= xMax and yMin <= yIntersect <= yMax. That said, what advantages does your approach offer over this? – Gaullism

@cp.engr: When you write it out in detail, you'll find that your version of the computation amounts to the same thing as mine — this is not surprising, since it's computing the same result in the same way. My approach has the advantage of working with vectors rather than x and y coordinates — this means that each step has one operation on a vector rather than two operations on x and y values. – Shellishellie

I don't understand case 1. What do you mean by this sentence? "express the endpoints of the second segment (q and q + s) in terms of the equation of the first line segment (p + t r)". Are t0 and t1 the endpoints of the the second segment? How did you derive those equations? ...Thanks for your responses on this old post of yours. :) – Gaullism

I tried to implement this method but hit a problem with the following segments: p = (399.827942f, 174.388153f); p+r = (385.034454f, 188.630402f); q = (375.360474f, 197.943909f); q+s = (365.106628f, 207.815674f). Using 32-bit FP arithmetic, I get rxs = -3.05175781e-5, (q-p)xs = -1.52587891e-5, t=0.5, u=0, indicating that the two segments intersect, which is not true. Above Gareth stated that there's no need to use an epsilon and we are just trying to avoid division by zero. However, this case seems to show that we do need to use an epsilon. – Anaclinal

(1, 1) is (anti)parallel to (−1, −1), not perpendicular. – Shellishellie

TypeScript implementation that not only reports whether the line segments intersect, but also where: gist.github.com/tp/75cb619a7e40e6ad008ef2a6837bbdb2 – Wadsworth

The pedants don't want signed, they want the z-component. "signed magnitude" is not a mathematical concept :) – Forefoot

Is there any way of adding a tolerance to this, to try and catch lines that don't quite intersect (but nearly do)? – Floatstone

I'd like to rewrite this with a better naming choice, u and v for vectors, A and B for points, α and β for scalars (but there's a figure too, can't do it) – Cufic

@caub: I like the names the way they are. – Shellishellie

Answer 2

FWIW, the following function (in C) both detects line intersections and determines the intersection point. It is based on an algorithm in Andre LeMothe's "Tricks of the Windows Game Programming Gurus". It's not dissimilar to some of the algorithm's in other answers (e.g. Gareth's). LeMothe then uses Cramer's Rule (don't ask me) to solve the equations themselves.

I can attest that it works in my feeble asteroids clone, and seems to deal correctly with the edge cases described in other answers by Elemental, Dan and Wodzu. It's also probably faster than the code posted by KingNestor because it's all multiplication and division, no square roots!

I guess there's some potential for divide by zero in there, though it hasn't been an issue in my case. Easy enough to modify to avoid the crash anyway.

// Returns 1 if the lines intersect, otherwise 0. In addition, if the lines 
// intersect the intersection point may be stored in the floats i_x and i_y.
char get_line_intersection(float p0_x, float p0_y, float p1_x, float p1_y, 
    float p2_x, float p2_y, float p3_x, float p3_y, float *i_x, float *i_y)
{
    float s1_x, s1_y, s2_x, s2_y;
    s1_x = p1_x - p0_x;     s1_y = p1_y - p0_y;
    s2_x = p3_x - p2_x;     s2_y = p3_y - p2_y;

    float s, t;
    s = (-s1_y * (p0_x - p2_x) + s1_x * (p0_y - p2_y)) / (-s2_x * s1_y + s1_x * s2_y);
    t = ( s2_x * (p0_y - p2_y) - s2_y * (p0_x - p2_x)) / (-s2_x * s1_y + s1_x * s2_y);

    if (s >= 0 && s <= 1 && t >= 0 && t <= 1)
    {
        // Collision detected
        if (i_x != NULL)
            *i_x = p0_x + (t * s1_x);
        if (i_y != NULL)
            *i_y = p0_y + (t * s1_y);
        return 1;
    }

    return 0; // No collision
}

BTW, I must say that in LeMothe's book, though he apparently gets the algorithm right, the concrete example he shows plugs in the wrong numbers and does calculations wrong. For example:

(4 * (4 - 1) + 12 * (7 - 1)) / (17 * 4 + 12 * 10)

= 844/0.88

= 0.44

That confused me for hours. :(

Armenian

Thanks Gavin - the solution you mention is the one that works best for me too. – Borgerhout

This is a fantastic implementation. Here is a conversion of this function to JavaScript. It returns a list of the intersection point if it exists, and null if the segments do not intersect. – Songful

function getLineIntersection(p0_x, p0_y, p1_x, p1_y, p2_x, p2_y, p3_x, p3_y) { var s1_x, s1_y, s2_x, s2_y; s1_x = p1_x - p0_x; s1_y = p1_y - p0_y; s2_x = p3_x - p2_x; s2_y = p3_y - p2_y; var s, t; s = (-s1_y * (p0_x - p2_x) + s1_x * (p0_y - p2_y)) / (-s2_x * s1_y + s1_x * s2_y); t = ( s2_x * (p0_y - p2_y) - s2_y * (p0_x - p2_x)) / (-s2_x * s1_y + s1_x * s2_y); – Songful

if (s >= 0 && s <= 1 && t >= 0 && t <= 1) { // Collision detected var intX = p0_x + (t * s1_x); var intY = p0_y + (t * s1_y); return [intX, intY]; } return null; // No collision } – Songful

Great answer! But crap. I just spent like 15 minutes porting it to JavaScript only to find that someone already did it. D: – Xerophilous

good algorithm, however fyi it doesn't handle cases where the determinant is 0. (the -s2_x * s1_y + s1_x * s2_y above). If it's 0 (or near 0) the lines are parallel or collinear. If it's collinear then the intersection may be another line segment. – Bio

Calculate the determinant just once, store it, and check, if it's lesser than EPSILON, return 0 right ahead. If not, continue the algorithm with the already calculated determinant. This should handle the parallel/collinear case for most people. – Montpelier

The two division operations can be avoided for speed (division costs more than multiplication); if the lines intersect you need one division, if they do not intersect you need zero. One should first calculate the denominator and stop early if it is zero (possibly adding code to detect colinearity.) Next, instead of calculating s and t directly, test the relationship between the two numerators and the denominator. Only if the lines are confirmed to intersect do you actually need to calculate the value of t (but not s). – Slue

Translated to F#, works fine. pastebin.com/nf56MHP7 – Preece

I did performance testing on all algorithms posted here, and this one is at least twice as fast as any of the others. Thanks for posting! – Irritable

@SargeBorsch In your code sdiv and tdiv always have the same value, and further you (correctly) don't use tdiv. So this does not look fine, remove tdiv. – Marchand

Python: def lineIntersection(self, A, B, C, D): Bx_Ax = B[0] - A[0] By_Ay = B[1] - A[1] Dx_Cx = D[0] - C[0] Dy_Cy = D[1] - C[1] determinant = (-Dx_Cx * By_Ay + Bx_Ax * Dy_Cy) if abs(determinant) < 1e-20: return None s = (-By_Ay * (A[0] - C[0]) + Bx_Ax * (A[1] - C[1])) / determinant t = ( Dx_Cx * (A[1] - C[1]) - Dy_Cy * (A[0] - C[0])) / determinant if s >= 0 and s <= 1 and t >= 0 and t <= 1: return (A[0] + (t * Bx_Ax), A[1] + (t * By_Ay)) return None – Pellegrino

If the line segments are collinear (like (-1, 2) to (1, 2), (1, 2) to (2, 2)), then s and t are NaN and 0 is returned. – Coleen

To catch colinear line, drop the equal sign in the comparison statement. This means that two parallel lines that have point on each other are not classified as intersecting. – Portis

Since there is some ambiguity in the parameters: For get_line_intersection(Ax, Ay, Bx, By, Cx, Cy, Dx, Dy, &Ox, &Oy), the segments tested are A–B vs C–D. – Accustom

@Qwertie, what if determinant < 0? Then you need to test d < s < 0, and if it's > 0, then 0 < s < d. Much simpler to make division. – Strait

For unambiguous completeness, are points 0 and 1 on one line? – Tiddly

I'm going to mess with this code a bit looks nice, but I was curious before I try to figure it out myself. If being able to keep s1, s2, s3, and s4 positive make any difference in reducing the code. Since I know the direction the line is going so I know how to keep the differences positive without additional checks. – Croup

Answer 3

The problem reduces to this question: Do two lines from A to B and from C to D intersect? Then you can ask it four times (between the line and each of the four sides of the rectangle).

Here's the vector math for doing it. I'm assuming the line from A to B is the line in question and the line from C to D is one of the rectangle lines. My notation is that Ax is the "x-coordinate of A" and Cy is the "y-coordinate of C." And "*" means dot-product, so e.g. A*B = Ax*Bx + Ay*By.

E = B-A = ( Bx-Ax, By-Ay )
F = D-C = ( Dx-Cx, Dy-Cy ) 
P = ( -Ey, Ex )
h = ( (A-C) * P ) / ( F * P )

This h number is the key. If h is between 0 and 1, the lines intersect, otherwise they don't. If F*P is zero, of course you cannot make the calculation, but in this case the lines are parallel and therefore only intersect in the obvious cases.

The exact point of intersection is C + F*h.

More Fun:

If h is exactly 0 or 1 the lines touch at an end-point. You can consider this an "intersection" or not as you see fit.

Specifically, h is how much you have to multiply the length of the line in order to exactly touch the other line.

Therefore, If h<0, it means the rectangle line is "behind" the given line (with "direction" being "from A to B"), and if h>1 the rectangle line is "in front" of the given line.

Derivation:

A and C are vectors that point to the start of the line; E and F are the vectors from the ends of A and C that form the line.

For any two non-parallel lines in the plane, there must be exactly one pair of scalar g and h such that this equation holds:

A + E*g = C + F*h

Why? Because two non-parallel lines must intersect, which means you can scale both lines by some amount each and touch each other.

(At first this looks like a single equation with two unknowns! But it isn't when you consider that this is a 2D vector equation, which means this is really a pair of equations in x and y.)

We have to eliminate one of these variables. An easy way is to make the E term zero. To do that, take the dot-product of both sides of the equation using a vector that will dot to zero with E. That vector I called P above, and I did the obvious transformation of E.

You now have:

A*P = C*P + F*P*h
(A-C)*P = (F*P)*h
( (A-C)*P ) / (F*P) = h

Zebrass

So, I see how that could be used to detect if the lines cross or not, but what about determining at what x/y it crosses? – Scour

No problem! They cross at: C + Dh So, in coordinates, at: x-intersect at Cx + Dxh y-intersect at Cy + Dy*h – Zebrass

@Jason Cohen, thank you very much. I would edit your answer and add in your additional information but I can not with my rep level. Do you think you could update your final answer to include the solving for x/y of the intersection so I can mark it as accepted? – Scour

This algorithm is nice. But there is a hole in it as pointed to by Dan @ #563698 & Elemental @ #563698 It would be cool if you would update your answer for future reference. Thanks. – Sweet

Is this algorithm numerically stable? I've tried a similliar aproach and it turned out to give weird results when working on floats. – Oratorian

There seems to be another problem with this algorithm. When it's fed the points A={1, 0} B={2, 0} C={0, 0} D={1,0}, although the line segments clearly touch at an end, FP (and also EQ, in line with the user below's fix) are both 0, thus causing division by 0 to find h and g. Still working on the solution for this one, but I thought the problem was worth pointing out. – Attenuant

Chantz comment must be taken into account in this answer! – Horselaugh

This answer is simply incorrect. Try A={0,0}, B={0,1}, C={0,2} D={2,0} – Gressorial

A + E*g = C + F*h The two lines intersect if and only if the solution to that equation (assuming they are not parallel) has both, g and h between 0 and 1 (in- or exclusive, depending on whether you count touching at an end point). – Suppress

Has this been fixed? – Hypothermia

Answer 4

I have tried to implement the algorithm so elegantly described by Jason above; unfortunately while working though the mathematics in the debugging I found many cases for which it doesn't work.

For example consider the points A(10,10) B(20,20) C(10,1) D(1,10) gives h=.5 and yet it is clear by examination that these segments are no-where near each other.

Graphing this makes it clear that 0 < h < 1 criteria only indicates that the intercept point would lie on CD if it existed but tells one nothing of whether that point lies on AB. To ensure that there is a cross point you must do the symmetrical calculation for the variable g and the requirement for interception is: 0 < g < 1 AND 0 < h < 1

Kalisz

I've been pulling my hair out trying to figure out why the accepted answer wasn't working for me. Thanks so much! – Theresatherese

Also notable that the boundary conditions work in this case (i.e for h=0 or h=1 or g=0 or g=1 the lines 'just' touch – Kalisz

For the people having trouble visualizing the result, I have made an implementation of this in Javascript: jsfiddle.net/ferrybig/eokwL9mp – Largo

Answer 5

Here's an improvement to Gavin's answer. marcp's solution is similar also, but neither postpone the division.

This actually turns out to be a practical application of Gareth Rees' answer as well, because the cross-product's equivalent in 2D is the perp-dot-product, which is what this code uses three of. Switching to 3D and using the cross-product, interpolating both s and t at the end, results in the two closest points between the lines in 3D. Anyway, the 2D solution:

int get_line_intersection(float p0_x, float p0_y, float p1_x, float p1_y, 
    float p2_x, float p2_y, float p3_x, float p3_y, float *i_x, float *i_y)
{
    float s02_x, s02_y, s10_x, s10_y, s32_x, s32_y, s_numer, t_numer, denom, t;
    s10_x = p1_x - p0_x;
    s10_y = p1_y - p0_y;
    s32_x = p3_x - p2_x;
    s32_y = p3_y - p2_y;

    denom = s10_x * s32_y - s32_x * s10_y;
    if (denom == 0)
        return 0; // Collinear
    bool denomPositive = denom > 0;

    s02_x = p0_x - p2_x;
    s02_y = p0_y - p2_y;
    s_numer = s10_x * s02_y - s10_y * s02_x;
    if ((s_numer < 0) == denomPositive)
        return 0; // No collision

    t_numer = s32_x * s02_y - s32_y * s02_x;
    if ((t_numer < 0) == denomPositive)
        return 0; // No collision

    if (((s_numer > denom) == denomPositive) || ((t_numer > denom) == denomPositive))
        return 0; // No collision
    // Collision detected
    t = t_numer / denom;
    if (i_x != NULL)
        *i_x = p0_x + (t * s10_x);
    if (i_y != NULL)
        *i_y = p0_y + (t * s10_y);

    return 1;
}

Basically it postpones the division until the last moment, and moves most of the tests until before certain calculations are done, thereby adding early-outs. Finally, it also avoids the division by zero case which occurs when the lines are parallel.

You also might want to consider using an epsilon test rather than comparison against zero. Lines that are extremely close to parallel can produce results that are slightly off. This is not a bug, it is a limitation with floating point math.

Trina

Fails if some of the points have a value of 0.. that should not happen right? – Escribe

I've made a correction for a bug introduced when deferring the divide. t could be positive when the numer and denom were both negative. – Trina

Does not work if p0-p1 is vertical and p2-p3 is horizontal and the two segments cross. (the first return is executed) – Corporate

The coolinear case has two possibilites: non onverlapping and overlapping. The first shoul return false the second true. In your code this is not tested. it always returns false as most answers here. It's a shame that no solution really seems to work. – Marchand

@fabio really? do you have a test case? – Thierry

It is dangerous to test floats for equality to zero - it avoids division by zero, but division by very small numbers is just as bad. – Walkway

There are several similar calculations within the code, applied to different pieces of data. I wonder whether it's possible to vectorize them for SIMD calculations? All 4 subtractions s10_x = p1_x - p0_x; .... s32_y = p3_y - p2_y can be done in 1 SIMD instruction; a couple of SIMD instructions for simultaneous s_numer = s10_x * s02_y - s10_y * s02_x and t_numer = s32_x * s02_y - s32_y * s02_x; even things like (s_numer < 0) == denomPositive and (t_numer < 0) == denomPositive can probably computed one shot, etc. – Pinprick

This only checks for intersections between the points listed, right? If so, any way to expand this to the entire line defined by the segment passed in? – Vincenz

Can you enlighten me why all these use such vague variable names as s32_y instead of something that describes what it is like point2YDifference? – Coleen

If s_numer == 0, then the test (s_numer < 0) == denomPositive will be false if denom > 0 but true if denom < 0. By symmetry, this does not make sense! – Unilingual

Answer 6

Question C: How do you detect whether or not two line segments intersect?

I have searched for the same topic, and I wasn't happy with the answers. So I have written an article that explains very detailed how to check if two line segments intersect with a lot of images. There is complete (and tested) Java-code.

Here is the article, cropped to the most important parts:

The algorithm, that checks if line segment a intersects with line segment b, looks like this:

What are bounding boxes? Here are two bounding boxes of two line segments:

If both bounding boxes have an intersection, you move line segment a so that one point is at (0|0). Now you have a line through the origin defined by a. Now move line segment b the same way and check if the new points of line segment b are on different sides of line a. If this is the case, check it the other way around. If this is also the case, the line segments intersect. If not, they don't intersect.

Question A: Where do two line segments intersect?

You know that two line segments a and b intersect. If you don't know that, check it with the tools I gave you in "Question C".

Now you can go through some cases and get the solution with 7th grade math (see code and interactive example).

Question B: How do you detect whether or not two lines intersect?

Let's say your point A = (x1, y1), point B = (x2, y2), C = (x_3, y_3), D = (x_4, y_4). Your first line is defined by AB (with A != B), and your second one by CD (with C != D).

function doLinesIntersect(AB, CD) {
    if (x1 == x2) {
        return !(x3 == x4 && x1 != x3);
    } else if (x3 == x4) {
        return true;
    } else {
        // Both lines are not parallel to the y-axis
        m1 = (y1-y2)/(x1-x2);
        m2 = (y3-y4)/(x3-x4);
        return m1 != m2;
    }
}

Question D: Where do two lines intersect?

Check with Question B if they intersect at all.

The lines a and b are defined by two points for each line. You can basically apply the same logic was used in Question A.

Amorino

To be clear, the Question B in this answer is truly about two lines intersecting, not line segments. I'm not complaining; it's not incorrect. Just don't want anyone to be misled. – Marlenmarlena

There' no "question C". And Question D only bounces back onto Question A. – Winglet

Answer 7

The answer once accepted here is incorrect (it has since been unaccepted, so hooray!). It does not correctly eliminate all non-intersections. Trivially it may appear to work but it can fail, especially in the case that 0 and 1 are considered valid for h.

Consider the following case:

Lines at (4,1)-(5,1) and (0,0)-(0,2)

These are perpendicular lines which clearly do not overlap.

A=(4,1)
B=(5,1)
C=(0,0)
D=(0,2)
E=(5,1)-(4,1)=(-1,0)
F=(0,2)-(0,0)=(0,-2)
P=(0,1)
h=((4,1)-(0,0)) dot (0,1) / ((0,-2) dot (0,1)) = 0

According to the above answer, these two line segments meet at an endpoint (values of 0 and 1). That endpoint would be:

(0,0)+(0,-2)*0=(0,0)

So, apparently the two line segments meet at (0,0), which is on line CD, but not on line AB. So what is going wrong? The answer is that the values of 0 and 1 are not valid and only sometimes HAPPEN to correctly predict endpoint intersection. When the extension of one line (but not the other) would meet the line segment, the algorithm predicts an intersection of line segments, but this is not correct. I imagine that by testing starting with AB vs CD and then also testing with CD vs AB, this problem would be eliminated. Only if both fall between 0 and 1 inclusively can they be said to intersect.

I recommend using the vector cross product method if you must predict end-points.

-Dan

Internist

The "accepted" answer can change, so you should call it something else. (In fact, I think it has changed since your comment) – Merissa

Answer 8

Python version of iMalc's answer:

def find_intersection( p0, p1, p2, p3 ) :

    s10_x = p1[0] - p0[0]
    s10_y = p1[1] - p0[1]
    s32_x = p3[0] - p2[0]
    s32_y = p3[1] - p2[1]

    denom = s10_x * s32_y - s32_x * s10_y

    if denom == 0 : return None # collinear

    denom_is_positive = denom > 0

    s02_x = p0[0] - p2[0]
    s02_y = p0[1] - p2[1]

    s_numer = s10_x * s02_y - s10_y * s02_x

    if (s_numer < 0) == denom_is_positive : return None # no collision

    t_numer = s32_x * s02_y - s32_y * s02_x

    if (t_numer < 0) == denom_is_positive : return None # no collision

    if (s_numer > denom) == denom_is_positive or (t_numer > denom) == denom_is_positive : return None # no collision


    # collision detected

    t = t_numer / denom

    intersection_point = [ p0[0] + (t * s10_x), p0[1] + (t * s10_y) ]


    return intersection_point

Caddoan

Remember that you need to make your numbers floats or change line 8 to use denom = float(...) – Pappano

Answer 9

Finding the correct intersection of two line segments is a non-trivial task with lots of edge cases. Here's a well documented, working and tested solution in Java.

In essence, there are three things that can happen when finding the intersection of two line segments:

1. The segments do not intersect

2. There is a unique intersection point

3. The intersection is another segment

NOTE: In the code, I assume that a line segment (x1, y1), (x2, y2) with x1 = x2 and y1 = y2 is a valid line segment. Mathematically speaking, a line segment consists of distinct points, but I am allowing segments to be points in this implementation for completeness.

Code is taken from my github repo

/**
 * This snippet finds the intersection of two line segments.
 * The intersection may either be empty, a single point or the
 * intersection is a subsegment there's an overlap.
 */

import static java.lang.Math.abs;
import static java.lang.Math.max;
import static java.lang.Math.min;

import java.util.ArrayList;
import java.util.List;

public class LineSegmentLineSegmentIntersection {

  // Small epsilon used for double value comparison.
  private static final double EPS = 1e-5;

  // 2D Point class.
  public static class Pt {
    double x, y;
    public Pt(double x, double y) {
      this.x = x; 
      this.y = y;
    }
    public boolean equals(Pt pt) {
      return abs(x - pt.x) < EPS && abs(y - pt.y) < EPS;
    }
  }

  // Finds the orientation of point 'c' relative to the line segment (a, b)
  // Returns  0 if all three points are collinear.
  // Returns -1 if 'c' is clockwise to segment (a, b), i.e right of line formed by the segment.
  // Returns +1 if 'c' is counter clockwise to segment (a, b), i.e left of line
  // formed by the segment.
  public static int orientation(Pt a, Pt b, Pt c) {
    double value = (b.y - a.y) * (c.x - b.x) - 
                   (b.x - a.x) * (c.y - b.y);
    if (abs(value) < EPS) return 0;
    return (value > 0) ? -1 : +1;
  }

  // Tests whether point 'c' is on the line segment (a, b).
  // Ensure first that point c is collinear to segment (a, b) and
  // then check whether c is within the rectangle formed by (a, b)
  public static boolean pointOnLine(Pt a, Pt b, Pt c) {
    return orientation(a, b, c) == 0 && 
           min(a.x, b.x) <= c.x && c.x <= max(a.x, b.x) && 
           min(a.y, b.y) <= c.y && c.y <= max(a.y, b.y);
  }

  // Determines whether two segments intersect.
  public static boolean segmentsIntersect(Pt p1, Pt p2, Pt p3, Pt p4) {

    // Get the orientation of points p3 and p4 in relation
    // to the line segment (p1, p2)
    int o1 = orientation(p1, p2, p3);
    int o2 = orientation(p1, p2, p4);
    int o3 = orientation(p3, p4, p1);
    int o4 = orientation(p3, p4, p2);

    // If the points p1, p2 are on opposite sides of the infinite
    // line formed by (p3, p4) and conversly p3, p4 are on opposite
    // sides of the infinite line formed by (p1, p2) then there is
    // an intersection.
    if (o1 != o2 && o3 != o4) return true;

    // Collinear special cases (perhaps these if checks can be simplified?)
    if (o1 == 0 && pointOnLine(p1, p2, p3)) return true;
    if (o2 == 0 && pointOnLine(p1, p2, p4)) return true;
    if (o3 == 0 && pointOnLine(p3, p4, p1)) return true;
    if (o4 == 0 && pointOnLine(p3, p4, p2)) return true;

    return false;
  }

  public static List<Pt> getCommonEndpoints(Pt p1, Pt p2, Pt p3, Pt p4) {

    List<Pt> points = new ArrayList<>();

    if (p1.equals(p3)) {
      points.add(p1);
      if (p2.equals(p4)) points.add(p2);

    } else if (p1.equals(p4)) {
      points.add(p1);
      if (p2.equals(p3)) points.add(p2);

    } else if (p2.equals(p3)) {
      points.add(p2);
      if (p1.equals(p4)) points.add(p1);

    } else if (p2.equals(p4)) {
      points.add(p2);
      if (p1.equals(p3)) points.add(p1);
    }

    return points;
  }

  // Finds the intersection point(s) of two line segments. Unlike regular line 
  // segments, segments which are points (x1 = x2 and y1 = y2) are allowed.
  public static Pt[] lineSegmentLineSegmentIntersection(Pt p1, Pt p2, Pt p3, Pt p4) {

    // No intersection.
    if (!segmentsIntersect(p1, p2, p3, p4)) return new Pt[]{};

    // Both segments are a single point.
    if (p1.equals(p2) && p2.equals(p3) && p3.equals(p4))
      return new Pt[]{p1};

    List<Pt> endpoints = getCommonEndpoints(p1, p2, p3, p4);
    int n = endpoints.size();

    // One of the line segments is an intersecting single point.
    // NOTE: checking only n == 1 is insufficient to return early
    // because the solution might be a sub segment.
    boolean singleton = p1.equals(p2) || p3.equals(p4);
    if (n == 1 && singleton) return new Pt[]{endpoints.get(0)};

    // Segments are equal.
    if (n == 2) return new Pt[]{endpoints.get(0), endpoints.get(1)};

    boolean collinearSegments = (orientation(p1, p2, p3) == 0) && 
                                (orientation(p1, p2, p4) == 0);

    // The intersection will be a sub-segment of the two
    // segments since they overlap each other.
    if (collinearSegments) {

      // Segment #2 is enclosed in segment #1
      if (pointOnLine(p1, p2, p3) && pointOnLine(p1, p2, p4))
        return new Pt[]{p3, p4};

      // Segment #1 is enclosed in segment #2
      if (pointOnLine(p3, p4, p1) && pointOnLine(p3, p4, p2))
        return new Pt[]{p1, p2};

      // The subsegment is part of segment #1 and part of segment #2.
      // Find the middle points which correspond to this segment.
      Pt midPoint1 = pointOnLine(p1, p2, p3) ? p3 : p4;
      Pt midPoint2 = pointOnLine(p3, p4, p1) ? p1 : p2;

      // There is actually only one middle point!
      if (midPoint1.equals(midPoint2)) return new Pt[]{midPoint1};

      return new Pt[]{midPoint1, midPoint2};
    }

    /* Beyond this point there is a unique intersection point. */

    // Segment #1 is a vertical line.
    if (abs(p1.x - p2.x) < EPS) {
      double m = (p4.y - p3.y) / (p4.x - p3.x);
      double b = p3.y - m * p3.x;
      return new Pt[]{new Pt(p1.x, m * p1.x + b)};
    }

    // Segment #2 is a vertical line.
    if (abs(p3.x - p4.x) < EPS) {
      double m = (p2.y - p1.y) / (p2.x - p1.x);
      double b = p1.y - m * p1.x;
      return new Pt[]{new Pt(p3.x, m * p3.x + b)};
    }

    double m1 = (p2.y - p1.y) / (p2.x - p1.x);
    double m2 = (p4.y - p3.y) / (p4.x - p3.x);
    double b1 = p1.y - m1 * p1.x;
    double b2 = p3.y - m2 * p3.x;
    double x = (b2 - b1) / (m1 - m2);
    double y = (m1 * b2 - m2 * b1) / (m1 - m2);

    return new Pt[]{new Pt(x, y)};
  }

}

Here is a simple usage example:

  public static void main(String[] args) {

    // Segment #1 is (p1, p2), segment #2 is (p3, p4)
    Pt p1, p2, p3, p4;

    p1 = new Pt(-2, 4); p2 = new Pt(3, 3);
    p3 = new Pt(0, 0);  p4 = new Pt(2, 4);
    Pt[] points = lineSegmentLineSegmentIntersection(p1, p2, p3, p4);
    Pt point = points[0];

    // Prints: (1.636, 3.273)
    System.out.printf("(%.3f, %.3f)\n", point.x, point.y);

    p1 = new Pt(-10, 0); p2 = new Pt(+10, 0);
    p3 = new Pt(-5, 0);  p4 = new Pt(+5, 0);
    points = lineSegmentLineSegmentIntersection(p1, p2, p3, p4);
    Pt point1 = points[0], point2 = points[1];

    // Prints: (-5.000, 0.000) (5.000, 0.000)
    System.out.printf("(%.3f, %.3f) (%.3f, %.3f)\n", point1.x, point1.y, point2.x, point2.y);
  }

Curettage

it worked for my geo-coordinates system! thanks! But it is for infinite lines intersection, and I am more looking for finite lines intersection. – Longdistance

Answer 10

Just wanted to mention that a good explanation and explicit solution can be found in the Numeric Recipes series. I've got the 3rd edition and the answer is on page 1117, section 21.4. Another solution with a different nomenclature can be found in a paper by Marina Gavrilova Reliable Line Section Intersection Testing. Her solution is, to my mind, a little simpler.

My implementation is below:

bool NuGeometry::IsBetween(const double& x0, const double& x, const double& x1){
   return (x >= x0) && (x <= x1);
}

bool NuGeometry::FindIntersection(const double& x0, const double& y0, 
     const double& x1, const double& y1,
     const double& a0, const double& b0, 
     const double& a1, const double& b1, 
     double& xy, double& ab) {
   // four endpoints are x0, y0 & x1,y1 & a0,b0 & a1,b1
   // returned values xy and ab are the fractional distance along xy and ab
   // and are only defined when the result is true

   bool partial = false;
   double denom = (b0 - b1) * (x0 - x1) - (y0 - y1) * (a0 - a1);
   if (denom == 0) {
      xy = -1;
      ab = -1;
   } else {
      xy = (a0 * (y1 - b1) + a1 * (b0 - y1) + x1 * (b1 - b0)) / denom;
      partial = NuGeometry::IsBetween(0, xy, 1);
      if (partial) {
         // no point calculating this unless xy is between 0 & 1
         ab = (y1 * (x0 - a1) + b1 * (x1 - x0) + y0 * (a1 - x1)) / denom; 
      }
   }
   if ( partial && NuGeometry::IsBetween(0, ab, 1)) {
      ab = 1-ab;
      xy = 1-xy;
      return true;
   }  else return false;
}

Idio

Does not work for p1=(0,0), p2=(10,0), p3=(9,0), p4=(20,0) – Coles

Depends on your definition of "does not work" I guess. Denom is 0 so it will return false which seems correct to me as they don't intersect. Colinear is not the same as intersecting. – Idio

Answer 11

C and Objective-C

Based on Gareth Rees' answer

const AGKLine AGKLineZero = (AGKLine){(CGPoint){0.0, 0.0}, (CGPoint){0.0, 0.0}};

AGKLine AGKLineMake(CGPoint start, CGPoint end)
{
    return (AGKLine){start, end};
}

double AGKLineLength(AGKLine l)
{
    return CGPointLengthBetween_AGK(l.start, l.end);
}

BOOL AGKLineIntersection(AGKLine l1, AGKLine l2, CGPoint *out_pointOfIntersection)
{
    // https://mcmap.net/q/41388/-how-do-you-detect-where-two-line-segments-intersect-closed

    CGPoint p = l1.start;
    CGPoint q = l2.start;
    CGPoint r = CGPointSubtract_AGK(l1.end, l1.start);
    CGPoint s = CGPointSubtract_AGK(l2.end, l2.start);
    
    double s_r_crossProduct = CGPointCrossProductZComponent_AGK(r, s);
    double t = CGPointCrossProductZComponent_AGK(CGPointSubtract_AGK(q, p), s) / s_r_crossProduct;
    double u = CGPointCrossProductZComponent_AGK(CGPointSubtract_AGK(q, p), r) / s_r_crossProduct;
    
    if(t < 0 || t > 1.0 || u < 0 || u > 1.0)
    {
        if(out_pointOfIntersection != NULL)
        {
            *out_pointOfIntersection = CGPointZero;
        }
        return NO;
    }
    else
    {
        if(out_pointOfIntersection != NULL)
        {
            CGPoint i = CGPointAdd_AGK(p, CGPointMultiply_AGK(r, t));
            *out_pointOfIntersection = i;
        }
        return YES;
    }
}

CGFloat CGPointCrossProductZComponent_AGK(CGPoint v1, CGPoint v2)
{
    return v1.x * v2.y - v1.y * v2.x;
}

CGPoint CGPointSubtract_AGK(CGPoint p1, CGPoint p2)
{
    return (CGPoint){p1.x - p2.x, p1.y - p2.y};
}

CGPoint CGPointAdd_AGK(CGPoint p1, CGPoint p2)
{
    return (CGPoint){p1.x + p2.x, p1.y + p2.y};
}

CGFloat CGPointCrossProductZComponent_AGK(CGPoint v1, CGPoint v2)
{
    return v1.x * v2.y - v1.y * v2.x;
}

CGPoint CGPointMultiply_AGK(CGPoint p1, CGFloat factor)
{
    return (CGPoint){p1.x * factor, p1.y * factor};
}

Many of the functions and structs are private, but you should pretty easy be able to know what's going on. This is public on this repo https://github.com/hfossli/AGGeometryKit/

Escribe

Where is AGPointZero coming from in this code? – Sunwise

@Sunwise updated example to use CGPoint instead – Escribe

Answer 12

Plenty of solutions are available above, but I think below solution is pretty simple and easy to understand.

Two segments Vector AB and Vector CD intersect if and only if

1. The endpoints a and b are on opposite sides of the segment CD.
2. The endpoints c and d are on opposite side of the segment AB.

More specifically a and b are on opposite side of segment CD if and only if exactly one of the two triples a,c,d and b,c,d is in counterclockwise order.

Intersect(a, b, c, d)
 if CCW(a, c, d) == CCW(b, c, d)
    return false;
 else if CCW(a, b, c) == CCW(a, b, d)
    return false;
 else
    return true;

Here CCW represent counterclockwise which returns true/false based on the orientation of the points.

Source : http://compgeom.cs.uiuc.edu/~jeffe/teaching/373/notes/x06-sweepline.pdf Page 2

Prather

I think you should be a bit more specific: how is the CCW test defined? With the sign of the outer product? – Cardioid

Thanks; this pseudo-code allowed for a very straightforward implementation in Scratch; see this project: scratch.mit.edu/projects/129319027 – Subreption

Answer 13

This is working well for me. Taken from here.

 // calculates intersection and checks for parallel lines.  
 // also checks that the intersection point is actually on  
 // the line segment p1-p2  
 Point findIntersection(Point p1,Point p2,  
   Point p3,Point p4) {  
   float xD1,yD1,xD2,yD2,xD3,yD3;  
   float dot,deg,len1,len2;  
   float segmentLen1,segmentLen2;  
   float ua,ub,div;  

   // calculate differences  
   xD1=p2.x-p1.x;  
   xD2=p4.x-p3.x;  
   yD1=p2.y-p1.y;  
   yD2=p4.y-p3.y;  
   xD3=p1.x-p3.x;  
   yD3=p1.y-p3.y;    

   // calculate the lengths of the two lines  
   len1=sqrt(xD1*xD1+yD1*yD1);  
   len2=sqrt(xD2*xD2+yD2*yD2);  

   // calculate angle between the two lines.  
   dot=(xD1*xD2+yD1*yD2); // dot product  
   deg=dot/(len1*len2);  

   // if abs(angle)==1 then the lines are parallell,  
   // so no intersection is possible  
   if(abs(deg)==1) return null;  

   // find intersection Pt between two lines  
   Point pt=new Point(0,0);  
   div=yD2*xD1-xD2*yD1;  
   ua=(xD2*yD3-yD2*xD3)/div;  
   ub=(xD1*yD3-yD1*xD3)/div;  
   pt.x=p1.x+ua*xD1;  
   pt.y=p1.y+ua*yD1;  

   // calculate the combined length of the two segments  
   // between Pt-p1 and Pt-p2  
   xD1=pt.x-p1.x;  
   xD2=pt.x-p2.x;  
   yD1=pt.y-p1.y;  
   yD2=pt.y-p2.y;  
   segmentLen1=sqrt(xD1*xD1+yD1*yD1)+sqrt(xD2*xD2+yD2*yD2);  

   // calculate the combined length of the two segments  
   // between Pt-p3 and Pt-p4  
   xD1=pt.x-p3.x;  
   xD2=pt.x-p4.x;  
   yD1=pt.y-p3.y;  
   yD2=pt.y-p4.y;  
   segmentLen2=sqrt(xD1*xD1+yD1*yD1)+sqrt(xD2*xD2+yD2*yD2);  

   // if the lengths of both sets of segments are the same as  
   // the lenghts of the two lines the point is actually  
   // on the line segment.  

   // if the point isn’t on the line, return null  
   if(abs(len1-segmentLen1)>0.01 || abs(len2-segmentLen2)>0.01)  
     return null;  

   // return the valid intersection  
   return pt;  
 }  

 class Point{  
   float x,y;  
   Point(float x, float y){  
     this.x = x;  
     this.y = y;  
   }  

   void set(float x, float y){  
     this.x = x;  
     this.y = y;  
   }  
 }  

Scour

There are several problems with this code. It can raise an exception due to division by zero; it's slow because it takes square roots; and it sometimes returns false positives because it uses a fudge factor. You can do better than this! – Shellishellie

Okay as a solution but that given by Jason is definitely computationally quicker and avoids a lot of the problems with this solution – Kalisz

Answer 14

I tried some of these answers, but they didnt work for me (sorry guys); after some more net searching I found this.

With a little modification to his code I now have this function that will return the point of intersection or if no intersection is found it will return -1,-1.

    Public Function intercetion(ByVal ax As Integer, ByVal ay As Integer, ByVal bx As Integer, ByVal by As Integer, ByVal cx As Integer, ByVal cy As Integer, ByVal dx As Integer, ByVal dy As Integer) As Point
    '//  Determines the intersection point of the line segment defined by points A and B
    '//  with the line segment defined by points C and D.
    '//
    '//  Returns YES if the intersection point was found, and stores that point in X,Y.
    '//  Returns NO if there is no determinable intersection point, in which case X,Y will
    '//  be unmodified.

    Dim distAB, theCos, theSin, newX, ABpos As Double

    '//  Fail if either line segment is zero-length.
    If ax = bx And ay = by Or cx = dx And cy = dy Then Return New Point(-1, -1)

    '//  Fail if the segments share an end-point.
    If ax = cx And ay = cy Or bx = cx And by = cy Or ax = dx And ay = dy Or bx = dx And by = dy Then Return New Point(-1, -1)

    '//  (1) Translate the system so that point A is on the origin.
    bx -= ax
    by -= ay
    cx -= ax
    cy -= ay
    dx -= ax
    dy -= ay

    '//  Discover the length of segment A-B.
    distAB = Math.Sqrt(bx * bx + by * by)

    '//  (2) Rotate the system so that point B is on the positive X axis.
    theCos = bx / distAB
    theSin = by / distAB
    newX = cx * theCos + cy * theSin
    cy = cy * theCos - cx * theSin
    cx = newX
    newX = dx * theCos + dy * theSin
    dy = dy * theCos - dx * theSin
    dx = newX

    '//  Fail if segment C-D doesn't cross line A-B.
    If cy < 0 And dy < 0 Or cy >= 0 And dy >= 0 Then Return New Point(-1, -1)

    '//  (3) Discover the position of the intersection point along line A-B.
    ABpos = dx + (cx - dx) * dy / (dy - cy)

    '//  Fail if segment C-D crosses line A-B outside of segment A-B.
    If ABpos < 0 Or ABpos > distAB Then Return New Point(-1, -1)

    '//  (4) Apply the discovered position to line A-B in the original coordinate system.
    '*X=Ax+ABpos*theCos
    '*Y=Ay+ABpos*theSin

    '//  Success.
    Return New Point(ax + ABpos * theCos, ay + ABpos * theSin)
End Function

Fugue

Answer 15

There seems to be some interest in Gavin's answer for which cortijon proposed a javascript version in the comments and iMalc provided a version with slightly fewer computations. Some have pointed out shortcomings with various code proposals and others have commented on the efficiency of some code proposals.

The algorithm provided by iMalc via Gavin's answer is the one that I am currently using in a javascript project and I just wanted to provide a cleaned up version here if it may help anyone.

// Some variables for reuse, others may do this differently
var p0x, p1x, p2x, p3x, ix,
    p0y, p1y, p2y, p3y, iy,
    collisionDetected;

// do stuff, call other functions, set endpoints...

// note: for my purpose I use |t| < |d| as opposed to
// |t| <= |d| which is equivalent to 0 <= t < 1 rather than
// 0 <= t <= 1 as in Gavin's answer - results may vary

var lineSegmentIntersection = function(){
    var d, dx1, dx2, dx3, dy1, dy2, dy3, s, t;

    dx1 = p1x - p0x;      dy1 = p1y - p0y;
    dx2 = p3x - p2x;      dy2 = p3y - p2y;
    dx3 = p0x - p2x;      dy3 = p0y - p2y;

    collisionDetected = 0;

    d = dx1 * dy2 - dx2 * dy1;

    if(d !== 0){
        s = dx1 * dy3 - dx3 * dy1;
        if((s <= 0 && d < 0 && s >= d) || (s >= 0 && d > 0 && s <= d)){
            t = dx2 * dy3 - dx3 * dy2;
            if((t <= 0 && d < 0 && t > d) || (t >= 0 && d > 0 && t < d)){
                t = t / d;
                collisionDetected = 1;
                ix = p0x + t * dx1;
                iy = p0y + t * dy1;
            }
        }
    }
};

Birck

I don't understand how you can understand what's going on with lines like t = dx2 * dy3 - dx3 * dy2;... – Coleen

@Supuhstar It has to do with vector math and the definition of dot product and cross product. For example the code you posted represents a cross product operation. It's a way of projecting one line segment onto another to determine where it falls on the other line segment, before the starting point somewhere in the middle or after the line. So t is a normalized value. If it's between 0 and 1 then the two segments intersect. If it's less than 0 or greater than one then they do not. – Birck

@Supuhstar Also note that in order for the projection to find the actual point the result must be scaled. That's where t/d comes in. – Birck

I mean how do you understand what's going on at a glance with variable names like that? Why not something like crossProduct = (line1XDifference * line2YDifference) - (line2XDifference * line1YDifference) and scaledResult = crossProduct / dotProduct? – Coleen

@Supuhstar Ah, I see what you mean. Erm, well I suppose there is really no good reason to speak of beyond obsessing over efficiency, but that is not a very good reason in itself because compilers do a pretty good job of taking most any code you give them and making it as efficient as possible while not changing what is should compute. On the other hand, the names p1x, p1y etc. are meant to describe points by their x and y values, so p1x is an abbreviation for point1x, likewise d1x, in my mind is an abbreviation for the greek letter deltaX or you could say differenceInX. (more) – Birck

@Supuhstar As for t that comes from the convention in algebra to call the parameter t in parametric equations - which is what this method is doing. It represents all of the x's and y's in the equations of both line segments by a single parameter - t. By going through the transformations, using differences and cross products, to put those line segments into a single form, that represented by only t, it allows the method to inspect only t to determine the answer to the question "do they intersect". Similarly d is an abbreviation for determinant. (more) – Birck

@Supuhstar So when you put those insights, steeped in the actual math that is being applied, it makes it clearer what all the variables mean and what they do - which by your quite valid critique, merely allows me to be lazy and make my code very abbreviated without it actually being obfuscated. I hope that the explanation helps, and in light of it you can forgive my terrible variable naming conventions. :-) – Birck

@Supuhstar On a side note, I must say FiM++ looks interesting. Any interest in using deep learning to extract facts, patterns, algorithms, etc by parsing conversational speech, i.e. not strictly structured language? In other words, it would be interesting to try to convert natural language into databases full of facts, patterns or filters, and sets of instructions. – Birck

that'd be neat! Hit me up on email: [email protected] – Coleen

I found the code pretty easy to follow. Thx @Birck – Psalm

Answer 16

I think there is a much much simpler solution for this problem. I came up with another idea today and it seems to work just fine (at least in 2D for now). All you have to do, is to calculate the intersection between two lines, then check if the calculated intersection point is within the boundig boxes of both line segments. If it is, the line segments intersect. That's it.

EDIT:

This is how I calculate the intersection (I don't know anymore where I found this code snippet)

Point3D

comes from

System.Windows.Media.Media3D

public static Point3D? Intersection(Point3D start1, Point3D end1, Point3D start2, Point3D end2) {

        double a1 = end1.Y - start1.Y;
        double b1 = start1.X - end1.X;
        double c1 = a1 * start1.X + b1 * start1.Y;

        double a2 = end2.Y - start2.Y;
        double b2 = start2.X - end2.X;
        double c2 = a2 * start2.X + b2 * start2.Y;

        double det = a1 * b2 - a2 * b1;
        if (det == 0) { // lines are parallel
            return null;
        }

        double x = (b2 * c1 - b1 * c2) / det;
        double y = (a1 * c2 - a2 * c1) / det;

        return new Point3D(x, y, 0.0);
    }

and this is my (simplified for the purpose of the answer) BoundingBox class:

public class BoundingBox {
    private Point3D min = new Point3D();
    private Point3D max = new Point3D();

    public BoundingBox(Point3D point) {
        min = point;
        max = point;
    }

    public Point3D Min {
        get { return min; }
        set { min = value; }
    }

    public Point3D Max {
        get { return max; }
        set { max = value; }
    }

    public bool Contains(BoundingBox box) {
        bool contains =
            min.X <= box.min.X && max.X >= box.max.X &&
            min.Y <= box.min.Y && max.Y >= box.max.Y &&
            min.Z <= box.min.Z && max.Z >= box.max.Z;
        return contains;
    }

    public bool Contains(Point3D point) {
        return Contains(new BoundingBox(point));
    }

}

Relly

Answer 17

This solution may help

public static float GetLineYIntesept(PointF p, float slope)
    {
        return p.Y - slope * p.X;
    }

    public static PointF FindIntersection(PointF line1Start, PointF line1End, PointF line2Start, PointF line2End)
    {

        float slope1 = (line1End.Y - line1Start.Y) / (line1End.X - line1Start.X);
        float slope2 = (line2End.Y - line2Start.Y) / (line2End.X - line2Start.X);

        float yinter1 = GetLineYIntesept(line1Start, slope1);
        float yinter2 = GetLineYIntesept(line2Start, slope2);

        if (slope1 == slope2 && yinter1 != yinter2)
            return PointF.Empty;

        float x = (yinter2 - yinter1) / (slope1 - slope2);

        float y = slope1 * x + yinter1;

        return new PointF(x, y);
    }

Longe

Answer 18

I ported Kris's above answer to JavaScript. After trying numerous different answers, his provided the correct points. I thought I was going crazy that I wasn't getting the points I needed.

function getLineLineCollision(p0, p1, p2, p3) {
    var s1, s2;
    s1 = {x: p1.x - p0.x, y: p1.y - p0.y};
    s2 = {x: p3.x - p2.x, y: p3.y - p2.y};

    var s10_x = p1.x - p0.x;
    var s10_y = p1.y - p0.y;
    var s32_x = p3.x - p2.x;
    var s32_y = p3.y - p2.y;

    var denom = s10_x * s32_y - s32_x * s10_y;

    if(denom == 0) {
        return false;
    }

    var denom_positive = denom > 0;

    var s02_x = p0.x - p2.x;
    var s02_y = p0.y - p2.y;

    var s_numer = s10_x * s02_y - s10_y * s02_x;

    if((s_numer < 0) == denom_positive) {
        return false;
    }

    var t_numer = s32_x * s02_y - s32_y * s02_x;

    if((t_numer < 0) == denom_positive) {
        return false;
    }

    if((s_numer > denom) == denom_positive || (t_numer > denom) == denom_positive) {
        return false;
    }

    var t = t_numer / denom;

    var p = {x: p0.x + (t * s10_x), y: p0.y + (t * s10_y)};
    return p;
}

Ablebodied

Answer 19

I tried lot of ways and then I decided to write my own. So here it is:

bool IsBetween (float x, float b1, float b2)
{
   return ( ((x >= (b1 - 0.1f)) && 
        (x <= (b2 + 0.1f))) || 
        ((x >= (b2 - 0.1f)) &&
        (x <= (b1 + 0.1f))));
}

bool IsSegmentsColliding(   POINTFLOAT lineA,
                POINTFLOAT lineB,
                POINTFLOAT line2A,
                POINTFLOAT line2B)
{
    float deltaX1 = lineB.x - lineA.x;
    float deltaX2 = line2B.x - line2A.x;
    float deltaY1 = lineB.y - lineA.y;
    float deltaY2 = line2B.y - line2A.y;

    if (abs(deltaX1) < 0.01f && 
        abs(deltaX2) < 0.01f) // Both are vertical lines
        return false;
    if (abs((deltaY1 / deltaX1) -
        (deltaY2 / deltaX2)) < 0.001f) // Two parallel line
        return false;

    float xCol = (  (   (deltaX1 * deltaX2) * 
                        (line2A.y - lineA.y)) - 
                    (line2A.x * deltaY2 * deltaX1) + 
                    (lineA.x * deltaY1 * deltaX2)) / 
                 ((deltaY1 * deltaX2) - (deltaY2 * deltaX1));
    float yCol = 0;
    if (deltaX1 < 0.01f) // L1 is a vertical line
        yCol = ((xCol * deltaY2) + 
                (line2A.y * deltaX2) - 
                (line2A.x * deltaY2)) / deltaX2;
    else // L1 is acceptable
        yCol = ((xCol * deltaY1) +
                (lineA.y * deltaX1) -
                (lineA.x * deltaY1)) / deltaX1;

    bool isCol =    IsBetween(xCol, lineA.x, lineB.x) &&
            IsBetween(yCol, lineA.y, lineB.y) &&
            IsBetween(xCol, line2A.x, line2B.x) &&
            IsBetween(yCol, line2A.y, line2B.y);
    return isCol;
}

Based on these two formulas: (I simplified them from equation of lines and other formulas)

Cantonese

Works but try to input this coordinate (if it's colinear/overlaping then it will return false result): PointA1=(0,0) PointA2=(0,2) and PointB1=(0,1) PointB2=(0,5) – Hypethral

@Hypethral Well, that's because the code returns false for parallel lines. I see the problem, however, I still don't know what the function should return as there is an infinite number of answers to it. – Cantonese

Answer 20

This based on Gareth Ree's answer. It also returns the overlap of the line segments if they do. Coded in C++, V is a simple vector class. Where the cross product of two vectors in 2D returns a single scalar. It was tested and passed by my schools automatic testing system.

//Required input point must be colinear with the line
bool on_segment(const V& p, const LineSegment& l)
{
    //If a point is on the line, the sum of the vectors formed by the point to the line endpoints must be equal
    V va = p - l.pa;
    V vb = p - l.pb;
    R ma = va.magnitude();
    R mb = vb.magnitude();
    R ml = (l.pb - l.pa).magnitude();
    R s = ma + mb;
    bool r = s <= ml + epsilon;
    return r;
}

//Compute using vector math
// Returns 0 points if the lines do not intersect or overlap
// Returns 1 point if the lines intersect
//  Returns 2 points if the lines overlap, contain the points where overlapping start starts and stop
std::vector<V> intersect(const LineSegment& la, const LineSegment& lb)
{
    std::vector<V> r;

    //https://mcmap.net/q/41388/-how-do-you-detect-where-two-line-segments-intersect-closed
    V oa, ob, da, db; //Origin and direction vectors
    R sa, sb; //Scalar values
    oa = la.pa;
    da = la.pb - la.pa;
    ob = lb.pa;
    db = lb.pb - lb.pa;

    if (da.cross(db) == 0 && (ob - oa).cross(da) == 0) //If colinear
    {
        if (on_segment(lb.pa, la) && on_segment(lb.pb, la))
        {
            r.push_back(lb.pa);
            r.push_back(lb.pb);
            dprintf("colinear, overlapping\n");
            return r;
        }

        if (on_segment(la.pa, lb) && on_segment(la.pb, lb))
        {
            r.push_back(la.pa);
            r.push_back(la.pb);
            dprintf("colinear, overlapping\n");
            return r;
        }

        if (on_segment(la.pa, lb))
            r.push_back(la.pa);

        if (on_segment(la.pb, lb))
            r.push_back(la.pb);

        if (on_segment(lb.pa, la))
            r.push_back(lb.pa);

        if (on_segment(lb.pb, la))
            r.push_back(lb.pb);

        if (r.size() == 0)
            dprintf("colinear, non-overlapping\n");
        else
            dprintf("colinear, overlapping\n");

        return r;
    }

    if (da.cross(db) == 0 && (ob - oa).cross(da) != 0)
    {
        dprintf("parallel non-intersecting\n");
        return r;
    }

    //Math trick db cross db == 0, which is a single scalar in 2D.
    //Crossing both sides with vector db gives:
    sa = (ob - oa).cross(db) / da.cross(db);

    //Crossing both sides with vector da gives
    sb = (oa - ob).cross(da) / db.cross(da);

    if (0 <= sa && sa <= 1 && 0 <= sb && sb <= 1)
    {
        dprintf("intersecting\n");
        r.push_back(oa + da * sa);
        return r;
    }

    dprintf("non-intersecting, non-parallel, non-colinear, non-overlapping\n");
    return r;
}

Iiette

Answer 21

Here's a basic implementation of a line segment in C#, with corresponding intersection detection code. It requires a 2D vector/point struct called Vector2f, though you can replace this with any other type that has X/Y properties. You could also replace float with double if that suits your needs better.

This code is used in my .NET physics library, Boing.

public struct LineSegment2f
{
    public Vector2f From { get; }
    public Vector2f To { get; }

    public LineSegment2f(Vector2f @from, Vector2f to)
    {
        From = @from;
        To = to;
    }

    public Vector2f Delta => new Vector2f(To.X - From.X, To.Y - From.Y);

    /// <summary>
    /// Attempt to intersect two line segments.
    /// </summary>
    /// <remarks>
    /// Even if the line segments do not intersect, <paramref name="t"/> and <paramref name="u"/> will be set.
    /// If the lines are parallel, <paramref name="t"/> and <paramref name="u"/> are set to <see cref="float.NaN"/>.
    /// </remarks>
    /// <param name="other">The line to attempt intersection of this line with.</param>
    /// <param name="intersectionPoint">The point of intersection if within the line segments, or empty..</param>
    /// <param name="t">The distance along this line at which intersection would occur, or NaN if lines are collinear/parallel.</param>
    /// <param name="u">The distance along the other line at which intersection would occur, or NaN if lines are collinear/parallel.</param>
    /// <returns><c>true</c> if the line segments intersect, otherwise <c>false</c>.</returns>
    public bool TryIntersect(LineSegment2f other, out Vector2f intersectionPoint, out float t, out float u)
    {
        var p = From;
        var q = other.From;
        var r = Delta;
        var s = other.Delta;

        // t = (q − p) × s / (r × s)
        // u = (q − p) × r / (r × s)

        var denom = Fake2DCross(r, s);

        if (denom == 0)
        {
            // lines are collinear or parallel
            t = float.NaN;
            u = float.NaN;
            intersectionPoint = default(Vector2f);
            return false;
        }

        var tNumer = Fake2DCross(q - p, s);
        var uNumer = Fake2DCross(q - p, r);

        t = tNumer / denom;
        u = uNumer / denom;

        if (t < 0 || t > 1 || u < 0 || u > 1)
        {
            // line segments do not intersect within their ranges
            intersectionPoint = default(Vector2f);
            return false;
        }

        intersectionPoint = p + r * t;
        return true;
    }

    private static float Fake2DCross(Vector2f a, Vector2f b)
    {
        return a.X * b.Y - a.Y * b.X;
    }
}

Saltigrade

Answer 22

A C++ program to check if two given line segments intersect

#include <iostream>
using namespace std;

struct Point
{
    int x;
    int y;
};

// Given three colinear points p, q, r, the function checks if
// point q lies on line segment 'pr'
bool onSegment(Point p, Point q, Point r)
{
    if (q.x <= max(p.x, r.x) && q.x >= min(p.x, r.x) &&
        q.y <= max(p.y, r.y) && q.y >= min(p.y, r.y))
       return true;

    return false;
}

// To find orientation of ordered triplet (p, q, r).
// The function returns following values
// 0 --> p, q and r are colinear
// 1 --> Clockwise
// 2 --> Counterclockwise
int orientation(Point p, Point q, Point r)
{
    // See 10th slides from following link for derivation of the formula
    // http://www.dcs.gla.ac.uk/~pat/52233/slides/Geometry1x1.pdf
    int val = (q.y - p.y) * (r.x - q.x) -
              (q.x - p.x) * (r.y - q.y);

    if (val == 0) return 0;  // colinear

    return (val > 0)? 1: 2; // clock or counterclock wise
}

// The main function that returns true if line segment 'p1q1'
// and 'p2q2' intersect.
bool doIntersect(Point p1, Point q1, Point p2, Point q2)
{
    // Find the four orientations needed for general and
    // special cases
    int o1 = orientation(p1, q1, p2);
    int o2 = orientation(p1, q1, q2);
    int o3 = orientation(p2, q2, p1);
    int o4 = orientation(p2, q2, q1);

    // General case
    if (o1 != o2 && o3 != o4)
        return true;

    // Special Cases
    // p1, q1 and p2 are colinear and p2 lies on segment p1q1
    if (o1 == 0 && onSegment(p1, p2, q1)) return true;

    // p1, q1 and p2 are colinear and q2 lies on segment p1q1
    if (o2 == 0 && onSegment(p1, q2, q1)) return true;

    // p2, q2 and p1 are colinear and p1 lies on segment p2q2
    if (o3 == 0 && onSegment(p2, p1, q2)) return true;

     // p2, q2 and q1 are colinear and q1 lies on segment p2q2
    if (o4 == 0 && onSegment(p2, q1, q2)) return true;

    return false; // Doesn't fall in any of the above cases
}

// Driver program to test above functions
int main()
{
    struct Point p1 = {1, 1}, q1 = {10, 1};
    struct Point p2 = {1, 2}, q2 = {10, 2};

    doIntersect(p1, q1, p2, q2)? cout << "Yes\n": cout << "No\n";

    p1 = {10, 0}, q1 = {0, 10};
    p2 = {0, 0}, q2 = {10, 10};
    doIntersect(p1, q1, p2, q2)? cout << "Yes\n": cout << "No\n";

    p1 = {-5, -5}, q1 = {0, 0};
    p2 = {1, 1}, q2 = {10, 10};
    doIntersect(p1, q1, p2, q2)? cout << "Yes\n": cout << "No\n";

    return 0;
}

Rosanarosane

Answer 23

Based on @Gareth Rees answer, version for Python:

import numpy as np

def np_perp( a ) :
    b = np.empty_like(a)
    b[0] = a[1]
    b[1] = -a[0]
    return b

def np_cross_product(a, b):
    return np.dot(a, np_perp(b))

def np_seg_intersect(a, b, considerCollinearOverlapAsIntersect = False):
    # https://mcmap.net/q/41388/-how-do-you-detect-where-two-line-segments-intersect-closed/565282#565282
    # http://www.codeproject.com/Tips/862988/Find-the-intersection-point-of-two-line-segments
    r = a[1] - a[0]
    s = b[1] - b[0]
    v = b[0] - a[0]
    num = np_cross_product(v, r)
    denom = np_cross_product(r, s)
    # If r x s = 0 and (q - p) x r = 0, then the two lines are collinear.
    if np.isclose(denom, 0) and np.isclose(num, 0):
        # 1. If either  0 <= (q - p) * r <= r * r or 0 <= (p - q) * s <= * s
        # then the two lines are overlapping,
        if(considerCollinearOverlapAsIntersect):
            vDotR = np.dot(v, r)
            aDotS = np.dot(-v, s)
            if (0 <= vDotR  and vDotR <= np.dot(r,r)) or (0 <= aDotS  and aDotS <= np.dot(s,s)):
                return True
        # 2. If neither 0 <= (q - p) * r = r * r nor 0 <= (p - q) * s <= s * s
        # then the two lines are collinear but disjoint.
        # No need to implement this expression, as it follows from the expression above.
        return None
    if np.isclose(denom, 0) and not np.isclose(num, 0):
        # Parallel and non intersecting
        return None
    u = num / denom
    t = np_cross_product(v, s) / denom
    if u >= 0 and u <= 1 and t >= 0 and t <= 1:
        res = b[0] + (s*u)
        return res
    # Otherwise, the two line segments are not parallel but do not intersect.
    return None

Yellowlegs

Answer 24

If each side of the rectangle is a line segment, and the user drawn portion is a line segment, then you need to just check the user drawn segment for intersection with the four side line segments. This should be a fairly simple exercise given the start and end points of each segment.

Cyzicus

Note that this was a reasonable answer to the question as originally framed but now that the question has been edited heavily it doesn't make so much sense. – Wines

Answer 25

Based on t3chb0t's answer:

int intersezione_linee(int x1, int y1, int x2, int y2, int x3, int y3, int x4, int y4, int& p_x, int& p_y)
{
   //L1: estremi (x1,y1)(x2,y2) L2: estremi (x3,y3)(x3,y3)
   int d;
   d = (x1-x2)*(y3-y4) - (y1-y2)*(x3-x4);
   if(!d)
       return 0;
   p_x = ((x1*y2-y1*x2)*(x3-x4) - (x1-x2)*(x3*y4-y3*x4))/d;
   p_y = ((x1*y2-y1*x2)*(y3-y4) - (y1-y2)*(x3*y4-y3*x4))/d;
   return 1;
}

int in_bounding_box(int x1, int y1, int x2, int y2, int p_x, int p_y)
{
    return p_x>=x1 && p_x<=x2 && p_y>=y1 && p_y<=y2;

}

int intersezione_segmenti(int x1, int y1, int x2, int y2, int x3, int y3, int x4, int y4, int& p_x, int& p_y)
{
    if (!intersezione_linee(x1,y1,x2,y2,x3,y3,x4,y4,p_x,p_y))
        return 0;

    return in_bounding_box(x1,y1,x2,y2,p_x,p_y) && in_bounding_box(x3,y3,x4,y4,p_x,p_y);
}

Trihedron

Answer 26

I read these algorithm from the book "multiple view geometry"

following text using

' as transpose sign

* as dot product

x as cross product, when using as operator

1. line definition

a point x_vec = (x, y)' lies on the line ax + by + c = 0

we denote L = (a, b, c)', the point as (x, y, 1)' as homogeneous coordinates

the line equation can be written as

(x, y, 1)(a, b, c)' = 0 or x' * L = 0

2. intersection of lines

we have two lines L1=(a1, b1, c1)', L2=(a2, b2, c2)'

assume x is a point, a vector, and x = L1 x L2 (L1 cross product L2).

be careful, x is always a 2D point, please read homogeneous coordinates if you are confused about (L1xL2) is a three elements vector, and x is a 2D coordinates.

according to triple product, we know that

L1 * ( L1 x L2 ) = 0, and L2 * (L1 x L2) = 0, because of L1,L2 co-plane

we substitute (L1xL2) with vector x, then we have L1*x=0, L2*x=0, which means x lie on both L1 and L2, x is the intersection point.

be careful, here x is homogeneous coordinates, if the last element of x is zero, it means L1 and L2 are parallel.

Evieevil

Answer 27

Many answers have wrapped up all the calculations into a single function. If you need to calculate the line slopes, y-intercepts, or x-intercepts for use elsewhere in your code, you'll be making those calculations redundantly. I have separated out the respective functions, used obvious variable names, and commented my code to make it easier to follow. I needed to know if lines intersect infinitely beyond their endpoints, so in JavaScript:

http://jsfiddle.net/skibulk/evmqq00u/

var point_a = {x:0, y:10},
    point_b = {x:12, y:12},
    point_c = {x:10, y:0},
    point_d = {x:0, y:0},
    slope_ab = slope(point_a, point_b),
    slope_bc = slope(point_b, point_c),
    slope_cd = slope(point_c, point_d),
    slope_da = slope(point_d, point_a),
    yint_ab = y_intercept(point_a, slope_ab),
    yint_bc = y_intercept(point_b, slope_bc),
    yint_cd = y_intercept(point_c, slope_cd),
    yint_da = y_intercept(point_d, slope_da),
    xint_ab = x_intercept(point_a, slope_ab, yint_ab),
    xint_bc = x_intercept(point_b, slope_bc, yint_bc),
    xint_cd = x_intercept(point_c, slope_cd, yint_cd),
    xint_da = x_intercept(point_d, slope_da, yint_da),
    point_aa = intersect(slope_da, yint_da, xint_da, slope_ab, yint_ab, xint_ab),
    point_bb = intersect(slope_ab, yint_ab, xint_ab, slope_bc, yint_bc, xint_bc),
    point_cc = intersect(slope_bc, yint_bc, xint_bc, slope_cd, yint_cd, xint_cd),
    point_dd = intersect(slope_cd, yint_cd, xint_cd, slope_da, yint_da, xint_da);

console.log(point_a, point_b, point_c, point_d);
console.log(slope_ab, slope_bc, slope_cd, slope_da);
console.log(yint_ab, yint_bc, yint_cd, yint_da);
console.log(xint_ab, xint_bc, xint_cd, xint_da);
console.log(point_aa, point_bb, point_cc, point_dd);

function slope(point_a, point_b) {
  var i = (point_b.y - point_a.y) / (point_b.x - point_a.x);
  if (i === -Infinity) return Infinity;
  if (i === -0) return 0;
  return i;
}

function y_intercept(point, slope) {
    // Horizontal Line
    if (slope == 0) return point.y;
  // Vertical Line
    if (slope == Infinity)
  {
    // THE Y-Axis
    if (point.x == 0) return Infinity;
    // No Intercept
    return null;
  }
  // Angled Line
  return point.y - (slope * point.x);
}

function x_intercept(point, slope, yint) {
    // Vertical Line
    if (slope == Infinity) return point.x;
  // Horizontal Line
    if (slope == 0)
  {
    // THE X-Axis
    if (point.y == 0) return Infinity;
    // No Intercept
    return null;
  }
  // Angled Line
  return -yint / slope;
}

// Intersection of two infinite lines
function intersect(slope_a, yint_a, xint_a, slope_b, yint_b, xint_b) {
  if (slope_a == slope_b)
  {
    // Equal Lines
    if (yint_a == yint_b && xint_a == xint_b) return Infinity;
    // Parallel Lines
    return null;
  }
  // First Line Vertical
    if (slope_a == Infinity)
  {
    return {
        x: xint_a,
      y: (slope_b * xint_a) + yint_b
    };
  }
  // Second Line Vertical
    if (slope_b == Infinity)
  {
    return {
        x: xint_b,
      y: (slope_a * xint_b) + yint_a
    };
  }
  // Not Equal, Not Parallel, Not Vertical
  var i = (yint_b - yint_a) / (slope_a - slope_b);
  return {
    x: i,
    y: (slope_a * i) + yint_a
  };
}

Subconscious