Meshes with Python & Blender: Icospheres

In the third part of the series we get into mak­ing Icospheres. In oth­er words gen­er­at­ing icosa­he­dra, sub­di­vid­ing and refin­ing them to spheres. To top it off, we’ll also look at two ways of set­ting a mesh’s shad­ing to smooth.

Tutorial Series

What is a icosahedron

An icosa­he­dron is a poly­he­dron with 20 faces. There are sev­er­al kinds of icosa­he­dra. However, to make icos­pheres we’ll only be look­ing at con­vex reg­u­lar icosa­he­dra (also the most famous kind). You can find more about them and their prop­er­ties on Wikipedia

So, why icos­pheres? Icospheres have a more even dis­tri­b­u­tion of geom­e­try than UV Spheres. Deforming UV Spheres often gives strange results near the poles due to the high­er den­si­ty of geom­e­try, while icos­pheres give a more even and organ­ic look­ing result. On top of that icos­pheres are asym­met­ri­cal which helps sell an organ­ic deformation.

This tuto­r­i­al is based on the orig­i­nal icosa­he­dron code from Andreas Kahler adapt­ed to Python 3 and Blender.

Setting it up

I bet you already know this by now. Let’s start import­ing and then move on to our usu­al scaffolding.

import bpy
from math import sqrt

# -----------------------------------------------------------------------------
# Settings
name = 'Icosomething'
scale = 1
subdiv = 5

# -----------------------------------------------------------------------------
# Add Object to Scene

mesh = bpy.data.meshes.new(name)
mesh.from_pydata(verts, [], faces)

obj = bpy.data.objects.new(name, mesh)
bpy.context.scene.collection.objects.link(obj)

bpy.context.view_layer.objects.active = obj
obj.select = True

In the set­tings, subdiv will con­trol how many times to sub­di­vide the mesh and scale will be a sim­ple uni­form scale para­me­ter much like the ones in the pre­vi­ous tuto­r­i­al. Setting subdiv to zero will cre­ate a icosa­he­dron (instead of an icos­phere). Note that a subdiv val­ue of 9 will result in a mesh with over 5 mil­lion faces. You may want to stay below that val­ue depend­ing on your hard­ware (or how bad­ly you want to crash Blender!).

Putting the sphere in ico-sphere

Simply sub­di­vid­ing the icosa­he­dron will only get us a refined icosa­he­dron. We need to make sure the ver­tices come togeth­er in a way that resem­bles a sphere.

To make this hap­pen we need to make sure the ver­tices we add lie on the unit sphere. The unit sphere is an “imag­i­nary” sphere of radius 1. We can deter­mine the posi­tion of each point (ver­tex) on the unit sphere with a sim­ple for­mu­la and then fix the coor­di­nates to it. More on Wikipedia

To do this we’ll have a vertex() func­tion that fix to the unit sphere (and does the scal­ing too).

def vertex(x, y, z):
    """ Return vertex coordinates fixed to the unit sphere """

    length = sqrt(x**2 + y**2 + z**2)

    return [(i * scale) / length for i in (x,y,z)]

Make the base Icosahedron

Now that we know the verts are falling on the unit sphere we can go ahead and cre­ate the base icosa­he­dron. Like the cube before, the eas­i­est way is to input the ver­tices and faces manually.

One of the ways to build an icosa­he­dron is to con­sid­er its ver­tices as the cor­ners of three orthog­o­nal gold­en planes. They are gold­en planes because they fol­low the gold­en ratio. The ver­tices of these planes lie on the coor­di­nates (0, ±1, ±φ), (±φ, 0, ±1) and (±1, ±φ, 0). Note that the let­ter φ (phi) rep­re­sents the gold­en ratio val­ue, while ± means “neg­a­tive or positive”.

These com­bi­na­tions result in 12 ver­tices, which cre­ate 20 equi­lat­er­al tri­an­gles with 5 tri­an­gles meet­ing at each ver­tex. Check out the fol­low­ing dia­gram (points are numer­at­ed for the faces list below).

Damn, that’s a lot of math! However it’s all fair­ly straight for­ward (and repet­i­tive) when you think about it. A math jock would prob­a­bly find this bor­ing. Here’s the code to make the icosahedron.

# --------------------------------------------------------------
# Make the base icosahedron

# Golden ratio
PHI = (1 + sqrt(5)) / 2

verts = [
          vertex(-1,  PHI, 0),
          vertex( 1,  PHI, 0),
          vertex(-1, -PHI, 0),
          vertex( 1, -PHI, 0),

          vertex(0, -1, PHI),
          vertex(0,  1, PHI),
          vertex(0, -1, -PHI),
          vertex(0,  1, -PHI),

          vertex( PHI, 0, -1),
          vertex( PHI, 0,  1),
          vertex(-PHI, 0, -1),
          vertex(-PHI, 0,  1),
        ]


faces = [
         # 5 faces around point 0
         [0, 11, 5],
         [0, 5, 1],
         [0, 1, 7],
         [0, 7, 10],
         [0, 10, 11],

         # Adjacent faces
         [1, 5, 9],
         [5, 11, 4],
         [11, 10, 2],
         [10, 7, 6],
         [7, 1, 8],

         # 5 faces around 3
         [3, 9, 4],
         [3, 4, 2],
         [3, 2, 6],
         [3, 6, 8],
         [3, 8, 9],

         # Adjacent faces
         [4, 9, 5],
         [2, 4, 11],
         [6, 2, 10],
         [8, 6, 7],
         [9, 8, 1],
]

Strategy for subdivision

We can grab a tri­an­gle and split each edge cre­at­ing three tri­an­gles in its place. Basically turn­ing tri­an­gles into lit­tle tri­forces. Note that when I say split I’m not talk­ing about actu­al­ly run­ning an oper­a­tor and cut­ting the edge, but rather plac­ing a new ver­tex in the mid­dle of the oth­er orig­i­nal two and mak­ing three new faces.

However if we went around split­ting all edges we would quick­ly run into the same edges we have already split. This would result in a lot of dupli­cat­ed verts and a headache when try­ing to build the faces.

To pre­vent that, let’s keep a list of the edges we have already split (a cache), and check it before split­ting. This cache will be a dic­tio­nary. The keys will be the index of the ver­tices, ordered from small­er to greater. That way the key will remain the same, no mat­ter how we pass the edge’s vertices.

middle_point_cache = {}

def middle_point(point_1, point_2):
    """ Find a middle point and project to the unit sphere """

    # We check if we have already cut this edge first
    # to avoid duplicated verts
    smaller_index = min(point_1, point_2)
    greater_index = max(point_1, point_2)

    key = '{0}-{1}'.format(smaller_index, greater_index)

    if key in middle_point_cache:
        return middle_point_cache[key]

    # If it's not in cache, then we can cut it
    vert_1 = verts[point_1]
    vert_2 = verts[point_2]
    middle = [sum(i)/2 for i in zip(vert_1, vert_2)]

    verts.append(vertex(*middle))

    index = len(verts) - 1
    middle_point_cache[key] = index

    return index

The mid­dle ver­tex is cal­cu­lat­ed by adding the coor­di­nates from both ver­tices and divid­ing them by 2. Finally we add it to the cache and return the index to make the faces list.

Subdividing

With the middle_point() func­tion done we can move over to loop­ing and mak­ing the subdivisions.

At each sub­di­vi­sions lev­el we cre­ate a new emp­ty list of faces, and at the end we replace the orig­i­nal faces list with the new one. Then we go through each face, find­ing the mid­dle point for it’s three edges, stor­ing the indices and cre­at­ing 4 new faces from them (remem­ber the dia­gram above).

# Subdivisions
# --------------------------------------------------------------

for i in range(subdiv):
    faces_subdiv = []

    for tri in faces:
        v1 = middle_point(tri[0], tri[1])
        v2 = middle_point(tri[1], tri[2])
        v3 = middle_point(tri[2], tri[0])

        faces_subdiv.append([tri[0], v1, v3])
        faces_subdiv.append([tri[1], v2, v1])
        faces_subdiv.append([tri[2], v3, v2])
        faces_subdiv.append([v1, v2, v3])

    faces = faces_subdiv

Make it smooth

We now have a mesh that approx­i­mates a sphere quite well, but still looks faceted. Time to smooth it out.

Smooth shad­ing is an attribute of faces. So, to set a whole mesh as smooth we need to enable smooth shad­ing on all of it’s faces. You can do this using the same oper­a­tor we use when we click the Set smooth but­ton in the UI:

bpy.ops.object.shade_smooth()

This will work just fine for this script because the con­text is right for the oper­a­tor. In oth­er cas­es, you might find the oper­a­tor refus­es to run due to an “incor­rect con­text”. context in Blender is a sort of god-vari­able that con­tains infor­ma­tion about the cur­rent state of the appli­ca­tion. This includes things like mouse cur­sor posi­tion, the active area and more. You can over­ride the con­text when call­ing an oper­a­tor but cur­rent­ly there’s no easy way of know­ing what each oper­a­tor wants to see in the context.

Luckily there’s anoth­er way to do this the “low-lev­el” way by set­ting each face in a loop.

for face in mesh.polygons:
    face.use_smooth = True

In the world of Blender script­ing “Low Level” refers to skip­ping oper­a­tors and going straight to the meth­ods and attrib­ut­es of objects. from_pydata() is anoth­er exam­ple of low-lev­el work.

The ben­e­fits of low-lev­el is that it doesn’t depend on con­text, it’s often more flex­i­ble and saves you the over­head of going through the oper­a­tors sys­tem. In this case you could also decide to apply smooth only on some faces.

Final Code

import bpy
from math import sqrt

# -----------------------------------------------------------------------------
# Settings

scale = 1
subdiv = 5
name = 'Icosomething'


# -----------------------------------------------------------------------------
# Functions

middle_point_cache = {}


def vertex(x, y, z):
    """ Return vertex coordinates fixed to the unit sphere """

    length = sqrt(x**2 + y**2 + z**2)

    return [(i * scale) / length for i in (x,y,z)]


def middle_point(point_1, point_2):
    """ Find a middle point and project to the unit sphere """

    # We check if we have already cut this edge first
    # to avoid duplicated verts
    smaller_index = min(point_1, point_2)
    greater_index = max(point_1, point_2)

    key = '{0}-{1}'.format(smaller_index, greater_index)

    if key in middle_point_cache:
        return middle_point_cache[key]

    # If it's not in cache, then we can cut it
    vert_1 = verts[point_1]
    vert_2 = verts[point_2]
    middle = [sum(i)/2 for i in zip(vert_1, vert_2)]

    verts.append(vertex(*middle))

    index = len(verts) - 1
    middle_point_cache[key] = index

    return index


# -----------------------------------------------------------------------------
# Make the base icosahedron

# Golden ratio
PHI = (1 + sqrt(5)) / 2

verts = [
          vertex(-1,  PHI, 0),
          vertex( 1,  PHI, 0),
          vertex(-1, -PHI, 0),
          vertex( 1, -PHI, 0),

          vertex(0, -1, PHI),
          vertex(0,  1, PHI),
          vertex(0, -1, -PHI),
          vertex(0,  1, -PHI),

          vertex( PHI, 0, -1),
          vertex( PHI, 0,  1),
          vertex(-PHI, 0, -1),
          vertex(-PHI, 0,  1),
        ]


faces = [
         # 5 faces around point 0
         [0, 11, 5],
         [0, 5, 1],
         [0, 1, 7],
         [0, 7, 10],
         [0, 10, 11],

         # Adjacent faces
         [1, 5, 9],
         [5, 11, 4],
         [11, 10, 2],
         [10, 7, 6],
         [7, 1, 8],

         # 5 faces around 3
         [3, 9, 4],
         [3, 4, 2],
         [3, 2, 6],
         [3, 6, 8],
         [3, 8, 9],

         # Adjacent faces
         [4, 9, 5],
         [2, 4, 11],
         [6, 2, 10],
         [8, 6, 7],
         [9, 8, 1],
        ]


# -----------------------------------------------------------------------------
# Subdivisions

for i in range(subdiv):
    faces_subdiv = []

    for tri in faces:
        v1 = middle_point(tri[0], tri[1])
        v2 = middle_point(tri[1], tri[2])
        v3 = middle_point(tri[2], tri[0])

        faces_subdiv.append([tri[0], v1, v3])
        faces_subdiv.append([tri[1], v2, v1])
        faces_subdiv.append([tri[2], v3, v2])
        faces_subdiv.append([v1, v2, v3])

    faces = faces_subdiv


# -----------------------------------------------------------------------------
# Add Object to Scene

mesh = bpy.data.meshes.new(name)
mesh.from_pydata(verts, [], faces)

obj = bpy.data.objects.new(name, mesh)
bpy.context.scene.collection.objects.link(obj)

bpy.context.view_layer.objects.active = obj
obj.select = True


# -----------------------------------------------------------------------------
# Smoothing

#bpy.ops.object.shade_smooth()

for face in mesh.polygons:
    face.use_smooth = True

Wrap up

That com­pletes the third part of this series. If you are inter­est­ed in mak­ing spheres or apply points/objects in a spher­i­cal way con­sid­er read­ing more about the unit sphere and how to apply it along normals.

Things you can do for yourself:

  • Optimize it (Hint: you don’t need to store the key as a string)
  • Bring the rota­tion and trans­la­tion matri­ces from pre­vi­ous part
  • Remove the scale set­ting and use a matrix instead

Next time we’ll go back to cubes and learn how to make a round cube, apply mod­i­fiers and more. Do you have any ques­tions, or sug­ges­tions for the next tuto­ri­als? Leave a com­ment below!

All the posts you can read
TutorialsBlender, Meshes, Python19.10.2020