Formula (bool): True indicates that the area of the circumscribing circle
    should be calculated as pi*r^2. This is a more perfect circle than any
    circle composed of squares can be. As a result, the circumscribing circle
    loses the property of representing the maximum number of cells that could
    possibly be in the patch - this can produce negative numbers, particularly
    for small patches. Setting this value to False will calculate a
    circumscribing circle by directly counting grid cells that would be
    included in the circle. This retains the property of never being less than
    patch area and so will never return a value less than 0. However, it
    produces some strange artifacts for small patches and less precisely
    approximates the values reported in the original paper introducing this
    metric (Baker and Cai, 1992). It will also be slightly slower.

def related_circumscribing_circle(patch, formula=True, world_size=(60, 60)):
    """
    Formula (bool): True indicates that the area of the circumscribing circle
    should be calculated as pi*r^2. This is a more perfect circle than any
    circle composed of squares can be. As a result, the circumscribing circle
    loses the property of representing the maximum number of cells that could
    possibly be in the patch - this can produce negative numbers, particularly
    for small patches. Setting this value to False will calculate a
    circumscribing circle by directly counting grid cells that would be
    included in the circle. This retains the property of never being less than
    patch area and so will never return a value less than 0. However, it
    produces some strange artifacts for small patches and less precisely
    approximates the values reported in the original paper introducing this
    metric (Baker and Cai, 1992). It will also be slightly slower.
    """
    patch_area = float(area(patch))
    max_dist = 0.0
    cell_pair = (None, None)

    try:
        hull = ConvexHull(patch)
        edge = list(np.array(patch)[hull.vertices])
    except:
        edge = patch

    for i, cell1 in enumerate(edge):
        for j, cell2 in enumerate(edge[i+1:]):
            squared_dist = squared_toroidal_dist(cell1, cell2, world_size)
            if squared_dist > max_dist:
                max_dist = squared_dist
                cell_pair = (cell1, cell2)

    radius = sqrt(max_dist)/2.0  # only take sqrt once

    if radius == 0:
        # This is a 1-cell patch - manually return 0
        return 0

    if formula:
        return 1-(patch_area/((radius**2)*pi))

    center = ((cell_pair[0][0]+cell_pair[1][0])/2.0,
              ((cell_pair[0][1]+cell_pair[1][1])/2.0))

    # Calculating area of circumscrbing circle
    # by brute force. Turns out that this is the
    # Gauss circle problem, which is solved by an
    # infinite sum, so brute force will be more
    # precise

    circle_area = 0.0
    x_floor = int(floor(center[0]-radius))
    x_ceil = int(ceil(center[0]+radius)+1)

    y_floor = int(floor(center[1]-radius))
    y_ceil = int(ceil(center[1]+radius)+1)

    for x in range(x_floor, x_ceil):
        for y in range(y_floor, y_ceil):
            if dist((x, y), center) <= radius:
                circle_area += 1

    return 1 - (patch_area/circle_area)