An efficient way of merging two convex hulls

Before reading this article, I recommend you to visit following two articles. Many concepts and codes are referred from there.

1. Determining if two consecutive segments turn left or right
2. Convex Hull Algorithms: Jarvis’s March (Introduction Part)

Introduction

In this article, I am going to talk about the linear time algorithm for merging two convex hulls. Given two convex hull as shown in the figure below.

The algorithm should produce the final merged convex hull as shown in the figure below.

Note: In order for this algorithm to work correctly, two convex hulls must be in the distinct left and right position and should not be overlapped. That means the x-coordinates of all the points of the left convex hull must be less than the x-coordinates of all the points of the right convex hull.

Algorithm

1. Find the rightmost point $(p)$ of the left convex hull and leftmost $(q)$ for the right convex hull.
2. Make two copies of $p$ and $q$. Now we have two $p$s and two $q$s.
3. Raise the first copy of $p$ and $q$ to the make the upper tangent.
4. Lower the second copy of $p$ and $q$ to make the lower tangent.

The details of steps 3, 4 are given in the next section with the help of an example.

Example

Suppose we want to merge the following two convex hulls.

First, we find the rightmost point ($p$) of the left convex hull ($q$) and leftmost point of the right convex hull and make two copies of the points $p$ and $q$.

Next, we raise the one copies of $pq$ (shown by a green line segment) to the upper tangent. This is how we do it. We keep the point $p$ fix and rotate the point $q$ in a clockwise direction as long as the sequence $p$, $q$ and the clockwise neighbor of $q$ makes the left turn (anti-clockwise turn). In this case, it makes a left turn so we move the point $q$ as shown in the figure below.

Again keep the $p$ fix in its previous position, we try to raise $q$. The sequence $p$, $q$ and the clockwise neighbor of $q$ constitutes a right turn and we stop. That means we can not raise $q$ further. Now, we keep the point $q$ fix and try to raise the point $p$ using the opposite logic. That is, we continue to raise $p$ as long as $q$, $p$ and anti-clockwise neighbor of $p$ makes a right turn. Since they make a right turn, we move $p$ as shown in the figure below.

Again, sequence $q$, $p$ and anti-clockwise neighbor of $p$ makes a right turn. Therefore, we move $p$.

Now, we try to raise $p$ but we can not do that because the sequence makes left turn. We stop here and repeat the process keeping $p$ fix. The sequence $p$, $q$ and clockwise neighbor of $q$ makes right turn. That means we can not raise $q$ further. Segment $pq$ is now the upper tangent.

We use the similar approach to lower another copies of $pq$ (shown by the red line). First keeping $p$ fix and lowering $q$ as long as the sequence $p$, $q$ and anti-clockwise neighbor of $q$ makes the right turn. Second we fix $q$ and lower $p$ as long as the sequence $q$, $p$ and clockwise neighbor of $p$ constitutes left turn. We stop when we cannot lower both $p$ and $q$ further.

In the above example, we lower the $q$ sine the sequence makes a right turn.

We can not lower $q$ further. We try to lower $p$ keeping $q$ fixed. Since the sequence makes a left turn, we lower it.

The sequence makes right turn. We stop and try to lower $q$ again keeping $p$ fixed. As the sequence makes a right turn, we lower $q$.

$pq$ is now the lower tangent as we cannot lower both of them any further. The algorithm terminates here.

Implementation

To implement the algorithm in linear time, a suitable data structure is needed. Here I am using a circular linked list where each point holds the references to the next nearest points in a clockwise and anti-clockwise direction. A node of the linked list can be implemented as below.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

class Point:
    def __init__(self, x, y):
        self.x = x
        self.y = y
        self.cw_next = None
        self.ccw_next = None

    def subtract(self, p):
        return Point(self.x - p.x, self.y - p.y)

    def __str__(self):
        return '(' + str(self.x) + ', ' + str(self.y) + ')'

    def __eq__(self, p):
        return self.x == p.x and self.y == p.y

Variable cw_next holds the adjacent point in the clockwise direction and variable ccw_next holds the adjacent point in the counter clockwise direction. The merge method is straight forward as described in the algorithm section.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64

def merge(chull1, chull2):
    # get the rightmost point of left convex hull
    p = max(chull1, key = lambda point: point.x)

    # get the leftmost poitn of right convex hull
    q = min(chull2, key = lambda point: point.x)

    # make copies of p and q
    cp_p = p
    cp_q = q

    # raise the bridge pq to the uper tangent
    prev_p = None
    prev_q = None
    while (True):
        prev_p = p
        prev_q = q
        if q.cw_next:
            # move p clockwise as long as it makes left turn
            while direction(p, q, q.cw_next) < 0:
                q = q.cw_next
        if p.ccw_next:
            # move p as long as it makes right turn
            while direction(q, p, p.ccw_next) > 0:
                p = p.ccw_next

        if p == prev_p and q == prev_q:
            break
    
    # lower the bridge cp_p cp_q to the lower tangent
    prev_p = None
    prev_q = None
    while (True):
        prev_p = cp_p
        prev_q = cp_q
        if cp_q.ccw_next:
            # move q as long as it makes right turn
            while direction(cp_p, cp_q, cp_q.ccw_next) > 0:
                cp_q = cp_q.ccw_next
        if cp_p.cw_next:
            # move p as long as it makes left turn
            while direction(cp_q, cp_p, cp_p.cw_next) < 0:
                cp_p = cp_p.cw_next
        if cp_p == prev_p and cp_q == prev_q:
            break

    # remove all other points
    p.cw_next = q
    q.ccw_next = p

    cp_p.ccw_next = cp_q
    cp_q.cw_next = cp_p

    # final result
    result = []
    start = p 
    while (True):
        result.append(p)
        p = p.ccw_next

        if p == start:
            break

    return result

The following code can be tested with the test case given below.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42

if __name__ == '__main__':
    p1 = Point(2, 2)
    p2 = Point(3, 4)
    p3 = Point(7, 3)
    p4 = Point(5, 2)
    p5 = Point(7, 5)
    p6 = Point(5, 6)

    p1.cw_next = p2 
    p1.ccw_next = p4
    p2.cw_next = p6
    p2.ccw_next = p1
    p3.cw_next = p4
    p3.ccw_next = p5
    p4.cw_next = p1 
    p4.ccw_next = p3 
    p5.cw_next = p3
    p5.ccw_next = p6
    p6.cw_next = p5
    p6.ccw_next = p2

    p7 = Point(11, 5)
    p8 = Point(15, 7)
    p9 = Point(12, 3)
    p10 = Point(15, 3)
    p11 = Point(16, 5)
    p12 = Point(12, 7)

    p7.cw_next = p12
    p7.ccw_next = p9
    p8.cw_next = p11
    p8.ccw_next = p12
    p9.cw_next = p7
    p9.ccw_next = p10
    p10.cw_next = p9
    p10.ccw_next = p11
    p11.cw_next = p10
    p11.ccw_next = p8
    p12.cw_next = p8
    p12.ccw_next = p7

    p = merge([Point(1, 1)], [Point(0, 0)]

Complexity

If you observe the code, you notice there is a nested while loop. You might be tempted to think this has a quadratic time complexity. But if you carefully observe, each point is covered at most twice. So the overall complexity is $O(n)$.

References

1. Mount, D. M. (n.d.). CMSC 754 Computational Geometry. Lecture. Retrieved August 23, 2018, from https://www.cs.umd.edu/class/spring2012/cmsc754/Lects/cmsc754-lects.pdf