Find out if 2 lines intersect [duplicate]

Question

Possible Duplicate:
How do you detect where two line segments intersect?
Determining if two line segments intersect?

Given are two lines l1=((A0, B0), (A1, B1)) and l2=((A2, B2), (A3, B3)); Ax, Bx are integers and (Ax, Bx) specify the starts and ends of the lines.

Is there a algorithm using only integer arithmetic that determines if l1 and l2 intersect? (Only a Boolean answer is needed.)

My own approach was to compute a point near the intersection point with fixed-point arithmetic. The solution (a, b) was then substituted in the following equations:

I: abs((A0 + a * (A1-A0)) - (A2 + b * (A3-A2))) < tolerance
II: abs((B0 + a * (B1-B0)) - (B2 + b * (B3-B2))) < tolerance

My method should return true if both I and II evaluate to true.

My C++-Code:
vec.h:

#ifndef __MY_VECTOR__
#define __MY_VECTOR__
#include <stdarg.h>
template<typename VType, unsigned int dim>
class vec {
private:
    VType data[dim];
public:
    vec(){}
    vec(VType v0, ...){
            data[0] = v0;
            va_list l;
            va_start(l, v0);
            for(unsigned int i=1; i<dim; ++i){
                    data[i] = va_arg(l, VType);
            }
            va_end(l);
    }
    ~vec(){}
    VType& operator[](unsigned int i){
            return data[i];
    }
    VType operator[](unsigned int i) const {
            return data[i];
    }};
    template<typename VType, unsigned int dim, bool doDiv>
    vec<VType, dim> helpArith1(const vec<VType, dim>& A, long delta){
            vec<VType, dim> r(A);
            for(unsigned int i=0; i<dim; ++i){
                    r[i] = doDiv ? (r[i] / delta) : (r[i] * delta);
            }
            return r;
    }
    template<typename VType, unsigned int dim>
    vec<VType, dim> operator*(const vec<VType, dim>& v, long delta) {
        return helpArith1<VType, dim, false>(A, delta);
    }
    template<typename VType, unsigned int dim>
    vec<VType, dim> operator*(long delta, const vec<VType, dim>& v){
        return v * delta;
    }
    template<typename VType,unsigned int dim>
    vec<VType, dim> operator/(const vec<VType, dim>& A, long delta) {
        return helpArith1<VType, dim, true>(A, delta);
    }
    template<typename VType, unsigned int dim, bool doSub>
    vec<VType, dim> helpArith2(const vec<VType, dim>& A, const vec<VType, dim>& B){
        vec<VType, dim> r;
        for(unsigned int i=0; i<dim; ++i){
            r[i] = doSub ? (A[i]-B[i]):(A[i]+B[i]);
        }
        return r;
    }
    template<typename VType, unsigned int dim>
    vec<VType, dim> operator+(const vec<VType, dim>& A, const vec<VType, dim>& B){
        return helpArith2<VType, dim, false>(A, B);
    }
    template<typename VType, unsigned int dim>
    vec<VType, dim> operator-(const vec<VType, dim>& A, const vec<VType, dim>& B){
        return helpArith2<VType, dim, true>(A, B);
    }
    template<typename VType, unsigned int dim>
    bool operator==(const vec<VType, dim>& A, const vec<VType, dim>& B) {
            for(unsigned int i==0; i<dim; ++i){
                if(A[i]!=B[i]){
                            return false;
                    }
            }
            return true;
    }
    template<typename VType, unsigned int dim>
    bool operator!=(const vec<VType, dim>& A, const vec<VType, dim>& B) {
            return !(A==B);
    }
    #endif

line.h:

#ifndef __MY_LINE__
#define __MY_LINE__
#include "vec.h"
unsigned long int ggt(unsigned long int A, unsigned long int B) {
    if(A==0) {
        if(B==0) {
            return 1;
        }
        return B;
    }
    while(B!=0) {
        unsigned long int temp = A % B;
        A = B;
        B = temp;
    }
    return A;
}
#define ABS(n) ( ((n)<0) ? (-n) : (n) )
struct line {
    vec<long int, 2> A, B;
    explicit line(long int iA_0, long int iA_1, long int iB_0, long int iB_1) :
        A(vec<long int, 2>(iA_0<<8, iA_1<<8)),
        B(vec<long int, 2>(iB_0<<8, iB_1<<8)){}
    vec<long int, 2> slope() const{
        vec<long int, 2> temp = A-B;
        if(temp[0]<0) {
            temp[0] = -1 * temp[0];
            temp[1] = -1 * temp[1];
        }
        return temp/ggt(ABS(temp[0]), ABS(temp[1]));
    }
};
bool intersect(line l1, line l2) {
    const long int epsilon = 1<<4;
    vec<long int, 2> sl1 = l1.slope(), sl2 = l2.slope();
    // l2.A + b*sl2 = l1.A + a*sl1
    // <=> l2.A - l1.A = a*sl1 - b*sl2  // = (I, II)^T
    // I': sl2[1] * I; II': sl2[0] * II
    vec<long int, 2> L = l2.A - l1.A, R = sl1;
    L[0] = L[0] * sl2[1];        R[0] = R[0] * sl2[1];
    L[1] = L[1] * sl2[0];        R[1] = R[1] * sl2[0];
    // I' - II'
    long int L_SUB = L[0] - L[1], R_SUB = R[0] - R[1];
    if(ABS(R_SUB) == 0) {
        return ABS(L_SUB) == 0;
    }
    long int temp = ggt(ABS(L_SUB), ABS(R_SUB));
    L_SUB /= temp; R_SUB /= temp;
    // R_SUB * a = L_SUB
    long int a = L_SUB/R_SUB, b = ((l1.A[0] - l2.A[0])*R_SUB + L_SUB * sl1[0])/R_SUB;
    // if the given lines intersect, then {a, b} must be the solution of
    // l2.A - l1.A = a*sl1 - b*sl2
    L = l2.A - l1.A;
    long x = ABS((L[0]- (a*sl1[0]-b*sl2[0]))), y = ABS((L[1]- (a*sl1[1]-b*sl2[1])));
    return x<epsilon && y < epsilon;
}
#endif

main.cpp:

#include "line.h"
int main(){
    line A(0, 0, 6, 0), B(3, 3, 4, -3);
    bool temp = intersect(A, B);
    return 0;
}

(I am not sure if my intersect function works for all lines, but with the test data I used so far it returned with the correct result.)

Answer 1

This is possible. We want to check whether both endpoins of l1 are on different sides of l2, and both endpoints of l2 are on opposite sides of l1.

To check on which side of l1=((A0, B0), (A1, B1)) a point (A, B) lies, we take:

• an arbitrary normal vector N perpendicular to the line; one such vector is (B1-B0, A1-A0)
• the vector P from the start of the line to the point (A, B), which is (A-A0, B-B0)

We then compute the dot product:

N · P = (A-A0, B-B0) · (B1-B0, A1-A0) = (A-A0) * (B1-B0) + (B-B0) * (A1-A0)

We're only interested in the sign: if it's positive, the point is on one side of the line; if it's negative, it's on the other. As you see, no floating point arithmetic required.

We can take advantage of the fact that numbers with opposite signs, when multiplied, always yield negative. So the full expression to determine whether two line segments ((A0, B0), (A1, B1)) and ((A2, B2), (A3, B3)) intersect is:

((A2-A0)*(B1-B0) - (B2-B0)*(A1-A0)) * ((A3-A0)*(B1-B0) - (B3-B0)*(A1-A0)) < 0
&&
((A0-A2)*(B3-B2) - (B0-B2)*(A3-A2)) * ((A1-A2)*(B3-B2) - (B1-B2)*(A3-A2)) < 0

Test Code

Some C++ code to test the above calculation:

#include <iostream>
#include <cstdlib>

struct Point {
    int x,y;
};

bool isIntersecting(Point& p1, Point& p2, Point& q1, Point& q2) {
    return (((q1.x-p1.x)*(p2.y-p1.y) - (q1.y-p1.y)*(p2.x-p1.x))
            * ((q2.x-p1.x)*(p2.y-p1.y) - (q2.y-p1.y)*(p2.x-p1.x)) < 0)
            &&
           (((p1.x-q1.x)*(q2.y-q1.y) - (p1.y-q1.y)*(q2.x-q1.x))
            * ((p2.x-q1.x)*(q2.y-q1.y) - (p2.y-q1.y)*(q2.x-q1.x)) < 0);
}

int main(int argc, char* argv[]) {
    if(argc != 9) {
        std::cout << "Call as " << argv[0] << " <p1.x> <p1.y> <p2.x> "
                  << "<p2.y> <q1.x> <q1.y> <q2.x> <q2.y>" << std::endl;
        return -1;
    }

    Point p1 = {.x = atoi(argv[1]), .y = atoi(argv[2])};
    Point p2 = {.x = atoi(argv[3]), .y = atoi(argv[4])};
    Point q1 = {.x = atoi(argv[5]), .y = atoi(argv[6])};
    Point q2 = {.x = atoi(argv[7]), .y = atoi(argv[8])};

    if(isIntersecting(p1,p2,q1,q2)) {
        std::cout << "Segments intersect" << std::endl;
        return 1;
    }
    else {
        std::cout << "Segments do not intersect" << std::endl;
        return 0;
    }
}

Results:

$ ./intersection_test 0 0 10 10 0 10 10 0 # example from the comments
Segments intersect
$ ./intersection_test 0 1 2 1 1 2 1 0
Segments intersect
$ ./intersection_test 0 0 0 1 1 1 1 0
Segments do not intersect
$ ./intersection_test 1 1 5 3 3 4 7 2 # q touches but not intersects at p2
Segments do not intersect                             
$ ./intersection_test 1 1 5 3 3 4 6 2
Segments intersect