Orientation and intersection

Do two line pieces meet? The question sounds like it needs line equations, slopes and division. It does not: four signs of cross products answer it, exactly, with whole numbers. The hard part is the special cases, where a point lies exactly on a line.

Orientation

The orientation of three points A, B, C is the sign of cross(A, B, C) from lesson 1:

• positive: A → B → C turns left (counterclockwise), so C is on the left of the line through A and B, looking from A towards B;
• negative: the turn goes right (clockwise);
• zero: the three points lie on one line (they are collinear).

A line goes on forever in both directions; a line piece (a segment) AB stops at A and B. Orientation is about the line; whether a point is on the piece needs one more check.

Do two pieces meet?

Take pieces AB and CD. If C and D are on opposite sides of the line AB, then the piece CD crosses that line somewhere. That alone is not enough: the crossing point may lie beyond A or B. But if also A and B are on opposite sides of the line CD, the crossing point lies on both pieces.

In signs: cross(A, B, C) and cross(A, B, D) have opposite signs, and so do cross(C, D, A) and cross(C, D, B). Step through the first two variants of the figure: in “Miss” only the first test passes.

When an orientation is 0

If one orientation is 0, say cross(A, B, C) = 0, then C is on the line AB. It is on the piece AB exactly when it lies in the box with corners A and B: its x between A.x and B.x, and its y between A.y and B.y. Checking only x is wrong for an upright piece, where every x is the same.

def cross(o, a, b):
    return (a[0] - o[0]) * (b[1] - o[1]) - (a[1] - o[1]) * (b[0] - o[0])

def within(p, q, r):            # r is in the box with corners p and q
    return (min(p[0], q[0]) <= r[0] <= max(p[0], q[0]) and
            min(p[1], q[1]) <= r[1] <= max(p[1], q[1]))

def pieces_meet(a, b, c, d):
    d1, d2 = cross(a, b, c), cross(a, b, d)
    d3, d4 = cross(c, d, a), cross(c, d, b)
    if ((d1 > 0 > d2) or (d1 < 0 < d2)) and ((d3 > 0 > d4) or (d3 < 0 < d4)):
        return True                         # a proper crossing
    return ((d1 == 0 and within(a, b, c)) or (d2 == 0 and within(a, b, d)) or
            (d3 == 0 and within(c, d, a)) or (d4 == 0 and within(c, d, b)))

Why the second part is enough. If the pieces meet but do not cross properly, the meeting point is an end of one of them: it is A, B, C or D, or the pieces lie on one line and overlap, and then an end of one lies on the other. In each case one orientation is 0 and that end is in the other piece’s box. The other way round, a 0 orientation with the end in the box means that end lies on the other piece.

Cases to test, with the variants “Touch” and “One line” of the figure:

• an end of one piece on the middle of the other (a T shape);
• two pieces that share only an end;
• two pieces on one line: they may be apart, touch at one end, or overlap along a stretch. All four orientations are 0, so only the boxes decide;
• a “piece” whose two ends are equal: a single point;
• upright pieces, where all x are equal.

To tell one common point from many (the practice problem): if the four points are not all on one line, the pieces share at most one point. If they are all on one line, measure each piece along the line (by x, or by y for an upright line) and compare the overlap of the two ranges: empty, a single value, or longer.

Where they cross

Write the pieces as A + t·r and C + u·s, with r = B − A, s = D − C, and t, u between 0 and 1. Taking the cross product of both sides of A + t·r = C + u·s with s removes u (because cross(s, s) = 0):

t = cross(C − A, s) / cross(r, s).

If cross(r, s) = 0, the pieces are parallel and there is no single crossing point. Otherwise the point is A + t·r. With fractions.Fraction it is exact:

from fractions import Fraction

def crossing_point(a, b, c, d):
    r = (b[0] - a[0], b[1] - a[1])
    s = (d[0] - c[0], d[1] - c[1])
    den = r[0] * s[1] - r[1] * s[0]
    if den == 0:
        return None                          # parallel lines
    t = Fraction((c[0] - a[0]) * s[1] - (c[1] - a[1]) * s[0], den)
    return (a[0] + t * r[0], a[1] + t * r[1])

crossing_point((0, 0), (4, 4), (0, 4), (4, 0))   # (Fraction(2, 1), Fraction(2, 1))

Only the yes-or-no question needs no division at all, which is why the orientation test is preferred.

Is a point inside a polygon?

Send a ray (a half-line) from the point P to the right, and count how many edges of the polygon it crosses. Each crossing moves between inside and outside, and far to the right is outside, so P is inside exactly when the count is odd. This is the ray casting or even-odd rule.

• First check the edges: if P lies on an edge (orientation 0 and inside the box), answer “on the edge”.
• Corners on the ray: count an edge only when exactly one of its ends is strictly above P ((y1 > py) != (y2 > py)). A corner the ray passes through is then counted once or not at all, correctly, and flat edges are never counted.
• No division: for an edge going up from (x₁, y₁) to (x₂, y₂), the ray meets it to the right of P exactly when P is on its left: cross(P₁, P₂, P) > 0. For an edge going down, the sign flips.

def where(poly, p):
    n = len(poly)
    inside = False
    for k in range(n):
        a, b = poly[k], poly[(k + 1) % n]
        if cross(a, b, p) == 0 and within(a, b, p):
            return "edge"
        if (a[1] > p[1]) != (b[1] > p[1]):       # the edge spans p's height
            c = cross(a, b, p)
            if (c > 0) == (b[1] > a[1]):        # it crosses to the right of p
                inside = not inside
    return "inside" if inside else "outside"

Each test is O(n) for n corners. For a convex polygon, binary search over the corners answers in O(log n).

Now try GEOM-P2, Do two tracks meet?.