Search code examples
pythonalgorithmgeometry

Gift Wrapping Algorithm in Python never terminating


I've been trying to implement the gift wrapping algorithm in Python, I currently have the following code:

def createIslandPolygon(particleCoords):

    startPoint = min(particleCoords.iteritems(),key = lambda x: x[1][1])[1]

    check = 1

    islandPolygon = []

    particleList = []

    for key in particleCoords:

        particleList.append(particleCoords[key])

    currentPoint = startPoint

    while(currentPoint != startPoint or check == 1):

        islandPolygon.append(currentPoint)

        check = 0

        angleDict = {}
        angleList = []

        for point in particleList:

            if point != currentPoint:

                angleDict[(angleBetweenTwoPoints(currentPoint, point))] = point
                angleList.append(angleBetweenTwoPoints(currentPoint, point))

        smallestAngle = min(angleList)

        currentPoint = angleDict[smallestAngle]

    return islandPolygon

and for calculating the polar coordinates:

def angleBetweenTwoPoints(p1, p2):

    p3 = (p1[0], p1[1] + 2)

    a = (p1[0] - p2[0], p1[1] - p2[1])
    b = (p1[0] - p3[0], p1[1] - p3[1])

    theta = ((a[0]*b[0]) + (a[1]*b[1]))
    theta = theta / (sqrt((a[0]*a[0]) + (a[1]*a[1])) * sqrt((b[0]*b[0]) + (b[1]*b[1])))
    theta = math.acos(theta)

    return theta

The problem is that the code never seems to leave the while loop, and I'm not sure why. Does anyone have any idea?

Thanks.

(yeah, the codes pretty shoddy, I just threw it together quickly)

Edit: Printing out the coordinates seems to show them jumping between just two coordinates.


Solution

  • According to http://en.wikipedia.org/wiki/Gift_wrapping_algorithm you need to do this:

       pointOnHull = leftmost point in S
       i = 0
       repeat
          P[i] = pointOnHull
          endpoint = S[0]         // initial endpoint for a candidate edge on the hull
          for j from 1 to |S|-1
             if (endpoint == pointOnHull) or (S[j] is on left of line from P[i] to endpoint)
                endpoint = S[j]   // found greater left turn, update endpoint
          i = i+1
          pointOnHull = endpoint
       until endpoint == P[0]      // wrapped around to first hull point
    

    You have this correct:

       pointOnHull = leftmost point in S
    

    And this:

      P[i] = pointOnHull
    

    But here's the part I'm not sure about:

      (S[j] is on left of line from P[i] to endpoint)
    

    Instead you find the smallest angle out of all the angles it makes with all other points. But according to wikipedia what you want is the leftmost angle out of all the angles it makes with all other points. I have some code for working with angles:

    def normalizeangle(radians):
        return divmod(radians, math.pi*2)[1]
    
    
    
    def arclength(radians1, radians2 = 0):
        radians1, radians2 = normalizeangle(radians1), normalizeangle(radians2)
        return min(normalizeangle(radians1 - radians2), normalizeangle(radians2 - radians1))
    
    
    
    def arcdir(radians1, radians2 = 0):
        radians1, radians2 = normalizeangle(radians1), normalizeangle(radians2)
        return cmp(normalizeangle(radians1 - radians2), normalizeangle(radians2 - radians1))
    

    arcdir will tell you if an angle is on the left or on the right of another, so you can use it to see if an angle is more left and so should be used.

    If you move along the points always picking the leftmost angle to the next point, you'll go around the perimeter of your polygon and reach the beginning again (since you picked the leftmost point, you know it has to be on the perimeter and will be reached again)