iso4.py

# Nicolas Seriot
# 2021-09-17
# 2021-09-18
# https://gist.github.com/nst/032a61feb1fc60a74e13d2ea994763c0
# Thread: https://twitter.com/nst021/status/1437889678947110912
# Typical output: https://seriot.ch/visualization/iso4.png

import cairo
import random
import numpy as np

import cProfile, pstats, io
from pstats import SortKey

MARGIN = 50
X_MAX, Y_MAX, Z_MAX = 60, 40, 80 # space size
DW, DH = 32, 16 # diamond size
DW2, DH2 = int(DW/2), int(DH/2)

def draw_surface(c, points, color, vertices = []):
    
    c.save()

    c.set_source_rgb(*color)
    
    [c.line_to(*p) for p in points]

    c.fill()
    
    c.set_source_rgb(0,0,0)
    c.set_line_width(1)
    
    for pts in vertices:
        p1, p2 = pts
        c.move_to(*p1)
        c.line_to(*p2)
        c.stroke()
    
    c.restore()

def draw_grid(c):

    for x in range(X_MAX):
        for y in range(Y_MAX):        
            draw_cube(c, x, y, 0, flat=True)

def model_to_canvas(x, y):

    o_x = DW2 * Y_MAX
    o_y = DH2
    
    return (o_x + DW2*(x-y), o_y + DH2*(x+y))

def draw_cube(c, x, y, z,
              nx = False, nx_ = False, ny = False, ny_ = False, nz = False, nz_ = False,
              nxy_=False, nxz = False, nx_z_ = False, nyz = False, ny_z_ = False,
              flat=False):
    
    x, y = model_to_canvas(x, y)
    
    c.save()

    c.translate(x, y)

    """
      5
    2   6
      3      - DH
    1   7
      4
      | z_offset
    """
    
    z_offset = DH * z

    p1 = (0,   z_offset + DH2)
    p2 = (0,   z_offset + DH2 + DH)
    p3 = (DW2, z_offset + DH)
    p4 = (DW2, z_offset + 0)
    p5 = (DW2, z_offset + 2*DH)
    p6 = (DW,  z_offset + DH2 + DH)
    p7 = (DW,  z_offset + DH2)
    
    #COLOR_RIGHT = (0.5,0.5,0.5)
    #COLOR_TOP = (1,1,1)
    #COLOR_LEFT = (0,0,0)
    
    z_ratio = 0.4 + z*0.6/Z_MAX
    
    COLOR_TOP   = (1 * z_ratio,   0.5 * z_ratio, 0.5 * z_ratio)
    COLOR_LEFT  = (0.6 * z_ratio, 0,             0)
    COLOR_RIGHT = (1 * z_ratio,   0,             0)

    vertices = []
    
    if flat:
        vertices.append((p1,p3))
        vertices.append((p3,p7))
        vertices.append((p7,p4))
        vertices.append((p4,p1))
        draw_surface(c, [], (1,1,1), vertices)
    else:
    
        if not nx:#
            vertices.append((p5, p6))
            vertices.append((p6, p7))
        if not ny:
            vertices.append((p2, p5))
            vertices.append((p2, p1))
        if not nz_:
            vertices.append((p4, p7))
            vertices.append((p1, p4))
        if not nx_ and not ny_:
            vertices.append((p3, p4))
        if not nz and not ny_:
            vertices.append((p6, p3))
        if not nz and not nx_:
            vertices.append((p2, p3))
        if nxy_:
            vertices.append((p6, p7))
        if nx_z_:
            vertices.append((p1, p4))
        if nxz:
            vertices.append((p5, p6))
        if nyz:
            vertices.append((p2, p5))
        if ny_z_:
            vertices.append((p4, p7))
        
        draw_surface(c, [p2, p5, p6, p3], COLOR_TOP, vertices)
        draw_surface(c, [p1, p2, p3, p4], COLOR_LEFT, vertices)
        draw_surface(c, [p3, p6, p7, p4], COLOR_RIGHT, vertices)
    
    c.restore()

def fill_model(m, density):

    for index, x in np.ndenumerate(m):
        m[index] = random.uniform(0, 1) < density

def visibility_matrix(m):

    X,Y,Z = m.shape

    # nothing is visible except the three visible faces of the space

    v = np.full(m.shape, False)
    
    for x in range(X):
        for y in range(Y):
            v[x][y][Z-1] = True

    for x in range(X):
        for z in range(Z):
            v[x][0][z] = True

    for y in range(Y):
        for z in range(Z):
            v[0][y][z] = True

    # iterating from user's standpoint
    
    for x in range(X):
        for y in range(Y):
            for z in range(Z)[::-1]:

                # if not visible
                # no need to update visibility of "back" cube
                # continue to the next cube
                if not v[x][y][z]:
                    continue

                # if m is empty
                # "back" cube becomes visible
                if not m[x][y][z]:
                    if x < (X-1) and y < (Y-1) and z > 0:
                        v[x+1][y+1][z-1] = True
        
    return v

def draw_model(c, m):
    
    pr = cProfile.Profile()
    pr.enable()




    X,Y,Z = m.shape
    
    v = visibility_matrix(m)
        
    for x in range(X)[::-1]:
        for y in range(Y)[::-1]:
            for z in range(Z_MAX):
                if m[x][y][z]:
                    
                    if not v[x][y][z]:
                        #print("-- skip", x, y, z)
                        continue
                    
                    nx  = m[x+1][y][z] if (x+1) < X else False
                    nx_ = m[x-1][y][z] if (x-1) >= 0 else False

                    ny  = m[x][y+1][z] if (y+1) < Y else False
                    ny_ = m[x][y-1][z] if (y-1) >= 0 else False

                    nz  = m[x][y][z+1] if (z+1) < Z else False
                    nz_ = m[x][y][z-1] if (z-1) >= 0 else False
                    
                    nxy_  = m[x+1][y-1][z] if (x+1) < X  and (y-1) >= 0 else False

                    nxz   = m[x+1][y][z+1] if (x+1) < X  and (z+1) < Z  else False
                    nx_z_ = m[x-1][y][z-1] if (x-1) >= 0 and (z-1) >= 0 else False

                    nyz   = m[x][y+1][z+1] if (y+1) < Y  and (z+1) < Z  else False
                    ny_z_ = m[x][y-1][z-1] if (y-1) >= 0 and (z-1) >= 0 else False
                    
                    draw_cube(c,x,y,z,nx,nx_,ny,ny_,nz,nz_,nxy_,nxz,nx_z_,nyz,ny_z_)





    pr.disable()
    s = io.StringIO()
    sortby = SortKey.CUMULATIVE
    ps = pstats.Stats(pr, stream=s).sort_stats(sortby)
    ps.print_stats()
    print(s.getvalue())

def fill_shape_1(m, o):

    for x in range(-2, 3):    
        m[o+x][o][o] = True

    for y in range(-2, 3):    
        m[o][o+y][o] = True

    for z in range(-2, 3):    
        m[o][o][o+z] = True

def fill_shape_2(m, o):

    fill_shape_1(m, o)

    m[o-1][o-1][o] = True
    m[o-1][o+1][o] = True
    m[o+1][o-1][o] = True

def fill_shape_3(m, o):
    for x in range(-1, 2):    
        for y in range(-1, 2):    
            for z in range(-1, 2):
                m[o+x][o+y][o+z] = True

def main():

    #random.seed(0)

    m = np.full((X_MAX, Y_MAX, Z_MAX), False)
    
    fill_model(m, 0.5)
    
    #fill_shape_1(m, 2)
    #fill_shape_2(m, 6)
    #fill_shape_3(m, 5)
        
    surface = cairo.ImageSurface(cairo.FORMAT_ARGB32,
                                 DW + (X_MAX+Y_MAX)*DW2 + MARGIN*2,
                                 DH + (X_MAX+Y_MAX)*DH2 + Z_MAX*DH + MARGIN*2)
    c = cairo.Context(surface)
    
    cm = cairo.Matrix(yy=-1, y0=surface.get_height())
    c.transform(cm)
    
    c.translate(MARGIN, MARGIN)
    
    c.set_antialias(cairo.ANTIALIAS_NONE)
    
    #background
    c.set_source_rgb(1,1,1)
    c.paint()
    
    #pen
    c.set_source_rgb(0,0,0)
    c.set_line_width(1)
    
    draw_grid(c)
    draw_model(c, m)
    
    #draw_cube(c, 1, 0, 0)
    
    surface.write_to_png("iso4.png")

if __name__ == "__main__":
        
    main()