Color Object Based on 3 Nearest Objects

Here is a working version of your code. The code is self explanatory but I added comments where I felt clarification was needed. It is also definitely not an optimized algorithm, since it loops through all the petals for every triangle that must be colored.

import matplotlib.pyplot as plt
from shapely.geometry import Point, LineString
from shapely.ops import unary_union, polygonize
from matplotlib.pyplot import cm
import numpy as np
from colorir import *


# Color format used by matplotlib for input and output of colors
mpl_cfmt = ColorFormat(sRGB, max_rgb=1, include_a=True, round_to=-1)
# Color format used for the linear interpolation, aka mixing
# Try changing HCLuv for CIELuv, HSV or other color spaces and see which you like the best
mix_cfmt = ColorFormat(HCLuv)
def mix_colors(c1, c2, c3):
    # Interpret input colors with matplotlab's color format
    c1, c2, c3 = mpl_cfmt.format(c1), mpl_cfmt.format(c2), mpl_cfmt.format(c3)
    # Converts colors to color format used for mixing
    c1, c2, c3 = np.array(mix_cfmt.format(c1)), np.array(mix_cfmt.format(c2)), np.array(mix_cfmt.format(c3))
    # Mix the three colors
    # Because we havent implemented a "shortest pathway" algorithm, mixing in hue based color systems
    # may prove problematic. But I just thought that HCLuv looked the best and it hasn't caused issues
    # so I stuck with it
    final_c = (c1 + c2 + c3) / 3
    # Convert back to matplotlib's color format
    return mpl_cfmt.format(mix_cfmt.format(final_c))


# Given a reference polygon and a set of petals, select the 3 closest ones and get its colors
def find_closest_poly_colors(reference, poly_color):
    ordered = sorted(poly_color, key=lambda tup: reference.centroid.distance(tup[0].centroid))
    return ordered[0][1], ordered[1][1], ordered[2][1]


def plot_coords(coords, color):
    pts = list(coords)
    x, y = zip(*pts)
    plt.plot(x,y, color='k', linewidth=1)
    plt.fill_between(x, y, facecolor=color)


def plot_polys(polys, color):
    for poly, color in zip(polys, color):
        plot_coords(poly.exterior.coords, color)

x = 0
y = 0
h = 1.73205080757

points = [Point(x, y),
          Point((x + 2), y),
          Point((x - 2), y),
          Point((x + 1), (y + h)),
          Point((x - 1), (y + h)),
          Point((x + 1), (y - h)),
          Point((x - 1), (y - h)),
          Point((x + 3), h),
          Point((x - 3), h),
          Point((x + 3), -h),
          Point((x - 3), -h),
          Point((x + 2), (h + h)),
          Point((x - 2), (h + h)),
          Point((x + 2), (-h + -h)),
          Point((x - 2), (-h + -h)),
          Point((x + 4), y),
          Point((x - 4), y),
          Point(x, (h + h)),
          Point(x, (-h + -h)),
          Point((x + 4), (h + h)),
          Point((x - 4), (h + h)),
          Point((x + 4), (-h + -h)),
          Point((x - 4), (-h + -h)),
          Point((x + 1), (h + h + h)),
          Point((x - 1), (h + h + h)),
          Point((x + 1), (-h + -h + -h)),
          Point((x - 1), (-h + -h + -h)),
          Point((x + 5), h),
          Point((x - 5), h),
          Point((x + 5), -h),
          Point((x - 5), -h)]


circles = []
for point in points:
    circles.append(point.buffer(2))



rings = [LineString(list(pol.exterior.coords)) for pol in circles]
union = unary_union(rings)
result_polys = [geom for geom in polygonize(union)]

threshold = 0.01
filtered_polys = [polygon for polygon in result_polys if polygon.area > threshold]

complete_polys = [polygon for polygon in filtered_polys if (polygon.centroid.x**2 + polygon.centroid.y**2 < 4**2)]

print("total polygons = " + str(len(result_polys)))
print("filtered polygons = " + str(len(filtered_polys)))
print("complete polygons = " + str(len(complete_polys)))


fig = plt.figure()
ax = fig.add_subplot()
fig.subplots_adjust(top=0.85)

limit = 0.66
petals = [polygon for polygon in complete_polys if polygon.area < limit]

triangles = [polygon for polygon in complete_polys if polygon.area > limit]

colors1 = cm.viridis(np.linspace(0, 1, len(triangles)))
poly_color = list(zip(triangles, colors1))
colors2 = []
for poly in petals:
    closest = find_closest_poly_colors(poly, poly_color)
    colors2.append(mix_colors(*closest))


plot_polys(petals, color=colors2)
plot_polys(triangles, color=colors1)

ax.set_aspect('equal')
plt.show()

I think that maybe in your exemple you mistyped:

petals = [polygon for polygon in complete_polys if polygon.area > limit]

For:

# Note the < sign
petals = [polygon for polygon in complete_polys if polygon.area < limit]

But I left that unchanged since it wasn't part of your question.

Lastly, you can browse a list of color systems that can be used for interpolation in colorir here. I experimented with a few of them and HCLuv seemed the best for me, but feel free to try them out yourself.

Best regards.