CG · Maths · Python

Convex Hull – Maths can be fun

I always strive to break the mantle of sophistication that is so attached with mathematics. In fact in all of the sciences, if the fundamentals are understood clearly, in simple terms, it is easy for anybody to understand the advanced concepts. But we have to have a solid fundamental base to start with. I always maintain that if we can learn the 26 alphabets of English, we can learn almost anything !!

Recently, I am working on creating a new algorithm for finding the convex hull from a given set of 2d points in a plane. In this post, we will discover the  concept of the convex hull. It will be a fun ride. In fact, the method that I devised for finding the convex hull is so easy that even a six-year-old can do it given a piece of paper and pencil.

So first thing first, what is a Convex Hull ?

Convex Hull

If we keep ourselves restricted to the computational geometry and 2 dimensions:

Convex Hull can be defined as the boundary polygon for a given set of points.

It is clearer from the illustration below.

convexHull

Recently, I am working on creating a new algorithm for finding the convex hull from a given set of 2d points in a plane. In this post, we will discover the  concept of the convex hull. It will be a fun ride, do believe me. In fact, the method that I devised for finding the convex hull is so easy that even a six-year-old can do it given a piece of paper and pencil.

 

So first thing first, what is a Convex Hull ?

 

Convex Hull

If we keep ourselves restricted to the computational geometry and 2 dimensions:

 

Convex Hull can be defined as the boundary polygon for a given set of points.

 

It is clearer from the illustration below.

convexHull

You can also imagine the convex hull by the following analogy. Consider this points as nails on a piece of wood. If you were to wrap a rubber band around these nails, it will snap to the boundary nails. The shape of the rubber band represents the convex hull for the given nail positions.

The convex hull has numerous application in CG. Namely, better bounding boxes for meshes, BSP calculations for render times, finding projection areas to drive asset resolutions and numerous other. Creativity has no bounds, you can use convex hulls to do many things, use your imagination.

Let’s  start,

Suppose we have a set of points in a given 2D plane. We want to find out a way to find all the points that define the boundary of all the given points.

As I said, it is going to be fun. Let us start with a few random points in a 2d plane. We assign a name to each of this points.

convexHullNamed1

Looking at the figure above, we prepare two lists of these points. First, going to from left to right, we take the left most point, and then the next left most point and so on, then we make another list going from top to bottom for the same points.

These are the two list that we will have:

 

Left to Right      : [f, j, d, a, b, h, l, m, c, g, i, e, k]

Top to Bottom  : [h, m, e, j, a, c, b, k, f, i, g, d, l]

The cool thing now is we have all the information we need to process the convex hull. Just grab a paper and pencil and play with these lists according to the rules of the game I define below. In fact, what we will be doing from now is the gist of my algorithm. Remember, we are working with alphabets, so no calculations are needed, not even an addition. It is so simple. Even the alphabets do not represent any order. They are just a bunch of symbols to distinguish between the points,  you can use smileys if you want instead !! and it will still work. We also don’t need to look at the figure even. Just the lists.

 

OK, so let’s start the game, here’s how we play:

 

We start from the Left-Right list. We take the first element of this list. In this ‘f’. This definitely is a border point. So our first point of the convex hull is ‘f’.

 

Convex Hull = [f, ]

 

Next, we look at the point next to ‘f’ in the Left-Right list, it is ‘j’. Now we should see whether ‘j’ also comes after ‘f’ in the Top-Bottom list. No. So we move to the next point in our Left-Right list, it is ‘d’. Does ‘d’ comes after ‘f’ in the Top-Bottom list. Yes, it does. So check-in this point next in our convexHull,

 

Convex Hull = [f, d, ]

As soon as we check-in a point in the Convex Hull list. We should remove the previous ‘f’ from the Left-Right and Top-Bottom lists. So these now become.

 

Left to Right      : [ j, d, a, b, h, l, m, c, g, i, e, k]

Top to Bottom  : [h, m, e, j, a, c, b, k, i, g, d, l]

 

So our focus now moves to ‘d’. What is the next point after ‘d’ in the Left-Right list –  ‘a’. Does ‘a’ comes after ‘d in the Top-Bottom list. No. So we move on, till we reach ‘l’. This comes after ‘d’ in both Left-Right and Top-Bottom list and so this is a convex hull point.

Convex Hull = [f, d, l, ]

 

again we remove previous point ‘d’ from both the lists,

 

Left to Right      : [ j, a, b, h, l, m, c, g, i, e, k]

Top to Bottom  : [h, m, e, j, a, c, b, k, i, g, l]

 

Now as we have reached the last point of the Top-Bottom list. We can no longer check whether there is any point after it. Whenever we reach the last point in any of the list. We should switch the list and change direction. Since we have reached ‘l’, which is the last point of the Top-Bottom list, we should now start looking for points in the Top-Bottom list but in the reverse direction. So what point comes before ‘l’ in the Top-Bottom list. ‘g’, does ‘g’ also comes after ‘l’ in the Left-Right list, yes. So we check in ‘g’ to the Convex Hull.

 

Convex Hull = [f, d, l, g]

 

removing ‘l’ we get,

Left to Right      : [ j, a, b, h, m, c, g, i, e, k]

Top to Bottom  : [h, m, e, j, a, c, b, k, i, g]

 

You get the idea, again, continuing in the same fashion we will have the convex hull as :

Convex Hull = [f, d, l, g, i, k]

 

and the lists as :

Left to Right      : [ j, a, b, h, m, c, e, k]

Top to Bottom  : [h, m, e, j, a, c, b, k]

when we reach ‘k’, this is the last element in the Left-Right list. So according to our rules, we switch the lists and start to look in the reverse in the Left-Right list.  At this moment, we are looking in both the lists in the reverse directions remember. So the point before ‘k’ in Left-Right list is ‘e’, does it also come before ‘k’ in the Top-Bottom list, the answer is yes and so we check in ‘e’ in the convex hull and remove the previous point ‘k’ from the list,

Convex Hull = [f, d, l, g, i, k, e]

Left to Right      : [ j, a, b, h, m, c, e]

Top to Bottom  : [h, m, e, j, a, c, b]

 

moving in the same fashion in both the lists we have,

 

Convex Hull = [f, d, l, g, i, k, e, m, h]

Left to Right      : [ j, a, b, h, c]

Top to Bottom  : [h,  j, a, c, b]

 

‘h’ being the first point of the Top-Bottom list, we again switch the lists and look in the Top-Bottom list in the forward direction, while reverse in the Left-Right list. And we keep on going to find the same patterns of occurrence, we are almost there,

 

in the Top-Bottom list, the element next to ‘h’ is ‘j’. Does ‘j’ comes before ‘h’ in the Left-Right list. Yes, so we check in ‘j’ and remove ‘h’ from the lists. At this point, we stop. Why ? because we have reached where we started, the first element of the Left-Right list. So our convex hull at this point look like.

 

Convex Hull = [f, d, l, g, i, k, e, m, h, j]

 

Let’s draw our hull and check what we have,

convexHullRough

Not bad, without any calculation. But we have a few concave points still inside. These are ‘g’ , ‘i’ and ‘m’. These points might or might not  appear in the hull depending upon the uniqueness of our data set. But we should always as a final step remove these out of our way. To do this we take three points of the hull at a time and check for convexity. Please refer to my post for more details finding convex points.

So we take three point sets from the hull in the order they appear, like ‘fdl’, ‘dlg’, ‘lgi’ and so on and so forth and remove all the concave points, we call this a single iteration. We do a few iterations of the same logic till we have no more convex points left.

This is how these points are removed in the iteration :

First Iteration : ‘i’ is removed from ‘gik’, ‘m’ is removed from ‘emh’.

convexHullFirstIteration

Second Iteration : ‘g’ is removed from ‘lgk’.

convexHullSecondIteration2

And lo! we finally have our convex hull as:

 

Final Convex Hull = [  f, d, l, k, e,  h, j  ]

convexHullFinal

Let’s use our algorithm inside Softimage to find the convex hull for meshes. To use the algorithm, we use a python implementation of the algorithm and feed to it the position array of a mesh. We create a curve based on the convex hull. Just select any mesh and run the below code. It will create a curve representing the convex hull.

 

Note that as we are on a 2D plane, we get the hull as we see the points in front view. Even if you rotate, translate or scale the mesh, the hull will still be created at the current position array. To create the hull for a different view of the mesh, select the points of the mesh and then rotate, scale or transform them. After applying your transformations to the mesh points, freeze the mesh.

 

Here are a few screenshots.

meshHulls

Please note that for the implementation of the algorithm in python, I am not using any optimizing techniques. This code creates a convex hull for 15000 mesh points in about 6 seconds. To deploy this as a production ready solution, I would use numpy, cython or build my own extensions written in c/c++ for performance enhancement and better memory management to handle a few million points. Since that is beyond the scope of this post, so let’s not discuss it.

This  is the  implementation:

from win32com.client import constants as c
import random

def printOnOff(inStr):
    printIt = False
    if printIt:
        print inStr
    else:
        pass

    return

def makeNulls(xCoords, yCoords, sortedPoints,curve):

    nbPoints = len(xCoords)
    mdl = Application.ActiveSceneRoot.AddModel(None,'randomPoints')
    mdl.Properties('Visibility').Parameters('viewvis').Value = False

    points = mdl.AddNull('points')
    points.Properties('Visibility').Parameters('viewvis').Value = False

    if curve:
        points.AddChild(curve)    

    for i in range(nbPoints):
        if i < 10:
            zeros = '000'
        elif i < 100:
            zeros = '00'
        elif i < 1000: zeros = '0' else: zeros = '' null = points.AddNull('point_%s%s'%(zeros, i)) null.Kinematics.Global.Parameters('PosX').Value = xCoords[i] null.Kinematics.Global.Parameters('PosY').Value = yCoords[i] null.Kinematics.Global.Parameters('PosZ').Value = 0 null.Parameters('primary_icon').Value = 0 null.Parameters('shadow_icon').Value = 5 null.Parameters('shadow_colour_custom').Value = True if i in sortedPoints: r = 0.0 else: r = 1.0 null.Parameters('R').Value = r null.Parameters('G').Value = 1.0 null.Parameters('B').Value = 0.0 return def createRandomPoints(nbPoints=20): xCoords = [] yCoords = [] for i in range(nbPoints): xRand = nbPoints * random.random() #* random.randint(1,3) yRand = nbPoints * random.random() #* random.randint(1,3) x = random.uniform(-1 * xRand, xRand) y = random.uniform(-1 * yRand, yRand) xCoords.append(x) yCoords.append(y) return xCoords, yCoords def makeOrderedTuples(xCoords, yCoords): data = None nbPoints = len(xCoords) xCoordsTuples = [] yCoordsTuples = [] for i in range(nbPoints): xRand = nbPoints * random.random() #* random.randint(1,3) yRand = nbPoints * random.random() #* random.randint(1,3) x = xCoords[i] y = yCoords[i] xCoordsTuples.append((x, i)) yCoordsTuples.append((y, i)) sortedXIndex = [] sortedXIndex = [ xTuple[1] for xTuple in sorted(xCoordsTuples, key=lambda xCoord: xCoord[0])] sortedYIndex = [] sortedYIndex = list(reversed([ yTuple[1] for yTuple in sorted(yCoordsTuples, key=lambda yCoord: yCoord[0])])) data = xCoordsTuples, yCoordsTuples, sortedXIndex, sortedYIndex return data def makeCoordString(coords): xCoords, yCoords,zCoords = coords pStr = '' nbPoints = len(xCoords) for i in range(nbPoints): pStr += '(' pStr += str(xCoords[i]) pStr += ',' pStr += str(yCoords[i]) pStr += ',' pStr += str(zCoords[i]) pStr += ')' if i != nbPoints - 1: pStr += ',' return pStr def makeCurve(pStr): curve = Application.CreateCurve( 1, 0, pStr)('Value') curve.Name = 'convexHull' Application.ApplyTopoOp('CrvOpenClose', curve, 3, c.siPersistentOperation) Application.FreezeObj(curve) #curve.Properties('Visibility').Parameters('selectability').Value = False return curve def isConvex(pointA, pointP, pointB): x1, y1 = pointA x2, y2 = pointB x3, y3 = pointP quadB = 0 if x2 > x1 and y2 > y1:
        quadB = 0
    elif x2 < x1 and y2 > y1:
        quadB = 1
    elif x2 < x1 and y2 < y1: quadB = 2 elif x2 > x1 and y2 < y1: quadB = 3 y = y1 + ((x3 - x1) * (y2 - y1) / (x2 - x1)) if quadB == 0 or quadB == 3: if y > y3:
            return 1
        else:
            return 0
    elif quadB == 1 or quadB == 2:
        if y > y3:
            return 0
        else:
            return 1

def getConvexHull(xCoords, yCoords):
    xCoords, yCoords, sortedX, sortedY = makeOrderedTuples(xCoords, yCoords)
    startPoint = sortedX[0]
    firstPoint = startPoint
    sortedPoints = [firstPoint]
    inflection = 0
    normalInflection = False
    inflectionPoints = []
    nextPoint = 0
    nbPoints = len(sortedX)
    for i in range(nbPoints):
        yIndexOfFirstPoint = sortedY.index(firstPoint)
        xIndexOfFirstPoint = sortedX.index(firstPoint)
        for xIndex in range(xIndexOfFirstPoint + 1, nbPoints):
            pointToCheck = sortedX[xIndex]
            yIndexOfPointToCheck = sortedY.index(pointToCheck)
            if yIndexOfPointToCheck > yIndexOfFirstPoint:
                nextPoint = pointToCheck
                sortedPoints.append(nextPoint)
                if yIndexOfPointToCheck + 1 == nbPoints:
                    printOnOff('First Inflection at %s !!'%nextPoint)
                    inflectionPoints.append(nextPoint)
                    inflection += 1
                    normalInflection = True
                break
        if inflection == 1:
            break

        firstPoint = nextPoint

    if not normalInflection:
        printOnOff('First Inflection at %s - Abnormal!!'%firstPoint)
        inflection += 1

    normalInflection = False
    firstPoint = nextPoint
    for i in range(nbPoints):
        xIndexOfFirstPoint = sortedX.index(firstPoint)
        yIndexOfFirstPoint = sortedY.index(firstPoint)
        for yIndex in range(yIndexOfFirstPoint -1, -1,-1):
            pointToCheck = sortedY[yIndex]
            xIndexOfPointToCheck = sortedX.index(pointToCheck)
            if xIndexOfPointToCheck > xIndexOfFirstPoint:
                nextPoint = pointToCheck
                sortedPoints.append(nextPoint)
                if xIndexOfPointToCheck + 1 == nbPoints:
                    printOnOff('Second Inflection at %s !!'%nextPoint)
                    inflectionPoints.append(nextPoint)
                    inflection += 1
                    normalInflection = True
                break
        if inflection == 2:
            break

        firstPoint = nextPoint

    if not normalInflection:
        printOnOff('Second Inflection at %s - Abnormal!!'%firstPoint)
        inflection += 1

    normalInflection = False
    firstPoint = nextPoint
    for i in range(nbPoints):
        yIndexOfFirstPoint = sortedY.index(firstPoint)
        xIndexOfFirstPoint = sortedX.index(firstPoint)
        for xIndex in range(xIndexOfFirstPoint - 1, -1,-1):
            pointToCheck = sortedX[xIndex]
            yIndexOfPointToCheck = sortedY.index(pointToCheck)
            if yIndexOfPointToCheck < yIndexOfFirstPoint:
                nextPoint = pointToCheck
                sortedPoints.append(nextPoint)
                if yIndexOfPointToCheck == 0:
                    printOnOff('Third Inflection at %s !!'%nextPoint)
                    inflectionPoints.append(nextPoint)
                    inflection += 1
                    normalInflection = True
                break
        if inflection == 3:
            break

        firstPoint = nextPoint

    if not normalInflection:
        printOnOff('Third Inflection at %s - Abnormal!!'%nextPoint)
        inflection += 1

    printOnOff('\n\n\n')
    normalInflection = False
    firstPoint = nextPoint
    for i in range(nbPoints):
        xIndexOfFirstPoint = sortedX.index(firstPoint)
        yIndexOfFirstPoint = sortedY.index(firstPoint)
        for yIndex in range(yIndexOfFirstPoint + 1, nbPoints):
            pointToCheck = sortedY[yIndex]
            xIndexOfPointToCheck = sortedX.index(pointToCheck)

            if xIndexOfPointToCheck < xIndexOfFirstPoint:
                if pointToCheck in sortedPoints:
                    printOnOff('Final Point %s added !!'%firstPoint)
                    inflection += 1
                    normalInflection = True
                    break
                #printOnOff('breaking out . . . ')
                nextPoint = pointToCheck
                sortedPoints.append(nextPoint)
                break
        if inflection == 4:
            break

        firstPoint = nextPoint

    if not normalInflection:
        printOnOff('Fourth Inflection at %s - Abnormal!!'%nextPoint)
        inflection += 1

    printOnOff('\n\n\n')
    printOnOff('No. of Inflections : %s'%inflection)

    sortedXCoords = [ xTuple[0] for xTuple in sorted(xCoords,
                                                  key=lambda xCoord: xCoord[1])]

    sortedYCoords = [ yTuple[0] for yTuple in sorted(yCoords,
                                                  key=lambda yCoord: yCoord[1])]

    printOnOff('Inflection Points : %s'%inflectionPoints)

    filterConcave = True
    if filterConcave:
        noConcaveLeft = False
        filteredPoints = sortedPoints
        while not noConcaveLeft:
            nbFilteredPoints = len(filteredPoints)
            concavePoints = []
            for i in range(nbFilteredPoints):
                if i == nbFilteredPoints - 1:
                    break

                j = i + 1
                if i == nbFilteredPoints - 2:
                    k = 0
                else:
                    k = i + 2

                indexStartPoint = filteredPoints[i]
                indexMiddlePoint = filteredPoints[j]
                indexEndPoint = filteredPoints[k]

                startPoint = (sortedXCoords[indexStartPoint], sortedYCoords[indexStartPoint])
                middlePoint = (sortedXCoords[indexMiddlePoint],
                                       sortedYCoords[indexMiddlePoint])
                endPoint = (sortedXCoords[indexEndPoint], sortedYCoords[indexEndPoint])

                if not isConvex(startPoint, middlePoint, endPoint):
                    concavePoints.append(indexMiddlePoint)

            if not len(concavePoints):
                noConcaveLeft = True

            for point in filteredPoints:
                if point in concavePoints:
                    filteredPoints.remove(point) 

        sortedPoints = filteredPoints

    makeConvexHullCurve = True
    curve = None
    if makeConvexHullCurve:
        nbBorderPoints = len(sortedPoints)

        xSortedPoints = [sortedXCoords[i] for i in sortedPoints]
        ySortedPoints = [sortedYCoords[i] for i in sortedPoints]
        zSortedPoints = [0.0 for i in range(nbBorderPoints)]

        coordsForString = (xSortedPoints, ySortedPoints, zSortedPoints)

        pStr = makeCoordString(coordsForString)
        curve = makeCurve(pStr)

    return sortedPoints, curve

def createHullFromSelectedMesh():
    mesh = Application.Selection(0)
    if not mesh or mesh.Type != 'polymsh':
        Application.LogMessage('Invalid Selection !!!\nPlease select a polymesh.', 2)
        return
    geom = mesh.ActivePrimitive.Geometry
    points = geom.Points
    posArr = points.PositionArray

    xCoords = list(posArr[0])
    yCoords = list(posArr[1])
    sortedPoints, curve = getConvexHull(xCoords, yCoords)
    curve.Name = 'convexHull%s'%mesh.FullName
    mesh.AddChild(curve)

    printOnOff('Convex Hull Created Successfully')

    return

def createHullFromRandomPoints(nbPoints = 100):
    xCoords, yCoords = createRandomPoints(nbPoints)
    sortedPoints, curve = getConvexHull(xCoords, yCoords)
    makeNulls(xCoords, yCoords, sortedPoints, curve)
    printOnOff('Convex Hull Created Successfully')

    return

#--------------------------------------------------------------------------
# If you wish test the algorithm on a random points then set the flag
#
# hullFromMesh = False
#
# else the hull will be created from the selected mesh
#
# This will create a bunch of random point and claculate the hull on these
# points. You can set the total number of points with the variable
# nbPoints.
#
# If you care for data logged out of algorithm then set the flag
# printIt = True in the top most function printOnOff. It is off by default.
#-------------------------------------------------------------------------

hullFromMesh = True

if hullFromMesh:
    createHullFromSelectedMesh()
else:
    createHullFromRandomPoints(nbPoints=100)

Enjoy!

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s