Find the intersection of two lines

Find the intersection of two lines is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Finding the intersection of two lines that are in the same plane is an important topic in collision detection.

Task

Find the point of intersection of two lines in 2D. The first line passes though (4.0,0.0) and (6.0,10.0). The second line passes though (0.0,3.0) and (10.0,7.0).

ALGOL 68

Using "school maths". <lang algol68>BEGIN

   # mode to hold a point #
   MODE POINT = STRUCT( REAL x, y );
   # mode to hold a line expressed as y = mx + c #
   MODE LINE  = STRUCT( REAL m, c );
   # returns the line that passes through p1 and p2 #
   PROC find line = ( POINT p1, p2 )LINE:
        IF x OF p1 = x OF p2 THEN
            # the line is vertical                 #
            LINE( 0, x OF p1 )
        ELSE
            # the line is not vertical             #
            REAL m = ( y OF p1 - y OF p2 ) / ( x OF p1 - x OF p2 );
            LINE( m, y OF p1 - ( m * x OF p1 ) )
        FI # find line # ;

   # returns the intersection of two lines - the lines must be distinct and not parallel #
   PRIO INTERSECTION = 5;
   OP   INTERSECTION = ( LINE l1, l2 )POINT:
        BEGIN
            REAL x = ( c OF l2 - c OF l1 ) / ( m OF l1 - m OF l2 );
            POINT( x, ( m OF l1 * x ) + c OF l1 )
        END # INTERSECTION # ;

   # find the intersection of the lines as per the task #
   POINT i = find line( POINT( 4.0, 0.0 ), POINT( 6.0, 10.0 ) )
             INTERSECTION find line( ( 0.0, 3.0 ), ( 10.0, 7.0 ) );
   print( ( fixed( x OF i, -8, 4 ), fixed( y OF i, -8, 4 ), newline ) )

END</lang>

Output:

  5.0000  5.0000

C#

<lang csharp>using System; using System.Drawing; public class Program {

   static PointF FindIntersection(PointF s1, PointF e1, PointF s2, PointF e2) {
       float d = (s1.X - e1.X) * (s2.Y - e2.Y) - (s1.Y - e1.Y) * (s2.X - e2.X);
       //If lines are parallel, the result will be (NaN, NaN).
       if (d == 0) return new PointF(float.NaN, float.NaN);
       float a = s1.X * e1.Y - s1.Y * e1.X;
       float b = s2.X * e2.Y - s2.Y * e2.X;
       float x = (a * (s2.X - e2.X) - (s1.X - e1.X) * b) / d;
       float y = (a * (s2.Y - e2.Y) - (s1.Y - e1.Y) * b) / d;
       return new PointF(x, y);
   }

   static void Main() {
       Func<float, float, PointF> p = (x, y) => new PointF(x, y);
       Console.WriteLine(FindIntersection(p(4f, 0f), p(6f, 10f), p(0f, 3f), p(10f, 7f)));
       Console.WriteLine(FindIntersection(p(0f, 0f), p(1f, 1f), p(1f, 2f), p(4f, 5f)));
   }

}</lang>

Output:

{X=5, Y=5}
{X=NaN, Y=NaN}

C++

<lang cpp>#include <iostream>

1. include <cmath>
2. include <assert.h>

using namespace std;

/** Calculate determinant of matrix: [a b] [c d]

• /

inline double Det(double a, double b, double c, double d) { return a*d - b*c; }

///Calculate intersection of two lines. ///\return true if found, false if not found or error bool LineLineIntersect(double x1, double y1, //Line 1 start double x2, double y2, //Line 1 end double x3, double y3, //Line 2 start double x4, double y4, //Line 2 end double &ixOut, double &iyOut) //Output { double detL1 = Det(x1, y1, x2, y2); double detL2 = Det(x3, y3, x4, y4); double x1mx2 = x1 - x2; double x3mx4 = x3 - x4; double y1my2 = y1 - y2; double y3my4 = y3 - y4;

double xnom = Det(detL1, x1mx2, detL2, x3mx4); double ynom = Det(detL1, y1my2, detL2, y3my4); double denom = Det(x1mx2, y1my2, x3mx4, y3my4); if(denom == 0.0)//Lines don't seem to cross { ixOut = NAN; iyOut = NAN; return false; }

ixOut = xnom / denom; iyOut = ynom / denom; if(!isfinite(ixOut) || !isfinite(iyOut)) //Probably a numerical issue return false;

return true; //All OK }

int main() { // **Simple crossing diagonal lines**

//Line 1 double x1=4.0, y1=0.0; double x2=6.0, y2=10.0;

//Line 2 double x3=0.0, y3=3.0; double x4=10.0, y4=7.0;

double ix = -1.0, iy = -1.0; bool result = LineLineIntersect(x1, y1, x2, y2, x3, y3, x4, y4, ix, iy); cout << "result " << result << "," << ix << "," << iy << endl;

double eps = 1e-6; assert(result == true); assert(fabs(ix - 5.0) < eps); assert(fabs(iy - 5.0) < eps);

}</lang>

Output:

result 1,5,5

Haskell

<lang haskell>type Line = (Point, Point)

type Point = (Float, Float)

-- INTERSECTION OF TWO LINES -------------------------------------------------- maybeIntersection :: Line -> Line -> Maybe Point maybeIntersection ((ax, ay), (bx, by)) ((px, py), (qx, qy)) =

 let (pqDX, abDX) = (px - qx, ax - bx)
     (pqDY, abDY) = (py - qy, ay - by)
     determinant = abDX * pqDY - abDY * pqDX
     f pq ab =
       ((((ax * by) - (ay * bx)) * pq) - 
        (((px * qy) - (py * qx)) * ab)) /
       determinant
 in case determinant of
      0 -> Nothing
      _ -> Just (f pqDX abDX, f pqDY abDY)

-- TEST ----------------------------------------------------------------------- ab :: Line ab = ((4.0, 0.0), (6.0, 10.0))

pq :: Line pq = ((0.0, 3.0), (10.0, 7.0))

interSection :: Maybe Point interSection = maybeIntersection ab pq

main :: IO () main =

 putStrLn $
 case interSection of
   Nothing -> "(parallel lines – no intersection)"
   Just x -> show x</lang>

Output:

(5.0,5.0)

JavaScript

Translation of: Haskell

ES6

<lang JavaScript>(() => {

   // INTERSECTION OF TWO LINES ----------------------------------------------

   // maybeIntersection :: Line -> Line -> Maybe Point
   const maybeIntersection = ([
       [ax, ay],
       [bx, by]
   ], [
       [px, py],
       [qx, qy]
   ]) => {
       const
           abDX = ax - bx,
           pqDX = px - qx,
           abDY = ay - by,
           pqDY = py - qy,
           determinant = abDX * pqDY - abDY * pqDX;

       return determinant !== 0 ? ({
           nothing: false,
           just: (() => {
               const
                   abD = ax * by - ay * bx,
                   pqD = px * qy - py * qx;
               return ap(
                   [([pq, ab]) =>
                       (abD * pq - ab * pqD) / determinant
                   ], [
                       [pqDX, abDX],
                       [pqDY, abDY]
                   ]
               );
           })()
       }) : {
           nothing: true
       };
   };

   // GENERIC FUNCTIONS ------------------------------------------------------

   // A list of functions applied to a list of arguments
   // <*> :: [(a -> b)] -> [a] -> [b]
   const ap = (fs, xs) => //
       [].concat.apply([], fs.map(f => //
           [].concat.apply([], xs.map(x => [f(x)]))));

   // show :: a -> String
   const show = x => JSON.stringify(x); //, null, 2);

   // TEST -------------------------------------------------------------------

   // ab :: Line
   const ab = [
       [4.0, 0.0],
       [6.0, 10.0]
   ];

   // pq :: Line
   const pq = [
       [0.0, 3.0],
       [10.0, 7.0]
   ];

   // intersection :: Maybe Point
   const intersection = maybeIntersection(ab, pq);

   return intersection.nothing ? (
       '(Parallel lines – no intersection)'
   ) : show(intersection.just);

})();</lang>

Output:

[5,5]

Kotlin

Translation of: C#

<lang scala>// version 1.1.2

class PointF(val x: Float, val y: Float) {

   override fun toString() = "{$x, $y}"

}

class LineF(val s: PointF, val e: PointF)

fun findIntersection(l1: LineF, l2: LineF): PointF {

   val a1 = l1.e.y - l1.s.y
   val b1 = l1.s.x - l1.e.x
   val c1 = a1 * l1.s.x + b1 * l1.s.y

   val a2 = l2.e.y - l2.s.y
   val b2 = l2.s.x - l2.e.x
   val c2 = a2 * l2.s.x + b2 * l2.s.y

   val delta = a1 * b2 - a2 * b1
   // If lines are parallel, intersection point will contain infinite values
   return PointF((b2 * c1 - b1 * c2) / delta, (a1 * c2 - a2 * c1) / delta)

}

fun main(args: Array<String>) {

   var l1: LineF
   var l2: LineF
   l1 = LineF(PointF(4f, 0f), PointF(6f, 10f))
   l2 = LineF(PointF(0f, 3f), PointF(10f, 7f))
   println(findIntersection(l1, l2))
   l1 = LineF(PointF(0f, 0f), PointF(1f, 1f))
   l2 = LineF(PointF(1f, 2f), PointF(4f, 5f))
   println(findIntersection(l1, l2))

}</lang>

Output:

{5.0, 5.0}
{-Infinity, -Infinity}

Lua

Translation of: C#

<lang lua>function intersection (s1, e1, s2, e2)

 local d = (s1.x - e1.x) * (s2.y - e2.y) - (s1.y - e1.y) * (s2.x - e2.x)
 local a = s1.x * e1.y - s1.y * e1.x
 local b = s2.x * e2.y - s2.y * e2.x
 local x = (a * (s2.x - e2.x) - (s1.x - e1.x) * b) / d
 local y = (a * (s2.y - e2.y) - (s1.y - e1.y) * b) / d
 return x, y

end

local line1start, line1end = {x = 4, y = 0}, {x = 6, y = 10} local line2start, line2end = {x = 0, y = 3}, {x = 10, y = 7} print(intersection(line1start, line1end, line2start, line2end))</lang>

Output:

5       5

Perl

Translation of: C#

If warning are enabled the second print will issue a warning since we are trying to print out an undef

<lang perl> sub intersect {

 my ($x1, $y1, $x2, $y2, $x3, $y3, $x4, $y4) = @_;
 my $a1 = $y2 - $y1;
 my $b1 = $x1 - $x2;
 my $c1 = $a1 * $x1 + $b1 * $y1;
 my $a2 = $y4 - $y3;
 my $b2 = $x3 - $x4;
 my $c2 = $a2 * $x3 + $b2 * $y3;
 my $delta = $a1 * $b2 - $a2 * $b1;
 return (undef, undef) if $delta == 0;
 # If delta is 0, i.e. lines are parallel then the below will fail
 my $ix = ($b2 * $c1 - $b1 * $c2) / $delta;
 my $iy = ($a1 * $c2 - $a2 * $c1) / $delta;
 return ($ix, $iy);

}

my ($ix, $iy) = intersect(4, 0, 6, 10, 0, 3, 10, 7); print "$ix $iy\n"; ($ix, $iy) = intersect(0, 0, 1, 1, 1, 2, 4, 5); print "$ix $iy\n"; </lang>

Perl 6

Works with: Rakudo version 2016.11

Translation of: zkl

<lang perl6>sub intersection (Real $ax, Real $ay, Real $bx, Real $by,

                 Real $cx, Real $cy, Real $dx, Real $dy ) {

   sub term:<|AB|> { determinate($ax, $ay, $bx, $by) }
   sub term:<|CD|> { determinate($cx, $cy, $dx, $dy) }

   my $ΔxAB = $ax - $bx;
   my $ΔyAB = $ay - $by;
   my $ΔxCD = $cx - $dx;
   my $ΔyCD = $cy - $dy;

   my $x-numerator = determinate( |AB|, $ΔxAB, |CD|, $ΔxCD );
   my $y-numerator = determinate( |AB|, $ΔyAB, |CD|, $ΔyCD );
   my $denominator = determinate( $ΔxAB, $ΔyAB, $ΔxCD, $ΔyCD );

   return 'Lines are parallel' if $denominator == 0;

   ($x-numerator/$denominator, $y-numerator/$denominator);

}

sub determinate ( Real $a, Real $b, Real $c, Real $d ) { $a * $d - $b * $c }

1. TESTING

say 'Intersection point: ', intersection( 4,0, 6,10, 0,3, 10,7 ); say 'Intersection point: ', intersection( 4,0, 6,10, 0,3, 10,7.1 ); say 'Intersection point: ', intersection( 0,0, 1,1, 1,2, 4,5 );</lang>

Output:

Intersection point: (5 5)
Intersection point: (5.010893 5.054466)
Intersection point: Lines are parallel

Python

<lang python>from __future__ import print_function from shapely.geometry import LineString

if __name__=="__main__": line1 = LineString([(4.0,0.0), (6.0,10.0)]) line2 = LineString([(0.0,3.0), (10.0,7.0)]) print (line1.intersection(line2))</lang>

Output:

POINT (5 5)

Racket

Translation of: C++

<lang racket>#lang racket/base (define (det a b c d) (- (* a d) (* b c))) ; determinant

(define (line-intersect ax ay bx by cx cy dx dy) ; --> (values x y)

 (let* ((det.ab (det ax ay bx by))
        (det.cd (det cx cy dx dy))
        (abΔx (- ax bx))
        (cdΔx (- cx dx))
        (abΔy (- ay by))
        (cdΔy (- cy dy))
        (xnom (det det.ab abΔx det.cd cdΔx))
        (ynom (det det.ab abΔy det.cd cdΔy))
        (denom (det abΔx abΔy cdΔx cdΔy)))
   (when (zero? denom)
     (error 'line-intersect "parallel lines"))
   (values (/ xnom denom) (/ ynom denom))))

(module+ test (line-intersect 4 0 6 10

                             0 3 10 7))</lang>

Output:

5
5

REXX

version 1

Naive implementation. To be improved for parallel lines and degenerate lines such as y=5 or x=8. <lang rexx>/* REXX */ Parse Value '(4.0,0.0)' With '(' xa ',' ya ')' Parse Value '(6.0,10.0)' With '(' xb ',' yb ')' Parse Value '(0.0,3.0)' With '(' xc ',' yc ')' Parse Value '(10.0,7.0)' With '(' xd ',' yd ')'

Say 'The two lines are:' Say 'yab='ya-xa*((yb-ya)/(xb-xa))'+x*'||((yb-ya)/(xb-xa)) Say 'ycd='yc-xc*((yd-yc)/(xd-xc))'+x*'||((yd-yc)/(xd-xc))

x=((yc-xc*((yd-yc)/(xd-xc)))-(ya-xa*((yb-ya)/(xb-xa))))/,

                        (((yb-ya)/(xb-xa))-((yd-yc)/(xd-xc)))

Say 'x='||x

       y=ya-xa*((yb-ya)/(xb-xa))+x*((yb-ya)/(xb-xa))

Say 'yab='y Say 'ycd='yc-xc*((yd-yc)/(xd-xc))+x*((yd-yc)/(xd-xc)) Say 'Intersection: ('||x','y')'</lang>

Output:

The two lines are:
yab=-20.0+x*5
ycd=3.0+x*0.4
x=5
yab=5.0
ycd=5.0
Intersection: (5,5.0)

version 2

complete implementation taking care of all possibilities.
Variables are named after the Austrian notation for a straight line: y=k*x+d <lang>say ggx1('4.0 0.0 6.0 10.0 0.0 3.0 10.0 7.0') say ggx1('0.0 0.0 0.0 10.0 0.0 3.0 10.0 7.0') say ggx1('0.0 0.0 0.0 10.0 0.0 3.0 10.0 7.0') say ggx1('0.0 0.0 0.0 1.0 1.0 0.0 1.0 7.0') say ggx1('0.0 0.0 0.0 0.0 0.0 3.0 10.0 7.0') say ggx1('0.0 0.0 3.0 3.0 0.0 0.0 6.0 6.0') say ggx1('0.0 0.0 3.0 3.0 0.0 1.0 6.0 7.0') Exit

ggx1: Procedure /*---------------------------------------------------------------------

• find the intersection of the lines AB and CD
• --------------------------------------------------------------------*/

Parse Arg xa ya xb yb xc yc xd yd Say 'A=('xa'/'ya') B=('||xb'/'yb') C=('||xc'/'yc') D=('||xd'/'yd')' res= If xa=xb Then Do /* AB is a vertical line */

 k1='*'                            /* slope is infinite             */
 x1=xa                             /* intersection's x is xa        */
 If ya=yb Then                     /* coordinates are equal         */
   res='Points A and B are identical' /* special case               */
 End

Else Do /* AB is not a vertical line */

 k1=(yb-ya)/(xb-xa)                /* compute the slope of AB       */
 d1=ya-k1*xa                /* and its intersection with the y-axis */
 End

If xc=xd Then Do /* CD is a vertical line */

 k2='*'                            /* slope is infinite             */
 x2=xc                             /* intersection's x is xc        */
 If yc=yd Then                     /* coordinates are equal         */
   res='Points C and D are identical' /* special case               */
 End

Else Do /* CD is not a vertical line */

 k2=(yd-yc)/(xd-xc)                /* compute the slope of CD       */
 d2=yc-k2*xc                /* and its intersection with the y-axis */
 End

If res= Then Do /* no special case so far */

 If k1='*' Then Do                 /* AB is vertical                */
   If k2='*' Then Do               /* CD is vertical                */
     If x1=x2 Then                 /* and they are identical        */
       res='Lines AB and CD are identical'
     Else                          /* not identical                 */
       res='Lines AB and CD are parallel'
     End
   Else Do
     x=x1                          /* x is taken from AB            */
     y=k2*x+d2                     /* y is computed from CD         */
     End
   End
 Else Do                           /* AB is not verical             */
   If k2='*' Then Do               /* CD is vertical                */
     x=x2                          /* x is taken from CD            */
     y=k1*x+d1                     /* y is computed from AB         */
     End
   Else Do                         /* AB and CD are not vertical    */
     If k1=k2 Then Do              /* identical slope               */
       If d1=d2 Then               /* same intersection with x=0    */
         res='Lines AB and CD are identical'
       Else                        /* otherwise                     */
         res='Lines AB and CD are parallel'
       End
     Else Do                       /* finally the normal case       */
       x=(d2-d1)/(k1-k2)           /* compute x                     */
       y=k1*x+d1                   /* and y                         */
       End
     End
   End
 End
 If res= Then                    /* not any special case          */
   res='Intersection is ('||x'/'y')'  /* that's the result          */
 Return '  -->' res</lang>

Output:

A=(4.0/0.0) B=(6.0/10.0) C=(0.0/3.0) D=(10.0/7.0)
  --> Intersection is (5/5.0)
A=(0.0/0.0) B=(0.0/10.0) C=(0.0/3.0) D=(10.0/7.0)
  --> Intersection is (0.0/3.0)
A=(0.0/0.0) B=(0.0/10.0) C=(0.0/3.0) D=(10.0/7.0)
  --> Intersection is (0.0/3.0)
A=(0.0/0.0) B=(0.0/1.0) C=(1.0/0.0) D=(1.0/7.0)
  --> Lines AB and CD are parallel
A=(0.0/0.0) B=(0.0/0.0) C=(0.0/3.0) D=(10.0/7.0)
  --> Points A and B are identical
A=(0.0/0.0) B=(3.0/3.0) C=(0.0/0.0) D=(6.0/6.0)
  --> Lines AB and CD are identical
A=(0.0/0.0) B=(3.0/3.0) C=(0.0/1.0) D=(6.0/7.0)
  --> Lines AB and CD are parallel

Sidef

Translation of: Perl 6

<lang ruby>func det(a, b, c, d) { a*d - b*c }

func intersection(ax, ay, bx, by,

                 cx, cy, dx, dy) {

   var detAB = det(ax,ay, bx,by)
   var detCD = det(cx,cy, dx,dy)

   var ΔxAB = (ax - bx)
   var ΔyAB = (ay - by)
   var ΔxCD = (cx - dx)
   var ΔyCD = (cy - dy)

   var x_numerator = det(detAB, ΔxAB, detCD, ΔxCD)
   var y_numerator = det(detAB, ΔyAB, detCD, ΔyCD)
   var denominator = det( ΔxAB, ΔyAB,  ΔxCD, ΔyCD)

   denominator == 0 && return 'lines are parallel'
   [x_numerator / denominator, y_numerator / denominator]

}

say ('Intersection point: ', intersection(4,0, 6,10, 0,3, 10,7)) say ('Intersection point: ', intersection(4,0, 6,10, 0,3, 10,7.1)) say ('Intersection point: ', intersection(0,0, 1,1, 1,2, 4,5))</lang>

Output:

Intersection point: [5, 5]
Intersection point: [2300/459, 2320/459]
Intersection point: lines are parallel

zkl

Translation of: C++

<lang zkl>fcn lineIntersect(ax,ay, bx,by, cx,cy, dx,dy){ // --> (x,y)

  detAB,detCD := det(ax,ay, bx,by), det(cx,cy, dx,dy);
  abDx,cdDx := ax - bx, cx - dx;	// delta x
  abDy,cdDy := ay - by, cy - dy;	// delta y

  xnom,ynom := det(detAB,abDx, detCD,cdDx), det(detAB,abDy, detCD,cdDy);
  denom     := det(abDx,abDy, cdDx,cdDy);
  if(denom.closeTo(0.0, 0.0001))
     throw(Exception.MathError("lineIntersect: Parallel lines"));

  return(xnom/denom, ynom/denom);

} fcn det(a,b,c,d){ a*d - b*c } // determinant</lang> <lang zkl>lineIntersect(4.0,0.0, 6.0,10.0, 0.0,3.0, 10.0,7.0).println();</lang>

Output:

L(5,5)